Anomaly Detection Methods for Sensor Networks #

An anomaly detector that cannot tell a fouled inlet from a wildfire plume is worse than no detector at all — it either buries operators in false alarms or silently drops the extreme values that a monitoring network exists to capture. Environmental sensor networks make this harder than textbook outlier detection: baselines drift with the seasons, every channel carries a strong diurnal cycle, hardware fails in ways that look statistically identical to real events, and the interesting signal is often the tail you are tempted to clip. Choosing a method is therefore less about picking the “best” algorithm and more about layering cheap univariate tripwires, a multivariate model, and a spatial gate so that each stage compensates for the blind spots of the others. This is the detection layer of the broader Automated Calibration, Validation & Anomaly Detection pipeline, and it assumes calibrated, quality-flagged input rather than raw telemetry.

The sections below focus on method selection — when a rolling statistic is enough, when you need Isolation Forest, and how PyOD ensembles stabilize a decision boundary that any single detector gets wrong — and on the spatial coherence step that turns raw anomaly scores into an actionable fault-versus-event verdict.

Prerequisites #

Anomaly detection is the last analytical stage before alerting, so it inherits every upstream assumption. Confirm these before you tune a single threshold.

• Python 3.11 with pandas==2.2.2, numpy==1.26.4, scikit-learn==1.4.2, and pyod==1.1.3. Pin PyOD explicitly — its default detector list and scoring conventions have changed across minor releases, and an unpinned upgrade will silently shift your operating point.
• Calibrated, quality-flagged input. Detectors must score corrected values, not raw ADC counts. Run sensor drift correction algorithms first so that a slow baseline migration is not repeatedly flagged as a stream of point anomalies.
• Harmonized units and scale across devices. A multivariate model reads every feature on its own numeric scale; mixing µg/m³ and counts-per-cubic-foot for the same pollutant poisons the joint distribution. Apply cross-device normalization techniques before fitting.
• Minimum feature schema per observation: timestamp (UTC), sensor_id, one or more calibrated numeric channels, latitude, longitude, and the upstream qc_flag (0–3). Detectors should train only on qc_flag == 0 rows so that known-invalid data does not define “normal”.
• A deseasonalization reference. Either a modeled diurnal profile per channel or enough history (14+ days at native cadence) to compute a rolling climatology. Detecting on the raw signal instead of the residual is the single most common cause of diurnal false positives.

Step-by-Step Workflow #

Step 1 — Remove the Predictable Signal #

Almost every environmental channel is dominated by a daily cycle. If you score the raw value, the detector spends its sensitivity budget re-discovering that ozone peaks in the afternoon. Subtract a per-hour climatological profile so the detector sees only the residual — the part the physics did not predict.

# python 3.11 · pandas==2.2.2 · numpy==1.26.4 · scikit-learn==1.4.2 · pyod==1.1.3
from __future__ import annotations
import pandas as pd
import numpy as np


def deseasonalize(df: pd.DataFrame, value_col: str, time_col: str = "timestamp") -> pd.DataFrame:
    """
    Subtract an hour-of-day climatology so detectors score the residual.

    Adds a ``<value_col>_resid`` column: raw minus the median value observed
    in the same hour-of-day bucket across the whole history.
    """
    df = df.copy()
    df[time_col] = pd.to_datetime(df[time_col], utc=True)
    hod = df[time_col].dt.hour
    climatology = df.groupby(hod)[value_col].transform("median")
    df[f"{value_col}_resid"] = df[value_col] - climatology
    return df

Parameters and complexity. Hour-of-day is a coarse but robust profile; for parameters with a weekly traffic signature (NO2 near roads) add day-of-week to the grouping key. Cost is O(n) time and O(24) extra space for the profile. Use median rather than mean so a past event does not inflate the baseline it will later be compared against.

Step 2 — Robust Univariate Baseline (Z-Score and IQR/MAD) #

The rolling Z-score is the cheapest tripwire in the toolbox, but the classic mean/standard-deviation form has a fatal weakness on sensor data: the very outlier you want to catch inflates the window’s own standard deviation and masks itself. The robust fix is the modified Z-score built on the median and the median absolute deviation (MAD), which is not dragged around by a single spike. Pair it with an IQR fence for a second, distribution-free opinion.

def robust_univariate_flags(
    resid: pd.Series,
    window: int = 24,
    mad_threshold: float = 3.5,
    iqr_k: float = 3.0,
) -> pd.DataFrame:
    """
    Flag univariate anomalies on a residual series using a rolling modified
    Z-score (median + MAD) and a rolling IQR fence. Returns a frame with both
    boolean flags so a downstream step can require agreement or either signal.

    Complexity: O(n * window) with pandas rolling; O(window) space.
    """
    med = resid.rolling(window, min_periods=window // 2).median()
    mad = (resid - med).abs().rolling(window, min_periods=window // 2).median()
    # 0.6745 scales MAD to be comparable to a standard deviation for normal data
    mod_z = 0.6745 * (resid - med) / mad.replace(0, np.nan)

    q1 = resid.rolling(window, min_periods=window // 2).quantile(0.25)
    q3 = resid.rolling(window, min_periods=window // 2).quantile(0.75)
    iqr = q3 - q1
    iqr_out = (resid < q1 - iqr_k * iqr) | (resid > q3 + iqr_k * iqr)

    return pd.DataFrame({
        "mad_flag": (mod_z.abs() > mad_threshold).fillna(False),
        "iqr_flag": iqr_out.fillna(False),
    })

Parameters and complexity. Set window to roughly one diurnal cycle in samples (24 for hourly data) even after deseasonalizing, to absorb residual autocorrelation. mad_threshold of 3.5 is the Iglewicz–Hoaglin convention; drop to 3.0 for stable indoor sensors, raise to 4.0 for noisy low-cost gas modules. Time is O(n·window); on multi-million-row backfills push this into a vectorized rolling call rather than a Python loop.

Step 3 — Multivariate Detection with Isolation Forest #

A univariate test can never see that 95 %RH at 38 °C is physically implausible when neither value alone is out of range. Isolation Forest scores the joint feature vector: it randomly partitions the space and measures how few splits it takes to isolate a point, so combinations that sit off the main data manifold get short paths and high scores. It scales near-linearly and needs no distance matrix, which is why it is the default multivariate workhorse here rather than kNN or LOF.

from sklearn.ensemble import IsolationForest


def fit_isolation_forest(
    train: pd.DataFrame,
    feature_cols: list[str],
    contamination: float = 0.01,
    n_estimators: int = 200,
    random_state: int = 42,
) -> IsolationForest:
    """
    Fit Isolation Forest on clean (qc_flag == 0) multivariate telemetry.

    ``contamination`` should match the true anomaly fraction; it sets the
    score threshold, not the model shape. ``max_samples='auto'`` (256) keeps
    trees shallow and training cheap.

    Complexity: ~O(n_estimators * max_samples * log(max_samples)) to fit.
    """
    X = train[feature_cols].dropna()
    clf = IsolationForest(
        n_estimators=n_estimators,
        contamination=contamination,
        max_samples="auto",
        random_state=random_state,
        n_jobs=-1,
    )
    clf.fit(X)
    return clf

Parameters and complexity. Unlike the univariate step, Isolation Forest’s decision boundary moves with contamination, so treat it as a real operating-point knob and set it from historical fault rates, not a round number. n_estimators=200 is comfortably past the point where scores stabilize for a handful of channels. Standardize features first only if you later mix in distance-based detectors — trees are scale-invariant, but the ensemble in the next step is not.

Step 4 — Stabilize the Boundary with a PyOD Ensemble #

Any single detector encodes an assumption: Isolation Forest assumes anomalies are isolable by axis-parallel cuts, a Gaussian test assumes light tails. Real sensor faults violate different assumptions on different days, so a single model’s boundary is brittle. PyOD’s value is a uniform interface across detectors plus a combination layer. Pair Isolation Forest with ECOD (empirical-CDF tail probability, parameter-free) and COPOD (copula-based), which read the tails directly and barely depend on contamination. Combine their normalized scores with average-of-maximum so one confident detector can still raise a flag without three having to agree.

from pyod.models.ecod import ECOD
from pyod.models.copod import COPOD
from pyod.models.iforest import IForest
from pyod.utils.utility import standardizer
from pyod.models.combination import average, maximization


def pyod_ensemble_scores(
    train: pd.DataFrame,
    score: pd.DataFrame,
    feature_cols: list[str],
    contamination: float = 0.01,
) -> np.ndarray:
    """
    Fit an ECOD + COPOD + IForest ensemble on clean data and return a single
    combined anomaly score per row of ``score`` (higher = more anomalous).

    Scores from each detector are standardized to zero mean / unit variance
    before combination so no detector's raw scale dominates.
    """
    # Standardize train and score together so both share one scale.
    X_train, X_score = standardizer(
        train[feature_cols].to_numpy(), score[feature_cols].to_numpy()
    )

    detectors = [ECOD(contamination=contamination),
                 COPOD(contamination=contamination),
                 IForest(contamination=contamination, n_estimators=200, random_state=42)]

    per_detector = np.zeros((X_score.shape[0], len(detectors)))
    for i, det in enumerate(detectors):
        det.fit(X_train)
        per_detector[:, i] = det.decision_function(X_score)

    per_detector = standardizer(per_detector)
    # Blend a robust average with a sensitive maximum so a lone strong
    # detector is not out-voted, but noise in one detector is smoothed.
    return 0.5 * average(per_detector) + 0.5 * maximization(per_detector)

Parameters and complexity. ECOD and COPOD are O(n·d) and parameter-free — ideal edge candidates. Isolation Forest anchors the ensemble on interaction effects the tail-based detectors miss. The 50/50 average–maximum blend is a starting point; shift toward maximization when missing a real event is costlier than a false alarm, and toward average when alert fatigue is the binding constraint.

Step 5 — Gate on Spatial Coherence #

Every method above produces a statistical verdict. None of them can answer the question operators actually care about: is this a broken sensor or a real event? That verdict is spatial. A genuine plume, front, or discharge propagates — it raises correlated readings at neighboring nodes within the wind-advection time — whereas a clogged inlet or a failing ADC spikes one node in isolation. Apply the gate after scoring so it never suppresses the raw detection, only reclassifies it.

def spatial_coherence_gate(
    anomaly: dict,
    neighbors: pd.DataFrame,
    advection_minutes: float,
    resid_col: str = "value_resid",
    z_neighbor: float = 2.0,
) -> str:
    """
    Reclassify a scored anomaly as 'event' or 'fault' by checking whether any
    neighbor's residual also departs from baseline within the advection window.

    ``neighbors`` holds recent residuals for nodes within the correlation
    radius; ``advection_minutes`` = node spacing / prevailing wind speed
    (with a 15-minute floor).
    """
    t0 = pd.Timestamp(anomaly["timestamp"])
    window = neighbors[
        (neighbors["timestamp"] >= t0 - pd.Timedelta(minutes=advection_minutes))
        & (neighbors["timestamp"] <= t0 + pd.Timedelta(minutes=advection_minutes))
    ]
    responded = (window[resid_col].abs() > z_neighbor * window[resid_col].std()).any()
    return "event" if responded else "fault"

Parameters and complexity. The gate is O(k) in the number of neighbors within the correlation radius. Set advection_minutes from node spacing divided by prevailing wind speed with a 15-minute floor; too tight a window classifies slow-moving events as faults. This step is the decisive difference between a detector that cries wolf and one operators trust.

Configuration and Tuning #

Detector choice and thresholds are sensor-specific. These are field-tested starting points, not universal constants.

Validation #

After wiring the stages together, confirm the output is structurally and statistically sane before it feeds alerting.

• Score distribution. Combined ensemble scores should be unimodal with a thin right tail. A bimodal distribution usually means an uncorrected drift is being scored as a persistent second cluster — return to drift correction.
• Flag rate. The fraction of observations flagged should sit within a factor of two of your contamination setting over a representative week. A flag rate an order of magnitude higher signals a non-stationary residual (deseasonalization failed) rather than a fault storm.
• Gate split. Track the event-versus-fault ratio from Step 5. In a healthy maintained network, the large majority of raw anomalies resolve to faults; a sudden surge in “event” classifications is either a real regional episode or a systematic coordinate error to investigate.
• Backtest against known incidents. Replay a labeled fault (a documented inlet clog) and a labeled event (a recorded burn) through the pipeline and confirm each lands in the correct class. This end-to-end check catches integration bugs that per-stage unit tests miss.

Failure Modes and Edge Cases #

Scoring the raw signal instead of the residual. Skipping Step 1 makes the afternoon of every warm day look anomalous. The symptom is a flag rate that oscillates with the diurnal cycle. Always detect on the deseasonalized residual.

Training on contaminated data. Fitting Isolation Forest or a PyOD detector on rows that include past faults teaches the model that faults are normal, collapsing its sensitivity. Train strictly on qc_flag == 0 history and refit after major maintenance.

Contamination set to a convenient default. Leaving contamination at a library default rather than the measured fault rate is the most common tuning error. Isolation Forest’s boundary in particular moves with it, so a wrong value shifts every downstream flag.

Masking by the window itself. A classic mean/standard-deviation Z-score lets a large spike inflate its own window variance and hide. The MAD-based modified Z-score in Step 2 exists specifically to avoid this; do not substitute the naive form.

Suppressing real events at the detector. Applying a spatial or plausibility filter inside the scoring step, rather than as a post-hoc gate, can delete genuine extremes before they are ever seen. Keep detection and the fault/event verdict as separate stages so the raw anomaly is always recorded.

Integration #

The detection layer sits between validation and storage. Its scored, gated output — anomaly_score, mad_flag, iqr_flag, and the event/fault class — travels as columns alongside the calibrated value and its qc_flag. Upstream, it depends on sensor drift correction algorithms so that slow baseline migration is not mistaken for a run of point anomalies, and on cross-device normalization techniques so the multivariate feature space is comparable across a heterogeneous fleet.

For the head-to-head reasoning behind picking a statistical or a model-based detector, see Isolation Forest vs Z-Score for Sensor Anomaly Detection. To see every stage — range checks, QC flags, drift correction, and PyOD scoring — assembled into one runnable recipe, follow Building a Full Air Quality QC Pipeline with pandas and PyOD. The flagged output is the input contract for the full Automated Calibration, Validation & Anomaly Detection pipeline and for downstream spatial storage and export.

FAQ #

Do I still need a Z-score if I run Isolation Forest?

Usually yes. A rolling Z-score catches a sharp single-channel spike within one observation and runs at the edge for pennies of compute, while Isolation Forest needs a fitted model and a full feature vector. Run the cheap univariate test as a first-pass tripwire and reserve the multivariate model for the joint-implausibility cases that no single-channel test can see.

What contamination value should I set for PyOD detectors?

Set contamination to your true expected anomaly fraction, not a convenient default. For a maintained air-quality network it is typically 0.5–2 percent. ECOD and COPOD estimate their own thresholds from the tails and are less sensitive to this value, which is one reason to include them in the ensemble alongside Isolation Forest, whose decision boundary shifts noticeably with contamination.

How do I stop the detector from flagging a real wildfire plume as a fault?

Never let a statistical or model score decide the fault-versus-event question alone. After scoring, apply a spatial coherence gate: a genuine plume raises PM2.5 at several downwind nodes within the advection time, whereas a clogged inlet spikes one node in isolation. Only anomalies with no neighbor response inside the advection window are reclassified as faults; correlated anomalies are surfaced as events.

Why deseasonalize before anomaly detection instead of using a wider window?

A wider rolling window smears the diurnal cycle into the baseline but also blunts the detector’s response to genuine short anomalies. Subtracting a modeled or climatological daily profile removes the predictable component directly, so the residual is closer to stationary and a tight window can stay sensitive. This matters most for parameters with strong daily structure such as temperature, ozone, and NO2.

Which PyOD detector is cheapest to run at the edge?

ECOD is the lightest — it is parameter-free, needs no distance matrix, and scales linearly with rows and columns, so it fits comfortably on a gateway-class device. COPOD is similarly cheap. Reserve kNN or LOF for the cloud tier because their pairwise-distance cost grows quadratically and dominates on long backfills.

Related #

• Automated Calibration, Validation & Anomaly Detection — the full quality pipeline this detection layer plugs into
• Isolation Forest vs Z-Score for Sensor Anomaly Detection — the decision guide for choosing a statistical or model-based detector
• Building a Full Air Quality QC Pipeline with pandas and PyOD — every stage assembled into one runnable PM2.5 recipe
• Sensor Drift Correction Algorithms — the upstream correction that must run before detection
• Cross-Device Normalization Techniques — harmonizing a heterogeneous fleet so the multivariate feature space is comparable