Multi-collinearity diagnostics
mcDiagnose.RdConducts a series of checks for multicollinearity.
Author
Paul E. Johnson pauljohn@ku.edu
Examples
library(rockchalk)
N <- 100
dat <- genCorrelatedData3(y~ 0 + 0.2*x1 + 0.2*x2, N=N, means=c(100,200),
sds=c(20,30), rho=0.4, stde=10)
dat$x3 <- rnorm(100, m=40, s=4)
m1 <- lm(y ~ x1 + x2 + x3, data=dat)
summary(m1)
#>
#> Call:
#> lm(formula = y ~ x1 + x2 + x3, data = dat)
#>
#> Residuals:
#> Min 1Q Median 3Q Max
#> -25.4507 -6.3771 0.4244 6.7919 27.8159
#>
#> Coefficients:
#> Estimate Std. Error t value Pr(>|t|)
#> (Intercept) -17.22483 14.47591 -1.190 0.23702
#> x1 0.19220 0.05764 3.335 0.00121 **
#> x2 0.22485 0.04012 5.605 1.99e-07 ***
#> x3 0.29984 0.27003 1.110 0.26960
#> ---
#> Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
#>
#> Residual standard error: 10.29 on 96 degrees of freedom
#> Multiple R-squared: 0.4107, Adjusted R-squared: 0.3923
#> F-statistic: 22.3 on 3 and 96 DF, p-value: 4.84e-11
#>
m1d <- mcDiagnose(m1)
#> The following auxiliary models are being estimated and returned in a list:
#> x1 ~ x2 + x3
#> x2 ~ x1 + x3
#> x3 ~ x1 + x2
#>
#> R_j Squares of auxiliary models
#> x1 x2 x3
#> 0.15876962 0.17577317 0.03519524
#> The Corresponding VIF, 1/(1-R_j^2)
#> x1 x2 x3
#> 1.188735 1.213258 1.036479
#> Bivariate Pearson Correlations for design matrix
#> x1 x2 x3
#> x1 1.00 0.40 -0.11
#> x2 0.40 1.00 -0.18
#> x3 -0.11 -0.18 1.00
m2 <- lm(y ~ x1 * x2 + x3, data=dat)
summary(m2)
#>
#> Call:
#> lm(formula = y ~ x1 * x2 + x3, data = dat)
#>
#> Residuals:
#> Min 1Q Median 3Q Max
#> -26.5835 -7.3616 0.4972 6.9620 27.5394
#>
#> Coefficients:
#> Estimate Std. Error t value Pr(>|t|)
#> (Intercept) 8.951069 36.626590 0.244 0.807
#> x1 -0.082785 0.357998 -0.231 0.818
#> x2 0.081736 0.188227 0.434 0.665
#> x3 0.326652 0.272770 1.198 0.234
#> x1:x2 0.001428 0.001835 0.778 0.438
#>
#> Residual standard error: 10.31 on 95 degrees of freedom
#> Multiple R-squared: 0.4145, Adjusted R-squared: 0.3898
#> F-statistic: 16.81 on 4 and 95 DF, p-value: 1.883e-10
#>
m2d <- mcDiagnose(m2)
#> The following auxiliary models are being estimated and returned in a list:
#> x1 ~ x2 + x3 + `x1:x2`
#> x2 ~ x1 + x3 + `x1:x2`
#> x3 ~ x1 + x2 + `x1:x2`
#> `x1:x2` ~ x1 + x2 + x3
#>
#> R_j Squares of auxiliary models
#> x1 x2 x3 x1:x2
#> 0.97810524 0.96240164 0.05057824 0.98978980
#> The Corresponding VIF, 1/(1-R_j^2)
#> x1 x2 x3 x1:x2
#> 45.673021 26.596903 1.053273 97.941235
#> Bivariate Pearson Correlations for design matrix
#> x1 x2 x3 x1:x2
#> x1 1.00 0.40 -0.11 0.88
#> x2 0.40 1.00 -0.18 0.78
#> x3 -0.11 -0.18 1.00 -0.18
#> x1:x2 0.88 0.78 -0.18 1.00
m3 <- lm(y ~ log(10+x1) + x3 + poly(x2,2), data=dat)
summary(m3)
#>
#> Call:
#> lm(formula = y ~ log(10 + x1) + x3 + poly(x2, 2), data = dat)
#>
#> Residuals:
#> Min 1Q Median 3Q Max
#> -25.1293 -6.1131 0.1664 6.8374 27.5673
#>
#> Coefficients:
#> Estimate Std. Error t value Pr(>|t|)
#> (Intercept) -43.9508 30.7476 -1.429 0.15617
#> log(10 + x1) 19.4623 6.0634 3.210 0.00181 **
#> x3 0.2783 0.2719 1.023 0.30869
#> poly(x2, 2)1 63.5698 11.4524 5.551 2.57e-07 ***
#> poly(x2, 2)2 -2.3473 10.4215 -0.225 0.82228
#> ---
#> Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
#>
#> Residual standard error: 10.36 on 95 degrees of freedom
#> Multiple R-squared: 0.4085, Adjusted R-squared: 0.3836
#> F-statistic: 16.41 on 4 and 95 DF, p-value: 2.993e-10
#>
m3d <- mcDiagnose(m3)
#> The following auxiliary models are being estimated and returned in a list:
#> `log(10 + x1)` ~ x3 + `poly(x2, 2)1` + `poly(x2, 2)2`
#> x3 ~ `log(10 + x1)` + `poly(x2, 2)1` + `poly(x2, 2)2`
#> `poly(x2, 2)1` ~ `log(10 + x1)` + x3 + `poly(x2, 2)2`
#> `poly(x2, 2)2` ~ `log(10 + x1)` + x3 + `poly(x2, 2)1`
#>
#> R_j Squares of auxiliary models
#> log(10 + x1) x3 poly(x2, 2)1 poly(x2, 2)2
#> 0.16946919 0.03517444 0.18171293 0.01181599
#> The Corresponding VIF, 1/(1-R_j^2)
#> log(10 + x1) x3 poly(x2, 2)1 poly(x2, 2)2
#> 1.204049 1.036457 1.222065 1.011957
#> Bivariate Pearson Correlations for design matrix
#> log(10 + x1) x3 poly(x2, 2)1 poly(x2, 2)2
#> log(10 + x1) 1.00 -0.10 0.40 -0.09
#> x3 -0.10 1.00 -0.18 -0.04
#> poly(x2, 2)1 0.40 -0.18 1.00 0.00
#> poly(x2, 2)2 -0.09 -0.04 0.00 1.00
N <- 100
x1 <- 50 + rnorm(N)
x2 <- log(rgamma(N, 2,1))
x3 <- rpois(N, lambda=17)
z1 <- gl(5, N/5)
dummies <- contrasts(z1)[ as.numeric(z1), ]
dimnames(dummies) <- NULL ## Avoids row name conflict in data.frame below
y3 <- x1 -.5 * x2 + 0.1 * x2^2 + dummies %*% c(0.1,-0.1,-0.2,0.2)+ 5 * rnorm(N)
dat <- data.frame(x1=x1, x2=x2, x3=x3, z1=z1, y3 = y3)
m3 <- lm(y3 ~ x1 + poly(x2,2) + log(x1) + z1, dat)
summary(m3)
#>
#> Call:
#> lm(formula = y3 ~ x1 + poly(x2, 2) + log(x1) + z1, data = dat)
#>
#> Residuals:
#> Min 1Q Median 3Q Max
#> -16.2968 -3.7781 0.1393 3.6930 13.8968
#>
#> Coefficients:
#> Estimate Std. Error t value Pr(>|t|)
#> (Intercept) 6366.1227 5897.4133 1.079 0.283
#> x1 44.8789 40.5480 1.107 0.271
#> poly(x2, 2)1 -9.1033 5.6107 -1.622 0.108
#> poly(x2, 2)2 12.0999 5.7158 2.117 0.037 *
#> log(x1) -2188.2967 2025.7917 -1.080 0.283
#> z12 0.9709 1.8133 0.535 0.594
#> z13 0.3067 1.7756 0.173 0.863
#> z14 -1.6760 1.7956 -0.933 0.353
#> z15 0.8750 1.7728 0.494 0.623
#> ---
#> Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
#>
#> Residual standard error: 5.579 on 91 degrees of freedom
#> Multiple R-squared: 0.1472, Adjusted R-squared: 0.07221
#> F-statistic: 1.963 on 8 and 91 DF, p-value: 0.06003
#>
mcDiagnose(m3)
#> The following auxiliary models are being estimated and returned in a list:
#> x1 ~ `poly(x2, 2)1` + `poly(x2, 2)2` + `log(x1)` + z12 + z13 +
#> z14 + z15
#> `poly(x2, 2)1` ~ x1 + `poly(x2, 2)2` + `log(x1)` + z12 + z13 +
#> z14 + z15
#> `poly(x2, 2)2` ~ x1 + `poly(x2, 2)1` + `log(x1)` + z12 + z13 +
#> z14 + z15
#> `log(x1)` ~ x1 + `poly(x2, 2)1` + `poly(x2, 2)2` + z12 + z13 +
#> z14 + z15
#> z12 ~ x1 + `poly(x2, 2)1` + `poly(x2, 2)2` + `log(x1)` + z13 +
#> z14 + z15
#> z13 ~ x1 + `poly(x2, 2)1` + `poly(x2, 2)2` + `log(x1)` + z12 +
#> z14 + z15
#> z14 ~ x1 + `poly(x2, 2)1` + `poly(x2, 2)2` + `log(x1)` + z12 +
#> z13 + z15
#> z15 ~ x1 + `poly(x2, 2)1` + `poly(x2, 2)2` + `log(x1)` + z12 +
#> z13 + z14
#>
#> R_j Squares of auxiliary models
#> x1 poly(x2, 2)1 poly(x2, 2)2 log(x1) z12 z13
#> 0.99980191 0.01125798 0.04728646 0.99980192 0.40835047 0.38298425
#> z14 z15
#> 0.39663613 0.38101947
#> The Corresponding VIF, 1/(1-R_j^2)
#> x1 poly(x2, 2)1 poly(x2, 2)2 log(x1) z12 z13
#> 5048.192161 1.011386 1.049633 5048.403009 1.690190 1.620704
#> z14 z15
#> 1.657375 1.615560
#> Bivariate Pearson Correlations for design matrix
#> x1 poly(x2, 2)1 poly(x2, 2)2 log(x1) z12 z13 z14 z15
#> x1 1.00 0.07 0.04 1.00 0.21 -0.05 -0.12 -0.04
#> poly(x2, 2)1 0.07 1.00 0.00 0.07 0.02 -0.06 0.02 0.05
#> poly(x2, 2)2 0.04 0.00 1.00 0.04 0.14 -0.01 0.08 -0.06
#> log(x1) 1.00 0.07 0.04 1.00 0.21 -0.05 -0.12 -0.04
#> z12 0.21 0.02 0.14 0.21 1.00 -0.25 -0.25 -0.25
#> z13 -0.05 -0.06 -0.01 -0.05 -0.25 1.00 -0.25 -0.25
#> z14 -0.12 0.02 0.08 -0.12 -0.25 -0.25 1.00 -0.25
#> z15 -0.04 0.05 -0.06 -0.04 -0.25 -0.25 -0.25 1.00