simulate.plcp: Perform a simulation study using a `study_parameters`-object

Description

Perform a simulation study using a study_parameters-object

Usage

# S3 method for plcp
simulate(object, nsim, seed = NULL, formula = NULL,
  satterthwaite = FALSE, CI = FALSE, cores = 1, progress = FALSE,
  batch_progress = TRUE, ...)
# S3 method for plcp_multi
simulate(object, nsim, seed = NULL,
  formula = NULL, satterthwaite = FALSE, CI = FALSE, cores = 1,
  progress = FALSE, batch_progress = TRUE, ...)

Arguments

object

An object created by study_parameters.

nsim

The number of simulations to run.

seed

Currently ignored.

formula

Model formula(s) used to analyze the data, see Details. Should be created using sim_formula. It is also possible to compare multiple models, e.g. a correct and a misspecified model, by combining the formulas using sim_formula_compare. See Examples. If NULL the formula is made automatically, using create_lmer_formula, which does not support objects with multiple simulation setups.

satterthwaite

Logical; if TRUE Satterthwaite's degrees of freedom approximation will be used when computing p-values. This is implemented using the lmerTest-package. See Details.

Logical; if TRUE coverage rates for 95 % confidence intervals will be calculated. See Details.

cores

Number of CPU cores to use. If called from a GUI environment (e.g. RStudio) or a computer running Microsoft Windows, PSOCK clusters will be used. If called from a non-interactive Unix environment forking is utilized.

progress

logical; will display progress if TRUE. Currently ignored on Windows. Package pbmclapply is used to display progress, which relies on forking. N.B using a progress bar will noticeably increase the simulation time, due to the added overhead.

batch_progress

logical; if TRUE progress will be shown for simulations with multiple setups.

...

Optional arguments, see Saving in Details section.

Details

See also vignette("simulations", package = "powerlmm") for a tutorial.

Model formula

If no data transformation is used, the available model terms are:

y the outcome vector, with potential missing data.
y_c the complete version of y, before dropout was simulated.
time the time vector.
treatment treatment indicator (0 = "control", 1 = "treatment").
subject subject-level id variable, from 1 to total number of subjects.
cluster for three-level models; the cluster-level id variable, from 1 to the total number of clusters.

See Examples and the simulation-vignette for formula examples. For objects that contain a single study setup, then the lmer formula can be created automatically using create_lmer_formula.

Satterthwaite's approximation, and CI coverage

To decrease the simulation time the default is to only calculate Satterthwaite's dfs and the CIs' coverage rates for the test of 'time:treatment'-interaction. This can be changed using the argument test in sim_formula.

Confidence intervals are both calculated using profile likelihood and by the Wald approximation, using a 95 % confidence level.

Saving intermediate results for multi-sims

Objects with multi-sims can be save after each batch is finished. This is highly recommended when many designs are simulated. The following additional arguments control saving behavior:

'save', logical, if TRUE each batch is saved as a RDS-file. Results are saved in your working directory, in the directory specified by save_folder.
'save_folder' a character indicating the folder name. Default is 'save'.
'save_folder_create', logical, if TRUE then save_folder will be created if it does not exist in your working directory.

Examples

Run this code

# NOT RUN {
# Two-level ---------------------------------------------------------------
p <- study_parameters(n1 = 11,
                      n2 = 25,
                      sigma_subject_intercept = 1.44,
                      sigma_subject_slope = 0.2,
                      sigma_error = 1.44,
                      effect_size = cohend(0.5))

f <- sim_formula("y ~ treatment * time + (1 + time | subject)")


res <- simulate(object = p,
                nsim = 1000,
                formula = f,
                satterthwaite = TRUE,
                progress = FALSE,
                cores = 1,
                save = FALSE)

summary(res)


# Three-level (nested) ------------------------------------------------------
p <- study_parameters(n1 = 10,
                      n2 = 20,
                      n3 = 4,
                      sigma_subject_intercept = 1.44,
                      icc_pre_cluster = 0,
                      sigma_subject_slope = 0.2,
                      icc_slope = 0.05,
                      sigma_error = 1.44,
                      effect_size = 0)

## compare correct and miss-specified model
f0 <- "y ~ treatment * time + (1 + time | subject)"
f1 <- "y ~ treatment * time + (1 + time | subject) + (0 + time | cluster)"
f <- sim_formula_compare("correct" = f1,
                         "wrong" = f0)

res <- simulate(object = p,
                nsim = 1000,
                formula = f,
                satterthwaite = TRUE,
                progress = FALSE,
                cores = 1,
                save = FALSE)

summary(res)

## Compare random effects using LRT,
## summarise based on best model from each sim
summary(res,
        model_selection = "FW",
        LRT_alpha = 0.1,
        para = "treatment:time")

# Partially nested design ---------------------------------------------------
p <- study_parameters(n1 = 11,
                      n2 = 10,
                      n3 = 4,
                      sigma_subject_intercept = 1.44,
                      icc_pre_cluster = 0,
                      sigma_subject_slope = 0.2,
                      cor_subject = -0.5,
                      icc_slope = 0.05,
                      sigma_error = 1.44,
                      partially_nested = TRUE,
                      effect_size = cohend(-0.5))

f <- sim_formula("y ~ treatment * time + (1 + time | subject) +
                  (0 + treatment:time | cluster)")

res <- simulate(object = p,
                nsim = 1000,
                formula = f,
                satterthwaite = TRUE,
                progress = FALSE,
                cores = 4,
                save = FALSE)

summary(res)

# Run multiple designs  -----------------------------------------------------
p <- study_parameters(n1 = 10,
                      n2 = 20,
                      n3 = c(2, 4, 6),
                      sigma_subject_intercept = 1.44,
                      icc_pre_cluster = 0,
                      sigma_subject_slope = 0.2,
                      icc_slope = 0.05,
                      sigma_error = 1.44,
                      effect_size = cohend(0.5))

f0 <- "y ~ treatment * time + (1 + time | subject)"
f1 <- "y ~ treatment * time + (1 + time | subject) + (0 + time | cluster)"
f <- sim_formula_compare("correct" = f1,
                         "wrong" = f0)

res <- simulate(object = p,
                nsim = 1000,
                formula = f,
                satterthwaite = TRUE,
                progress = FALSE,
                cores = 1,
                save = FALSE)

# Summarize 'time:treatment' results for n3 = c(2, 4, 6) for 'correct' model
summary(res, para = "time:treatment", model = "correct")

# Summarize cluster-level random slope  for n3 = c(2, 4, 6) for 'correct' model
summary(res, para = "cluster_slope", model = "correct")
# }

Run the code above in your browser using DataLab

Description

Usage

Arguments

Details

See Also

Examples