Numeric Models¶

Provides an API to define epidemiological models.

class epipack.numeric_epi_models.ConstantBirthRate(rate)[source]¶

Bases: object

Will be used as a function of time t and state y, returning a rate value.

Parameters: rate (float) -- Constant rate value

class epipack.numeric_epi_models.ConstantLinearRate(rate, comp0)[source]¶

Bases: object

Will be used as a function of time t and state y, returning a rate value.

Parameters

rate (float) -- Constant rate value
comp0 (int) -- Index of the corresponding reacting component. The incidence of this component will be multiplied with the value of rate.

class epipack.numeric_epi_models.ConstantQuadraticRate(rate, comp0, comp1)[source]¶

Bases: object

Will be used as a function of time t and state y, returning a rate value.

Parameters

rate (float) -- Constant rate value
comp0 (int) -- Index of one of the reacting components. The incidence of this component will be multiplied with the value of rate.
comp1 (int) -- Index of the other reacting component. The incidence of this component will be multiplied with the value of rate.

class epipack.numeric_epi_models.DynamicBirthRate(rate)[source]¶

Bases: object

Will be used as a function of time t and state y, returning a rate value.

Parameters: rate (function) -- Function of time t and state y

class epipack.numeric_epi_models.DynamicLinearRate(rate, comp0)[source]¶

Bases: object

Will be used as a function of time t and state y, returning a rate value.

Parameters

rate (function) -- Function of time t and state y
comp0 (int) -- Index of the corresponding reacting component. The incidence of this component will be multiplied with the value of rate.

class epipack.numeric_epi_models.DynamicQuadraticRate(rate, comp0, comp1)[source]¶

Bases: object

Will be used as a function of time t and state y, returning a rate value.

Parameters

rate (function) -- Function of time t and state y
comp0 (int) -- Index of one of the reacting components. The incidence of this component will be multiplied with the value of rate.
comp1 (int) -- Index of the other reacting component. The incidence of this component will be multiplied with the value of rate.

class epipack.numeric_epi_models.EpiModel(compartments, initial_population_size=1, correct_for_dynamical_population_size=False, integral_solver='solve_ivp')[source]¶

Bases: epipack.integrators.IntegrationMixin

A general class to define a standard mean-field compartmental epidemiological model, based on reaction events.

Parameters

compartments (list of string) -- A list containing compartment strings.
initial_population_size (float, default = 1.0) -- The population size at \(t = 0\).
correct_for_dynamical_population_size (bool, default = False) -- If True, the quadratic coupling terms will be divided by the population size.
integral_solver (str, default = 'solve_ivp') -- Whether or not to use the initial-value solver solve_ivp. to determine a time leap for time-varying rates. If not 'solve_ivp', a Newton-Raphson method will be used on the upper bound of a quad-integrator.

compartments¶

A list containing strings or hashable types that describe each compartment, (e.g. "S", "I", etc.).

Type: list

compartment_ids¶

Maps compartments to their indices in self.compartments.

Type: dict

N_comp¶

Number of compartments (including population number)

Type: int

initial_population_size¶

The population size at \(t = 0\).

Type: float

correct_for_dynamical_population_size¶

If True, the quadratic coupling terms will be divided by the sum of all compartments, otherwise they will be divided by the initial population size.

Type: bool

birth_rate_functions¶

A list of functions that return rate values based on time t and state vector y. Each entry corresponds to an event update in self.birth_event_updates.

Type: list of ConstantBirthRate or DynamicBirthRate

birth_event_updates¶

A list of vectors. Each entry corresponds to a rate in birth_rate_functions and quantifies the change in individual counts in the compartments.

Type: list of numpy.ndarray

linear_rate_functions¶

A list of functions that return rate values based on time t and state vector y. Each entry corresponds to an event update in self.linear_event_updates.

Type: list of ConstantLinearRate or DynamicLinearRate

linear_event_updates¶

A list of vectors. Each entry corresponds to a rate in linear_rate_functions and quantifies the change in individual counts in the compartments.

Type: list of numpy.ndarray

quadratic_rate_functions¶

A list of functions that return rate values based on time t and state vector y. Each entry corresponds to an event update in self.quadratic_event_updates.

Type: list of ConstantQuadraticRate or DynamicQuadraticRate

quadratic_event_updates¶

A list of vectors. Each entry corresponds to a rate in quadratic_rate_functions and quantifies the change in individual counts in the compartments.

Type: list of numpy.ndarray

y0¶

The initial conditions.

Type: numpy.ndarray

rates_have_explicit_time_dependence¶

Internal switch that's flipped when a non-constant rate is passed to the model.

Type: bool

use_ivp_solver¶

Whether or not to use the initial-value solver to determine a time leap for time-varying rates. If False, a Newton-Raphson method will be used on the upper bound of a quad-integrator.

Type: bool

Example

>>> epi = EpiModel(["S","I","R"])
>>> print(epi.compartments)
[ "S", "I", "R" ]

add_fission_processes(process_list)[source]¶

Define linear fission processes between compartments.

Parameters

process_list (list of tuple) --

A list of tuples that contains fission rates in the following format:

[
    ("source_compartment", rate, "target_compartment_0", "target_compartment_1" ),
    ...
]

Example

For pure exponential growth of compartment B.

epi.add_fission_processes([
    ("B", growth_event, "B", "B" ),
])

add_fusion_processes(process_list)[source]¶

Define fusion processes between compartments.

Parameters

process_list (list of tuple) --

A list of tuples that contains fission rates in the following format:

[
    ("coupling_compartment_0", "coupling_compartment_1", rate, "target_compartment_0" ),
    ...
]

Example

Fusion of reactants "A", and "B" to form "C".

epi.add_fusion_processes([
    ("A", "B", reaction_rate, "C" ),
])

add_linear_events(event_list, allow_nonzero_column_sums=False)[source]¶: Add linear events without resetting the existing event terms. See epipack.numeric_epi_models.EpiModel.set_linear_events() for docstring.

add_quadratic_events(event_list, allow_nonzero_column_sums=False)[source]¶: Add quadratic events without resetting the existing event terms. See epipack.numeric_epi_models.EpiModel.set_quadratic_events() for docstring.

add_transition_processes(process_list)[source]¶

Define the linear transition processes between compartments.

Parameters

process_list (list of tuple) --

A list of tuples that contains transitions events in the following format:

[
    ( source_compartment, rate, target_compartment ),
    ...
]

Example

For an SEIR model.

epi.add_transition_processes([
    ("E", symptomatic_rate, "I" ),
    ("I", recovery_rate, "R" ),
])

add_transmission_processes(process_list)[source]¶

A wrapper to define quadratic process rates through transmission reaction equations. Note that in stochastic network/agent simulations, the transmission rate is equal to a rate per link. For the mean-field ODEs, the rates provided to this function will just be equal to the prefactor of the respective quadratic terms.

For instance, if you analyze an SIR system and simulate on a network of mean degree \(k_0\), a basic reproduction number \(R_0\), and a recovery rate \(\mu\), you would define the single link transmission process as

("I", "S", R_0/k_0 * mu, "I", "I")

For the mean-field system here, the corresponding reaction equation would read

("I", "S", R_0 * mu, "I", "I")

Parameters

process_list (list of tuple) --

A list of tuples that contains transitions rates in the following format:

[
    ("source_compartment",
     "target_compartment_initial",
     rate
     "source_compartment",
     "target_compartment_final",
     ),
    ...
]

Example

For an SEIR model.

epi.add_transmission_processes([
    ("I", "S", +1, "I", "E" ),
])

dydt(t, y)[source]¶

Compute the current momenta of the epidemiological model.

Parameters

t (float) -- Current time
y (numpy.ndarray) -- The entries correspond to the compartment frequencies (or counts, depending on population size).

get_compartment(iC)[source]¶: Get the compartment, given an integer ID iC

get_compartment_changes(rates)[source]¶

Sample a state change vector with probability proportional to its rate in rates.

Needed for stochastic simulations.

Parameters: rates (numpy.ndarray) -- A non-zero list of rates. Expects rates to be sorted according to self.birth_event_updates + self.linear_event_updates + self.quadratic_event_updates.
Returns: dy -- A state change vector.
Return type: numpy.ndarray

get_compartment_id(C)[source]¶: Get the integer ID of a compartment C

get_event_rates(t, y)[source]¶

Get a list of rate values corresponding to the previously set events.

Parameters

t (float) -- Current time
y (numpy.ndarray) -- Current state vector

Returns

rates -- A list of rate values corresponding to rates. Ordered as birth_rate_functions + linear_rate_functions + quadratic_rate_functions.

Return type

list

get_numerical_dydt()[source]¶: Return a function that obtains t and y as an input and returns dydt of this system

get_numerical_event_and_rate_functions()[source]¶

This function is needed to generalize stochastic simulations for child classes.

Returns

get_event_rates (func) -- A function that takes the current time t and state vector y and returns numerical event rate lists.
get_compartment_changes (funx) -- A function that takes a numerical list of event rates and returns a random event state change vector with probability proportional to its entry in rates.

get_time_leap_and_proposed_compartment_changes(t, current_event_rates=None, get_event_rates=None, get_compartment_changes=None)[source]¶

For the current event rates, obtain a proposed time leap and concurrent state change vector.

This method is needed for stochastic simulations.

Parameters

t (float) -- current time
current_event_rates (list, default = None) -- A list of constant rate values. Will be ignored if self.rates_have_explicit_time_dependence is True, which is why None is a valid value.
get_event_rates (function, default = None) -- A function that takes time t and current state y as input and computes the rates of all possible events. If None, will attempt to set this to self.get_event_rates().
get_compartment_changes (function, default = None) -- A function that takes computed event rates and returns a random state change with probability proportional to its rate. If None, will attempt to set this to self.get_compartment_changes().

Returns

tau (float) -- A time leap.
dy (numpy.ndarray) -- A state change vector.

set_linear_events(event_list, allow_nonzero_column_sums=False, reset_events=True)[source]¶

Define the linear transition events between compartments.

Parameters

event_list (list of tuple) --

A list of tuples that contains transition events in the following format:

[
    (
        ("affected_compartment_0",),
        rate,
        [
            ("affected_compartment_0", dN0),
            ("affected_compartment_1", dN1),
            ...
        ],
     ),
    ...
]

allow_nonzero_column_sums (bool, default : False) -- Traditionally, epidemiological models preserve the total population size. If that's not the case, switch off testing for this.
reset_events (bool, default : True) -- Whether to reset all linear events to zero before converting those.

Example

For an SEIR model with infectious period tau and incubation period theta.

epi.set_linear_events([
    ( ("E",),
      1/theta,
      [ ("E", -1), ("I", +1) ]
    ),
    ( ("I",),
      1/tau,
      [ ("I", -1), ("R", +1) ]
    ),
])

Read as "compartment E reacts with rate \(1/\theta\) which leads to the decay of one E particle to one I particle."

set_processes(process_list, allow_nonzero_column_sums=False, reset_events=True, ignore_rate_position_checks=False)[source]¶

Converts a list of reaction process tuples to event tuples and sets the rates for this model.

Parameters

process_list (list of tuple) --

A list containing reaction processes in terms of tuples.

[
    # transition process
    ( source_compartment, rate, target_compartment),

    # transmission process
    ( coupling_compartment_0, coupling_compartment_1, rate, target_compartment_0, target_ccompartment_1),

    # fission process
    ( source_compartment, rate, target_compartment_0, target_ccompartment_1),

    # fusion process
    ( source_compartment_0, source_compartment_1, rate, target_compartment),

    # death process
    ( source_compartment, rate, None),

    # birth process
    ( None, rate, target_compartment),
]

allow_nonzero_column_sums (bool, default : False) -- Traditionally, epidemiological models preserve the total population size. If that's not the case, switch off testing for this.
reset_events (bool, default : True) -- If this is True, reset all events to zero before setting the new ones.
ignore_rate_position_checks (bool, default = False) -- This function usually checks whether the rate of a reaction is positioned correctly. You can turn this behavior off for transition, birth, death, and transmission processes. (Useful if you want to define symbolic transmission processes that are compartment-dependent).

set_quadratic_events(event_list, allow_nonzero_column_sums=False, reset_events=True)[source]¶

Define quadratic transition events between compartments.

Parameters

event_list (list of tuple) --

A list of tuples that contains transmission events in the following format:

[
    (
        ("coupling_compartment_0", "coupling_compartment_1"),
        rate,
        [
            ("affected_compartment_0", dN0),
            ("affected_compartment_1", dN1),
            ...
        ],
     ),
    ...
]

allow_nonzero_column_sums (bool, default : False) -- Traditionally, epidemiological models preserve the total population size. If that's not the case, switch off testing for this.
reset_events (bool, default : True) -- Whether to reset all linear events to zero before converting those.

Example

For an SEIR model with infection rate eta.

epi.set_quadratic_events([
    ( ("S", "I"),
      eta,
      [ ("S", -1), ("E", +1) ]
    ),
])

Read as

"Coupling of S and I leads to the decay of one S particle to one E particle with rate \(\eta\).".

simulate(tmax, return_compartments=None, sampling_dt=None, sampling_callback=None, adopt_final_state=False)[source]¶

Returns values of the given compartments at the demanded time points (as a numpy.ndarray of shape (return_compartments), len(time_points).

If return_compartments is None, all compartments will be returned.

Parameters

tmax (float) -- maximum length of the simulation
return_compartments (list of compartments, default = None:) -- The compartments for which to return time series. If None, all compartments will be returned.
sampling_dt (float, default = None) -- Temporal distance between samples of the compartment counts. If None, every change will be returned.
sampling_callback (funtion, default = None) -- A function that's called when a sample is taken

Returns

t (numpy.ndarray) -- times at which compartment counts have been sampled
result (dict) -- Dictionary mapping a compartment to a time series of its count.

class epipack.numeric_epi_models.SEIRModel(infection_rate, recovery_rate, symptomatic_rate, initial_population_size=1.0)[source]¶

Bases: epipack.numeric_epi_models.EpiModel

An SEIR model derived from epipack.numeric_epi_models.EpiModel.

class epipack.numeric_epi_models.SIModel(infection_rate, initial_population_size=1.0)[source]¶

Bases: epipack.numeric_epi_models.EpiModel

An SI model derived from epipack.numeric_epi_models.EpiModel.

class epipack.numeric_epi_models.SIRModel(infection_rate, recovery_rate, initial_population_size=1.0)[source]¶

Bases: epipack.numeric_epi_models.EpiModel

An SIR model derived from epipack.numeric_epi_models.EpiModel.

class epipack.numeric_epi_models.SIRSModel(infection_rate, recovery_rate, waning_immunity_rate, initial_population_size=1.0)[source]¶

Bases: epipack.numeric_epi_models.EpiModel

An SIRS model derived from epipack.numeric_epi_models.EpiModel.

class epipack.numeric_epi_models.SISModel(infection_rate, recovery_rate, initial_population_size=1.0)[source]¶

Bases: epipack.numeric_epi_models.EpiModel

An SIS model derived from epipack.numeric_epi_models.EpiModel.

Parameters

R0 (float) -- The basic reproduction number. From this, the infection rate is computed as infection_rate = R0 * recovery_rate
recovery_rate (float) -- Recovery rate
population_size (float, default = 1.0) -- Number of people in the population.

epipack.numeric_epi_models.custom_choice(p)[source]¶: Return an index of normalized probability vector p with probability equal to the corresponding entry in p.