graph_tool.dynamics.SIRSState#

class graph_tool.dynamics.SIRSState(g, beta=1, gamma=0.1, mu=0.1, r=0, exposed=False, epsilon=0.1, v0=None, s=None, constant_beta=True)[source]#

Bases: DiscreteStateBase

SIRS compartmental epidemic model.

Parameters:
gGraph

Graph to be used for the dynamics

betafloat or EdgePropertyMap (optional, default: 1.)

Transmission probability.

gammafloat or VertexPropertyMap (optional, default: .1)

I to R recovery probability.

mufloat or VertexPropertyMap (optional, default: .1)

R to S recovery probability.

rfloat or VertexPropertyMap (optional, default: 0.)

Spontaneous infection probability.

exposedboolean (optional, default: False)

If True, an SEIRS model is simulated, with an additional “exposed” state.

epsilonfloat or VertexPropertyMap (optional, default: .1)

Susceptible to exposed transition probability. This only has an effect if exposed=True.

v0int or Vertex (optional, default: None)

Initial infectious vertex. If not provided, and if the global state is also not provided via paramter s, a random vertex will be chosen.

sVertexPropertyMap (optional, default: None)

Initial global state. If not provided, all vertices will be initialized to the susceptible state.

constant_betaboolean (optional, default: True)

If True, and beta is an edge property map, it will be assumed that the beta values do not change, such that the probability values can be pre-computed for efficiency. If beta is a float, this option has no effect.

Notes

This implements an SIRS epidemic process [pastor-satorras-epidemic-2015], where nodes in the susceptible state (value 0) are infected by neighbours in the infectious state (value 1), which can then eventually recover to a recovered state (value 2), and finally back to the susceptible state.

If a node \(i\) is updated at time \(t\), the transition probabilities from state \(s_i(t)\) to state \(s_i(t+1)\) are given as follows:

  1. If \(s_i(t) = 0\), we have \(s_i(t+1) = 1\) with probability
    \[(1-r_i)\left[1-\prod_j(1-\beta_{ij})^{A_{ij}\delta_{s_j(t),1}}\right] + r_i,\]

    otherwise \(s_i(t+1) = 0\).

  2. If \(s_i(t) = 1\), we have \(s_i(t+1) = 2\) with probability \(\gamma_i\), or \(s_i(t+1) = 1\) with probability \(1-\gamma_i\).

  3. If \(s_i(t) = 2\), we have \(s_i(t+1) = 1\) with probability \(\mu_i\), or \(s_i(t+1) = 2\) with probability \(1-\mu_i\).

If the option exposed == True is given, then the states transit first from 0 to -1 (exposed) with probability given by 1. above, and then finally from -1 to 1 with probability \(\epsilon_i\).

References

[pastor-satorras-epidemic-2015]

Romualdo Pastor-Satorras, Claudio Castellano, Piet Van Mieghem, and Alessandro Vespignani, “Epidemic processes in complex networks”, Rev. Mod. Phys. 87, 925 (2015) DOI: 10.1103/RevModPhys.87.925 [sci-hub, @tor], arXiv: 1408.2701

Examples

>>> g = gt.collection.data["pgp-strong-2009"]
>>> state = gt.SIRSState(g, beta=0.2, gamma=0.025, mu=0.02)
>>> S, X, R = [], [], []
>>> for t in range(2000):
...     ret = state.iterate_sync()
...     s = state.get_state().fa
...     S.append((s == 0).sum())
...     X.append((s == 1).sum())
...     R.append((s == 2).sum())

>>> figure(figsize=(6, 4))
<...>
>>> plot(S, label="Susceptible")
[...]
>>> plot(X, label="Infectious")
[...]
>>> plot(R, label="Recovered")
[...]
>>> xlabel(r"Time")
Text(...)
>>> ylabel(r"Number of nodes")
Text(...)
>>> legend(loc="best")
<...>
>>> tight_layout()
>>> savefig("SIRS.svg")
../_images/SIRS.svg

Number of susceptible, infectious, and recovered nodes vs. time for an SIRS dynamics.#

Methods

copy()

Return a copy of the state.

get_active()

Returns list of "active" nodes, for states where this concept is used.

get_state()

Returns the internal VertexPropertyMap with the current state.

iterate_async([niter])

Updates nodes asynchronously (i.e. single vertex chosen randomly), niter number of times.

iterate_sync([niter])

Updates nodes synchronously (i.e. a full "sweep" of all nodes in parallel), niter number of times.

reset_active()

Resets list of "active" nodes, for states where this concept is used.

copy()#

Return a copy of the state.

get_active()#

Returns list of “active” nodes, for states where this concept is used.

get_state()#

Returns the internal VertexPropertyMap with the current state.

iterate_async(niter=1)#

Updates nodes asynchronously (i.e. single vertex chosen randomly), niter number of times. This function returns the number of nodes that changed state.

iterate_sync(niter=1)#

Updates nodes synchronously (i.e. a full “sweep” of all nodes in parallel), niter number of times. This function returns the number of nodes that changed state.

If enabled during compilation, this algorithm runs in parallel (i.e. using more than one thread.)

reset_active()#

Resets list of “active” nodes, for states where this concept is used.