基于回归方法分析波士顿房价数据间的相关关系

Analysis of the Correlation between Housing Price Data in Boston Based on the Regression Method

摘要: 根据波士顿房价数据集中的变量使用R软件对波士顿房价建立线性回归模型，对回归方程和回归系数进行显著性检验，针对违背基本假设的情况使用Box-Cox变换后再建立模型。为适当精简方程使用Lasso回归，但其建立的模型回归系数很小，原因是此数据中的变量并没有多重共线性，与使用R软件判断结果一致。最后，数据中的响应变量与其相关系数的绝对值大于0.5的自变量建立线性回归方程，并对房价进行预测。由于波士顿房价的分布范围会随着影响因素的变化而发生变化，且中位数具有一定的稳健性，因而我们对房价的中位数建立回归模型，即分位数回归模型。

Abstract: According to the variables in the Boston housing price data set, a linear regression model was es-tablished for the Boston housing price by using R software. The significance test of the regression equation and regression coefficient was carried out. The model was established after the Box-Cox transformation was used for the case that the basic assumptions were violated. Lasso regression was used to simplify the equation appropriately, but the regression coefficient of the model estab-lished by lasso regression was small, because the variables in this data were not multicollinearity, which was consistent with the judgment results of R software. Finally, the response variable in the data and the independent variable whose absolute value of its correlation coefficient is greater than 0.5 establish a linear regression equation and predict the housing price. Because the distribution range of housing price in Boston will change with the change of influencing factors, and the median has certain robustness, we establish a regression model for the median of housing price, namely quantile regression model.

1. 引言

波士顿房价数据集是统计的20世纪70年代中期波士顿郊区房价的中位数，统计了当时教区部分的犯罪率、房产税等共计13个指标，统计出房价，试图能找到指标与房价的关系并进行预测。

2. 材料与方法

2.1. 变量名称与建模目的

2.1.1. 变量名称简介

分析波士顿房价数据集(Boston House Price Dataset)可知影响响应变量MEDV的因素可能有13个，以下为各个属性的介绍，见表1。

Table 1. Introduction of related variables

表1. 相关变量的介绍

本例是属于回归模型的案例，在数据集中包含506组数据。通过对波士顿房地产数据进行初步的观察并分析找出影响房价中位数的因素，希望建立一个能够预测房屋价值的多元线性回归模型。

2.1.2. 多元线性回归模型的一般形式

设随机变量y与一般变量 $x_1,x_2,\cdots ,x_p$ 的线性回归模型为

$y={\beta }_0+{\beta }_1{x}_1+{\beta }_2{x}_2+\cdots +{\beta }_p{x}_p+\epsilon$

式中， ${\beta }_0,{\beta }_1,\cdots ,{\beta }_p$ 是 $p+1$ 个未知参数， ${\beta }_0$ 称为回归常数， ${\beta }_1,\cdots ,{\beta }_p$ 称为回归系数。y称为解释变量(因变量)， $x_1,x_2,\cdots ,x_p$ 是p个可以精确测量并控制的一般变量，称为解释变量(自变量)。

$\epsilon$ 是随机误差，并且假定

$\{\begin{array}{l}E\left(\epsilon \right)=0\hfill \\ \mathrm{var}\left(\epsilon \right)={\sigma }^2\hfill \end{array}$

2.2. 问题解决方法与知识依托

在本例中我们使用R软件解决相应问题。部分代码见参考文献 [1]。

2.2.1. 预处理

首先将数据导入R软件中，为了消除量纲不同和数量级差异带来的影响，就需要对数据进行标准化处理，然后用最小二乘法估计未知参数，求得标准化回归系数。

样本数据的标准化公式为：

$x_{ij}^{*}=\frac{x_{ij}-{\bar{x}}_j}{\sqrt{L_{jj}/n}},i=1,2,\cdots ,n;j=1,2,\cdots ,p$

$y_i^{*}=\frac{y_i-\bar{y}}{\sqrt{L_{yy}/n}},i=1,2,\cdots ,n$

式中

$L_{jj}={\displaystyle \underset{i=1}{\overset{n}{\sum }}{\left(x_{ij}-{\bar{x}}_j\right)}^2}$

是自变量 $x_j\left(j=1,2,\cdots ,p\right)$ 的离差平方和。用最小二乘法求出标准化的样本数据的经验回归方程，记为：

${\hat{y}}^{*}={\hat{\beta }}_1^{*}{x}_1^{*}+{\hat{\beta }}_2^{*}{x}_2^{*}+\cdots +{\hat{\beta }}_p^{*}{x}_p^{*}$

式中， ${\hat{\beta }}_1^{*},{\hat{\beta }}_2^{*},\cdots ,{\hat{\beta }}_p^{*}$ 是y对自变量 $x_1,x_2,\cdots ,x_p$ 的标准化回归系数。

2.2.2. 回归参数的普通最小二乘估计

即寻找参数 ${\beta }_0,{\beta }_1,\cdots ,{\beta }_p$ 的估计值 ${\hat{\beta }}_1,{\hat{\beta }}_2,\cdots ,{\hat{\beta }}_p$，使离差平方和 $Q\left({\beta }_0,{\beta }_1,\cdots ,{\beta }_p\right)={\displaystyle \underset{i=1}{\overset{n}{\sum }}{\left(y_i-{\beta }_0-{\beta }_1{x}_{i1}-\cdots -{\beta }_p{x}_{ip}\right)}^2}$ 达到极小。

当 ${\left(X^{\prime }X\right)}^{-1}$ 存在时，即得回归参数的最小二乘估计为：

$\hat{\beta }={\left(X^{\prime }X\right)}^{-1}{X}^{\prime }y$

2.2.3. 回归方程、回归系数的检验

1) F检验

对多元线性回归方程的显著性检验就是要看自变量 $x_1,x_2,\cdots ,x_p$ 从整体上对随机变量y是否有明显的影响。

原假设 $H_0:{\beta }_1={\beta }_2=\cdots ={\beta }_p=0$

构造F检验统计量如下：

$F=\frac{SSR/p}{SSE/\left(n-p-1\right)}$

当原假设成立时，F服从自由度为 $\left(p,n-p-1\right)$ 的F分布。

当 $F>F_{\alpha }\left(p,n-p-1\right)$ 时，拒绝原假设 $H_0$，否则认为在显著性水平 $\alpha$ 下，y与 $x_1,x_2,\cdots ,x_p$ 有显著的线性关系，即回归方程是显著的。

2) t检验

检验 $x_j$ 是否显著等价于检验

$H_{0j}:{\beta }_j=0,j=1,2,\cdots ,p$

如果接受原假设，则 $x_j$ 不显著；如果拒绝原假设，则 $x_j$ 是显著的。

据此可以构造t统计量

$t_j=\frac{\hat{\beta }}{\sqrt{c_{jj}}\hat{\sigma }}$

式中

$\hat{\sigma }=\sqrt{\frac{1}{n-p-1}{\displaystyle \underset{i=1}{\overset{n}{\sum }}{e}_i^2}}$

2.2.4. 违背基本假设情况的检验

1) 异方差性

违背了回归模型的基本假定，即

$\mathrm{var}\left({\epsilon }_i\right)\ne \mathrm{var}\left({\epsilon }_j\right)$，当 $i\ne j$ 时

诊断方法：绘制残差图 等级相关系数法

解决方法：多元加权最小二乘估计BOX-COX变换

2) 自相关性

违背基本假设，即

$\mathrm{cov}\left({\epsilon }_i,{\epsilon }_j\right)\ne 0$，当 $i\ne j$ 时

诊断方法：图示检验法 自相关系数法DW检验

解决方法：迭代法 差分法BOX-COX变换

2.2.5. 多重共线性

1) 共线性诊断

① 方差扩大因子法

$c_{jj}=\frac{1}{1-R_j^2}$ 作为方差扩大因子的定义，证明见参考文献 [2]，当 $VI{F}_j\ge 10$ 时，说明自变量 $x_j$ 与其余自变量之间有严重的多重共线性。(注意：有些教材认为 $vif>4$ 即存在多重共线性。详见参考文献 [3]。)

② 条件数

记 $X^{\prime }X$ 的最大特征根为 ${\lambda }_m$，称

$k_i=\sqrt{\frac{{\lambda }_m}{{\lambda }_i}},i=0,1,\cdots ,p$

为特征根 ${\lambda }_i$ 的条件数。

通常认为 $0<k<10$ 时，设计矩阵X没有多重共线性； $10\le k<100$ 时，存在较强的多重共线性； $k\ge 100$ 时，存在严重的多重共线性。

2) 解决方法

剔除不重要的解释变量，在此例中我们将看到不显著的回归系数，当回归系数不显著时，剔除变量。

由于样本量足够大，因而增大样本量已经无法解决问题。

岭回归，详见参考文献 [4]。

主成分回归与偏最小二乘估计。

2.2.6. Lasso回归

Lasso回归又称为套索回归，并提供了从零开始到最小二乘拟合的系数和拟合的整个序列。Lasso是一种收缩估计方法，其基本思想是在回归系数的绝对值之和小于一个常数的约束条件下，使残差平方和最小化，从而能够产生某些严格等于0的回归系数，进一步得到可以解释的模型。R语言中有多个包可以实现Lasso回归，这里使用lars包实现。

3. 结果与分析

3.1. 回归方程的建立

3.1.1. 回归方程的初步建立

由于数据为多元的，因而无法用一元回归分析的方法绘制散点图。为探究各个属性与响应变量的关系，我们先对其建立线性回归模型，讨论模型的合理性。回归系数及p值见表2。

Table 2. Coefficients of regression equation and their p values

表2. 回归方程的系数及其p值

统计量 $F=108.1$， $p<\text{2}\text{.2e}-\text{16}$，给定显著性水平 $\alpha =0.05$，则 $p<\alpha$，因而拒绝原假设，认为回归方程是显著的。但是，根据上表可知部分回归系数不显著。此时残差的标准差为0.516， $R^2=0.7406$，调整的 $R^2=0.7338$，拟合效果一般。

3.1.2. 回归方程的进一步分析

由于回归方程中部分回归系数不显著，因而剔除不显著的变量。首先剔除变量中p值最大的，进行回归分析，然后在剩下的变量中剔除最大的，进行分析，依次进行，直至回归方程中所有的回归系数都显著为止。回归方程与逐步回归选择变量结果相同。

建立的回归方程为：

$\begin{array}{c}\hat{M}=-0.101C\hat{R}+0.116\hat{Z}+0.075C\hat{H}-0.219\hat{N}+0.290RM\\ -0.342\hat{D}+0.284RA-0.216\hat{T}-0.223\hat{P}+0.092\hat{B}-0.406\hat{L}\end{array}$

说明：由于变量较多，为适当精简方程，将变量的首字母作为变量，同时回归系数保留三位小数。

此时，统计量 $F=128.5$， $p<\text{2}\text{.2e}-\text{16}$，回归方程显著，调整的 $R^2=0.7348$，并且回归方程的各个回归系数都已显著。

3.2. 违背基本假设情况的检验与解决

回归方程的残差为

$\hat{\epsilon }=Y-\hat{Y}=\left(I-H\right)Y$

其中

$H=X{\left(X^{\prime }X\right)}^{-1}{X}^{\prime }$

称H为帽子矩阵。

在得到回归方程后，计算残差，可以对残差进行正态性检验。

检验结果 $p<\alpha$，因而认为残差不满足正态性假设。另外由残差的QQ图(见图1)也可以看出残差不满足正态性假设。

左上图为残差与拟合图，用来检验线性，若散点集中分布在一条直线附近，则表示线性关系良好；

右上图为QQ图，用来检验正态性，若散点集中分布在Q-Q图中的直线上，则表示残差正态性良好；

左下图为位置尺度图，用来检验同方差性，若点在曲线周围随机分布，则表示同方差性成立；

右下图为残差与杠杆图，可以观测出离群点、高杠杆点和强影响点。独立性是无法从图中分辨出来的。

从图上可以看出该模型的残差并未随机分布，而是呈现异方差的问题。

同样，未标准化的模型也具有异方差性，下面对未标准化的数据进行变换。

对模型进行Box-Cox正态变换，求得 $EstPower=0.1158$，这里取0.116。

记变换后的MEDV为y，回归系数见表3，则回归方程为

$\begin{array}{c}y=5.0748-0.01354\text{CRIM}+1.7157\text{e}-03\text{ZN}+0.1516\text{CHAS}-1.0424\text{NOX}\\ +0.1422\text{RM}-0.0761\text{DIS}+0.0191\text{RAD}-7.896\text{e}-04\text{TAX}\\ -0.0542\text{PTRATIO}+5.917\text{e}-04\text{B}-0.0397\text{LSTAT}\end{array}$

回归方程、回归系数皆通过显著性检验，且 $\hat{\sigma }=0.2688$，比之前模型有所降低， $R^2=0.7884$。比之前有所提高，绘制残差图发现残差也有所改善，因而模型较之前有所改善。变化前后的残差比较见图2。

Figure 1. Residual and fitting diagram, QQ diagram, position scale diagram, residual and lever diagram are drawn

图1. 绘制残差与拟合图、QQ图、位置尺度图、残差与杠杆图

Figure 2. Comparison of two residuals before and after transformation

图2. 变换前后两残差的比较

Table 3. Regression coefficient after Box-Cox transformation

表3. Box-Cox变换后的回归系数

由表3可以看出，部分回归系数较小，且自变量较多。

3.3. 多重共线性的诊断

使用函数vif(myfit)，可以求出各个自变量的方差扩大因子。方差扩大因子均小于10，不存在多重共线性。

使用条件数求得 $k=\text{62}\text{.47931}$， $\sqrt{k}=\text{7}\text{.904386}$，设计矩阵X没有多重共线性，同样由岭迹图也可以看出。

3.4. 降维

虽然以上求得的模型通过了检验，但是自变量数量较多，尝试使用达到降维的目的。

通过Lasso回归得到的 $R^2=0.788$，截距项为4.332224。

Lasso回归后不为零的回归系数见表4。

Table 4. Non-zero Lasso regression coefficient

表4. 不为零的Lasso回归系数

从结果可以看到，ZN项的系数为0，TAX系数的绝对值是剩下的所有项中值最大的，这里也可以看出来，其他项虽然系数都非常小但不为0，这是因为这些项之间的关系是非线性的，无法用线性组合互相表示。

由图3可以看到图中的竖线对应于Lasso中迭代的次数，对应的系数值不为0的自变量即为选入的，竖线的标号与step相对应。

在进行Lasso回归后，自变量的数量只减少一个，且各回归系数取值较小，不好处理。

Figure 3. Shows the order in which the independent variables are selected

图3. 展示自变量被选入的顺序

4. 讨论

4.1. 响应变量与部分自变量的回归模型

1) 通过计算自变量与响应变量的相关系数，可以发现与响应变量有较大相关关系的有RM、PTRATIO、LSTAT三个变量，因此对其建立线性回归模型。

方法与上同，经过Box-Cox变换后的回归方程 $R^2=\text{0}\text{.7167}$， $\hat{\sigma }=0.3975$。

回归方程为：

$yy=4.889430+0.228333\text{RM}-0.072931\text{PTRATIO}-0.061199\text{LSTAT}$

2) 回归系数的解释

RM增加，MEDV也会增加。因为随着房屋数量的增加，相对房屋价格应该会减小。

LSTAT增加，MEDV会减小。因为低收入者多的地方，他们居住的地区房屋价格会低一些。

PTRATIO增加，MEDV会减小。因为师生数量比表明了一个地方教育发展状况，比值越大，说明该地区缺老师，教育状况较差，因此该地区房价也会低。

4.2. 利用回归模型对自有住房的中位数MEDV进行预测

假设王某是一个在波士顿地区的房屋经纪人，使用此模型对客户进行评估他们想要出售的房屋的中位数报价，王某会建议每位客户的房屋销售价格大约为多少？客户的信息见表5，建议见表6。

Table 5. Information collected by three customers

表5. 三个客户收集到的信息

Table 6. Suggests that the mean of the median house price is

表6. 建议房价中位数的均值为

4.3. 模型分析

虽然模型的误差标准差较小，但是模型的拟合优度一般，可能是线性回归模型不合适，也有可能采集的数据不能充分解释响应变量的值，应该尝试建立非线性模型或其他模型提高拟合优度；数据中可以看到异常值，既不能盲目删除，也不应该置之不理，应具体分析，具体情况具体讨论。

4.4. 使用性探讨

1978年采集的数据，在考虑通货膨胀的前提下，由于相关的政策发生了变化，因此在今天不适用；

以上数据所采集的变量，不能够完全描述一个房屋，房屋价值还受房屋外观、新旧程度等因素的影响；

像波士顿这样的大城市，回归模型仅适用于它本身，不能适用于其他乡镇。