mgcv::gam 无法正确分解平滑的线性分量

问题描述 投票:0回答:1

我有兴趣在 GAMM 中单独测试线性趋势和非线性分量。

我遵循了 mgcv 文档中的示例以及 stackexchange 上相关问题的答案,以使用

m=c(2,0)
和默认的
tp
平滑。该示例符合我的预期:(1) x 的参数系数估计的斜率非常接近简单线性回归的估计; (2) 平滑看起来去趋势且相关(非常小)大小。

但是,当我尝试用我的数据执行此操作时,结果对我来说毫无意义。看起来线性分量是任意的(有时甚至没有正确的符号,如下例所示)并且平滑没有去趋势(它显然包含强线性分量)。

我想知道我是否误解了这应该做什么,或者我是否以某种方式错误地编码了。由于文档示例工作正常,我不认为这是 mgcv 中的错误。

我还能如何在 GAMM 中独立于非线性分量来测试线性趋势?

这是 mgcv 文档中的相关示例

require(mgcv)
n <- 100
set.seed(2)
x <- runif(n)
y <- x + x^2*.2 + rnorm(n) *.1
ed <- data.frame(y,x)
mg <- gam(y~s(x,m=c(2,0))+x,data=ed,method="REML")

据我了解,这从平滑中删除了线性分量,以便可以单独估计。

这是 gam 和 lm 的斜率估计值,非常接近:

> summary(mg)

Family: gaussian 
Link function: identity 

Formula:
y ~ s(x, m = c(2, 0)) + x

Parametric coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept) -0.02133    0.05315  -0.401    0.689    
x            1.18249    0.10564  11.193   <2e-16 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Approximate significance of smooth terms:
        edf Ref.df     F p-value  
s(x) 0.9334      8 0.304   0.076 .
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

R-sq.(adj) =   0.91   Deviance explained = 91.1%
-REML = -70.567  Scale est. = 0.012767  n = 100

> summary(lm(y~x,ed))


Call:
lm(formula = y ~ x, data = ed)

Residuals:
     Min       1Q   Median       3Q      Max 
-0.24082 -0.07942 -0.01136  0.09430  0.22836 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept) -0.03019    0.02212  -1.365    0.175    
x            1.20052    0.03851  31.175   <2e-16 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 0.1144 on 98 degrees of freedom
Multiple R-squared:  0.9084,    Adjusted R-squared:  0.9075 
F-statistic: 971.9 on 1 and 98 DF,  p-value: < 2.2e-16

这是带有 lm(红色)和 gam(蓝色)的叠加预测的数据图:

plot(y~x,ed)
abline(lm(y~x,ed),col="red")
points(ed$x,predict(mg,ed),col="blue")

绘制 lm 和 gam 的数据和模型预测

这是平滑的,它看起来完全去趋势化并且非常小,与图表一致:

plot.gam(mg)

x 平滑图

这是使用我的数据进行的相同分析。

read.table("ev.txt") -> ev
me <- gam(logart~s(item,m=c(2,0))+item,data=ev,method="REML")

由 gam 估计的项目的线性斜率是负值 (-0.11),尽管实际斜率是正值 (+0.32):

> summary(me)

Family: gaussian 
Link function: identity 

Formula:
logart ~ s(item, m = c(2, 0)) + item

Parametric coefficients:
             Estimate Std. Error  t value Pr(>|t|)    
(Intercept)  6.305740   0.003587 1758.100   <2e-16 ***
item        -0.112921   0.049174   -2.296   0.0217 *  
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Approximate significance of smooth terms:
          edf Ref.df     F  p-value    
s(item) 6.125      8 3.583 3.15e-05 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

R-sq.(adj) =  0.0298   Deviance explained = 3.14%
-REML = -842.54  Scale est. = 0.039768  n = 4453


> summary(lm(logart~item,ev))

Call:
lm(formula = logart ~ item, data = ev)

Residuals:
     Min       1Q   Median       3Q      Max 
-0.65009 -0.12830  0.00627  0.13301  0.61885 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept) 6.311575   0.003001 2103.47   <2e-16 ***
item        0.031751   0.003048   10.42   <2e-16 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 0.2001 on 4451 degrees of freedom
Multiple R-squared:  0.0238,    Adjusted R-squared:  0.02359 
F-statistic: 108.5 on 1 and 4451 DF,  p-value: < 2.2e-16

这是数据图以及 lm 和 gam 的预测,两者看起来完全合理:

plot(logart~item,ev,col="#20202020")
abline(lm(logart~item,ev),col="red")
points(ev$item,predict(me,ev),col="blue")

我的数据与 lm 和 gam 预测的关系图

因此,整体模型看起来不错,只是平滑实际上并未去趋势,并且似乎完成了大部分工作(注意 y 轴范围):

plot.gam(me)

平滑的情节

(结果与

bam
相同,并且明确指定
bs='tp'
没有区别,正如预期的那样。)

这可能吗?就与 lm/lmer 估计相匹配的 gam/bam (mis) 线性分量而言,不同的数据集如何导致如此不同的结果?我是否误解了单独测试非线性的要点(或过程)? 我还应该/还能怎么做?

mgcv
1个回答
0
投票

以下答案由 Simon Wood 通过电子邮件发送:

平滑基础与线性趋势不正交,因此您 不会期望斜率估计值与 lm 的斜率估计值相匹配。重点 这种结构的真正目的是为了测试是否有任何东西 是否需要超越线性趋势,即是否需要平滑 而不是简单的线性效应。

我不能说我理解这个模型中的线性斜率到底代表什么,但听起来该模型不应该像我错误地假设的那样产生线性和去趋势非线性分量的分离。为此,伍德建议可以

克隆平滑构造函数并添加代码以正交化基础

© www.soinside.com 2019 - 2024. All rights reserved.