Skip to course content
Free Python statistics course

Statistical Data Analytics with Python

Unit 01.04: A test, or a framework you can extend

scipy runs the test. statsmodels builds the model. Knowing which you need saves a great deal of fighting with the wrong library.

A test, or a framework you can extend

scipy.stats provides individual tests. ttest_ind, chi2_contingency, pearsonr โ€” each takes arrays and returns a compact result object. It is the right tool when the question is 'is there a difference'.

statsmodels provides models. Its formula API takes a DataFrame and an R-style formula and returns coefficients, standard errors, confidence intervals and diagnostics. It is the right tool the moment you need to adjust for something, because a covariate is one more term in the formula.

For a two-group comparison they answer the same question, and the statsmodels version is the one you can grow.

This block runs the same comparison both ways.

import numpy as np
from scipy import stats
import statsmodels.api as sm
import statsmodels.formula.api as smf
import pandas as pd

rng = np.random.default_rng(11)
d = pd.DataFrame({
    "site": ["comparison"] * 40 + ["pilot"] * 40,
    "wellbeing": np.concatenate([rng.normal(58, 9, 40), rng.normal(64, 9, 40)]),
})

# scipy.stats: a test, returning a small result object.
a = d.loc[d.site == "comparison", "wellbeing"]
b = d.loc[d.site == "pilot", "wellbeing"]
t_res = stats.ttest_ind(b, a, equal_var=False)
print("scipy.stats.ttest_ind")
print(f"  t = {t_res.statistic:.3f}   p = {t_res.pvalue:.5f}")
ci = t_res.confidence_interval()
print(f"  95% CI for the difference: [{ci.low:.2f}, {ci.high:.2f}]")
print()

# statsmodels: a model, returning coefficients, intervals and diagnostics.
model = smf.ols("wellbeing ~ site", data=d).fit()
print("statsmodels.formula.api.ols")
print(model.summary2().tables[1].round(4).to_string())
print()
print("The same comparison, two tools. scipy answers 'is there a difference';")
print("statsmodels puts it in a framework you can add covariates to.")
print(f"R-squared: {model.rsquared:.4f}   n = {int(model.nobs)}")

scipy gives t = 4.449, p = 0.00003 and a 95% interval for the difference of [4.30, 11.26]. statsmodels reports the same comparison as a coefficient: site[T.pilot] at 7.7779 with the identical interval of [4.2974, 11.2583] and the same t statistic. R-squared is 0.2024 on 80 observations. Same arithmetic, two presentations โ€” and only the second accepts + department as the next term.

The mistake this prevents

The mistake is reaching for scipy when the analysis will need covariates. You end up hand-rolling adjustments that statsmodels does in one formula.

Takeaway

Use scipy for a single self-contained test and statsmodels when the analysis is a model or will grow covariates. Both give intervals โ€” read them before the p-value.