R smooth.spline 擬合平滑樣條曲線

說明

將三次平滑樣條擬合到提供的數據。

用法

smooth.spline(x, y = NULL, w = NULL, df, spar = NULL, lambda = NULL, cv = FALSE,
              all.knots = FALSE, nknots = .nknots.smspl,
              keep.data = TRUE, df.offset = 0, penalty = 1,
              control.spar = list(), tol = 1e-6 * IQR(x), keep.stuff = FALSE)

.nknots.smspl(n)

參數

x	給出預測變量值的向量，或指定 x 和 y 的列表或兩列矩陣。
y	回應。如果 y 丟失或 NULL ，則假定響應由 x 指定，其中 x 為索引向量。
w	與 x 長度相同的可選權重向量；默認全部為1。
df	所需的等效自由度數(平滑矩陣的跡)。必須在 (1,n_x] 、 n_x 中唯一 x 值的數量，請參見下文。
spar	平滑參數，通常(但不一定)位於 (0,1] 中。當指定 spar 時，擬合(懲罰對數似然)準則中的二階導數平方積分的係數 \lambda 是 spar 的單調函數，請參閱下麵的詳細信息。或者，可以指定 lambda 來代替無標度 spar = s 。
lambda	如果需要，可以指定內部(取決於設計)平滑參數 \lambda 而不是 spar 。這對於交叉驗證或引導程序等重采樣算法可能是理想的。
cv	當 FALSE 時，普通留一法 ( TRUE ) 或 ‘generalized’ 交叉驗證 (GCV) ；僅當spar和df均未指定時才用於平滑參數計算；然而，它用於確定結果中的cv.crit。將其設置為 NA 以獲得加速會跳過杠杆和任何分數的評估。
all.knots	如果 TRUE ，則 x 中的所有不同點都用作結。如果 FALSE (默認)，則使用 x[] 的子集，特別是 x[j]，其中 nknots 索引在 1:n 中均勻分布，另請參閱下一個參數 nknots 。 或者，使用嚴格遞增的numeric向量，指定要使用的“all the knots”；必須重新縮放為 [0, 1]，使其對應於返回的 ans $ fit$knots 序列，而不重複邊界結。
nknots	整數或 function 給出 all.knots = FALSE 時要使用的結數。如果是函數(默認情況下)，則結數為 nknots(nx) 。默認情況下使用 .nknots.smspl() ，對於 n_x > 49 這小於 n_x ，即唯一 x 值的數量，請參閱注釋。
keep.data	邏輯指定輸入數據是否應保留在結果中。如果TRUE(默認)，則可從結果中獲得擬合值和殘差。
df.offset	允許 GCV 準則中的自由度增加df.offset。
penalty	GCV 準則中自由度的懲罰係數。
control.spar	計算平滑參數 spar 時控製根查找的命名組件的可選列表，即缺少或 NULL ，請參見下文。 請注意，這部分是實驗性的，可能會隨著一般晶石計算的改進而改變！ 低的： 下界為spar;默認為 -1.5(用於隱式默認為 0R1.4 之前的版本)。 高的： spar 的上限；默認為+1.5。 托爾： 使用的絕對精度(公差)；默認為 1e-4(以前為 1e-3)。 每股收益： 使用的相對精度；默認為 2e-8(以前為 0.00244)。 痕跡： 邏輯指示是否應該跟蹤迭代。 馬克西特： 給出最大迭代次數的整數；默認為 500。 請注意，僅在區間 [low, high] 中搜索 spar 。
tol	same-ness 的容差或 x 值的唯一性。這些值被分箱到大小為 tol 的箱中，並且落入同一箱中的值被視為相同。必須嚴格為正(且有限)。
keep.stuff	實驗性的 logical 指示結果是否應保留內部計算的額外內容。應該允許重建 X 矩陣等。
n	對於.nknots.smspl；通常是唯一 x 值(又名 n_x )的數量。

細節

x 和 y 都不允許包含缺失值或無限值。

x 向量應至少包含四個不同的值。這裏的‘Distinct’由tol控製：被視為相同的值被它們的第一個值替換，並且相應的y和w被相應地池化。

除非指定 lambda 而不是 spar ，否則使用的計算 \lambda (作為 s=spar 的函數)為 \lambda = r \cdot 256^{3 s - 1} ，其中 r = tr(X' W X) / tr(\Sigma) 、 \Sigma 是由 \Sigma_{ij} = \int B_i''(t) B_j''(t) dt 給出的矩陣， X 由 X_{ij} = B_j(x_i) 給出，W 是權重的對角矩陣(縮放後其軌跡為 n，即原始觀測數)，B_k(.) 是 k -th B-spline 。

請注意，使用這些定義 f_i = f(x_i) 和 B-spline 基礎表示 f = X c (即 c 是樣條係數的向量)，懲罰對數似然為 L = (y - f)' W (y - f) + \lambda c' \Sigma c ，因此 c 是(嶺回歸)(X' W X + \lambda \Sigma) c = X' W y 的解決方案。

如果 spar 和 lambda 缺失或 NULL ，則使用 df 的值來確定平滑程度。如果 df 也丟失，則使用留一交叉驗證(普通或 ‘generalized’ 由 cv 確定)來確定 \lambda 。

請注意，根據上述關係，spar是s = s0 + 0.0601 \cdot \log\lambda，這是故意的不同的來自 S-PLUS 的實現smooth.spline(在哪裏spar正比於\lambda)。在R的(\log \lambda)規模，改變更有意義spar線性地。

但請注意，目前對於小於 -1 或 -2 的 spar 值，結果可能變得非常不可靠。對於大於 2 左右的值也可能會發生同樣的情況。不要考慮將 spar 或控件 low 和 high 設置在這樣的安全範圍之外，除非您知道自己在做什麽！同樣，指定 lambda 而不是 spar 也很微妙，特別是 lambda 的 “safe” 值的範圍不是比例不變的，因此完全依賴於數據。

當 x 中存在重複點時，‘generalized’ 交叉驗證方法 GCV 將正常工作。然而，留一交叉驗證對於重複點的含義並不明確，並且內部代碼使用涉及省略重複點組的近似值。在這種情況下最好避免cv = TRUE。

值

帶有組件的 "smooth.spline" 類的對象

x	不同的 x 值按升序排列，請參閱上麵的“詳細信息”。
y	對應於 x 的擬合值。
w	在 x 的唯一值處使用的權重。
yin	用於唯一 y 值的 y 值。
tol	tol 參數(其默認值取決於 x )。
data	僅當keep.data = TRUE：本身是一個list，其組件x、y和w具有相同的長度。這些是原始的 (x_i,y_i,w_i), i = 1, \dots, n 值，其中 data$x 可能具有重複值，因此比上述 x 組件長；查看具體信息。
n	一個整數； (原始)樣本量。
lev	(當 cv 不是 NA 時)利用平滑矩陣的對角線值。
cv	使用的 cv 參數；即 FALSE 、 TRUE 或 NA 。
cv.crit	交叉驗證分數，‘generalized’ 或 true，具體取決於 cv 。 CV 分數通常稱為 “PRESS”(並標記為 print() )，表示“預測平方和”。請注意，這與擬合期間最小化的(CV 或 GCV)分數不同(並在 crit 中返回)，例如，在 nx < n 的情況下(其中 nx =n_x 是唯一的 x 值)。
pen.crit	懲罰標準，一個非負數；隻是(加權)殘差平方和 (RSS)， sum(.$w * residuals(.)^2) 。
crit	底層標準值最小化.Fortran常規 'sslvrg’。什麽時候df已指定，標準為3 + (tr(S_\lambda) - df)^2, 其中3 +是出於數字(和曆史)原因嗎？
df	使用等效自由度。請注意，(目前)當真實的 df 介於 1 和 2 之間時，該值可能會變得相當不精確。
spar	計算或給出的 spar 的值，除非在此處設置為 NA 時將其指定為 c(lambda = *) 。
ratio	(當上麵的 spar 不是 NA 時)，值 r ，兩個矩陣跡線的比率。
lambda	\lambda 對應 spar 的值，請參見上麵的詳細信息。
iparms	命名的整數(3)向量，其中..$ipars["iter"]給出了使用的spar計算迭代次數。
auxMat	實驗性的；當 keep.stuff 為 true 時，包含部分內部計算的 “flat” 數值向量。
fit	供 predict.smooth.spline 使用的列表，包含組件 結： 結序列(包括重複的邊界結)，縮放為 [0, 1] (通過 min 和 range )。 NK: 係數數或 ‘proper’ 結數加 2。 係數： 使用的樣條基礎的係數。 最小值，範圍： 給出 x 相應數量的數字。
call	匹配的調用。

method(class = "smooth.spline") 顯示了基於上麵lev 向量的hatvalues() 方法。

注意

唯一 x 值的數量 \mathbf{nx} = n_x 由 tol 參數確定，相當於

    nx <- length(x) - sum(duplicated( round((x - mean(x)) / tol) ))
  

默認all.knots = FALSE和nknots = .nknots.smspl，隻需要使用O({n_x}^{0.2})結而不是n_x為了n_x > 49。這降低了速度和內存需求，但不再大幅降低，因為R版本 1.5.1 僅適用於O(n_k) + O(n)其中n_k是節數。

在這種情況下，並非所有唯一的 x 值都用作結，結果是回歸樣條而不是嚴格意義上的平滑樣條，但非常接近，除非使用較小的平滑參數(或較大的 df )。

例子

require(graphics)
plot(dist ~ speed, data = cars, main = "data(cars)  &  smoothing splines")
cars.spl <- with(cars, smooth.spline(speed, dist))
cars.spl
## This example has duplicate points, so avoid cv = TRUE

lines(cars.spl, col = "blue")
ss10 <- smooth.spline(cars[,"speed"], cars[,"dist"], df = 10)
lines(ss10, lty = 2, col = "red")
legend(5,120,c(paste("default [C.V.] => df =",round(cars.spl$df,1)),
               "s( * , df = 10)"), col = c("blue","red"), lty = 1:2,
       bg = 'bisque')


## Residual (Tukey Anscombe) plot:
plot(residuals(cars.spl) ~ fitted(cars.spl))
abline(h = 0, col = "gray")

## consistency check:
stopifnot(all.equal(cars$dist,
                    fitted(cars.spl) + residuals(cars.spl)))
## The chosen inner knots in original x-scale :
with(cars.spl$fit, min + range * knot[-c(1:3, nk+1 +1:3)]) # == unique(cars$speed)

## Visualize the behavior of  .nknots.smspl()
nKnots <- Vectorize(.nknots.smspl) ; c.. <- adjustcolor("gray20",.5)
curve(nKnots, 1, 250, n=250)
abline(0,1, lty=2, col=c..); text(90,90,"y = x", col=c.., adj=-.25)
abline(h=100,lty=2); abline(v=200, lty=2)

n <- c(1:799, seq(800, 3490, by=10), seq(3500, 10000, by = 50))
plot(n, nKnots(n), type="l", main = "Vectorize(.nknots.smspl) (n)")
abline(0,1, lty=2, col=c..); text(180,180,"y = x", col=c..)
n0 <- c(50, 200, 800, 3200); c0 <- adjustcolor("blue3", .5)
lines(n0, nKnots(n0), type="h", col=c0)
axis(1, at=n0, line=-2, col.ticks=c0, col=NA, col.axis=c0)
axis(4, at=.nknots.smspl(10000), line=-.5, col=c..,col.axis=c.., las=1)

##-- artificial example
y18 <- c(1:3, 5, 4, 7:3, 2*(2:5), rep(10, 4))
xx  <- seq(1, length(y18), length.out = 201)
(s2   <- smooth.spline(y18)) # GCV
(s02  <- smooth.spline(y18, spar = 0.2))
(s02. <- smooth.spline(y18, spar = 0.2, cv = NA))
plot(y18, main = deparse(s2$call), col.main = 2)
lines(s2, col = "gray"); lines(predict(s2, xx), col = 2)
lines(predict(s02, xx), col = 3); mtext(deparse(s02$call), col = 3)

## Specifying 'lambda' instead of usual spar :
(s2. <- smooth.spline(y18, lambda = s2$lambda, tol = s2$tol))



## The following shows the problematic behavior of 'spar' searching:
(s2  <- smooth.spline(y18, control =
                      list(trace = TRUE, tol = 1e-6, low = -1.5)))
(s2m <- smooth.spline(y18, cv = TRUE, control =
                      list(trace = TRUE, tol = 1e-6, low = -1.5)))
## both above do quite similarly (Df = 8.5 +- 0.2)

作者

R implementation by B. D. Ripley and Martin Maechler (spar/lambda, etc).

來源

該函數基於 T. Hastie 和 R. Tibshirani 編寫的 GAMFIT Fortran 程序中的代碼(最初取自 http://lib.stat.cmu.edu/general/gamfit )，該程序使用 Finbarr O'Sullivan 的樣條代碼。其設計與 Chambers & Hastie (1992) 的 smooth.spline 函數類似。

參考

Chambers, J. M. and Hastie, T. J. (1992) Statistical Models in S, Wadsworth & Brooks/Cole.

Green, P. J. and Silverman, B. W. (1994) Nonparametric Regression and Generalized Linear Models: A Roughness Penalty Approach. Chapman and Hall.

Hastie, T. J. and Tibshirani, R. J. (1990) Generalized Additive Models. Chapman and Hall.

也可以看看

predict.smooth.spline 用於評估樣條曲線及其導數。