Taller R para Ciencias Sociales - Sesión 4

Flacso - México

Author

Elvira-Torres & Gomez-Santiago

Published

October 15, 2023

Sesión 4 - Covarianza, Correlación y Regresión

En esta sesión se introducirá el siguiente contenido:

Covarianza Cov(x,y).
Correlación de Pearson \rho(x,y).
Regresión lineal simple.
Mínimos cuadrados ordinarios (OLS).
Intervalos de confianza para los coeficientes, la línea de regresión y predicción de valores.

Cargamos la base de volatilidad electoral en las elecciones de Jalisco para el periodo entre 2006 y 2018 y se realizan las transformaciones necesarias a las variables de tipo categóricas -de numeric a factor-.

A continuación se presentan una serie de medidas de asociación para variables cuantitativas continuas que permiten explorar la relación entre dos variables aleatorias x,y a partir de su variación conjunta. Previo a presentar dichas medidas y su corrrspondiente código de aplicación en R, presentemos el caso de como podrían relacionarse dos variables de tipo numeric.

La base de datos de volatilidad electoral en las elecciones de Jalisco para el periodo entre 2006 y 2018 contiene la medida de Volatilidad Total de Pedersen (Vt=\frac{\sum_{i=1}^{n}(p_{t-1}-p_t)}{2}) registrada para todas las elecciones municipales de la entidad en el periodo. Este índice indica la cantidad mínima total del cambio de preferencias electorales entre elecciones inmediatas de manera porcentual. Este será la variable de interés a explicar, o variable de salida de nuestros modelos. La segunda variable de interés por su relación con la primera, o variable predictora, será el Número efectivo de partidos de Laakso y Tageepera (NEP=\frac{1}{\sum_{i=1}^{n}p_i^2}). Este ínidce permite hacer un conteo de partidos relevantes de acuerdo a su porcentaje de votación recibida en la elección.

Ambas variables guardan una relación estudiada a nivel nacional en análisis realizados a partir de estudios longitudinales extensos. Desde esta evidencia surge la posibilidad de relacionar dichas variables a nivel subnacional. De esta manera evaluamos la relación en un primer momento a partir de un gráfico de dispersión utilizando la librería de GGplot2 con la función de geom_point().

A partir del gráfico anterior podemos concluir ciertas afirmaciones:
1. Existe una relación positiva aparente entre ambas variables.
2. La relación parece guardar una forma lineal.

A partir de esta primera exploración, se presentan medidas de asociación y prueba de hipótesis para respaldar estas conclusiones aparentes obtendias del análisis visual del gráfico.

Covarianza

La Covarianza nos indica el grado en que dos variables varian conjuntamente, y permite medir -en términos de la unidad de medida de las variables- en que sentido opera esta variación. La fórmula para obtener la covarianza es la siguiente: Cov(x,y)=\frac{1}{n-1}\sum_{i=1}^{n}(x_i-\bar{x})(y_i-\bar{y})

Para obtener la covarianza de dos variables de tipo numeric en R se utiliza el siguiente comando de la función cov():

#Covarianza

cov(x = base$nep,y = base$vt,use = "everything",method = "pearson")

[1] 4.132081

A pesar de que la dirección positiva de la relación indicada por el análisis de covarianza confirma en cierto grado la relación observada en el gráfico de dispersión, esta medida de asociación no nos permite obtener más conclsuiones debido al problema de la escala de medición (hecha en términos de la multiplicación de ambas escalas de medición, a priori no estandarizadas).

Correlación

EL méotodo de correlación de Pearson permite medir el grado de asociación lineal entre dos variables en una escala de entre 0 y 1, donde valores cercanos a la unidad correponden a una asociación más fuerte entre dos variables aleatorias. La fórmula para obtener la covarianza es la siguiente: r_{xy}=(\frac{1}{n-1})\sum_{i=1}^{n}(x_i-\bar{x})(y_i-\bar{y})(\frac{1}{S_xS_y})=\frac{\sum_{i=1}^{n}(x_i-\bar{x})(y_i-\bar{y})}{\sqrt{\sum_{i=1}^{n}(x_i-\bar{x})^2}\sqrt{\sum_{i=1}^{n}(y_i-\bar{y})^2}}

Para obtener el estadístico de correlación de dos variables de tipo numeric en R se utiliza el siguiente comando de la función cor():

#Correlación

cor(x = base$nep,y = base$vt,use = "everything",method = "pearson")

[1] 0.2910677

La interpretación del estadístico de correlación de Pearson generalmente define que una relación de 0.3 es una relación lineal moderada.

También se puede eleaborar un test de correlación con la función de cor.test() (para más detalles sobre el cálculo ver ?cor.test):

#Test de correlación de Pearson

cor.test(x = base$nep,y = base$vt,alternative = "two.sided",method = "pearson")


    Pearson's product-moment correlation

data:  base$nep and base$vt
t = 6.7894, df = 498, p-value = 3.214e-11
alternative hypothesis: true correlation is not equal to 0
95 percent confidence interval:
 0.2087041 0.3693315
sample estimates:
      cor 
0.2910677

Regresión lineal simple

La Regresión lineal es una herramienta del análisis estadístico que permite relacionar dos variables aleatorias con valores numéricos a partir del método de los Mínimos Cuadrados Ordinarios (OLS). Este método de estimación se realiza por medio de la Suma de los errores cuadrados (SSE) SSE=\epsilon^2_1+\epsilon^2_2+...+\epsilon^2_n=\sum_{i=1}^{n}\epsilon^2_i derivada de los valores de los parámetros \alpha (intercepto) y \beta (pendiente) que minimicen dicha suma de errores cuadrados.

Dada una relación simple de carácter y=f(x) entre dos variable aleatorias x,y donde la variable de respuesta (y) es explicada por la variable explicativa (x) se tiene la siguiente ecuación de la estimación de la regresión: y_i=\beta_0+\beta_1x_i+\epsilon_i

De esta ecuación -y únicamente para el caso de la regresión lineal simple- los parámetros que minimizan dicha suma de errores cuadrados (\beta_0 y \beta_1) se pueden calcular de la siguiente manera por el método de OLS. El cálculo para obtener la pendiente es: \hat{\beta_1}=\frac{\sum_{i=1}^n(x_i-\bar{x})(y_i-\bar{y})}{\sum_{i=1}^n(x_i-\bar{x})^2}=\frac{S_{xy}}{S_{xx}}

Donde S_{xy} es la Suma de los productos cruzados de x,y y S_{xx} corresponde a la Suma de cuadrados de x. La estimación para el intercepto se calcula a partir de la media de la variable y. La fórmula es la siguiente: \hat{\beta_0}=\bar{y}-\hat{\beta_1}\bar{x}

Prueba de hipótesis con regresión lineal

La Variabilidad Total (STT) de la variable y puede ser dividida de la siguiente manera: S_{yy}=STT donde S_{yy} representa la Suma de cuadrados en y o la Variabilidad Total en y.

La Variabilidad Total (STT) puede ser dividida entre la Variabilidad explicada en y (SSR) SSR=\sum_{i=1}^{n}(\hat{y_i}-\bar{y}) y la Variabilidad no explicada en y (SSE) SSE=\sum_{i=1}^{n}(y_i-\hat{y})^2=\sum_{i=1}^{n}\epsilon^2_i, de forma tal que: STT=SSR+SSE=\sum_{i=1}^{n}(\hat{y_i}-\bar{y}) + \sum_{i=1}^{n}(y_i-\hat{y})

Asímismo, el Error cuadrado medio (promedio de la Variabilidad no explicada en y (SSE)) denominado como \hat{\sigma^2} se obtiene de la siguiente manera: \hat{\sigma^2}=MSE=(\frac{1}{n-2})SSE

Gracias a la posibilidad de dividir la varianza es que existe una medida que indica la proporción de la varianza explicada por el modelo de regresión respecto de la varianza atribuida al error a partir de la ecuación de forma. Esta medida es la R^2 que toma valores entre 0 y 1 representando la dicha proporción de varianza atendida por la ecuación de regresión con los parámtros estimados \hat{y_i}=\beta_0+\beta_1x_i: R^2=\frac{SSR}{SST}=1-\frac{SSE}{SST}=\frac{\sum_{i=1}^n(\hat{y_i}-\bar{y})^2}{\sum_{i=1}^n(y_i-\bar{y})^2}

Una vez que se conoce el método de estimación de Mínimos Cuadrados Ordinarios (OLS) es posible plantear la relación entre las variables numéricas establecidas a manera de hipótesis falseables, donde: H_0\colon\beta_1=0 \text{ }\text{ ;}\text{ }H_1\colon\beta_1\neq 0

De forma que H_1 equivale a establcer que el efecto de x sobre la variable y es distitnto de cero. A partir de una prueba t se peude probar esta hipótesis. La fórmula de la prueba es la siguiente: t=\frac{\beta_1-0}{se_{\beta_{1}}}

donde se_{\beta_1} equivale a: se_{\beta_1}=\sqrt{\frac{MSE}{S_{xx}}}

Una vez que tenemos los elementos para tratar la variabilidad obtenida de la estimación de la regresión, procedemos a ejemplificar una regresión lineal simple a partir de la relación de dos variables de tipo numeric en R.

Regresión en R

Una vez más se busca analizar la relación entre las variables de tipo numeric de Volatilidad electoral y Número Efectivo de Partidos de las elecciones municiaples en la entidad de Jalisco para el periodo 2006-2018. De esta manera se podrá expresar que la fragmentación del sistema de partidos en la elección tiene una relación significativa en la volatilidad electoral medida en la forma de una hipótesis falseable de carácter H_1: \beta_1\neq0.

La ecuación de regresión a estimar tiene la siguiente forma: Vt_i=\beta_0+\beta1*NEP+\epsilon_i donde el valor de la Volatilidad Electoral medida depende del valor del intercepto y del efecto de la fragmentación, así como de un término aleatorio de error. La función que permite realizar una regresión en R es lm() del paquete stats():

#Regresión lineal simple

#estimación de coeficientes
lm(formula = vt ~ nep, data = base)


Call:
lm(formula = vt ~ nep, data = base)

Coefficients:
(Intercept)          nep  
     13.309        3.539

Esta función obtiene como resultado los parámetros de la regresión \beta_0 y \beta_1 que minimizan la recta de mínimos cuadrados (Least squares line) con intercepto y pendiente correspondiente a los parámetros descritos que describe a la ecuación de la estimación: \hat{y_i}=\hat{Vt}=\beta_0+\beta_1NEP

A partir de la sustitución de los valores obtenidos en la salida de la función lm() en la ecuación de la estimación se puede obtener un valor estimado para y. De esta ecuación es posible concluir que el Efecto Marginal (el efecto de un cambio unitario en x sobre el valor de y) de la fragmentación del sistema de partidos medida por el Número Efectivo de Partidos sobre el nivel de Volatilidad Electoral observado es de 3.5. Dada la forma particular de la ecuación de regresión con un solo regresor, el valor del coeficiente estimado para \beta_1 equivale al efecto marginal de x sobre y.

Sin embargo, de la salida de código previa con las estimaciones de los parámetros no se tiene la información completa sobre el proceso de regresión, ni tampoco los elementos necesarios para la prueba de hipótesis. Para obtener más información sobre la regresión guardamos la regresión en un objeto y utilizamos la función summary() para desplegar el resultado. Adicionalmente, la función confint() despliega los intervalos de confianza del coeficiente estimado a partir de la siguiente fórmula: \beta_1 \pm t_{\alpha/2,n-2}(\frac{\sqrt{MSE}}{\sqrt{(S_{xx})}})

# Regresión lineal simple 2

rls <- lm(formula = vt ~ nep, data = base)
summary(rls)


Call:
lm(formula = vt ~ nep, data = base)

Residuals:
    Min      1Q  Median      3Q     Max 
-21.590  -9.508  -1.197   6.961  61.363 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept)  13.3091     1.9456   6.841 2.32e-11 ***
nep           3.5386     0.5212   6.789 3.21e-11 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 12.58 on 498 degrees of freedom
Multiple R-squared:  0.08472,   Adjusted R-squared:  0.08288 
F-statistic:  46.1 on 1 and 498 DF,  p-value: 3.214e-11

# Intervalos de confianza al 95%

confint(rls)

               2.5 %    97.5 %
(Intercept) 9.486511 17.131693
nep         2.514574  4.562582

De esta salida obtenemos mayor información:

La distribución de los residuales (errores de estimación) de la regresión por sus medidas de distribución.
Los coeficientes estimados (parámetros) con su respectivo error estándar (SE).
La prueba t correspondiente a cada variable y su respectiva probabilidad acumulada en la forma de p-valor.
El error estándar medio (MSE) y los grados de libertad (df) de la regresión.
La medida de R^2 y un valor de R^2 ajustado al número de regresores.
Una prueba F para determinar si el modelo es mejor explicando la variación que un modelo nulo que solo considere al intercepto.
Los intervalos de confianza para los coeficientes.

En este caso, a pesar del bajo valor de R^2, la estimación del parámetro de \beta_1 para la relación entre fragmentación y se sostiene con un valor t de 6.78 y un p-valor menor a un nivel de confianza del 95% (\alpha=.05). A continuación se muestra de manera gráfica la prueba t:

A continuación se presenta una forma gráfica de visualizar a la Línea de Mínimos Cuadrados -con su respectivo intervalo de confianza para la pendiente- estimada a partir de la regresión simple entre dos variables con un diagrama de dispersión a partir de la extensión geom_smooth() de ggplot2. El siguiente gráfico de dispersión de los valores reales de ambas variables ilustra la relación de correspondencia entre la predicción del valor estimado de Volatilidad electoral a partir de la ecuación de regresión para cada valor particular del Número Efectivo de Partidos y la línea de mínimos cuadrados.

Es necesario destacar que los valores predichos de y para cada valor de x no están dados de manera automática. Por ser el resultado de la ecuación de un modelo de regresión, existe una posibilidad de incertidumbre acerca de los valores particulares de la pendiente de la línea de mínimos cuadrados estimada. Esta incertidumbre se refleja en la sección oscurecida de la línea de ajuste, que ilustra el rango de valores que podría tomar el efecto de x sobre y para cada valor en la recta de x. La fórmula para obtener los intervalos para cada valor de x_i se compone del producto entre la raíz del error de la línea de regresión y la raíz del error estándar de la estimación para x_k. Su estructura es la siguiente: \hat{y_h}\pm t_{\alpha/2,n-2}\sqrt{MSE}\sqrt{\frac{1}{n}+\frac{(x-\bar{x})^2}{S_{xx}}} donde \hat{y_h} es la estimación para el valor de y cuando x=h y \sqrt{MSE}\sqrt{\frac{1}{n}+\frac{(x-\bar{x})^2}{S_{xx}}} es la predicción para el error estándar cuando x=h.

Además de poder establecer intervalos de confianza para la pendiente de la recta de regresión, es posible establecer bajo un nivel de confianza un rango de valores para las estimaciones de y a partir de la recta de regresión. La f´romula para su cálculo es la siguiente: \hat{y_h}\pm t_{\alpha/2,n-2}\sqrt{MSE(1+\frac{1}{n}+\frac{(x-\bar{x})^2}{S_{xx}})}

A continuación se presenta una forma de visualizar los intervalos a manera de bandas por debajo y por encima de la línea de regresión:

El gráfico de dispersión anterior contiene toda la información obtenida a partir del modelo de regresión lineal simple sobre la relación de las variable x y y en la forma de f(x)=y.

Supuestos de la regresión lineal

Por último, es necesario destacar que la regresión lineal debe realizarse sobre datos que cumplan las siguientes condiciones:

Linealidad - la relación estimada por la regresión lineal es de carácter lineal, por lo que la relación entre las variables debe cumplir esta forma funcional. Este supuesto puede verificarse de manera gráfica con un diagrama de dispersión o con las medidas de asociación de correlación.

Independencia - los errores de las variables aleatorias deben ser independientes. Si se trabaja con datos de varios periodos temporales puede verificarse si los errores muestran tendencias lineales en el tiempo.

Normalidad de los errores - los errores deben distribuirse de manera normal al rededor de cero. Este supuesto permite realizar las pruebas de hipótesis. Puede comprobarse a partir de la distribución de los errores indicada en la salida del comando lm(), donde la mediana debe acercarse a cero, o a partir de un Gráfico de distribución normal con la función plot(), la cual requiere del objeto creado con la función de lm() seguido de un parámetro wich=2:

Homocedasticidad - la variabilidad de las observaciones al rededor de la Línea de Mínimos Cuadrados debe ser constante. Puede verificarse a partir de un gráfico de residuales frente a valores predecidos. A continuación se presenta un gráfico de residuales utilizando el código plot()la cual requiere del objeto creado con la función de lm() seguido de un parámetro wich=1:

Multicolinearidad - o dependencia lineal cercana entre los regresores. No debe existir relación lineal entre las variables regresoras. Pueden realizarse pruebas con el Variance Inflation Factor o simplemente descartarse variables que tengan una correlación fuerte.

--- title: "Taller R para Ciencias Sociales - Sesión 4" subtitle: "Flacso - México" author: "Elvira-Torres & Gomez-Santiago" date: "October 15, 2023" toc: true format: html: css: styles.css html-math-method: katex code-tools: true self-contained: true execute: warning: false --- ## Sesión 4 - Covarianza, Correlación y Regresión En esta sesión se introducirá el siguiente contenido: 1. Covarianza $Cov(x,y)$. 2. Correlación de Pearson $\rho(x,y)$. 3. Regresión lineal simple. 4. Mínimos cuadrados ordinarios ($OLS$). 5. Intervalos de confianza para los coeficientes, la línea de regresión y predicción de valores. Cargamos la base de volatilidad electoral en las elecciones de Jalisco para el periodo entre 2006 y 2018 y se realizan las transformaciones necesarias a las variables de tipo categóricas -de `numeric` a `factor`-. ```{r,include=FALSE} rm(list = ls()) ``` ```{r, include=FALSE} library(tidyverse) library(readxl) library(glue) library(ggtext) library(marginaleffects) ``` ```{r,warning=FALSE,include=FALSE} # Nueva base con variables de volatilidad nuevas. # leer base panel base <- read_excel("C:/Users/diego/Documents/FLACSO/Tesis Maestría/Jalisco municipios V.xlsx",4) names(base) # factor base$alt <- as.factor(base$alt) levels(base$alt) <- c("No", "Sí") base$urbano <- as.factor(base$urbano) levels(base$urbano) <- c("Interior", "Metropolitano") base$altac <- as.factor(base$altac) ``` A continuación se presentan una serie de medidas de asociación para variables cuantitativas continuas que permiten explorar la relación entre dos variables aleatorias $x,y$ a partir de su variación conjunta. Previo a presentar dichas medidas y su corrrspondiente código de aplicación en `R`, presentemos el caso de como podrían relacionarse dos variables de tipo `numeric`. La base de datos de volatilidad electoral en las elecciones de Jalisco para el periodo entre 2006 y 2018 contiene la medida de **Volatilidad Total de Pedersen** ($Vt=\frac{\sum_{i=1}^{n}(p_{t-1}-p_t)}{2}$) registrada para todas las elecciones municipales de la entidad en el periodo. Este índice indica la cantidad mínima total del cambio de preferencias electorales entre elecciones inmediatas de manera porcentual. Este será la variable de interés a explicar, o *variable de salida* de nuestros modelos. La segunda variable de interés por su relación con la primera, o *variable predictora*, será el **Número efectivo de partidos de Laakso y Tageepera** ($NEP=\frac{1}{\sum_{i=1}^{n}p_i^2}$). Este ínidce permite hacer un conteo de partidos relevantes de acuerdo a su porcentaje de votación recibida en la elección. Ambas variables guardan una relación estudiada a nivel nacional en análisis realizados a partir de estudios longitudinales extensos. Desde esta evidencia surge la posibilidad de relacionar dichas variables a nivel subnacional. De esta manera evaluamos la relación en un primer momento a partir de un **gráfico de dispersión** utilizando la librería de `GGplot2` con la función de `geom_point()`. ```{r,echo=FALSE} #Gráfico de dispersión básico base %>% ggplot(aes(x = nep, y = vt)) + geom_point(size=0) + geom_text(label= base$alfa,size=3)+ labs(x = "NEP",y = "Volatilidad", title = "Volatilidad electoral ~ Fragmentación")+ theme_bw() ``` A partir del gráfico anterior podemos concluir ciertas afirmaciones:\ 1. Existe una relación positiva aparente entre ambas variables.\ 2. La relación parece guardar una forma lineal.\ A partir de esta primera exploración, se presentan medidas de asociación y prueba de hipótesis para respaldar estas conclusiones aparentes obtendias del análisis visual del gráfico. ## Covarianza La *Covarianza* nos indica el grado en que dos variables varian conjuntamente, y permite medir -en términos de la unidad de medida de las variables- en que sentido opera esta variación. La fórmula para obtener la covarianza es la siguiente: $$Cov(x,y)=\frac{1}{n-1}\sum_{i=1}^{n}(x_i-\bar{x})(y_i-\bar{y})$$ Para obtener la covarianza de dos variables de tipo `numeric` en `R` se utiliza el siguiente comando de la función `cov()`: ```{r} #Covarianza cov(x = base$nep,y = base$vt,use = "everything",method = "pearson") ``` A pesar de que la dirección positiva de la relación indicada por el análisis de covarianza confirma en cierto grado la relación observada en el gráfico de dispersión, esta medida de asociación no nos permite obtener más conclsuiones debido al problema de la escala de medición (hecha en términos de la multiplicación de ambas escalas de medición, *a priori* no estandarizadas). ## Correlación EL méotodo de correlación de Pearson permite medir el grado de asociación lineal entre dos variables en una escala de entre $0$ y $1$, donde valores cercanos a la unidad correponden a una asociación más fuerte entre dos variables aleatorias. La fórmula para obtener la covarianza es la siguiente: $$r_{xy}=(\frac{1}{n-1})\sum_{i=1}^{n}(x_i-\bar{x})(y_i-\bar{y})(\frac{1}{S_xS_y})=\frac{\sum_{i=1}^{n}(x_i-\bar{x})(y_i-\bar{y})}{\sqrt{\sum_{i=1}^{n}(x_i-\bar{x})^2}\sqrt{\sum_{i=1}^{n}(y_i-\bar{y})^2}}$$ Para obtener el estadístico de correlación de dos variables de tipo `numeric` en `R` se utiliza el siguiente comando de la función `cor()`: ```{r} #Correlación cor(x = base$nep,y = base$vt,use = "everything",method = "pearson") ``` La interpretación del estadístico de *correlación de Pearson* generalmente define que una relación de $0.3$ es una relación lineal moderada. También se puede eleaborar un test de correlación con la función de `cor.test()` (para más detalles sobre el cálculo ver `?cor.test`): ```{r} #Test de correlación de Pearson cor.test(x = base$nep,y = base$vt,alternative = "two.sided",method = "pearson") ``` ## Regresión lineal simple La **Regresión lineal** es una herramienta del análisis estadístico que permite relacionar dos variables aleatorias con valores numéricos a partir del método de los **Mínimos Cuadrados Ordinarios (OLS)**. Este método de estimación se realiza por medio de la *Suma de los errores cuadrados (SSE)* $SSE=\epsilon^2_1+\epsilon^2_2+...+\epsilon^2_n=\sum_{i=1}^{n}\epsilon^2_i$ derivada de los valores de los parámetros $\alpha$ (intercepto) y $\beta$ (pendiente) que minimicen dicha suma de errores cuadrados. Dada una relación simple de carácter $y=f(x)$ entre dos variable aleatorias $x,y$ donde la *variable de respuesta (y)* es explicada por la *variable explicativa (x)* se tiene la siguiente ecuación de la estimación de la regresión: $$y_i=\beta_0+\beta_1x_i+\epsilon_i$$ De esta ecuación -y únicamente para el caso de la regresión lineal simple- los parámetros que minimizan dicha suma de errores cuadrados ($\beta_0$ y $\beta_1$) se pueden calcular de la siguiente manera por el método de **OLS**. El cálculo para obtener la pendiente es: $$\hat{\beta_1}=\frac{\sum_{i=1}^n(x_i-\bar{x})(y_i-\bar{y})}{\sum_{i=1}^n(x_i-\bar{x})^2}=\frac{S_{xy}}{S_{xx}}$$ Donde $S_{xy}$ es la *Suma de los productos cruzados de x,y* y $S_{xx}$ corresponde a la *Suma de cuadrados de x*. La estimación para el *intercepto* se calcula a partir de la media de la variable $y$. La fórmula es la siguiente: $$\hat{\beta_0}=\bar{y}-\hat{\beta_1}\bar{x}$$ ## Prueba de hipótesis con regresión lineal La *Variabilidad Total (STT)* de la variable $y$ puede ser dividida de la siguiente manera: $$S_{yy}=STT$$ donde $S_{yy}$ representa la *Suma de cuadrados en y* o la *Variabilidad Total en y*. La *Variabilidad Total (STT)* puede ser dividida entre la *Variabilidad explicada en y (SSR)* $SSR=\sum_{i=1}^{n}(\hat{y_i}-\bar{y})$ y la *Variabilidad no explicada en y (SSE)* $SSE=\sum_{i=1}^{n}(y_i-\hat{y})^2=\sum_{i=1}^{n}\epsilon^2_i$, de forma tal que: $$STT=SSR+SSE=\sum_{i=1}^{n}(\hat{y_i}-\bar{y}) + \sum_{i=1}^{n}(y_i-\hat{y})$$ Asímismo, el *Error cuadrado medio* (promedio de la *Variabilidad no explicada en y (SSE)*) denominado como $\hat{\sigma^2}$ se obtiene de la siguiente manera: $$\hat{\sigma^2}=MSE=(\frac{1}{n-2})SSE$$ Gracias a la posibilidad de dividir la varianza es que existe una medida que indica la proporción de la varianza explicada por el modelo de regresión respecto de la varianza atribuida al error a partir de la ecuación de forma. Esta medida es la $R^2$ que toma valores entre $0$ y $1$ representando la dicha proporción de varianza atendida por la ecuación de regresión con los parámtros estimados $\hat{y_i}=\beta_0+\beta_1x_i$: $$R^2=\frac{SSR}{SST}=1-\frac{SSE}{SST}=\frac{\sum_{i=1}^n(\hat{y_i}-\bar{y})^2}{\sum_{i=1}^n(y_i-\bar{y})^2}$$ Una vez que se conoce el método de estimación de **Mínimos Cuadrados Ordinarios (OLS)** es posible plantear la relación entre las variables numéricas establecidas a manera de hipótesis falseables, donde: $$H_0\colon\beta_1=0 \text{ }\text{ ;}\text{ }H_1\colon\beta_1\neq 0$$ De forma que $H_1$ equivale a establcer que el efecto de $x$ sobre la variable $y$ es distitnto de cero. A partir de una prueba *t* se peude probar esta hipótesis. La fórmula de la prueba es la siguiente: $$t=\frac{\beta_1-0}{se_{\beta_{1}}}$$ donde $se_{\beta_1}$ equivale a: $$se_{\beta_1}=\sqrt{\frac{MSE}{S_{xx}}}$$ Una vez que tenemos los elementos para tratar la variabilidad obtenida de la estimación de la regresión, procedemos a ejemplificar una regresión lineal simple a partir de la relación de dos variables de tipo `numeric` en `R`. ## Regresión en R Una vez más se busca analizar la relación entre las variables de tipo `numeric` de *Volatilidad electoral* y *Número Efectivo de Partidos* de las elecciones municiaples en la entidad de Jalisco para el periodo 2006-2018. De esta manera se podrá expresar que la fragmentación del sistema de partidos en la elección tiene una relación significativa en la volatilidad electoral medida en la forma de una hipótesis falseable de carácter H_1: $\beta_1\neq0$. La ecuación de regresión a estimar tiene la siguiente forma: $$Vt_i=\beta_0+\beta1*NEP+\epsilon_i$$ donde el valor de la Volatilidad Electoral medida depende del valor del intercepto y del efecto de la fragmentación, así como de un término aleatorio de error. La función que permite realizar una regresión en `R` es `lm()` del paquete `stats()`: ```{r} #Regresión lineal simple #estimación de coeficientes lm(formula = vt ~ nep, data = base) ``` Esta función obtiene como resultado los *parámetros* de la regresión $\beta_0$ y $\beta_1$ que minimizan la recta de mínimos cuadrados (*Least squares line*) con intercepto y pendiente correspondiente a los parámetros descritos que describe a la ecuación de la estimación: $$\hat{y_i}=\hat{Vt}=\beta_0+\beta_1NEP$$ A partir de la sustitución de los valores obtenidos en la salida de la función `lm()` en la ecuación de la estimación se puede obtener un valor estimado para $y$. De esta ecuación es posible concluir que el *Efecto Marginal* (el efecto de un cambio unitario en $x$ sobre el valor de $y$) de la fragmentación del sistema de partidos medida por el *Número Efectivo de Partidos* sobre el nivel de *Volatilidad Electoral* observado es de $3.5$. Dada la forma particular de la ecuación de regresión con un solo regresor, el valor del coeficiente estimado para $\beta_1$ equivale al efecto marginal de $x$ sobre $y$. Sin embargo, de la salida de código previa con las estimaciones de los parámetros no se tiene la información completa sobre el proceso de regresión, ni tampoco los elementos necesarios para la prueba de hipótesis. Para obtener más información sobre la regresión guardamos la regresión en un objeto y utilizamos la función `summary()` para desplegar el resultado. Adicionalmente, la función `confint()` despliega los intervalos de confianza del coeficiente estimado a partir de la siguiente fórmula: $$\beta_1 \pm t_{\alpha/2,n-2}(\frac{\sqrt{MSE}}{\sqrt{(S_{xx})}})$$ ```{r} # Regresión lineal simple 2 rls <- lm(formula = vt ~ nep, data = base) summary(rls) # Intervalos de confianza al 95% confint(rls) ``` De esta salida obtenemos mayor información: 1. La distribución de los residuales (errores de estimación) de la regresión por sus medidas de distribución. 2. Los coeficientes estimados (parámetros) con su respectivo error estándar (*SE*). 3. La prueba $t$ correspondiente a cada variable y su respectiva probabilidad acumulada en la forma de $p-valor$. 4. El error estándar medio (*MSE*) y los grados de libertad (*df*) de la regresión. 5. La medida de $R^2$ y un valor de $R^2$ ajustado al número de regresores. 6. Una prueba $F$ para determinar si el modelo es mejor explicando la variación que un modelo nulo que solo considere al intercepto. 7. Los intervalos de confianza para los coeficientes. En este caso, a pesar del bajo valor de $R^2$, la estimación del parámetro de $\beta_1$ para la relación entre fragmentación y se sostiene con un valor $t$ de $6.78$ y un *p-valor* menor a un nivel de confianza del 95% ($\alpha=.05$). A continuación se muestra de manera gráfica la prueba $t$: ```{r,echo=FALSE} curve(dt(x,df = 498),xlim = c(-7, 7)) abline(v = 0,col="blue") text(0, 0.15, "H0", pos = 1, col = "blue") text(0, 0.13, "B1 = 0", pos = 1, col = "blue") abline(v = qnorm(.025),col="red") text(qnorm(.025), 0.15, "Valor", pos = 1, col = "red") text(qnorm(.025), 0.13, "Crítico", pos = 1, col = "red") abline(v = qnorm(.975),col="red") text(qnorm(.975), 0.15, "Valor", pos = 1, col = "red") text(qnorm(.975), 0.13, "Crítico", pos = 1, col = "red") abline(v = 6.7,col="purple") text(5.4, 0.15, "Coeficiente", pos = 1, col = "purple") text(5.4, 0.13, "Estimado", pos = 1, col = "purple") ``` A continuación se presenta una forma gráfica de visualizar a la *Línea de Mínimos Cuadrados* -con su respectivo intervalo de confianza para la pendiente- estimada a partir de la regresión simple entre dos variables con un diagrama de dispersión a partir de la extensión `geom_smooth()` de `ggplot2`. El siguiente gráfico de dispersión de los valores reales de ambas variables ilustra la relación de correspondencia entre la predicción del valor estimado de *Volatilidad electoral* a partir de la ecuación de regresión para cada valor particular del *Número Efectivo de Partidos* y la línea de mínimos cuadrados. ```{r,echo=FALSE,message=FALSE} # Asignamos la predicción del modelo para la variable "y" a la base original base$fit <- rls$fitted.values # Extraemos los coeficientes del modelo de regresión rls en un dataframe coefsrls <- tibble(coef(object = rls)) coef_nep <- coefsrls$`coef(object = rls)`[2] intercept <- coefsrls$`coef(object = rls)`[1] # Creamos un dataframe con los valores de "x" para estimar "y" a <- tibble(valor_nep = c(2,4,6)) # Creamos una función para calcular los valores ajustados de "x" con la ecuación de regresión nep_slope <- function(x) (coef_nep*x)+intercept # Creamos un dataframe con los valores de "y" estimados a partir de valores de "x" de interés b <- tibble(fit_rls=nep_slope(c(2,4,6))) # Creamos un solo dtaframe con los valores para "x" y los valores ajustados de "y" z <- cbind(a,b) # Creamos un vector de etiquetas para graficar z <- z %>% mutate(nice_label = glue("NEP: {valor_nep} ", "Vt ajustada: {fit_rls} ", "Slope: **{3.54}**")) # Gráfico de dispersión con linea de mínimos cuadrados e intervalo de confianza para la pendiente base %>% ggplot(aes(x = nep, y = vt)) + geom_point() + geom_smooth(method="lm",formula = y ~ x,color = "cyan",linewidth = .5,alpha=.1,se=.1)+ geom_point(data=base,aes(x=nep,y=fit),color="red",alpha=.1)+ labs(x = "NEP",y = "Valor de ajuste de Volatilidad", title = "Ajuste de Volatilidad electoral ~ Fragmentación")+ theme_bw()+ geom_point(data=z,aes(x=valor_nep,y=fit_rls), size = 4, shape = 18, color = "red")+ geom_richtext(data = z, aes(x = valor_nep, y = fit_rls, label = nice_label), nudge_y = -5,nudge_x = +.3,size=2,alpha=.7) + labs(x = "Número Efectivo de Partidos", y = "Volatilidad electoral ajustada", subtitle = "Valores ajustados - Línea de Mínimos Cuadrados")+ theme(plot.subtitle = element_markdown())+ scale_y_continuous(limits = c(0, 90), breaks = seq(0, 90, by = 10)) + scale_x_continuous(limits = c(1, 10.5), breaks = seq(0, 10, by = 1)) ``` Es necesario destacar que los valores predichos de $y$ para cada valor de $x$ no están dados de manera automática. Por ser el resultado de la ecuación de un modelo de regresión, existe una posibilidad de incertidumbre acerca de los valores particulares de la pendiente de la línea de mínimos cuadrados estimada. Esta incertidumbre se refleja en la sección oscurecida de la línea de ajuste, que ilustra el rango de valores que podría tomar el efecto de $x$ sobre $y$ para cada valor en la recta de $x$. La fórmula para obtener los intervalos para cada valor de $x_i$ se compone del producto entre la raíz del error de la línea de regresión y la raíz del error estándar de la estimación para $x_k$. Su estructura es la siguiente: $$\hat{y_h}\pm t_{\alpha/2,n-2}\sqrt{MSE}\sqrt{\frac{1}{n}+\frac{(x-\bar{x})^2}{S_{xx}}}$$ donde $\hat{y_h}$ es la estimación para el valor de $y$ cuando $x=h$ y $\sqrt{MSE}\sqrt{\frac{1}{n}+\frac{(x-\bar{x})^2}{S_{xx}}}$ es la predicción para el error estándar cuando $x=h$. ```{r,echo=FALSE} #Creamos una nueva base base2 <- base %>% select(vt,nep,fit) # Calcular valores de SSE, MSE y Sxx sse <- sum((base2$vt - base2$fit)^2) # suma de los residuales cuadrados SEE mse <- sse / (length(base$vt) - 2) #Mean Squared Error sxx <- sum((base$nep-mean(base$nep))^2) #variabiliadd de x sxi <- (base$nep-mean(base$nep))^2 #variabción de i cuando xh #Valor de prueba t n-2 df t.val <- qt(0.975, length(base$vt) - 2) #Fórmula del se para cada valor de xh base2$se <- (sqrt(mse)) * (sqrt(((1 / length(base$vt))) + ((sxi) / (sxx)))) #Columnas con los intervalos para cada observación base2$slope.upper <- (base2$fit + (t.val * base2$se)) base2$slope.lower <- (base2$fit - (t.val * base2$se)) # Gráfico con ggplot a partir de los valores calculados base2 %>% ggplot(aes(x = nep, y = vt)) + geom_point(alpha=.5,color="grey") + geom_smooth(method="lm",formula = y ~ x,color = "cyan",linewidth = .5,alpha=FALSE)+ geom_ribbon(aes(ymin = slope.lower, ymax = slope.upper), alpha = .1,fill="black")+ geom_point(aes(x=nep,y=fit),color="red",alpha=.3)+ geom_errorbar(mapping = aes(ymin = slope.lower, ymax = slope.upper),color="purple")+ labs(x = "NEP",y = "Valor de ajuste de Volatilidad", title = "Ajuste de Volatilidad electoral ~ Fragmentación")+ #theme_bw()+ #geom_point(data=z,aes(x=valor_nep,y=fit_rls), size = 4, shape = 18, color = "red")+ #geom_richtext(data = z, aes(x = valor_nep, y = fit_rls, label = nice_label), # nudge_y = -5,nudge_x = +.3,size=2,alpha=.7) + labs(x = "Número Efectivo de Partidos", y = "Volatilidad electoral ajustada", subtitle = "Valores ajustados Int. Conf. 95% Línea de Mínimos Cuadrados")+ theme(plot.subtitle = element_markdown())+ scale_y_continuous(limits = c(0, 90), breaks = seq(0, 90, by = 10)) + scale_x_continuous(limits = c(1, 10.5), breaks = seq(0, 10, by = 1)) ``` ```{r,include=FALSE} #Cálculo paralelo con `marginaleffects` mareff <- slopes(rls) #conf.low y conf.high son los valores que puede tomar B1 en cada obs. #si se quiere calcular el listón se debe calcular el valor predicho sustituyendo B1 en la ecuación con cl/ch #Cálculo a partir de objeto con `slpoes()` mareff$intl <- (mean(mareff$predicted)-(mareff$conf.low*mean(mareff$nep))) mareff$intu <- (mean(mareff$predicted)-(mareff$conf.high*mean(mareff$nep))) mareff$mcl <- (mareff$conf.low*mareff$nep)+mareff$intl mareff$mcu <- (mareff$conf.high*mareff$nep)+mareff$intu #Gráfico simple de los valores predichos de los coeficientes para cada valor de x plot_slopes(rls, variables = "nep", condition = "nep")+ scale_y_continuous(limits = c(0, 10), breaks = seq(0, 10, by = 1)) + scale_x_continuous(limits = c(1, 10.5), breaks = seq(0, 10, by = 1)) #Gráfico con ggplot para comprobar el cálculo de las predicciones para cada valor de xh con sus ic mareff %>% ggplot(aes(x = nep, y = vt)) + geom_point(alpha=.5,color="grey") + geom_smooth(method="lm",formula = y ~ x,color = "cyan",linewidth = .5,alpha=.1)+ geom_ribbon(aes(ymin = mcl, ymax = mcu), alpha = .1,fill="black") # Gráfico con ggplot base2 %>% ggplot(aes(x = nep, y = vt)) + geom_point(alpha=.5,color="grey") + geom_smooth(method="lm",formula = y ~ x,color = "cyan",linewidth = .5,alpha=FALSE)+ geom_ribbon(aes(ymin = slope.lower, ymax = slope.upper), alpha = .1,fill="black")+ geom_point(aes(x=nep,y=fit),color="red",alpha=.3)+ geom_errorbar(mapping = aes(ymin = slope.lower, ymax = slope.upper),color="purple")+ labs(x = "NEP",y = "Valor de ajuste de Volatilidad", title = "Ajuste de Volatilidad electoral ~ Fragmentación")+ #theme_bw()+ #geom_point(data=z,aes(x=valor_nep,y=fit_rls), size = 4, shape = 18, color = "red")+ #geom_richtext(data = z, aes(x = valor_nep, y = fit_rls, label = nice_label), # nudge_y = -5,nudge_x = +.3,size=2,alpha=.7) + labs(x = "Número Efectivo de Partidos", y = "Volatilidad electoral ajustada", subtitle = "Valores ajustados Int. Conf. 95% Línea de Mínimos Cuadrados")+ theme(plot.subtitle = element_markdown())+ scale_y_continuous(limits = c(0, 90), breaks = seq(0, 90, by = 10)) + scale_x_continuous(limits = c(1, 10.5), breaks = seq(0, 10, by = 1)) ``` Además de poder establecer intervalos de confianza para la pendiente de la recta de regresión, es posible establecer bajo un nivel de confianza un rango de valores para las estimaciones de $y$ a partir de la recta de regresión. La f´romula para su cálculo es la siguiente: $$\hat{y_h}\pm t_{\alpha/2,n-2}\sqrt{MSE(1+\frac{1}{n}+\frac{(x-\bar{x})^2}{S_{xx}})}$$ A continuación se presenta una forma de visualizar los intervalos a manera de bandas por debajo y por encima de la línea de regresión: ```{r,echo=FALSE,warning=FALSE} ## INTERVALOS PARA LA LÍNEA DE REGRESIÓN #Calculo del error estándar para predicciones a partir de un df con datos de interés # Creamos un dataframe con los valores de "x" para estimar "y" c <- tibble(nep = seq(1.5,10,by=.5)) # Creamos un dataframe con los valores de "y" estimados a partir de valores de "x" de interés d <- tibble(vt=nep_slope(seq(1.5,10,by=.5))) # Creamos un solo dtaframe con los valores para "x" y los valores ajustados de "y" w <- cbind(c,d) #Nuevos valores individuales para las predicciones sxi <- (w$nep-mean(base$nep))^2 #Cálculo de se para cada valor de xh w$se <- (sqrt(mse)) * (sqrt(((1 / length(base$vt))) + ((sxi) / (sxx)))) #Valores de intervalos de confianza para linea de regresión w$slope.upper <- (w$vt + (t.val * w$se)) w$slope.lower <- (w$vt - (t.val * w$se)) ## ETIQUETAS CON VALORES DE INTRVALOS DE CONFIANZA DE PREDICCIÓN z$nep <- z$valor_nep z$vt <- z$fit_rls #Nuevos valores individuales para las predicciones sxi <- (z$nep-mean(base$nep))^2 #Cálculo de se para predicciones z$seb <- sqrt(mse*(1+(1/length(base$vt))+(sxi/sxx))) #Valores de intervalos para las predicciones z$confl <- z$vt-t.val*z$seb z$confu <- z$vt+t.val*z$seb #Etiquetas para el gráfico z <- z %>% mutate(nice_label = glue("NEP: {valor_nep} ", "Vt ajustada: {fit_rls} ", "Lim. Sup.: {round(confu,1)} ", "Lim. Inf.: {round(confl,1)} ", "Beta1: {round(confl,1)}")) ## BANDAS DE CONFIANZA PARA VALORES DE PREDICCIÓN #Nuevos valores individuales para las predicciones sxi <- (base2$nep-mean(base2$nep))^2 #Cálculo de se para predicciones base2$seb <- sqrt(mse*(1+(1/length(base2$vt))+(sxi/sxx))) #Valores de intervalos para las predicciones base2$confl <- base2$fit-t.val*base2$seb base2$confu <- base2$fit+t.val*base2$seb #Gráfico de efectos marginales base2 %>% ggplot(aes(x = nep, y = vt)) + geom_point(alpha=.5,color="grey") + geom_smooth(method="lm",formula = y ~ x,color = "cyan",linewidth = .5,alpha=FALSE)+ geom_ribbon(aes(ymin = slope.lower, ymax = slope.upper), alpha = .1,fill="black")+ geom_point(data=z,aes(x=valor_nep,y=fit_rls), size = 4, shape = 18, color = "red")+ geom_errorbar(data = w,mapping = aes(ymin = slope.lower,ymax = slope.upper),color="purple3",width=.1)+ geom_line(aes(y=confu), color = "purple4", linetype = "dashed",linewidth=.8)+ geom_line(aes(y=confl), color = "purple4", linetype = "dashed",linewidth=.8)+ labs(x = "NEP",y = "Valor de ajuste de Volatilidad", title = "Ajuste de Volatilidad electoral ~ Fragmentación")+ geom_richtext(data = z, aes(x = valor_nep, y = fit_rls, label = nice_label), nudge_y = -7,nudge_x = .5,size=2,alpha=.7) + labs(x = "Número Efectivo de Partidos", y = "Volatilidad electoral ajustada", subtitle = "Valores ajustados Int. Conf. 95% Línea de Mínimos Cuadrados")+ theme(plot.subtitle = element_markdown())+ scale_y_continuous(limits = c(0, 90), breaks = seq(0, 90, by = 10)) + scale_x_continuous(limits = c(1, 10.5), breaks = seq(0, 10, by = 1)) ``` El gráfico de dispersión anterior contiene toda la información obtenida a partir del modelo de regresión lineal simple sobre la relación de las variable $x$ y $y$ en la forma de $f(x)=y$. ```{r,include=FALSE} #Forma alternativa con predict predicc <- predict(rls, interval = "prediction") base5 <- cbind(base2,predicc) ``` ## Supuestos de la regresión lineal Por último, es necesario destacar que la regresión lineal debe realizarse sobre datos que cumplan las siguientes condiciones: **Linealidad** - la relación estimada por la regresión lineal es de carácter lineal, por lo que la relación entre las variables debe cumplir esta forma funcional. Este supuesto puede verificarse de manera gráfica con un diagrama de dispersión o con las medidas de asociación de correlación. **Independencia** - los errores de las variables aleatorias deben ser independientes. Si se trabaja con datos de varios periodos temporales puede verificarse si los errores muestran tendencias lineales en el tiempo. **Normalidad de los errores** - los errores deben distribuirse de manera normal al rededor de cero. Este supuesto permite realizar las pruebas de hipótesis. Puede comprobarse a partir de la distribución de los errores indicada en la salida del comando `lm()`, donde la mediana debe acercarse a cero, o a partir de un *Gráfico de distribución normal* con la función `plot()`, la cual requiere del objeto creado con la función de `lm()` seguido de un parámetro $wich=2$: ```{r,echo=FALSE} #Gráfico de distribución normal de los errores plot(rls,which=2) ``` **Homocedasticidad** - la variabilidad de las observaciones al rededor de la *Línea de Mínimos Cuadrados* debe ser constante. Puede verificarse a partir de un *gráfico de residuales* frente a valores predecidos. A continuación se presenta un gráfico de residuales utilizando el código `plot()`la cual requiere del objeto creado con la función de `lm()` seguido de un parámetro $wich=1$: ```{r,echo=FALSE} #Gráfico de residuales ~ valores predichos plot(rls, which = 1) ``` **Multicolinearidad** - o dependencia lineal cercana entre los regresores. No debe existir relación lineal entre las variables regresoras. Pueden realizarse pruebas con el *Variance Inflation Factor* o simplemente descartarse variables que tengan una correlación fuerte.