深入理解线性模型（二）---基于似然函数的估计

更新时间：2019.10.31

目录

• 1. 引言
• 2. 关于$varepsilon$假设
• 3. 基于似然函数的估计
  • 3.1 基于假设1
  • 3.2 基于假设2
  • 3.3. 基于假设3
• 4. 估计的优良性
• 5. 假设的场景

1. 引言

  在上一篇中，我们从损失函数的角度出发讨论了(eta)和(sigma)的估计。在本篇将换一种极具统计味道的角度，从似然函数出发来讨论了(eta)和(sigma)的估计。从中我们也将看见，在不同的假设中，损失函数将会发生不同的变化。

2. 关于(varepsilon)假设

  在上一篇(基于损失函数的估计)中，我们提到，对于线性模型，我们常常使用Guass-Markov假设，即：

1. (E(varepsilon) = 0)
2. (cov(varepsilon) = sigma^2 I_n)

  但是，实际上我们同方差的假设是总是不满足的，完整来说，对(varepsilon)的假设应该有三种：

1. 同方差，且各个随机误差变量不相关：(cov(varepsilon) = sigma^2 I_n)
2. 异常差，但各个随机误差变量不相关，(cov(varepsilon) = diag(sigma_1^2, sigma_2^2, cdots, sigma_n^2))
3. 异方差，且各个随机误差变量是相关的，

[cov(varepsilon) = egin{pmatrix} sigma_{11}^2 & cov(varepsilon_1, varepsilon_2) & cdots & cov(varepsilon_1, varepsilon_n)\ cov(varepsilon_2, varepsilon_1) & sigma_{22}^2 & cdots & cov(varepsilon_2, varepsilon_n)\ vdots & vdots & & vdots\ cov(varepsilon_n, varepsilon_1) & cov(varepsilon_n, varepsilon_2) & cdots & sigma_{nn}^2 end{pmatrix} ]

  此时，记(cov(varepsilon) = Sigma)

3. 基于似然函数的估计

  之前是从损失函数的角度进行参数的估计，但是实际上每个损失函数都应该对应着一个分布，并使得分布的似然函数达到最大
  我们知道在X给定的情况下，似然函数(L( heta;Y,X) = P_{ heta}(Y_1 = y_1, Y_2 = y_2, cdots, Y_n = y_n))。假设(Y_1, Y_2, cdots, Y_n)是独立的，有(L( heta;Y,X) = prod_{i=1}^nP(Y = y_i))。当是离散情况的时候，可以进一步化为：(L( heta;Y,X) = prod_{i=1}^nP_i( heta))。当是连续情况的时候，则可以化为：(L( heta;Y,X) = prod_{i=1}^n f(y_i; heta))

3.1 基于假设1

  如果满足假设1，(cov(varepsilon) = sigma^2 I_n)， 并加上一个正态性的假设，即有(varepsilon_i sim N(0, sigma^2))，那么，(y_i = x_ieta + varepsilon_i sim N(x_ieta, sigma^2))，那么有似然函数：
egin{equation}
egin{split}
L(eta, sigma^2, Y, X) & = prod_{i=1}^n f(y_i)\
& = prod_{i=1}^n frac{1}{sqrt{2pi}sigma} e^{- frac{(y_i - x_i eta)^2}{2 sigma^2}}\
& = (frac{1}{sqrt{2pi}sigma})^n e^{- frac{1}{2 sigma^2} displaystyle sum_{i=1}^n(y_i - x_ieta)^2}
end{split}
end{equation}

  可以看到，似然函数中含有的(sum_{i=1}^n(y_i - x_ieta)^2)部分正是我们之前讨论的二次损失形式。那么我们便了解到，基于假设1时，确实是应该采用我们之前所使用的二次损失形式
  通常为了简便计算，我们都会将似然函数对数化

egin{equation}
egin{split}
lnL(eta, sigma^2, Y, X) & = -nln(sqrt{2pi}sigma)- frac{1}{2 sigma^2} sum_{i=1}^n(y_i - x_ieta)^2
end{split}
end{equation}

  记(G(eta, sigma^2) = nln(sqrt{2pi}sigma) + frac{1}{2 sigma^2} sum_{i=1}^n(y_i - x_ieta)^2)，令似然函数最大化，即是求(min hspace{1mm}G(eta, sigma^2))

  对(G(eta, sigma^2))求关于(eta)的偏导有

egin{equation}
egin{split}
frac {partial G(eta, sigma^2)}{partial eta}
&= 0 + frac{1}{2 sigma^2}2 displaystyle sum_{i=1}^n (y_i - x_i eta)x_i\
& = frac{1}{2 sigma^2} displaystyle sum_{i=1}^n 2(x_i y_i - x_i^2 eta) = 0
end{split}
\
=> displaystyle sum_{i=1}^n (x_i y_i - x_i^2 eta) = 0 => displaystyle sum_{i=1}^n x_iy_i = displaystyle sum_{i=1}^n x_i^2 eta\
=> X^T Y = X^T X eta => hat eta = (X^T X)^{-1} X^T Y
end{equation}

  对(G(eta, sigma^2))求关于(sigma)的偏导有

egin{equation}
egin{split}
frac {partial G(eta, sigma^2)}{partial sigma}
&= nfrac{1}{sqrt{2pi} sigma} sqrt{2pi} - frac{2}{2 sigma^3} sum_{i=1}^n(y_i - x_i eta)^2 \
& = frac{n}{sigma} + frac{1}{sigma^3} sum_{i=1}^n(y_i - x_i eta)^2 = 0
end{split}
\
=> frac{1}{sigma^3} sum_{i=1}^n(y_i - x_i eta)^2 = frac{n}{sigma}
=> hat sigma^2 = frac{displaystyle sum_{i=1}^n(y_i - x_ieta)^2}{n}
end{equation}

  从这里便可以看出，通过似然函数，一次就搞定了参数(eta)和(sigma)的估计，而基于损失函数的估计只是估计出了(eta)，而(sigma)是另外造一套理论估计的

• tips：但是基于似然函数的(sigma)估计有一个小问题，它所得到的不是一个无偏估计（和(hat sigma^2 = frac{SSE}{n-p})略显不同）。因此，有的人也采用限制似然估计（REML）来进行代替。

3.2 基于假设2

  如果满足假设2，(cov(varepsilon) = cov(varepsilon) = diag(sigma_1^2, sigma_2^2, cdots, sigma_n^2))， 并加上一个正态性的假设，即有(varepsilon_i sim N(0, sigma^2_{ii}))，那么，(y_i = x_ieta + varepsilon_i sim N(x_ieta, sigma^2_{ii}))，那么有似然函数：

egin{equation}
egin{split}
L(eta, sigma^2, Y, X) & = prod_{i=1}^n f(y_i)\
& = prod_{i=1}^n frac{1}{sqrt{2pi} sigma_{ii}} e^{- frac{(y_i - x_i eta)^2}{2 sigma^2_{ii}}}\
& = (frac{1}{sqrt{2pi}})^n prod_{i=1}^n(frac{1}{sigma_{ii}}) e^{- frac{1}{2} displaystyle sum_{i=1}^n(frac {y_i - x_i eta}{sigma_{ii}})^2}
end{split}
end{equation}

  我们可以发现基于假设2下，似然函数的核心部分发生了变化，不再是(sum_{i=1}^n(y_i - x_ieta)^2)。因此，根据之前的经验，基于假设2，所采用的损失函数也应该发生变化。此时采用的损失函数应该是标准化的二次损失(displaystyle sum_{i=1}^n(frac {y_i - x_ieta}{sigma_{ii}})^2)，我们也把这称为加权最小二乘估计。
  将似然函数对数化：
egin{equation}
egin{split}
lnL(eta, sigma^2, Y, X) = -nln(sqrt{2pi})- sum_{i=1}^nlnsigma_{ii} - frac{1}{2} displaystyle sum_{i=1}^n(frac {y_i - x_ieta}{sigma_{ii}})^2
end{split}
end{equation}

  记(G(eta, sigma_{ii}^2) = nln(sqrt{2pi}) + sum_{i=1}^nlnsigma_{ii} + frac{1}{2} displaystyle sum_{i=1}^n(frac {y_i - x_ieta}{sigma_{ii}})^2)，令似然函数最大化，即是求(min hspace{1mm}G(eta, sigma_{ii}^2))
  对(G(eta, sigma_{ii}^2))求关于(eta)的偏导有

egin{equation}
egin{split}
frac {partial G(eta, sigma_{ii}^2)}{partial sigma_{ii}}
&= 0 + 0 - frac{1}{2}2 displaystyle sum_{i=1}^n (frac {y_i - x_ieta}{sigma_{ii}})frac{x_i}{sigma_{ii}}\
& = - displaystyle sum_{i=1}^n (frac {x_iy_i - x_i^2 eta}{sigma_{ii}^2}) = 0
end{split}
\
=> displaystyle sum_{i=1}^n (frac {x_i y_i}{sigma_{ii}^2}) = displaystyle sum_{i=1}^n (frac {x_i^2 eta}{sigma_{ii}^2}) \
=> X_c^T Y_c = X_c^T X_ceta => hat eta = (X_c^T X_c)^{-1} X_c^T Y_c
end{equation}

  记(X_c = (frac{x_1}{sigma_{11}}, frac{x_2}{sigma_{22}}, cdots, frac{x_n}{sigma_{nn}})^T, Y_c = (frac{y_1}{sigma_{11}}, frac{y_2}{sigma_{22}}, cdots, frac{y_n}{sigma_{nn}})^T)
  对(G(eta, sigma_{ii}^2))求关于(sigma_{ii})的偏导有，以(sigma_{11})为例

egin{equation}
egin{split}
frac {partial G(eta, sigma_{ii}^2)}{partial sigma_{11}}
&= 0 + frac{1}{sigma_{11}} - frac{1}{2} 2 frac{(y_1 - x_1 eta)^2}{ sigma_{11}^3} \
& = frac{1}{sigma_{11}} - frac{(y_1 - x_1 eta)^2}{ sigma_{11}^3} = 0
end{split}
\
=> frac{1}{sigma_{11}} = frac{(y_1 - x_1 eta)^2}{ sigma_{11}^3}
=> hat sigma_{11}^2 = (y_1 - x_1 eta)^2
end{equation}

  类似地，也就有(hat sigma_{ii}^2 = (y_i - x_ieta)^2)

3.3. 基于假设3

  如果满足假设3，(cov(varepsilon) = Sigma)， 并加上一个正态性的假设，即有(varepsilon)满足多维正态分布，(varepsilon sim N_n(0, sigma^2_{ii}))，那么，(Y = Xeta + varepsilon sim N_n(Xeta, Sigma))，那么有似然函数

egin{equation}
egin{split}
L(eta, Sigma Y, X) & =P(Y_1 = y_1, Y_2 = y_2, cdots, Y_n = y_n) = P(Y=y)\
& = frac{1}{(sqrt{2pi})^n |Sigma|^{frac{1}{2}}}e ^{- frac{1}{2}(Y - Xeta)^T sum^{-1} (Y - Xeta)}
end{split}
end{equation}

  其中，(|Sigma|)是(Sigma)的行列式
  我们可以发现基于假设3下，似然函数的核同样也发生了变化。那么，基于这种假设，此时采用的损失函数应该是((y - xeta)^T Sigma^{-1} (y - xeta))。将似然函数对数化：

[lnL(eta, Sigma, Y, X) = -nln(sqrt{2pi})- frac{1}{2} ln|Sigma| - frac{1}{2} (Y - Xeta)^T (Sigma)^{-1} (Y - Xeta) ]

  记(G(eta, Sigma) = nln(sqrt{2pi}) + frac{1}{2}ln|Sigma| + frac{1}{2} (Y - Xeta)^T Sigma^{-1} (Y - Xeta))，令似然函数最大化，即是求(min hspace{1mm}G(eta, Sigma))
  对(G(eta, Sigma))求关于(eta)的偏导有

egin{equation}
egin{split}
frac {partial G(eta, Sigma)}{partial eta}
&= 0 + 0 - frac{1}{2}2 X^T Sigma^{-1} (Y - Xeta)\
& = X^T Sigma^{-1}(Xeta - Y) = 0
end{split}
\
=> X^T Sigma^{-1}Xeta = X^T Sigma^{-1}Y \
=> hat eta = (X^T Sigma^{-1} X)^{-1} X^T Sigma^{-1} Y
end{equation}

  对(G(eta, Sigma))求关于(Sigma)的偏导有

egin{equation}
egin{split}
mathrm{d}G & = frac{1}{2} |Sigma|^{-1} d |Sigma| - frac{1}{2}(Y - X eta)^T Sigma^{-1}d Sigma Sigma^{-1}(Y-X eta)\
& = frac{1}{2}tr(Sigma^{-1} d Sigma) - tr(frac{1}{2}(Y - X eta)^T Sigma^{-1} d Sigma Sigma^{-1}(Y-X eta))\
& = frac{1}{2}tr(Sigma^{-1} d Sigma) - tr(frac{1}{2}Sigma^{-1}(Y-X eta)(Y - Xeta)^T Sigma^{-1} d Sigma)\
& = tr(frac{1}{2}((Sigma^{-1} - Sigma^{-1}(Y-X eta)(Y - X eta)^T Sigma^{-1}))d Sigma)
end{split}
\
=> frac{partial G}{partial Sigma} = frac{1}{2}(Sigma^{-1} - Sigma^{-1} (Y-X eta)(Y - X eta)^T Sigma^{-1} )^T = 0\
=> Sigma^{-1}(Y-X eta)(Y - X eta)^T Sigma^{-1} = Sigma^{-1} \
=> hat Sigma = (Y-X eta)(Y - X eta)^T
end{equation}

4. 估计的优良性

  在基于损失函数的估计中，我们讨论了估计的优良性，那么当换了假设和损失函数后，我们的估计是否还是具有优良的性质呢
  对于假设3中，有
egin{equation}
egin{split}
L_3(eta) & = (Y - Xeta)^T Sigma^{-1} (Y - Xeta) \
& = (Y - Xeta)^T Sigma^{- frac{1}{2}} Sigma^{-frac{1}{2}} (Y - Xeta)\
& = (Sigma^{-frac{1}{2}}Y - Sigma^{-frac{1}{2}}X eta)^T( Sigma^{- frac{1}{2}}Y - Sigma^{-frac{1}{2}}X eta)\
& = (Y^* - X^* eta)^T (Y^* - X^* eta)
end{split}
end{equation}

  其中，记(Sigma^{-frac{1}{2}}Y - Sigma^{-frac{1}{2}}Xeta)为(Y^* - X^* eta)，由于(L_1(eta) = (Y-Xeta)^T(Y - Xeta))具有优良的性质，那么(L_3(eta) = (Y^* - X^* eta)^T(Y^* - X^* eta))的估计也应该具有优良的性质。

5. 假设的场景

  为什么总假设线性模型符合假设1呢？实际上当我们基于假设2时，要估计的参数有n+p个(n个不同的(sigma_{ii})，和p个(eta_i))，而我们只有n个样本，这样就出现自由度不足的情况；而当我们基于假设3时，要估计的参数就更多了（有(frac{n^2 + n}{2}+p)个）。这样基本很难做估计，即使是做出出来了，估计也不一定唯一。

  面对这种情况，通常我们都要加大样本量，像可以一个个体测m次，得到mn个数据，当然这时模型也变成了混合模型。因此，对于假设2和假设3，更加适合一些纵向数据（经济上的面板数据、心理学上的重复测量数据、社会学上的多水平数据）