3.3.2　期望最大化填补法

3.3.1节详细介绍了极大似然估计法，期望最大化填补法是基于极大似然估计法的缺失值填补方法。在该方法中，缺失值作为变量参与到参数估计过程中，以迭代的方式交替更新缺失值与待估计参数，以此达到缺失值填补的目标。

针对数据集X中第c类样本构成的子集X^(c)，假设其中的所有现有值构成集合Xp^(c)，而所有的缺失值构成集合X_m^(c)。对参数向量β^(c)进行估计时，需最大化式（3-85）所示的对数似然函数：

由于X_m^(c)为缺失值，故无法参照3.3.1节中的描述估计参数向量β^(c)。此时，可采用期望最大化算法以迭代的方式在估计参数的同时进行缺失值填补。该算法每轮迭代包括两个步骤，即期望步和最大化步，以下简称E步和M步。在E步中，根据参数向量β^(c)和现有值计算缺失值；在M步中，根据填补结果对参数向量β^(c)进行极大似然估计。为了更加具体地说明填补过程，下面给出一个例子。

现在有3个编号分别为A、B、C的箱子，每个箱子中都有黑白两种球，各箱子中黑球所占的比例分别为p_A、p_B、p_C。实验步骤如下：从箱子A中取出一个球，若为黑球，则从箱子B中取出一个球；若为白球，则从箱子C中取出一个球。箱子B或箱子C所抽取球的颜色即为实验结果。对于抽取球的颜色，如果是黑色则记为1，如果是白色则记为0。独立重复实验N次，假设每次实验只能观测到最终实验结果，不能观测抽取过程，则对于此实验而言，箱子A的抽取结果可视作缺失值，待估计参数p_A、p_B、p_C组成参数向量β=[p_A,p_B,p_C]。下面通过EM算法对参数向量β进行极大似然估计，并同时填补箱子A的抽取结果。

如3.3.1节所述，极大似然估计法的目的是：利用已知的样本结果，反推最大概率导致这种结果的参数值。本节加入EM算法的目的则是在估计参数值的同时进行缺失值填补。在上述黑白球问题中，每次独立重复实验的结果是完整的，N次实验的结果被记为Y=[y₁,y₂,…,y_N]；箱子A的抽取结果是缺失的，该抽取结果记为Z=[z₁,z₂,…,z_N]，其仅有两种可能的取值，即1或0。

根据式（3-75），参数向量β=[p_A,p_B,p_C]对于Y=[y₁,y₂,…,y_N]和Z=[z₁,z₂,…,z_N]的似然函数如式（3-86）所示：

式中，P(y_i,z_i|β)为在参数向量β条件下y_i、z_i的概率分布。基于式（3-86）所得的对数似然函数如式（3-87）所示：

极大似然估计的目标为寻找使得该对数似然函数取最大值的参数向量，其目标函数如式（3-88）所示：

由于Z=[z₁,z₂,…,z_N]为缺失数据，故采用EM算法在填补缺失值的基础上进行参数估计。在EM算法开始之前，需随机初始化参数向量β=[p_A,p_B,p_C]，记该初始参数向量为β⁽⁰⁾。EM算法采用迭代的方式交替更新填补值和待估计参数，每轮迭代分为E步和M步。假设迭代次数的上限为L，则第l次迭代中，E步可分为两个阶段，第一阶段基于现有数据Y=[y₁,y₂,…,y_N]和当前参数向量β^(l-1)=[p_A^(l-1),p_B^(l-1),p_C^(l-1)]填补缺失值，第二阶段基于填补结果求对数似然函数H(β)的期望。在第一阶段中，由于Z为离散变量，故在缺失值填补时将该变量取各离散值的概率作为填补结果。对于z_i∈Z，其取值范围为{0,1}，根据概率学的相关知识求得其取值为1的概率，如式（3-89）所示：

式（3-89）中，P(y_i,z_i=1|β^(l-1))的计算过程如式（3-90）所示：

同理可得，P(y_i,z_i=0|β^(l-1))如式（3-91）所示：

将式（3-90）和式（3-91）带入式（3-89），所得结果如式（3-92）所示：

由式（3-92）可见，z_i=1的概率完全可由y_i和β^(l-1)表示，即根据现有值y_i和参数向量β^(l-1)计算出离散变量z_i的填补值。为方便描述，记μ_i^(l)=P(z_i=1|y_i,β^(l-1))为第l次迭代中z_i取值为1的概率，则z_i=0的概率如式（3-93）所示：

在第二阶段中，基于填补结果组成的向量Z求对数似然函数H(β)的期望，计算方法如式（3-94）所示：

式（3-94）中，P(y_i,z_i=1|β^(l-1))和P(y_i,z_i=0|β^(l-1))分别如式（3-90）和式（3-91）所示，P(z_i=1|y_i,β^(l-1))和P(z_i=0|y_i,β^(l-1))分别如式（3-92）和式（3-93）所示。因此，第二阶段基于填补结果所得对数似然函数H(β)的期望如式（3-95）所示：

在M步中，对各参数求偏导得到使式（3-95）期望最大化的参数值。对于第l轮迭代(l=1,2,…,δ)，M步所求参数组成的参数向量为β^(l)=[p_A^(l),p_B^(l),p_C^(l)]。

该期望对参数p_A的偏导如式（3-96）所示：

对式（3-96）简化得式（3-97）：

同理，可求得参数p_B和p_C的更新规则，分别如式（3-98）和式（3-99）所示：

每一轮迭代分E步和M步进行，在E步，首先根据式（3-89）和式（3-93）更新填补值，随后基于填补结果根据式（3-95）计算对数似然函数的期望；在M步，根据式（3-97）、式（3-98）和式（3-99）进行参数估计。重复上述迭代直到两次迭代之间的参数更新幅度小于设定的阈值，或迭代次数达到上限。最后一轮迭代的填补值为最终填补结果，该轮迭代中M步估计的参数为极大似然估计结果。

基于上述黑白球实例对期望最大化填补法进行提炼总结。假设数据集X中第c类样本构成子集X^(c)，其中现有值的集合为X_p^(c)，缺失值的集合为X_m^(c)，其中样本分布由参数向量β^(c)唯一确定。以下为期望最大化算法填补缺失值的流程。

步骤1：初始化参数向量β^(c)，EM算法的阈值为ε，迭代次数上限为L；

步骤2：构建参数向量β^(c)对于子集X^(c)的对数似然函数，如式（3-100）所示：

在式（3-100）的基础上设定如式（3-101）所示的优化目标：

式中，^(c)表示使H(β^(c))达到最大值的参数向量β^(c)取值；

步骤3：根据当前参数和现有值进行缺失值填补，即E步的第一阶段。第l轮迭代(l=1,2,…,L)获得的填补值可表示为_m^(c)(l)；

步骤4：基于填补结果_m^(c)(l)计算极大似然函数的期望，即E步的第二阶段，记为；

步骤5：对各参数求偏导计算使该期望取最大值的参数值，即M步。第l轮迭代（l=1,2,…,δ）获得的参数向量记为β^(c)(l)；

步骤6：如果两轮迭代参数的变化幅度小于阈值ε，或迭代次数达到上限L，则迭代结束。将最后一轮的填补值作为最终填补结果，记为_m^(c)，最后一轮求得的参数向量为极大似然估计结果，记为^(c)；否则，返回步骤3。

期望最大化填补法以极大似然估计法作为理论基础，采用迭代的方式填补缺失值，对数据集中完整数据的利用较为充分，通常能够获得较为精确的填补结果。但该方法的精度与缺失率相关，当缺失率太大时，上述迭代优化过程容易陷入局部最优解，在影响填补精度的同时还会导致方法的收敛速度显著降低。