PowRPriori (version 0.1.2)

power_sim: Perform a Power Analysis for (Generalized) Linear Mixed-Effects Models via Data Simulation

Description

This is the main function of the PowRPriori package. It iteratively simulates datasets for increasing sample sizes to determine the required sample size to achieve a desired level of statistical power for specific model parameters.

Usage

power_sim(
  formula,
  design,
  test_parameter = NULL,
  fixed_effects,
  random_effects = NULL,
  icc_specs = NULL,
  overall_variance = NULL,
  family = "gaussian",
  power_crit = 0.8,
  n_start,
  n_increment,
  max_simulation_steps = 100,
  n_issue_stop_prop = 0.2,
  n_is_total = TRUE,
  n_sims = 2000,
  alpha = 0.05,
  parallel_plan = "multisession"
)

Value

An object of class PowRPriori, which is a list containing the power table, a sample dataset, all simulation parameters, and detailed results from all runs (coefficients and random effect estimates).

Arguments

formula

An lme4-style model formula (e.g. outcome ~ predictor1 * predictor2 + (1 | id)).

design

A PowRPriori_design object created by define_design().

test_parameter

A character vector of the variable names to test for power. If NULL (default), power is calculated for all fixed effects except the intercept. Note: The parameter names need to comply with the names expected by the model. Correctly naming of the variables is aided by the output of the get_fixed_effects_structure() helper function.

fixed_effects

A named list of the fixed-effects coefficients. It is highly recommended to generate this using get_fixed_effects_structure() or fixed_effects_from_average_outcome().

random_effects

A named, nested list specifying the standard deviations (SDs) and (if applicable) correlations of the random effects. It is highly recommended to generate this using get_random_effects_structure(). If this parameter is not used, icc_specs and overall_variance need to be supplied.

icc_specs

Optional. A named list of Intraclass Correlation Coefficients for defining simple random-intercept models. Must be used with overall_variance.

overall_variance

The total variance of the outcome, required when icc_specs is used.

family

The model family: "gaussian" (for LMMs), "binomial" (for logistic GLMMs), or "poisson" (for poisson GLMMs).

power_crit

The desired statistical power level (e.g., 0.80 for 80%).

n_start

The starting sample size for the simulation.

n_increment

The step size for increasing the sample size in each iteration.

max_simulation_steps

A hard stop for the simulation, limiting the number of sample size steps to prevent infinite loops. Defaults to 100 steps.

n_issue_stop_prop

The proportion of model issues (e.g., singular fits, non-convergence) at which the simulation will be automatically canceled. Defaults to a proportion of 20%.

n_is_total

Boolean that controls how sample sizes are interpreted. If TRUE (default), n_start refers to the total sample size. If FALSE, it refers to the sample size per cell (see define_design() for details on nested designs).

n_sims

The number of simulations to run for each sample size step. Defaults to 2000.

alpha

The significance level (alpha) for the power calculation. Defaults to 0.05.

parallel_plan

A string specifying the future plan for parallel processing. Defaults to "multisession" to enable parallel computing. Use "sequential" for debugging.

Details

The function supports parallel computation using future. Simple linear models (i.e. regression models) can also be analyzed using this function. In this case, no specification of the random_effects or icc_specs parameter is necessary.icc_specs should only be used when simulating a model containing only random intercepts and no random slopes. Refer to the vignette for a more detailed description of the complete workflow for using this function.

Examples

  design <- define_design(
    id = "subject",
    between = list(group = c("Control", "Treatment")),
    within = list(time = c("pre", "post"))
  )

  fixed_effects <- list(
    `(Intercept)` = 10,
    groupTreatment = 2,
    timepost = 1,
    `groupTreatment:timepost` = 1.5
  )

  random_effects <- list(
    subject = list(`(Intercept)` = 3),
    sd_resid = 5
  )
# \donttest{
  power_results <- power_sim(
    formula = y ~ group * time + (1|subject),
    design = design,
    fixed_effects = fixed_effects,
    random_effects = random_effects,
    test_parameter = "groupTreatment:timepost",
    n_start = 20,
    n_increment = 5,
    n_sims = 100, # Use low n_sims for quick examples
    parallel_plan = "multisession"
  )

  summary(power_results)
  plot_sim_model(power_results)
# }