我想将响应变量 ( y
) 建模为两个解释变量 (x
和z
) 的函数。是采用世界公认的方法(方法)y
对人体生理参数的浓度测量。是一种新颖的方法(方法),可能是替代品,因为它便宜得多。然而,人们认为该方法或多或少是准确的,具体取决于使用所需设备的时间( )。所以,我的目标是测试和准确测量. 我有 6 个人 ( : , , , , , )。下面我展示和之间的关系:A
x
B
B
z
x
z
y
ID
A
B
C
D
E
F
x
y
*注意:这里我分类Z
是为了便于说明,但是z
, as x
, 是数字,范围从 1 (hr) 到 7 (hr)。
鉴于我的响应变量具有明显的非正态分布,方差随着增加而x
增加,并且我有几个人,但我对个体之间的差异不感兴趣,我虽然在 GLMM 上使用 GAMMA 分布来测试显着效果x
并z
解释y
。
我使用 gamma 分布并使用日志链接运行 GLMM。
## Setting random structure ####
mod1<-glmer(Y~ 1 + (1|ID),data = df, family=Gamma(link=log))
mod2<-glmer(Y~ 1 + (X|ID),data = df, family=Gamma(link=log))
AIC(mod1,mod2)
## Setting fixed structure ####
mod1<-glmer(Y~ 1 + (X|ID),data = df, family=Gamma(link=log), control = glmerControl(optimizer ="Nelder_Mead"))
mod2<-glmer(Y~X + (X|ID),data = df, family=Gamma(link=log), control = glmerControl(optimizer ="Nelder_Mead"))
mod3<-glmer(Y~X + Z + (X|ID),data = df, family=Gamma(link=log), control = glmerControl(optimizer ="Nelder_Mead"))
mod4<-glmer(Y~X + X:Z + (X|ID),data = df, family=Gamma(link=log), control = glmerControl(optimizer ="Nelder_Mead"))
r.squaredGLMM(mod4)[1,c(1,2)]
mod.list <- list(mod1,mod2,mod3,mod4)
model.sel<-model.sel(mod.list, rank="AIC",extra=c(r2=function(x) round(r.squaredGLMM(x)[1,c(1,2)],2)))
为了测试是否显着x
,z
我比较了 AIC 的模型。
model.sel
Model selection table
(Intr) ns(VDB.V13,3) n.V13 cnd((Int)) dsp((Int)) cnd(ns(VDB.V13,3)) cnd(n.V13:VDB.V13) r2.R2m r2.R2c class control ziformula dispformula random df logLik AIC delta weight
3 -2.567 + -0.04178 0.66 0.71 glmerMod gC(Nl_Md) VD.V1|I 9 2017.580 -4017.2 0.00 1
2 -2.645 + 0.65 0.70 glmerMod gC(Nl_Md) VD.V1|I 8 2006.875 -3997.7 19.41 0
4 -2.661 + + + 0.66 0.76 glmmTMB ~0 ~1 c(VD.V1|I) 9 2001.622 -3985.2 31.92 0
1 -1.559 0.00 0.36 glmerMod gC(Nl_Md) VD.V1|I 5 1682.428 -3354.9 662.31 0
Abbreviations:
control: gC(Nl_Md) = ‘glmerControl(Nelder_Mead)’
Models ranked by AIC(x)
Random terms:
VD.V1|I = ‘VeDBA.V13AP | ID’
c(VD.V1|I) = ‘cond(VeDBA.V13AP | ID)’
当我观察残差模式与预测值时,我的问题就出现了,因为我认为有明确的模式。在这里,我没有显示残差的分布,因为如果我理解正确,则不需要它们的正态分布。
有谁知道我能做什么?有什么建议吗?我没有分享数据,因为太长(n=2027)。
提前致谢
我的数据框负责人:
ID Y X Z
1 A 0.34136077 1.55682000 2
2 A 0.05124066 0.05766000 2
3 A 0.05901189 0.05125333 3
4 A 0.05213855 0.05766000 2
5 A 0.05437708 0.05125333 3
6 A 0.08433229 0.05766000 3
7 A 0.03618396 0.04484667 3
8 A 0.03622474 0.05766000 1
9 A 0.18244336 0.05125333 3
10 A 0.03625487 0.03844000 2
11 A 0.03840890 0.04484667 3
12 A 0.04235018 0.03844000 3
13 A 0.03862926 0.03844000 3
14 A 0.03749647 0.02883000 2
15 A 0.04395015 0.03844000 2
16 A 0.04040225 0.04805000 2
17 A 0.04419507 0.05766000 3
18 A 0.33186947 2.53704000 1
19 A 0.31986092 0.74958000 1
20 A 0.08127853 0.05766000 1
")
根据 JTH 的提议遵循的程序
mod1<-glmer(Y~ 1 + (X|ID),data = df, family=Gamma(link=log), control = glmerControl(optimizer ="Nelder_Mead"))
mod2<- glmer(Y~ ns(X, 4) + (X|ID),data = df, family=Gamma(link=log), control = glmerControl(optimizer ="Nelder_Mead"))
mod3<-glmer(Y~ ns(X, 4) + Z + (X|ID),data = df, family=Gamma(link=log), control = glmerControl(optimizer ="Nelder_Mead"))
mod4<-glmer(Y~ X:Z + ns(X, 4) + (X|ID),data = df, family=Gamma(link=log), control = glmerControl(optimizer ="Nelder_Mead"))
mod5<-glmer(Y~ ns(X, 4):Z + (X|ID),data = df, family=Gamma(link=log), control = glmerControl(optimizer ="Nelder_Mead"))
AIC(mod1,mod2,mod3,mod4,mod5)
结果如下:
r.squaredGLMM(mod4)[1,c(1,2)]
mod.list <- list(mod1,mod2,mod3,mod4)
model.sel<-model.sel(mod.list, rank="AIC",extra=c(r2=function(x) round(r.squaredGLMM(x)[1,c(1,2)],2)))
model.sel
Int ns(X) Z Z:X Z:ns(X) r2m r2c df logLik AIC delta weight
4 -2.827029 + 0.1871558 0.65 0.69 10 320.9342 -621.8684 0.00000 1.000000e+00
2 -2.660326 + 0.48 0.54 9 289.8442 -561.6884 60.17996 8.552394e-14
3 -2.660204 + 0.0001000814 0.48 0.54 10 289.8443 -559.6886 62.17976 3.146559e-14
5 -2.220731 + 0.04 0.72 9 259.3610 -500.7219 121.14643 4.936132e-27
1 -1.377842 0.00 0.27 5 246.9520 -483.9039 137.96446 1.100015e-30
我的残差与预测值的图现在是这样的:
尽管模式已急剧减少,但我怀疑仍然存在一些异方差性。但是,我无法删除它。