- Email: [email protected]
Stochastic spatial structured model for vertically and horizontally transmitted infection
Physica A (
)
–
Contents lists available at ScienceDirect
Physica A journal homepage: www.elsevier.com/locate/physa
Stochastic spatial structured model for vertically and horizontally transmitted infection Ana T.C. Silva a , Vladimir R.V. Assis a , Suani T.R. Pinho b , Tânia Tomé c,∗ , Mário J. de Oliveira c a
Departamento de Física, Universidade Estadual de Feira de Santana, 44036-900 Feira de Santana, BA, Brazil
b
Instituto de Física, Universidade Federal da Bahia, 40210-340 Salvador, BA, Brazil Instituto de Física, Universidade de São Paulo, Rua do Matão, 1371, 05508-090 São Paulo, SP, Brazil
c
highlights • • • •
An epidemic process with horizontal and vertical transmission is analyzed. The epidemic stochastic lattice model is studied by Monte Carlo simulations. The epidemic stochastic lattice model is studied also by a pair mean-field approach. The phase diagram may display a healthy and infected coexisting state.
article
info
Article history: Received 6 July 2016 Received in revised form 1 October 2016 Available online xxxx Keywords: Stochastic lattice model Model for spreading of disease Model for vertical and horizontal transmitted infection
abstract We study a space structured stochastic model for vertical and horizontal transmitted infection. By means of simple and pair mean-field approximation as well as Monte Carlo simulations, we construct the phase diagram, which displays four states: healthy (H), infected (I), extinct (E), and coexistent (C). In state H only healthy hosts are present, whereas in state I only infected hosts are present. The state E is characterized by the extinction of the hosts whereas in state C there is a coexistence of infected and healthy hosts. In addition to the usual scenario with continuous transition between the I, C and H phases, we found a different scenario with the suppression of the C phase and a discontinuous phase transition between I and H phases. © 2016 Elsevier B.V. All rights reserved.
1. Introduction The transmission of disease between hosts of the same species occurs either vertically, when the pathogen is transmitted directly from the parents to the child, or horizontally, when the hosts are not in a parent–child relationship. Several diseases are transmitted in both horizontal and vertical ways [1–11]. A question that arises in the study of horizontally and vertically transmitted diseases is the determination of the amount of horizontal transmission needed to maintain a vertical transmission. This problem was analyzed by Lipsitch et al. [7] who proposed that in the stationary state, the prevalence of infection can be higher when there are both types of transmissions as compared to the case when only one type is present. This result implies that for the case of virulent pathogens, the high number of infected cases coming from vertical transmission is a consequence of a high rate of horizontal transmission.
∗
Corresponding author. E-mail address: [email protected] (T. Tomé).
http://dx.doi.org/10.1016/j.physa.2016.10.048 0378-4371/© 2016 Elsevier B.V. All rights reserved.
2
A.T.C. Silva et al. / Physica A (
)
–
Fig. 1. Illustration of the reactions of the stochastic lattice model. The states 0, 1 and 2 represent, respectively, an empty site, a site occupied by a healthy host and by an infected host. The quantities u1 , b1 , u2 , b2 , and β are the transition rates.
One approach in the study of population biology of infectious diseases consists in setting up differential equations for the densities of the various types of population [12], similar to equations used in chemical kinetics in which the law of mass action is employed. The model of Lipsitch et al. [7] is of this type and consists of two differential equations for the densities of infected and uninfected populations. It is possible, on the other hand, to set up stochastic population models with space structure to describe the spreading of diseases and, more generally, population dynamics [13–21]. In fact, a space structured stochastic model with the same reactions as the model used by Lipsitch et al. has been introduced by Schinazi [17]. For a particular case, he showed that if the vertical transmission is sufficiently high then infected hosts can be eliminated even if the horizontal transmission is high. Here, we analyze a space structured stochastic model which is a slight modification of the Schinazi’s model. The model is investigated by means of two mean-field approaches, the simple mean-field approximation and the pair mean-field approximation, and by numerical simulations on a square lattice. The evolution equations at the level of simple mean-field approximation are found to be equivalent to the differential equations used by Lipsitch et al. [7]. The analysis of the mean-field equations and the numerical simulations on a square lattice reveal that the model may display four states: an extinction (E) state, where the hosts have disappeared, a healthy (H) state, with healthy hosts and the absence of infected hosts, an infected (I) state, with infected hosts and the absence of healthy hosts, and a coexistence (C) state, where both infected and healthy hosts are present. For some set of values of the parameters, the H and I states are not contiguous and are separated by the C state. The transition lines HC and IC of the phase diagram correspond to continuous phase transitions. However, for other sets of values of the parameters, it is possible that the H and I states be contiguous. In this case the transition line HI corresponds to a discontinuous phase transition. The phase diagram obtained by the simple and pair mean-field approximation as well as the numerical simulations are qualitatively similar. However, there are important quantitative changes. For instance, we have found that the H, I and C regions of the phase diagram obtained by numerical simulations are smaller than their counterpart obtained by pair mean-field approximations which in turn are smaller than those obtained by simple mean-field approximation. 2. Model The stochastic lattice model is defined as follows. Each site of a regular lattice is either empty or occupied by a host that can be healthy or infected. A healthy host dies spontaneously with rate u1 and is born by a catalytic process with rate strength b1 . Analogously, an infected host dies spontaneously with rate u2 and is born by a catalytic process with rate strength b2 , a reaction related to vertical transmission. A healthy host becomes infected catalytically with rate strength β , a reaction understood as horizontal transmission. These reactions are illustrated in Fig. 1. Denoting the state of site i by the stochastic variable ηi , then ηi takes the values 0, 1, or 2, according to whether the site is empty, occupied by a healthy host or occupied by an infected host. The possible transitions of the state of a given site i with the respective rates are as follows Reaction Rate
0→1 b 1 ni
1→0 u1
0→2 b2 m i
2→0 u2
1→2 β mi
where ni and mi denote the number of healthy hosts and the number of infected hosts in the neighborhood of site i, respectively. The neighborhood of a site in a regular lattice is chosen to be its nearest neighbor sites. The model has five parameters: b1 , u1 , b2 , u2 , and β . Sometimes it will be convenient to use the parameters r1 and r2 defined by r1 =
u1 b1
,
and r2 =
u2 b2
.
(1)
When u1 = u2 , the present stochastic model reduces to a particular case of the stochastic model introduced by Schinazi [17]. The transition rates above define a stochastic process characterized by a master equation, which gives the evolution of the probability distribution P (η) associated to the microscopic state of the system η = {ηi }. From the master equation it is possible to write down the time evolution of the several marginal probability distributions such as P (ηi ), the probability
A.T.C. Silva et al. / Physica A (
)
–
3
of site i being in state ηi , and P (ηi , ηj ), the probability of neighboring sites i and j being in states ηi and ηj , respectively. Particularly, the time evolution equations for the one-site probabilities are d dt d
P (1) = b1 P (10) − u1 P (1) − β P (12),
P (2) = b2 P (20) − u2 P (2) + β P (12). dt Since P (0) + P (1) + P (2) = 1, the evolution equation for P (0) can be obtained from Eqs. (2) and (3).
(2) (3)
3. Simple mean-field approximation We solve Eqs. (2) and (3) by the use of a simple mean-field approach which amounts to employ the approximation P (ηi , ηj ) = P (ηi )P (ηj ) [21], which means to say that the correlations are neglected in this approximation. Using the notation x = P (1), which is the density of healthy hosts, y = P (2), which is the density of infected hosts, and z = P (0) = 1 − x − y, we may write P (10) = xz, P (20) = yz, and P (12) = xy, so that the set of Eqs. (2) and (3) reduces to the following equations dx dt dy
= b1 xz − u1 x − β xy,
(4)
= b2 yz − u2 y + β xy. (5) dt These equations are equivalent to the equations analyzed by Lipsitch et al. [7]. To see that this is indeed the case, it suffices to perform the rescaling x → x/K and y → y/K , where K is the carrying capacity. Here, the property that the host population is density limited is taken into account by considering that hosts proliferate only on empty sites of the lattice. This procedure has been employed for population biology systems described by stochastic lattice models [13]. The simple mean-field equations predict four types of stationary solutions, or states, whose stability can be found by analyzing the eigenvalues λ1 and λ2 of the Hessian
b 1 ( z − x ) − u1 − β y −(b2 − β)y
−(b1 + β)x . b2 (z − y) − u2 + β x
(6)
If the real part of the eigenvalues is negative, the corresponding state is stable. One solution is the trivial solution x = 0, y = 0, and z = 1, characterized by the extinction of hosts, which we call the extinction (E) state. This trivial solution is stable when u1 > b1 and u2 > b2 so that the lines of stability are given by r1 = 1,
(r2 ≥ 1), (r1 ≥ 1),
r2 = 1,
(7) (8)
which are, respectively, the vertical and horizontal lines shown in Fig. 2. When we cross the line defined by Eq. (7), there appears a state with the absence of infection, where all hosts are healthy, which we call healthy (H) state, characterized by x ̸= 0 and y = 0, that is, x = 1 − r1
and
y = 0.
(9)
From the stability condition, we found that this state is stable in the region limited by the line given by Eq. (7) and by the line given by b2 (r2 − r1 ) + β(r1 − 1) = 0,
(10)
as shown in Fig. 2. The line given by (10) is the upper inclined line in Fig. 2(a) and the lower inclined dashed line in Fig. 2(b). When we cross the line defined by Eq. (8), there emerges a state with the prevalence of infection, with no healthy hosts, characterized by x = 0 and y ̸= 0, that is, x=0
and
y = 1 − r2 .
(11)
From the stability condition, we found that this state is stable in the region limited by the line given by Eq. (8) and by the line given by b1 (r2 − r1 ) + β(r2 − 1) = 0,
(12)
as shown in Fig. 2. The line given by (12) is the lower inclined line in Fig. 2(a) and the upper inclined dashed line in Fig. 2(b). In addition to these three states, there is an additional state where the populations of healthy and infected hosts coexist, which we call coexistence (C) state, characterized by x ̸= 0 and y ̸= 0. From the stationary solution of Eqs. (4) and (5), the densities of healthy and infected hosts are x= y=
b1 b2 (r2 − r1 ) + β b2 (r2 − 1)
β(b1 − b2 + β) b1 b2 (r1 − r2 ) − β b1 (r1 − 1)
β(b1 − b2 + β)
,
(13)
.
(14)
4
A.T.C. Silva et al. / Physica A (
)
–
Fig. 2. Phase diagram in the variables r2 versus r1 obtained by the simple mean-field approximation. (a) b1 = 1, b2 = 2/3 and β = 1. (b) b1 = 1/3, b2 = 3 and β = 1. The states are: extinct (E), healthy (H), infected (I) and a state where there is a coexistence (C) of infected and healthy hosts. Solid lines represent continuous phase transitions. The upper and lower dashed lines are the spinodal lines, that is, the limits of stability of states I and H, respectively.
These solutions are the same as found in Ref. [7]. Taking the limits x → 0 and y → 0, we may find the boundaries of the state C in the phase diagram. We see that these two boundaries coincide with the lines given by Eqs. (12) and (10), respectively. We remark that, for some values of the parameters, the state C may not exist. If
β ≤ b2 − b1 ,
(15)
then the line of stability of state H, given by Eq. (10), will remain below the line of stability of state I, given by Eq. (12), as shown in Fig. 2(b), and state C does not exist anymore, given rise to a co-stability of states H and I. It is worth mentioning that the C state should not be confused with the co-stability of states H and I.1 In the y–x plane, the former corresponds to a single stable fixed point with x ̸= 0 and y ̸= 0, whereas the later corresponds to two stable fixed points, one with x = 0 and y ̸= 0 and the other with x ̸= 0 and y = 0. It is convenient to represent the possible states in a phase diagram in the variables h and v , defined by
v= h=
r1 r2
,
(16)
β(1 − r1 ) b2 r2
,
(17)
as introduced by Lipsitch et al. [7]. In these variables the two lines given by Eqs. (12) and (10) are h=
b1 r1 (1 − r1 )(1 − v) b2 r2 (v − r1 )
,
h = 1 − v,
(18) (19)
respectively. The simple mean-field results using these variables are shown in Fig. 3. It should be remarked that two types of phase diagrams are obtained. If u1 /u2 ≥ 1 than the line given by Eq. (18) is always above the line given by Eq. (19), as shown in Fig. 3(a). However, if u1 /u2 < 1 then these two lines may cross each other at a nonzero value of h and 1 − v , as shown in Fig. 3(b). In this case the two states I and H become contiguous and separated by a discontinuous phase transition up to a bicritical point. 4. Pair approximation Now we take into account that the states of neighboring sites might be correlated. To this end we derive from the master equation the time evolution equations for the marginal probability distributions P (ηi , ηj ) of two nearest neighbor sites i and j being in states ηi and ηj , respectively. Thus, in addition to Eqs. (2) and (3), we write down the evolution equations for these marginal probability distributions. Using a regular lattice, the evolution equations for the pair probabilities P (ηi , ηj ) are d dt
P (10) = u2 P (12) + u1 P (11) − u1 P (10) − b1 gP (10)
+ b1 fP (100) − b1 fP (101) − b2 fP (102) − β fP (210),
(20)
1 Using the thermodynamic language, the co-stability of state H and I corresponds to a two-phase coexistence whereas the coexisting C state corresponds to a single phase.
A.T.C. Silva et al. / Physica A (
)
–
5
Fig. 3. Phase diagram, obtained by the simple mean-field approximation, in the variables h = (β/b2 r2 )(1 − r1 ) and 1 − v where v = r1 /r2 , for the following values of the parameters: (a) b1 = 0.6, u1 = 0.1 and u2 = 0.08; (b) b1 = 0.6, u1 = 0.1 and u2 = 0.12. The states are: healthy (H), infected (I) and a state where there is a coexistence (C) of infected and healthy hosts. Solid lines represent continuous phase transitions. The upper and lower dashed lines are the spinodal lines, that is, the limits of stability of states H and I, respectively. The full circle is a bicritical point.
d dt
P (20) = u1 P (21) + u2 P (22) − u2 P (20) − b2 gP (20)
+ b2 fP (200) − b2 fP (202) − b1 fP (201) + β fP (210),
(21)
d
P (21) = −u1 P (21) − u2 P (21) + b2 fP (201) + b1 fP (201) + β fP (211) − β fP (212) − β gP (21), (22) dt where f = (ζ − 1)/ζ and g = 1 − f , with ζ being the coordination number of the lattice, that is, the number of neighboring sites of a given site of the regular lattice. In these equations, P (ηi , ηj , ηk ) denotes the three-site probability distribution related to the sites i, j, and k, where i and k are neighboring sites of a central site j, and i ̸= k. It should be noted that there are several types of three-site probabilities which should be distinguished by a different notation. However, we are not using distinct notations in view of the approximation we are about to use. In this approximation the several types of three-site probabilities become equal to each other. The set of Eqs. (20), (21), (22) cannot be solved by considering only the two-site correlations P (01), P (02), and P (12) because these equations involve correlations of greater order. A truncation scheme [21] will then be used, which amounts to employ the following approximation for the three-site correlations P (ηi , ηj , ηk ) =
P (ηi , ηj )P (ηj , ηk ) P (ηj )
.
(23)
This is the pair approximation: correlations of three or more sites are written in terms of one-site and two-site correlations Using the notation s = P (10), v = P (20), and w = P (12), so that P (11) = x − s −w , P (22) = y −v−w , and P (00) = z − s −v , Eqs. (2) and (3), and Eqs. (20), (21), and (22) are written as dx dt dy dt ds dt dv dt dw
= b1 s − u1 x − βw,
(24)
= b2 v − u2 y + βw,
(25)
= (u2 − u1 )w + u1 x − (2u1 + b1 g )s + b1 f
s(z − s − v)
= (u1 − u2 )w + u2 y − (2u2 + b2 g )v + b2 f = −(u1 + u2 + β g )w + (b1 + b2 )f
sv
z
− b1 f
v(z − s − v)
+ βf
z
s2
− b2 f
w(x − s − w)
− b2 f
z
v2 z
− βf
sv z
− b1 f
w2
− βf
sv z
ws x
+ βf
,
ws x
(26)
,
(27)
. (28) dt z x x This set of five equations makes up a closed set of equations for the variables x, y, s, v and w . We have solved these equations numerically and found that they predict the same four types of solutions, or states, found in Section 3 by the use of the simple mean-field approximation. The state E, characterized by x = 0, y = 0, s = 0, v = 0, and w = 0, is stable as long as r1 > f and r2 > f . When this condition is not met, other states emerge. When the state E becomes unstable, it gives rise to the state H or to state I, as can be seen in Fig. 4. The state C lies between these two states from which it is separated by continuous transition lines. If the ratio β/b2 ≪ 1, the state C disappears and the H and I states become contiguous in the phase diagram, and separated by a discontinuous phase transition.
6
A.T.C. Silva et al. / Physica A (
)
–
Fig. 4. Phase diagram in the variables r2 versus r1 obtained by the pair mean-field approximation (f = 3/4). (a) b1 = 1, b2 = 2/3 and β = 1. (b) b1 = 1/3, b2 = 3 and β = 1. The initial conditions were: x = 0.4, y = 0.4, s = 0.07, v = 0.07, and w = 0.15. A continuous line represents a continuous phase transition whereas a dashed line represents a discontinuous phase transition. The states are: extinct (E), healthy (H), infected (I) and a state where there is a coexistence (C) of infected and healthy hosts.
Fig. 5. Phase diagram, obtained by the pair mean-field approximation (f = 3/4), in the variables h = (β/b2 r2 )(f − r1 ) and 1 − v where v = r1 /r2 , for the following values of the parameters: (a) b1 = 0.6, u1 = 0.1, u2 = 0.08; (b) b1 = 0.6, u1 = 0.1, u2 = 0.2. The states are: healthy (H), infected (I) and a state where there is a coexistence (C) of infected and healthy hosts. A continuous line represents a continuous phase transition whereas a dashed line represents a discontinuous phase transition. The full circle is a bicritical point.
Again it is convenient to represent the possible states in a phase diagram in the variable h and v , which are now defined by
v= h=
r1 r2
,
β b2 r2
(29)
(f − r1 ).
(30)
In these variables, the phase diagram is shown in Fig. 5. The state C is found to lie between the states I and H, separated from them by continuous transition lines, as shown in Fig. 5(a). However, for some values of the parameters, the state I and H might be contiguous and separated by a discontinuous phase transition line which ends on a bicritical point, as can be seen in Fig. 5(b). 5. Simulations We performed Monte Carlo simulations of the stochastic lattice model by the use of the following algorithm. At each time step, a site of the lattice is chosen at random, say site i. (a) If site i is empty, then one of its four neighbor sites is chosen at random, say site j. If the neighboring site j is occupied by a healthy host, then site i becomes occupied by a healthy host with probability p01 ; if j is occupied by an infected host then site i becomes occupied by an infected host with probability p02 ; otherwise site i remains unchanged. (b) If site i is occupied by a healthy host then with probability p10 site i becomes empty and with probability p12 a neighboring site is chosen at random; if the neighboring site is occupied by an infected
A.T.C. Silva et al. / Physica A (
)
–
7
Fig. 6. Phase diagram in the variables r2 versus r1 obtained by numerical simulations on a square lattice. (a) b1 = 1, b2 = 2/3 and β = 1. (b) b1 = 1/3, b2 = 3 and β = 1. A continuous line represents a continuous phase transition whereas a dashed line represents a discontinuous phase transition. The states are: extinct (E), healthy (H), infected (I) and a state where there is a coexistence (C) of infected and healthy hosts.
host, site i becomes infected. Otherwise, the state of site i remains unchanged. (c) If site i is occupied by an infected host, it becomes empty with probability p20 . Otherwise, the site remains in the same state. Notice that the parameters p01 , p02 , p10 , p20 , and p12 are, proportional to the rates b1 , b2 , u1 , u2 , and β , respectively. We have performed numerical simulations on a square lattice with N = 100 × 100 sites with periodic boundary conditions. We started with a random configuration in which each site was empty, occupied by healthy host or by an infected host, with equal probabilities. For a fixed set of values of the parameters, we used up to 106 Monte Carlo steps to measure the quantities of interest. By observing these quantities as functions of time, it is possible to know when these quantities has reached a steady state. We observed that, after 10% of the total time, the system had already reached the steady state. As a rule, we have discarded the first 10% of the iterations. We have determined the averages ⟨N1 ⟩ and ⟨N2 ⟩ of the number of healthy and infected hosts and the densities of healthy and infected hosts, x = ⟨N1 ⟩/N and y = ⟨N2 ⟩/N. From these results we have determined the phase diagram shown in Fig. 6. As before, the state E is defined by x = y = 0, the state H by x ̸= 0 and y = 0, the state I by x = 0 and y ̸= 0, and C by x ̸= 0 and y ̸= 0. When one crosses the critical line EH of the phase diagram from the region E to the region H, the infected hosts are absent in both sides of the line. Because the infected hosts are absent, the reactions involve only the empty sites and the healthy hosts. The model reduces thus to the contact model with the two reactions: a catalytic infection reaction 0 → 1, with rate b1 and a spontaneous recovery reaction 1 → 0, with rate u1 . It is well known that the contact process undergoes a phase transition when the ratio r between the recovery rate and the infection rate equals rc = 0.608 [21]. Therefore, the line EH is defined by r1 = u1 /b1 = 0.608 as shown in Fig. 6. By a similar reasoning, we conclude that the transition line IE occurs when r2 = u2 /b2 = 0.608 as shown in Fig. 6. 6. Discussion and conclusion remarks We have analyzed a stochastic lattice model for vertical and horizontal transmitted infection by means of mean-field approximations and numerical simulations. The results for simple mean-field approximation are summarized in Figs. 2 and 3. In Fig. 2, we analyze the contribution of horizontal transmission, fixing the value of parameter β , in two cases: β > b2 − b1 (see Fig. 2(a)) and β ≤ b2 − b1 (see Fig. 2(b)). In Fig. 2(a), we observe a continuous phase transition based on the critical lines given by Eqs. (10) and (12), that correspond, respectively, to healthy–coexistence (HC) transition and coexistence–infected (CI) transition. In Fig. 2(b), for which the action of horizontal transmission is smaller than in previous case, we observe the bi-stability of phases H and I, limited by the two spinodal lines, with no coexistence (C). The particular case β = 0, for which there is no horizontal transmission, there is no coexistence of infected and uninfected hosts, but there is a co-stability. This case may be relevant for the infections caused by parasites, like endophytic fungi, for which only vertical transmission occurs [22]. In Fig. 3, assuming the change of variables previously proposed by Lipsitch et al. [7], the phase diagram takes into account variations of both vertical and horizontal transmissions. In Fig. 3(a), for which u2 < u1 , we recover the results previously observed in Ref. [7]. However, for a different region of parameter space, for which u2 > u1 , Fig. 3(b) reveals a discontinuous phase transition between phases I and H, for low values of horizontal transmission h and high values of vertical transmission v . It may be relevant since, for some parasites like the intestinal microbiota of pups, the vertical transmission is dominant [23]. In the literature, discontinuous phase transitions are reported in different situations in the context of theoretical epidemiology. For instance, in the analysis of reinfection based on the SIRI model [24], when primary infection is less frequent than secondary ones, the epidemic transition turns from continuous to discontinuous. It corresponds to a different region of parameter space since the parameter of reinfection changes from small to large values, a typical situation in many social
8
A.T.C. Silva et al. / Physica A (
)
–
contagious processes but not so common for transmitted diseases. Recently, there are also some theoretical results indicating first-order phase transitions in outbreaks of co-infectious diseases [25]; nowadays co-infections are very frequent for both directly transmitted diseases like Tuberculosis–AIDS, and vector-transmitted diseases like Dengue–Zika–Chikungunya. In Figs. 4 and 5, we have done similar study for pair mean-field approximation. Comparing Figs. 2(a) and 4(a), a better approximation leads to a smaller region of C-phase. The same conclusion can be said when one compares Figs. 4(a) and 6(a), obtained by numerical simulations. Notice, however, that the decrease is due mainly because the line HC is shifted down whereas the line CI remains almost in the same place. Comparing Figs. 2(b) and 4(b), a better approximation leads to a decreasing of the bi-stability region. In numerical simulations, as seen in Fig. 6(b), this region has disappeared giving rise to a discontinuous line between the H and I phases. In summary, in addition to the usual scenario with continuous transition between the I, C and H phases, a different scenario with the suppression of the C phase and a discontinuous phase transition between I and H phases is presented for a model with both horizontal and vertical transmissions. In future work we intend to improve the numerical simulations with the purpose of analyzing the critical behavior of the model. We also intend to study the model on other regular and irregular networks in order to analyze the role of the horizontal and vertical transmissions in terms of the spatial–temporal patterns. Acknowledgments The authors thank FAPESB (contract number PNX 006/2009) for supporting the research. STRP, TT and MJO thank CNPq, grant numbers 306458/2015, 309747/2014 and 306210/2011, respectively. SRTP also thanks Intituto Nacional de Ciência e Tecnologia: Ciência, Inovação, Tecnologia em Saúde (contract number 57386/2008-9). TT and MJO acknowledge Instituto Nacional de Ciência e Tecnologia: Fluidos Complexos (contract number 573560/2008-0). References [1] J.E. Smith, A.M. Dunn, Transovarial transmission, Parasitol. Today 7 (1991) 146–148. [2] E. Zchori-Fein, C.J. Geden, D.A. Rutz, Microsporidioses of Muscidifurax raptor (Hymenoptera: Pteromalidae) and other pteromalid parasitoids of muscoid flies, J. Invertebr. Pathol. 60 (1992) 292–298. [3] W.L. Shoop, Vertical transmission of helmnths: hypobiosis and amphiparatensis, Parasitol. Today 7 (1991) 51–54. [4] A.H. Purcell, K.G. Suslow, M. Klein, Transmission via plants of an insect pathogenic bacterium that does not multiply or move in plants, Microb. Ecol. 27 (1994) 19–26. [5] C.A. Mims, Vertical transmission of viruses, Microbiol. Rev. 45 (1981) 267–286. [6] M. Narita, M. Shibata, T. Togashi, Y. Koga, Vertical transmission of human T-cell leukemia virus type I, J. Infect. Dis. 163 (1991) 204. [7] M. Lipsitch, M. Nowak, D. Ebert, R. May, The population dynamics of vertically and horizontally transmitted parasites, Proc. R. Soc. B 260 (1995) 321–327. [8] M. Lipsitch, S. Siller, M.A. Nowak, The evolution of virulence in pathogens with vertically and horizontally transmission, Evolution 50 (1996) 1729–1741. [9] A.P. Galvani, Epidemiology meets evolutionary ecology, Trends Ecol. Evol. 18 (2003) 132–139. [10] M. Besnard, S. Lastère, A. Teissier, V.M. Cao-Lormeau, D. Musso, Evidence of perinatal transmission of Zika Virus, French Polynesia, December 2013 and February 2014, Eurosurveillance 19 (2014) 20751. [11] G. Cavet, R.S. Aguiar, A.S.O. Melo, S.A. Sampaio, I. de Filippis, A. Fabri, E.S.M. Araújo, P.C. de Siqueira, M.C.L. de Mendonça, L. de Oliveira, D.A. Tschoeke, C.G. Schrago, F.L. Thompson, P. Brasil, F.B. dos Santos, R.M.R. Nogueira, A. Tanuri, A.M.B. de Filippis, Detection and sequencing of Zika virus from amniotic fluid of fetuses with microcephaly in Brazil: a case study, Lancet Infect. Dis. 16 (2016) 653–660. [12] R.M. Anderson, R.M. May, Population biology of infectious diseases: Part I, Nature 280 (1979) 361–367. [13] J.E. Satulovsky, T. Tomé, Stochastic lattice gas model for a predator–prey system, Phys. Rev. E 49 (1994) 5073–5079. [14] K. Satô, H. Matsuda, A. Sasaki, Pathogen invasion and host extinction in lattice sctructured populations, J. Math. Biol. 32 (1994) 251–268. [15] Y. Haraguchi, A. Sasaki, The evolution of parasite virulence and transmisstion rate in a spatially structured population, J. Theoret. Biol. 203 (2000) 85–96. [16] M. Boots, A. Sasaki, The evolutionary dynamics of local infection and global reproduction in host-parasite interactions, Ecol. Lett. 3 (2000) 181–185. [17] R.B. Schinazi, Horizontal versus vertical transmission of parasites in a stochastic spatial model, Math. Biosci. 168 (2000) 1–8. [18] M.J. Keeling, P. Rohani, Modeling Infectious Diseases in Human and Animals, Princeton University Press, Princeton, 2008. [19] A.L. Rodrigues, T. Tomé, Reactions-diffusion stochastic lattice model for a predator–prey system, Braz. J. Phys. 38 (2008) 87–93. [20] D.R. de Souza, T. Tomé, Stochastic lattice gas model describing the dynamics of the SIRS epidemic process, Physica A 389 (2010) 1142–1150. [21] T. Tomé, M.J. de Oliveira, Stochastic Dynamics and Irreversibility, Springer, Heidelberg, 2015. [22] S.H. Faeth, K.P. Hadeler, H.R. Thieme, An apparent paradox of horizontal and vertical disease transmission, J. Biol. Dyn. 1 (2007) 45–62. [23] R. Inoue, K. Ushida, Vertical and horizontal transmission of intestinal commensal bacteria in the rat model, FEMS Microbiol. Ecol. 46 (2003) 213–219. [24] J. Gomez-Gardenes, A.S. de Barros, S.T.R. Pinho, R.F.S. Andrade, Abrupt transitions from reinfections in social contagions, Europhys. Lett. 110 (2015) 58006. [25] H.-K. Janssen, O. Stenull, First-order phase transitions in outbreaks of co-infectious diseases an the extended general epidemic process, Europhys. Lett. 113 (2016) 26005.
)
–
Contents lists available at ScienceDirect
Physica A journal homepage: www.elsevier.com/locate/physa
Stochastic spatial structured model for vertically and horizontally transmitted infection Ana T.C. Silva a , Vladimir R.V. Assis a , Suani T.R. Pinho b , Tânia Tomé c,∗ , Mário J. de Oliveira c a
Departamento de Física, Universidade Estadual de Feira de Santana, 44036-900 Feira de Santana, BA, Brazil
b
Instituto de Física, Universidade Federal da Bahia, 40210-340 Salvador, BA, Brazil Instituto de Física, Universidade de São Paulo, Rua do Matão, 1371, 05508-090 São Paulo, SP, Brazil
c
highlights • • • •
An epidemic process with horizontal and vertical transmission is analyzed. The epidemic stochastic lattice model is studied by Monte Carlo simulations. The epidemic stochastic lattice model is studied also by a pair mean-field approach. The phase diagram may display a healthy and infected coexisting state.
article
info
Article history: Received 6 July 2016 Received in revised form 1 October 2016 Available online xxxx Keywords: Stochastic lattice model Model for spreading of disease Model for vertical and horizontal transmitted infection
abstract We study a space structured stochastic model for vertical and horizontal transmitted infection. By means of simple and pair mean-field approximation as well as Monte Carlo simulations, we construct the phase diagram, which displays four states: healthy (H), infected (I), extinct (E), and coexistent (C). In state H only healthy hosts are present, whereas in state I only infected hosts are present. The state E is characterized by the extinction of the hosts whereas in state C there is a coexistence of infected and healthy hosts. In addition to the usual scenario with continuous transition between the I, C and H phases, we found a different scenario with the suppression of the C phase and a discontinuous phase transition between I and H phases. © 2016 Elsevier B.V. All rights reserved.
1. Introduction The transmission of disease between hosts of the same species occurs either vertically, when the pathogen is transmitted directly from the parents to the child, or horizontally, when the hosts are not in a parent–child relationship. Several diseases are transmitted in both horizontal and vertical ways [1–11]. A question that arises in the study of horizontally and vertically transmitted diseases is the determination of the amount of horizontal transmission needed to maintain a vertical transmission. This problem was analyzed by Lipsitch et al. [7] who proposed that in the stationary state, the prevalence of infection can be higher when there are both types of transmissions as compared to the case when only one type is present. This result implies that for the case of virulent pathogens, the high number of infected cases coming from vertical transmission is a consequence of a high rate of horizontal transmission.
∗
Corresponding author. E-mail address: [email protected] (T. Tomé).
http://dx.doi.org/10.1016/j.physa.2016.10.048 0378-4371/© 2016 Elsevier B.V. All rights reserved.
2
A.T.C. Silva et al. / Physica A (
)
–
Fig. 1. Illustration of the reactions of the stochastic lattice model. The states 0, 1 and 2 represent, respectively, an empty site, a site occupied by a healthy host and by an infected host. The quantities u1 , b1 , u2 , b2 , and β are the transition rates.
One approach in the study of population biology of infectious diseases consists in setting up differential equations for the densities of the various types of population [12], similar to equations used in chemical kinetics in which the law of mass action is employed. The model of Lipsitch et al. [7] is of this type and consists of two differential equations for the densities of infected and uninfected populations. It is possible, on the other hand, to set up stochastic population models with space structure to describe the spreading of diseases and, more generally, population dynamics [13–21]. In fact, a space structured stochastic model with the same reactions as the model used by Lipsitch et al. has been introduced by Schinazi [17]. For a particular case, he showed that if the vertical transmission is sufficiently high then infected hosts can be eliminated even if the horizontal transmission is high. Here, we analyze a space structured stochastic model which is a slight modification of the Schinazi’s model. The model is investigated by means of two mean-field approaches, the simple mean-field approximation and the pair mean-field approximation, and by numerical simulations on a square lattice. The evolution equations at the level of simple mean-field approximation are found to be equivalent to the differential equations used by Lipsitch et al. [7]. The analysis of the mean-field equations and the numerical simulations on a square lattice reveal that the model may display four states: an extinction (E) state, where the hosts have disappeared, a healthy (H) state, with healthy hosts and the absence of infected hosts, an infected (I) state, with infected hosts and the absence of healthy hosts, and a coexistence (C) state, where both infected and healthy hosts are present. For some set of values of the parameters, the H and I states are not contiguous and are separated by the C state. The transition lines HC and IC of the phase diagram correspond to continuous phase transitions. However, for other sets of values of the parameters, it is possible that the H and I states be contiguous. In this case the transition line HI corresponds to a discontinuous phase transition. The phase diagram obtained by the simple and pair mean-field approximation as well as the numerical simulations are qualitatively similar. However, there are important quantitative changes. For instance, we have found that the H, I and C regions of the phase diagram obtained by numerical simulations are smaller than their counterpart obtained by pair mean-field approximations which in turn are smaller than those obtained by simple mean-field approximation. 2. Model The stochastic lattice model is defined as follows. Each site of a regular lattice is either empty or occupied by a host that can be healthy or infected. A healthy host dies spontaneously with rate u1 and is born by a catalytic process with rate strength b1 . Analogously, an infected host dies spontaneously with rate u2 and is born by a catalytic process with rate strength b2 , a reaction related to vertical transmission. A healthy host becomes infected catalytically with rate strength β , a reaction understood as horizontal transmission. These reactions are illustrated in Fig. 1. Denoting the state of site i by the stochastic variable ηi , then ηi takes the values 0, 1, or 2, according to whether the site is empty, occupied by a healthy host or occupied by an infected host. The possible transitions of the state of a given site i with the respective rates are as follows Reaction Rate
0→1 b 1 ni
1→0 u1
0→2 b2 m i
2→0 u2
1→2 β mi
where ni and mi denote the number of healthy hosts and the number of infected hosts in the neighborhood of site i, respectively. The neighborhood of a site in a regular lattice is chosen to be its nearest neighbor sites. The model has five parameters: b1 , u1 , b2 , u2 , and β . Sometimes it will be convenient to use the parameters r1 and r2 defined by r1 =
u1 b1
,
and r2 =
u2 b2
.
(1)
When u1 = u2 , the present stochastic model reduces to a particular case of the stochastic model introduced by Schinazi [17]. The transition rates above define a stochastic process characterized by a master equation, which gives the evolution of the probability distribution P (η) associated to the microscopic state of the system η = {ηi }. From the master equation it is possible to write down the time evolution of the several marginal probability distributions such as P (ηi ), the probability
A.T.C. Silva et al. / Physica A (
)
–
3
of site i being in state ηi , and P (ηi , ηj ), the probability of neighboring sites i and j being in states ηi and ηj , respectively. Particularly, the time evolution equations for the one-site probabilities are d dt d
P (1) = b1 P (10) − u1 P (1) − β P (12),
P (2) = b2 P (20) − u2 P (2) + β P (12). dt Since P (0) + P (1) + P (2) = 1, the evolution equation for P (0) can be obtained from Eqs. (2) and (3).
(2) (3)
3. Simple mean-field approximation We solve Eqs. (2) and (3) by the use of a simple mean-field approach which amounts to employ the approximation P (ηi , ηj ) = P (ηi )P (ηj ) [21], which means to say that the correlations are neglected in this approximation. Using the notation x = P (1), which is the density of healthy hosts, y = P (2), which is the density of infected hosts, and z = P (0) = 1 − x − y, we may write P (10) = xz, P (20) = yz, and P (12) = xy, so that the set of Eqs. (2) and (3) reduces to the following equations dx dt dy
= b1 xz − u1 x − β xy,
(4)
= b2 yz − u2 y + β xy. (5) dt These equations are equivalent to the equations analyzed by Lipsitch et al. [7]. To see that this is indeed the case, it suffices to perform the rescaling x → x/K and y → y/K , where K is the carrying capacity. Here, the property that the host population is density limited is taken into account by considering that hosts proliferate only on empty sites of the lattice. This procedure has been employed for population biology systems described by stochastic lattice models [13]. The simple mean-field equations predict four types of stationary solutions, or states, whose stability can be found by analyzing the eigenvalues λ1 and λ2 of the Hessian
b 1 ( z − x ) − u1 − β y −(b2 − β)y
−(b1 + β)x . b2 (z − y) − u2 + β x
(6)
If the real part of the eigenvalues is negative, the corresponding state is stable. One solution is the trivial solution x = 0, y = 0, and z = 1, characterized by the extinction of hosts, which we call the extinction (E) state. This trivial solution is stable when u1 > b1 and u2 > b2 so that the lines of stability are given by r1 = 1,
(r2 ≥ 1), (r1 ≥ 1),
r2 = 1,
(7) (8)
which are, respectively, the vertical and horizontal lines shown in Fig. 2. When we cross the line defined by Eq. (7), there appears a state with the absence of infection, where all hosts are healthy, which we call healthy (H) state, characterized by x ̸= 0 and y = 0, that is, x = 1 − r1
and
y = 0.
(9)
From the stability condition, we found that this state is stable in the region limited by the line given by Eq. (7) and by the line given by b2 (r2 − r1 ) + β(r1 − 1) = 0,
(10)
as shown in Fig. 2. The line given by (10) is the upper inclined line in Fig. 2(a) and the lower inclined dashed line in Fig. 2(b). When we cross the line defined by Eq. (8), there emerges a state with the prevalence of infection, with no healthy hosts, characterized by x = 0 and y ̸= 0, that is, x=0
and
y = 1 − r2 .
(11)
From the stability condition, we found that this state is stable in the region limited by the line given by Eq. (8) and by the line given by b1 (r2 − r1 ) + β(r2 − 1) = 0,
(12)
as shown in Fig. 2. The line given by (12) is the lower inclined line in Fig. 2(a) and the upper inclined dashed line in Fig. 2(b). In addition to these three states, there is an additional state where the populations of healthy and infected hosts coexist, which we call coexistence (C) state, characterized by x ̸= 0 and y ̸= 0. From the stationary solution of Eqs. (4) and (5), the densities of healthy and infected hosts are x= y=
b1 b2 (r2 − r1 ) + β b2 (r2 − 1)
β(b1 − b2 + β) b1 b2 (r1 − r2 ) − β b1 (r1 − 1)
β(b1 − b2 + β)
,
(13)
.
(14)
4
A.T.C. Silva et al. / Physica A (
)
–
Fig. 2. Phase diagram in the variables r2 versus r1 obtained by the simple mean-field approximation. (a) b1 = 1, b2 = 2/3 and β = 1. (b) b1 = 1/3, b2 = 3 and β = 1. The states are: extinct (E), healthy (H), infected (I) and a state where there is a coexistence (C) of infected and healthy hosts. Solid lines represent continuous phase transitions. The upper and lower dashed lines are the spinodal lines, that is, the limits of stability of states I and H, respectively.
These solutions are the same as found in Ref. [7]. Taking the limits x → 0 and y → 0, we may find the boundaries of the state C in the phase diagram. We see that these two boundaries coincide with the lines given by Eqs. (12) and (10), respectively. We remark that, for some values of the parameters, the state C may not exist. If
β ≤ b2 − b1 ,
(15)
then the line of stability of state H, given by Eq. (10), will remain below the line of stability of state I, given by Eq. (12), as shown in Fig. 2(b), and state C does not exist anymore, given rise to a co-stability of states H and I. It is worth mentioning that the C state should not be confused with the co-stability of states H and I.1 In the y–x plane, the former corresponds to a single stable fixed point with x ̸= 0 and y ̸= 0, whereas the later corresponds to two stable fixed points, one with x = 0 and y ̸= 0 and the other with x ̸= 0 and y = 0. It is convenient to represent the possible states in a phase diagram in the variables h and v , defined by
v= h=
r1 r2
,
(16)
β(1 − r1 ) b2 r2
,
(17)
as introduced by Lipsitch et al. [7]. In these variables the two lines given by Eqs. (12) and (10) are h=
b1 r1 (1 − r1 )(1 − v) b2 r2 (v − r1 )
,
h = 1 − v,
(18) (19)
respectively. The simple mean-field results using these variables are shown in Fig. 3. It should be remarked that two types of phase diagrams are obtained. If u1 /u2 ≥ 1 than the line given by Eq. (18) is always above the line given by Eq. (19), as shown in Fig. 3(a). However, if u1 /u2 < 1 then these two lines may cross each other at a nonzero value of h and 1 − v , as shown in Fig. 3(b). In this case the two states I and H become contiguous and separated by a discontinuous phase transition up to a bicritical point. 4. Pair approximation Now we take into account that the states of neighboring sites might be correlated. To this end we derive from the master equation the time evolution equations for the marginal probability distributions P (ηi , ηj ) of two nearest neighbor sites i and j being in states ηi and ηj , respectively. Thus, in addition to Eqs. (2) and (3), we write down the evolution equations for these marginal probability distributions. Using a regular lattice, the evolution equations for the pair probabilities P (ηi , ηj ) are d dt
P (10) = u2 P (12) + u1 P (11) − u1 P (10) − b1 gP (10)
+ b1 fP (100) − b1 fP (101) − b2 fP (102) − β fP (210),
(20)
1 Using the thermodynamic language, the co-stability of state H and I corresponds to a two-phase coexistence whereas the coexisting C state corresponds to a single phase.
A.T.C. Silva et al. / Physica A (
)
–
5
Fig. 3. Phase diagram, obtained by the simple mean-field approximation, in the variables h = (β/b2 r2 )(1 − r1 ) and 1 − v where v = r1 /r2 , for the following values of the parameters: (a) b1 = 0.6, u1 = 0.1 and u2 = 0.08; (b) b1 = 0.6, u1 = 0.1 and u2 = 0.12. The states are: healthy (H), infected (I) and a state where there is a coexistence (C) of infected and healthy hosts. Solid lines represent continuous phase transitions. The upper and lower dashed lines are the spinodal lines, that is, the limits of stability of states H and I, respectively. The full circle is a bicritical point.
d dt
P (20) = u1 P (21) + u2 P (22) − u2 P (20) − b2 gP (20)
+ b2 fP (200) − b2 fP (202) − b1 fP (201) + β fP (210),
(21)
d
P (21) = −u1 P (21) − u2 P (21) + b2 fP (201) + b1 fP (201) + β fP (211) − β fP (212) − β gP (21), (22) dt where f = (ζ − 1)/ζ and g = 1 − f , with ζ being the coordination number of the lattice, that is, the number of neighboring sites of a given site of the regular lattice. In these equations, P (ηi , ηj , ηk ) denotes the three-site probability distribution related to the sites i, j, and k, where i and k are neighboring sites of a central site j, and i ̸= k. It should be noted that there are several types of three-site probabilities which should be distinguished by a different notation. However, we are not using distinct notations in view of the approximation we are about to use. In this approximation the several types of three-site probabilities become equal to each other. The set of Eqs. (20), (21), (22) cannot be solved by considering only the two-site correlations P (01), P (02), and P (12) because these equations involve correlations of greater order. A truncation scheme [21] will then be used, which amounts to employ the following approximation for the three-site correlations P (ηi , ηj , ηk ) =
P (ηi , ηj )P (ηj , ηk ) P (ηj )
.
(23)
This is the pair approximation: correlations of three or more sites are written in terms of one-site and two-site correlations Using the notation s = P (10), v = P (20), and w = P (12), so that P (11) = x − s −w , P (22) = y −v−w , and P (00) = z − s −v , Eqs. (2) and (3), and Eqs. (20), (21), and (22) are written as dx dt dy dt ds dt dv dt dw
= b1 s − u1 x − βw,
(24)
= b2 v − u2 y + βw,
(25)
= (u2 − u1 )w + u1 x − (2u1 + b1 g )s + b1 f
s(z − s − v)
= (u1 − u2 )w + u2 y − (2u2 + b2 g )v + b2 f = −(u1 + u2 + β g )w + (b1 + b2 )f
sv
z
− b1 f
v(z − s − v)
+ βf
z
s2
− b2 f
w(x − s − w)
− b2 f
z
v2 z
− βf
sv z
− b1 f
w2
− βf
sv z
ws x
+ βf
,
ws x
(26)
,
(27)
. (28) dt z x x This set of five equations makes up a closed set of equations for the variables x, y, s, v and w . We have solved these equations numerically and found that they predict the same four types of solutions, or states, found in Section 3 by the use of the simple mean-field approximation. The state E, characterized by x = 0, y = 0, s = 0, v = 0, and w = 0, is stable as long as r1 > f and r2 > f . When this condition is not met, other states emerge. When the state E becomes unstable, it gives rise to the state H or to state I, as can be seen in Fig. 4. The state C lies between these two states from which it is separated by continuous transition lines. If the ratio β/b2 ≪ 1, the state C disappears and the H and I states become contiguous in the phase diagram, and separated by a discontinuous phase transition.
6
A.T.C. Silva et al. / Physica A (
)
–
Fig. 4. Phase diagram in the variables r2 versus r1 obtained by the pair mean-field approximation (f = 3/4). (a) b1 = 1, b2 = 2/3 and β = 1. (b) b1 = 1/3, b2 = 3 and β = 1. The initial conditions were: x = 0.4, y = 0.4, s = 0.07, v = 0.07, and w = 0.15. A continuous line represents a continuous phase transition whereas a dashed line represents a discontinuous phase transition. The states are: extinct (E), healthy (H), infected (I) and a state where there is a coexistence (C) of infected and healthy hosts.
Fig. 5. Phase diagram, obtained by the pair mean-field approximation (f = 3/4), in the variables h = (β/b2 r2 )(f − r1 ) and 1 − v where v = r1 /r2 , for the following values of the parameters: (a) b1 = 0.6, u1 = 0.1, u2 = 0.08; (b) b1 = 0.6, u1 = 0.1, u2 = 0.2. The states are: healthy (H), infected (I) and a state where there is a coexistence (C) of infected and healthy hosts. A continuous line represents a continuous phase transition whereas a dashed line represents a discontinuous phase transition. The full circle is a bicritical point.
Again it is convenient to represent the possible states in a phase diagram in the variable h and v , which are now defined by
v= h=
r1 r2
,
β b2 r2
(29)
(f − r1 ).
(30)
In these variables, the phase diagram is shown in Fig. 5. The state C is found to lie between the states I and H, separated from them by continuous transition lines, as shown in Fig. 5(a). However, for some values of the parameters, the state I and H might be contiguous and separated by a discontinuous phase transition line which ends on a bicritical point, as can be seen in Fig. 5(b). 5. Simulations We performed Monte Carlo simulations of the stochastic lattice model by the use of the following algorithm. At each time step, a site of the lattice is chosen at random, say site i. (a) If site i is empty, then one of its four neighbor sites is chosen at random, say site j. If the neighboring site j is occupied by a healthy host, then site i becomes occupied by a healthy host with probability p01 ; if j is occupied by an infected host then site i becomes occupied by an infected host with probability p02 ; otherwise site i remains unchanged. (b) If site i is occupied by a healthy host then with probability p10 site i becomes empty and with probability p12 a neighboring site is chosen at random; if the neighboring site is occupied by an infected
A.T.C. Silva et al. / Physica A (
)
–
7
Fig. 6. Phase diagram in the variables r2 versus r1 obtained by numerical simulations on a square lattice. (a) b1 = 1, b2 = 2/3 and β = 1. (b) b1 = 1/3, b2 = 3 and β = 1. A continuous line represents a continuous phase transition whereas a dashed line represents a discontinuous phase transition. The states are: extinct (E), healthy (H), infected (I) and a state where there is a coexistence (C) of infected and healthy hosts.
host, site i becomes infected. Otherwise, the state of site i remains unchanged. (c) If site i is occupied by an infected host, it becomes empty with probability p20 . Otherwise, the site remains in the same state. Notice that the parameters p01 , p02 , p10 , p20 , and p12 are, proportional to the rates b1 , b2 , u1 , u2 , and β , respectively. We have performed numerical simulations on a square lattice with N = 100 × 100 sites with periodic boundary conditions. We started with a random configuration in which each site was empty, occupied by healthy host or by an infected host, with equal probabilities. For a fixed set of values of the parameters, we used up to 106 Monte Carlo steps to measure the quantities of interest. By observing these quantities as functions of time, it is possible to know when these quantities has reached a steady state. We observed that, after 10% of the total time, the system had already reached the steady state. As a rule, we have discarded the first 10% of the iterations. We have determined the averages ⟨N1 ⟩ and ⟨N2 ⟩ of the number of healthy and infected hosts and the densities of healthy and infected hosts, x = ⟨N1 ⟩/N and y = ⟨N2 ⟩/N. From these results we have determined the phase diagram shown in Fig. 6. As before, the state E is defined by x = y = 0, the state H by x ̸= 0 and y = 0, the state I by x = 0 and y ̸= 0, and C by x ̸= 0 and y ̸= 0. When one crosses the critical line EH of the phase diagram from the region E to the region H, the infected hosts are absent in both sides of the line. Because the infected hosts are absent, the reactions involve only the empty sites and the healthy hosts. The model reduces thus to the contact model with the two reactions: a catalytic infection reaction 0 → 1, with rate b1 and a spontaneous recovery reaction 1 → 0, with rate u1 . It is well known that the contact process undergoes a phase transition when the ratio r between the recovery rate and the infection rate equals rc = 0.608 [21]. Therefore, the line EH is defined by r1 = u1 /b1 = 0.608 as shown in Fig. 6. By a similar reasoning, we conclude that the transition line IE occurs when r2 = u2 /b2 = 0.608 as shown in Fig. 6. 6. Discussion and conclusion remarks We have analyzed a stochastic lattice model for vertical and horizontal transmitted infection by means of mean-field approximations and numerical simulations. The results for simple mean-field approximation are summarized in Figs. 2 and 3. In Fig. 2, we analyze the contribution of horizontal transmission, fixing the value of parameter β , in two cases: β > b2 − b1 (see Fig. 2(a)) and β ≤ b2 − b1 (see Fig. 2(b)). In Fig. 2(a), we observe a continuous phase transition based on the critical lines given by Eqs. (10) and (12), that correspond, respectively, to healthy–coexistence (HC) transition and coexistence–infected (CI) transition. In Fig. 2(b), for which the action of horizontal transmission is smaller than in previous case, we observe the bi-stability of phases H and I, limited by the two spinodal lines, with no coexistence (C). The particular case β = 0, for which there is no horizontal transmission, there is no coexistence of infected and uninfected hosts, but there is a co-stability. This case may be relevant for the infections caused by parasites, like endophytic fungi, for which only vertical transmission occurs [22]. In Fig. 3, assuming the change of variables previously proposed by Lipsitch et al. [7], the phase diagram takes into account variations of both vertical and horizontal transmissions. In Fig. 3(a), for which u2 < u1 , we recover the results previously observed in Ref. [7]. However, for a different region of parameter space, for which u2 > u1 , Fig. 3(b) reveals a discontinuous phase transition between phases I and H, for low values of horizontal transmission h and high values of vertical transmission v . It may be relevant since, for some parasites like the intestinal microbiota of pups, the vertical transmission is dominant [23]. In the literature, discontinuous phase transitions are reported in different situations in the context of theoretical epidemiology. For instance, in the analysis of reinfection based on the SIRI model [24], when primary infection is less frequent than secondary ones, the epidemic transition turns from continuous to discontinuous. It corresponds to a different region of parameter space since the parameter of reinfection changes from small to large values, a typical situation in many social
8
A.T.C. Silva et al. / Physica A (
)
–
contagious processes but not so common for transmitted diseases. Recently, there are also some theoretical results indicating first-order phase transitions in outbreaks of co-infectious diseases [25]; nowadays co-infections are very frequent for both directly transmitted diseases like Tuberculosis–AIDS, and vector-transmitted diseases like Dengue–Zika–Chikungunya. In Figs. 4 and 5, we have done similar study for pair mean-field approximation. Comparing Figs. 2(a) and 4(a), a better approximation leads to a smaller region of C-phase. The same conclusion can be said when one compares Figs. 4(a) and 6(a), obtained by numerical simulations. Notice, however, that the decrease is due mainly because the line HC is shifted down whereas the line CI remains almost in the same place. Comparing Figs. 2(b) and 4(b), a better approximation leads to a decreasing of the bi-stability region. In numerical simulations, as seen in Fig. 6(b), this region has disappeared giving rise to a discontinuous line between the H and I phases. In summary, in addition to the usual scenario with continuous transition between the I, C and H phases, a different scenario with the suppression of the C phase and a discontinuous phase transition between I and H phases is presented for a model with both horizontal and vertical transmissions. In future work we intend to improve the numerical simulations with the purpose of analyzing the critical behavior of the model. We also intend to study the model on other regular and irregular networks in order to analyze the role of the horizontal and vertical transmissions in terms of the spatial–temporal patterns. Acknowledgments The authors thank FAPESB (contract number PNX 006/2009) for supporting the research. STRP, TT and MJO thank CNPq, grant numbers 306458/2015, 309747/2014 and 306210/2011, respectively. SRTP also thanks Intituto Nacional de Ciência e Tecnologia: Ciência, Inovação, Tecnologia em Saúde (contract number 57386/2008-9). TT and MJO acknowledge Instituto Nacional de Ciência e Tecnologia: Fluidos Complexos (contract number 573560/2008-0). References [1] J.E. Smith, A.M. Dunn, Transovarial transmission, Parasitol. Today 7 (1991) 146–148. [2] E. Zchori-Fein, C.J. Geden, D.A. Rutz, Microsporidioses of Muscidifurax raptor (Hymenoptera: Pteromalidae) and other pteromalid parasitoids of muscoid flies, J. Invertebr. Pathol. 60 (1992) 292–298. [3] W.L. Shoop, Vertical transmission of helmnths: hypobiosis and amphiparatensis, Parasitol. Today 7 (1991) 51–54. [4] A.H. Purcell, K.G. Suslow, M. Klein, Transmission via plants of an insect pathogenic bacterium that does not multiply or move in plants, Microb. Ecol. 27 (1994) 19–26. [5] C.A. Mims, Vertical transmission of viruses, Microbiol. Rev. 45 (1981) 267–286. [6] M. Narita, M. Shibata, T. Togashi, Y. Koga, Vertical transmission of human T-cell leukemia virus type I, J. Infect. Dis. 163 (1991) 204. [7] M. Lipsitch, M. Nowak, D. Ebert, R. May, The population dynamics of vertically and horizontally transmitted parasites, Proc. R. Soc. B 260 (1995) 321–327. [8] M. Lipsitch, S. Siller, M.A. Nowak, The evolution of virulence in pathogens with vertically and horizontally transmission, Evolution 50 (1996) 1729–1741. [9] A.P. Galvani, Epidemiology meets evolutionary ecology, Trends Ecol. Evol. 18 (2003) 132–139. [10] M. Besnard, S. Lastère, A. Teissier, V.M. Cao-Lormeau, D. Musso, Evidence of perinatal transmission of Zika Virus, French Polynesia, December 2013 and February 2014, Eurosurveillance 19 (2014) 20751. [11] G. Cavet, R.S. Aguiar, A.S.O. Melo, S.A. Sampaio, I. de Filippis, A. Fabri, E.S.M. Araújo, P.C. de Siqueira, M.C.L. de Mendonça, L. de Oliveira, D.A. Tschoeke, C.G. Schrago, F.L. Thompson, P. Brasil, F.B. dos Santos, R.M.R. Nogueira, A. Tanuri, A.M.B. de Filippis, Detection and sequencing of Zika virus from amniotic fluid of fetuses with microcephaly in Brazil: a case study, Lancet Infect. Dis. 16 (2016) 653–660. [12] R.M. Anderson, R.M. May, Population biology of infectious diseases: Part I, Nature 280 (1979) 361–367. [13] J.E. Satulovsky, T. Tomé, Stochastic lattice gas model for a predator–prey system, Phys. Rev. E 49 (1994) 5073–5079. [14] K. Satô, H. Matsuda, A. Sasaki, Pathogen invasion and host extinction in lattice sctructured populations, J. Math. Biol. 32 (1994) 251–268. [15] Y. Haraguchi, A. Sasaki, The evolution of parasite virulence and transmisstion rate in a spatially structured population, J. Theoret. Biol. 203 (2000) 85–96. [16] M. Boots, A. Sasaki, The evolutionary dynamics of local infection and global reproduction in host-parasite interactions, Ecol. Lett. 3 (2000) 181–185. [17] R.B. Schinazi, Horizontal versus vertical transmission of parasites in a stochastic spatial model, Math. Biosci. 168 (2000) 1–8. [18] M.J. Keeling, P. Rohani, Modeling Infectious Diseases in Human and Animals, Princeton University Press, Princeton, 2008. [19] A.L. Rodrigues, T. Tomé, Reactions-diffusion stochastic lattice model for a predator–prey system, Braz. J. Phys. 38 (2008) 87–93. [20] D.R. de Souza, T. Tomé, Stochastic lattice gas model describing the dynamics of the SIRS epidemic process, Physica A 389 (2010) 1142–1150. [21] T. Tomé, M.J. de Oliveira, Stochastic Dynamics and Irreversibility, Springer, Heidelberg, 2015. [22] S.H. Faeth, K.P. Hadeler, H.R. Thieme, An apparent paradox of horizontal and vertical disease transmission, J. Biol. Dyn. 1 (2007) 45–62. [23] R. Inoue, K. Ushida, Vertical and horizontal transmission of intestinal commensal bacteria in the rat model, FEMS Microbiol. Ecol. 46 (2003) 213–219. [24] J. Gomez-Gardenes, A.S. de Barros, S.T.R. Pinho, R.F.S. Andrade, Abrupt transitions from reinfections in social contagions, Europhys. Lett. 110 (2015) 58006. [25] H.-K. Janssen, O. Stenull, First-order phase transitions in outbreaks of co-infectious diseases an the extended general epidemic process, Europhys. Lett. 113 (2016) 26005.