ABSTRACT

From the mathematical point of view, the modeling of epidemics—in other words, the spread of infectious diseases transmitted from individual to individual—is very similar to the modeling of the magnetic systems studied by statistical physics. In this work, we use this analogy between mathematical epidemiology and statistical physics to study the classical mathematical model of epidemiology SI (Susceptible-Infected) approached through the Ising-Glauber model, in which individuals would be represented by atoms with spins −1 (susceptible) and 1 (infected). A Monte Carlo computational simulation was also performed for the Ising-Glauber model in a square network, where each network point represents an individual and the down and up spins represent susceptible and infected individuals.

Keywords:

Ising-Glauber Model, Epidemiology, Computer Simulation

1. Introduction

Mathematical models and computational simulations of infectious diseases are widely used tools in epidemiology. They make possible the understanding of the basic mechanisms whereby epidemic outbreaks can settle in a given population. In this way, it is possible to develop interventions that control or prevent the epidemic [1] . The most used mathematical models in epidemiology to describe the spread of infectious diseases divide the population into compartments related to the state in which individuals are found in the development of the disease [2] .

One of the best-known models is the SI (Susceptible-Infected) model, which divides the population into two compartments: susceptible and infected [3] . However, this compartmental model can be studied in a deterministic and stochastic way [4] . The deterministic models present a better outcome for large populations and are usually translated into systems of ordinary nonlinear differential equations whose resolution requires complex analytical and numerical methods [5] . On the other hand, because of their conceptual simplicity, they are more useful for identifying the most important variables and parameters in each process.

Stochastic models have a general application and emphasize the randomic aspect of the occurrence of a phenomenon. They involve more details about the properties of the system and can be elaborated based on the characteristics of the individuals and their structure of contacts [5] . However, the classic stochastic models of epidemiology are very complex and difficult to address mathematically.

In the field of Material Physics, one of the models used to study the transition properties of the ferro-paramagnetic phase of ferromagnetic materials is the Ising model. This model can also be used in the area of epidemiology to study the properties that are responsible for the spread of a disease since its mathematical structure resembles those used in epidemiological models. In this way, the present article aims to make an analogy between the SI model and the Ising model according to Glauber’s dynamics that has been studied extensively by physicists for more than sixty years and makes it easier to analyze similar epidemiological models.

The methodology for conducting the comparative study of the SI epidemiological model and the Ising-Glauber physical model took the following path: 1) initially, a deterministic and stochastic study of the SI model was performed. 2) A study of the Ising model was made in order to show that the phase transition occurs only in one dimension. 3) A study of the model proposed by Glauber for the Ising model was carried out. 4) An analytical comparison of the Ising-Glauber model was made for a single spin and also in the absence of a magnetic field, with the SI model (approached in a deterministic and stochastic way). 5) Finally, an analytical comparison and a computational simulation, via Monte Carlo, with N spin and in the absence of the magnetic field were carried out, with the SI model (approached deterministically and stochastically).

2. Epidemiological Model SI

One of the simplest models in epidemiology separates individuals into two states: susceptible (S) and infected (I). The SI model consists in analyzing a population of N susceptible people (S) in which an infected individual (I) is inserted within this population, so when the susceptible individual is infected, it remains in that state (infected) forever [5] .

2.1. Deterministic SI Model

The SI model is described by two non-linear ordinary differential equations.

$\frac{\text{d}S}{\text{d}t}=-\beta S\left(t\right)I\left(t\right)$ (1)

$\frac{\text{d}I}{\text{d}t}=\beta S\left(t\right)I\left(t\right)$ (2)

In these equations, $\beta$ is the rate of infective contacts that is a proportionality constant in density-dependent or frequency-dependent transmission [6] . Solving the differential equations above, by the method of separation of variables and considering the initial conditions, $S_0=N$ e $I_0=1$ is obtained.

$S\left(t\right)=\frac{N\left(1+N\right)}{N+{\text{e}}^{\beta \left(N+1\right)t}}$ (3)

$I\left(t\right)=\frac{\left(N+1\right){\text{e}}^{\beta \left(N+1\right)t}}{N+{\text{e}}^{\beta \left(N+1\right)t}}$ (4)

It can be noted that $S\left(t\right)$ is a continuous and decrescent function in which when time tends to infinity ( $t\to \infty$ ) the number of susceptible tends to zero (Equation (3)). The speed of the decrease depends solely on the value of β. On the other hand, $I\left(t\right)$ is a continuous and crescent function in which when time tends to infinity ( $t\to \infty$ ) all individuals have ended up infected ( $I\left(t\right)\to N+1$ ).

2.2. Stochastic SI Model

The deterministic SI model is, from the mathematical point of view, relatively simple, because it only has the rate of producing infected. However, if stochastically analyzed, it becomes somewhat more complex. In a stochastic SI model with infection force λ = constant, the possible states and the transition probabilities between them are represented in the following Table 1.

Thus, the probability of having x susceptible individuals and y infected individuals—as a function of the number of susceptible and infected individuals at time t and the transition probabilities—at time t + Δt is:

$P_{x;y}\left(t+\Delta t\right)=P_{x;y}\left(t\right)T{P}_{\left(x;y\to x;y\right)}+P_{x+1;y-1}\left(t\right)T{P}_{\left(x+1;y-1\to x;y\right)}$ (5)

Substituting the transition probabilities $T{P}_{\left(x;y\to x;y\right)}$ and $T{P}_{\left(x+1;y-1\to x;y\right)}$, for those indicated in the table below, we have:

$\frac{P_{x;y}\left(t+\Delta t\right)-P_{x;y}\left(t\right)}{\Delta t}=-\lambda x{P}_{x;y}\left(t\right)+\lambda \left(x+1\right)P_{x+1;y-1}\left(t\right)$ (6)

And, taking the limit $\Delta t\to 0$, we obtain the master equation for the SI model with constant λ.

$\frac{\text{d}{P}_{x;y}}{\text{d}t}=-\lambda x{P}_{x;y}+\lambda \left(x+1\right)P_{x+1;y-1}$ (7)

In some cases, the master equation can be solved analytically in a relatively simple way by the method of the Generating Function. A generating function is a

function that “generates” a probability distribution and can always be written as a summation [7] .

$G\left(u,v,t\right)={\displaystyle \underset{x,y}{\sum }{u}^x{v}^y{P}_x\left(t\right)P_y\left(t\right)}$ (8)

Thus, the expected mean value of the susceptible $\langle x\rangle$ is obtained from the generative function.

$\langle x\rangle ={\displaystyle \underset{x}{\sum }x\cdot P_x}$ (9)

so that

$\langle x\rangle =N{\text{e}}^{-\lambda t}$ (10)

It can be noted that the function found above has the same behavior as the exponential distribution in the radioactive decay or process of a population age

evolution. Thus, $\lambda$ is the decay rate (infection) and, naturally, $\frac{1}{\lambda }$ is the mean

decay time or mean time to infection, as expected.

3. Ising Model

Ising’s model was introduced by Wilhelm Lenz in 1920 to his doctoral student Ernest Ising. The main objective of this model is to study the ferromagnetism in a one-dimensional structure, and it is defined by the following Hamiltonian:

$\mathcal{H}=-J{\displaystyle \underset{\left(ij\right)}{\overset{N}{\sum }}{\sigma }_i{\sigma }_j}-H{\displaystyle \underset{i}{\overset{N}{\sum }}{\sigma }_i}$ (11)

The Hamiltonian ( $\mathcal{H}$ ) of the Ising Model relates the exchange interaction J to the sum over the nearest neighbor pairs ( $s_i{s}_{i+1}$ ). In this way, the spin variable $s_i$ that takes +1 or −1 values corresponding to the states up or down for each site $i=1,2,\cdots ,N$ of a d-dimensional structure. Everything in nature or systems of nature have a tendency to achieve the lowest energy state, so this model produces an ordered ferromagnetic state when the term J is positive, in which the parallel alignment of spins corresponds to an energy −J and the antiparallel alignment of spins (spins in opposite directions) corresponds to an energy +J. On the other hand, the term H is related to the external magnetic field in the system that is associated to each of the spins ( $s_i$ ).

In 1925 [8] , this model was solved exactly in one dimension by Ising and has reached the following the expression of magnetization.

$m\left(T,H\right)=\frac{\text{senh}\left(\frac{H}{k_{B}T}\right)}{{\left[{\text{senh}}^2\left(\frac{H}{k_{B}T}+\exp \left(\frac{-4J}{k_{B}T}\right)\right)\right]}^{1/2}}$ (12)

where $k_B$ is the Boltzmann constant and T is the temperature. Notice that the magnetization is zero when H = 0. When T > 0, the magnetization also tends to zero in which there is no spontaneous magnetization in one dimension, so the model does not explain ferromagnetism. These results led Ernest Ising to erroneously conclude that the same would happen in two and three dimensions.

Later in 1944, Lars Onsager obtained the exact expression of the partition function and of the free energy of Helmholtz in two dimensions for a square network at zero external field [9] . The result shows phase transition between the ferromagnetic and paramagnetic phases at a nonzero temperature. Since then, many derivations have been obtained for the Onsager solution in other arrangements and through numerical techniques, such as the Metropolis algorithm for Monte Carlo simulation.

4. Glauber Model

The Glauber model is a dynamic for the Ising model because it describes an Ising system in contact with a heat reservoir at a given temperature. The variables of spins ( ${\sigma }_i=\pm 1$ ) associated to each site i are random time variables also known as variables stochastic [10] . This model is governed by a master equation that in equilibrium reduces to the Ising model by the transition rate (transition probabilities per unit of time), which obeys the condition of detailed balance. This system consists of N sites, so the total number of configurations or states of the system is calculated by 2^N. A possible system configuration will be a vector $\sigma =\left({\sigma }_1,{\sigma }_2,{\sigma }_3,\cdots ,{\sigma }_{N-1},{\sigma }_N\right)$ whose i-th component is the variable associated with site i [7] .

The time evolution of the average of the system is governed by the master equation.

$\frac{\text{d}}{\text{d}t}P\left(\sigma ,t\right)=\underset{{\sigma }^{\prime }\left(\ne \sigma \right)}{{\displaystyle \sum }}\left\{W\left(\sigma ,{\sigma }^{\prime }\right)P\left({\sigma }^{\prime },t\right)-W\left({\sigma }^{\prime },\sigma \right)P\left(\sigma ,t\right)\right\}$ (13)

where $P\left(\sigma ,t\right)$ is the probability of occurrence of configuration $\sigma$ at time t and $W\left(\sigma ,{\sigma }^{\prime }\right)$ is the transition rate from ${\sigma }^{\prime }$ to $\sigma$. The transitions occur between configurations that differ only by the state of one site, so the transition rate is given by:

$W\left(\sigma ,{\sigma }^{\prime }\right)={\displaystyle \underset{i}{\sum }\delta \left({{\sigma }^{\prime }}_1,{\sigma }_1\right)\delta \left({{\sigma }^{\prime }}_2,{\sigma }_2\right)\cdots \delta \left({{\sigma }^{\prime }}_i,{\sigma }_i\right)\cdots \delta \left({{\sigma }^{\prime }}_N,{\sigma }_N\right)w_i\left(\sigma \right)}$ (14)

where

$\delta \left({{\sigma }^{\prime }}_j,{\sigma }_j\right)=\{\begin{array}{l}0se{{\sigma }^{\prime }}_j\ne {\sigma }_j\\ 1se{{\sigma }^{\prime }}_j={\sigma }_j\end{array}$ (15)

And $w_i\left(\sigma \right)$ is the rate of inversion of the state of the i-th site, from ${\sigma }_i$ to $-{\sigma }_i$. In this way, the master equation is given by:

$\frac{\text{d}}{\text{d}t}P\left(\sigma ,t\right)=\underset{i=1}{\overset{N}{{\displaystyle \sum }}}\left\{w_i\left({\sigma }^i\right)P\left({\sigma }^i,t\right)-w_i\left(\sigma \right)P\left(\sigma ,t\right)\right\}$ (16)

where the configuration $\sigma =\left({\sigma }_1,{\sigma }_2,\cdots ,{\sigma }_i,\cdots ,{\sigma }_N\right)$ is different from the configuration ${\sigma }_i=\left({\sigma }_1,{\sigma }_2,\cdots ,-{\sigma }_i,\cdots ,{\sigma }_N\right)$ only by the exchange of ${\sigma }_i$ to $-{\sigma }_i$. Thus, the temporal evolution of the average of the variable $\langle {\sigma }_i\rangle$ is

$\frac{\text{d}}{\text{d}t}\langle {\sigma }_i\rangle =-2\langle {\sigma }_i{w}_i\left(\sigma \right)\rangle$ (17)

The correspondence of the Glauber model with the Ising model is given through the choosing of transition rate, in a way that it depends on the values of the first neighboring spin. The rate chosen by Glauber [11] is given by

$w_i\left({\sigma }_i\right)=\frac{1}{2}\alpha \left[1-\frac{1}{2}\gamma {\sigma }_i\left({\sigma }_{i-1}+{\sigma }_{i+1}\right)\right]$ (18)

It is noticed that if γ is positive it has a tendency to a parallel configuration (that is the case of the ferromagnetic state) and if γ is negative it has a tendency to antiparallel configurations (that is the case of the paramagnetic state). On the other hand, the parameter α describes the time scale on which all transitions take place. In addition, the transition rate at equilibrium ( $t\to \infty$ ) has the same properties of the Ising model.

The partition function is given by:

$Z=\underset{\left\{\sigma \right\}}{{\displaystyle \sum }}\exp \left[-\frac{\mathcal{H}}{k_{B}T}\right]$ (19)

where {σ} denotes that the sum must go through all possible system configurations. The stationary probability is given by:

$P\left(\sigma \right)=\frac{1}{Z}\exp \left\{K{\displaystyle \underset{\left(ij\right)}{\overset{N}{\sum }}{\sigma }_i{\sigma }_j}\right\}$ (20)

where the sum extends over the pairs of sites that are nearest neighbors and

$K=\frac{J}{k_{B}T}$.

The transition rate satisfies the detailed balancing,

$\frac{P\left(-{\sigma }_i\right)}{P\left({\sigma }_i\right)}=\frac{w_i\left({\sigma }_i\right)}{w_i\left(-{\sigma }_i\right)}$ (21)

Substituting in the equation above the Glauber’s transition rate and stationary probability is. We found,

$\gamma =\text{tanh}\left(2K\right)$ (22)

Rewriting the transition rate for the d-dimensional case, we have:

$w_i\left({\sigma }_i\right)=\frac{1}{2}\alpha \left[1-{\sigma }_i\tanh \left(K{\displaystyle \underset{\delta }{\sum }{\sigma }_{i+\delta }}\right)\right]$ (23)

where the sum in $\delta$ extends over the first neighbors of ${\sigma }_i$. Glauber’s model has only one exact solution in one dimension, so

$\tanh \left[K\left({\sigma }_{i-1}+{\sigma }_{i+1}\right)\right]=\frac{\left({\sigma }_{i-1}+{\sigma }_{i+1}\right)}{2}\tanh \left(2K\right)$ is only true in the one-dimensional case [12] .

5. Results

The Ising-Glauber model can also be seen as a model for an epidemic spread, so the following is a comparison between the dynamic Ising model and epidemiological models in literature.

5.1. Analogy between the Glauber Model for a Single Spin in the Absence of External Magnetic Field and the SI Model

If we consider an “energy of infection” by analyzing the structure of compartments of the SI epidemiological model with the Hamiltonian elements of the Ising model (in Equation (11)), it is possible to draw an initial comparison. In this way, we can see that the Ising ferromagnetic model with external field equal to zero (H = 0) can be used to model a constant population where each individual would be represented by a spin that assumes the values of down (−1) and up (+1). The susceptible individual would be represented by a spin down and the infected individual would be represented by a spin up. The term interaction between the magnetic dipoles J would be related to the force of infection of the epidemiological model SI.

The physical magnitude present in Ising’s model that would be interesting to study and calculate is magnetization, since it calculates the sum of all spins configurations. That is, it measures how much the spins are aligned, which would be equivalent to calculating in the epidemiological model SI the number of susceptible or infected individuals. If the magnetization is negative, it would represent that the majority of the population is susceptible. If the magnetization is positive, it would represent that the majority of the population is infected and if the magnetization is zero, it would represent that half of the population is infected and the other half is susceptible.

In order to make the analogy between the physical model and the epidemiological model, the Ising model with the Glauber dynamics and the SI stochastic epidemiological model were specifically used. Firstly, the Ising-Glauber model was compared in the absence of magnetic field for a single spin ½ (spin that can assume the value of +1 or −1) with the stochastic model SI.

Initially, Glauber [11] considers in his work a system constituted solely by a single particle in the absence of an external magnetic field. This system is in contact with a heat reservoir, causing its spin to flip between the values +1 and −1 randomly. The transition rate per unit of time at which the particle makes transitions from either state to the opposite one is given by α/2 and let $P\left(\sigma ,t\right)$ be the probability of finding the spin in the state of $\sigma =\pm 1$ at time t, then we can write the master equation.

$\left(\text{d}/\text{d}t\right)P\left(\sigma ,t\right)=-\frac{1}{2}\alpha P\left(\sigma ,t\right)+\frac{1}{2}\alpha P\left(-\sigma ,t\right)$ (24)

In equilibrium when time tends to infinity the probabilities for the two states are equal. Taking into account the normalization condition, we have:

$P\left(1,t\right)+P\left(-1,t\right)=1$ (25)

Thus, we could write the expected value of the spin as a function of time,

$m\left(t\right)=P\left(1,t\right)-P\left(-1,t\right)=\underset{\sigma =\pm 1}{{\displaystyle \sum }}\sigma P\left(\sigma ,t\right)$ (26)

Therefore, the differential equation is obtained,

$\frac{\text{d}}{\text{d}t}m\left(t\right)=-\alpha m\left(t\right)$ (27)

Whose solution is,

$m\left(t\right)=m\left(0\right){\text{e}}^{-\alpha t}$ (28)

One can also write the individual probability,

$P\left(\sigma ,t\right)=\frac{1}{2}\left[1+\sigma m\left(t\right)\right]$ (29)

It can be seen that the average spin decays exponentially with a relaxation time 1/α towards the equilibrium value (null) [7] [13] .

The result of the average magnetization of equation 29 is similar to the result of the mean of susceptible of equation 10 of the stochastic model SI. The two equations decay exponentially with time and the mathematics of the two equations has the same behavior, which lead us to think in a similarity between the parameters $\alpha$ and $\lambda$, since both act as decay constants. The higher the $\alpha$ and $\lambda$, the faster the magnetization decreases and the susceptible ones respectively.

It is important to emphasize that the deterministic SI model considers interactions between individuals through the contact rate $\beta$ that is linked to the force of infection $\lambda$, whereas the Glauber model for a single spin that disregards interactions. However, these models are equivalent because even Glauber’s model for a single spin, not considering interactions, is in contact with a heat reservoir, causing the spin to change randomly between the values of +1 and −1 which would correspond to the susceptible and infected states.

5.2. Analogy between the Glauber Model for N Spin in the Absence of External Magnetic Field and the SI Model

A closed population can be represented by a one-dimensional chain of N individuals or atoms (Figure 1), with the boundary conditions ( $m_k\left(t\right)=m_{k+N}\left(t\right)$ ), taking into account that each atom has its spin and that each individual or spin can assume values of +1 and −1.

This one-dimensional model proposed by Glauber [11] can be solved in an exact way, that is, all correlations can be obtained in closed form. Substituting

the transition rate given by Equation (19) in the time evolution of the average given by Equation (17), we have

$\frac{\text{d}}{\text{d}t}{m}_k\left(t\right)=-\alpha m_k\left(t\right)+\frac{\alpha }{2}\tanh \left(2\beta J\right)\left[m_{k-1}\left(t\right)+m_{k+1}\left(t\right)\right]$ (30)

onde β = 1/K_BT.

Using the Fourier transform (FT) to solve the above equation, we write:

${\hat{m}}_k\left(t\right)=\frac{1}{N}{\displaystyle \underset{q}{\sum }{\hat{m}}_q\left(t\right)\exp \left(iqk\right)}$ (31)

where q belongs to the first zone of Brillouin

$q=0,\pm \frac{2\text{π}}{N},2\frac{2\text{π}}{N},\cdots ,\text{π}$ (32)

for N pair. Since magnetizations are real, you should have:

${\hat{m}}_q\left(t\right)={\hat{m}}_q^{\ast }\left(t\right)$ (33)

Thus, returning Equation (31)

$\frac{\text{d}}{\text{d}t}{\hat{m}}_q\left(t\right)=-\alpha {\hat{m}}_q\left(t\right)+\left(\alpha \gamma \cos q\right){\hat{m}}_q\left(t\right)$ (34)

where $\gamma =\text{tanh}\left(2\beta J\right)$.

Therefore

${\hat{m}}_q\left(t\right)=-A_q\exp \left\{-\alpha \left[1-\gamma \cos q\right]t\right\}$ (35)

where $A_q$ is a constant and $A_q=A_{-q}^{\ast }$.

Inserting this in Equation (32), we obtain the result

$m_k\left(t\right)=\frac{1}{N}{\displaystyle \underset{q}{\sum }{A}_q\exp \left\{-\alpha \left[1-\gamma \cos q\right]t+iqk\right\}}$ (36)

In particular

$m_0\left(t\right)=\exp \left\{-\alpha \left[1-\gamma \cos q\right]t\right\}$ (37)

Decreasing to the equilibrium null with the time constant ${\tau }^{-1}=\alpha \left(1-\gamma \right)$. The magnetization per site will be given by:

$m\left(t\right)=\frac{1}{N}{\displaystyle \underset{k=1}{\overset{N}{\sum }}{m}_k\left(t\right)}=\frac{1}{N}\exp \left\{-\alpha \left(1-\gamma \right)t\right\}$ (38)

It is noticed that the magnetization vanishes when time tends to infinity, as long as γ ≠ 1 (T ≠ 0), and that the relaxation time ( $\tau =\alpha {\left(1-\gamma \right)}^{-1}$ ) diverges when γ→1 (T→0). This is, when the temperature is zero, the magnetization is nonzero, which corresponds to the result obtained by Ising for its one-dimensional model that only presents phase transition (disordered phase to ordered phase) at zero temperature [12] . It is also noticed that the constant γ is linked to the exchange interaction term J and the temperature as shown in Equation (23), so when there is interaction between spins, the constant γ can be analogous to the force of infection λ of the epidemiological model SI.

These results indicate that when time tends to infinity, magnetization tends to decrease to zero, so 50% of the population will be infected. The temperature concept of the Ising model, according to Crisostomo et al. [14] , may be related in epidemiology to the level of aggression of a virus or associated to external parameters such as the ambient temperature and humidity. In addition, it can be associated to the effects of cultural and socio-economic risk factors.

5.3. Monte Carlo Simulation of the 2D Glauber-Ising Model

In this simulation, R codes for the Glauber-Ising model were developed from existing codes of the 2D static Ising model [15] . Then, the Monte Carlo method of the Glauber-Ising model will be presented in a square network with N = L × L sites. For each site of the network, a variable σ_i is associated, which assumes the values of +1 and −1. Periodic boundary conditions are also used to eliminate edge effects, so all spins will have the same number of neighbors.

Settings are generated according to the probability below.

$p=\frac{1}{2}\left[1-{\sigma }_i\tanh \left(\frac{J}{k_{B}T}{\displaystyle \underset{\delta }{\sum }{\sigma }_{i+\delta }}\right)\right]$ (39)

As the external magnetic field is zero, the transition rate has an inversion symmetry, $w_i\left(-{\sigma }_i\right)=w_i\left({\sigma }_i\right)$ [7] . To configure an algorithm for numerical simulation of the Glauber-Ising model, the following steps are followed. Initially, an arbitrary initial $\left({\sigma }_i={\sigma }_1,\cdots ,{\sigma }_i,\cdots ,{\sigma }_N\right)$ configuration is generated. After an i-th site of this network is chosen randomly, then the probability of inversion given by Equation (40) is calculated and a random number $\xi$ is generated in the interval of {0,1}. If $\xi \le p$, then the variable ${\sigma }_i$ is inverted $\left({\sigma }_i\to -{\sigma }_i\right)$, so a new system is generated $\left({\sigma }_i={\sigma }_1,\cdots ,-{\sigma }_i,\cdots ,{\sigma }_N\right)$. But, if $\xi >p$, then the variable ${\sigma }_i$ is not inverted and the system continues in the same state $\left({\sigma }_i={\sigma }_1,\cdots ,{\sigma }_i,\cdots ,{\sigma }_N\right)$ [16] .

A Monte Carlo step is completed when this process, described above, is done N (number of network sites) times. In this simulation, the computationally intensive part of this process comes from having to repeat these steps enough times that equilibrium should have been achieved. In this way, after several repetitions, the statistical sampling is obtained. A small network N = 12 × 12 has been chosen so that the equilibrium requires about 10⁶ steps Monte Carlo [15] .

The average magnetization (Figure 2) is calculated directly by the sum of all the spins over the entire network.

$M=\frac{1}{N}\left({\displaystyle \underset{i=1}{\overset{N}{\sum }}{\sigma }_i}\right)$ (40)

Note that the magnetization is finite for any value of temperature (T) and decreases sharply near the critical temperature. That is, in the limit where $\to \infty$, the magnetization ( $m\to 0$ ) is zero for $T\ge T_c$. The mean of spins up and down were also calculated (Figure 3 and Figure 4).

Observing the average behavior of the spins, it is noticed that above the critical temperature the average of the spins down decreases and the average of the spins up increases to a certain value. Analyzed these results from an epidemiological point of view, in which the spins down represents the susceptible and the spins up represents the infected, it is perceived that the disease was able to propagate when it reached the critical temperature.

6. Conclusions

The deterministic SI and SIS mathematical models seek, in a simplified way, to generate good predictions, aiding the understanding of the epidemics and the study of control mechanisms. Nonetheless, these models have limitations, since they assume either that all individuals in a population are intermingling freely with one another. In other words, it assumes that individuals are distributed approximately evenly across the population, so that all susceptibles are subject to roughly the same probability of a contact being with an infectious individual. In reality, however, transmission typically occurs locally, between nearby individuals [17] . The Ising model, however, relates the interaction factor or the contact factor only among the first neighbors.

As previously seen, the Ising model was described by Glauber’s dynamics for a single spin and, in the absence of a magnetic field, it has the same behavior as the stochastic SI model. When analyzed in one dimension, this model has the behavior of Equation (39), where the temperature is an important parameter and the system only undergoes a phase transition when the temperature is zero. On the other hand, when analyzed in two or more dimensions, the model presents a phase transition when it reaches the critical temperature, in which the number of infected individuals corresponds to half of the population.

The formalism of the Ising model is completely different from that of the epidemic models found in the literature. To make an analogy between the variables and the relations established between them is not trivial, even in the simplest case of magnetization (m). Although the Ising model considers only interactions between close individuals that make it more like reality, where individuals are considered static, which does not occur in populations of animals or people. Nevertheless, the study of epidemic systems through the Ising model shows another point of view that perhaps brings to light properties not described in the epidemiological models.

For example, Anderson [18] discusses a discrepancy between the contact rate and the basic number of reproduction of the infection. The author shows a case in which an individual interviewed contracted a particular sexually transmitted disease and reported having had 1.46 partners for a given time interval. It has been observed that this value is very low, since in the endemic equilibrium (R₀ = 1)¹ the theory suggests that, on average, each infected person must have two effective contacts during the course of an infection, i.e., each individual must be infected by another infected individual and infect another susceptible individual.

This example described by Andreson [18] can be linked to the fact that the one-dimensional Ising model does not undergo phase transition and its magnetization is characterized only by the random formation of agglomerates. This shows that some infections will not evolve in the population if each individual has no more than two partners to transmit it, even though the basic number of infection in this population, R₀, is greater than 1. This hypothesis may be very interesting for the study of some sexually transmitted diseases. Another point is that the Ising model in two or more dimensions presents a phase transition for a given set of values of its parameters. It is known that epidemiological models also undergo a kind of phase transition when R₀ = 1.

Finally, this work sought to investigate and model epidemic systems, SI, using the stochastic formalism of statistical mechanics, based on the study of the Ising-Glauber model. However, there are many other research possibilities between the epidemic compartmental models and the Ising model. For example, a comparison of the SIS (Suspected-Infected-Susceptible) model with Ising’s model described by Glauber’s dynamics with external magnetic field can be studied. We can also study a comparison of the SIR model with the Ising model in the absence of magnetic field and with spin variables (−1.0 and +1)—where the spins down configuration would represent the susceptible individuals, the configuration of spins up would represent the infected individuals and the configuration of spins zero would represent the individuals recovered. Another possibility is to represent the SIRS model using the Ising model with spin variables (−1, 0 and +1) in a magnetic field oriented to the down configuration. In this sense, it is evident that the Ising model allows many possibilities of rearrangements of its terms, in a way that allows us to create analogs to the epidemiological models found in the literature.

Funding

This work was supported by the foundation CAPES, entity of the Brazilian Government focused on the training of human resources.

Conflicts of Interest

The authors declare no conflicts of interest regarding the publication of this paper.

Cite this paper

Lara, C.G.A., Massad, E., Lopez, L.F. and Amaku, M. (2019) Analogy between the Formulation of Ising-Glauber Model and Si Epidemiological Model. Journal of Applied Mathematics and Physics, 7, 1052-1066. https://doi.org/10.4236/jamp.2019.75071

References

NOTES

¹R₀: measures the number of infections that result from the introduction of an infected individual into a fully susceptible population (McCallum; Barlow; Hone, 2001).