Skip to main content
4 of 4
added 124 characters in body
JElder
  • 1.3k
  • 8
  • 17

PowerSim function (simr package) indicates 100% power for small sample dataset, regardless of effect size, alpha level, or if fixed/random slopes

I am using the simr package to do power analyses for lmer multilevel models I have run, to determine the power of a pilot dataset for future research.

The dataset consists of 46 subjects with approximately 148 trials within subjects at level 1, and I've tried the power simulations at the original sample, at an extended sample of 800, with fixed slopes, random slopes, and at alpha of .001 and .05. The effect size from the original dataset was .10, but I have also expanded it to 1 and 4.

Reliably, I am getting a power of 100% for every variant of this simulation I've run. I did a power curve analysis and plotted it as well. It suggests that I have above 80% power by as low as 10 subjects for the random slopes model and as low as 2-3 subjects for the fixed slopes model. Frankly, this seems very odd to me and doesn't make much sense. The only post hoc justification I can make here is that my level 1 sample is decently large (148 per level 2 sample).

Primarily, I am confused how/why I am getting such high power and if it is because I executed something incorrectly or there is something off with my data or execution. If it is not something I did wrong, why would a sample of 2-3 subjects (with 148 trials within) at fixed slopes and 10 subjects at random slopes give you power of over 80%? That does not make sense to me and seems very off. It also seems off that I reliably get a power of 100% regardless of what parameters I change.

Any insights or thoughts would be great.

Further details are below...

Random slopes model appears as:

Linear mixed model fit by REML ['lmerMod']
Formula: selfRespT2 ~ feedback + selfRespT1 + valEstF + (feedback + selfRespT1 +      valEstF | subID)
   Data: fullDf1
Control: lmerControl(optimizer = "bobyqa", optCtrl = list(maxfun = 2e+05))

REML criterion at convergence: 21206.9

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-4.5158 -0.6078  0.0341  0.6287  3.6362 

Random effects:
 Groups   Name        Variance Std.Dev. Corr             
 subID    (Intercept) 1.450505 1.20437                   
          feedback    0.005708 0.07555  -0.55            
          selfRespT1  0.030467 0.17455  -0.78  0.15      
          valEstF     0.007138 0.08448  -0.73  0.53  0.19
 Residual             1.449100 1.20379                   
Number of obs: 6525, groups:  subID, 46

Fixed effects:
            Estimate Std. Error t value
(Intercept)  0.98865    0.20401   4.846
feedback     0.10127    0.01313   7.712
selfRespT1   0.51978    0.02795  18.599
valEstF      0.16395    0.02160   7.590

Correlation of Fixed Effects:
           (Intr) fedbck slfRT1
feedback   -0.458              
selfRespT1 -0.662  0.074       
valEstF    -0.698  0.268  0.016

Then simulation appears as:

sim2 <- powerSim(MLM.1.1, fixed("feedback", "z"), seed = 2, nsim = 800, alpha = .05)

Power for predictor 'feedback', (95% confidence interval):
      100.0% (99.54, 100.0)

Test: z-test
      Effect size for feedback is 0.10

Based on 800 simulations, (4 warnings, 0 errors)
alpha = 0.05, nrow = 6525

Time elapsed: 0 h 37 m 36 s

nb: result might be an observed power calculation

Fixed slopes model appears as:

Linear mixed model fit by REML ['lmerMod']
Formula: selfRespT2 ~ feedback + selfRespT1 + valEstF + (1 | subID)
   Data: fullDf1
Control: lmerControl(optimizer = "bobyqa", optCtrl = list(maxfun = 2e+05))

REML criterion at convergence: 21547.8

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-4.3314 -0.6033  0.0533  0.6255  3.6128 

Random effects:
 Groups   Name        Variance Std.Dev.
 subID    (Intercept) 0.1058   0.3252  
 Residual             1.5597   1.2489  
Number of obs: 6525, groups:  subID, 46

Fixed effects:
            Estimate Std. Error t value
(Intercept) 0.941441   0.107714    8.74
feedback    0.097369   0.007106   13.70
selfRespT1  0.516865   0.010486   49.29
valEstF     0.187282   0.018181   10.30

Correlation of Fixed Effects:
           (Intr) fedbck slfRT1
feedback   -0.179              
selfRespT1 -0.120 -0.207       
valEstF    -0.758 -0.006 -0.301

Simulation for fixed slopes model appears as:

sim4 <- powerSim(MLM.1.2, fixed("feedback", "z"), seed = 2, nsim = 800, alpha = .05)

Power for predictor 'feedback', (95% confidence interval):
      100.0% (99.54, 100.0)

Test: z-test
      Effect size for feedback is 0.10

Based on 800 simulations, (4 warnings, 0 errors)
alpha = 0.05, nrow = 6525

Time elapsed: 0 h 36 m 55 s

nb: result might be an observed power calculation

I receive the same 100% power when I manually change the effect size of .1 to 1 or to 4, as well as when I change the alpha to .001 instead of .05.

Power curve analysis for random slopes model appears as:

enter image description here

Power curve analysis for fixed slopes model appears as:

enter image description here

UPDATE: I manipulated the effect size and here too it rapidly gets to 100% power.

Fixed slopes, changing effect size

enter image description here

Random slopes, changing effect size

enter image description here

Visualized changing effect sizes:

enter image description here

JElder
  • 1.3k
  • 8
  • 17