我有一个纵向研究的数据集,其中有10个结果变量和16个自变量。但为了简化我的问题,我编辑了数据集,现在它包含:
- 4列作为结果变量
- 1列作为自变量(TimePointNum2)
- 6例受试者(受试者)
- 每例受试者随访5次(时间点编号2)
数据集可以通过here下载。如果您不满意,我已经使用dput()
为dataframe m2创建了对象(参见下面的数据)
library(tidyverse)
library(broom)
library(lme4)
library(nlme)
m2 <- read_csv('question.csv')
我的目标是:
- 针对 Dataframe 的每个结果变量(MUFA、SFA、体重、体脂)运行(迭代或循环)回归模型,其中每个结果的自变量是TimePointNum2。
所以我运行了很多次lmer
(这是违反DRY规则的)
mod_SFA_num <- lmer(SFA ~ TimePointNum2 + (1| Subject), data = m2)
mod_MUFA_num <- lmer(MUFA ~ TimePointNum2 + (1| Subject), data = m2)
mod_BodyFat_num <- lmer(BodyFat ~ TimePointNum2 + (1| Subject), data = m2)
mod_weight_num <- lmer(Weight ~ TimePointNum2 + (1| Subject), data = m2)
我看到了一些使用purrr::map
的建议,最相似的方法是
mtcars %>%
names() %>%
paste('mpg ~', .) %>%
map(as.formula) %>%
map(lm, data = mtcars) %>%
map_df(~broom::tidy(.))
受上面例子的启发,所以我运行这些代码:
m2 %>%
names() %>%
paste(., '~ TimePoint + (1| Subject)') %>%
map(as.formula) %>%
map(lmer, data = m2) %>%
map(~broom::tidy(.))
但是我收到了一个错误(很可能是第5行的lmer
)
Error in `map()`:
ℹ In index: 1.
Caused by error in `mkRespMod()`:
! response must be numeric
Run `rlang::last_trace()` to see where the error occurred.
我的期望输出类似于下表:
| 模型||价值|Std.Error| DF| t值|p值|
| - -----|- -----|- -----|- -----|- -----|- -----|- -----|
| | | | | | | |
| 身体脂肪|(截距)|28.4583| 0.68122|三百七十五|41.7758| 0|
| | 时间点编号2| 0.06894| 0.02536|三百七十五|2.71861| 0.0069|
| 重量|(截距)|16.4157| 0.42985|三百七十五|38.1895| 0|
| | 时间点编号2| 0.08064| 0.01945|三百七十五|4.14603| 0|
| SFA|(截距)|67.1564| 1.79525|三百七十四|37.4079| 0|
| | 时间点编号2| 0.73767| 0.48301|三百七十四|1.52724| 0.1275|
| PUFA|(截距)|12.8161| 0.85069|三百六十五|十五点零六五六|0|
| | 时间点编号2|-0.5978| 0.23454|三百六十五|-2.5491| 0.0112|
我的问题是:
1.我们怎样才能纠正这个错误呢?
1.如何使用purrr::map
对每个感兴趣的列进行循环回归?
1.如何生成所需的表(如上所示,可能使用tidy()
)?
注:tidyverse友好方法优先
数据
m2 <- structure(
list(
Subject = c(
"B001",
"B001",
"B001",
"B001",
"B001",
"B002",
"B002",
"B002",
"B002",
"B002",
"B003",
"B003",
"B003",
"B003",
"B003",
"B004",
"B004",
"B004",
"B004",
"B004",
"B005",
"B005",
"B005",
"B005",
"B005",
"B006",
"B006",
"B006",
"B006",
"B006"
),
Weight = c(
68.5,
69.1,
70.2,
70,
70.1,
56.7,
56.1,
55.8,
56.5,
55.2,
62.9,
62.5,
62,
62.7,
63.3,
63.9,
62.6,
62.1,
61.5,
62.6,
53.5,
52,
51.8,
52.4,
52.9,
77,
77.1,
77.8,
79.2,
78.7
),
BodyFat = c(
21.8,
21.1,
21.7,
22.2,
22,
34,
33.3,
33,
32.4,
32.8,
36.9,
36.6,
36.1,
36.8,
36.9,
33.6,
32.5,
32.7,
31.9,
32.2,
12.2,
10.6,
11,
11.5,
11.6,
23,
22.5,
22.6,
23.1,
23.9
),
SFA = c(
62.895,
84.705,
83.49,
66.64,
72.19,
56.195,
93.945,
92.635,
88.51,
83.505,
92.67,
90.81,
83.37,
90.205,
84.195,
88.065,
53.69,
93.14,
52.57,
4.505,
63.505,
92.59,
80.87,
89.125,
67.305,
71.41,
89.1,
91.175,
80.635,
89.425
),
MUFA = c(
14.455,
11.71,
12.58,
21.135,
26.175,
13.285,
4.91,
6.08,
6.735,
12.745,
7.325,
7.605,
15.735,
7.985,
12.32,
7.92,
42.045,
4.57,
24.305,
0,
18.585,
5.955,
14.815,
8.775,
20.295,
11.485,
6.745,
6.32,
10.805,
6.8
),
TimePointNum2 = c(
1,
2,
3,
4,
5,
1,
2,
3,
4,
5,
1,
2,
3,
4,
5,
1,
2,
3,
4,
5,
1,
2,
3,
4,
5,
1,
2,
3,
4,
5
)
),
row.names = c(NA, -30L),
class = "data.frame"
)
2条答案
按热度按时间nimxete21#
要从混合效应模型获得整洁的输出,您需要使用
tidy()
from{broom.mixed}
。您可以按照nest-unnest工作流在每列上Map回归模型,并从中获得整洁的结果。cygmwpex2#