在混合效应模型中绘制二次函数
Drawing a quadratic function in mixed effect model
这个问题有点长,感谢您的耐心等待。
这是我的数据
https://www.dropbox.com/s/jo22d68a8vxwg63/data.csv?dl=0
我构建了一个混合效应模型
library(lme4)
mod <- lmer(sqrt(y) ~ x1 + I(x1^2) + x2 + I(x2^2) + x3 + I(x3^2) + x4 + I(x4^2) + x5 + I(x5^2) +
x6 + I(x6^2) + x7 + I(x7^2) + x8 + I(x8^2) + (1|loc) + (1|year), data = data)
所有预测变量都是标准化的,我很想知道 y
如何随着 x5
的变化而变化,同时将其他变量保持在它们的平均值(等于 0,因为所有变量都是标准化的).
我就是这样做的。
# make all predictors except x5 equal to zero
data$x1<-0
data$x2<-0
data$x3<-0
data$x4<-0
data$x6<-0
data$x7<-0
data$x8<-0
# Use the predict function
library(merTools)
fitted <- predictInterval(merMod = mod, newdata = data, level = 0.95, n.sims = 1000,stat = "median",include.resid.var = TRUE)
现在我想将拟合绘制为 x5
的二次函数。我这样做:
i<-order(data$x5)
plot(data$x5[i],fitted$fit[i],type="l")
我预计这会产生一个 y
的图作为 x5
的二次函数。但如您所见,我得到了以下没有任何二次曲线的图。谁能告诉我我在这里做错了什么?
我不确定 predictInterval
的来源,但您可以使用 predict
来做到这一点。诀窍就是确保将随机效果设置为 0。以下是如何做到这一点
newdata <- data
newdata[,paste0("x", setdiff(1:8,5))] <- 0
y <- predict(mod, newdata=newdata, re.form=NA)
plot(data$x5, y)
re.form=NA
部分去掉随机效应
这个问题有点长,感谢您的耐心等待。
这是我的数据
https://www.dropbox.com/s/jo22d68a8vxwg63/data.csv?dl=0
我构建了一个混合效应模型
library(lme4)
mod <- lmer(sqrt(y) ~ x1 + I(x1^2) + x2 + I(x2^2) + x3 + I(x3^2) + x4 + I(x4^2) + x5 + I(x5^2) +
x6 + I(x6^2) + x7 + I(x7^2) + x8 + I(x8^2) + (1|loc) + (1|year), data = data)
所有预测变量都是标准化的,我很想知道 y
如何随着 x5
的变化而变化,同时将其他变量保持在它们的平均值(等于 0,因为所有变量都是标准化的).
我就是这样做的。
# make all predictors except x5 equal to zero
data$x1<-0
data$x2<-0
data$x3<-0
data$x4<-0
data$x6<-0
data$x7<-0
data$x8<-0
# Use the predict function
library(merTools)
fitted <- predictInterval(merMod = mod, newdata = data, level = 0.95, n.sims = 1000,stat = "median",include.resid.var = TRUE)
现在我想将拟合绘制为 x5
的二次函数。我这样做:
i<-order(data$x5)
plot(data$x5[i],fitted$fit[i],type="l")
我预计这会产生一个 y
的图作为 x5
的二次函数。但如您所见,我得到了以下没有任何二次曲线的图。谁能告诉我我在这里做错了什么?
我不确定 predictInterval
的来源,但您可以使用 predict
来做到这一点。诀窍就是确保将随机效果设置为 0。以下是如何做到这一点
newdata <- data
newdata[,paste0("x", setdiff(1:8,5))] <- 0
y <- predict(mod, newdata=newdata, re.form=NA)
plot(data$x5, y)
re.form=NA
部分去掉随机效应