🧐 lme4 | 多层线性模型小彩蛋解答
1写在前面
本期回答一下上期中小彩蛋部分的问题，如下：
Q: 不同的
department的base不同，raise也不同，我们得出不同的α和β。 可否等价为，先按照department分组，然后分别计算α和β。 A: 不等价！
2用到的包
rm(list = ls())
library(tidyverse)
library(lme4)
library(modelr)
library(broom)
library(ggsci)
library(broom.mixed)
3示例数据
数据描述的是不同部门(department)的老师的收入(salary)情况。
具体可见上期。
4多层线性模型
4.1 建模
m1 <- lmer(salary ~ experience + (1 + experience | department), data = df)
m1
broom.mixed::tidy(m1, effects = "ran_vals")
4.2 加上predictions
df1 <- df %>%
add_predictions(m1)
df1
4.3 可视化
p1 <- df1 %>%
ggplot(aes(
x = experience, y = salary, group = department,
colour = department
)) +
geom_point() +
geom_line(aes(x = experience, y = pred)) +
labs(x = "Experience", y = "Predicted Salary") +
ggtitle("Varying Intercept and Slopes Salary Prediction") +
scale_color_npg()
p1
5分组建立简单线性回归模型
5.1 批量建模
这里我们使用nest函数容纳一下建模的大量数据。这里我就不做具体展示了，大家可以自己试一下。
m2 <- df %>%
group_by(department) %>%
nest() %>%
mutate(mdl = map(data, ~ lm(salary ~ 1 + experience, data=.))) %>%
mutate(fit = map(mdl, ~ .$fitted.values))
m2
5.2 加上predictions
df2 <- m2 %>%
mutate(., data = map2(data, mdl, add_predictions)) %>%
select(., -mdl, -fit) %>%
unnest()
df2
5.3 可视化
p2 <- df2 %>%
ggplot(aes(
x = experience, y = salary, group = department,
colour = department
)) +
geom_point() +
geom_line(aes(x = experience, y = pred)) +
labs(x = "Experience", y = "Predicted Salary") +
ggtitle("Varying Intercept and Slopes Salary Prediction") +
scale_color_npg()
p2
6比较两种模型
6.1 比较pred值
差异不是很大，但还是不同的。🤣
df3 <- df1 %>%
dplyr::select(.,ids,pred) %>%
left_join(.,df2[,c(2,7)],by = "ids")
df3
6.2 可视化
library(patchwork)
p1 + p2
Note! 大家不要认为差异不大就随便选用建模方法，当我们纳入更多变量的时候，可能不同建模方式的差异就会显现啦！🥰
最后祝大家早日不卷!~
本文参与 腾讯云自媒体同步曝光计划，分享自微信公众号。
原始发表：2022-09-10，如有侵权请联系 cloudcommunity@tencent.com 删除
评论
登录后参与评论
推荐阅读
目录
腾讯云开发者
