袋裝AIC值來自混合效應線性迴歸模型R

Question

我已經用1000次數據集的替換重新採樣，現在想要爲這1000個數據集中的每一個擬合三個模型幷包裝它們的AIC分數。該程序的最終目標是獲得所有模型中每個模型的平均AIC評分及其95％置信區間。下面的代碼是錯誤的，我不知道我犯了什麼錯誤。發生的是，最終矩陣僅包含來自前幾次迭代的AIC得分向量（即，不是全部1000）。在每次迭代初始化主矩陣或向量時，是否有錯誤？或者，也許我的行添加過程有缺陷？或者，如果代碼是正確的，它是否可以與正在輸入此代碼的數據集一起使用？如果後者是這種情況，那麼爲什麼當代碼讀入這些數據集時只是跳過它們而沒有得到錯誤？我一直在努力掙扎幾天，感到非常困惑，所以任何幫助，將不勝感激。

require(lme4) 
require(lmerTest) 

# initializing an empty matrix for storing each vector of AIC scores from each iteration 
# the matrix has width 3 because three models are fitted at each iteration 
AIC.scores = data.frame(matrix(, nrow = 0, ncol = 3)) 

#fit regression models to each of 1000 datasets 
for(iter in 1:1000){ 

    #retrieving the data set, named accordingly, for the current iteration 
    data = read.csv(paste("data_set_", iter,".csv", sep=""), header=TRUE) 

    #initializing vector of AICs from models in current iteration 
    AIC.score = vector(mode="numeric", length=3) 

    mod1 = lmer(RT.log ~ crit.var1.log.std + 
         (1|Subject) + 
         (1|Item), 
          data = data, 
          REML=FALSE) 

    AIC.score[1] = summary(mod1)$AIC[1] 

    mod2 = lmer(RT.log ~ crit.var2.log.std + 
         (1|Subject) + 
         (1|Item), 
          data = data, 
          REML=FALSE) 

    AIC.score[2] = summary(mod2)$AIC[1] 

    mod3 = lmer(RT.log ~ crit.var3.log.std + 
         (1|Subject) + 
         (1|Item), 
          data = data, 
          REML=FALSE) 

    AIC.score[3] = summary(mod3)$AIC[1] 

    #adding vector of AICs scores from current iteration to main matrix 
    AIC.scores = rbind(AIC.scores, t(AIC.score)) 

    cat("bagging iteration", iter, "completed!\n") 

} 

#renaming column names in AIC score matrix 
colnames(AIC.scores) = c("model1", "model2", "model3") 

# function for calculating mean AIC and 95% C.I.s for each model across all iterations 
norm.interval = function(data, z=1.96) { 
    mean = mean(data) 
    variance = var(data) 
    sd = sqrt(variance/length(data)) 
    c(mean, mean - z * sd, mean + z * sd) 
} 

for (i in 1:3) { 

    cat("The mean, lCI, uCI for model", i, "are:", norm.interval(AIC.scores[,i]), "\n") 

} 

Answer 1

在黑暗中拍攝而不知道你的lmer模型是什麼或你的數據是什麼。

閱讀您的所有data.frame作爲一個單獨的列表：

all.data <- lapply(paste0("data_set_", 1:1000, ".csv"), read.csv, header=TRUE) 

生產出的文本形式的公式11聚物，基於上述數據的列表：

all.form <- lapply(paste0("all.data[[", 1:1000, "]]"), function(x) list(
    mod1 = paste0("lmer(RT.log ~ crit.var1.log.std + (1|Subject) + (1|Item), REML=FALSE, data =", x, ")"), 
    mod2 = paste0("lmer(RT.log ~ crit.var2.log.std + (1|Subject) + (1|Item), REML=FALSE, data =", x, ")"), 
    mod3 = paste0("lmer(RT.log ~ crit.var3.log.std + (1|Subject) + (1|Item), REML=FALSE, data =", x, ")") 
)) 

執行11聚物的公式並將模型寫入單個列表中：

all.lmer.mod <- lapply(all.form. function(x) lapply(x, function(y) eval(parse(text=y)))) 

提取lmer型號的所有AIC值：

all.AIC <- lapply(all.lmer.mod, function(x) lapply(x, AIC)) 

請注意，如果您的數據量較大，並且您有很多主題和項目級別，則過程將花費很長時間。通過首先更改1:1000進行測試，以便首先說明1:2。