Datos del modelo ejemplo

library(haven)
hprice1<-read_dta("C:/Users/Jacqueline Vanessa/Desktop/UES/Ciclo I - 2022/EMA118/TAREAS/UNIDAD 2/Prueba clase - sesion sincronica/hprice1.dta")
head(hprice1,n=6)

## # A tibble: 6 x 10
##   price assess bdrms lotsize sqrft colonial lprice lassess llotsize lsqrft
##   <dbl>  <dbl> <dbl>   <dbl> <dbl>    <dbl>  <dbl>   <dbl>    <dbl>  <dbl>
## 1  300    349.     4    6126  2438        1   5.70    5.86     8.72   7.80
## 2  370    352.     3    9903  2076        1   5.91    5.86     9.20   7.64
## 3  191    218.     3    5200  1374        0   5.25    5.38     8.56   7.23
## 4  195    232.     3    4600  1448        1   5.27    5.45     8.43   7.28
## 5  373    319.     4    6095  2514        1   5.92    5.77     8.72   7.83
## 6  466.   414.     5    8566  2754        1   6.14    6.03     9.06   7.92

Modelo estimado

library(stargazer)
modelo_estimado<-lm(price~assess+bdrms+lotsize+colonial+llotsize,data = hprice1)
stargazer(modelo_estimado,type = "html",title = "modelo estimado")

**modelo estimado**

	Dependent variable:

	price

assess	0.940^***
	(0.072)

bdrms	8.620
	(6.791)

lotsize	0.001
	(0.001)

colonial	10.031
	(10.580)

llotsize	-13.357
	(17.813)

Constant	68.090
	(146.133)


Observations	88
R²	0.832
Adjusted R²	0.822
Residual Std. Error	43.364 (df = 82)
F Statistic	81.224^*** (df = 5; 82)

Note:	p<0.1; p<0.05; p<0.01

Indice de Condición

El número de condición mide la sensibilidad de las estimaciones mínimo cuadráticas ante pequeños cambios en los datos.

El número de condición,κ(x) , es igual a la raíz cuadrada de la razón entre la raíz característica más grande (λmax) y la raíz característica más pequeña (λmin) de la matriz XtX, normalizada.

Interpretación

* Sí κ(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

Cálculo “manual”

Se basa la matriz XtX:

library(stargazer)
X_mat<-model.matrix(modelo_estimado)
stargazer(head(X_mat,n=6),type="html")


	(Intercept)	assess	bdrms	lotsize	colonial	llotsize

1	1	349.100	4	6,126	1	8.720
2	1	351.500	3	9,903	1	9.201
3	1	217.700	3	5,200	0	8.556
4	1	231.800	3	4,600	1	8.434
5	1	319.100	4	6,095	1	8.715
6	1	414.500	5	8,566	1	9.056

XX_matrix<-t(X_mat)%*%X_mat
stargazer(XX_matrix,type = "text")

## 
## ============================================================================================
##             (Intercept)     assess         bdrms        lotsize      colonial    llotsize   
## --------------------------------------------------------------------------------------------
## (Intercept)     88        27,784.800        314         793,748         61        783.649   
## assess      27,784.800   9,563,053.000  102,507.500 278,300,049.000 19,578.900  250,005.600 
## bdrms           314       102,507.500      1,182       2,933,767       228       2,802.953  
## lotsize       793,748   278,300,049.000  2,933,767  16,165,159,010   555,967   7,457,452.000
## colonial        61        19,578.900        228         555,967         61        544.060   
## llotsize      783.649     250,005.600    2,802.953   7,457,452.000   544.060     7,004.230  
## --------------------------------------------------------------------------------------------

La matriz XtX debe ser normalizada, para evitar sesgo de la escala de las variables.

library(stargazer)
options(scipen = 999)
Sn<-solve(diag(sqrt(diag(XX_matrix))))
stargazer(Sn,type = "text")

## 
## ======================================
## 0.107   0      0      0      0     0  
## 0     0.0003   0      0      0     0  
## 0       0    0.029    0      0     0  
## 0       0      0   0.00001   0     0  
## 0       0      0      0    0.128   0  
## 0       0      0      0      0   0.012
## --------------------------------------

XtX normalizada:

library(stargazer)
XX_norm<-(Sn%*%XX_matrix)%*%Sn
stargazer(XX_norm,type = "text",digits = 4)

## 
## =========================================
## 1      0.9578 0.9736 0.6655 0.8326 0.9982
## 0.9578   1    0.9642 0.7078 0.8106 0.9660
## 0.9736 0.9642   1    0.6712 0.8491 0.9742
## 0.6655 0.7078 0.6712   1    0.5599 0.7008
## 0.8326 0.8106 0.8491 0.5599   1    0.8323
## 0.9982 0.9660 0.9742 0.7008 0.8323   1   
## -----------------------------------------

Autovalores de XtX Normalizada:

library(stargazer)
#autovalores
lambdas<-eigen(XX_norm,symmetric = TRUE)
stargazer(lambdas$values,type = "text")

## 
## ====================================
## 5.193 0.492 0.237 0.049 0.028 0.0005
## ------------------------------------

Cálculo de k(x):

K<-sqrt(max(lambdas$values)/min(lambdas$values))
print(K)

## [1] 106.4903

Como κ(x)>30 se considera que la multicolinealidad es severa.

Cálculo del Indice de Condición usando librería “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.1420         0
## Farrar Chi-Square:       164.9525         1
## Red Indicator:             0.3832         0
## Sum of Lambda Inverse:    12.3289         0
## Theil's Method:           -0.8940         0
## Condition Number:        106.4903         1
## 
## 1 --> COLLINEARITY is detected by the test 
## 0 --> COLLINEARITY is not detected by the test

Cálculo del Indice de Condición usando librería “olsrr”

library(olsrr)
ols_eigen_cindex(model = modelo_estimado)

##     Eigenvalue Condition Index     intercept       assess       bdrms
## 1 5.1933823834        1.000000 0.00003563203 0.0013514449 0.001229095
## 2 0.4923113019        3.247919 0.00004447199 0.0001049883 0.001469387
## 3 0.2365049515        4.686030 0.00028607748 0.0112620908 0.004101462
## 4 0.0490596083       10.288762 0.00483247587 0.4576465408 0.014764598
## 5 0.0282837924       13.550531 0.00119864705 0.2138228412 0.907363692
## 6 0.0004579625      106.490334 0.99360269558 0.3158120940 0.071071766
##       lotsize    colonial      llotsize
## 1 0.003714445 0.007882944 0.00003054565
## 2 0.296269130 0.061464755 0.00001170881
## 3 0.032237901 0.854490636 0.00020426106
## 4 0.010083947 0.002127499 0.00265275719
## 5 0.009254094 0.069640172 0.00167335230
## 6 0.648440482 0.004393993 0.99542737500

Prueba de Farrar-Glaubar

(implementación de la prueba de Bartlett) Esta prueba identifica si a nivel poblacional, los regresores del modelo presentan independencia estadistica (son ortogonales), a través de la matriz de correlación muestral R, y se verifica si a nivel poblacional dicha matriz de correlación corresponde a una matriz identidad.

- Si no se rechaza H0, no hay evidencia de multicolinealidad, caso contrario - Si se rechaza H0 hay evidencia de multicolinealidad.

Cálculo “manual”

Calculo de |R|

library(stargazer)
Zn<-scale(X_mat[,-1])
stargazer(head(Zn,n=6),type = "html")


	assess	bdrms	lotsize	colonial	llotsize

1	0.350	0.513	-0.284	0.662	-0.340
2	0.375	-0.675	0.087	0.662	0.543
3	-1.029	-0.675	-0.375	-1.495	-0.641
4	-0.881	-0.675	-0.434	0.662	-0.866
5	0.035	0.513	-0.287	0.662	-0.349
6	1.036	1.702	-0.045	0.662	0.277

* Calcular la matriz R:

library(stargazer)
n<-nrow(Zn)
R<-(t(Zn)%*%Zn)*(1/(n-1))
#También se puede calcular R a través de cor(X_mat[,-1])
stargazer(R,type = "html",digits = 4)


	assess	bdrms	lotsize	colonial	llotsize

assess	1	0.4825	0.3281	0.0829	0.5717
bdrms	0.4825	1	0.1363	0.3046	0.1695
lotsize	0.3281	0.1363	1	0.0140	0.8079
colonial	0.0829	0.3046	0.0140	1	0.0386
llotsize	0.5717	0.1695	0.8079	0.0386	1

* Calcular |R|:

determinante_R<-det(R)
print(determinante_R)

## [1] 0.1419755

* Aplicando la prueba de Farrer Glaubar (Bartlett):

Estadistico χ2FG

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] 164.9525

* Valor Critico:

gl<-m*(m-1)/2
VC<-qchisq(p = 0.95,df = gl)
print(VC)

## [1] 18.30704

* Regla de desición:

Como χ2FG≥V.C. se rechaza H0, por lo tanto hay evidencia de colinealidad en los regresores.

Cálculo de FG usando “mactest”

library(mctest)
mctest::omcdiag(mod = modelo_estimado)

## 
## Call:
## mctest::omcdiag(mod = modelo_estimado)
## 
## 
## Overall Multicollinearity Diagnostics
## 
##                        MC Results detection
## Determinant |X'X|:         0.1420         0
## Farrar Chi-Square:       164.9525         1
## Red Indicator:             0.3832         0
## Sum of Lambda Inverse:    12.3289         0
## Theil's Method:           -0.8940         0
## Condition Number:        106.4903         1
## 
## 1 --> COLLINEARITY is detected by the test 
## 0 --> COLLINEARITY is not detected by the test

Cálculo de FG usando la “psych”

library(psych)
FG_test<-cortest.bartlett(X_mat[,-1])
print(FG_test)

## $chisq
## [1] 164.9525
## 
## $p.value
## [1] 0.000000000000000000000000000003072151
## 
## $df
## [1] 10

Factores Inflacionarios de la Varianza (FIV)

Los denominados, variance inflation factor (VIF), por sus siglas en inglés, determinan el tamaño relativo de la varianza del j-ésimo parámetro estimado, respecto a la varianza esperada del estimador en ausencia de colinealidad.

Rererencia 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  10

Cálculo manual

Matriz de Correlación de los regresores del modelo (Como se obtuvo con anterioridad):

print(R)

##              assess     bdrms    lotsize   colonial  llotsize
## assess   1.00000000 0.4824739 0.32814633 0.08293582 0.5716654
## bdrms    0.48247394 1.0000000 0.13632563 0.30457549 0.1694902
## lotsize  0.32814633 0.1363256 1.00000000 0.01401865 0.8078552
## colonial 0.08293582 0.3045755 0.01401865 1.00000000 0.0386421
## llotsize 0.57166539 0.1694902 0.80785523 0.03864210 1.0000000

Inversa de la matriz de correlación R−1:

inversa_R<-solve(R)
print(inversa_R)

##              assess      bdrms    lotsize   colonial   llotsize
## assess    2.1535576 -0.9010888  0.8347216  0.1520968 -1.7586001
## bdrms    -0.9010888  1.5104833 -0.3640744 -0.4021983  0.5687703
## lotsize   0.8347216 -0.3640744  3.2049651  0.1130042 -3.0089889
## colonial  0.1520968 -0.4021983  0.1130042  1.1142184 -0.1531266
## llotsize -1.7586001  0.5687703 -3.0089889 -0.1531266  4.3456744

VIF’s para el modelo estimado:

VIFs<-diag(inversa_R)
print(VIFs)

##   assess    bdrms  lotsize colonial llotsize 
## 2.153558 1.510483 3.204965 1.114218 4.345674

Cálculo de los VIF’s usando “performance”

library(performance)
VIFs<-multicollinearity(x = modelo_estimado,verbose = FALSE)
VIFs

## # Check for Multicollinearity
## 
## Low Correlation
## 
##      Term  VIF Increased SE Tolerance
##    assess 2.15         1.47      0.46
##     bdrms 1.51         1.23      0.66
##   lotsize 3.20         1.79      0.31
##  colonial 1.11         1.06      0.90
##  llotsize 4.35         2.08      0.23

plot(VIFs)

Cálculo de los VIF’s usando “car”

library(car)
VIFs_car<-vif(modelo_estimado)
print(VIFs_car)

##   assess    bdrms  lotsize colonial llotsize 
## 2.153558 1.510483 3.204965 1.114218 4.345674

Cálculo de los VIF’s usando “mactest”

library(mctest)
mc.plot(mod = modelo_estimado,vif = 2)

Práctica de clase - multicolinealidad

Jacqueline Vanessa Claros Ruiz

2022-04-28

Datos del modelo ejemplo

Modelo estimado

Indice de Condición

Interpretación

Cálculo “manual”

Cálculo del Indice de Condición usando librería “mctest”

Cálculo del Indice de Condición usando librería “olsrr”

Prueba de Farrar-Glaubar

Cálculo “manual”

* Calcular la matriz R:

* Calcular |R|:

* Aplicando la prueba de Farrer Glaubar (Bartlett):

* Valor Critico:

* Regla de desición:

Cálculo de FG usando “mactest”

Cálculo de FG usando la “psych”

Factores Inflacionarios de la Varianza (FIV)

Rererencia entre R2j

Cálculo manual

Cálculo de los VIF’s usando “performance”

Cálculo de los VIF’s usando “car”

Cálculo de los VIF’s usando “mactest”