key: cord-0497737-pi3rmrnq
authors: Sakaguchi, Hidetsugu; Nakao, Yuta
title: Slow decay of infection in the inhomogeneous SIR model
date: 2020-12-24
journal: nan
DOI: nan
sha: e2bd6905d50212cc9ce55ccd6cd273aece62e38f
doc_id: 497737
cord_uid: pi3rmrnq

The SIR model with spatially inhomogeneous infection rate is studied with numerical simulations in one, two, and three dimensions, considering the case that the infection spreads inhomogeneously in densely populated regions or hot spots. We find that the total population of infection decays very slowly in the inhomogeneous systems in some cases, in contrast to the exponential decay of the infected population I(t) in the SIR model of the ordinary differential equation. The slow decay of the infected population suggests that the infection is locally maintained for long and it is difficult for the disease to disappear completely.

the first-order phase transition, where both the diffusion and coarsening processes are important [30, 31] . We think that the diffusion, coarsening, and quenched randomness are important for the slow decay in our model, however, the theoretical understanding is not sufficient yet and the details are left to future study.

We make a very brief review of the dynamics of the Kermack-McKendrick model as an ordinary differential equation [3, 4] . The Kermack-McKendrick model is expressed as

where S and I denote susceptible and infected populations. The parameters β and γ denote infection and recovering rates. The recovered population R is calculated from

Since the three variables S, I, and R are used, the Kermack-McKendrick model is called the SIR model. If the population E of exposed persons before the appearance of symptoms is included, the SIR model is generalized to the SEIR model:

where δ denotes the incidence rate. In the SIR model, recovered persons are assumed not to be infected again owing to the immunity. However, there are diseases such as Malaria, for which there is a possibility that recovered persons are infected again. For such diseases, the SIS model where dS/dt = −βSI + γI is used instead of Eq. (1). In the SIR and SEIR models, I decays to zero finally, but I does not decay to zero in the SIS model. The contact process is a stochastic version of SIS model. In this paper, we consider mainly the SIR model, however, slow dynamics is observed even in the SEIR model as shown in Sec. IV. There is a stationary solution S = S 0 and I = 0 to Eqs. (1) and (2) . The stationary state is unstable for S 0 β > γ. Figure 1 (a) shows trajectories in (S, I) space starting from the initial condition S(0) = 4, 2.5, and 1.3 at β = γ = 1. The initial condition for I is fixed to be 0.00001. Figure 1 (b) shows time evolutions of I(t) for S(0) = 4 and 2.5 at β = γ = 1 in the semi-logarithmic scale. The infection I(t) spreads initially and then decays to zero exponentially since the susceptible population S(t) becomes below γ/β. This is a state that the herd immunity is attained. The final state is (S, I) = (S ∞ , 0) where S ∞ < γ/β depends on the initial value S(0). The uninfected population S ∞ takes a smaller value for a larger initial value S(0) of S(t). S ∞ can be calculated for the conserved quantity Q of this equation:

If I(0) is sufficiently small, S ∞ is a solution of

Figure 1(c) shows S ∞ as a function of S(0) at γ = β = 1. The ratio S ∞ /S 0 takes any value between 0 and 1, depending on the initial value S 0 .

Hereafter, we consider spatially extended systems. The spread of infection occurs in densely populated regions. In this paper, we consider mainly SIR models with inhomogeneous infection rate on one-, two-, ang three-dimensional lattices. In one dimension, the model equation is written:

where D S and D I are diffusion constants. Periodic boundary conditions are imposed and the system size is N . The initial conditions were set to be S i = 1 and I i = 0.00001 for numerical simulations. Fig. 2 (a). The infection occurs locally near i = N/2, S i diffuses into the infection region, and is infected near i = N/2. Outside of the infection region, the profile of S i has a nearly constant slope, that is, Figure 4 (c) shows the relationship between β N/2 and S N/2 . The localized structure disappears for β N/2 < 1.5 and S i = S N/2 = S B = 1 is satisfied. Next, we make an analysis of the localized state and the power-law decay of exponent 1/2. If the continuum approximation is taken, Eqs. (5) and (6) are rewritten as

where ∇ 2 = ∂ 2 /∂x 2 in one dimension, however, the model will be later extended to two and three dimensions. For the stationary state in one dimension, S(x) and I(x) satisfy

If I(x) and β(x) are assumed to be

Equation (10) yields

at x = N/2, which leads to

Equations (11) and (13) yield

The dashed line in Fig. 4 (c) shows the relationship between β N/2 = β 0 + β 1 and S(N/2) by Eq. (14) . Fairly good agreement with direct numerical results is seen. If the infection region is sufficiently small, Eq. (9) gives

From Eq. (15), I(N/2) is approximated at

If the boundary conditions are not fixed to S B , the power-law decay is observed. Similar power-law decay is observed even for D I = 0. In the continuum approximation, the solution satisfies

for x = N/2, since I(x, t) is almost zero for x = N/2, and 2D S ∂S/∂x = β N/2 S(N/2)I(N/2) at x = (N/2) + . Because S(N/2) is fixed to be γ/β N/2 in case of D I = 0, I(N/2) is determined to be 2D S ∂S/∂x/γ. That is, the population of infection is determined by the diffusion process of susceptible population for large t. The slope ∂S/∂x at x = (N/2) + can be calculated from the solution of S(x, t) to the diffusion equation. The diffusion equation can be solved by the Fourier transform

I(N/2) is evaluated as

Since the summation is negligible for n satisfying 4D S π 2 n 2 t/N 2 >> 1 or n >> N/(4D S π 2 t) 1/2 , I(N/2) is approximated as

if the contribution (1 − S(N/2))/(N/2) in the second term of Eq. (18) is neglected for large N/2. This is a reason of the power law of exponent 1/2. The dashed line in Fig. 5 (a) denotes this relation, which is good approximation to the direct numerical simulation. Two-and three dimensional models are expressed with Eqs. (7) and (8) if ∇ 2 is rewritten as ∂ 2 /∂x 2 + ∂ 2 /∂y 2 in two dimensions and ∂ 2 /∂x 2 + ∂ 2 /∂y 2 + ∂ 2 /∂z 2 in three dimensions. If S, I and β depends only the radius r, Eqs. (7) and (8) are expressed as

where d denotes the dimension 2 or 3. We have performed numerical simulation of Eqs. (19) and (20) as a onedimensional discrete system similar to Eqs. (5) and (6) 

In this section, we study random SIR models in one, two and three dimensions. The one-dimensional model is again expressed as

where β i is a uniform random number between 0 and β m . Fig. 7(a) . Figure  7 (c) shows four snapshot profiles of S i at t = 200, 300, 2000, and 10000 in the same frame. The infection occurs at many hot spots at t = 200, and the spatial profile is intermittent as shown in Fig. 7(b) . The number of hot spots decrease with time. The cusp points in Fig, 7 (c) correspond the hot spots. S i decreases with time by the diffusion to the hot spots and infection at the hot spots. Figure 7 (d) shows the local average of β i (red dotted line), that is,

j=i−10 β i , 100I i (blue dashed line) at t = 1000 and the maximum value of I i (t) (green solid line) for 0 < t < 10000 in the range of 2000 < i < 3000. In most cases, hot spots appear near the points where the locally averaged infection rateβ i is large. Infection is stamped out at some hot spots, however, strong hot spots survive for long. The lifetime of a hot spot with an interval L is estimated as O(L 2 ), because the width of the diffusion field increases as t 1/2 as shown in Fig. 5 (b) and the power-law decay changes to an exponential decay when the width reaches the size of interval. If the intervals between strong hot spots become longer after the burnout of some hot spots, the lifetime of the survived hot spot becomes even longer. The power-law decay of t −1/2 by the diffusion effect and the increase of lifetime by the coarsening process might be the origin of the slow decay in one dimension.

We study the SEIR model to check the generality of the slow dynamics: The model equation is

where S i , E i and I i denote respectively the susceptible, exposed, and infected populations, and δ i denotes the incidence rate. The parameters are set to be D S = D E = 1, D I = 0.5, γ = 1, and δ = 2. Figure 8 The two-dimensional random SIR model is expressed as The dotted line is SI ∝ 1/t 0.85 . The system size is 600 × 600. The initial condition is S i,j = 1 and I i,j = 0.00001. The total population of infection seems to decay in a power law in these numerical simulations. Figure 9 (b) shows some snapshots of I i,j at a section of j = N/2 when β i,j takes a uniform random value between 0 and 3.2 at t = 50, 100, · · · , 500. Localized clusters of infection survive for long. The number of hot spots decreases with time, or a coarsening occurs also in two dimensions.

The three-dimensional model is expressed as

dI i,j,k dt = β i,j,k S i,j,k I i,j,k − γI i,j,k + D I (I i+1,j,k + I i−1,j,k + I i,j+1,k + I i,j−1,k + I i,j,k+1 + I i,j,k−1 − 6I i,j,k ), (24) to the case shown in Fig. 6(d) . I N/2+1,N/2+1,N/2 is almost constant after t > 50, that is, the localized infection is maintained for very long. Figure 10 

We have found slow decay of infection in the Kermack-McKendrick model with spatially inhomogeneous infection rate in some parameter range. First, we have studied the Kermack-McKendrick model with one hot spot where the infection rate is locally higher than the surrounding region. We have shown theoretically a power-law decay of exponent 1/t 1/2 in the one-dimensional system with a spatially localized hot shot. The slow decay occurs as 1/ log t in the two-dimensional system with one hot spot, and the decay seems to be even slower in three dimensions. Next, we have studied the random Kermack-McKendrick model, and found the power-law type slow decay in one, two, and three dimensions. We found that the infection occurs locally at hot spots. The uninfected persons in the surrounding area around the hot spots diffuse into the hot spots and are infected at the hot spots. The number of hot spots decreases in time and the lifetime of the survived hot spots become even longer because the surrounding areas of the survived hot spots become larger. The mechanism of the slow decay in our system is unique in that the diffusion, coarsening, and quenched randomness are important, although there is some similarity to the slows dynamics in the Griffiths phase in the contact process [18, 19] and the phase transition dynamics [30, 31] . However, the slower decay

our two-and three-dimensional SIR models with one hot spot, and the exponent of the power-law decay in one-, two-, and three-dimensional random SIR models are not well understood. They are left to future study. Our finding of the slow decay of the infected population suggests that the infection is locally maintained for long time and the diseases are hardly stamped out in some cases

Mathematical Epidemiology of Infectious Diseases: Model Building, Analysis, and Interpretation

Non-equilibrium phase transitions

Phase Transition Dynamics