Ejercicio de Multicolinealidad
VANESSA IVETH LOPEZ GONZALEZ
2024-05-13
Ciudad Universitaria, 13 de Mayo de 2024
Utilizando los datos del dataframe hprice1: disponible en el paquete wooldridge use el siguiente código para generar el dataframe:
library(wooldridge)
data(hprice1)
head(force(hprice1), n=5) ## 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.8296301. Estime el siguiente modelo:
\[price=\hat{\alpha}+\hat{\alpha}_1\times\text{lotsize}+\hat{\alpha}_2 \times\text{sqrft} +\hat{\alpha}_3 \times\text{bdrms}+\epsilon\]
library(stargazer)
modelo_precio<-lm(formula = price~lotsize+sqrft+bdrms,data = hprice1)
stargazer(modelo_precio,title = "Modelo Precio",type = "text",digits = 4)##
## Modelo Precio
## ===============================================
## Dependent variable:
## ---------------------------
## price
## -----------------------------------------------
## lotsize 0.0021***
## (0.0006)
##
## sqrft 0.1228***
## (0.0132)
##
## bdrms 13.8525
## (9.0101)
##
## Constant -21.7703
## (29.4750)
##
## -----------------------------------------------
## Observations 88
## R2 0.6724
## Adjusted R2 0.6607
## Residual Std. Error 59.8335 (df = 84)
## F Statistic 57.4602*** (df = 3; 84)
## ===============================================
## Note: *p<0.1; **p<0.05; ***p<0.012. Verifique si hay evidencia de la independencia de los regresores (no colinealidad), a través de:
A) Indice de condición y prueba de FG, presente sus resultados de manera tabular en ambos casos y para la prueba de FG presente también sus resultados de forma gráfica usando la librería fastGraph
Indice de condición de forma manual
Basado en la matriz XtX :
Matriz X
matriz_x <- model.matrix(modelo_precio)
stargazer(head(matriz_x,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
## ---------------------------------Matriz XX
matriz_xx<-t(matriz_x)%*%matriz_x
stargazer(head(matriz_xx),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
## --------------------------------------------------------------Cálculo de la matriz de normalización:
options(scipen = 999)
sn<-solve(diag(sqrt(diag(matriz_xx))))
stargazer(head(sn),type = "text")##
## ==========================
## 0.107 0 0 0
## 0 0.00001 0 0
## 0 0 0.0001 0
## 0 0 0 0.029
## --------------------------XtX Normalizada
matrizxx_norm<-(sn%*%matriz_xx)%*%sn
stargazer(matrizxx_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 XtX Normalizada
lambdas<-eigen(matrizxx_norm, symmetric = TRUE)
stargazer(lambdas$values,type = "text")##
## =======================
## 3.482 0.455 0.039 0.025
## -----------------------Cálculo de κ(x)=√λmax/λmin
k<-sqrt(max(lambdas$values)/min(lambdas$values))
print(k)## [1] 11.86778Como k(x)< 20, la multicolinealidad se considera que no es un problema, ya que es leve
Indice de condición usando la librería “mctest”
library(mctest)
mctest(mod=modelo_precio)##
## 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 testIndice de condición usando la librería “olsrr”
library(olsrr)
ols_eigen_cindex(model=modelo_precio)## 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.975025872Prueba de Farrar-Glaubar de forma manual
Calculo de |R|
zn<-scale(matriz_x[,-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
n<-nrow(zn)
matriz_r<-(t(zn)%*%zn)*(1/(n-1))
stargazer(matriz_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
## -----------------------------R a través de “cor”
matrizr_cor<-cor(matriz_x[,-1])
stargazer(matrizr_cor,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
## -----------------------------Determinante |R|
determinante_r<-det(matriz_r)
print(determinante_r)## [1] 0.6917931Aplicando la prueba de Farrer Glaubar (Bartlett)
Estadístico χ2FG
m<-ncol(matriz_x[,-1])
n<-nrow(matriz_x[,-1])
chi_FG<--(n-1-(2*m+5)/6)*log(determinante_r)
print(chi_FG)## [1] 31.38122Valor crítico
gl<-m*(m-1)/2
vc<-qchisq(p = 0.95,df = gl)
print(vc)## [1] 7.814728Regla de desición:
Como el Estadistico ≥ Valor Critico se rechaza H0, por lo tanto hay evidencia de colinealidad en los regresores
Cálculo de FG usando “mctest”
library(mctest)
mctest::omcdiag(mod = modelo_precio)##
## Call:
## mctest::omcdiag(mod = modelo_precio)
##
##
## 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 testCálculo de FG usando la “psych”
library(psych)
FG_test<-cortest.bartlett(matriz_x[,-1])
print(FG_test)## $chisq
## [1] 31.38122
##
## $p.value
## [1] 0.0000007065806
##
## $df
## [1] 3Gráfica
chi_FG<--(n-1-(2*m+5)/6)*log(determinante_r)
gl<-m*(m-1)/2
vc<-qchisq(p = 0.95,df = gl)
alpha_sig <- 0.05
library(fastGraph)
shadeDist(chi_FG,ddist = "dchisq",parm1 = gl,lower.tail = FALSE,xmin = 0,xlab = "Valor de chi-cuadrado",main = "Prueba de Bartlett")
abline(v = vc,col = "red",lty = 2)
axis(1,at = vc,labels = paste("vc:",round(vc, 2)),col.axis = "black",las = 1)
if (chi_FG > vc) {
text(x = vc + 0, y = 0.22, labels = "Rechazar H0", col = "blue", cex = 0.8)
} else {
text(x = vc + 0, y = 0.22, labels = "No rechazar H0", col = "blue", cex = 0.8)
}
text(vc, 0, expression(alpha == 0.05), pos = 4, col = "black", cex = 0.8)
B) Factores inflacionarios de la varianza, presente sus resultados de forma tabular y de forma gráfica.
Referencia entre R2j:
library(dplyr)
r.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(matriz_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:
matrizr_inversa<-solve(matriz_r)
print(matrizr_inversa)## 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(matrizr_inversa)
print(VIFs)## lotsize sqrft bdrms
## 1.037211 1.418654 1.396663Cálculo de los VIF’s usando “performance”
library(performance)
VIFs<-multicollinearity(x = modelo_precio,verbose = FALSE)
print(VIFs)## # Check for Multicollinearity
##
## Low Correlation
##
## Term VIF VIF 95% CI Increased SE 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]library(see)
plot(VIFs)## Variable `Component` is not in your data frame :/
Cálculo de los VIF’s usando “car”
library(car)## Loading required package: carData##
## Attaching package: 'car'## The following object is masked from 'package:dplyr':
##
## recode## The following object is masked from 'package:psych':
##
## logitVIFs_car<-vif(modelo_precio)
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_precio ,vif = 2)