Я пытаюсь провести анализ мощности для проекта, где мы ожидаем 40 000 наблюдений (по крайней мере) на 40-60 "единицах".
Когда я делаю этот анализ мощности, используя гораздо меньшее число наблюдения в единицах (скажем, 100 на единицу), моделирование работает нормально. Однако требования к памяти, кажется, взрываются, когда я увеличиваю количество наблюдений на единицу до 300 или 500. Мне не удалось запустить симуляцию даже на компьютере с пределом памяти 70 ГБ.
Что происходит на? Почему требования к памяти так высоки и как я могу сделать их более управляемыми?
Ниже приведен общий скрипт имитации мощности c для справки.
#psacr001_power.R: Run a job using specifications determined using a shell script
#Get the arguments passed from the shell script
args <- commandArgs(trailingOnly=TRUE)
job <- as.numeric(args[[1]])
sample_size <- as.numeric(args[[2]])
number_of_labs <- as.numeric(args[[3]])
icc_size <- as.numeric(args[[4]])
slope_size <- as.numeric(args[[5]])
effect <- as.numeric(args[[6]])
#Required R libraries
library(lme4)
library(simr)
library(pbkrtest)
library(tibble)
#create dataset
simulated.df <- tibble(
x = sample(rep(0:1, sample_size * number_of_labs / 2)), #group assignment
g = rep(1:number_of_labs, sample_size), #number of clusters
y = rnorm(sample_size * number_of_labs) #outcome
)
#Create simple model
model.of.interest <- lmer(y ~ x + (x|g), data=simulated.df)
#Create a small effect of interest
fixef(model.of.interest)['x'] <- effect
#create various ICC intercept
VarCorr(model.of.interest)["g"][["g"]][1] <- icc_size
#create varioous slopes
VarCorr(model.of.interest)["g"][["g"]][4] <- slope_size
#try not to be singular
attr(VarCorr(model.of.interest)[["g"]], "correlation")[2:3] <- .3
power_summary <- tibble(
job = numeric(),
successes = numeric(),
trials = numeric(),
mean = numeric(),
lower = numeric(),
upper = numeric()
)
#simulate those models
temp <- powerSim(fit = model.of.interest,
nsim = 200,
progress=FALSE
)
power_summary[1 , 1] <- job
power_summary[1 , 2:6] <- summary(temp)
write.csv(power_summary, paste(c("res_", job, ".csv"), collapse=""), row.names=FALSE)