时间戳

首页 > 解决方案 > How can I loop a list of models to get slope estimate

r - How can I loop a list of models to get slope estimate

问题描述

I have a list of models as specified by the following code:

varlist <- list("PRS_Kunkle", "PRS_Kunkle_e07", 
   "PRS_Kunkle_e06","PRS_Kunkle_e05", "PRS_Kunkle_e04", 
   "PRS_Kunkle_e03", "PRS_Kunkle_e02", "PRS_Kunkle_e01",  
   "PRS_Kunkle_e00", "PRS_Jansen", "PRS_deroja_KANSL")

PRS_age_pacc3 <- lapply(varlist, function(x) {    
    lmer(substitute(z_pacc3_ds ~ i*AgeAtVisit  + i*I(AgeAtVisit^2)  +  
             APOE_score + gender + EdYears_Coded_Max20 +  
             VisNo + famhist + X1 + X2 + X3 + X4 + X5 +
             (1 |family/DBID),
       list(i=as.name(x))), data = WRAP_all, REML = FALSE)
})

I want to obtain the slope of PRS at different age points in each of the models. How can I write code to achieve this goal? Without loop, the raw code should be:

test_stat1 <- simple_slopes(PRS_age_pacc3[[1]], levels=list(AgeAtVisit=c(55,60,65,70,75,80)))
test_stat2 <- simple_slopes(PRS_age_pacc3[[2]], levels=list(AgeAtVisit=c(55,60,65,70,75,80)))
test_stat3 <- simple_slopes(PRS_age_pacc3[[3]], levels=list(AgeAtVisit=c(55,60,65,70,75,80)))
test_stat4 <- simple_slopes(PRS_age_pacc3[[4]], levels=list(AgeAtVisit=c(55,60,65,70,75,80)))
test_stat5 <- simple_slopes(PRS_age_pacc3[[5]], levels=list(AgeAtVisit=c(55,60,65,70,75,80)))
test_stat6 <- simple_slopes(PRS_age_pacc3[[6]], levels=list(AgeAtVisit=c(55,60,65,70,75,80)))
test_stat7 <- simple_slopes(PRS_age_pacc3[[7]], levels=list(AgeAtVisit=c(55,60,65,70,75,80)))
test_stat8 <- simple_slopes(PRS_age_pacc3[[8]], levels=list(AgeAtVisit=c(55,60,65,70,75,80)))
test_stat9 <- simple_slopes(PRS_age_pacc3[[9]], levels=list(AgeAtVisit=c(55,60,65,70,75,80)))
test_stat10 <- simple_slopes(PRS_age_pacc3[[10]], levels=list(AgeAtVisit=c(55,60,65,70,75,80)))
test_stat11 <- simple_slopes(PRS_age_pacc3[[11]], levels=list(AgeAtVisit=c(55,60,65,70,75,80)))

标签： rloopslme4mixed-models

解决方案

library(lme4)
library(reghelper)
set.seed(101)
## add an additional factor variable so we can use it for an interaction
sleepstudy$foo <- factor(sample(LETTERS[1:3], size = nrow(sleepstudy),
                         replace = TRUE))
m1 <- lmer(Reaction ~ Days*foo + I(Days^2)*foo + (1|Subject), data = sleepstudy)
s1 <- simple_slopes(m1, levels=list(Days = c(5, 10, 15)))

查看这些结果，s1是一个具有 6 行（级别数foo×Days指定值数）和 5 列（Days、foo、估计、标准误差、t 值）的数据框。

最简单的方法是：

res <- list()
for (i in seq_along(varlist)) {
   res[[i]] <- simple_slopes(model_list[[i]], ...)  ## add appropriate args here
}
res <- do.call("rbind", res)  ## collapse elements to a single data frame
## add an identifier column
res_final <- data.frame(model = rep(varlist, each = nrow(res[[1]])), res)

如果您想变得更漂亮，可以将for循环替换为适当的lapply. 如果你想比这更漂亮：

library(tidyverse)
(model_list
  %>% setNames(varlist)
  ## map_dfr runs the function on each element, collapses results to
  ##  a single data frame. `.id="model"` adds the names of the list elements
  ##  (set in the previous step) as a `model` column
  %>% purrr::map_dfr(simple_slopes, ... <extra args here>, .id = "model")
)

顺便说一句，simple_slopes当您在模型中也有二次项时，我会非常小心。计算的斜率将（可能）仅适用于模型中任何其他连续变量为零的情况。您可能希望将变量集中在 Schielzeth 2010生态学和进化方法中（“简单意味着改进......”）

推荐阅读