DejanDraschkow / mixedpower

GNU General Public License v3.0
16 stars 13 forks source link

Run Time Ongoing #2

Closed mpipoly closed 1 year ago

mpipoly commented 2 years ago

Hello,

I am attempting to perform a post-hoc power analysis on a complex linear mixed effects model from real data. We want to see how powered we are with each effect of interest. The model has a total of 6 predictor terms and two random effects. The full model assessed a three-way interaction and thus contains 15 total estimates (including the constant).

I treated the following variables as:

model <- lme4 model object
data <-  exact data frame used to create lme4 model object (was filtered to only consider accurate responses)
fixed_effects <- vector of individual variables terms, as strings used in the model : example c("variable 1", "variable 2", etc..)
simvar <- Subject dummy coded, as used in lme4 

# SIMULATION PARAMETERS
steps <- c(30, 40, 50, 60, 70, 90, 120, 160, 180) # sample sizes that we want to estimate!
critical_value <- 2 # which t/z value do we want to use to test for significance?
n_sim <- 1000 # how many single simulations should be used to estimate power?

# ------------------------------------------ #
# INCLUDE SESOI SIMULATION

SESOI <- vector of 15 effect sizes that correspond to the 15 betas produced by the lme4 model: example c(45.3, .076, etc..)

# ------------------------------------------ #
# RUN SIMULATION
level_sim <- mixedpower(model = model, data = subset(data, stimulus.ACC == 1),
                              fixed_effects = fixed_effects,
                              simvar = "dum_sum", steps = steps,
                              critical_value = 2.0, n_sim = n_sim,
                              SESOI = SESOI, databased = T)

My issue is that for the past few days I have been stuck on "step 30" even though my r-script is running on a high-performance computer with 32 cores and 500 gigabytes of RAM. It may be that this is expected considering what I asked mixed power to do but I am unsure of why that is the case and how it is occurring under the hood. Could this be solved with a verbose output? Is this a flag I am missing that I should know about?

Any advice and direction here is greatly appreciated!

lkumle commented 2 years ago

Hi,

this response might come to late, but I would recommend reducing n_sim to e.g., n_sim = 5 and steps to c(30, 40) to get an idea of how long the simulations take first. If the model/data is very large, mixedpower does indeed take quite some time. If the smaller (n_sim = 5) simulations run without any issues, the runtime might just be a product of your complex model. However, it should usually not take "a few days".

mpipoly commented 2 years ago

Thanks for the reply. I did try something like (but not that low of sims) that using all of my model terms and covariates. It just occurred to me now I might want to try the simulation with just the interaction of interest alongside the random effects to see if reducing terms helps reduce time (to figure out a baseline to run).

I will ping back what I find out. Any other advice that comes to mind is greatly appreciated!

lkumle commented 2 years ago

You could also check out this Julia package, which should offer faster simulations!