Mathematical description

In order to solve the time evolution of the compartmental model, we require some translation into mathematics. The different solving methods of PyGOM is a subject for a later chapter, for now we introduce the two ways such systems may be viewed. Recall the model of the previous section:

Ordinary Differential Equations

One way to frame these mathematically is via Ordinary Differential Equations (ODEs): a relationship between quantities and their rates of change with respect to one independent variable - in this case, time. In general, \(n\) variables, \(\mathbf{y}(t) = \left(y_{1}(t),\ldots,y_{n}(t)\right)\), describe the state of the system at time, \(t\). The rates of change of these variables are governed by the current state of the system, time and parameters, \(\boldsymbol{\theta} = \left(\theta_{1},\ldots, \theta_{p}\right)\):

\[\begin{aligned} \frac{\mathrm{d} \mathbf{y}}{\mathrm{d} t} = f(\mathbf{y},\boldsymbol{\theta}, t) \end{aligned}\]

where \(f\) is a vector function giving the rates of change of each dependent variable. For the example in question, the ODEs are:

\[\begin{split}\begin{aligned} \frac{\mathrm{d} S}{\mathrm{d} t} &= -\frac{\beta S I}{N} + \mu N - \mu S\\ \frac{\mathrm{d} I}{\mathrm{d} t} &= \frac{\beta S I}{N} - \gamma I - \mu I\\ \frac{\mathrm{d} R}{\mathrm{d} t} &= \gamma I - \mu R\\ \frac{\mathrm{d} I_{\mathrm{tot}}}{\mathrm{d} t} &= \frac{\beta S I}{N}\\ \frac{\mathrm{d} M}{\mathrm{d} t} &= \frac{c \beta S I}{N} \end{aligned}\end{split}\]

In framing the problem in this way we are making two approximations. Firstly, we are allowing the population numbers to be treated as continuous variables. Having non integer numbers of people is clearly not realistic, but can be acceptable when considering large populations. Secondly, in taking the net effect of the underlying processes, we confine ourselves to study only the deterministic evolution of the system.

Continuous Time Markov Chains

In a Continuous Time Markov Chain (CTMC), we again have the vector of variables, \(\mathbf{y}(t)\), describing the state of the system. There are \(N_E\) different types of event and the effect of an event occurring is to change the state by some amount. Mathematically speaking, if event \(j\) occurs then the variable \(y_i\) changes by \(\mathcal{D}_{ij}\), which, as with ODEs, may depend on the current state, time and parameters:

\[\begin{aligned} \Delta y_i \text{ due to event of transition type, j } = \mathcal{D}_{ij}(\mathbf{y}, \boldsymbol{\theta}, t) \end{aligned}\]

where \(\mathcal{D}\) is a \(n \times N_E\) matrix. These events occur randomly in time and independently of each other at rates (also possibly dependent on the state of the system):

\[\begin{aligned} \boldsymbol{\lambda}(\mathbf{y}, \boldsymbol{\theta}, t) = \left(\lambda_{1}(\mathbf{y}, \boldsymbol{\theta}, t),\ldots, \lambda_{N_E}(\mathbf{y}, \boldsymbol{\theta}, t)\right) \end{aligned}\]

For our case:

\[\begin{split}\begin{aligned} \mathbf{y}= \begin{pmatrix} S\\ I\\ R\\ I_{\mathrm{tot}}\\ M \end{pmatrix} \text{, } \boldsymbol{\lambda} = \begin{pmatrix} \frac{\beta S I}{N}, \gamma I, \mu N, \mu S, \mu I, \mu R \end{pmatrix} \text{, } \mathcal{D} = \begin{pmatrix} -1 & 0 & 1 & -1 & 0 & 0\\ 1 & -1 & 0 & 0 & -1 & 0\\ 0 & 1 & 0 & 0 & 0 & -1\\ 1 & 0 & 0 & 0 & 0 & 0\\ c & 0 & 0 & 0 & 0 & 0 \end{pmatrix} \end{aligned}\end{split}\]

where the column and row vector formatting of \(\mathbf{y}\) and \(\boldsymbol{\lambda}\) respectively is to better illustrate their context with the matrix, \(\mathcal{D}\). How we solve a mathematical problem framed in this way will be covered later on, for now we note that it is possible to derive exact, stochastic solutions, unlike in the ODE case.

Warning

Not all ODEs are CTMCs

It is worth including a warning at this point. We’ve seen that our SIR example model is amenable to both an ODE and a CTMC approach, however, this is not the case for all systems. Borrowing an example from physics, Hooke’s Law describes how the position, \(x\), and velocity, \(v\), of a mass, \(m\), attached to a spring (with spring constant, \(k\)) evolves with time. The set of ODEs are as follows:

\[\begin{split}\begin{aligned} \frac{\mathrm{d}x}{\mathrm{d}t}&=v\\ \frac{\mathrm{d}v}{\mathrm{d}t}&= -\frac{k}{m}x \end{aligned}\end{split}\]

Whilst it is possible to introduce stochasticity to this system, perhaps due to random fluctuations in the spring properties, this will not take the form of a CTMC. There is actually nothing to mathematically forbid one from modelling the system as a CTMC with \(\boldsymbol{\lambda}= \begin{pmatrix} 1\\ 1 \end{pmatrix}\), \(\mathcal{D}=\begin{pmatrix} v & 0 \\ 0 & -\frac{k}{m} x \end{pmatrix}\), but the resulting dynamics will not have a basis in reality since the underlying system is not a compartmental model.

It is therefore useful to remember that whilst PyGOM can be, in principle, used for any system of ODEs, it should only be used to provide stochastic results for true compartmental models.