Message Boards Message Boards

0
|
5175 Views
|
3 Replies
|
3 Total Likes
View groups...
Share
Share this post:

GammaDistribution different in R and Mathematica?

Posted 11 years ago
I have been doing some Bayesian data analysis using the Gamma Distribution as conjuage prior to the Poisson Distribution for a particular problem. The book that I am using as a reference refers to the Gamma Distribution with parameters Gamma(3,5). When I go through the exercise in the book using R I get a mean value of 0.6, which is what I expect and a plot of the distribution also shows that values from about 5 onwards are close to zero. This reflects what the book also says. I like to use Mathematica and I often use it to complement to R. However, in this case when I plot GammaDistribution[3,5] in Mathematica I get quite a different curve to what I have in R.

Any ideas why this might be?
POSTED BY: David Post
3 Replies
Are you sure that you are giving/using the parameters correctly?  

In[1]:= Mean[GammaDistribution[a, b]]


Out[1]= a b
 
So that one would expect your mean to be 15 for 
GammaDistribution[3, 5]
POSTED BY: David Reiss
In these cases it's a good idea to look up the precise definition and usage of these functions in the respective documentation pages.  Often there are differences between system, as even in math books there's no one true convention that everyone follows.

The Mathematica documentation is very good in that it will give a precise definition of every mathematical functin of distribution (usually in a very clear and understandable manner).  In my experience the documentation of R packages is usually not nearly as good: it will often just say "this is a gamma distribution" or "this is an incidence matrix" without clarifying which meaning of those words it's referring to.

So you can look up that Mathematica's GammaDistribution[a,b] represents a dsitrbution with PDF proportional to Exp[-x/b] x^(a-1).

It seems R has several packages and several functions for gamma distritbuions.  The set of R functions named as *gamma() have two settable parameters: the shape and scale.  "Shape" corresponds to 'a' in the above formula and "scale" corresponds to 'b'.

As David said, Mean@GammaDistribution[a,b] is a*b.

I'm not very fluent in R, but if I do
x <- rgamma(1000, shape=2, scale=3)
mean(x)

I get a value close to 6=2*3, i.e. consistent with Mathematica.

Conclusion: you'll need to look up the precise meaning of the parameters the books uses and see how they transform to the parameters Mathematica uses.  Most likely the second parameter in the book refers to "1/scale", which R calls "rate".
POSTED BY: Szabolcs Horvát
Posted 11 years ago

Szabolocs' advice is excellent in general and specifically correct here. R allows two ways to specify the gamma distributions and one obtains what David obtained when the shape and scale keywords are not used explicitly:

mean(rgamma(1000, 3, 5))

The difference is that the above is equivalent to

mean(rgamma(1000, 3, scale=1/5))
POSTED BY: Jim Baldwin
Reply to this discussion
Community posts can be styled and formatted using the Markdown syntax.
Reply Preview
Attachments
Remove
or Discard

Group Abstract Group Abstract