多元回归理论 多元回归的模型

文章目录

• 模型

• 回归模型
• 估计的多元回归方程

• 最小二乘估计

• β ^ \hat{\beta} β^的性质

• 拟合优度检验

• 多重判定系数
• 估计标准误差

• 显著性检验

• 线性关系的检验
• 回归系数的检验

• 多重共线性
• 预测
• 变量选择与逐步回归
• 实战

模型

一元回归模型在二维空间中是一条直线，而多元回归模型，由于自变量的增加，它在更高维的空间中显示，因此很难将它画出来，相比一元回归模型，用归纳的方法将多元回归模型表述出来。

回归模型

多元线性回归的模型为$y=\beta_0+\beta_1x_1++\beta_2x_2+\cdots+\beta_kx_k+\varepsilon$，$\beta_0,\beta_1,\beta_2,\cdots,\beta_k$是回归系数，$\varepsilon$为误差项。

于是得到多元回归方程：
$E(y)=\beta_0+\beta_1x_1++\beta_2x_2+\cdots+\beta_kx_k$
它描述了$y$与$x_1,x_2,\cdots,x_k$之间的关系。

现有$n$个独立观测数据$(y_i , x_{i1},\cdots, x_{ik} )$ ，$i = 1,\cdots,n,n>k$，有
$y_i=\beta_0+\beta_1x_{i1}++\beta_2x_{i2}+\cdots+\beta_kx_{ik}+\varepsilon_i$

为了将模型简化，用矩阵的形式表示该模型，

记
$Y =\begin{bmatrix} y_1 \\ y_2 \\ \vdots \\ y_n \\ \end{bmatrix} , X= \begin{bmatrix} 1 & x_{11} & \cdots\ &x_{1k}\\ 1 & x_{21} & \cdots\ &x_{2k}\\ \vdots & \vdots & \ddots & \vdots \\ 1 & x_{n1} & \cdots\ & x_{nk}\\ \end{bmatrix} , \varepsilon =\begin{bmatrix} \varepsilon_1 \\ \varepsilon_2 \\ \vdots \\ \varepsilon_n \\ \end{bmatrix} , \beta =\begin{bmatrix} \beta_0 \\ \beta_1 \\ \vdots \\ \beta_k \\ \end{bmatrix}$

则有
$Y=X\beta+\varepsilon$

和一元线性回归类似，假设：

• $r(X)=k+1$，即该矩阵列满秩
• 满足$G-M$条件
• $\varepsilon\sim N(0,\sigma^2I_n)$，其中$I_n$为$n$阶单位阵

根据假设，则有
$Y \sim N(X\beta,\sigma^2I_n )$

估计的多元回归方程

估计样本的多元回归方程，总体回归参数$\beta_0,\beta_1,\beta_2,\cdots,\beta_k$是未知的，因此我们需要用样本去估计它们，其估计的回归方程为：
$\hat{y}=\hat{\beta_0}+\hat{\beta_1}x_1+\hat{\beta_2}x_2+\cdots+\hat{\beta_k}x_k$

最小二乘估计

和一元的一样，我们要估计未知的参数$\beta_0,\beta_1,\beta_2,\cdots,\beta_k$，采用最小二乘法。
在一元中，我们知道要使得拟合误差(即残差)平方和最小，即$\sum e_i^2$最小。
同样的，记
$Q(\beta_0,\beta_1,\cdots,\beta_k)=\sum_{i=1}^{n}(y_i-\beta_0-\beta_1x_{i1}-\beta_2x_{i2}- \cdots - \beta_kx_{ik})^2$
则
$\begin{aligned} Q(\hat{\beta_0},\hat{\beta_1},\cdots,\hat{\beta_k}) &=\mathop{min}\limits_{\beta_0,\beta_1,\cdots,\beta_k}{Q(\beta_0,\beta_1,\cdots,\beta_k)}\\ &=\sum_{i=1}^{n}(y_i-\hat{\beta_0}-\hat{\beta_1}x_{i1}-\cdots-\hat{\beta_k}x_{ik})^2 \end{aligned}$
这时对每个估计参数求偏导，得到以下偏导方程组，
$\begin{aligned} \begin{cases} &\frac{\partial Q}{\partial \beta_0} =-2\sum_{i=1}^{n}(y_i-\hat{\beta_0}-\hat{\beta_1}x_{i1}-\cdots-\hat{\beta_k}x_{ik})=0 \\ \\ &\frac{\partial Q}{\partial \beta_1} =-2\sum_{i=1}^{n}x_{i1}(y_i-\hat{\beta_0}-\hat{\beta_1}x_{i1}-\cdots-\hat{\beta_k}x_{ik})=0\\ &\vdots \\ &\frac{\partial Q}{\partial \beta_k} =-2\sum_{i=1}^{n}x_{ik}(y_i-\hat{\beta_0}-\hat{\beta_1}x_{i1}-\cdots-\hat{\beta_k}x_{ik})=0\\ \end{cases} \end{aligned}$

整理得到，
$\begin{cases} \hat{\beta_0}n+ \hat{\beta_1}\sum x_{i1}+\hat{\beta_2}\sum x_{i2}+\cdots+\hat{\beta_k}\sum x_{ik}=\sum y_i\\ \\ \hat{\beta_0}\sum x_{i1}+\hat{\beta_1}\sum x_{i1}^2+\hat{\beta_2}\sum x_{i1}x_{i2}+\cdots+\hat{\beta_k}\sum x_{i1}x_{ik}=\sum x_{i1}y_i\\ \vdots \\ \hat{\beta_0}\sum x_{ik}+\hat{\beta_1}\sum x_{ik}x_{i1}+\hat{\beta_2}\sum x_{ik}x_{i2}+\cdots+\hat{\beta_k}\sum x_{ik}^2=\sum x_{ik}y_k\\ \end{cases}$
由于式子太多，为了方便，用矩阵表示，即
$X^TX\hat{\beta}=X^TY$

前面我们假设$X$矩阵列满秩，则$X^TX$为可逆矩阵，那么就可以得到
$\hat{\beta}=(X^TX)^{-1}X^TY$
最后得到估计的方程
$\hat{Y}=X\hat{\beta}=X(X^TX)^{-1}X^TY$

并且由方程组我们可以知道，
$\begin{aligned} \begin{cases} \sum_{i=1}^{n} e_i=0\\ \\ \sum_{i=1}^{n}\sum_{j=1}^{k} e_ix_{ij}=0 \end{cases} \end{aligned}$

将上面的等式写成矩阵的形式，即
$e^TX=0$
其中$e^T=[e_1,e_2,\cdots,e_n]$，$X$为一开始表示的矩阵。

$\hat{\beta}$的性质

$\hat{\beta}$是$Y$的线性函数，这是显然的。并且$\hat{\beta}$是$\beta$的无偏估计，下面给出证明：
$\begin{aligned} E(\hat{\beta}) &=E((X^TX)^{-1}X^TY)\\ &=E((X^TX)^{-1}X^T(X\beta+\varepsilon))\\ &=E((X^TX)^{-1}X^TX\beta)\\ &=\beta \end{aligned}$
关于协方差矩阵的理解$\begin{aligned} Var(\hat{\beta}) &=Cov(\hat{\beta},\hat{\beta})\\ &=Cov((X^TX)^{-1}X^TY,(X^TX)^{-1}X^TY)\\ &=(X^TX)^{-1}X^TCov(Y,Y)[(X^TX)^{-1}X^T]^T\\ &=(X^TX)^{-1}X^T\sigma^2I_nX(X^TX)^{-1}\\ &=\sigma^2(X^TX)^{-1} \end{aligned}$

于是我们估计的$\hat{\beta}$
$\hat{\beta} \sim N(\beta,\sigma^2(X^TX)^{-1})$

拟合优度检验

多重判定系数

与一元回归类似，对于多元回归，需要用多重判定系数评价其拟合程度。
同样有 $SST=SSR+SSE$ 成立。
为了度量多元回归方程拟合程度，我们定义多重判定系数 ,
$R^2=\frac{SSR}{SST}=1-\frac{SSE}{SST}$

和一元一样，只要回归平方和占比越大，其拟合效果越好。

对于多重判定系数有一点特别重要的需要说明：自变量个数的增加将影响到因变量中被估计的回归方程所解释的比例。当增加自变量时，会使预测误差变得较小，从而减小残差平方和$SSE$，那么$SSR$就会变大。从而使$R^2$变大。
这就会引发一个问题，如果模型中增加一个自变量，即使这个自变量在统计上并不显著， $R^2$也会变大。为了避免这个问题。统计学家提出了调整的多重判定系数，记为$R^2_\alpha$ 。
$R^2_\alpha=1-(1-R^2)(\frac{n-1}{n-k-1})$

$R_a^2$同时考虑了样本量$n$和模型中自变量的个数$k$的影响，这就使得$R_a^2$的值永远小于$R^2$，而且$R_a^2$的值不会因为模型中自变量的个数增多而逐渐接近于$1$。
在多元回归分析中，通常采用调整的多重判定回归系数。

估计标准误差

和一元回归一样，$S_e$是总体回归模型中$\sigma^2=Var(\varepsilon_i)$的无偏估计量。

$S_e=\sqrt{MSE}=\sqrt{\frac{1}{n-k-1}\sum_{i=1}^{n}(y_i-\hat{y_i})^2}$

多元回归对$S_e$的解释和一元类似，这里就不再叙述了。

显著性检验

我们仍需要对线性关系$（F检验）$和回归系数$（t检验）$的检验。

线性关系的检验

为检验 $y$ 与 $k$ 个自变量之间的关系是否显著，进行$F$检验。

其中
$\begin{aligned} &\frac{SSR}{\sigma^2}\sim \chi^2(k)\\ &\frac{SSE}{\sigma^2}\sim \chi^2(n-k-1) \end{aligned}$

构造统计量
$F=\frac{MSR}{MSE} \sim F(k,n-k-1)$

检验步骤如下，
第一步： 提出假设。
    $H_0：\beta_1=\beta_2=\cdots=\beta_k=0，H_1：\beta_1,\beta_2,\cdots,\beta_k不全为0$
第二步：检验统计量。
    $F=\frac{MSR}{MSE}$
第三步：作出决策。
    给定显著性水平$\alpha$，查表，
    若$F>F_\alpha$，拒绝$H_0$，表明$k$个自变量和$y$的线性关系是显著的，
    若$F<F_\alpha$，接受$H_0$，没有证据表明$k$个自变量和$y$的线性关系是显著的。

回归系数的检验

进行完线性关系的检验，就可以对回归系数$\beta_i$进行检验，即每个自变量对$y$的影响进行 $t$

前面我们有
$\hat{\beta} \sim N(\beta,\sigma^2(X^TX)^{-1})$
对于每个$\hat{\beta_i}，i=1,2,\cdots,k$，有
$\hat{\beta_i} \sim N(\beta_i,\sigma^2c_{ii})$
其中$c_{ii}$是$(X^TX)^{-1}$中主对角线第$(i)$个元素。

构造统计量
$t=\frac{\hat{\beta_i}}{S_e/\sqrt{c_{ii}}} \sim t(n-k-1)$

检验步骤如下，
第一步： 提出假设。
    $H_0：\beta_i=0，H_1：\beta_i\ne0$
第二步：检验统计量。
    $t=\frac{\hat{\beta_i}}{S_e/\sqrt{c_{ii}}}$
第三步：作出决策。
    给定显著性水平$\alpha$，并根据自由度$df=n-k-1$查表，
    若$\mid t \mid > t_{\alpha/2}$，拒绝$H_0$，表明$k个x$对$y$的影响是显著的，
    若$\mid t \mid < t_{\alpha/2}$，接受$H_0$，没有证据表明$k个x$对$y$的影响是显著的。

多重共线性

当回归模型中两个或两个以上的变量彼此相关时，则回归模型中存在多重共线性。

多重共线性的判定：

• 模型中中各对自变量之间显著相关
• 当模型的线性关系检验$(F$检验$)$显著时，几乎所有的回归系数$β_i$的$t$检验却不显著。
• 回归系数的正负号与预期的相反。
• 容忍度$(tolerance)$ 与 方差扩大因子$(VIF)$。

容忍度：
某个变量的容忍度等于 $1$ 减去该自变量为因变量而其他 $k−1$ 个自变量为预测变量时所得到的线性回归模型的判定系数。即 $1−R^2_i$。 容忍度越小，多重共线性越严重。通常认为容忍度小于 $0.1$

方差扩大因子：
容忍度的倒数。 因此，$VIF$ 越大，多重共线性越严重，一般认为 $VIF$ 的值大于 $10$

多重共线性问题的处理：

• 将一个或多个相关的自变量从模型中剔除，使保留的自变量尽可能不相关。
• 如果要在模型中保留所有的自变量，那么应该：
  (1)避免根据 $t$ 统计量对单个参数 $\beta$
  (2)对因变量$y$值的推断（预测和估计）限定在自变量样本值的范围内。

预测

当回归模型和系数通过检验后，可由给定的$x_0 = (x_{01},\cdots, x_{0k} )$预测 $y_0$，$y_0$是随机的，显然其预测值为
$\hat{y_0}=\hat{\beta_0}+\hat{\beta_1}x_{01}+\cdots+\hat{\beta_k}x_{0k}$

给定$\alpha$，即可得到$y_0$的预测区间，结果较复杂，

若$n$较大时，且 $x_{0i}$ 接近 $\overline{x}_i$，$y_0$的预测区间可以简化为

$[\hat{y_0}-z_{\alpha/2}S_e,\hat{y_0}+z_{\alpha/2}S_e]$

变量选择与逐步回归

在建立回归模型时，希望用最少的变量来建立模型，有效地引进自变量，使得模型更有操作性。

若新变量引入模型后，之前某个变量的$t检验不显著，则这个新变量将被模型剔除，在这种情况下，回归分析就很难存在多重共线性的影响，因此，要逐步回归搜索自变量，避免造成多重共线性。

逐步回归：逐步回归（Stepwise Regression）是一种常用的消除多重共线性、选取“最优”回归方程的方法。其做法是将逐个引入自变量，引入的条件是该自变量经F检验是显著的，每引入一个自变量后，对已选入的变量进行逐个检验，如果原来引入的变量由于后面变量的引入而变得不再显著，那么就将其剔除。引入一个变量或从回归方程中剔除一个变量，为逐步回归的一步，每一步都要进行F 检验，以确保每次引入新变量之前回归方程中只包含显著的变量。这个过程反复进行，直到既没有不显著的自变量选入回归方程，也没有显著自变量从回归方程中剔除为止。

多元回归分析大概就介绍这么多吧，中心化，标准化，样本相关阵等知识可以看这篇。

实战

下面来实战吧~~~

import pandas as pd
import matplotlib.pyplot as plt
from sklearn.linear_model import LinearRegression
from sklearn import model_selection
import seaborn as sns

# 读取文件
data = pd.read_excel(r'Predict to Profit.xlsx')
examDf = pd.DataFrame(data)
print(examDf.head())

数据文件为Predict to Profit.xlsx，运行结果前五行如下：

该数据有三个自变量影响利润的大小，下面进行回归分析。

# 数据清洗,比如第一列有可能是日期，这样的话我们就只需要从第二列开始的数据，
# 这个情况下，把下面中括号中的0改为1就好，要哪些列取哪些列
new_examDf = examDf.ix[:, 0:]

# 检验数据
print(new_examDf.describe())  # 数据描述，会显示最值，平均数等信息，可以简单判断数据中是否有异常值
print(new_examDf[new_examDf.isnull() == True].count())  # 检验缺失值，若输出为0，说明该列没有缺失值

# 输出相关系数，判断是否值得做线性回归模型
print(new_examDf.corr())  # 0-0.3弱相关；0.3-0.6中相关；0.6-1强相关；

# 通过seaborn添加一条最佳拟合直线和95%的置信带，直观判断相关关系
sns.pairplot(data, aspect=0.8, kind='reg')
plt.show()

运行结果：

训练线性回归模型

# 拆分训练集和测试集
X_train, X_test, Y_train, Y_test = model_selection.train_test_split(new_examDf.ix[:, :3], new_examDf.Profit,
                                                                    train_size=0.8)
print("自变量---源数据:", new_examDf.ix[:, :3].shape, "；  训练集:", X_train.shape, "；  测试集:", X_test.shape)
print("因变量---源数据:", examDf.Profit.shape, "；  训练集:", Y_train.shape, "；  测试集:", Y_test.shape)

# 调用线性规划包
model = LinearRegression()

model.fit(X_train, Y_train)  # 线性回归训练

a = model.intercept_  # 截距
b = model.coef_  # 回归系数
print("拟合参数:截距", a, ",回归系数：", b)

# 显示线性方程，并限制参数的小数位为两位
print("最佳拟合线: Y = ", round(a, 2), "+", round(b[0], 2), "* X1 + ", round(b[1], 2), "* X2 + ", round(b[2], 2), "* X3")

Y_pred = model.predict(X_test)  # 对测试集数据，用predict函数预测

plt.plot(range(len(Y_pred)), Y_pred, 'red', linewidth=2.5, label="predict data")
plt.plot(range(len(Y_test)), Y_test, 'green', label="test data")
plt.legend(loc=3)
plt.show()  # 显示预测值与测试值曲线

运行结果：

对比spss结果，如下图

两种方法，得到两个拟合方程
$y=51130.77+0.77x_1-0.01x_2+0.03x_3$

$y=52382.44+0.78x_1-0.22x_2+0.025x_3$