基于Python的多元线性回归及其应用 多元线性回归介绍、建模与Python实现，使用sklearn库
基于Python的多元线性回归及其应用
1. 多元线性回归模型
在之前文章基于Python的一元线性回归中，我们对回归分析及一元线性回归模型做了介绍，但实际问题中，影响因素往往不止一个，当因变量
Y
Y
Y与多个自变量
X
1
,
X
2
,
.
.
.
,
X
k
X_1,X_2,...,X_k
X1,X2,...,Xk有线性关系时，就需要用到多元线性回归模型描述：
Y
=
β
0
+
β
1
X
1
+
β
2
X
2
+
.
.
.
+
β
k
X
k
+
ϵ
ϵ
i
∼
N
(
0
,
σ
2
)
Y = \beta_0 + \beta_1 X_1 + \beta_2 X_2 + ... + \beta_k X_k +\epsilon \epsilon_i \sim N(0,\sigma^2)
Y=β0+β1X1+β2X2+...+βkXk+ϵϵi∼N(0,σ2)
上式中，
β
j
(
j
=
0
,
1
,
.
.
.
,
k
)
\beta_j(j=0,1,...,k)
βj(j=0,1,...,k)称为回归系数或模型待估参数。
β
j
\beta_j
βj与
σ
2
\sigma^2
σ2均未知，同一元回归分析一样，需要利用样本数据进行估计。
上式称为总体回归函数的随机表达式，它的非随机表达式为:
E
(
Y
∣
X
1
,
X
2
,
.
.
.
,
X
k
)
=
β
0
+
β
1
X
1
+
β
2
X
2
+
.
.
.
+
β
k
X
k
ϵ
i
∼
N
(
0
,
σ
2
)
E(Y|X_1,X_2,...,X_k) = \beta_0 + \beta_1 X_1 + \beta_2 X_2 + ... + \beta_k X_k \epsilon_i \sim N(0,\sigma^2)
E(Y∣X1,X2,...,Xk)=β0+β1X1+β2X2+...+βkXkϵi∼N(0,σ2)
设
(
X
i
1
,
X
i
2
,
.
.
.
,
X
i
k
,
Y
i
)
:
(
i
=
1
,
2
,
.
.
.
,
n
)
{(X_{i1},X_{i2},...,X_{ik},Y_{i}):(i=1,2,...,n)}
(Xi1,Xi2,...,Xik,Yi):(i=1,2,...,n)为一组观测值，则总体回归模型为：
{
Y
=
β
0
+
β
1
X
i
1
+
β
2
X
i
2
+
.
.
.
+
β
k
X
i
k
+
ϵ
i
,
i
=
1
,
2
,
.
.
.
,
n
ϵ
i
∼
N
(
0
,
σ
2
)
C
o
v
(
ϵ
i
,
ϵ
j
)
=
0
,
i
≠
j
,
i
,
j
=
1
,
2
,
.
.
.
,
n
\left\{ \begin{aligned} & Y = \beta_0 + \beta_1 X_{i1} + \beta_2 X_{i2} + ... + \beta_k X_{ik} +\epsilon_i,i=1,2,...,n\\ & \epsilon_i \sim N(0,\sigma^2) \\ & Cov(\epsilon_i,\epsilon_j)=0,i\not=j,i,j=1,2,...,n \end{aligned} \right.
⎩⎨⎧Y=β0+β1Xi1+β2Xi2+...+βkXik+ϵi,i=1,2,...,nϵi∼N(0,σ2)Cov(ϵi,ϵj)=0,i=j,i,j=1,2,...,n
多元线性回归与一元回归模型类似，为保证普通最小二乘估计量的优良性质，多元回归模型也需满足经典假设条件（约束条件），即随机误差项服从正态分布，且具有0均值、同方差、序列不相关性。此外，由于有多个自变量，自变量之间也不能存在完全多重共线性。
2. 基于Python的多元线性回归
随机创建 X 1 , X 2 , X 3 , Y X1,X2,X3,Y X1,X2,X3,Y四个数组，使 Y = 4 ∗ X 1 − 3 ∗ X 2 + 2 X 3 + 5 Y=4*X1-3*X2+2X3+5 Y=4∗X1−3∗X2+2X3+5。然后加入一些干扰噪声，再尝试做线性回归。
from sklearn.linear_model import LinearRegression
import numpy as np
import random
if __name__ == '__main__':
# 随机创建X1,X2,X3,Y。使Y=4*X1-3*X2+2X3+5
X1 = [random.randint(0,100) for i in range(0, 50)]
X2 = [random.randint(0,60) for i in range(0, 50)]
X3 = [random.randint(0, 35) for i in range(0, 50)]
Y =[4*x1-3*x2+2*x3+5 for x1,x2,x3 in zip(X1,X2,X3)]
# 组合X1，X2成n行2列数据
X_train = np.array(X1+X2+X3).reshape((len(X1), 3), order="F")
Y_train = np.array(Y).reshape((len(Y), 1))
# 加入噪声干扰
noise = np.random.randn(50, 1)
noise = noise - np.mean(noise)
Y_train = Y_train+noise
#新建一个线性回归模型，并把数据放进去对模型进行训练
lineModel = LinearRegression()
lineModel.fit(X_train, Y_train)
#用训练后的模型，进行预测
Y_predict = lineModel.predict(X_train)
f=""
#coef_是系数，intercept_是截距
a_arr = lineModel.coef_[0]
b = lineModel.intercept_[0]
for i in range(0,len(a_arr)):
ai=a_arr[i]
if ai>=0:
ai = "+%.4f" %(ai)
else:
ai = "%.4f" % (ai)
f = f+"%s*x%s"%(ai, str(i+1))
f="y=%s+%.4f" % (f[1:],b)
print("拟合方程",f)
#对回归模型进行评分，这里简单使用训练集进行评分，实际很多时候用其他的测试集进行评分
print("得分", lineModel.score(X_train, Y_train))拟合方程 y=4.0056*x1-3.0077*x2+1.9824*x3+5.2960
得分 0.99994863924353833. 多元线性回归的应用
3.1 点预测
- 由多元线性回归模型
Y
^
=
X
β
^
\hat{Y}=X\hat{\beta}
Y^=Xβ^，对于给定的样本以外的自变量的具体值
X
0
=
(
1
,
X
01
,
X
02
,
.
.
.
,
X
0
k
)
X_0=(1,X_{01},X_{02},...,X_{0k})
X0=(1,X01,X02,...,X0k)，可以得到解释变量的预测值的估计值：
Y 0 ^ = X 0 β ^ \hat{Y_0}=X_0\hat{\beta} Y0^=X0β^ - 将回归值 Y 0 ^ \hat{Y_0} Y0^作为 Y 0 Y_0 Y0的估计值，称 Y 0 ^ \hat{Y_0} Y0^为 Y 0 Y_0 Y0的点预测值。
3.2 总体均值的区间预测
- 因为回归值 Y 0 ^ \hat{Y_0} Y0^服从正态分布,通过构造 t t t统计量，并根据给定得置信水平可进行总体均值的区间预测。
3.3 总体个别值的区间预测
- 因为 Y 0 ^ \hat{Y_0} Y0^和 Y 0 Y_0 Y0均服从正态分布，所以预测误差也服从正态分布。
3.4 多元线性回归的应用案例
- 为了全面反映中国“人口自然增长率”的全貌，选择人口增长率作为被解释变量，来反映中国人口的增长；选择“国民收入（亿元）” 与“人均GDP”作为经济增长的指标；并以“居民消费价格指数增长率”作为居民消费水平的指标。试用回归模型分析它们的关系，若2017年的居民消费增长率为2.00%，国民总收入为760000亿元，并且人均GDP为55000元，那么预测人口自然增长率是多少？从《中国统计局》收集到1989年至2016年的相关数据，如下表所示。
设定的线性回归模型为：
Y
t
=
β
0
+
β
1
X
1
t
+
β
2
X
2
t
+
.
.
.
+
β
3
X
3
t
+
ϵ
t
Y_t = \beta_0 + \beta_1 X_{1t} + \beta_2 X_{2t} + ... + \beta_3 X_{3t} +\epsilon _t
Yt=β0+β1X1t+β2X2t+...+β3X3t+ϵt
上式中，人口自然增长率(‰)为被解释变量
Y
t
Y_t
Yt，居民消费价格指数增长率
X
1
t
X_{1t}
X1t、国民收入（亿元）
X
2
t
X_{2t}
X2t与人均GDP
X
3
t
X_{3t}
X3t为三个解释变量，
ϵ
t
\epsilon_{t}
ϵt为误差项。
from sklearn.linear_model import LinearRegression
import numpy as np
#人口自然增长率
y=[15.04,14.39,12.98,11.6,11.45,11.21,10.55,
10.42,10.06,9.14,8.18,7.58,6.95,6.45,6.01,
5.87,5.89,5.28,5.17,5.08,4.87,4.79,4.79,
4.95,4.92,5.21,4.96,5.86]
#居民消费增长率
x1=[18,3.1,3.4,6.4,14.7,24.1,17.1,8.3,2.8,
-0.79,-1.41,0.42,0.69,-0.8,1.2,3.9,1.8,
1.51,4.8,5.9,-0.71,3.29,5.39,2.6,2.6,2,1.4,2]
# 国民总收入(亿元)
x2=[17001,18718,21826,26937,35260,48108,59811,
70142,78802.9,83817.6,89366.5,99066.1,109276.2,
120480.4,136576.3,161415.4,185998.9,219028.5,
270844,321500.5,348498.5,411265.2,484753.2,
539116.5,590422.4,644791.1,686449.6,740598.7]
#人均GDP（元）
x3=[1519,1644,1893,2311,2998,4044,5046,5846,
6481,6860,7229,7942,8717,9506,10666,12487,
14368,16738,20505,24121,26222,30876,36403,
40007,43852,47203,50251,53980]
#年份
year=["1989年","1990年","1991年","1992年","1993年",
"1994年","1995年","1996年","1997年","1998年",
"1999年","2000年","2001年","2002年","2003年",
"2004年","2005年","2006年","2007年","2008年",
"2009年","2010年","2011年","2012年","2013年",
"2014年","2015年","2016年"]
# 组合X1,X2,x3为3列数据
X_train = np.array(x1+x2+x3).reshape((len(x1), 3), order="F")
Y_train = np.array(y).reshape((len(y), 1))
lineModel = LinearRegression()
lineModel.fit(X_train, Y_train)
f=""
a_arr = lineModel.coef_[0]
b = lineModel.intercept_[0]
for i in range(0,len(a_arr)):
ai=a_arr[i]
if ai>=0:
ai = "+%.4f" %(ai)
else:
ai = "%.4f" % (ai)
f = f+"%s*x%s"%(ai, str(i+1))
f="y=%s+%.4f" % (f[1:],b)
print("拟合方程",f)
#对回归模型进行评分，这里简单使用训练集进行评分，实际很多时候用其他的测试集进行评分
print("得分", lineModel.score(X_train, Y_train))拟合方程 y=0.1162*x1+0.0004*x2-0.0061*x3+13.3553
得分 0.8355549431684318结果表明，在其他变量不变的情况下，当居民消费价格指数增长率上升1%，人口增长率上升0.1162%；在其他变量不变的情况下，当年国民总收入每增长1亿元，人口增长率上升0.0004%；在其他变量不变的情况下，人均GDP上升1元，人口增长率就会降低0.0061%。
设2017年居民消费价格指数增长率2.00,国民总收入760000,人均GDP为55000时，代入回归模型，得：
print(lineModel.predict([[2,760000,55000]]))[[8.46908693]]即预测得2017年人口增长率为8.469‰