Symbolic Models

Provides an API to define epidemiological models in terms of sympy symbolic expressions.

class epipack.symbolic_epi_models.SymbolicEpiModel(compartments, initial_population_size=1, correct_for_dynamical_population_size=False)[source]

Bases: SymbolicMixin, EpiModel

Define a model based on the analytical framework offered by Sympy.

This class uses the event-based framework where state-change vectors are associated with event rates.

Parameters

compartments (list) -- A list of sympy.Symbol instances that symbolize compartments.
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 sum of all compartments, otherwise they will be divided by the initial population size.

compartments

A list of sympy.Symbol instances that symbolize compartments.

Type: list

N_comp

Number of compartments (including population number)

Type: int

parameter_values

Maps parameter symbols to numerical values,

Type: dict

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 symbolic expressions

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 sympy.Matrix

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 symbolic expressions

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 sympy.Matrix

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 symbolic expressions

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 sympy.Matrix

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

dydt()[source]

Compute the momenta of the epidemiological model as symbolic expressions.

Parameters

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

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_quadratic_events(event_list, allow_nonzero_column_sums=False, reset_events=True)[source]

Define the 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\).".

class epipack.symbolic_epi_models.SymbolicMixin[source]

Bases: object

Provides methods that are useful to both epipack.symbolic_epi_models.SymbolicEpiModel and epipack.symbolic_matrix_epi_models.SymbolicMatrixEpiModel

ODEs()[source]: Obtain the equations of motion for this model in form of equations.

ODEs_jupyter()[source]: Pretty-print the equations of motion for this model in a Jupyter notebook.

find_fixed_points()[source]: Find states in which this model is stationary (fixed points).

get_eigenvalues_at_disease_free_state(disease_free_state=None)[source]

Obtain the Jacobian's eigenvalues at the disease free state.

Parameters

disease_free_state (dict, default = None) --

A dictionary where a compartment symbol maps to an expression (the value of this compartment in the fixed point). If compartments are missing, it is implicitly assumed that this compartment has a value of zero.

If None, the disease_free_state is assumed to be at disease_free_state = { S: 1 }.

Returns

eigenvalues -- Each entry maps an eigenvalue expression to its multiplicity.

Return type

dict

get_eigenvalues_at_fixed_point(fixed_point_dict)[source]

Obtain the Jacobian's eigenvalues at a given fixed point.

Parameters: fixed_point_dict (dict) -- A dictionary where a compartment symbol maps to an expression (the value of this compartment in the fixed point). If compartments are missing, it is implicitly assumed that this compartment has a value of zero.
Returns: eigenvalues -- Each entry maps an eigenvalue expression to its multiplicity.
Return type: dict

get_jacobian_at_fixed_point(fixed_point_dict, simplify=True)[source]

Obtain the Jacobian at a given fixed point.

Parameters

fixed_point_dict (dict) -- A dictionary where a compartment symbol maps to an expression (the value of this compartment in the fixed point). If compartments are missing, it is implicitly assumed that this compartment has a value of zero.
simplify (bool) -- whether or not to let sympy try to simplify the expressions

Returns

J -- The Jacobian matrix at the given fixed point.

Return type

sympy.Matrix

get_numerical_dydt(lambdify_modules='numpy')[source]

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

Parameters

time_points (np.ndarray) -- An array of time points at which the compartment values should be evaluated and returned.
return_compartments (list, default = None) -- A list of compartments for which the result should be returned. If return_compartments is None, all compartments will be returned.
integrator (str, default = 'dopri5') -- Which method to use for integration. Currently supported are 'euler' and 'dopri5'. If 'euler' is chosen, \(\delta t\) will be determined by the difference of consecutive entries in time_points.
adopt_final_state (bool, default = False) -- Whether or not to adopt the final state of the integration

Returns

dydt -- A function dydt(t, y, *args, **kwargs) that returns the numerical momenta of this system at time t and state vector y.

Return type

func

get_numerical_event_and_rate_functions()[source]

Converts the symbolic event lists and corresponding symbolic rates to functions that return numeric event lists and numeric rates based on the current time and state vector.

This function is needed in the epipack.numeric_epi_models.EpiModel base class for stochastic simulations.

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.

jacobian(simplify=True)[source]

Obtain the Jacobian for this model.

Parameters: simplify (bool, default = True) -- If True, epipack will try to simplify the evaluated Jacobian. This might not be desirable in some cases due to its long evaluation time, which is why it can be turned off.

set_parameter_values(parameter_values)[source]

Set numerical values for the parameters of this model. This might include free symbols that are part of symbolized rate functions.

Parameters: parameter_values (dict) -- A dictionary mapping symbols to numerical values.

class epipack.symbolic_epi_models.SymbolicODEModel(ODEs)[source]

Bases: SymbolicEpiModel

Define a model purely based on a list of ODEs.

Parameters

ODEs (list) --

A list of symbolic ODEs in format

sympy.Eq(sympy.Derivative(Y, t), expr)

add_fission_processes(*args, **kwargs)[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(*args, **kwargs)[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(*args, **kwargs)[source]: Add linear events without resetting the existing event terms. See epipack.numeric_epi_models.EpiModel.set_linear_events() for docstring.

add_quadratic_events(*args, **kwargs)[source]: Add quadratic events without resetting the existing event terms. See epipack.numeric_epi_models.EpiModel.set_quadratic_events() for docstring.

add_transition_processes(*args, **kwargs)[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(*args, **kwargs)[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.

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()[source]: Return the momenta of the epidemiological model as symbolic expressions.

set_linear_events(*args, **kwargs)[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(*args, **kwargs)[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(*args, **kwargs)[source]

Define the 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(*args, **kwargs)[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
use_ivp_solver (bool, default = None) -- Wether or not to use an initial value problem solver to obtain a time leap in explicitly time-dependent problems. If None, will use the value of the class attribute self.use_ivp_solver.
rates_have_explicit_time_dependence (bool, default = None) -- Wether or not the problem is explicitly time-dependent. If None, will use the value of the class attribute self.rates_have_explicit_time_dependence.
ignore_warnings (bool, default = False) -- wether or not to raise warnings about unset explicit time.

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.symbolic_epi_models.SymbolicSEIRModel(infection_rate, latency_rate, recovery_rate, initial_population_size=1)[source]

Bases: SymbolicEpiModel

An SEIR model derived from epipack.symbolic_epi_models.SymbolicEpiModel.

class epipack.symbolic_epi_models.SymbolicSIModel(infection_rate, initial_population_size=1)[source]

Bases: SymbolicEpiModel

An SI model derived from epipack.symbolic_epi_models.SymbolicEpiModel.

class epipack.symbolic_epi_models.SymbolicSIRModel(infection_rate, recovery_rate, initial_population_size=1)[source]

Bases: SymbolicEpiModel

An SIR model derived from epipack.symbolic_epi_models.SymbolicEpiModel.

class epipack.symbolic_epi_models.SymbolicSIRSModel(infection_rate, recovery_rate, waning_immunity_rate, initial_population_size=1)[source]

Bases: SymbolicEpiModel

An SIRS model derived from epipack.symbolic_epi_models.SymbolicEpiModel.

class epipack.symbolic_epi_models.SymbolicSISModel(infection_rate, recovery_rate, initial_population_size=1)[source]

Bases: SymbolicEpiModel

An SIS model derived from epipack.symbolic_epi_models.SymbolicEpiModel.

epipack.symbolic_epi_models.get_temporal_interpolation(time_data, value_data, interpolation_degree=1)[source]

Obtain a symbolic piecewise function that interpolates between values given in value_data for the intervals defined in time_data, based on a spline interpolation of degree interpolation_degree. If interpolation_degree == 0, the function changes according to step functions. In this case time_data needs to have one value more than value_data.

The values in time_data and value_data can be symbols or numeric values.

Parameters

time_data (list) -- Sorted list of time values.
value_data (list) -- List of values corresponding to the times given in time_data.
interpolation_degree (int) -- The degree of the polynomial that interpolates between values.