Hello, is there an example on how to code up a custom RNG for a non-standard discrete distribution?

The examples in the manual are only for continuous distribution?

I am interested in the Discrete-Weibull distribution.

But any example is welcome and I will try to adapt it for my need. Thank you.

Iām assuming you are interested in discrete user-defined probability functions (and not in the _rng functions)?

Here is a model that fits a discrete weibull distribution (parameterized by \alpha and \beta as in the wikipedia article) to some data `y`

:

```
functions {
real discrete_weibull_lpmf(int x, real alpha, real beta) {
return log_diff_exp(-pow(x * 1.0 / alpha, beta), -pow(((x + 1.0) / alpha), beta));
}
}
data {
int N;
int<lower=0> y[N];
}
parameters {
real<lower=0> alpha;
real<lower=0> beta;
}
model {
alpha ~ normal(0, 5);
beta ~ normal(0, 5);
for (i in 1:N) {
y[i] ~ discrete_weibull(alpha, beta);
}
}
```

log_diff_exp(a, b) is the same as log(exp(a) - exp(b)), just a bit more numerically stable.

In case you are actually interested in a discrete_weibull_rng function, a common trick is to draw from a standard uniform distribution and transform these by the inverse cdf.

Thanks. I actually AM interested in a discrete_weibull_rng function!

Once I solve for *x* using the inverse cdf, do I use the floor function to ensure it is an integer?

That would be my first attempt. You can always produce samples using your scheme and then simulate a large number of random draws and compare the frequencies with the true probabilities.

Note that all rounding functions return a real, so your random draws are always represented as `2.0, 3.0`

etc. There is no good workaround to convert them to integers right now. Depending on your use case the easiest solution is probably to do the sampling outside of Stan.

2 Likes