Ingreso vs Educación

Descripción

Para este ejemplo consideramos los datos sobre desarrollo humano de las entidades del país en los años 2008, 2010 y 2014. El archivo se puede descargar desde aquí. La fuente original de los datos esta acá.

Primero hacermos una exploración del conjunto de datos

file_url <- 'http://sigma.iimas.unam.mx/jsantibanez/Cursos/Regresion/2017_2/dh.csv'
idh.ent <- read.csv(file_url)
str(idh.ent)

## 'data.frame':    96 obs. of  10 variables:
##  $ Ent : Factor w/ 32 levels "Ags","BC","BCS",..: 1 2 3 4 8 9 7 6 5 10 ...
##  $ Year: int  2008 2008 2008 2008 2008 2008 2008 2008 2008 2008 ...
##  $ Esp : num  75.1 73.1 75.3 74.4 75.1 ...
##  $ APE : num  8.7 8.89 8.99 7.78 8.97 ...
##  $ AEE : num  11.8 12.3 13.4 11.9 12.4 ...
##  $ ING : num  17507 21346 25716 16127 17266 ...
##  $ IS  : num  0.847 0.817 0.851 0.837 0.848 ...
##  $ IE  : num  0.618 0.638 0.672 0.588 0.643 ...
##  $ II  : num  0.78 0.81 0.838 0.768 0.778 ...
##  $ IDH : num  0.742 0.75 0.782 0.723 0.751 ...

El conjunto de datos contiene 10 variables, una alfanumérica y nueve numéricas. La descripción de las variables es la siguiente:

• Ent: entidad federativa
• Year: Año de referencia
• Esp: Esperanza de vida al nacer
• APE: Años promedio de escolaridad
• AEE: Años esperados de escolaridad
• ING: Ingreso nacional bruto (per cápita y en USD PPC)
• IS: Índice de salud
• IE: Índice de educación
• II: Índice de ingreso
• IDH: Índice de desarrollo humano

Nos interesa explorar las relaciones entre ingreso y eduación, en específico, nos interesa saber cómo el nivel educativo (APE) es explicado por el ingreso (ING). En este caso usaremos el lograritmo del ingreso.

library(dplyr)
idh10 <- filter(idh.ent, Year == 2010)
x <- log(idh10$ING, 2)
y <- idh10$APE
n <- length(x)
range(x)

## [1] 12.91585 14.83605

Primero graficamos las variables de interés:

Recordamos que los estimadores de minímos cuadrados son:

\[\begin{align*} \hat{\beta}_0 &= \bar{y}_n - \frac{S_{xy}}{S_{xx}} \bar{x}_n \\ \hat{\beta}_1 &= \frac{S_{xy}}{S_{xx}} \\ \hat{\sigma}^2 &= \frac{1}{n-2}\sum_{i=1}^n (y_i - \hat{\beta}_0 - \hat{\beta}_1x_i)^2 \end{align*}\]

sxx <- sum((x-mean(x))**2)
sxy <- sum((x-mean(x)) * (y - mean(y)))

b0.h <- mean(y) - sxy/sxx * mean(x); b0.h

## [1] -21.04291

b1.h <- sxy / sxx; b1.h

## [1] 2.092595

s2.h <- sum((y - b0.h - b1.h * x)**2) / (32 - 2); s2.h

## [1] 0.1861986

Entonces \(\hat{\beta}_0 =\) -21.043, \(\hat{\beta}_1 =\) 2.093 y \(\hat{\sigma}^2 =\) 0.186. Podemos graficar la recta ajustada con la nube de puntos.

También podemos estimar las varianzas de nuestros estimadores. Recordemos:

\[\begin{align*} V(\hat{\beta}_0) &= \sigma^2 \left(\frac{1}{n} + \frac{\bar{x}_n^2}{S_{xx}} \right) \\ V(\hat{\beta}_1) &= \frac{\sigma^2}{S_{xx}} \end{align*}\]

Podemos reemplazar \(\sigma^2\) por \(\hat{\sigma}^2\), entonces:

v_b0.h <- (1/n + mean(x)**2/sxx)*s2.h; v_b0.h

## [1] 6.013039

v_b1.h <- s2.h / sxx; v_b1.h

## [1] 0.03086201

De lo anterior se sigue que \(V(\hat{\beta}_0) =\) 6.013 y \(V(\hat{\beta}_1) =\) 0.031.

Intepretación de los resultados