sampleSize.RdComputes the required sample size to achieve a target power in studies with multiple endpoints and treatment arms. The function employs modified root-finding algorithms to estimate sample size while accounting for correlation structures, variance assumptions, and equivalence bounds across endpoints. It is particularly useful for bioequivalence trials and multi-arm studies with complex endpoint structures.
sampleSize(
distribution = c("norm", "lnorm", "pois", "nbinom"),
mu_list = NULL,
varcov_list = NA,
sigma_list = NA,
cor_mat = NA,
sigmaB = NA,
rate_list = NULL,
exposure = 1,
dispersion = 0.1,
Eper = c(0, 0),
Eco = c(0, 0),
rho = 0,
TAR = rep(1, length(mu_list)),
arm_names = NA,
ynames_list = NA,
type_y = NA,
list_comparator = NA,
list_y_comparator = NA,
power = 0.8,
alpha = 0.05,
lequi.tol = NA,
uequi.tol = NA,
list_lequi.tol = NA,
list_uequi.tol = NA,
dtype = "parallel",
ctype = "ROM",
vareq = TRUE,
k = NA,
adjust = "no",
dropout = NA,
nsim = 5000,
seed = 1234,
ncores = 1,
optimization_method = "fast",
lower = 2,
upper = 500,
step.power = 6,
step.up = TRUE,
pos.side = FALSE,
maxiter = 1000,
verbose = FALSE,
keep_sim_data = FALSE,
.warn_redundant_bon = TRUE
)Outcome distribution. Choose the R distribution names
"norm", "lnorm", "pois", or "nbinom". Matching is case-insensitive.
The longer labels such as "normal", "lognormal", and "poisson" are
accepted for compatibility; results store the R-style labels above.
Named list of arithmetic means per treatment arm. Each element is a vector representing expected outcomes for all endpoints in that arm.
List of variance-covariance matrices, where each element corresponds to a comparator. Each matrix has dimensions: number of endpoints × number of endpoints.
List of standard deviation vectors, where each element corresponds to a comparator and contains one standard deviation per endpoint.
Matrix specifying the correlation structure between
endpoints. For continuous outcomes it is used with sigma_list to
calculate varcov_list; for count outcomes it is used by the joint count
engine.
Numeric. Between-subject standard deviation parameter for the continuous 2×2 design; for count outcomes it is the log-rate standard deviation in the 2×2 kernel.
Named list of equal-length endpoint-rate vectors, one per count-outcome arm.
Exposure per subject for count outcomes. Supply a scalar or endpoint vector shared by arms, or a named list of arm-specific values.
Positive negative-binomial dispersion parameter. The
per-subject negative-binomial size is 1 / dispersion; parallel-arm
totals use size n / dispersion. Supply a scalar or endpoint vector shared
by arms, or a named list of arm-specific values.
Optional numeric vector of length 2 specifying period effects. For count outcomes these are log-rate effects applied to periods 1 and 2.
Optional numeric vector of length 2 specifying carry-over effects in the order reference carry-over and treatment carry-over. For count outcomes these are log-rate effects in period 2.
Numeric. Correlation parameter applied uniformly across all endpoint pairs. Used with sigma_list to compute varcov_list when cor_mat or varcov_list are not provided.
Numeric vector specifying treatment allocation rates per arm. The order must match arm_names. Defaults to equal allocation across arms if not provided.
Optional character vector of treatment names. If not supplied, names are derived from mu_list.
Optional list of vectors specifying endpoint names per arm. If names are missing, arbitrary names are assigned based on order.
Integer vector indicating endpoint types: 1 for co-primary endpoints, 2 for secondary endpoints.
List of comparators. Each element must be a vector
of length 2 in the form c(test, reference). The first arm is treated as
the treatment arm and the second arm as the reference arm throughout the
simulation and estimand calculations.
List of endpoint sets per comparator. Each element is a vector containing endpoint names to compare. If not provided, all endpoints common to both comparator arms are used.
For count outcomes, the selected endpoints define the count multiplicity
and effective k; joint count analyses require the same endpoint set for
every comparison.
Numeric. Target power (default = 0.8).
Numeric. Significance level (default = 0.05).
Numeric. Lower equivalence bounds (e.g., -0.5) applied uniformly across all endpoints and comparators.
Numeric. Upper equivalence bounds (e.g., 0.5) applied uniformly across all endpoints and comparators.
List of numeric vectors specifying lower equivalence bounds per comparator.
List of numeric vectors specifying upper equivalence bounds per comparator.
Character. Trial design: "parallel" (default) for parallel-group or "2x2" for crossover (only for 2-arm studies).
Character. Continuous-outcome test type: "DOM" (Difference of Means) or "ROM" (Ratio of Means). For Poisson and negative-binomial outcomes, "RR" (event-rate ratio) is used; an unavailable value triggers a warning and the applicable default is used.
Logical. Assumes equal variances across arms if TRUE (default = FALSE).
Integer vector. Minimum number of successful endpoints required for global bioequivalence per comparator. Defaults to all endpoints per comparator.
Character. Alpha adjustment method: "k" (K-fold),
"bon" (Bonferroni across all selected endpoints), "sid"
(Sidak), "t" (Mielke's strong \(k\)-out-of-\(m\) adjustment using
alpha / (m - k + 1); legacy "pc" aliases are accepted),
"no" (none), or
"seq" (sequential).
Numeric vector specifying dropout proportion per arm.
Integer. Number of simulated studies (default = 5000).
Integer. Seed for reproducibility.
Integer. Number of processing cores for parallel computation. Defaults to 1. Set to NA for automatic detection (ncores - 1). For count outcomes, Monte Carlo trials are split into independent seeded chunks and evaluated by the count C++ kernel on each worker.
Character. Sample size optimization method: "fast" (default, root-finding algorithm) or "step-by-step".
Integer. Minimum sample size for search range (default = 2).
Integer. Maximum sample size for the search range (default = 500). For count outcomes, this is the maximum number of subjects per arm; the plotted and returned total sample size is this value multiplied by the number of trial arms.
Numeric. Initial step size for sample size search, defined as 2^step.power. Used when optimization_method = "fast".
Logical. If TRUE (default), search increments upward from lower; if FALSE, decrements downward from upper. Used when optimization_method = "fast".
Logical. If TRUE, finds the smallest integer i closest to the root such that f(i) > 0. Used when optimization_method = "fast".
Integer. Maximum iterations allowed for sample size estimation (default = 1000). Used when optimization_method = "fast".
Logical. If TRUE, prints progress and messages during execution (default = FALSE).
Logical. If TRUE, retain model-scale observations
for every simulated trial in sim_data for distribution diagnostics.
Defaults to FALSE because retained data can be large.
Logical. If TRUE, warn when a requested
multiplicity adjustment is redundant or uncalibrated for the selected
endpoint decision.
A list containing:
responseArray summarizing simulation results, including estimated sample sizes, achieved power, and confidence intervals.
table.iterData frame showing estimated sample sizes and calculated power at each iteration. For count outcomes, one row is retained for every evaluated candidate.
table.testData frame containing test results for all simulated trials. For count outcomes, this contains complete-trial, comparator, and endpoint decision indicators for each simulated trial and candidate; the count kernel returns aggregate decision counts rather than raw endpoint-level test statistics.
param.uOriginal input parameters.
paramFinal adjusted parameters used in sample size calculation.
param.dTrial design parameters used in the simulation.
sim_dataOptional long-format simulated observations,
returned when keep_sim_data = TRUE.
The common planning arguments are power, alpha,
list_comparator, list_lequi.tol, list_uequi.tol,
k, adjust, dtype, dropout, nsim,
seed, lower, and upper. Use the following
distribution-specific arguments in addition to those common arguments:
Supply mu_list, sigma_list
or varcov_list; use cor_mat or rho for endpoint
dependence. The ctype argument selects DOM or ROM testing.
Supply rate_list,
list_comparator, and comparator-specific list_lequi.tol
and list_uequi.tol. Use exposure and, for negative-binomial
outcomes, dispersion; both may be scalar, endpoint-specific, or
named arm-specific lists. Continuous-outcome arguments are ignored.
For count outcomes, optimization_method = "fast" brackets the first
sample size whose simulated power reaches the target and refines the
bracket by integer bisection. The "step-by-step" option remains available
when a complete candidate-by-candidate power table is preferred. The fast
method assumes the usual approximately monotone power curve and uses the
same seed at each candidate to reduce simulation noise.
The effective endpoint count is comparator-specific: when
list_y_comparator is omitted, only endpoints present in both arms are
tested; when it is supplied, only the listed endpoints are tested. k is
validated against that comparator-specific count and oversized values are
capped with a warning. Formal endpoint-wise adjustment is unnecessary when
all selected endpoints are required (k = m), although requested
Bonferroni or Sidak adjustment remains available with a warning. For
k < m, adjust = "no" is explicitly reported as an uncalibrated choice.
For a k-of-m decision, adjust = "t" applies Mielke's strong
\(k\)-out-of-\(m\) calibration alpha / (m - k + 1). The legacy
adjust = "pc" label is accepted as an alias.
For continuous and count outcomes, type_y is used with
adjust = "seq"; named endpoint vectors are aligned to the selected
comparator endpoints. Count analyses use the same primary-gate and
secondary-family decision rule as the continuous kernels.
The unified function returns primary class simss for all outcome
distributions. Count results retain countss as a secondary
compatibility class. Use summary() and plot() to inspect the
result.
Schuirmann, D. J. (1987). A comparison of the Two One-Sided Tests procedure and the Power approach for assessing the equivalence of average bioavailability. Journal of Pharmacokinetics and Biopharmaceutics, 15(6), 657-680. doi:10.1007/BF01068419
Mielke, J., Jones, B., Jilma, B., & König, F. (2018). Sample size for multiple hypothesis testing in biosimilar development. Statistics in Biopharmaceutical Research, 10(1), 39-49. doi:10.1080/19466315.2017.1371071
Berger, R. L., & Hsu, J. C. (1996). Bioequivalence trials, intersection-union tests, and equivalence confidence sets. Statistical Science, 283-302.
Sozu, T., Sugimoto, T., Hamasaki, T., & Evans, S. R. (2015). "Sample Size Determination in Clinical Trials with Multiple Endpoints." SpringerBriefs in Statistics. doi:10.1007/978-3-319-22005-5
mu_list <- list(SB2 = c(AUCinf = 38703, AUClast = 36862, Cmax = 127.0),
EUREF = c(AUCinf = 39360, AUClast = 37022, Cmax = 126.2),
USREF = c(AUCinf = 39270, AUClast = 37368, Cmax = 129.2))
sigma_list <- list(SB2 = c(AUCinf = 11114, AUClast = 9133, Cmax = 16.9),
EUREF = c(AUCinf = 12332, AUClast = 9398, Cmax = 17.9),
USREF = c(AUCinf = 10064, AUClast = 8332, Cmax = 18.8))
# Equivalent boundaries
lequi.tol <- c(AUCinf = 0.8, AUClast = 0.8, Cmax = 0.8)
uequi.tol <- c(AUCinf = 1.25, AUClast = 1.25, Cmax = 1.25)
# Arms to be compared
list_comparator <- list(EMA = c("SB2", "EUREF"),
FDA = c("SB2", "USREF"))
# Endpoints to be compared
list_y_comparator <- list(EMA = c("AUCinf", "Cmax"),
FDA = c("AUClast", "Cmax"))
# Equivalence boundaries for each comparison
lequi_lower <- c(AUCinf = 0.80, AUClast = 0.80, Cmax = 0.80)
lequi_upper <- c(AUCinf = 1.25, AUClast = 1.25, Cmax = 1.25)
# Run the simulation
sampleSize(power = 0.9, alpha = 0.05, mu_list = mu_list,
sigma_list = sigma_list, list_comparator = list_comparator,
list_y_comparator = list_y_comparator,
list_lequi.tol = list("EMA" = lequi_lower, "FDA" = lequi_lower),
list_uequi.tol = list("EMA" = lequi_upper, "FDA" = lequi_upper),
adjust = "no", dtype = "parallel", ctype = "ROM", vareq = FALSE,
distribution = "lnorm", ncores = 1, nsim = 50, seed = 1234)
#> Sample Size Calculation Results
#> -------------------------------------------------------------
#> Study Design: parallel trial targeting 90% power with a 5% type-I error.
#>
#> Comparisons:
#> SB2 vs. EUREF
#> - Endpoints Tested: AUCinf, Cmax
#> (multiple co-primary endpoints, m = 2 )
#> SB2 vs. USREF
#> - Endpoints Tested: AUClast, Cmax
#> (multiple co-primary endpoints, m = 2 )
#> -------------------------------------------------------------
#> Parameter Value
#> Total Sample Size 129
#> Achieved Power 92
#> Power Confidence Interval 79.9 - 97.4
#> -------------------------------------------------------------
# The same entry point for a two-arm Poisson count-rate calculation:
sampleSize(power = 0.80, distribution = "Poisson",
rate_list = list(TEST = 0.21, REF = 0.20),
list_comparator = list(TEST_vs_REF = c("TEST", "REF")),
list_lequi.tol = list(TEST_vs_REF = 0.80),
list_uequi.tol = list(TEST_vs_REF = 1.25),
exposure = 10, lower = 20, upper = 500,
nsim = 100, seed = 1234)
#> Count-rate equivalence sample size
#> Subjects per arm: 186
#> Total subjects: 372
#> Endpoints: 1 (required: 1 )
#> Alpha adjustment: none
#> Achieved power: 0.8100 [0.7167, 0.8789]