2016-04-05 6 views
0

の印刷をキャッチ。この関数にはeigchkと呼ばれる別の関数が含まれ、係数行列が特異行列かどうかがチェックされます。ここで は、私は特定の機能<code>fRegress</code>パッケージ<code>fda</code>を使用しています機能

はそれを書いたパッケージの所有者(J. O.ラムゼイ、ジャイルズフッカー、およびスペンサーグラーヴ)としての機能です。

eigchk <- function(Cmat) { 

    # check Cmat for singularity 

    eigval <- eigen(Cmat)$values 
    ncoef <- length(eigval) 
    if (eigval[ncoef] < 0) { 
     neig <- min(length(eigval),10) 
     cat("\nSmallest eigenvalues:\n") 
     print(eigval[(ncoef-neig+1):ncoef]) 
     cat("\nLargest eigenvalues:\n") 
     print(eigval[1:neig]) 
     stop("Negative eigenvalue of coefficient matrix.") 
    } 
    if (eigval[ncoef] == 0) stop("Zero eigenvalue of coefficient matrix.") 
    logcondition <- log10(eigval[1]) - log10(eigval[ncoef]) 
    if (logcondition > 12) { 
     warning("Near singularity in coefficient matrix.") 
     cat(paste("\nLog10 Eigenvalues range from\n", 
       log10(eigval[ncoef])," to ",log10(eigval[1]),"\n")) 
    } 
} 

logconditionは12とプリント固有値のその後の範囲よりも大きい場合は、もし最後の条件チェックを見ることができるように。

次のコードは、粗さpennaltyと正則化のuseageを実現します。このコードは、「RとMatlabを使用した機能的データ分析」という本から引用されています。

annualprec = log10(apply(daily$precav,2,sum)) 
tempbasis =create.fourier.basis(c(0,365),65) 
tempSmooth=smooth.basis(day.5,daily$tempav,tempbasis) 
tempfd =tempSmooth$fd 

templist = vector("list",2) 
templist[[1]] = rep(1,35) 
templist[[2]] = tempfd 

conbasis = create.constant.basis(c(0,365)) 
betalist = vector("list",2) 
betalist[[1]] = conbasis 


SSE = sum((annualprec - mean(annualprec))^2) 
Lcoef = c(0,(2*pi/365)^2,0) 
harmaccelLfd = vec2Lfd(Lcoef, c(0,365)) 
betabasis = create.fourier.basis(c(0, 365), 35) 
lambda = 10^12.5 
betafdPar = fdPar(betabasis, harmaccelLfd, lambda) 
betalist[[2]] = betafdPar 

annPrecTemp = fRegress(annualprec, templist, betalist) 
betaestlist2 = annPrecTemp$betaestlist 
annualprechat2 = annPrecTemp$yhatfdobj 

SSE1.2 = sum((annualprec-annualprechat2)^2) 
RSQ2 = (SSE - SSE1.2)/SSE 
Fratio2 = ((SSE-SSE1.2)/3.7)/(SSE1/30.3) 

resid = annualprec - annualprechat2 
SigmaE. = sum(resid^2)/(35-annPrecTemp$df) 
SigmaE = SigmaE.*diag(rep(1,35)) 
y2cMap = tempSmooth$y2cMap 

stderrList = fRegress.stderr(annPrecTemp, y2cMap, SigmaE) 

betafdPar = betaestlist2[[2]] 
betafd = betafdPar$fd 
betastderrList = stderrList$betastderrlist 
betastderrfd = betastderrList[[2]] 

ペナルティ係数として、著者は特定のlambdaを使用します。 次のコードは、適切な `lambda。 loglamの値とクロスバリデーションを変更することにより

loglam = seq(5,15,0.5) 
nlam = length(loglam) 
SSE.CV = matrix(0,nlam,1) 
for (ilam in 1:nlam) { 
    lambda = 10ˆloglam[ilam] 
    betalisti = betalist 
    betafdPar2 = betalisti[[2]] 
    betafdPar2$lambda = lambda 
    betalisti[[2]] = betafdPar2 
    fRegi = fRegress.CV(annualprec, templist, 
       betalisti) 
    SSE.CV[ilam] = fRegi$SSE.CV 
} 

私はlambda最善をequaireすると仮定し、まだloglamの長さが大きいか、その値にある場合singulrityに係数行列を導きます。私はすでに上述したように機能eigchkによって作成された

Log10 Eigenvalues range from 
-5.44495317739048 to 6.78194912518214 

:私は、次のメッセージを受け取ります。 今、私の質問は、このいわゆる警告をキャッチする方法はありますか?キャッチとは、これが起こったときに私に警告する機能や方法を意味し、loglamの値を調整することができます。メッセージのこのプリントの横にある機能に実際の警告定義がないので、アイデアがなくなりました。

は、ご提案のためにあなたにすべてをたくさんありがとうございます。

+0

メッセージの場合、 'tryCatch'で捕まえることができますが、その後に何をしたいのか分かりません – rawr

+0

@rawr私の目標は、 'loglam'のベクトル値を調整する関数ですこの印刷物にしかし問題は、このプリントが何らかの追加の通知であり、実際の警告ではないということです。 – user3185925

答えて

0

loglamに潜在的な問題があることを警告する場合は、tryおよびtryCatchの機能を参照することをお勧めします。次に、警告条件が満たされた場合に実装する動作を定義できます。

警告の出力を保存したいだけであれば(質問のタイトルから推測されるかもしれませんが、必要でないかもしれません)、capture.outputを調べてみてください。

+0

私は既にtryCatch関数を見てきました。私の理解が正しければ、私の関数が出力としてwarnungを持っていたとしても、このプリントは評価を止めない何らかの追加の通知です。 – user3185925

関連する問題