0
$\begingroup$

Given an easily-computable probability density function $f(x)$, what algorithm can we use to numerically approximate percentiles?

For instance, we might be looking for $x$ such that given $X \sim f(x)$, $p(X \le x) = 95\%$.


Here is a use case that motivates this:

We're estimating the shape $\alpha$ and inverse scale $\theta$ of a gamma distribution (from a sequence of samples). Via Wikipedia we find a conjugate prior of the following form, originally developed by Miller (1980):

$$ f(\alpha, \theta) \propto \frac{\theta^{v \alpha - 1} p^{\alpha - 1} e^{-s\theta}}{\Gamma(\alpha)^n} $$

This somewhat unwieldy function represents our belief about likely values for $\alpha$ and $\theta$. We are only interested in $\theta$, so we find its pdf by integration:

$$ g(\theta) \propto \int f(\alpha, \theta) d\alpha $$

Miller (1980) suggests that we need numerical integration to compute this function.

Now we can easily plot this $g(\theta$) function by repeatedly approximating it for many $\theta$ values. From this density plot we can visually glean the most probably region for $\theta$. So far so good.

But we also want to compute the $\theta$ values for the 5th and 95th percentiles, $G^{-1}(0.05)$ and $G^{-1}(0.95)$, where $G^{-1}(y)$ is the inverse cdf. So my question is, how do we approximate those given $g(\theta)$?

There must surely be a way to numerically approximate this by harnessing some well-known algorithm, but I cannot find a way to express $G^{-1}$ in terms of $g$. Am I missing some intermediate step to make it work? The best I can come up with is getting a lot of samples for $g$ and summing them.

References:

Miller, Robert B. "Bayesian Analysis of the Two-Parameter Gamma Distribution" Technometrics, 1980, 22(1), 65-69

I haven't done math in a while, so it's all rusty -- if there's mistakes or misphrasings in this question, please edit away!

$\endgroup$
2
$\begingroup$

You should calculate the cdf G as the integral of g and solve G(x) = 0.95. If you word it: G(x)-0.95 = 0, you can use the newton method.

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.