Help for package pmsims

Title:

Simulation-Based Sample Size Tools for Prediction Models

Version:

1.1.0

Description:

Provides a flexible, simulation-based toolkit for exploring how much data are needed to develop reliable prediction models. It works by repeatedly generating data, fitting models, and evaluating performance to show how sample size affects predictive accuracy, calibration, and overfitting. The package supports continuous, binary, and time-to-event outcomes and can be used with both regression-based modelling approaches and machine-learning methods. It is designed to help researchers plan studies, assess feasibility, and build more robust and generalisable models. The methods are described in Olaniran et al. (2026) <doi:10.1186/s12874-026-02935-9> and Shamsutdinova et al. (2026) <doi:10.48550/arXiv.2602.23507>.

License:

GPL (≥ 3)

URL:

https://docs.pmsims.com/, https://github.com/pmsims-package/pmsims

BugReports:

https://github.com/pmsims-package/pmsims/issues

Depends:

R (≥ 4.1.0)

Imports:

cli, ggplot2, lifecycle, mlpwr, parallel, pROC, stats, survival (≥ 3.5-0), utils

Suggests:

covr, DescTools, glmnet (≥ 5.0), knitr, mlbench, mlr, randomForestSRC, ranger, rmarkdown, synthpop, testthat (≥ 3.0.0), timeROC, tuneRanger, xgboost

Encoding:

UTF-8

LazyData:

true

LazyDataCompression:

VignetteBuilder:

knitr

Config/testthat/edition:

Config/Needs/website:

glmmTMB

Config/roxygen2/version:

8.0.0

NeedsCompilation:

Packaged:

2026-11-10 21:24:10 UTC; k1644956

Author:

Ewan Carr

[aut, cre], Gordon Forbes

[aut], Ridwan Olaniran

[aut], Diana Shamsutdinova

[aut], Daniel Stahl

[aut], Sarah Markham

[aut], Felix Zimmer

[aut]

Maintainer:

Ewan Carr <ewan.carr@kcl.ac.uk>

Repository:

CRAN

Date/Publication:

2026-11-10 21:50:02 UTC

pmsims: Simulation-Based Sample Size Tools for Prediction Models

Description

Provides a flexible, simulation-based toolkit for exploring how much data are needed to develop reliable prediction models. It works by repeatedly generating data, fitting models, and evaluating performance to show how sample size affects predictive accuracy, calibration, and overfitting. The package supports continuous, binary, and time-to-event outcomes and can be used with both regression-based modelling approaches and machine-learning methods. It is designed to help researchers plan studies, assess feasibility, and build more robust and generalisable models. The methods are described in Olaniran et al. (2026) doi:10.1186/s12874-026-02935-9 and Shamsutdinova et al. (2026) doi:10.48550/arXiv.2602.23507.

Author(s)

Maintainer: Ewan Carr ewan.carr@kcl.ac.uk (ORCID)

Authors:

Ewan Carr ewan.carr@kcl.ac.uk (ORCID)
Gordon Forbes (ORCID)
Ridwan Olaniran (ORCID)
Diana Shamsutdinova (ORCID)
Daniel Stahl (ORCID)
Sarah Markham (ORCID)
Felix Zimmer (ORCID)

default_model_generators Generate appropriate model based on input arguments

Description

default_model_generators Generate appropriate model based on input arguments

Usage

.xgb_cv_nrounds(
  dtrain,
  params,
  nrounds_max = 500L,
  nfold = 5L,
  early_stopping_rounds = 20L
)

Format

A named list of default model generator functions grouped by outcome type.

Value

default_models is a list containing built-in model generators for binary, continuous, and survival outcomes.

adaptive_startvalues Derive adaptive sample sizes

Description

adaptive_startvalues Derive adaptive sample sizes

Usage

adaptive_startvalues(
  output,
  aggregate_fun,
  var_bootstrap,
  target,
  ci_q = 0.975
)

Arguments

output

List-like object containing track_bisection, produced by calculate_bisection().

aggregate_fun

Function used to summarise replicate performance values (for example, mean or a quantile function).

var_bootstrap

Function returning the bootstrap variance of the aggregated performance.

target

Numeric target performance threshold.

ci_q

Numeric quantile for confidence-interval construction (default 0.975 gives a two-sided 95% interval).

Tuning function for a binary outcome model

Description

Finds the beta_signal and log-odds intercept (mu_lp) that jointly produce the target AUC and prevalence under the exact data-generating settings (complexity, nonlinear_strength, correlation, distribution, predictor_type).

Usage

binary_tuning(
  target_prevalence,
  target_performance,
  candidate_features,
  proportion_noise_features,
  complexity = 1,
  nonlinear_strength = NULL,
  correlation = 0.3,
  distribution = "normal",
  predictor_type = "continuous",
  binary_prevalence = 0,
  n_sim = 3e+05,
  n_validate = NULL,
  beta_interval = c(1e-04, 20),
  beta_tol = 1e-04,
  mu_tol = 1e-06,
  tolerance = 0.02
)

Arguments

target_prevalence

Target outcome prevalence. Must be in (0, 1).

target_performance

Target AUC (c-statistic). Must be in (0.5, 1).

candidate_features

Total number of predictors (signal + noise).

proportion_noise_features

Proportion of candidate_features that are noise predictors. Must be in [0, 1).

complexity

Integer 1-4. Default = 1.

nonlinear_strength

Fraction of signal variance carried by the nonlinear component (C2/C3 only), in [0, 1). When NULL (default), the complexity-level default is used: C1 = 0, C2 = 0.2, C3 = 0.3, C4 = 0.

correlation

Common pairwise predictor correlation. Defaults to the SAME value as generate_binary_data() so the two never silently disagree. Whatever you pass here must match the generator call.

distribution

Global continuous predictor distribution family. Defaults to match generate_binary_data().

predictor_type

"continuous" (default) or "binary". Defaults to match generate_binary_data().

binary_prevalence

Bernoulli probability for binary predictors. Defaults to match generate_binary_data().

n_sim

Sample size for the internal LP simulation. Larger values give more stable AUC estimates within the bisection. Default = 300 000.

n_validate

Sample size for the independent end-to-end validation. Defaults to n_sim.

beta_interval

Search interval c(lo, hi) for beta_signal. Default = c(1e-4, 20). The upper bound is doubled automatically if the AUC at hi is still below target_performance.

beta_tol

Absolute convergence tolerance for the beta bisection. Default = 1e-4.

mu_tol

Absolute convergence tolerance for the mu_lp root-finding. Default = 1e-6.

tolerance

Acceptable absolute deviation of the validated prevalence and AUC from their targets before a warning is raised. Default = 0.02.

Value

A named numeric vector:

mu_lp: Tuned log-odds intercept.
beta_signal: Tuned effect size.
prevalence_achieved: Prevalence from an independent generate_binary_data() validation draw.
auc_achieved: AUC of the oracle LP from the same draw.
var_lp_unit: Estimated Var(LP) at beta_signal = 1.

The Bisection Engine

Description

Runs a bisection search over sample size using repeated simulations and summaries of the chosen performance metric.

Usage

calculate_bisection(
  data_function = data_function,
  model_function = model_function,
  metric_function = metric_function,
  value_on_error = value_on_error,
  min_sample_size = min_sample_size,
  max_sample_size = max_sample_size,
  test_n = test_n,
  n_reps_total = n_reps_total,
  n_reps_per = n_reps_per,
  target_performance = target_performance,
  c_statistic,
  mean_or_assurance = mean_or_assurance,
  tol = 0.001,
  parallel = FALSE,
  cores = 20,
  verbose = FALSE,
  budget = TRUE,
  evaluator = NULL,
  max_n = 1e+06
)

Arguments

data_function

Function taking a single argument, n, giving the training sample size, and returning a dataset that can be passed to model_function.

model_function

Function that fits a model to the dataset returned by data_function. It must take the generated dataset as its only argument and return a fitted model object.

metric_function

Function that evaluates predictive performance on test data. It must take three positional arguments in the order ⁠(test_data, fitted_model, model_name)⁠ and return a single numeric value. Optionally, users may set attr(metric_function, "value_on_error") to a single numeric fallback value to be returned if model fitting or metric evaluation fails during a simulation run.

value_on_error

Numeric fallback returned when a simulation run fails.

min_sample_size

Optional integer lower bound for the sample-size search. If supplied, max_sample_size must also be supplied.

max_sample_size

Optional integer upper bound for the sample-size search. If supplied, min_sample_size must also be supplied. Supplying both bounds defines the search space directly, so the adaptive starting-value search is skipped. Because the runtime estimate is extrapolated from that stage, no long-run warning is issued either.

test_n

Integer size of the test dataset used to evaluate model performance. This should usually be large enough that test-set variability is negligible relative to the training-sample search.

n_reps_total

Integer total number of simulation replications allocated to the search. The search evaluates approximately n_reps_total / n_reps_per candidate sample sizes.

n_reps_per

Integer number of simulation replications performed at each candidate sample size.

target_performance

Numeric target value for the chosen performance metric. The search aims to find the smallest sample size n for which the selected criterion is met relative to this threshold.

c_statistic

Optional numeric value used only by the internal start-value heuristics for some outcome and metric combinations. In most custom workflows this should be left as NULL.

mean_or_assurance

Character string specifying the criterion used to define the minimum sample size. Must be either "mean" or "assurance".

tol

Numeric tolerance controlling when the bisection loop stops.

parallel

Logical; if TRUE the replicates at each sample size run in parallel (forked processes; serial on Windows). Results do not depend on the number of cores.

cores

Integer number of cores to use when parallel = TRUE.

verbose

Logical flag passed to mlpwr; when TRUE verbose output is printed.

budget

Logical; if TRUE the algorithm halts once the evaluation budget is exhausted instead of using tol.

evaluator

Optional evaluator from new_evaluator(), shared by every stage of the search. Created from the data, model and metric functions when NULL.

max_n

Largest sample size the adaptive start-value search may try.

Value

A list containing the simulation results, performance summaries, optional tracking history, and the track_bisection records.

Learning-curve engine

Description

Learning-curve engine

Usage

calculate_curve(
  test_n,
  n_reps_total,
  n_reps_per,
  se_final = NULL,
  min_sample_size,
  max_sample_size,
  target_performance,
  c_statistic,
  mean_or_assurance,
  progress = TRUE,
  data_function,
  model_function,
  metric_function,
  value_on_error,
  evaluator = NULL,
  max_n = 1e+06,
  parallel = FALSE,
  cores = 1L,
  min_n_floor = NULL,
  boot_reps = 200L,
  start_n = NULL,
  live_plot = FALSE,
  ...
)

Arguments

test_n

Integer size of the test dataset used to evaluate model performance. This should usually be large enough that test-set variability is negligible relative to the training-sample search.

n_reps_total

Integer total number of simulation replications allocated to the search. The search evaluates approximately n_reps_total / n_reps_per candidate sample sizes.

n_reps_per

Integer number of simulation replications performed at each candidate sample size.

se_final

Not used: this engine stops when the replicate budget is spent, and reports a bootstrap interval for its answer in diagnostics$curve$n_ci.

min_sample_size

Optional integer lower bound for the sample-size search. If supplied, max_sample_size must also be supplied.

max_sample_size

target_performance

Numeric target value for the chosen performance metric. The search aims to find the smallest sample size n for which the selected criterion is met relative to this threshold.

c_statistic

Optional numeric value used only by the internal start-value heuristics for some outcome and metric combinations. In most custom workflows this should be left as NULL.

mean_or_assurance

Character string specifying the criterion used to define the minimum sample size. Must be either "mean" or "assurance".

progress

Logical; show a progress bar over the replicate budget.

data_function

Function taking a single argument, n, giving the training sample size, and returning a dataset that can be passed to model_function.

model_function

Function that fits a model to the dataset returned by data_function. It must take the generated dataset as its only argument and return a fitted model object.

metric_function

value_on_error

Numeric value recorded for a replicate whose fit or metric fails.

evaluator

Optional evaluator from new_evaluator(), shared by every stage of the search. Created from the data, model and metric functions when NULL.

max_n

Largest sample size the search may try when no min_sample_size and max_sample_size are given.

parallel, cores

Run the replicates of each batch in parallel with parallel::mclapply(); used only when evaluator is NULL.

min_n_floor

Smallest sample size the search may evaluate. Defaults to max(10, p + 5) for p predictors.

boot_reps

Number of bootstrap refits for the interval of the answer.

start_n

Optional first sample size of the pilot; defaults to the heuristic start value.

live_plot

Redraw the learning curve after each batch of replicates (interactive sessions with a graphics device only).

...

Unused.

Calculate performance metrics Calculates the performance metrics for a model given a sample size n.

Description

Calculate performance metrics Calculates the performance metrics for a model given a sample size n.

Usage

calculate_metrics_perf(
  n,
  data_function,
  model_function,
  metric_function,
  value_on_error
)

Arguments

n

Integer sample size.

value_on_error

Numeric fallback returned if the metric cannot be computed.

Value

The calculated performance metric

mlpwr engine

Description

mlpwr engine

Usage

calculate_mlpwr(
  test_n,
  n_reps_total,
  n_reps_per,
  se_final,
  min_sample_size,
  max_sample_size,
  target_performance,
  c_statistic,
  mean_or_assurance,
  n_init,
  progress = TRUE,
  verbose,
  data_function,
  model_function,
  metric_function,
  value_on_error,
  evaluator = NULL,
  max_n = 1e+06,
  parallel = FALSE,
  cores = 1L,
  ...
)

Arguments

test_n

Integer size of the test dataset used to evaluate model performance. This should usually be large enough that test-set variability is negligible relative to the training-sample search.

n_reps_total

Integer total number of simulation replications allocated to the search. The search evaluates approximately n_reps_total / n_reps_per candidate sample sizes.

n_reps_per

Integer number of simulation replications performed at each candidate sample size.

min_sample_size

Optional integer lower bound for the sample-size search. If supplied, max_sample_size must also be supplied.

max_sample_size

target_performance

Numeric target value for the chosen performance metric. The search aims to find the smallest sample size n for which the selected criterion is met relative to this threshold.

c_statistic

Optional numeric value used only by the internal start-value heuristics for some outcome and metric combinations. In most custom workflows this should be left as NULL.

mean_or_assurance

Character string specifying the criterion used to define the minimum sample size. Must be either "mean" or "assurance".

n_init

Integer number of initial sample sizes simulated before the Gaussian process search begins.

progress

Logical flag controlling whether the mlpwr progress bar is shown.

verbose

Logical flag passed to mlpwr; when TRUE verbose output is printed.

data_function

Function taking a single argument, n, giving the training sample size, and returning a dataset that can be passed to model_function.

model_function

Function that fits a model to the dataset returned by data_function. It must take the generated dataset as its only argument and return a fitted model object.

metric_function

value_on_error

Numeric fallback value used if model fitting or metric calculation fails.

evaluator

Optional evaluator from new_evaluator(), shared by every stage of the search. Created from the data, model and metric functions when NULL.

max_n

Largest sample size the adaptive start-value search may try.

parallel, cores

Run the replicates of each batch in parallel (see calculate_bisection()); with mlpwr, replicates are requested one at a time, so this mainly affects the adaptive stage and verification.

...

Additional options passed to mlpwr::find.design().

mlpwr-bs Hybrid engine using bisection to determine initial range and mlpwr for search

Description

mlpwr-bs Hybrid engine using bisection to determine initial range and mlpwr for search

Usage

calculate_mlpwr_bs(
  test_n,
  n_reps_total,
  n_reps_per,
  se_final,
  min_sample_size,
  max_sample_size,
  target_performance,
  c_statistic,
  mean_or_assurance,
  progress = TRUE,
  verbose,
  data_function,
  model_function,
  metric_function,
  value_on_error,
  evaluator = NULL,
  max_n = 1e+06,
  parallel = FALSE,
  cores = 1L,
  ...
)

Arguments

test_n

Integer size of the test dataset used to evaluate model performance. This should usually be large enough that test-set variability is negligible relative to the training-sample search.

n_reps_total

Integer total number of simulation replications allocated to the search. The search evaluates approximately n_reps_total / n_reps_per candidate sample sizes.

n_reps_per

Integer number of simulation replications performed at each candidate sample size.

min_sample_size

Optional integer lower bound for the sample-size search. If supplied, max_sample_size must also be supplied.

max_sample_size

target_performance

Numeric target value for the chosen performance metric. The search aims to find the smallest sample size n for which the selected criterion is met relative to this threshold.

c_statistic

Optional numeric value used only by the internal start-value heuristics for some outcome and metric combinations. In most custom workflows this should be left as NULL.

mean_or_assurance

Character string specifying the criterion used to define the minimum sample size. Must be either "mean" or "assurance".

progress

Logical flag controlling whether the mlpwr progress bar is shown.

verbose

Logical flag passed to mlpwr; when TRUE verbose output is printed.

data_function

Function taking a single argument, n, giving the training sample size, and returning a dataset that can be passed to model_function.

model_function

Function that fits a model to the dataset returned by data_function. It must take the generated dataset as its only argument and return a fitted model object.

metric_function

value_on_error

Numeric fallback value used if model fitting or metric calculation fails.

evaluator

Optional evaluator from new_evaluator(), shared by every stage of the search. Created from the data, model and metric functions when NULL.

max_n

Largest sample size the adaptive start-value search may try.

parallel, cores

Run the replicates of each batch in parallel (see calculate_bisection()); with mlpwr, replicates are requested one at a time, so this mainly affects the adaptive stage and verification.

...

Additional options passed to mlpwr::find.design().

Value

List containing the combined bisection and mlpwr results (results, summaries, min_n, perf_n, and mlpwr_ds).

Tuning function for a continuous outcome model

Description

Finds the beta_signal value that produces a target large-sample R^2 under the exact data-generating settings (complexity, nonlinear_strength, correlation, distribution, predictor_type).

Usage

continuous_tuning(
  r2,
  candidate_features,
  proportion_noise_features,
  complexity = 1,
  nonlinear_strength = NULL,
  correlation = 0,
  distribution = "normal",
  predictor_type = "continuous",
  binary_prevalence = 0,
  n_sim = 1e+05
)

Arguments

r2

Target large-sample R^2 (proportion of variance explained). Must be in (0, 1).

candidate_features

Total number of predictors (signal + noise).

proportion_noise_features

Proportion of candidate_features that are noise predictors (zero coefficient). Must be in [0, 1).

complexity

Integer 1-4 controlling the functional form of the linear predictor passed to the data generator. Default = 1.

nonlinear_strength

Fraction of signal variance carried by the nonlinear component (C2/C3 only), in [0, 1). When NULL (default), the complexity-level default is used: C1 = 0, C2 = 0.2, C3 = 0.3, C4 = 0.

correlation

Common pairwise predictor correlation. Default = 0.

distribution

Global continuous predictor distribution family passed to the data generator. Default = "normal".

predictor_type

"continuous" (default) or "binary".

binary_prevalence

Bernoulli probability for binary predictors. Required (and used) only when predictor_type = "binary".

n_sim

Sample size used to estimate Var(LP). Larger values give more stable estimates. Default = 100 000.

Value

A named numeric vector:

beta_signal: Tuned effect size.
r2_achieved: Empirical R^2 verified in a large simulation.
var_lp_unit: Estimated Var(LP) at beta_signal = 1 (the scaling constant).

Create default data generating functions

Description

Create default data generating functions

Usage

default_data_generators(opts)

Arguments

opts

A list with two elements:

type: Outcome type: "continuous", "binary", or "survival".
args: Named list of arguments to pre-set on the corresponding generator function.

Value

A partially-applied generator function whose formals have been set to the values in opts$args.

Simulate binary outcome data

Description

Simulate binary outcome data

Usage

generate_binary_data(
  n,
  n_signal_parameters,
  noise_parameters,
  beta_signal,
  complexity = 1,
  nonlinear_strength = NULL,
  predictor_type = "continuous",
  binary_prevalence = 0,
  correlation = 0.3,
  distribution = "normal",
  mu_lp = 0,
  baseline_prob = 0.5
)

Arguments

n

Sample size.

n_signal_parameters

Number of signal predictors. These occupy the first n_signal_parameters columns (x1 ... x_S).

noise_parameters

Number of noise predictors (zero coefficient).

beta_signal

Base effect size / overall scale of the signal.

complexity

Integer 1-4 specifying the functional form of the linear predictor:

Linear — lp = \alpha + \beta\sum_j x_j.
Quadratic — linear + a quadratic nonlinear component.
Quadratic + Interaction — linear + quadratic and pairwise interaction nonlinear component.
Friedman — canonical Friedman (1991) benchmark.

nonlinear_strength

Fraction of signal variance carried by the nonlinear component (C2/C3 only), in [0, 1). When NULL (default), the complexity-level default is used: C1 = 0, C2 = 0.2, C3 = 0.3, C4 = 0. Ignored for C1 (pure linear) and C4 (canonical Friedman).

predictor_type

Type of predictors: "continuous" (default) or "binary". When "binary", all predictors are drawn as Bernoulli(binary_prevalence); binary_prevalence must be in (0, 1] and distribution is ignored.

binary_prevalence

Scalar in (0, 1]. Bernoulli probability applied to all predictors when predictor_type = "binary". Default = 0.

correlation

Scalar in [-1, 1]. Common pairwise correlation applied via a Gaussian copula: correlation is the correlation of the underlying normal variables, which equals the Pearson correlation of the predictors only when they are normal. Default = 0.3. Set to 0 for independence.

distribution

Distribution family for all continuous predictors. Default = "normal". For complexity 4, if left at "normal" the framework uses "uniform" (Friedman canonical).

mu_lp

Intercept on the log-odds scale. Default = 0.

baseline_prob

Nominal baseline event probability (documentation only; the realised probability is determined by mu_lp).

Value

A data frame with columns y (0/1), x1, x2, ...

Simulate continuous outcome data

Description

Simulate continuous outcome data

Usage

generate_continuous_data(
  n,
  n_signal_parameters,
  noise_parameters,
  beta_signal,
  complexity = 1,
  nonlinear_strength = NULL,
  predictor_type = "continuous",
  binary_prevalence = 0,
  correlation = 0.3,
  distribution = "normal",
  intercept = 0
)

Arguments

n

Sample size.

n_signal_parameters

Number of signal predictors. These occupy the first n_signal_parameters columns (x1 ... x_S).

noise_parameters

Number of noise predictors (zero coefficient).

beta_signal

Base effect size / overall scale of the signal.

complexity

Integer 1-4 specifying the functional form of the linear predictor:

Linear — lp = \alpha + \beta\sum_j x_j.
Quadratic — linear + a quadratic nonlinear component.
Quadratic + Interaction — linear + quadratic and pairwise interaction nonlinear component.
Friedman — canonical Friedman (1991) benchmark.

nonlinear_strength

predictor_type

binary_prevalence

Scalar in (0, 1]. Bernoulli probability applied to all predictors when predictor_type = "binary". Default = 0.

correlation

distribution

Distribution family for all continuous predictors. Default = "normal". For complexity 4, if left at "normal" the framework uses "uniform" (Friedman canonical).

intercept

Scalar intercept added to the linear predictor. Default = 0.

Value

A data frame with columns y, x1, x2, ...

References

Friedman, J. H. (1991). Multivariate adaptive regression splines. The Annals of Statistics, 19(1), 1-67. doi:10.1214/aos/1176347963

Construct the linear predictor

Description

Construct the linear predictor

Usage

generate_linear_predictor(
  X,
  n_signal_parameters,
  noise_parameters,
  intercept,
  beta_signal,
  complexity,
  nonlinear_strength = NULL,
  predictor_type,
  binary_prevalence,
  correlation,
  distribution
)

Arguments

X

n x p predictor matrix (colnames: x1, x2, ...).

n_signal_parameters

Signal predictor count.

noise_parameters

Noise predictor count.

intercept

Scalar intercept (log-odds intercept for binary, log-hazard intercept for survival, mean intercept for continuous).

beta_signal

Base effect size (overall scale). The tuning routines scale this to hit the target oracle metric.

complexity

Integer 1-4.

nonlinear_strength

Fraction of signal variance carried by the nonlinear component (C2/C3 only), in [0, 1). NULL uses the complexity-level default in COMPLEXITY_NONLINEAR_STRENGTH_DEFAULTS.

predictor_type, binary_prevalence, correlation, distribution

The settings X was drawn with (see generate_predictors). Required for complexity 2 and 3, where they select the reference estimates for the nonlinear component; they must match how X was drawn, or the component is no longer orthogonal to the linear terms.

Details

Why variance fraction rather than an R^2 split. For a continuous outcome whose target is R^2, splitting the signal by latent-scale variance is exact: R^2 is a variance ratio, so a correctly specified linear model recovers exactly the linear share. For binary (AUC) and survival (Harrell's C) outcomes that logic is only approximate, because AUC and C are rank/link-mediated. What carries over is the key property: the nonlinear component is built to be inaccessible to any linear-in-X model, so the discrimination it carries genuinely requires a nonlinear learner and more sample size.

Construction. The linear predictor is

lp = \alpha + \beta\,\Bigl(\,\underbrace{\textstyle\sum_j x_j}_{L} \; + \; \kappa\,\mathrm{sd}(L)\, \underbrace{\tilde N}_{N\,\mathrm{std}} \Bigr)

where L = \sum_j x_j is the linear score (identical to complexity 1) and \tilde N is the nonlinear aggregate (\sum_j x_j^2 for C2; plus \sum_{j<k} x_j x_k for C3) residualised against the full design [1, x_1, \ldots, x_S] and standardised to unit SD. The projection coefficients and both SDs are estimated once per predictor configuration from a 500,000-row reference draw (see nonlinear_reference()), so every row is generated by the same model whatever the number of rows drawn. Because L and \tilde N are orthogonal in the population (up to the error in these estimates), the nonlinear variance fraction is

f = \frac{\kappa^2}{1 + \kappa^2},\qquad\text{equivalently}\qquad \kappa = \sqrt{\frac{f}{1 - f}}.

The user supplies nonlinear_strength = f (e.g. 0.2 = 20% of the signal variance is nonlinear) and \kappa is computed internally. The population fraction is f up to the error in the reference draw: within 0.002 for normal, uniform, Laplace and exponential predictors, but only within about 0.02 for t and lognormal predictors, whose squares are heavy-tailed.

Because the bracket is fixed given X, lp is linear in beta_signal, so the tuning routines (which scale beta_signal to hit the target oracle AUC / C-index) work unchanged.

Calibrating difficulty. f sets the difficulty order but not, by itself, a precise AUC/C-index gap: the metric reduction a given f produces also depends on the tuned scale, the prevalence/censoring, and the predictor family. To target a specific linear-accessible metric, calibrate f per complexity by simulation (raise it until a fitted GLM/Cox on the linear terms reaches the desired value while beta_signal keeps the oracle metric on target).

Complexity 1 — Linear: lp = \alpha + \beta \sum_j x_j.

Complexity 4 — Friedman (1991) (canonical, nonlinear_strength ignored): lp = \alpha + \beta[10\sin(\pi x_1 x_2) + 20(x_3-0.5)^2 + 10x_4 + 5x_5]. Friedman contains linear-accessible terms (10x_4 + 5x_5); to place C4 on the same controlled ladder as C2/C3, residualise the Friedman vector against Xs and apply the same kappa split.

Value

Numeric vector of length n.

References

Friedman, J. H. (1991). Multivariate adaptive regression splines. The Annals of Statistics, 19(1), 1-67.

Generate the n x p predictor matrix

Description

Draws all predictors from a single global distribution family. With a non-zero correlation, they are drawn through a Gaussian copula with a common pairwise correlation.

Usage

generate_predictors(
  n,
  n_signal_parameters,
  noise_parameters,
  complexity = 1,
  predictor_type = "continuous",
  binary_prevalence = 0,
  correlation = 0.3,
  distribution = "normal"
)

Arguments

n

Sample size.

n_signal_parameters

Number of signal predictors.

noise_parameters

Number of noise predictors.

complexity

Integer 1-4 (used to resolve the C4 distribution default).

predictor_type

"continuous" (default) or "binary".

binary_prevalence

Bernoulli probability; used when predictor_type = "binary".

correlation

Scalar common pairwise correlation; 0 = independent.

distribution

Global continuous distribution family; used when predictor_type = "continuous".

Value

Named n x p numeric matrix (column names: x1, x2, ...).

Simulate survival outcome data

Description

Simulate survival outcome data

Usage

generate_survival_data(
  n,
  n_signal_parameters,
  noise_parameters,
  beta_signal,
  baseline_hazard,
  censoring_rate,
  complexity = 1,
  nonlinear_strength = NULL,
  predictor_type = "continuous",
  binary_prevalence = 0,
  correlation = 0.3,
  distribution = "normal",
  intercept = 0
)

Arguments

n

Sample size.

n_signal_parameters

Number of signal predictors. These occupy the first n_signal_parameters columns (x1 ... x_S).

noise_parameters

Number of noise predictors (zero coefficient).

beta_signal

Base effect size / overall scale of the signal.

baseline_hazard

Baseline hazard rate (exponential survival model).

censoring_rate

Administrative censoring proportion in (0, 1).

complexity

Integer 1-4 specifying the functional form of the linear predictor:

Linear — lp = \alpha + \beta\sum_j x_j.
Quadratic — linear + a quadratic nonlinear component.
Quadratic + Interaction — linear + quadratic and pairwise interaction nonlinear component.
Friedman — canonical Friedman (1991) benchmark.

nonlinear_strength

predictor_type

binary_prevalence

Scalar in (0, 1]. Bernoulli probability applied to all predictors when predictor_type = "binary". Default = 0.

correlation

distribution

Distribution family for all continuous predictors. Default = "normal". For complexity 4, if left at "normal" the framework uses "uniform" (Friedman canonical).

intercept

Scalar intercept added to the linear predictor. Default = 0.

Value

A data frame with columns time, event (0 = censored, 1 = event), x1, x2, ...

get_min_sample_size: Heuristic starting-n for binary/continuous/survival prediction

Description

get_min_sample_size: Heuristic starting-n for binary/continuous/survival prediction

Usage

get_min_sample_size(
  npar,
  prevalence = NULL,
  c_stat = NULL,
  calibration_slope = NULL,
  epv_value = NULL,
  outcome_type = c("binary", "survival", "continuous")
)

Arguments

npar

Integer; number of predictors in the model.

prevalence

Numeric in ⁠[0, 1]⁠; optional event rate or case fraction used for EPV calculations.

c_stat

Numeric in (0.5, 1]; anticipated discrimination (C-statistic). Lower values inflate the heuristic.

calibration_slope

Numeric; anticipated calibration slope. Values below 1 trigger a modest inflation.

epv_value

Numeric; target events-per-variable (EPV) value applied when prevalence is supplied.

outcome_type

Character string; must be one of "binary", "survival", or "continuous".

Value

Integer recommended starting value from which to calculate the minimum sample size.

Objective function for pmsims

Description

Computes the objective value

\left| M(n) - M^* - \lambda\,C(n) \right|,

where M(n) is performance at sample size n, M^* is the target performance, C(n) is a cost term increasing with n, and \lambda (set via penalty_weight) controls the penalty on large samples.

Usage

objective_function(
  n,
  penalty_weight,
  target_performance,
  min_sample_size,
  max_sample_size,
  value_on_error,
  data_function = NULL,
  model_function = NULL,
  metric_function = NULL
)

Arguments

n

Integer sample size.

penalty_weight

Penalty weight \lambda on the cost term C(n); 0 minimises only the absolute performance gap.

target_performance

Target performance M^* (numeric).

min_sample_size

Minimum sample size considered (lower bound of search).

max_sample_size

Maximum sample size considered (upper bound of search).

value_on_error

Value to return if the objective cannot be evaluated.

data_function

Data generator passed to calculate_metrics_perf().

model_function

Model-fitting function passed to calculate_metrics_perf().

metric_function

Metric function passed to calculate_metrics_perf().

Value

A single numeric: the objective value at n.

Examples

# Example usage (assuming helper functions exist):
# objective_function(
#   n = 500,
#   penalty_weight = 0.1,
#   target_performance = 0.75,
#   min_sample_size = 100,
#   max_sample_size = 5000,
#   value_on_error = Inf
# )

Parse and validate input specifications

Description

This function validates the provided data, model, and metric specifications, and returns corresponding generator functions for each. It ensures that all required inputs are provided and correctly configured.

Usage

parse_inputs(data_spec, metric, model)

Arguments

data_spec

A list containing two elements:

type: A character string indicating the outcome type.
args: A list of arguments to be passed to the data-generating function.

metric

A character vector specifying one or more metrics to be used. Currently, only the first element is used.

model

A character string specifying the model to be used.

Details

This function calls default_data_generators(), default_model_generators(), and default_metric_generator() to construct the appropriate functions based on the supplied inputs.

Value

A list containing three elements:

data_function: The data-generating function.
model_function: The model-generating function.
metric_function: The metric function corresponding to the chosen metric.

Plot sample-size learning curves for `pmsims` outputs

Description

Produces a ggplot showing the simulated points and the fitted learning curve (the learning-curve fit for method = "curve", the Gaussian-process surrogate for the mlpwr engines) stored inside a pmsims object. Optionally returns the underlying data instead of drawing the plot.

Usage

## S3 method for class 'pmsims'
plot(x, metric_label = NULL, plot = TRUE, ...)

Arguments

x

A pmsims object returned by simulate_binary(), simulate_continuous(), simulate_survival(), or simulate_custom().

metric_label

Optional string used for the y-axis label when the object does not already record the metric name.

plot

Logical; if TRUE (default) the function prints the plot. If FALSE, the data used to build the plot are returned instead of drawing anything.

...

Currently unused.

Value

Invisibly returns the ggplot object when plot = TRUE. When plot = FALSE, returns a list with two data frames: for the mlpwr engines, observed_data (simulated points) and predicted_data (Gaussian-process predictions); for method = "curve", observed_data (the criterion at each simulated sample size, with replicates, standard error and a +/- 2 SE interval) and fitted_curve (the fitted learning curve, or NULL when none was fitted).

Precomputed sample-size simulations

Description

Results from four completed sample-size searches, one for each of the package's entry points. A search repeatedly simulates datasets, fits a prediction model to each and evaluates its performance across a range of candidate sample sizes, so a realistic run takes minutes to hours. These objects were computed once, with the calls shown under Source, so that the examples and the vignette can demonstrate the output without recomputing it.

Usage

binary_example

continuous_example

survival_example

custom_example

Format

A list of class "pmsims". The most useful components are:

min_n: The estimated minimum sample size.
perf_n: Expected performance at min_n.
target_performance: The performance target the search aimed for.
metric: The performance metric used.
model: The model fitted to each simulated dataset.
status: "ok" when the search found min_n and a check at min_n confirmed it; see simulate_custom() for other values.
data: The simulated performance values at each sample size, used by plot().
diagnostics: Details of the search, including the fitted learning curve and the interval for min_n (diagnostics$curve).
summaries: Aggregated performance summaries across replications.
simulation_time: Elapsed time of the original run, in seconds.

Objects produced by the wrapper functions additionally record the data-generating configuration (for example outcome_prevalence, correlation and complexity); see simulate_binary() for the full set.

An object of class pmsims of length 43.

An object of class pmsims of length 42.

An object of class pmsims of length 44.

An object of class pmsims of length 31.

Details

Each object is a "pmsims" object, as returned by simulate_binary(), simulate_continuous(), simulate_survival() and simulate_custom(), and can be used with print() and plot() in the usual way.

Source

Generated by data-raw/precomputed-examples.R.

binary_example:

set.seed(123)
simulate_binary(
  signal_parameters = 20,
  noise_parameters = 0,
  complexity = 1,
  data_control = list(correlation = 0.3),
  outcome_prevalence = 0.30,
  maximum_achievable_cstatistic = 0.80,
  model = "glm",
  metric = "calibration_slope",
  target_performance = 0.85,
  n_reps_total = 1000,
  mean_or_assurance = "assurance"
)

continuous_example:

set.seed(123)
simulate_continuous(
  signal_parameters = 15,
  noise_parameters = 0,
  complexity = 1,
  data_control = list(correlation = 0.3),
  maximum_achievable_rsquared = 0.50,
  model = "lm",
  metric = "calibration_slope",
  target_performance = 0.95,
  n_reps_total = 1000,
  mean_or_assurance = "assurance"
)

survival_example:

set.seed(123)
simulate_survival(
  signal_parameters = 15,
  noise_parameters = 0,
  complexity = 1,
  data_control = list(correlation = 0.3),
  maximum_achievable_cindex = 0.70,
  baseline_hazard = 0.01,
  censoring_rate = 0.30,
  model = "coxph",
  metric = "calibration_slope",
  target_performance = 0.90,
  n_reps_total = 1000,
  mean_or_assurance = "assurance"
)

custom_example: a linear model with three independent predictors and a population R^2 of 0.5, targeting a calibration slope of 0.9 with 80% assurance. The search uses bounds of 25 to 1,000 participants, 1,000 replications and 30,000 test observations. See simulate_custom() for the generating functions and executable call (seed 123).

Examples

binary_example
binary_example$min_n

Print a `pmsims` result

Description

Print a pmsims result

Usage

## S3 method for class 'pmsims'
print(x, ..., max_width = 80, verbose = FALSE)

Arguments

x

A pmsims object.

...

Currently unused.

max_width

Maximum console width used for the summary.

verbose

Logical. If TRUE, add implementation-level detail: internal metric identifiers, the engine settings used for the search, and any quantities recorded on an internal search scale. summary.pmsims() is the same display with verbose = TRUE by default.

Value

x, invisibly.

Minimum sample size for binary-outcome prediction models

Description

Compute the minimum sample size required to develop a prediction model with a binary outcome. The function wraps a simulation-based search that, by default, fits a learning curve to every simulated replicate and simulates further where the curve crosses the target (method = "curve"; see simulate_custom()); method = "mlpwr" selects the Gaussian-process search of earlier versions. From user inputs (outcome prevalence, maximum achievable performance, target performance, etc.) it constructs a data-generating function, a model-fitting function, and a metric function, then searches for the smallest n that meets the chosen performance criterion.

Usage

simulate_binary(
  signal_parameters,
  noise_parameters = 0,
  complexity = 1,
  data_control = NULL,
  outcome_prevalence,
  maximum_achievable_cstatistic,
  model = c("glm", "lasso", "ridge", "rf", "xgboost"),
  metric = "calibration_slope",
  target_performance,
  n_reps_total = 1000,
  mean_or_assurance = "assurance",
  auto_n_reps_total = TRUE,
  ...
)

Arguments

signal_parameters

Integer. Number of candidate predictors associated with the outcome (i.e., true signal features).

noise_parameters

Integer. Number of candidate predictors not associated with the outcome (noise features). Default is 0.

complexity

Integer in 1:4 selecting the data-generating signal structure (see Data control). Default 1. Complexity 4 requires signal_parameters = 5.

data_control

Optional named list controlling the predictors (see Data control). Default NULL (generator defaults).

outcome_prevalence

Numeric in (0, 1). Target prevalence of the binary outcome in the intended modelling context.

maximum_achievable_cstatistic

Numeric in (0, 1). Maximum achievable C-statistic with effectively unlimited data. This is used to calibrate the data-generating mechanism and is not the minimum acceptable threshold.

model

Character string specifying the modelling algorithm. One of "glm" (logistic regression), "lasso", "ridge", "rf" (random forest), or "xgboost" (gradient-boosted trees).

metric

Character string naming the performance metric used to assess the sample size; defaults to "calibration_slope". Metric identifiers use one canonical form throughout the package: "calibration_slope", "auc", "r2", "cindex" and "csse". Metrics where smaller is better, such as "calibration_in_the_large", cannot be used as search targets.

"calibration_slope" is the slope from regressing the observed outcome on the model's linear predictor in held-out data; 1 indicates perfect calibration, and values below 1 indicate overfitting. Note that for the machine-learning models ("lasso", "ridge", "rf", "xgboost") this is converted internally to the calibration slope squared error for optimisation and translated back before results are returned; you don't need to do anything, and target_performance is still given on the calibration slope scale. Results derived this way carry a footnote marker in the printed output.

"csse" is the calibration slope squared error, -(1 - s)^2 for a calibration slope s, so that larger is better and 0 is perfect calibration. It can be requested directly, which is mainly useful for advanced use and for comparison against the internal conversion described above. When requesting it directly you are responsible for supplying target_performance on the CSSE scale: a calibration slope target of 0.9 corresponds to a CSSE target of -0.01. No adjustment is applied on your behalf, and results are reported on the CSSE scale.

target_performance

Numeric. Minimum acceptable value of the selected performance metric M^*; the algorithm searches for the smallest n meeting the chosen criterion with respect to this threshold.

n_reps_total

Integer. Total number of simulation replications used by the engine across the search.

mean_or_assurance

Character string, either "mean" or "assurance". Controls whether the minimum n is defined by the mean-based criterion or the assurance-based criterion (with the assurance level \delta controlled by the engine's defaults or additional arguments in ...).

auto_n_reps_total

Logical. When TRUE (default) and method = "mlpwr", n_reps_total is raised to 2000 for designs with at most 10 signal parameters fitted with "lm", "glm" or "coxph", where the Gaussian-process search is limited by simulation budget rather than by the design. It has no effect with the default method = "curve". The rule only ever increases the budget, and reports the change with a message. Set FALSE to use n_reps_total exactly as supplied.

...

Additional options passed to simulate_custom() (e.g., assurance level \delta, per-iteration settings).

Value

An object of class "pmsims" containing the estimated minimum sample size and simulation diagnostics (inputs, fitted GP curve, intermediate evaluations, and summary metrics).

Criteria

Two formulations are supported.

Mean-based: find the smallest n such that the expected model performance exceeds the target M^*, i.e.

\min_n \; \mathbb{E}_{D_n}\{ M \mid D_n \} \ge M^*.
Assurance-based: find the smallest n such that the probability the performance exceeds M^* is at least \delta (e.g. 0.80), i.e.

\min_n \; \mathbb{P}_{D_n}\!\left( M \mid D_n \ge M^* \right) \ge \delta.

Here, M is the chosen performance metric and the probability/expectation is over repeated samples of training data of size n. The assurance criterion explicitly accounts for variability across training sets; models with higher variance typically require larger n to satisfy it.

Data control

complexity selects the signal structure of the data-generating mechanism: 1 purely linear, 2 linear + quadratic, 3 linear + quadratic + interaction, 4 the Friedman function. The Friedman function uses five predictors, so complexity 4 requires exactly five signal predictors (signal_parameters = 5) and stops with an error otherwise; any number of noise predictors can be added. data_control is an optional list fine-tuning the predictors:

nonlinear_strength: Numeric in ⁠[0, 1)⁠. Fraction of the signal variance carried by the nonlinear, linearly-inaccessible component. Applies to complexity 2 and 3 only; ignored (with a warning) for 1 and 4. If omitted, the generator's per-complexity default is used.
correlation: Numeric in [-1, 1]. Pairwise correlation among the candidate predictors, induced through a Gaussian copula. Default 0.3. This is the correlation of the underlying normal variables: for "normal" predictors it is also the observed correlation, but for other distributions the observed correlation is lower. For example, binary predictors with prevalence 0.3 correlate at about 0.18 when correlation = 0.3.
predictor_distribution: One of "normal", "uniform", "binary", "exponential", "lognormal", "t", "laplace". "binary" selects 0/1 predictors and requires binary_predictor_prevalence; any other value selects continuous predictors from that family. Default "normal".
binary_predictor_prevalence: Numeric in ⁠(0, 1)⁠. Prevalence of the binary predictors; required when predictor_distribution = "binary", ignored (with a warning) otherwise. Note: binary predictors are incompatible with complexity 2/3 because squaring a 0/1 variable returns itself.

Examples


set.seed(123)
est <- simulate_binary(
  signal_parameters = 3,
  noise_parameters = 0,
  complexity = 1,
  data_control = list(correlation = 0),
  outcome_prevalence = 0.50,
  maximum_achievable_cstatistic = 0.80,
  model = "glm",
  metric = "calibration_slope",
  target_performance = 0.9,
  mean_or_assurance = "assurance",
  min_sample_size = 50,
  max_sample_size = 1000,
  n_reps_total = 1000,
  auto_n_reps_total = FALSE,
  test_n = 30000,
  progress = FALSE
)
est
est$min_n
plot(est)

Minimum sample size for continuous‐outcome prediction models

Description

Compute the minimum sample size required to develop a prediction model with a continuous outcome. This wraps the same simulation engine as simulate_binary(): by default a learning curve fitted to every simulated replicate (method = "curve"), or the Gaussian-process search with method = "mlpwr". From user inputs (maximum achievable performance, target performance, etc.) it constructs a data-generating function, model-fitting function, and metric function, then searches for the smallest n meeting the chosen criterion.

Usage

simulate_continuous(
  signal_parameters,
  noise_parameters = 0,
  complexity = 1,
  data_control = NULL,
  maximum_achievable_rsquared,
  model = c("lm", "lasso", "ridge", "rf", "xgboost"),
  metric = "calibration_slope",
  target_performance,
  n_reps_total = 1000,
  mean_or_assurance = "assurance",
  auto_n_reps_total = TRUE,
  ...
)

Arguments

signal_parameters

Integer. Number of candidate predictors associated with the outcome (i.e., true signal features).

noise_parameters

Integer. Number of candidate predictors not associated with the outcome (noise features). Default is 0.

complexity

Integer in 1:4 selecting the data-generating signal structure (see Data control). Default 1. Complexity 4 requires signal_parameters = 5.

data_control

Optional named list controlling the predictors (see Data control). Default NULL (generator defaults).

maximum_achievable_rsquared

Numeric in (0, 1). Maximum achievable R^2 with effectively unlimited data. This is used to calibrate the data-generating mechanism and is not the minimum acceptable threshold.

model

Character string specifying the modelling algorithm. One of "lm" (linear regression), "lasso", "ridge", "rf" (random forest), or "xgboost" (gradient-boosted trees).

metric

target_performance

Numeric. Minimum acceptable value of the selected performance metric M^*; the algorithm searches for the smallest n meeting the chosen criterion with respect to this threshold.

n_reps_total

Integer. Total number of simulation replications used by the engine across the search.

mean_or_assurance

auto_n_reps_total

...

Additional options passed to simulate_custom() (e.g., assurance level \delta, per-iteration settings).

Value

An object of class "pmsims" containing the estimated minimum sample size and simulation diagnostics (inputs, fitted GP curve, intermediate evaluations, and summary metrics).

Criteria

Two formulations are supported.

Mean-based: find the smallest n such that the expected model performance exceeds the target M^*, i.e.

\min_n \; \mathbb{E}_{D_n}\{ M \mid D_n \} \ge M^*.
Assurance-based: find the smallest n such that the probability the performance exceeds M^* is at least \delta (e.g. 0.80), i.e.

\min_n \; \mathbb{P}_{D_n}\!\left( M \mid D_n \ge M^* \right) \ge \delta.

Data control

nonlinear_strength: Numeric in ⁠[0, 1)⁠. Fraction of the signal variance carried by the nonlinear, linearly-inaccessible component. Applies to complexity 2 and 3 only; ignored (with a warning) for 1 and 4. If omitted, the generator's per-complexity default is used.
correlation: Numeric in [-1, 1]. Pairwise correlation among the candidate predictors, induced through a Gaussian copula. Default 0.3. This is the correlation of the underlying normal variables: for "normal" predictors it is also the observed correlation, but for other distributions the observed correlation is lower. For example, binary predictors with prevalence 0.3 correlate at about 0.18 when correlation = 0.3.
predictor_distribution: One of "normal", "uniform", "binary", "exponential", "lognormal", "t", "laplace". "binary" selects 0/1 predictors and requires binary_predictor_prevalence; any other value selects continuous predictors from that family. Default "normal".
binary_predictor_prevalence: Numeric in ⁠(0, 1)⁠. Prevalence of the binary predictors; required when predictor_distribution = "binary", ignored (with a warning) otherwise. Note: binary predictors are incompatible with complexity 2/3 because squaring a 0/1 variable returns itself.

Examples


set.seed(123)
est <- simulate_continuous(
  signal_parameters = 3,
  noise_parameters = 0,
  complexity = 1,
  data_control = list(correlation = 0),
  maximum_achievable_rsquared = 0.50,
  model = "lm",
  metric = "calibration_slope",
  target_performance = 0.9,
  mean_or_assurance = "assurance",
  min_sample_size = 50,
  max_sample_size = 1000,
  n_reps_total = 1000,
  auto_n_reps_total = FALSE,
  test_n = 30000,
  progress = FALSE
)
est
est$min_n
plot(est)

Minimum sample size for custom simulation workflows

Description

Compute the minimum sample size required to achieve a target level of predictive performance using user-defined simulation components. simulate_custom() is the low-level interface in pmsims: users supply a data-generating function, a model-fitting function, and a metric function, and the chosen search engine estimates the smallest n meeting the selected performance criterion.

Usage

simulate_custom(
  data_function,
  model_function,
  metric_function,
  target_performance,
  c_statistic = NULL,
  mean_or_assurance = "assurance",
  test_n = 30000,
  min_sample_size = NULL,
  max_sample_size = NULL,
  n_reps_total = 1000,
  n_reps_per = 20,
  method = "curve",
  progress = TRUE,
  verbose = FALSE,
  verify_reps = 100,
  max_n = NULL,
  cores = 1L,
  live_plot = FALSE,
  ...
)

Arguments

data_function

Function taking a single argument, n, giving the training sample size, and returning a dataset that can be passed to model_function.

model_function

Function that fits a model to the dataset returned by data_function. It must take the generated dataset as its only argument and return a fitted model object.

metric_function

target_performance

Numeric target value for the chosen performance metric. The search aims to find the smallest sample size n for which the selected criterion is met relative to this threshold.

c_statistic

Optional numeric value used only by the internal start-value heuristics for some outcome and metric combinations. In most custom workflows this should be left as NULL.

mean_or_assurance

Character string specifying the criterion used to define the minimum sample size. Must be either "mean" or "assurance".

test_n

Integer size of the test dataset used to evaluate model performance. This should usually be large enough that test-set variability is negligible relative to the training-sample search.

min_sample_size

Optional integer lower bound for the sample-size search. If supplied, max_sample_size must also be supplied.

max_sample_size

n_reps_total

Integer total number of simulation replications allocated to the search. The search evaluates approximately n_reps_total / n_reps_per candidate sample sizes.

n_reps_per

Integer number of simulation replications performed at each candidate sample size.

method

Character string specifying the search engine: "curve" (default; fits a learning curve to all replicates and searches where it crosses the target, see Details), "mlpwr" (adaptive start values, then mlpwr's Gaussian-process search), "bisection" or "mlpwr-bs".

progress

Logical flag controlling whether a progress bar is shown: the curve engine's own, or mlpwr's for the mlpwr-based methods.

verbose

Logical flag controlling engine-specific diagnostic output when supported. For the bisection engine, setting verbose = TRUE stores the iteration history on the returned object.

verify_reps

Integer number of fresh replicates simulated at the returned sample size to check that it meets the target (see status). Set to 0 to skip the check.

max_n

Largest sample size the search may try. If performance there is still clearly below the target, the search stops with status "not_bracketed". Defaults to 200,000 for random forests and xgboost (one batch of replicates at a million rows would take hours) and 1,000,000 otherwise.

cores

Number of cores for running the replicates of each batch in parallel (Unix-like systems only). Results are the same for any number of cores, because each replicate has its own random-number stream.

live_plot

Logical; with method = "curve", redraw the learning curve after each batch of replicates so the search can be watched as it runs (interactive sessions with a graphics device only; it does not change the results).

...

Additional arguments passed to the selected search engine.

Details

With method = "curve", every replicate is kept and a monotone learning curve C(n) = a - b n^{-c} is fitted to the criterion at each sample size evaluated (weighted by a smoothed model of the replicate spread; the ceiling a is capped at the metric's maximum). A pilot doubles (or halves) from a heuristic start until the criterion is seen on both sides of the target; each further batch of n_reps_per replicates is placed where the fitted curve crosses the target. The search range is therefore never fixed. The answer is where the final fitted curve crosses the target. A target is declared unreachable if a bootstrap upper bound for the curve's ceiling a is below it; the search also stops (status "not_bracketed") when the crossing is very likely beyond max_n. These stops extrapolate the fitted curve and should be read as heuristics.

The answer is a median-unbiased estimate of where the criterion reaches the target. The interval in diagnostics$curve$n_ci reflects Monte Carlo error given the learning-curve shape, not uncertainty about the shape; a shape-free check (diagnostics$crosscheck_n) is reported alongside, and flagged when it differs by more than 10%. The verification of the answer has limited power for skewed metrics such as CSSE: passing it does not show that the target is met.

Value

An object of class "pmsims" containing the estimated minimum sample size min_n (numeric; NA when no sample size could be found) and a status:

"ok": The search found min_n, and simulating it again confirmed the target is met (or the check was skipped).
"not_verified": Simulating min_n again gave performance clearly below the target, so min_n is likely too small.
"not_bracketed": No sample size up to max_n met the target, or the search returned no sample size.
"unreachable": With method = "curve": the learning curve levels off below the target, so no sample size reaches it.
"replicates_failed": At least half the simulation replicates failed to fit or score the model, so no sample size was estimated. With method = "curve", also when 20% or more failed at every sample size tried, so no learning curve could be fitted.

status_message explains a status other than "ok", verification holds the check at min_n, and diagnostics records the search bounds, whether min_n lies on one of them, mlpwr's warnings, Gaussian-process restarts and failed replicates.

Random numbers

Each simulation replicate runs on its own random-number stream, derived from a single draw from the session's RNG. Results are therefore reproducible with set.seed(), and a replicate's data do not depend on how many random numbers other parts of the search consumed.

Examples

# Three independent predictors with a population R-squared of 0.5.
data_fun <- function(n) {
  x1 <- rnorm(n)
  x2 <- rnorm(n)
  x3 <- rnorm(n)
  y <- (x1 + x2 + x3) / sqrt(3) + rnorm(n)
  data.frame(y = y, x1 = x1, x2 = x2, x3 = x3)
}

model_fun <- function(dat) {
  stats::lm(y ~ ., data = dat)
}

# Calibration slope evaluated on independent test data.
metric_fun <- function(test_data, fit, model) {
  preds <- stats::predict(fit, newdata = test_data)
  unname(stats::coef(stats::lm(test_data$y ~ preds))[2])
}
attr(metric_fun, "metric") <- "calibration_slope"


set.seed(123)
est <- simulate_custom(
  data_function = data_fun,
  model_function = model_fun,
  metric_function = metric_fun,
  target_performance = 0.9,
  mean_or_assurance = "assurance",
  min_sample_size = 25,
  max_sample_size = 1000,
  n_reps_total = 1000,
  test_n = 30000,
  progress = FALSE
)
est
est$min_n
plot(est)

Minimum sample size for survival‐outcome prediction models

Description

Compute the minimum sample size required to develop a prediction model with a time-to-event (survival) outcome. As with the other wrappers, this uses a simulation-based learning curve, fitted by default to every simulated replicate (method = "curve"; method = "mlpwr" uses a Gaussian-process surrogate instead), to locate the smallest n meeting the chosen performance criterion.

Usage

simulate_survival(
  signal_parameters,
  noise_parameters = 0,
  complexity = 1,
  data_control = NULL,
  maximum_achievable_cindex,
  baseline_hazard = 1,
  censoring_rate,
  model = c("coxph", "lasso", "ridge", "rf", "xgboost"),
  metric = "calibration_slope",
  target_performance,
  n_reps_total = 1000,
  mean_or_assurance = "assurance",
  auto_n_reps_total = TRUE,
  ...
)

Arguments

signal_parameters

Integer. Number of candidate predictors associated with the outcome (i.e., true signal features).

noise_parameters

Integer. Number of candidate predictors not associated with the outcome (noise features). Default is 0.

complexity

Integer in 1:4 selecting the data-generating signal structure (see Data control). Default 1. Complexity 4 requires signal_parameters = 5.

data_control

Optional named list controlling the predictors (see Data control). Default NULL (generator defaults).

maximum_achievable_cindex

Numeric in (0, 1). Maximum achievable C-index with effectively unlimited data. This is used to calibrate the data-generating mechanism and is not the minimum acceptable threshold.

baseline_hazard

Numeric greater than 0. Baseline hazard level used by the data-generating mechanism (e.g., the constant hazard in an exponential baseline). Larger values imply shorter event times, all else equal.

censoring_rate

Numeric in [0, 1). Proportion of individuals expected to be censored in the simulated datasets (administrative or random censoring). Higher values imply fewer observed events for a fixed n.

model

Character string specifying the modelling algorithm. One of "coxph" (Cox proportional hazards), "lasso", "ridge", "rf" (random survival forest), or "xgboost" (gradient boosting with a Cox objective).

metric

target_performance

Numeric. Minimum acceptable value of the selected performance metric M^*; the algorithm searches for the smallest n meeting the chosen criterion with respect to this threshold.

n_reps_total

Integer. Total number of simulation replications used by the engine across the search.

mean_or_assurance

auto_n_reps_total

...

Additional options passed to simulate_custom() (e.g., assurance level \delta, per-iteration settings).

Value

An object of class "pmsims" containing the estimated minimum sample size and simulation diagnostics (inputs, fitted GP curve, intermediate evaluations, and summary metrics).

Criteria

Two formulations are supported.

Mean-based: find the smallest n such that the expected model performance exceeds the target M^*, i.e.

\min_n \; \mathbb{E}_{D_n}\{ M \mid D_n \} \ge M^*.
Assurance-based: find the smallest n such that the probability the performance exceeds M^* is at least \delta (e.g. 0.80), i.e.

\min_n \; \mathbb{P}_{D_n}\!\left( M \mid D_n \ge M^* \right) \ge \delta.

Data control

nonlinear_strength: Numeric in ⁠[0, 1)⁠. Fraction of the signal variance carried by the nonlinear, linearly-inaccessible component. Applies to complexity 2 and 3 only; ignored (with a warning) for 1 and 4. If omitted, the generator's per-complexity default is used.
correlation: Numeric in [-1, 1]. Pairwise correlation among the candidate predictors, induced through a Gaussian copula. Default 0.3. This is the correlation of the underlying normal variables: for "normal" predictors it is also the observed correlation, but for other distributions the observed correlation is lower. For example, binary predictors with prevalence 0.3 correlate at about 0.18 when correlation = 0.3.
predictor_distribution: One of "normal", "uniform", "binary", "exponential", "lognormal", "t", "laplace". "binary" selects 0/1 predictors and requires binary_predictor_prevalence; any other value selects continuous predictors from that family. Default "normal".
binary_predictor_prevalence: Numeric in ⁠(0, 1)⁠. Prevalence of the binary predictors; required when predictor_distribution = "binary", ignored (with a warning) otherwise. Note: binary predictors are incompatible with complexity 2/3 because squaring a 0/1 variable returns itself.

Examples


set.seed(123)
est <- simulate_survival(
  signal_parameters = 1,
  noise_parameters = 0,
  complexity = 1,
  data_control = list(correlation = 0),
  maximum_achievable_cindex = 0.70,
  baseline_hazard = 0.01,
  censoring_rate = 0.30,
  model = "coxph",
  metric = "calibration_slope",
  target_performance = 0.9,
  mean_or_assurance = "assurance",
  min_sample_size = 25,
  max_sample_size = 500,
  n_reps_total = 1000,
  auto_n_reps_total = FALSE,
  test_n = 30000,
  progress = FALSE
)
est
est$min_n
plot(est)

Detailed summary of a `pmsims` result

Description

The same display as print.pmsims(), with the implementation-level detail that the default print method hides: internal metric identifiers, the engine settings used for the search, and any quantities recorded on an internal search scale.

Usage

## S3 method for class 'pmsims'
summary(object, ..., max_width = 80, verbose = TRUE)

Arguments

object

A pmsims object.

...

Currently unused.

max_width

Maximum console width used for the summary.

verbose

Logical. Set to FALSE for the default print display.

Value

object, invisibly.

Tuning function for a survival outcome model

Description

Finds the beta_signal that produces the target C-index under the exact data-generating settings (complexity, nonlinear_strength, correlation, distribution, predictor_type). The event rate is set exactly by quantile censoring, and the baseline hazard only fixes the time scale.

Usage

survival_tuning(
  target_prevalence,
  target_performance,
  candidate_features,
  proportion_noise_features,
  complexity = 1,
  nonlinear_strength = NULL,
  correlation = 0.3,
  distribution = "normal",
  predictor_type = "continuous",
  binary_prevalence = 0,
  n_sim = 50000,
  n_validate = 1e+05,
  beta_interval = c(1e-04, 20),
  beta_tol = 1e-04,
  tolerance = 0.02
)

Arguments

target_prevalence

Target event rate (proportion of events) in (0, 1). Equivalently 1 - censoring_rate.

target_performance

Target C-index in (0.5, 1).

candidate_features

Total number of predictors (signal + noise).

proportion_noise_features

Proportion of candidate_features that are noise predictors. Must be in [0, 1).

complexity

Integer 1-4. Default = 1.

nonlinear_strength

Fraction of signal variance carried by the nonlinear component (C2/C3 only), in [0, 1). When NULL (default), the complexity-level default is used: C1 = 0, C2 = 0.2, C3 = 0.3, C4 = 0.

correlation

Common pairwise predictor correlation. Default = 0.3 (to match generate_survival_data()).

distribution

Global continuous predictor distribution family. Default = "normal".

predictor_type

"continuous" (default) or "binary".

binary_prevalence

Bernoulli probability for binary predictors.

n_sim

Sample size for LP simulation and the bisection. Default 50000.

n_validate

Sample size for the independent validation. Default 100000.

beta_interval

Search interval c(lo, hi) for beta_signal. The upper bound is doubled automatically if the C-index at hi is below target.

beta_tol

Absolute convergence tolerance for the beta bisection.

tolerance

Acceptable absolute deviation of the validated event rate and C-index from their targets before a warning is raised. Default 0.02.

Value

A named numeric vector: lambda_opt (baseline hazard / time scale), beta_signal, event_rate, cindex, var_lp_unit.

Package {pmsims}

pmsims: Simulation-Based Sample Size Tools for Prediction Models

Description

Author(s)

See Also

default_model_generators Generate appropriate model based on input arguments

Description

Usage

Format

Value

adaptive_startvalues Derive adaptive sample sizes

Description

Usage

Arguments

Tuning function for a binary outcome model

Description

Usage

Arguments

Value

The Bisection Engine

Description

Usage

Arguments

Value

Learning-curve engine

Description

Usage

Arguments

Calculate performance metrics Calculates the performance metrics for a model given a sample size n.

Description

Usage

Arguments

Value

mlpwr engine

Description

Usage

Arguments

mlpwr-bs Hybrid engine using bisection to determine initial range and mlpwr for search

Description

Usage

Arguments

Value

Tuning function for a continuous outcome model

Description

Usage

Arguments

Value

Create default data generating functions

Description

Usage

Arguments

Value

Simulate binary outcome data

Description

Usage

Arguments

Value

Simulate continuous outcome data

Description

Usage

Arguments

Value

References

Construct the linear predictor

Description

Usage

Arguments

Details

Value

References

Generate the n x p predictor matrix

Description

Usage

Arguments

Value

Simulate survival outcome data

Description

Usage

Arguments

Value

Plot sample-size learning curves for `pmsims` outputs

Print a `pmsims` result