线性回归系数的几个性质

摘要

我们知道一元和多元线性回归系数都有解析解，本文将简要介绍总结线性回归系数的几个常见的性质。

线性回归问题的描述

我们回忆一下，单变量线性回归问题是指，给定了 $n$ 个观察量 $(x_1,y_1),(x_2,y_2),\cdots ,(x_n,y_n)$。我们希望用一个线性的关系 $y={\beta }_{1}x+{\beta }_0+ϵ$ 来描述这些观察量的规律。

这里，我们把方程
$$y={\beta }_{1}x+{\beta }_0+ϵ$$

称为 总体回归模型 (population regression model)。

而当给定的 $n$ 个观察量 $(x_1,y_1),(x_2,y_2),\cdots ,(x_n,y_n)$，我们称方程
$$y_i={\beta }_1{x}_i+{\beta }_0+{ϵ}_i$$
为 样本回归模型 （sample regression model）。
其中 ${ϵ}_i$ 为误差项 (error)，每个 ${ϵ}_i$ 均是一个随机变量，独立且均服从 一个均值为0, 方差为 ${\sigma }^2$ 的概率分布。其中 ${\sigma }^2$ 为误差项的方差，在我们所考虑的线性模型中，我们已知 ${\sigma }^2$ 是固定的，但是我们不知道其具体的数值。

所以我们去做线性回归“拟合”模型的参数 ${\beta }_1$ 和 ${\beta }_0$ 时，实际上是根据 $(x_1,y_1),(x_2,y_2),\cdots ,(x_n,y_n)$ 去做点估计 (point estimator)。

这里值得注意的是，我们把 $x_1,x_2,\cdots ,x_n$ 当作给定的值，即可以认为是常量，而每一个 $y_i,i=1,2,\cdots n$ 均是一个随机变量。

所以， $y_i$ 作为随机变量，可以用下图来体现。

上图中，$x$ 轴上的绿色的圆点表示固定的 $x_i$ ，而橙色的直线表示真实的线性关系，而蓝色的钟型曲线代表随机变量 $y_i$ 服从的概率密度函数。

单变量线性回归系数的公式

我们回顾单变量线性回归问题的公式，有
$$\{\begin{array}{ll}& \widehat{{\beta }_1}=\frac{\sum (x_i-\bar{x})(y_i-\bar{y})}{\sum (x_i-\bar{x}{)}^2}\\ & \widehat{{\beta }_0}=\bar{y}-\widehat{{\beta }_1}\bar{x}\end{array}$$

即我们对参数 ${\beta }_1,{\beta }_0$ 的估计是 $\widehat{{\beta }_1}=\frac{\sum (x_i-\bar{x})(y_i-\bar{y})}{\sum (x_i-\bar{x}{)}^2}$，对参数 ${\beta }_0$ 的估计是 $\widehat{{\beta }_0}=\bar{y}-\widehat{{\beta }_1}\bar{x}$。

为方便起见，我们记

$S_{xx}=\sum _{i=1}^n(x_i-\bar{x}{)}^2$,
$S_{xy}=\sum _{i=1}^n(x_i-\bar{x})(y_i-\bar{y})$。

我们可以进一步简化 $S_{xy}$ 为 $S_{xy}=\sum _{i=1}^n(x_i-\bar{x})y_i$。

那么我们可以将 $\widehat{{\beta }_1}$ 与 $\widehat{{\beta }_0}$ 写成：
$\widehat{{\beta }_1}=\sum _{i=1}^n\frac{(x_i-\bar{x})}{S_{xx}}\cdot y_i$,
$\widehat{{\beta }_0}=\bar{y}-\widehat{{\beta }_1}\bar{x}$。

回忆起之前我们提到的，每一个 $y_i,i=1,2,\cdots n$ 均是一个随机变量，这里我们分别对 ${\beta }_1$ 和 ${\beta }_0$ 的估计 $\widehat{{\beta }_1}$ 与 $\widehat{{\beta }_0}$ 就是 $y_i,i=1,2,\cdots n$ 的函数。

无偏估计

首先，我们证明上述估计 $\widehat{{\beta }_1}$ 与 $\widehat{{\beta }_0}$ 是无偏估计 (unbiased)。

要证明我们的估计是无偏估计，我们须要证明我们估计的期望恒等于所估计的参数，即我们须要证明：
$$\{\begin{array}{ll}& \mathbb{E}[\widehat{{\beta }_1}]={\beta }_1\\ & \\ & \mathbb{E}[\widehat{{\beta }_0}]={\beta }_0\end{array}$$

证明过程十分直接，我们直接将上一节的表达式代入。

$\mathbb{E}[\widehat{{\beta }_1}]=\mathbb{E}[\sum _{i=1}^n\frac{(x_i-\bar{x})}{S_{xx}}\cdot y_i]$，注意到 $y_i={\beta }_1{x}_i+{\beta }_0+{ϵ}_i$，我们有 $\mathbb{E}(y_i)={\beta }_1{x}_i+{\beta }_0$。这里我们用到了 $\mathbb{E}[{ϵ}_i]=0$。

代入，我们有

$\begin{array}{rl}\mathbb{E}[\widehat{{\beta }_1}]& =\mathbb{E}[\sum _{i=1}^n\frac{(x_i-\bar{x})}{S_{xx}}\cdot y_i]\\ & =\sum _{i=1}^n\frac{(x_i-\bar{x})}{S_{xx}}\mathbb{E}[y_i]\\ & =\sum _{i=1}^n\frac{(x_i-\bar{x})}{S_{xx}}\cdot ({\beta }_1{x}_i+{\beta }_0)\\ & ={\beta }_1\cdot [\sum _{i=1}^n\frac{(x_i-\bar{x})}{S_{xx}}\cdot x_i]+{\beta }_0\cdot \sum _{i=1}^n\frac{(x_i-\bar{x})}{S_{xx}}\\ & ={\beta }_1\end{array}$。

注意，上式中，我们用到了 $\sum _{i=1}^n\frac{(x_i-\bar{x})}{S_{xx}}=0$，以及 $\sum _{i=1}^n(x_i-\bar{x})\cdot x_i=S_{xx}$。

同样的，我们有，
$\begin{array}{rl}\mathbb{E}(\widehat{{\beta }_0})& =\mathbb{E}[\bar{y}-\bar{x}\cdot \widehat{{\beta }_1}]\\ & =\frac{1}{n}\sum _{i=1}^n({\beta }_0+{\beta }_1\cdot x_i)-\bar{x}\cdot {\beta }_1\\ & ={\beta }_0\end{array}$。

所以，$\widehat{{\beta }_0}$ 也是 ${\beta }_0$ 的无偏估计。

线性回归系数的方差

下面我们来看线性回归系数 $\widehat{{\beta }_1}$ 和 $\widehat{{\beta }_0}$ 的方差。直观上看，当误差项的方差 ${\sigma }^2$ 越大时，$\mathrm{Var}(\widehat{{\beta }_1})$ 与 $\mathrm{Var}(\widehat{{\beta }_0})$ 也应该越大。因为误差项的方差越大，随机变量 $y_i$ 的取值就会越不确定，从而使得线性拟合的直线不确定性越大。

具体地，我们有
$\mathrm{Var}(\widehat{{\beta }_1})=\frac{{\sigma }^2}{S_{xx}}$,
$\mathrm{Var}(\widehat{{\beta }_0})={\sigma }^2\left(\frac{1}{n}+\frac{{\bar{x}}^2}{S_{xx}}\right)$。

对于 $\widehat{{\beta }_1}$，我们可以直接计算如下，
因为 $\widehat{{\beta }_1}=\sum _{i=1}^n\frac{(x_i-\bar{x})}{S_{xx}}\cdot y_i$, 我们记 $c_i=\frac{(x_i-\bar{x})}{S_{xx}}$。于是，

$\begin{array}{rl}\mathrm{Var}(\widehat{{\beta }_1})& =\sum _{i=1}^n\left(c_i^2\cdot \mathrm{Var}(y_i)\right)=\sum _{i=1}^n\frac{(x_i-\bar{x}{)}^2}{S_{xx}^2}\cdot \mathrm{Var}(y_i)\\ & =\frac{{\sigma }^2}{S_{xx}^2}\sum _{i=1}^n(x_i-\bar{x}{)}^2=\frac{{\sigma }^2}{S_{xx}}。\end{array}$

下面我们来看 $\mathrm{Var}(\widehat{{\beta }_0})$。

因为 $\widehat{{\beta }_0}=\bar{y}-\bar{x}\cdot \widehat{{\beta }_1}$，我们有 $\mathrm{Var}(\widehat{{\beta }_0})=\mathrm{Var}(\bar{y})+\mathrm{Var}(\bar{x}\cdot \widehat{{\beta }_1})-2\mathrm{Cov}(\bar{y},\bar{x}\cdot \widehat{{\beta }_1})$。

我们知道，$\mathrm{Var}(\bar{y})=\frac{1}{n^2}\sum _{i=1}^n\mathrm{Var}(y_i)=\frac{1}{n^2}\sum _{i=1}^n\mathrm{Var}({ϵ}_i)=\frac{{\sigma }^2}{n}$。
而 $\mathrm{Var}(\bar{x}\cdot \widehat{{\beta }_1})={\bar{x}}^2\cdot \mathrm{Var}(\widehat{{\beta }_1})={\bar{x}}^2\cdot \frac{{\sigma }^2}{S_{xx}}$。
下面我们来看协方差项 $\mathrm{Cov}(\bar{y},\bar{x}\cdot \widehat{{\beta }_1})=\bar{x}\cdot \mathrm{Cov}(\bar{y},\widehat{{\beta }_1})$。

$\begin{array}{rl}\mathrm{Cov}(\bar{y},\widehat{{\beta }_1})& =\mathrm{Cov}(\frac{1}{n}\sum _{i=1}^n{y}_i,\sum _{i=1}^n\frac{x_i-\bar{x}}{S_{xx}}{y}_i)\\ & =\frac{1}{n}\sum _{i=1}^n\mathrm{Cov}(y_i,\frac{x_i-\bar{x}}{S_{xx}}{y}_i)\\ & =\frac{1}{n}\sum _{i=1}^n\frac{x_i-\bar{x}}{S_{xx}}\cdot {\sigma }^2\\ & =\frac{{\sigma }^2}{n}\frac{x_i-\bar{x}}{S_{xx}}=0。\end{array}$
于是，我们有，

$\mathrm{Var}(\widehat{{\beta }_0})=\mathrm{Var}(\bar{y})+\mathrm{Var}(\bar{x}\cdot \widehat{{\beta }_1})-2\mathrm{Cov}(\bar{y},\bar{x}\cdot \widehat{{\beta }_1})=\frac{{\sigma }^2}{n}+{\bar{x}}^2\cdot \frac{{\sigma }^2}{S_{xx}}={\sigma }^2\left(\frac{1}{n}+\frac{{\bar{x}}^2}{S_{xx}}\right)$。

其余的几个性质

残差项之和为0

我们记残差项为 $e_i$，$e_i=y_i-\widehat{y_i}$。即 $\sum _{i=1}^n(y_i-\widehat{y_i})=0$。

这里我们可以用求 $\widehat{{\beta }_1}$ 和 $\widehat{{\beta }_0}$ 的公式来证明。

在文章 单变量线性回归的最小二乘法公式 中，我们提到在用偏导数求 $\widehat{{\beta }_1}$ 和 $\widehat{{\beta }_0}$ 时，我们有$\widehat{{\beta }_1}$ 和 $\widehat{{\beta }_0}$ 的表达式如下：

$$\{\begin{array}{ll}& \frac{\partial \text{ RSS}}{\partial {\beta }_0}=2n{\beta }_0+2\sum _{i=1}^n{x}_i{\beta }_1-2\sum _{i=1}^n{y}_i=0\\ & \\ & \frac{\partial \text{ RSS}}{\partial {\beta }_1}=2\sum _{i=1}^n{x}_i^2{\beta }_1+2\sum _{i=1}^n{x}_i{\beta }_0-2\sum _{i=1}^n{x}_i{y}_i=0\end{array}$$

根据第一个式子，我们把 $\widehat{{\beta }_1}$ 和 $\widehat{{\beta }_0}$ 代入，我们有

$n\widehat{{\beta }_0}+\sum _{i=1}^n{x}_i\widehat{{\beta }_1}-\sum _{i=1}^n{y}_i=0$。
即，$\sum _{i=1}^n\left(\widehat{{\beta }_0}+x_i\widehat{{\beta }_1}\right)-\sum _{i=1}^n{y}_i=0$，亦 $\sum _{i=1}^n(y_i-\widehat{y_i})=0$。

线性拟合直线总会经过 $(\bar{x},\bar{y})$ 这个点

拟合直线为 $y=\widehat{{\beta }_0}+\widehat{{\beta }_1}\cdot x$。而我们有
$\widehat{{\beta }_0}=\bar{y}-\widehat{{\beta }_1}\bar{x}$，所以 线性拟合直线总会经过 $(\bar{x},\bar{y})$ 这个点。

在 $x_i$ 权重下，残差和为0，即$\sum _{i=1}^n{x}_i{e}_i=0$

我们可以直接把 $e_i=y_i-\widehat{y_i}$ 代入。注意到 $\widehat{y_i}=\widehat{{\beta }_0}+\widehat{{\beta }_1}\cdot x_i$，我们有

$\begin{array}{rl}\sum _{i=1}^n{x}_i{e}_i& =\sum _{i=1}^n{x}_i\cdot (y_i-\widehat{y_i})=\sum _{i=1}^n{x}_i\cdot (y_i-\widehat{{\beta }_0}-\widehat{{\beta }_1}\cdot x_i)\\ & =\sum _{i=1}^n{x}_i{y}_i-\widehat{{\beta }_0}\sum _{i=1}^n{x}_i-\widehat{{\beta }_1}\sum _{i=1}^n{x}_i^2\\ & =\sum _{i=1}^n{x}_i{y}_i-\widehat{{\beta }_1}\sum _{i=1}^n{x}_i^2-(\bar{y}-\widehat{{\beta }_1}\cdot \bar{x})\cdot \sum _{i=1}^n{x}_i\\ & =\sum _{i=1}^n{x}_i{y}_i-\widehat{{\beta }_1}\sum _{i=1}^n{x}_i^2-n\bar{x}\bar{y}+n\widehat{{\beta }_1}{\bar{x}}^2\\ & =\left(\sum _{i=1}^n{x}_i{y}_i-n\bar{x}\bar{y}\right)-\left(\sum _{i=1}^n{x}_i^2-n{\bar{x}}^2\right)\widehat{{\beta }_1}\\ & =S_{xy}-S_{xx}\widehat{{\beta }_1}\\ & =0\end{array}$

故得证。

在 $\widehat{y_i}$ 权重下，残差和为0，即$\sum _{i=1}^n\widehat{y_i}{e}_i=0$

这里我们只须要利用上一个公式，即 $\sum _{i=1}^n{x}_i{e}_i=0$ 即可。

因为我们有 $\sum _{i=1}^n\widehat{y_i}{e}_i=\sum _{i=1}^n(\widehat{{\beta }_0}+\widehat{{\beta }_1}\cdot x_i)e_i=\widehat{{\beta }_0}\sum _{i=1}^n{e}_i+\widehat{{\beta }_1}\sum _{i=1}^n{x}_i{e}_i=0$。

模拟

最后我们用 python 程序来模拟“证明” $\widehat{{\beta }_1}$ 和 $\widehat{{\beta }_0}$ 分别时 ${\beta }_1$ 和 ${\beta }_0$ 的无偏估计。

class unbiased_beta:
    
    def __init__(self, arr_x: np.array, beta1: float, beta0: float, epsilon: float):
        #self.N = N
        self.X = arr_x
        self.beta1 = beta1
        self.beta0 = beta0
        self.epsilon = epsilon
        
        self.Sxx = ((self.X - self.X.mean()) ** 2).sum()
        self.X_bar = self.X.mean()
        
    def estimate_beta(self, N: int) -> tuple:
        res_beta1, res_beta0 = [], []
        for i in range(N):
            #print(i)
            cur_error = np.random.normal(0, self.epsilon, arr_x.shape)
            cur_y = self.beta0 + self.beta1 * self.X + cur_error
            cur_y_bar = cur_y.mean()
            Sxy = ((self.X - self.X.mean()) * (cur_y - cur_y_bar)).sum()
            cur_beta1 = Sxy / self.Sxx
            cur_beta0 = cur_y_bar - cur_beta1 * self.X_bar
            res_beta1.append(cur_beta1)
            res_beta0.append(cur_beta0)
        return np.mean(res_beta1), np.mean(res_beta0)

arr_x = np.array(range(1, 11))
a = unbiased_beta(arr_x, 2, 3, 1)
res = a.estimate_beta(10 ** 5)

res

(1.9988026861047237, 3.0029188805679303)

可以发现，在经过多次的实验之后，我们得到的 $\widehat{{\beta }_1}$ 和 $\widehat{{\beta }_0}$ 的平均值是非常接近真实值 ${\beta }_1$ 和 ${\beta }_0$ 的。

plt.figure(figsize=(8, 6), dpi=100)
plt.hist(res[0], bins=50, density=True);
line_vert = [[2, c] for c in np.linspace(0, 4, 100)]
plt.plot([c[0] for c in line_vert], [c[1] for c in line_vert], '-', linewidth=4)
plt.xlabel("estimated beta1 value", fontsize=20)
plt.ylabel("count", fontsize=20)
plt.xticks(fontsize=20)
plt.yticks(fontsize=20)

多元线性回归的系数也是无偏估计

上文中我们介绍了单变量线性回归中，$\widehat{{\beta }_1}$ 与 $\widehat{{\beta }_0}$ 是 ${\beta }_1$ 与 ${\beta }_0$ 的无偏估计 (unbiased) 。下面我们来看多变量线性回归中系数估计的无偏性质。

根据多变量线性回归的 normal equation，我们有

$\hat{\beta }=(X^{T}X{)}^{-1}{X}^T\mathbf{y}$。

注意到，$\mathbf{y}=X\beta +ϵ$。
代入 $\hat{\beta }$ 的表达式，我们有

$\begin{array}{rl}\mathbb{E}[\hat{\beta }]& =\mathbb{E}[(X^{T}X{)}^{-1}{X}^T\mathbf{y}]=\mathbb{E}[(X^{T}X{)}^{-1}{X}^T(X\beta +ϵ)]\\ & =\mathbb{E}[(X^{T}X{)}^{-1}{X}^{T}X\beta ]+\mathbb{E}[(X^{T}X{)}^{-1}{X}^Tϵ]\\ & =\beta \end{array}$。

于是，我们知道，对于多元线性回归，$\hat{\beta }$ 也是 $\beta$ 的无偏估计。

我们用程序验证如下。

class unbiased_beta_mv:
    
    def __init__(self, arr_x: np.array, beta: np.array, epsilon: float):
        #self.N = N
        self.X = arr_x
        self.beta = beta
        self.epsilon = epsilon
        
    def estimate_beta(self, N: int) -> tuple:
        res_beta = []
        for i in range(N):
            #print(i)
            cur_error = np.random.normal(0, self.epsilon, arr_x.shape[0])
            cur_y = self.X @ self.beta + cur_error.reshape((-1, 1))
            cur_beta_hat = np.linalg.inv(self.X.T @ self.X) @ self.X.T @ cur_y
            res_beta.append(cur_beta_hat)
        return np.hstack(res_beta)

beta = np.array([[1], [2], [3], [4]])
arr_x = np.random.normal(0, 1, (10, 3))
arr_x = np.hstack((np.ones((10, 1)), arr_x))

a = unbiased_beta_mv(arr_x, beta, 1)

res = a.estimate_beta(10 ** 5)

np.mean(res, axis=1)

array([1.00255519, 1.99960968, 2.99919333, 4.00008326])

参考文献

[1] Introduction to linear regression analysis, Douglas C. Montgomery, Elizabeth A. Peck, G. Geoffrey Vining, John Wiley & Sons (2021)