当内存装不下数据时
在机器学习特征工程或数据分析中，我们经常遇到这样的场景：手头有成百上千个独立的特征文件（CSV、Parquet 或 Numpy 格式），总量达到了几百 GB 甚至 TB 级别。现在需要计算这些特征的全局统计量（平均值、方差、标准差）来进行归一化（Standardization）。然而，开发机内存只有 16GB。如果尝试简单的 pandas.read_csv() 或 numpy.concatenate() 把所有数据读入内存，程序会瞬间 OOM（Out of Memory）崩溃。面对据量 >> 内存的场景，我们需要一种流式（Streaming）或增量（Incremental）的计算方案——Welford 算法。
教科书公式的陷阱
在统计学教科书中，我们学到的方差计算公式通常是这样的：
$$\sigma^2 = \frac{\sum_{i=1}^n x_i^2 - \frac{(\sum_{i=1}^n x_i)^2}{n}}{n}$$
或者它的变形：
$$\sigma^2 = \frac{1}{n} \sum_{i=1}^n x_i^2 - \bar{x}^2$$
这种公式非常适合手算，但在计算机工程实现中，它有两个致命缺陷：
- 内存不友好：你需要维护两个巨大的和（$\sum x$ 和 $\sum x^2$），虽然这可以通过累加解决，但无法避免第二个问题。
- 数值稳定性极差（Catastrophic Cancellation）：这是最严重的问题。假设你的数据数值很大（例如 $10^9$），那么 $\sum x^2$ 会变成一个天文数字。当公式中计算 $\sum x^2 - n\bar{x}^2$ 时，实际上是在做两个非常巨大的数字相减。在浮点数运算中，大数相减会丢失大量的有效位数，导致结果精度极低，甚至算出负的方差（这在数学上是不可能的，但在计算机里常发生）。
我们需要一种算法，它既不需要存储历史数据，又能在计算过程中保持数值较小，这就是 Welford 在线算法。
Welford 算法的核心思想
B. P. Welford 于 1962 年提出了一种迭代算法。它的核心思想是：每读入一个新的数据点，就利用旧的统计量，修正出新的统计量。
我们只需要维护三个变量：
- $n$：当前的样本计数。
- $\mu$：当前的平均值（Mean）。
- $M_2$：当前的平方差聚合值（用于计算方差）。
迭代公式每当流式数据中进来一个新的数值 $x$，更新步骤如下：
- 计数加一：$n \leftarrow n + 1$
- 更新均值：$\delta = x - \mu_{\text{old}}, \quad \mu_{\text{new}} \leftarrow \mu_{\text{old}} + \frac{\delta}{n}$
- 更新 $M_2$（最精妙的一步）：$\delta_2 = x - \mu_{\text{new}}, \quad M_2 \leftarrow M_2 + \delta \times \delta_2$
- 最终的标准差计算：$\text{std} = \sqrt{\frac{M_2}{n-1}}$
为什么它更优秀？
在这个公式中，我们始终在计算 $x - \mu$（数据点与均值的差值）。无论原始数据 $x$ 有多大，这个差值通常都比较小。计算小数字的平方和，永远比计算大数字的平方和要精确得多。
代码实践
import numpy as np
class WelfordStats:
"""Welford's Online Algorithm for calculating Mean and Variance.
适合流式数据、大数据集的增量计算。
"""
def __init__(self):
self.n = 0
self.mean = 0.0
self.M2 = 0.0
def update(self, x):
"""处理单个数值"""
self.n += 1
delta = x - self.mean
self.mean += delta / self.n
delta2 = x - self.mean
self.M2 += delta * delta2
def update_batch(self, chunk):
""" 处理一批数据 (Numpy Array)
在实际文件读取中，分块读取效率远高于逐行读取。
"""
chunk = np.asarray(chunk)
# 针对这一批数据进行循环更新
for x in chunk:
self.update(x)
@property
def variance(self):
"""样本方差 (Sample Variance, 分母 n-1)"""
if self.n < 2:
return 0.0
return self.M2 / (self.n - 1)
@property
def std(self):
"""样本标准差"""
return np.sqrt(self.variance)
def __str__(self):
return f"Count: {self.n}, Mean: {self.mean:.6f}, Std: {self.std:.6f}"
回到文章开头的需求，假设有一堆 CSV 文件，我们可以这样计算全局分布：
import pandas as pd
import glob
# 1. 初始化统计器
stats = WelfordStats()
# 2. 获取文件列表
file_list = glob.glob("features/*.csv")
# 3. 流式读取并计算
for file_path in file_list:
print(f"Processing {file_path}...")
# 分块读取 (chunksize)，避免一次性读入整个文件
# 假设特征列名为 'feature_value'
for chunk in pd.read_csv(file_path, chunksize=10000):
data_chunk = chunk['feature_value'].values
# 喂给 Welford 算法
stats.update_batch(data_chunk)
# 4. 输出最终结果
print("-"*30)
print("全局统计结果：")
print(stats)
推导过程
公式 $M_{2,n} = M_{2,n-1} + (x_n - \mu_{n-1})(x_n - \mu_n)$ 是成立的核心在于利用均值的递推关系进行代数变换。我们要计算的是前 $n$ 个数据的离差平方和：
$M_{2,n} = \sum_{i=1}^n (x_i - \mu_n)^2$
将其拆解为前 $n-1$ 项和第 $n$ 项，并引入旧均值 $\mu_{n-1}$ 进行配方，经过一系列消项和化简（利用 $\sum(x_i - \mu)=0$ 的性质），最终可以证明：
$\sum_{i=1}^n (x_i - \mu_n)^2 = \sum_{i=1}^{n-1} (x_i - \mu_{n-1})^2 + (x_n - \mu_{n-1})(x_n - \mu_n)$
这正是 Welford 算法中 $M_2$ 的更新逻辑：新的 $M_2$ = 旧的 $M_2$ + (当前值 - 旧均值) × (当前值 - 新均值)
总结
Welford 算法是处理大规模数据统计的'最佳实践'。
- 空间复杂度：$O(1)$。无论数据多少，只占 3 个变量的内存。
- 时间复杂度：$O(N)$。只需遍历一次数据。
- 数值稳定性：极高，避免了大数吃小数的问题。
