1. Tarea Binomial

Intuición antes del código

Dada una \(Bin(20, \theta)\), una distribución Beta(α, β) suele usarse como a priori. Su interpretación es muy intuitiva:

  • Si α y β son pequeños (ej: α=β=1, a priori uniforme), la posterior queda fuertemente influenciada por la verosimilitud.
  • Si α y β son grandes (ej: α=β=50), la prior está muy concentrada y domina incluso con un número moderado de observaciones → la posterior se mueve poco. (a priori arrastra a la posterior)
  • Si α > β, la prior “empuja” a la posterior hacia valores grandes de θ (cercanos a 1).
  • Si β > α, la prior favorece valores bajos de θ.

Lo que tenemos:

  • Datos: \(X \sim \text{Binomial}(n, \theta)\), donde \(\theta\) es la probabilidad de éxito.

  • Prior: \(\theta \sim \text{Beta}(\alpha, \beta)\).

  • Posterior: \(\theta | X \sim \text{Beta}(\alpha + x, \, \beta + n - x)\)

Código en R

library(ggplot2)
library(dplyr)
library(patchwork) #bloques de graficas

theta <- seq(0.001, 0.999, length.out = 500)

# Datos observados (ejemplo1)
n <- 20     # número de ensayos
x <- 14     # éxitos observados

graficar_priorios <- function(alpha, beta, n=20, x){
  ## a priori ##
  priori <- dbeta(theta, alpha, beta)
  ## Verosimilitud ## L(θ∣x)∝θ^{x}(1−θ)^{n−x}.
  veros <- dbinom(x, size = n, prob = theta)
  veros <- veros/max(veros)*max(priori)
  ## a posteriori ##
  post <- dbeta(theta, alpha + x, beta + n - x)
  # Dataframe unificado ## 
  df <- data.frame(
    theta = rep(theta, 3),
    densidad = c(priori, veros, post),
    tipo = rep(c("Prior", "Verosimilitud",
                 "Posterior"), 
               each=length(theta))
  )
    ggplot(df, aes(x=theta, y=densidad, color=tipo)) +
    geom_line(size=1.2) +
    labs(title = paste0("α=", alpha, ", β=", beta, 
                        " | n=", n, ", x=", x),
         y = "Densidad", x = expression(theta)) +
    theme_minimal() +
    scale_color_manual(values=c("blue", "red",
                                "darkgreen"))
}

Ejemplos comparativos

Prior uniforme (\(α=1\), \(β=1\))

graficar_priorios(1, 1, n, x)

-> Aquí la posterior queda guiada casi únicamente por los datos.

Prior sesgada a la izquierda (α=2, β=10)

graficar_priorios(2, 10, n, x)

-> La posterior favorece valores bajos de θ.

Prior sesgada a la derecha (α=10, β=2)

graficar_priorios(10, 2, n, x)

-> La posterior favorece valores altos de θ.

Prior muy informativa (α=50, β=50)

graficar_priorios(50, 50, n, x)

-> Muy informativa (α y β grandes): la prior domina y la posterior apenas se mueve.

Intuición al variar los tamaños de muestra

  • ¿Con \(n\) pequeño qué pasa?: la prior influye mucho.
  • ¿Con \(n\) grande? la posterior se concentra fuertemente alrededor de la proporción observada \(x/n\), y la prior importa menos.

Ejemplo 1: Prior uniforme (α=1, β=1) y distintos n

comparar_n(alpha=1, beta=1, tasas=0.56, n_vals=c(10, 50, 200)) + patchwork::plot_layout(ncol = 2)

Amedida que \(n\) crece, la distribución posterior se concentra alrededor de la proporción observada \(\frac{x}{n} \approx 0.56\).

Ejemplo 2: Prior sesgada (α=10, β=2)

comparar_n(alpha=10, beta=2, tasas=0.56, n_vals=c(10, 50, 200)) + patchwork::plot_layout(ncol = 2)

Con \(n=10\), la prior sesgada domina; con \(n=200\), los datos imponen la forma de la posterior.

Ejemplo 3: Prior muy informativa (α=50, β=50)

comparar_n(alpha=50, beta=50, tasas=0.56, n_vals=c(10, 50, 200)) + patchwork::plot_layout(ncol = 2)

Incluso con \(n=50\), la prior fuerte todavía frena el movimiento de la posterior; con \(n=200\), los datos ya empujan más fuerte.

Prior de Jeffreys

En el modelo binomial:

\[ X \sim \text{Binomial}(n, \theta) \]

La prior de Jeffreys es:

\[ \pi(\theta) \propto \frac{1}{\sqrt{\theta (1-\theta)}} \;\;\;\;\;\Rightarrow\;\;\;\;\; \theta \sim \text{Beta}\left(\tfrac{1}{2}, \tfrac{1}{2}\right) \]

Es una prior no informativa, y suele preferirse frente a la prior uniforme:

\[ \theta \sim \text{Beta}(1,1) \]

La posterior resultante es:

\[ \theta \mid X \sim \text{Beta}\left(\tfrac{1}{2} + x, \; \tfrac{1}{2} + n - x\right) \]

Prior no informativa de Jeffreys (α=0.5, β=0.5)

comparar_n(alpha = 0.5, beta = 0.5, tasas = 0.56, n_vals = c(10, 50, 200))

Con \(n=10\), la forma de la prior de Jeffreys todavía se percibe en la posterior.

Con \(n = 50\), ya se observa que los datos comienzan a dominar, la posterior se concentra alrededor de la proporción observada \(x/n \approx 0.56\), aunque la curvatura en los bordes de la prior todavía influye ligeramente.

Con \(n=200\), la influencia de la prior prácticamente desaparece y los datos mandan.

La prior de Jeffreys es menos “plana” que la uniforme (α=1, β=1), asigna más peso a valores extremos (cercanos a 0 y 1). Con \(n=50\) y \(x=28\), la posterior se concentra igual que en el caso uniforme, pero se nota que la forma inicial de la prior es diferente.

Ejemplo explicativo: un dado cargado

Imaginemos que queremos estimar la probabilidad \(\theta\) de que un dado cargado saque “6”.

  • Sin información previa: usamos una prior uniforme \(\theta \sim \text{Beta}(1,1)\).
  • Creencia de que casi nunca sale 6: prior sesgada \(\theta \sim \text{Beta}(2,10)\).
  • Creencia de que sale 6 muy seguido: prior sesgada \(\theta \sim \text{Beta}(10,2)\).
  • Prior muy informativa: \(\theta \sim \text{Beta}(50,50)\), como si ya tuviéramos muchas tiradas previas en mente.

Ahora tiramos el dado \(n=20\) veces y observamos \(x=14\) veces el valor “6”.

graficar_priorios(1, 1, n=20, x=14)   # Uniforme

graficar_priorios(2, 10, n=20, x=14)  # Sesgada a valores bajos

graficar_priorios(10, 2, n=20, x=14)  # Sesgada a valores altos

graficar_priorios(50, 50, n=20, x=14) # Muy informativa

Por tanto, se ilustra cómo la posterior combina la creencia previa y los datos observados:
- Con prior uniforme, los datos dominan.
- Con prior sesgada, los datos van corrigiendo la creencia previa.
- Con prior informativa, la posterior apenas se mueve al inicio.

\((X \mid \theta, \phi) \sim \mathcal{N}(\theta, \phi^{-1})\)

2. Tarea Normal-Gamma-Normal

Normal-Normal

Inicamos con la parte a) la cual es asignar una distribución normal para la apriori. De lo cual tenemos:

Distribución Datos: \(X \sim \mathcal{N}(\theta, \phi^{-1})\) tal que \(\phi=\frac{1}{\sigma^2}\)

Verosimilitud: \(f(x)=\frac{\phi^{1/2}}{\sqrt(2\pi)} \exp\left[-\tfrac{\phi}{2}(x - \theta)^2\right]\)

Verosimilitud: \(f(x)=P(X|\theta,\phi)\varpropto (\phi^{1/2}) \exp\left[-\tfrac{\phi}{2}(x - \theta)^2\right]\)

La priori conjunta \(P(\theta, \phi)\) se construye como el producto de \(P(\theta \mid \phi)\) y \(P(\phi)\), es decir: \(P(\theta, \phi) = P(\theta \mid \phi)\, P(\phi)\).

De esta manera desarrollamos que: \(P(\theta | \phi)=\phi^{1/2}\exp\left[-\tfrac{\phi}{2}(\theta - \mu_0)^2\right]\), y además, \(\phi \sim \mathcal{G}(\alpha, \beta)\) con parametros conocidos. Por tanto, obtenemos la segunda parte \(P(\phi)=\phi^{\alpha-1}\exp \left[ -\beta\phi \right]\) donde \(\phi>0\) y \(\alpha, \beta\) positivos.

Es decir que al calcular nuestra funcion a posteriori dada por \(P(\theta,\phi|X)\varpropto P(X|\theta,\phi)*P(\theta | \phi)*P(\phi)\), resulta en la siguiente expresión:

\[ P(\theta |\phi,X)\varpropto (\phi^{1/2}) \exp\left[-\tfrac{\phi}{2}(x - \theta)^2\right] * \phi^{1/2}\exp\left[-\tfrac{\phi}{2}(\theta - \mu_0)^2\right] * \phi^{\alpha-1}\exp \left[ -\beta\phi \right] \]

Cuyo desarrollo se realizó en clase y se llegó a que la funcion a posteriori se distribuye normal como se ve a continuación:

\[ P(\theta | \phi,X) \varpropto \exp \left[ -\frac{2\phi}{2}*[\theta - \frac{x+\mu_0}{2}]^2\right] \]

en otras palabras; \((\theta | \phi,X)\sim \mathcal{N}(\frac{x+\mu_0}{2}, \frac{1}{2\sigma^2})\).

Normal-Gamma

Análogamente y razonando de con las mismas distribuciones descritas, pero, manteniendo a \(\phi\) en la verosimilitud y la a priori como nuestro parámetro de interes y además, no tenemos en cuenta el parámetro \(\theta\). Por tanto, se llega a la conlcusión de que: \[ P(\phi | \theta,X) \varpropto \phi^{(\alpha+1)-1} \exp - \left[ \frac{\phi}{2}*[(x-\theta)^2 + (\theta-\mu_0)^2]+ \phi\beta\right] \] De esta manera obtenemos que \((\phi | \theta,X)\sim \mathcal{G}(\alpha+1, \tau)\) donde \(\tau=- \left[ \frac{\phi}{2}*[(x-\theta)^2 + (\theta-\mu_0)^2]+ \phi\beta\right]\).

set.seed(357)
# Dato observado
x <- 5   # observación única
# Parámetros asumimos conocidos
mu0 <- 0
alpha <- 2
beta <- 2

Itera <- 5000
theta <- numeric(Itera)
phi <- numeric(Itera)
# Valores iniciales # Inicio G
theta[1] <- mean(x)
phi[1] <- 1  # inicializamos arbitrariamente
# Gibbs Sampling
for (t in 2:Itera) {
        # Muestrear [θ|φ,X] ~  N((x + mu0) / 2 , 1 / (2 * phi) )
   mu_post <- (x+mu0)/2
   var_post <- 1/(2*phi[t-1])
   theta[t] <- rnorm(1, mean = mu_post, sd = sqrt(var_post))
        # Muestrear [φ|θ,X]  ~ Gamma(α + 1, β + 0.5 * [(x - θ)^2 + (θ - μ0)^2])
   alpha_post <- (alpha+1)
   beta_post <- (beta + 0.5*((x-theta[t])^2+(theta[t]-mu0)^2))
   phi[t] <- rgamma(1, shape = alpha_post, rate = beta_post)
}
# Convertimos a σ^2
sigma2 <- 1/phi
par(mfrow = c(2,2))

# [θ|φ,X] ~  N((x + mu0) / 2 , 1 / (2 * phi) )
plot(theta, type = "l", col = "blue", main = expression(theta), xlab="Iteraciones" )
hist(theta, breaks = 30, col = "skyblue", main = expression(theta),xlab="Distribucion θ" )
# [φ|θ,X]  ~ Gamma(α + 1, β + 0.5 * [(x - θ)^2 + (θ - μ0)^2])
plot(sigma2, type = "l", col = "red", main = expression(sigma^2),xlab="Iteraciones" )
hist(sigma2, breaks = 30, col = "pink", main = expression(sigma^2),xlab="Distribucion σ2")

# Estimación de los parámetros desconocidos
mean(theta)
## [1] 2.482873
mean(sigma2)
## [1] 5.526018
Componente Distribución condicional Parte del modelo
\((X \mid \theta, \phi)\) \(N(\theta, \phi^{-1})\) Datos
\((\theta \mid \phi)\) \(N(\mu_0,\phi^{-1})\) Prior condicional
\((\phi)\) \(G(\alpha, \beta)\) Prior marginal
\((\theta \mid \phi, X)\) \(N\left(\frac{x+\mu_0}{2}, \frac{1}{2\phi}\right)\) Condicional completa
\((\phi \mid \theta, X)\) \(G\left(\alpha+1,; \beta+0.5[(x-\theta)^2+(\theta-\mu_0)^2]\right)\) Condicional completa

Modelo Normal-Gamma con n > 1

# Datos observados
n <- 11
x <- rnorm(n, mean = 5, sd = 2)   # ejemplo de datos simulados
# H. Parámetros de la priori 
mu0 <- 0
alpha <- 2
beta <- 2
# Configuración del Gibbs
Itera <- 5000
theta <- numeric(Itera)
phi <- numeric(Itera)

## Estimación posterior media de θ: 4.955864
## Estimación posterior media de σ²: 1.672762

Comportamiento de las posteriori según cambios en n

Para \((\theta|\phi, X)\) cuando variamos los valores de \(G(\alpha,\beta)\); la distribución priori. Obtenemos que si (\(\alpha, \beta\) grandes) entonces hablarémos de una A priori informativa, en caso contrario con (\(\alpha, \beta\) pequeños) estaríamos frente a una A priori difusa. De esta manera, se muestra como cambian las posteriori según la apriori (informativa o difusa) y según el tamaño de muestra \(n\).

\(X \sim N(\mu=5,\sigma^2=4)\)

A medida que se incrementa el tamaño de la muestra \(n\), las distribuciones a posteriori´s se vuelven más precisas.

Con una previa informativa, la posterior se ve sesgada hacia la creencia previa \(\mu_0\). Con una previa difusa, la posterior se aproxima al promedio muestral, reflejando la influencia dominante de los datos.

Para \((\phi|\theta, X)\), tenemos de manera análoga lo siguiente:

El tamaño muestral \(n\) y la elección de la previa determinan el equilibrio entre la creencia inicial y la evidencia empírica:

Con pocos datos (\(n\) pequeño) y una previa informativa, nuestra estimación posterior se mantiene cerca de la creencia inicial. Sin embargo, al recolectar más datos (\(n\) grande) o al usar una previa poco informativa, el peso de la evidencia empírica domina, atrayendo la posterior hacia el valor real del parámetro.

La forma de la posterior se vuelve más angosta (menor varianza) conforme crece \(n\). Más certeza

Observar Verosimilitud, a priori y posteriori

Otros gráficos

\(\theta\)

\(\sigma^2\)
\(\sigma^2\) y \(\theta\)