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][1]][1]

Power curve analysis for fixed slopes model appears as:

[![enter image description here][2]][2]


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

Fixed slopes, changing effect size

[![enter image description here][3]][3]

Random slopes, changing effect size

[![enter image description here][4]][4]

Visualized changing effect sizes:

[![enter image description here][5]][5]


  [1]: https://i.sstatic.net/A6gYd.png
  [2]: https://i.sstatic.net/GCTPf.png
  [3]: https://i.sstatic.net/Mq9JP.png
  [4]: https://i.sstatic.net/E7Wej.png
  [5]: https://i.sstatic.net/p6neQ.png