Procurando uma distribuição, talvez incomum, consistente com dois pontos de dados e restrições de especialistas?

8

Estou tentando declarar uma distribuição anterior para uma meta-análise bayesiana.

Eu tenho as seguintes informações sobre uma variável aleatória:

  1. Duas observações: 3.0, 3.6
  2. um cientista que estuda a variável me disse que , e que valores tão altos quanto 6 têm probabilidade diferente de zero.P(X<2)=P(X>8)=0

Eu usei a seguinte abordagem para otimização (o modo de log-N = :eμσ2)

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

texto alternativo

Na figura, você pode ver a distribuição que isso retorna, mas eu gostaria de encontrar algo mais parecido com as linhas vermelhas que desenhei.

Isso fornece a mesma distribuição modelada usando o lognormal, gama ou normal, e resulta em uma distribuição com e , ou seja:P ( X = 6 ) < 0,01P(X=5)<0.05P(X=6)<0.01

 plnorm(c(5,6), v$estimate[1],v$estimate[2])

Alguém pode sugerir alternativas? Eu preferiria ficar com uma distribuição única do que com uma mistura.

Obrigado!

David LeBauer
fonte
1
Na sua pergunta, você tem e . Você quer dizer: e ?? P(X=2)=0P(X=8)=0P(X<2)=0P(X>8)=0
M. Tibbits 12/12
1
@ David Depois de ler esta pergunta muitas e muitas vezes ainda não a entendo. Você está usando as duas observações e a opinião do cientista para estimar um anterior para uma análise bayesiana? O seu prévio será baseado apenas na opinião do cientista e você deseja atualizá-lo com as observações? Qual é a distinção entre "não surpreendente" e ter probabilidade zero?
whuber
1
@ David Eu aprecio o seu cuidado em obter as informações. Mas atribuir valores específicos de 0,05 e 0,01 a essas probabilidades é questionável. Isso não é culpa sua; é assim que as coisas são. Não podemos esperar que as pessoas identifiquem tão bem as probabilidades (é para isso que servem os dados). Talvez você gostaria de representar aqueles com hiperpriors :-)?
whuber
1
@whuber Obrigado por acompanhar. Agora entendo o conceito da solução de entropia máxima que você sugeriu anteriormente, embora isso exija algum aprendizado antes que eu queira aplicá-la. Sua alternativa, maximizando a variação, parece suficientemente consistente com o meu objetivo e direta para eu implementar. Obrigado novamente.
David LeBauer
1
@ David Esteja ciente de que um problema com a maximização da variação é que a solução será uma distribuição discreta. Espero que concentre 0,01 de probabilidade em X = 8, 0,04 em X = 6, 157/300 em X = 2 e o restante (32/75) em X = 5. (Essa variação é igual a 2,59.)
whuber

Respostas:

5

Se, dada uma resposta ao meu comentário acima, você deseja limitar o intervalo da distribuição, por que não ajustar simplesmente uma distribuição Beta na qual você redimensiona o intervalo de unidade? Em outras palavras, se você souber que o parâmetro de interesse deve estar entre , por que não definir . Onde eu centralizei o intervalo em zero, dividido pela largura para que Y tenha um intervalo de 1 e, em seguida, adicionei volta para que o intervalo de Y seja . (Você pode pensar de qualquer maneira: diretamente de ou de[2,8]Y=X56+12=X2612[0,1][2,8][0,1][2,8][12,12][0,1], mas achei que o último poderia ser mais fácil a princípio).

Então, com dois pontos de dados, você poderia ajustar um beta posterior com um beta uniforme antes?

M. Tibbits
fonte
2

E a distribuição Kumaraswamy , que possui o seguinte pdf:

a > 0 b > 0 0 < x < 1

f(x;a,b)=abxa1(1xa)b1
para , , . Essa distribuição pode ser redimensionada para ter o suporte necessário.a>0b>00<x<1
csgillespie
fonte
2

Como a distribuição log-normal possui dois parâmetros, você não pode ajustá-la satisfatoriamente a três restrições que não se ajustam naturalmente a ela. Com quantis extremos de 2,5 e 7,5, o modo é ~ 4 e não há muito que você possa fazer sobre isso. Como a escala dos erros para ae bé muito menor que para c, um deles será praticamente ignorado durante a otimização.

Para um melhor ajuste, você pode escolher uma distribuição de três parâmetros, por exemplo, a distribuição gama generalizada (implementada no VGAMpacote) ou adicionar um parâmetro shift à distribuição lognormal (ou gama, ...).

Como última observação, como a distribuição que você está procurando claramente não é simétrica, a média das duas observações fornecidas não é o valor correto para o modo. Eu maximizaria a soma das densidades em 3,0 e 3,6, mantendo os quantis extremos em 2,5 e 7,5 - isso é possível se você tiver três parâmetros.

Aniko
fonte
1

Você também pode tentar a distribuição triangular. Para ajustar isso, você basicamente especifica um limite inferior (esse seria X = 2), um limite superior (esse seria X = 8) e um valor "mais provável". A página da wikepedia http://en.wikipedia.org/wiki/Triangular_distribution possui mais informações sobre esta distribuição. Se não houver muita fé no valor "mais provável" (como parece ser, antes da observação de qualquer dado), pode ser uma boa ideia colocar uma distribuição prévia não informativa e usar os dois dados pontos para estimar esse valor. Um bom é o prior de jeffrey, que para esse problema seria p (c) = 1 / (pi * sqrt ((c-2) * (c-8))), onde "c" é o "valor mais provável "(consistente com a notação da Wikipedia).

Dado isso anteriormente, é possível calcular a distribuição posterior de c analiticamente ou via simulação. A forma analítica da probabilidade não é particularmente agradável, portanto a simulação parece ser mais atraente. Este exemplo é particularmente adequado para amostragem por rejeição (consulte a página da wiki para obter uma descrição geral da amostragem por rejeição), porque a probabilidade máxima é de 1/3 ^ n, independentemente do valor de c, que fornece o "limite superior". Então, você gera um "candidato" a partir do prior de jeffrey (chame de c_i) e avalia a probabilidade desse candidato L (x1, .., xn | c_i) e divide pela probabilidade maximizada para dar (3 ^ n ) * L (x1, .., xn | c_i). Você gera uma variável aleatória U (0,1) e, se u for menor que (3 ^ n) * L (x1, .., xn | c_i), aceite c_i como um valor amostrado posterior, caso contrário, jogue fora c_i e comece de novo. Repita esse processo até ter amostras aceitas suficientes (100, 500, 1.000 ou mais, dependendo da precisão que você deseja). Em seguida, faça a média da amostra de qualquer função de c em que você esteja interessado (a probabilidade de uma nova observação é um candidato óbvio para sua aplicação).

Uma alternativa para aceitar-rejeitar é usar o valor da probabilidade como um peso (e não gerar u) e, em seguida, continuar com as médias ponderadas usando todos os candidatos, em vez de médias não ponderadas com os candidatos aceitos.

probabilityislogic
fonte