Hello,

considering that in RStan ~ beta_binomial should not be used with the current release (Very different sampling with ~beta_binomial than with += beta_binomial_lpmf - #4 by jonah) I have implemented the truncated version as shown in the manual

```
y ~ poisson(3.7);
if (y < 2 || y > 10)
target += negative_infinity();
else
target += -log_sum_exp(poisson_lpmf(2 | 3.7),
log_diff_exp(poisson_lcdf(10 | 3.7),
poisson_lcdf(2 | 3.7)));
```

My function is this

```
real beta_binomial_truncated_lpmf(int p, int exposure, real alpha, real beta, int lower, int upper){
real lp = 0;
lp += beta_binomial_lpmf(p | exposure, alpha, beta );
if (p < lower || p > upper)
lp += negative_infinity();
else
lp += -log_sum_exp(
beta_binomial_lpmf(lower | exposure, alpha, beta),
log_diff_exp(beta_binomial_lcdf(upper | exposure, alpha, beta),
beta_binomial_lcdf(lower | exposure, alpha, beta))
);
return (lp);
}
```

However, a simple model takes an unreasonable amount of time

```
1000 transitions using 10 leapfrog steps per transition would take 37900 seconds.
```

Does anyone know about any bug in any of those densities calculations? Or Am I implementing it wrong?