赞
踩
系列文章:
写的比较仓促,代码中如有错误欢迎指正!
上一次对多元线性回归的估计以及参数和方程的显著性进行了python实现。但是这些都是建立多元线性回归的几个假设基础之上的:
若果模型违反了相应的假设就会犯对应的错误,我们在计量经济分析中的检验就是检验出是否可能犯了某一类错误,若果极有可能犯了一种错误时,我们应该怎么修正它,才能保证分析的结果是有效的。
先来考虑一个现实问题,我们在预测一个活动的KPI时,利用多个变量作为解释变量,历史KPI做被解释变量,如果解释变量中同时包含PV和UV,这样的回归是好还是坏?
1.1.完全多重共线性与不完全多重共线性
上面对于多元线性回归的假设中有一条是要求 X X X满秩,即: R a n k ( X ) = k Rank(X)=k Rank(X)=k,当X中存在一个变量完全可以由另一个变量表示的时候,变量间就出现了完全多重共线,即:
λ 1 x 1 + λ 2 x 2 + . . . + λ k x k = 0 λ_1x_1+λ_2x_2+...+λ_kx_k=0 λ1x1+λ2x2+...+λkxk=0
由秩的定义可知,此时 X X X不再满秩,即 R a n k ( X ) < k Rank(X)<k Rank(X)<k,这违反了多元线性回归的基本假定。
实际上当 x i x_i xi与 x j x_j xj之间存在很想的线性相关性时,对模型的估计(区间估计)和检验也是非常不利的。出现这种情况,我们称为不完全多重共线性
1.2.多重共线性的原因
由于经济和现实的生产过程中,许多变量往往是强相关和有相同变化趋势的。这就导致了多重共线性有事在所难免。实际上,多重共线性大多不是模型问题,而是数据问题。模型中我们大量的采用滞后变量或者变量选取不当,也会导致共线性的产生。
1.3.完全多重共线性的后果
当存在完全多重共线性时,多元线性回归可能存在一下后果:
1.4.不完全多重共线性的后果
1.4.多重共线性的检验
简单相关系数矩阵法、
我们可以计算
X
X
X中每两个变量组合之间的相关系数,来确定变变量间是否存在很强的线性相关性,利用dataframe的df.corr()
计算即可完成。
条件数
计算
(
X
′
X
)
(X'X)
(X′X)的最大特征根与最小特征根比值的平方根,当做条件数:
条件数=
m
a
x
(
λ
i
)
/
m
i
n
(
λ
i
)
\sqrt{max(λ_i)/min(λ_i)}
max(λi)/min(λi)
如果条件数超过20,即认为存在比较严重的多重共线性。(回想一下PCA的原理,和特征值
λ
i
λ_i
λi的含义,其实条件数还是很直观的)
在下面实际的python测试中,条件数达到了4573个!但是根据其他法则判断,共线性又不是十分的明显,个人感觉最好不要单一的用条件数去判断多重共线性。
方差膨胀因子
偷懒一下,直接粘贴网上的定义吧
一般
V
I
F
i
VIF_i
VIFi大于10,我们即认为变量
x
i
x_i
xi存在多重共线性。
V
I
F
VIF
VIF的优点是可以定位到存在共线性的变量,缺点是计算比较复杂。我们可以计算出方程的条件数,判断是否存在多重共线性,然后再利用X的相关系数矩阵定位存在多重共线性的变量(计算相关系数矩阵也是比较耗时的)。
1.5.多重共线性的补救
增加样本量
增加样本量可以达到减小参数估计方差的目的,在存在共线性时提高估计的可靠性
逐步回归法
每种变量组合都进行一次回归,选出最好的组合,这种做法实在是太麻烦和耗时了
PCA
我们可以对变量数据
X
X
X做一次PCA,得到新的矩阵
Z
Z
Z,Z中的主成分时正交的,完全没有共线性。但是这样缺点也很明显,回归是没有经济意义的,我偏离了我们做计量经济分析的初衷!
岭回归
在机器学习领域我们知道有两种添加正则化的回归方式:岭回归和Lasso回归两种方式,对应着L2正则和L1正则。岭回归让参数变得很小,Lasso回归让变量变得系数,可用于特征选择。。。。
在解决多重共线性中的岭回归办法和机器学习中添加L1正则的岭回归其实就是一个东西!。岭回归是最小二乘的一种改进方法,用岭回归时参数的估计变为:
β
=
(
X
′
X
+
k
I
)
−
1
X
′
y
β=(X'X+kI)^{-1}X'y
β=(X′X+kI)−1X′y
其中
k
k
k为岭回归参数(在机器学习中回应惩罚的力度),需要注意的是,岭回归的参数估计是有偏的,参数
β
β
β的方差估计也需要跟着一起变。此时
β
β
β的方差估计为:
$
V
a
r
(
β
∣
k
)
=
σ
2
(
X
′
X
+
k
I
)
−
1
∗
(
X
′
X
)
∗
(
X
′
X
+
k
I
)
Var(β|k) = σ^2(X'X+kI)^{-1}*(X'X)*(X'X+kI)
Var(β∣k)=σ2(X′X+kI)−1∗(X′X)∗(X′X+kI) $
对应的t检验的统计量也会跟着一起改变。
我们任沿用上一节的多元线性回归代码,再其基础上进行添加,我们需要实现的有:计算被解释变量的相关系数矩阵、没个变量的反差膨胀因子VIF,方程拟合后的条件数、对方程惊醒岭回归。先把框架定义好:
class mul_linear_model(object): #用于计算数据的相关系数矩阵 def data_relation(self): pass #用于计算每个变量的VIF def VIF_cul(self, ind): pass #所有变量的分行差膨胀因子检验 def VIF_test(self): pass #计算拟合的条件数 def condition_num(self): pass #岭回归 def Ridge_fit(self, k = 0.1): pass
下面就是完整的多重共线性的检验和补救方法。
# -*- coding: utf-8 -*- """ Created on Fri Apr 3 15:40:38 2020 @author: nbszg """ import numpy as np import pandas as pd import scipy.stats as st class mul_linear_model(object): #intercept定义是否有截距项,在计量经济分析中,基本都要有截距项,否则计算出的拟合优度是没有意义的 def __init__(self, y, X, intercept=True): ... #拟合模型 def fit(self): ... return self.b #预测函数 def predict(self, X_p): ... #参数T检验 def T_test(self, alpha): ... #拟合优度计算 def R_square(self): ... #整体F检验 def F_test(self, alpha): ... #用于计算数据的相关系数矩阵 def data_relation(self): #需要区分有截距项和没有截距项的情况 if self.intercept: return self.data_x.loc[:,self.data_x.columns!='intercept'].corr() else: return self.data_x.corr() #用于计算每个变量的VIF def VIF_cul(self, ind): #使用最小二乘法(X'X)的逆成X'y yk = np.array(self.X[:, ind].T) #有截距项时: if self.intercept: Xk = self.X[:, [i for i in range(1,ind)]+[j for j in range(ind+1, self.K)]] #无截距项时: else: Xk = self.X[:, [i for i in range(ind)]+[j for j in range(ind+1, self.K)]] #直接调用自身的类方法记性计算其他Xj对Xk回归的可决系数 vif_model = mul_linear_model(yk, pd.DataFrame(Xk), intercept=False) vif_model.fit(output=False) Rk_2 = vif_model.R_square() #计算方差膨胀因子 VIF_k = 1/(1-Rk_2) return VIF_k #所有变量的分行差膨胀因子检验 def VIF_test(self): #循环计算每个变量的方差膨胀因子 for i in range(1, self.K): VIF_k = self.VIF_cul(i) print("variable {0}'s VIF estimate is : {1}".format(self.columns[i], VIF_k)) #计算拟合的条件数 def condition_num(self): #有截距项时 if self.intercept: M = np.dot(self.XT[[i for i in range(1, self.K)],:], self.X[:, [j for j in range(1, self.K)]]) else: M = np.dot(self.XT, self.X) #计算X'X的特征值与特征向量 eigenvalue, eigenvector = np.linalg.eig(M) #计算条件数 cond_num = np.sqrt(eigenvalue[0]/eigenvalue[-1]) print("data X's condition number is: {0}".format(cond_num)) #岭回归 def Ridge_fit(self, k = 0.1, output=True): #使用岭回归(X'X+kI)的逆成X'y self.b = np.array(np.dot(np.dot(np.linalg.inv(np.dot(self.XT,self.X)+np.mat(np.eye(self.K))),self.XT),self.y).T) if output: for i in range(self.K): print("variable {0}'s cofe estimate is: {1}".format(self.columns[i], self.b[i][0])) #计算t检验的因子需要跟着一起变动!!! self.S = np.array(np.dot(np.dot(np.linalg.inv(np.dot(self.XT,self.X)+np.mat(np.eye(self.K))), np.dot(self.XT,self.X)), np.linalg.inv(np.dot(self.XT,self.X)+np.mat(np.eye(self.K))))) return self.b
利用上一节的数据和例子继续进行检验:
print('data columns related coefficient is: \n{0}'.format(linear_model.data_relation()))
print('\n')
linear_model.VIF_test()
print('\n')
linear_model.condition_num()
print('\n')
linear_model.Ridge_fit()
得到的检验结果和岭回归系数如下,可以看到条件数达到了4573个!但是根据其他法则判断,共线性又不是十分的明显,个人感觉最好不要单一的用条件数去判断多重共线性。
Copyright © 2003-2013 www.wpsshop.cn 版权所有,并保留所有权利。