# サンプル作成
Rで確率・統計:多重共線性:分散拡大係数(variance inflation factor, VIF)
Rでデータサイエンス
多重共線性:分散拡大係数(variance inflation factor, VIF)
# 重回帰分析のサンプルは次の通り。
# サンプルサイズは30、定数項と説明変数(x1,x2,x3)の係数は下記。
<- 30
n <- 5
const <- c(1.5, -2, 4) # a1,a2,a3
coeff # x1を作成。標準偏差の数値に意味はなし(以降、同様)
<- rnorm(n = n, sd = 3)
x1 # x2を作成。x1との線形結合。係数と定数項の数値に意味はなし。
<- 4 * x1 + rnorm(n = n, sd = 2)
x2 # x3を作成。
<- rnorm(n = n, sd = 2)
x3 # 説明変数を行列に。
<- cbind(x1, x2, x3)
X # yを作成。
<- (X %*% coeff) %>%
y0 as.vector() #;y0
<- y0 + const + rnorm(n = n) #;y
y # (y0-(x1*coeff[1]+x2*coeff[2]+x3*coeff[3])) %>% unique()
# 重相関分析
<- y ~ x1 + x2 + x3
formula <- lm(formula = formula)
lmresult %>%
lmresult ::summ(confint = T) jtools
Observations | 30 |
Dependent variable | y |
Type | OLS linear regression |
F(3,26) | 6140.68 |
R² | 1.00 |
Adj. R² | 1.00 |
Est. | 2.5% | 97.5% | t val. | p | |
---|---|---|---|---|---|
(Intercept) | 4.91 | 4.51 | 5.32 | 24.91 | 0.00 |
x1 | 1.51 | 0.56 | 2.45 | 3.28 | 0.00 |
x2 | -2.01 | -2.25 | -1.78 | -17.31 | 0.00 |
x3 | 4.05 | 3.87 | 4.23 | 46.47 | 0.00 |
Standard errors: OLS |
# 但し、x1とx2には有意な線形相関があり。
<- lm(x2 ~ x1)
lmresultx2x1 %>%
lmresultx2x1 ::summ(confint = T)
jtoolspar(mfrow = c(1, 2))
plot(x1, x2)
abline(lmresultx2x1, col = "red")
$residuals %>%
lmresultx2x1pacf()
Observations | 30 |
Dependent variable | x2 |
Type | OLS linear regression |
F(1,28) | 1952.25 |
R² | 0.99 |
Adj. R² | 0.99 |
Est. | 2.5% | 97.5% | t val. | p | |
---|---|---|---|---|---|
(Intercept) | -0.23 | -0.89 | 0.43 | -0.72 | 0.48 |
x1 | 3.93 | 3.75 | 4.11 | 44.18 | 0.00 |
Standard errors: OLS |
# 分散拡大係数をパッケージcarの関数vifで確認。
# 引数mod(method)にVIFを求める重相関分析の結果、ここではlmresultを指定。
# x1とx2には多重共線性が現れるはず。
::vif(mod = lmresult) car
x1 x2 x3
72.395449 72.792073 1.047477
# x1、x2は共に10を超えている。
# VIFの基準として10以上を問題にしている。https://www.heisei-u.ac.jp/ba/fukui/tips/tip006.pdf
# つまりx1、x2の多重共線性が示唆される。
# VIFを求める関数
<- function(formula) {
fun_vif # 対象とする説明変数を「その他の説明変数」で線形回帰し、その決定係数を求める。
<- lm(formula = formula) %>%
lmresult summary()
<- 1/(1 - lmresult$r.squared)
vif # 同時に重相関係数も求める。
<- lmresult$r.squared %>%
r sqrt()
return(list(vif = vif, r = r))
}
# x1のVIF。x1を目的変数、x2、x3を説明変数とし線形回帰を取る。
<- x1 ~ x2 + x3
formula fun_vif(formula = formula)
# x1のVIFは関数vif{car}の結果と同じ。
$vif
[1] 72.39545
$r
[1] 0.9930695
# x2のVIF
<- x2 ~ x1 + x3
formula fun_vif(formula = formula)
# x2のVIFは関数vif{car}の結果と同じ。
$vif
[1] 72.79207
$r
[1] 0.9931074
# 最後はx3のVIF
<- x3 ~ x1 + x2
formula fun_vif(formula = formula)
# x3のVIFは関数vif{car}の結果と同じ。
$vif
[1] 1.047477
$r
[1] 0.2128971
# VIFが10の場合の相関係数を最後に確認。
1 - 1/10)^0.5 (
[1] 0.9486833
参考引用資料
最終更新
Sys.time()
[1] "2024-04-23 06:25:36 JST"
R、Quarto、Package
R.Version()$version.string
[1] "R version 4.3.3 (2024-02-29 ucrt)"
::quarto_version() quarto
[1] '1.4.553'
packageVersion(pkg = "tidyverse")
[1] '2.0.0'