EAMM: Simulation function for exploratory power analysis for random... In pamm: Power Analysis for Random Effects in Mixed Models

Description

Given a specific sample size, fixed number of group and replicates per group, the function simulate different variance-covariance structure and assess p-values and power of random intercept and random slope

Usage

 ```1 2 3``` ```EAMM(numsim, group, repl, fixed = c(0, 1, 0), VI = seq(0.05, 0.95,0.05), VS = seq(0.05, 0.5, 0.05), CoIS = 0, relIS = "cor", n.X, autocorr.X, X.dist, intercept = 0, heteroscedasticity = c("null"), mer.sim, mer.model) ```

Arguments

 `numsim` number of simulation for each step `group` number of group `repl` number of replicates per group `fixed` vector of lenght 3 with mean, variance and estimate of fixed effect to simulate. Default: `c(0,1,0)` `VI` variance component of intercept. Could be specified as a vector. Default: `seq(0.05,0.95,0.05)` `VS` variance component of slope. Could be specified as a vector. Default: `seq(0.05,0.5,0.05)` `CoIS` value of correlation or covariance between random intercept and random slope. Default: 0 `relIS` "cor" or "cov" set the type of relation give in CoIS. By default the relation is set to correlation `n.X` number of different values to simulate for the fixed effect (covariate). If `NA`, all values of X are independent between groups. If the value specified is equivalent to the number of replicates per group, `repl`, then all groups are observed for the same values of the covariate. Default: `NA` `autocorr.X` correlation between two successive covariate value for a group. Default: `0` `X.dist` specify the distribution of the fixed effect. Only "gaussian" (normal distribution) and "unif" (uniform distribution) are accepted actually. Default: `"gaussian"` `intercept` a numeric value giving the expected intercept value. Default: 0 `heteroscedasticity` a vector specifying heterogeneity in residual variance across X. If `c("null")` residual variance is homogeneous across X. If `c("power",t1,t2)` models heterogeneity with a constant plus power variance function. Letting v denote the variance covariate and s2(v) denote the variance function evaluated at v, the constant plus power variance function is defined as s2(v) = (t1 + |v|^t2)^2, where t1, t2 are the variance function coefficients. If `c("exp",t)`,models heterogeneity with an exponential variance function. Letting v denote the variance covariate and s2(v) denote the variance function evaluated at v, the exponential variance function is defined as s2(v) = exp(2* t * v), where t is the variance function coefficient. Default:"Null" `mer.sim` Use the simluate.merMod function to simulate the data. Potentially faster for large dataset but more restricted in terms of options `mer.model` Simulate the data based on a existing data and model structure from a lmer object. Should be specified as a list of 3 components: a mer object fitted via lmer, an environmental covariate for which to test the random slope, a random effect (e.g. `list(fm1,"Days","Subject"`)

Details

P-values for random effects are estimated using a log-likelihood ratio test between two models with and without the effect. Power represent the percentage of simulations providing a significant p-value for a given random structure. Residual variance, e, is calculted as 1-VI.

Value

data frame reporting estimated P-values and power with CI for random intercept and random slope

Warning

the simulation is based on a balanced data set with unrelated group

Julien Martin

References

Martin, Nussey, Wilson and Reale Submitted Measuring between-individual variation in reaction norms in field and experimental studies: a power analysis of random regression models. Methods in Ecology and Evolution.

`PAMM`, `SSF`, `plot.EAMM`
 ``` 1 2 3 4 5 6 7 8 9 10 11``` ```## Not run: ours <- EAMM(numsim=10,group=10,repl=4,fixed=c(0,1,1),VI=seq(0.1,0.3,0.05), VS=seq(0.05,0.2,0.05) ) plot(ours, "both") (fm1 <- lmer(Reaction ~ Days + (Days | Subject), sleepstudy)) ours2 <- EAMM(numsim=10, mer.model=list(model=fm1,env="Days",random="Subject"), VI=seq(0.3,0.5,0.1), VS=seq(0.05,0.2,0.05) ) plot(ours2, "both") ## End(Not run) ```