Phase III Clinical Trial Design
Build a Python toolkit to calculate sample sizes and generate randomization schedules.
Introduction
30 Second Summary
Every new medicine that reaches pharmacy shelves has passed through a Phase III clinical trial. Getting the trial design wrong means either enrolling thousands more patients than necessary, or missing a real treatment effect because the study was too small.
In this project, you will build a Python toolkit that calculates how many patients a trial needs, generates a stratified randomization schedule assigning 300 subjects to treatment arms, and validates the balance of that schedule with statistical tests.
What You'll Build
A complete Phase III parallel-group trial design toolkit consisting of three Python modules and a generated CSV dataset.
By the end of this project, you'll have:
- A working sample size calculator that computes patient counts for both continuous and time-to-event endpoints under standard Phase III assumptions.
- A reproducible stratified block randomization schedule that assigns 300 subjects across six strata defined by age group and disease severity.
- An automated validation report with a chi-squared balance test and a clear PASS/FAIL verdict.
- Secret Mission: a power curve visualization showing how required sample size changes across a range of effect sizes.
Are there any prerequisites?
This project assumes familiarity with hypothesis testing fundamentals (alpha, power, Type I and Type II errors) and basic Python (functions, loops, and data structures). You do not need prior experience with clinical trial design or biostatistics libraries.
Set Up Your Python Environment
Clinical trial design relies on precise statistical calculations that must produce the same results every time you run them. A dedicated Python environment keeps your scientific libraries isolated from other projects, ensuring reproducibility.
In this step, you'll create a project folder, set up a virtual environment, and install the four scientific computing libraries your trial design simulator depends on.
In this step, get ready to:
- Verify Python 3.11+ is available on your machine.
- Create a project folder with an isolated virtual environment.
- Install the scientific computing libraries (NumPy, SciPy, pandas, statsmodels).
Verify Python 3.11+ is installed
The scientific libraries you'll install require Python 3.11 or higher. Let's confirm your installation is ready.
- Press Cmd+Space (macOS) or the Windows key (Windows) to open your search bar.
- Type Terminal and press Enter to open it.
- Check your Python version by running this command:
macOS
python3 --version
Windows
python --version
You should see a version number of 3.11 or higher in the output.
Not seeing 3.11 or higher?
If you see a version below 3.11, download the latest Python from python.org and run the installer.
If the command is not found at all, visit python.org to install Python for your operating system.
help me install Python 3.11 on my operating system
Create your project folder and virtual environment
A virtual environment is an isolated space where Python packages live for one specific project. This prevents version conflicts between projects.
- Move to your Desktop by running this command:
cd ~/Desktop
- Create your project folder and move into it by running these commands:
mkdir trial-design-simulator
cd trial-design-simulator
- Confirm the folder exists by running:
ls
You should see an empty directory since you just created it.
- Create a virtual environment inside your project folder by running:
macOS
python3 -m venv .venv
Windows
python -m venv .venv
What does this command do?
python -m venv .venv tells Python to create a virtual environment in a hidden folder called .venv. All packages you install next will live inside this folder rather than in your system Python.
- Activate the virtual environment by running:
macOS
source .venv/bin/activate
Windows
.venv\Scripts\activate
You should see (.venv) appear at the beginning of your terminal prompt. This confirms the virtual environment is active.
Don't see (.venv) in your prompt?
On Windows, you may need to allow script execution first. Run Set-ExecutionPolicy -Scope CurrentUser RemoteSigned in PowerShell, then try the activate command again.
On macOS or Linux, make sure you typed source before the path. Without it, the command fails silently.
help me activate my Python virtual environment
Install the scientific computing dependencies
Your trial design simulator uses four Python libraries for statistical calculations, data manipulation, and randomization.
- Install all four packages by running:
pip install numpy==2.2.5 scipy==1.15.3 pandas==2.2.3 statsmodels==0.14.6
What are these libraries?
- numpy provides fast numerical arrays and random number generation for your randomization engine.
- scipy supplies statistical distributions and tests, including the chi-squared test for balance validation.
- pandas handles tabular data. You'll use it to build and export the randomization schedule as a CSV.
- statsmodels includes power analysis tools. It calculates how many patients each treatment arm needs.
The installation may take a minute. When it finishes, you should see a success message for each package.
- Verify the packages installed correctly by running:
pip list
You should see numpy, scipy, pandas, and statsmodels listed with their version numbers.
Seeing installation errors?
Make sure (.venv) appears in your prompt before running pip install. If not, re-run the activation command from the previous substep.
If you see a permissions error, you may have accidentally used your system Python. Deactivate with deactivate, re-activate the virtual environment, and try again.
help me fix my pip install errors
Your environment is ready. Next up, you'll build the sample size calculator that determines how many patients each arm of your trial needs.
Build the Sample Size Calculator
Before a Phase III trial enrolls a single patient, the study team must answer a critical question: how many subjects do we need?
Enrolling too few patients risks missing a real treatment effect. Enrolling too many wastes resources and exposes extra patients to experimental therapy unnecessarily. In this step, you will build a calculator that answers that question for two common endpoint types.
In this step, get ready to:
- Calculate sample size for a continuous endpoint using Cohen's d and TTestIndPower.
- Calculate sample size for a survival endpoint using the Schoenfeld formula.
- Print a formatted summary table showing all assumptions and results.
Create the continuous-endpoint calculator
A continuous endpoint measures a numeric outcome like blood pressure or cholesterol level. To size a trial for this type of endpoint, you convert the expected treatment difference into a standardized effect size (Cohen's d), then solve for the number of subjects that gives you the desired statistical power.
- Create a new file called sample_size.py inside your trial-design-simulator folder.
- Copy the following code into the file:
"""\nPhase III Clinical Trial - Sample Size & Power Calculator\nCovers: (1) continuous endpoint via two-sample t-test\n (2) time-to-event endpoint via Schoenfeld formula\n"""
import math
from scipy.stats import norm
from statsmodels.stats.power import TTestIndPower
def sample_size_continuous(
mean_control: float,
mean_treatment: float,
std_dev: float,
alpha: float = 0.05,
power: float = 0.80,
ratio: float = 1.0,
alternative: str = "two-sided",
) -> dict:
"""
Calculate sample size for a two-sample t-test (continuous endpoint).
Parameters
----------
mean_control : float
Expected mean in the control arm.
mean_treatment : float
Expected mean in the treatment arm.
std_dev : float
Common standard deviation (assumed equal across arms).
alpha : float
Significance level (two-sided).
power : float
Desired statistical power.
ratio : float
Allocation ratio n2/n1 (default 1:1).
alternative : str
'two-sided', 'larger', or 'smaller'.
Returns
-------
dict with effect_size, n_per_arm, n_total, and parameters.
"""
# Convert raw mean difference to standardized effect size
effect_size = abs(mean_treatment - mean_control) / std_dev
# Use statsmodels power analysis to solve for n
analysis = TTestIndPower()
n1 = analysis.solve_power(
effect_size=effect_size,
nobs1=None,
alpha=alpha,
power=power,
ratio=ratio,
alternative=alternative,
)
# Round up to whole subjects
n_per_arm = math.ceil(n1)
n_total = n_per_arm + math.ceil(n_per_arm * ratio)
return {
"effect_size_cohens_d": round(effect_size, 4),
"n_per_arm": n_per_arm,
"n_total": n_total,
"alpha": alpha,
"power": power,
"ratio": ratio,
"alternative": alternative,
}
What does this code do?
The function sample_size_continuous() takes the expected means for both arms and a common standard deviation, then computes Cohen's d (the standardized difference).
It passes that effect size to TTestIndPower.solve_power() from statsmodels, which uses the exact t-distribution to find how many subjects per arm you need. Finally, it rounds up to whole subjects because you cannot enroll half a patient.
What is Cohen's d?
Cohen's d expresses the treatment effect in units of standard deviation. A d of 0.84 means the treatment arm's mean is 0.84 standard deviations away from the control arm's mean.
Standardizing lets you compare effect sizes across studies that measure different things (blood pressure in mmHg, cholesterol in mg/dL, pain scores, etc.).
- Save sample_size.py.
- You will not run it yet because there is no main() block. The next chunk adds the survival calculator.
Seeing a red underline on the imports?
If your editor shows an import error for scipy or statsmodels, make sure your virtual environment is activated. In your terminal, run .venv\Scripts\activate (Windows) or source .venv/bin/activate (macOS/Linux).
help me fix my import errors in sample_size.py
Add the survival-endpoint calculator
Many Phase III trials measure time-to-event outcomes (time to disease progression, overall survival). The Schoenfeld formula calculates how many events you need to observe to detect a given hazard ratio with the log-rank test.
- In sample_size.py, add the following function below the sample_size_continuous() function:
def sample_size_survival(
hazard_ratio: float,
alpha: float = 0.05,
power: float = 0.80,
allocation_ratio: float = 0.5,
event_probability: float = 0.70,
) -> dict:
"""
Calculate sample size for a two-arm survival trial using the
Schoenfeld (1981) formula for the log-rank test.
Required number of events:
d = (z_{alpha/2} + z_{beta})^2 / (p * (1-p) * (log(HR))^2)
Parameters
----------
hazard_ratio : float
Postulated hazard ratio (treatment / control). < 1 means benefit.
alpha : float
Two-sided significance level.
allocation_ratio : float
Proportion allocated to treatment arm (0.5 = equal allocation).
event_probability : float
Expected overall probability of observing the event during follow-up.
Returns
-------
dict with required events, n_per_arm, n_total, and parameters.
"""
# Critical values from the standard normal distribution
z_alpha = norm.ppf(1 - alpha / 2)
z_beta = norm.ppf(power)
log_hr = math.log(hazard_ratio)
p = allocation_ratio
# Schoenfeld formula: required number of events
d = (z_alpha + z_beta) ** 2 / (p * (1 - p) * log_hr**2)
d = math.ceil(d)
# Convert events to total subjects using expected event probability
n_total = math.ceil(d / event_probability)
n_treatment = math.ceil(n_total * p)
n_control = n_total - n_treatment
return {
"hazard_ratio": hazard_ratio,
"required_events": d,
"n_treatment": n_treatment,
"n_control": n_control,
"n_total": n_total,
"alpha": alpha,
"power": power,
"allocation_ratio": allocation_ratio,
"event_probability": event_probability,
}
What does this code do?
The function sample_size_survival() implements the Schoenfeld (1981) formula for the log-rank test. It first computes the required number of events (deaths or progressions), then divides by the expected event probability to get total subjects.
The formula uses norm.ppf() (the inverse of the standard normal CDF) to convert alpha and power into z-scores. A hazard ratio of 0.75 means the treatment arm has 25% fewer events per unit time than the control.
Why implement the formula manually instead of using a library?
The Schoenfeld formula is only five lines of math. Implementing it directly lets you see each assumption (proportional hazards, equal allocation, expected event rate) rather than hiding them inside a library function.
For a regulatory submission, you need to justify every parameter. Seeing the formula in your code makes that justification transparent.
- Save sample_size.py.
- You still cannot run it yet. The next chunk adds the main() block that calls both functions and prints results.
Syntax error on the type hint?
If you see an error about list[str] or dict as a type hint, confirm your Python version is 3.11 or higher. Lowercase generic type hints (dict instead of Dict) require Python 3.9+.
help me fix type hint errors in my sample size calculator
Add the main block and run the calculator
The main() block ties everything together. It defines the clinical scenario (systolic blood pressure trial with a 15% reduction target) and prints a clear summary table with all assumptions labeled.
- At the bottom of sample_size.py, add the following code below the sample_size_survival() function:
def main():
print("=" * 65)
print(" PHASE III SAMPLE SIZE CALCULATIONS")
print("=" * 65)
# --- Continuous Endpoint ---
print("\n[1] CONTINUOUS ENDPOINT (Two-Sample t-Test)")
print("-" * 50)
# Example: systolic blood pressure reduction
# Control mean = 140 mmHg, Treatment achieves 15% reduction = 119 mmHg
# SD = 25 mmHg (typical for SBP)
control_mean = 140.0
treatment_mean = control_mean * (1 - 0.15) # 15% reduction
sd = 25.0
result_cont = sample_size_continuous(
mean_control=control_mean,
mean_treatment=treatment_mean,
std_dev=sd,
)
print(f" Control mean: {control_mean:.1f}")
print(f" Treatment mean: {treatment_mean:.1f} (15% reduction)")
print(f" Common SD: {sd:.1f}")
print(f" Cohen's d: {result_cont['effect_size_cohens_d']}")
print(f" Alpha (two-sided): {result_cont['alpha']}")
print(f" Power: {result_cont['power']}")
print(f" Allocation ratio: 1:1")
print(f" ---------------------------------")
print(f" N per arm: {result_cont['n_per_arm']}")
print(f" N total: {result_cont['n_total']}")
# --- Survival Endpoint ---
print("\n[2] TIME-TO-EVENT ENDPOINT (Log-Rank Test, Schoenfeld Formula)")
print("-" * 50)
result_surv = sample_size_survival(
hazard_ratio=0.75,
event_probability=0.70,
)
print(f" Hazard Ratio: {result_surv['hazard_ratio']}")
print(f" Alpha (two-sided): {result_surv['alpha']}")
print(f" Power: {result_surv['power']}")
print(f" Allocation: 1:1 (p = {result_surv['allocation_ratio']})")
print(f" Event probability: {result_surv['event_probability']}")
print(f" ---------------------------------")
print(f" Required events: {result_surv['required_events']}")
print(f" N per arm (approx): {result_surv['n_treatment']}")
print(f" N total: {result_surv['n_total']}")
print("\n" + "=" * 65)
if __name__ == "__main__":
main()
What does this code do?
The main() function sets up a realistic clinical scenario: a systolic blood pressure trial where the treatment arm targets a 15% reduction from a control mean of 140 mmHg (SD = 25 mmHg).
It calls both calculator functions and prints a formatted report. Every assumption (alpha, power, allocation ratio, event probability) is labeled so a reviewer can trace the calculation.
The if __name__ == "__main__" guard ensures the report only prints when you run the file directly, not when you import the functions into another module.
- Save sample_size.py.
- Run the calculator by executing this command in your terminal:
python sample_size.py
- You should see a formatted report with two sections: the continuous endpoint showing N per arm and N total, and the survival endpoint showing required events and total sample size.
Seeing ModuleNotFoundError?
If you see ModuleNotFoundError: No module named 'statsmodels', your virtual environment is not active. Activate it first, then re-run:
.venv\Scripts\activate
python sample_size.py
On macOS or Linux, use source .venv/bin/activate instead.
help me fix ModuleNotFoundError when running sample_size.py
✔️ Awesome, I've got everything!
Great. Double check you have saved sample_size.py and that running it prints both the continuous and survival calculations without errors.
ⓧ I'd like to double check the full code
Here is the complete sample_size.py file. Compare it against yours:
"""
Phase III Clinical Trial - Sample Size & Power Calculator
Covers: (1) continuous endpoint via two-sample t-test
(2) time-to-event endpoint via Schoenfeld formula
"""
import math
from scipy.stats import norm
from statsmodels.stats.power import TTestIndPower
def sample_size_continuous(
mean_control: float,
mean_treatment: float,
std_dev: float,
alpha: float = 0.05,
power: float = 0.80,
ratio: float = 1.0,
alternative: str = "two-sided",
) -> dict:
"""
Calculate sample size for a two-sample t-test (continuous endpoint).
Parameters
----------
mean_control : float
Expected mean in the control arm.
mean_treatment : float
Expected mean in the treatment arm.
std_dev : float
Common standard deviation (assumed equal across arms).
alpha : float
Significance level (two-sided).
power : float
Desired statistical power.
ratio : float
Allocation ratio n2/n1 (default 1:1).
alternative : str
'two-sided', 'larger', or 'smaller'.
Returns
-------
dict with effect_size, n_per_arm, n_total, and parameters.
"""
effect_size = abs(mean_treatment - mean_control) / std_dev
analysis = TTestIndPower()
n1 = analysis.solve_power(
effect_size=effect_size,
nobs1=None,
alpha=alpha,
power=power,
ratio=ratio,
alternative=alternative,
)
n_per_arm = math.ceil(n1)
n_total = n_per_arm + math.ceil(n_per_arm * ratio)
return {
"effect_size_cohens_d": round(effect_size, 4),
"n_per_arm": n_per_arm,
"n_total": n_total,
"alpha": alpha,
"power": power,
"ratio": ratio,
"alternative": alternative,
}
def sample_size_survival(
hazard_ratio: float,
alpha: float = 0.05,
power: float = 0.80,
allocation_ratio: float = 0.5,
event_probability: float = 0.70,
) -> dict:
"""
Calculate sample size for a two-arm survival trial using the
Schoenfeld (1981) formula for the log-rank test.
Required number of events:
d = (z_{alpha/2} + z_{beta})^2 / (p * (1-p) * (log(HR))^2)
Parameters
----------
hazard_ratio : float
Postulated hazard ratio (treatment / control). < 1 means benefit.
alpha : float
Two-sided significance level.
allocation_ratio : float
Proportion allocated to treatment arm (0.5 = equal allocation).
event_probability : float
Expected overall probability of observing the event during follow-up.
Returns
-------
dict with required events, n_per_arm, n_total, and parameters.
"""
z_alpha = norm.ppf(1 - alpha / 2)
z_beta = norm.ppf(power)
log_hr = math.log(hazard_ratio)
p = allocation_ratio
d = (z_alpha + z_beta) ** 2 / (p * (1 - p) * log_hr**2)
d = math.ceil(d)
n_total = math.ceil(d / event_probability)
n_treatment = math.ceil(n_total * p)
n_control = n_total - n_treatment
return {
"hazard_ratio": hazard_ratio,
"required_events": d,
"n_treatment": n_treatment,
"n_control": n_control,
"n_total": n_total,
"alpha": alpha,
"power": power,
"allocation_ratio": allocation_ratio,
"event_probability": event_probability,
}
def main():
print("=" * 65)
print(" PHASE III SAMPLE SIZE CALCULATIONS")
print("=" * 65)
# --- Continuous Endpoint ---
print("\n[1] CONTINUOUS ENDPOINT (Two-Sample t-Test)")
print("-" * 50)
# Example: systolic blood pressure reduction
# Control mean = 140 mmHg, Treatment achieves 15% reduction = 119 mmHg
# SD = 25 mmHg (typical for SBP)
control_mean = 140.0
treatment_mean = control_mean * (1 - 0.15) # 15% reduction
sd = 25.0
result_cont = sample_size_continuous(
mean_control=control_mean,
mean_treatment=treatment_mean,
std_dev=sd,
)
print(f" Control mean: {control_mean:.1f}")
print(f" Treatment mean: {treatment_mean:.1f} (15% reduction)")
print(f" Common SD: {sd:.1f}")
print(f" Cohen's d: {result_cont['effect_size_cohens_d']}")
print(f" Alpha (two-sided): {result_cont['alpha']}")
print(f" Power: {result_cont['power']}")
print(f" Allocation ratio: 1:1")
print(f" ---------------------------------")
print(f" N per arm: {result_cont['n_per_arm']}")
print(f" N total: {result_cont['n_total']}")
# --- Survival Endpoint ---
print("\n[2] TIME-TO-EVENT ENDPOINT (Log-Rank Test, Schoenfeld Formula)")
print("-" * 50)
result_surv = sample_size_survival(
hazard_ratio=0.75,
event_probability=0.70,
)
print(f" Hazard Ratio: {result_surv['hazard_ratio']}")
print(f" Alpha (two-sided): {result_surv['alpha']}")
print(f" Power: {result_surv['power']}")
print(f" Allocation: 1:1 (p = {result_surv['allocation_ratio']})")
print(f" Event probability: {result_surv['event_probability']}")
print(f" ---------------------------------")
print(f" Required events: {result_surv['required_events']}")
print(f" N per arm (approx): {result_surv['n_treatment']}")
print(f" N total: {result_surv['n_total']}")
print("\n" + "=" * 65)
if __name__ == "__main__":
main()
Your sample size calculator is working. Next up, you will build the stratified block randomization engine that assigns those subjects to treatment arms.
Build the Stratified Block Randomization Engine
Your sample size calculator now tells you how many patients each arm needs. The next question is equally critical: how do you assign those patients to Treatment or Control without introducing bias?
In this step, you will build a stratified block randomization engine. It generates a complete allocation schedule that balances treatment arms within prognostically important subgroups, using permuted blocks of varying sizes to prevent investigators from guessing the next assignment.
In this step, get ready to:
- Implement a permuted-block generation function that shuffles treatment labels within fixed-size blocks.
- Build a stratified randomization function that fills 6 strata with permuted blocks until each stratum reaches its target enrollment.
- Generate a 300-subject randomization schedule and export it to CSV.
Create the randomization module and permuted-block function
A permuted block is the smallest unit of randomization. For a block of size 4 with 1:1 allocation, the block contains exactly 2 Treatment and 2 Control labels, shuffled into a random order.
What is a permuted block?
Imagine putting 2 blue balls (Treatment) and 2 red balls (Control) into a bag and drawing them one at a time. Each complete draw of 4 guarantees perfect 1:1 balance at that point. By varying the bag size (4 or 6), you make the sequence harder to predict.
- Create a new file called randomization.py in your trial-design-simulator folder.
- Add the module docstring, imports, and the permuted-block function by pasting this code:
"""
Phase III Clinical Trial - Stratified Block Randomization Engine
Generates a reproducible randomization schedule stratified by:
- Age Group: <65, >=65
- Disease Severity: Mild, Moderate, Severe
Uses permuted blocks of varying sizes (4 and 6) with a fixed seed.
"""
import numpy as np
import pandas as pd
def generate_permuted_block(
block_size: int, rng: np.random.Generator
) -> list[str]:
"""
Generate a single permuted block of treatment assignments.
For a 1:1 allocation and block_size B, the block contains
B/2 'Treatment' and B/2 'Control' labels in random order.
Parameters
----------
block_size : int
Must be an even number (multiple of 2 for 1:1 allocation).
rng : numpy.random.Generator
Random number generator for reproducibility.
Returns
-------
list of str
Shuffled treatment assignments for the block.
"""
half = block_size // 2
block = ["Treatment"] * half + ["Control"] * half
rng.shuffle(block)
return block
What does this code do?
The function accepts a block size and a NumPy random Generator object. It builds a list with equal Treatment and Control labels, shuffles them in place, and returns the shuffled list.
Passing the rng object (rather than calling numpy.random.shuffle globally) ensures every call draws from the same reproducible sequence.
Why build randomization from scratch?
Libraries like stochatreat or randomizr can generate schedules for you. However, the core algorithm is only a few lines of code. Implementing it yourself exposes the mechanism a regulatory reviewer will ask about: how exactly are allocations generated and why are they balanced?
- Save randomization.py.
Seeing a syntax error on the type hint?
The list[str] return annotation requires Python 3.9+. If you see an error here, confirm you are running Python 3.11 or later inside your virtual environment.
help me fix a type hint syntax error in my randomization.py
Implement the stratified block randomization function
The stratification function iterates over each stratum, calculates how many subjects belong there, and fills that stratum with permuted blocks until the target count is reached.
- In randomization.py, add the following function below generate_permuted_block():
def stratified_block_randomize(
n_total: int,
strata_proportions: dict[str, float],
block_sizes: list[int],
seed: int = 42,
) -> pd.DataFrame:
"""
Generate a full stratified block randomization schedule.
Parameters
----------
n_total : int
Total number of subjects to randomize.
strata_proportions : dict
Mapping of stratum name to its expected proportion of enrollment.
Values must sum to 1.0.
block_sizes : list of int
Available block sizes (randomly chosen per block).
seed : int
Random seed for full reproducibility.
Returns
-------
pd.DataFrame
Randomization schedule with columns:
SubjectID, AgeGroup, Severity, Stratum, BlockID, Treatment
"""
rng = np.random.default_rng(seed)
records = []
subject_counter = 1
for stratum_name, proportion in strata_proportions.items():
n_stratum = round(n_total * proportion)
age_group, severity = stratum_name.split(" | ")
assigned = 0
block_counter = 1
while assigned < n_stratum:
# Randomly pick a block size
bs = rng.choice(block_sizes)
block = generate_permuted_block(bs, rng)
for assignment in block:
if assigned >= n_stratum:
break
records.append(
{
"SubjectID": f"SUBJ-{subject_counter:04d}",
"AgeGroup": age_group,
"Severity": severity,
"Stratum": stratum_name,
"BlockID": f"{stratum_name}_B{block_counter:03d}",
"Treatment": assignment,
}
)
subject_counter += 1
assigned += 1
block_counter += 1
return pd.DataFrame(records)
What does this code do?
- np.random.default_rng(seed) creates a NumPy Generator seeded at 42, so the schedule is identical every run.
- The outer loop iterates over each stratum. The stratum name (e.g. <65 | Mild) is split into its age group and severity components.
- The inner while loop picks a random block size (4 or 6), generates a permuted block, and assigns subjects one at a time until the stratum is full.
- The if assigned >= n_stratum: break guard handles the case where the last block exceeds the remaining slots.
- Save randomization.py.
Seeing an indentation error?
Make sure there is exactly one blank line between the end of the docstring in the previous function and the start of this function definition. Also check that the nested while/for loops use consistent 4-space indentation throughout.
help me fix an indentation error in randomization.py
Generate the 300-subject schedule and export to CSV
The main block defines the 6 strata (2 age groups crossed with 3 severity levels), calls the randomization function, and exports the result to a CSV file that downstream systems can consume.
- Add the following main() function at the bottom of randomization.py:
def main():
# Define strata: 2 age groups x 3 severity levels = 6 strata
# Proportions reflect typical enrollment patterns
strata_proportions = {
"<65 | Mild": 0.20,
"<65 | Moderate": 0.20,
"<65 | Severe": 0.10,
">=65 | Mild": 0.15,
">=65 | Moderate": 0.20,
">=65 | Severe": 0.15,
}
n_total = 300
block_sizes = [4, 6]
seed = 42
print("=" * 65)
print(" STRATIFIED BLOCK RANDOMIZATION")
print("=" * 65)
print(f"\n Total subjects: {n_total}")
print(f" Block sizes: {block_sizes}")
print(f" Random seed: {seed}")
print(f" Strata: {len(strata_proportions)}")
print(f" Stratification: Age Group (<65, >=65) x Severity (Mild, Moderate, Severe)")
schedule = stratified_block_randomize(
n_total=n_total,
strata_proportions=strata_proportions,
block_sizes=block_sizes,
seed=seed,
)
# Export to CSV
output_file = "randomization_schedule.csv"
schedule.to_csv(output_file, index=False)
print(f"\n Schedule exported: {output_file}")
print(f" Total records: {len(schedule)}")
# Summary
print("\n Allocation Summary:")
print(" " + "-" * 40)
summary = schedule.groupby("Stratum")["Treatment"].value_counts().unstack(fill_value=0)
print(summary.to_string(index=True))
print("\n" + "=" * 65)
if __name__ == "__main__":
main()
What does this code do?
- strata_proportions defines 6 strata with their expected enrollment fractions. The proportions sum to 1.0 and reflect that younger/milder patients enroll more readily.
- schedule.to_csv(output_file, index=False) writes the Pandas DataFrame to randomization_schedule.csv without the row index.
- The allocation summary uses groupby().value_counts().unstack() to produce a cross-tabulation of Treatment vs Control counts per stratum.
- Save randomization.py.
- Run the randomization engine by running this command:
python randomization.py
You should see a formatted report showing 300 total subjects, the 6 strata, and an allocation summary table with Treatment and Control counts per stratum.
- Confirm the CSV was created by running:
ls
You should see randomization_schedule.csv alongside sample_size.py and randomization.py in the output.
Seeing a ModuleNotFoundError for numpy or pandas?
Make sure your virtual environment is activated. Your terminal prompt should show (.venv) at the beginning. If not, activate it again with .venv\Scripts\activate (Windows) or source .venv/bin/activate (macOS/Linux), then re-run the script.
help me fix a ModuleNotFoundError when running randomization.py
✔️ Awesome, I've got everything!
Great. Double check you have saved randomization.py and that randomization_schedule.csv exists in your project folder.
ⓧ I'd like to double check the full code
"""
Phase III Clinical Trial - Stratified Block Randomization Engine
Generates a reproducible randomization schedule stratified by:
- Age Group: <65, >=65
- Disease Severity: Mild, Moderate, Severe
Uses permuted blocks of varying sizes (4 and 6) with a fixed seed.
"""
import numpy as np
import pandas as pd
def generate_permuted_block(
block_size: int, rng: np.random.Generator
) -> list[str]:
"""
Generate a single permuted block of treatment assignments.
For a 1:1 allocation and block_size B, the block contains
B/2 'Treatment' and B/2 'Control' labels in random order.
Parameters
----------
block_size : int
Must be an even number (multiple of 2 for 1:1 allocation).
rng : numpy.random.Generator
Random number generator for reproducibility.
Returns
-------
list of str
Shuffled treatment assignments for the block.
"""
half = block_size // 2
block = ["Treatment"] * half + ["Control"] * half
rng.shuffle(block)
return block
def stratified_block_randomize(
n_total: int,
strata_proportions: dict[str, float],
block_sizes: list[int],
seed: int = 42,
) -> pd.DataFrame:
"""
Generate a full stratified block randomization schedule.
Parameters
----------
n_total : int
Total number of subjects to randomize.
strata_proportions : dict
Mapping of stratum name to its expected proportion of enrollment.
Values must sum to 1.0.
block_sizes : list of int
Available block sizes (randomly chosen per block).
seed : int
Random seed for full reproducibility.
Returns
-------
pd.DataFrame
Randomization schedule with columns:
SubjectID, AgeGroup, Severity, Stratum, BlockID, Treatment
"""
rng = np.random.default_rng(seed)
records = []
subject_counter = 1
for stratum_name, proportion in strata_proportions.items():
n_stratum = round(n_total * proportion)
age_group, severity = stratum_name.split(" | ")
assigned = 0
block_counter = 1
while assigned < n_stratum:
# Randomly pick a block size
bs = rng.choice(block_sizes)
block = generate_permuted_block(bs, rng)
for assignment in block:
if assigned >= n_stratum:
break
records.append(
{
"SubjectID": f"SUBJ-{subject_counter:04d}",
"AgeGroup": age_group,
"Severity": severity,
"Stratum": stratum_name,
"BlockID": f"{stratum_name}_B{block_counter:03d}",
"Treatment": assignment,
}
)
subject_counter += 1
assigned += 1
block_counter += 1
return pd.DataFrame(records)
def main():
# Define strata: 2 age groups x 3 severity levels = 6 strata
# Proportions reflect typical enrollment patterns
strata_proportions = {
"<65 | Mild": 0.20,
"<65 | Moderate": 0.20,
"<65 | Severe": 0.10,
">=65 | Mild": 0.15,
">=65 | Moderate": 0.20,
">=65 | Severe": 0.15,
}
n_total = 300
block_sizes = [4, 6]
seed = 42
print("=" * 65)
print(" STRATIFIED BLOCK RANDOMIZATION")
print("=" * 65)
print(f"\n Total subjects: {n_total}")
print(f" Block sizes: {block_sizes}")
print(f" Random seed: {seed}")
print(f" Strata: {len(strata_proportions)}")
print(f" Stratification: Age Group (<65, >=65) x Severity (Mild, Moderate, Severe)")
schedule = stratified_block_randomize(
n_total=n_total,
strata_proportions=strata_proportions,
block_sizes=block_sizes,
seed=seed,
)
# Export to CSV
output_file = "randomization_schedule.csv"
schedule.to_csv(output_file, index=False)
print(f"\n Schedule exported: {output_file}")
print(f" Total records: {len(schedule)}")
# Summary
print("\n Allocation Summary:")
print(" " + "-" * 40)
summary = schedule.groupby("Stratum")["Treatment"].value_counts().unstack(fill_value=0)
print(summary.to_string(index=True))
print("\n" + "=" * 65)
if __name__ == "__main__":
main()
Your randomization engine is working. You now have a reproducible 300-subject schedule exported as a dataset. Next up, you will build a validation report that statistically confirms the schedule is balanced across all strata.
Build the Balance Validation Report
You now have a 300-subject randomization schedule sitting in a CSV file. A regulatory reviewer will ask one critical question: does the randomization actually produce balanced treatment arms across all strata?
In this step, you will build an automated validation report that reads the schedule, runs a statistical test for balance, and outputs a clear PASS or FAIL verdict. You will also create a README.md documenting the methodology behind all three modules.
In this step, get ready to:
- Load and analyze the randomization schedule for treatment balance.
- Test for systematic imbalance using a chi-squared test of independence.
- Generate a formatted PASS/FAIL validation report and document your project.
Create the validation script and load the schedule
The validation module needs to read the CSV that randomization.py generated and prepare the data for statistical testing.
- Create a new file called validation.py inside your trial-design-simulator folder and add the following code:
"""
Phase III Clinical Trial - Randomization Balance Validation Report
Reads the generated randomization schedule and checks:
1. Overall treatment balance
2. Per-stratum treatment balance
3. Chi-squared test of independence (Treatment vs Stratum)
4. Maximum absolute imbalance within any stratum
"""
import pandas as pd
from scipy.stats import chi2_contingency
def load_schedule(filepath: str = "randomization_schedule.csv") -> pd.DataFrame:
"""Load the randomization schedule from CSV."""
return pd.read_csv(filepath)
What does this code do?
The module docstring outlines the four checks the report will perform. The imports bring in pandas for data loading and chi2_contingency from scipy for the statistical test.
The load_schedule() function defaults to reading randomization_schedule.csv from the current directory. This keeps the module flexible for testing with different schedule files.
- Save validation.py.
You won't see output yet. The next chunk adds the analysis logic that this loader feeds into.
Compute balance metrics with chi-squared testing
A chi-squared test of independence checks whether treatment allocation is related to stratum membership. If your randomization is working correctly, the two should be independent.
- In validation.py, add the compute_balance_metrics() function below load_schedule():
def compute_balance_metrics(df: pd.DataFrame) -> dict:
"""
Compute balance metrics across strata.
Returns
-------
dict with overall_counts, stratum_counts, contingency_table,
chi2_result, and max_imbalance.
"""
# Overall counts
overall = df["Treatment"].value_counts().to_dict()
# Per-stratum counts
contingency = pd.crosstab(df["Stratum"], df["Treatment"])
stratum_counts = contingency.to_dict(orient="index")
# Chi-squared test of independence
chi2, p_value, dof, expected = chi2_contingency(contingency.values)
# Maximum absolute imbalance within any stratum
if "Treatment" in contingency.columns and "Control" in contingency.columns:
imbalance = abs(contingency["Treatment"] - contingency["Control"])
else:
cols = contingency.columns.tolist()
imbalance = abs(contingency[cols[0]] - contingency[cols[1]])
max_imbalance = int(imbalance.max())
worst_stratum = imbalance.idxmax()
return {
"overall_counts": overall,
"stratum_counts": stratum_counts,
"contingency_table": contingency,
"chi2": chi2,
"p_value": p_value,
"dof": dof,
"max_imbalance": max_imbalance,
"worst_stratum": worst_stratum,
}
What does this code do?
- overall counts how many subjects ended up in Treatment vs Control across the entire schedule.
- pd.crosstab() builds a contingency table crossing each stratum with each treatment arm. This is the input the chi-squared test needs.
- chi2_contingency() returns the test statistic, p-value, degrees of freedom, and expected frequencies. A high p-value (above 0.05) means no evidence of systematic imbalance.
- The imbalance calculation finds the worst stratum (the one with the largest difference between Treatment and Control counts). This is an operational check separate from the statistical test.
- Save validation.py.
You won't see a change yet. The next chunk adds the report printer that displays these metrics.
Seeing a NameError on chi2_contingency?
Make sure your import at the top reads from scipy.stats import chi2_contingency. If you imported scipy.stats without the specific function, the name won't be available directly.
Help me debug a NameError in my validation script.
Print the formatted validation report and run it
The report needs to present results in a format a study statistician can quickly scan. It ends with a definitive PASS or FAIL verdict based on two criteria: the chi-squared p-value and the maximum within-stratum imbalance.
- In validation.py, add the print_report() function below compute_balance_metrics():
def print_report(metrics: dict, max_block_size: int = 6) -> None:
"""Print a formatted validation report."""
print("=" * 65)
print(" RANDOMIZATION BALANCE VALIDATION REPORT")
print("=" * 65)
# Overall balance
print("\n [1] OVERALL TREATMENT ALLOCATION")
print(" " + "-" * 40)
for arm, count in metrics["overall_counts"].items():
print(f" {arm:12s}: {count}")
total = sum(metrics["overall_counts"].values())
print(f" {'Total':12s}: {total}")
# Per-stratum balance
print("\n [2] PER-STRATUM BALANCE")
print(" " + "-" * 40)
ct = metrics["contingency_table"]
print(ct.to_string(index=True))
# Chi-squared test
print("\n [3] CHI-SQUARED TEST OF INDEPENDENCE")
print(" " + "-" * 40)
print(f" H0: Treatment allocation is independent of stratum")
print(f" Chi-squared statistic: {metrics['chi2']:.4f}")
print(f" Degrees of freedom: {metrics['dof']}")
print(f" p-value: {metrics['p_value']:.4f}")
# Max imbalance
print("\n [4] MAXIMUM WITHIN-STRATUM IMBALANCE")
print(" " + "-" * 40)
print(f" Max |Treatment - Control|: {metrics['max_imbalance']}")
print(f" Worst stratum: {metrics['worst_stratum']}")
# Verdict
print("\n [5] VERDICT")
print(" " + "-" * 40)
chi2_pass = metrics["p_value"] > 0.05
imbalance_pass = metrics["max_imbalance"] <= max_block_size
if chi2_pass and imbalance_pass:
verdict = "PASS"
detail = "No significant imbalance detected."
elif not chi2_pass:
verdict = "FAIL"
detail = f"Chi-squared test significant (p={metrics['p_value']:.4f} <= 0.05)."
else:
verdict = "FAIL"
detail = f"Max imbalance ({metrics['max_imbalance']}) exceeds block size ({max_block_size})."
print(f" Chi-squared p > 0.05: {'PASS' if chi2_pass else 'FAIL'}")
print(f" Max imbalance <= {max_block_size}: {'PASS' if imbalance_pass else 'FAIL'}")
print(f" ---------------------------------")
print(f" OVERALL: ** {verdict} **")
print(f" {detail}")
print("\n" + "=" * 65)
What does the verdict logic check?
The report issues a PASS only when both criteria are met:
- Chi-squared p-value > 0.05 confirms no statistically significant relationship between treatment and stratum.
- Max imbalance <= block size confirms no single stratum has a treatment count difference exceeding the largest permuted block (6). This is an operational guardrail on top of the statistical test.
- Now add the main() block at the bottom of validation.py to tie everything together:
def main():
df = load_schedule()
metrics = compute_balance_metrics(df)
print_report(metrics)
if __name__ == "__main__":
main()
What does this code do?
The main() function chains the three steps: load the CSV, compute all metrics, and print the formatted report. The if __name__ guard means the report only prints when you run the file directly.
- Save validation.py.
- Run the validation report by executing this command:
python validation.py
You should see a multi-section report ending with a VERDICT section showing ** PASS ** and the message "No significant imbalance detected."
Seeing a FileNotFoundError?
This means randomization_schedule.csv is not in the same directory as validation.py. Make sure you ran python randomization.py in the previous step and that your terminal is in the trial-design-simulator folder.
Help me fix a FileNotFoundError in my validation script.
✔️ Awesome, I've got everything!
Great. Double check you have saved validation.py and confirmed the PASS verdict in your terminal.
ⓧ I'd like to double check the full code
"""
Phase III Clinical Trial - Randomization Balance Validation Report
Reads the generated randomization schedule and checks:
1. Overall treatment balance
2. Per-stratum treatment balance
3. Chi-squared test of independence (Treatment vs Stratum)
4. Maximum absolute imbalance within any stratum
"""
import pandas as pd
from scipy.stats import chi2_contingency
def load_schedule(filepath: str = "randomization_schedule.csv") -> pd.DataFrame:
"""Load the randomization schedule from CSV."""
return pd.read_csv(filepath)
def compute_balance_metrics(df: pd.DataFrame) -> dict:
"""
Compute balance metrics across strata.
Returns
-------
dict with overall_counts, stratum_counts, contingency_table,
chi2_result, and max_imbalance.
"""
# Overall counts
overall = df["Treatment"].value_counts().to_dict()
# Per-stratum counts
contingency = pd.crosstab(df["Stratum"], df["Treatment"])
stratum_counts = contingency.to_dict(orient="index")
# Chi-squared test of independence
chi2, p_value, dof, expected = chi2_contingency(contingency.values)
# Maximum absolute imbalance within any stratum
if "Treatment" in contingency.columns and "Control" in contingency.columns:
imbalance = abs(contingency["Treatment"] - contingency["Control"])
else:
cols = contingency.columns.tolist()
imbalance = abs(contingency[cols[0]] - contingency[cols[1]])
max_imbalance = int(imbalance.max())
worst_stratum = imbalance.idxmax()
return {
"overall_counts": overall,
"stratum_counts": stratum_counts,
"contingency_table": contingency,
"chi2": chi2,
"p_value": p_value,
"dof": dof,
"max_imbalance": max_imbalance,
"worst_stratum": worst_stratum,
}
def print_report(metrics: dict, max_block_size: int = 6) -> None:
"""Print a formatted validation report."""
print("=" * 65)
print(" RANDOMIZATION BALANCE VALIDATION REPORT")
print("=" * 65)
# Overall balance
print("\n [1] OVERALL TREATMENT ALLOCATION")
print(" " + "-" * 40)
for arm, count in metrics["overall_counts"].items():
print(f" {arm:12s}: {count}")
total = sum(metrics["overall_counts"].values())
print(f" {'Total':12s}: {total}")
# Per-stratum balance
print("\n [2] PER-STRATUM BALANCE")
print(" " + "-" * 40)
ct = metrics["contingency_table"]
print(ct.to_string(index=True))
# Chi-squared test
print("\n [3] CHI-SQUARED TEST OF INDEPENDENCE")
print(" " + "-" * 40)
print(f" H0: Treatment allocation is independent of stratum")
print(f" Chi-squared statistic: {metrics['chi2']:.4f}")
print(f" Degrees of freedom: {metrics['dof']}")
print(f" p-value: {metrics['p_value']:.4f}")
# Max imbalance
print("\n [4] MAXIMUM WITHIN-STRATUM IMBALANCE")
print(" " + "-" * 40)
print(f" Max |Treatment - Control|: {metrics['max_imbalance']}")
print(f" Worst stratum: {metrics['worst_stratum']}")
# Verdict
print("\n [5] VERDICT")
print(" " + "-" * 40)
chi2_pass = metrics["p_value"] > 0.05
imbalance_pass = metrics["max_imbalance"] <= max_block_size
if chi2_pass and imbalance_pass:
verdict = "PASS"
detail = "No significant imbalance detected."
elif not chi2_pass:
verdict = "FAIL"
detail = f"Chi-squared test significant (p={metrics['p_value']:.4f} <= 0.05)."
else:
verdict = "FAIL"
detail = f"Max imbalance ({metrics['max_imbalance']}) exceeds block size ({max_block_size})."
print(f" Chi-squared p > 0.05: {'PASS' if chi2_pass else 'FAIL'}")
print(f" Max imbalance <= {max_block_size}: {'PASS' if imbalance_pass else 'FAIL'}")
print(f" ---------------------------------")
print(f" OVERALL: ** {verdict} **")
print(f" {detail}")
print("\n" + "=" * 65)
def main():
df = load_schedule()
metrics = compute_balance_metrics(df)
print_report(metrics)
if __name__ == "__main__":
main()
Now that all three modules are complete, the project needs a README.md that documents the statistical methodology a reviewer would need to understand.
- Create a new file called README.md inside your trial-design-simulator folder and add the following content:
# Phase III Clinical Trial Design Simulator
A Python toolkit for the statistical design phase of a Phase III parallel-group randomized controlled trial.
## Components
| Module | Purpose |
|--------|----------|
| `sample_size.py` | Sample size & power calculations |
| `randomization.py` | Stratified block randomization schedule |
| `validation.py` | Balance validation report |
## Statistical Methodology
### Sample Size: Continuous Endpoint
The required sample size per arm for a two-sample t-test is derived from Cohen's d and implemented via `statsmodels.stats.power.TTestIndPower.solve_power`.
### Sample Size: Time-to-Event Endpoint
Schoenfeld (1981) formula for the log-rank test: d = (z_{alpha/2} + z_{beta})^2 / (p*(1-p)*(log HR)^2)
### Randomization: Stratified Permuted Blocks
Within each stratum, subjects are allocated using permuted blocks of randomly varying size (4 or 6). Strata: Age Group (<65, >=65) x Disease Severity (Mild, Moderate, Severe) = 6 strata total. Fixed seed=42 for reproducibility.
### Validation
Chi-squared test of independence (Treatment x Stratum) and maximum within-stratum imbalance check.
## Quick Start
```bash
python sample_size.py
python randomization.py
python validation.py
```
## Requirements
- Python >= 3.11
- numpy, scipy, pandas, statsmodels
Why include a README?
A regulatory statistician reviewing trial design documents needs to quickly understand the methodology without reading source code. The README serves as the entry point, summarizing each module's purpose and the statistical formulas behind the calculations.
- Save README.md.
- Confirm all project files exist by running:
ls
You should see sample_size.py, randomization.py, validation.py, randomization_schedule.csv, and README.md listed in the output.
Your Phase III trial design toolkit is complete. You have a sample size calculator, a stratified randomization engine, a validated schedule, and a balance report that confirms everything is working as expected.
Secret mission
Power Curve Visualization
Your sample size calculator gives a single number, but protocol teams need to see how sensitive that number is to their assumptions. Build a power curve that visualizes how required sample size changes across a range of effect sizes, helping the team judge whether the study is robust to optimistic assumptions.
Clean Up Your Resources
Clean Up Your Resources
Decide whether to keep your resources running, pause them to come back later, or delete them entirely. This project runs entirely locally, so there are no ongoing costs.
Resources you used:
- Project folder trial-design-simulator (includes virtual environment, Python source files, generated CSV, and power curve PNG).
Keep everything running
No action needed. Choose this if you're still actively building or want to extend the toolkit further.
- Your project folder stays on your Desktop with all files intact.
- The virtual environment remains ready to activate whenever you return.
- Python and pip remain installed on your system.
Pause - I'll come back to this later
Deactivate the virtual environment to free up your terminal session. Your files stay on disk for when you return.
- Deactivate the virtual environment by running:
deactivate
- When you come back, reactivate with:
cd ~/Desktop/trial-design-simulator
.venv\Scripts\activate
All your source files, the generated CSV, and the power curve PNG will still be there.
Delete - I don't want to use this again
Remove the entire project folder and all generated files. This cannot be undone.
- First, deactivate the virtual environment if it is still active:
deactivate
- Delete the project folder by running this command (macOS/Linux):
cd ~/Desktop
rm -rf trial-design-simulator
- On Windows, run this instead:
cd ~\Desktop
rmdir /s /q trial-design-simulator
This removes all Python source files, the virtual environment with installed packages, the randomization schedule CSV, the power curve PNG, and the README.
Nice Work!
Nice Work!
Nice work! You have built a complete Phase III clinical trial design simulator from scratch using Python.
This toolkit covers every stage a biostatistician needs before a study start-up meeting: sample size justification, randomization methodology, and validation.
You've learned how to:
- Build a dual-endpoint sample size calculator using Cohen's d for continuous outcomes and the Schoenfeld formula for time-to-event outcomes.
- Implement stratified permuted-block randomization from scratch with reproducible seeds and generate a 300-subject schedule stratified by age group and disease severity.
- Create an automated balance validation report using a chi-squared test of independence with a PASS/FAIL verdict.
- Produce a power curve visualization demonstrating how sample size responds to effect size assumptions.
Ready to quiz yourself?