7.6 Estimación y predicción
Según Neter, Wasserman, y Kutner (1996), los modelos de regresión lineales son utilizado esencialmente con dos fines. Uno es tratar de explicar la variable de interés en términos de covariables que pueden encontrarse en la encuesta o en registros administrativos, censos, etc. Adicionalmente, también son usados para predecir valores de la variable en estudio ya sea dentro del intervalo de valores recogidos en la muestra o por fuera de dicho intervalo. Lo primero se ha abordado a lo largo de todo el capítulo y lo segundo se obtiene de la siguiente manera:
\[ \hat{E}(y_{i}\mid\boldsymbol{x}_{obs,i})=\boldsymbol{x}_{obs,i}\hat{\boldsymbol{\beta}} \]
De manera explícita, en el modelo que ejemplifica este capítulo, la expresión para las predicciones sería la siguiente:
\[ \hat{E}(y_{i}\mid\boldsymbol{x}_{obs,i})=\hat{\beta}_{0}+\hat{\beta}_{1}x_{1i} \]
La varianza de la estimación se calcula de la siguiente manera:
\[ var(t(\hat{E}(t(y_{i}\mid x_{obs,i})) = x'_{obs,i}cov(t(\hat{\beta})x{}_{obs,i} \]
A continuación, se presenta cómo se realiza la estimación del valor esperado utilizando R. Por ende, primero se estiman los parámetros del modelo:
| term | estimate | std.error | statistic | p.value |
|---|---|---|---|---|
| (Intercept) | 72.6993 | 68.1486 | 1.067 | 0.2883 |
| Expenditure | 1.1886 | 0.2158 | 5.508 | 0.0000 |
| ZoneUrban | 70.9745 | 42.1315 | 1.685 | 0.0948 |
| SexMale | 20.9344 | 15.9909 | 1.309 | 0.1931 |
| Age2 | 0.0082 | 0.0056 | 1.460 | 0.1469 |
Por lo anterior, la estimación del valor esperado o predicción queda:
\[ \hat{E}(y_{i}\mid\boldsymbol{x}_{obs,i})=103.136+1.263x_{1i} \] Para calcular la varianza de la estimación, primero se deben obtener las varianzas de la estimación de los parámetros:
| (Intercept) | Expenditure | ZoneUrban | SexMale | Age2 | |
|---|---|---|---|---|---|
| (Intercept) | 4644.2350 | -13.5405 | 579.8621 | 310.2784 | -0.2086 |
| Expenditure | -13.5405 | 0.0466 | -4.5712 | -0.4802 | 0.0006 |
| ZoneUrban | 579.8621 | -4.5712 | 1775.0665 | -148.0881 | -0.0647 |
| SexMale | 310.2784 | -0.4802 | -148.0881 | 255.7098 | -0.0034 |
| Age2 | -0.2086 | 0.0006 | -0.0647 | -0.0034 | 0.0000 |
Ahora bien, se procede a realizar los cálculos como lo indica la expresión mostrada anteriormente:
xobs <- model.matrix(fit_svy) %>%
data.frame() %>% slice(1) %>% as.matrix()
cov_beta <- vcov(fit_svy) %>% as.matrix()
as.numeric(xobs %*% cov_beta %*% t(xobs))## [1] 1999Si el objetivo ahora es calcular el intervalo de confianza para la predicción se utiliza la siguiente ecuación:
\[ \boldsymbol{x}_{obs,i}\hat{\beta}\pm t_{(t(1-\frac{\alpha}{2},n-p)}\sqrt{var(t(\hat{E}(t(y_{i}\mid\boldsymbol{x}_{obs,i}))} \]
Para realizar los cálculos en R, se utiliza la función confint y predict como sigue:
pred <- data.frame(predict(fit_svy, type = "link"))
pred_IC <- data.frame(confint(predict(fit_svy, type = "link")))
colnames(pred_IC) <- c("Lim_Inf", "Lim_Sup")
pred_ICAhora, de manera gráfica las estimaciones y los intervalos de confianza para algunas observaciones se verían de la siguiente forma:
pred <- bind_cols(pred, pred_IC)
pred$Expenditure <- encuesta$Expenditure
pred %>% slice(1:6L)
pd <- position_dodge(width = 0.2)
ggplot(pred %>% slice(1:100L),
aes(x = Expenditure , y = link)) +
geom_errorbar(aes(ymin = Lim_Inf,
ymax = Lim_Sup),
width = .1,
linetype = 1) +
geom_point(size = 2, position = pd) +
theme_bw()Por último, si el interés es hacer una predicción fuera del rango de valores que fue capturado en la muestra, se debe tener en cuenta que la expresión de la varianza cambia ligeramente. En ese caso, la varianza para la predicción se hace siguiendo la siguiente ecuación:
\[ var(t(\hat{E}(t(y_{i}\mid\boldsymbol{x}_{obs,i}))=\boldsymbol{x}_{obs,i}^{t}cov(t(\boldsymbol{\beta})\boldsymbol{x}_{obs,i} + \hat{\sigma}^2_{yx} \]
Para ejemplificarlos en R, supongamos que se desea predecir el valor del ingreso para un gasto de 2600 (valor que no fue observado en la muestra). El primer paso es crear un conjunto de datos nuevos, así:
Luego, se construye la matriz de observaciones y se calcula la varianza como sigue:
x_noObs = matrix(c(1, 1600, 1, 1, 40 ^ 2), nrow = 1)
as.numeric(sqrt(x_noObs %*% cov_beta %*% t(x_noObs)))## [1] 264.7Por último, el intervalo de confianza sigue la siguiente ecuación:
\[ \boldsymbol{x}_{obs,i}\hat{\beta}\pm t_{(t(1-\frac{\alpha}{2},n-p)}\sqrt{var(t(\hat{E}(t(y_{i}\mid\boldsymbol{x}_{obs,i}))+\hat{\sigma}_{yx}^{2}} \]
En R la predicción se realiza usando la función predict sobre el objeto fit_svy de la siguiente manera. Esta salida provee la predicción puntual y su correspondiente error estándar, el cual es idéntico al calculado anteriormente.
## link SE
## 1 2080 265Finalmente, el intervalo de confianza se calcula usando la función confint, así:
| 2.5 % | 97.5 % |
|---|---|
| 1561 | 2598 |