
Boost Users : 
Subject: Re: [Boostusers] [distributions]: Inverse Gamma
From: Thomas Mang (Thomas.Mang_at_[hidden])
Date: 20100804 11:01:44
On 04.08.2010 10:37, Paul A. Bristow wrote:
>
>
>> Original Message
>> From: boostusersbounces_at_[hidden] [mailto:boostusersbounces_at_[hidden]] On Behalf Of Thomas Mang
>> Sent: Tuesday, August 03, 2010 3:07 PM
>> To: boostusers_at_[hidden]
>> Subject: Re: [Boostusers] [distributions]: Inverse Gamma
>>
> >>
>>>> any plans of implementing the inverse gamma distribution as part of the distributions library ?
>>>
>>> This looks possible  but I'm curious about applications  you obviously have one, but Wikipedia doesn't mention any
>>>
>>> http://en.wikipedia.org/wiki/Inversegamma_distribution
>>>
>>> But you obviously have one ;)
>>
>> Yes I truly have one ;) The inverse gamma distribution and its special
>> case, the scaled inverse chisquare distribution, is the conjugate prior
>> to the normal distribution variance parameter in Bayesian statistics.
>> Pretty much as uncommon and unheard of as it is outside Bayes world [to
>> the best of my knowledge], it's very much central to Bayesian stats and
>> appears in every textbook right after the introduction chapter ;)
>>
>> http://en.wikipedia.org/wiki/Scaled_inverse_chisquare_distribution
>> http://en.wikipedia.org/wiki/Conjugate_prior
>>
>> Hence I wonder it has not been requested so far  but being a Bayesian
>> C++ / booster I definitely want / need it :).
>>
>> @John: Yes it is a transformation deviate of the gamma, and an easy so.
>> And it should be fairly easy to implement IMHO.
>
> The pdf and pdf etc looks fairly straightforward (only uses exp, pow and gamma?) so I might be persuaded to do these.
> But the inverses (qhantiles) may prove more troublesome if have to be done by brute force numerically
>
> R library does it numerically http://rss.acs.unt.edu/Rdoc/library/pscl/html/igamma.html
>
>  or are there analytic expressions for these?
>
> Or can it use the inverse of the gamma distribution?
Yes pdf just uses exp, pow and gamma function, for the cdf /
complementary cdf also the incomplete gamma function is needed  all in
the library available.
For quantiles yes you can use quantiles of the gamma distribution.
Indeed all of pdf, cdf , complementary cdf and quantiles can be derived
using calls to the corresponding functions of the gamma distribution:
If we define both the gamma and inverse gamma distribution with the
shape and scale parameter (as is done for the boost implementation of
the gamma) and in that constructor order, then the relations are [should
be] as follows:
typedef gamma_distribution<...> gamma_dist;
typedef inverse_gamma_distribution<...> inv_gamma_dist;
pdf(inv_gamma_dist(a, b), x) == pdf(gamma_dist(a, 1/b), 1/x) * 1/(x*x);
cdf(inv_gamma_dist(a, b), x) == 1  cdf(gamma_dist(a, 1/b), 1/x);
quantile(inv_gamma(a, b), p) == 1 / quantile(gamma_dist(a, 1/b), 1/x);
where for the pdf the extra term 1/(x*x) is the density transformation
Jacobian term while for the probabilities no equivalent term is needed.
Note that I just wrote "1" as nominator to save space and avoid a
newsgrouplinebreak; it should be a floating point type in real code of
course. And my "[should be]" refers to that I just worked that out and
tested it in R, but not in C++. But it will be doublechecked anyway
after it's implemented. Why R evaluates the quantile numerically is
beyond my understanding.
Whether you implement the pdf / cdf directly or go for above relations
is probably an implementation issue, in terms of performance it
shouldn't make a big difference.
If you go for implementation that would be great and very much
appreciated with lots of thanks. Could you then also implement the
inverse chisquare and scaled inverse chisquare distribution as well
(they are just special inverse gamma distributions, but their names
might ring more bells and I'd find it appealing to offer them straight
away, similar to Chisquare and Gamma distribution cases) ?
>
>> Is contribution on my side expected (can be done just notice I am a
>> [heavy !] user of the stats library only, not familiar with code /
>> numerical stability issues).
>
> I'm not mathematician enough to deal with this  but I can deal with the obfuscated code (by templating and policies) if
> you can provide the equations.
>
> (And there is the question of testing  some parameters and value combinations (preferably exact) are needed for sanity
> and accuracy checks).
>
> Are there use cases for the inverses?
What is "usually" (with my subjective definition of "usually") needed
is pdf and, fairly important actually for various applications, random
numbers. Out of memory not a single application for cdf or quantiles
comes to mind. Of much greater importance are spherical properties such
as Mahalanobis distance and area for confidence ellipses etc.
To the best of my knowledge more advanced matrix algebra is required
however (inverses, decompositions). I have to check with the uBLAS
library what they offer as base, and what is needed.
best,
Thomas
Boostusers list run by williamkempf at hotmail.com, kalb at libertysoft.com, bjorn.karlsson at readsoft.com, gregod at cs.rpi.edu, wekempf at cox.net