I have a system of ordinary differential equations model describing the dynamics of transmission of malaria fitted to malaria incidence data.
After the fitting process, I got the posterior distributions of my parameters but I want to get the set of parameter values that generated the posterior predictions that I got such that when I use them for numerical simulation of the system of ODEs model my outcome will be the same as the posterior predictions I got from my fitting process.
Hi, @Idowu. Sorry this has taken so long to answer.
The posterior draws from the sampler are the parameter values that generate posterior predictions. In general, if you have a joint density p(\theta, y) for data y and parameters \theta and you want to do posterior predictive inference for new data \widetilde{y}, then it looks like
where \theta^{(m)} \sim p(\theta \mid y) are your posterior draws.
I wasn’t quite sure what you meant here. There aren’t any posterior predictions from fitting in general, just draws from the posterior. You can put posterior predictions in the generated quantities block. The third part of the User’s Guide explains how to do posterior predictive inference.
Thank you for your response.
I just realised that I did not pay attention to some details in the third part of the users guide on how to do posterior predictive inference. I had to go through it again, together with the work of Alex Pavlakis on “Making Predictions from Stan models in R”.
I was able to get a hint that I am now trying to implement.
I did not know I could use the draws from the posterior distributions to do posterior prediction. I was initially wrongfully seeking how to obtain point estimates from the draws that I could use to do simulation in R which would yield the exact fitted result that I got.
Good day @Bob_Carpenter and all
Please, one more thing to clarify.
I am trying to adapt the concept of posterior predictive inference to my model of ODEs using stand alone generated quantity but I am finding it difficult to do, please is there any resource where posterior predictive inference with this type of model is exemplified?
Any help to get around this will be appreciated.
I found @Bob_Carpenter’s case study on the Lotka-Volterra ODE very helpful to get my head around ODEs in Stan. In that model, the posterior predictions for the system state at each observed times point are treated as parameters. I’ve built on Bob’s example here to make posterior plots of the system state and of a derived quantity (the oscillation period).
If you want state estimates at unobserved time points (e.g., interpolating or fore/backcasting), you have several options. There’s a user guide entry here for calculating these directly in Stan, which is perhaps the most elegant approach. Alternatively, you can solve the ODE with R. Define your ODE function for use with a solver like deSolve. Then, extract the joint posterior draws of the ODE parameters from your model fit with posterior::as_draws_df(mymodel). Finally, you can supply the solver with the R function describing the ODE, the time points you want to solve at, and row of parameter draws. Looping through each row of joint posterior ODE parameter draws will give you the posterior distribution at the desired time point.