自己造轮子——计算VIF_基于statsmodels

我知道用stata和R的更多，此处仅基于兴趣摸索，本文数据和方法均来自于朱顺全老师的教程

核心代码

书上老师使用的方法是一个一个计算，可维护性有点低，所以我就优化了一下

X_with_const = sm.add_constant(X)
vif_values = [variance_inflation_factor(X_with_const.values, i) for i in range(X_with_const.shape[1])]

使用了statsmodels.stats.outliers_influence 模块里面的 variance_inflation_factor 函数

首先，在X矩阵中添加常数项以符合 variance_inflation_factor 函数的计算要求

然后，使用 X_with_const.value 来返回 variance_inflation_factor 函数计算需要的 NumPy 数组

然后，使用 range(X_with_const.shape[1]) 来生成一个从 0 到 X_with_const 矩阵的列数（即特征数加一，因为有常数项）减一的整数序列。因为 shape[1] 返回的是矩阵的列数，所以这个序列包含了从 0 到列数减一的所有索引。

最后，使用 variance_inflation_factor 函数来计算给定 X_with_const.values 中第 i 个特征的 VIF

完整代码

import pandas as pd
import numpy as np
import statsmodels.api as sm
from statsmodels.stats.outliers_influence import variance_inflation_factor
from patsy import dmatrices

#### 数据导入
data=pd.DataFrame()
data=pd.DataFrame(pd.read_excel('C:\\Users\\cyx94a\\Desktop\\Python\\data\\al7-1.xls'))
vars=['run','z1','z2','z3']
df=data[vars]

#### 生成设计矩阵
y,X=dmatrices('run~z1+z2+z3',data=df,return_type='dataframe')

#### 描述性统计和相关分析
print(data.describe())
print(data.corr())

#### OLS
model=sm.OLS(y,X)
fit=model.fit()
print(fit.summary())

#### 计算VIF
X_with_const = sm.add_constant(X) # 为矩阵X增加常数项

vif_values = [variance_inflation_factor(X_with_const.values, i) for i in range(X_with_const.shape[1])] # 计算每个特征的 VIF

print("VIF values:", vif_values) # 从列表的第二个数字及之后就是各个自变量的VIF值