**Thiago G. Martins » R**, and kindly contributed to R-bloggers]. (You can report issue about the content on this page here)

Want to share your content on R-bloggers? click here if you have a blog, or here if you don't.

Recently I had to define a R function that would compute the -th quantile of a continuous random variable based on an user-defined density function. Since the main objective is to have a general function that computes the quantiles for any user-defined density function it needs be done numerically.

** Problem statement **

Suppose we are interested in numerically computing the -th quantile of a continuous random variable . Assume for now that is defined on . So, the problem is to find such that , where is the cumulative distribution function of ,

and is its density function.

** Optimization **

This problem can be cast into an optimization problem,

That said, all we need is a numerical optimization routine to find that minimize and possibly a numerical integration routine to compute in case it is not available in closed form.

Below is one possible implementation of this strategy in R, using the `integrate`

function for numerical integration and `nlminb`

for the numerical optimization.

CDF <- function(x, dist){ integrate(f=dist, lower=-Inf, upper=x)$value } objective <- function(x, quantile, dist){ (CDF(x, dist) - quantile)^2 } find_quantile <- function(dist, quantile){ result = nlminb(start=0, objective=objective, quantile = quantile, dist = dist)$par return(result) }

and we can test if everything is working as it should using the following

simple task of computing the 95th-quantile of the standard Gaussian distribution.

std_gaussian <- function(x){ dnorm(x, mean = 0, sd = 1) } find_quantile(dist = std_gaussian, quantile = 0.95) [1] 1.644854 qnorm(0.95) [1] 1.644854

**Random variables with bounded domain **

If has a bounded domain, the optimization problem above becomes more tricky and unstable. To facilitate the numerical computations we can transform using a monotonic and well behaved function , so that is defined on the real line. Then we apply the numerical solution given above to compute the -th quantile of , say , and then transform the quantile back to the original scale . Finally, is the -th quantile of that we were looking for.

For example, if we could use .

**leave a comment**for the author, please follow the link and comment on their blog:

**Thiago G. Martins » R**.

R-bloggers.com offers

**daily e-mail updates**about R news and tutorials about learning R and many other topics. Click here if you're looking to post or find an R/data-science job.

Want to share your content on R-bloggers? click here if you have a blog, or here if you don't.