Hypothesis Testing in Python

A hypothesis test turns “this looks different” into a number you can defend. This guide covers the framework, z-scores and standardisation, proportion tests, t-tests, ANOVA, chi-square tests, power, the non-parametric fallbacks, and the decision tree for picking the right test in the first place.

Most analysis ends with a comparison. This month is up on last month. The treatment group converted better than the control. The Air customers pay more than the Road customers. Every one of those observations arrives with the same question attached: is this real, or did the numbers just land that way?

Hypothesis testing is the machinery for answering that. It does not tell you whether an effect exists, and it does not tell you whether an effect matters. What it tells you is narrower and more useful than either: how surprising your data would be in a world where nothing was going on. If the answer is “extremely”, you have grounds to argue that something is going on.

The mechanics below are mostly a page of code each. The part worth slowing down on is what the output means, because the p-value is the most widely misread number in applied statistics.

The data

Two fictional tables carry every example. deliveries holds courier records: delayed as “Yes” or “No”, unit_cost in pounds per parcel, service_level as Air, Road or Rail, contract_type and fuel_band as further categoricals, and parcel_weight_kg. depots holds paired measurements from a warehouse layout change, pick_time_before and pick_time_after, in minutes per pick at the same depot.

import numpy as np
import pandas as pd
import pingouin
from scipy.stats import norm, chisquare
from statsmodels.stats.proportion import proportions_ztest, proportion_confint

pingouin provides the t-tests, ANOVA, chi-square and non-parametric tests and is installed with pip install pingouin. normconverts test statistics into p-values, chisquare runs the goodness of fit test, and the two statsmodels imports handle proportions.

1. The framework

Every test in this guide follows the same four steps.

1. State the hypotheses:  H0 (null) against H1 (alternative)
2. Calculate a test statistic from the data
3. Convert the test statistic into a p-value
4. Compare the p-value to the significance level, alpha
   p <= alpha  ->  reject H0, the result is statistically significant
   p >  alpha  ->  fail to reject H0

Hypothesis testing is built like a courtroom. The null hypothesis is the presumption of innocence: there is no effect, no difference, nothing to see. It holds until the evidence against it becomes implausible. The alternative hypothesis is the case you are trying to make.

The test statistic measures how far your data sits from what the null predicted, scaled by the variation you would expect anyway. The p-value converts that distance into a probability: how often data at least this extreme would show up in a world where H0 is true. And alpha is the false positive rate you are willing to accept, conventionally 0.05, agreed before you look.

Two phrasings matter. “Fail to reject H0” is not “H0 is true”: a court that acquits has not proved innocence, it has found the evidence insufficient. And alpha belongs at the start of the analysis, not the end. Deciding on 0.05 and then adjusting to 0.1 because your p-value came out at 0.07 is not a threshold, it is a rationalisation.

2. What a p-value actually is

The p-value is the probability of observing data at least as extreme as yours, given that the null hypothesis is true.

Read that conditional carefully, because almost every misuse comes from reversing it. Three things a p-value is not.

It is not the probability that the null hypothesis is true. The calculation assumes H0 is true and works forwards from there. It cannot then turn round and tell you how likely that assumption was. Getting from a p-value to the probability of a hypothesis requires a prior, which is Bayesian territory and a different calculation entirely.

It is not the probability your result is due to chance. A p-value of 0.03 does not mean there is a 3% chance you got lucky. It means that if nothing were going on, data this extreme would show up 3% of the time.

It is not a measure of how large or important the effect is. With a large enough sample, a difference of no practical consequence returns a tiny p-value. Section 13 covers what to report alongside it.

Two kinds of mistake are possible and they trade off against each other. Rejecting a null hypothesis that was actually true is a Type I error, a false positive, and it happens at exactly the rate you set with alpha. Failing to reject a null that was actually false is a Type II error, a false negative, and its rate is called beta.

Setting alpha at 0.05 means accepting that 5% of the time, when nothing is going on, you will claim something is. Lower it to 0.01 and you make fewer false positives while missing more real effects. You cannot reduce both rates at once by adjusting alpha, only by collecting more data.

The complement of beta, one minus the false negative rate, is the test’s power: the probability of detecting an effect that genuinely exists. That is the number the next section is about, and it is the one almost nobody calculates.

3. Power, and how many observations you need

Everything above assumes the data already exists. Before it does, the useful question is how much of it you need, and that is a power calculation.

Four quantities are locked together: the significance level, the power you want, the size of the effect you care about detecting, and the sample size. Fix any three and the fourth is determined. In practice you fix alpha at 0.05, power at 0.8, and the smallest effect that would change a decision, then solve for n.

from statsmodels.stats.power import TTestIndPower
analysis = TTestIndPower()
n = analysis.solve_power(effect_size=0.5, alpha=0.05, power=0.8)
print(f"{n:.0f} observations per group")

solve_power solves for whichever argument you leave out, so omitting nobs1 returns the required sample size per group. effect_size here is Cohen’s d, the difference in means expressed in standard deviations, where 0.2 is small, 0.5 is medium and 0.8 is large. A medium effect at 80% power needs roughly 64 per group, and halving the effect you want to detect roughly quadruples the sample you need.

For proportions the effect size is Cohen’s h rather than d, and statsmodels will compute it for you:

from statsmodels.stats.proportion import proportion_effectsize
from statsmodels.stats.power import NormalIndPower
effect = proportion_effectsize(0.10, 0.08)      # detecting 10% against 8%
n = NormalIndPower().solve_power(effect_size=effect, alpha=0.05, power=0.8)

Run that and the answer is a few thousand per arm, which is why detecting small conversion differences takes far more traffic than people expect.

The reason this matters is that an underpowered study is worse than no study. At 30% power you will miss most real effects, so a null result tells you nothing. And the effects you do detect are systematically overstated, because only the samples that happened to land unusually far from zero cleared the threshold. That is why small studies produce large, exciting, unreplicable findings.

Running the calculation afterwards on the effect you observed, sometimes called post hoc power, is not informative and is a known statistical error: it is just a transformation of the p-value you already have. Power analysis belongs before data collection, and afterwards the confidence interval answers the question you were actually asking.

4. Z-scores and standardisation

Almost every test in this guide is the same calculation in different clothes: take the gap between what you observed and what H0 predicted, then divide by how much that quantity naturally varies.

/math z = \frac{\hat{\theta} – \theta_0}{\text{SE}}

z_score = (sample_statistic - hypothesised_value) / standard_error
p_value = norm.cdf(z_score)                   # H1: less than
p_value = 1 - norm.cdf(z_score)               # H1: greater than
p_value = 2 * (1 - norm.cdf(abs(z_score)))    # H1: different, either direction

A z-score is a distance measured in standard errors. A z of 2 says the sample landed two standard errors from the claim, which is far enough to be interesting. A z of 0.3 says it landed roughly where the claim expected.

Which tail you take has to match the alternative hypothesis. norm.cdf(z) gives the probability of a standard normal value at or below z, which is what you want when H1 says “less than”. One minus that gives the right tail for “greater than”. And when H1 just says “different”, both directions count as evidence, so you take the tail beyond the absolute value and double it. For a two-sided test at alpha 0.05 the critical value is 1.96, so any absolute z above that is significant.

Pick the direction before you see the data. Running a two-sided test, finding p equal to 0.08, and switching to a one-sided test to halve it is not a different analysis, it is the same analysis with the threshold moved after the fact.

5. One-sample proportion test

Testing a claim about a single rate: are more than 6% of parcels delayed?

p_0 = 0.06                                          # the claim under H0
p_hat = (deliveries["delayed"] == "Yes").mean()     # the observed proportion
n = len(deliveries)
standard_error = np.sqrt(p_0 * (1 - p_0) / n)
z_score = (p_hat - p_0) / standard_error
p_value = 1 - norm.cdf(z_score)                     # H1: the rate is above 6%

The boolean trick in line two is worth internalising. Python treats True as 1 and False as 0, so averaging a comparison gives you the proportion that satisfied it. It works for any condition, not just equality.

The standard error comes from the binomial distribution, which describes how much a proportion bounces around from sample to sample:

Note that it uses the hypothesised proportion rather than the observed one. The whole calculation is performed inside the null hypothesis’s view of the world, so every quantity in it should come from that view.

This test relies on the normal approximation to the binomial, which needs a reasonable number of outcomes on both sides. The usual check is at least 10 successes and 10 failures. With 5 delayed parcels out of 900 the approximation breaks down and you want an exact binomial test instead, from scipy.stats.binomtest.

6. Two-sample proportion test

Comparing rates between two groups.

success_counts = np.array([45, 16])    # delayed parcels in each group
group_sizes = np.array([545, 455])     # total parcels in each group
stat, p_value = proportions_ztest(
    count=success_counts,
    nobs=group_sizes,
    alternative="larger",              # or "smaller", or "two-sided"
)

count holds the successes per group and nobs the totals. alternative="larger" tests whether the first group’s proportion exceeds the second, "smaller" reverses it, and "two-sided" tests only for a difference. The function returns the statistic and the p-value as a tuple.

The manual version shows what it is doing:

p_hat_pooled = success_counts.sum() / group_sizes.sum()
standard_error = np.sqrt(
    p_hat_pooled * (1 - p_hat_pooled) * (1 / group_sizes[0] + 1 / group_sizes[1])
)
p_hat_a, p_hat_b = success_counts / group_sizes
z_score = (p_hat_a - p_hat_b) / standard_error

Under H0 the two groups share a single underlying rate, so the best estimate of it pools all the successes over all the observations:

That pooled rate then feeds the standard error of the difference:

Using the pooled rate rather than the two separate ones is the part worth noticing: the standard error is computed under the assumption H0 is true, because that is the world the p-value describes.

Write it out once for understanding, then use the function. Hand-rolled statistics are a good way to learn and a bad way to ship.

7. T-tests for comparing means

When the outcome is a quantity rather than a rate, you compare means, and the reference distribution changes from normal to t. The t distribution has heavier tails, which accounts for the fact that you are estimating the standard deviation from the same small sample you are testing. As the sample grows, t converges on normal and the distinction stops mattering.

# One-sample: is the mean difference different from zero?
pingouin.ttest(x=depots["pick_time_before"] - depots["pick_time_after"],
               y=0, alternative="two-sided")
# Paired: the same depots measured twice
pingouin.ttest(x=depots["pick_time_before"], y=depots["pick_time_after"],
               paired=True, alternative="greater")

Passing y=0 tests a single sample against a fixed value. Passing paired=True tells the test the two columns are linked measurements on the same units. alternative takes "two-sided", "greater" or "less", matching your H1. Pingouin returns a DataFrame with the t-statistic, degrees of freedom, p-value, confidence interval and Cohen’s d, which is more than most alternatives give you.

Pairing is not a formatting detail, it changes the test. When the same depot is measured before and after, each depot acts as its own control, so the depot-to-depot variation that would otherwise swamp the signal cancels out. Analysing paired data as if it were two independent groups throws that away and can easily turn a clear effect into a null result. If the rows are linked, say so.

For two separate groups:

air = deliveries.loc[deliveries["service_level"] == "Air", "unit_cost"]
road = deliveries.loc[deliveries["service_level"] == "Road", "unit_cost"]
pingouin.ttest(x=air, y=road, alternative="two-sided", correction=True)

.loc[condition, column] pulls each group out as a Series. Without paired=True this runs the independent two-sample test, and correction=True applies Welch’s correction, which does not assume the two groups have equal variance.

The manual version, and a correction worth knowing:

numerator = air.mean() - road.mean()
denominator = np.sqrt(air.var(ddof=1) / len(air) + road.var(ddof=1) / len(road))
t_stat = numerator / denominator

That denominator keeps the two variances separate, which makes it the Welch form. You will often see it paired with degrees_of_freedom = n_a + n_b - 2, and those two do not go together: that formula is the degrees of freedom for the pooled test, which assumes both groups share a variance. Welch’s version uses the Welch-Satterthwaite approximation instead, which is messier and which every library computes for you. Calculating by hand, either pool the variances and use n_a + n_b - 2, or use the separate-variance denominator and let pingouin.ttest(..., correction=True) sort out the degrees of freedom.

For everyday work, use Welch. It costs almost nothing in power when the variances happen to be equal and it stays correct when they are not, which is why it is the better default. Note that scipy.stats.ttest_ind defaults to equal_var=True, so you have to ask for Welch explicitly there.

8. ANOVA for three or more groups

With three groups you might be tempted to run three t-tests. With five groups that becomes ten tests, and the false positive problem takes over.

anova_results = pingouin.anova(
    data=deliveries,
    dv="unit_cost",             # dependent variable, numeric
    between="service_level",    # grouping variable, categorical
)
pairwise_results = pingouin.pairwise_tests(
    data=deliveries,
    dv="unit_cost",
    between="service_level",
    padjust="holm",
)

ANOVA is an omnibus test answering one question: is there any difference among these group means? It reports an F-statistic and a p-value, and a significant result tells you the groups are not all the same and nothing about which ones differ. pairwise_tests then runs a t-test for every pair, with padjust applying a correction for the number of comparisons.

The reason the correction is necessary is arithmetic. Each test at alpha 0.05 has a 5% false positive rate, so the chance of at least one false positive across m independent tests is:

Run twenty tests at 0.05 and that comes to roughly 0.64. You are more likely than not to find something significant even when nothing is going on.

Three corrections are worth knowing. padjust="bonf" multiplies each p-value by the number of comparisons, which is simple, reliable and heavy-handed: it controls false positives well and costs you real effects in the process. padjust="holm" is uniformly more powerful than Bonferroni while controlling the same error rate, so there is little reason to prefer plain Bonferroni once you know it exists. And padjust="fdr_bh" controls the false discovery rate instead, which suits exploratory work with many comparisons where you can tolerate some false positives among your hits.

For correcting a set of p-values that did not come from pairwise_tests, statsmodels has the general version:

from statsmodels.stats.multitest import multipletests
reject, p_adjusted, _, _ = multipletests(p_values, alpha=0.05, method="holm")

That takes any list of p-values, which is what you want when you have run a dozen unrelated tests across a report and need to be honest about it.

On older pingouin versions pairwise_tests is called pairwise_ttests. Same arguments, renamed.

9. Chi-square tests for categorical data

When both variables are categorical, means have no meaning and t-tests do not apply. Chi-square compares counts.

expected, observed, stats = pingouin.chi2_independence(
    data=deliveries, x="contract_type", y="fuel_band"
)
print(stats[stats["test"] == "pearson"])

The test builds a cross-tabulation, works out what counts you would expect in each cell if the two variables were completely unrelated, and measures the gap. Large gaps mean the variables carry information about each other, so a small p-value means they are associated rather than independent.

It returns three objects: the expected counts under independence, the observed counts, and a table of test statistics. Pingouin computes several variants, and filtering to "pearson" gives you the standard chi-square that people mean by the name.

Look at expected before trusting the result, because that is the assumption most often skipped. Chi-square approximates a discrete distribution with a continuous one, and the approximation degrades when expected cell counts are small. The conventional check is that every expected count should be at least 5, or at minimum that no more than 20% of cells fall below it. When cells are sparse, collapse categories together or use Fisher’s exact test.

Goodness of fit asks a different question: does one column match a hypothesised distribution?

observed_counts = deliveries["service_level"].value_counts().sort_index()
hypothesised_props = pd.Series({"Air": 0.20, "Rail": 0.30, "Road": 0.50})
expected_counts = hypothesised_props.sort_index() * len(deliveries)
gof_test = chisquare(f_obs=observed_counts, f_exp=expected_counts)

value_counts() counts the rows in each category and .sort_index() puts them in a known order. The hypothesised proportions get multiplied by the total to become expected counts, and a small p-value means the observed distribution does not match.

Two requirements on those arrays cause more failures than the statistics do. chisquare expects counts rather than proportions, and it raises if the two arrays do not sum to within a rounding tolerance of each other. And they must be in the same category order, which is what sorting both by index guarantees. Skip that and you get the quieter version of the bug, where the totals match but the categories are lined up wrongly and the test silently answers a question you did not ask.

10. Checking assumptions before you test

Every parametric test above rests on assumptions, and checking them takes one line.

counts = deliveries["service_level"].value_counts()
print(counts)
print((counts >= 30).all())

.all() returns True only if every group clears the threshold. Print the counts themselves as well, because a bare False does not tell you which group is short.

The thirty per group figure is a rule of thumb rather than a law. It exists because parametric tests rely on the Central Limit Theorem to treat sample means as approximately normal, and thirty is roughly where that becomes safe for moderately behaved data. Heavily skewed data needs considerably more and near-symmetric data is fine with less.

Proportion tests need a different check entirely. There the requirement is on outcomes rather than rows: at least 10 successes and 10 failures in each group. A sample of 900 with only 4 successes fails badly even though 900 looks like a comfortable number.

The rest, briefly. The t-tests want roughly 30 per group, independent observations and no extreme outliers, and the paired version applies that to the differences rather than the raw values. ANOVA wants the same plus reasonably similar variances across groups. Chi-square wants expected counts of at least 5 in essentially every cell.

When a check fails, the next section is the fallback.

11. Non-parametric alternatives

Non-parametric tests throw away the raw values and work with ranks instead. Sort every observation, replace each with its position in the order, and run the test on those positions. The largest value becomes “rank n” whether it exceeded the next one by 2 or by 2,000, which is exactly why outliers stop dominating and why nothing needs to be normally distributed.

Each parametric test has a counterpart. The one-sample and paired t-tests map to the Wilcoxon signed-rank test, the two-sample t-test maps to Mann-Whitney U, and ANOVA maps to Kruskal-Wallis.

# Wilcoxon signed-rank: paired data
pingouin.wilcoxon(x=depots["pick_time_before"],
                  y=depots["pick_time_after"],
                  alternative="two-sided")
# Mann-Whitney U: two independent groups
on_time = deliveries.loc[deliveries["delayed"] == "No", "parcel_weight_kg"]
delayed = deliveries.loc[deliveries["delayed"] == "Yes", "parcel_weight_kg"]
pingouin.mwu(x=on_time, y=delayed, alternative="two-sided")
# Kruskal-Wallis: three or more groups
pingouin.kruskal(data=deliveries, dv="parcel_weight_kg", between="service_level")

All three take the same alternative argument as their parametric counterparts, and pulling each group out with .loc[condition, column] is cleaner than reshaping the frame to wide format first.

The cost of switching is a modest loss of power when the parametric assumptions did hold. The gain is that the result stays valid when they did not, which is usually the better trade when you have any doubt.

One honest caveat. These tests are often described as testing medians, and that is a simplification. Mann-Whitney tests whether values drawn from one group tend to exceed values drawn from the other, which becomes a statement about medians only when the two distributions have a similar shape. Two groups with identical medians but very different spreads can still return a significant result. Plot the two distributions before you describe what the test found.

12. Confidence intervals

A test returns a verdict. An interval returns a range, and it is almost always the more informative output.

For a difference in means, pingouin already gives it to you:

result = pingouin.ttest(x=air, y=road, correction=True)
print(result[["T", "p-val", "CI95%", "cohen-d"]])

The CI95% column bounds the true difference between the two group means, in pounds per parcel. That is the number to quote.

For proportions, statsmodels has the equivalent:

proportion_confint(count=45, nobs=545, alpha=0.05, method="wilson")

Use method="wilson" rather than the default. The normal approximation interval misbehaves badly when the proportion is near 0 or 1 or the sample is small, sometimes producing bounds below zero, and the Wilson score interval stays sensible throughout.

Interval and test are two views of the same calculation. A 95% interval that excludes the null value corresponds to a two-sided test rejecting at alpha 0.05, and one that includes it corresponds to failing to reject. The interval carries strictly more information, because it also tells you how large the effect might plausibly be.

That extra information changes what a null result means. “Not significant” is ambiguous on its own: it could mean the effect is near zero, or it could mean your sample was too small to see anything. An interval from minus 0.2 to plus 0.3 says the first. An interval from minus 40 to plus 45 says the second. The p-value cannot distinguish them and the interval does it immediately.

This is also the honest answer to the post hoc power question in section 3. If you want to know what a null result rules out, read the interval.

13. Significance is not importance

A p-value answers “is this distinguishable from nothing?” It does not answer “is this worth acting on?”

Those questions come apart as the sample grows. Significance depends on both the size of the effect and the size of the sample, so with enough data any difference that is not exactly zero eventually crosses 0.05. Run an A/B test on ten million sessions and a conversion difference of 0.01 percentage points comes back significant, correct, and completely irrelevant to the business.

The fix is to report an effect size alongside the p-value. Pingouin gives you Cohen’s d in the t-test output without being asked, expressing the difference between two means in standard deviations, so it is comparable across studies and unaffected by sample size. The conventional reading is that 0.2 is small, 0.5 is medium and 0.8 is large, and those labels are rough guides rather than thresholds.

Better still, report the difference in the units people actually use. “Air costs 4.10 pounds more per parcel, 95% interval 3.20 to 5.00” is a sentence someone can make a decision from. “p less than 0.001” is not, however impressive it looks.

14. Common pitfalls

Reading the p-value backwards. It is the probability of the data given the null, never the probability of the null given the data.

Treating “not significant” as “no effect”. It means the evidence was insufficient, which can just as easily mean the sample was too small. Check the confidence interval before concluding anything.

Choosing alpha after seeing the p-value. Setting the threshold to whatever your result cleared is not a decision rule.

Switching to a one-sided test to get under 0.05. The direction of H1 belongs in place before the data arrives.

Skipping the power calculation. An underpowered study cannot produce an informative null, and the effects it does detect are systematically overstated.

Computing power after the fact. Post hoc power is a transformation of the p-value and tells you nothing new. Report the interval instead.

Running many tests without correction. Twenty tests at alpha 0.05 give roughly a 64% chance of at least one false positive. Use padjust="holm", or multipletests for a set of tests from anywhere.

Running pairwise t-tests instead of ANOVA. With more than two groups the omnibus test comes first and the pairwise comparisons follow it with a correction applied.

Analysing paired data as independent. It discards the pairing and throws away most of your power.

Assuming equal variances. scipy.stats.ttest_ind defaults to equal_var=True. Pass equal_var=False, or use pingouin.ttest(..., correction=True).

Mixing the Welch standard error with pooled degrees of freedom. n_a + n_b - 2 belongs with the pooled variance formula, not the separate-variance one.

Skipping the chi-square expected count check. Sparse cells break the approximation. Look at the expected table first.

Passing proportions to chisquare. It wants counts, both arrays must sum to the same total, and they must be in the same category order.

Using the normal approximation for a proportion interval. Near 0 or 1 it produces impossible bounds. Use method="wilson".

Reporting significance without an effect size. With a large enough sample everything is significant and almost nothing matters.

16. The decision tree

BEFORE YOU COLLECT ANYTHING
  how many observations?         ->  solve_power, fix alpha and power at
                                     0.05 and 0.8, then name the smallest
                                     effect worth detecting
WHAT ARE YOU COMPARING?
One group
  a proportion against a claim   ->  one-proportion z-test
  a mean against a claim         ->  one-sample t-test
                                     (Wilcoxon if assumptions fail)
Two groups
  proportions                    ->  two-proportion z-test
  means, same units measured
    twice (paired)               ->  paired t-test
                                     (Wilcoxon signed-rank if assumptions fail)
  means, separate groups         ->  two-sample t-test, Welch
                                     (Mann-Whitney U if assumptions fail)
Three or more groups
  means                          ->  ANOVA, then pairwise tests with padjust
                                     (Kruskal-Wallis if assumptions fail)
Two categorical variables        ->  chi-square test of independence
One categorical variable
  against a known distribution   ->  chi-square goodness of fit
CHECK BEFORE YOU TRUST ANY OF THEM
  counts per group >= 30              for t-tests and ANOVA
  >= 10 successes and failures        for proportion tests
  expected cell counts >= 5           for chi-square
REPORT AFTERWARDS
  the effect size, in the units people use
  the confidence interval, not just the verdict

Three things carry the whole subject. The p-value answers one narrow question, how surprising the data would be if nothing were going on, and every common misreading comes from asking it to answer a different one. Choosing the right test is mostly a matter of what you are comparing and whether its assumptions hold, which the decision tree above settles in a few seconds. And the number worth reporting is almost never the p-value on its own, because an effect size and an interval say how big the thing is and how sure you are, which is what anyone reading it actually needs.

State the hypothesis first, set alpha first, work out the sample size first, then look. Everything after that is arithmetic.

Thanks, Andrei.