Amostragem eficiente de uma distribuição Beta com limite

Como devo coletar amostras com eficiência da seguinte distribuição?

x \sim B (α, β), x > k

$x \sim B(\alpha, \beta),\space x > k$

Se não for muito grande, a amostragem por rejeição pode ser a melhor abordagem, mas não sei como proceder quando for grande. Talvez exista alguma aproximação assintótica que possa ser aplicada? $k$ $k$

random-generation beta-distribution truncation user1502040
fonte

Não é inequivocamente claro o que você pretende lá com " ". Você quer dizer uma distribuição beta truncada (truncada à esquerda em )?

x \sim B (α, β), x > k

$x \sim B(\alpha, \beta),\space x > k$

k

$k$

Glen_b -Reinstala Monica

@Glen_b exatamente.

precisa saber é o seguinte

Para ambos os parâmetros de forma maiores que 1, a distribuição beta é côncava em log; portanto, envelopes exponenciais podem ser usados para amostragem por rejeição. Para gerar variáveis beta não truncadas, você já está amostrando a partir de distribuições exponenciais truncadas (o que é fácil de fazer), deve ser fácil adaptar esse método.

Scortchi - Reinstate Monica

Respostas:

A maneira mais simples e mais geral que se aplica a qualquer distribuição truncada (também pode ser generalizada para truncamento nos dois lados) é usar a amostragem por transformação inversa . Se é a distribuição cumulativa de juros, defina e tome $F$ $p_0 = F(k)$

U \sim U (p_{0}, 1) X = F^{- 1} (U)

$U \sim \mathcal{U}(p_0, 1) \\ X = F^{-1}(U)$

onde é uma amostra de truncada à esquerda em . A função quantil vai mapear probabilidades para amostras de . Como tomamos valores de somente da "área" que corresponde aos valores da distribuição beta da região não truncada, você amostrará apenas esses valores. $X$ $F$ $k$ $F^{-1}$ $F$ $U$

Este método é ilustrado na imagem abaixo, onde a área truncada é marcada por um retângulo cinza, os pontos em vermelho são desenhados da distribuição e depois transformados em amostras. $\mathcal{U}(p_0, 1)$ $\mathcal{B}(2, 8)$

Tim
fonte

(+1) Vale ressaltar que a função quantil não é tão facilmente avaliada.

Scortchi - Restabelece Monica

@ Scortchi Se a ou b são 1 ou, pelo menos, um número inteiro, existe uma forma não tão ruim (consulte a Wikipedia ). E em Python existe scipy.special.betainco inverso e em R existe pbeta.

Graipher

@ Graipher: Eu deveria ter dito "barato, em geral" - seria melhor evitar Newton-Raphson ou outras soluções iterativas, se possível. (BTW é qbetapara a função quantílica em R.)

Scortchi - Reinstate Monica

@ Scortchi, você está certo, mas na maioria dos casos, para computadores modernos, isso não deve ser um grande problema. Eu também recomendo essa abordagem, uma vez que está diretamente disponível na maioria dos softwares e pode ser generalizada para qualquer distribuição truncada, apenas se alguém tiver acesso à função quantil.

Tim

Sem dúvida, é bom ter um método geralmente aplicável e de fácil implementação, cujo tempo de execução não cresça com ; e para distribuições com funções quantílicas de forma fechada, por exemplo, o Weibull, deve ser o melhor possível. No entanto, desconfio que terá que ser configurado para cortar uma parte bastante grande da distribuição beta antes que ela supere os algoritmos eficientes de amostragem por rejeição que também estão disponíveis na maioria dos softwares e que dependem apenas do cálculo da densidade de probabilidade do beta.

k

$k$

k

$k$

Scortchi - Restabelece Monica

@ A resposta de Tim mostra como a amostragem por transformação inversa pode ser adaptada para distribuições truncadas, liberando o tempo de execução da dependência do limite . Eficiências adicionais podem ser obtidas evitando a avaliação numérica dispendiosa da função beta quantil e usando a amostragem por transformada inversa como parte da amostragem por rejeição. $k$

A função densidade de uma distribuição beta com parâmetros de forma & duplamente truncada em (para um pouco mais de generalidade) é $\alpha$ $\beta$ $k_1<k_2$

f (x) = \frac{x^{(α - 1)} (1 - x)^{(β - 1)}}{B (k 2, α, β) - B (k_{1}, α, β)}

$f(x) = \frac{x^{(\alpha-1)}(1-x)^{(\beta-1)}}{\operatorname{B}(k2, \alpha, \beta) - \operatorname{B}(k_1, \alpha, \beta)}$

Pegue qualquer parte monotonicamente crescente da densidade entre e : para é log-côncava, para que você possa envolvê-lo com uma função exponencial desenhada tangente a qualquer ponto ao longo dele: $x_\mathrm{L}$ $x_\mathrm{U}$ $\alpha,\beta>1$

g (x) = c \cdot λ e^{- λ (x - x_{L})}

$g(x) = c \cdot \lambda \mathrm{e}^{-\lambda (x-x_\mathrm{L})}$

Encontre definindo os gradientes das densidades de log iguais $\lambda$

- λ = \frac{a - 1}{x} - \frac{b - 1}{1 - x}

$-\lambda = \frac{a-1}{x} - \frac{b-1}{1-x}$ e encontre calculando quanto a densidade exponencial precisa ser dimensionada para atender à densidade naquele ponto

c

$c$

c = \frac{f (x)}{λ e^{- λ (x - x_{L})}}

$c = \frac{f(x)}{\lambda\mathrm{e}^{-\lambda(x-x_\mathrm{L})}}$

O melhor envelope para amostragem de rejeição é aquele com a menor área abaixo dela. A área é Substituindo expressões em por & , e largando fatores constantes dá

A = c \cdot (1 - e^{- λ (x_{U} - x_{L})})

$A = c \cdot (1 - \mathrm{e}^{-\lambda(x_\mathrm{U}-x_\mathrm{L})})$

x

$x$

λ

$\lambda$

c

$c$

\begin{aligned} Q (x) = & \frac{x^{a} (1 - x)^{b}}{(a + b - 2) x - a + 1} \cdot \\ [\exp (\frac{(b - 1) (x - x_{L})}{1 - x} + \frac{x_{L} (a - 1)}{x} - (a - 1)) - \\ \exp (\frac{(b - 1) (x - x_{U})}{1 - x} + \frac{x_{U} (a - 1)}{x} - (a - 1))] \end{aligned}

$\begin{align} Q(x)= & \frac{x^a (1-x)^b}{(a+b-2)x - a+1} \cdot\\ & \left[\exp\left(\frac{(b-1)(x-x_L)}{1-x} + \frac{x_L (a-1)}{x} - (a-1)\right) - \right.\\ & \left. \exp\left(\frac{(b-1)(x-x_U)}{1-x} + \frac{x_U(a-1)}{x} - (a-1)\right)\right]\\ \end{align}$

Escrever a derivada é deixado como um exercício para os leitores ou seus computadores. Qualquer algoritmo de busca de raiz pode então ser usado para encontrar o para o qual . $\frac{\mathrm{d} Q}{\mathrm{d} x}$ $x$ $\frac{\mathrm{d} Q}{\mathrm{d} x} = 0$

O mesmo vale, mutatis mutandis, para partes monotonicamente decrescentes da densidade. Portanto, a distribuição beta truncada pode ser muito bem envolvida por duas funções exponenciais, digamos, uma de para o modo e outra do modo para , pronta para amostragem por rejeição. (Para uma variável aleatória uniforme truncada , tem uma distribuição exponencial truncada com o parâmetro de taxa .) $k_1$ $k_2$ $U$ $\frac{- \log(1-U)}{\lambda}$ $\lambda$

A beleza dessa abordagem é que todo o trabalho duro é feito. Uma vez definida a função de envelope, calculada a constante de normalização para a densidade beta truncada, tudo o que resta é gerar variáveis aleatórias uniformes e executar algumas operações aritméticas, logs e potências e comparações simples. Apertar a função de envelope - com linhas horizontais ou mais curvas exponenciais - pode, obviamente, reduzir o número de rejeições.

Scortchi - Restabelecer Monica
fonte

+1 Boa ideia. Como o Beta é aproximadamente normal para valores modestos a grandes de seus parâmetros, dependendo da proximidade entre si, o uso de um envelope gaussiano pode ser um pouco mais eficiente ainda.

whuber

@whuber: eu queria uma distribuição com uma distribuição quantílica de forma fechada para o envelope, mas suponho que você pudesse gerar as variáveis gaussianas truncadas eficientemente usando uma daquelas boas aproximações à função quantílica gaussiana. Ainda estou interessado no que você faria quando ou .

α < 1

$\alpha<1$

β < 1

$\beta<1$

Scortchi - Reinstate Monica

Para esses pequenos e , convém mudar para uma cauda exponencial. Não sei ao certo o que você quer dizer com "forma fechada": quando você olha atentamente para implementações computacionais de exponenciais, gaussianas, funções hipergeométricas etc. - ou seja, todas as funções não algébricas - você descobre que nenhuma das eles têm uma "forma fechada" geral: são calculados através de aproximações sucessivas, como séries de Taylor, frações parciais ou expansões assintóticas. Não há muita diferença entre calcular um exponencial e calcular um quantil gaussiano.

α

$\alpha$

β

$\beta$

whuber

@ whuber: (1) A abordagem adotada aqui para a construção de envelopes não funcionaria porque as densidades não são côncavas em log. (2) (a) eu quis dizer certamente funções algébricas + logs e potências, trig. funções se me perguntassem, e talvez até funções gama - confesso que não tinha uma noção precisa. (b) Ponto tomado - as rápidas avaliações de funções não se limitam àquelas com formulários fechados.

Scortchi - Restabelecer Monica

Bom argumento sobre falha de concavidade de log. Eu suspeito que uma distribuição de lei de energia deva criar um bom envelope para ou .

α < 1

$\alpha\lt 1$

β < 1

$\beta \lt 1$

whuber