R - использовать simsem на карте

У меня есть список из более чем 100 моделей SEM, рассчитанных за lavaan. В этом примере я буду использовать только две модели:

fit.model1.cfa <- '
    ind60 =~ x1 + x2 + x3
    dem60 =~ y1 + y2 + y3 + y4
    dem65 =~ y5 + y6 + y7 + y8
    ind60 ~~ ind60
    dem60 ~~ dem60
    dem65 ~~ dem65' %>%
       sem(PoliticalDemocracy)

fit.model2.cfa <- '
    ind60 =~ x1 + x2
    dem60 =~ y1 + y2 + y3
    dem65 =~ y5 + y6 + y7
    ind60 ~~ ind60
    dem60 ~~ dem60
    dem65 ~~ dem65' %>%
       sem(PoliticalDemocracy)

Не так уж важно, что это за модели на самом деле. Теперь мы объединим их в список:

models <- list(fit.model1.cfa, fit.model2.cfa)

Я хочу провести апостериорный анализ мощности на своих моделях, и для этой цели я использую следующий цикл for:

pwr <- list()

for(i in 1:2){
   sim(nRep = 1000, model = models[[i]], n = models[[i]]@SampleStats@ntotal, generate = models[[i]]) %>%
       getPower() -> pwr[[i]]
}

Однако я бы хотел избежать цикла for и вместо этого использовать функцию map. Я пробовал следующее:

models %>%
   map(., sim, nRep = 1000, model = ., n = pluck(pluck(., "SampleStats"), "ntotal"), generate = .) %>%
   map(getPower)

Однако это не работает. Так что я делаю не так?


person J. Doe    schedule 24.02.2020    source источник
comment
Из какого пакета поступает функция sim?   -  person rdornas    schedule 24.02.2020
comment
@rdornas из пакета simsem   -  person J. Doe    schedule 24.02.2020


Ответы (1)


Я впервые слышу о lavaan, simsem и их моделях, но я пытаюсь помочь вам с помощью purrr map. Я сделал простое решение с nRep = 10 просто для тестирования, и вы можете воспроизвести его с помощью nRep = 1000. Также было возвращено много предупреждений и прогресса, но я удалил их в приведенном ниже примере для чистоты.

library(tidyverse, verbose = F)
library(lavaan)
library(simsem)

fit.model1.cfa <- '
    ind60 =~ x1 + x2 + x3
    dem60 =~ y1 + y2 + y3 + y4
    dem65 =~ y5 + y6 + y7 + y8
    ind60 ~~ ind60
    dem60 ~~ dem60
    dem65 ~~ dem65' %>%
  sem(PoliticalDemocracy)

fit.model2.cfa <- '
    ind60 =~ x1 + x2
    dem60 =~ y1 + y2 + y3
    dem65 =~ y5 + y6 + y7
    ind60 ~~ ind60
    dem60 ~~ dem60
    dem65 ~~ dem65' %>%
  sem(PoliticalDemocracy)

models <- list(model1 = fit.model1.cfa, model2 = fit.model2.cfa)


power <- map(models, function(x){

  n <- x@SampleStats@ntotal

  z <- sim(model = x, nRep = 10, n = n, generate = x)

  getPower(z)

})

power

#> $model1
#>    ind60=~x2    ind60=~x3    dem60=~y2    dem60=~y3    dem60=~y4    dem65=~y6 
#>          1.0          1.0          1.0          1.0          1.0          1.0 
#>    dem65=~y7    dem65=~y8 ind60~~ind60 dem60~~dem60 dem65~~dem65       x1~~x1 
#>          1.0          1.0          1.0          1.0          1.0          1.0 
#>       x2~~x2       x3~~x3       y1~~y1       y2~~y2       y3~~y3       y4~~y4 
#>          0.5          1.0          1.0          1.0          1.0          1.0 
#>       y5~~y5       y6~~y6       y7~~y7       y8~~y8 ind60~~dem60 ind60~~dem65 
#>          1.0          1.0          1.0          1.0          1.0          1.0 
#> dem60~~dem65 
#>          1.0 
#> 
#> $model2
#>    ind60=~x2    dem60=~y2    dem60=~y3    dem65=~y6    dem65=~y7 ind60~~ind60 
#>          1.0          1.0          1.0          1.0          1.0          1.0 
#> dem60~~dem60 dem65~~dem65       x1~~x1       x2~~x2       y1~~y1       y2~~y2 
#>          1.0          1.0          0.4          0.4          0.8          1.0 
#>       y3~~y3       y5~~y5       y6~~y6       y7~~y7 ind60~~dem60 ind60~~dem65 
#>          1.0          1.0          1.0          1.0          0.9          1.0 
#> dem60~~dem65 
#>          1.0
person rdornas    schedule 24.02.2020