深入理解线性模型(二)---基于似然函数的估计
更新时间:2019.10.31
1. 引言
在上一篇中,我们从损失函数的角度出发讨论了\(\beta\)和\(\sigma\)的估计。在本篇将换一种极具统计味道的角度,从似然函数出发来讨论了\(\beta\)和\(\sigma\)的估计。从中我们也将看见,在不同的假设中,损失函数将会发生不同的变化。
2. 关于\(\varepsilon\)假设
在上一篇(基于损失函数的估计)中,我们提到,对于线性模型,我们常常使用Guass-Markov假设,即:
- \(E(\varepsilon) = 0\)
- \(cov(\varepsilon) = \sigma^2 I_n\)
但是,实际上我们同方差的假设是总是不满足的,完整来说,对\(\varepsilon\)的假设应该有三种:
- 同方差,且各个随机误差变量不相关:\(cov(\varepsilon) = \sigma^2 I_n\)
- 异常差,但各个随机误差变量不相关,\(cov(\varepsilon) = diag(\sigma_1^2, \sigma_2^2, \cdots, \sigma_n^2)\)
- 异方差,且各个随机误差变量是相关的,
此时,记\(cov(\varepsilon) = \Sigma\)
3. 基于似然函数的估计
之前是从损失函数的角度进行参数的估计,但是实际上每个损失函数都应该对应着一个分布,并使得分布的似然函数达到最大
我们知道在X给定的情况下,似然函数\(L(\theta;Y,X) = P_{\theta}(Y_1 = y_1, Y_2 = y_2, \cdots, Y_n = y_n)\)。假设\(Y_1, Y_2, \cdots, Y_n\)是独立的,有\(L(\theta;Y,X) = \prod_{i=1}^nP(Y = y_i)\)。当是离散情况的时候,可以进一步化为:\(L(\theta;Y,X) = \prod_{i=1}^nP_i(\theta)\)。当是连续情况的时候,则可以化为:\(L(\theta;Y,X) = \prod_{i=1}^n f(y_i;\theta)\)
3.1 基于假设1
如果满足假设1,\(cov(\varepsilon) = \sigma^2 I_n\), 并加上一个正态性的假设,即有\(\varepsilon_i \sim N(0, \sigma^2)\),那么,\(y_i = x_i\beta + \varepsilon_i \sim N(x_i\beta, \sigma^2)\),那么有似然函数:
\begin{equation}
\begin{split}
L(\beta, \sigma^2, Y, X) & = \prod_{i=1}^n f(y_i)\\
& = \prod_{i=1}^n \frac{1}{\sqrt{2\pi}\sigma} e^{- \frac{(y_i - x_i \beta)^2}{2 \sigma^2}}\\
& = (\frac{1}{\sqrt{2\pi}\sigma})^n e^{- \frac{1}{2 \sigma^2} \displaystyle \sum_{i=1}^n(y_i - x_i\beta)^2}
\end{split}
\end{equation}
可以看到,似然函数中含有的\(\sum_{i=1}^n(y_i - x_i\beta)^2\)部分正是我们之前讨论的二次损失形式。那么我们便了解到,基于假设1时,确实是应该采用我们之前所使用的二次损失形式
通常为了简便计算,我们都会将似然函数对数化
\begin{equation}
\begin{split}
lnL(\beta, \sigma^2, Y, X) & = -nln(\sqrt{2\pi}\sigma)- \frac{1}{2 \sigma^2} \sum_{i=1}^n(y_i - x_i\beta)^2
\end{split}
\end{equation}
记\(G(\beta, \sigma^2) = nln(\sqrt{2\pi}\sigma) + \frac{1}{2 \sigma^2} \sum_{i=1}^n(y_i - x_i\beta)^2\),令似然函数最大化,即是求\(min \hspace{1mm}G(\beta, \sigma^2)\)
对\(G(\beta, \sigma^2)\)求关于\(\beta\)的偏导有
\begin{equation}
\begin{split}
\frac {\partial G(\beta, \sigma^2)}{\partial \beta}
&= 0 + \frac{1}{2 \sigma^2}2 \displaystyle \sum_{i=1}^n (y_i - x_i \beta)x_i\\
& = \frac{1}{2 \sigma^2} \displaystyle \sum_{i=1}^n 2(x_i y_i - x_i^2 \beta) = 0
\end{split}
\\
=> \displaystyle \sum_{i=1}^n (x_i y_i - x_i^2 \beta) = 0 => \displaystyle \sum_{i=1}^n x_iy_i = \displaystyle \sum_{i=1}^n x_i^2 \beta\\
=> X^T Y = X^T X \beta => \hat \beta = (X^T X)^{-1} X^T Y
\end{equation}
对\(G(\beta, \sigma^2)\)求关于\(\sigma\)的偏导有
\begin{equation}
\begin{split}
\frac {\partial G(\beta, \sigma^2)}{\partial \sigma}
&= n\frac{1}{\sqrt{2\pi} \sigma} \sqrt{2\pi} - \frac{2}{2 \sigma^3} \sum_{i=1}^n(y_i - x_i \beta)^2 \\
& = \frac{n}{\sigma} + \frac{1}{\sigma^3} \sum_{i=1}^n(y_i - x_i \beta)^2 = 0
\end{split}
\\
=> \frac{1}{\sigma^3} \sum_{i=1}^n(y_i - x_i \beta)^2 = \frac{n}{\sigma}
=> \hat \sigma^2 = \frac{\displaystyle \sum_{i=1}^n(y_i - x_i\beta)^2}{n}
\end{equation}
从这里便可以看出,通过似然函数,一次就搞定了参数\(\beta\)和\(\sigma\)的估计,而基于损失函数的估计只是估计出了\(\beta\),而\(\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_i\beta + \varepsilon_i \sim N(x_i\beta, \sigma^2_{ii})\),那么有似然函数:
\begin{equation}
\begin{split}
L(\beta, \sigma^2, Y, X) & = \prod_{i=1}^n f(y_i)\\
& = \prod_{i=1}^n \frac{1}{\sqrt{2\pi} \sigma_{ii}} e^{- \frac{(y_i - x_i \beta)^2}{2 \sigma^2_{ii}}}\\
& = (\frac{1}{\sqrt{2\pi}})^n \prod_{i=1}^n(\frac{1}{\sigma_{ii}}) e^{- \frac{1}{2} \displaystyle \sum_{i=1}^n(\frac {y_i - x_i \beta}{\sigma_{ii}})^2}
\end{split}
\end{equation}
我们可以发现基于假设2下,似然函数的核心部分发生了变化,不再是\(\sum_{i=1}^n(y_i - x_i\beta)^2\)。因此,根据之前的经验,基于假设2,所采用的损失函数也应该发生变化。此时采用的损失函数应该是标准化的二次损失\(\displaystyle \sum_{i=1}^n(\frac {y_i - x_i\beta}{\sigma_{ii}})^2\),我们也把这称为加权最小二乘估计。
将似然函数对数化:
\begin{equation}
\begin{split}
lnL(\beta, \sigma^2, Y, X) = -nln(\sqrt{2\pi})- \sum_{i=1}^nln\sigma_{ii} - \frac{1}{2} \displaystyle \sum_{i=1}^n(\frac {y_i - x_i\beta}{\sigma_{ii}})^2
\end{split}
\end{equation}
记\(G(\beta, \sigma_{ii}^2) = nln(\sqrt{2\pi}) + \sum_{i=1}^nln\sigma_{ii} + \frac{1}{2} \displaystyle \sum_{i=1}^n(\frac {y_i - x_i\beta}{\sigma_{ii}})^2\),令似然函数最大化,即是求\(min \hspace{1mm}G(\beta, \sigma_{ii}^2)\)
对\(G(\beta, \sigma_{ii}^2)\)求关于\(\beta\)的偏导有
\begin{equation}
\begin{split}
\frac {\partial G(\beta, \sigma_{ii}^2)}{\partial \sigma_{ii}}
&= 0 + 0 - \frac{1}{2}2 \displaystyle \sum_{i=1}^n (\frac {y_i - x_i\beta}{\sigma_{ii}})\frac{x_i}{\sigma_{ii}}\\
& = - \displaystyle \sum_{i=1}^n (\frac {x_iy_i - x_i^2 \beta}{\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 \beta}{\sigma_{ii}^2}) \\
=> X_c^T Y_c = X_c^T X_c\beta => \hat \beta = (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(\beta, \sigma_{ii}^2)\)求关于\(\sigma_{ii}\)的偏导有,以\(\sigma_{11}\)为例
\begin{equation}
\begin{split}
\frac {\partial G(\beta, \sigma_{ii}^2)}{\partial \sigma_{11}}
&= 0 + \frac{1}{\sigma_{11}} - \frac{1}{2} 2 \frac{(y_1 - x_1 \beta)^2}{ \sigma_{11}^3} \\
& = \frac{1}{\sigma_{11}} - \frac{(y_1 - x_1 \beta)^2}{ \sigma_{11}^3} = 0
\end{split}
\\
=> \frac{1}{\sigma_{11}} = \frac{(y_1 - x_1 \beta)^2}{ \sigma_{11}^3}
=> \hat \sigma_{11}^2 = (y_1 - x_1 \beta)^2
\end{equation}
类似地,也就有\(\hat \sigma_{ii}^2 = (y_i - x_i\beta)^2\)
3.3. 基于假设3
如果满足假设3,\(cov(\varepsilon) = \Sigma\), 并加上一个正态性的假设,即有\(\varepsilon\)满足多维正态分布,\(\varepsilon \sim N_n(0, \sigma^2_{ii})\),那么,\(Y = X\beta + \varepsilon \sim N_n(X\beta, \Sigma)\),那么有似然函数
\begin{equation}
\begin{split}
L(\beta, \Sigma Y, X) & =P(Y_1 = y_1, Y_2 = y_2, \cdots, Y_n = y_n) = P(Y=y)\\
& = \frac{1}{(\sqrt{2\pi})^n |\Sigma|^{\frac{1}{2}}}e ^{- \frac{1}{2}(Y - X\beta)^T \sum^{-1} (Y - X\beta)}
\end{split}
\end{equation}
其中,\(|\Sigma|\)是\(\Sigma\)的行列式
我们可以发现基于假设3下,似然函数的核同样也发生了变化。那么,基于这种假设,此时采用的损失函数应该是\((y - x\beta)^T \Sigma^{-1} (y - x\beta)\)。将似然函数对数化:
记\(G(\beta, \Sigma) = nln(\sqrt{2\pi}) + \frac{1}{2}ln|\Sigma| + \frac{1}{2} (Y - X\beta)^T \Sigma^{-1} (Y - X\beta)\),令似然函数最大化,即是求\(min \hspace{1mm}G(\beta, \Sigma)\)
对\(G(\beta, \Sigma)\)求关于\(\beta\)的偏导有
\begin{equation}
\begin{split}
\frac {\partial G(\beta, \Sigma)}{\partial \beta}
&= 0 + 0 - \frac{1}{2}2 X^T \Sigma^{-1} (Y - X\beta)\\
& = X^T \Sigma^{-1}(X\beta - Y) = 0
\end{split}
\\
=> X^T \Sigma^{-1}X\beta = X^T \Sigma^{-1}Y \\
=> \hat \beta = (X^T \Sigma^{-1} X)^{-1} X^T \Sigma^{-1} Y
\end{equation}
对\(G(\beta, \Sigma)\)求关于\(\Sigma\)的偏导有
\begin{equation}
\begin{split}
\mathrm{d}G & = \frac{1}{2} |\Sigma|^{-1} d |\Sigma| - \frac{1}{2}(Y - X \beta)^T \Sigma^{-1}d \Sigma \Sigma^{-1}(Y-X \beta)\\
& = \frac{1}{2}tr(\Sigma^{-1} d \Sigma) - tr(\frac{1}{2}(Y - X \beta)^T \Sigma^{-1} d \Sigma \Sigma^{-1}(Y-X \beta))\\
& = \frac{1}{2}tr(\Sigma^{-1} d \Sigma) - tr(\frac{1}{2}\Sigma^{-1}(Y-X \beta)(Y - X\beta)^T \Sigma^{-1} d \Sigma)\\
& = tr(\frac{1}{2}((\Sigma^{-1} - \Sigma^{-1}(Y-X \beta)(Y - X \beta)^T \Sigma^{-1}))d \Sigma)
\end{split}
\\
=> \frac{\partial G}{\partial \Sigma} = \frac{1}{2}(\Sigma^{-1} - \Sigma^{-1} (Y-X \beta)(Y - X \beta)^T \Sigma^{-1} )^T = 0\\
=> \Sigma^{-1}(Y-X \beta)(Y - X \beta)^T \Sigma^{-1} = \Sigma^{-1} \\
=> \hat \Sigma = (Y-X \beta)(Y - X \beta)^T
\end{equation}
4. 估计的优良性
在基于损失函数的估计中,我们讨论了估计的优良性,那么当换了假设和损失函数后,我们的估计是否还是具有优良的性质呢
对于假设3中,有
\begin{equation}
\begin{split}
L_3(\beta) & = (Y - X\beta)^T \Sigma^{-1} (Y - X\beta) \\
& = (Y - X\beta)^T \Sigma^{- \frac{1}{2}} \Sigma^{-\frac{1}{2}} (Y - X\beta)\\
& = (\Sigma^{-\frac{1}{2}}Y - \Sigma^{-\frac{1}{2}}X \beta)^T( \Sigma^{- \frac{1}{2}}Y - \Sigma^{-\frac{1}{2}}X \beta)\\
& = (Y^* - X^* \beta)^T (Y^* - X^* \beta)
\end{split}
\end{equation}
其中,记\(\Sigma^{-\frac{1}{2}}Y - \Sigma^{-\frac{1}{2}}X\beta\)为\(Y^* - X^* \beta\),由于\(L_1(\beta) = (Y-X\beta)^T(Y - X\beta)\)具有优良的性质,那么\(L_3(\beta) = (Y^* - X^* \beta)^T(Y^* - X^* \beta)\)的估计也应该具有优良的性质。
5. 假设的场景
为什么总假设线性模型符合假设1呢?实际上当我们基于假设2时,要估计的参数有n+p个(n个不同的\(\sigma_{ii}\),和p个\(\beta_i\)),而我们只有n个样本,这样就出现自由度不足的情况;而当我们基于假设3时,要估计的参数就更多了(有\(\frac{n^2 + n}{2}+p\)个)。这样基本很难做估计,即使是做出出来了,估计也不一定唯一。
面对这种情况,通常我们都要加大样本量,像可以一个个体测m次,得到mn个数据,当然这时模型也变成了混合模型。因此,对于假设2和假设3,更加适合一些纵向数据(经济上的面板数据、心理学上的重复测量数据、社会学上的多水平数据)