Могут ли результаты моделей, пригодных для умножения вмененных наборов данных, быть извлечены в кадр данных? - PullRequest
Новая
       101

Могут ли результаты моделей, пригодных для умножения вмененных наборов данных, быть извлечены в кадр данных?

0 голосов
/ 06 сентября 2018

Можно ли извлечь в кадр данных объединенные оценки из нескольких моделей, пригодных для умножения вмененных данных?

Ниже показано, как я выполняю это для полного фрейма данных случая (т. Е. Без отсутствующих данных) - я хотел бы выполнить аналогичный процесс для извлечения аналогичных результатов нескольких моделей, подходящих для вмененных данных: 

library(tidyverse)
library(broom)
library(mice)

data <- nhanes
sapply(data, function(x) sum(is.na(x))) #check missing data
data <- data %>% filter(bmi !="NA" & hyp != "NA" & chl != "NA") # remove missing data

out <-c("bmi")
exp <- c("chl","age","factor(hyp)")

#run models and extract to tidy data frame
models <- expand.grid(out, exp) %>%
group_by(Var1) %>% rowwise() %>%
summarise(frm = paste0(Var1, "~", Var2)) %>%
group_by(model_id = row_number(),frm) %>%
do(tidy(lm(.$frm, data = data))) %>%
mutate(lci = estimate-(1.96*std.error),
     uci = estimate+(1.96*std.error))

Ниже приведен пример вычисления отсутствующих данных с использованием mice и подгонки только одной регрессионной модели: 

# Impute missing data using mice
data <- nhanes
imp <- mice(data, print = F)

#Fit single model
fit <- with(imp, lm(bmi ~ chl))

#Get pooled estimates
a <- pool(fit)

summary(a)

1 Ответ

0 голосов
/ 06 сентября 2018

Ключевым моментом здесь является начало с complete(imp, "long"), поскольку оно предоставляет все вмененные наборы данных. После этого вам придется немного поиграть с некоторыми функциями tidyverse и broom, особенно nest() и tidy(), которые здесь очень полезны. Попробуйте это: 

library(tidyverse)
library(broom)
library(mice)
data <- nhanes # data
imp <- mice(data, print = F) # imputation
# complete data
data.complete <- complete(imp, "long")  
glimpse(data.complete) # all the 5 imputations are here

data.complete %>% 
        select(-.id) %>% 
        nest(-.imp) %>%
        mutate(model = map(data, ~lm(bmi ~ chl, data = .)),
               tidied = map(model, tidy)) %>%
        unnest(tidied) %>%
        filter(term == "chl") %>%
        mutate(adjusted = p.adjust(p.value),
               lci = estimate-(1.96*std.error),
               uci = estimate+(1.96*std.error))
# output
  .imp term   estimate  std.error statistic     p.value   adjusted           lci        uci
1    1  chl 0.01972747 0.01755024  1.124057 0.272584078 0.45430932 -0.0146709916 0.05412594
2    2  chl 0.02133664 0.01719462  1.240891 0.227154661 0.45430932 -0.0123648105 0.05503808
3    3  chl 0.03070542 0.01534959  2.000407 0.057397701 0.22512674  0.0006202261 0.06079062
4    4  chl 0.04109955 0.02044568  2.010183 0.056281686 0.22512674  0.0010260220 0.08117308
5    5  chl 0.05448964 0.01585764  3.436175 0.002251967 0.01125984  0.0234086522 0.08557062
...