powerSim for a non existing fixed effect

[maelle]

Nice package, thanks!

I've just noticed this

library("simr")
#> Loading required package: lme4
#> Loading required package: Matrix
#> 
#> Attaching package: 'lme4'
#> The following object is masked from 'package:stats':
#> 
#>     sigma
#> Warning: replacing previous import 'lme4::sigma' by 'stats::sigma' when
#> loading 'pbkrtest'
#> Warning: replacing previous import 'lme4::sigma' by 'stats::sigma' when
#> loading 'simr'
fm1 <- lmer(y ~ x + (1|g), data=simdata)
powerSim(fm1, nsim=10, test = fixed("notavariable"))
#> Simulating: |                                                             |Simulating: |======                                                       |Simulating: |============                                                 |Simulating: |==================                                           |Simulating: |========================                                     |Simulating: |==============================                               |Simulating: |====================================                         |Simulating: |==========================================                   |Simulating: |================================================             |Simulating: |======================================================       |Simulating: |=============================================================|
#> Warning in observedPowerWarning(sim): This appears to be an "observed
#> power" calculation
#> Power for predictor 'notavariable', (95% confidence interval):
#>        0.00% ( 0.00, 30.85)
#> 
#> Test: Kenward Roger (package pbkrtest)
#> 
#> Based on 10 simulations, (10 warnings, 10 errors)
#> alpha = 0.05, nrow = 30
#> 
#> Time elapsed: 0 h 0 m 0 s
#> 
#> nb: result might be an observed power calculation
lastResult()$pval
#>  [1] NA NA NA NA NA NA NA NA NA NA

Session info

``` r devtools::session_info() #> Session info -------------------------------------------------------------- #> setting value #> version R version 3.3.1 (2016-06-21) #> system x86_64, mingw32 #> ui RTerm #> language (EN) #> collate Spanish_Spain.1252 #> tz Europe/Paris #> date 2017-04-13 #> Packages ------------------------------------------------------------------ #> package * version date source #> backports 1.0.5 2017-01-18 CRAN (R 3.3.2) #> binom 1.1-1 2014-01-02 CRAN (R 3.3.3) #> car 2.1-4 2016-12-02 CRAN (R 3.3.2) #> devtools 1.12.0 2016-06-24 CRAN (R 3.3.2) #> digest 0.6.12 2017-01-27 CRAN (R 3.3.2) #> evaluate 0.10 2016-10-11 CRAN (R 3.3.1) #> htmltools 0.3.5 2016-03-21 CRAN (R 3.2.4) #> iterators 1.0.8 2015-10-13 CRAN (R 3.2.3) #> knitr 1.15.1 2016-11-22 CRAN (R 3.3.2) #> lattice 0.20-34 2016-09-06 CRAN (R 3.3.2) #> lme4 * 1.1-12 2016-04-16 CRAN (R 3.2.5) #> magrittr 1.5 2014-11-22 CRAN (R 3.2.2) #> MASS 7.3-45 2015-11-10 CRAN (R 3.2.2) #> Matrix * 1.2-7 2016-08-28 CRAN (R 3.3.1) #> MatrixModels 0.4-1 2015-08-22 CRAN (R 3.2.5) #> memoise 1.0.0 2016-01-29 CRAN (R 3.2.3) #> mgcv 1.8-17 2017-02-08 CRAN (R 3.3.2) #> minqa 1.2.4 2014-10-09 CRAN (R 3.2.2) #> nlme 3.1-131 2017-02-06 CRAN (R 3.3.2) #> nloptr 1.0.4 2014-08-04 CRAN (R 3.2.2) #> nnet 7.3-12 2016-02-02 CRAN (R 3.3.1) #> pbkrtest 0.4-7 2017-03-15 CRAN (R 3.3.3) #> plotrix 3.6-4 2016-12-30 CRAN (R 3.3.2) #> plyr 1.8.4 2016-06-08 CRAN (R 3.2.5) #> quantreg 5.29 2016-09-04 CRAN (R 3.3.2) #> Rcpp 0.12.10 2017-03-19 CRAN (R 3.3.3) #> RLRsim 3.1-3 2016-11-04 CRAN (R 3.3.3) #> rmarkdown 1.4 2017-03-24 CRAN (R 3.3.3) #> rprojroot 1.2 2017-01-16 CRAN (R 3.3.2) #> simr * 1.0.2-1 2017-04-13 Github (pitakakariki/simr@a58d9b8) #> SparseM 1.76 2017-03-09 CRAN (R 3.3.3) #> stringi 1.1.3 2017-03-21 CRAN (R 3.3.3) #> stringr 1.2.0 2017-02-18 CRAN (R 3.3.3) #> withr 1.0.2 2016-06-20 CRAN (R 3.2.5) #> yaml 2.1.14 2016-11-12 CRAN (R 3.3.2) ```

Should the function output a warning in that case?

[pitakakariki]

If you look at lastResult()$errors you'll see error messages due to the incorrect test variable. The error message could certainly be more informative though.

I might add a warning if there are any errors caught in the simulations, it's easy to miss in the printout.

[maelle]

Oh true I had forgotten about that. Maybe I'd expect the function to actually fail when the user inputs a wrong fixed effect name?

Thanks for your extremely quick answer! :rocket:

[pitakakariki]

You might want to reinstall lme4 to get rid of those masking warnings. The version you have was built before they moved sigma from nlme to stats, tends to cause all sorts of subtle problems.

[pitakakariki]

I'd like it to fail in that case, but I still have to work out a good way to distinguish errors that should be logged from errors that should stop the simulation.

Also fixed doesn't always know which variables are legal, b/c

• no access to fit, which might be a function instead of a merMod anyway.
• categorical variables might mean different fixed effect names.

I'll leave the issue open until I've worked those bits out.

[maelle]

Thanks, I'm currently re-installing lme4 (from Github, no newer version on CRAN)

[maelle]

I thought I had another issue yesterday but I think it's now gone after the re-installation, no idea why, but thank you so much! :grin:

[liznecka]

I'm getting the same issue with powerSim output 0.00% (0.00, 30.85) even when I use a variable that is in my model, and the lastResult()$errors output seems to be an issue coming from simulate.MerMod... Are there clearer messages to explain the problem?

`

test<- lmer(vas~x + (1|subj), data = dfdata) power<- powerSim(test, test = fixed("x", method = "z"), nsim = 10) power Power for predictor 'x', (95% confidence interval): 0.00% ( 0.00, 30.85)

Test: z-test Effect size for x is 3.8

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

Time elapsed: 0 h 0 m 2 s lastResult()$errors stage index message 1 Simulating 1 new levels detected in newdata 2 Simulating 2 new levels detected in newdata 3 Simulating 3 new levels detected in newdata 4 Simulating 4 new levels detected in newdata 5 Simulating 5 new levels detected in newdata 6 Simulating 6 new levels detected in newdata 7 Simulating 7 new levels detected in newdata 8 Simulating 8 new levels detected in newdata 9 Simulating 9 new levels detected in newdata 10 Simulating 10 new levels detected in newdata `

[pitakakariki]

I think that error usually happens when your data has missing values. Could you please try it again but with an na.omit first?

dfdata2 <- na.omit(dfdata)

[liznecka]

Yes, that fixed it! Thank you so much! Great package!

[pitakakariki]

I've seen the error before but I've never been able to reproduce it. Would you mind sharing your dataset so I can track it down?

(Either here or by email, simr.peter at gmail.com)