Overview
multigrain (Multiple Testing using Graphical Approaches; Improve power with Numerical Optimisation) is an R package that finds the optimal graphical multiple testing procedure for your clinical trial.
Confirmatory trials routinely test multiple hypotheses — across endpoints, doses, or subpopulations — and must strongly control the family-wise error rate (FWER). Graphical approaches are a popular framework for this, but choosing the hypothesis weights and transition matrix that make up the graph is typically done by hand.
multigrain automates this: given simulated p-values and a user-defined trial success measure (a function that scores each possible pattern of hypothesis rejections), it searches over hypothesis weights and transition weights to find the graph that maximises expected trial success. The result is a valid graphical test that strongly controls FWER.
Key features
| Step | Function | Purpose |
|---|---|---|
| Simulate p-values | simulate_pvalues() |
Draw p-values from a multivariate normal test-statistic model |
| Define success | trial_success() |
Specify what “trial success” means (compiled to C++ for speed) |
| Constrain the graph | graph_constraint() |
Fix weights, edges, or testing hierarchies |
| Optimise | graph_optimise() |
Find the graph that maximises expected trial success |
| Evaluate | calc_power_pvals() |
Compute local power, disjunctive/conjunctive power, and custom metrics |
Installation
multigrain is not on CRAN yet. Install the development version from GitHub with:
# install.packages("pak")
pak::pak("GSK-Biostatistics/multigrain")Quick start
Below is an end-to-end example with four hypotheses: two primary (H1, H2) and two secondary (H3, H4). Secondaries are gated by their corresponding primaries (H3 can only be tested after H1 is rejected; H4 after H2).
1. Simulate p-values
Assume one-sided test statistics are jointly normal. Provide the per-hypothesis nominal power (unadjusted for multiplicity) and a correlation matrix:
library(multigrain)
power_nominal <- c(0.95, 0.90, 0.85, 0.80)
corr_mat <- matrix(
c(
1.0, 0.5, 0.4, 0.2,
0.5, 1.0, 0.2, 0.4,
0.4, 0.2, 1.0, 0.5,
0.2, 0.4, 0.5, 1.0
),
nrow = 4, byrow = TRUE
)
set.seed(1)
pvals <- simulate_pvalues(
power_nominal = power_nominal,
corr_matrix = corr_mat,
nsim = 1e5
)2. Define trial success measures
Use trial_success() to express what “success” means for the trial. The variables r1, r2, r3, r4 are binary indicators — ri = 1 if hypothesis Hi is rejected, 0 otherwise. You can combine them with arithmetic (+, *) and logical (&&, ||) operators.
# Average power: mean proportion of hypotheses rejected
avg_power <- trial_success(0.25 * (r1 + r2 + r3 + r4))
# Disjunctive power: at least one rejection
disjunctive <- trial_success(r1 || r2 || r3 || r4)
# Custom: require a primary rejection, with extra credit for the
# gated secondary (H3 adds value only if H1 is also rejected, etc.)
custom_success <- trial_success(
0.25 * (2 * (r1 || r2) + r1 * r3 + r2 * r4)
)3. Constrain the graph (optional)
Use graph_constraint() to encode structural rules — for example, that only primaries receive initial alpha, and that secondaries can only receive alpha recycled from their gating primary:
gc <- graph_constraint(
hyp_constraint = c(NA, NA, 0, 0),
trans_constraint = matrix(
c(
0, NA, NA, 0,
NA, 0, 0, NA,
0, 1, 0, 0,
1, 0, 0, 0
),
nrow = 4, byrow = TRUE
)
)Here NA means “free to optimise” and a fixed number means “held constant”.
4. Optimise
g_opt <- graph_optimise(
pvals = pvals,
trial_success = custom_success,
graph_constraint = gc
)
summary(g_opt)
plot(g_opt)5. Evaluate
Compute multiplicity-adjusted power and any custom trial success measures for the optimised graph:
calc_power_pvals(
pvals,
hyp_weight = g_opt$hyp_weight,
trans_matrix = g_opt$trans_matrix,
custom_power = list(
AvgPower = avg_power,
Disjunctive = disjunctive,
CustomSuccess = custom_success
)
)Documentation
- Getting Started — full walkthrough of the asthma example above
- Graph Constraints — specifying and visualising constraints
- Function Reference