`smooth.spline`が長い（周期的な）時系列に深刻に劣る

Question

非常に長くてノイズの多いデータをRで滑らかにしたいと思います。しかし、非常に周期的なデータの場合、すぐに使用できるsmooth.spline()はすぐに壊れます平滑化されたデータが鳴り始める。

（またはノイズなし）コサイン時系列を考えてみましょう

t <- seq(0,100*2*pi,length.out=3000) 
y <- cos(t)# + rnorm(length(t), 0,0.05) 

y100_s <- smooth.spline(y)$y 

plot(y~t, type="l") 
lines(y100_s~t, col="blue") 

我々はフィット感が（鳴っていても、より低い周波数で

# rms increases as points are added to smooth.spline 
rms <- sapply(seq(250,3000,by=250), function(i) 
    sqrt(mean((y[1:i] - smooth.spline(y[1:i])$y)^2))) 

plot(rms) 

、 smooth.spline()に複数の値を追加することの影響を調べることができますオプション）。

t <- seq(0,50*2*pi,length.out=3000) 
y <- cos(t)# + rnorm(length(t), 0,0.05) 

y50_s <- smooth.spline(y)$y 

require(pracma) 

peaks <- list(findpeaks(y50_s),findpeaks(-y50_s)) 

plot(y~t, type="l") 
lines(y50_s~t, col="red") 

lines(peaks[[1]][,1]~t[peaks[[1]][,2]], type="l") 
lines(-peaks[[2]][,1]~t[peaks[[2]][,2]], type="l") 

ビットのために探索した後、この動作はスパー引数の関数であるように見えますが、私は効果を排除するために十分小さな値に設定することはできません。これはスプラインフィットの明白な結果であり、すぐに使える方法に頼っていないという欠点がありますが、私はいくつかの洞察に感謝します。 smooth.spline()で指定できるコントロール、または平滑化のための代替の推奨事項/戦略がありますか？

Answer 1

定期的な信号を常にフィッティングしているかどうかわかりません。その場合は、mgcv::gamの周期的スプラインを使用する方がはるかに優れています。しかし、今はこの問題を忘れてみましょう。

データが頻繁に振動している場合は、十分な数のノット、つまり適切なノット密度を選択する必要があります。そうしないとオーバースムージング（アンダーフィット）が発生します。

t <- seq(0, 100 * 2 * pi, length.out = 3000) 
y <- cos(t) # + rnorm(length(t), 0, 0.05) 
fit <- smooth.spline(t, y) 

あなたはn = 3000データポイントを持っている：

はあなたの例を見てください。デフォルトではn > 49の場合、smooth.splineはデータよりもはるかに少ないノット数を使用します。正確にはサービスルーチン.nknots.smsplによって選択されます。しかし、これには最適性の正当性はない。これが妥当かどうかは、あなた次第です。

これは194ノットのみを使用し、モデルはペナルティ効果なしで194自由度で終了します。私が先に言ったように、あなただけのアンダーフィッティングで終わる：

plot(t, y, type = "l", col = "gray") 
lines(fit, col = 2) 

は、理想的には、罰せ回帰は、ノットの数よりも実質的に小さい自由度で終わります。当初のペナルティなしの回帰から生じる過剰フィット問題を修正するためにペナルティが使用されることは、しばしば忘れられる。ペナルティ効果が見られない場合でも、オリジナルのペナルティなしモデルはアンダーフィッティングデータですので、オーバーフィット状態になるまでノット数を増やしてください。これについて考えるのが怠け者なら、all.knots = TRUEと設定してください。単変量平滑化スプラインは、計算コストがO(n)で非常に安いです。すべてのデータをノットとして使用しても、効率の問題にはなりません。

fit <- smooth.spline(t, y, all.knots = TRUE) 

length(fit$fit$nk) - 2L 
# [1] 3000 

fit$df 
# [1] 3000 

ああ、私たちはまだペナルティの影響を見ませんでした、なぜですか？私たちは騒々しいデータを持っていないからです。 yにノイズを加えなかったので、すべてのノットを使って補間を行います。yに雑音を加えて、私がペナルティについて何を説明しているかを本当に理解してください。

set.seed(0) 
t <- seq(0, 100 * 2 * pi, length.out = 3000) 
y <- cos(t) + rnorm(length(t), 0, 0.05) 

fit <- smooth.spline(t, y, all.knots = TRUE) 

length(fit$fit$nk) 
# [1] 3002 

fit$df 
# [1] 705.0414 

3000と比べると705の方が小さいことに注意してください。スプラインを見てみましたか？

plot(t, y, type = "l", col = "gray") 
lines(fit, col = 2) 

アンダーフィットもオーバーフィットもありません。ペナルティングはバイアスと分散の最適なトレードオフをもたらします。