search

pitakakariki / simr

Power Analysis of Generalised Linear Mixed Models by Simulation
69 stars 19 forks source link

nrow=Na and Error in all simulations #144

Open Hayafh opened 5 years ago

Hayafh commented 5 years ago

Hi, I am trying to do a power analysis for the following model: BPD_FINAL ~ SEX + AGE + SES + ethnicity

I am trying to predict offspring BPD (BPD_FINAL) with parent substance use (DAD_DRUGS_SX) and seeing if the adoptive status (idabd) of the offspring moderate the association (DAD_DRUGS_SX*idabd). The data has sibling pairs, nested within families (IDYrFamf). I am trying to use SIMR to figure out what effect size I’m powered to detect. Here's what I have tried, and what I get:

library(simr) library(lme4) model1<- glmer(BPD_FINAL ~ SEX + AGE + SES + ethnicity

       + idabd + AGE_M + AGE_D + 
         DAD_DRUGS_SX + DAD_DRUGS_SX*idabd +(1|IDYrFamf), data=zinb, na.action = na.exclude)

Warning message: In glmer(BPD_FINAL ~ SEX + AGE + SES + ethnicity + idabd + AGE_M + : calling glmer() with family=gaussian (identity link) as a shortcut to lmer() is deprecated; please call lmer() directly

summary(model1) Linear mixed model fit by REML ['lmerMod'] Formula: BPD_FINAL ~ SEX + AGE + SES + ethnicity + idabd + AGE_M + AGE_D + DAD_DRUGS_SX + DAD_DRUGS_SX * idabd + (1 | IDYrFamf) Data: zinb

REML criterion at convergence: 7403.1

Scaled residuals: Min 1Q Median 3Q Max -2.44186 -0.65810 -0.02018 0.59392 3.00327

Random effects: Groups Name Variance Std.Dev. IDYrFamf (Intercept) 13.11 3.620 Residual 72.40 8.509 Number of obs: 1019, groups: IDYrFamf, 520

Fixed effects: Estimate Std. Error t value (Intercept) 44.20132 4.56210 9.689 SEX -2.26402 0.59287 -3.819 AGE 0.02689 0.15163 0.177 SES 0.01861 0.72896 0.026 ethnicityKorean -1.82528 0.96996 -1.882 ethnicityOther -2.93594 1.71184 -1.715 idabd -3.24498 0.99987 -3.245 AGE_M -0.05144 0.11631 -0.442 AGE_D 0.04405 0.10230 0.431 DAD_DRUGS_SX 0.08398 0.50120 0.168 idabd:DAD_DRUGS_SX 1.26421 0.57741 2.189

powerSim(model1) Power for predictor 'SEX', (95% confidence interval):================================================================| 0.00% ( 0.00, 0.37)

Test: Kenward Roger (package pbkrtest) Effect size for SEX is -2.3

Based on 1000 simulations, (0 warnings, 1000 errors) alpha = 0.05, nrow = NA

Time elapsed: 0 h 1 m 12 s

I have tried na.omit, and specifying a certain variable for fixed effects (fixef(model1)["SEX"] fixef(model1)["SEX"]<- -2.26402) but I sitll get the same error, and power of 0. What am I missing? :)

pitakakariki commented 5 years ago

Try fitting your model with lmer, as per the warning above and see if that helps.

If you're still getting errors, have a look at lastResult()$errors to see what kind of errors you're getting.

munoztd0 commented 3 years ago

Same issue here when fitted with lmer lastResult()$errors outputs: Simulating object is not a matrix