##Data frame
library(wooldridge)
## Warning: package 'wooldridge' was built under R version 4.5.3
data(hprice1)
head(force(hprice1),n=6) #mostrar las primeras 5 observaciones
## price assess bdrms lotsize sqrft colonial lprice lassess llotsize
## 1 300.000 349.1 4 6126 2438 1 5.703783 5.855359 8.720297
## 2 370.000 351.5 3 9903 2076 1 5.913503 5.862210 9.200593
## 3 191.000 217.7 3 5200 1374 0 5.252274 5.383118 8.556414
## 4 195.000 231.8 3 4600 1448 1 5.273000 5.445875 8.433811
## 5 373.000 319.1 4 6095 2514 1 5.921578 5.765504 8.715224
## 6 466.275 414.5 5 8566 2754 1 6.144775 6.027073 9.055556
## lsqrft
## 1 7.798934
## 2 7.638198
## 3 7.225482
## 4 7.277938
## 5 7.829630
## 6 7.920810
#Modelo estimado
library(stargazer)
##
## Please cite as:
## Hlavac, Marek (2022). stargazer: Well-Formatted Regression and Summary Statistics Tables.
## R package version 5.2.3. https://CRAN.R-project.org/package=stargazer
modelo_estimado<-lm(price~assess+bdrms+lotsize+colonial+llotsize,data = hprice1)
stargazer(modelo_estimado,type = "text",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
## R2 0.832
## Adjusted R2 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 ##Calculo manual de la matriz X^(t)*X
library(stargazer)
X_mat<-model.matrix(modelo_estimado)
stargazer(head(X_mat,n=6),type="text")
##
## =====================================================
## (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
## -----------------------------------------------------
##Calculo de la matriz XX
XX_matrix<-t(X_mat)%*%X_mat
stargazer(head(XX_matrix,n=6), 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
## --------------------------------------------------------------------------------------------
##Calculo de la matriz de normalización
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
## --------------------------------------
##X(t)*X 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 X(t)*X 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
## ------------------------------------
##Calculo del K(X)
K<-sqrt(max(lambdas$values)/min(lambdas$values))
print(K)
## [1] 106.4903
#Conclusión: Como κ(x)>30 se considera que la multicolinealidad es severa
#CALCULO DEL INDICE DE CONDICION USANDO LA LIBRERIA “mctest”
library(mctest)
## Warning: package 'mctest' was built under R version 4.5.2
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)
## Warning: package 'olsrr' was built under R version 4.5.3
##
## Adjuntando el paquete: 'olsrr'
## The following object is masked from 'package:wooldridge':
##
## cement
## The following object is masked from 'package:datasets':
##
## rivers
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
#Calculo de la prueba Farrar Gable de forma manual ##Calculo de |R|
library(stargazer)
Zn<-scale(X_mat[,-1])# Zn es una matriz
stargazer(head(Zn,n=6),type = "text")
##
## =========================================
## 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)# para saber la cantidad de filas de Zn
R<-(t(Zn)%*%Zn)*(1/(n-1))# formula de R
#También se puede calcular R a través de cor(X_mat[,-1])
stargazer(R,type = "text",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
#PARTE ESTADISTICA #APLICANDO la prueba de FARRER GLAUBAR(Bartlet) #CALCULO DEL ESTADISTICO X^2FG
m<-ncol(X_mat[,-1]) #obtengo el numero de variable.
n<-nrow(X_mat[,-1]) # numero de casos, dado por el numero de filas.
chi_FG<--(n-1-(2*m+5)/6)*log(determinante_R)#formula del estadistico.
print(chi_FG)
## [1] 164.9525
#CALCULO DEL VALOR CRITICO
gl<-m*(m-1)/2
VC<-qchisq(p = 0.95,df = gl)# queremos el cuantil,rchisq sirve para calcular valores criticos.
print(VC)
## [1] 18.30704
##REGLA DE DESCISIÓN ##Dado que X^(2)FG>=v.c entonces se rechaza la hipótesis nula. Por tanto hay evidencia de colinealidad en los regresores.
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
##Calculo de FG usando la “psych”
library(psych)
## Warning: package 'psych' was built under R version 4.5.3
FG_test<-cortest.bartlett(X_mat[,-1])
## R was not square, finding R from data
print(FG_test)
## $chisq
## [1] 164.9525
##
## $p.value
## [1] 0.000000000000000000000000000003072151
##
## $df
## [1] 10
#Factores Inflacionarios de la Varianza(FIV) ##CALCULO MANUAL ###inversa de R
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
###CALCULO DEL VIF
VIFs<-diag(inversa_R)
print(VIFs)
## assess bdrms lotsize colonial llotsize
## 2.153558 1.510483 3.204965 1.114218 4.345674
#SI en rango es de 2 entonces, lotsize y llotsize, tiene colinealidad.
##CALCULO DEL VIFs UTILIZANADO LA LIBRERIA PERFORMANCE
library(performance)
## Warning: package 'performance' was built under R version 4.5.3
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
## assess 2.15 [1.67, 2.98] 1.47 0.46 [0.34, 0.60]
## bdrms 1.51 [1.24, 2.08] 1.23 0.66 [0.48, 0.81]
## lotsize 3.20 [2.39, 4.51] 1.79 0.31 [0.22, 0.42]
## colonial 1.11 [1.01, 1.87] 1.06 0.90 [0.53, 0.99]
## llotsize 4.35 [3.17, 6.17] 2.08 0.23 [0.16, 0.32]
##CALCULO DEL VIFs UTILIZANDO LA LIBRERIA CAR
library(car)
## Cargando paquete requerido: carData
##
## Adjuntando el paquete: 'car'
## The following object is masked from 'package:psych':
##
## logit
VIFs_car<-vif(modelo_estimado)
print(VIFs_car)
## assess bdrms lotsize colonial llotsize
## 2.153558 1.510483 3.204965 1.114218 4.345674
##calculo DE LOS VIFS USANDO MCTEST
library(mctest)
mc.plot(mod = modelo_estimado,vif = 2)#ESTO ME modifica la frontera