functions { real lgammaJacLinf(vector x, real alpha, real lambda); } data { int N; vector[N] y; } parameters { real alpha; real lambda; } model { y ~ gamma_lpdf(alpha,1/lambda); target+=lgammaJacLinf(y,alpha,lambda); }