SEX-BIASED PREVALENCE IN INFECTIONS WITH HETEROSEXUAL, DIRECT, AND VECTOR-MEDIATED TRANSMISSION: A THEORETICAL ANALYSIS

. Three deterministic Kermack-McKendrick-type models for study-ing the transmission dynamics of an infection in a two-sex closed population are analyzed here. In each model it is assumed that infection can be transmitted through heterosexual contacts, and that there is a higher probability of transmission from one sex to the other than vice versa. The study is focused on understanding whether and how this bias in transmission reﬂects in sex diﬀerences in ﬁnal attack ratios (i.e. the fraction of individuals of each sex that eventually gets infected). In the ﬁrst model, where the other two transmis- sion modes are not considered, the attack ratios (fractions of the population of each sex that will eventually be infected) can be obtained as solutions of a system of two nonlinear equations, that has a unique solution if the net reproduction number exceeds unity. It is also shown that the ratio of attack ratios depends solely on the ratio of gender-speciﬁc susceptibilities and on the basic reproductive number of the epidemic R 0 , and that the gender-speciﬁc ﬁnal attack-ratio is biased in the same direction as the gender-speciﬁc susceptibilities. The second model allows also for infection transmission through direct, non-sexual, contacts. In this case too, an analytical expression is derived from which the attack ratios can be obtained. The qualitative results are similar to those obtained for the previous model, but another important parameter for determining the value of the ratio between the attack ratios in the two sexes is obtained, the relative weight of direct vs. heterosexual transmission (namely, ρ ). Quantitatively, the ratio of ﬁnal attack ratios generally will not exceed 1.5, if non-sexual transmission accounts for most transmission events ( ρ ≥ 0 . 6) and the ratio of gender-speciﬁc susceptibilities is not too large (say, 5 at most). The third model considers vector-borne, instead of direct transmission. In this case, we were not able to ﬁnd an analytical expression for the ﬁnal attack ratios, but used instead numerical simulations. The results on ﬁnal attack ratios are actually quite similar to those obtained with the second model. It is interesting to note that transient patterns can diﬀer from ﬁnal attack ratios, as new cases will tend to occur more often in the more susceptible sex, while later depletion of susceptibles may bias the ratio in the opposite direction. The analysis of these simple models, despite their lack of realism, can help in providing insight into, and assessment of, the potential role of gender-speciﬁc transmission in infections with multiple modes of transmission, such as Zika virus (ZIKV), by gauging what can be expected to be seen from epidemiological reports of new cases, disease incidence and seroprevalence surveys.


1.
Introduction. Infectious diseases of humans may be transmitted via several modes, such as direct (e.g., influenza, measles, HIV/AIDS), indirect (e.g., gastrointestinal diseases such as cholera, rotavirus infection, cryptosporidiosis), vertical (e.g., HIV/AIDS, HSV-2), vector-borne (e.g., malaria, Zika virus (ZIKV), dengue fever). Multiple modes of transmission have been ascertained for some of those: for instance Toxoplasma gondii can be transmitted through the environment or sexual contacts [15]; hepatitis A can be transmitted through contaminated food or water, or direct contacts including sexual contacts [4]. While the actual occurrence of a particular transmission mode can be ascertained by isolating active virions (or more generally pathogen agents) in a specific site, or occasionally by case reports that unequivocally establish the path through which an individual has become infected, it may be difficult to recognize the relative importance of different transmission modes in sustaining an epidemic.
The transmission of the flavivirus Zika (ZIKV) was initially described as occurring solely through the bite of infected adult female mosquitoes of the genus Aedes. From its discovery in 1952 and the first confirmed human case in 1954 [19] until 2007, confirmed cases of ZIKV infection from Africa and Southeast Asia have been rare [12]. The first widespread epidemic was in 2007 in the Micronesian island of Yap, and there was no evidence of direct human-to-human transmission during this outbreak [8] .
However, data collected during the 2015 ZIKV outbreaks in South America (the largest in history) shows the presence of Zika virions in the semen of some infected males, even several weeks after the likely time of disease exposure [17]. Moreover, several cases of ZIKV infections were documented outside the epidemic area, especially in the United States, that could be explained only through sexual contacts with partners that had travelled to Zika-infected areas [10]. Male-to-female sexual transmission of ZIKV has been clearly established, even when the ZIKV-infected male is asymptomatic [3], and there was a (single) well-documented case of femaleto-male sexual transmission in the USA [6].
Still, the general understanding in the epidemiological community is that transmission through infected mosquitoes is dominant in sustaining the Zika epidemic [1,20]. Correspondingly, mathematical models being proposed for the Zika epidemics (see, for instance, [9,14,22]) generally do not incorporate sexual transmission, mostly for mathematical tractability and lack of reliable data (see however [18]).
However, recently, Coelho et al. [5] suggested instead that sexual transmission may have a much greater role than previously envisaged. Their argument is based on the much larger number of cases of ZIKV reported among women than among men, even after correcting for the bias due to the systematic screening of pregnant women. They argue that the disproportionately large number of infected women can be explained by the fact that male-to-female ZIKV transmission is much larger than the female-to-male transmission in sexual contacts (thereby implying that sexual transmission has, all along, been an important factor in the spread of the Zika flavivirus).
It is quite possible that the significant gender disparity in reported cases is due to a higher attention (of the public health agencies, and of the general public) to ZIKV infection in women, including those that are not pregnant but of childbearing age. We do not intend to discuss this specific issue, as we have no data or direct information source to assess it. Rather, we ask the question of what can be theoretically inferred about the relevance of different transmission modes on unequal prevalence amongst the sexes and how such relevance can be determined from epidemic data, assuming that they have been collected without bias (or that existing biases have been removed).
In the next section, we examine a simple Kermack-McKendrick type [13] deterministic model for ZIKV dynamics that only accounts for heterosexual transmission. In this simple context, we are able to derive a relation between sex-biased transmission of the infection and sex-biased incidence. This model is further extended to also allow for direct transmission of ZIKV (as a proxy for vector-borne transmission), to study how sex-biased incidence depends both on sex-biased infectiousness and on the relative importance of heterosexual transmission. Finally, we consider a model with vector transmission. Although this final model represents a major simplification of the reality of vector-based transmission of ZIKV, it still poses major challenges in deriving analytical expressions for theoretically measuring the impact of such transmission mode. Thus, we perform some numerical simulations to see whether the conclusions obtained in the case of direct transmission also hold in that case.
2. The model with heterosexual transmission only. The simplest model for heterosexual transmission of a disease takes the form of an SIR (susceptible-exposedinfected) Kermack-McKendrick formulation [13]. As we are interested in the epidemic development over a short time-course, we neglect human demography.
It is assumed, for the sake of simplicity, that males and females have the same population size N f = N m = N , and we neglect heterogeneities in sexual behavior. Thus, any individual has, on average, c new sexual partners per unit time. On the other hand, we assume that the probability of ZIKV transmission from an infected female to a susceptible male, denoted by ψ m , may be different from ψ f , the probability of ZIKV transmission from an infected male to a susceptible female. It should be noted that for many infectious diseases, it has been estimated that ψ m < ψ f [2].
Furthermore, we assume that the average duration of the infection may differ between males and females. In other words, the recovery rates -γ m for males and γ f for females -are not necessarily equal. This leads to the following model, where S f and I f denote the density of susceptible and infected females (and similarly for males), and the dot above a state variable denotes its derivative with respect to time. It is worth mentioning that the equations for removed individuals, R f and R m , are not written explicitly because S f + I f + R f ≡ N ≡ S m + I m + R m is assumed constant: where , and β f = ψ f c and β m = ψ m c are, respectively, the male-to-female and female-to-male transmission rates. To analyze the model (1), it is convenient to use a minor modification of the conservation method used for the standard Kermack-McKendrick model, as follows It can be shown, following the same argument as for the standard model (see, e.g., [7]), that lim S ∞ m = lim t→∞ S m (t) (it is obvious that these limits exist and are positive), and applying the conservation equations (2), it follows that We now analyze Equation (3) under the assumption that no individuals were already removed at time t = 0, and that only a small fraction were already infected. That denote the final attack ratios for males and females, respectively, Equation (3) gives By eliminating one unknown at a time in Equation (4), it can be shown that the system has a unique solution z = (z f , z m ) ∈ (0, 1) 2 for each value of ε m and ε f in (0, 1). If we consider the limit of z as ε m and ε f tend to 0, we see that it has to be solution of H(z) = H 1 (z), H 2 (z) = 0, where Concerning Equation (5) Proof. Consider Equation (5); clearly, (0, 0) is a solution of H(z) = 0 for any parameter values. Further, any other solution must lie at the intersection of the curves z 2 = g 1 (z 1 ) and . Both g 1 and g 2 are increasing functions, with g 1 (x) < 0 < g 2 (x) for x ∈ (0, 1), and g 1 (1) < 1, while lim z→1 − g 2 (z) = +∞. Hence, there exists a unique solution z = (z 1 ,z 2 ) of H(z) = 0 if, and only if, g 1 (0) > g 2 (0). Computing g 1 (x) and g 2 (x), one sees this is equivalent to R s 0 > 1. The threshold quantity R s 0 is the basic reproduction number for this model, and can be obtained using the method of next-generation matrix [7]. The superscript s on R s 0 is used to indicate that it corresponds to (hetero)sexual transmission. We shall use other superscripts later to indicate other modes of transmission.
The parameter z f [z m ] represents the fraction of females [males] that has been infected by the end of the epidemic; in epidemiology, this quantity is often termed attack rate (see. e.g., [11]), but we shall use the name attack ratio, as it is a ratio rather than a rate. The focus of our work will actually be z f /z m , i.e. the sex ratio of attack ratios.
Proof. Obviously, we need only to prove the first assertion. Suppose, by contradiction, thatz 1 ≤z 2 . Then, Hence, contradicting the assumption that H(z) = 0.
In order to quantify the relative sizes of the attack ratiosz 1 andz 2 , we find numerical solutions of H(z) = 0. We fix a value for R s 0 and vary the ratio of susceptibilities of females to males, (β f γ f )/(β m γ m ), between 1 and 100 to see the effect onz 1 andz 2 . The results are shown in Fig. 1 for R 0 = 1.25 on the left and for R s 0 = 1.5 on the right. We note that, indeed, increasing the susceptibility of females relative to that of males increases the attack ratio in females and decreases the rate in males, with an overall modest decrease. βmγm (which we called relative susceptibility of female to male, even though the recovery rates are also included), and the y-axis displays the corresponding values z = (z 1 ,z 2 ) solving H(z) = 0, as well as their average (z 1 +z 2 )/2.

3.
Model with heterosexual and direct contact transmission. In this section, we model the scenario where the disease can be transmitted both heterosexually and through other direct modes. Although this is probably not realistic for ZIKV, it may be of interest in other contexts, and is also useful as an approximation for models with vector-borne transmission (see next Section). We further assume that such contacts are independent of the individuals' sex, and that the probability of getting infected per contact is the same for all individuals. Specifically, we assume that each susceptible individual can be infected through such contacts at a rateβ d (I f + I m )/(2N ), whereβ d is an aggregate parameter depending on the number of such contacts per unit time and on the probability of getting infected per contact. The denominator is equal to 2N as this is the total population density, including females and males. Based on the above assumptions, we obtain the following system: It follows, by modifying the conservation equations (2), that having introduced, for ease of notation, the constant β d =β d /2.
Integrating (8) for t from 0 to ∞, and using the same notation as in the previous section, gives (9) By taking the limit as ε m and ε f tend to 0, it follows that of z = (z f , z m ) has to be solution of Clearly, z = (0, 0) is a solution of H(z) = 0. In order to look for other solutions, the first equation of (10) is solved for z 1 to get while solving the second equation for z 2 gives Both g 1 and h 2 are convex functions on [0, 1) satisfying g 1 (0) = h 2 (0) = 0 and Thus, if g 1 (0) ≥ 0, g 1 is invertible as a function from [0, 1) onto [0, +∞). Moreover, if g 1 (0) < 0, then there exists z m > 0 such that g 1 (z m ) = 0 and g 1 is invertible as a function from (z m , +∞) onto (0, +∞). In both cases, we can define an inverse function of g 1 , h 1 : (0, ∞) → (0, 1) say, which is increasing, concave and satisfies lim x→∞ h 1 (x) = 1. Furthermore, if g 1 (0) ≥ 0, then lim x→0 + h 1 (x) = 0 and, conversely, if g 1 (0) < 0, then lim x→0 + h 1 (x) = z m > 0. Note that positive solutions of H(z) = 0 correspond to points (z 1 , z 2 ) such that It follows, from the properties of the functions h 1 and h 2 , that (13) has a unique solution, unless lim that h 1 (0) = 1/g 1 (0), these conditions are equivalent to g 1 (0) ≥ 0 and h 2 (0)g 1 (0) ≥ 1, that can be rewritten as Inequality (15) can be further simplified to give Since it is obvious that inequality (16) implies those in (14), the latter are superfluous and we only need to avoid (16) in order to have a unique solution. Therefore, we have proved the following result.
Theorem 3.1. The system H(z) = 0, defined by (10), admits a unique solution if, and only if, We shall see that inequality (17) is equivalent to R sn 0 > 1, where R sn 0 , defined below, is the reproduction number corresponding to both sexual and non-sexual modes of transmission. Proof. We shall use the approach of the next-generation matrix. The next-generation matrix associated with model (7) is given by Indeed, one can easily see that k i1 = (K) i1 represents the expected number of secondary infections produced in class i = 1 (females) or i = 2 (males) by an average infected female over her infectious period, assuming that the whole population is suceptible, and similarly for k i2 . Then, It follows that It should now be noted that the inequalities in (14) imply the first inequality in (20). Expanding the second inequality in (20), it can be seen that it is equivalent to (16). This implies (14), as we have already seen. Hence, (20) is equivalent to (16), so that we can say that a unique solution of H(z) = 0 (corresponding to an epidemic outbreak) exists if, and only if, R sn 0 > 1. Finally, as in the case of Theorem 2.2, the following result can be proven. Proof. We prove the first assertion. Assume, by contradiction, thatz 1 ≤z 2 . Then, Hence, contradicting the assumption that H(z) = 0.
This Theorem states that if there exists a bias in heterosexual transmission, such bias is always reflected in the expression for the final attack ratio, even when there is another form of direct transmission that is independent of sex. However, we wish to understand how the ratio of final attack ratios z f /z m changes, depending on the relative importance of non-sexual (direct) vs. (hetero)sexual disease transmission. The relative weight of direct transmission is henceforth defined as It can be seen that R sn 0 > R n 0 + R s 0 . Thus, ρ cannot be exactly considered as the contribution of direct transmission to the overall R sn 0 , but is simply based on the proportions of reproduction numbers relative to the two modes. Fig. 2 shows the effect of the differences in male-female susceptibility in heterosexual transmission on the ratio of the sex-specific attack ratios, depending on whether direct transmission also occurs. It can be seen from this figure that with only heterosexual transmission, the final attack ratios can be rather different for the two sexes, although not as much as the differences in susceptibility. For instance, if females were 10 times more susceptible than males, the final attack ratio in females would less than 3 times higher than the final attack ratio in males. If the direct transmission route is responsible for a sizeable portion of all transmission, then the difference becomes quite smaller. For instance, if direct and heterosexual   (1) including heterosexual transmission only and R 0 is given by (6); the red curves using the model (7) that includes both types of transmission, where R 0 is given by (19) and ρ defined in (21) equal to 50%. Here, for the sake of simplicity, we have set γ m = γ f . transmission are equal in transmission potential, the ratio of final attack ratios is always smaller than 3, even assuming a 50-fold difference in susceptibility.
An alternative way to assess the joint impact of ρ (relative importance of direct transmission) and β f /β m (ratio of susceptibilities) can be seen in Fig. 3, where a contour plot of the ratioz f /z m , as a function of ρ and β f /β m , is depicted.
In summary, the analyses in this section show that the final attack ratios can be obtained as the unique solution of H(z) = 0, when the reproduction number R sn 0 , defined in (19), exceeds unity. The ratio between female and male attack ratios always exceeds 1 when male-to-female transmission is higher than female-to-male transmission. Its numerical value depends, as illustrated in Figures 2 and 3, on three parameters: q (the ratio between female and male susceptiblity), R 0 (the basic reproduction number) and ρ (the relative importance of direct vs. heterosexual transmission). In principle, this allows for a retrospective estimate of ρ from the sex-ratio of reported cases, if q and ρ are known, assuming that the model is approximately correct, and that there is no bias in reporting.

The model including vector-borne transmission.
We assume now that infection can also be transmitted by effective contact with the disease vector (adult female Aedes mosquitoes in the case of ZIKV). We modify the system (1) by adding terms corresponding to infection from mosquito bites, as well as an equation for the dynamics of the infected mosquito population whose size we shall denote by V I = V I (t). For mathematical tractability, we shall consider a constant total mosquito population of size V ; we are implicitly assuming a source of mosquitoes to balance the loss due to natural mortality. We shall not model that source explicitly, since it is not needed here. The mortality rate of adult mosquitos is denoted by µ V ; the effective transmission rate of infection per human from mosquito to human, β V , is the product of the average number of bites a human receives per unit time and the probability that the bite from an infected mosquito results in the infection of the susceptible human. Similarly, the effective transmission rate of infection per mosquito from an infected human, β H , is the product of the average number of bites a mosquito makes per unit time and the probability that the bite from an uninfected mosquito to an infected human actually infects the mosquito. For the model (22), it is not possible to obtain analytical expressions for the final attack ratios (similar to those obtained for the previous models). Hence, we resort to estimating the attack ratios numerically. The following parameter values, reported in the ZIKV study in Latin America by Ferguson et al. [9], are used: β H = 0.035, β V = 0.05, γ f = γ m = 1/6 and µ V = 0.2 using 1 day as the time unit. It is worth mentioning that Ferguson et al. [9] also give estimates for R 0 varying in time, but with an average slightly above 2. Furthermore, their model does not consider sexual transmission, so that R 0 can be obtained (as in Ross-Mac Donald [16,21]) as Model (22) includes infection transmission both through mosquito bites and heterosexual contacts. To obtain the associated reproduction number, the infected compartments in the Equation (22) (that is, the second, fourth and fifth equations), will be used. At the disease-free equilibrium, (0, 0, 0), the Jacobian of the aforementioned system, containing just the infected compartments, is where we have used S f = S m = N . The associated next generation matrix is given by , also having assumed, for the sake of simplicity, γ m = γ f .
Hence, it follows that, the basic reproduction number, R sv 0 , is the spectral radius of K. Computing the determinant of K − λI, one sees that An analytical solution to this equation is very cumbersome, so that R sv 0 has to be computed numerically. Moreover, the relation between R 0 and the final attack ratios is different for models (22) and (7).
Thus, in order to compare the results obtained using model (7) (and illustrated in Figures 2 and 3) to those obtained with model (22), we concentrated on final attack ratios rather than on R 0 (given by either (19) or (24)). Model (22)  female/male attack ratio q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q direct transmission; β f /β m =5 direct transmission; β f /β m =10 direct transmission; β f /β m =20 vector transmission; β f /β m =5 vector transmission; β f /β m =10 vector transmission; β f /β m =20 Figure 4. Ratio of sex-specific final attack ratios,z f /z m , as a function of average attack ratio (z f +z m )/2 for models (7) (lines) and (22) for different values of q (see legend). Here ρ = 1/2 and parameters are varied to keep q and ρ at these values.
More specifically, we varied β m and adjusted the other parameters so as to keep q and ρ at the set values. This is compared to what is obtained by solving (10) for the same values of q and ρ (given by (21) in that case). We remark that the ratio of final attack ratios obtained by simulating (22) depend on the values chosen for b = β f /γ, c = β m /γ and ad = (R v 0 ) 2 , but do not depend on individual parameters while keeping those quantities fixed. It should be noted that the results from model (22) show a very similar pattern to those obtained from model (7) for each value of q, both qualitatively and quantitatively. Thus, the results obtained in the analysis of the model with direct transmission are quite useful also in providing insight into the model with vector-based transmission.
Finally, we note that examining the transient phase of the epidemic may yield a different picture. Indeed, initially, the sex-ratio in new cases is closer to the value q representing the ratio of susceptibilities, but then drops because fewer susceptibles are left in that class. In Fig. 5, we show one simulation in time, illustrating this phenomenon. Even if the sex-ratio of final attack ratios is 1.25, the sex-ratio in new cases during the growing phase of the epidemic is around 2.

5.
Conclusions. Three Kermack-Mckendrick type deterministic models have been developed and used to gain insight into the transmission dynamics of a disease within a two-sex closed population. The first model, which considers heterosexual contact as the sole mode of transmission, enables us to theoretically assess whether the frequently observed disparity in transmission rates per contact between maleto-female and female-to-male (usually the former being quite larger than the latter) necessarily leads to similar differences in the final attack ratios between the genders. The attack ratios can be obtained as solutions of a system of two nonlinear equations; we proved that the system has a unique solution if the net reproduction number exceeds unity, extending the classical result by Kermack and McKendrick [13] for an SIR model in a closed homogeneous popualtion. It is further shown, as intuitively expected, that the gender-specific final attack-ratios are biased in the same direction as the gender-specific susceptibilities. We have also shown that the ratio of attack ratios depends solely on the ratio of gender-specific susceptibilities and on the basic reproductive number of the epidemic, R 0 .
The second model extends the first by also allowing for infection transmission through direct contact in a gender-independent way. Here, too, we derived analytical expressions to obtain those attack ratios and proved that the attack ratios are non-zero if the net reproduction number exceeds unity; furthermore, the gender specific final attack ratios in that case are biased in the same direction as the gender-specific susceptibilities. Through some numerical examples, we show how the ratio of the gender-specific final attack ratios varies, depending on the relative importance of the two transmission pathways and on R 0 (Figure 2). In principle, the relative importance ρ of the two transmission pathways could be ascertained retrospectively, if we knew the ratio of gender-specific susceptibilities, the ratio of final attack ratios and the value of R 0 (see Figure 3) and, of course, if this simple model were correct. Qualitatively, one can say that the ratio of final attack ratios will not exceed 1.5 if non-sexual transmission accounts for most transmission events (ρ ≥ 0.6) and the ratio of gender-specific susceptibilities is not too large (say, 5 at most).
The third model, which also has two modes of transmission (namely, heterosexual contact and through a vector (mosquito bites)), can be considered a model for ZIKV that is very simplified in several respects; for instance, we assumed that the total vector population is constant, which implicitly requires a source for new unifected mosquitos to replace the ones that die. We derived, following standard methods [7], an expression for the net reproduction number as solution of a cubic equation, but were not able to find a simple equation that allows to obtain the final attack ratios. We showed, via numerical simulations, that results are quite similar to those obtained with second model, that allows for heterosexual and direct transmission. Thus, it seems that Figure 3 can be used to estimate the sex-ratio of final attack ratios for a given value of q (ratio of male-to-female vs. female-to-male transmission probability) and ρ (relative weight of vector-borne vs. heterosexual transmission).
It should, finally, be remarked that transient patterns can differ from final attack ratios. Indeed, it is quite intuitive that, if one sex succumbs to most of the infections during the early stages (as perhaps [5] was the case with females during the 2015 ZIKV outbreaks in the Americas) due to their larger susceptibility, there will be fewer susceptibles of that sex later in the epidemic and this will bias new infections towards the opposite sex. Clearly, this effect is partiuclarly strong when final attack ratios are relatively large, so that depletion of susceptibles is relevant. An example of such a pattern is provided by the simulation shown in Fig.5. Thus, data showing no gender-bias in final attack ratios (as in the analysis of Zika epidemics in Zap island [8]) are not necessarily in contradiction with gender-bias in reported new cases during the expanding phase of an epidemic (as suspected by some for Brazil [5]), even without assuming virus evolution, or behavioral differences in the populations.
We think that these models, while certainly overly simplistic, can help in assessing the potential role of gender-specific transmission in infections with multiple modes of transmission, by gauging what can be expected to be seen from reports of new cases, or seroprevalence surveys. Certainly, investigations of specific cases of sexual transmission remain essential for obtaining empirical information about the role of sexual transmission for Zika [1].