Profesor:Dr. José Manuel MAGALLANES REYES, Ph.D

  • Profesor Principal del Departamento de Ciencias Sociales, Sección de Ciencia Política y Gobierno.

  • Oficina 223 - Edificio CISEPA / ECONOMIA / CCSS

  • Telefono: (51) 1 - 6262000 anexo 4302

  • Correo Electrónico: jmagallanes@pucp.edu.pe


Sesión 2: Regresión Gaussiana


DOI

1 Introducción

La regresión es una técnica donde hay que definir una variable dependiente y una o más independientes. Las independientes pueden tener rol predictor, dependiendo del diseño de investigación, pero cumple siempre un rol asociativo; así, la regresión quiere informar cuánto la variabilidad de la variable (independiente) puede asociarse a la variabilidad de la dependiente, controlando el efecto de otras variables, de ahí que es una técnica para probar hipótesis direccionales o asimétricas (las correlaciones tiene hipótesis no direccionales o simétricas).

La regresión Gaussiana busca proponer un modelo relacional entre variable, y es aplicable sólo cuando la variable dependiente (la Y) es numérica, continua, y no acotada.

Caso de Estudio: Pavimentando con votos

Profesores de la Universidad de los Andes (Mejia Guinand, Botero, and Rodriguez Raga 2008) decidieron estudiar cómo la distribución de fondos públicos fue afectada por factores políticos durante el primer periodo del Presidente Uribe (2002-2006). Las hipótesis que se plantean son:

  • H1: la asignación presupuestal en infraestructura vial en Colombia responde a los criterios técnicos y económicos determinados en el Plan Nacional de Desarrollo y otros documentos de carácter técnico elaborados por el gobierno.

  • H2: la asignación presupuestal en infraestructura vial en Colombia responde a negociaciones bilaterales entre el ejecutivo y el legislativo basadas en necesidades políticas y electorales de corto plazo.

  • H3: la asignación presupuestal en infraestructura vial en Colombia responde al esfuerzo del gobierno por fortalecer su base social de apoyo local a través de los Consejos Comunales de Gobierno.

Para ello, organizaron estos datos:

1.1 Preparación de los datos

Una vez se tienen en claro las hipótesis, podemos pasar a organizar los datos. Veamos los datos del artículo en cuestíon:

rm(list = ls()) # limpiar el working environment

#carga de data
linkToData='https://docs.google.com/spreadsheets/d/e/2PACX-1vSRpCC8gKIxMxpK0wjgLcl-GQWdw6sAeB16Sixkq6kZXeTfMS8_n70twbEbQ2tMcGp8tm-6x8qf8ieo/pub?gid=234740779&single=true&output=csv'
pavi=read.csv(linkToData)
str(pavi)
## 'data.frame':    1096 obs. of  9 variables:
##  $ municipio       : chr  "m1" "m2" "m3" "m4" ...
##  $ apropiaciondolar: num  102.17 4.19 1.59 0 0 ...
##  $ consejocomunal  : int  0 0 0 0 0 0 0 0 0 0 ...
##  $ ejecucion       : int  0 0 0 0 0 0 0 0 0 0 ...
##  $ uribista        : int  0 1 1 1 1 1 1 1 1 NA ...
##  $ pctopo          : num  14.82 14.51 15.08 6.15 47.31 ...
##  $ priorizado      : int  0 1 0 0 0 0 0 0 1 0 ...
##  $ poblacioncienmil: num  20.9155 0.2406 0.0398 0.0558 0.2723 ...
##  $ nbi             : num  12.2 33.8 28.5 33.1 27.1 ...

Recuerda que es mejor usar datos sin valores perdidos en los análisis multivariados 1.

pavi=pavi[complete.cases(pavi),]

Los datos debemos darles el formato adecuado. En este caso hay que cambiar a categóricos algunas variables:

seleccion=c("consejocomunal","ejecucion","uribista","priorizado")
pavi[,seleccion]=lapply(pavi[,seleccion],as.factor)

1.2 Explorando las variables

paviStats=summary(pavi[,-1])
paviStats
##  apropiaciondolar  consejocomunal ejecucion uribista     pctopo      
##  Min.   :  0.000   0:827          0:844     0:325    Min.   : 0.000  
##  1st Qu.:  0.000   1: 49          1: 32     1:551    1st Qu.: 6.285  
##  Median :  0.000                                     Median :19.607  
##  Mean   :  8.926                                     Mean   :27.306  
##  3rd Qu.: 11.436                                     3rd Qu.:44.455  
##  Max.   :132.643                                     Max.   :99.419  
##  priorizado poblacioncienmil        nbi       
##  0:639      Min.   : 0.00629   Min.   : 6.84  
##  1:237      1st Qu.: 0.07460   1st Qu.:27.33  
##             Median : 0.14593   Median :40.26  
##             Mean   : 0.45639   Mean   :41.72  
##             3rd Qu.: 0.27433   3rd Qu.:53.78  
##             Max.   :69.26836   Max.   :97.32

¿Qué puedes decir de estos resultados?

Luego, es importante saber qué relación se revela entre la variable dependiente y las demás. Se puede comenzar con las que son numéricas, es decir, haciendo un análisis de correlación. Eso lo apreciamos en la Figura 1.1.

library(ggcorrplot)
colNums=names(pavi)[c(2,6,8,9)]
numXs=pavi[,colNums]
ggcorrplot(cor(numXs),lab = T,show.diag = F)
Correlación entre la VD y sus VIs

Figure 1.1: Correlación entre la VD y sus VIs

Por otro lado, podemos revisar la relación con las categóricas. Rapidamente podemos ver la diferencia estadística, parametrica y no paramétrica, por grupos en relación a nuestra dependiente en la Tabla 1.1.

library(magrittr)
library(kableExtra)

colCats=setdiff(names(pavi), colNums)[-1]
diffPara=c()
diffNoPara=c()

for (col in colCats){
    diffPara=c(diffPara,t.test(pavi[,"apropiaciondolar"]~pavi[,col])['p.value']<=0.05)
    diffNoPara=c(diffNoPara,wilcox.test(pavi[,"apropiaciondolar"]~pavi[,col])['p.value']<=0.05)
}
data.frame(cbind(colCats,diffPara,diffNoPara),
           row.names = 1:length(colCats))%>%
           kable(caption = "Diferencia de 'VD:Apropiacion' por Grupo")%>%
            kableExtra::kable_styling(full_width = FALSE)
Table 1.1: Table 1.2: Diferencia de ‘VD:Apropiacion’ por Grupo
colCats diffPara diffNoPara
consejocomunal TRUE TRUE
ejecucion FALSE FALSE
uribista TRUE FALSE
priorizado FALSE TRUE

Lo mostrado en la Tabla 1.1 puede visualizarse en la Figura 1.2

par(mfrow = c(2, 2))  

for (col in colCats) { 
  boxplot(pavi$apropiaciondolar~pavi[,col],
       main = col,xlab="",ylab = "")
}
Apropiacion por Grupos

Figure 1.2: Apropiacion por Grupos

2 Análisis de Regresión

Veamos una hipotesis:

El Beneficio recibido en un municipio ha sido afectado por el porcentaje de votos recibidos por los candidatos de la oposición a Uribe a la camara de representantes, controlando por tamaño de población.

# hipotesis en R
modelo1=formula(apropiaciondolar~pctopo+poblacioncienmil)

Si deseamos probar esa hipotesis en R, podemos hacer lo siguiente:

reg1=lm(modelo1,data=pavi)
summary(reg1)
## 
## Call:
## lm(formula = modelo1, data = pavi)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -53.809  -8.671  -7.448   2.602  90.426 
## 
## Coefficients:
##                  Estimate Std. Error t value Pr(>|t|)    
## (Intercept)       8.80319    0.79976  11.007   <2e-16 ***
## pctopo           -0.03146    0.02150  -1.463    0.144    
## poblacioncienmil  2.15125    0.20073  10.717   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 15.85 on 873 degrees of freedom
## Multiple R-squared:  0.1181, Adjusted R-squared:  0.1161 
## F-statistic: 58.45 on 2 and 873 DF,  p-value: < 2.2e-16

El resultado anterior puede ser presentado de mejor manera usando modelsummary package, como se puede ver en la Tabla 2.1 (el error típico está entre parentesis).

library(modelsummary)
model1=list('apropiacion (I)'=reg1)
modelsummary(model1, title = "Regresion: modelo 1",
             stars = TRUE,
             output = "kableExtra")
Table 2.1: Regresion: modelo 1
&nbsp;apropiacion (I)
(Intercept) 8.803***
(0.800)
pctopo −0.031
(0.021)
poblacioncienmil 2.151***
(0.201)
Num.Obs. 876
R2 0.118
R2 Adj. 0.116
AIC 7332.6
BIC 7351.7
Log.Lik. −3662.298
F 58.447
RMSE 15.83
+ p < 0.1, * p < 0.05, ** p < 0.01, *** p < 0.001


Al probar esta hipótesis vemos…

  1. que pctopo tiene signo negativo (relación inversa con la VD),
  2. que la magnitud de ese efecto es -0.031, lo que indica cuanto varía apropiaciondolar en promedio cuando pctopo se incremente en una unidad, controlando por poblacioncienmil.
  3. que pctopo NO tiene efecto significativo.

Esto es información suficiente para representar esa relación con la Ecuación (2.1).

\[\begin{equation} apropiaciondolar = 8.8031897 + -0.0314583 \cdot pctopo + 2.1512469\cdot poblacioncienmil + \epsilon \tag{2.1} \end{equation}\]

Justamente el R cuadrado ajustado (0.1180881) nos brinda un porcentaje (multiplicalo por 100) que da una pista de nuestra cercanía a una situación perfecta (cuando vale 1).

¿Y sí queremos ver el efecto de consejo comunal (consejocomunal)?

# modelo 1 mas una nueva independiente
modelo2=formula(apropiaciondolar~pctopo+consejocomunal+poblacioncienmil)

Esta nueva hipótesis desea evaluar si la visita de Uribe a un Consejo Comunal influye en la asignación de presupuesto. El resultado de este proceso lo vemos a continuación.

reg2=lm(modelo2,data=pavi)
summary(reg2)
## 
## Call:
## lm(formula = modelo2, data = pavi)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -62.378  -7.933  -6.834   1.573  91.247 
## 
## Coefficients:
##                  Estimate Std. Error t value Pr(>|t|)    
## (Intercept)       8.04119    0.79129  10.162  < 2e-16 ***
## pctopo           -0.02987    0.02103  -1.420    0.156    
## consejocomunal1  14.79601    2.32127   6.374 2.98e-10 ***
## poblacioncienmil  1.91236    0.19987   9.568  < 2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 15.51 on 872 degrees of freedom
## Multiple R-squared:  0.1573, Adjusted R-squared:  0.1545 
## F-statistic: 54.28 on 3 and 872 DF,  p-value: < 2.2e-16

Visualizando de manera alternativa en la Tabla 2.2:

model2=list('apropiacion (II)'=reg2)
modelsummary(model2, title = "Regresion: modelo 2",
             stars = TRUE,
             output = "kableExtra")
Table 2.2: Regresion: modelo 2
&nbsp;apropiacion (II)
(Intercept) 8.041***
(0.791)
pctopo −0.030
(0.021)
consejocomunal1 14.796***
(2.321)
poblacioncienmil 1.912***
(0.200)
Num.Obs. 876
R2 0.157
R2 Adj. 0.154
AIC 7294.7
BIC 7318.6
Log.Lik. −3642.351
F 54.277
RMSE 15.47
+ p < 0.1, * p < 0.05, ** p < 0.01, *** p < 0.001

Al probar esta hipótesis vemos que…

  • pctopo tiene signo negativo; NO tiene efecto significativo; y la magnitud de ese efecto es -0.03, lo que indica cuanto varíaría apropiaciondolar en promedio cuando pctopo se incremente en una unidad, controlando por las demás variables.
  • consejocomunal SÍ tiene efecto significativo al 0.001; ese efecto es directo, pues el coeficiente calculado es positivo; y la magnitud de ese efecto es 14.796, lo que indica cuanto varía apropiaciondolar en promedio cuando consejocomunal es 1 y no 0, también controlando por las demás variables.

Nótese la lectura del efecto cuando la variable independiente es categórica (o factor). Primero, nota que R indica el valor 1 con el nombre de la variable: eso indica que el valor 0 (cuando el consejo comunal de ese municipio NO fue visitado) es la categoría de referencia; es decir, el coeficiente nos indica cuanto se modifica el valor promedio de la dependiente cuando se pasa de 0 a 1. Si la variable independiente es politómica (no ordinal), aparecerá cada categoría menos la de referencia, y el efecto siempre se debe interpretar como el efecto de la variable mostrada versus la de referencia. Con esta información podemos proponer la Ecuación (2.2)

\[\begin{equation} apropiaciondolar = 8.0411893 + -0.0298691 \cdot pctopo + 14.7960149 \cdot consejocomunal + 1.9123628\cdot poblacioncienmil + \epsilon \tag{2.2} \end{equation}\]

Incluyamos ahora el predictor uribista al reciente modelo (Tabla 2.2), obtenemos estos resultados:

# modelo 2 mas 'uribista'
modelo3=formula(apropiaciondolar~pctopo+consejocomunal+uribista+poblacioncienmil)

reg3=lm(modelo3,data=pavi)
summary(reg3)
## 
## Call:
## lm(formula = modelo3, data = pavi)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -61.004  -7.546  -6.482   1.847  91.958 
## 
## Coefficients:
##                  Estimate Std. Error t value Pr(>|t|)    
## (Intercept)       9.80542    1.11705   8.778  < 2e-16 ***
## pctopo           -0.03753    0.02126  -1.765   0.0779 .  
## consejocomunal1  14.73795    2.31613   6.363 3.19e-10 ***
## uribista1        -2.45126    1.09801  -2.232   0.0258 *  
## poblacioncienmil  1.89052    0.19965   9.469  < 2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 15.47 on 871 degrees of freedom
## Multiple R-squared:  0.1621, Adjusted R-squared:  0.1583 
## F-statistic: 42.14 on 4 and 871 DF,  p-value: < 2.2e-16

De manera más presentable en la Tabla 2.3:

model3=list('apropiacion (III)'=reg3)
modelsummary(model3, title = "Regresion: modelo 3",
             stars = TRUE,
             output = "kableExtra")
Table 2.3: Regresion: modelo 3
&nbsp;apropiacion (III)
(Intercept) 9.805***
(1.117)
pctopo −0.038+
(0.021)
consejocomunal1 14.738***
(2.316)
uribista1 −2.451*
(1.098)
poblacioncienmil 1.891***
(0.200)
Num.Obs. 876
R2 0.162
R2 Adj. 0.158
AIC 7291.7
BIC 7320.4
Log.Lik. −3639.852
F 42.140
RMSE 15.43
+ p < 0.1, * p < 0.05, ** p < 0.01, *** p < 0.001

Al probar esta hipótesis vemos que..

  • pctopo tiene signo negativo; y ahora SÍ tiene efecto significativo; y la magnitud de ese efecto es -0.038, lo que indica cuanto varíaría apropiaciondolar en promedio cuando pctopo se incremente en una unidad, controlando por las demás variables.

  • consejocomunal MANTIENE tiene efecto significativo al 0.001; ese efecto es directo, pues el coeficiente calculado es positivo; y la magnitud de ese efecto es 14.738, lo que indica cuanto varía apropiaciondolar en promedio cuando consejocomunal es 1 y no 0, también controlando por las demás variables.

  • uribista tiene efecto significativo al 0.05; ese efecto es inverso, pues el coeficiente calculado es negativo; y la magnitud de ese efecto es -2.451, lo que indica cuanto varía apropiaciondolar en promedio cuando uribista es 1 y no 0, también controlando por las demás variables.

Con esta información podemos proponer la Ecuación (2.3)

\[\begin{equation} apropiaciondolar = 9.8054173 + -0.0375299 \cdot pctopo + 14.7379501 \cdot consejocomunal + -2.4512573 \cdot uribista+ 1.8905189\cdot poblacioncienmil + \epsilon \tag{2.3} \end{equation}\]

2.1 Estandarización de Coeficientes

Del resultado de la Tabla 2.3 NO podemos directamente decir que consejo comunal tiene más efecto que los demás por el solo hecho que el valor estimado sea mayor a los demás. Para saber cuál tiene más efecto, cuando los predictores tienen, como en este caso unidades diferentes, estandarizamos los datos y volvemos a correr la regresión. Veamos la Tabla 2.4.

modelo3_st=formula(scale(apropiaciondolar)~scale(pctopo)+scale(as.numeric(consejocomunal))+scale(as.numeric(uribista))+scale(poblacioncienmil))

modelo3_st=lm(modelo3_st,data=pavi)

modelo3_st=list('apropiacion (III_st)'=modelo3_st)
modelsummary(modelo3_st, title = "Regresion: modelo 3 con \ncoeficientes estandarizados",
             stars = TRUE,
             output = "kableExtra")
Table 2.4: Regresion: modelo 3 con coeficientes estandarizados
&nbsp;apropiacion (III_st)
(Intercept) 0.000
(0.031)
scale(pctopo) −0.055+
(0.031)
scale(as.numeric(consejocomunal)) 0.201***
(0.032)
scale(as.numeric(uribista)) −0.070*
(0.031)
scale(poblacioncienmil) 0.299***
(0.032)
Num.Obs. 876
R2 0.162
R2 Adj. 0.158
AIC 2342.0
BIC 2370.7
Log.Lik. −1165.004
F 42.140
RMSE 0.91
+ p < 0.1, * p < 0.05, ** p < 0.01, *** p < 0.001

De los resultados de la Tabla 2.4 se despeja que, dejando de lado la variable de control (poblacion por cien mil), la que tiene mayor efecto es consejo comunal, algo que no era evidente en la Tabla 2.3. Esto nos indica cuantas desviaciones estándar varía la variable dependiente (apropiacion dolar) cuando la dependiente varía en una (1) desviación estándar 2. Nota además que hemos tenido que estandarizar variables categóricas, lo cuál NO es posible. Ello no debe preocupar pues esta etapa es sólo para comparar tamaño de efecto entre variables.

Podemos simplificar los pasos si utilizamos funciones del paquete lm.beta. La función que lleva el mismo nombre nos da el resultado rápidamente, como se ve en la Tabla 2.5.

library(lm.beta)
## Warning: package 'lm.beta' was built under R version 4.4.1
model3beta=list('apropiacion (III)'=lm.beta(reg3))
modelsummary(model3beta, title = "Regresion: modelo 3 con \ncoeficientes estandarizados usando lm.beta()",
             stars = TRUE,
             output = "kableExtra")
Table 2.5: Regresion: modelo 3 con coeficientes estandarizados usando lm.beta()
&nbsp;apropiacion (III)
pctopo −0.055
(0.021)
consejocomunal1 0.201*
(2.316)
uribista1 −0.070
(1.098)
poblacioncienmil 0.299
(0.200)
Num.Obs. 876
R2 0.162
R2 Adj. 0.158
AIC 7291.7
BIC 7320.4
Log.Lik. −3639.852
RMSE 15.43
+ p < 0.1, * p < 0.05, ** p < 0.01, *** p < 0.001

Puedes complementar esta información revisando esta nota en Wikipedia (2022).

3 Selección de modelos

Veamos todos los modelos anteriores en la tabla 3.1:

models=list('apropiacion (I)'=reg1,
            'apropiacion (II)'=reg2,
            'apropiacion (III)'=reg3)
modelsummary(models, title = "Resultados de todos los modelos",
             stars = TRUE,
             output = "kableExtra")
Table 3.1: Resultados de todos los modelos
&nbsp;apropiacion (I) &nbsp;apropiacion (II) &nbsp;apropiacion (III)
(Intercept) 8.803*** 8.041*** 9.805***
(0.800) (0.791) (1.117)
pctopo −0.031 −0.030 −0.038+
(0.021) (0.021) (0.021)
poblacioncienmil 2.151*** 1.912*** 1.891***
(0.201) (0.200) (0.200)
consejocomunal1 14.796*** 14.738***
(2.321) (2.316)
uribista1 −2.451*
(1.098)
Num.Obs. 876 876 876
R2 0.118 0.157 0.162
R2 Adj. 0.116 0.154 0.158
AIC 7332.6 7294.7 7291.7
BIC 7351.7 7318.6 7320.4
Log.Lik. −3662.298 −3642.351 −3639.852
F 58.447 54.277 42.140
RMSE 15.83 15.47 15.43
+ p < 0.1, * p < 0.05, ** p < 0.01, *** p < 0.001

Hagamos un cambio a la tabla anterior, para que en vez de los errores típicos, se muestre el intervalo de confianza del coeficiente estimado:

models=list('apropiacion (I)'=reg1,
            'apropiacion (II)'=reg2,
            'apropiacion (III)'=reg3)
modelsummary(models, title = "Resultados de todos los modelos",statistic = "conf.int",
             stars = TRUE,
             output = "kableExtra")
Table 3.2: Resultados de todos los modelos
&nbsp;apropiacion (I) &nbsp;apropiacion (II) &nbsp;apropiacion (III)
(Intercept) 8.803*** 8.041*** 9.805***
[7.234, 10.373] [6.488, 9.594] [7.613, 11.998]
pctopo −0.031 −0.030 −0.038+
[−0.074, 0.011] [−0.071, 0.011] [−0.079, 0.004]
poblacioncienmil 2.151*** 1.912*** 1.891***
[1.757, 2.545] [1.520, 2.305] [1.499, 2.282]
consejocomunal1 14.796*** 14.738***
[10.240, 19.352] [10.192, 19.284]
uribista1 −2.451*
[−4.606, −0.296]
Num.Obs. 876 876 876
R2 0.118 0.157 0.162
R2 Adj. 0.116 0.154 0.158
AIC 7332.6 7294.7 7291.7
BIC 7351.7 7318.6 7320.4
Log.Lik. −3662.298 −3642.351 −3639.852
F 58.447 54.277 42.140
RMSE 15.83 15.47 15.43
+ p < 0.1, * p < 0.05, ** p < 0.01, *** p < 0.001

La tabla 3.2 muestra que los valores no significativos incluyen al cero en el intervalo de confianza (salvo en el caso de pctopo para el modelo 3, con significancia al 0.1). Gráficamente:

library(ggplot2)
library(sjPlot)
## Warning: package 'sjPlot' was built under R version 4.4.1
plot_models(reg1,reg2,reg3,vline.color = "black",m.labels=c("Modelo 1","Modelo 2","Modelo 3"),dot.size = 1,line.size = 0.6)

En todos los modelos, el \(\epsilon\) no tiene coeficiente, representamos su variación usando el error típico de los residuos o residual standard error (RSE). Nótese que éste ha variado de un modelo ha otro tanto en las Tablas 3.1 y 3.2, y en el último modelo se tiene un R2Adj mayor (el menor RSE). Aquí vale la pena preguntarse si esta disminución del error es significativa. Los resultados los vemos en la Tabla 3.3.

library(magrittr)
library(knitr)
tanova=anova(reg1,reg2,reg3)

kable(tanova,
      caption = "Tabla ANOVA para comparar modelos")%>%kableExtra::kable_styling(full_width = FALSE)
Table 3.3: Table 3.4: Tabla ANOVA para comparar modelos
Res.Df RSS Df Sum of Sq F Pr(>F)
873 219454.2 NA NA NA NA
872 209684.3 1 9769.883 40.814973 0.0000000
871 208491.3 1 1192.989 4.983868 0.0258382


La comparación de modelos usando la tabla de análisis de varianza (anova) propone como hipótesis nula que los modelos no difieren (no se ha reducido el error al pasar de un modelo a otro). Cuando la comparación es significativa (vea el Pr(>F)), rechazamos igualdad de modelos: el modelo 2 presenta menos error que el modelo 1, pero el modelo 2 sí presenta menor error al modelo 2, por lo que el modelo 3 debe ser el elegido (hasta aquí).

4 Diagnósticos de la Regresión

Con el apoyo de las computadoras y el software estadístico es relativamente facil calcular una regresión. Sin embargo, hay que analizar los resultados obtenidos para poder rener una mejor conclusión. Revisemos el resultado de la segunda regresion en las siguientes sub secciones.

4.1 Linealidad

Se asume relación lineal entre la variable dependiente (Y) y las independientes y Xs. Para ello analizamos la relación entre los residuos y los valores que predice el modelo de regresión. Ver Figura 4.1.

# linea roja debe tender a horizontal
plot(reg2, 1)
Evaluando Linealidad

Figure 4.1: Evaluando Linealidad

La falta de linearidad provocaría que el modelo no sirva para explicar las mismas variables con datos diferentes en otros estudios.

4.2 Homocedasticidad

Se asume que la dispersión de los errores de la estimación (\(\hat{apropiaciondolar}\)) mantiene una variación homogenea. Si analizamos las raices de los errores estandarizados versus los valores estimados, la distancia de los puntos a la linea de referencia debe ser similar. Ver Figura 4.2.

# linea roja debe tender a horizontal
plot(reg2, 3)
Evaluando Homocedasticidad

Figure 4.2: Evaluando Homocedasticidad

También podemos utilizar el test de Breusch-Pagan, como se ve en la Tabla 4.1.

library(lmtest)
library(kableExtra)
# null: modelo homocedastico
resBP=bptest(reg2)
data.frame(list('BP'=resBP$statistic,
             'df'=resBP$parameter,
             "p-value"=resBP$p.value))%>%
    kable(caption = resBP$method)%>%kable_styling(full_width = F)
Table 4.1: Table 4.2: studentized Breusch-Pagan test
BP df p.value
BP 34.83587 3 1e-07

Es estadístico de BP sale 34.8358715 con un p-valor de 1.3195008^{-7}; de ahí que se rechaza que el modelo muestre homocedasticidad.

La presencia de heterocedasticidad afecta el cálculo de los p-valores, lo que afectará la validez de todo el modelo.

4.3 Normalidad de los residuos

Los residuos, la diferencia entre apropiaciondolar y \(\hat{apropiaciondolar}\), deben distribuirse de manera normal. Un qq-plot nos permite revelar si hay normalidad. Ver Figura 4.3.

# puntos cerca a la diagonal?
plot(reg2, 2)
Evaluando Normalidad de residuos

Figure 4.3: Evaluando Normalidad de residuos

Otra opción, es aplicar el test de Shapiro a los residuos, como se ve en la Tabla 4.3.

#NULL: Datos se distribuyen de manera normal
resSW=shapiro.test(reg2$residuals)
data.frame(list('SW'=resSW$statistic,
             "p-value"=resSW$p.value))%>%
    kable(caption = resSW$method)%>%kable_styling(full_width = F)
Table 4.3: Table 4.4: Shapiro-Wilk normality test
SW p.value
W 0.6820371 0

La falta de normalidad limita la capacidad de hacer inferencias a partir de lo encontrado.

4.4 No multicolinelidad

Si los predictores tienen una correlación muy alta entre sí, hay multicolinealidad. Ver Tabla 4.5.

library(DescTools)
## Warning: package 'DescTools' was built under R version 4.4.1
# > 5 es problematico
VIF(reg2) %>%kable(col.names = "VIF",caption ="Evaluando Multicolinealidad usando VIF (Variance Inflation Factors)" )%>%kable_styling(full_width = F)
Table 4.5: Table 4.6: Evaluando Multicolinealidad usando VIF (Variance Inflation Factors)
VIF
pctopo 1.000152
consejocomunal 1.036568
poblacioncienmil 1.036454

La presencia de la multicolinealidad no perjudica tanto el calculo de (\(\hat{apropiaciondolar}\)), pero evita calcular bien el efecto de cada regresor.

4.5 Valores influyentes

Hay casos particulares, que tienen la capacidad de trastocar lo que el modelo representa. A veces detectándolos y suprimiéndolos, podemos ver un mejor modelo. La Figura 4.4 muestra la posible existencia de los valores influyentes.

plot(reg2, 5)
Existencia de Valores Influyentes

Figure 4.4: Existencia de Valores Influyentes

Si queremos ver más de cerca los posible casos influyentes, le prestamos atencion al indice de Cook y a los valores predecidos (los hat values) de la Tabla 4.7

checkReg2=as.data.frame(influence.measures(reg2)$is.inf)
checkReg2[checkReg2$cook.d & checkReg2$hat,c('cook.d','hat')]%>%kable(caption = "Valores Influyentes criticos")%>%kable_styling(full_width = F)
Table 4.7: Table 4.8: Valores Influyentes criticos
cook.d hat
146 TRUE TRUE

Si alguna fila aparece en el resultado, ese caso está afectando los cálculos de la regresión (sin él habría otro resultado).


Bibliografía

Mejia Guinand, Luis Bernardo, Felipe Botero, and Juan Carlos Rodriguez Raga. 2008. “¿Pavimentando con votos? Apropiación presupuestal para proyectos de infraestructura vial en Colombia, 2002-2006.” Colombia Internacional, no. 68: 14–42. http://www.redalyc.org/articulo.oa?id=81211204002.
Wikipedia. 2022. “Standardized Coefficient.” https://en.wikipedia.org/w/index.php?title=Standardized_coefficient&oldid=1084836823.



al INICIO


  1. Esto puede significar que encuentres resultados distintos a los del trabajo original.↩︎

  2. Recuerda que el signo del coeficiente no importa para determinar el tamaño del efecto comparativo, importan su valor absoluto↩︎