<jupyter_start><jupyter_text>Model 4.3: SEIR in an open population in Python/PyGOM*Author*: <NAME> @twomagpi*Date*: 2018-10-03<jupyter_code># housekeeping
import numpy
from pygom import DeterministicOde, TransitionType, SimulateOde, Transition<jupyter_output><empty_output><jupyter_text>Build the model system<jupyter_code># first set up the states
states = ['Susceptible', 
         'Preinfectious', 
         'Infectious', 
         'Immune'
        ]
# now set up the parameters
parameters = ['beta',
              'infectious_rate',
              'rec_rate',
              'total_popn',
              'b_rate',
              'm_rate']

# setup the ODEs
odes = [
    Transition(origin='Susceptible', 
               equation='-beta*Susceptible*Infectious + total_popn*b_rate - Susceptible*m_rate',
               transition_type=TransitionType.ODE),
    Transition(origin='Preinfectious', 
               equation='beta*Susceptible*Infectious - Preinfectious*infectious_rate - Preinfectious*m_rate',
               transition_type=TransitionType.ODE), 
    Transition(origin='Infectious', 
               equation='Preinfectious*infectious_rate - Infectious*rec_rate - Infectious*m_rate',
               transition_type=TransitionType.ODE), 
    Transition(origin='Immune', 
               equation='Infectious*rec_rate - Immune*m_rate',
               transition_type=TransitionType.ODE)
    ]
# now create the system
ode_obj = DeterministicOde(states, 
                           parameters,
                           ode=odes)
ode_obj.print_ode()<jupyter_output>⎡ dSusceptible/dt=          -Infectious⋅Susceptible⋅β - Susceptible⋅mᵣₐₜₑ + bᵣ
⎢                                                                             
⎢dPreinfectious/dt=  Infectious⋅Susceptible⋅β - Preinfectious⋅infectiousᵣₐₜₑ -
⎢                                                                             
⎢ dInfectious/dt=        -Infectious⋅mᵣₐₜₑ - Infectious⋅recᵣₐₜₑ + Preinfectiou
⎢                                                                             
⎣   dImmune/dt=                           -Immune⋅mᵣₐₜₑ + Infectious⋅recᵣₐₜₑ  

ₐₜₑ⋅totalₚₒₚₙ       ⎤
                    ⎥
 Preinfectious⋅mᵣₐₜₑ⎥
                    ⎥
s⋅infectiousᵣₐₜₑ    ⎥
                    ⎥
                    ⎦<jupyter_text>Now set up the parameters<jupyter_code>## Demography ##
# total population size
total_popn = 100000
# life expectancy in years
life_expectancy_yrs = 70
# per capita mortality rate
m_rate = 1/(life_expectancy_yrs*365)
# per capita birth rate
b_rate = m_rate

## Transmission and Infection parameters ##
# average pre-infectious period
preinfectious_period = 8
# average infectious period
infectious_period = 7
# R0
R0 = 13
# rate at which individuals become infectious
infectious_rate = 1/preinfectious_period
#rate at which individuals recover to become immune
rec_rate = 1/infectious_period
# rate at which two specific individuals come into effective contact per unit time
beta = R0/(total_popn*infectious_period)
prop_immune_0 = 0.0    

## Inital Values ##
#initial number of immune individuals
Immune_0 = total_popn*prop_immune_0
# initial number of infectious individuals
Infectious_0 = 1
# initial number of susceptible individuals
Sus_0 = total_popn-Infectious_0 -Immune_0
# initial number of pre-infectious individuals
Preinfectious_0 = 0<jupyter_output><empty_output><jupyter_text>Solve the model system<jupyter_code># Create the timeline
t = numpy.linspace(0, 36500, 1001)
# Set the inital values for the states
ode_obj.initial_values = ([Sus_0, 
                           Preinfectious_0,
                           Infectious_0,
                           Immune_0], 
                          t[0])

#set the values for the parameters
ode_obj.parameters = {'beta': beta,
                      'infectious_rate': infectious_rate,
                      'rec_rate': rec_rate,
                      'total_popn': total_popn,
                      'b_rate': b_rate,
                      'm_rate': m_rate
                     }

solution,output = ode_obj.integrate(t[1::], full_output=True)
ode_obj.plot()
# Now for a pretty graph
from bokeh.plotting import figure, show, output_notebook
from bokeh.palettes import Accent4

output_notebook()
colors = Accent4

p = figure(title="Measles Example",
          y_axis_type="log"
          )

for i in range(solution.shape[1]):
    p.line(t, 
           solution[:,i], 
           line_color=colors[i],
           legend=states[i],
           line_width=2)

show(p)<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.08371888081075127, 'fraction_numerical': 0.013879709187045605, 'mean_word_length': 2.2244318181818183, 'pattern_counts': {'":': 0, '<': 17, '<?xml version=': 0, '>': 17, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 9, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '4', 'dir': 'jupyter-structured-clean-dedup', 'id': '46338', 'license': 'MIT', 'n_tokens_mistral': 1561, 'n_tokens_neox': 1400, 'n_words': 296, 'path': 'notebooks/vw10/4_3/python_pygom.ipynb', 'repo_name': 'epirecipes/epiRecipes', 'size': '4546'}