lmer模型 – 访问模型的元素

huangapple go评论61阅读模式
英文:

lmer model - accessing elements of the model

问题

我是R初学者,对lmer和混合效应建模也不太了解。

如果有一个lmer模型:

model <- lmer(C ~ Group * Time + (1|ID), REML = FALSE, data=data)

你可以使用以下方法访问模型的细节以进行绘图:

plot(model)

或者

plot(predict(model), resid(model))

还有其他访问模型元素/细节的方法吗?可以访问哪些模型元素?我知道可以使用summary()函数查看随机效应和固定效应,但如何绘制残差与模型的其他元素,例如残差与固定效应或随机效应?

英文:

I am an R beginner and even newer to lmer and mixed effects modelling.

If I have an lmer model of:

model &lt;- lmer(C ~ Group * Time + (1|ID), REML = FALSE, data=data)

How do I access details of the model in order to plot them?

I know I can plot the residuals against the fitted values using either:

plot(model)

OR

plot(predict(model), resid(model))

But is there a way to access other elements/details of the model? What other model elements can be accessed? I am aware that I can use the function 'summary()' to look up the random and fixed effect. But what about plotting the residuals against other elements of the model e.g. residuals against fixed effects or random effect?

答案1

得分: 3

当你有这样的疑问时,你应该查看包的文档,比如在这种情况下,可以使用 vignette(lme4)。仅通过简单的搜索,我找到了这个表格,也许对你有用:

方法 描述
anova 固定效应贡献的分解或模型比较。
as.function 返回轮廓化的偏差或REML准则的函数。
coef 每个水平的随机效应和固定效应之和。
confint 线性混合模型参数的置信区间。
deviance 减去两倍的最大对数似然。 (使用 REMLcrit 获取REML准则。)
df.residual 残差自由度。
drop1 从模型中去掉可允许的单一项。
extractAIC 广义Akaike信息准则。
fitted 在给定条件模式下的拟合值。
fixef 固定效应系数的估计。
formula 拟合模型的混合模型公式。
logLik 最大对数似然。
model.frame 适合模型所需的数据。
model.matrix 固定效应模型矩阵。
ngrps 每个分组因子中的水平数。
nobs 观测数。
plot 混合模型拟合的诊断图。
predict 各种类型的预测值。
print 混合模型对象的基本打印输出。
profile 不同模型参数的轮廓化似然。
ranef 随机效应的条件模式。
refit 重新适应到响应变量的新观测的模型。
refitML 最大似然重新适应的模型。
residuals 各种类型的残差值。
sigma 残差标准偏差。
simulate 从适合的混合模型生成的模拟数据。
summary 混合模型的摘要。
terms 混合模型的术语表示。
update 使用修订后的公式或其他参数的更新模型。
VarCorr 估计的随机效应方差、标准差和相关性。
vcov 固定效应估计值的协方差矩阵。
weights 用于模型拟合的先验权重。
英文:

When you have doubts like this you shoud check the Vignette of the package, in this case vignette(lme4). Only with a superficial search I found this table that maybe is of use:

Method Description
anova Decomposition of fixed-effects contributions or model comparison.
as.function Function returning profiled deviance or REML criterion.
coef Sum of the random and fixed effects for each level.
confint Confidence intervals on linear mixed-model parameters.
deviance Minus twice maximum log-likelihood. (Use REMLcrit for the REML criterion.)
df.residual Residual degrees of freedom.
drop1 Drop allowable single terms from the model.
extractAIC Generalized Akaike information criterion
fitted Fitted values given conditional modes.
fixef Estimates of the fixed-effects coefficients
formula Mixed-model formula of fitted model.
logLik Maximum log-likelihood.
model.frame Data required to fit the model.
model.matrix Fixed-effects model matrix
ngrps Number of levels in each grouping factor.
nobs Number of observations.
plot Diagnostic plots for mixed-model fits.
predict Various types of predicted values.
print Basic printout of mixed-model objects.
profile Profiled likelihood over various model parameters.
ranef Conditional modes of the random effects.
refit A model (re)fitted to a new set of observations of the response variable.
refitML A model (re)fitted by maximum likelihood.
residuals Various types of residual values.
sigma Residual standard deviation.
simulate Simulated data from a fitted mixed model.
summary Summary of a mixed model.
terms Terms representation of a mixed model.
update An updated model using a revised formula or other arguments.
VarCorr Estimated random-effects variances, standard deviations, and correlations.
vcov Covariance matrix of the fixed-effect estimates.
weights Prior weights used in model fitting.

答案2

得分: 2

这是一个绘制残差与时间因子的方法之一:

library(broom.mixed)
aa <- augment(model)
boxplot(.resid ~ Time, data = aa)

你还可以使用内置的绘图方法(参见help("plot.merMod")):

plot(model, resid(.) ~ Time)
英文:

Here's one way to plot residuals vs your time factor:

library(broom.mixed)
aa &lt;- augment(model)
boxplot(.resid ~ Time, data = aa)

You could also use the built-in plot method (see help(&quot;plot.merMod&quot;)):

plot(model, resid(.) ~ Time)

答案3

得分: 0

是的,您可以使用summary()函数,就像对待任何其他lmglm对象一样,以获得对实际模型的评估。

同样,您可以将这些对象视为数据框的一种。

唯一的区别是,您应该使用@符号,而不是使用$

例如:

lm <- lmer(Petal.Length ~ Sepal.Length * Sepal.Width + (1|Species), data = iris)

lm@frame %>% head()

 Petal.Length Sepal.Length Sepal.Width Species
1          1.4          5.1         3.5  setosa
2          1.4          4.9         3.0  setosa
3          1.3          4.7         3.2  setosa
4          1.5          4.6         3.1  setosa
5          1.4          5.0         3.6  setosa
6          1.7          5.4         3.9  setosa

从那一点开始,您可以深入探索对象。

英文:

Yes you can use summary() function just like with any other lm or glm object to get an assessment of the actual model.

Similarly you can treat these objects like you would a dataframe.

The only difference is instead of using the $ you should use the @ symbol.

For example:

lm &lt;- lmer(Petal.Length ~ Sepal.Length * Sepal.Width + (1|Species), data = iris)

lm@frame %&gt;% head()

 Petal.Length Sepal.Length Sepal.Width Species
1          1.4          5.1         3.5  setosa
2          1.4          4.9         3.0  setosa
3          1.3          4.7         3.2  setosa
4          1.5          4.6         3.1  setosa
5          1.4          5.0         3.6  setosa
6          1.7          5.4         3.9  setosa

You can from that point explore deeper in the object.

huangapple
  • 本文由 发表于 2023年4月17日 21:29:40
  • 转载请务必保留本文链接:https://go.coder-hub.com/76035709.html
匿名

发表评论

匿名网友

:?: :razz: :sad: :evil: :!: :smile: :oops: :grin: :eek: :shock: :???: :cool: :lol: :mad: :twisted: :roll: :wink: :idea: :arrow: :neutral: :cry: :mrgreen:

确定