Epidemiology

Epidemiology: compartmental, network, metapopulation, and agent-based systems of infection and recovery.

ActiveHostElement(label, visual=None, *, i=0.0, a=0.1, lam=0.3, mu=0.1)

Bases: Element

Node of an activity-driven network, susceptible at i = 0 or infected at i = 1, in discrete time.

A TransitionInteraction recovers an infected node with probability mu in a step, and another infects a susceptible node with probability lam by every infected neighbour along the edges of the previous step. A TopologyInteraction joins the node to other nodes with the probability a, its activity.

Parameters:

Name	Type	Description	Default
label	str	The label of the node or of the group of nodes.	required
visual	VisualSpec | None	How the element is drawn, or None for default_visual_.	None
i	Annotated[float, StateSpec(domain=closed_open(0, oo))]	Whether the node is infected.	0.0
a	Annotated[float, ParameterSpec(domain=closed_open(0, oo))]	The activity, the probability of firing at a step.	0.1
lam	Annotated[float, ParameterSpec(domain=closed_open(0, oo))]	The probability of infection along one edge.	0.3
mu	Annotated[float, ParameterSpec(domain=closed_open(0, oo))]	The probability of recovery at a step.	0.1

Examples:

>>> from skcomplex.systems.biology.epidemiology import ActiveHostElement
>>> ActiveHostElement('node').mu
0.1

Source code in src/skcomplex/systems/biology/epidemiology/_adaptive.py

def __init__(
    self,
    label: str,
    visual: VisualSpec | None = None,
    *,
    i: Annotated[float, StateSpec(domain=RealSet.closed_open(0, oo))] = 0.0,
    a: Annotated[float, ParameterSpec(domain=RealSet.closed_open(0, oo))] = 0.1,
    lam: Annotated[float, ParameterSpec(domain=RealSet.closed_open(0, oo))] = 0.3,
    mu: Annotated[float, ParameterSpec(domain=RealSet.closed_open(0, oo))] = 0.1,
) -> None:
    super().__init__(label, visual=visual)
    self.i = i
    self.a = a
    self.lam = lam
    self.mu = mu

CarrierElement(label, visual=None, *, p=0.05, lam=0.2)

Bases: Element

Node of a network whose state is the probability that it is infected, in quenched mean field.

The node recovers at rate 1 and every infected neighbour infects it at rate lam, with the probabilities of the neighbours taken as independent, so dp_i / dt = -p_i + lam (1 - p_i) sum_j A_ij p_j over the neighbours j.

Parameters:

Name	Type	Description	Default
label	str	The label of the node or of the group of nodes.	required
visual	VisualSpec | None	How the element is drawn, or None for default_visual_.	None
p	Annotated[float, StateSpec(domain=closed_open(0, oo))]	The initial probability of being infected.	0.05
lam	Annotated[float, ParameterSpec(domain=closed_open(0, oo))]	The spreading rate, the ratio of the infection rate per edge and the recovery rate.	0.2

Examples:

>>> from skcomplex.systems.biology.epidemiology import CarrierElement
>>> CarrierElement('node').lam
0.2

Source code in src/skcomplex/systems/biology/epidemiology/_adaptive.py

def __init__(
    self,
    label: str,
    visual: VisualSpec | None = None,
    *,
    p: Annotated[float, StateSpec(domain=RealSet.closed_open(0, oo))] = 0.05,
    lam: Annotated[float, ParameterSpec(domain=RealSet.closed_open(0, oo))] = 0.2,
) -> None:
    super().__init__(label, visual=visual)
    self.p = p
    self.lam = lam

ColocationInteraction(label, *, people, p=0.05, gamma=0.1, closed=False)

Bases: Interaction

Daily movement of people between their homes and workplaces, and infection between people at one location.

A step is half a day: in even steps every person is at home and in odd steps at the workplace, unless workplaces are closed, the time-labelled visits of the bipartite graph of people and locations of EpiSims. At every step a susceptible person is infected with the probability 1 - (1 - p)^k by the k infected people at the same location, and an infected person recovers with the probability gamma.

Parameters:

Name	Type	Description	Default
label	str	The label of the interaction.	required
people	Annotated[str, RoleSpec()]	The path of the agent group of people.	required
p	float	The probability of transmission per infected contact and step.	0.05
gamma	float	The probability of recovery per step.	0.1
closed	bool	Whether the workplaces are closed, so people stay at home.	False

Examples:

>>> from skcomplex.systems.biology.epidemiology import ColocationInteraction
>>> ColocationInteraction('contacts', people='people', p=0.02).p
0.02

Source code in src/skcomplex/systems/biology/epidemiology/_cities.py

def __init__(
    self,
    label: str,
    *,
    people: Annotated[str, RoleSpec()],
    p: float = 0.05,
    gamma: float = 0.1,
    closed: bool = False,
) -> None:
    super().__init__(label)
    self.people = people
    self.p = p
    self.gamma = gamma
    self.closed = closed

visit(system)

Move every person to the location of the half day, then infect and recover.

Parameters:

Name	Type	Description	Default
system	SystemView	The system.	required

Source code in src/skcomplex/systems/biology/epidemiology/_cities.py

def visit(self, system: SystemView) -> None:
    """Move every person to the location of the half day, then infect and recover.

    Args:
        system:
            The system.
    """
    people = system.agents(self.people, CommuterView)
    at_work = round(system.time) % 2 == 1 and not self.closed
    places = np.array([int(person.work if at_work else person.home) for person in people])
    health = np.array([person.health for person in people])
    infected = np.bincount(places[health == INFECTED], minlength=places.max() + 1)[places]
    draws = system.rng.random(len(people))
    infection = (health == SUSCEPTIBLE) & (draws < 1 - (1 - self.p) ** infected)
    recovery = (health == INFECTED) & (draws < self.gamma)
    health = np.where(infection, INFECTED, np.where(recovery, RECOVERED, health))
    for person, place, state in zip(people, places, health, strict=True):
        person.move_to(int(place))
        person.health = float(state)

CommuterElement(label, visual=VisualSpec('agents', {'color_by': 'health', 'color_range': (0.0, 2.0)}), *, home=0.0, work=0.0, health=0.0)

Bases: Element

Person of a synthetic city who sleeps at home and spends the day at a workplace, with an SIR health state.

The home and the workplace are the locations the person visits, the edges of the person in the bipartite graph of people and locations. Health is 0 susceptible, 1 infected, or 2 recovered.

Parameters:

Name	Type	Description	Default
label	str	The label of the person or of the group of people.	required
visual	VisualSpec | None	How the element is drawn, or None for default_visual_.	VisualSpec('agents', {'color_by': 'health', 'color_range': (0.0, 2.0)})
home	Annotated[float, StateSpec(domain=NN)]	The location of the home.	0.0
work	Annotated[float, StateSpec(domain=NN)]	The location of the workplace or school.	0.0
health	Annotated[float, StateSpec(domain=NN)]	The health state.	0.0

Examples:

>>> from skcomplex.systems.biology.epidemiology import CommuterElement
>>> CommuterElement('person', home=3.0, work=12.0).work
12.0

Source code in src/skcomplex/systems/biology/epidemiology/_cities.py

def __init__(
    self,
    label: str,
    visual: VisualSpec | None = VisualSpec('agents', {'color_by': 'health', 'color_range': (0.0, 2.0)}),
    *,
    home: Annotated[float, StateSpec(domain=NN)] = 0.0,
    work: Annotated[float, StateSpec(domain=NN)] = 0.0,
    health: Annotated[float, StateSpec(domain=NN)] = 0.0,
) -> None:
    super().__init__(label, visual=visual)
    self.home = home
    self.work = work
    self.health = health

DegreeClassElement(label, visual=None, *, rho=0.01, lam=0.2, k=1.0, P=1.0)

Bases: Element

Infected density of the nodes of one degree of a network, in heterogeneous mean field.

A node of degree k recovers at rate 1 and is infected at rate lam k Theta, where Theta = sum_k k P(k) rho_k / <k> is the probability that an edge points to an infected node, a total over the degree classes, so d rho_k / dt = -rho_k + lam k (1 - rho_k) Theta. The classes are the members of a group, and the recovery and the infection are interactions of the group, which heterogeneous_mean_field_sis adds.

Parameters:

Name	Type	Description	Default
label	str	The label of the class or the group of classes.	required
visual	VisualSpec | None	How the element is drawn, or None for default_visual_.	None
rho	Annotated[float, StateSpec(domain=closed_open(0, oo))]	The initial infected density.	0.01
lam	Annotated[float, ParameterSpec(domain=closed_open(0, oo))]	The spreading rate, the ratio of the infection and recovery rates.	0.2
k	Annotated[float, ParameterSpec(domain=open(0, oo))]	The degree of the class.	1.0
P	Annotated[float, ParameterSpec(domain=closed_open(0, oo))]	The probability P(k) of the degree.	1.0

Examples:

>>> from skcomplex.core.components.elements import GroupElement
>>> from skcomplex.systems.biology.epidemiology import DegreeClassElement
>>> GroupElement(
...     'network', element=DegreeClassElement, size=2, values={'k': [1.0, 2.0], 'P': [0.5, 0.5]}
... ).prepare().size_
2

Source code in src/skcomplex/systems/biology/epidemiology/_networks.py

def __init__(
    self,
    label: str,
    visual: VisualSpec | None = None,
    *,
    rho: Annotated[float, StateSpec(domain=RealSet.closed_open(0, oo))] = 0.01,
    lam: Annotated[float, ParameterSpec(domain=RealSet.closed_open(0, oo))] = 0.2,
    k: Annotated[float, ParameterSpec(domain=RealSet.open(0, oo))] = 1.0,
    P: Annotated[float, ParameterSpec(domain=RealSet.closed_open(0, oo))] = 1.0,
) -> None:
    super().__init__(label, visual=visual)
    self.rho = rho
    self.lam = lam
    self.k = k
    self.P = P

EpidemicCenterElement(label, visual=None, *, x=0.0, a=2.2, recovery=6, e=0.005)

Bases: Element

City in a discrete-time SIRS epidemic, where an infection lasts one step and immunity recovery steps.

The state keeps the infected fractions of the last 12 steps, so the recovery time is at most 12.

Parameters:

Name	Type	Description	Default
label	str	The label of the center.	required
visual	VisualSpec | None	How the element is drawn, or None for default_visual_.	None
x	Annotated[float, StateSpec(domain=closed_open(0, oo), lags=MAX_RECOVERY, initial=(0.0,) * MAX_RECOVERY)]	The initial fraction of infected individuals.	0.0
a	Annotated[float, ParameterSpec(domain=closed_open(0, oo))]	The rate of transmission per contact.	2.2
recovery	Annotated[int, ParameterSpec(domain=NN)]	The number of steps a recovered individual stays immune, at most 12.	6
e	Annotated[float, ParameterSpec(domain=closed_open(0, oo))]	The fraction of the infected of the center who visit another center.	0.005

Examples:

>>> from skcomplex.systems.biology.epidemiology import EpidemicCenterElement
>>> EpidemicCenterElement('city', x=0.01).recovery
6

Source code in src/skcomplex/systems/biology/epidemiology/_recurrent.py

def __init__(
    self,
    label: str,
    visual: VisualSpec | None = None,
    *,
    x: Annotated[
        float,
        StateSpec(domain=RealSet.closed_open(0, oo), lags=MAX_RECOVERY, initial=(0.0,) * MAX_RECOVERY),
    ] = 0.0,
    a: Annotated[float, ParameterSpec(domain=RealSet.closed_open(0, oo))] = 2.2,
    recovery: Annotated[int, ParameterSpec(domain=NN)] = 6,
    e: Annotated[float, ParameterSpec(domain=RealSet.closed_open(0, oo))] = 0.005,
) -> None:
    super().__init__(label, visual=visual)
    self.x = x
    self.a = a
    self.recovery = recovery
    self.e = e

HostElement(label, visual=None, *, i=0.0, w=1.0)

Bases: Element

Node of a network that is susceptible at i = 0 or infected at i = 1, in continuous time.

A SpreadingInteraction infects and recovers it, and a TopologyInteraction rewires its edges at its rate w.

Parameters:

Name	Type	Description	Default
label	str	The label of the node or of the group of nodes.	required
visual	VisualSpec | None	How the element is drawn, or None for default_visual_.	None
i	Annotated[float, StateSpec(domain=closed_open(0, oo))]	Whether the node is infected.	0.0
w	Annotated[float, ParameterSpec(domain=closed_open(0, oo))]	The rate at which the node, while susceptible, rewires an edge to an infected neighbour.	1.0

Examples:

>>> from skcomplex.systems.biology.epidemiology import HostElement
>>> HostElement('node').w
1.0

Source code in src/skcomplex/systems/biology/common/networks/_host.py

def __init__(
    self,
    label: str,
    visual: VisualSpec | None = None,
    *,
    i: Annotated[float, StateSpec(domain=RealSet.closed_open(0, oo))] = 0.0,
    w: Annotated[float, ParameterSpec(domain=RealSet.closed_open(0, oo))] = 1.0,
) -> None:
    super().__init__(label, visual=visual)
    self.i = i
    self.w = w

ImmunityLossInteraction(label, *, recovered, susceptible, xi=0.05, stochastic=False)

Bases: ReactionInteraction

Loss of immunity, the reaction R -> S at the rate xi R from recovered back to susceptible.

It is a ReactionInteraction between the states n of the two compartments its roles bind, a rate equation, or in stochastic mode a random event that moves one individual.

Parameters:

Name	Type	Description	Default
label	str	The label of the interaction.	required
recovered	Annotated[str, RoleSpec()]	The path of the recovered compartment.	required
susceptible	Annotated[str, RoleSpec()]	The path of the susceptible compartment.	required
xi	Annotated[float, ParameterSpec(domain=closed_open(0, oo))]	The rate at which immunity wanes.	0.05
stochastic	bool	Whether losses of immunity are random events.	False

Examples:

>>> from skcomplex.systems.biology.epidemiology import ImmunityLossInteraction
>>> ImmunityLossInteraction('waning', recovered='R', susceptible='S', xi=0.05).xi
0.05

Source code in src/skcomplex/systems/biology/epidemiology/_factories.py

def __init__(
    self,
    label: str,
    *,
    recovered: Annotated[str, RoleSpec()],
    susceptible: Annotated[str, RoleSpec()],
    xi: Annotated[float, ParameterSpec(domain=RealSet.closed_open(0, oo))] = 0.05,
    stochastic: bool = False,
) -> None:
    super().__init__(label)
    self.recovered = recovered
    self.susceptible = susceptible
    self.xi = xi
    self.stochastic = stochastic

ImmunizingDegreeClassElement(label, visual=None, *, s=0.999, i=0.001, r=0.0, beta=0.1, mu=1.0, k=1.0, P=1.0)

Bases: Element

Susceptible, infected, and recovered densities of the nodes of one degree, in heterogeneous mean field.

A susceptible node of degree k is infected at rate beta k theta and an infected node recovers at rate mu, so d i_k / dt = -mu i_k + beta k s_k theta and d r_k / dt = mu i_k. On an uncorrelated network, the probability that a contact of a susceptible node is infected is theta = sum_k (k - 1) P(k) i_k / <k>, a total over the degree classes, since an infected neighbour of degree k has at most k - 1 other contacts. The classes are the members of a group, and the infections and the recoveries are interactions of the group, which heterogeneous_mean_field_sir adds.

Parameters:

Name	Type	Description	Default
label	str	The label of the class or the group of classes.	required
visual	VisualSpec | None	How the element is drawn, or None for default_visual_.	None
s	Annotated[float, StateSpec(domain=closed_open(0, oo))]	The initial susceptible density.	0.999
i	Annotated[float, StateSpec(domain=closed_open(0, oo))]	The initial infected density.	0.001
r	Annotated[float, StateSpec(domain=closed_open(0, oo))]	The initial recovered density.	0.0
beta	Annotated[float, ParameterSpec(domain=closed_open(0, oo))]	The rate of transmission along one contact.	0.1
mu	Annotated[float, ParameterSpec(domain=open(0, oo))]	The recovery rate.	1.0
k	Annotated[float, ParameterSpec(domain=open(0, oo))]	The degree of the class.	1.0
P	Annotated[float, ParameterSpec(domain=closed_open(0, oo))]	The probability P(k) of the degree.	1.0

Examples:

>>> from skcomplex.systems.biology.epidemiology import ImmunizingDegreeClassElement
>>> GroupElement(
...     'network',
...     element=ImmunizingDegreeClassElement,
...     size=3,
...     values={'k': [1.0, 2.0, 3.0], 'P': [0.5, 0.3, 0.2]},
... ).size
3

Source code in src/skcomplex/systems/biology/epidemiology/_outbreaks.py

def __init__(
    self,
    label: str,
    visual: VisualSpec | None = None,
    *,
    s: Annotated[float, StateSpec(domain=RealSet.closed_open(0, oo))] = 0.999,
    i: Annotated[float, StateSpec(domain=RealSet.closed_open(0, oo))] = 0.001,
    r: Annotated[float, StateSpec(domain=RealSet.closed_open(0, oo))] = 0.0,
    beta: Annotated[float, ParameterSpec(domain=RealSet.closed_open(0, oo))] = 0.1,
    mu: Annotated[float, ParameterSpec(domain=RealSet.open(0, oo))] = 1.0,
    k: Annotated[float, ParameterSpec(domain=RealSet.open(0, oo))] = 1.0,
    P: Annotated[float, ParameterSpec(domain=RealSet.closed_open(0, oo))] = 1.0,
) -> None:
    super().__init__(label, visual=visual)
    self.s = s
    self.i = i
    self.r = r
    self.beta = beta
    self.mu = mu
    self.k = k
    self.P = P

InfectionInteraction(label, *, susceptible, infected, beta=0.001, stochastic=False)

Bases: ReactionInteraction

Mass-action infection, the reaction S + I -> 2 I at the rate beta S I from susceptible to infected.

It is a ReactionInteraction between the states n of the two compartments its roles bind, a rate equation, or in stochastic mode a random event that moves one individual.

Parameters:

Name	Type	Description	Default
label	str	The label of the interaction.	required
susceptible	Annotated[str, RoleSpec()]	The path of the susceptible compartment.	required
infected	Annotated[str, RoleSpec()]	The path of the infected compartment.	required
beta	Annotated[float, ParameterSpec(domain=closed_open(0, oo))]	The transmission rate.	0.001
stochastic	bool	Whether infections are random events.	False

Examples:

>>> from skcomplex.core.components.elements import CompartmentElement
>>> from skcomplex.core.components.system import System
>>> from skcomplex.core.engine.runtime import Simulation
>>> from skcomplex.systems.biology.epidemiology import InfectionInteraction
>>> system = System(
...     elements=[CompartmentElement('S', n=99.0), CompartmentElement('I', n=1.0)],
...     interactions=[InfectionInteraction('infection', susceptible='S', infected='I', beta=0.01)],
... )
>>> trajectory = Simulation(system, dt=0.1).prepare().run(until=5.0).trajectory_
>>> round(float(trajectory.values('S.n')[-1] + trajectory.values('I.n')[-1]), 6)
100.0

Source code in src/skcomplex/systems/biology/epidemiology/_factories.py

def __init__(
    self,
    label: str,
    *,
    susceptible: Annotated[str, RoleSpec()],
    infected: Annotated[str, RoleSpec()],
    beta: Annotated[float, ParameterSpec(domain=RealSet.closed_open(0, oo))] = 0.001,
    stochastic: bool = False,
) -> None:
    super().__init__(label)
    self.susceptible = susceptible
    self.infected = infected
    self.beta = beta
    self.stochastic = stochastic

PartnershipEpidemicElement(label, visual=None, *, S=0.99, E=0.0, I=0.01, R=0.0, J=0.0, c=2.0, beta=0.2, tau=1.0, T=0.5, nu=0.1)

Bases: Element

Population of Watts and May whose infections pass along partnerships that last, with a fixed exposed period.

Susceptible people S acquire new partners at the rate c, each infecting with the probability beta over the partnership, and are infected at the rate pi = c beta (I / A + J) with the active population A = S + E + I. The first term is the infection by a new partner who is infectious, and the memory J is the infection by a current partner who turned infectious after the partnership formed, dJ/dt = pi(t - T) S(t - T) / A - J / tau, where partnerships last a mean time tau. The exposed E turn infectious after the fixed time T, and the infectious I develop the disease, which ends their sexual activity, at the rate nu, into R. The infection starts at time zero, so the infectious and the memory are zero before it.

Parameters:

Name	Type	Description	Default
label	str	The label of the population.	required
visual	VisualSpec | None	How the element is drawn, or None for default_visual_.	None
S	Annotated[float, StateSpec(domain=FRACTION)]	The initial fraction of susceptible people.	0.99
E	Annotated[float, StateSpec(domain=FRACTION)]	The initial fraction of exposed people.	0.0
I	Annotated[float, StateSpec(domain=FRACTION, history=0.0)]	The initial fraction of infectious people.	0.01
R	Annotated[float, StateSpec(domain=FRACTION)]	The initial fraction of people who developed the disease.	0.0
J	Annotated[float, StateSpec(domain=closed_open(0, oo), history=0.0)]	The initial infection memory of current partnerships.	0.0
c	Annotated[float, ParameterSpec(domain=closed_open(0, oo))]	The rate of acquiring new partners.	2.0
beta	Annotated[float, ParameterSpec(domain=FRACTION)]	The probability of infection by an infectious partner over a partnership.	0.2
tau	Annotated[float, ParameterSpec(domain=open(0, oo))]	The mean duration of a partnership.	1.0
T	Annotated[float, ParameterSpec(domain=open(0, oo))]	The duration of the exposed period.	0.5
nu	Annotated[float, ParameterSpec(domain=closed_open(0, oo))]	The rate at which the infectious develop the disease.	0.1

Examples:

>>> from skcomplex.systems.biology.epidemiology import PartnershipEpidemicElement
>>> PartnershipEpidemicElement('population').tau
1.0

Source code in src/skcomplex/systems/biology/epidemiology/_concurrency.py

def __init__(
    self,
    label: str,
    visual: VisualSpec | None = None,
    *,
    S: Annotated[float, StateSpec(domain=FRACTION)] = 0.99,
    E: Annotated[float, StateSpec(domain=FRACTION)] = 0.0,
    I: Annotated[float, StateSpec(domain=FRACTION, history=0.0)] = 0.01,
    R: Annotated[float, StateSpec(domain=FRACTION)] = 0.0,
    J: Annotated[float, StateSpec(domain=RealSet.closed_open(0, oo), history=0.0)] = 0.0,
    c: Annotated[float, ParameterSpec(domain=RealSet.closed_open(0, oo))] = 2.0,
    beta: Annotated[float, ParameterSpec(domain=FRACTION)] = 0.2,
    tau: Annotated[float, ParameterSpec(domain=RealSet.open(0, oo))] = 1.0,
    T: Annotated[float, ParameterSpec(domain=RealSet.open(0, oo))] = 0.5,
    nu: Annotated[float, ParameterSpec(domain=RealSet.closed_open(0, oo))] = 0.1,
) -> None:
    super().__init__(label, visual=visual)
    self.S = S
    self.E = E
    self.I = I
    self.R = R
    self.J = J
    self.c = c
    self.beta = beta
    self.tau = tau
    self.T = T
    self.nu = nu

PathogenElement(label, visual=None, *, eta=0.0)

Bases: Element

Pathogen whose strength in every season is drawn anew, uniformly from zero to one.

Parameters:

Name	Type	Description	Default
label	str	The label of the pathogen.	required
visual	VisualSpec | None	How the element is drawn, or None for default_visual_.	None
eta	Annotated[float, StateSpec(domain=closed_open(0, 1))]	The strength of the strain of the current season.	0.0

Examples:

>>> from skcomplex.systems.biology.epidemiology import PathogenElement
>>> PathogenElement('pathogen').eta
0.0

Source code in src/skcomplex/systems/biology/epidemiology/_records.py

def __init__(
    self,
    label: str,
    visual: VisualSpec | None = None,
    *,
    eta: Annotated[float, StateSpec(domain=RealSet.closed_open(0, 1))] = 0.0,
) -> None:
    super().__init__(label, visual=visual)
    self.eta = eta

PersonElement(label, visual=VisualSpec('agents', {'color_by': 'health', 'color_range': (0.0, 4.0)}), *, health=SUSCEPTIBLE, d=0.8, beta=0.6, gamma=0.01, mu=0.0, alpha=0.0, lockdown_start=0.0, lockdown_length=0.0, lockdown_d=0.1)

Bases: Element

Person of the agent-based SIR system, with incubation, mortality, waning immunity, and lockdown.

A MovementInteraction walks a living person to a neighbouring cell with probability d, or lockdown_d during the lockdown, which lasts lockdown_length steps from lockdown_start. Transitions of the agent form, which agent_based_sir adds, move the health of the person. An infectious person infects, with probability beta, every susceptible person in its cell, who becomes exposed. An infectious person recovers with probability gamma or dies with probability mu, and a recovered person loses immunity with probability alpha. Health is 0 susceptible, 1 exposed, 2 infected, 3 recovered, or 4 dead.

The incubation follows the convention of the residence time after of TransitionInteraction: the step of the exposure counts as the first step in the exposed state, so a person exposed in step n becomes infectious at the end of step n + incubation - 1. The clock that the person kept before counted the step of the exposure only when the person acted after its infector in that step, so the onset now comes half a step earlier on average.

Parameters:

Name	Type	Description	Default
label	str	The label of the agent or group.	required
visual	VisualSpec | None	How the element is drawn, or None for default_visual_.	VisualSpec('agents', {'color_by': 'health', 'color_range': (0.0, 4.0)})
health	Annotated[float, StateSpec(domain=NN)]	The health state.	SUSCEPTIBLE
d	Annotated[float, ParameterSpec(domain=closed_open(0, oo))]	The probability of moving at every step.	0.8
beta	Annotated[float, ParameterSpec(domain=closed_open(0, oo))]	The probability that an infectious person infects the susceptible persons of its cell.	0.6
gamma	Annotated[float, ParameterSpec(domain=closed_open(0, oo))]	The probability of recovering at every step.	0.01
mu	Annotated[float, ParameterSpec(domain=closed_open(0, oo))]	The probability of dying at every step of the infection.	0.0
alpha	Annotated[float, ParameterSpec(domain=closed_open(0, oo))]	The probability of losing immunity at every step.	0.0
lockdown_start	Annotated[float, ParameterSpec(domain=closed_open(0, oo))]	The step at which the lockdown starts.	0.0
lockdown_length	Annotated[float, ParameterSpec(domain=closed_open(0, oo))]	The number of steps of the lockdown, 0 for none.	0.0
lockdown_d	Annotated[float, ParameterSpec(domain=closed_open(0, oo))]	The probability of moving during the lockdown.	0.1

Examples:

>>> from skcomplex.core.components.elements import GroupElement
>>> from skcomplex.systems.biology.epidemiology import PersonElement
>>> GroupElement('people', element=PersonElement, size=100).size
100

Source code in src/skcomplex/systems/biology/epidemiology/_agents.py

def __init__(
    self,
    label: str,
    visual: VisualSpec | None = VisualSpec('agents', {'color_by': 'health', 'color_range': (0.0, 4.0)}),
    *,
    health: Annotated[float, StateSpec(domain=NN)] = SUSCEPTIBLE,
    d: Annotated[float, ParameterSpec(domain=RealSet.closed_open(0, oo))] = 0.8,
    beta: Annotated[float, ParameterSpec(domain=RealSet.closed_open(0, oo))] = 0.6,
    gamma: Annotated[float, ParameterSpec(domain=RealSet.closed_open(0, oo))] = 0.01,
    mu: Annotated[float, ParameterSpec(domain=RealSet.closed_open(0, oo))] = 0.0,
    alpha: Annotated[float, ParameterSpec(domain=RealSet.closed_open(0, oo))] = 0.0,
    lockdown_start: Annotated[float, ParameterSpec(domain=RealSet.closed_open(0, oo))] = 0.0,
    lockdown_length: Annotated[float, ParameterSpec(domain=RealSet.closed_open(0, oo))] = 0.0,
    lockdown_d: Annotated[float, ParameterSpec(domain=RealSet.closed_open(0, oo))] = 0.1,
) -> None:
    super().__init__(label, visual=visual)
    self.health = health
    self.d = d
    self.beta = beta
    self.gamma = gamma
    self.mu = mu
    self.alpha = alpha
    self.lockdown_start = lockdown_start
    self.lockdown_length = lockdown_length
    self.lockdown_d = lockdown_d

RecoveryInteraction(label, *, infected, recovered, gamma=0.1, stochastic=False)

Bases: ReactionInteraction

Recovery, the reaction I -> R at the rate gamma I from infected to recovered.

It is a ReactionInteraction between the states n of the two compartments its roles bind, a rate equation, or in stochastic mode a random event that moves one individual.

Parameters:

Name	Type	Description	Default
label	str	The label of the interaction.	required
infected	Annotated[str, RoleSpec()]	The path of the infected compartment.	required
recovered	Annotated[str, RoleSpec()]	The path of the recovered compartment.	required
gamma	Annotated[float, ParameterSpec(domain=closed_open(0, oo))]	The recovery rate.	0.1
stochastic	bool	Whether recoveries are random events.	False

Examples:

>>> from skcomplex.core.components.elements import CompartmentElement
>>> from skcomplex.core.components.system import System
>>> from skcomplex.core.engine.runtime import Simulation
>>> from skcomplex.systems.biology.epidemiology import RecoveryInteraction
>>> system = System(
...     elements=[CompartmentElement('I', n=10.0), CompartmentElement('R')],
...     interactions=[RecoveryInteraction('recovery', infected='I', recovered='R', gamma=0.5)],
... )
>>> trajectory = Simulation(system, dt=0.01).prepare().run(until=2.0).trajectory_
>>> round(float(trajectory.values('I.n')[-1]), 3)
3.679

Source code in src/skcomplex/systems/biology/epidemiology/_factories.py

def __init__(
    self,
    label: str,
    *,
    infected: Annotated[str, RoleSpec()],
    recovered: Annotated[str, RoleSpec()],
    gamma: Annotated[float, ParameterSpec(domain=RealSet.closed_open(0, oo))] = 0.1,
    stochastic: bool = False,
) -> None:
    super().__init__(label)
    self.infected = infected
    self.recovered = recovered
    self.gamma = gamma
    self.stochastic = stochastic

ResidentElement(label, visual=None, *, R=0.5, infected=0.0)

Bases: Element

Resident of a population with a fixed resistance to a disease, immune for good once infected.

Parameters:

Name	Type	Description	Default
label	str	The label of the resident or of the group of residents.	required
visual	VisualSpec | None	How the element is drawn, or None for default_visual_.	None
R	Annotated[float, StateSpec(domain=closed(0, 1))]	The resistance, between 0 and 1, where 1 is immunity.	0.5
infected	Annotated[float, StateSpec(domain=FiniteEnumeratedSet([0, 1]))]	Whether the resident was infected in the last outbreak.	0.0

Examples:

>>> from skcomplex.systems.biology.epidemiology import ResidentElement
>>> ResidentElement('resident', R=0.3).R
0.3

Source code in src/skcomplex/systems/biology/epidemiology/_records.py

def __init__(
    self,
    label: str,
    visual: VisualSpec | None = None,
    *,
    R: Annotated[float, StateSpec(domain=RealSet.closed(0, 1))] = 0.5,
    infected: Annotated[float, StateSpec(domain=FiniteEnumeratedSet([0, 1]))] = 0.0,
) -> None:
    super().__init__(label, visual=visual)
    self.R = R
    self.infected = infected

ResistanceInfectionInteraction(label, *, residents, pathogen)

Bases: Interaction

Infection of every resident whose resistance is below the strength of the season's strain, who becomes immune.

A resident with R < eta is infected in the season and its resistance becomes 1, so it is never infected again.

Parameters:

Name	Type	Description	Default
label	str	The label of the interaction.	required
residents	Annotated[str, RoleSpec()]	The path of the group of residents.	required
pathogen	Annotated[str, RoleSpec()]	The path of the pathogen.	required

Examples:

>>> from skcomplex.systems.biology.epidemiology import ResistanceInfectionInteraction
>>> ResistanceInfectionInteraction('infection', residents='people', pathogen='flu').pathogen
'flu'

Source code in src/skcomplex/systems/biology/epidemiology/_records.py

def __init__(
    self,
    label: str,
    *,
    residents: Annotated[str, RoleSpec()],
    pathogen: Annotated[str, RoleSpec()],
) -> None:
    super().__init__(label)
    self.residents = residents
    self.pathogen = pathogen

SpreadingInteraction(label, *, target, state='state', beta=1.0, gamma=0.0, omega=0.0, inactive=0.0, active=1.0, removed=None, form='number', layer=None)

Bases: TransitionInteraction

Spreading of the active value of a discrete state along the edges of a network or the bonds of a lattice.

An inactive member becomes active at the rate beta times the number of its active neighbours, or times their fraction of its neighbours. An active member recovers at the rate gamma, to the inactive value, or to the removed value when there is one, and a removed member returns to the inactive value at the rate omega. Every move is a transition at a rate, so the state changes by random events in continuous time.

Parameters:

Name	Type	Description	Default
label	str	The label of the interaction, unique among the parts of its system.	required
target	Annotated[str, RoleSpec()]	The path of the group on the nodes or sites.	required
state	str	The name of the discrete state of the members.	'state'
beta	Annotated[float, ParameterSpec(domain=closed_open(0, oo))]	The rate of spreading from one active neighbour, or from all neighbours when they are all active.	1.0
gamma	Annotated[float, ParameterSpec(domain=closed_open(0, oo))]	The rate of recovery of an active member.	0.0
omega	Annotated[float, ParameterSpec(domain=closed_open(0, oo))]	The rate at which a removed member returns to the inactive value.	0.0
inactive	float	The value of the state of an inactive member.	0.0
active	float	The value of the state of an active member.	1.0
removed	float | None	The value of the state of a recovered member, or None for members that recover to the inactive value.	None
form	str	'number' for a rate proportional to the number of active neighbours, or 'fraction' for one proportional to their fraction of the neighbours.	'number'
layer	str | None	The layer of the network whose edges the spreading runs on, or None for the graph of the network or the neighbourhood of the lattice.	None

Raises:

Type	Description
ValidationError	The form is unknown, or the values of the state are not distinct.

Examples:

>>> from skcomplex.systems.biology.epidemiology import SpreadingInteraction
>>> spreading = SpreadingInteraction('spreading', target='sites', beta=2.0, gamma=1.0, form='fraction')
>>> spreading.form
'fraction'

Source code in src/skcomplex/systems/biology/common/networks/_spreading.py

def __init__(
    self,
    label: str,
    *,
    target: Annotated[str, RoleSpec()],
    state: str = 'state',
    beta: Annotated[float, ParameterSpec(domain=RealSet.closed_open(0, oo))] = 1.0,
    gamma: Annotated[float, ParameterSpec(domain=RealSet.closed_open(0, oo))] = 0.0,
    omega: Annotated[float, ParameterSpec(domain=RealSet.closed_open(0, oo))] = 0.0,
    inactive: float = 0.0,
    active: float = 1.0,
    removed: float | None = None,
    form: str = 'number',
    layer: str | None = None,
) -> None:
    super().__init__(label, target=target, state=state)
    self.beta = beta
    self.gamma = gamma
    self.omega = omega
    self.inactive = inactive
    self.active = active
    self.removed = removed
    self.form = form
    self.layer = layer

SubpopulationElement(label, visual=None, *, S=1000.0, I=0.0, R=0.0, beta=0.3, gamma=0.1, D=0.01)

Bases: Element

Well-mixed SIR population at a node of a network, whose individuals hop along each edge at the rate D.

Inside the node, beta S I / N individuals are infected and gamma I recover per unit time, with N = S + I + R. Hopping along every edge adds D times the graph Laplacian to every compartment.

Parameters:

Name	Type	Description	Default
label	str	The label of the subpopulation or the field.	required
visual	VisualSpec | None	How the element is drawn, or None for default_visual_.	None
S	Annotated[float, StateSpec(domain=closed_open(0, oo))]	The initial number of susceptible individuals.	1000.0
I	Annotated[float, StateSpec(domain=closed_open(0, oo))]	The initial number of infected individuals.	0.0
R	Annotated[float, StateSpec(domain=closed_open(0, oo))]	The initial number of recovered individuals.	0.0
beta	Annotated[float, ParameterSpec(domain=closed_open(0, oo))]	The infection rate.	0.3
gamma	Annotated[float, ParameterSpec(domain=closed_open(0, oo))]	The recovery rate.	0.1
D	Annotated[float, ParameterSpec(domain=closed_open(0, oo))]	The rate at which an individual hops along one edge.	0.01

Examples:

>>> from skcomplex.systems.biology.epidemiology import SubpopulationElement
>>> SubpopulationElement('city', S=1000.0).beta
0.3

Source code in src/skcomplex/systems/biology/epidemiology/_networks.py

def __init__(
    self,
    label: str,
    visual: VisualSpec | None = None,
    *,
    S: Annotated[float, StateSpec(domain=RealSet.closed_open(0, oo))] = 1000.0,
    I: Annotated[float, StateSpec(domain=RealSet.closed_open(0, oo))] = 0.0,
    R: Annotated[float, StateSpec(domain=RealSet.closed_open(0, oo))] = 0.0,
    beta: Annotated[float, ParameterSpec(domain=RealSet.closed_open(0, oo))] = 0.3,
    gamma: Annotated[float, ParameterSpec(domain=RealSet.closed_open(0, oo))] = 0.1,
    D: Annotated[float, ParameterSpec(domain=RealSet.closed_open(0, oo))] = 0.01,
) -> None:
    super().__init__(label, visual=visual)
    self.S = S
    self.I = I
    self.R = R
    self.beta = beta
    self.gamma = gamma
    self.D = D

TravellerElement(label, visual=VisualSpec('agents', {'color_by': 'health', 'color_range': (0.0, 3.0)}), *, health=SUSCEPTIBLE, tau=0.05, m=0.05)

Bases: Element

Resident of a city at a network node, who infects fellow residents and travels to neighbouring cities.

In its own activation, an infected traveller infects every susceptible traveller of its node with probability tau, by a TransitionInteraction whose members act on their cell, and another transition makes it recover after a residence of a number of steps. A MovementInteraction moves it to a random neighbouring node with probability m. Health is 0 susceptible, 2 infected, or 3 recovered.

Parameters:

Name	Type	Description	Default
label	str	The label of the agent or group.	required
visual	VisualSpec | None	How the element is drawn, or None for default_visual_.	VisualSpec('agents', {'color_by': 'health', 'color_range': (0.0, 3.0)})
health	Annotated[float, StateSpec(domain=NN)]	The health state.	SUSCEPTIBLE
tau	Annotated[float, ParameterSpec(domain=closed_open(0, oo))]	The probability of infecting a susceptible traveller of the same node at every step.	0.05
m	Annotated[float, ParameterSpec(domain=closed_open(0, oo))]	The probability of moving to a neighbouring node at every step.	0.05

Examples:

>>> from skcomplex.core.components.elements import GroupElement
>>> from skcomplex.systems.biology.epidemiology import TravellerElement
>>> GroupElement('travellers', element=TravellerElement, size=10).size
10

Source code in src/skcomplex/systems/biology/epidemiology/_agents.py

def __init__(
    self,
    label: str,
    visual: VisualSpec | None = VisualSpec('agents', {'color_by': 'health', 'color_range': (0.0, 3.0)}),
    *,
    health: Annotated[float, StateSpec(domain=NN)] = SUSCEPTIBLE,
    tau: Annotated[float, ParameterSpec(domain=RealSet.closed_open(0, oo))] = 0.05,
    m: Annotated[float, ParameterSpec(domain=RealSet.closed_open(0, oo))] = 0.05,
) -> None:
    super().__init__(label, visual=visual)
    self.health = health
    self.tau = tau
    self.m = m

VaccinatedPopulationElement(label, visual=None, *, S=0.99, I=0.01, R=0.0, u=0.0, beta=0.5, gamma=0.1, a=1.0, b=1.0)

Bases: Element

Susceptible, infected, and recovered fractions of a population vaccinated at a bounded rate.

The fractions follow S' = -beta S I - u S, I' = beta S I - gamma I, and R' = gamma I + u S, where the vaccination rate u is a control input between 0 and its largest value. The planner minimises the burden of the infections and the cost of the vaccines, ∫ (a I + b u^2 / 2) dt, after Lenhart and Workman.

Parameters:

Name	Type	Description	Default
label	str	The label of the population.	required
visual	VisualSpec | None	How the element is drawn, or None for default_visual_.	None
S	Annotated[float, StateSpec(domain=closed_open(0, oo))]	The initial susceptible fraction.	0.99
I	Annotated[float, StateSpec(domain=closed_open(0, oo))]	The initial infected fraction.	0.01
R	Annotated[float, StateSpec(domain=closed_open(0, oo))]	The initial recovered fraction.	0.0
u	Annotated[float, ParameterSpec(adaptive=True, bounds=(0.0, MAX_RATE))]	The vaccination rate, the control input.	0.0
beta	Annotated[float, ParameterSpec(domain=closed_open(0, oo))]	The transmission rate.	0.5
gamma	Annotated[float, ParameterSpec(domain=closed_open(0, oo))]	The recovery rate.	0.1
a	Annotated[float, ParameterSpec(domain=closed_open(0, oo))]	The weight of the infected fraction in the cost.	1.0
b	Annotated[float, ParameterSpec(domain=open(0, oo))]	The weight of the squared vaccination rate in the cost.	1.0

Examples:

>>> from skcomplex.systems.biology.epidemiology import VaccinatedPopulationElement
>>> VaccinatedPopulationElement('population', I=0.01).I
0.01

Source code in src/skcomplex/systems/biology/epidemiology/_vaccination.py

def __init__(
    self,
    label: str,
    visual: VisualSpec | None = None,
    *,
    S: Annotated[float, StateSpec(domain=RealSet.closed_open(0, oo))] = 0.99,
    I: Annotated[float, StateSpec(domain=RealSet.closed_open(0, oo))] = 0.01,
    R: Annotated[float, StateSpec(domain=RealSet.closed_open(0, oo))] = 0.0,
    u: Annotated[float, ParameterSpec(adaptive=True, bounds=(0.0, MAX_RATE))] = 0.0,
    beta: Annotated[float, ParameterSpec(domain=RealSet.closed_open(0, oo))] = 0.5,
    gamma: Annotated[float, ParameterSpec(domain=RealSet.closed_open(0, oo))] = 0.1,
    a: Annotated[float, ParameterSpec(domain=RealSet.closed_open(0, oo))] = 1.0,
    b: Annotated[float, ParameterSpec(domain=RealSet.open(0, oo))] = 1.0,
) -> None:
    super().__init__(label, visual=visual)
    self.S = S
    self.I = I
    self.R = R
    self.u = u
    self.beta = beta
    self.gamma = gamma
    self.a = a
    self.b = b

WalkerElement(label, visual=VisualSpec('agents', {'color_by': 'health', 'color_range': (0.0, 2.0)}), *, health=HEALTHY, infection=1.0, travel=1.0)

Bases: Element

Random walker that carries a deadly disease, of Charbonneau's epidemic system.

In its own activation, a living walker moves to a neighbouring cell or node with probability travel by a MovementInteraction, and a TransitionInteraction whose members act on their cell makes a sick walker infect every healthy walker of its cell with probability infection. Another transition kills a sick walker after a residence of a number of steps. Health is 0 healthy, 1 sick, or 2 dead.

Parameters:

Name	Type	Description	Default
label	str	The label of the agent or group.	required
visual	VisualSpec | None	How the element is drawn, or None for default_visual_.	VisualSpec('agents', {'color_by': 'health', 'color_range': (0.0, 2.0)})
health	Annotated[float, StateSpec(domain=NN)]	The health state.	HEALTHY
infection	Annotated[float, ParameterSpec(domain=closed_open(0, oo))]	The probability that a sick walker infects a healthy walker of its cell at every step.	1.0
travel	Annotated[float, ParameterSpec(domain=closed_open(0, oo))]	The probability of moving at every step.	1.0

Examples:

>>> from skcomplex.core.components.elements import GroupElement
>>> from skcomplex.systems.biology.epidemiology import WalkerElement
>>> GroupElement('walkers', element=WalkerElement, size=10).size
10

Source code in src/skcomplex/systems/biology/epidemiology/_agents.py

def __init__(
    self,
    label: str,
    visual: VisualSpec | None = VisualSpec('agents', {'color_by': 'health', 'color_range': (0.0, 2.0)}),
    *,
    health: Annotated[float, StateSpec(domain=NN)] = HEALTHY,
    infection: Annotated[float, ParameterSpec(domain=RealSet.closed_open(0, oo))] = 1.0,
    travel: Annotated[float, ParameterSpec(domain=RealSet.closed_open(0, oo))] = 1.0,
) -> None:
    super().__init__(label, visual=visual)
    self.health = health
    self.infection = infection
    self.travel = travel

activity_driven_epidemics(n=500, *, m=2, exponent=2.2, least=0.01, lam=0.9, mu=0.1, infected=0.05, seed=1)

Return SIS spreading on an activity-driven temporal network, whose edges last one step.

At every step every node fires with probability its activity, drawn from a^-exponent above least, and joins m random nodes. The infection travels along the edges of a step, which disappear at the next step, so it spreads only when lam / mu exceeds a threshold set by the moments of the activity, 2 <a> / (m (<a> + sqrt(<a^2>))) in the limit of small probabilities, and not by the degrees of the network aggregated over time.

Parameters:

Name	Type	Description	Default
n	int	The number of nodes.	500
m	int	The number of edges an active node creates.	2
exponent	float	The exponent of the power-law distribution of activities.	2.2
least	float	The smallest activity.	0.01
lam	float	The probability of infection along one edge.	0.9
mu	float	The probability of recovery at a step.	0.1
infected	float	The initial share of infected nodes.	0.05
seed	int | None	The seed of the activities and the initial infections.	1

Returns:

Name	Type	Description
system	System	The system, with the states hosts.i and activities hosts.a of the group hosts on the nodes, and the variable group contacts on the edges of every step.

Examples:

>>> from skcomplex.systems.biology.epidemiology import activity_driven_epidemics
>>> float(activity_driven_epidemics(n=10, seed=0).elements[0].values['a'].min()) >= 0.01
True

Source code in src/skcomplex/systems/biology/epidemiology/_adaptive.py

@register_system(
    field='biology/epidemiology',
    tags=('temporal network', 'activity-driven', 'stochastic map', 'epidemic threshold'),
    references=('Perra et al. 2012', 'Pastor-Satorras et al. 2015'),
    defaults={'until': 150, 'dt': 1.0},
)
def activity_driven_epidemics(
    n: int = 500,
    *,
    m: int = 2,
    exponent: float = 2.2,
    least: float = 0.01,
    lam: float = 0.9,
    mu: float = 0.1,
    infected: float = 0.05,
    seed: int | None = 1,
) -> System:
    """Return SIS spreading on an activity-driven temporal network, whose edges last one step.

    At every step every node fires with probability its activity, drawn from `a^-exponent` above `least`, and joins
    `m` random nodes. The infection travels along the edges of a step, which disappear at the next step, so it spreads
    only when `lam / mu` exceeds a threshold set by the moments of the activity, `2 <a> / (m (<a> + sqrt(<a^2>)))` in
    the limit of small probabilities, and not by the degrees of the network aggregated over time.

    Args:
        n:
            The number of nodes.

        m:
            The number of edges an active node creates.

        exponent:
            The exponent of the power-law distribution of activities.

        least:
            The smallest activity.

        lam:
            The probability of infection along one edge.

        mu:
            The probability of recovery at a step.

        infected:
            The initial share of infected nodes.

        seed:
            The seed of the activities and the initial infections.

    Returns:
        system:
            The system, with the states `hosts.i` and activities `hosts.a` of the group `hosts` on the nodes, and the
            variable group `contacts` on the edges of every step.

    Examples:
        >>> from skcomplex.systems.biology.epidemiology import activity_driven_epidemics
        >>> float(activity_driven_epidemics(n=10, seed=0).elements[0].values['a'].min()) >= 0.01
        True
    """
    rng = np.random.default_rng(seed)
    activities = np.minimum(least * (1 - rng.random(n)) ** (-1 / (exponent - 1)), 1.0)
    initial = (rng.random(n) < infected).astype(float)
    space = NetworkSpace(nx.empty_graph(n))
    hosts = GroupElement(
        'hosts',
        element=ActiveHostElement,
        space=space,
        on='sites',
        variable=True,
        capacity=n,
        values={'i': initial, 'a': activities, 'lam': lam, 'mu': mu},
        visual=VisualSpec('sites', {'states': ['i'], 'low': 0, 'high': 1}),
    )
    contacts = GroupElement(
        'contacts',
        element=Element,
        space=space,
        on='edges',
        variable=True,
        capacity=2 * m * n,
        visual=HIDDEN,
    )
    return System(
        elements=[hosts, contacts],
        interactions=[
            TopologyInteraction('firing', target='hosts', action='add-edge', probability=_get_activity, edges=m),
            TopologyInteraction('expiry', target='hosts', action='remove-edge'),
            TransitionInteraction(
                'infection',
                target='hosts',
                state='i',
                source=0,
                destination=1,
                probability=_get_infection,
            ),
            TransitionInteraction(
                'recovery',
                target='hosts',
                state='i',
                source=1,
                destination=0,
                probability=_get_recovery,
            ),
        ],
        time=DiscreteTime(),
    )

adaptive_sis(n=200, k=6.0, *, beta=0.4, gamma=1.0, w=1.0, infected=0.2, seed=1)

Return the SIS epidemic on a network whose susceptible nodes cut their edges to infected neighbours.

Every edge between a susceptible and an infected node is rewired at rate w: the susceptible end keeps it and joins it to another susceptible node. Rewiring isolates the infected nodes and raises the epidemic threshold above its value on the fixed network, while the susceptible nodes cluster.

Parameters:

Name	Type	Description	Default
n	int	The number of nodes.	200
k	float	The mean degree of the initial Erdős-Rényi network.	6.0
beta	float	The infection rate per edge to an infected neighbour.	0.4
gamma	float	The recovery rate.	1.0
w	float	The rewiring rate of an edge between a susceptible and an infected node.	1.0
infected	float	The initial share of infected nodes.	0.2
seed	int | None	The seed of the network and of the initial infections.	1

Returns:

Name	Type	Description
system	System	The system, with the states hosts.i of the variable group hosts on the nodes and the variable group contacts on the edges.

Examples:

>>> from skcomplex.systems.biology.epidemiology import adaptive_sis
>>> [component.label for component in adaptive_sis(n=20).elements]
['hosts', 'contacts']

Source code in src/skcomplex/systems/biology/epidemiology/_adaptive.py

@register_system(
    field='biology/epidemiology',
    tags=('network', 'adaptive network', 'rewiring', 'jumps', 'epidemic threshold'),
    references=("Gross, Dommar D'Lima, and Blasius 2006", 'Thurner et al. 2018', 'Pastor-Satorras et al. 2015'),
    defaults={'until': 30.0, 'dt': 0.25},
)
def adaptive_sis(
    n: int = 200,
    k: float = 6.0,
    *,
    beta: float = 0.4,
    gamma: float = 1.0,
    w: float = 1.0,
    infected: float = 0.2,
    seed: int | None = 1,
) -> System:
    """Return the SIS epidemic on a network whose susceptible nodes cut their edges to infected neighbours.

    Every edge between a susceptible and an infected node is rewired at rate `w`: the susceptible end keeps it and
    joins it to another susceptible node. Rewiring isolates the infected nodes and raises the epidemic threshold above
    its value on the fixed network, while the susceptible nodes cluster.

    Args:
        n:
            The number of nodes.

        k:
            The mean degree of the initial Erdős-Rényi network.

        beta:
            The infection rate per edge to an infected neighbour.

        gamma:
            The recovery rate.

        w:
            The rewiring rate of an edge between a susceptible and an infected node.

        infected:
            The initial share of infected nodes.

        seed:
            The seed of the network and of the initial infections.

    Returns:
        system:
            The system, with the states `hosts.i` of the variable group `hosts` on the nodes and the variable group
            `contacts` on the edges.

    Examples:
        >>> from skcomplex.systems.biology.epidemiology import adaptive_sis
        >>> [component.label for component in adaptive_sis(n=20).elements]
        ['hosts', 'contacts']
    """
    space = _random_network(n, k, seed)
    initial = (np.random.default_rng(seed).random(n) < infected).astype(float)
    hosts = GroupElement(
        'hosts',
        element=HostElement,
        space=space,
        on='sites',
        variable=True,
        capacity=n,
        values={'i': initial, 'w': w},
        visual=VisualSpec('sites', {'states': ['i'], 'low': 0, 'high': 1}),
    )
    edges = space.graph.number_of_edges()
    contacts = GroupElement(
        'contacts',
        element=Element,
        space=space,
        on='edges',
        variable=True,
        capacity=edges,
        visual=HIDDEN,
    )
    spreading = SpreadingInteraction('spreading', target='hosts', state='i', beta=beta, gamma=gamma)
    return System(
        elements=[hosts, contacts],
        interactions=[spreading, _get_rewiring('source'), _get_rewiring('target')],
    )

agent_based_sir(width=50, height=50, *, people=500, infected=5, d=0.8, beta=0.6, gamma=0.02, mu=0.0, alpha=0.0, incubation=0.0, lockdown_start=0.0, lockdown_length=0.0, lockdown_d=0.1, random_state=0)

Return people who walk on a grid and infect those in their cell, with optional extensions of the SIR system.

The basic system infects with probability beta and recovers with probability gamma. Its final size depends on beta and gamma separately, not only on their ratio, because walking limits the contacts. Mortality mu, waning immunity alpha, an incubation of incubation steps, and a lockdown that lowers the probability of moving to lockdown_d extend it.

Parameters:

Name	Type	Description	Default
width	int	The number of columns of the grid.	50
height	int	The number of rows of the grid.	50
people	int	The number of people.	500
infected	int	The initial number of infected people.	5
d	float	The probability of moving at every step.	0.8
beta	float	The probability that an infectious person infects the susceptible persons of its cell.	0.6
gamma	float	The probability of recovering at every step.	0.02
mu	float	The probability of dying at every step of the infection.	0.0
alpha	float	The probability of losing immunity at every step.	0.0
incubation	float	The number of steps an infected person stays exposed before it becomes infectious, counting the step of the exposure, as the residence time of TransitionInteraction counts them.	0.0
lockdown_start	float	The step at which the lockdown starts.	0.0
lockdown_length	float	The number of steps of the lockdown, 0 for none.	0.0
lockdown_d	float	The probability of moving during the lockdown.	0.1
random_state	int | Generator | None	The seed or generator of the initial cells and infections.	0

Returns:

Name	Type	Description
system	System	The system, with the group people.

Examples:

>>> from skcomplex.systems.biology.epidemiology import agent_based_sir
>>> int(agent_based_sir(people=100, infected=3, random_state=0).elements[0].values['health'].sum())
6

Source code in src/skcomplex/systems/biology/epidemiology/_agents.py

@register_system(
    field='biology/epidemiology',
    tags=('compartments', 'agents', 'grid', 'lockdown', 'mortality'),
    references=('Argun et al. 2021', 'Kermack and McKendrick 1927'),
    defaults={'until': 300},
)
def agent_based_sir(
    width: int = 50,
    height: int = 50,
    *,
    people: int = 500,
    infected: int = 5,
    d: float = 0.8,
    beta: float = 0.6,
    gamma: float = 0.02,
    mu: float = 0.0,
    alpha: float = 0.0,
    incubation: float = 0.0,
    lockdown_start: float = 0.0,
    lockdown_length: float = 0.0,
    lockdown_d: float = 0.1,
    random_state: int | np.random.Generator | None = 0,
) -> System:
    """Return people who walk on a grid and infect those in their cell, with optional extensions of the SIR system.

    The basic system infects with probability `beta` and recovers with probability `gamma`. Its final size depends on
    `beta` and `gamma` separately, not only on their ratio, because walking limits the contacts. Mortality `mu`, waning
    immunity `alpha`, an incubation of `incubation` steps, and a lockdown that lowers the probability of moving to
    `lockdown_d` extend it.

    Args:
        width:
            The number of columns of the grid.

        height:
            The number of rows of the grid.

        people:
            The number of people.

        infected:
            The initial number of infected people.

        d:
            The probability of moving at every step.

        beta:
            The probability that an infectious person infects the susceptible persons of its cell.

        gamma:
            The probability of recovering at every step.

        mu:
            The probability of dying at every step of the infection.

        alpha:
            The probability of losing immunity at every step.

        incubation:
            The number of steps an infected person stays exposed before it becomes infectious, counting the step of
            the exposure, as the residence time of `TransitionInteraction` counts them.

        lockdown_start:
            The step at which the lockdown starts.

        lockdown_length:
            The number of steps of the lockdown, 0 for none.

        lockdown_d:
            The probability of moving during the lockdown.

        random_state:
            The seed or generator of the initial cells and infections.

    Returns:
        system:
            The system, with the group `people`.

    Examples:
        >>> from skcomplex.systems.biology.epidemiology import agent_based_sir
        >>> int(agent_based_sir(people=100, infected=3, random_state=0).elements[0].values['health'].sum())
        6
    """
    rng = np.random.default_rng(random_state)
    cells = rng.integers(0, [width, height], size=(people, 2))
    group = GroupElement(
        'people',
        element=PersonElement,
        size=people,
        values={
            'health': _initial_health(people, infected, rng, INFECTED),
            'd': d,
            'beta': beta,
            'gamma': gamma,
            'mu': mu,
            'alpha': alpha,
            'lockdown_start': lockdown_start,
            'lockdown_length': lockdown_length,
            'lockdown_d': lockdown_d,
        },
        location=[tuple(int(item) for item in cell) for cell in cells],
    )
    walk = MovementInteraction('walk', target='people', condition=_is_living_person, accept=_get_person_step)
    space = LatticeSpace((width, height), boundary='fixed', neighborhood='von_neumann')
    return System(elements=[group], interactions=[walk, *_get_person_transitions(incubation)], space=space)

delayed_sirs(first=0.01, second=0.0, *, a=2.0, recovery=6, e=0.005)

Return two cities with recurrent outbreaks of a disease of fixed infection and immunity periods.

With one infected step and recovery immune steps, each city follows x' = a x (1 - x - x(t-1) - ... - x(t-recovery)), which reduces to the logistic map without immunity. Outbreaks recur like relaxation oscillations, and the visits of a fraction e of the infected couple the two cities, whose outbreaks lock in or out of phase.

Parameters:

Name	Type	Description	Default
first	float	The initial infected fraction of the first city.	0.01
second	float	The initial infected fraction of the second city.	0.0
a	float	The rate of transmission per contact.	2.0
recovery	int	The number of steps a recovered individual stays immune, at most 12.	6
e	float	The fraction of the infected who visit the other city.	0.005

Returns:

Name	Type	Description
system	System	The system, with the centers first and second.

Examples:

>>> from skcomplex.systems.biology.epidemiology import delayed_sirs
>>> delayed_sirs().elements[0].x
0.01

Source code in src/skcomplex/systems/biology/epidemiology/_recurrent.py

@register_system(
    field='biology/epidemiology',
    tags=('compartments', 'discrete map', 'delay', 'synchronization', 'relaxation oscillation'),
    references=('Gros 2015', 'He and Stone 2003'),
    defaults={'until': 200, 'dt': 1},
)
def delayed_sirs(
    first: float = 0.01,
    second: float = 0.0,
    *,
    a: float = 2.0,
    recovery: int = 6,
    e: float = 0.005,
) -> System:
    """Return two cities with recurrent outbreaks of a disease of fixed infection and immunity periods.

    With one infected step and `recovery` immune steps, each city follows `x' = a x (1 - x - x(t-1) - ... -
    x(t-recovery))`, which reduces to the logistic map without immunity. Outbreaks recur like relaxation oscillations,
    and the visits of a fraction `e` of the infected couple the two cities, whose outbreaks lock in or out of phase.

    Args:
        first:
            The initial infected fraction of the first city.

        second:
            The initial infected fraction of the second city.

        a:
            The rate of transmission per contact.

        recovery:
            The number of steps a recovered individual stays immune, at most 12.

        e:
            The fraction of the infected who visit the other city.

    Returns:
        system:
            The system, with the centers `first` and `second`.

    Examples:
        >>> from skcomplex.systems.biology.epidemiology import delayed_sirs
        >>> delayed_sirs().elements[0].x
        0.01
    """
    return System(
        elements=[
            EpidemicCenterElement('first', x=first, a=a, recovery=recovery, e=e),
            EpidemicCenterElement('second', x=second, a=a, recovery=recovery, e=e),
        ],
        interactions=[
            MapInteraction('visits_first', target='first', state='x', function=_get_next_infected, sources=('second',)),
            MapInteraction(
                'visits_second',
                target='second',
                state='x',
                function=_get_next_infected,
                sources=('first',),
            ),
        ],
    )

episims_city(people=1000, *, household=3.0, workplace=20.0, p=0.05, gamma=0.1, infected=5, closed=False, seed=0)

Return a synthetic city in the spirit of EpiSims, where people infect those they meet at home and at work.

EpiSims models a city as a bipartite graph of people and locations whose edges are visits labelled with their times, and spreads disease between people at the same location. Here every person belongs to a household of mean size household and a workplace or school of mean size workplace, the people sleep at home and work by day, and the disease spreads by SIR contacts at every location. The locations form a network whose edges join a home to the workplaces of its members, the projection of the bipartite graph onto the locations. Households spread the disease within families and workplaces between them, so closing the workplaces, one of the interventions EpiSims was built to compare, confines outbreaks to households.

Parameters:

Name	Type	Description	Default
people	int	The number of people.	1000
household	float	The mean size of a household.	3.0
workplace	float	The mean size of a workplace or school.	20.0
p	float	The probability of transmission per infected contact and half day.	0.05
gamma	float	The probability of recovery per half day.	0.1
infected	int	The initial number of infected people.	5
closed	bool	Whether the workplaces are closed.	False
seed	int | None	The seed of the households, the workplaces, and the initial infections.	0

Returns:

Name	Type	Description
system	System	The system, with the homes, workplaces, and health states of the agent group people on the locations.

Examples:

>>> from skcomplex.systems.biology.epidemiology import episims_city
>>> episims_city(100).elements[0].size
100

Source code in src/skcomplex/systems/biology/epidemiology/_cities.py

@register_system(
    field='biology/epidemiology',
    tags=('agents', 'contact network', 'bipartite', 'commuting', 'episims'),
    references=('Eubank et al. 2004', 'Raman 2021'),
    defaults={'until': 200},
)
def episims_city(
    people: int = 1000,
    *,
    household: float = 3.0,
    workplace: float = 20.0,
    p: float = 0.05,
    gamma: float = 0.1,
    infected: int = 5,
    closed: bool = False,
    seed: int | None = 0,
) -> System:
    """Return a synthetic city in the spirit of EpiSims, where people infect those they meet at home and at work.

    EpiSims models a city as a bipartite graph of people and locations whose edges are visits labelled with their
    times, and spreads disease between people at the same location. Here every person belongs to a household of mean
    size `household` and a workplace or school of mean size `workplace`, the people sleep at home and work by day,
    and the disease spreads by SIR contacts at every location. The locations form a network whose edges join a home to
    the workplaces of its members, the projection of the bipartite graph onto the locations. Households spread the
    disease within families and workplaces between them, so closing the workplaces, one of the interventions EpiSims
    was built to compare, confines outbreaks to households.

    Args:
        people:
            The number of people.

        household:
            The mean size of a household.

        workplace:
            The mean size of a workplace or school.

        p:
            The probability of transmission per infected contact and half day.

        gamma:
            The probability of recovery per half day.

        infected:
            The initial number of infected people.

        closed:
            Whether the workplaces are closed.

        seed:
            The seed of the households, the workplaces, and the initial infections.

    Returns:
        system:
            The system, with the homes, workplaces, and health states of the agent group `people` on the locations.

    Examples:
        >>> from skcomplex.systems.biology.epidemiology import episims_city
        >>> episims_city(100).elements[0].size
        100
    """
    rng = np.random.default_rng(seed)
    homes = rng.integers(0, max(1, round(people / household)), people)
    n_homes = int(homes.max()) + 1
    works = n_homes + rng.integers(0, max(1, round(people / workplace)), people)
    graph = nx.Graph()
    graph.add_nodes_from(range(int(works.max()) + 1))
    graph.add_edges_from(zip(homes.tolist(), works.tolist(), strict=True))
    health = np.zeros(people)
    health[rng.choice(people, size=infected, replace=False)] = 1.0
    group = GroupElement(
        'people',
        element=CommuterElement,
        size=people,
        values={'home': homes.astype(float), 'work': works.astype(float), 'health': health},
        location=homes.tolist(),
    )
    contacts = ColocationInteraction('contacts', people='people', p=p, gamma=gamma, closed=closed)
    return System(elements=[group], interactions=[contacts], space=NetworkSpace(graph), time=DiscreteTime())

heterogeneous_mean_field_sir(exponent=2.5, k_min=3, k_max=100, *, beta=0.1, mu=1.0, infected=0.0001)

Return the SIR epidemic on an uncorrelated network of power-law degrees, in heterogeneous mean field.

The degrees run from k_min to k_max with P(k) ~ k^-exponent. An outbreak grows when beta / mu > <k> / (<k^2> - <k>), a threshold close to <k> / <k^2> on heterogeneous networks, which vanishes on scale-free networks as k_max grows. Below it the epidemic dies out after infecting a share of the order of the initial one, and above it a finite share of the population ends recovered. The hubs are infected first.

Parameters:

Name	Type	Description	Default
exponent	float	The exponent of the degree distribution.	2.5
k_min	int	The smallest degree.	3
k_max	int	The largest degree.	100
beta	float	The rate of transmission along one contact.	0.1
mu	float	The recovery rate.	1.0
infected	float	The initial infected density of every class.	0.0001

Returns:

Name	Type	Description
system	System	The system, with the densities network.s, network.i, and network.r of the degree classes from k_min to k_max, members of the group network, and the shares infected and recovered of the whole population, the last one the final epidemic size at the end of the outbreak.

Examples:

>>> from skcomplex.systems.biology.epidemiology import heterogeneous_mean_field_sir
>>> heterogeneous_mean_field_sir(k_min=2, k_max=5).elements[0].prepare().size_
4

Source code in src/skcomplex/systems/biology/epidemiology/_outbreaks.py

@register_system(
    field='biology/epidemiology',
    tags=('network', 'mean field', 'ode', 'epidemic threshold', 'scale-free', 'final size'),
    references=(
        'Bertin 2021',
        'Moreno, Pastor-Satorras, and Vespignani 2002',
        'Barrat, Barthelemy, and Vespignani 2008',
    ),
    defaults={'until': 100.0, 'dt': 0.5},
    observables=lambda system: [
        ObservableSpec('infected', _infected_share, group='network', reduce='sum'),
        ObservableSpec('recovered', _recovered_share, group='network', reduce='sum'),
    ],
)
def heterogeneous_mean_field_sir(
    exponent: float = 2.5,
    k_min: int = 3,
    k_max: int = 100,
    *,
    beta: float = 0.1,
    mu: float = 1.0,
    infected: float = 1e-4,
) -> System:
    """Return the SIR epidemic on an uncorrelated network of power-law degrees, in heterogeneous mean field.

    The degrees run from `k_min` to `k_max` with `P(k) ~ k^-exponent`. An outbreak grows when
    `beta / mu > <k> / (<k^2> - <k>)`, a threshold close to `<k> / <k^2>` on heterogeneous networks, which vanishes on
    scale-free networks as `k_max` grows. Below it the epidemic dies out after infecting a share of the order of the
    initial one, and above it a finite share of the population ends recovered. The hubs are infected first.

    Args:
        exponent:
            The exponent of the degree distribution.

        k_min:
            The smallest degree.

        k_max:
            The largest degree.

        beta:
            The rate of transmission along one contact.

        mu:
            The recovery rate.

        infected:
            The initial infected density of every class.

    Returns:
        system:
            The system, with the densities `network.s`, `network.i`, and `network.r` of the degree classes from `k_min`
            to `k_max`, members of the group `network`, and the shares `infected` and `recovered` of the whole
            population, the last one the final epidemic size at the end of the outbreak.

    Examples:
        >>> from skcomplex.systems.biology.epidemiology import heterogeneous_mean_field_sir
        >>> heterogeneous_mean_field_sir(k_min=2, k_max=5).elements[0].prepare().size_
        4
    """
    degrees = np.arange(k_min, k_max + 1, dtype=float)
    weights = degrees**-exponent / np.sum(degrees**-exponent)
    classes = GroupElement(
        'network',
        element=ImmunizingDegreeClassElement,
        size=len(degrees),
        values={'s': 1.0 - infected, 'i': infected, 'r': 0.0, 'beta': beta, 'mu': mu, 'k': degrees, 'P': weights},
    )
    interactions = [
        CouplingInteraction(
            'infection',
            target='network',
            read='i',
            write='i',
            function=_infection,
            reduce='normalized-mean',
            mode='group',
            weights=_edge_weight,
            include_self=True,
        ),
        CouplingInteraction(
            'contagion',
            target='network',
            read='s',
            write='s',
            function=_susceptible_loss,
            reduce='normalized-mean',
            mode='group',
            weights=_edge_weight,
            include_self=True,
        ),
        TransferInteraction('recovery', source='network', target='network', state='i', into='r', amount=_recovery),
    ]
    return System(elements=[classes], interactions=interactions)

heterogeneous_mean_field_sis(exponent=2.5, k_min=3, k_max=30, *, lam=0.2, rho=0.01)

Return the SIS epidemic on a network with a power-law degree distribution, in heterogeneous mean field.

The degrees run from k_min to k_max with P(k) ~ k^-exponent. The infection persists above the epidemic threshold lam_c = <k> / <k^2>, which vanishes on scale-free networks as k_max grows, and dies out below it.

Parameters:

Name	Type	Description	Default
exponent	float	The exponent of the degree distribution.	2.5
k_min	int	The smallest degree.	3
k_max	int	The largest degree.	30
lam	float	The spreading rate, the ratio of the infection and recovery rates.	0.2
rho	float	The initial infected density of every class.	0.01

Returns:

Name	Type	Description
system	System	The system, with the densities network.rho of the degree classes from k_min to k_max, members of the group network.

Examples:

>>> from skcomplex.systems.biology.epidemiology import heterogeneous_mean_field_sis
>>> heterogeneous_mean_field_sis(k_min=2, k_max=5).elements[0].prepare().size_
4

Source code in src/skcomplex/systems/biology/epidemiology/_networks.py

@register_system(
    field='biology/epidemiology',
    tags=('network', 'mean field', 'ode', 'epidemic threshold', 'scale-free'),
    references=('Pastor-Satorras and Vespignani 2001', 'Pastor-Satorras et al. 2015'),
    defaults={'until': 50.0, 'dt': 0.5},
)
def heterogeneous_mean_field_sis(
    exponent: float = 2.5,
    k_min: int = 3,
    k_max: int = 30,
    *,
    lam: float = 0.2,
    rho: float = 0.01,
) -> System:
    """Return the SIS epidemic on a network with a power-law degree distribution, in heterogeneous mean field.

    The degrees run from `k_min` to `k_max` with `P(k) ~ k^-exponent`. The infection persists above the epidemic
    threshold `lam_c = <k> / <k^2>`, which vanishes on scale-free networks as `k_max` grows, and dies out below it.

    Args:
        exponent:
            The exponent of the degree distribution.

        k_min:
            The smallest degree.

        k_max:
            The largest degree.

        lam:
            The spreading rate, the ratio of the infection and recovery rates.

        rho:
            The initial infected density of every class.

    Returns:
        system:
            The system, with the densities `network.rho` of the degree classes from `k_min` to `k_max`, members of the
            group `network`.

    Examples:
        >>> from skcomplex.systems.biology.epidemiology import heterogeneous_mean_field_sis
        >>> heterogeneous_mean_field_sis(k_min=2, k_max=5).elements[0].prepare().size_
        4
    """
    degrees = np.arange(k_min, k_max + 1, dtype=float)
    weights = degrees**-exponent / np.sum(degrees**-exponent)
    classes = GroupElement(
        'network',
        element=DegreeClassElement,
        size=len(degrees),
        values={'rho': rho, 'lam': lam, 'k': degrees, 'P': weights},
    )
    interactions = [
        RateInteraction('recovery', target='network', state='rho', rate=_recovery),
        CouplingInteraction(
            'infection',
            target='network',
            read='rho',
            write='rho',
            function=_infection,
            reduce='normalized-mean',
            mode='group',
            weights=_edge_weight,
            include_self=True,
        ),
    ]
    return System(elements=[classes], interactions=interactions)

metapopulation(graph=None, *, population=1000.0, infected=10.0, beta=0.3, gamma=0.1, D=0.01)

Return SIR epidemics in subpopulations at the nodes of a network, linked by individuals who travel.

The infection starts in the first node and spreads to the others with the travellers. The total number of individuals stays constant.

Parameters:

Name	Type	Description	Default
graph	object	The NetworkX graph, or None for a Barabási-Albert graph of 30 nodes.	None
population	float	The initial number of individuals of every node.	1000.0
infected	float	The initial number of infected individuals in the first node.	10.0
beta	float	The infection rate.	0.3
gamma	float	The recovery rate.	0.1
D	float	The rate at which an individual hops along one edge.	0.01

Returns:

Name	Type	Description
system	System	The system, with the compartments cities.S, cities.I, and cities.R, one value per node.

Examples:

>>> import networkx as nx
>>> from skcomplex.systems.biology.epidemiology import metapopulation
>>> metapopulation(nx.path_graph(4)).elements[0].prepare().size_
4

Source code in src/skcomplex/systems/biology/epidemiology/_networks.py

@register_system(
    field='biology/epidemiology',
    tags=('network', 'metapopulation', 'reaction-diffusion', 'ode', 'field'),
    references=('Colizza, Pastor-Satorras, and Vespignani 2007', 'Pastor-Satorras et al. 2015'),
    defaults={'until': 150.0, 'dt': 1.0},
)
def metapopulation(
    graph: object = None,
    *,
    population: float = 1000.0,
    infected: float = 10.0,
    beta: float = 0.3,
    gamma: float = 0.1,
    D: float = 0.01,
) -> System:
    """Return SIR epidemics in subpopulations at the nodes of a network, linked by individuals who travel.

    The infection starts in the first node and spreads to the others with the travellers. The total number of
    individuals stays constant.

    Args:
        graph:
            The NetworkX graph, or `None` for a Barabási-Albert graph of 30 nodes.

        population:
            The initial number of individuals of every node.

        infected:
            The initial number of infected individuals in the first node.

        beta:
            The infection rate.

        gamma:
            The recovery rate.

        D:
            The rate at which an individual hops along one edge.

    Returns:
        system:
            The system, with the compartments `cities.S`, `cities.I`, and `cities.R`, one value per node.

    Examples:
        >>> import networkx as nx
        >>> from skcomplex.systems.biology.epidemiology import metapopulation
        >>> metapopulation(nx.path_graph(4)).elements[0].prepare().size_
        4
    """
    graph = nx.barabasi_albert_graph(30, 2, seed=0) if graph is None else graph
    size = nx.number_of_nodes(graph)
    sick = np.zeros(size)
    sick[0] = infected
    values = {
        'S': np.full(size, population) - sick,
        'I': sick,
        'R': np.zeros(size),
        'beta': beta,
        'gamma': gamma,
        'D': D,
    }
    return System(
        elements=[FieldElement('cities', element=SubpopulationElement, space=NetworkSpace(graph), values=values)],
    )

network_sir(graph=None, *, residents=20, tau=0.05, days=14.0, m=0.05, random_state=0)

Return SIR epidemics in cities at the nodes of a network, carried between cities by travellers.

Every node starts with the same number of susceptible residents, and one resident of a random node is infected.

Parameters:

Name	Type	Description	Default
graph	object	The NetworkX graph of the cities, or None for a Barabási-Albert graph of 20 nodes.	None
residents	int	The initial number of travellers per node.	20
tau	float	The probability of infecting a susceptible traveller of the same node at every step.	0.05
days	float	The number of steps an infection lasts.	14.0
m	float	The probability of moving to a neighbouring node at every step.	0.05
random_state	int | Generator | None	The seed or generator of the first infection.	0

Returns:

Name	Type	Description
system	System	The system, with the group travellers.

Examples:

>>> import networkx as nx
>>> from skcomplex.systems.biology.epidemiology import network_sir
>>> network_sir(nx.cycle_graph(5), residents=3).elements[0].prepare().size_
15

Source code in src/skcomplex/systems/biology/epidemiology/_agents.py

@register_system(
    field='biology/epidemiology',
    tags=('compartments', 'agents', 'network', 'transport'),
    references=('Datseris and Parlitz 2022', 'Thurner et al. 2018'),
    defaults={'until': 150},
)
def network_sir(
    graph: object = None,
    *,
    residents: int = 20,
    tau: float = 0.05,
    days: float = 14.0,
    m: float = 0.05,
    random_state: int | np.random.Generator | None = 0,
) -> System:
    """Return SIR epidemics in cities at the nodes of a network, carried between cities by travellers.

    Every node starts with the same number of susceptible residents, and one resident of a random node is infected.

    Args:
        graph:
            The NetworkX graph of the cities, or `None` for a Barabási-Albert graph of 20 nodes.

        residents:
            The initial number of travellers per node.

        tau:
            The probability of infecting a susceptible traveller of the same node at every step.

        days:
            The number of steps an infection lasts.

        m:
            The probability of moving to a neighbouring node at every step.

        random_state:
            The seed or generator of the first infection.

    Returns:
        system:
            The system, with the group `travellers`.

    Examples:
        >>> import networkx as nx
        >>> from skcomplex.systems.biology.epidemiology import network_sir
        >>> network_sir(nx.cycle_graph(5), residents=3).elements[0].prepare().size_
        15
    """
    rng = np.random.default_rng(random_state)
    graph = nx.barabasi_albert_graph(20, 2, seed=0) if graph is None else graph
    nodes = list(nx.nodes(graph))
    size = residents * len(nodes)
    group = GroupElement(
        'travellers',
        element=TravellerElement,
        size=size,
        values={'health': _initial_health(size, 1, rng, INFECTED), 'tau': tau, 'm': m},
        location=[node for node in nodes for _ in range(residents)],
    )
    infection = TransitionInteraction(
        'infection',
        target='travellers',
        state='health',
        source=SUSCEPTIBLE,
        destination=INFECTED,
        probability=_get_traveller_infection,
        over='cell',
        schedule='member',
        form='agent',
    )
    recovery = TransitionInteraction(
        'recovery',
        target='travellers',
        state='health',
        source=INFECTED,
        destination=RECOVERED,
        condition=_is_due,
        after=math.floor(days) + 1,
        schedule='member',
        form='agent',
    )
    travel = MovementInteraction('travel', target='travellers', accept=_get_traveller_step)
    return System(elements=[group], interactions=[infection, recovery, travel], space=NetworkSpace(graph))

pair_formation_sis(n=400, k=2, *, beta=None, gamma=0.5, sigma=None, pressure=1.5, turnover=2.0, infected=0.05, seed=1)

Return the SIS epidemic on partnerships that dissolve and are replaced at once, with k partners each.

Every person keeps k partners at any time, on a random regular network. A partnership dissolves at the rate sigma, and the person who keeps it finds a new partner at once, at random. An infectious person infects a susceptible partner at the rate beta and recovers at the rate gamma. By default the infection pressure k beta and the rate k sigma of new partners per person are held fixed at pressure and turnover, so k = 1 is serial monogamy and a larger k is concurrency at the same contact rate. Concurrency lets the infection pass on before the partnership that brought it ends, so the epidemic grows faster and reaches more people.

Parameters:

Name	Type	Description	Default
n	int	The number of people, even when k is odd.	400
k	int	The number of concurrent partners of every person.	2
beta	float | None	The infection rate per partnership, or None for pressure / k.	None
gamma	float	The recovery rate.	0.5
sigma	float | None	The dissolution rate of a partnership, or None for turnover / k.	None
pressure	float	The infection pressure of a person whose partners are all infectious, used when beta is None.	1.5
turnover	float	The rate of new partners per person, used when sigma is None.	2.0
infected	float	The initial share of infectious people.	0.05
seed	int | None	The seed of the network and of the initial infections.	1

Returns:

Name	Type	Description
system	System	The system, with the states people.i of the variable group people on the nodes and the variable group partnerships on the edges.

Examples:

>>> from skcomplex.systems.biology.epidemiology import pair_formation_sis
>>> pair_formation_sis(n=20, k=4).interactions[0].beta
0.375

Source code in src/skcomplex/systems/biology/epidemiology/_concurrency.py

@register_system(
    field='biology/epidemiology',
    tags=('concurrency', 'partnerships', 'temporal network', 'sis', 'jumps'),
    references=('Morris and Kretzschmar 1997', 'Masuda, Miller, and Holme 2021'),
    defaults={'until': 40.0, 'dt': 0.5},
)
def pair_formation_sis(
    n: int = 400,
    k: int = 2,
    *,
    beta: float | None = None,
    gamma: float = 0.5,
    sigma: float | None = None,
    pressure: float = 1.5,
    turnover: float = 2.0,
    infected: float = 0.05,
    seed: int | None = 1,
) -> System:
    """Return the SIS epidemic on partnerships that dissolve and are replaced at once, with `k` partners each.

    Every person keeps `k` partners at any time, on a random regular network. A partnership dissolves at the rate
    `sigma`, and the person who keeps it finds a new partner at once, at random. An infectious person infects a
    susceptible partner at the rate `beta` and recovers at the rate `gamma`. By default the infection pressure `k beta`
    and the rate `k sigma` of new partners per person are held fixed at `pressure` and `turnover`, so `k = 1` is serial
    monogamy and a larger `k` is concurrency at the same contact rate. Concurrency lets the infection pass on before
    the partnership that brought it ends, so the epidemic grows faster and reaches more people.

    Args:
        n:
            The number of people, even when `k` is odd.

        k:
            The number of concurrent partners of every person.

        beta:
            The infection rate per partnership, or `None` for `pressure / k`.

        gamma:
            The recovery rate.

        sigma:
            The dissolution rate of a partnership, or `None` for `turnover / k`.

        pressure:
            The infection pressure of a person whose partners are all infectious, used when `beta` is `None`.

        turnover:
            The rate of new partners per person, used when `sigma` is `None`.

        infected:
            The initial share of infectious people.

        seed:
            The seed of the network and of the initial infections.

    Returns:
        system:
            The system, with the states `people.i` of the variable group `people` on the nodes and the variable group
            `partnerships` on the edges.

    Examples:
        >>> from skcomplex.systems.biology.epidemiology import pair_formation_sis
        >>> pair_formation_sis(n=20, k=4).interactions[0].beta
        0.375
    """
    beta = pressure / k if beta is None else beta
    sigma = turnover / k if sigma is None else sigma
    space = NetworkSpace(nx.random_regular_graph(k, n, seed=seed))
    initial = (np.random.default_rng(seed).random(n) < infected).astype(float)
    people = GroupElement(
        'people',
        element=HostElement,
        space=space,
        on='sites',
        variable=True,
        capacity=n,
        values={'i': initial},
        visual=VisualSpec('sites', {'states': ['i'], 'low': 0, 'high': 1}),
    )
    partnerships = GroupElement(
        'partnerships',
        element=Element,
        space=space,
        on='edges',
        variable=True,
        capacity=space.graph.number_of_edges(),
        visual=VisualSpec('sites', {'states': []}),
    )
    return System(
        elements=[people, partnerships],
        interactions=[
            SpreadingInteraction('spreading', target='people', state='i', beta=beta, gamma=gamma),
            TopologyInteraction('dissolution', target='people', action='rewire', rate=sigma, keep='source'),
        ],
    )

quenched_mean_field_sis(graph=None, *, n=300, m=3, lam=0.2, p=0.05, seed=1)

Return the SIS epidemic on a fixed network in quenched mean field, one infection probability per node.

The infection persists above the threshold lam_c = 1 / Lambda, the inverse of the largest eigenvalue of the adjacency matrix, and dies out below it. On scale-free networks the hubs keep it alive at small lam.

Parameters:

Name	Type	Description	Default
graph	Graph | None	The network, or None for a Barabási-Albert network of n nodes with m edges per new node.	None
n	int	The number of nodes of the default network.	300
m	int	The number of edges of every new node of the default network.	3
lam	float	The spreading rate, the ratio of the infection rate per edge and the recovery rate.	0.2
p	float	The initial probability that every node is infected.	0.05
seed	int | None	The seed of the default network.	1

Returns:

Name	Type	Description
system	System	The system, with the probabilities nodes.p of the group nodes on the nodes of the network.

Examples:

>>> import networkx as nx
>>> from skcomplex.systems.biology.epidemiology import quenched_mean_field_sis
>>> quenched_mean_field_sis(nx.cycle_graph(10)).elements[0].prepare().size_
10

Source code in src/skcomplex/systems/biology/epidemiology/_adaptive.py

@register_system(
    field='biology/epidemiology',
    tags=('network', 'quenched mean field', 'ode', 'epidemic threshold', 'spectral radius'),
    references=('Pastor-Satorras et al. 2015', 'Wang et al. 2003'),
    defaults={'until': 40.0, 'dt': 0.2},
)
def quenched_mean_field_sis(
    graph: nx.Graph | None = None,
    *,
    n: int = 300,
    m: int = 3,
    lam: float = 0.2,
    p: float = 0.05,
    seed: int | None = 1,
) -> System:
    """Return the SIS epidemic on a fixed network in quenched mean field, one infection probability per node.

    The infection persists above the threshold `lam_c = 1 / Lambda`, the inverse of the largest eigenvalue of the
    adjacency matrix, and dies out below it. On scale-free networks the hubs keep it alive at small `lam`.

    Args:
        graph:
            The network, or `None` for a Barabási-Albert network of `n` nodes with `m` edges per new node.

        n:
            The number of nodes of the default network.

        m:
            The number of edges of every new node of the default network.

        lam:
            The spreading rate, the ratio of the infection rate per edge and the recovery rate.

        p:
            The initial probability that every node is infected.

        seed:
            The seed of the default network.

    Returns:
        system:
            The system, with the probabilities `nodes.p` of the group `nodes` on the nodes of the network.

    Examples:
        >>> import networkx as nx
        >>> from skcomplex.systems.biology.epidemiology import quenched_mean_field_sis
        >>> quenched_mean_field_sis(nx.cycle_graph(10)).elements[0].prepare().size_
        10
    """
    space = NetworkSpace(graph) if graph is not None else NetworkSpace(nx.barabasi_albert_graph(n, m, seed=seed))
    return System(
        elements=[LatticeElement('nodes', element=CarrierElement, space=space, values={'p': p, 'lam': lam})],
    )

random_walker_epidemics(width=64, height=64, *, walkers=1000, lifetime=20.0, infection=1.0, travel=1.0, graph=None, random_state=0)

Return random walkers among whom one sick walker starts a deadly epidemic, on a grid or on a network.

Sick walkers infect the walkers they meet in a cell and die lifetime steps later, so the epidemic spreads as a front like a forest fire and burns out. Above a critical density of walkers it percolates across the grid. On a network, such as a small world of cities, walkers travel along edges and long links carry the disease far.

Parameters:

Name	Type	Description	Default
width	int	The number of columns of the grid, whose edges walkers cannot cross.	64
height	int	The number of rows of the grid.	64
walkers	int	The number of walkers.	1000
lifetime	float	The number of steps a walker lives after its infection.	20.0
infection	float	The probability that a sick walker infects a healthy walker of its cell at every step.	1.0
travel	float	The probability of moving at every step.	1.0
graph	object	A NetworkX graph whose nodes the walkers occupy instead of the grid, or None for the grid.	None
random_state	int | Generator | None	The seed or generator of the initial places and the first sick walker.	0

Returns:

Name	Type	Description
system	System	The system, with the group walkers.

Examples:

>>> from skcomplex.systems.biology.epidemiology import random_walker_epidemics
>>> float(random_walker_epidemics(walkers=50).elements[0].values['health'].sum())
1.0

Source code in src/skcomplex/systems/biology/epidemiology/_agents.py

@register_system(
    field='biology/epidemiology',
    tags=('agents', 'grid', 'network', 'random walk', 'percolation'),
    references=('Charbonneau 2017',),
    defaults={'until': 300},
)
def random_walker_epidemics(
    width: int = 64,
    height: int = 64,
    *,
    walkers: int = 1000,
    lifetime: float = 20.0,
    infection: float = 1.0,
    travel: float = 1.0,
    graph: object = None,
    random_state: int | np.random.Generator | None = 0,
) -> System:
    """Return random walkers among whom one sick walker starts a deadly epidemic, on a grid or on a network.

    Sick walkers infect the walkers they meet in a cell and die `lifetime` steps later, so the epidemic spreads as a
    front like a forest fire and burns out. Above a critical density of walkers it percolates across the grid. On a
    network, such as a small world of cities, walkers travel along edges and long links carry the disease far.

    Args:
        width:
            The number of columns of the grid, whose edges walkers cannot cross.

        height:
            The number of rows of the grid.

        walkers:
            The number of walkers.

        lifetime:
            The number of steps a walker lives after its infection.

        infection:
            The probability that a sick walker infects a healthy walker of its cell at every step.

        travel:
            The probability of moving at every step.

        graph:
            A NetworkX graph whose nodes the walkers occupy instead of the grid, or `None` for the grid.

        random_state:
            The seed or generator of the initial places and the first sick walker.

    Returns:
        system:
            The system, with the group `walkers`.

    Examples:
        >>> from skcomplex.systems.biology.epidemiology import random_walker_epidemics
        >>> float(random_walker_epidemics(walkers=50).elements[0].values['health'].sum())
        1.0
    """
    rng = np.random.default_rng(random_state)
    if graph is None:
        cells = rng.integers(0, [width, height], size=(walkers, 2))
        locations: list[object] = [tuple(int(item) for item in cell) for cell in cells]
        space: LatticeSpace | NetworkSpace = LatticeSpace((width, height), boundary='fixed', neighborhood='von_neumann')
    else:
        nodes = list(nx.nodes(graph))
        locations = [nodes[int(pick)] for pick in rng.integers(0, len(nodes), size=walkers)]
        space = NetworkSpace(graph)
    health = _initial_health(walkers, 1, rng, SICK)
    group = GroupElement(
        'walkers',
        element=WalkerElement,
        size=walkers,
        values={'health': health, 'infection': infection, 'travel': travel},
        location=locations,
    )
    interactions = [
        MovementInteraction(
            'walk',
            target='walkers',
            condition=_is_living_walker,
            accept=_get_walker_step,
            schedule='member',
        ),
        TransitionInteraction(
            'infection',
            target='walkers',
            state='health',
            source=HEALTHY,
            destination=SICK,
            probability=_get_walker_infection,
            over='cell',
            schedule='member',
            form='agent',
        ),
        TransitionInteraction(
            'death',
            target='walkers',
            state='health',
            source=SICK,
            destination=DECEASED,
            condition=_is_due,
            after=math.ceil(lifetime),
            schedule='member',
            form='agent',
        ),
    ]
    return System(elements=[group], interactions=interactions, space=space)

record_dynamics_epidemic(n=10000, *, seed=0)

Return yearly strains that infect every resident less resistant than their strength, who become immune.

Every resident keeps a resistance R drawn uniformly from [0, 1]. Once a year a strain of strength eta, drawn uniformly from [0, 1), infects the residents with R < eta, whose resistance becomes 1. In a large population an outbreak happens only in a year whose strain is stronger than every earlier one, a record, so outbreaks thin out as time goes on: the mean number of outbreaks in t years is close to the harmonic number H_t, about ln t + 0.577, a Poisson process in the logarithm of time.

Parameters:

Name	Type	Description	Default
n	int	The number of residents.	10000
seed	int | None	The seed of the resistances.	0

Returns:

Name	Type	Description
system	System	The system, with the resistances residents.R, the infections of the last year residents.infected, and the strength pathogen.eta of the strain of the coming year.

Examples:

>>> from skcomplex.systems.biology.epidemiology import record_dynamics_epidemic
>>> record_dynamics_epidemic(100).elements[0].prepare().size_
100

Source code in src/skcomplex/systems/biology/epidemiology/_records.py

@register_system(
    field='biology/epidemiology',
    tags=('record dynamics', 'intermittency', 'immunity', 'stochastic map'),
    references=('Jensen 2023',),
    defaults={'until': 1000, 'dt': 1},
)
def record_dynamics_epidemic(n: int = 10_000, *, seed: int | None = 0) -> System:
    """Return yearly strains that infect every resident less resistant than their strength, who become immune.

    Every resident keeps a resistance `R` drawn uniformly from `[0, 1]`. Once a year a strain of strength `eta`, drawn
    uniformly from `[0, 1)`, infects the residents with `R < eta`, whose resistance becomes 1. In a large population an
    outbreak happens only in a year whose strain is stronger than every earlier one, a record, so outbreaks thin out
    as time goes on: the mean number of outbreaks in `t` years is close to the harmonic number `H_t`, about
    `ln t + 0.577`, a Poisson process in the logarithm of time.

    Args:
        n:
            The number of residents.

        seed:
            The seed of the resistances.

    Returns:
        system:
            The system, with the resistances `residents.R`, the infections of the last year `residents.infected`, and
            the strength `pathogen.eta` of the strain of the coming year.

    Examples:
        >>> from skcomplex.systems.biology.epidemiology import record_dynamics_epidemic
        >>> record_dynamics_epidemic(100).elements[0].prepare().size_
        100
    """
    resistance = np.random.default_rng(seed).random(n)
    residents = GroupElement('residents', element=ResidentElement, size=n, values={'R': resistance})
    return System(
        elements=[residents, PathogenElement('pathogen')],
        interactions=[ResistanceInfectionInteraction('infection', residents='residents', pathogen='pathogen')],
        time=DiscreteTime(),
    )

sir(susceptible=990.0, infected=10.0, recovered=0.0, *, beta=0.0003, gamma=0.1, stochastic=False)

Return an epidemic that moves individuals from susceptible to infected to recovered.

Parameters:

Name	Type	Description	Default
susceptible	float	The initial number of susceptible individuals.	990.0
infected	float	The initial number of infected individuals.	10.0
recovered	float	The initial number of recovered individuals.	0.0
beta	float	The infection rate per susceptible-infected pair.	0.0003
gamma	float	The recovery rate.	0.1
stochastic	bool	Whether infections and recoveries are random events of single individuals.	False

Returns:

Name	Type	Description
system	System	The system, with the compartments S, I, and R.

Examples:

>>> from skcomplex.systems.biology.epidemiology import sir
>>> [element.label for element in sir().elements]
['S', 'I', 'R']

Source code in src/skcomplex/systems/biology/epidemiology/_factories.py

@register_system(
    field='biology/epidemiology',
    tags=('compartments', 'ode', 'jump'),
    references=('Kermack and McKendrick 1927',),
    defaults={'until': 160.0, 'dt': 1.0},
)
def sir(
    susceptible: float = 990.0,
    infected: float = 10.0,
    recovered: float = 0.0,
    *,
    beta: float = 0.0003,
    gamma: float = 0.1,
    stochastic: bool = False,
) -> System:
    """Return an epidemic that moves individuals from susceptible to infected to recovered.

    Args:
        susceptible:
            The initial number of susceptible individuals.

        infected:
            The initial number of infected individuals.

        recovered:
            The initial number of recovered individuals.

        beta:
            The infection rate per susceptible-infected pair.

        gamma:
            The recovery rate.

        stochastic:
            Whether infections and recoveries are random events of single individuals.

    Returns:
        system:
            The system, with the compartments `S`, `I`, and `R`.

    Examples:
        >>> from skcomplex.systems.biology.epidemiology import sir
        >>> [element.label for element in sir().elements]
        ['S', 'I', 'R']
    """
    return System(
        elements=[
            CompartmentElement('S', n=susceptible),
            CompartmentElement('I', n=infected),
            CompartmentElement('R', n=recovered),
        ],
        interactions=[
            InfectionInteraction('infection', susceptible='S', infected='I', beta=beta, stochastic=stochastic),
            RecoveryInteraction('recovery', infected='I', recovered='R', gamma=gamma, stochastic=stochastic),
        ],
    )

sir_vaccination(S=0.99, I=0.01, *, beta=0.5, gamma=0.1, a=1.0, b=1.0)

Return an SIR epidemic with the vaccination rate that minimises the burden of infections and vaccines.

Run with an optimal control solver, vaccination is strongest early, while many are susceptible, and the epidemic peaks lower than without vaccination. With no weight a on the infections the optimal rate is 0, the epidemic without control.

Parameters:

Name	Type	Description	Default
S	float	The initial susceptible fraction.	0.99
I	float	The initial infected fraction.	0.01
beta	float	The transmission rate.	0.5
gamma	float	The recovery rate.	0.1
a	float	The weight of the infected fraction in the cost.	1.0
b	float	The weight of the squared vaccination rate in the cost.	1.0

Returns:

Name	Type	Description
system	System	The system, with the fractions population.S, population.I, and population.R, the vaccination rate population.u, and the observables vaccinated and cost, the running cost.

Examples:

>>> from skcomplex.systems.biology.epidemiology import sir_vaccination
>>> sir_vaccination(beta=0.3).elements[0].beta
0.3

Source code in src/skcomplex/systems/biology/epidemiology/_vaccination.py

@register_system(
    field='biology/epidemiology',
    tags=('optimal control', 'vaccination', 'sir', 'pontryagin', 'ode'),
    references=('Lenhart and Workman 2007',),
    defaults={'until': 30.0, 'dt': 0.2, 'solver': PontryaginSolver(objective=_cost)},
    observables=lambda system: [ObservableSpec('vaccinated', _vaccinated), ObservableSpec('cost', _cost)],
)
def sir_vaccination(
    S: float = 0.99,
    I: float = 0.01,
    *,
    beta: float = 0.5,
    gamma: float = 0.1,
    a: float = 1.0,
    b: float = 1.0,
) -> System:
    """Return an SIR epidemic with the vaccination rate that minimises the burden of infections and vaccines.

    Run with an optimal control solver, vaccination is strongest early, while many are susceptible, and the epidemic
    peaks lower than without vaccination. With no weight `a` on the infections the optimal rate is 0, the epidemic
    without control.

    Args:
        S:
            The initial susceptible fraction.

        I:
            The initial infected fraction.

        beta:
            The transmission rate.

        gamma:
            The recovery rate.

        a:
            The weight of the infected fraction in the cost.

        b:
            The weight of the squared vaccination rate in the cost.

    Returns:
        system:
            The system, with the fractions `population.S`, `population.I`, and `population.R`, the vaccination rate
            `population.u`, and the observables `vaccinated` and `cost`, the running cost.

    Examples:
        >>> from skcomplex.systems.biology.epidemiology import sir_vaccination
        >>> sir_vaccination(beta=0.3).elements[0].beta
        0.3
    """
    population = VaccinatedPopulationElement(
        'population',
        S=S,
        I=I,
        R=max(1 - S - I, 0.0),
        beta=beta,
        gamma=gamma,
        a=a,
        b=b,
    )
    return System(elements=[population])

sirs(susceptible=990.0, infected=10.0, recovered=0.0, *, beta=0.0003, gamma=0.1, xi=0.02, stochastic=False)

Return an epidemic whose immunity wanes, which returns in damped waves to an endemic equilibrium.

The equilibrium has S = gamma / beta susceptible and I = xi (N - S) / (gamma + xi) infected.

Parameters:

Name	Type	Description	Default
susceptible	float	The initial number of susceptible individuals.	990.0
infected	float	The initial number of infected individuals.	10.0
recovered	float	The initial number of recovered individuals.	0.0
beta	float	The infection rate per susceptible-infected pair.	0.0003
gamma	float	The recovery rate.	0.1
xi	float	The rate at which immunity wanes.	0.02
stochastic	bool	Whether infections, recoveries, and losses of immunity are random events of single individuals.	False

Returns:

Name	Type	Description
system	System	The system, with the compartments S, I, and R.

Examples:

>>> from skcomplex.systems.biology.epidemiology import sirs
>>> sirs(xi=0.1).interactions[-1].xi
0.1

Source code in src/skcomplex/systems/biology/epidemiology/_factories.py

@register_system(
    field='biology/epidemiology',
    tags=('compartments', 'ode', 'jump', 'endemic', 'waning immunity'),
    references=('Kermack and McKendrick 1932', 'Argun et al. 2021'),
    defaults={'until': 300.0, 'dt': 1.0},
)
def sirs(
    susceptible: float = 990.0,
    infected: float = 10.0,
    recovered: float = 0.0,
    *,
    beta: float = 0.0003,
    gamma: float = 0.1,
    xi: float = 0.02,
    stochastic: bool = False,
) -> System:
    """Return an epidemic whose immunity wanes, which returns in damped waves to an endemic equilibrium.

    The equilibrium has `S = gamma / beta` susceptible and `I = xi (N - S) / (gamma + xi)` infected.

    Args:
        susceptible:
            The initial number of susceptible individuals.

        infected:
            The initial number of infected individuals.

        recovered:
            The initial number of recovered individuals.

        beta:
            The infection rate per susceptible-infected pair.

        gamma:
            The recovery rate.

        xi:
            The rate at which immunity wanes.

        stochastic:
            Whether infections, recoveries, and losses of immunity are random events of single individuals.

    Returns:
        system:
            The system, with the compartments `S`, `I`, and `R`.

    Examples:
        >>> from skcomplex.systems.biology.epidemiology import sirs
        >>> sirs(xi=0.1).interactions[-1].xi
        0.1
    """
    return System(
        elements=[
            CompartmentElement('S', n=susceptible),
            CompartmentElement('I', n=infected),
            CompartmentElement('R', n=recovered),
        ],
        interactions=[
            InfectionInteraction('infection', susceptible='S', infected='I', beta=beta, stochastic=stochastic),
            RecoveryInteraction('recovery', infected='I', recovered='R', gamma=gamma, stochastic=stochastic),
            ImmunityLossInteraction('waning', recovered='R', susceptible='S', xi=xi, stochastic=stochastic),
        ],
    )

sis(susceptible=990.0, infected=10.0, *, beta=0.0003, gamma=0.1, stochastic=False)

Return an infection without immunity, where the recovered become susceptible again and the disease is endemic.

Above the threshold beta N > gamma, the infected settle at I = N - gamma / beta, and below it they die out.

Parameters:

Name	Type	Description	Default
susceptible	float	The initial number of susceptible individuals.	990.0
infected	float	The initial number of infected individuals.	10.0
beta	float	The infection rate per susceptible-infected pair.	0.0003
gamma	float	The recovery rate.	0.1
stochastic	bool	Whether infections and recoveries are random events of single individuals.	False

Returns:

Name	Type	Description
system	System	The system, with the compartments S and I.

Examples:

>>> from skcomplex.systems.biology.epidemiology import sis
>>> [element.label for element in sis().elements]
['S', 'I']

Source code in src/skcomplex/systems/biology/epidemiology/_factories.py

@register_system(
    field='biology/epidemiology',
    tags=('compartments', 'ode', 'jump', 'endemic'),
    references=('Kermack and McKendrick 1932', 'Argun et al. 2021', 'Datseris and Parlitz 2022'),
    defaults={'until': 100.0, 'dt': 1.0},
)
def sis(
    susceptible: float = 990.0,
    infected: float = 10.0,
    *,
    beta: float = 0.0003,
    gamma: float = 0.1,
    stochastic: bool = False,
) -> System:
    """Return an infection without immunity, where the recovered become susceptible again and the disease is endemic.

    Above the threshold `beta N > gamma`, the infected settle at `I = N - gamma / beta`, and below it they die out.

    Args:
        susceptible:
            The initial number of susceptible individuals.

        infected:
            The initial number of infected individuals.

        beta:
            The infection rate per susceptible-infected pair.

        gamma:
            The recovery rate.

        stochastic:
            Whether infections and recoveries are random events of single individuals.

    Returns:
        system:
            The system, with the compartments `S` and `I`.

    Examples:
        >>> from skcomplex.systems.biology.epidemiology import sis
        >>> [element.label for element in sis().elements]
        ['S', 'I']
    """
    return System(
        elements=[CompartmentElement('S', n=susceptible), CompartmentElement('I', n=infected)],
        interactions=[
            InfectionInteraction('infection', susceptible='S', infected='I', beta=beta, stochastic=stochastic),
            RecoveryInteraction('recovery', infected='I', recovered='S', gamma=gamma, stochastic=stochastic),
        ],
    )

watts_may_concurrency(infected=0.01, *, c=2.0, beta=0.2, tau=1.0, T=0.5, nu=0.1)

Return the SEIR epidemic of Watts and May, where partnerships that last longer spread the infection faster.

A susceptible person is infected at the rate pi(t) = c beta [I / (S + E + I) + int exp(-(t - t') / tau) pi(t' - T) S(t' - T) / (S + E + I)(t') dt'], by a new partner who is infectious or by a current partner who turned infectious during the partnership. The mean number of current partners is c tau, so a longer partnership raises the mean degree of the momentary network, and the infection rate grows with tau at fixed c and beta.

Parameters:

Name	Type	Description	Default
infected	float	The initial fraction of infectious people.	0.01
c	float	The rate of acquiring new partners.	2.0
beta	float	The probability of infection by an infectious partner over a partnership.	0.2
tau	float	The mean duration of a partnership.	1.0
T	float	The duration of the exposed period.	0.5
nu	float	The rate at which the infectious develop the disease.	0.1

Returns:

Name	Type	Description
system	System	The system, with the fractions population.S, population.E, population.I, population.R, and the infection memory population.J.

Examples:

>>> from skcomplex.systems.biology.epidemiology import watts_may_concurrency
>>> watts_may_concurrency(tau=3.0).elements[0].tau
3.0

Source code in src/skcomplex/systems/biology/epidemiology/_concurrency.py

@register_system(
    field='biology/epidemiology',
    tags=('concurrency', 'partnerships', 'seir', 'delay', 'hiv'),
    references=('Watts and May 1992', 'Masuda, Miller, and Holme 2021'),
    defaults={'until': 40.0, 'dt': 0.25},
)
def watts_may_concurrency(
    infected: float = 0.01,
    *,
    c: float = 2.0,
    beta: float = 0.2,
    tau: float = 1.0,
    T: float = 0.5,
    nu: float = 0.1,
) -> System:
    """Return the SEIR epidemic of Watts and May, where partnerships that last longer spread the infection faster.

    A susceptible person is infected at the rate `pi(t) = c beta [I / (S + E + I) + int exp(-(t - t') / tau) pi(t' - T)
    S(t' - T) / (S + E + I)(t') dt']`, by a new partner who is infectious or by a current partner who turned infectious
    during the partnership. The mean number of current partners is `c tau`, so a longer partnership raises the mean
    degree of the momentary network, and the infection rate grows with `tau` at fixed `c` and `beta`.

    Args:
        infected:
            The initial fraction of infectious people.

        c:
            The rate of acquiring new partners.

        beta:
            The probability of infection by an infectious partner over a partnership.

        tau:
            The mean duration of a partnership.

        T:
            The duration of the exposed period.

        nu:
            The rate at which the infectious develop the disease.

    Returns:
        system:
            The system, with the fractions `population.S`, `population.E`, `population.I`, `population.R`, and the
            infection memory `population.J`.

    Examples:
        >>> from skcomplex.systems.biology.epidemiology import watts_may_concurrency
        >>> watts_may_concurrency(tau=3.0).elements[0].tau
        3.0
    """
    population = PartnershipEpidemicElement(
        'population',
        S=1 - infected,
        I=infected,
        c=c,
        beta=beta,
        tau=tau,
        T=T,
        nu=nu,
    )
    return System(elements=[population])

zombie_outbreak(susceptible=500.0, zombies=1.0, removed=0.0, *, birth=0.0, beta=0.0095, alpha=0.005, delta=0.0001, zeta=0.0001, stochastic=False)

Return the outbreak of zombies of Munz et al., where the dead rise again and humans are eventually overrun.

Humans S are born at the rate birth and die naturally at the rate delta. A zombie Z bites a human at the rate beta S Z, and humans destroy zombies at the rate alpha S Z. The dead and destroyed R rise as zombies at the rate zeta. Without births the total S + Z + R is conserved, and with resurrection the only stable state is the doomsday state without humans, whatever the initial numbers, whenever beta > alpha.

Parameters:

Name	Type	Description	Default
susceptible	float	The initial number of humans.	500.0
zombies	float	The initial number of zombies.	1.0
removed	float	The initial number of dead and destroyed.	0.0
birth	float	The birth rate of humans.	0.0
beta	float	The rate at which zombies turn humans.	0.0095
alpha	float	The rate at which humans destroy zombies.	0.005
delta	float	The rate of natural death of humans.	0.0001
zeta	float	The rate at which the removed rise as zombies.	0.0001
stochastic	bool	Whether every event moves one individual at random times.	False

Returns:

Name	Type	Description
system	System	The system, with the compartments S, Z, and R.

Examples:

>>> from skcomplex.systems.biology.epidemiology import zombie_outbreak
>>> [element.label for element in zombie_outbreak().elements]
['S', 'Z', 'R']

Source code in src/skcomplex/systems/biology/epidemiology/_factories.py

@register_system(
    field='biology/epidemiology',
    tags=('compartments', 'ode', 'jump', 'outbreak', 'zombies'),
    references=('Munz et al. 2009', 'Raman 2021'),
    defaults={'until': 10.0, 'dt': 0.05},
)
def zombie_outbreak(
    susceptible: float = 500.0,
    zombies: float = 1.0,
    removed: float = 0.0,
    *,
    birth: float = 0.0,
    beta: float = 0.0095,
    alpha: float = 0.005,
    delta: float = 0.0001,
    zeta: float = 0.0001,
    stochastic: bool = False,
) -> System:
    """Return the outbreak of zombies of Munz et al., where the dead rise again and humans are eventually overrun.

    Humans `S` are born at the rate `birth` and die naturally at the rate `delta`. A zombie `Z` bites a human at the
    rate `beta S Z`, and humans destroy zombies at the rate `alpha S Z`. The dead and destroyed `R` rise as zombies at
    the rate `zeta`. Without births the total `S + Z + R` is conserved, and with resurrection the only stable state
    is the doomsday state without humans, whatever the initial numbers, whenever `beta > alpha`.

    Args:
        susceptible:
            The initial number of humans.

        zombies:
            The initial number of zombies.

        removed:
            The initial number of dead and destroyed.

        birth:
            The birth rate of humans.

        beta:
            The rate at which zombies turn humans.

        alpha:
            The rate at which humans destroy zombies.

        delta:
            The rate of natural death of humans.

        zeta:
            The rate at which the removed rise as zombies.

        stochastic:
            Whether every event moves one individual at random times.

    Returns:
        system:
            The system, with the compartments `S`, `Z`, and `R`.

    Examples:
        >>> from skcomplex.systems.biology.epidemiology import zombie_outbreak
        >>> [element.label for element in zombie_outbreak().elements]
        ['S', 'Z', 'R']
    """
    mode = 'jump' if stochastic else 'rate'
    return System(
        elements=[
            CompartmentElement('S', n=susceptible),
            CompartmentElement('Z', n=zombies),
            CompartmentElement('R', n=removed),
        ],
        interactions=[
            InfectionInteraction('bites', susceptible='S', infected='Z', beta=beta, stochastic=stochastic),
            ReactionInteraction(
                'destruction',
                reactants={('Z', 'n'): 1},
                products={('R', 'n'): 1},
                modifiers={('S', 'n'): 1},
                k=alpha,
                mode=mode,
            ),
            ReactionInteraction('death', reactants={('S', 'n'): 1}, products={('R', 'n'): 1}, k=delta, mode=mode),
            ReactionInteraction('resurrection', reactants={('R', 'n'): 1}, products={('Z', 'n'): 1}, k=zeta, mode=mode),
            ReactionInteraction('birth', products={('S', 'n'): 1}, k=birth, mode=mode),
        ],
    )