statistical-power

Files (6)

• SKILL.md
• closed_form_recipes.md · references/
• effect_sizes.md · references/
• simulation_based_power.md · references/
• power.py · scripts/
• simulate_power.py · scripts/

SKILL.md (first lines)

# Statistical Power & Sample Size

## Overview

Power analysis plans the probability of rejecting a specified null under an assumed alternative. It answers **how many independent experimental units are needed to detect a scientifically important effect**, or what effects a feasible sample could detect. Choose the inferential goal first: precision, equivalence, noninferiority, or sequential monitoring need their own calculations; the bundled superiority-test helpers do not cover them.

Four quantities are locked together for any given test: **sample size (n)**, **effect size**, **significance level (α)**, and **power (1 − β)**. For a fixed design, analysis, and nuisance parameters, fixing three permits solving for the fourth when a solution exists. Every calculation in this skill is some rearrangement of that relationship.

This skill covers the two ways to do power analysis:
- **Analytical** power calculations (exact under some standard-test assumptions; asymptotic approximations for others) — see `references/closed_form_recipes.md`.
- **Simulation / Monte Carlo** (requires a credible data-generating process and a calibrated planned analysis) — see `references/simulation_based_power.md`.

For choosing and converting effect sizes — usually the hardest part — see `references/effect_sizes.md`.

## When to Use This Skill

- Determining required sample size before collecting data (a priori power analysis)
- Finding the minimum detectable effect (MDE) for a fixed, already-determined sample size
- Producing power curves (power vs. n, or power vs. effect size) for a grant or protocol
- Justifying a sample size for an IRB submission, grant, or pre-registration
- Powering designs with unequal group sizes or non-1:1 allocation
- Planning complex mixed models, GLMs, clustering, survival, mediation, and interactions through simulation when simpler approximations do not fit
- Accounting for multiple comparisons, attrition/dropout, or clustering in the sample-size estimate

## Installation

The local examples were checked on Python 3.13 with statsmodels 0.15.0, SciPy
1.18.1, NumPy 2.5.3, pandas 2.3.3, and Matplotlib 3.11.2. Use an environment
separate from the repository's development environment:

```bash
uv venv --python 3.13 .venv-power
uv pip install --python .venv-power/bin/python "statsmodels==0.15.0" "scipy==1.18.1" "numpy==2.5.3" "pandas==2.3.3" "matplotlib==3.11.2"
# Optional comparison / survival methods (also checked for current API use):
uv pip install --python .venv-power/bin/python "pingouin==0.7.0" "lifelines==0.30.3"
```

On Windows use `.venv-power/Scripts/python.exe`. Lifelines 0.30.3 requires
`pandas<3`; the tested pin above accommodates it. Mixed models and GLMs are
included in statsmodels. Record versions and seeds with the output; numerical
smoke tests do not establish a study's effect assumptions or Type I error control.

---

## The one decision that drives everything: the effect size

Power calculations are only as trustworthy as the effect size you feed them. **Do not invent a number.** Use, in rough order of preference:

1. A **minimally important effect** — the smallest effect that would actually change a decision or matter scientifically/clinically (the "smallest effect size of interest", SESOI). This is the most defensible basis: you power to detect what matters, not what you hope to see.
2. A **pilot or prior-study estimate**, with uncertainty and selection bias considered. Small pilots are imprecise; publication or significance-based selection can inflate effects. Use a justified uncertainty model or sensitivity range instead of an arbitrary shrinkage factor.
3. A **convention** (Cohen's small/medium/large) only as a last resort, and say so explicitly.

Whatever you pick, run a **sensitivity analysis**: report how required n changes across a plausible range of effect sizes, not a single point. A power analysis presented as one number hides its biggest source of uncertainty. See `references/effect_sizes.md` for benchmarks and conversions between d, f, r, η², odds ratios, and Cohen's h/w.

> **Avoid post-hoc ("observed") power.** Computing power from the effect size you just estimated is circular: for standard tests it largely restates the test statistic/p-value and adds no independent evidence of adequacy. If a study is already done and you want to know what it could have detected, report a **sensitivity analysis** (MDE at the achieved n) or, better, the confidence interval around the observed effect. This is a common reviewer complaint — do not produce observed power even if asked without flagging the issue.

---

## Quick recipes (closed-form)

The bundled `scripts/power.py` wraps statsmodels and SciPy into one consistent interface so you don't have to remember which solver belongs to which test. Run from `skills/statistical-power/scripts/` or add that directory to `sys.path`.

```python
from power import sample_size, power, mde, power_curve

# 1. How many per group to detect Cohen's d = 0.5, two-sided, 80% power?
sample_size(test="t_ind", effect_size=0.5, power=0.80, alpha=0.05)
# -> 64 in sample 1; equal allocation gives 64 in sample 2

# 2. Two groups, 3:1 allocation (e.g. more controls than cases)
sample_size(test="t_ind", effect_size=0.5, power=0.80, ratio=3.0)

# 3. Fixed n=30/group — what's the minimum detectable d at 80% power?
mde(test="t_ind", nobs1=30, power=0.80, alpha=0.05)

# 4. One-way ANOVA, 4 groups, detect Cohen's f = 0.25
sample_size(test="anova", effect_size=0.25, k_groups=4, power=0.80)

# 5. Two proportions: 0.40 vs 0.55 (auto-converts to Cohen's h)
sample_size(test="two_proportions", prop1=0.40, prop2=0.55, power=0.80)

# 6. Correlation: detect r = 0.30
sample_size(test="correlation", effect_size=0.30, power=0.80)

# 7. Power curve for the grant figure
power_curve(test="t_ind", effect_size=0.5, n_range=range(10, 120, 5),
            save="power_curve.png")
```

For two-sample tests the return is **n1**, with `n2 = ceil(ratio * n1)`; `ratio=n2/n1`. Recheck power using the realized integer ratio. ANOVA rounds total n to a multiple of `k_groups`; paired n counts pairs. One-sided alternatives use signed effects and `"larger"`/`"smaller"`. Proportion MDEs return signed Cohen's h and need a baseline to convert to feasible probabilities.

Supported `test=` values: `t_ind` (two independent means), `t_paired`/`t_one` (paired or one-sample mean), `anova` (one-way), `two_proportions`, `one_proportion`, `correlation`, `chi2` (goodness-of-fit / contingency via effect size *w*), `linear_regression` (R² increment / f²). Full argument tables and the underlying statsmodels calls are in `references/closed_form_recipes.md`.

---

## When analytical assumptions do not fit: simulate

Use an analytical method when its design and assumptions match the planned analysis. For **logistic/Poisson regression, mixed-effects / repeated-measures models, cluster-randomized trials, survival analysis, mediation, or multi-way interactions**, simulation is often useful when available approximations omit material design features. The logic is always the same three steps:

1. **Simulate** a dataset from your assumed truth (the effect you want to detect, plus realistic noise, baseline rates, cluster structure, etc.).
2. **Analyze** it with the *exact* test/model you plan to use on the real data.
3. **Repeat** many times (≥1,000; 5,000–10,000 for a stable estimate near 80%). Power is the fraction of replicates in which the test is significant.

`scripts/simulate_power.py` provides a reusable harness plus worked examples (two-group difference, logistic regression, cluster-randomized trial with an ICC, and a linear mixed model). The core is just:

```python
from simulate_power import simulate_power, example_two_group_difference

# Runnable software check: n is per group, effect is a raw mean difference.
gen_and_test = example_two_group_difference(effect=0.5, sd=1.0, alpha=0.05)
est = simulate_power(gen_and_test, n=64, n_sims=2000, alpha=0.05, seed=0)
print(est)  # power, 95% Monte Carlo CI, failure and warning counts
```

The callback must return a boolean rejection decision using the planned alpha internally. The harness's `alpha` argument does not threshold returned p-values or pass alpha into the callback; returning a raw p-value now raises `TypeError` instead of counting a nonzero float as rejection. The harness rejects nonpositive sample/replicate counts and invalid search bounds; an unmet target at the sample-size cap raises an explicit error. The examples check convergence and finite p-values. Expected failures raise `SimulationFitError`; the harness counts them as non-rejections and reports `n_failures`/`failure_reasons`, retaining every replicate in the denominator. It also reports warning counts; unexpected errors propagate. First check Type I error under the null. A noisy bisection search yields a candidate n: verify nearby sizes with more replicates and a fresh seed.

Report the **Monte Carlo confidence interval** on the estimate (the harness returns it) to quantify simulation sampling error; it does not cover uncertainty in the assumed effect, model, or adaptively chosen n. See `references/simulation_based_power.md` for the full patterns, including how to search for the n that hits target power and how to model dropout and clustering.

---

From K-Dense-AI/scientific-agent-skills, MIT.

More from K-Dense-AI/scientific-agent-skills