# How to fit epidemic model in Stan

I would like to forecast R_t in an epidemic (the expected number of new infectious produced per single infection at time t) and have the following (simplified) model that I would like to fit in Stan:

I_t\sim \text{Poisson}(R_t\sum_{s=1}^t I_{t-s} w_s),

where I_t is the true incidence at time t and w_s is a weighting variable that accounts for the infectious period etc. of the disease which is strongly determined by prior epidemiological knowledge. We are dealing with daily data and so, for the diseases we are considering, s is typically around 50 (days).

Additionally, we don’t actually observe the true incidence at any given time; we see an observed incidence assuming a fixed observation process:

I^{obs}_t\sim \text{binomial}(I_t, p)

where p is the probability that we detect a given case which is heavily informed by another analysis.

I haven’t a clue how to marginalise out all the I_t parameters and whether this is possible in Stan.
Currently, we are using a Gibbs sampler to do the MCMC instead but wanted to ask if anyone has ideas as to how to estimate this model in Stan?

I can’t see an easy way to marginalise them out. How important is it that values of I_{t} are actually discrete? If the I_{t} values are all sufficiently large then you might be able to pretend that it is continuous, or at least evaluate the binomial density function using its continuous analog.

Thanks for your comment. No, unfortunately, I_t is small and very stochastic (we’re looking at small regions).

True “I” can be marginalized out. Try

I_t^{observation}\sim \text{Poisson}(p \times R_t\sum_{s=1}^t I_{t-s} w_s)

You just need I_t as a transformed parameter. Using I_t^{observation}/p^t is a common if not great choice. If$I$ is not huge you can explicitly marginalized it. There’s a variety of solutions but they depend on the context.