R choose.k 平滑的基本尺寸选择

说明

使用惩罚回归平滑器时，选择基本维度并检查选择。

惩罚回归平滑器通过使用相对适中大小的基础 k 进行定义，从而获得计算效率。在 mgcv 包中设置模型时，在模型公式中使用 s 或 te 项，必须选择 k：默认值本质上是任意的。

实际上， k-1 (或 k )设置与 s 平滑相关的自由度上限(平滑上的可识别性约束通常会丢失 1 个自由度)。对于 te 平滑，自由度的上限由为约束的每个边平滑提供的 k 值的乘积给出。然而，实际有效自由度由拟合期间选择的惩罚程度(GCV、AIC、REML 或任何指定的值)控制。例外情况是如果使用 fx=TRUE 选项指定平滑，在这种情况下它不会受到惩罚。

因此，k 的精确选择通常并不重要：它应该选择足够大，以便您有足够的自由度来相当好地表示底层 ‘truth’，但又足够小以保持合理的计算效率。显然‘large’和‘small’取决于正在解决的特定问题。

与所有模型假设一样，能够非正式地检查 k 的选择非常有用。如果模型项的有效自由度估计远小于k-1，那么这不太值得，但随着 EDF 方法k-1，检查可能很重要。有用的通用方法如下：(i)拟合模型并提取偏差残差； (ii) 对于模型中的每个平滑项，使用大幅增加的 k 来拟合等效的、单一的、平滑的残差，以查看残差中是否存在可以通过增加 k 来解释的模式。下面提供了示例。

显而易见但成本更高的替代方案是简单地增加可疑的 k 并重新安装原始模型。如果这样做没有导致统计上重要的变化，那么 k 就足够大了。 (当 k 增加时，平滑度选择标准和/或有效自由度的变化提供了明显的数值测量，表明拟合是否发生了实质性变化。)

gam.check 对基本尺寸运行基于模拟的简单检查，这有助于标记 k 过低的术语。从 plot.gam 可用的部分残差中也可以看到严重太小的 k 。

可能导致混淆的一种情况是：模型使用 k=10 拟合以获得平滑项，并且该项的 EDF 估计为 7.6，低于最大值 9。然后使用 k=20 重新拟合该模型EDF 增加到 8.7 - 发生了什么 - 为什么 EDF 第一次不是 8.7？解释是，k=20 的函数空间包含 EDF 8.7 的函数子空间比 k=10 的函数空间更大：这个较大子空间中的函数之一比 EDF 8.7 中的任何函数都更适合数据。更小的子空间。然而，这些微妙之处很少对从模型拟合中得出的统计结论产生太大影响。

例子

## Simulate some data ....
library(mgcv)
set.seed(1) 
dat <- gamSim(1,n=400,scale=2)

## fit a GAM with quite low `k'
b<-gam(y~s(x0,k=6)+s(x1,k=6)+s(x2,k=6)+s(x3,k=6),data=dat)
plot(b,pages=1,residuals=TRUE) ## hint of a problem in s(x2)

## the following suggests a problem with s(x2)
gam.check(b)

## Another approach (see below for more obvious method)....
## check for residual pattern, removeable by increasing `k'
## typically `k', below, chould be substantially larger than 
## the original, `k' but certainly less than n/2.
## Note use of cheap "cs" shrinkage smoothers, and gamma=1.4
## to reduce chance of overfitting...
rsd <- residuals(b)
gam(rsd~s(x0,k=40,bs="cs"),gamma=1.4,data=dat) ## fine
gam(rsd~s(x1,k=40,bs="cs"),gamma=1.4,data=dat) ## fine
gam(rsd~s(x2,k=40,bs="cs"),gamma=1.4,data=dat) ## `k' too low
gam(rsd~s(x3,k=40,bs="cs"),gamma=1.4,data=dat) ## fine

## refit...
b <- gam(y~s(x0,k=6)+s(x1,k=6)+s(x2,k=20)+s(x3,k=6),data=dat)
gam.check(b) ## better

## similar example with multi-dimensional smooth
b1 <- gam(y~s(x0)+s(x1,x2,k=15)+s(x3),data=dat)
rsd <- residuals(b1)
gam(rsd~s(x0,k=40,bs="cs"),gamma=1.4,data=dat) ## fine
gam(rsd~s(x1,x2,k=100,bs="ts"),gamma=1.4,data=dat) ## `k' too low
gam(rsd~s(x3,k=40,bs="cs"),gamma=1.4,data=dat) ## fine

gam.check(b1) ## shows same problem

## and a `te' example
b2 <- gam(y~s(x0)+te(x1,x2,k=4)+s(x3),data=dat)
rsd <- residuals(b2)
gam(rsd~s(x0,k=40,bs="cs"),gamma=1.4,data=dat) ## fine
gam(rsd~te(x1,x2,k=10,bs="cs"),gamma=1.4,data=dat) ## `k' too low
gam(rsd~s(x3,k=40,bs="cs"),gamma=1.4,data=dat) ## fine

gam.check(b2) ## shows same problem

## same approach works with other families in the original model
dat <- gamSim(1,n=400,scale=.25,dist="poisson")
bp<-gam(y~s(x0,k=5)+s(x1,k=5)+s(x2,k=5)+s(x3,k=5),
        family=poisson,data=dat,method="ML")

gam.check(bp)

rsd <- residuals(bp)
gam(rsd~s(x0,k=40,bs="cs"),gamma=1.4,data=dat) ## fine
gam(rsd~s(x1,k=40,bs="cs"),gamma=1.4,data=dat) ## fine
gam(rsd~s(x2,k=40,bs="cs"),gamma=1.4,data=dat) ## `k' too low
gam(rsd~s(x3,k=40,bs="cs"),gamma=1.4,data=dat) ## fine

rm(dat) 

## More obvious, but more expensive tactic... Just increase 
## suspicious k until fit is stable.

set.seed(0) 
dat <- gamSim(1,n=400,scale=2)
## fit a GAM with quite low `k'
b <- gam(y~s(x0,k=6)+s(x1,k=6)+s(x2,k=6)+s(x3,k=6),
         data=dat,method="REML")
b 
## edf for 3rd smooth is highest as proportion of k -- increase k
b <- gam(y~s(x0,k=6)+s(x1,k=6)+s(x2,k=12)+s(x3,k=6),
         data=dat,method="REML")
b 
## edf substantially up, -ve REML substantially down
b <- gam(y~s(x0,k=6)+s(x1,k=6)+s(x2,k=24)+s(x3,k=6),
         data=dat,method="REML")
b 
## slight edf increase and -ve REML change
b <- gam(y~s(x0,k=6)+s(x1,k=6)+s(x2,k=40)+s(x3,k=6),
         data=dat,method="REML")
b 
## defintely stabilized (but really k around 20 would have been fine)

作者

Simon N. Wood simon.wood@r-project.org

参考

Wood, S.N. (2017) Generalized Additive Models: An Introduction with R (2nd edition). CRC/Taylor & Francis.

https://www.maths.ed.ac.uk/~swood34/