← Back to Autonomy

Multivariate monitoring · system health management

Reading the Whole Machine

A field guide to the methods that watch many sensors at once: is the system still in its normal pattern, which channel broke it, and which fault is it.

Single-signal limits catch the obvious. The hard faults hide in the relationships between sensors: every channel stays inside its band while the pattern that ties them together quietly breaks. Multivariate monitoring is the family of methods that watch that joint pattern. They split cleanly by what you are given and what you are asking:

Question 1
Still normal?

Does a new sample fit the healthy multi-sensor pattern, and if not, which channel is responsible?

Question 2
Related how?

How does one block of variables, the commands, relate to another, the responses?

Question 3
Which fault?

Given labeled history, separate or classify the known fault modes.

Each card gives the idea in plain words, a system-health example, a figure, and runnable Python. Two live explorers let you inject a fault, and a worked case study runs the whole stack on one record. A decision map keys every method back to the question it answers.

Interactive explorers

The figures above are frozen. These two move. Inject a fault and watch the monitoring statistics respond, or reshape a covariance and watch the outlier metric react.

A

PCA monitoring: inject a fault

T² / Q live

A PCA model is fit on healthy data. Drift one sensor and watch T² and Q cross their limits, then read the contribution bar to see which channel is blamed.

B

Mahalanobis vs. Euclidean

covariance-aware

Reshape the healthy correlation and drag the test point. Euclidean distance ignores the cloud's shape; Mahalanobis does not.

Euclidean0.0
Mahalanobis0.0

Subspace monitoring, the workhorse

Build one model of healthy behaviour, then ask of every new sample: is it still inside the normal pattern, and if not, which sensor broke it. This is what runs continuously on a fleet.

01

PCA with Hotelling's T² and the SPE / Q statistic

unsupervised

The classic process-monitoring pair. T² watches variation inside the model subspace; Q watches the residual subspace. Contribution plots then name the offending variable.

In plain words

Principal component analysis learns the few directions along which a healthy multi-sensor system actually varies. Two numbers then summarize any new sample. Hotelling's T² measures how far it sits inside that learned subspace, so it catches unusual but in-pattern operation. The squared prediction error, called SPE or Q, measures how far it falls outside the subspace, so it catches a brand-new pattern the model never saw, such as two sensors that used to track each other drifting apart. When either crosses its control limit, a contribution plot breaks the statistic back down onto the raw sensors and points at the culprit.

Health-management example

A converter has a dozen correlated temperatures, currents, and voltages. PCA on a month of healthy data keeps three components. A developing solder crack slowly decouples one thermistor from the rest: T² stays calm because each reading is individually normal, but Q climbs past its limit because the joint correlation is broken. The Q contribution plot puts almost all the blame on that one thermistor channel.

2026-06-19T20:05:46.239655 image/svg+xml Matplotlib v3.10.8, https://matplotlib.org/
import numpy as np
from scipy import stats

class PCAMonitor:
    def fit(self, X, var=0.9, alpha=0.99):
        self.mu = X.mean(0); self.sd = X.std(0, ddof=1); self.sd[self.sd==0]=1
        Z = (X - self.mu) / self.sd
        n, m = Z.shape
        _, s, Vt = np.linalg.svd(Z, full_matrices=False)
        eig = s**2 / (n - 1)
        k = int(np.searchsorted(np.cumsum(eig)/eig.sum(), var) + 1)
        self.P, self.eig, self.res = Vt[:k].T, eig[:k], eig[k:]
        self.k = k
        self.t2_lim = k*(n-1)*(n+1)/(n*(n-k))*stats.f.ppf(alpha, k, n-k)
        t1,t2,t3 = self.res.sum(), (self.res**2).sum(), (self.res**3).sum()
        h0 = 1 - 2*t1*t3/(3*t2**2); ca = stats.norm.ppf(alpha)
        self.q_lim = t1*(ca*np.sqrt(2*t2*h0**2)/t1 + 1 + t2*h0*(h0-1)/t1**2)**(1/h0)
        return self
    def _z(self, X): return (np.atleast_2d(X) - self.mu) / self.sd
    def t2(self, X): T = self._z(X) @ self.P; return (T**2 / self.eig).sum(1)
    def q(self, X):
        Z = self._z(X); E = Z - Z @ self.P @ self.P.T; return (E**2).sum(1)
    def q_contrib(self, x):
        z = self._z(x); e = z - z @ self.P @ self.P.T; return (e[0]**2)

mon = PCAMonitor().fit(X_healthy)
flag = (mon.t2(X_new) > mon.t2_lim) | (mon.q(X_new) > mon.q_lim)
culprit = np.argmax(mon.q_contrib(X_new[flag][-1]))  # which sensor
Practical note

Always scale to unit variance first, or the largest-range sensor dominates the components. Choose the rank from a scree plot or cumulative variance, not reflexively at 95%. Recompute limits on a held-out healthy set, since in-sample limits are optimistic. Watch T² and Q together: T² alone misses new patterns, Q alone misses extreme but in-pattern excursions. Contribution plots indicate, they do not prove, since a fault smears across correlated sensors.

Jackson & Mudholkar (1979), Technometrics 21; Wise & Gallagher (1996), J. Process Control 6; Qin (2003), J. Chemometrics 17.

02

Mahalanobis distance

unsupervised

A covariance-aware outlier metric: the multivariate generalization of the z-score. It rescales by the healthy covariance, so it flags points that violate the joint correlation even when each axis looks normal.

In plain words

A plain z-score asks how many standard deviations a value sits from the mean, one variable at a time. That misses faults that keep every single sensor in its normal range while breaking the relationship between them. Mahalanobis distance fixes this by dividing out the full covariance of the healthy data. Geometrically it stretches and rotates space so the healthy cloud becomes a sphere; distance is then measured in that whitened space. A point sitting off the main correlation axis is far in Mahalanobis terms even if it is close in ordinary Euclidean terms. Under a Gaussian assumption the squared distance follows a chi-square law, which gives a clean control limit.

Health-management example

Two coolant sensors normally rise and fall together. After a fault one reads high while the other reads low, both still within their individual operating bands. Euclidean distance from the mean barely moves; Mahalanobis distance jumps well past its chi-square limit because the pair has left the diagonal ridge of normal joint behaviour.

2026-06-19T20:05:46.928776 image/svg+xml Matplotlib v3.10.8, https://matplotlib.org/
import numpy as np
from scipy import stats

def mahalanobis(X_ref, x):
    mu = X_ref.mean(0)
    S = np.cov(X_ref, rowvar=False)
    Sinv = np.linalg.pinv(S)
    d = np.atleast_2d(x) - mu
    dist = np.sqrt(np.einsum("ij,jk,ik->i", d, Sinv, d))
    limit = np.sqrt(stats.chi2.ppf(0.99, df=X_ref.shape[1]))
    return dist, limit

dist, limit = mahalanobis(X_healthy, X_new)
alarms = dist > limit
Practical note

The covariance must be estimated on genuinely healthy, representative data, and you need roughly ten or more samples per variable or the inverse is unstable; regularize, shrink, or use a robust covariance (MCD) for heavy tails and outliers in the reference set. The chi-square limit assumes approximate multivariate normality, so verify it or set the threshold empirically from a clean validation set.

Mahalanobis (1936), Proc. Natl. Inst. Sci. India 2; Rousseeuw & Van Zomeren (1990), JASA 85 (robust covariance).

03

Independent component analysis (ICA)

source separation

Separates statistically independent sources from their mixtures. The right tool when a fault arrives mixed together with other signals at a shared sensor.

In plain words

PCA finds uncorrelated, orthogonal directions ranked by variance. ICA instead finds directions that are statistically independent, which need not be orthogonal, by pushing the recovered signals as far from Gaussian as possible. The payoff is source separation: when one sensor records a sum of several physical processes, ICA can pull them apart into their underlying generators. In monitoring it is used both to unmix overlapping fault signatures and as a feature step before the same T² and Q machinery, since a fault often lives in one independent component rather than one raw channel.

Health-management example

A single accelerometer on a gearbox hears two tones at once: a healthy meshing vibration and a faint bearing defect tone. They overlap in the raw spectrum. ICA across a few accelerometers separates the two independent sources, and the bearing component, once isolated, shows the rising amplitude that the mixed signal hid.

2026-06-19T20:05:46.832246 image/svg+xml Matplotlib v3.10.8, https://matplotlib.org/
import numpy as np
from sklearn.decomposition import FastICA

ica = FastICA(n_components=3, whiten="unit-variance",
              max_iter=1000, random_state=0).fit(X_obs)
sources = ica.transform(X_obs)   # recovered independent signals
mixing  = ica.mixing_            # how sources combine into sensors
# monitor the most fault-sensitive source as its own signal
track = sources[:, 0]
Practical note

ICA needs non-Gaussian sources; it cannot separate Gaussian mixtures, and it fixes neither the order nor the sign or scale of the recovered sources, so label them by inspection or by correlation to known references. Whiten first (FastICA does this internally), and never feed it more components than sensors. It is sensitive to strong noise, so denoise or band-limit before unmixing.

Hyvarinen & Oja (2000), Neural Networks 13 (FastICA); Lee et al. (2004), J. Process Control 14 (ICA monitoring).

Relating two blocks of variables

Sometimes the question is not whether one cloud is normal, but how one set of variables relates to another: commands to responses, inputs to outputs.

04

Canonical correlation analysis (CCA)

two-block

Finds the linear combination of one variable block that correlates most strongly with a combination of another. Natural for relating an actuator or command block to a response block.

In plain words

Given two groups of variables measured on the same samples, CCA looks for one weighted mixture from each group such that the two mixtures correlate as highly as possible, then a second orthogonal pair, and so on. Each pair is a canonical mode, and its canonical correlation says how tightly that input combination is tied to that output combination. It generalizes ordinary correlation from two scalars to two vectors, and it answers a question single correlations cannot: which combined drive mode explains which combined response mode.

Health-management example

On the command side you have throttle, gear, and steering; on the response side, three accelerations. CCA reveals that one weighted blend of commands explains almost all of one blend of responses, with canonical correlation near one. A new operating regime where that canonical correlation collapses signals that the command-to-response coupling has changed, often an early sign of an actuator or linkage fault.

2026-06-19T20:05:46.692869 image/svg+xml Matplotlib v3.10.8, https://matplotlib.org/
import numpy as np
from sklearn.cross_decomposition import CCA

cca = CCA(n_components=2, scale=True).fit(X_cmd, Y_resp)
U, V = cca.transform(X_cmd, Y_resp)        # canonical variates
corrs = [np.corrcoef(U[:, i], V[:, i])[0, 1] for i in range(2)]
# a drop in the leading canonical correlation flags decoupling
print("canonical correlations:", np.round(corrs, 2))
Practical note

Canonical correlations are optimistically high when either block has many variables relative to the sample count, so regularize (regularized or sparse CCA) or reduce dimension first, and validate the correlations on held-out data. Scale both blocks. Interpret the weight vectors with care, since correlated variables within a block share and trade weight.

Hotelling (1936), Biometrika 28; Hardoon, Szedmak & Shawe-Taylor (2004), Neural Computation 16.

Supervised, with a fault label in hand

When historical data is labeled by quality outcome or fault class, supervised projections beat unsupervised ones: they aim the model at the distinction you actually care about.

05

PLS and discriminant PLS (PLS-DA)

supervised

The supervised counterpart to PCA. It builds latent directions that covary with a quality or fault label, so it separates classes the unsupervised subspace may blur.

In plain words

PCA chooses directions by variance alone, indifferent to any label. Partial least squares instead chooses latent directions that maximize covariance between the inputs and a target. With a continuous target it is a regression; with a class label coded as an indicator it becomes discriminant PLS, a classifier. Because it leans toward the label, PLS often needs fewer components than PCA to expose a fault, and its VIP score, variable importance in projection, ranks how much each sensor contributes, which doubles as fault localization.

Health-management example

You have runs labeled healthy or degraded. The degradation is spread thinly across several sensors, so it hides below the dominant load variation that PCA's first components capture. PLS-DA finds the latent direction that tracks the label and cleanly separates the two groups, and its VIP scores flag the three sensors carrying the fault signature.

2026-06-19T20:05:46.550291 image/svg+xml Matplotlib v3.10.8, https://matplotlib.org/
import numpy as np
from sklearn.cross_decomposition import PLSRegression

Y = (y == 1).astype(float)[:, None]                 # fault = 1
pls = PLSRegression(n_components=2, scale=True).fit(X, Y)
T = pls.x_scores_                                   # latent scores
pred = (pls.predict(X)[:, 0] > 0.5).astype(int)
# VIP score per variable (importance in the projection)
W, Q = pls.x_weights_, pls.y_loadings_
ss = (T**2).sum(0) * (Q**2).sum(0); p = X.shape[1]
vip = np.sqrt(p * ((ss * (W/np.linalg.norm(W,axis=0))**2).sum(1)) / ss.sum())
Practical note

Pick the number of latent variables by cross-validation, not by eye, and beware that PLS-DA on a label will manufacture apparent separation even from noise if components are overfit, so always validate on held-out runs and, with imbalanced classes, on a permutation test. VIP above one is the usual relevance cutoff but is a heuristic. Scale the inputs.

Wold, Sjostrom & Eriksson (2001), Chemom. Intell. Lab. Syst. 58; Barker & Rayens (2003), J. Chemometrics 17 (PLS-DA).

06

Linear and quadratic discriminant analysis (LDA / QDA)

supervised

Probabilistic classifiers for known, labeled fault classes. LDA assumes one shared covariance and draws straight boundaries; QDA fits a covariance per class and draws curved ones.

In plain words

When the fault modes are known and you have labeled examples of each, discriminant analysis models every class as a Gaussian and assigns a new sample to the most probable one. LDA ties all classes to a single shared covariance, which yields linear decision boundaries and very few parameters, so it is stable on small data. QDA gives each class its own covariance, which yields quadratic boundaries that fit classes of different shape or spread, at the cost of estimating far more parameters and a real risk of overfitting when classes are small. LDA also doubles as a supervised dimensionality reduction that maximizes between-class separation.

Health-management example

A pump has three labeled states: healthy, cavitation, and bearing wear. Their feature clouds have different shapes, cavitation being far more spread out. LDA's shared-covariance boundary clips part of the cavitation class; QDA's per-class covariance wraps each cloud correctly and recovers the misclassified runs, because here the extra parameters are justified by enough labeled data.

2026-06-19T20:05:47.018072 image/svg+xml Matplotlib v3.10.8, https://matplotlib.org/
from sklearn.discriminant_analysis import (
    LinearDiscriminantAnalysis, QuadraticDiscriminantAnalysis)

lda = LinearDiscriminantAnalysis().fit(X, y)
qda = QuadraticDiscriminantAnalysis().fit(X, y)
print("LDA acc", lda.score(X, y), " QDA acc", qda.score(X, y))
# LDA as supervised projection that maximizes class separation
X_lda = LinearDiscriminantAnalysis(n_components=2).fit_transform(X, y)
Practical note

LDA's shared-covariance assumption is also its strength: it is stable when classes are small or unequal. Switch to QDA only when you have enough labeled samples per class to estimate a separate covariance each, otherwise QDA overfits. Both assume Gaussian classes; strongly non-Gaussian or multi-modal classes need a different classifier. Always report cross-validated, not resubstitution, accuracy.

Fisher (1936), Annals of Eugenics 7; Hastie, Tibshirani & Friedman (2009), Elements of Statistical Learning, ch. 4.

Decision map

Start from the question you actually have. The badge marks whether the method needs labels, separates sources, relates two blocks, or just monitors.

If you are asking...Reach forKind
Is any new sample still inside the normal multi-sensor pattern?PCA T² + Qunsupervised
Which sensor broke the pattern?Q / T² contribution plotunsupervised
Is this point an outlier given the joint correlation?Mahalanobis distanceunsupervised
A fault shows up only as a mixed signalICAsource separation
How does a command block relate to a response block?CCAtwo-block
I have a quality or fault label to separate againstPLS-DAsupervised
Which sensors carry the fault signature (with a label)?PLS VIP scoressupervised
Fault classes are known and labeled; classify new runsLDA / QDAsupervised

Decision tree

new multi-sensor sample have labels? no (monitor) yes (classify) one block or two? one two PCA T² / Q + Mahalanobis CCA subtle fault? yes known modes PLS-DA + VIP LDA / QDA Then adapt to reality: fault arrives mixed at one sensor → ICA to unmix, then monitor the component healthy manifold is curved → kernel PCA / isolation forest in place of linear Q autocorrelated or multi-mode → dynamic PCA, EWMA limits, mode-aware models

A working order of escalation

  1. Scale every sensor to unit variance, and set aside a clean healthy reference period.
  2. Fit PCA on the healthy data; monitor new samples with T² and Q against their limits.
  3. On an alarm, read the Q or T² contribution plot to localize the suspect channel.
  4. Cross-check with Mahalanobis distance, which shares none of PCA's rank choices.
  5. If a fault arrives mixed, unmix with ICA; to relate command and response blocks, use CCA.
  6. Once history is labeled, switch to PLS-DA for subtle faults and LDA or QDA to classify known modes.
# common stack
pip install numpy scipy scikit-learn
# PCAMonitor, pls_da, cca_fit, ica_sources, mahalanobis, lda_qda all ship in the kit engine

Rules of thumb

Heuristics that hold across most multivariate monitoring work, before the method-specific detail.

  1. Scale before you decompose. Unit-variance every sensor first, or PCA, PLS, CCA and Mahalanobis all just track whichever channel has the biggest numbers.
  2. Fit on healthy, score on everything. The model and its control limits come from a clean reference period; new data is only ever scored, never refit into the baseline.
  3. Use T² and Q as a pair, never alone. T² catches in-pattern excursions, Q catches new patterns. Most real faults trip Q first.
  4. Reach for a label when you have one. PLS-DA and LDA beat PCA and Mahalanobis at separating a known fault, because they optimize for the distinction instead of for variance.
  5. Localize, then confirm. Contribution and VIP plots point at suspect sensors but smear blame across correlated channels; treat them as leads, not verdicts.
  6. Mind the sample-to-variable ratio. Covariance-based methods need many more samples than variables; otherwise regularize, shrink, or reduce dimension first.
  7. Validate out of sample. In-sample limits and accuracies are optimistic; hold out clean data, cross-validate latent-variable counts, and permutation-test PLS-DA separation.
  8. Detrend slow operating drift. A baseline that wanders across operating modes inflates every statistic; model per mode or remove the trend before monitoring.
Pick by scenario
ScenarioPick & whyKind
Continuous fleet monitoring, many correlated sensors, no labelsPCA T² + Qthe standard always-on pairunsupervised
Alarm fired, which sensor is responsibleQ / T² contributionsdecompose the statistic onto variablesunsupervised
Point looks normal per channel but feels wrongMahalanobis distanceit scores the joint correlationunsupervised
Overlapping vibration or current sources at one sensorICAunmix independent generatorssource separation
Relate a command block to a response blockCCAcorrelation between two variable setstwo-block
Labeled healthy vs degraded runs, fault is subtlePLS-DAsupervised, aims at the labelsupervised
Known, labeled fault modes, classify new runsLDA / QDAGaussian class modelssupervised
Few labeled samples per classLDAshared covariance is stable on small datasupervised
Plenty of labels, classes differ in shapeQDAper-class covariance fits curved boundariessupervised

The six at a glance

One table to keep the family straight: what each method needs, what it returns, and how its alarm threshold is set.

MethodNeeds labelsNonlinearityPrimary outputLimit / thresholdCost
PCA T² / Qnolinear (curved via kernel PCA)two statistics plus contributionsF-distribution, Jackson-Mudholkarlow
Mahalanobisnono, covariance onlyone distancechi-squarelow
ICAnonon-Gaussian sourcesindependent componentsnone native, monitor a componentmedium
CCAno, two blockslinearcanonical correlationspermutation or F-testlow
PLS-DAyesmostly linearlatent scores plus VIPcross-validatedmedium
LDA / QDAyesLDA linear, QDA curvedclass and posteriorGaussian posteriorlow
Which statistic catches which fault
Fault typeSignatureCaught by
Sensor bias (constant offset)one channel shifts off the othersQ + contribution; Mahalanobis
Slow driftsmall persistent climbEWMA / CUSUM of Q, not a per-sample limit
Stuck-at / flatlinevariance collapses, correlation breaksQ + contribution
Gain changeamplitude scales if in-pattern, Q if it breaks correlation
Two-sensor decouplingeach in range, joint pattern brokenQ + contribution; Mahalanobis
Multivariate mean shiftoperating point moves inside the subspace
New healthy operating modewhole cloud relocates, no faultmode-aware monitoring (a global model false-alarms)

Worked case study: a battery-pack sensor array

One synthetic but structured fault, a cell-group cooling fault that slowly decouples one temperature sensor, taken through the monitoring stack. The pack reports 8 channels; every number below is computed by the engine on the same 300-sample record.

Step 1

Fit the healthy model

PCA on the clean reference period keeps 2 components, capturing 95% of the variance. That fixes the control limits: T² at 9.3 and Q at 1.50. Everything afterward is only ever scored against this model, never folded back into the baseline.

2026-06-19T19:38:54.217637 image/svg+xml Matplotlib v3.10.8, https://matplotlib.org/
Step 2

Detect

On the healthy stretch the Q statistic exceeds its limit 0% of the time, the expected false-alarm rate at the 99% limit. Once the cooling fault develops, Q breaches its limit on 73% of samples while T² stays largely calm: the readings are each still individually plausible, but their joint correlation has broken. That T²-quiet, Q-loud signature is the fingerprint of a new pattern rather than a larger version of an old one.

Step 3

Localize

The Q contribution plot on the worst fault sample puts 82% of the residual on a single channel, T_cell3, the decoupled temperature sensor. The neighbouring cell temperatures pick up a little blame because they are correlated with it, which is exactly why contributions are read as leads, not proof.

2026-06-19T19:38:54.338218 image/svg+xml Matplotlib v3.10.8, https://matplotlib.org/
Step 4

Cross-check and classify

Mahalanobis distance on the same worst sample reads 22 against a limit of 4.5, an independent confirmation that does not rely on the PCA rank choice. With the fault period now labeled, PLS-DA separates healthy from faulted at 88% accuracy and its top VIP variable is again T_cell3, agreeing with the contribution plot. LDA on the labeled features classifies at 88%.

Verdict

What the stack delivered

An unsupervised PCA monitor caught a fault that no single-channel limit would have flagged, Q localized it to the right sensor, Mahalanobis confirmed it without sharing PCA's assumptions, and once labeled the supervised tools both separated the fault and re-identified the culprit channel. No single statistic did all of this; the workflow did. The reflex order is: scale, fit on healthy, watch T² and Q, read contributions to localize, then confirm and classify.

Q contribution on the worst fault sample, top channels:

channelQ contribution
T_cell320.21
T_cell11.52
T_cell41.45
T_cell20.87
T_outlet0.37

The time dimension: sequential detection

A control limit asks one question per sample, in isolation. But a developing fault is a sequence, and a small persistent shift can sit just under the limit for a long time. The fix is to accumulate evidence across samples.

EWMA and CUSUM

Two charts turn a per-sample statistic into a memory of recent behaviour. The exponentially weighted moving average smooths the statistic so a small sustained shift compounds into a clear alarm; the cumulative sum adds up small deviations until they cross a bound. Both detect small drifts far sooner than a raw limit, at the cost of being slower on large sudden jumps. The honest figure of merit is not a single false-positive rate but the average run length: how long, on average, until an alarm, in control and after a shift.

2026-06-19T20:05:47.151073 image/svg+xml Matplotlib v3.10.8, https://matplotlib.org/
2026-06-19T20:05:48.475089 image/svg+xml Matplotlib v3.10.8, https://matplotlib.org/
Practical note

Set the EWMA and CUSUM reference from a clean healthy phase, not from the stream you are watching, or a developing fault inflates the limit and hides itself. Tune the memory to the drift you fear: small lambda for slow faults, larger for abrupt ones. Judge alarms by run length, not single crossings, and report both in-control and post-shift run length when you set a limit. The Tennessee Eastman Process is the canonical public benchmark for comparing these schemes.

Page (1954), Biometrika 41 (CUSUM); Roberts (1959), Technometrics 1 (EWMA); Downs & Vogel (1993), Comput. Chem. Eng. 17 (Tennessee Eastman).

C

Per-sample limit vs. EWMA

run length

A small step appears at the shaded line. Watch the per-sample chart trip late and intermittently while the EWMA accumulates evidence and alarms cleanly.

Operating modes and autocorrelation

The two assumptions most quietly violated on a real fleet: that there is one operating regime, and that successive samples are independent. Both wreck a naive monitor.

One model per mode

A single PCA model assumes one mean and one covariance. The moment the system shifts to a healthy but different regime, cold start versus fast charge versus cruise, every sample looks anomalous and the alarm floods. Mode-aware monitoring clusters the healthy data into regimes, fits a model per regime, and scores each sample against its own mode, so a healthy mode change is no longer a fault. The alternative is one adaptive model that updates slowly, but that risks adapting to a real fault.

2026-06-19T20:05:48.637311 image/svg+xml Matplotlib v3.10.8, https://matplotlib.org/
Practical note

Detect the mode from operating variables you trust, load, speed, temperature setpoint, not from the same sensors you are monitoring. Keep a minimum dwell time before switching models to avoid thrashing at transitions, and treat the transitions themselves as their own regime or exclude them. Recursive or moving-window PCA handles slow baseline drift, but cap the adaptation rate so it cannot quietly learn a fault as normal.

Hwang & Han (1999), Control Eng. Practice 7; Zhao, Wang & Jia (2004), Chemom. Intell. Lab. Syst. 73.

D

One global model across two modes

non-stationary

Mode B is perfectly healthy, just a different operating regime. Separate the modes and watch a global model raise false alarms; switch to mode-aware and they vanish.

Autocorrelation breaks the i.i.d. limit

The control limits assume independent samples. Real sensor streams are autocorrelated, so the monitoring statistic is too: excursions arrive in clusters rather than one at a time, and the true alarm behaviour drifts from the nominal even when the average rate looks right. Dynamic PCA absorbs the temporal structure by stacking each sample with its recent past before decomposing, which restores the independence the limits rely on; alternatively, set limits empirically to respect the autocorrelation and judge alarms by run length.

2026-06-19T20:05:48.788477 image/svg+xml Matplotlib v3.10.8, https://matplotlib.org/
Practical note

Choose the number of lags from the autocorrelation of the residual, adding lags until it whitens, usually one or two. Time-lag embedding multiplies the variable count, so guard the sample-to-variable ratio. When in doubt, set the control limit from a held-out healthy stream rather than the closed-form formula.

Ku, Storer & Georgakis (1995), Chemom. Intell. Lab. Syst. 30.

Modern relatives of the Q statistic

Linear PCA measures distance to a flat plane. When the healthy region is curved, that plane is the wrong reference and Q misfires. The learned one-class detectors are the nonlinear descendants of the same idea.

When the healthy manifold curves

Kernel PCA computes the same reconstruction error in a nonlinear feature space, so the healthy support may bend. One-class SVM learns a tight boundary around the healthy region directly, and isolation forest scores how easily a point is separated from the bulk. All three play the role of Q for data whose normal operating envelope is not a straight subspace, and all three pay for it with more data, more tuning, and less transparent contributions. The classical T-squared and Q remain the first choice when the healthy region is roughly linear, precisely because they are cheap and their contributions localize a fault.

2026-06-19T20:05:49.108266 image/svg+xml Matplotlib v3.10.8, https://matplotlib.org/
Practical note

Reach for these only when linear PCA visibly fails, for example when healthy Q is high on data you know is healthy. They need enough data to trace the manifold and careful bandwidth or contamination settings, and they do not hand you a clean contribution plot, so localization needs a separate step such as feature ablation. An autoencoder's reconstruction error is the same construction at larger scale.

Scholkopf et al. (2001), Neural Computation 13 (one-class SVM); Liu, Ting & Zhou (2008), ICDM (isolation forest); Lee et al. (2004), Chem. Eng. Sci. 59 (kernel PCA monitoring).

From detection to diagnosis

Detection says something is wrong. Diagnosis says which sensor, and keeps the monitor running when one drops out. Two refinements of the contribution idea do most of this work.

Reconstruction-based contributions

A plain contribution plot smears blame onto every sensor correlated with the faulty one, because they all share the broken pattern. Reconstruction-based contribution instead asks, for each variable, how much the statistic would fall if that variable alone were corrected back onto the model. It guarantees the true single-variable fault gets the largest score and sharply reduces the smearing, which matters when neighbours are tightly coupled.

2026-06-19T20:05:49.210103 image/svg+xml Matplotlib v3.10.8, https://matplotlib.org/
Practical note

RBC isolates a single faulty variable cleanly but can still be fooled by multiple simultaneous faults; for those, combine it with a structured fault library. It is a diagnosis tool, not a detection one, so apply it only after a statistic has already alarmed.

Alcala & Qin (2009), Automatica 45.

Running through a sensor dropout

The same healthy correlation structure that detects a fault can also fill a gap. When a sensor fails or its reading is rejected, the PCA model reconstructs it from the others by projecting the remaining measurements onto the model subspace, so monitoring continues on the recovered value instead of going blind. The reconstruction degrades gracefully as more sensors drop, and the residual on the reconstructed channel is itself a check on the assumption.

2026-06-19T20:05:49.298716 image/svg+xml Matplotlib v3.10.8, https://matplotlib.org/
Practical note

Reconstruction is trustworthy only while the missing variable is well predicted by the retained components; check the model's reconstruction error per variable on healthy data first. Never reconstruct a sensor and then use that reconstruction to declare the same sensor faulty, which is circular.

Nelson, Taylor & MacGregor (1996), Chemom. Intell. Lab. Syst. 35.

From detection to prognosis

The monitoring statistic is not only an alarm. Tracked over a component's life, a rising Q or its smoothed version is a health index, and a health index can be projected forward to a remaining-useful-life estimate.

The Q statistic as a degradation feature

A fault rarely appears all at once; it grows. If the EWMA of Q trends upward as a bearing wears or a cell ages, that trend is a degradation signal. Fitting a model to the trend and extrapolating it to the failure threshold turns detection into prediction, the bridge from condition monitoring into prognostics and remaining-useful-life estimation, where survival models and Gaussian-process regression take over with their own uncertainty bands.

2026-06-19T20:05:49.403490 image/svg+xml Matplotlib v3.10.8, https://matplotlib.org/
Practical note

A monitoring statistic makes a usable health index only when it trends monotonically with damage; verify that on run-to-failure data before trusting an extrapolation. Report a confidence interval on the remaining life, not a point, and re-estimate as new data arrives, since early-life trends are weak predictors. This is the handoff point to a dedicated prognostics model.

Lee et al. (2014), Mech. Syst. Signal Process. 42 (PHM review); Si et al. (2011), Eur. J. Oper. Res. 213 (RUL estimation).

References

  1. Jackson & Mudholkar (1979), Technometrics 21; Wise & Gallagher (1996), J. Process Control 6; Qin (2003), J. Chemometrics 17.
  2. Mahalanobis (1936), Proc. Natl. Inst. Sci. India 2; Rousseeuw & Van Zomeren (1990), JASA 85 (robust covariance).
  3. Hyvarinen & Oja (2000), Neural Networks 13 (FastICA); Lee et al. (2004), J. Process Control 14 (ICA monitoring).
  4. Hotelling (1936), Biometrika 28; Hardoon, Szedmak & Shawe-Taylor (2004), Neural Computation 16.
  5. Wold, Sjostrom & Eriksson (2001), Chemom. Intell. Lab. Syst. 58; Barker & Rayens (2003), J. Chemometrics 17 (PLS-DA).
  6. Fisher (1936), Annals of Eugenics 7; Hastie, Tibshirani & Friedman (2009), Elements of Statistical Learning, ch. 4.

Take it into your editor

The same content is packaged as a VS Code Copilot skill in the progressive-disclosure layout: a routing SKILL.md, per-topic reference docs, a tested multivariate_monitoring.py engine that implements every method with numpy, scipy and scikit-learn only, and Copilot integration files (a monitoring-analyst chat mode, a /monitor prompt, and path-scoped instructions). Drop the folder into your kit, or point Copilot at the skill directly.

multivariate-monitoring-kit.zip
SKILL.md · references · multivariate_monitoring.py · Copilot agent + prompt

Link resolves when the zip sits beside this file. Run python scripts/multivariate_monitoring.py for the self-test on synthetic fault data.