ODE Event/Forcing Examples

[imadmali]

People have asked me about how to implement events and forcing functions in Stan and there's no straightforward way to implement it like in the deSolve R package. And there doesn't seem to be any centralized information about how to implement them in Stan (just various postings on google groups and discourse). The purpose of these examples is show how to implement events and forcing functions in Stan using a stylized example. I'm not sure but maybe this would be more helpful for users if it was located in the Stan User's Guide, in the ODE chapter.

[bob-carpenter]

It doesn't have to be either-or. It'd be great if you could contribute a patch to the user's guide ODE chapter.

[bob-carpenter]

I'm not the best one to review this, by the way---pinging @wds15 (who's very busy) and @charlesm93, who may also know others who could review. If nobody else wants to look at it, I can take a run from a beginner's perspective (I don't even know what a forcing function is).

[imadmali]

Ok cool. This isn't supposed to be complicated/advanced in anyway. I'm writing it from a beginners perspective myself. It's just a basic outline of how to implement the event/forcing function methodology. But regardless I'm happy to take feedback from Sebastian and Charles if they have the time.

[bbbales2]

Yo I had a look at this. I've got comments but I think it'll be easier to do them inline on paper and scan em' back.

[bbbales2]

Aka tomorrow!

[charlesm93]

Sorry I missed this. I'll also take a look today.

[charlesm93]

I think this works well as an introduction! Could you clarify what your goal and target audience is? The case study should reference previous work written on the subject, and also recommend further reading for folks who want to learn more.

A few notes:

• the definition of the forcing function is loose. Are you defining a precise technical term or a generic notion used in some sciences?
• "An event happens when the dependent variable suddenly changes." You haven't defined what the dependent variable is, and the term "suddenly" is imprecise. Also, this definition is restrictive. An event can be the start of an infusion, in which case your drug mass still changes continuously.
• "This causes a discontinuity in the variable's sequence." What kind of discontinuity (value, first derivative,...)? What is the variable's sequence?
• I'm not sure "intestine" is the right term. In the PK literature, I've usually seen "gut"
• The ODE system you give corresponds to a "one compartment model with a first-order absorption from the gut", which in the PK literature is different than the two compartment model.
• The ODE you give also has an analytical solution, so it's best not to solve it numerically. Even if you don't work out the analytical solution, you should use a matrix exponential. And even if you do decide to do it numerically, you better explain why you're using the BDF integrator instead of rk45 -- is your system stiff?

Here are some references I would include:

• First, there was this discussion on discourse. You may want to look at the case study @wds15 started writing.
• "Differential Equation Based Models in Stan" (Margossian & Gillespie, 2017, Stan Con 2017), which contains PK and PK/PD examples, and shows Stan code to handle forcing functions.
• I would mention Torsten, an extension of Stan for pharmacometrics, which has specialized routines to handle events from clinical trials. Here are some conference publications on the subject: (Margossian & Gillespie, 2016, Journal of Pharmacokinetics and Pharmacodynamics (JPP), volume 43) and (Gillespie & Zhang, 2018, JPP, volume 45) and (Zhang & Gillespie, 2019, PAGE 2019). @wds15 also wrote some code to facilitate the use of an event schedule, though I can't find the code right now.

Obviously there's much more that can be written on the subject: you can talk about PK/PD models, doing autodiff on ODEs, reference papers that use Stan and ODEs, doing PDEs in Stan, computing steady states, etc. I drafted a technical appendix for Torsten a while ago that covers some of these subjects, which I'm hoping to finish when Torsten v1.0 is released. It's a draft, very far from finished.

Tagging @yizhang-yiz.

[yizhang-yiz]

Would be nice to point out that events are static, in the sense that their temporal appearance doesn’t depend on the ODE unknowns. CVODES we have in Stan supports dynamic event and I simply don’t have time to implement it.

-YZ

On Oct 25, 2019 at 1:32 PM, <Charles Margossian (mailto:[email protected])> wrote:

I think this works well as an introduction, though I'm not sure what your goal and target audience is. Why did you write this? The case study should reference previous work written on the subject, and also recommend further reading for folks who want to learn more.

A few notes:

the definition of the forcing function is loose. Are you defining a precise technical term or a generic notion used in some sciences?

"An event happens when the dependent variable suddenly changes." You haven't defined what the dependent variable is, and the term "suddenly" is imprecise. Also, this definition is restrictive. An event can be the start of an infusion, in which case your drug mass still changes continuously.

"This causes a discontinuity in the variable's sequence." What kind of discontinuity (value, first derivative,...)? What is the variable's sequence?

I'm not sure "intestine" is the right term. In the PK literature, I've usually seen "gut"

The ODE system you give corresponds to a "one compartment model with a first-order absorption from the gut", which in the PK literature is different than the two compartment model.

The ODE you give also has an analytical solution, so it's best not to solve it numerically. Even if you don't work out the analytical solution, you should use a matrix exponential. And even if you do decide to do it numerically, you better explain why you're using the BDF integrator instead of rk45 -- is your system stiff?

Here are some references I would include:

First, there was this discussion on discourse (https://discourse.mc-stan.org/t/forced-odes-a-start-for-a-case-study). You may want to look at the case study @wds15 (https://github.com/wds15) started writing.

"Differential Equation Based Models in Stan" (Margossian & Gillespie, 2017 (https://mc-stan.org/events/stancon2017-notebooks/stancon2017-margossian-gillespie-ode.html), Stan Con 2017), which contains PK and PK/PD examples, and shows Stan code to handle forcing functions.

I would mention Torsten (https://github.com/metrumresearchgroup/Torsten), an extension of Stan for pharmacometrics, which has specialized routines to handle events from clinical trials. Here are some conference publications on the subject: (Margossian & Gillespie (https://www.researchgate.net/publication/323834461_Stan_functions_for_pharmacometrics_modeling), 2016, Journal of Pharmacokinetics and Pharmacodynamics (JPP), volume 43) and (Gillespie & Zhang (https://link.springer.com/article/10.1007/s10928-018-9606-9), 2018, JPP, volume 45) and (Zhang & Gillespie (https://metrumrg.com/wp-content/uploads/Pubs/torsten_page_2019_poster.pdf), 2019, PAGE 2019). @wds15 (https://github.com/wds15) also wrote some code to facilitate the use of an event schedule, though I can't find the code right now.

Obviously there's much more that can be written on the subject: you can talk about PK/PD models, doing autodiff on ODEs, reference papers that use Stan and ODEs, doing PDEs in Stan, computing steady states, etc. I drafted a technical appendix for Torsten (https://github.com/charlesm93/presentations/blob/master/TorstenAppendix/Torsten_appendix.pdf) a while ago that covers some of these subjects, which I'm hoping to finish when Torsten v1.0 is released. It's a draft, very far from finished.

Tagging @yizhang-yiz (https://github.com/yizhang-yiz).

— You are receiving this because you were mentioned. Reply to this email directly, view it on GitHub (https://github.com/stan-dev/example-models/pull/159?email_source=notifications&email_token=AH4SDMACJWYE4AFVFSH65C3QQNJWBA5CNFSM4JBGB6U2YY3PNVWWK3TUL52HS4DFVREXG43VMVBW63LNMVXHJKTDN5WW2ZLOORPWSZGOECJPQIQ#issuecomment-546502690), or unsubscribe (https://github.com/notifications/unsubscribe-auth/AH4SDMCJDQ44MGZWN65RTP3QQNJWBANCNFSM4JBGB6UQ).

[imadmali]

Thanks for the feedback! I'll incorporate some of the suggestions.

@charlesm93 the main audience for this was anyone who uses ODEs and events/forcing functions and would like to understand how to implement them in Stan; not just the people in the pharma community. Sorry for the confusion about that. I made it seem like I was building a case study around a pharma application when in fact that was the easiest way for me to present the method of implementing events/forcing functions in Stan. That's why I omitted a lot of detail (e.g. choice of integrator, compartment definition, etc).

Originally I planned for the PR to have just R scripts and Stan files like the examples in /misc. However, because the implementation of events and forcing functions is not straightforward I thought it would be best to frame this with an overly simplified PK/PD example. I'm thinking it might be better to revert back to the option of putting code in /misc and putting more detail in the documentation as @bob-carpenter mentioned. That way I can do a better job of abstracting away the pharma example.