1. EM算法与高斯混合模型基础
EM算法（Expectation-Maximization）是处理含隐变量概率模型参数估计的经典方法。我第一次接触这个概念是在处理用户行为数据时，发现传统聚类方法无法解释某些异常点分布。EM算法的精妙之处在于它通过"猜测-验证"的迭代方式，逐步逼近真实参数。
高斯混合模型（GMM）是EM算法最典型的应用场景之一。想象你面前有一堆不同品种的葡萄，有的大而圆，有的小而椭圆。如果只用单一高斯分布描述，相当于假设所有葡萄来自同一品种，这显然不合理。GMM通过多个高斯分布的线性组合，可以更准确地描述复杂数据分布。
在数学表达上，一个K成分的GMM概率密度函数为：
p(x) = Σ[π_k * N(x|μ_k, Σ_k)] # k=1到K其中π_k是混合系数，满足Σπ_k=1。这个公式就像用多个不同形状的钟形曲线叠加来描述数据。
2. EM算法核心原理拆解
EM算法的迭代过程让我想起教孩子认动物的过程。E步相当于让孩子根据现有认知猜测图片中的动物（建立隐变量分布），M步则是根据孩子的猜测更新认知模型（更新参数）。这种类比帮助我理解了算法的本质。
算法流程具体分为：
- 初始化：随机设置各高斯成分的参数
- E步骤：计算后验概率γ(z_nk)，表示数据点n属于第k个成分的概率
# E步骤示例代码 gamma = pi * multivariate_normal.pdf(X, mu, sigma) # 未归一化 gamma /= gamma.sum(axis=1, keepdims=True)- M步骤：根据当前γ重新估计参数
# M步骤示例代码 Nk = gamma.sum(axis=0) mu = (gamma.T @ X) / Nk sigma = np.array([(gamma[:,k]*(X-mu[k]).T)@(X-mu[k])/Nk[k] for k in range(K)]) pi = Nk / len(X)实际项目中我遇到的一个坑是初始化敏感性问题。有次用随机初始化导致算法收敛到糟糕的局部最优，后来改用k-means初始化才解决。这也印证了EM算法对初始值依赖的特性。
3. 时空数据缺失值填补实战
处理PM2.5监测数据时，我发现缺失情况主要分两类：
- 部分缺失：某站点少数时间点数据缺失
- 完全缺失：某站点全天数据缺失
对于部分缺失，传统插值方法可能足够。但遇到某城市监测站全天宕机的情况时，就需要结合空间信息。我的解决方案是构建时空GMM模型：
- 数据准备：将每个站点的每小时数据视为一个特征向量，加入经纬度信息
- 模型构建：定义包含空间权重的协方差矩阵
# 空间权重计算（以经纬度为例） def spatial_kernel(lon1, lat1, lon2, lat2): return np.exp(-0.5*((lon1-lon2)**2 + (lat1-lat2)**2)/h)- EM优化：在E步计算时，对缺失数据只使用已知维度计算似然
在具体实现时，需要注意：
- 对完全缺失的站点，先赋予邻近站点的均值
- 设置合理的空间衰减系数h
- 加入时间趋势项处理周期性变化
4. 工程实践中的优化技巧
经过多个环保监测项目实践，我总结了以下经验：
参数初始化：
- 使用k-means++初始化均值
- 协方差矩阵初始化为各维度方差的单位矩阵
- 混合系数初始化为均匀分布
加速收敛：
# 加入动量项的参数更新 mu_new = (1-momentum)*mu_mle + momentum*mu_old处理奇异矩阵：
sigma[k] += 1e-6*np.eye(n_features) # 加入微小正则项一个有趣的发现是：当处理长三角地区数据时，加入风速风向作为辅助变量后，填补精度提升了约15%。这说明领域知识的融入能显著提升模型效果。
对于超参数选择，我通常采用：
- 通过交叉验证选择最佳成分数K
- 用BIC准则平衡模型复杂度
- 空间衰减系数h根据站点平均距离设定
最终完整的时空填补流程包括：
- 数据标准化
- 构建时空特征矩阵
- EM迭代优化
- 后处理平滑
- 结果验证
这种方案在某省级环保项目中，将日均填补误差从12.5μg/m³降低到7.8μg/m³，验证了其有效性。
