How to calculate a Wilson interval for a binomial proportion in Statsmodels (with Code)

Editorial · Published 31 October 2026

The Short Answer (TL;DR)

Call proportion_confint with the observed success count, total trial count and method="wilson" to calculate a Wilson score interval for one binomial proportion. Set alpha to one minus the confidence level and state that the trials are independent with a common success probability.

from statsmodels.stats.proportion import proportion_confint
successes, trials, alpha = 12, 40, 0.05
low, high = proportion_confint(successes,trials,alpha=alpha,method="wilson",alternative="two-sided"); print(successes/trials,low,high)

Detailed Implementation

This example uses synthetic data and does not need a network connection. Install the recorded package versions in a separate Python environment before running it.

python -m pip install numpy==2.5.4 scipy==1.18.1 pandas==3.0.6 statsmodels==0.15.0

from pathlib import Path
import numpy as np
from statsmodels.stats.proportion import proportion_confint

def wilson_interval(successes, trials, confidence=.95):
    # This example accepts actual binomial counts, not percentages or weights.
    for value in (successes, trials):
        if isinstance(value, (bool, np.bool_)) or not isinstance(value, (int, np.integer)):
            raise ValueError("Successes and trials must be integer counts.")
    if trials <= 0 or not 0 <= successes <= trials:
        raise ValueError("Require trials > 0 and 0 <= successes <= trials.")
    if not np.isscalar(confidence) or not np.isfinite(confidence) or not 0 < confidence < 1:
        raise ValueError("Confidence must be a finite number strictly between zero and one.")
    return proportion_confint(
        count=successes, nobs=trials, alpha=1-confidence,
        method="wilson", alternative="two-sided",
    )

# Synthetic binary outcomes: 12 successes among 40 independent trials.
successes, trials, confidence = 12, 40, .95
low, high = wilson_interval(successes, trials, confidence)
p_hat = successes/trials
report = (
    f"successes={successes}; trials={trials}; sample proportion={p_hat:.6f}\n"
    f"{100*confidence:.0f}% Wilson score CI=[{low:.9f}, {high:.9f}]\n"
    f"percentage form=[{100*low:.6f}%, {100*high:.6f}%]\n"
)
Path("wilson-binomial.txt").write_text(report, encoding="utf-8")
print(report)

The sample proportion is 0.300000, with 95% Wilson interval [0.180748452, 0.454300188], or [18.074845%, 45.430019%].

Which values go into the function?

The example supplies 12 successes and 40 trials, giving an observed proportion of 0.30. The function receives the count 12, not the percentage 30 or a rounded rate. Setting alpha=0.05 requests 95% confidence; it does not request a 5% confidence interval. The code explicitly selects a two-sided Wilson interval rather than relying on the API defaults.

What does the resulting interval describe?

The estimated population success probability lies between about 0.180748 and 0.454300 under this interval procedure. Report it as a fraction, or multiply both endpoints by 100 to obtain about 18.074845% to 45.430019%. The confidence statement concerns repeated-sampling coverage of a method, not a 95% posterior probability for this particular fixed parameter. It is not a prediction interval for the next binary outcome.

Why do zero or all successes still leave uncertainty?

With zero successes in 40 trials, the Wilson upper endpoint is positive; with 40 successes, the lower endpoint is below one. An observed boundary proportion does not prove the underlying probability equals that boundary. Wilson endpoints remain in [0,1] and need not be symmetric around the observed proportion. The checks cover every possible success count for six trial sizes, including n=1, and separately reconstruct the score-equation roots.

Which assumptions and coverage limits matter?

A binomial model treats the trial count as fixed and the outcomes as independent Bernoulli observations sharing a common success probability. Repeated observations from the same person, household clusters, survey weights and sampling without replacement from a small finite population need design-aware inference. This synthetic count example does not establish those conditions for a real dataset. Wilson inverts a normal-score test; its finite-sample coverage is approximate and can fall below the nominal confidence level. The Clopper-Pearson method is a different exact binomial construction, available as method="beta", with its own conservatism and discrete behavior. An interval for one group does not by itself provide an interval or test for the difference between two groups.

Mathematical Context

Write p_hat=k/n and let z be the standard-normal quantile at 1-alpha/2. Inverting the score inequality n(p_hat-p)^2 <= z^2 p(1-p) produces the quadratic (n+z^2)p^2-(2k+z^2)p+k^2/n <= 0. Its two roots are the Wilson endpoints. Equivalently, divide p_hat+z^2/(2n), plus or minus z times sqrt(p_hat(1-p_hat)/n+z^2/(4n^2)), by 1+z^2/n. For k=12 and n=40, the 95% endpoints are 0.180748452 and 0.454300188. The probability, its estimate and the endpoints are dimensionless. More trials at the same observed rate narrow this example interval, but they do not repair dependence or sampling bias.

\[ p_{L,U}=\frac{\widehat p+\frac{z^2}{2n}\mp z\sqrt{\frac{\widehat p(1-\widehat p)}{n}+\frac{z^2}{4n^2}}}{1+\frac{z^2}{n}} \]

Reproducibility

Executed with Python 3.12.14; numpy 2.5.4, scipy 1.18.1, pandas 3.0.6, statsmodels 0.15.0. Automated checks confirm the stated output properties.