Estoy tratando de establecer una distribución previa para un metanálisis bayesiano.
Tengo la siguiente información sobre una variable aleatoria:
- Dos observaciones: 3.0, 3.6
- Un científico que estudia la variable me ha dicho que , y que valores tan altos como 6 tienen una probabilidad distinta de cero.
He utilizado el siguiente enfoque para la optimización (el modo de log-N = :
prior <- function(parms, x, alpha) {
a <- abs(plnorm(x[1], parms[1], parms[2]) - (alpha/2))
b <- abs(plnorm(x[2], parms[1], parms[2]) - (1-alpha/2))
mode <- exp(parms[1] - parms[2]^2)
c <- abs(mode-3.3)
return(a + b + c)
}
v = nlm(prior,c(log(3.3),0.14),alpha=0.05,x=c(2.5,7.5))
x <- seq(1,10,0.1)
plot(x, dlnorm(x, v$estimate[1], v$estimate[2]))
abline(v=c(2.5,7.5), lty=2) #95%CI
En la figura, puede ver la distribución que esto devuelve, pero me gustaría encontrar algo más parecido a las líneas rojas que he dibujado.
Esto proporciona la misma distribución conformada utilizando lognormal, gamma o normal, y da como resultado una distribución con y , es decir:
plnorm(c(5,6), v$estimate[1],v$estimate[2])
¿Alguien puede sugerir alternativas? Preferiría seguir con una sola distribución en lugar de una mezcla.
¡Gracias!
r
distributions
probability
bayesian
optimization
David LeBauer
fuente
fuente
Respuestas:
Si, dada una respuesta a mi comentario anterior, desea vincular el rango de la distribución, ¿por qué no simplemente ajustar una distribución Beta donde reescala al intervalo de la unidad? En otras palabras, si sabe que el parámetro de interés debe estar entre , entonces ¿por qué no definir . Donde primero centré el intervalo en cero, dividido por el ancho para que Y tenga un rango de 1, y luego agregué para que el rango de Y sea . (Puede pensarlo de cualquier manera: directamente desde o desde[2,8] Y=X−56+12=X−26 12 [0,1] [2,8]→[0,1] [2,8]→[−12,12]→[0,1] , pero pensé que esto último podría ser más fácil al principio).
Entonces, con dos puntos de datos, ¿podría colocar una beta posterior con una beta uniforme anterior?
fuente
¿Qué pasa con la distribución Kumaraswamy , que tiene el siguiente pdf:
a > 0 b > 0 0 < x < 1
fuente
Dado que la distribución logarítmica normal tiene dos parámetros, no puede ajustarla satisfactoriamente a tres restricciones que no la ajustan naturalmente. Con cuantiles extremos de 2.5 y 7.5, el modo es ~ 4, y no hay mucho que pueda hacer al respecto. Dado que la escala de los errores para
a
yb
es mucho menor que parac
, uno de ellos será ignorado durante la optimización.Para un mejor ajuste, puede elegir una distribución de tres parámetros, por ejemplo, la distribución gamma generalizada (implementada en el
VGAM
paquete), o agregar un parámetro de desplazamiento a la distribución lognormal (o gamma, ...).Como última nota, dado que la distribución que está buscando claramente no es simétrica, el promedio de las dos observaciones dadas no es el valor correcto para el modo. Maximizaría la suma de las densidades en 3.0 y 3.6 mientras mantenía los cuantiles extremos en 2.5 y 7.5; esto es posible si tiene tres parámetros.
fuente
También puedes probar la distribución triangular. Para ajustar esto, básicamente especifica un límite inferior (sería X = 2), un límite superior (sería X = 8) y un valor "más probable". La página de wikepedia http://en.wikipedia.org/wiki/Triangular_distribution tiene más información sobre esta distribución. Si no hay mucha fe en el valor "más probable" (como parece ser, antes de observar cualquier dato), puede ser una buena idea colocar una distribución previa no informativa y luego usar los dos datos puntos para estimar este valor. Una buena es la anterior de Jeffrey, que para este problema sería p (c) = 1 / (pi * sqrt ((c-2) * (c-8))), donde "c" es el "valor más probable "(consistente con la notación de wikipedia).
Dado esto antes, puede calcular la distribución posterior de c analíticamente o mediante simulación. La forma analítica de la probabilidad no es particularmente agradable, por lo que la simulación parece ser más atractiva. Este ejemplo es particularmente adecuado para el muestreo de rechazo (consulte la página wiki para obtener una descripción general del muestreo de rechazo), porque la probabilidad maximizada es 1/3 ^ n independientemente del valor de c, que proporciona el "límite superior". Entonces, genera un "candidato" a partir del previo de Jeffrey (llámelo c_i), y luego evalúa la probabilidad en este candidato L (x1, .., xn | c_i), y divide por la probabilidad maximizada, para dar (3 ^ n ) * L (x1, .., xn | c_i). Luego genera una variable aleatoria U (0,1), y si u es menor que (3 ^ n) * L (x1, .., xn | c_i), entonces acepta c_i como un valor de muestra posterior, de lo contrario deseche c_i y empezar de nuevo Repita este proceso hasta que tenga suficientes muestras aceptadas (100, 500, 1,000 o más, según la precisión que desee). Luego, simplemente tome el promedio de la muestra de cualquier función de c que le interese (la probabilidad de una nueva observación es un candidato obvio para su aplicación).
Una alternativa para aceptar-rechazar es usar el valor de la probabilidad como una ponderación (y no generar la u), y luego proceder a tomar promedios ponderados con todos los candidatos, en lugar de promedios no ponderados con los candidatos aceptados
fuente