如何从Tidymodels混合效应模型中提取random intercepts（随机截距）

我尝试使用lme4和tidymodels生态下的multilevelmod提取随机截距（random intercepts），用lme4可以通过以下方式实现：

使用R和lme4实现： #

library("tidyverse")
library("lme4")

# 构建模型
mod <- lmer(Reaction ~ Days + (1|Subject),data=sleepstudy)

# 生成扩展数据集
expanded_df <- with(sleepstudy,
                    data.frame(
                      expand.grid(Subject=levels(Subject),
                                  Days=seq(min(Days),max(Days),length=51))))

# 生成包含随机截距的预测结果
predicted_df <- data.frame(expanded_df,resp=predict(mod,newdata=expanded_df))

predicted_df 

# 绘制截距对比
ggplot(predicted_df,aes(x=Days,y=resp,colour=Subject))+
       geom_line() 

使用tidymodels实现： #

# 示例来自multilevelmod官方仓库
library("multilevelmod")
library("tidymodels")
library("tidyverse")
library("lme4")
#> Loading required package: parsnip
data(sleepstudy, package = "lme4")


# 设置模型引擎为lme4
mixed_model_spec <- linear_reg() %>% set_engine("lmer")


# 拟合模型
mixed_model_fit_tidy <- 
  mixed_model_spec %>% 
  fit(Reaction ~ Days + (1 | Subject), data = sleepstudy)


expanded_df_tidy <- with(sleepstudy,
                    data.frame(
                      expand.grid(Subject=levels(Subject),
                                  Days=seq(min(Days),max(Days),length=51))))



predicted_df_tidy <- data.frame(expanded_df_tidy,resp=predict(mixed_model_fit_tidy,new_data=expanded_df_tidy))


ggplot(predicted_df_tidy,aes(x=Days,y=.pred,colour=Subject))+
       geom_line()

我发现用predict()函数似乎仅返回固定效应的预测结果。有没有办法从tidymodels和multilevelmod中提取random intercepts？ 我知道这个包还在开发阶段，所以不确定现在有没有支持这个功能。

predict默认仅返回固定效应结果，是因为multilevelmod调用底层lme4的预测方法时，默认设置了仅计算固定效应，可通过以下两种方式获取随机截距相关结果：

方法1：预测时包含随机效应 #

在调用predict时传入re.form = NULL参数，该参数会直接传递给lme4的预测方法，代表纳入所有随机效应计算预测值，得到的结果就和原生lme4的输出完全一致：

# 沿用你之前生成的拟合对象和扩展数据集
predicted_df_tidy_random <- data.frame(
  expanded_df_tidy,
  predict(mixed_model_fit_tidy, new_data = expanded_df_tidy, re.form = NULL)
)

# 此时绘图就会显示每个Subject的独立截距
ggplot(predicted_df_tidy_random, aes(x = Days, y = .pred, colour = Subject)) +
  geom_line()

方法2：直接提取随机截距数值 #

如果不需要预测值，只想单独获取每个分组的随机截距估计值，可以直接调用拟合对象内部原生lme4模型的ranef()方法：

# 提取所有随机效应（当前场景下就是每个Subject的随机截距）
random_intercepts <- ranef(mixed_model_fit_tidy$fit)$Subject

# 若需要得到每个Subject的实际总截距，可以和固定效应截距相加
fixed_intercept <- fixef(mixed_model_fit_tidy$fit)[["(Intercept)"]]
random_intercepts$total_intercept <- fixed_intercept + random_intercepts[["(Intercept)"]]

# 查看结果
random_intercepts

内容的提问来源于stack exchange，提问作者daszlosek