I’m looking for a way to do the same for the beta distribution: to sample from the beta distribution with an upper bound. But it seems that stan does not have a built-in inverse cdf for the beta distribution. Any help is much appreciated!
real beta_ub_rng(real alpha, real beta, real ub) {
real p_ub = beta_cdf(ub, alpha, beta); // cdf for upper bound
real u = uniform_rng(0, p_ub);
return ???? // inverse cdf for value
}
Thanks! I was also wondering for gamma: How would I sample from gamma with lower or upper bound?
This is what I have but I’m too new to this to be sure that I’m doing it right:
real gamma_lb_rng(real alpha, real beta, real lb) {
real p_lb = gamma_cdf(lb, alpha, beta);
real u = uniform_rng(p_lb, 1);
real y = gamma_cdf(u, alpha, beta);
return y;
}
Unfortunately, that’s not correct, you need real y = gamma_inverse_cdf(u, alpha, beta);
Where gamma_inverse_cdf is the inverse CDF (also called quantile functions and PPF) of the gamma function, which AFAIK is not available in Stan. But if you are willing to invoke C++, it should be easy to implement it using boost (e.g. as discussed at Quantile functions in boost (C++) - Stack Overflow).
Also note that there is a possibility for confusion as there is also “inverse gamma” distribution, so that “CDF of (inverse gamma distribution)” is something completely different than “(inverse CDF) of gamma distribution” (which is the one you want).
Thanks @martinmodrak! so inv_gamma_cdf() which exists in stan is the cdf of inverse gamma distribution, which is NOT what I want, correct? Thanks for pointing this out!
Correct. Also, instead of dipping into C++ it might just be easier to prepare the necessary data in Stan but call the inverse CDF from R or whatever language you invoke Stan from (in R it is called qgamma)