Prueba de Multicolinealidad
KENNEDY_ULISES_VASQUEZ_PEREZ
Carga de Datos
library(wooldridge)
data(hprice1)
head(force(hprice1),n=5) #mostrar las primeras 5 observaciones## price assess bdrms lotsize sqrft colonial lprice lassess llotsize lsqrft
## 1 300 349.1 4 6126 2438 1 5.703783 5.855359 8.720297 7.798934
## 2 370 351.5 3 9903 2076 1 5.913503 5.862210 9.200593 7.638198
## 3 191 217.7 3 5200 1374 0 5.252274 5.383118 8.556414 7.225482
## 4 195 231.8 3 4600 1448 1 5.273000 5.445875 8.433811 7.277938
## 5 373 319.1 4 6095 2514 1 5.921578 5.765504 8.715224 7.829630- Estimacion del Modelo
library(stargazer)
modelo_estimado<-lm(price~lotsize+sqrft+bdrms, data = hprice1)
stargazer(modelo_estimado, title = "Prueba de Multicolinealidad", type = "text")##
## Prueba de Multicolinealidad
## ===============================================
## Dependent variable:
## ---------------------------
## price
## -----------------------------------------------
## lotsize 0.002***
## (0.001)
##
## sqrft 0.123***
## (0.013)
##
## bdrms 13.853
## (9.010)
##
## Constant -21.770
## (29.475)
##
## -----------------------------------------------
## Observations 88
## R2 0.672
## Adjusted R2 0.661
## Residual Std. Error 59.833 (df = 84)
## F Statistic 57.460*** (df = 3; 84)
## ===============================================
## Note: *p<0.1; **p<0.05; ***p<0.012. Verifique si hay evidencia de la indenpencia de los regrosores (no colinealidad), a traves de: a) Indice de condicion y prueba de FG, presente sus resultadod de manera tabular en ambos casos y para la prueba FG presente tambien sus resultados de forma grafica usando la libreria fastaGrahp
Calculo Manual
library(stargazer)
x_mat<-model.matrix(modelo_estimado)
stargazer(head(x_mat, n=6), type = "text")##
## =================================
## (Intercept) lotsize sqrft bdrms
## ---------------------------------
## 1 1 6,126 2,438 4
## 2 1 9,903 2,076 3
## 3 1 5,200 1,374 3
## 4 1 4,600 1,448 3
## 5 1 6,095 2,514 4
## 6 1 8,566 2,754 5
## ---------------------------------xx_matrix<-t(x_mat) %*% x_mat
stargazer(xx_matrix, type = "text")##
## ==============================================================
## (Intercept) lotsize sqrft bdrms
## --------------------------------------------------------------
## (Intercept) 88 793,748 177,205 314
## lotsize 793,748 16,165,159,010 1,692,290,257 2,933,767
## sqrft 177,205 1,692,290,257 385,820,561 654,755
## bdrms 314 2,933,767 654,755 1,182
## --------------------------------------------------------------library(stargazer)
options(scipen = 9999)
Sn<-solve(diag(sqrt(diag(xx_matrix))))
stargazer(Sn, type = "text")##
## ==========================
## 0.107 0 0 0
## 0 0.00001 0 0
## 0 0 0.0001 0
## 0 0 0 0.029
## --------------------------X^tX Normalizada: Calcular la Matriz Normalizada
library(stargazer)
xx_norm<-(Sn%*%xx_matrix)%*%Sn
stargazer(xx_norm, type = "text", digits = 4)##
## ===========================
## 1 0.6655 0.9617 0.9736
## 0.6655 1 0.6776 0.6712
## 0.9617 0.6776 1 0.9696
## 0.9736 0.6712 0.9696 1
## ---------------------------Autovalores de X^tX
library(stargazer)
#autovalores
lamdas<-eigen(xx_norm, symmetric = TRUE)
stargazer(lamdas$values, type = "text")##
## =======================
## 3.482 0.455 0.039 0.025
## -----------------------CALCULO DE K(X)=√(λmax/λmin)
k<-sqrt(max(lamdas$values)/min(lamdas$values))
print(k)## [1] 11.86778Sí κ(x) es inferior o igual a 20, la multicolinealidad es leve, no se considera un problema. Para 20<κ(x)<30, la multicolinealidad se considera moderada. En el caso de que κ(x)≥30 la multicolinealidad es severa.
INTERPRETACION * Como k(x)<=20 con un valor de 11.86 la multicolinealidad no se considera un problema*
Calculo del indice de condicion utilizando la libreria mctest
library(mctest)
x_mat<-model.matrix(modelo_estimado)
mctest(mod = modelo_estimado)##
## Call:
## omcdiag(mod = mod, Inter = TRUE, detr = detr, red = red, conf = conf,
## theil = theil, cn = cn)
##
##
## Overall Multicollinearity Diagnostics
##
## MC Results detection
## Determinant |X'X|: 0.6918 0
## Farrar Chi-Square: 31.3812 1
## Red Indicator: 0.3341 0
## Sum of Lambda Inverse: 3.8525 0
## Theil's Method: -0.7297 0
## Condition Number: 11.8678 0
##
## 1 --> COLLINEARITY is detected by the test
## 0 --> COLLINEARITY is not detected by the testComo se puede observar donde dice “condition Number” es el indice y podemos concluir que: * Como k(x)<=20 con un valor de 11.86 la multicolinealidad no se considera un problema.
Calculo del Indice de Condicion usando la libreria “Otsrr”.
library(olsrr)
ols_eigen_cindex(model=modelo_estimado)## Eigenvalue Condition Index intercept lotsize sqrft bdrms
## 1 3.48158596 1.000000 0.003663034 0.0277802824 0.004156293 0.002939554
## 2 0.45518380 2.765637 0.006800735 0.9670803174 0.006067321 0.005096396
## 3 0.03851083 9.508174 0.472581427 0.0051085488 0.816079307 0.016938178
## 4 0.02471941 11.867781 0.516954804 0.0000308514 0.173697079 0.975025872CALCULO MANUAL PARA REALIZAR LA PRUEBA DE FARRAR-GLAUBAR Caculo lRl
library(stargazer)
Zn<-scale(x_mat[, -1])
stargazer(head(Zn, n=6), type = "text")##
## =======================
## lotsize sqrft bdrms
## -----------------------
## 1 -0.284 0.735 0.513
## 2 0.087 0.108 -0.675
## 3 -0.375 -1.108 -0.675
## 4 -0.434 -0.980 -0.675
## 5 -0.287 0.867 0.513
## 6 -0.045 1.283 1.702
## -----------------------Calcular la Matriz R
library(stargazer)
n<-nrow(Zn)
R<-(t(Zn)%*%Zn)*(1/(n-1))
stargazer(R, type = "text", digits = 4)##
## =============================
## lotsize sqrft bdrms
## -----------------------------
## lotsize 1 0.1838 0.1363
## sqrft 0.1838 1 0.5315
## bdrms 0.1363 0.5315 1
## -----------------------------Calcular lRl
determinante_R<-det(R)
print(determinante_R)## [1] 0.6917931Aplicando la prueba de Farrer Glaubar (Barlett)
m<-ncol(x_mat[, -1])
n<-nrow(x_mat[, -1])
chi_FG<--(n-1-(2*m+5)/6)*log(determinante_R)
print(chi_FG)## [1] 31.38122Encontrar el valor critico V.C.
gl<-m*(m-1)/2
VC<-qchisq(p=0.95, df=gl)
print(VC)## [1] 7.814728Regla de desicion: Como X^2FG es de 31.38122 >= que el valor critico de 7.814728 se rechaza H0, por lo tanto hay evidencia de colinealidad en los regresores.
** Calculando del estadistico utilizando la libreria “mctest and pysich”
Calculo de FG usando mctest
library(mctest)
mctest::omcdiag(mod = modelo_estimado)##
## Call:
## mctest::omcdiag(mod = modelo_estimado)
##
##
## Overall Multicollinearity Diagnostics
##
## MC Results detection
## Determinant |X'X|: 0.6918 0
## Farrar Chi-Square: 31.3812 1
## Red Indicator: 0.3341 0
## Sum of Lambda Inverse: 3.8525 0
## Theil's Method: -0.7297 0
## Condition Number: 11.8678 0
##
## 1 --> COLLINEARITY is detected by the test
## 0 --> COLLINEARITY is not detected by the testy nuevamente tenemos que el Chi-Square es de 31.3812 superior al valor critico de 7.8 encontrado por lo tanto se detecta como bien lo dice el diagnostico colinealidad en el modelo.
Calculo de FG usando la pysych
library(psych)
FG_test<-cortest.bartlett(x_mat[, -1])
print(FG_test)## $chisq
## [1] 31.38122
##
## $p.value
## [1] 0.0000007065806
##
## $df
## [1] 3Resultado de forma grafica
library(fastGraph)## Warning: package 'fastGraph' was built under R version 4.5.3# Usamos tus variables calculadas: chi_FG (31.38) y gl (3)
shadeDist(xshade = chi_FG,
ddist = "dchisq",
parm1 = gl,
lower.tail = FALSE,
main = "Prueba de Farrar-Glauber (Bartlett)",
sub = paste("Chi-cuadrado =", round(chi_FG, 4)))
Nuevamente se confirma que es de 31.38122 el estadistico FG
b) Factores Inflacionarios de la varianza, presente sus resultados de forma tabular y de forma grafica
Referencia entre RJ^2
library(dplyr)## Warning: package 'dplyr' was built under R version 4.5.3##
## Adjuntando el paquete: 'dplyr'## The following objects are masked from 'package:stats':
##
## filter, lag## The following objects are masked from 'package:base':
##
## intersect, setdiff, setequal, unionR.cuadrado.regresores<-c(0,0.5,.8,.9)
as.data.frame(R.cuadrado.regresores) %>% mutate(VIF=1/(1-R.cuadrado.regresores))## R.cuadrado.regresores VIF
## 1 0.0 1
## 2 0.5 2
## 3 0.8 5
## 4 0.9 10Cálculo manual: Matriz de Correlación de los regresores del modelo (Como se obtuvo con anterioridad):
print(R)## lotsize sqrft bdrms
## lotsize 1.0000000 0.1838422 0.1363256
## sqrft 0.1838422 1.0000000 0.5314736
## bdrms 0.1363256 0.5314736 1.0000000Inversa de la matriz de correlación R−1:
inversa_R<-solve(R)
print(inversa_R)## lotsize sqrft bdrms
## lotsize 1.03721145 -0.1610145 -0.05582352
## sqrft -0.16101454 1.4186543 -0.73202696
## bdrms -0.05582352 -0.7320270 1.39666321VIF’s para el modelo estimado:
VIFs<-diag(inversa_R)
print(VIFs)## lotsize sqrft bdrms
## 1.037211 1.418654 1.396663Cálculo de los VIF’s usando “performance”
library(performance)
VIFs<-multicollinearity(x = modelo_estimado,verbose = FALSE)
VIFs## # Check for Multicollinearity
##
## Low Correlation
##
## Term VIF VIF 95% CI adj. VIF Tolerance Tolerance 95% CI
## lotsize 1.04 [1.00, 11.02] 1.02 0.96 [0.09, 1.00]
## sqrft 1.42 [1.18, 1.98] 1.19 0.70 [0.51, 0.85]
## bdrms 1.40 [1.17, 1.95] 1.18 0.72 [0.51, 0.86]# Convertimos el objeto de 'performance' a un vector numérico para barplot
vif_vector <- as.numeric(VIFs$VIF)
names(vif_vector) <- VIFs$Term # Le ponemos los nombres de las variables
# Ahora graficamos usando el vector
color_palette <- ifelse(vif_vector > 5, "red", "skyblue")
barplot(vif_vector,
main = "Factores Inflacionarios de la Varianza (VIF)",
col = color_palette,
las = 1,
ylim = c(0, max(vif_vector) + 1),
ylab = "Valor del VIF")
# Líneas de referencia
abline(h = 5, col = "orange", lwd = 2, lty = 2)
abline(h = 10, col = "red", lwd = 2, lty = 2)
# 1. Cargar la librería necesaria
library(performance)
library(see) # Esta librería es la que permite que el gráfico se vea así de bien## Warning: package 'see' was built under R version 4.5.3# 2. Generar el análisis de colinealidad del modelo que ya estimaste
colinealidad <- check_collinearity(modelo_estimado)
# 3. Graficar
plot(colinealidad)
Cálculo de los VIF’s usando “car”
library(car)
VIFs_car<-vif(modelo_estimado)
print(VIFs_car)## lotsize sqrft bdrms
## 1.037211 1.418654 1.396663Cálculo de los VIF’s usando “mctest”
library(mctest)
mc.plot(mod = modelo_estimado,vif = 2)