setwd("~/ESTADISTICA1011")

Teorema de Bayes

El teorema de Bayes, en la teorĆ­a de la probabilidad, es una proposiciĆ³n planteada por el matemĆ”tico inglĆ©s Thomas Bayes (1702-1761)1 y publicada pĆ³stumamente en 1763, que expresa:

En tƩrminos mƔs generales y menos matemƔticos, el teorema de Bayes es de enorme relevancia puesto que vincula la probabilidad de A dado B con la probabilidad de B dado A.

Es decir, por ejemplo, que sabiendo la probabilidad de tener un dolor de cabeza dado que se tiene gripe, se podrĆ­a saber (si se tiene algĆŗn dato mĆ”s), la probabilidad de tener gripe si se tiene un dolor de cabeza. Muestra este sencillo ejemplo la alta relevancia del teorema en cuestiĆ³n para la ciencia en todas sus ramas, puesto que tiene vinculaciĆ³n Ć­ntima con la comprensiĆ³n de la probabilidad de aspectos causales dados los efectos observados.

la probabilidad condicional de un evento aleatorio A dado B en terminos de la distribucion de probabilidad condicional del evento B dado A y la distribucion de probabilidad marginal de solo A.

Formula de Bayes

Con base en la definicion de probabilidad condicionada se obtiene la Formula de Bayes, tambien conocida como Regla de Bayes:

ESTADISTICA1011

ESTADISTICA1011

Analisis Bayesianos

Los analisis Bayesianos son similares a los que vimos en en el sentido que dependen explicitamente de modelos probabilisticos para los datos. Es decir, tenemos que definir un modelo de datos. La gran diferencia es que con Bayes, podemos obtener distribuciones de probabilidades para todas las cantidades no observadas, incluyendo parƃĀ”metros, valores perdidos o nuevas observaciones que todavia no hemos hecho. De esta manera, los analisis Bayesianos nos permiten cuantificar incertidumbre y armar modelos realistas que tienen en cuenta por ejemplo observaciones imperfectas.

Como vimos en la teorica, la regla de Bayes planteada en tƃĀ©rminos de datos y parametros es:

\[ p(\boldsymbol{\theta} \lvert \boldsymbol{y}) = \frac{p(\boldsymbol{y} \lvert \boldsymbol{\theta}) p(\boldsymbol{\theta)}}{\int p(\boldsymbol{y} \lvert \boldsymbol{\theta)} p(\boldsymbol{\theta)} d \boldsymbol{\theta} } \] Likehood : verosimilitud

Es decir que la probabilidad posterior de los parƃĀ”metros ƎĀø dado que observamos los datos y es igual al likelihood multiplicado por las previas y dividido por la probabilidad total de los datos. La funcion de likelihood nos da la probabilidad de observar los datos condicional al valor de los parƃĀ”metros p(y|ƎĀø). La previa de los parametros p(ƎĀø)

refleja los posibles valores de los parametros de acuerdo con nuestras creencias previas, o los resultados de estudios anteriores, o lo que nos parece que tiene sentido para el sistema de estudio (en definitiva, en base a informacion previa). Finalmente, la probabilidad total de los datos se obtiene integrando la funcion de lilkelihood sobre los posibles valores de los parƃĀ”metros que define la previa. Como veremos mas adelante, los analisis Bayesianos combinados con metodos numericos permiten analizar modelos con muchos parametros, niveles de variabilidad y variables pero primero vamos a empezar por casos simples donde podemos calcular las posteriores directamente.

  • Ejercicio 1: Analisis bayesiano

Imaginen que queremos estudiar la remociĆ³n de frutos en 30 plantas. En cada una de las plantas marcamos 20 frutos y contamos cuĆ”ntos son removidos por dispersores luego de un tiempo fijo. Si suponemos que un buen modelo para este tipo de datos es una distribuciĆ³n Binomial con una probabilidad de Ć©xito fija hacemos:

set.seed(1234)
nobs <- 30 #NĆŗmero de observaciones (plantas)
frutos <- rep(20,nobs)
p_rem <- 0.2 #probabilidad de remociĆ³n por fruto
removidos <- rbinom(nobs, size=frutos,prob=p_rem)
removidos
##  [1] 2 4 4 4 6 5 0 3 5 4 5 4 3 7 3 6 3 3 2 3 3 3 2 1 3 6 4 7 6 1

El modelo de datos (cuĆ”ntos frutos son removidos) es una Binomial con nĆŗmero de pruebas (la cantidad de frutos disponibles) conocido. Para hacer un anĆ”lisis Bayesiano de estos datos, tenemos que definir una previa para la probabilidad de Ć©xito (p_rem).

Esa previa tiene que tomar valores continuos entre 0 y 1. Una opciĆ³n serĆ­a una distribuciĆ³n uniforme con esos lĆ­mites, pero si usamos una distribuciĆ³n Beta, es posible obtener un resultado analĆ­tico para la posterior.

En este caso, la posterior es otra distribuciĆ³n Beta pero con sus parĆ”metros actualizados en base a las observaciones. Se dice entonces que la distribuciĆ³n Beta es la conjugada de la Binomial. Si la previa de la tasa de remociĆ³n por fruto es una distribuciĆ³n Beta con parĆ”metros Ī± y Ī², actualizamos los valores de Ī± y Ī² en base a la cantidad de Ć©xitos y fracasos obervados.

La posterior de la tasa de remociĆ³n por fruto es entonces una Beta con Ī±=āˆ‘y, Ī²=āˆ‘(nāˆ’y) donde y representa a los frutos removidos de los n disponibles. Veamos como hacer esto en R.

#previas
alpha <- 1
beta <- 1
alpha_p <- alpha + sum (removidos)
beta_p <- alpha + sum(frutos-removidos)

Para obtener el valor esperado de una distribucion Beta hacermos

alpha_p / (alpha_p + beta_p)
## [1] 0.1877076

Eso nos da un estimador puntual de la probabilidad de remociĆ³n por fruto p_rem. Para tener una medida de incertidumbre alrededor de este valor, podemos ver los cuantiles de la posterior:

qbeta(c(0.025, 0.975),alpha_p,beta_p     )
## [1] 0.1575462 0.2198340

TambiĆ©n podemos graficar la distribuciĆ³n posterior y compararla con la previa para ver cuĆ”nto aprendimos haciendo el anĆ”lisis.

op <- par(cex.lab = 1.5, font.lab = 2, cex.axis = 1.3, las = 1, bty = "n")
curve(dbeta(x, alpha + sum(removidos), beta + sum(frutos - removidos)), lwd = 2, 
    ylab = "Densidad de probabilidad", xlab = "Probabilidad de remociĆ³n")
curve(dbeta(x, 1, 1), lwd = 2, col = "gray", add = TRUE)
text(0.6, 2.5, "previa")
text(0.35, 12, "posterior")

par(op)

ĀæCĆ³mo sabemos si estas estimaciones tienen sentido para nuestros datos? En este caso, la pregunta es trivial porque conocemos cĆ³mo se generaron los datos, pero cuando trabajamos con datos de verdad, el modelo de datos es un supuesto y tenemos que ver si ese supuesto tiene sentido.

Una opciĆ³n para contestar esa pregunta es hacer simulaciones a partir de la posterior y compararlas con los datos.

nreps <- 10000
vals <- 0:20 #posibles valores de remociĆ³n 
res <- matrix(NA, nreps, length(vals)-1) #Matriz para resultados
p_sim <- rbeta(nreps,alpha_p,beta_p) #muestra aleatoria de la posterior

for (i in 1:nreps) {
    tmp <- rbinom(nobs, frutos, p_sim[i])
    res[i, ] <- hist(tmp, right = FALSE, breaks = vals, plot = FALSE)$density
}

plot(table(removidos)/nobs, xlim = c(0, 10), ylim = c(0, 0.6), ylab = "frecuencia", 
    type = "p", pch = 19)
library(coda)
## Warning: package 'coda' was built under R version 4.0.3

podemos ver que las observaciones (puntos), estƔn dentro del intervalo de credibilidad para muestras de una Binomial con probabilidad de Ʃxito dada por la posterior de p_rem.

Veamos ahora un anƔlisis con datos del mundo real.

url <- "https://sites.google.com/site/modelosydatos/quintral.txt"
quintral <- read.table(url, header = TRUE)

# previas
alpha <- 1
beta <- 1

alpha_p <- alpha + sum(quintral$Removidos)
beta_p <- beta + sum(quintral$Frutos - quintral$Removidos)

op <- par(cex.lab = 1.5, font.lab = 2, cex.axis = 1.3, las = 1, bty = "n")
curve(dbeta(x, alpha_p, beta_p), lwd = 2, ylab = "Densidad de probabilidad", 
    xlab = "Probabilidad de remociĆ³n")
curve(dbeta(x, 1, 1), lwd = 2, col = "gray", add = TRUE)
text(0.2, 2.5, "previa")
text(0.65, 20, "posterior")

Vemos que la posterior es bastante acotada, implicando que tendrĆ­amos bastante seguridad respecto a que la probabilidad de remociĆ³n de un fruto de quintral es cercana al 60%.

Veamos quƩ pasa si simulamos datos de este modelo y contrastamos con las observaciones.

Asignacion:

-> Utilice los principios de bayes obtenidos hasta ahora para realizar un ejercicio de aplicacion de estos a un caso de su area del conocimiento

** En la academia de MateMovil, la probabilidad de que a un alumno seleccionado al azar le guste el helado es del 60 %, mientras que la probabilidad de que a un alumno le guste la torta es del 36 %. AdemƔs, se sabe que la probabilidad de que a un alumno le guste la torta dado que le gusta el helado es del 40 %. Calcular la probabilidad de que a un alumno le guste el helado, dado que le gusta la torta.**

Primero definimos los 2 eventos con los que vamos a trabajar:

  • h: que a un alumno le guste el helado.
  • t: que a un alumno le guste la torta.

Tenemos los siguientes datos:

  • P(h) = 0,6.
  • P(t) = 0,36.
  • P(t|h) = 0,4.

Nos piden calcular P(h|t).

ESTADISTICA1011

ESTADISTICA1011