Si es Bayes, es bueno…

October 23, 2011
By

This post was kindly contributed by Farmacovigilancia nerd2!=nerd - go there to comment and to read the full post.


Los métodos de inferencia vistos hasta el momento formanparte de lo que se da en llamar “inferencia frecuentista”. El modelo subyacentees que, verificando ciertos supuestos, es posible ajustar una distribución(normal u otra) a mis datos para inferir conclusiones acerca de la población.
La inferencia bayesiana tiene otra visión. En general partede que la distribución “real” de la población es algo imposible de conocer. Enlugar de tomar una distribución como la verdadera, se toma lo que conoce “apriori” el investigador para elegir una distribución que se ajuste a ello ytrabajar inicialmente bajo ese supuesto. Luego se toma una muestra de lapoblación y con los datos obtenidos de ella se realiza un “ajuste” de ladistribución empírica, llevando a una distribución que podría ser diferente ala anterior. Es lo que se conoce como distribución “a posteriori”.
O sea que, proponemos un modelo de acuerdo a lo queconocemos, tomamos evidencia para contrastarla contra el modelo y corregimos elmodelo de acuerdo a lo aprendido.

Para tener un ejemplo intentaré estimar cuál es laproporción de estudiantes de una universidad que duerme al menos 8 horas.
Este ejemplo está tomado por completo del libro “Bayesian computation with R” de Jim Albert.

Supongamos que un investigador está interesado en conoceresa proporción. El verdadero valor de la proporción le es desconocido.
De acuerdo a la información de otras publicaciones yexperiencia previa, piensa 2 cosas:
·        La mejor predicción que puede haceres de 30%
·        Es muy probable que se halle pordebajo del 50%

Construimos un rango de valores posibles para la proporción.Llamaremos “p” a la proporción. En R escribimos:

>p=seq(0.05,0.95,by=0.1)
> p
 [1] 0.05 0.15 0.25 0.35 0.45 0.55 0.65 0.750.85 0.95

Esta instrucción construye un rango entre 0.05 (5%) y 0.95(95%) que aumenta de a 0.1 (10%). Es un rango de valores posibles de p que nosservirá para construir una distribución a priori.
Ahora a cada uno de estos valores le otorgamos un puntajeentre 0 y 10 de acuerdo a lo que pensamos (esta escala es arbitraria). Estospuntajes son llamados “pesos”.

> pesos = c(1, 5.2, 8, 7.2, 4.6, 2.1, 0.7, 0.1, 0, 0)
> pesos
[1] 1.0 5.2 8.0 7.2 4.6 2.1 0.7 0.1 0.0 0.0

Y convertimos esos pesos a un formato de probabilidad. Parahacerlo, dividimos cada peso por la suma total de los mismos. A estasprobabilidades se las llama “priors” o probabilidades “a priori”.

> prior = pesos/sum(pesos)
> prior

Esto nos permite crear una distribución que indica unaprobabilidad para cada valor de p. Para visualizarla podemos usar un gráfico:

> plot(p, prior, type = "h")

El gráfico obtenido nos muestra nuestra distribución “apriori”. Ahora la tarea será obtener una muestra de estudiantes. Lo hacemos yobtenemos datos de 27 estudiantes. Luego de revisar las encuestas notamos quehay 11 estudiantes que duermen más de 8 horas. Esto nos da un número de“éxitos” (11) y “fracasos” (16). Colocamos esos valores dentro de un objeto:

> data=c(11,16)
> data
[1] 11 16

El uso de la instrucción “c” es para concatenar variosvalores y en este caso, asignarlos a un solo objeto.
Ahora hay que cargar un paquete de R llamado “LearnBayes”.Por supuesto, tiene que estar previamente instalado.

> library("LearnBayes")

Si no estuviese instalado nos daría un error diciendo que noencuentra el paquete. Procedemos como en el post que explica sobre instalación de paquetes de R y lo instalamos. Luego volvemos a tipear la instrucciónanterior.

Luego de cargado el paquete escribimos el siguiente comando:

> post = pdisc(p, prior, data)
> post

Esto nos da las probabilidades a posteriori. Es bastanteclaro que son diferentes a las probabilidades a priori.
Para visualizarlas en forma comparativa escribimos una seriede instrucciones que no explicaré en detalle, pero nos permiten ver una gráficasobre la otra. En caso de que a alguien le interese la explicación detallada melo puede pedir. No obstante, se puede obtener ayuda de R escribiendo “?”seguido de la instrucción, o help(“instrucción”).

> library(lattice)
> PRIOR=data.frame("prior",p,prior)
> POST=data.frame("posterior",p,post)
> names(PRIOR)=c("Type","P","Probability")
> names(POST)=c("Type","P","Probability")
> data=rbind(PRIOR,POST)
>xyplot(Probability~P|Type,data=data,layout=c(1,2),type="h",lwd=3,col="black")

Esto nos muestra que las 2 distribuciones, la a prioriconstruida con las suposiciones del investigador y la a posteriori, calculadausando la a priori y corrigiendo de acuerdo a las observaciones, sondiferentes.

Existen funciones matemáticas para calcular lasdistribuciones a posteriori a partir de las a priori. Estas funciones sondiferentes en cada situación y mi objetivo en este post no es ser extenso enéste tema.
Baste decir, que con esta metodología es posible estimardistribuciones pero también parámetros, como medias, sd y otros.

Con estos conocimientos podemos ir, ahora sí, a continuar con el aprendizaje del paquete PhViD para detección de señales.



Tags: , ,

Comments are closed.