3

我使用包中的lmer()函数构建了一个混合效果模型lme4。出于某种良好的哲学原因，该lme4软件包不输出系数的 p 值。但是，我仍然需要 p 值在我的出版物中报告。我知道有多种方法可以使用由创建的模型计算 p 值lmer()，例如这里。

我的问题是：我想使用包中的tidy()函数提取 p 值broom。在这里，我真的很想坚持下去，tidy()因为我想维护以下管道： 

data_frame %>% group_by(grouping variables) %>% do(tidy(fitted_model))

一种选择是创建一个自定义函数并将其附加到管道中。但是，broom包的手册页（http://rpackages.ianhowson.com/cran/broom/man/lme4_tidiers.html）说：

"p.value  P-value computed from t-statistic (may be missing/NA)".

通过这个，我假设一个从 lmer 模型给出的 t 值计算 p 值的函数已经在 broom 中实现。所以，我不愿意重新发明轮子。

问题是我根本没有得到名为 p.value 的列。我期待一个名为 p.value 的列，其中 NA 是最坏的情况。

代码：

library(lme4)
library(broom)
lme <- lmer(Reaction ~ Days + (1 + Days | Subject), sleepstudy)
tidy(lme)
tidy(lme, effects = "fixed")

输出：

> tidy(lme)
                         term     estimate std.error statistic    group
1                  (Intercept) 251.40510485  6.824557 36.838306    fixed
2                         Days  10.46728596  1.545789  6.771485    fixed
3       sd_(Intercept).Subject  24.74045195        NA        NA  Subject
4              sd_Days.Subject   5.92213312        NA        NA  Subject
5 cor_(Intercept).Days.Subject   0.06555113        NA        NA  Subject
6      sd_Observation.Residual  25.59181564        NA        NA Residual
> tidy(lme, effects = "fixed")
         term  estimate std.error statistic
1 (Intercept) 251.40510  6.824557 36.838306
2        Days  10.46729  1.545789  6.771485
4

2 回答 2

5

您将需要包lmerTest来获取 p 值。tidy不适用于该lme对象，您需要将其附加到您的格式中。

attach(mtcars)
lme <- lmer(mpg ~ cyl + (1 + cyl | carb), mtcars)
summary(lme)
于 2016-05-13T03:43:31.050 回答
1

您可以使用sjstats::p_value()从许多不同模型（包括（广义）混合模型）中获取 p 值。

于 2018-08-10T13:38:47.750 回答