Statistics and signal processing

scipy.stats for distributions, fitting and hypothesis tests; scipy.signal for filter design, zero-phase filtering and peak detection.

scipy.stats: distributions and tests

Every distribution in scipy.stats exposes the same small interface — pdf, cdf, ppf (the inverse CDF), rvs for random variates, and fit — so learning one is learning all of them.

import numpy as np
from scipy import stats

rng = np.random.default_rng(42)
sample = stats.norm.rvs(loc=5.0, scale=2.0, size=200, random_state=rng)

print(sample.mean(), sample.std(ddof=1))
print(stats.norm.cdf(6.0, loc=5.0, scale=2.0))     # P(X <= 6)
print(stats.norm.ppf([0.025, 0.975], 5.0, 2.0))    # central 95% interval

mu, sigma = stats.norm.fit(sample)                 # maximum likelihood
print(mu, sigma)

other = stats.norm.rvs(loc=5.4, scale=2.0, size=200, random_state=rng)
print(stats.ttest_ind(sample, other, equal_var=False))   # Welch's t-test
print(stats.mannwhitneyu(sample, other))                 # non-parametric alternative

Freeze the parameters once: d = stats.norm(5.0, 2.0), then call d.cdf(6.0) instead of repeating them at every call site.
Use stats.ttest_ind(..., equal_var=False) unless you have evidence the variances are equal. Student's t-test is the default, and it is the wrong default for this kind of comparison.
Check assumptions before the test: normality for a t-test, independence always, and similar distribution shapes before reading Mann-Whitney as a difference in medians.
Vectorise: pass a 2-D array with axis=0 and one call tests thousands of features at once.
A p-value is not an effect size. Report the difference in means and a confidence interval next to it.

scipy.signal: filtering and peaks

import numpy as np
from scipy import signal

fs = 500.0
t = np.arange(0.0, 2.0, 1.0 / fs)
x = np.sin(2 * np.pi * 5 * t) + 0.5 * np.sin(2 * np.pi * 60 * t)

sos = signal.butter(4, 20.0, btype="low", fs=fs, output="sos")
y = signal.sosfiltfilt(sos, x)         # zero-phase: no shift, offline only

peaks, props = signal.find_peaks(y, height=0.5, distance=int(fs / 10))
print(t[peaks][:5])
print(props["peak_heights"][:5])

Task	Function	Note
Design a filter	signal.butter, cheby1, ellip	Pass output="sos" rather than the legacy b, a coefficients
Apply to data	signal.sosfilt, sosfiltfilt	filtfilt is zero-phase but doubles the order and needs the whole signal
Find peaks	signal.find_peaks	Set height, distance or prominence or you will find noise
Spectral estimate	signal.welch, signal.periodogram	welch averages segments for a smoother, lower-variance estimate
Convolution or correlation	signal.fftconvolve, correlate	fftconvolve is far faster for long signals
Resample	signal.resample_poly	Polyphase filtering beats naive interpolation

⚠️

Designing a high-order IIR filter with output="ba" produces coefficients that are numerically unstable — round-off moves the poles and the filter rings or diverges. Ask for output="sos" (second-order sections) and apply it with sosfilt or sosfiltfilt.

FAQ

Which filter should I choose?

Butterworth when a flat passband matters, Chebyshev when steeper roll-off matters more than ripple, and a rolling mean or median when the signal is short and the requirement is simple. Design in SOS form either way.

Why does my filtered signal start late?

Any causal filter introduces group delay. sosfiltfilt removes it by filtering forwards and then backwards, at the cost of needing the entire signal and slightly extending transients at both ends.

NumPy arrays Machine learning in one page

Last refreshed 2026-09-18.

Statistics and signal processing

scipy.stats: distributions and tests

scipy.signal: filtering and peaks

FAQ

Related