【Python计量】多重共线性的检验和解决

此文章首发于公众号：Python for Finance

链接：https://mp.weixin.qq.com/s/uvYV8qNLg0lWUdWIL-1vmA

我们以伍德里奇《计量经济学导论：现代方法》的”第4章 多元回归分析：推断”的案例4.1为例，利用WAGE1中526个工人的观测数据，研究其工资和受教育程度、工作经验和任职时长的关系。构建多元线性回归模型如下：
l o g ( w a g e ) = β 0 + β 1 e d u c + β 2 e x p e r + β 3 t e n u r e + u log(wage)=\beta_0+\beta_1educ+\beta_2exper+\beta_3tenure+u l o g (w a g e )=β0 ​+β1 ​e d u c +β2 ​e x p er +β3 ​t e n u re +u

• log(wage)：工资的对数，也记作lwage
• educ：受教育程度
• exper：工作经验
• tenure：任职时长

代码如下：

import wooldridge as woo
import statsmodels.formula.api as smf

wage1 = woo.dataWoo('wage1')

reg = smf.ols(formula='lwage ~ educ + exper + tenure', data = wage1)
results = reg.fit()
print(results.summary())

模型回归结果如下：
l o g ( w a g e ) ^ = 0.2844 + 0.0920 e d u c + 0.0041 e x p e r + 0.0221 t e n u r e + u \widehat{log(wage)}=0.2844+0.0920educ+0.0041exper+0.0221tenure+u l o g (w a g e )​=0.2844 +0.0920 e d u c +0.0041 e x p er +0.0221 t e n u re +u
本例研究自变量educ、exper、tenure之间是否存在多重共线性。

一般地，产生多重共线性的主要原因有以下几个方面：

• 模型设定错误：这个不需要解释，模型设定错误的时候什么情况都有可能发生，一定要克服。
• 数据采集方法不当：如果在总体中的一个较小的范围内抽样，自变量没有显著的波动会导致自变量和截距项之间产生多重共线性的现象。
• 经济变量之间具有共同变化趋势：如时间序列数据中，GDP、就业人口、消费等数据的变化常常会具有相同的时间趋势。
• 模型中包含滞后变量：在经济计量模型中，往往需要引入滞后经济变量来反映真实的经济关系。例如，当期收入和前期收入间有较强的线性相关性。
• 多项式的引入：如模型中包括X X X,X 2 X^2 X 2,X 3 X^3 X 3作为解释变量，当 X变化不大时会呈现出严重的多重共线性。

1、简单相关系数法

简单相关系数检验法是利用解释变量之间的线性相关程度去判断是否存在严重多重共线性的一种简便方法。一般而言，如果每两个解释变量的简单相关系数比较高，例如 ∣ r ∣ > 0.8 ，则可认为存在着较严重的多重共线性。

但我们需要注意的是，较高的简单相关系数只是多重共线性存在的充分条件，而不是必要条件。特别是在多于两个解释变量的回归模型中，有时较低的相关系数也可能存在多重共线性。因此并不能简单地依据相关系数进行多重共线性的准确判断。

data = wage1[['educ','exper','tenure']]
print(data.corr())

结果为：

            educ     exper    tenure
educ    1.000000 -0.299542 -0.056173
exper  -0.299542  1.000000  0.499291
tenure -0.056173  0.499291  1.000000

本例中，自变量的相关系数均小于0.5。

2、方差膨胀因子法(variance inflation factor，VIF)

V I F i = 1 1 − R i 2 VIF_i=\frac{1}{1-R_i^2}V I F i ​=1 −R i 2 ​1 ​

VIF值越大则共线性问题越明显，通常以10作为判断边界。当VIF

方法一：自行计算VIF

reg1 = smf.ols(formula=' educ ~ exper + tenure', data = wage1)
results1 = reg1.fit()
vif1 = (1 - results1.rsquared)**(-1)

reg2 = smf.ols(formula=' exper ~ educ + tenure', data = wage1)
results2 = reg2.fit()
vif2 = (1 - results2.rsquared)**(-1)

reg3 = smf.ols(formula=' tenure ~ educ + exper', data = wage1)
results3 = reg3.fit()
vif3 = (1 - results3.rsquared)**(-1)

print(vif1,vif2,vif3)

方法二：运用statsmodel库中的variance_inflation_factor,可以简便的计算所有变量的vif

statsmodels.stats.outliers_influence.variance_inflationfact(exog,exogidx)

exog即为自变量矩阵，exog_idx为自变量所在的列

from statsmodels.stats.outliers_influence import variance_inflation_factor

exog = reg.exog
exog_names = reg.exog_names

vif = []
for i in range(exog.shape[1]-1):
    vif.append(variance_inflation_factor(exog,i+1))
    print(f'{exog_names[i+1]} 的方差膨胀因子VIF = {vif[i]}')

运行结果为：

educ 的方差膨胀因子VIF = 1.1127707502838835
exper 的方差膨胀因子VIF = 1.4776177726317783
tenure 的方差膨胀因子VIF = 1.3492955605611776

我们也可以制作VIF的表格，编制代码如下：

from statsmodels.stats.outliers_influence import variance_inflation_factor

exog = reg.exog
exog_names = reg.exog_names

vif_info = pd.DataFrame()
vif_info['exog'] = exog_names[1:]
vif_info['VIF'] = [variance_inflation_factor(exog, i+1) for i in range(exog.shape[1]-1)]
vif_info = vif_info.sort_values('VIF', ascending=False)
print(vif_info)

结果为：

如果在模型中引入解释变量的平方项，则容易引起多重共线性。考虑在上述回归中加入工作经验exper的平方，在进行多重共线性检验。

#建立回归模型:
reg = smf.ols(formula='lwage ~ educ + exper + tenure + I(exper**2)', data = wage1)
results = reg.fit()

exog = reg.exog #获取自变量矩阵（第一列为截距项）
exog_names = reg.exog_names #获取自变量名称（第一个为截距项）

vif_info = pd.DataFrame()
vif_info['exog'] = exog_names[1:]
vif_info['VIF'] = [variance_inflation_factor(exog, i+1) for i in range(exog.shape[1]-1)]
vif_info = vif_info.sort_values('VIF', ascending=False)
print(vif_info)

结果为：

            exog        VIF
1          exper  13.866460
3  I(exper ** 2)  13.492686
2         tenure   1.355685
0           educ   1.140644

从上表可知，变量e x p e r exper e x p er和e x p e r 2 exper^2 e x p e r 2的VIF均大于10，故存在多重共线性。进一步，将e x p e r 2 exper^2 e x p e r 2对e x p e r exper e x p er进行回归。

reg = smf.ols(formula=' I(exper**2) ~  exper ', data = wage1)
results = reg.fit()
print(results.summary())

结果为：

                            OLS Regression Results
==============================================================================
Dep. Variable:          I(exper ** 2)   R-squared:                       0.923
Model:                            OLS   Adj. R-squared:                  0.923
Method:                 Least Squares   F-statistic:                     6323.

Date:                Sun, 07 Aug 2022   Prob (F-statistic):          1.35e-294
Time:                        16:23:54   Log-Likelihood:                -3448.6
No. Observations:                 526   AIC:                             6901.

Df Residuals:                     524   BIC:                             6910.

Df Model:                           1
Covariance Type:            nonrobust
==============================================================================
                 coef    std err          t      P>|t|      [0.025      0.975]
Intercept      0.9981      0.040     24.718      0.000       0.919       1.077
exper_s        0.7062      0.040     17.472      0.000       0.627       0.786
==============================================================================
Omnibus:                       44.401   Durbin-Watson:                   1.867
Prob(Omnibus):                  0.000   Jarque-Bera (JB):               53.524
Skew:                           0.770   Prob(JB):                     2.39e-12
Kurtosis:                       3.264   Cond. No.                         1.00
==============================================================================

Notes:
[1] Standard Errors assume that the covariance matrix of the errors is correctly specified.

上表显示，此回归的R 2 R^2 R 2为0.268，相对于e x p e r 2 exper^2 e x p e r 2对e x p e r exper e x p er回归R 2 R^2 R 2(0.923)大大下降。

Original: https://blog.csdn.net/mfsdmlove/article/details/126213172
Author: Python for Finance
Title: 【Python计量】多重共线性的检验和解决