Skip to contents

State Updates

Source: R/formula_list_generators.R
mp_euler.Rd

Use these functions to update a model spec so that the state variables are updated according to different numerical methods.

Usage

mp_euler(model)

mp_rk4(model)

mp_euler_multinomial(model)

mp_hazard(model)

Arguments

model

Object with quantities that have been explicitly marked as state variables.

Details

The default update method for model specifications produced using mp_tmb_model_spec is mp_euler. This update method yields a difference-equation model where the state is updated once per time-step using the absolute flow rate as the difference between steps.

These state update functions are used to modify a model specification to use a particular kind of state update. To see these modified models for a particular example one may use the mp_expand function (see examples).

Functions

  • mp_rk4(): ODE solver using Runge-Kutta 4. Any formulas that appear before model flows in the during list will only be updated with RK4 if they do contain functions in getOption("macpan2_non_iterable_funcs") and if they do not make any state variable assignments (i.e., the left-hand-side does not contain state variables). Each formula that does not meet these conditions will be evaluated only once at each time-step before the other three RK4 iterations are taken. By default, the time_var function and functions that generate random numbers (e.g., rbinom) are not iterable. Functions that generate random numbers will only be called once with state update methods that do not repeat expressions more than once per time-step (e.g., mp_euler), and so repeating these functions with RK4 could make it difficult to compare methods. If you really do want to regenerate random numbers at each RK4 iteration, you can do so by setting the above option appropriately. The time_var function assumes that it will only be called once per time-step, and so it should never be removed from the list of non-iterable functions. Although in principle it could make sense to update state variables manually, it currently causes us to be confused. We therefore require that all state variables updates are set explicitly (e.g., with mp_per_capita_flow) if any are explicit.

  • mp_euler_multinomial(): Update state with process error given by the Euler-multinomial distribution.

  • mp_hazard(): Update state with hazard steps, which is equivalent to taking the step given by the expected value of the Euler-multinomial distribution.

Examples

sir = mp_tmb_library("starter_models", "sir", package = "macpan2")
sir
#> ---------------------
#> Default values:
#> ---------------------
#>  matrix row col value
#>    beta           0.2
#>   gamma           0.1
#>       N         100.0
#>       I           1.0
#>       R           0.0
#> 
#> ---------------------
#> Before the simulation loop (t = 0):
#> ---------------------
#> 1: S ~ N - I - R
#> 
#> ---------------------
#> At every iteration of the simulation loop (t = 1 to T):
#> ---------------------
#> 1: mp_per_capita_flow(from = "S", to = "I", rate = infection ~ I * 
#>      beta/N)
#> 2: mp_per_capita_flow(from = "I", to = "R", rate = recovery ~ gamma)
#> 
sir |> mp_euler()             |> mp_expand()
#> ---------------------
#> Default values:
#> ---------------------
#>  matrix row col value
#>    beta           0.2
#>   gamma           0.1
#>       N         100.0
#>       I           1.0
#>       R           0.0
#> 
#> ---------------------
#> Before the simulation loop (t = 0):
#> ---------------------
#> 1: S ~ N - I - R
#> 
#> ---------------------
#> At every iteration of the simulation loop (t = 1 to T):
#> ---------------------
#> 1: infection ~ S * (I * beta/N)
#> 2: recovery ~ I * (gamma)
#> 3: S ~ S - infection
#> 4: I ~ I + infection - recovery
#> 5: R ~ R + recovery
#> 
sir |> mp_rk4()               |> mp_expand()
#> ---------------------
#> Default values:
#> ---------------------
#>  matrix row col value
#>    beta           0.2
#>   gamma           0.1
#>       N         100.0
#>       I           1.0
#>       R           0.0
#> 
#> ---------------------
#> Before the simulation loop (t = 0):
#> ---------------------
#> 1: S ~ N - I - R
#> 
#> ---------------------
#> At every iteration of the simulation loop (t = 1 to T):
#> ---------------------
#>  1: infection ~ S * (I * beta/N)
#>  2: recovery ~ I * (gamma)
#>  3: k1_S ~ -infection
#>  4: k1_I ~ +infection - recovery
#>  5: k1_R ~ +recovery
#>  6: infection ~ (S + (k1_S/2)) * ((I + (k1_I/2)) * beta/N)
#>  7: recovery ~ (I + (k1_I/2)) * (gamma)
#>  8: k2_S ~ -infection
#>  9: k2_I ~ +infection - recovery
#> 10: k2_R ~ +recovery
#> 11: infection ~ (S + (k2_S/2)) * ((I + (k2_I/2)) * beta/N)
#> 12: recovery ~ (I + (k2_I/2)) * (gamma)
#> 13: k3_S ~ -infection
#> 14: k3_I ~ +infection - recovery
#> 15: k3_R ~ +recovery
#> 16: infection ~ (S + k3_S) * ((I + k3_I) * beta/N)
#> 17: recovery ~ (I + k3_I) * (gamma)
#> 18: k4_S ~ -infection
#> 19: k4_I ~ +infection - recovery
#> 20: k4_R ~ +recovery
#> 21: S ~ S + (k1_S + 2 * k2_S + 2 * k3_S + k4_S)/6
#> 22: I ~ I + (k1_I + 2 * k2_I + 2 * k3_I + k4_I)/6
#> 23: R ~ R + (k1_R + 2 * k2_R + 2 * k3_R + k4_R)/6
#> 
sir |> mp_euler_multinomial() |> mp_expand()
#> ---------------------
#> Default values:
#> ---------------------
#>  matrix row col value
#>    beta           0.2
#>   gamma           0.1
#>       N         100.0
#>       I           1.0
#>       R           0.0
#> 
#> ---------------------
#> Before the simulation loop (t = 0):
#> ---------------------
#> 1: S ~ N - I - R
#> 
#> ---------------------
#> At every iteration of the simulation loop (t = 1 to T):
#> ---------------------
#> 1: infection ~ reulermultinom(S, (I * beta/N))
#> 2: recovery ~ reulermultinom(I, gamma)
#> 3: S ~ S - infection
#> 4: I ~ I + infection - recovery
#> 5: R ~ R + recovery
#>