R:smooth.spline LOOCV-错误取决于数据点的顺序?
R: smooth.spline LOOCV-error depends on order of datapoints?
我想对一些数据进行平滑样条拟合,我注意到内部计算的 LOOCV 误差似乎取决于数据是否无序。具体来说,我只是在数据排序的时候才得到预期的结果。
我不明白为什么会出现这种情况?有帮助吗?
set.seed(0)
x <- seq(1:10)
y <- x^2 + rnorm(10,0,2)
fit.ss <- smooth.spline(x=x, y=y, cv=TRUE)
cat("CV ordered: ",format(fit.ss$cv.crit))
# CV ordered: 13.46173
xu <- sample(x)
yu <- y[xu]
fit.ss.u <- smooth.spline(x=xu, y=yu, cv=TRUE)
cat("CV unorderd: ",format(fit.ss.u$cv.crit))
# CV unorderd: 65552.74
spar.opt <- fit.ss$spar
preds <- rep(NA, 10)
for (i in 1:10){
ss <- smooth.spline(x=x[-i], y=y[-i], cv=TRUE, spar=spar.opt)
preds[i] <- predict(ss,x=x[i])$y
}
cat("CV manual: ",format(mean((preds - y)**2)))
# CV manual: 13.49424
CV 订购和 CV 手册(几乎)相同并且符合预期,而未订购的版本完全不同。
请注意,这是 https://stats.stackexchange.com/q/561802/213798 的副本,我似乎没有得到任何输入。
看起来像是 smooth.spline
中的错误。当它在内部计算 cv.crit
时,它会将原始顺序的观察结果与 x
顺序的预测进行比较。 (我不确定确切的区别是什么,但大概是某种“留一”计算。)
代码如下:
cv.crit <-
if(is.na(cv)) NA
else {
r <- y - fit$ty[ox]
if(cv) {
ww <- wbar
ww[ww == 0] <- 1
r <- r / (1 - (lev[ox] * w)/ww[ox])
if(no.wgts) mean(r^2) else weighted.mean(r^2, w)
} else
(if(no.wgts) mean(r^2) else weighted.mean(r^2, w)) /
(1 - (df.offset + penalty * df)/n)^2
}
第 4 行看起来不对劲。在这一点上,你的未排序数据,我看到
Browse[2]> y
[1] 47.142866 80.988466 104.809307 25.829283 63.410559 3.525909 32.920100 3.347533 18.544859 11.659599
和
Browse[2]> fit$ty[ox]
[1] 2.458502 5.274807 11.019719 17.995820 25.281214 34.165585 46.918576 63.054358 82.093996 103.915902
所以看起来 fit$ty[ox]
是基于有序的 x
值,而 y
是原始顺序。
不幸的是,更正并不明显:此时 ox
是 TRUE
,所以它没有做任何事情。他们真正需要做的是以与 fit$ty
排序相同的方式对 y
进行排序。但是其他地方可能还有其他问题,因为当我尝试这样做时,这还不足以解决问题。
我想对一些数据进行平滑样条拟合,我注意到内部计算的 LOOCV 误差似乎取决于数据是否无序。具体来说,我只是在数据排序的时候才得到预期的结果。
我不明白为什么会出现这种情况?有帮助吗?
set.seed(0)
x <- seq(1:10)
y <- x^2 + rnorm(10,0,2)
fit.ss <- smooth.spline(x=x, y=y, cv=TRUE)
cat("CV ordered: ",format(fit.ss$cv.crit))
# CV ordered: 13.46173
xu <- sample(x)
yu <- y[xu]
fit.ss.u <- smooth.spline(x=xu, y=yu, cv=TRUE)
cat("CV unorderd: ",format(fit.ss.u$cv.crit))
# CV unorderd: 65552.74
spar.opt <- fit.ss$spar
preds <- rep(NA, 10)
for (i in 1:10){
ss <- smooth.spline(x=x[-i], y=y[-i], cv=TRUE, spar=spar.opt)
preds[i] <- predict(ss,x=x[i])$y
}
cat("CV manual: ",format(mean((preds - y)**2)))
# CV manual: 13.49424
CV 订购和 CV 手册(几乎)相同并且符合预期,而未订购的版本完全不同。
请注意,这是 https://stats.stackexchange.com/q/561802/213798 的副本,我似乎没有得到任何输入。
看起来像是 smooth.spline
中的错误。当它在内部计算 cv.crit
时,它会将原始顺序的观察结果与 x
顺序的预测进行比较。 (我不确定确切的区别是什么,但大概是某种“留一”计算。)
代码如下:
cv.crit <-
if(is.na(cv)) NA
else {
r <- y - fit$ty[ox]
if(cv) {
ww <- wbar
ww[ww == 0] <- 1
r <- r / (1 - (lev[ox] * w)/ww[ox])
if(no.wgts) mean(r^2) else weighted.mean(r^2, w)
} else
(if(no.wgts) mean(r^2) else weighted.mean(r^2, w)) /
(1 - (df.offset + penalty * df)/n)^2
}
第 4 行看起来不对劲。在这一点上,你的未排序数据,我看到
Browse[2]> y
[1] 47.142866 80.988466 104.809307 25.829283 63.410559 3.525909 32.920100 3.347533 18.544859 11.659599
和
Browse[2]> fit$ty[ox]
[1] 2.458502 5.274807 11.019719 17.995820 25.281214 34.165585 46.918576 63.054358 82.093996 103.915902
所以看起来 fit$ty[ox]
是基于有序的 x
值,而 y
是原始顺序。
不幸的是,更正并不明显:此时 ox
是 TRUE
,所以它没有做任何事情。他们真正需要做的是以与 fit$ty
排序相同的方式对 y
进行排序。但是其他地方可能还有其他问题,因为当我尝试这样做时,这还不足以解决问题。