Skip to content

SciPy

The scientific toolbox built on NumPy: probability and statistics, optimization, sparse matrices, interpolation, integration, signal filtering, and more — everything a data scientist needs beyond arrays and tables.

Table of Contents

Installation / Setup

pip install scipy          # or: uv add scipy

Requires NumPy. Import the submodules you need:

import scipy.stats as stats
import scipy.optimize as opt
import scipy.linalg as la
import scipy.sparse as sp
import scipy.interpolate as interp
import scipy.integrate as integrate
import scipy.signal as signal
import scipy.ndimage as ndimage
import scipy.cluster as cluster
import scipy.spatial as spatial

Core Syntax

Statistics (scipy.stats)

Distributions (frozen objects)

norm = stats.norm(loc=0.0, scale=1.0)
norm.pdf(0.0)        # probability density
norm.cdf(1.96)       # cumulative: P(X <= 1.96)
norm.ppf(0.975)      # quantile / inverse CDF (z for 95% CI)
norm.rvs(size=100, random_state=42)   # sample

# Common distributions: norm, t, chi2, f, uniform, poisson, binom,
# exponential, lognorm, gamma, beta, weibull_min, pareto

Hypothesis tests

Test Function Use
Two-sample t-test stats.ttest_ind(a, b) Compare two groups
Paired t-test stats.ttest_rel(a, b) Before/after
One-sample t-test stats.ttest_1samp(a, mu) Mean vs constant
Mann–Whitney U stats.mannwhitneyu(a, b) Non-parametric groups
Wilcoxon stats.wilcoxon(a, b) Non-parametric paired
Pearson / Spearman stats.pearsonr(a, b), stats.spearmanr(a, b) Correlation
Chi-squared stats.chi2_contingency(table) Independence in a table
ANOVA stats.f_oneway(a, b, c) 3+ group means
Normality stats.shapiro(x), stats.normaltest(x) Check Gaussian
KS test stats.kstest(x, "norm") Fit to a distribution 
stat, p = stats.ttest_ind(control, treatment)
print(f"t={stat:.3f}, p={p:.4f}")     # p < 0.05 -> significant difference

Descriptive stats

stats.describe(x)                 # count, mean, min, max, variance, skew, kurtosis
stats.zscore(x)                   # standardize
stats.sem(x)                      # standard error of the mean
stats.linregress(x, y)            # slope, intercept, rvalue, pvalue, stderr

Optimization (scipy.optimize)

# Minimize a scalar/multivariate function
res = opt.minimize(lambda v: (v[0]-3)**2 + (v[1]+2)**2, x0=[0, 0])
res.x                       # array([ 3., -2.])
res.fun                     # objective at optimum
res.success, res.message

# Fit a curve to data
def decay(t, a, k):
    return a * np.exp(-k * t)

params, cov = opt.curve_fit(decay, t, y_obs, p0=[1, 1])
perr = np.sqrt(np.diag(cov))            # parameter errors

# Root finding & least squares
root = opt.root(lambda x: x**2 - 4, x0=1.5)
res = opt.least_squares(lambda p: model(p, x) - y, x0=[1, 1])

For linear least squares prefer np.linalg.lstsq (see NumPy).

Linear algebra (scipy.linalg)

la.inv(A); la.det(A); la.norm(v)
la.svd(A, full_matrices=False)     # U, s, Vt
la.eig(A); la.eigvals(A)
la.cholesky(A); la.lu(A); la.qr(A)
la.solve(A, b)
la.expm(A)                         # matrix exponential (Markov chains!)
la.block_diag(A, B); la.kron(A, B)

SciPy's linalg wraps LAPACK more completely than NumPy's (e.g. expm, cholesky flavors, lu).

Sparse matrices (scipy.sparse)

# Formats: csr (row), csc (column), coo (coordinates), lil, dok
M = sp.csr_matrix((data, (rows, cols)), shape=(1000, 1000))
M = sp.eye(1000, format="csr")
M.toarray()                          # dense (only for small!)
x = sp.linalg.spsolve(M, b)          # efficient solve
M @ v                                # sparse matmul stays sparse/fast

Sparse is the right tool for one-hot encodings, document-term matrices, graphs, and embeddings where >90% of entries are zero.

Interpolation (scipy.interpolate)

f = interp.interp1d(x, y, kind="linear")   # or kind="cubic"
f(2.5)                                     # interpolated value
cs = interp.CubicSpline(x, y)
grid = interp.griddata(points, values, xi, method="linear")  # scattered -> grid

Integration (scipy.integrate)

integrate.quad(np.sin, 0, np.pi)                # -> (2.0, error_estimate)
integrate.simpson(y, x=x)                       # numeric integral of samples
sol = integrate.solve_ivp(ode_fn, [t0, t1], y0) # ODE initial-value problems

Signal processing (scipy.signal)

y_smooth = signal.savgol_filter(y, window_length=11, polyorder=3)  # smooth curve
y_smooth = signal.medfilt(y, kernel_size=5)                        # denoise spikes
resampled = signal.resample(y, num=1000)                           # change length
f, t, Sxx = signal.spectrogram(x, fs=100)                          # time-frequency
b, a = signal.butter(4, 0.1, btype="low")        # design filter
y_filt = signal.filtfilt(b, a, y)                # zero-phase filtering

Image processing (scipy.ndimage)

blur = ndimage.gaussian_filter(img, sigma=2)
zoomed = ndimage.zoom(img, 2)                    # 2x upscale
rotated = ndimage.rotate(img, 45)
labeled, n = ndimage.label(mask)                 # connected components
centers = ndimage.center_of_mass(mask, labeled, range(1, n+1))

Clustering & spatial

centroids, labels = cluster.vq.kmeans2(X, k=5)   # k-means (fast, simple)
Z = cluster.hierarchy.linkage(X, method="ward")  # for dendrograms
dendro = cluster.hierarchy.dendrogram(Z)

spatial.distance.euclidean(a, b)
D = spatial.distance.pdist(X)                    # condensed pairwise distances
Dm = spatial.distance.cdist(A, B)                # cross distances (n x m)
tree = spatial.KDTree(X)
dist, idx = tree.query(point, k=5)               # fast nearest neighbors

Practical Examples

1. A/B test decision

import numpy as np
import scipy.stats as stats

rng = np.random.default_rng(7)
control = rng.normal(50, 10, 200)
variant = rng.normal(53, 10, 200)

stat, p = stats.ttest_ind(control, variant)
ci = stats.norm.interval(0.95,
        loc=(variant.mean() - control.mean()),
        scale=np.sqrt(control.var()/200 + variant.var()/200))
print(f"delta mean = {variant.mean()-control.mean():.2f}, p = {p:.4f}, 95% CI {ci}")

2. Fit an exponential decay

t = np.linspace(0, 5, 50)
y_true = 3.0 * np.exp(-1.2 * t)
y_obs = y_true + rng.normal(0, 0.1, 50)

(a, k), _ = opt.curve_fit(decay, t, y_obs, p0=[1, 1])
# a ~ 3.0, k ~ 1.2

3. Sparse one-hot matrix for modeling

from sklearn.preprocessing import OneHotEncoder
import scipy.sparse as sp

X_ohe = OneHotEncoder().fit_transform(categories_df)   # already sparse
# X_ohe: (n_samples, n_categories) in csr format — feed directly to solvers

4. Nearest-neighbor lookup with KDTree

refs = rng.normal(size=(10_000, 3))
queries = rng.normal(size=(100, 3))
tree = spatial.KDTree(refs)
d, i = tree.query(queries, k=1)    # distance + index of nearest ref

Common Patterns

# z-score significance helper
def report(p, alpha=0.05):
    return "significant" if p < alpha else "not significant"

# 95% CI of the mean
ci = stats.t.interval(0.95, df=len(x)-1, loc=x.mean(), scale=stats.sem(x))

# Smooth-then-differentiate curves
y_s = signal.savgol_filter(y, 11, 3)
grad = np.gradient(y_s)

# Fast pairwise distances
D = spatial.distance.squareform(spatial.distance.pdist(X))   # full n x n matrix

Best Practices

  • Frozen distributions (norm(loc, scale)) are cleaner than passing args to every call.
  • Report effect sizes and CIs along with p-values — a tiny p-value can hide a trivial effect.
  • Check test assumptions (normality for t-tests; use Mann–Whitney/Wilcoxon otherwise).
  • Keep matrices sparse (csr/csc) whenever >90% zeros; never call .toarray() on large ones.
  • Use savgol_filter or filtfilt for smooth, phase-preserving signal cleaning.
  • Always pass random_state/seed to rvs for reproducibility.
  • Prefer optimize.curve_fit with explicit p0 and inspect the covariance for fit quality.

Common Mistakes

Mistake Consequence Correct approach
Testing many pairs without correction Inflated false positives Use statsmodels multiple-comparison methods or Bonferroni
chi2_contingency on raw counts vs proportions Wrong p-value Pass the contingency count table
.toarray() on huge sparse matrix MemoryError Keep sparse; use spsolve
Interpolating beyond data range Wild extrapolation Use bounds_error=True or clip
Applying t-test to non-normal data Invalid inference mannwhitneyu / wilcoxon
Ignoring res.success from minimize Trusting a failed optimum Check res.success and res.message

Troubleshooting

Symptom Cause Fix
OptimizeWarning: Covariance of the parameters could not be estimated Poor fit / p0 Better p0; more data; check model
ValueError: operands could not be broadcast Wrong shapes for test inputs Same length; drop NaN first
NaN p-values NaN in data x = x[~np.isnan(x)] or nan_policy="omit"
SparseEfficiencyWarning Operation better in another format Convert: M.tocsr() / .tocsc()
solve_ivp not converging Stiff ODE method="LSODA" or method="Radau"
ImportError for a submodule Wrong name scipy.signal, not scipy.signals; check dir(scipy)

Official Documentation

Last Updated

2026-09-06 — Initial version.