VIF(분산팽창요인), 결정계수

이전 포스팅에서 다중공선성에 대해서 다루었습니다. 그렇다면 다중공선성을 어떻게 진단을 할 수 있을까요? 완전한 다중공선성(Perfect Multicollinearity)라면, 계수가 추정이 안되거나 경고문구를 보고 확인할 수 있었지만, Imperfect Multicollinearity라면 우리는 사전에 진단할 수 있는 방법이 필요합니다.

[VIF(분산팽창요인)]

이러한 다중공선성을 진단하기 위해 우리는 VIF라는 개념을 많이 활용합니다. VIF는 독립변수가 여러개 있을 때, 특정 독립변수를 종속변수로하고 나머지 독립변수를 독립변수로 하여 회귀분석을 수행하여 변수간에 관계성을 측정합니다. 

예를들어 독립변수 $X_1$,$X_2$,$X_3$이 있다고 가정해보겠습니다. 그렇다면 우리는 아래와 같이 3번의 회귀분석을 수행할 수 있씁니다.

$X_1$ ~ $X_2$+$X_3$ ---(1)

$X_2$ ~ $X_1$+$X_3$ ---(2)

$X_3$ ~ $X_1$+$X_2$ ---(3)

3번의 회귀분석을 수행하여, 각 회귀식에서 결정계수를 추정하고, 결정계수를 활용하여 VIF를 결정합니다. 첫 번째식의 결정계수를 $R_1^2$라고 한다면, 첫번째 독립변수의 VIF는 아래와 같이 정의 됩니다.

$VIF_1$=1/(1-$R_1^2$)

위의 식을 일반화 하면, j번째 변수의 VIF 값은 아래와 같습니다.

$VIF_j$=1/(1-$R_j^2$)

다중공선성을 정리했을 때(direction-f.tistory.com/45), 독립변수간에 상관관계가 있다면 추정한 계수들의 분산이 커진다는 것을 우리는 확인했습니다. $VIF_j$값은 얼마나 커지는지에 대한 값을 대표한다고 볼 수 있습니다. 다시 말하면, 독립변수간 상관관계가 있을때 $VIF_j$값만큼 $Var(\beta_j)$가 커지게됩니다.

독립변수간 상관관계가 있을 때 j번째 변수의 분산은 아래와 같이 정의 됩니다.

https://en.wikipedia.org/wiki/Variance_inflation_factor

$R_j^2$값은 항상 1보다 작기때문에 독립변수간 상관관계가 있을 때 $Var(\beta_j)$는 항상 커지게 됩니다.

최종적으로 정리하게 되면, 독립변수간 상관관계가 있으면 추정된 계수들의 분산은 무조건 커지게 되고, 얼마나 커지는지를 보는 것이 VIF값입니다. 따라서 VIF가 크다는 것은 다른 변수들에 의해서 설명이 많이 된다는 것이고 회귀분석시 독립변수에서 제외되어도 해당 회귀모형은 종속변수를 설명하는데 충분한 설명력을 가질 수 있다라는 것을 뜻합니다.

그렇다면, Python에서 VIF를 실습하기 전에, 결정계수가 무엇인지 정의하고 실습을 진행하겠습니다.

[결정계수]

결정계수는 아래와 같이 정의 됩니다.

해당 결정계수 ($R^2$)는 독립변수가 늘어날 수록 커지게 됩니다. 따라서 이러한 단점을 보완하기 위해 Adjusted-$R^2$을 활용하며, Adjusted-$R^2$는 아래와 같이 표현됩니다.

이 때, $n$은 데이터의 수이며, $k$는 독립변수의 수입니다. 즉 $k$가 커질수록 $\widehat{R^2}$가 작아지도록 설계되어 있습니다.

Python을 통해 먼저 결정계수를 추정해보고 VIF를 계산해보겠습니다.

import statsmodels.api as sm
import statsmodels.formula.api as smf
import matplotlib.pyplot as plt
import seaborn as sns
import numpy as np
import pandas as pd
from math import sqrt
import random
import scipy.stats as stats

## load the data
CASchools = pd.read_csv("CASchools.csv")

## define variable
CASchools["STR"]=CASchools["students"]/CASchools["teachers"]
CASchools["score"]=(CASchools["read"]+CASchools["math"])/2
CASchools["FracEL"]=CASchools["english"]/100

## adjusted R-square by hand
multi_mod_ = smf.ols("score ~ STR+english",data=CASchools).fit()
multi_mod_.summary()

n=len(CASchools)
k=2

y_mean = CASchools["score"].mean()

SSR = multi_mod_.resid.pow(2).sum()
TSS = (CASchools["score"]-y_mean).pow(2).sum()
ESS = (multi_mod_.fittedvalues-y_mean).pow(2).sum()

Rsq = 1-(SSR/TSS)
adf_Rsq = 1-(n-1)/(n-k-1)*SSR/TSS

summary()를 통해 확인한 $R^2$, $\widehat{R^2}$값과 직접 계산한 값과 동일함을 확인 할 수 있습니다.

추가적으로 VIF를 계산해보겠습니다. statsmodels 패키지를 활용하였으며, dmatrix를 활용하여 Input/output 메트릭스를 구성하였습니다. dmatrix를 활용하면 손쉽게 우리가 활용할 Input Matrix와 Output Matrix를 구현할 수 있습니다.

from patsy import dmatrices
from statsmodels.stats.outliers_influence import variance_inflation_factor

y, X = dmatrices('score ~ STR + english+read', data=CASchools, return_type = 'dataframe')

vif = pd.DataFrame()
vif["VIF Factor"] = [variance_inflation_factor(X.values, i) for i in range(X.shape[1])]
vif["features"] = X.columns 
vif

'''
    VIF Factor   features
0  2440.389705  Intercept
1     1.065403        STR
2     1.911380    english
3     1.963476       read
'''

위의 결과를 보면 우리가 활용한 변수가 모두 VIF가 2미만으로 다중공선성 문제는 일으키지 않을 것으로 판단 할 수 있습니다. 일반적으로 VIF값이 10이상이면 문제가 있을 수 있다고 판단하고, 제외합니다.

참조:

www.econometrics-with-r.org/index.html

bkshin.tistory.com/entry/DATA-20-%EB%8B%A4%EC%A4%91%EA%B3%B5%EC%84%A0%EC%84%B1%EA%B3%BC-VIF