On the effect of postponing pregnancy in a Zika transmission model

We construct a Zika transmission model to investigate the effect of postponing pregnancy on the infection intensity. We perform analytical and numerical investigations for deterministic and stochastic analysis to obtain the basic reproductive ratio, endemic state, probability of disease extinction, and the probability of outbreak. The results indicate that by reducing the pregnancy rate the mosquito-to-human ratio increases, and, consequently, the basic reproductive ratio increases. Simultaneously, the probability of disease extinction decreases, and the probability of disease outbreak increases. On the other hand, the endemic state of infected infants initially increases with the decrease of the pregnancy recruitment rate, up to a certain level, and decreases as the recruitment rate of pregnancy tends to zero. This work highlights that postponing pregnancy that gives the individual temporary protection for unexpected infected newborns may increase the population infectivity.


Introduction
Zika virus disease (ZIKV) has been known as an infectious disease that is transmitted through the bite of Aedes type mosquitoes, including Aedes aegypti. The virus was first identified in Uganda in 1947 in monkeys from the Zika forest. Further transmission in a human was then reported in 1952 [1]. A large outbreak of Zika disease was first reported in Micronesia in 2007 and transmitted to French Polinesia in 2012-2013 and to South America in 2015 [2]. Since then, the disease has been identified in 34 countries [3]. The incubation period for Zika is approximately 2-12 days [4], and most infected people do not indicate symptoms. Simple symptoms that may appear are fever, itchy maculopapular rash, nonpurulent conjunctivitis, and arthralgia [5]. These simple symptoms are indistinguishable from those of other arbovirus diseases. The disease got worldwide attention after the report from Brazil in 2015, confirming 141 suspected microcephaly cases among newborns, and the number is still growing [2]. This microcephaly disease has created huge concern even in unaffected countries, especially among pregnant women who plan to travel to high-risk countries [6]. WHO issued an interim report in September 2016 stating that "Couples or women planning a pregnancy and returning from Zika's affected area are advised to wait at least six months before trying to conceive to ensure that possible Zika virus infection has cleared" [7]. This advice raises some questions regarding the role of postponing pregnancy in the control and prevention of Zika transmission.
Simple mathematical models for Zika transmission have been developed recently. Due to the similarity of the infection process to dengue, the authors in [8,9] adopted SIR-SI and SEIR-SEI models, which are commonly used in dengue cases for Zika. Another model, given in [10], takes into account the infection process human-mosquito through mosquito bites and human-human through sexual transmission route. Intervention strategy such as reducing contacts between humans and mosquitoes can be well accommodated in the model as shown in [11]. In those models, the main goal is to obtain an expression for the basic reproductive ratio R 0 .
For practical application, the interest is how to estimate R 0 . In the case of dengue, in which daily incidence data are available in most affecting countries, estimation of R 0 can be well estimated from the data during the early infection period [12,13], in which the exponential growth still took place. With limited Zika data, estimation of R 0 in endemic countries significantly varies. Using dengue information, the authors in [14] estimated that R 0 for Zika in French Polynesia was comparable with dengue and chikungunya. Using a similar dengue model and surveillance data for Zika, the authors in [15] and [16] show different estimates for R 0 . These differences are understandable due to various assumptions and limited and less accurate data being used in the estimation. The basic reproductive ratio R 0 is important and so far is "the only" existing biological indicator for identifying the potential occurrence of an endemic state at the early stages of an outbreak. However, more indicators are needed to give a comprehensive understanding of the intensity of disease transmission. In the present work, we construct a dynamical model for Zika disease, in which compartments for pregnant women and infected infants are considered.

The model
Simple models for vector-borne diseases have been developed from the simple dengue transmission model SIR-SI [17]. In this model the basic reproductive ratio R 0 is obtained as an endemic indicator representing the average number of secondary infections resulting from a "single" infection entering an all-susceptible population [18]. It is well known that when R 0 ≤ 1, the disease-free equilibrium is stable and becomes unstable if R 0 ≥ 1, followed by the appearance of a stable endemic state. More elaborate structures involving different compartments such as with multistrain dengue virus [19] and involving antibodydependent enhancement (ADE) [20,21], age structure [22] or gender division [10] can be constructed. In the case of Zika, this SIR-SI model is not sufficient to describe the complicated dynamics of disease transmission and the effects of infection on newborns. In the following subsection, we construct a deterministic model for Zika transmission.

Deterministic model
We consider Zika transmission in a closed human-mosquito population in which pregnant women and infected infants are grouped into separate compartments, and newborns are delivered only from pregnant woman compartments. Let S h (t), I h (t), and R h (t) be susceptible, infectious, and recovered (nonpregnant) human subpopulations at time t, respectively. Let S p (t), I p (t), and R p (t) be the numbers of susceptible, infectious, and recovered pregnant women at time t, respectively, and let B(t) denote the number of infected infant subpopulation at time t. We assume that infected infants remain infected throughout their life. In this context, infection of an infant is considered to be microcephaly and never recovered. With this assumption, the dynamic of B(t) is decoupled from the rest of the states. Let S v (t) and I v (t) denote the susceptible and infectious vector subpopulations at time t, respectively. The total host and vector populations are tively. The total population of pregnant women is given by N p (t) = S p (t) + I p (t) + R p (t).
The recruitment rate A v > 0 for mosquitoes (vectors) is assumed to be constant, and the mosquito death rate is μ v > 0. The average lifespan for humans is 1/μ h , so that μ h > 0 is the natural mortality rate for humans. There is no available data on the average lifetime of microcephaly persons. We assume that the average lifetime of microcephaly person is less than the average lifetime of nonmicrocephaly, that is, 1 Let the pregnant women give birth at rate α > 0, whereas pregnant women who are infectious or recovered (I p or R p ) have the possibility of giving birth to an infected infant with probability ρ ∈ [0, 1]. The recruitment rate of pregnancy is assumed to be a constant W > 0, and the rate at which pregnant women revert to a state of nonpregnancy is δ > 0.
It is natural to take δ = α as it is done in the simulation. Throughout this manuscript, we assume that μ h < α, so that the total host population is positive. A description of the state variables and parameters, along with their dimensions, is given in Table 1.
The deterministic model is formulated by the following ordinary differential equations: (1) The total populations satisfy For simplification of the dynamical analysis, we assume that infected infants remain in infected status for the rest of their life and that the population sizes N p and N v are constant in time. Note that in our simulation, we take α = δ = 1/270. Under this assumption, we have

Dynamical analysis
We write the deterministic model (1) in the form Note that with the biological restriction (3), we conclude that F(X(t), e i ≤ 0, which implies that is positively invariant.
Equations (4) have a unique disease-free equilibrium (DFE) given by wherẽ With the constant subpopulations (3), we can write the first equation in (2) as We conclude that the total human population remains bounded and not larger thanN h . As a consequence, the solutions of (1) are bounded in . Noting that infected infants do not produce any secondary infections, there are only three infectious compartments I h , I p , and I v .
Linearizing the system at DFE for these infectious variables, we obtain the Jacobian matrix Rewrite the Jacobian matrix in the form The next-generation matrix NGM for Equations (1) is The basic reproduction ratio R 0 is obtained from the spectral radius of the next-generation matrix NGM: In a more concise form, we can rewrite where π = N v /N h is the mosquito-human ratio before the infection takes place. This threshold is exactly the basic reproductive ratio of host-vector dengue type [17]. Note that R 0 increases as the rate of pregnancy W decreases. This is a direct consequence of reducing the total human population, which increases the mosquito-human ratio π .
To analyze the existence of endemic equilibrium, we rewrite Equation (4) in the form of a linear perturbation problem: where e is the unit vector with e B = 1, e j = 0, j = B, and 0 < 1. It is natural that DFE becomes unstable and a unique endemic equilibrium EE = (S h ,Ī h ,R h ,S p ,Ī p ,R p ,B,S v ,Ī v ) exists for R 0 > 1. The proof of this statement can be referred to [23].
We first consider the unperturbed case of Equation (15), that is, for = 0. We further will analyze the relation between the basic reproductive ratio R 0 and the endemic state of infected infantsB for the case of = 0. By eliminating the equilibrium state (1) we obtain the endemic equilibrium in terms ofB: whereB satisfies the equation with Note that from the equilibriumB in Equation (1) the existence of the endemic stateB has to satisfȳ Substituting B 0 into Equation (17), we have G(B 0 ) < 0. Hence from (17) and (19) we obtain the unique endemic equilibrium The graphs ofB for different values of ρ are given in Fig. 1. Figure 1 shows that as the recruitment rate of pregnant women decreases, the basic reproductive ratio R 0 increases, and the endemic state of infected infantsB increases but eventually decreases as W → 0. This phenomenon shows that the basic reproductive ratio as an endemic threshold only provides an early measurement at the start of the infection process and is not a good representation for measurement of the size of the endemic state.
The existence of a unique endemic equilibrium for the unperturbed system (15) is shown by using the implicit function theorem. With the restriction within two manifolds N p and N v in (3), it suffices to reduce the dynamics of (S h , I h , R h , S p , I p , I v ) in the six-dimensional space. Note that the dynamic of B can be decoupled from the system. In this reduced form the determinant of the Jacobian of (15) at the endemic equilibriumĒE is given by where, under condition (19), The equation for D 2 (B) is given by with where From (25) and (26) we have that det(S(D 2 (B), G(B))) > 0. We conclude that the Sylvester matrix S(D 2 (B), G(B)) has a full rank and {D 2 (B), G(B)} has no common denominator [24]. Using the implicit function theorem, we conclude that the perturbed system (15) has a unique endemic equilibrium for relatively small . Here we use numerical simulation of the unperturbed system to represent the dynamical behavior for 0 < 1. In the next subsection, we develop and analyze a stochastic model to compute the probability of disease extinction prior to a relatively large outbreak.

Continuous-time Markov chain
Consider a discrete-valued random vector where the components S h (t), I h (t), R h (t), . . . , I v (t) are discrete-valued random variables for the numbers of susceptible humans, infectious humans, recovered (nonpregnant) humans, . . . , and infectious vectors at time t, respectively. We construct a continuous-time Markov chain (CTMC) model of Zika transmission with the transitions in the stochastic process {X(t)|t ∈ [0, ∞)} during a small time period t. The possible transitions and their corresponding rates are given in Table 2.

Branching process approximation
As in [25,26], the nonlinear CTMC dynamics is approximated near the DFE using a multitype branching process. The only sources of infection for our model are the states I h , I p , and I v . Therefore we apply the branching process approximation only to these infectious states, and the numbers of susceptible humans, pregnant women, and mosquitoes are assumed to be near the DFE: S h (0) ≈S h , S p (0) ≈S p , and S v (0) ≈S v . We begin by constructing an offspring probability generating function (p.g.f ) f i : [0, 1] 3 → [0, 1], i = 1, 2, 3, for each of the infectious states. Here we use the variables u 1 , u 2 , and u 3 as "dummy" variables representing the three types of infectious states. If I h (0) = 1, I p (0) = 0, and I v (0) = 0, then the offspring p.g.f. for I h is given by represents the infection probability per contact between an infectious host I h and a susceptible vector S v . The term (μ h + γ h )/(β v + μ h + γ h ) represents the probability that an infectious host dies or recovers before causing any secondary infections.
If I h (0) = 0, I p (0) = 1, and I v (0) = 0, then the offspring p.g.f. for I p is given by The term α/(2α +β v + μ h + γ h ) represents the probability that an infectious pregnant woman loses her baby and is considered an infectious host. The termβ v /(2α +β v + μ h + γ h ) represents the infection probability that an infectious pregnant woman infects a susceptible vector. The term (μ h + γ h )/(2α +β v + μ h + γ h ) represents the probability that an infectious pregnant woman dies or recovers before causing any secondary infections.
If I h (0) = 0, I p (0) = 0, and I v (0) = 1, then the offspring p.g.f. for I v is given by represents the probability that an infectious vector infects a susceptible host. The termβ h /(β h +β h + μ v ) represents the probability that an infectious vector infects a susceptible pregnant woman. The term μ v /(β h +β h + μ v ) represents the probability that an infectious vector dies before causing any secondary infections. The expectation matrix for the branching process is the nonnegative 3 × 3 matrix M = [m ij ], where m ij = ∂f j u i , and the partial derivatives are evaluated at (u 1 , u 2 , u 3 ) = (1, 1, 1), that is, the m ij entry represents the average number of type i infectious individuals produced by one infectious individual of type j. In particular, If ρ(M) > 1, then the offspring p.g.f.s have a unique fixed point (q 1 , q 2 , q 3 ) ∈ (0, 1) 3 [26]. The term q 1 represents the probability of disease extinction within the nonpregnant human population, q 2 represents the probability of disease extinction within the pregnant human population, and q 3 represents the probability of disease extinction within the vector population. Thus, for ρ(M) > 1, we have the probability of ultimate disease extinction and the probability of an outbreak can be defined as 1 -P 0 , where an outbreak refers to anything other than disease extinction. Although the spectral radius of M cannot be obtained explicitly, there is a relationship between ρ(M) and the basic reproductive ratio R 0 given by the threshold theorem of Allen and van den Driessche [26]. Since M is irreducible, the matrix F in (11) is nonnegative, the matrix V in (11) is a nonsingular M-matrix, and ρ(M) > 1 if and only if R 0 > 1. The coordinates of the fixed point (q 1 , q 2 , q 3 ) ∈ (0, 1) 3 of the offspring p.g.f.s are given by Note that q 1 = q 2 , which implies that introducing one infectious nonpregnant human into the population will result in the same probability of disease extinction as introducing one infectious pregnant woman.

Numerical simulations
Basic reproductive ratio and probability of disease extinction are important indicators for a comprehensive measure in the field applications. We perform numerical simulations for the deterministic and stochastic systems with the parameter values taken from Table 3 and compare the results. We use Euler's method and Gillespie's algorithm [27] to obtain a numerical solution of the deterministic model and to simulate sample paths for the CTMC model, respectively.

Probability of disease extinction
Using the values of parameters in Table 3, we obtain the basic reproductive ratio R 0 ≈ 3.0.
Solutions of the deterministic model exhibit a major outbreak in both host and vector populations. It is shown that there is either a major outbreak followed by disease extinction or a rapid disease extinction prior to a major outbreak for sample paths of the CTMC model. In the case of an outbreak the sample path agrees well with the ODE solution. The ODEs solution and one sample path of the CTMC model are given in Fig. 2. The probability of extinction P 0 is calculated for various initial conditions using the values of parameters in Table 3. The results are compared to the numerical approximation (Approx.) calculated from the proportion of 10,000 sample paths of the CTMC model with I h (t) + I p (t) + I v (t) = 0 (disease extinction) before reaching an outbreak size of 20. If I h (t) + I p (t) + I v (t) ≥ 20, then it is considered an outbreak. The results are summarized in Table 4.

Comparison of deterministic and stochastic models
In Table 5, it is shown that reducing the pregnancy rate W directly increases the vectorhuman population ratio π , which then increases the basic reproductive ratio R 0 and the probability of an outbreak 1 -P 0 . On the other hand, the endemic level of infected infants B increases up to a certain level and decreases as W → 0.

Conclusions
In this paper, we presented a mathematical model for Zika disease transmission to investigate the effect of postponing pregnancy on the transmission intensity. The model considers separate compartments for pregnant women and infected infants as well as the recruitment rate of pregnancy as a control parameter. By exploiting the deterministic and stochastic models we are able to perform analysis of important indicators, that is, the basic reproductive ratio, endemic states of infected infants, probability of disease extinction,  Table 3 for given initial conditions S h (0) = 994, I h (0) = 1, R h (0) = 0, S p (0) = 4, I p (0) = 1, R p (0) = 0, B(0) = 0, S v (0) = 749, and I v (0) = 1. The probability of disease extinction is P 0 = 0.0368. Table 4 Probability of disease extinction P 0 and results of numerical approximation (Approx.) using 10,000 sample paths of the CTMC model for parameter values taken from Table 3.  and probability of disease outbreak. It is clear that postponing pregnancy only is not always an ideal control for the total population. Although this model is restricted to several assumptions, we believe that this study contributes an important insight to understand the complicated Zika transmission process.