Introducción a la Regresión Lineal Bayesiana

Parte I: Teorema de Bayes, Prior y Likelihood

library(tidyverse)
library(gridExtra)
theme_set(theme_bw())

1 Objetivos de la clase

Al finalizar esta sesión el estudiante será capaz de:

• Comprender la diferencia entre la inferencia clásica y la inferencia bayesiana.
• Interpretar el significado del Prior.
• Entender qué representa el Likelihood.
• Aplicar el Teorema de Bayes a un problema de regresión.
• Visualizar cómo cambia nuestra incertidumbre cuando llegan nuevos datos.

---

2 ¿Por qué Bayes?

En econometría tradicional normalmente estimamos un modelo

\[ Y=\beta_0+\beta_1X+\varepsilon \]

y obtenemos un único estimador

\[ \hat{\beta}_1=0.28 \]

En el enfoque clásico este número resume toda la información.

Sin embargo, surge una pregunta importante.

¿Qué tan seguros estamos de que el verdadero parámetro vale exactamente 0.28?

La estadística bayesiana responde de manera distinta.

No busca un único valor para el parámetro.

Busca una distribución de probabilidad.

En lugar de preguntar

¿Cuál es el valor de \(\beta\)?

preguntamos

\[ P(\beta|Datos) \]

es decir,

¿Cuál es la distribución del parámetro después de observar los datos?

Esta diferencia conceptual será el eje de toda la clase.

---

3 Caso de estudio

Supongamos que trabajamos para un banco.

Queremos estudiar cómo el ingreso mensual de un cliente influye sobre el saldo promedio de su tarjeta de crédito.

La variable respuesta será

\[ Y=\text{Saldo de la tarjeta} \]

La variable explicativa será

\[ X=\text{Ingreso mensual} \]

---

4 Importación de datos

Utilizaremos datos reales de clientes bancarios.

library(readxl)

datos_raw <- read_excel(
  "C:/Cesar/Cesar Lectures/Bayesian_econometrics/Bayesian_regress/Datos_Bayes_Banking.xlsx"
)

datos <- datos_raw |>
  dplyr::select(Ingreso = Ingreso_Mensual,
                Saldo   = Saldo_Tarjeta)

head(datos)

## # A tibble: 6 × 2
##   Ingreso Saldo
##     <dbl> <dbl>
## 1   3067.  737.
## 2   2524.  822.
## 3   4202. 1170.
## 4   3456. 1078.
## 5   2420.  711.
## 6   2619.  670.

---

Veamos la relación entre ambas variables.

ggplot(datos,aes(Ingreso,Saldo))+

geom_point(size=2,
color="steelblue")+

geom_smooth(method="lm",
se=FALSE,
color="red")+

labs(title="Ingreso vs Saldo",
subtitle="Datos Bayes Banking")

Observamos una relación positiva.

En econometría clásica ajustaríamos inmediatamente una regresión.

En Bayes aún no.

Antes debemos hablar de nuestras creencias.

---

5 El pensamiento bayesiano

Supongamos que antes de observar cualquier cliente un analista del banco afirma

“Creo que el efecto del ingreso sobre el saldo es cercano a 0.20.”

Naturalmente esta afirmación tiene incertidumbre.

No sabemos si realmente será 0.20.

Podría ser

• 0.18
• 0.22
• 0.25

Por tanto representamos esa incertidumbre mediante una distribución.

---

6 El Prior

Definimos el prior para el parámetro

\[ \beta_1 \]

como

\[ \beta_1\sim N(0.20,0.08^2) \]

Esto significa

• media = 0.20

• desviación estándar = 0.08

---

7 Visualizando el Prior

beta <- seq(-0.1,0.5,length=500)

prior <- dnorm(beta,
mean=.20,
sd=.08)

prior_df <- data.frame(beta,prior)

ggplot(prior_df,
aes(beta,prior))+

geom_line(size=1.2,
color="blue")+

labs(title="Distribución Prior",
subtitle="Nuestra creencia antes de observar datos",
x=expression(beta[1]),
y="Densidad")

Observe que todavía no hemos utilizado ninguna observación.

Toda esta distribución proviene únicamente del conocimiento previo.

---

8 Interpretación

La figura anterior nos dice que

• valores cercanos a 0.20 son muy probables

• valores cercanos a 0.40 son poco probables

• valores negativos prácticamente no son creíbles

Esta es nuestra información inicial.

---

9 ¿Qué sucede cuando llegan datos?

Supongamos que llegan los primeros clientes.

head(datos,8)

## # A tibble: 8 × 2
##   Ingreso Saldo
##     <dbl> <dbl>
## 1   3067.  737.
## 2   2524.  822.
## 3   4202. 1170.
## 4   3456. 1078.
## 5   2420.  711.
## 6   2619.  670.
## 7   2803.  853.
## 8   1846.  682.

Ahora aparece un concepto nuevo.

El Likelihood.

---

10 ¿Qué es el Likelihood?

En una regresión suponemos

\[ Y_i\sim N(\mu_i,\sigma^2) \]

donde

\[ \mu_i=\beta_0+\beta_1X_i \]

Si proponemos un determinado valor de

\[ \beta_1 \]

podemos calcular

¿Qué tan compatibles son los datos con ese valor?

Eso es exactamente el Likelihood.

---

11 Intuición

Supongamos que proponemos

\[ \beta_1=0.05 \]

La recta probablemente ajuste muy mal.

Por tanto

Likelihood pequeño.

Ahora proponemos

\[ \beta_1=0.25 \]

La recta ajusta bastante bien.

Likelihood grande.

En consecuencia,

el likelihood mide

Qué tan bien explica un parámetro los datos observados.

---

12 La función de verosimilitud

Matemáticamente

\[ L(\beta)= \prod_{i=1}^{n} f(y_i|\beta) \]

No debemos interpretar esta expresión como una probabilidad del parámetro.

Los parámetros son constantes desconocidas.

Lo que hace el likelihood es responder

Si \(\beta\) tomara este valor, ¿qué tan probables serían los datos observados?

---

13 Comparando distintos valores de beta

Supongamos tres posibles pendientes

par(mfrow=c(1,3))

plot(datos$Ingreso,
datos$Saldo,
main=expression(beta==0.10),
pch=19)

abline(200,.10,col=2,lwd=2)

plot(datos$Ingreso,
datos$Saldo,
main=expression(beta==0.25),
pch=19)

abline(200,.25,col=2,lwd=2)

plot(datos$Ingreso,
datos$Saldo,
main=expression(beta==0.40),
pch=19)

abline(200,.40,col=2,lwd=2)

par(mfrow=c(1,1))

Visualmente podemos notar que una pendiente cercana a 0.25 explica mejor los datos.

Esto significa que su likelihood será mayor.

---

14 Idea fundamental

Hasta este momento tenemos dos piezas.

14.1 Prior

Lo que creíamos antes de observar datos.

14.2 Likelihood

Lo que dicen los datos.

La siguiente pregunta será

¿Cómo combinamos ambas fuentes de información?

La respuesta será mediante el Teorema de Bayes.

---

15 Resumen

En esta primera parte aprendimos que

• La inferencia bayesiana trabaja con distribuciones.
• El prior representa nuestras creencias iniciales.
• El likelihood resume la información contenida en los datos.
• Ambos elementos serán combinados mediante el Teorema de Bayes para obtener la distribución posterior.

En la siguiente parte construiremos paso a paso la distribución posterior y veremos cómo cada nuevo conjunto de datos actualiza nuestras creencias sobre el parámetro de interés.

16 Parte II. Del Likelihood a la Distribución Posterior

Hasta este momento tenemos dos ingredientes fundamentales.

1. Una creencia inicial sobre el parámetro (Prior).
2. Los datos observados.

La pregunta natural es:

¿Cómo combinamos ambas fuentes de información?

La respuesta la proporciona el Teorema de Bayes.

---

17 El Teorema de Bayes

En términos generales,

\[ P(\theta|Datos)= \frac{P(Datos|\theta)P(\theta)} {P(Datos)} \]

En regresión lineal escribimos

\[ P(\beta|Y)= \frac{P(Y|\beta)P(\beta)} {P(Y)} \]

Cada término tiene un significado diferente.

Componente	Interpretación
\(P(\beta)\)	Prior
\(P(Y|\beta)\)	Likelihood
\(P(\beta|Y)\)	Posterior
\(P(Y)\)	Constante de normalización

Observe que

\[ Posterior = Prior \times Likelihood \]

salvo por una constante que garantiza que la distribución resultante integre uno.

---

18 ¿Qué representa el denominador?

Muchos estudiantes piensan que

\[ P(Y) \]

es una parte “misteriosa”.

En realidad solamente sirve para normalizar la distribución.

Es decir,

si multiplicamos

Prior × Likelihood

obtendremos una curva que todavía no es una distribución de probabilidad.

Necesitamos dividir entre

\[ P(Y) \]

para que el área bajo la curva sea exactamente igual a uno.

Por ello normalmente escribimos

\[ P(\beta|Y) \propto P(Y|\beta)P(\beta) \]

donde el símbolo

\[ \propto \]

se lee

“proporcional a”.

---

19 Construyendo el Likelihood

Supongamos que conocemos

\[ \beta_0=200 \]

y queremos evaluar distintos valores para

\[ \beta_1. \]

Nuestro modelo es

\[ Y_i=200+\beta_1X_i+\varepsilon_i \]

donde

\[ \varepsilon_i\sim N(0,\sigma^2). \]

Supongamos además que

\[ \sigma=180. \]

---

20 Probabilidad de una observación

Tomemos únicamente el primer cliente.

datos[1,]

## # A tibble: 1 × 2
##   Ingreso Saldo
##     <dbl> <dbl>
## 1   3067.  737.

Supongamos que proponemos

\[ \beta_1=0.20 \]

El modelo predice

beta <- 0.20

mu <- 200 + beta*datos$Ingreso[1]

mu

## [1] 813.436

La diferencia entre el valor observado y el esperado es el error.

error <- datos$Saldo[1]-mu

error

## [1] -76.486

Si los errores siguen una distribución Normal,

la probabilidad de observar ese dato viene dada por

densidad <- dnorm(datos$Saldo[1],
                  mean=mu,
                  sd=180)

densidad

## [1] 0.002025022

Esta cantidad representa

la densidad de probabilidad asociada al primer cliente.

---

21 Visualizando una distribución Normal

x <- seq(mu-600,mu+600,length=500)

normal_df <- data.frame(
x=x,
y=dnorm(x,mu,180)
)

ggplot(normal_df,
aes(x,y))+

geom_line(color="steelblue",
linewidth=1.2)+

geom_vline(xintercept=datos$Saldo[1],
linetype=2,
color="red")+

labs(title="Distribución del primer cliente",
subtitle="La línea roja corresponde al valor observado",
x="Saldo",
y="Densidad")

Mientras más cerca esté el dato de la media,

mayor será su densidad.

---

22 Del primer cliente a todos los clientes

Con una sola observación el likelihood sería

\[ L(\beta)=f(y_1|\beta) \]

Con dos observaciones

\[ L(\beta)= f(y_1|\beta) \times f(y_2|\beta) \]

Con

100 observaciones

\[ L(\beta)= \prod_{i=1}^{100} f(y_i|\beta) \]

Esta multiplicación incorpora la información contenida en todos los clientes.

---

23 Calculando el Likelihood para distintos valores de beta

Construiremos una pequeña función.

likelihood <- function(beta){

mu <- 200 + beta*datos$Ingreso

sum(dnorm(datos$Saldo,
          mean=mu,
          sd=180,
          log=TRUE))

}

¿Por qué utilizamos logaritmos?

Porque multiplicar cientos de probabilidades produce números extremadamente pequeños.

En su lugar calculamos

\[ \log(L)= \sum \log(f(y_i)) \]

Esto evita problemas numéricos.

---

24 Evaluando muchos posibles valores

beta_grid <- seq(0.05,
                 0.40,
                 length=250)

logLik <- sapply(beta_grid,
                 likelihood)

lik_df <- data.frame(beta_grid,
                     logLik)

---

25 El Log-Likelihood

ggplot(lik_df,
aes(beta_grid,
logLik))+

geom_line(linewidth=1.2,
color="darkgreen")+

labs(title="Log-Likelihood",
x=expression(beta[1]),
y="Log-Likelihood")

Observe que

existe un valor de

\[ \beta \]

que hace máximo el likelihood.

Ese valor coincide aproximadamente con la pendiente que obtendríamos mediante OLS.

Hasta aquí todavía no estamos haciendo Bayes.

Simplemente estamos evaluando

qué tan compatibles son los datos con distintos valores del parámetro.

---

26 Recuperando el Likelihood

Como el gráfico anterior está en escala logarítmica,

podemos recuperar una versión proporcional del likelihood.

lik_df$Likelihood <-
exp(logLik-max(logLik))

Restamos el máximo únicamente para evitar problemas numéricos.

---

27 Visualizando el Likelihood

ggplot(lik_df,
aes(beta_grid,
Likelihood))+

geom_line(color="firebrick",
linewidth=1.2)+

labs(title="Likelihood",
x=expression(beta[1]),
y="Likelihood")

Ahora sí observamos claramente

qué valores del parámetro son más compatibles con los datos.

---

28 Comparando Prior y Likelihood

prior <- dnorm(beta_grid,
               mean=.20,
               sd=.08)

comparacion <- data.frame(

beta=beta_grid,

Prior=prior/

max(prior),

Likelihood=lik_df$Likelihood/

max(lik_df$Likelihood)

)

ggplot(comparacion)+

geom_line(aes(beta,Prior,
color="Prior"),
linewidth=1.2)+

geom_line(aes(beta,Likelihood,
color="Likelihood"),
linewidth=1.2)+

scale_color_manual(values=c("blue","red"))+

labs(title="Prior vs Likelihood",
y="Escala relativa",
color="")

El Prior representa la información previa.

El Likelihood representa la evidencia contenida en los datos.

La Posterior combinará ambas fuentes de información.

---

29 Construcción intuitiva de la Posterior

Recordemos

\[ Posterior \propto Prior \times Likelihood. \]

Podemos construir una versión proporcional de la posterior multiplicando ambas curvas.

posterior <- prior*lik_df$Likelihood

posterior <- posterior/sum(posterior)

posterior_df <- data.frame(
beta=beta_grid,
Posterior=posterior)

---

30 La Distribución Posterior

ggplot()+

geom_line(data=comparacion,
aes(beta,Prior/max(Prior),
color="Prior"),
linewidth=1)+

geom_line(data=comparacion,
aes(beta,Likelihood/max(Likelihood),
color="Likelihood"),
linewidth=1)+

geom_line(data=posterior_df,
aes(beta,
Posterior/max(Posterior),
color="Posterior"),
linewidth=1.4)+

scale_color_manual(values=c(
Prior="blue",
Likelihood="red",
Posterior="darkgreen"))+

labs(title="Actualización Bayesiana",
x=expression(beta[1]),
y="Escala relativa",
color="")

Observe cómo la distribución posterior se concentra alrededor de los valores que:

• eran razonables según el conocimiento previo; y
• explican bien los datos observados.

Esta es la esencia del enfoque bayesiano.

---

31 Idea clave

En econometría clásica obtenemos

\[ \hat{\beta}=0.24 \]

En Bayes obtenemos una distribución completa.

Esa distribución responde preguntas como

• ¿Cuál es el valor más probable?
• ¿Qué incertidumbre existe?
• ¿Qué tan probable es que el parámetro sea mayor que 0.20?
• ¿Qué intervalo contiene el 95% de la probabilidad?

Estas preguntas serán el tema de la siguiente parte.

---

32 Resumen

En esta sección aprendimos que:

• El Likelihood mide qué tan compatibles son los datos con cada valor posible del parámetro.
• El Prior resume nuestras creencias antes de observar información.
• El Teorema de Bayes combina ambas fuentes de información.
• La Posterior representa nuestra nueva creencia después de observar los datos.
• En el enfoque bayesiano, el objetivo ya no es estimar un único número, sino una distribución completa para los parámetros del modelo.

---

33 Interpretando la Distribución Posterior

Hasta ahora hemos obtenido una distribución posterior para el parámetro \(\beta_1\).

Sin embargo, surge una pregunta importante:

¿Qué información podemos extraer de esta distribución?

A diferencia de la econometría clásica, donde obtenemos un único estimador, la inferencia bayesiana nos permite resumir la distribución mediante varias medidas.

Las más importantes son:

• La media posterior.
• La moda posterior.
• El intervalo de credibilidad del 95%.

---

34 Media Posterior

La media posterior corresponde al valor esperado del parámetro después de haber observado los datos.

Matemáticamente,

\[ E(\beta|Y)= \int \beta P(\beta|Y)d\beta \]

En nuestro ejemplo aproximaremos esta cantidad numéricamente.

posterior_mean <-
sum(
posterior_df$beta *
posterior_df$Posterior
)

posterior_mean

## [1] 0.2165078

Observe que este valor resume el “centro de gravedad” de la distribución posterior.

---

35 Moda Posterior (MAP)

Otra medida importante es el valor más probable del parámetro.

Este estimador recibe el nombre de

Maximum A Posteriori (MAP).

Matemáticamente,

\[ \hat{\beta}_{MAP} = \arg\max_{\beta} P(\beta|Y) \]

En R,

posterior_mode <-
posterior_df$beta[
which.max(posterior_df$Posterior)
]

posterior_mode

## [1] 0.2158635

Observe que el MAP corresponde al punto donde la distribución posterior alcanza su máxima altura.

---

36 Intervalo de Credibilidad

Una gran ventaja del enfoque bayesiano es que podemos responder preguntas como

¿En qué rango se encuentra el parámetro con un 95% de probabilidad?

Para aproximarlo generaremos muestras de la distribución posterior.

set.seed(123)

beta_sim <- sample(
posterior_df$beta,
size = 10000,
replace = TRUE,
prob = posterior_df$Posterior
)

IC95 <- quantile(
beta_sim,
c(0.025,0.975)
)

IC95

##      2.5%     97.5% 
## 0.2060241 0.2271084

Este resultado corresponde al intervalo de credibilidad del 95%.

En términos bayesianos podemos afirmar que

Existe aproximadamente un 95% de probabilidad de que el verdadero parámetro se encuentre dentro de este intervalo.

Esta interpretación es mucho más intuitiva que la utilizada en la inferencia frecuentista.

---

37 Comparación con OLS

Veamos ahora cuál habría sido el resultado utilizando la regresión clásica.

modelo_ols <-
lm(Saldo ~ Ingreso,
data = datos)

summary(modelo_ols)

## 
## Call:
## lm(formula = Saldo ~ Ingreso, data = datos)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -408.60 -105.88   15.42  107.53  479.46 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)    
## (Intercept) 310.64774   86.82609   3.578  0.00054 ***
## Ingreso       0.18300    0.02694   6.794 8.53e-10 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 179.5 on 98 degrees of freedom
## Multiple R-squared:  0.3202, Adjusted R-squared:  0.3133 
## F-statistic: 46.16 on 1 and 98 DF,  p-value: 8.526e-10

La pendiente estimada es

coef(modelo_ols)[2]

##   Ingreso 
## 0.1830031

Observe que OLS produce un único valor.

Bayes, por el contrario, produce una distribución completa.

---

38 Comparación gráfica

ggplot() +

geom_line(
data = comparacion,
aes(beta,
Prior/max(Prior),
color="Prior"),
linewidth=1.1)+

geom_line(
data = comparacion,
aes(beta,
Likelihood/max(Likelihood),
color="Likelihood"),
linewidth=1.1)+

geom_line(
data = posterior_df,
aes(beta,
Posterior/max(Posterior),
color="Posterior"),
linewidth=1.4)+

geom_vline(
xintercept = coef(modelo_ols)[2],
linetype=2,
linewidth=1,
color="black")+

annotate(
"text",
x = coef(modelo_ols)[2],
y = 1.05,
label = "OLS",
angle = 90,
size = 4)+

scale_color_manual(
values=c(
Prior="blue",
Likelihood="red",
Posterior="darkgreen"))+

labs(
title="Comparación entre OLS y Bayes",
subtitle="La línea punteada corresponde al estimador OLS",
x=expression(beta[1]),
y="Escala relativa",
color="")

Observe cómo

• el Prior representa nuestro conocimiento inicial;
• el Likelihood resume la evidencia contenida en los datos;
• la Posterior combina ambas fuentes de información;
• el estimador OLS corresponde a un único punto dentro de esa distribución.

---

39 Conclusiones

En este ejemplo observamos que la inferencia bayesiana no produce únicamente un estimador puntual.

Produce una distribución completa del parámetro.

Esta distribución permite responder preguntas como

• ¿Cuál es el valor más probable del parámetro?
• ¿Qué incertidumbre existe?
• ¿Cuál es el intervalo de credibilidad del 95%?
• ¿Cómo cambia nuestra creencia cuando incorporamos nuevos datos?

Estas preguntas constituyen la principal diferencia entre la econometría clásica y la econometría bayesiana.