Using the gamma function of a fraction in a custom pdf

I have a custom pdf to model in Stan. This is problem 16.3 in Ben Lambert’s newish book. The pdf is:
p(X)=[∛b/(Γ (4/3)]*[ exp(-b(X^3 )]

We are asked to estimate b,but there is a gamma function of a fraction(4/3) in the pdf which Stan doesn’t like.
My code:

    real example_lpdf(real aX,real aMu){
return(log((1/(tgamma(4/3)))*(aMu^(1/3))*(exp(-aMu*(aX^3))))); //doesn't want to parse the gamma of a fraction

int N;
real X[N];

real b;

for (i in 1:N){
X[i] ~ example(b);
b~normal(0,2); //this looks like where most of it would be

It runs if you leave the gamma function out, and still runs with it but gives odd results.The Stan manual says “gammat” should work with a positive number,but it wants to round the fraction to an integer.

Maybe gamma(4/3) simplifies to a function of pi,but I couldn’t see this anywhere.
The datatset (stan_survival.csv)for this problem is on the website for the book “A Student’s Guide to Bayesian Statistics” by Ben Lambert(Sage Publications).

Thanks for any advice.


I guess 4/3 is 1 and gamma of 1 is 1 and log(1/1) = 0 and then the whole expression gets 0.
I think the most easy solution is to replace 4/3 by 4.0/3.0. Since this is a constant you can
omit the gamma part. I did not further check your statement, the function uses aX, aMu
the specification b and X. You may change this for clarity.

Thanks Andre.
It does seem easier to leave the constant/gamma term out but in some of the examples in this text book, the author leaves them in the model(pi,for example).

You can leave the constant in, just precalculate it as data.

1 Like

Thanks for your advice. This does make it more tractable.