赞
踩
摘要:本文基于波士顿房价的公开数据集,寻找影响波士顿房价的因素,统计数据包括城镇人均犯罪率、占地面积超过2.5万平方英尺的住宅用地比例、城镇非零售业务地区的比例以及查尔斯河虚拟变量等十三个因素。本文使用 R 语言,对各个影响因素的相关性进行了筛选分析,并采用逐步回归法得到了最优的多元线性回归模型。在讨论中,对线性回归模型的回归显著性以及拟合优度进行了检验,最后使用最优的多元线性回归模型对波士顿的房价进行预测。
关键字:多元线性回归,逐步回归,回归显著性,拟合优度,回归预测
影响房价的因素多种多样,涉及到包括经济、政治、文化、社会和生态的方方面面。而房价的预测,对于房产中介的定价,房地产企业的房展和人们的正常住房需求有着重要的作用。因此研究影响房价的因素,和房价之间的关系显得十分重要,建立有效的数学模型,根据各种因素对房价进行合理的预测。
本文使用波士顿的房价公开数据集,其包括以1000美元计算的自有住房价格的中位数,以及可能会影响房价的十三个因素,包括:城镇人均犯罪率,占地面积超过2.5万平方英尺的住宅用地比例,城镇非零售业务地区的比例,查尔斯河虚拟变量 (= 1 如果土地在河边;否则是0),一氧化氮浓度,平均每居民房数,在1940年之前建成的所有者占用单位的比例,与五个波士顿就业中心的加权距离,辐射状公路的可达性指数,每10,000美元的全额物业税率,城镇师生比例,城镇黑人的比例,人口中地位较低人群的百分数等因素。多个因素和房价的相关性高低不同,分析各个因素和房价的相关性,寻找主要的影响因素,并建立多元回归模型是一个有效的分析方法。
在现实问题中,变量与变量之间的关系大致可以分类两类:确定性关系和非确定关系。变量间的确定性关系是完全已知的,可用函数来描述;变量之间的非确定性关系称为相关关系,不能用某个确定函数来描述。回归分析就是这样一种应用广泛的研究变量与变量之间的相关关系的数理统计方法。回归分析的基本思想为:虽然变量间的相关关系无法用完全确定的函数关系来描述,但在平均意义下有一定的函数表达式。估计这种函数表达式是回归分析的主要任务,而根据观测数据建立的反映变量间相关关系的统计模型称为回归模型。
回归模型中涉及的相关研究主要分为两类:一类是研究随机变量间的相关关系。另一类是研究随机变量与普通变量间的相关关系。在本文的研究中,除了房价之外,都是普通变量,所以属于后一种类型的回归模型。并且在本文的研究中,涉及到一个响应变量和多个解释变量,属于多元回归模型。并且本文使用线性的统计模型进行分析,因此问题可以进一步归属为多元线性回归模型。
在多元线性回归模型中,问题的数学描述为:设随即变量 y 与 p(p≥2)个普通变量 x 1 , x 2 , … … , x p x_1,x_2,……,x_p x1,x2,……,xp有关,且满足关系式:
{ y = β 0 + β 1 x 1 + β 2 x 2 + … … + β p x p + ε E ( ε ) = 0 , V a r ( ε ) = σ 2 < + o o (1) \left\lbrace \right. \tag{1} {y=β0+β1x1+β2x2+……+βpxp+εE(ε)=0,Var(ε)=σ2<+oo(1)
其中 β 0 , β 1 , β 2 , … … , β p , σ 2 β_0,β_1,β_2,……,β_p,σ^2 β0,β1,β2,……,βp,σ2是与 x 1 , x 2 , … … , x p x_1,x_2,……,x_p x1,x2,……,xp无关的未知参数, ε ε ε是不可观测的随机变量。
在数据观测中,设有n组不全相同的样本观察值 ( x i 1 , x i 2 , … … , x i p ; y i ) ( i = 1 , 2 , … … , n ) (x_{i1},x_{i2},……,x_{ip};y_i)(i=1,2,……,n) (xi1,xi2,……,xip;yi)(i=1,2,……,n),由式(1)得:
{ y = β 0 + β 1 x i 1 + β 2 x i 2 + … … + β p x i p + ε i E ( ε ) = 0 , V a r ( ε ) = σ 2 < + o o i = 1 , 2 , … … , n (2) \left\lbrace \right. i=1,2,……,n \tag{2} {y=β0+β1xi1+β2xi2+……+βpxip+εiE(ε)=0,Var(ε)=σ2<+ooi=1,2,……,n(2)
其中 ε 1 , ε 2 , ε 3 , … … , ε n ε_1,ε_2,ε_3,……,ε_n ε1,ε2,ε3,……,εn相互独立。为了方便,使用矩阵表达式记为:
Y = ( y 1 y 2 ⋮ y n ) , X = ( 1 x 11 ⋯ x 1 p 1 x 21 ⋯ x 2 p ⋮ ⋮ ⋱ ⋮ 1 x n 1 ⋯ x n n ) , β = ( β 1 β 2 ⋮ β n ) , ε = ( ε 1 ε 2 ⋮ ε n ) (3) Y= \left( \right) , X= \left( \right) , β= \left( \right) , ε= \left( \right) \tag{3} Y=⎝⎜⎜⎜⎛y1y2⋮yn⎠⎟⎟⎟⎞,X=⎝⎜⎜⎜⎛11⋮1x11x21⋮xn1⋯⋯⋱⋯x1px2p⋮xnn⎠⎟⎟⎟⎞,β=⎝⎜⎜⎜⎛β1β2⋮βn⎠⎟⎟⎟⎞,ε=⎝⎜⎜⎜⎛ε1ε2⋮εn⎠⎟⎟⎟⎞(3)
则式(2)可表示成矩阵形式:
{ Y = X β + ε E ( ε ) = 0 , V a r ( ε ) = σ 2 I n (4) \left\lbrace \right. \tag{4} {Y=Xβ+εE(ε)=0,Var(ε)=σ2In(4)
其中 I n I_n In为 n 阶单位矩阵,Y 为随机变量的观测向量,β 为未知参数向量,X 为设计矩阵,ε 为 n 维随机误差向量。为了对参数做区间估计与假设检验,需要给出 Y 的确切分布,即给出 ε 的分布,常常假定 ε ~ N ( 0 , σ 2 ) ε ~ N(0, σ^2) ε~N(0,σ2),则式(4)为:
{ Y = X β + ε ε ~ N ( 0 , σ 2 I n ) (5) \left\lbrace \right. \tag{5} {Y=Xβ+εε~N(0,σ2In)(5)
使用最小二乘法估计回归系数β,就是最小化误差平方和:
Q ( β 0 , β 1 , β 2 , … … , β p ) = ∑ i = 1 n ( y i − β 0 − β 1 x i 1 − β 2 x i 2 − … … − β p x i p ) 2 (6) Q(β_0,β_1,β_2,……,β_p)=\sum_{i=1}^n(y_i-β_0-β_1x_{i1}-β_2x_{i2}-……-β_px_{ip})^2 \tag{6} Q(β0,β1,β2,……,βp)=i=1∑n(yi−β0−β1xi1−β2xi2−……−βpxip)2(6)
其误差平方和 Q ( β 0 , β 1 , β 2 , … … , β p ) Q(β_0,β_1,β_2,……,β_p) Q(β0,β1,β2,……,βp)的矩阵表示为:
Q ( β ) = ( Y − X β ) ′ ( Y − X β ) (7) Q(β)=(Y-Xβ)^{'}(Y-Xβ) \tag{7} Q(β)=(Y−Xβ)′(Y−Xβ)(7)
推导得:
∂ Q ( β ) ∂ β = − 2 X ′ Y + 2 X ′ X β (8) \frac{\partial Q(β)}{\partial β} = -2X^{'}Y + 2X^{'}Xβ \tag{8} ∂β∂Q(β)=−2X′Y+2X′Xβ(8)
令 ∂ Q ( β ) ∂ β = 0 \frac{\partial Q(β)}{\partial β}=0 ∂β∂Q(β)=0,可得:
β ^ = ( X ′ X ) − 1 X ′ Y (9) \hat β = (X^{'}X)^{-1}X^{'}Y \tag{9} β^=(X′X)−1X′Y(9)
则可得 p 元经验线性回归方程为:
y ^ = β ^ 0 + β ^ 1 x 1 + β ^ 2 x 2 + … … + β ^ p x p (10) \hat y=\hat β_0+\hat β_1x_1+\hat β_2x_2+……+\hat β_px_p \tag{10} y^=β^0+β^1x1+β^2x2+……+β^pxp(10)
在进行数据分析的时候,需要一个数量指标来衡量回归方程对样本点(x,y)本身的变化趋势描述好坏的程度,即需要一个数量指标来度量拟合好坏的程度。
数据的总波动可用离差平方和来描述:
L y y = ∑ i = 1 n ( y i − y ‾ ) 2 (11) L_{yy}=\sum_{i=1}^n(y_i-\overline y)^2 \tag{11} Lyy=i=1∑n(yi−y)2(11)
其可以分解为两部分:
L y y = Q + U , Q = ∑ i = 1 n ( y i − y ^ i ) 2 , Q = ∑ i = 1 n ( y ^ i − y ‾ ) 2 (12) L_{yy} = Q+U, Q=\sum_{i=1}^n(y_i-\hat y_i)^2, Q=\sum_{i=1}^n(\hat y_i-\overline y)^2 \tag{12} Lyy=Q+U,Q=i=1∑n(yi−y^i)2,Q=i=1∑n(y^i−y)2(12)
平方和的分解公式,式(12)表明, y 1 , y 2 , … … , y n y_1,y_2,……,y_n y1,y2,……,yn的分散程度,可以分解为两部分,一部分是由于x对y的线性关系引起的 y 的分散性,即回归平方和 U;另一部分是剩余部分引起的分散性,即残差平方和 Q,由随机误差引起的分散性包含在 Q 中。
在总和 L y y L_{yy} Lyy中,U 所占的比重越大,说明随机误差所占的比重越小,回归的效果就越显著,因此可以使用如下公式来度量线性回归效果的显著性:
R 2 = U L y y = 1 − Q L y y (13) R^2=\frac{U}{L_{yy}}=1-\frac{Q}{L_{yy}} \tag{13} R2=LyyU=1−LyyQ(13)
其中 R 2 R^2 R2描述了回归平方和 U 占离差平方和 L y y L_{yy} Lyy的比重,称 R 2 R^2 R2为决定系数,也称为拟合优度,表示用直线来拟合数据的好坏程度。而其开根号后的 r 称为相关系数,可以描述两个变量 y 与 x 之间的线性相关的密切程度。
由于 L y y ≥ U ≥ 0 L_{yy}≥U≥0 Lyy≥U≥0,故 1 ≥ R 2 ≥ 0 1≥R^2≥0 1≥R2≥0,即 |r| 越接近于 1,y 与 x 之间的线性相关程度就越大,即线性趋势越明显。特别地,当 |r| = 1 时,所有的点都在回归直线上;当 |r| = 0 时,回归直线平行于 y 轴,说明 y 与 x 线性无关。
综上所述,|r| 越小,则线性回归的效果越差,y 与 x 之间的线性相关性就不显著;而 |r| 越大,线性回归的效果越好,y 与 x 之间的线性相关性就越显著。因此,在本文中,使用 |r| 对线性相关关系是否显著进行检验。
前面的线性回归效果的显著性检验可以转化为如下的假设检验问题:
H 0 : b = 0 , H 1 : b ≠ 0 (14) H_0: b = 0, H_1: b ≠ 0 \tag{14} H0:b=0,H1:b=0(14)
这其实也是检验响应变量 y 是否显著依赖于解释变量 x 的假设检验问题。但对于多元线性回归,回归效果的显著性检验问题与回归系数的显著性检验问题是不同的。
回归效果的显著性检验主要是检验响应变量对所有解释变量的线性依赖性是否显著,而回归系数的显著性检验主要是检验响应变量对某个解释变量或某几个解释变量的线性依赖性是否显著。
在本文中,对多元线性回归的回归效果的显著性检验使用 F 检验法;对多元线性回归的回归系数的显著性检验使用 t 检验法。
线性回归方程的最优性包括在两个方面:
(1)该模型中包含所有对因变量 y 有显著影响的自变量。
(2)该模型中所包含的自变量个数尽可能少,不含有无意义的变量,而且还应该是这类方程中 R 2 R^2 R2达到最大者。
选择最优回归方程的方法有四种:全部比较法、向后法(只出不进法)、向前法(只进不出法)、逐步回归法(有进有出法)。在本文中,使用逐步回归法进行最优方程的选择。
逐步回归法的基本思想为:将变量一个个引入,引入的条件是该变量的偏 F 检验是显著的。同时,每引入一个新变量后又要对老变量逐个检验,将变得不显著的变量从回归模型中剔除。即逐步回归法实际上是对向前法和向后法的一种结合。重复上述步骤,直到所有模型外的变量都不能引入,模型内的变量都不能被剔除,即表示是找到了最优的回归方程。
在波士顿房价的公开数据集中,共包含506条观测数据。在本文中,选取前500条观测数据建立多元线性回归模型,并使用后6条数据进行多元线性回归模型的预测和分析。
在波士顿房价的公开数据集,各个因素的含义如下表所示:
序号 | 因素 | 含义 |
---|---|---|
1 | CRIM | 镇人均犯罪率 |
2 | ZN | 占地面积超过2.5万平方英尺的住宅用地比例 |
3 | INDUS | 城镇非零售业务地区的比例 |
4 | CHAS | 查尔斯河虚拟变量 (= 1 如果土地在河边;否则是0) |
5 | NOX | 一氧化氮浓度 |
6 | RM | 平均每居民房数 |
7 | AGE | 在1940年之前建成的所有者占用单位的比例 |
8 | DIS | 与五个波士顿就业中心的加权距离 |
9 | RAD | 辐射状公路的可达性指数 |
10 | TAX | 每10,000美元的全额物业税率 |
11 | PTRATIO | 城镇师生比例 |
12 | B | 城镇黑人的比例 |
13 | LSTAT | 人口中地位较低人群的百分数 |
14 | MEDV | 自住房屋的房价中位数,单位1000美元 |
首先读取 csv 格式的波士顿房价的公开数据集,其代码为:
> Boston_Housing_Data=read.csv("/Users/limengfan/Course_2020_2/数理统计/大作业:回归/Boston_Housing.csv")
>
上述代码将 csv 格式的数据集读取到了 R 语言中,存储格式为 R 语言的 data.frame 数据结构,显示前 3 行数据如下:
> Boston_Housing_Data[1:3,]
CRIM ZN INDUS CHAS NOX RM AGE DIS RAD TAX PTRATIO B LSTAT
1 0.00632 18 2.31 0 0.538 6.575 65.2 4.0900 1 296 15.3 396.90 4.98
2 0.02731 0 7.07 0 0.469 6.421 78.9 4.9671 2 242 17.8 396.90 9.14
3 0.02729 0 7.07 0 0.469 7.185 61.1 4.9671 2 242 17.8 392.83 4.03
MEDV
1 24.0
2 21.6
3 34.7
>
先针对各个影响因素,计算影响因素和房价中位数的相关系数,代码如下:
> cor(Boston_Housing_Data$CRIM, Boston_Housing_Data$MEDV) [1] -0.3883046 > cor(Boston_Housing_Data$ZN, Boston_Housing_Data$MEDV) [1] 0.3604453 > cor(Boston_Housing_Data$INDUS, Boston_Housing_Data$MEDV) [1] -0.4837252 > cor(Boston_Housing_Data$CHAS, Boston_Housing_Data$MEDV) [1] 0.1752602 > cor(Boston_Housing_Data$NOX, Boston_Housing_Data$MEDV) [1] -0.4273208 > cor(Boston_Housing_Data$RM, Boston_Housing_Data$MEDV) [1] 0.6953599 > cor(Boston_Housing_Data$AGE, Boston_Housing_Data$MEDV) [1] -0.3769546 > cor(Boston_Housing_Data$DIS, Boston_Housing_Data$MEDV) [1] 0.2499287 > cor(Boston_Housing_Data$RAD, Boston_Housing_Data$MEDV) [1] -0.3816262 > cor(Boston_Housing_Data$TAX, Boston_Housing_Data$MEDV) [1] -0.4685359 > cor(Boston_Housing_Data$PTRATIO, Boston_Housing_Data$MEDV) [1] -0.5077867 > cor(Boston_Housing_Data$B, Boston_Housing_Data$MEDV) [1] 0.3334608 > cor(Boston_Housing_Data$LSTAT, Boston_Housing_Data$MEDV) [1] -0.7376627 >
对各个单独的影响因素和房价,作散点图,代码如下:
> par(mfrow=c(4,4)) > plot( + x=Boston_Housing_Data$CRIM, xlab="CRIM", + y=Boston_Housing_Data$MEDV, ylab="MEDV", + pch='.' + ) > plot( + x=Boston_Housing_Data$ZN, xlab="ZN", + y=Boston_Housing_Data$MEDV, ylab="MEDV", + pch='.' + ) > plot( + x=Boston_Housing_Data$INDUS, xlab="INDUS", + y=Boston_Housing_Data$MEDV, ylab="MEDV", + pch='.' + ) > plot( + x=Boston_Housing_Data$CHAS, xlab="CHAS", + y=Boston_Housing_Data$MEDV, ylab="MEDV", + pch='.' + ) > plot( + x=Boston_Housing_Data$NOX, xlab="NOX", + y=Boston_Housing_Data$MEDV, ylab="MEDV", + pch='.' + ) > plot( + x=Boston_Housing_Data$RM, xlab="RM", + y=Boston_Housing_Data$MEDV, ylab="MEDV", + pch='.' + ) > plot( + x=Boston_Housing_Data$AGE, xlab="AGE", + y=Boston_Housing_Data$MEDV, ylab="MEDV", + pch='.' + ) > plot( + x=Boston_Housing_Data$DIS, xlab="DIS", + y=Boston_Housing_Data$MEDV, ylab="MEDV", + pch='.' + ) > plot( + x=Boston_Housing_Data$RAD, xlab="RAD", + y=Boston_Housing_Data$MEDV, ylab="MEDV", + pch='.' + ) > plot( + x=Boston_Housing_Data$TAX, xlab="TAX", + y=Boston_Housing_Data$MEDV, ylab="MEDV", + pch='.' + ) > plot( + x=Boston_Housing_Data$PTRATIO, xlab="PTRATIO", + y=Boston_Housing_Data$MEDV, ylab="MEDV", + pch='.' + ) > plot( + x=Boston_Housing_Data$B, xlab="B", + y=Boston_Housing_Data$MEDV, ylab="MEDV", + pch='.' + ) > plot( + x=Boston_Housing_Data$LSTAT, xlab="LSTAT", + y=Boston_Housing_Data$MEDV, ylab="MEDV", + pch='.' + ) >
散点图如下图所示:
可以看到,针对相关系数的符号,部分影响因素的与房价的相关关系表现为正相关,部分表现为负相关。针对相关系数的绝对值大小,部分影响因素对房价的影响显著,部分影响不显著。
针对所有的影响因素,建立多元线性回归方程:
M E D V = I n t e r c e p t + k 1 ⋅ C R I M + k 2 ⋅ Z N + k 3 ⋅ I N D U S + k 4 ⋅ C H A S + k 5 ⋅ N O X + k 6 ⋅ R M + k 7 ⋅ A G E + k 8 ⋅ D I S + k 9 ⋅ R A D + k 1 0 ⋅ T A X + k 1 1 ⋅ P T R A T I O + k 1 2 ⋅ B + k 1 3 ⋅ L S T A T (15) MEDV = Intercept+k_1·CRIM+k_2·ZN+k_3·INDUS+k_4·CHAS+k_5·NOX+k_6·RM \\ +k_7·AGE+k_8·DIS+k_9·RAD+k_10·TAX+k_11·PTRATIO+k_12·B+k_13·LSTAT \tag{15} MEDV=Intercept+k1⋅CRIM+k2⋅ZN+k3⋅INDUS+k4⋅CHAS+k5⋅NOX+k6⋅RM+k7⋅AGE+k8⋅DIS+k9⋅RAD+k10⋅TAX+k11⋅PTRATIO+k12⋅B+k13⋅LSTAT(15)
下面进行线性回归方程模型的建立,首先取前500条数据建立多元线性回归模型,并输出所求得的模型的信息:
> tlm <- lm(MEDV~CRIM+ZN+INDUS+CHAS+NOX+RM+AGE+DIS+RAD+TAX+PTRATIO+B+LSTAT, data=Boston_Housing_Data[1:500,]) > tlm Call: lm(formula = MEDV ~ CRIM + ZN + INDUS + CHAS + NOX + RM + AGE + DIS + RAD + TAX + PTRATIO + B + LSTAT, data = Boston_Housing_Data[1:500, ]) Coefficients: (Intercept) CRIM ZN INDUS CHAS NOX 35.406834 -0.105856 0.048752 0.018293 2.611199 -17.006242 RM AGE DIS RAD TAX PTRATIO 3.797362 0.002020 -1.514495 0.290883 -0.012565 -0.885637 B LSTAT 0.009296 -0.538782 >
可得,使用所有影响因素建立的初始多元线性回归方程为:
M E D V = 35.406834 + ( − 0.105856 ⋅ C R I M ) + ( 0.048752 ⋅ Z N ) + ( 0.018293 ⋅ I N D U S ) + ( 2.611199 ⋅ C H A S ) + ( − 17.006242 ⋅ N O X ) + ( 3.797362 ⋅ R M ) + ( 0.002020 ⋅ A G E ) + ( − 1.514495 ⋅ D I S ) + ( 0.290883 ⋅ R A D ) + ( − 0.012565 ⋅ T A X ) + ( − 0.885637 ⋅ P T R A T I O ) + ( 0.009296 ⋅ B ) + ( − 0.538782 ⋅ L S T A T ) (16) MEDV = 35.406834+(-0.105856·CRIM)+(0.048752·ZN)+(0.018293·INDUS) \\ +(2.611199·CHAS)+(-17.006242·NOX)+(3.797362·RM)+(0.002020·AGE) \\ +(-1.514495·DIS)+(0.290883·RAD)+(-0.012565·TAX)+(-0.885637·PTRATIO) \\ +(0.009296·B)+(-0.538782·LSTAT) \tag{16} MEDV=35.406834+(−0.105856⋅CRIM)+(0.048752⋅ZN)+(0.018293⋅INDUS)+(2.611199⋅CHAS)+(−17.006242⋅NOX)+(3.797362⋅RM)+(0.002020⋅AGE)+(−1.514495⋅DIS)+(0.290883⋅RAD)+(−0.012565⋅TAX)+(−0.885637⋅PTRATIO)+(0.009296⋅B)+(−0.538782⋅LSTAT)(16)
由上文的相关性分析可得,不同的影响因素,其对房价的影响显著性不同,我们应该选择显著性较高的影响因素。在自变量个数尽可能少的同时,使得 R 2 R^2 R2达到最大。
由此使用逐步回归法求得最优的多元线性回归模型,并输出所求得的模型的信息:
> tstep <- step(tlm) Start: AIC=1569.38 MEDV ~ CRIM + ZN + INDUS + CHAS + NOX + RM + AGE + DIS + RAD + TAX + PTRATIO + B + LSTAT Df Sum of Sq RSS AIC - AGE 1 0.52 10910 1567.4 - INDUS 1 1.99 10911 1567.5 <none> 10909 1569.4 - CHAS 1 206.51 11116 1576.8 - CRIM 1 233.27 11142 1578.0 - TAX 1 250.67 11160 1578.7 - B 1 269.66 11179 1579.6 - ZN 1 282.35 11192 1580.2 - RAD 1 428.01 11337 1586.6 - NOX 1 442.50 11352 1587.3 - PTRATIO 1 976.57 11886 1610.2 - DIS 1 1287.34 12196 1623.2 - RM 1 1851.73 12761 1645.8 - LSTAT 1 2506.89 13416 1670.8 Step: AIC=1567.4 MEDV ~ CRIM + ZN + INDUS + CHAS + NOX + RM + DIS + RAD + TAX + PTRATIO + B + LSTAT Df Sum of Sq RSS AIC - INDUS 1 1.99 10912 1565.5 <none> 10910 1567.4 - CHAS 1 208.16 11118 1574.8 - CRIM 1 233.22 11143 1576.0 - TAX 1 250.19 11160 1576.7 - B 1 272.31 11182 1577.7 - ZN 1 283.63 11193 1578.2 - RAD 1 428.09 11338 1584.7 - NOX 1 466.45 11376 1586.3 - PTRATIO 1 978.13 11888 1608.3 - DIS 1 1420.14 12330 1626.6 - RM 1 1947.00 12857 1647.5 - LSTAT 1 2813.68 13723 1680.1 Step: AIC=1565.49 MEDV ~ CRIM + ZN + CHAS + NOX + RM + DIS + RAD + TAX + PTRATIO + B + LSTAT Df Sum of Sq RSS AIC <none> 10912 1565.5 - CHAS 1 214.62 11126 1573.2 - CRIM 1 235.20 11147 1574.2 - B 1 271.08 11183 1575.8 - ZN 1 281.80 11194 1576.2 - TAX 1 286.11 11198 1576.4 - RAD 1 446.76 11358 1583.6 - NOX 1 486.15 11398 1585.3 - PTRATIO 1 982.85 11895 1606.6 - DIS 1 1518.18 12430 1628.6 - RM 1 1952.66 12864 1645.8 - LSTAT 1 2820.20 13732 1678.4 > > tstep Call: lm(formula = MEDV ~ CRIM + ZN + CHAS + NOX + RM + DIS + RAD + TAX + PTRATIO + B + LSTAT, data = Boston_Housing_Data[1:500, ]) Coefficients: (Intercept) CRIM ZN CHAS NOX RM 35.259036 -0.106217 0.048060 2.644310 -16.553028 3.799019 DIS RAD TAX PTRATIO B LSTAT -1.535918 0.284626 -0.012056 -0.879144 0.009297 -0.534946 >
可以看到,在逐步回归法中,进行了多次的迭代,最终确定的最优多元线性回归方程为:
M E D V = I n t e r c e p t + k 1 ⋅ C R I M + k 2 ⋅ Z N + k 4 ⋅ C H A S + k 5 ⋅ N O X + k 6 ⋅ R M + k 8 ⋅ D I S + k 9 ⋅ R A D + k 1 0 ⋅ T A X + k 1 1 ⋅ P T R A T I O + k 1 2 ⋅ B + k 1 3 ⋅ L S T A T (17) MEDV = Intercept+k_1·CRIM+k_2·ZN+k_4·CHAS+k_5·NOX+k_6·RM \\ +k_8·DIS+k_9·RAD+k_10·TAX+k_11·PTRATIO+k_12·B+k_13·LSTAT \tag{17} MEDV=Intercept+k1⋅CRIM+k2⋅ZN+k4⋅CHAS+k5⋅NOX+k6⋅RM+k8⋅DIS+k9⋅RAD+k10⋅TAX+k11⋅PTRATIO+k12⋅B+k13⋅LSTAT(17)
其选用了显著性水平较高的11个影响因素,同时保持了 R 2 R^2 R2达到最大:
M E D V = 35.259036 + ( − 0.106217 ⋅ C R I M ) + ( 0.048060 ⋅ Z N ) + ( 2.644310 ⋅ C H A S ) + ( − 16.553028 ⋅ N O X ) + ( 3.799019 ⋅ R M ) + ( − 1.535918 ⋅ D I S ) + ( 0.284626 ⋅ R A D ) + ( − 0.012056 ⋅ T A X ) + ( − 0.879144 ⋅ P T R A T I O ) + ( 0.009297 ⋅ B ) + ( − 0.534946 ⋅ L S T A T ) (18) MEDV = 35.259036+(-0.106217·CRIM)+(0.048060·ZN)+(2.644310·CHAS) \\ +(-16.553028·NOX)+(3.799019·RM)+(-1.535918·DIS)+(0.284626·RAD) \\ +(-0.012056·TAX)+(-0.879144·PTRATIO)+(0.009297·B)+(-0.534946·LSTAT) \tag{18} MEDV=35.259036+(−0.106217⋅CRIM)+(0.048060⋅ZN)+(2.644310⋅CHAS)+(−16.553028⋅NOX)+(3.799019⋅RM)+(−1.535918⋅DIS)+(0.284626⋅RAD)+(−0.012056⋅TAX)+(−0.879144⋅PTRATIO)+(0.009297⋅B)+(−0.534946⋅LSTAT)(18)
为了测试多元线性回归模型的正确性,使用波士顿房价的公开数据集中未用来进行模型建立的后六条观测数据进行模型预测分析。
数据 | 第一组 | 第二组 | 第三组 | 第四组 | 第五组 | 第六组 |
---|---|---|---|---|---|---|
CRIM | 0.22438 | 0.06263 | 0.04527 | 0.06076 | 0.10959 | 0.04741 |
ZN | 0 | 0 | 0 | 0 | 0 | 0 |
INDUS | 9.69 | 11.93 | 11.93 | 11.93 | 11.93 | 11.93 |
CHAS | 0 | 0 | 0 | 0 | 0 | 0 |
NOX | 0.585 | 0.573 | 0.573 | 0.573 | 0.573 | 0.573 |
RM | 6.027 | 6.593 | 6.12 | 6.976 | 6.794 | 6.03 |
AGE | 79.7 | 69.1 | 76.7 | 91 | 89.3 | 80.8 |
DIS | 2.4982 | 2.4786 | 2.2875 | 2.1675 | 2.3889 | 2.505 |
RAD | 6 | 1 | 1 | 1 | 1 | 1 |
TAX | 391 | 273 | 273 | 273 | 273 | 273 |
PTRATIO | 19.2 | 21 | 21 | 21 | 21 | 21 |
B | 396.9 | 391.99 | 396.9 | 396.9 | 393.45 | 396.9 |
LSTAT | 14.33 | 9.67 | 9.08 | 5.64 | 6.48 | 7.88 |
MEDV(真值) | 16.8 | 22.4 | 20.6 | 23.9 | 22 | 11.9 |
最初模型:式(16) | 20.67286194 | 23.99180352 | 22.86578569 | 28.17872333 | 26.65904323 | 22.84921432 |
最优模型:式(18) | 20.74983713 | 24.01022171 | 22.86990999 | 28.14474935 | 26.62665978 | 22.83564401 |
根据模型预测的值可得,最初模型和最优模型的差距并不大,且最优模型的参数更少,计算更简单,而且其值也更接近真值。
但是两个模型和真值都有一定的误差,说明使用此多元线性回归模型完整的建模此问题还有一定的差距。
在原始的 Excel 表中使用公式计算最初模型和最优模型的预测值,其公式如下:
最初模型:式(16)
=35.406834+(-0.105856*A2)+(0.048752*B2)+(0.018293*C2)+(2.611199*D2)
+(-17.006242*E2)+(3.797362*F2)+(0.00202*G2)+(-1.514495*H2)+(0.290883*I2)
+(-0.012565*J2)+(-0.885637*K2)+(0.009296*L2)+(-0.538782*M2)
最优模型:式(18)
= 35.259036+(-0.106217*A2)+(0.04806*B2)+(2.64431*D2)+(-16.553028*E2)
+(3.799019*F2)+(-1.535918*H2)+(0.284626*I2)+(-0.012056*J2)+(-0.879144*K2)
+(0.009297*L2)+(-0.534946*M2)
针对真值,最初模型的预测值和最优模型的预测值绘制折线图如下:
从上图可以得出,两个模型都在一定程度上起到了很好的预测效果,其中最优模型的方差较小,和真值更为接近。
针对逐步回归法求得的最优模型,对其回归显著性以及拟合优度进行检验。使用summary()
函数获取逐步回归的具体信息,代码及输入如下所示:
> summary(tstep) Call: lm(formula = MEDV ~ CRIM + ZN + CHAS + NOX + RM + DIS + RAD + TAX + PTRATIO + B + LSTAT, data = Boston_Housing_Data[1:500, ]) Residuals: Min 1Q Median 3Q Max -15.6886 -2.6695 -0.5701 1.8036 26.0814 Coefficients: Estimate Std. Error t value Pr(>|t|) (Intercept) 35.259036 5.097433 6.917 1.45e-11 *** CRIM -0.106217 0.032750 -3.243 0.001263 ** ZN 0.048060 0.013538 3.550 0.000422 *** CHAS 2.644310 0.853525 3.098 0.002060 ** NOX -16.553028 3.549999 -4.663 4.03e-06 *** RM 3.799019 0.406532 9.345 < 2e-16 *** DIS -1.535918 0.186399 -8.240 1.60e-15 *** RAD 0.284626 0.063676 4.470 9.74e-06 *** TAX -0.012056 0.003370 -3.577 0.000382 *** PTRATIO -0.879144 0.132603 -6.630 8.92e-11 *** B 0.009297 0.002670 3.482 0.000543 *** LSTAT -0.534946 0.047633 -11.231 < 2e-16 *** --- Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 Residual standard error: 4.729 on 488 degrees of freedom Multiple R-squared: 0.7436, Adjusted R-squared: 0.7379 F-statistic: 128.7 on 11 and 488 DF, p-value: < 2.2e-16 >
由上面的数据分析得到,逐步回归法求得的最优模型,其最终的公式为:MEDV ~ CRIM + ZN + CHAS + NOX + RM + DIS + RAD + TAX + PTRATIO + B + LSTAT
,即只使用了11个影响因素,并且各个影响因素的显著性都很高,通过对每个影响因素进行 t-检验获取的 t 值可以得到,公式中使用的各个影响因素的显著性水平都在0.005以上。
并且分析可得,残差的中位数为-0.5701
,即最优线性回归模型的残差在可接受的范围之内。其中相关系数为
R
2
=
0.7436
R^2=0.7436
R2=0.7436,表示 y 与 x 的线性相关性较强,且对回归的整体效果进行 F-检验,可得 p 值为 <2.2e-16,表示最优线性回归模型的显著性水平很高,最终回归的整体效果较为理想。
[1] 波士顿房价数据集官网. https://www.cs.toronto.edu/~delve/data/boston/bostonDetail.html
[2] 波士顿房价数据集. http://lib.stat.cmu.edu/datasets/boston
[3] 孙海燕,周梦,李卫国,冯伟. 数理统计[M]. 北京:北京航空航天大学出版社, 2016.
[4] 方小敏,齐德胜,张文霖,冯伟. 谁说菜鸟不会数据分析(R 语言篇)[M]. 电子工 业出版社, 2020.
Copyright © 2003-2013 www.wpsshop.cn 版权所有,并保留所有权利。