Ingreso vs Educación: Intervalos de confianza y pruebas de hipótesis

Primera parte

Considerar el conjunto de datos de desarrollo humano 2010 de las entidades del país. Construir intervalos de confianza $95\%$ para los parámetros del modelo RLS de logaritmo del ingreso contra años promedio de escolaridad.

Primero cargamos los datos

file_url <- 'http://sigma.iimas.unam.mx/jsantibanez/Cursos/Regresion/2017_2/dh.csv'
idh.ent <- read.csv(file_url)
library(dplyr)
idh10 <- filter(idh.ent, Year == 2010)
x <- log(idh10$ING)
y <- idh10$APE

Estimamos los parámetros del modelo

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

## [1] -21.04291

b1 <- sxy / sxx; b1

## [1] 3.018977

n <- length(x)
s2 <- sum((y - b0 - b1*x)**2) / (n - 2); s2

## [1] 0.1861986

Graficamos los datos y el modelo ajustado

par(mar = c(4.1, 4.1, 1, 1))
plot(x, y, xlab = 'log(Ingreso)', ylab = 'Educación',
     pch = 16, axes = F, xlim = c(8.75, 10.5), ylim = c(5.5, 10.75))
axis(1, seq(8.75, 10.5, by = 0.25))
axis(2, seq(5.5, 10.75, by = 1), las = 2)
box()
abline(a = b0, b = b1, col = 'red2')

Estimamos las varianzas de los estimadores

vb0 <- (1/n + mean(x)**2/sxx) * s2; vb0

## [1] 6.013039

vb1 <- s2 / sxx; vb1

## [1] 0.06423522

Calculamos los intervalos de confianza 95%

tt <- qt(0.025, n-2, lower.tail = F)
IC_b0 <- b0 + c(-1, 1) * tt * sqrt(vb0); IC_b0

## [1] -26.05087 -16.03495

IC_b1 <- b1 + c(-1, 1) * tt * sqrt(vb1); IC_b1

## [1] 2.501369 3.536584

cl <- qchisq(0.025, n-2)
cu <- qchisq(0.975, n-2)
IC_s2 <- (n-2) * s2 / c(cu, cl); IC_s2

## [1] 0.1189027 0.3326803

Segunda parte

¿Hay algún efecto del ingreso en el valor promedio de los años promedio de escolaridad?
¿Sería mejor ajustar un modelo sin intercepto?
¿La varianza del modelo es menor a 0.25?
Calcular un intervalo de confianza para la media de los años promedio de escolaridad para $40,000.00 de ingreso.
Calcular un intervalo de predicción de los años promedio de escolaridad para $40,000.00 de ingreso.

Efecto del ingreso en la educación

Lo que se pide es contrastar las hipótesis

\[\begin{equation*} H_0: \beta_1 = 0 \qquad \text{vs.} \qquad H_1: \beta_1 \neq 0 \end{equation*}\]

Calculamos el estadístico de prueba y lo comparamos con el valor crítico

tp <- abs(b1/sqrt(vb1)); tp

## [1] 11.91168

tt

## [1] 2.042272

Como el estadístico de prueba es mayor que el valor crítico, se rechaza $H_0$ con una significancia $\alpha = 0.05$. Los datos indican que hay algún efecto del ingreso en la educación.

Modelo con o sin intercepto

Lo que se pide es contrastar las hipótesis

\[\begin{equation*} H_0: \beta_0 = 0 \qquad \text{vs.} \qquad H_1: \beta_0 \neq 0 \end{equation*}\]

Calculamos el estadístico de prueba y lo comparamos con el valor crítico

tp <- abs(b0/sqrt(vb0)); tp

## [1] 8.581412

tt

## [1] 2.042272

Como el estadístico de prueba es mayor que el valor crítico, se rechaza $H_0$ con una significancia $\alpha = 0.05$. Los datos indican que el modelo de regresión debe incluir a $\beta_0$.

Hipótesis sobre la varianza del modelo

Lo que se pide es contrastar las hipótesis

\[\begin{equation*} H_0: \sigma^2 \leq 0.25 \text{vs.} \qquad H_1: \sigma^2 > 0.25 \end{equation*}\]

Calculamos el estadístico de prueba y lo comparamos con el valor crítico

tp <- (n - 2)*s2 / 0.25; tp

## [1] 22.34384

qchisq(0.95, n - 2)

## [1] 43.77297

El estadístico de prueba es menor que el valor critico, no se rechaza $H_0$ con una significancia $\alpha = 0.05$. Los datos indican que la varianza no es mayor a 0.25.

Intervalo de confianza de la respuesta media

Primero debemos obtener la estimación puntual y su varianza estimada

m40 <- b0 + b1 * log(40000); m40

## [1] 10.94808

vm40 <- (1/n + (log(40000) - mean(x))^2/sxx) * s2; vm40

## [1] 0.06091164

Ahora se calcula el intervalo de confianza

IC_m40 <- m40 + c(-1, 1) * tt * sqrt(vm40); IC_m40

## [1] 10.44405 11.45212

Intervalo de predicción

Ya tenemos la predicción puntual y hace falta la varianza de predicción

vp40 <- (1 + 1/n + (log(40000) - mean(x))^2/sxx) * s2; vp40

## [1] 0.2471103

Ahora se calcula el intervalo de confianza

IC_p40 <- m40 + c(-1, 1) * tt * sqrt(vp40); IC_p40

## [1]  9.932866 11.963301

ajuste <- lm(y ~ x)
summary(ajuste)

## 
## Call:
## lm(formula = y ~ x)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.87394 -0.23529 -0.03492  0.19134  1.27516 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept) -21.0429     2.4521  -8.581 1.42e-09 ***
## x             3.0190     0.2534  11.912 6.71e-13 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.4315 on 30 degrees of freedom
## Multiple R-squared:  0.8255, Adjusted R-squared:  0.8197 
## F-statistic: 141.9 on 1 and 30 DF,  p-value: 6.71e-13

Tercera parte

Calcular $SC_{TC}$, $SC_{reg}$ y $SC_{error}$ y verificicar la relación \[\begin{equation*} SC_{TC} = SC_{reg} + SC_{error} \end{equation*}\]
Calcular el coeficiente de determinación del modelo y verificar la relación con el coeficiente de correlación de Pearson.
Ajustar el modelo con la función lm y comparar los resultados con los anteriores.

Descomposición de la suma de cuadrados

Calculamos

SCTC  <- sum((y - mean(y))^2); SCTC

## [1] 32.00535

yh <- b0 + b1*x
SCREG <- sxx*b1^2; SCREG

## [1] 26.41939

sum((yh - mean(y))^2)

## [1] 26.41939

SCERR <- sum((y - yh)^2); SCERR

## [1] 5.585959

SCTC - (SCREG + SCERR)

## [1] 7.105427e-15

Coeficiente de determinación

Calculamos

R2 <- SCREG / SCTC; R2

## [1] 0.8254679

r  <- cor(x, y, method = 'pearson'); r

## [1] 0.9085527

R2 - r^2

## [1] -2.220446e-16

Ajuste del modelo con `lm`

Ajustamos el modelo y mostramos los resultados

ajuste <- lm(y ~ x)
summary(ajuste)

## 
## Call:
## lm(formula = y ~ x)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.87394 -0.23529 -0.03492  0.19134  1.27516 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept) -21.0429     2.4521  -8.581 1.42e-09 ***
## x             3.0190     0.2534  11.912 6.71e-13 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.4315 on 30 degrees of freedom
## Multiple R-squared:  0.8255, Adjusted R-squared:  0.8197 
## F-statistic: 141.9 on 1 and 30 DF,  p-value: 6.71e-13

Bandas de confianza para la respuesta media

ic_l <- function(x_n, alfa){
  (b0 + b1*x_n) - 2*qf(1 - alfa, 2, n-2)*sqrt(1/n + (x_n - mean(x))^2/sxx)*s2
}

ic_u <- function(x_n, alfa){
  (b0 + b1*x_n) + 2*qf(1 - alfa, 2, n-2)*sqrt(1/n + (x_n - mean(x))^2/sxx)*s2
}

par(mar = c(4.1, 4.1, 1, 1))
plot(x, y, xlab = 'log(Ingreso)', ylab = 'Educación',
     pch = 16, axes = F, xlim = c(8.75, 10.5), ylim = c(5.5, 10.75))
axis(1, seq(8.75, 10.5, by = 0.25))
axis(2, seq(5.5, 10.75, by = 1), las = 2)
box()
abline(a = b0, b = b1, col = 'red2')
curve(ic_l(x, 0.05), add = T, col = 'blue2')
curve(ic_u(x, 0.05),  add = T, col ='blue2')
curve(ic_l(x, 0.10), add = T, col = 'green4')
curve(ic_u(x, 0.10),  add = T, col ='green4')
points(mean(x), mean(y), pch = 16, cex = 2, col = 'steelblue3')
legend('topleft', c('0.05', '0.10'), lwd = 2, col = c('blue2', 'green4'))

Ingreso vs Educación: Intervalos de confianza y pruebas de hipótesis

Javier Santibáñez

9 de marzo de 2017

Primera parte

Primero cargamos los datos

Estimamos los parámetros del modelo

Graficamos los datos y el modelo ajustado

Estimamos las varianzas de los estimadores

Calculamos los intervalos de confianza 95%

Segunda parte

Efecto del ingreso en la educación

Modelo con o sin intercepto

Hipótesis sobre la varianza del modelo

Intervalo de confianza de la respuesta media

Intervalo de predicción

Tercera parte

Descomposición de la suma de cuadrados

Coeficiente de determinación

Ajuste del modelo con `lm`

Bandas de confianza para la respuesta media

Ingreso vs Educación: Intervalos de confianza y pruebas de hipótesis

Javier Santibáñez

9 de marzo de 2017

Primera parte

Primero cargamos los datos

Estimamos los parámetros del modelo

Graficamos los datos y el modelo ajustado

Estimamos las varianzas de los estimadores

Calculamos los intervalos de confianza 95%

Segunda parte

Efecto del ingreso en la educación

Modelo con o sin intercepto

Hipótesis sobre la varianza del modelo

Intervalo de confianza de la respuesta media

Intervalo de predicción

Tercera parte

Descomposición de la suma de cuadrados

Coeficiente de determinación

Ajuste del modelo con lm

Bandas de confianza para la respuesta media

Ajuste del modelo con `lm`