common_models

A set of commonly used models

pygom.model.common_models.FitzHugh(param=None)

The standard FitzHugh model without external input [FitzHugh1961]

\[\begin{split}\frac{dV}{dt} &= c ( V - \frac{V^{3}}{3} + R) \\ \frac{dR}{dt} &= -\frac{1}{c}(V - a + bR).\end{split}\]

Examples

>>> import numpy as np
>>> from pygom import common_models
>>> ode = common_models.FitzHugh({'a':0.2, 'b':0.2, 'c':3.0})
>>> t = np.linspace(0, 20, 101)
>>> x0 = [1.0, -1.0]
>>> ode.initial_values = (x0, t[0])
>>> solution = ode.integrate(t[1::])
>>> ode.plot()

pygom.model.common_models.Influenza_SLIARN(param=None)

A simple influenza model from [Brauer2008], page 323.

\[\begin{split}\frac{dS}{dt} &= -S \beta (I + \delta A) \\ \frac{dL}{dt} &= S \beta (I + \delta A) - \kappa L \\ \frac{dI}{dt} &= p \kappa L - \alpha I \\ \frac{dA}{dt} &= (1 - p) \kappa L - \eta A \\ \frac{dR}{dt} &= f \alpha I + \eta A \\ \frac{dN}{dt} &= -(1 - f) \alpha I\end{split}\]

pygom.model.common_models.Legrand_Ebola_SEIHFR(param=None)

The Legrand Ebola model [Legrand2007] with 6 compartments that includes the H = hospitalization and F = funeral state. Note that because this is an non-autonomous system, there are in fact a total of 7 states after conversion. The set of equations that describes the model are

\[\begin{split}\frac{dS}{dt} &= -(\beta_{I}SI + \beta_{H}SH + \beta_{F}SF) \\ \frac{dE}{dt} &= (\beta_{I}SI + \beta_{H}SH + \beta_{F}SF) - \alpha E \\ \frac{dI}{dt} &= \alpha E - (\gamma_{H} \theta_{1} + \gamma_{I}(1-\theta_{1})(1-\delta_{1}) + \gamma_{D}(1-\theta_{1})\delta_{1})I \\ \frac{dH}{dt} &= \gamma_{H}\theta_{1}I - (\gamma_{DH}\delta_{2} + \gamma_{IH}(1-\delta_{2}))H \\ \frac{dF}{dt} &= \gamma_{D}(1-\theta_{1})\delta_{1}I + \gamma_{DH}\delta_{2}H - \gamma_{F}F \\ \frac{dR}{dt} &= \gamma_{I}(1-\theta_{1})(1-\delta_{1})I + \gamma_{IH}(1-\delta_{2})H + \gamma_{F}F.\end{split}\]

Examples

>>> import numpy as np
>>> from pygom import common_models
>>> x0 = [1.0, 3.0/200000.0, 0.0, 0.0, 0.0, 0.0, 0.0]
>>> t = np.linspace(0, 25, 100)
>>> ode = common_models.Legrand_Ebola_SEIHFR([('beta_I',0.588),('beta_H',0.794),('beta_F',7.653),('omega_I',10.0/7.0),('omega_D',9.6/7.0),('omega_H',5.0/7.0),('omega_F',2.0/7.0),('alphaInv',7.0/7.0),('delta',0.81),('theta',0.80),('kappa',300.0),('interventionTime',7.0)])
>>> ode.initial_values = (x0, t[0])
>>> solution = ode.integrate(t[1::])
>>> ode.plot()

pygom.model.common_models.Lorenz(param=None)

Lorenz attractor define by three parameters, \(\beta,\sigma,\rho\) as per [Lorenz1963].

\[\begin{split}\frac{dx}{dt} &= \sigma (y-x) \\ \frac{dy}{dt} &= x (\rho - z) - y \\ \frac{dz}{dt} &= xy - \beta z\end{split}\]

Examples

>>> import matplotlib.pyplot as plt
>>> import numpy
>>> from pygom import common_models
>>> t = numpy.linspace(0, 20, 101)
>>> params = {'beta':8.0/3.0, 'sigma':10.0, 'rho':28.0}
>>> ode = common_models.Lorenz(params)
>>> ode.initial_values = ([1., 1., 1.], t[0])
>>> solution = ode.integrate(t[1::])
>>> plt.plot(solution[:,0], solution[:,2])
>>> plt.show()

pygom.model.common_models.Lotka_Volterra(param=None)

Standard Lotka-Volterra model with two states and four parameters [Lotka1920]

\[\begin{split}\frac{dx}{dt} &= \alpha x - cxy \\ \frac{dy}{dt} &= -\delta y + \gamma xy\end{split}\]

Examples

>>> import numpy as np
>>> from pygom import common_models
>>> params = {'alpha':1, 'delta':3, 'c':2, 'gamma':6}
>>> ode = common_models.Lotka_Volterra(params)
>>> ode.initial_values = ([2.0, 6.0], 0)
>>> t = np.linspace(0.1, 100, 10000)
>>> ode.integrate(t)
>>> ode.plot()

pygom.model.common_models.Lotka_Volterra_4State(param=None)

The four state Lotka-Volterra model [Lotka1920]. A common interpretation is that a = Grass, x = rabbits, y = foxes and b is the death of foxes.

\[\begin{split}\frac{da}{dt} &= k_{0} a x \\ \frac{dx}{dt} &= k_{0} a x - k_{1} x y \\ \frac{dy}{dt} &= k_{1} x y - k_{2} y \\ \frac{db}{dt} &= k_{2} y\end{split}\]

Examples

>>> import numpy as np
>>> from pygom import common_models
>>> x0 = [150.0, 10.0, 10.0, 0.0]
>>> t = np.linspace(0, 15, 100)
>>> params = [0.01, 0.1, 1.0]
>>> ode = common_models.Lotka_Volterra_4State(params)
>>> ode.initial_values = (x0, t[0])
>>> ode.integrate(t[1::])
>>> ode.plot()

pygom.model.common_models.Robertson(param=None)

The so called Robertson problem [Robertson1966], which is a standard example used to test stiff integrator.

\[\begin{split}\frac{dy_{1}}{dt} &= -0.04 y_{1} + 1 \cdot 10^{4} y_{2} y_{3} \\ \frac{dy_{2}}{dt} &= 0.04 y_{1} - 1 \cdot 10^{4} y_{2} y_{3} - 3 \cdot 10^{7} y_{2}^{2}\\ \frac{dy_{3}}{dt} &= 3 \cdot 10^{7} y_{2}^{2}\end{split}\]

Examples

>>> from pygom import common_models
>>> import numpy
>>> t = numpy.append(0, 4*numpy.logspace(-6, 6, 1000))
>>> ode = common_models.Robertson()
>>> ode.initial_values = ([1.0,0.0,0.0], t[0])
>>> solution = ode.integrate(t[1::])
>>> ode.plot() # note that this is not being plotted in the log scale

pygom.model.common_models.SEIR(param=None)

A standard SEIR model [Brauer2008], defined by the ode

\[\begin{split}\frac{dS}{dt} &= -\beta SI \\ \frac{dE}{dt} &= \beta SI - \alpha E \\ \frac{dI}{dt} &= \alpha E - \gamma I \\ \frac{dR}{dt} &= \gamma I\end{split}\]

See also

SEIR_Birth_Death()

Examples

>>> import numpy as np
>>> from pygom import common_models
>>> ode = common_models.SEIR({'beta':1800, 'gamma':100, 'alpha':35.84})
>>> t = np.linspace(0, 50, 1001)
>>> x0 = [0.0658, 0.0007, 0.0002, 0.0]
>>> ode.initial_values = (x0, t[0])
>>> solution,output = ode.integrate(t[1::], full_output=True)
>>> ode.plot()

pygom.model.common_models.SEIR_Birth_Death(param=None)

A standard SEIR model with birth and death [Aron1984], defined by the ode

\[\begin{split}\frac{dS}{dt} &= \mu - \beta SI - \mu S \\ \frac{dE}{dt} &= \beta SI - (\mu + \alpha) E \\ \frac{dI}{dt} &= \alpha E - (\mu + \gamma) I \\ \frac{dR}{dt} &= \gamma I\end{split}\]

See also

SEIR()

Examples

Uses the same set of parameters as the examples in SEIR() apart from \(\mu\) which is new.

>>> import numpy as np
>>> from pygom import common_models
>>> params = {'beta':1800, 'gamma':100, 'alpha':35.84, 'mu':0.02}
>>> ode = common_models.SEIR_Birth_Death(params)
>>> t = np.linspace(0, 50, 1001)
>>> x0 = [0.0658, 0.0007, 0.0002, 0.0]
>>> ode.initial_values = (x0, t[0])
>>> solution,output = ode.integrate(t[1::], full_output=True)
>>> ode.plot()

pygom.model.common_models.SEIR_Birth_Death_Periodic(param=None)

A SEIR birth death model with periodic contact [Aron1984], defined by the ode

\[\begin{split}\frac{dS}{dt} &= \mu - \beta(t)SI - \mu S \\ \frac{dE}{dt} &= \beta(t)SI - (\mu + \alpha) E \\ \frac{dI}{dt} &= \alpha E - (\mu + \gamma) I \\ \frac{dR}{dt} &= \gamma I\end{split}\]

where

\[\beta(t) = \beta_{0} (1 + \beta_{1} \cos(2 \pi t)).\]

An extension of an SEIR birth death model by varying the contact rate through time.

See also

SEIR(), SEIR_Birth_Death(), SIR_Periodic()

Examples

Uses the same set of parameters as the examples in SEIR_Birth_Death() but now we have two beta parameters instead of one.

>>> import matplotlib.pyplot as plt
>>> import numpy as np
>>> from pygom import common_models
>>> params = {'beta0':1800, 'beta1':0.2, 'gamma':100, 'alpha':35.84, 'mu':0.02}
>>> ode = common_models.SEIR_Birth_Death_Periodic(params)
>>> t = np.linspace(0, 50, 1001)
>>> x0 = [0.0658, 0.0007, 0.0002, 0.0]
>>> ode.initial_values = (x0, t[0])
>>> solution,output = ode.integrate(t[1::], full_output=True)
>>> ode.plot()
>>> plt.plot(np.log(solution[:,0]), np.log(solution[:,1]))
>>> plt.show()
>>> plt.plot(np.log(solution[:,0]), np.log(solution[:,2]))
>>> plt.show()

pygom.model.common_models.SEIR_Multiple(n=2, param=None)

An SEIR model that describe spatial heterogeneity [Brauer2008], page 180. The model originated from [Lloyd1996] and notations used here follows [Brauer2008].

\[\begin{split}\frac{dS_{i}}{dt} &= dN_{i} - dS_{i} - \lambda_{i} S_{i} \\ \frac{dE_{i}}{dt} &= \lambda_{i}S_{i} - (d + \epsilon) E_{i} \\ \frac{dI_{i}}{dt} &= \epsilon E_{i} - (d + \gamma) I_{i} \\ \frac{dR_{i}}{dt} &= \gamma I_{i} - dR_{i}\end{split}\]

where

\[\lambda_{i} = \sum_{j=1}^{n} \beta_{i,j} I_{j} (1\{i \neq j\} p)\]

with \(n\) being the number of patch and \(p\) the coupled factor.

Examples

Use the initial conditions that were derived from the stationary condition specified in [Brauer2008].

>>> import numpy as np
>>> from pygom import common_models
>>> paramEval = {'beta_00':0.0010107, 'beta_01':0.0010107,
>>>              'beta_10':0.0010107, 'beta_11':0.0010107,
>>>              'd':0.02,'epsilon':45.6, 'gamma':73.0,
>>>              'N_0':10**6,'N_1':10**6,'p':0.01}
>>> x0 = [36139.3224081278, 422.560577637822,
>>>       263.883351688369, 963174.233662546]
>>> ode = common_models.SEIR_Multiple()
>>> t = np.linspace(0, 40, 100)
>>> x01 = []
>>> for s in x0:
>>>     x01 += [s]
>>>     x01 += [s]
>>> ode.parameters = paramEval
>>> ode.initial_values = (x01, t[0])
>>> solution, output = ode.integrate(t[1::], full_output=True)
>>> ode.plot()

pygom.model.common_models.SIR(param=None)

A standard SIR model as per [Brauer2008]

\[\begin{split}\frac{dS}{dt} &= -\beta SI \\ \frac{dI}{dt} &= \beta SI - \gamma I \\ \frac{dR}{dt} &= \gamma I\end{split}\]

Examples

The model that produced top two graph in Figure 1.3 of the reference above. First, when everyone is susceptible and only one individual was infected.

>>> import numpy as np
>>> from pygom import common_models
>>> ode = common_models.SIR({'beta':3.6, 'gamma':0.2})
>>> t = np.linspace(0, 730, 1001)
>>> N = 7781984.0
>>> x0 = [1.0, 10/N, 0.0]
>>> ode.initial_values = (x0, t[0])
>>> solution = ode.integrate(t[1::])
>>> ode.plot()

Second model with a more realistic scenario

>>> import numpy as np
>>> from pygom import common_models
>>> ode = common_models.SIR({'beta':3.6, 'gamma':0.2})
>>> t = np.linspace(0, 730, 1001)
>>> N = 7781984.0
>>> x0 = [0.065, 123*(5.0/30.0)/N, 0.0]
>>> ode.initial_values = (x0, t[0])
>>> solution = ode.integrate(t[1::])
>>> ode.plot()

pygom.model.common_models.SIR_Birth_Death(param=None)

Extension of the standard SIR model [Brauer2008] to also include birth and death

\[\begin{split}\frac{dS}{dt} &= B -\beta SI - \mu S \\ \frac{dI}{dt} &= \beta SI - \gamma I - \mu I \\ \frac{dR}{dt} &= \gamma I\end{split}\]

See also

SIR()

Examples

The model that produced bottom graph in Figure 1.3 of the reference above.

>>> import numpy as np
>>> from pygom import common_models
>>> B = 126372.0/365.0
>>> N = 7781984.0
>>> params = {'beta':3.6, 'gamma':0.2, 'B':B/N, 'mu':B/N}
>>> ode = common_models.SIR_Birth_Death(params)
>>> t = np.linspace(0, 35*365, 10001)
>>> x0 = [0.065, 123.0*(5.0/30.0)/N, 0.0]
>>> ode.initial_values = (x0, t[0])
>>> solution,output = ode.integrate(t[1::], full_output=True)
>>> ode.plot()

pygom.model.common_models.SIR_N(param=None)

A standard SIR model [Brauer2008] with population N. This is the unnormalized version of the SIR model.

\[\begin{split}\frac{dS}{dt} &= -\beta SI/N \\ \frac{dI}{dt} &= \beta SI/N- \gamma I \\ \frac{dR}{dt} &= \gamma I\end{split}\]

Examples

The model that produced top two graph in Figure 1.3 of the reference above. First, when everyone is susceptible and only one individual was infected.

>>> import numpy as np
>>> from pygom import common_models
>>> ode = common_models.SIR({'beta':3.6, 'gamma':0.2})
>>> t = np.linspace(0, 730, 1001)
>>> N = 7781984.0
>>> x0 = [N, 1.0, 0.0]
>>> ode.initial_values = (x0, t[0])
>>> solution = ode.integrate(t[1::])
>>> ode.plot()

Second model with a more realistic scenario

>>> import numpy as np
>>> from pygom import common_models
>>> ode = common_models.SIR({'beta':3.6, 'gamma':0.2})
>>> t = np.linspace(0, 730, 1001)
>>> N = 7781984.0
>>> x0 = [int(0.065*N), 21.0, 0.0]
>>> ode.initial_values = (x0, t[0])
>>> solution = ode.integrate(t[1::])
>>> ode.plot()

pygom.model.common_models.SIS(param=None)

A standard SIS model

\[\begin{split}\frac{dS}{dt} &= -\beta SI + \gamma I \\ \frac{dI}{dt} &= \beta SI - \gamma I\end{split}\]

Examples

>>> import numpy as np
>>> from pygom import common_models
>>> ode = common_models.SIS({'beta':0.5, 'gamma':0.2})
>>> t = np.linspace(0, 20, 101)
>>> x0 = [1.0, 0.1]
>>> ode.initial_values = (x0, t[0])
>>> solution = ode.integrate(t[1::])
>>> ode.plot()

pygom.model.common_models.SIS_Periodic(param=None)

A SIS model with periodic contact, defined by the ode as per [Hethcote1973]

\[\frac{dI}{dt} = (\beta(t)N - \alpha) I - \beta(t)I^{2}\]

where

\[\beta(t) = 2 - 1.8 \cos(5t).\]

As the name suggests, it achieves a (stable) periodic solution.

Examples

>>> from pygom import common_models
>>> import numpy as np
>>> ode = common_models.SIS_Periodic({'alpha':1.0})
>>> t = np.linspace(0, 10, 101)
>>> x0 = [0.1, 0.0]
>>> ode.initial_values = (x0, t[0])
>>> solution = ode.integrate(t[1::])
>>> ode.plot()

pygom.model.common_models.vanDelPol(param=None)

The van der Pol equation [vanderpol1926], a second order ode

\[y^{\prime\prime} - \mu (1-y^{2}) y^{\prime} + y = 0\]

where \(\mu > 0\). This can be converted to a first order ode by equating \(x = y^{\prime}\)

\[x^{\prime} - \mu (1 - y^{2}) x + y = 0\]

which result in a coupled ode

\[\begin{split}x^{\prime} &= \mu (1 - y^{2}) x - y \\ y^{\prime} &= x\end{split}\]

and this can be solved via standard method.

Examples

>>> from pygom import common_models
>>> import numpy
>>> t = numpy.linspace(0, 20, 1000)
>>> ode = common_models.vanDelPol({'mu':1.0})
>>> ode.initial_values = ([2.0,0.0], t[0])
>>> solution = ode.integrate(t[1::])
>>> ode.plot()