Outlier Detection — Demo

Practical walkthrough of three streaming outlier detectors (IQR, MAD, FrugalMAD) from detector.py. Run cells in order.

Show the code
import sys, ossys.path.insert(0, os.path.dirname(os.path.abspath("__file__")))import numpy as npimport matplotlib.pyplot as pltfrom detector import OutlierDetector, OutlierDetectorMAD, OutlierDetectorFrugalMADrng = np.random.default_rng(42)plt.rcParams.update({"figure.figsize": (12, 3.5), "figure.dpi": 100})

1. Generate test dataA sine wave riding on a linear trend, plus Gaussian noise and 8 injected outlier spikes.

Show the code
N = 500t = np.arange(N)clean = 0.5 * np.sin(2 * np.pi * t / 80) + 0.002 * tnoise = 0.15 * rng.standard_normal(N)signal = clean + noise# Inject outliers at known positionsoutlier_idx = np.array([47, 123, 198, 257, 310, 355, 401, 462])outlier_sign = rng.choice([-1, 1], size=len(outlier_idx))signal[outlier_idx] += outlier_sign * rng.uniform(1.5, 3.0, size=len(outlier_idx))plt.plot(t, signal, lw=0.8, label="signal")plt.scatter(outlier_idx, signal[outlier_idx], c="red", zorder=5, label="injected outliers")plt.legend()plt.title("Test signal with injected outliers")plt.tight_layout()plt.show()

2. Compare all three detectorsRun IQR, MAD, and FrugalMAD on the same signal and mark which samples each flags.

Show the code
def run_detector(det, data):
    """Return boolean array of outlier flags."""
    return np.array([bool(det.process(x)) for x in data])

buf = 51
detectors = {
    "IQR (k=1.5)": OutlierDetector(buffer_size=buf, k=1.5),
    # k=4 on the MAD puts the fences in the same place as Tukey's k=1.5 on
    # the IQR (both ~2.7 sigma for Gaussian data); see the topic page.
    "MAD (k=4.0)": OutlierDetectorMAD(buffer_size=buf, k=4.0),
    "FrugalMAD (k=4.0)": OutlierDetectorFrugalMAD(k=4.0),
}

flags = {name: run_detector(det, signal) for name, det in detectors.items()}

# Evaluate each detector
true_set = set(outlier_idx)
print(f"{'Detector':<22s} {'TP':>3s} {'FP':>3s} {'FN':>3s}")
print("-" * 35)
for name, f in flags.items():
    detected = set(np.where(f)[0])
    tp = len(detected & true_set)
    fp = len(detected - true_set)
    fn = len(true_set - detected)
    print(f"{name:<22s} {tp:3d} {fp:3d} {fn:3d}")
Show the code
fig, axes = plt.subplots(3, 1, figsize=(12, 8), sharex=True)for ax, (name, f) in zip(axes, flags.items()):    ax.plot(t, signal, lw=0.6, color="steelblue")    detected = np.where(f)[0]    ax.scatter(detected, signal[detected], c="red", s=30, zorder=5, label="flagged")    ax.scatter(outlier_idx, signal[outlier_idx], facecolors="none",               edgecolors="green", s=80, lw=1.5, zorder=4, label="true outlier")    ax.set_ylabel("x")    ax.set_title(name)    ax.legend(loc="upper left", fontsize=8)axes[-1].set_xlabel("sample index")plt.tight_layout()plt.show()

3. Parameter sensitivity

The threshold sweep that used to live here (detection rate vs false positive rate as k varies) has been promoted to the topic page, where it is discussed properly alongside the fence-scale equivalence between IQR and MAD units, and re-run on every site build:

Which detector, at which threshold? The ROC curve

Keeping a second copy here would only let the two drift apart, so the cell below just reproduces the operating points the page reports, as a cross-check that this notebook and the page agree.

Show the code
# Cross-check: the operating points the topic page reports for this signal.
k_values = np.linspace(0.5, 5.0, 30)
n_out = len(outlier_idx)
n_norm = N - n_out

for name, Det, kw in [("IQR window", OutlierDetector, {"buffer_size": buf}),
                      ("MAD window", OutlierDetectorMAD, {"buffer_size": buf}),
                      ("FrugalMAD", OutlierDetectorFrugalMAD, {})]:
    tpr, fpr = [], []
    for k in k_values:
        det = Det(k=k, **kw)
        found = set(np.where([bool(det.process(x)) for x in signal])[0])
        tpr.append(len(found & true_set) / n_out)
        fpr.append(len(found - true_set) / n_norm)
    tpr, fpr = np.array(tpr), np.array(fpr)
    best = np.argmax(tpr - fpr)
    print(f"{name:>12s}: best at k={k_values[best]:.2f} -> "
          f"TPR {tpr[best]:.2f}, FPR {fpr[best]:.3f} (max TPR {tpr.max():.2f})")

4. Real-world scenario — sensor stream with interferenceSimulate a temperature sensor (slowly drifting around 22 C) with occasional electrical spikes. Show how FrugalMAD tracks the median online using only O(1) memory.

Show the code
# Simulate 1000-sample temperature stream
N_sensor = 1000
t_sensor = np.arange(N_sensor)
temp_base = 22.0 + 0.5 * np.sin(2 * np.pi * t_sensor / 400)  # slow drift
temp_noise = 0.2 * rng.standard_normal(N_sensor)
temp = temp_base + temp_noise

# Inject random spikes (electrical interference)
spike_idx = rng.choice(N_sensor, size=15, replace=False)
temp[spike_idx] += rng.choice([-1, 1], size=15) * rng.uniform(3, 8, size=15)

# Run FrugalMAD and record its internal state over time
det = OutlierDetectorFrugalMAD(k=4.0)
med_trace, mad_trace, flagged = [], [], []
for i, x in enumerate(temp):
    result = det.process(x)
    med_trace.append(det.med)
    # The detector seeds only its median on the first sample, so the scale
    # estimate is None until a deviation has been seen. NaN keeps the trace
    # numeric (and leaves a gap in the plot) instead of poisoning the array.
    mad_trace.append(det.mad if det.mad is not None else np.nan)
    if result:
        flagged.append(i)

med_trace = np.array(med_trace)
mad_trace = np.array(mad_trace)

fig, (ax1, ax2) = plt.subplots(2, 1, figsize=(12, 6), sharex=True)

# Top: signal + median estimate + fences
ax1.plot(t_sensor, temp, lw=0.5, color="steelblue", label="sensor reading")
ax1.plot(t_sensor, med_trace, color="orange", lw=1.5, label="frugal median")
ax1.fill_between(t_sensor, med_trace - 4 * mad_trace, med_trace + 4 * mad_trace,
                  alpha=0.15, color="orange", label="fences (k=4)")
ax1.scatter(flagged, temp[flagged], c="red", s=25, zorder=5, label="flagged")
ax1.set_ylabel("Temperature (C)")
ax1.set_title("FrugalMAD on a simulated temperature sensor stream")
ax1.legend(loc="upper right", fontsize=8)

# Bottom: MAD estimate convergence
ax2.plot(t_sensor, mad_trace, color="green", lw=1)
ax2.axhline(np.median(np.abs(temp_base + temp_noise - np.median(temp_base + temp_noise))),
            ls="--", color="gray", lw=0.8, label="true MAD (no spikes)")
ax2.set_ylabel("MAD estimate")
ax2.set_xlabel("Sample index")
ax2.set_title("MAD convergence over time")
ax2.legend(fontsize=8)

plt.tight_layout()
plt.show()

print(f"Spikes injected: {len(spike_idx)}, flagged: {len(flagged)}")

5. Key takeaways

  • IQR (sliding window): Simple and interpretable. Fences are computed from the buffer before the new sample is appended, so a sample never influences its own test. Blind until the buffer fills. Tukey’s k=1.5 is the standard starting point.
  • MAD (sliding window): A more robust spread estimate than the IQR for asymmetric or heavy-tailed noise, at the cost of noisier fences. Note that k is not on the same scale here: k=4 on the MAD puts the fences where k=1.5 on the IQR puts them.
  • FrugalMAD: O(1) memory, so the obvious choice for a microcontroller, and it has no buffer-fill blind spot. In exchange it lags fast-moving signals, needs a cold-start warm-up, and its scale estimate collapses on a flat stretch (it then goes quiet rather than flooding, and recovers over tens of samples). min_scale is a real, signal-dependent knob you must set for very small-amplitude inputs.
  • Threshold k is the main tuning knob, but compare detectors on a ROC curve rather than at equal k, since k multiplies a different spread measure in each. The topic page works this through and shows that the ranking here is a property of this test signal.