statspai.power

power

Power and sample size calculations for causal inference and epidemiological designs.

Supports RCT, DID, RD, IV, cluster RCT, and OLS — plus epidemiological study designs (two-proportion, log-rank/survival, case-control) — with power curves, minimum detectable effect (MDE), and sample-size solving.

RetrodesignResult dataclass

Bases: ResultProtocolMixin

Result of :func:retrodesign.

Attributes:

Name	Type	Description
table	DataFrame	One row per hypothesised effect: effect, se, power, type_s and type_m.
power	float	Probability that the estimate is statistically significant (first row).
type_s	float	Probability that a significant estimate has the wrong sign.
type_m	float	Exaggeration ratio: expected absolute value of a significant estimate divided by the true effect.
alpha, dof	float

Examples:

>>> import statspai as sp
>>> out = sp.retrodesign(0.1, 3.28)
>>> round(out.power, 3), round(out.type_s, 2), round(out.type_m)
(0.05, 0.46, 77)

PowerResult

Bases: ResultProtocolMixin

Container for power analysis results.

Attributes:

Name	Type	Description
power	float or ndarray	Computed power value(s).
n	int, float, or np.ndarray	Sample size(s) used.
effect_size	float or ndarray	Effect size(s) used.
design	str	Name of the research design.
params	dict	All parameters passed to the power function.

Examples:

>>> import statspai as sp
>>> res = sp.power_rct(n=500, effect_size=0.3)
>>> isinstance(res, sp.PowerResult)
True
>>> res.design
'rct'
>>> round(float(res.power), 4)
0.9184

summary

summary() -> str

Return a formatted summary string.

plot

plot(ax: Any = None, figsize: tuple[float, float] = (8, 5), **kwargs: Any) -> Any

Plot power curve.

Works when n or effect_size was supplied as an array / range.

Parameters:

Name	Type	Description	Default
ax	matplotlib Axes		None
figsize	tuple		(8, 5)
**kwargs	passed to ``ax.plot``		{}

Returns:

Type	Description
matplotlib Axes

retrodesign

retrodesign(effect: Any, se: Any, alpha: float = 0.05, dof: Optional[float] = None, method: str = 'exact') -> RetrodesignResult

Power, sign-error rate and exaggeration ratio of a design.

For a hypothesised true effect and the standard error a study delivers, answers three questions about its statistically significant results: how often there is one (power), how often it points the wrong way (type S error) and by what factor it overstates the truth on average (type M error, the exaggeration ratio). With low power, significance selects the estimates that happen to be far from zero, so a published significant estimate is biased away from zero and may well have the wrong sign.

Parameters:

Name	Type	Description	Default
effect	float or array - like	Hypothesised true effect, from outside information (earlier literature, a plausible range), never the estimate of the study under review. Several values give one row each.	required
se	float or array - like	Standard error of the estimate.	required
alpha	float	Two-sided significance level.	0.05
dof	float	Degrees of freedom of the standard error, when it is estimated and the test is a t-test. Omitted: the standard error is known and the test is a z-test.	None
method	('exact', 'shifted')	Only matters with dof. 'exact': the estimate is normal around the true effect, its estimated standard error is independent of it with dof degrees of freedom, and the test compares their ratio with the t critical value, so the test statistic is noncentral t. 'shifted': the estimate is effect + se * t_dof and is compared with se times the critical value, the approximation in the code published with Gelman and Carlin (2014).	'exact'

Returns:

Type	Description
RetrodesignResult

Notes

Nothing is simulated. With mu = effect / se and z the critical value, the known-variance case is power = P(Z > z - mu) + P(Z < -z - mu), type_s the smaller term over the sum, and the exaggeration ratio follows from E[Z 1(Z > c)] = pdf(c). With dof and method='exact' power and sign error are noncentral-t probabilities and the exaggeration ratio integrates the known-variance expression over the distribution of the estimated standard error. With method='shifted' the partial expectation of Student's t, dof / (dof - 1) (1 + c^2 / dof) pdf(c), gives a closed form.

The two methods differ noticeably at low power and few degrees of freedom (effect 0.5, se 1, 20 degrees of freedom: exaggeration 4.63 exact, 5.18 shifted). A simulation of the t-test itself reproduces the exact values.

Examples:

>>> import statspai as sp
>>> out = sp.retrodesign([0.5, 2.8], se=1.0)
>>> [round(v, 2) for v in out.table["power"]]
[0.08, 0.8]
>>> round(float(out.table["type_m"].iloc[1]), 2)
1.13

References

gelman2014beyond

mde

mde(design: str, *, n: Any = None, power_target: float = 0.8, **kwargs: Any) -> PowerResult

Compute the Minimum Detectable Effect (MDE) for a given design.

Inverts the power function to find the smallest effect size that achieves power_target at the given sample size.

Parameters:

Name	Type	Description	Default
design	str	Research design (see :func:power).	required
n	int	Sample size. For 'cluster_rct' this is n_clusters.	None
power_target	float	Desired power (default 0.80).	0.8
**kwargs	Any	Design-specific parameters.	{}

Returns:

Type	Description
PowerResult	With .effect_size set to the MDE.

Examples:

>>> import statspai as sp
>>> sp.mde("did", n=1000, n_periods=10, n_treated_periods=5)
PowerResult(design='did', power=0.8000, n=1000, effect_size=0.1524)

power_cluster_rct

power_cluster_rct(n_clusters: Any, cluster_size: float, effect_size: Any, icc: float, alpha: float = 0.05, sigma: float = 1.0) -> PowerResult

Power for a Cluster-Randomised Controlled Trial.

Parameters:

Name	Type	Description	Default
n_clusters	int or array - like	Total number of clusters (treatment + control).	required
cluster_size	int or float	Average number of individuals per cluster.	required
effect_size	float or array - like	Effect in the units of the outcome. With the default sigma=1 this is the standardised effect (effect / standard deviation).	required
icc	float	Intra-cluster correlation coefficient.	required
alpha	float	Significance level.	0.05
sigma	float	Individual-level outcome standard deviation.	1.0

Returns:

Type	Description
PowerResult

Examples:

>>> import statspai as sp
>>> res = sp.power_cluster_rct(n_clusters=40, cluster_size=30,
...                            effect_size=0.3, icc=0.05)
>>> round(float(res.power), 4)
0.913

power_did

power_did(n: Any, effect_size: Any, n_periods: int, n_treated_periods: int, rho: float = 0.5, alpha: float = 0.05, sigma: float = 1.0, prop_treat: float = 0.5) -> PowerResult

Power for Difference-in-Differences.

Two groups observed for n_periods periods, the last n_treated_periods of them after treatment, with errors that follow a stationary AR(1) process within unit: the serially correlated setting of Burlig, Preonas & Woerman (2020). The estimator is the difference in post-minus-pre mean changes between the groups, which is the two-way fixed-effects coefficient on a balanced panel. Its variance is exact under the stated error process:

.. math::

\mathrm{Var}(\hat\tau) = \frac{\sigma^2\, c' R c}{p(1-p)\, n},

where c puts 1 / n_treated_periods on each post period and -1 / (n_periods - n_treated_periods) on each pre period, and R is the AR(1) correlation matrix rho ** |s - t|.

Parameters:

Name	Type	Description	Default
n	int or array - like	Total number of units (treated + control).	required
effect_size	float or array - like	Effect in the units of the outcome. With the default sigma=1 this is the standardised effect (effect / standard deviation).	required
n_periods	int	Total number of time periods.	required
n_treated_periods	int	Number of post-treatment periods.	required
rho	float	First-order autocorrelation of the errors within a unit. Positive autocorrelation lowers the variance of a pre/post change when the panel is short and raises it as the panel lengthens.	0.5
alpha	float	Significance level (two-sided).	0.05
sigma	float	Standard deviation of the error in one unit-period.	1.0
prop_treat	float	Share of units that are treated.	0.5

Returns:

Type	Description
PowerResult

Notes

Until the 2026-10 fix the variance was sigma^2 (1 + (T - 1) rho) T / (n T_pre T_post): a cluster design effect, which does not apply to a within-unit change, and without the 1 / (p (1 - p)) factor of a two-group comparison. With rho = 0 it overstated power (the standard error was half its value); with a large rho it understated it. For 200 units, ten periods of which five treated, rho = 0.5 and an effect of 0.2 it returned 0.48 where simulation gives 0.38.

Examples:

>>> import statspai as sp
>>> res = sp.power_did(n=1000, effect_size=0.1, n_periods=10,
...                    n_treated_periods=5)
>>> round(float(res.power), 4)
0.4516

power_iv

power_iv(n: Any, effect_size: Any, first_stage_f: Optional[float] = None, r2_z: Optional[float] = None, alpha: float = 0.05, sigma: float = 1.0) -> PowerResult

Power for Instrumental Variables / 2SLS estimation.

The 2SLS estimator has asymptotic variance sigma^2 / (n R2_z) for a regressor with unit variance, where R2_z is the share of the endogenous regressor's variance the instruments explain: only that share identifies the effect. Power is therefore that of OLS on a sample of n R2_z observations.

Parameters:

Name	Type	Description	Default
n	int or array - like	Sample size.	required
effect_size	float or array - like	Effect of a one-standard-deviation change in the endogenous variable, in the units of the outcome (standardised when sigma=1).	required
first_stage_f	float or None	First-stage F-statistic of a single instrument, converted with R2_z = F / (F + n). With several instruments pass r2_z.	None
r2_z	float or None	Partial R-squared of the first-stage regression. Alternative to first_stage_f; if both are given, first_stage_f takes precedence.	None
alpha	float	Significance level.	0.05
sigma	float	Error standard deviation.	1.0

Returns:

Type	Description
PowerResult

Notes

With neither first_stage_f nor r2_z the result is the OLS benchmark. The normal approximation is unreliable for a weak first stage (F below about 10), where 2SLS is biased and its t test over-rejects.

Until the 2026-10 fix the result was the OLS power multiplied by F / (F + 1), which barely moves for any F above 10 however little of the regressor the instrument explains: 400 observations, an effect of 0.2 and a first-stage F of 37 returned 0.95 where simulation gives 0.21.

References

stock2005testing

Examples:

>>> import statspai as sp
>>> res = sp.power_iv(n=1000, effect_size=0.2, first_stage_f=20)
>>> round(float(res.power), 4)
0.1413

power_ols

power_ols(n: Any, effect_size: Any, n_covariates: int = 0, r2_other: float = 0.0, alpha: float = 0.05, sigma: float = 1.0) -> PowerResult

Power for OLS regression (single coefficient of interest).

Parameters:

Name	Type	Description	Default
n	int or array - like	Sample size.	required
effect_size	float or array - like	Effect of the variable of interest in the units of the outcome (standardised when sigma=1).	required
n_covariates	int	Number of other covariates in the model.	0
r2_other	float	R-squared attributable to other covariates (reduces residual variance and thus improves power).	0.0
alpha	float	Significance level.	0.05
sigma	float	Outcome standard deviation.	1.0

Returns:

Type	Description
PowerResult

Examples:

>>> import statspai as sp
>>> res = sp.power_ols(n=500, effect_size=0.2)
>>> round(float(res.power), 4)
0.9939

power_rct

power_rct(n: Any, effect_size: Any, alpha: float = 0.05, ratio: float = 1.0, sigma: float = 1.0) -> PowerResult

Power for a two-arm Randomised Controlled Trial.

Parameters:

Name	Type	Description	Default
n	int or array - like	Total sample size (treatment + control).	required
effect_size	float or array - like	Effect in the units of the outcome. With the default sigma=1 this is the standardised effect (effect / standard deviation).	required
alpha	float	Significance level (two-sided).	0.05
ratio	float	Treatment / control allocation ratio (1 = equal allocation).	1.0
sigma	float	Outcome standard deviation, in the units of effect_size. Power depends on effect_size / sigma.	1.0

Returns:

Type	Description
PowerResult

Examples:

>>> import statspai as sp
>>> res = sp.power_rct(n=500, effect_size=0.3)
>>> round(float(res.power), 4)
0.9184

power_rd

power_rd(n: Any, effect_size: Any, bandwidth: Optional[float] = None, kernel: str = 'triangular', density_at_cutoff: float = 1.0, alpha: float = 0.05, sigma: float = 1.0) -> PowerResult

Power for Regression Discontinuity designs.

Power of the conventional local-linear t test at a fixed bandwidth. Each side's boundary intercept has asymptotic variance sigma^2 C_K / (n h f), where n h f is the expected number of observations within h of the cutoff on that side and C_K is the kernel's boundary constant (4 for uniform, 4.8 for triangular, about 4.5 for Epanechnikov), so

.. math::

\mathrm{se}(\hat\tau) = \sigma \sqrt{\frac{2\, C_K}{n\, h\, f}}.

The robust bias-corrected test that sp.rdrobust reports by default is less powerful than this; for a design calculation anchored on pilot data use sp.rdpower (Cattaneo, Titiunik & Vazquez-Bare 2019).

Parameters:

Name	Type	Description	Default
n	int or array - like	Total sample size in the data.	required
effect_size	float or array - like	Effect at the cutoff in the units of the outcome. With the default sigma=1 this is the standardised effect.	required
bandwidth	float or None	Bandwidth on each side of the cutoff, in units of the running variable. If None, defaults to 0.5 (the whole support of a running variable on [0, 1] with the cutoff at its midpoint).	None
kernel	('triangular', 'uniform', 'epanechnikov')	Kernel used for local weighting.	'triangular'
density_at_cutoff	float	Density of the running variable at the cutoff (default 1.0 for a uniform running variable on [0,1]).	1.0
alpha	float	Significance level.	0.05
sigma	float	Conditional outcome std dev near the cutoff.	1.0

Returns:

Type	Description
PowerResult

Notes

Until the 2026-10 fix the standard error was that of a difference between two local means with an efficiency discount, about half the local-linear value: 2,000 observations and an effect of 0.25 returned 0.93 where simulation of the local-linear test gives 0.74.

References

cattaneo2019power

Examples:

>>> import statspai as sp
>>> res = sp.power_rd(n=2000, effect_size=0.25)
>>> round(float(res.power), 4)
0.7229

power_case_control

power_case_control(n_cases: Optional[ArrayLike] = None, odds_ratio: float = 2.0, exposure_prevalence: float = 0.3, *, ratio: float = 1.0, alpha: float = 0.05, alternative: str = 'two-sided', power_target: Optional[float] = None, test: str = 'chi2') -> PowerResult

Power (or number of cases) for an unmatched case-control study.

Parameters:

Name	Type	Description	Default
n_cases	int, array-like, or None	Number of cases. None + power_target solves for the number of cases.	None
odds_ratio	float	Exposure odds ratio to detect (must be > 0, != 1).	2.0
exposure_prevalence	float	Exposure prevalence among controls (the source-population exposure probability), in (0, 1).	0.3
ratio	float	Number of controls per case.	1.0
alpha	float	Significance level.	0.05
alternative	('two-sided', 'one-sided')		"two-sided"
power_target	float	Desired power; solve for the number of cases when n_cases=None.	None
test	('chi2', 'wald')	"wald": the unpooled two-proportion z test, Phi(|p1 - p0| / se1 - z) with se1 the unpooled standard error (both tails counted when two-sided, since 1.42). "chi2": the Pearson chi-squared test, whose null standard error uses the pooled proportion, with both tails counted when two-sided -- Stata power twoproportions (controls as group 1, cases as group 2). .. versionchanged:: 1.30.0 The default moved from "wald" (no package reference) to "chi2", the test of Stata power twoproportions and R power.prop.test. Pass test="wald" for the old power.	"chi2"

Returns:

Type	Description
PowerResult

Notes

The control exposure prevalence p0 and the odds ratio imply a case exposure prevalence p1 = (OR p0) / (1 + p0 (OR - 1)). Power is then a two-proportion comparison between cases (n_cases) and controls (ratio * n_cases).

Examples:

>>> import statspai as sp
>>> res = sp.power_case_control(n_cases=200, odds_ratio=2.0,
...                             exposure_prevalence=0.3)
>>> res.design
'case_control'
>>> round(float(res.power), 4)  # Pearson chi2 test, as Stata power twoproportions
0.9171
>>> # Number of cases for 80% power at an odds ratio of 2.0.
>>> int(sp.power_case_control(odds_ratio=2.0, exposure_prevalence=0.3,
...                           power_target=0.8).n)
141

power_logrank

power_logrank(n: Optional[ArrayLike] = None, hazard_ratio: float = 0.5, prob_event: float = 1.0, *, ratio: float = 1.0, alpha: float = 0.05, alternative: str = 'two-sided', power_target: Optional[float] = None) -> PowerResult

Power (or sample size) for a two-arm log-rank / survival comparison.

Implements the Schoenfeld (1983) formula: power depends on the number of observed events, D = n * prob_event, and the log hazard ratio.

Parameters:

Name	Type	Description	Default
n	int, array-like, or None	Total sample size. None + power_target solves for n.	None
hazard_ratio	float	Hazard ratio between the two arms (must be > 0, != 1).	0.5
prob_event	float	Probability that a randomly chosen subject is observed to have the event during follow-up (the overall event rate). D = n*prob_event.	1.0
ratio	float	Allocation ratio n2 / n1.	1.0
alpha	float	Significance level.	0.05
alternative	('two-sided', 'one-sided')		"two-sided"
power_target	float	Desired power; solve for n when supplied with n=None.	None

Returns:

Type	Description
PowerResult

Notes

With allocation share p = ratio/(1+ratio), the required number of events is D = (z_alpha + z_beta)^2 / (p(1-p) (ln HR)^2) and the power for a given D is Phi(sqrt(D p(1-p)) |ln HR| - z_alpha), plus Phi(-sqrt(D p(1-p)) |ln HR| - z_alpha) when two-sided (counted since 1.42).

Examples:

>>> import statspai as sp
>>> res = sp.power_logrank(n=300, hazard_ratio=0.6, prob_event=0.7)
>>> res.design
'logrank'
>>> round(float(res.power), 4)
0.9592
>>> # Sample size for 80% power at a hazard ratio of 0.6.
>>> int(sp.power_logrank(hazard_ratio=0.6, prob_event=0.7,
...                      power_target=0.8).n)
172

power_two_proportions

power_two_proportions(n: Optional[ArrayLike] = None, p1: float = 0.5, p2: float = 0.5, *, ratio: float = 1.0, alpha: float = 0.05, alternative: str = 'two-sided', power_target: Optional[float] = None, test: str = 'wald') -> PowerResult

Power (or sample size) to detect a difference between two proportions.

Parameters:

Name	Type	Description	Default
n	int, array-like, or None	Total sample size (both groups). Pass None together with power_target to solve for the smallest n achieving that power.	None
p1	float	Outcome probabilities in group 1 (reference) and group 2.	0.5
p2	float	Outcome probabilities in group 1 (reference) and group 2.	0.5
ratio	float	Allocation ratio n2 / n1 (1.0 = equal allocation).	1.0
alpha	float	Significance level.	0.05
alternative	('two-sided', 'one-sided')	Test sidedness (R's "two.sided" is accepted; other values are refused).	"two-sided"
power_target	float	Desired power; when supplied with n=None the function returns the required total sample size.	None
test	('wald', 'chi2')	"wald": the unpooled-variance z test, Phi(|p1 - p2| / se - z_alpha) with se = sqrt(p1(1-p1)/n1 + p2(1-p2)/n2). "chi2": the Pearson chi-squared test, whose null standard error uses the pooled proportion; this is Stata power twoproportions (its default test(chi2)) and R power.prop.test(strict = TRUE), and the default of :func:power_case_control. Two-sided power counts both tails under either test. Group sizes are n / (1 + ratio) and n ratio / (1 + ratio) without rounding, as Stata's nfractional option (Stata otherwise rounds each group to an integer). .. versionadded:: 1.42 1.42 keeps "wald" as the default; see MIGRATION.md for the planned change to "chi2".	"wald"

Returns:

Type	Description
PowerResult

Notes

Two-sided power is Phi(shift - z) + Phi(-shift - z), so a zero difference has power alpha. Before 1.42 only the first term was counted.

Examples:

>>> import statspai as sp
>>> res = sp.power_two_proportions(n=400, p1=0.5, p2=0.65)
>>> isinstance(res, sp.PowerResult)
True
>>> res.design
'two_proportions'
>>> round(float(res.power), 4)
0.8665
>>> # Solve for the total sample size achieving 80% power.
>>> int(sp.power_two_proportions(p1=0.5, p2=0.65, power_target=0.8).n)
334

power_ttest

power_ttest(n: Optional[float] = None, delta: Optional[float] = None, power: Optional[float] = None, *, sd: float = 1.0, alpha: float = 0.05, type: str = 'two-sample', alternative: str = 'two-sided', ratio: float = 1.0) -> PowerResult

Power, sample size or detectable difference of a t test.

Give two of n, delta and power; the third is computed. The calculation uses the noncentral t distribution, as Stata's power onemean / twomeans / pairedmeans and R's pwr.t.test do.

Parameters:

Name	Type	Description	Default
n	float	Sample size: the total over both groups for a two-sample test (Stata's n(); pwr.t.test takes the size of one group), the number of observations for a one-sample test, the number of pairs for a paired test.	None
delta	float	The difference to detect, in the units of the outcome: the difference between the two means, the distance of the mean from its null value, or the mean of the paired differences. With the default sd=1 it is the standardised effect (Cohen's d).	None
power	float	Target power, between alpha and 1.	None
sd	float	Standard deviation of the outcome, common to the two groups; of the paired differences for type='paired'.	1
alpha	float	Significance level.	0.05
type	('two-sample', 'one-sample', 'paired')		'two-sample'
alternative	('two-sided', 'greater', 'less')	'greater' and 'less' are one-sided tests in the direction of a positive and a negative delta.	'two-sided'
ratio	float	n2 / n1 for a two-sample test (Stata's nratio()).	1

Returns:

Type	Description
PowerResult	power, n and effect_size (delta / sd). When n is solved for it is rounded up to whole observations in each group, as Stata reports it, and power is the power actually achieved at that size; params['n_exact'] keeps the fractional solution that pwr.t.test prints. params also holds n1, n2, delta and df.

Notes

The two-sided power counts rejections in both tails. The detectable difference is returned as a positive number (negative for alternative='less').

Examples:

>>> import statspai as sp
>>> res = sp.power_ttest(n=128, delta=0.5)        # 64 per group
>>> round(res.power, 4)
0.8015
>>> sp.power_ttest(delta=0.5, power=0.8).n        # total, both groups
128
>>> round(sp.power_ttest(n=128, power=0.8).params["delta"], 4)
0.4991
>>> sp.power_ttest(delta=0.3, power=0.9, alpha=0.01, type="paired").n
169