Package {stepwedgepower}

stepwedgepower: Simulation-Based Power Evaluation for Stepped-Wedge Trials

Description

Provides design objects, binary-outcome data generation, mixed-model analysis, simulation-based power and type I error evaluation, sensitivity grids, and validation tools for stepped-wedge cluster randomized trials. The package supports asynchronous cumulative-intervention schedules and a general component-based engine for Control, A, B, and A+B. The component engine permits arbitrary withdrawal and reintroduction schedules, structurally unobserved cluster-periods, wash-in, user-defined carryover decay, an A-by-B interaction, estimability audits, multiplicity-aware contrasts, failure-aware operating characteristics, and explicit model-fitting outcome diagnostics. A batch-aware interface supports classic and batched stepped-wedge designs with delayed cluster initiation and calendar-time, time-on-trial, or separate batch-specific time effects.

Author(s)

Maintainer: Lin (Amanda) Li amandali14124277@gmail.com

Authors:

• Lin (Amanda) Li amandali14124277@gmail.com

• Florin Vaida (degree: PhD)

See Also

Useful links:

• https://github.com/AmandaLinLi/stepwedgepower

• Report bugs at https://github.com/AmandaLinLi/stepwedgepower/issues

Analyze aggregated grouped outcomes

Description

Generic form of the applied outcome analysis: fits GLM and GLMM rate models for one or more success/trial outcome definitions, across one or more links.

Usage

analyze_aggregated_outcomes(
  data,
  cluster = "cluster_id",
  group = "group",
  outcomes,
  links = c("logit", "identity"),
  nAGQ = 10
)

Arguments

Value

A nested list of fitted models and group-rate tables.

Reproduce the core Lp(a) outcome analyses

Description

Fits the main outcome models from the original script for both overall Lp(a) testing and Lp(a) testing among patients with elevated LDL.

Usage

analyze_lpa_outcomes(
  data,
  provider_var = "prov_id",
  specialty_var = "specialty",
  outcomes = list(
    overall = list(successes = "n_lpa_pat", trials = "n_total_pat"),
    high_ldl = list(successes = "n_ldl_lpa_pat", trials = "n_ldl_pat")
  ),
  links = c("logit", "identity"),
  nAGQ = 10
)

Arguments

Value

A nested list containing fitted models and specialty-rate tables.

Analyze a simulated stepped-wedge trial

Description

Convenience wrapper around [fit_stepwedge_model()] for data produced by [simulate_swcrt()].

Usage

analyze_swcrt(
  trial,
  ...
)

Arguments

Value

The list returned by [fit_stepwedge_model()].

Convert a Cumulative Design to a Component Design

Description

Maps the version 0.3.0 Control/A/A+B representation to the general four-state component representation.

Usage

as_component_design(design)

Arguments

Value

An sw_component_design object.

Examples

old <- sw_multistage_design(c(3, 3), a_start = c(2, 3),
                            b_start = c(4, 5), n_periods = 5)
as_component_design(old)

Audit a Classic or Batched Stepped-Wedge Design

Description

Checks batch initiation, observed-state support, fixed-effect rank, and formal estimability under the selected time model.

Usage

audit_batched_design(
  design,
  assumptions = NULL,
  time_model = c("calendar", "time_on_trial", "separate"),
  contrasts = NULL,
  include_interaction = FALSE,
  adjust_sequence = FALSE
)
## S3 method for class 'sw_batched_audit'
print(x, ...)

Arguments

Value

An object of class sw_batched_audit.

Audit a Component-Based Stepped-Wedge Design

Description

Checks observed treatment-state support, overlap by calendar period, resource counts, withdrawal and restart patterns, fixed-effect rank, and estimability of requested contrasts. Structurally unobserved cells are excluded from rank and overlap calculations but retained for latent exposure history.

Usage

audit_component_design(
  design,
  assumptions = NULL,
  contrasts = NULL,
  include_interaction = TRUE,
  adjust_sequence = FALSE,
  period_effect = c("categorical", "linear", "none")
)

## S3 method for class 'sw_component_audit'
print(x, ...)

Arguments

Value

An object of class sw_component_audit.

Examples

state <- rbind(S1 = c(0, 1, 3, 2), S2 = c(0, 2, 3, 3),
               S3 = c(0, 1, 1, 3), S4 = c(0, 2, 2, 3))
d <- sw_component_design(rep(3, 4), state = state)
audit_component_design(d)

Audit identifiability of an asynchronous cumulative stepped-wedge design

Description

Checks whether a Control -> A -> A+B schedule provides the concurrent comparisons needed to distinguish intervention effects from calendar time. It also flags two common threats to interpretation: adding B in one common calendar period, and adding B after exactly the same duration of A exposure in every sequence. The latter can make an incremental B effect difficult to distinguish from maturation of A when A itself may change with exposure time.

Usage

audit_multistage_design(
  design,
  delay_a = 0L,
  delay_b = 0L,
  adjust_sequence = FALSE
)

## S3 method for class 'sw_multistage_audit'
print(x, ...)

Arguments

Value

An object of class "sw_multistage_audit" containing overlap periods, rank diagnostics, schedule flags, state counts, and explanatory messages.

Examples

design <- sw_multistage_design(
  clusters_per_sequence = rep(4, 4),
  a_start = c(2, 3, 4, 5),
  b_start = c(5, 7, 6, 8),
  n_periods = 8
)
audit_multistage_design(design, delay_a = 1, delay_b = 1)

Convert a random-intercept standard deviation to an intraclass correlation

Description

Inverse of icc_to_cluster_sd.

Usage

cluster_sd_to_icc(cluster_sd)

Arguments

Value

The latent-scale intraclass correlation.

Examples

cluster_sd_to_icc(0.416)

Compare Batched Stepped-Wedge Designs and Time Models

Description

Combines operating characteristics, fit diagnostics, time parameterizations, and resource totals for design or time-model comparisons.

Usage

compare_batched_designs(
  designs, assumptions, time_models = "calendar", labels = NULL,
  nsim = 1000, alpha = 0.05, contrasts = NULL,
  include_interaction = FALSE,
  multiplicity = c("holm", "bonferroni", "none"),
  multiplicity_family = NULL, fit_link = c("logit", "identity"),
  nAGQ = 1, analysis_args = list(), n_cores = 1, seed = NULL,
  check_design = TRUE, warn_on_design = TRUE
)
compare_swd_bswd(
  swd, bswd, assumptions,
  swd_time_model = "calendar", bswd_time_model = "calendar", ...
)
compare_batched_time_models(
  design, assumptions,
  time_models = c("calendar", "time_on_trial", "separate"), ...
)
## S3 method for class 'sw_batched_comparison'
print(x, ...)

Arguments

Value

An object of class sw_batched_comparison.

Compare Component-Based Stepped-Wedge Designs

Description

Runs the component power engine for multiple schedules and reports operating characteristics and resource totals side by side. Calendar and observed cluster-period totals are reported separately so complete and batched stepped-wedge designs can be compared under equal observed resources.

Usage

compare_component_designs(
  designs,
  assumptions,
  labels = NULL,
  nsim = 1000,
  alpha = 0.05,
  contrasts = NULL,
  include_interaction = TRUE,
  multiplicity = c("holm", "bonferroni", "none"),
  multiplicity_family = NULL,
  fit_link = c("logit", "identity"),
  nAGQ = 1,
  analysis_args = list(),
  n_cores = 1,
  seed = NULL,
  check_design = TRUE,
  warn_on_design = TRUE
)

## S3 method for class 'sw_component_comparison'
print(x, ...)

Arguments

Value

An object of class sw_component_comparison. The fit_diagnostics element combines the five fit-category counts and rates across designs.

Examples

state <- rbind(S1 = c(0, 1, 3, 2), S2 = c(0, 2, 3, 3),
               S3 = c(0, 1, 1, 3), S4 = c(0, 2, 2, 3))
d <- sw_component_design(rep(3, 4), state = state)
a <- sw_component_assumptions(baseline_prob = 0.15,
  treatment_or_a = 1.5, treatment_or_b = 1.4,
  interaction_or = 1.1, icc = 0.05)
compare_component_designs(list(first = d, second = d), a,
                          nsim = 10, seed = 1,
                          warn_on_design = FALSE)

Compare two-state and cumulative three-state stepped-wedge designs

Description

Runs power_multistage_swcrt() for two candidate schedules, commonly a Control -> A design and an asynchronous Control -> A -> A+B design. The function reports the operating characteristics side by side and checks whether the designs use the same number of clusters, periods, and cluster-period observations.

Usage

compare_multistage_designs(
  design_two,
  design_three,
  assumptions_two,
  assumptions_three = assumptions_two,
  labels = c("Control -> A", "Control -> A -> A+B"),
  nsim = 1000,
  alpha = 0.05,
  multiplicity = c("holm", "bonferroni", "none"),
  fit_link = c("logit", "identity"),
  nAGQ = 1,
  analysis_args = list(),
  n_cores = 1,
  seed = NULL,
  check_design = TRUE,
  warn_on_design = TRUE
)

## S3 method for class 'sw_multistage_comparison'
print(x, ...)

Arguments

Details

When the same seed is used, both designs start from the same reproducible simulation streams. This can reduce simulation noise in a like-for-like comparison when their dimensions match.

Value

An object of class "sw_multistage_comparison" containing both power objects, a stacked comparison table, a combined fit-diagnostics table, and a resource check.

Examples

three <- sw_multistage_design(
  clusters_per_sequence = rep(4, 4),
  a_start = c(2, 3, 4, 5), b_start = c(5, 7, 6, 8),
  n_periods = 8
)
two <- sw_multistage_design(
  clusters_per_sequence = rep(4, 4),
  a_start = c(2, 3, 4, 5), b_start = rep(Inf, 4),
  n_periods = 8
)
a <- sw_multistage_assumptions(
  baseline_prob = 0.15, treatment_or_a = 1.4,
  incremental_or_b = 1.3, delay_a = 1, delay_b = 1,
  icc = 0.05, n_per_cluster_period = 30
)
compare_multistage_designs(two, three, a, nsim = 100, seed = 1)

Informational comparison with a closed-form power package

Description

Runs power_swcrt alongside the closed-form calculation in the swdpwr package. Offered for information, not as a validation target: closed-form stepped-wedge power is asymptotic in the number of clusters and can be markedly optimistic at the small cluster counts these trials typically use.

Usage

compare_with_swdpwr(
  clusters_per_sequence,
  n_per_cluster_period,
  baseline_prob,
  treatment_or,
  icc,
  alpha = 0.05,
  nsim = 500,
  seed = NULL,
  n_cores = 1,
  adjust_sequence = FALSE
)

Arguments

Value

A one-row data frame with the simulated power, its Monte Carlo interval, the closed-form power, and their difference.

References

Chen J, Zhou X, Li F, Spiegelman D (2022). "swdpwr: A SAS macro and an R package for power calculations in stepped wedge cluster randomized trials." Computer Methods and Programs in Biomedicine, 213, 106522. doi:10.1016/j.cmpb.2021.106522

Zhou X, Liao X, Kunz LM, Normand SLT, Wang M, Spiegelman D (2020). "A maximum likelihood approach to power calculations for stepped wedge designs of binary outcomes." Biostatistics, 21(1), 102-121. doi:10.1093/biostatistics/kxy031

Summarize Resources for a Component-Based Design

Description

Reports calendar and observed sequence-periods and cluster-periods. Structurally unobserved cells are omitted from observed-resource totals.

Usage

component_resource_summary(design, assumptions = NULL)

Arguments

Value

A one-row data frame containing calendar, observed, and missing resource counts and, when available, total individual observations.

Examples

design <- sw_incomplete_component_design(
  c(2, 2, 3, 3),
  rbind(
    c("0", "1", "1+2", "1+2", "-"),
    c("0", "0", "1", "1+2", "-"),
    c("-", "0", "1", "1+2", "1+2"),
    c("-", "0", "0", "1", "1+2")
  )
)
component_resource_summary(design)

Estimate group-specific probabilities from a fitted rate model

Description

Extracts group-level probabilities from a model produced by fit_grouped_rate_model.

Usage

estimate_group_rates(
  model,
  group_levels = NULL,
  group = "group",
  link = c("logit", "identity"),
  approximate_marginal = TRUE,
  logit_scale_factor = 0.346
)

Arguments

Value

A data frame of group-level linear predictors and probabilities.

Estimate power by repeated stepped-wedge simulation

Description

Repeats the stepped-wedge simulation-and-analysis workflow and estimates power.

Usage

estimate_power(
  n_simulations = 100,
  alpha = 0.05,
  treatment_or = NULL,
  n_clusters_per_sequence = NULL,
  sequence_names = NULL,
  cluster_sd = NULL,
  icc = NULL,
  baseline_probs = NULL,
  n_per_cluster_period = NULL,
  n_steps = NULL,
  fit_link = c("logit", "identity"),
  seed = NULL,
  nAGQ = 1,
  effect_size_or = NULL,
  n_providers_per_specialty = NULL,
  specialty_names = NULL,
  base_probs = NULL,
  tau_provider = NULL,
  pts_per_step = NULL
)

Arguments

Value

A list with the estimated power and the vector of p-values.

Estimate specialty-specific probabilities from a fitted model

Description

Extracts specialty-level probabilities from a fitted specialty model.

Usage

estimate_specialty_rates(
  model,
  specialty_levels = NULL,
  specialty_var = "specialty",
  link = c("logit", "identity"),
  approximate_marginal = TRUE,
  logit_scale_factor = 0.346
)

Arguments

Value

A data frame with specialty-level linear predictors and probabilities.

Estimate type I error by repeated stepped-wedge simulation

Description

A convenience wrapper around estimate_power that sets the treatment odds ratio to 1.

Usage

estimate_type1_error(
  n_simulations = 100,
  alpha = 0.05,
  n_clusters_per_sequence = NULL,
  sequence_names = NULL,
  cluster_sd = NULL,
  icc = NULL,
  baseline_probs = NULL,
  n_per_cluster_period = NULL,
  n_steps = NULL,
  fit_link = c("logit", "identity"),
  seed = NULL,
  nAGQ = 1,
  n_providers_per_specialty = NULL,
  specialty_names = NULL,
  base_probs = NULL,
  tau_provider = NULL,
  pts_per_step = NULL
)

Arguments

Value

A list like estimate_power with the type I error estimate added.

Fit a Classic or Batched Stepped-Wedge Analysis Model

Description

Fits the component GLMM with shared categorical calendar effects, shared categorical time-on-trial effects, or a categorical batch-by-local-time factor giving each batch a separate time profile.

Usage

fit_batched_model(
  data,
  time_model = c("calendar", "time_on_trial", "separate"),
  formula = NULL,
  outcome = c("events", "n"),
  contrasts = NULL,
  include_interaction = FALSE,
  multiplicity = c("holm", "bonferroni", "none"),
  multiplicity_family = NULL,
  adjust_sequence = FALSE,
  link = c("logit", "identity"),
  nAGQ = 1,
  alpha = 0.05
)

Arguments

Value

The result from fit_component_model with batched-time metadata.

Fit a Component-Based Stepped-Wedge Model

Description

Fits the planned mixed model and returns A versus Control, B versus Control, B versus A, combined-treatment, incremental, and interaction Wald contrasts with covariance-correct standard errors. If the data contain an observed column, rows marked false or zero are removed before outcome validation and fitting.

Usage

fit_component_model(
  data,
  formula = NULL,
  outcome = c("events", "n"),
  weight_a = "a_effect_weight",
  weight_b = "b_effect_weight",
  weight_ab = "ab_effect_weight",
  term_a = weight_a,
  term_b = weight_b,
  term_ab = weight_ab,
  include_interaction = TRUE,
  contrasts = NULL,
  multiplicity = c("holm", "bonferroni", "none"),
  multiplicity_family = NULL,
  period = "period",
  cluster = "cluster_id",
  sequence = "sequence",
  adjust_sequence = FALSE,
  period_effect = c("categorical", "linear", "none"),
  link = c("logit", "identity"),
  nAGQ = 1,
  alpha = 0.05
)

Arguments

Value

A list containing the fitted model, contrast table, multiplicity results, and fitting diagnostics. The fit_status field distinguishes hard_glmm_error, nonconverged_fit, converged_nonfinite_contrast, and successful_fit. Singularity is returned separately in singular; a singular fit is not automatically classified as unsuccessful. Error, optimizer, convergence-message, and non-finite-contrast details are retained.

Examples

state <- rbind(S1 = c(0, 1, 3, 2), S2 = c(0, 2, 3, 3),
               S3 = c(0, 1, 1, 3), S4 = c(0, 2, 2, 3))
d <- sw_component_design(rep(4, 4), state = state)
a <- sw_component_assumptions(baseline_prob = 0.15,
  treatment_or_a = 1.6, treatment_or_b = 1.5,
  interaction_or = 1.1, icc = 0.05, n_per_cluster_period = 40)
fit_component_model(simulate_component_swcrt(d, a, seed = 1))$contrasts

Fit a grouped rate model

Description

Fits a binomial GLM or cluster-random-intercept GLMM for an aggregated success/trial outcome grouped by an arbitrary grouping variable.

Usage

fit_grouped_rate_model(
  data,
  successes,
  trials,
  group = "group",
  cluster = NULL,
  link = c("logit", "identity"),
  random_intercept = !is.null(cluster),
  nAGQ = 10
)

Arguments

Value

A fitted glm or merMod object.

Fit a cumulative A then A+B stepped-wedge analysis model

Description

Fits a binomial cluster-random-intercept model with separate indicators for the delayed effect of A and the delayed incremental effect of B. The default model is

Usage

fit_multistage_model(
  data,
  formula = NULL,
  outcome = c("events", "n"),
  intervention_a = "a_effective",
  intervention_b = "b_effective",
  term_a = intervention_a,
  term_b = intervention_b,
  period = "period",
  cluster = "cluster_id",
  sequence = "sequence",
  adjust_sequence = FALSE,
  period_effect = c("categorical", "linear", "none"),
  link = c("logit", "identity"),
  nAGQ = 1,
  alpha = 0.05,
  multiplicity = c("holm", "bonferroni", "none")
)

Arguments

Details

cbind(events, n - events) ~ a_effective + b_effective + factor(period) + (1 | cluster_id).

The A coefficient compares A with control after the A wash-in period. The B coefficient compares A+B with A after the B wash-in period. The total A+B effect relative to control is the linear contrast beta_A + beta_B and its standard error includes the covariance between the two estimates.

Value

A list containing the fitted model, a contrast table, raw and adjusted p-values, detailed fit status, convergence and singularity diagnostics, captured warnings, and a joint-success indicator for the A and B component tests. Fit status distinguishes hard GLMM errors, nonconvergence, converged fits with non-finite contrasts, and successful fits; singularity is reported separately.

Examples

design <- sw_multistage_design(
  clusters_per_sequence = rep(5, 4),
  a_start = c(2, 3, 4, 5), b_start = c(5, 7, 6, 8),
  n_periods = 8
)
assumptions <- sw_multistage_assumptions(
  baseline_prob = 0.15, treatment_or_a = 1.6,
  incremental_or_b = 1.4, delay_a = 1, delay_b = 1,
  icc = 0.05, n_per_cluster_period = 40
)
fit <- fit_multistage_model(
  simulate_multistage_swcrt(design, assumptions, seed = 1)
)
fit$contrasts

Fit a specialty-level testing-rate model

Description

Fits either a binomial GLM or a provider-random-intercept binomial GLMM for aggregated success/trial data.

Usage

fit_specialty_rate_model(
  data,
  successes,
  trials,
  specialty_var = "specialty",
  provider_var = NULL,
  link = c("logit", "identity"),
  random_intercept = !is.null(provider_var),
  nAGQ = 10
)

Arguments

Value

A fitted glm or merMod object.

Fit a stepped-wedge analysis model

Description

Fits a binomial random-intercept model to aggregated cluster-period data, configurable through variable names or a user-supplied formula. Independent of any application-specific column names.

Usage

fit_stepwedge_model(
  data,
  formula = NULL,
  outcome = c("events", "n"),
  treatment = "intervention",
  period = "period",
  cluster = "cluster_id",
  sequence = "sequence",
  adjust_sequence = FALSE,
  period_effect = c("categorical", "linear", "none"),
  link = c("logit", "identity"),
  nAGQ = 1
)

Arguments

Value

A list with the fitted model, treatment coefficient row, estimate, standard error, p-value, and convergence/singularity flags.

Examples

d <- sw_design(c(6,6,6,6), crossover_period = c(2,3,4,5), n_periods = 5)
a <- sw_assumptions(baseline_prob = 0.1, treatment_or = 2, icc = 0.03)
fit_stepwedge_model(simulate_swcrt(d, a, seed = 1))$p_value

Convert an intraclass correlation to a random-intercept standard deviation

Description

For a logistic random-intercept model the latent-scale residual variance is \pi^2/3. This maps a latent-scale ICC to the cluster-level random-intercept standard deviation.

Usage

icc_to_cluster_sd(icc)

Arguments

Value

The random-intercept standard deviation.

Examples

icc_to_cluster_sd(0.05)

Power and Type I Error for Batched Stepped-Wedge Designs

Description

Estimates operating characteristics under a selected analysis time model. The generating time model remains separately specified in the assumptions object, allowing matched or misspecified time-model studies.

Usage

power_batched_swcrt(
  design, assumptions,
  analysis_time_model = c("calendar", "time_on_trial", "separate"),
  nsim = 1000, alpha = 0.05, contrasts = NULL,
  include_interaction = FALSE,
  multiplicity = c("holm", "bonferroni", "none"),
  multiplicity_family = NULL, fit_link = c("logit", "identity"),
  nAGQ = 1, analysis_args = list(), n_cores = 1, seed = NULL,
  check_design = TRUE, warn_on_design = TRUE
)

type1_batched_swcrt(
  design, assumptions,
  analysis_time_model = c("calendar", "time_on_trial", "separate"),
  nsim = 1000, alpha = 0.05, contrasts = NULL,
  include_interaction = FALSE,
  multiplicity = c("holm", "bonferroni", "none"),
  multiplicity_family = NULL, fit_link = c("logit", "identity"),
  nAGQ = 1, analysis_args = list(), n_cores = 1, seed = NULL,
  check_design = TRUE, warn_on_design = TRUE
)
## S3 method for class 'sw_batched_power'
print(x, ...)

Arguments

Value

An object of class sw_batched_power, inheriting from sw_component_power.

Estimate power for a component-based stepped-wedge design

Description

Repeatedly simulates a Control/A/B/A+B trial and fits the planned component mixed model. The output distinguishes conditional power among evaluable fits from failure-aware power, which counts unusable analyses as unsuccessful trials. It also reports estimation bias, standard-error calibration, confidence-interval coverage, convergence, singularity, and Monte Carlo uncertainty for every requested contrast. Fit diagnostics separately report hard GLMM errors, nonconverged fits, converged fits with non-finite requested contrasts, singular fits, and successful fits. Singularity is a non-exclusive flag and does not automatically make a fit unsuccessful.

Usage

power_component_swcrt(
  design,
  assumptions,
  nsim = 1000,
  alpha = 0.05,
  contrasts = NULL,
  include_interaction = TRUE,
  multiplicity = c("holm", "bonferroni", "none"),
  multiplicity_family = NULL,
  fit_link = c("logit", "identity"),
  nAGQ = 1,
  analysis_args = list(),
  n_cores = 1,
  seed = NULL,
  check_design = TRUE,
  warn_on_design = TRUE
)

Arguments

Value

An object of class "sw_component_power". In addition to the power and estimation tables, fit_diagnostics gives aggregate counts and rates for the five fit categories, and replicate_diagnostics gives one row per simulation with error and convergence messages.

Examples

state <- rbind(
  S1 = c(0, 1, 1, 3, 3, 2),
  S2 = c(0, 0, 2, 2, 3, 3),
  S3 = c(0, 0, 1, 1, 3, 3),
  S4 = c(0, 0, 2, 3, 3, 3)
)
d <- sw_component_design(rep(4, 4), state = state)
a <- sw_component_assumptions(
  baseline_prob = 0.15, treatment_or_a = 1.5,
  treatment_or_b = 1.4, interaction_or = 1.1,
  carryover_periods_a = 1, carryover_weights_a = 0.5,
  icc = 0.05, n_per_cluster_period = 30
)
p <- power_component_swcrt(d, a, nsim = 20, seed = 1,
                           warn_on_design = FALSE)
p

Power across a grid of design and assumption values

Description

Runs power_swcrt across the Cartesian product of one or more varying parameters, holding everything else fixed. Supports parallel execution, with results reproducible under a given seed regardless of the number of cores.

Usage

power_grid(
  design,
  assumptions,
  vary,
  nsim = 500,
  alpha = 0.05,
  n_cores = 1,
  seed = NULL,
  ...
)

Arguments

Value

An object of class "sw_power_grid": a data frame with one row per grid cell, carrying the full sw_power objects in an attribute.

Estimate power for an asynchronous Control to A to A+B design

Description

Repeatedly simulates a cumulative-intervention stepped-wedge trial and fits the prespecified mixed model. Separate operating characteristics are reported for A versus Control, incremental B (A+B versus A), total A+B versus Control, multiplicity-adjusted component tests, and joint success of both component tests.

Usage

power_multistage_swcrt(
  design,
  assumptions,
  nsim = 1000,
  alpha = 0.05,
  multiplicity = c("holm", "bonferroni", "none"),
  fit_link = c("logit", "identity"),
  nAGQ = 1,
  analysis_args = list(),
  n_cores = 1,
  seed = NULL,
  check_design = TRUE,
  warn_on_design = TRUE
)

Arguments

Details

Two power definitions are reported. Conditional power uses only replicates in which the planned analysis converged and produced the relevant test. Failure-aware power divides the same number of rejections by all simulated trials, treating an unusable analysis as an unsuccessful trial. Their difference reveals how much nominal performance depends on silently discarding failed fits.

Value

An object of class "sw_multistage_power". The fit_diagnostics element reports hard GLMM errors, nonconverged fits, converged fits with non-finite contrasts, singular fits, and successful fits. The replicate_diagnostics element retains one row per simulation.

Examples

design <- sw_multistage_design(
  clusters_per_sequence = rep(5, 5),
  a_start = c(2, 3, 4, 5, 6),
  b_start = c(6, 8, 7, 9, 10),
  n_periods = 10
)
assumptions <- sw_multistage_assumptions(
  baseline_prob = 0.15,
  treatment_or_a = 1.4,
  incremental_or_b = 1.3,
  delay_a = 1, delay_b = 1,
  icc = 0.05,
  n_per_cluster_period = 30
)
result <- power_multistage_swcrt(
  design, assumptions, nsim = 100, seed = 1
)
result

Estimate power for a stepped-wedge trial by simulation

Description

Repeatedly simulates from an [sw_design()] and [sw_assumptions()], fits the analysis model, and estimates rejection probability at level alpha. Returns an [sw_power] object with the Monte Carlo standard error, a confidence interval for power, per-simulation estimates and standard errors, convergence and singular-fit diagnostics, and empirical bias and coverage.

Usage

power_swcrt(
  design,
  assumptions,
  nsim = 1000,
  alpha = 0.05,
  fit_link = c("logit", "identity"),
  nAGQ = 1,
  analysis_args = list(),
  n_cores = 1,
  seed = NULL
)

Arguments

Value

An object of class "sw_power".

Examples

d <- sw_design(c(8,8,8,8), crossover_period = c(2,3,4,5), n_periods = 5)
a <- sw_assumptions(baseline_prob = 0.05, treatment_or = 2, icc = 0.03)
summary(power_swcrt(d, a, nsim = 50, seed = 1))

Prepare grouped cluster-level analysis data

Description

Filters a cluster-level dataset to the groups of interest, keeps clusters within a size window, and sorts the result.

Usage

prepare_cluster_data(
  data,
  groups = NULL,
  min_size = 0,
  max_size = Inf,
  group_var = "group",
  size_var = "n",
  order_var = NULL
)

Arguments

Value

A filtered and sorted data frame.

Prepare physician-level stepped-wedge analysis data

Description

Filters the input data to the specialties of interest, applies panel-size thresholds, removes extreme outliers, and sorts the output.

Usage

prepare_physician_data(
  data,
  specialties = c("CARDIOLOGY", "FAMILY MEDICINE", "INTERNAL MEDICINE", "NEUROLOGY"),
  min_patients = 100,
  max_patients = 10000,
  specialty_var = "specialty",
  patient_var = "n_total_pat",
  provider_name_var = "PROV_NAME"
)

Arguments

Value

A filtered and sorted data frame.

Read the bundled example physician data

Description

Reads a small synthetic physician-level example dataset bundled with the package.

Usage

read_example_physician_data()

Value

A data frame.

Fit the stepped-wedge analysis model to a simulated dataset

Description

Fits the mixed-effects stepped-wedge analysis model used in the original script.

Usage

run_stepwedge_analysis(
  sim_data,
  fit_link = c("logit", "identity"),
  nAGQ = 1
)

Arguments

Value

A list with the fitted model, coefficient table, and treatment p-value.

Simulate a Classic or Batched Stepped-Wedge Trial

Description

Simulates aggregated binary cluster-period outcomes and appends batch, calendar-period, and time-on-trial variables.

Usage

simulate_batched_swcrt(design, assumptions, seed = NULL, include_unobserved = FALSE)

Arguments

Value

A data frame of cluster-period outcomes and treatment-history variables.

Simulate a Component-Based Stepped-Wedge Trial

Description

Generates aggregated binary cluster-period outcomes under arbitrary Control/A/B/A+B schedules with wash-in, withdrawal, carryover, restart, secular trend, component effects, interaction, and structurally unobserved periods.

Usage

simulate_component_swcrt(
  design,
  assumptions,
  seed = NULL,
  include_unobserved = FALSE
)

Arguments

Value

A data frame containing assignments, exposure histories, effect weights, outcomes, and true probabilities.

Examples

state <- rbind(S1 = c(0, 1, 3, 2), S2 = c(0, 2, 3, 3))
d <- sw_component_design(c(3, 3), state = state)
a <- sw_component_assumptions(baseline_prob = 0.15,
  treatment_or_a = 1.4, treatment_or_b = 1.3,
  interaction_or = 1.1, icc = 0.05)
head(simulate_component_swcrt(d, a, seed = 1))

Simulate an asynchronous Control to A to A+B stepped-wedge trial

Description

Generates aggregated cluster-period binary outcomes for a cumulative intervention design. Sequences can start A and add B in different periods, and each intervention can have a prespecified wash-in delay.

Usage

simulate_multistage_swcrt(design, assumptions, seed = NULL)

Arguments

Value

A data frame with one row per cluster-period. In addition to the standard count columns, it contains the active treatment state, exposure times since A and B began, and the delayed-effect indicators used by the default analysis.

Examples

design <- sw_multistage_design(
  clusters_per_sequence = rep(4, 4),
  a_start = c(2, 3, 4, 5),
  b_start = c(5, 7, 6, 8),
  n_periods = 8
)
assumptions <- sw_multistage_assumptions(
  baseline_prob = 0.15, treatment_or_a = 1.4,
  incremental_or_b = 1.3, delay_a = 1, delay_b = 1,
  icc = 0.05, n_per_cluster_period = 30
)
head(simulate_multistage_swcrt(design, assumptions, seed = 1))

Simulate one stepped-wedge trial dataset

Description

Generates aggregated provider-by-step binomial data for a sequential stepped-wedge design.

Usage

simulate_stepwedge_trial(
  treatment_or = NULL,
  n_clusters_per_sequence = NULL,
  sequence_names = NULL,
  baseline_probs = NULL,
  cluster_sd = NULL,
  icc = NULL,
  n_per_cluster_period = NULL,
  n_steps = NULL,
  seed = NULL,
  effect_size_or = NULL,
  n_providers_per_specialty = NULL,
  specialty_names = NULL,
  base_probs = NULL,
  tau_provider = NULL,
  pts_per_step = NULL
)

Arguments

Value

A data frame with one row per provider-step combination.

Simulate one stepped-wedge cluster-randomized trial

Description

Generates an aggregated cluster-by-period dataset from an [sw_design()] and [sw_assumptions()] under a logistic random-intercept model with treatment, period, and baseline effects. Application-agnostic: clusters may be hospitals, clinics, schools, practices, physicians, or regions.

Usage

simulate_swcrt(
  design,
  assumptions,
  seed = NULL
)

Arguments

Value

A data frame with one row per cluster-period (columns cluster_id, sequence, sequence_idx, period, intervention, n, events, true_prob).

Examples

d <- sw_design(c(8,8,8,8), crossover_period = c(2,3,4,5), n_periods = 5)
a <- sw_assumptions(baseline_prob = 0.05, treatment_or = 1.5, icc = 0.05)
head(simulate_swcrt(d, a, seed = 1))

Summarize numeric variables by group

Description

Computes the sample size and summary statistics of one or more numeric variables within each level of a grouping variable.

Usage

summarize_by_group(
  data,
  group_var = "group",
  vars,
  na.rm = TRUE
)

Arguments

Value

A data frame with one row per group-variable combination.

Summarize physician counts by specialty

Description

Computes common summary statistics for one or more numeric variables within each specialty.

Usage

summarize_by_specialty(
  data,
  specialty_var = "specialty",
  vars = c("n_total_pat", "n_ldl_pat"),
  na.rm = TRUE
)

Arguments

Value

A data frame with one row per specialty-variable combination.

Specify the data-generating assumptions for a stepped-wedge trial

Description

Collects baseline risk, treatment effect, cluster random-effect magnitude (via icc or cluster_sd), secular period trend, and cluster-period sample size for simulation. Baseline and sample size may be scalar, per-sequence, or per sequence-period.

Usage

sw_assumptions(
  outcome = "binary",
  baseline_prob = NULL,
  baseline_logit = NULL,
  treatment_or = NULL,
  treatment_effect = NULL,
  icc = NULL,
  cluster_sd = NULL,
  period_effects = NULL,
  n_per_cluster_period = 20
)

Arguments

Value

An object of class "sw_assumptions".

Examples

sw_assumptions(baseline_prob = 0.05, treatment_or = 1.5, icc = 0.05,
               period_effects = c(0,0.05,0.1,0.15,0.2), n_per_cluster_period = 20)

Specify Assumptions for a Batched Stepped-Wedge Trial

Description

Extends component assumptions with one of the three time parameterizations used for batched stepped-wedge designs.

Usage

sw_batched_assumptions(
  ...,
  time_model = c("calendar", "time_on_trial", "separate"),
  time_effects = NULL
)

## S3 method for class 'sw_batched_assumptions'
print(x, ...)

Arguments

Value

An object of class sw_batched_assumptions, inheriting from sw_component_assumptions.

Construct a Classic or Batched Stepped-Wedge Design

Description

Creates a batch-aware component design. A batched stepped-wedge design groups sequences into batches that begin trial participation at staggered calendar periods. Dashes are structurally unobserved cluster-periods and are not Control observations.

Usage

sw_batched_design(
  clusters_per_sequence,
  state,
  batch = NULL,
  latent_state = NULL,
  sequence_names = NULL,
  missing = c("-", "--", "–", "—", "", "NA", ".")
)

as_batched_design(design, batch = NULL)

## S3 method for class 'sw_batched_design'
print(x, ...)

Arguments

Details

All sequences in a batch must share the same first observed calendar period. The returned object records calendar period, time on trial, batch starts, delays from the first batch, gaps from the preceding batch, cluster totals by batch, and the treatment schedule and observation mask.

Value

An object of class sw_batched_design, inheriting from sw_component_design.

Examples

bswd <- sw_batched_design(
  c(2, 2, 3, 3),
  rbind(
    c("0", "1", "1+2", "1+2", "-"),
    c("0", "0", "1", "1+2", "-"),
    c("-", "0", "1", "1+2", "1+2"),
    c("-", "0", "0", "1", "1+2")
  )
)
bswd

Specify Assumptions for a Component-Based Trial

Description

Defines component main effects, interaction, wash-in, withdrawal carryover, restart rules, correlation, secular trend, and sample size for the four-state engine.

Usage

sw_component_assumptions(
  outcome = "binary",
  baseline_prob = NULL,
  baseline_logit = NULL,
  treatment_or_a = NULL,
  treatment_effect_a = NULL,
  treatment_or_b = NULL,
  treatment_effect_b = NULL,
  interaction_or = NULL,
  interaction_effect = NULL,
  delay_a = 0L,
  delay_b = 0L,
  carryover_periods_a = 0L,
  carryover_periods_b = 0L,
  carryover_weights_a = NULL,
  carryover_weights_b = NULL,
  restart_rule_a = c("reset", "resume"),
  restart_rule_b = c("reset", "resume"),
  interaction_mode = c("effective_overlap", "assigned_overlap", "none"),
  icc = NULL,
  cluster_sd = NULL,
  period_effects = NULL,
  n_per_cluster_period = 20
)

## S3 method for class 'sw_component_assumptions'
print(x, ...)

Arguments

Value

An object of class sw_component_assumptions.

Examples

a <- sw_component_assumptions(
  baseline_prob = 0.15,
  treatment_or_a = 1.35,
  treatment_or_b = 1.25,
  interaction_or = 1.10,
  carryover_periods_a = 2,
  carryover_weights_a = c(0.5, 0.25),
  icc = 0.05
)
a

Construct a Component-Based Stepped-Wedge Design

Description

Creates an arbitrary asynchronous schedule containing Control, A, B, and A+B. The schedule may include withdrawal, reintroduction, and structurally unobserved cluster-periods, and need not be monotone. Use sw_incomplete_component_design() for schedules written directly with dashes or missing values.

Usage

sw_component_design(
  clusters_per_sequence,
  state = NULL,
  component_a = NULL,
  component_b = NULL,
  sequence_names = NULL,
  observed = NULL
)

## S3 method for class 'sw_component_design'
print(x, ...)

Arguments

Value

An object of class sw_component_design.

Examples

state <- rbind(
  S1 = c("Control", "A", "A", "A+B", "A+B", "B"),
  S2 = c("Control", "Control", "B", "B", "A+B", "A+B"),
  S3 = c("Control", "Control", "A", "A", "A+B", "A+B")
)
design <- sw_component_design(c(4, 4, 4), state = state)
design

Methods for component-based power objects

Description

Methods for component-based power objects

Usage

## S3 method for class 'sw_component_power'
print(x, ...)

## S3 method for class 'sw_component_power'
summary(object, ...)

## S3 method for class 'sw_component_power'
plot(x, type = c("conditional", "failure-aware"), ...)

Arguments

Construct a stepped-wedge design

Description

Builds a stepped-wedge design object with a 0/1 treatment matrix (sequences x periods), specified either by per-sequence crossover periods or by an explicit treatment matrix. Supports multiple baseline periods, unequal clusters per sequence, and custom or simultaneous crossovers.

Usage

sw_design(
  clusters_per_sequence,
  crossover_period = NULL,
  n_periods = NULL,
  treatment = NULL,
  sequence_names = NULL
)

Arguments

Value

An object of class "sw_design".

Examples

sw_design(clusters_per_sequence = c(10,10,10,10),
         crossover_period = c(2,3,4,5), n_periods = 5)

Construct a Component Design with Structurally Missing Periods

Description

Creates a component-based stepped-wedge design from a manuscript-style schedule in which dashes, blanks, or missing values mean that no outcome is observed. Leading missing cells inherit the first observed state and trailing missing cells inherit the last observed state. The inferred latent state is used only to preserve wash-in, withdrawal, carryover, and restart history; unobserved cells are excluded from simulation outcomes, model fitting, power, and resource totals.

Usage

sw_incomplete_component_design(
  clusters_per_sequence,
  state,
  latent_state = NULL,
  sequence_names = NULL,
  missing = c("-", "--", "–", "—", "", "NA", ".")
)

Arguments

Value

An object of class sw_component_design with an observation mask.

Examples

bswd <- sw_incomplete_component_design(
  c(2, 2, 3, 3),
  rbind(
    c("0", "1", "1+2", "1+2", "-"),
    c("0", "0", "1", "1+2", "-"),
    c("-", "0", "1", "1+2", "1+2"),
    c("-", "0", "0", "1", "1+2")
  )
)
bswd

Specify assumptions for a cumulative A then A+B stepped-wedge trial

Description

Collects the data-generating assumptions for an asynchronous cumulative intervention design created by sw_multistage_design(). The binary-outcome generating model is

Usage

sw_multistage_assumptions(
  outcome = "binary",
  baseline_prob = NULL,
  baseline_logit = NULL,
  treatment_or_a = NULL,
  treatment_effect_a = NULL,
  incremental_or_b = NULL,
  incremental_effect_b = NULL,
  delay_a = 0L,
  delay_b = 0L,
  icc = NULL,
  cluster_sd = NULL,
  period_effects = NULL,
  n_per_cluster_period = 20
)

## S3 method for class 'sw_multistage_assumptions'
print(x, ...)

Arguments

Details

\mathrm{logit}(p_{ijt}) = \beta_{0,jt} + \gamma_t + b_i + \beta_A A^*_{ijt} + \beta_B B^*_{ijt},

where ⁠A*⁠ becomes one after A has completed its prespecified delay and ⁠B*⁠ becomes one after B has completed its delay. Because B is never delivered alone, exp(beta_B) is the incremental effect of adding B among clusters already receiving A. The total A+B effect relative to control is exp(beta_A + beta_B).

A delay of zero means that the effect begins in the first active period. A delay of one treats the first active period as wash-in and begins the effect in the second exposure period.

Value

An object of class "sw_multistage_assumptions".

Examples

assumptions <- sw_multistage_assumptions(
  baseline_prob = 0.15,
  treatment_or_a = 1.40,
  incremental_or_b = 1.30,
  delay_a = 1,
  delay_b = 1,
  icc = 0.05,
  n_per_cluster_period = 30
)
assumptions

Construct an asynchronous cumulative-intervention stepped-wedge design

Description

Builds a design in which clusters move through the ordered treatment states Control -> A -> A+B. Sequences may start intervention A and add intervention B in different calendar periods. B is therefore interpreted as an add-on intervention delivered only after A has begun; the design contains no B-only condition.

Usage

sw_multistage_design(
  clusters_per_sequence,
  a_start = NULL,
  b_start = NULL,
  n_periods = NULL,
  state = NULL,
  sequence_names = NULL
)

## S3 method for class 'sw_multistage_design'
print(x, ...)

Arguments

Details

The design can be specified either by sequence-specific start periods or by an explicit state matrix. In a state matrix, rows are sequences, columns are periods, and entries are 0 (Control), 1 (A), or 2 (A+B). Every row must be non-decreasing and must contain at least one A-only period before A+B.

Value

An object of class "sw_multistage_design".

Examples

design <- sw_multistage_design(
  clusters_per_sequence = rep(6, 5),
  a_start = c(2, 3, 4, 5, 6),
  b_start = c(6, 7, 8, 9, 10),
  n_periods = 10
)
design

Methods for multistage stepped-wedge power objects

Description

Methods for multistage stepped-wedge power objects

Usage

## S3 method for class 'sw_multistage_power'
print(x, ...)

## S3 method for class 'sw_multistage_power'
summary(object, ...)

## S3 method for class 'sw_multistage_power'
plot(x, type = c("conditional", "failure-aware"), ...)

Arguments

Methods for sw_power objects

Description

print, summary, and plot methods for the object returned by power_swcrt. plot shows the distribution of estimated treatment effects with the true effect and mean estimate marked.

Usage

## S3 method for class 'sw_power'
print(x, ...)

## S3 method for class 'sw_power'
summary(object, ...)

## S3 method for class 'sw_power'
plot(x, ...)

Arguments

Value

Invisibly, the sw_power object. Called for their printing / plotting side effects.

Estimate global-null type I error for a component design

Description

Sets the A, B, and interaction effects to zero while preserving the schedule, carryover structure, ICC, secular trend, sample size, and planned analysis.

Usage

type1_component_swcrt(
  design,
  assumptions,
  nsim = 1000,
  alpha = 0.05,
  contrasts = NULL,
  include_interaction = TRUE,
  multiplicity = c("holm", "bonferroni", "none"),
  multiplicity_family = NULL,
  fit_link = c("logit", "identity"),
  nAGQ = 1,
  analysis_args = list(),
  n_cores = 1,
  seed = NULL,
  check_design = TRUE,
  warn_on_design = TRUE
)

Arguments

Value

An sw_component_power object evaluated under the global null.

Examples

state <- rbind(
  S1 = c(0, 1, 1, 3), S2 = c(0, 2, 2, 3),
  S3 = c(0, 0, 1, 3), S4 = c(0, 0, 2, 3)
)
d <- sw_component_design(rep(3, 4), state = state)
a <- sw_component_assumptions(
  baseline_prob = 0.15, treatment_or_a = 1.4,
  treatment_or_b = 1.3, interaction_or = 1.1, icc = 0.05
)
type1_component_swcrt(d, a, nsim = 20, seed = 1,
                      warn_on_design = FALSE)

Estimate global-null type I error for a multistage stepped-wedge design

Description

Convenience wrapper around power_multistage_swcrt() that sets both the A and incremental-B effects to zero while preserving every other assumption.

Usage

type1_multistage_swcrt(
  design,
  assumptions,
  nsim = 1000,
  alpha = 0.05,
  multiplicity = c("holm", "bonferroni", "none"),
  fit_link = c("logit", "identity"),
  nAGQ = 1,
  analysis_args = list(),
  n_cores = 1,
  seed = NULL,
  check_design = TRUE,
  warn_on_design = TRUE
)

Arguments

Value

An sw_multistage_power object evaluated under the global null.

Examples

design <- sw_multistage_design(
  clusters_per_sequence = rep(4, 4),
  a_start = c(2, 3, 4, 5), b_start = c(5, 7, 6, 8),
  n_periods = 8
)
a <- sw_multistage_assumptions(
  baseline_prob = 0.15, treatment_or_a = 1.4,
  incremental_or_b = 1.3, icc = 0.05
)
type1_multistage_swcrt(design, a, nsim = 100, seed = 1)

Validate the simulation engine against an independent reference

Description

Cross-checks power_swcrt against a reference implementation that shares none of its code, and reports whether the two agree within Monte Carlo error. This is a self-contained check and deliberately does not treat any external power package as ground truth.

Usage

validate_engine(
  design,
  assumptions,
  nsim = 500,
  alpha = 0.05,
  seed = NULL,
  tolerance = 0.05
)

Arguments

Value

A one-row data frame comparing the engine's power with the reference's.