Variance is one of the central dispersion measures in statistics. In practical systems, however, variance is often not computed once on a fixed table, but updated continuously while measurements arrive, are corrected, or are replaced in a sliding window.
The following develops numerically stable formulas for these scenarios. We begin with the standard sample-variance definition, explain why naive one-pass formulas are unstable in floating-point arithmetic, derive Welford's stable recurrence, and then derive a replacement update that avoids full recomputation.
Definitions and Notation
Let the sample be \(x_1, x_2, \dots, x_n\) with \(n \ge 2\). Define:
\[ \bar{x} = \frac{1}{n}\sum_{i=1}^{n}x_i, \qquad s^2 = \frac{1}{n-1}\sum_{i=1}^{n}(x_i-\bar{x})^2. \]
The quantity
\[ M_2 := \sum_{i=1}^{n}(x_i-\bar{x})^2 \]
is the corrected sum of squares, so \(s^2 = M_2/(n-1)\).
Classical Equivalent Form and Numerical Instability
Algebra yields the equivalent identity
\[ s^2 = \frac{1}{n-1}\left(\sum_{i=1}^{n}x_i^2 - n\bar{x}^2\right). \]
Although exact in real arithmetic, this form is often unstable in floating-point arithmetic because two large, nearly equal terms are subtracted. This catastrophic cancellation can cause severe relative error when true variance is small compared with the data magnitude.
Welford's Stable One-Pass Algorithm
Welford's method maintains the running mean \(M_k\) and running corrected sum \(S_k\) after \(k\) samples.
Initialization:
\[ M_1 = x_1, \qquad S_1 = 0. \]
For each new sample \(x_k\), \(k \ge 2\):
\[ M_k = M_{k-1} + \frac{x_k - M_{k-1}}{k}, \]
\[ S_k = S_{k-1} + (x_k - M_{k-1})(x_k - M_k). \]
Final sample variance after \(n\) samples:
\[ s^2 = \frac{S_n}{n-1}. \]
This recurrence is stable because it avoids directly subtracting two large accumulated quantities.
Variance Update Under Single-Element Replacement
Consider a fixed-size sample of size \(N\), with known mean \(\bar{x}\) and variance \(s^2\). Suppose one element \(x_{old}\) is replaced by \(x_{new}\).
Updated Mean
The updated mean is
\[ \bar{x}' = \bar{x} + \frac{x_{new} - x_{old}}{N}. \]
Updated Variance Through \(M_2\)
Let \(M_2 = (N-1)s^2\). A stable replacement update is:
\[ M_2' = M_2 + (x_{new}-x_{old})(x_{new}-\bar{x}' + x_{old}-\bar{x}). \]
Then
\[ s'^2 = \frac{M_2'}{N-1}. \]
This avoids recomputing \(M_2\) from all \(N\) values and is well suited to rolling-window pipelines.
Algorithmic Complexity
- Per new sample (Welford): \(O(1)\) time, \(O(1)\) memory.
- Per replacement update: \(O(1)\) time, \(O(1)\) memory.
- No full pass over \(N\) values is required for each update.
JavaScript Implementation
/**
* RollingStat.js
* Stable one-pass variance (Welford) + O(1) replacement update.
*/
class RollingStat {
constructor() {
this.n = 0;
this.mean = 0;
this.M2 = 0;
}
variance() {
return this.n > 1 ? this.M2 / (this.n - 1) : 0;
}
standardDeviation() {
return Math.sqrt(this.variance());
}
add(x) {
this.n++;
if (this.n === 1) {
this.mean = x;
this.M2 = 0;
return;
}
const delta = x - this.mean;
this.mean += delta / this.n;
this.M2 += delta * (x - this.mean);
}
replace(oldValue, newValue) {
if (this.n < 2) return;
const delta = newValue - oldValue;
const oldMean = this.mean;
const newMean = oldMean + delta / this.n;
this.mean = newMean;
this.M2 += delta * (newValue - newMean + oldValue - oldMean);
}
} Practical Notes
- Prefer storing \(M_2\) instead of variance directly to reduce repeated scaling noise.
- For very long streams, monitor for overflow if value magnitudes are extreme.
- After replacement-heavy workloads, occasional full recomputation can be used as an integrity check.
References
- Welford62Welford, B.P. (1962). Note on a method for calculating corrected sums of squares and products. Technometrics, 4(3), 419-420.
- Knuth97Knuth, D.E. (1997). The Art of Computer Programming, Vol. 2: Seminumerical Algorithms, 3rd ed., Addison-Wesley.
- Chan83Chan, T.F., Golub, G.H., LeVeque, R.J. (1983). Algorithms for computing the sample variance: analysis and recommendations. The American Statistician, 37(3), 242-247.
- Higham02Higham, N.J. (2002). Accuracy and Stability of Numerical Algorithms. SIAM.