Ich habe ein multiples imputiertes Datenset (mit 10 sets) und möchte nun eine LMM in R rechnen. Der erste Schritt wäre die verschiedenen Nullmodelle (mit verschiedenen Levels - Schüler, Klasse, Schule) aufzustellen, diese bezüglich -2log lik zu vergleichen, um zu entscheiden welche Model Variante die passendste ist, um weiterzurechnen. Damit ich die Modelle vergleichen kann möchte ich gerne gepoolte Outputs von meinen 10 Datensets (-2 log lik, AIC, BIC, Estimates etc.). Ich habe mal gestartet mit folgendem Code:
Code: Alles auswählen
##########################################################################
#2.0 Read Datafile
##########################################################################
## 2.1) .sav from SPSS:
mi_long <- read_sav("name_datafile.sav")
outcome <- "P_dur_day_total_min_MVPA_WD"
## 2.2) Delete original Data with missing values (Imputation_=0)
if ("Imputation_" %in% names(mi_long)) {
mi_long <- mi_long %>% filter(Imputation_ > 0)
}
## 2.3) Activate Mahalanobis-Filter - Multiple Ausreisser werden entfernt
mi_long <- mi_long %>%
filter(!is.na(filter_1)) %>%
filter(filter_1 == 1)
## 2.4) Control: Cases per Imputation after Filter Mahalo
mi_long %>% count(Imputation_) %>% print(n = Inf)
## 2.5) Split in List of 10 Datasets
imp_list <- split(mi_long, mi_long$Imputation_)
# Spaltennamen säubern und Typen setzen
need <- c(outcome, "Code", "Klasse_1", "Schule_1", "MZP_1")
imp_list <- lapply(imp_list, function(d){
names(d) <- trimws(names(d)) # evtl. Leerzeichen loswerden
stopifnot(all(need %in% names(d))) # brich ab, wenn was fehlt
d %>%
mutate(
Code = droplevels(factor(trimws(Code))),
Klasse_1 = droplevels(factor(trimws(Klasse_1))),
Schule_1 = droplevels(factor(trimws(Schule_1))),
MZP_1 = droplevels(factor(MZP_1)),
!!outcome := as.numeric(haven::zap_labels(haven::zap_missing(.data[[outcome]])))
)
})
## 2.7) mitml-compatible list
imp_mi <- as.mitml.list(imp_list)
##########################################################################
#3.0 Compare Nullmodel LEVEL 2, 3, 4
##########################################################################
## 3.1) Modellform (Nullmodel WITH REML)
nullmodell_1 <- lmer(P_dur_day_total_min_MVPA_WD ~ 1 + (1|Code),
data = imp_list[[1]], REML = TRUE)
summary(nullmodell_1)
nullmodell_2 <- lmer(P_dur_day_total_min_MVPA_WD ~ 1 + (1|Klasse_1/Code),
data = imp_list[[1]], REML = TRUE)
summary(nullmodell_2)
nullmodell_3 <- lmer(P_dur_day_total_min_MVPA_WD ~ 1 + (1|Schule_1/Klasse_1/Code),
data = imp_list[[1]], REML = TRUE)
summary(nullmodell_3)
Die Outputs sehen wie folgt aus:
> summary(nullmodell_1)
Linear mixed model fit by REML ['lmerMod']
Formula: P_dur_day_total_min_MVPA_WD ~ 1 + (1 | Code)
Data: imp_list[[1]]
REML criterion at convergence: 5856.1
Scaled residuals:
Min 1Q Median 3Q Max
-2.0224 -0.5515 -0.0530 0.4790 3.8756
Random effects:
Groups Name Variance Std.Dev.
Code (Intercept) 3.517 1.875
Residual 3.610 1.900
Number of obs: 1249, groups: Code, 736
Fixed effects:
Estimate Std. Error t value
(Intercept) 8.9330 0.0885 100.9
> nullmodell_2 <- lmer(P_dur_day_total_min_MVPA_WD ~ 1 + (1|Klasse_1/Code),
+ data = imp_list[[1]], REML = TRUE)
> summary(nullmodell_2)
Linear mixed model fit by REML ['lmerMod']
Formula: P_dur_day_total_min_MVPA_WD ~ 1 + (1 | Klasse_1/Code)
Data: imp_list[[1]]
REML criterion at convergence: 5838.4
Scaled residuals:
Min 1Q Median 3Q Max
-2.12781 -0.54841 -0.04246 0.47018 2.94882
Random effects:
Groups Name Variance Std.Dev.
Code:Klasse_1 (Intercept) 3.4458 1.8563
Klasse_1 (Intercept) 0.3647 0.6039
Residual 3.3648 1.8343
Number of obs: 1249, groups: Code:Klasse_1, 786; Klasse_1, 45
Fixed effects:
Estimate Std. Error t value
(Intercept) 8.958 0.125 71.64
> nullmodell_3 <- lmer(P_dur_day_total_min_MVPA_WD ~ 1 + (1|Schule_1/Klasse_1/Code),
+ data = imp_list[[1]], REML = TRUE)
> summary(nullmodell_3)
Linear mixed model fit by REML ['lmerMod']
Formula: P_dur_day_total_min_MVPA_WD ~ 1 + (1 | Schule_1/Klasse_1/Code)
Data: imp_list[[1]]
REML criterion at convergence: 5829.9
Scaled residuals:
Min 1Q Median 3Q Max
-2.1006 -0.5449 -0.0332 0.4707 2.9786
Random effects:
Groups Name Variance Std.Dev.
Code:Klasse_1:Schule_1 (Intercept) 3.4565 1.8592
Klasse_1:Schule_1 (Intercept) 0.0709 0.2663
Schule_1 (Intercept) 0.3171 0.5631
Residual 3.3615 1.8334
Number of obs: 1249, groups: Code:Klasse_1:Schule_1, 786; Klasse_1:Schule_1, 45; Schule_1, 12
Fixed effects:
Estimate Std. Error t value
(Intercept) 9.0247 0.1902 47.45!!Meines Erachtens werden nun aber nicht gepoolte Ergebnisse ausgegeben sondern die Ergebnisse vom ersten imputierten Datenset? Wie müsste ich den Code anpassen sodass ich gepoolte Ergebnisse bekomme sodass ich die Modelle anschliessend vergleichen kann?
Merci für die Hilfe!