Skip to content

Failing to reproduce a model #16

@sdwfrost

Description

@sdwfrost

I'm trying to reproduce a 'standard' SIR model in ReactiveDynamics, but I get very different results (my hand-coded discrete problem gives the expected result). Can you point me in the direction of where I'm going wrong in defining my model with ReactiveDynamics?

using ReactiveDynamics
using OrdinaryDiffEq
using Random
using Distributions
using Plots

function sir_markov(u,p,t)
    (S, I, R) = u
    (β, γ, N, δt) = p
    infrate = abs*I/N*S*δt)
    recrate = abs*I*δt)
    infection = rand(Poisson(infrate))
    recovery = rand(Poisson(recrate))
    S = S-infection
    I = I+infection-recovery
    R = R+recovery
    [S, I, R]
end

tspan = (0.0,40.0)
u0 = [990, 10, 0] # S, I, R
β = 0.5
γ = 0.25
N = 1000.0
δt = 0.1
p = [β, γ, N, δt] # β, γ, N
seed = 1234
Random.seed!(seed)

prob = DiscreteProblem(sir_markov, u0, tspan, p, dt=δt)
sol = solve(prob, FunctionMap())
plot(sol)

# model dynamics
sir_acs = @ReactionNetwork begin
        β*S*I, S+I --> 2I, name=>infection
        γ*I, I --> R, name=>recovery 
end

# simulation parameters
## initial values
@prob_init sir_acs S=990 I=10 R=0
## uncertainty in initial values (Gaussian)
@prob_uncertainty sir_acs S=0.0 I=0.0
## parameters
@prob_params sir_acs β=0.5 γ=0.25
## other arguments passed to the solver
@prob_meta sir_acs tspan=40 dt=0.01
# turn model into a discrete? problem
sir_prob = @problematize sir_acs
# solve the problem over multiple trajectories
sir_sol = @solve sir_prob
# plot the solution
@plot sir_sol

Metadata

Metadata

Assignees

No one assigned

    Labels

    No labels
    No labels

    Type

    No type

    Projects

    No projects

    Milestone

    No milestone

    Relationships

    None yet

    Development

    No branches or pull requests

    Issue actions