不一致的 mgcv gam.check 结果

机器算法验证 回归 广义加法模型 毫克CV
2022-04-15 10:14:27

我不确定这个问题是否适合交叉验证，但我不确定在哪里发布它。

我已经使用mgcv包构建了一个简单的模型。

a <- gam(x ~ s(y), method="REML", data=dat100k)

但是，当我运行时gam.check(a)，每次gam.check(a)使用相同的模型运行时都会得到不一致的输出。在第一次调用中 k 被认为是好的gam.check(a)（尽管 k 接近 edf 并且 p 不是很大），并且在第二次调用中根据 p 值认为 k 太低了。这是正常的吗？

这是运行 gam.check(a) 两次的输出：

gam.check(a)

Method: REML   Optimizer: outer newton
full convergence after 8 iterations.
Gradient range [-0.04365533,0.04277738]
(score 407913 & scale 214.1436).
Hessian positive definite, eigenvalue range [4.001815,49713.04].
Model rank =  10 / 10 

Basis dimension (k) checking results. Low p-value (k-index<1) may
indicate that k is too low, especially if edf is close to k'.

                    k'  edf k-index p-value
s(y) 9.00 8.98       1    0.66

> gam.check(a)

Method: REML   Optimizer: outer newton
full convergence after 8 iterations.
Gradient range [-0.04365533,0.04277738]
(score 407913 & scale 214.1436).
Hessian positive definite, eigenvalue range [4.001815,49713.04].
Model rank =  10 / 10 

Basis dimension (k) checking results. Low p-value (k-index<1) may
indicate that k is too low, especially if edf is close to k'.

                    k'  edf k-index p-value  
s(y) 9.00 8.98    0.98    0.04 *
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
1个回答

问题是由于gam.check()基于模型残差排列的基础维度测试。这些排列是使用伪随机数生成器计算的；通过设计，每次调用gam.check()（或直接调用k.check()自身）时，都会产生一组不同的排列，这会巧妙地改变所执行的排列测试的 p 值。

如果在调用gam.check()or之前设置了伪随机数生成器的种子k.check()，结果会变得一致：

library(mgcv)
set.seed(0)
dat <- gamSim(1,n=200)
b <- gam(y ~ s(x0) + s(x1) + s(x2) + s(x3), data = dat)

set.seed(1)
k.check(b)
set.seed(1)
k.check(b)

产生：

> set.seed(1)
> k.check(b)
      k'      edf   k-index p-value
s(x0)  9 2.318172 0.9959628  0.4950
s(x1)  9 2.305695 0.9693887  0.3225
s(x2)  9 7.654740 0.9605490  0.2825
s(x3)  9 1.232618 1.0372831  0.6850
> set.seed(1)
> k.check(b)
      k'      edf   k-index p-value
s(x0)  9 2.318172 0.9959628  0.4950
s(x1)  9 2.305695 0.9693887  0.3225
s(x2)  9 7.654740 0.9605490  0.2825
s(x3)  9 1.232618 1.0372831  0.6850

如果您执行的排列数量比默认值多n.rep = 400，那么所有其他等于 p 值的东西都应该开始稳定，但这会带来更大的计算成本。

我也会将此测试的输出（当然是 p 值）作为指导。在您的示例中，平滑的 EDF 几乎处于其最大可能值。无论简单的启发式测试的 p 值是多少，我都想增加（在公式中的 s(...)` 调用中将其k加倍）并查看估计的平滑度是否发生了变化；k = 20它是否使用比原始尺寸更大的 EDF？

其它你可能感兴趣的问题
上一篇AIC 或 BIC 值差异的实际意义是什么？ 下一篇神经网络提供了哪些传统非线性统计模型不提供的功能？