Using the gamma function of a fraction in a custom pdf

Modeling
1

Hello
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:

functions{
    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
}
}

data{
int N;
real X[N];
}

parameters{
real b;
}

model{
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.

Chris

2

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.

3

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).
Regards
Chris

4

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

1 Like
5

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