4. Statistical anomaly detection methods for equipment monitoring

Statistical anomaly detection treats a machine's healthy behavior as a probability distribution and flags observations that are unlikely under that distribution. For industrial equipment monitoring, this framing is powerful because it produces a continuous anomaly score rather than a binary alarm, and because the score has a clear interpretation: how many standard deviations, or how many units of Mahalanobis distance, separate the current observation from normal operation. This chapter covers the four methods that dominate practice — Gaussian models, Principal Component Analysis (PCA), robust covariance estimation, and threshold setting under non-stationary operating conditions — and shows how they connect into a single monitoring pipeline.

Gaussian Models and the Mahalanobis Distance

The simplest statistical detector assumes that a feature vector x (for example, RMS vibration on three bearings, motor current, and oil temperature) follows a multivariate normal distribution N(μ, Σ) during healthy operation. The anomaly score is the squared Mahalanobis distance:

d²(x) = (x − μ)ᵀ Σ⁻¹ (x − μ)

Under the Gaussian assumption, d² follows a chi-square distribution with p degrees of freedom, where p is the number of features. This gives a principled threshold: for a false-alarm rate of 1%, set the threshold at the 99th percentile of χ²(p). For p = 5, that value is approximately 15.09. An observation with d² > 15.09 is flagged.

The critical weakness is that Σ⁻¹ amplifies noise in directions where normal data has little variance. If two features are nearly collinear — say, two accelerometers mounted on the same bearing housing — Σ becomes ill-conditioned, and tiny measurement noise produces enormous d² values. Always check the condition number of Σ before trusting a Gaussian detector; values above roughly 10⁴ signal that regularization or dimensionality reduction is needed.

PCA-Based Detection

PCA addresses collinearity by projecting data onto a lower-dimensional subspace that captures the dominant variance. Given a centered training matrix X, compute the eigendecomposition of the covariance matrix and retain the top k eigenvectors as the principal subspace P. The remaining eigenvectors form the residual subspace P̃. Two statistics are then monitored:

  • T² statistic — variation within the principal subspace: T² = xᵀ P Λ⁻¹ Pᵀ x, where Λ holds the retained eigenvalues. This detects shifts that stay inside normal operating patterns but move too far along them.
  • SPE (Squared Prediction Error) or Q statistic — variation outside the principal subspace: SPE = ‖(I − PPᵀ)x‖². This detects new correlation structures that were absent during training, which is exactly what a developing fault looks like.

The choice of k controls the trade-off. Retaining too few components pushes normal variation into the residual and inflates false alarms; retaining too many hides faults inside the principal subspace. Cumulative variance thresholds (90–95%) are a common starting point, but cross-validation against known fault episodes is more reliable. A practical rule: if the SPE statistic fires constantly during known-healthy periods, increase k; if faults are missed, decrease it.

Robust Covariance Estimation

Both Gaussian and PCA detectors estimate μ and Σ from training data using sample means and sample covariance. A single contaminated training window — a brief overload, a sensor glitch, a maintenance event — shifts the mean and inflates the covariance, desensitizing the detector permanently. Robust estimators solve this by down-weighting outliers during fitting.

The Minimum Covariance Determinant (MCD) estimator is the standard choice. It searches for the subset of h observations (typically h ≈ 0.75n) whose covariance matrix has the smallest determinant, then rescales that covariance to be consistent for the full dataset. In Python, sklearn.covariance.MinCovDet implements this directly:

from sklearn.covariance import MinCovDet
mcd = MinCovDet(support_fraction=0.75, random_state=0).fit(X_train)
mu_robust = mcd.location_
sigma_robust = mcd.covariance_

Compare the condition number and the diagonal entries of sigma_robust against the sample covariance. If the robust estimate is dramatically smaller, the training set contained outliers that would have masked real faults. Robust covariance is not optional in industrial settings — it is the difference between a detector that works on the bench and one that survives a year on the plant floor.

Threshold Setting Under Non-Stationary Conditions

Industrial equipment rarely operates at a single setpoint. Load, ambient temperature, and speed vary continuously, so the "normal" distribution drifts. A fixed threshold calibrated on one operating condition will either miss faults at other conditions or alarm constantly. Three techniques handle this:

  1. Operating-condition segmentation. Partition historical data by regime (e.g., load bins, speed bins) and fit a separate model per regime. At runtime, classify the current regime and apply the corresponding threshold. This is simple and interpretable but requires enough data per regime.
  2. Exponentially weighted moving statistics. Update μ and Σ online with a forgetting factor λ (typically 0.95–0.99): μ_t = λμ_{t−1} + (1−λ)x_t. The detector adapts to slow drift while still responding to abrupt faults. The trade-off is that a slow-developing fault can be absorbed into the baseline; pair EWMA with a separate rate-of-change monitor.
  3. Quantile-based thresholds on the score. Instead of assuming a chi-square distribution, compute the empirical 99th percentile of the anomaly score over a rolling window of recent healthy data. This is distribution-free and naturally adapts to non-stationarity, at the cost of needing a clean reference window.
Do not combine EWMA adaptation with a fixed chi-square threshold. The chi-square threshold assumes a stationary distribution; if the mean is being updated online, the theoretical false-alarm rate no longer holds. Use empirical quantiles when adapting online.

Putting It Together: A Monitoring Pipeline

The figure below shows how these components connect. Raw sensor data flows through feature extraction and scaling, then into a robust covariance model that produces a Mahalanobis score. In parallel, a PCA model produces T² and SPE scores. All three scores pass through a regime-aware threshold stage that selects the appropriate cutoff based on the current operating condition, and the final decision is a logical OR of the individual flags.

Sensor Inputs Vibration RMS Current, Temp, Speed Feature Extraction Windowing, FFT bands Standardization Robust Covariance MCD → Mahalanobis d² PCA Model T² and SPE scores Anomaly Scores d², T², SPE (continuous) Regime-Aware Threshold Load/speed bin → cutoff Decision OR of flags Alarm / Normal Operating Regime Load, speed, ambient EWMA baseline update (λ = 0.97) Score interpretation: • d² > χ²₀.₉₉(p) → deviation from robust normal distribution • T² > threshold → abnormal magnitude within normal correlation structure • SPE > threshold → new correlation pattern (developing fault signature)

The pipeline is deliberately modular. If the robust covariance model produces too many false alarms during a particular load regime, you can adjust that regime's threshold without retraining the PCA model. If a new fault type appears that the PCA residual misses, you can add a dedicated feature without disturbing the Mahalanobis detector. This separation of concerns is what makes statistical monitoring maintainable over the equipment's service life.

When commissioning a new detector, run it in shadow mode for at least two weeks: compute scores and log decisions, but do not connect them to alarms. Compare the flagged timestamps against maintenance records. A detector that fires during known-good periods needs its threshold raised or its training data cleaned; one that stays silent through a known fault needs more features or a lower threshold.

The methods in this chapter share a common assumption: that healthy behavior can be characterized by a distribution estimated from historical data. When that assumption breaks — during a gearbox rebuild, a sensor replacement, or a process change — the model must be retrained. Building a retraining trigger into the monitoring system, based on sustained score elevation without a corresponding maintenance event, prevents the slow drift of a stale model into irrelevance.