如何解决比较gammmgcv模型来获得交互效果是否合适以及如何去做?
我想知道随着时间的推移,音量的发展(我指的是年龄,而不是波浪)是否在组之间有所不同。我也有一些协变量。 A做了一个模拟数据集:
library(simstudy)
def <- defData(varname = "volume",dist = "normal",formula = 1000,variance = 100)
def <- defData(def,varname = "group",formula = "1;2;3",dist = "categorical")
def <- defData(def,varname = "age",formula = 14,variance = 2)
def <- defData(def,varname = "scanner",formula = "1;2",dist = "categorical")
set.seed(1)
sim.data <- genData(220,def)
我想出了以下模型:
library(mgcv)
sim.data$group <- as.factor(sim.data$group)
m1 <- gamm(volume ~ group + s(age,by = group,k = 20,bs = "tp") + scanner,data = sim.data,random=list(id=~1))
summary(m1$gam)
结果如下
Family: gaussian
Link function: identity
Formula:
volume ~ group + s(age,bs = "tp") +
scanner
Parametric coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 997.3096 2.7303 365.280 <2e-16 ***
group2 0.5983 1.9557 0.306 0.760
group3 1.2819 1.7592 0.729 0.467
scanner 1.3806 1.3627 1.013 0.312
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Approximate significance of smooth terms:
edf Ref.df F p-value
s(age):group1 1 1 1.850 0.175
s(age):group2 1 1 0.001 0.969
s(age):group3 1 1 0.194 0.660
R-sq.(adj) = -0.0133
Scale est. = 2.1375e-07 n = 220
我希望我可以用我的 s(age,by=group) 术语解释各组之间发展的差异,但是它给了我每个组的平滑年龄参数。我怎样才能最好地评估交互作用(虽然我也在某处读到交互可能不是这个加性模型中的合适术语)?
我想过使用 anova(),比较一个模型和一个没有 by term 的模型
m2 <- gamm(volume ~ group + s(age,random=list(id=~1))
anova(m1$gam,m2$gam)
但是根据我将 lm 模型与 anova() 进行比较的经验,这并没有给我预期的输出
Family: gaussian
Link function: identity
Formula:
volume ~ group + s(age,bs = "tp") +
scanner
Parametric Terms:
df F p-value
group 2 0.297 0.743
scanner 1 1.026 0.312
Approximate significance of smooth terms:
edf Ref.df F p-value
s(age):group1 1 1 1.850 0.175
s(age):group2 1 1 0.001 0.969
s(age):group3 1 1 0.194 0.660
有谁知道比较模型是否是研究交互效果的合适方法以及如何获得实际效果?最好是一种我也可以用于成对分析的方法。
非常感谢所有帮助!
解决方法
有人离线向我提供了解决方案。您需要对分组变量进行排序,然后添加如下对比度:
sim.data$group<-as.factor(sim.data$group)
sim.data$group<-as.ordered(sim.data$group)
contrasts(sim.data$group) <- 'contr.treatment'
您可以像往常一样进行分析
m1 <- gamm(volume ~ group + s(age,by = group,k = 20,bs = "tp") + scanner,data = sim.data,random=list(id=~1))
summary(m1$gam)
希望这能帮助任何与我有同样问题的人。
版权声明:本文内容由互联网用户自发贡献,该文观点与技术仅代表作者本人。本站仅提供信息存储空间服务,不拥有所有权,不承担相关法律责任。如发现本站有涉嫌侵权/违法违规的内容, 请发送邮件至 dio@foxmail.com 举报,一经查实,本站将立刻删除。