Impact of Heterogeneity in Sexual Behavior on Effectiveness in Reducing HIV Transmission with Test-and-Treat Strategy

Abstract

The WHO's early-release guideline for antiretroviral treatment (ART) of HIV infection based on a recent trial conducted in 34 countries recommends starting treatment immediately upon an HIV diagnosis. Therefore, the test-and-treat strategy may become more widely used in an effort to scale up HIV treatment and curb further transmission. Here we examine behavioural determinants of HIV transmission and how heterogeneity in sexual behaviour influences the outcomes of this strategy. Using a deterministic model, we perform a systematic investigation into the effects of various mixing patterns in a population of men who have sex with men (MSM), stratified by partner change rates, on the elimination threshold and endemic HIV prevalence. We find that both the level of overdispersion in the distribution of the number of sexual partners and mixing between population subgroups have a large influence on endemic prevalence before introduction of ART and on possible long term effectiveness of ART. Increasing heterogeneity in risk behavior may lead to lower endemic prevalence levels, but requires higher coverage levels of ART for elimination. Elimination is only feasible for populations with a rather low degree of assortativeness of mixing and requires treatment coverage of almost 80% if rates of testing and treatment uptake by all population subgroups are equal. In this case, for fully assortative mixing and 80% coverage endemic prevalence is reduced by 57%. In the presence of heterogeneity in ART uptake, elimination is easier to achieve when the subpopulation with highest risk behavior is tested and treated more often than the rest of the population, and vice versa when it is less. The developed framework can be used to extract information on behavioral heterogeneity from existing data which is otherwise hard to determine from population surveys.

Abstract

The WHO's early-release guideline for antiretroviral treatment (ART) of HIV infection based on a recent trial conducted in 34 countries recommends starting treatment immediately upon an HIV diagnosis. Therefore, the test-and-treat strategy may become more widely used in an effort to scale up HIV treatment and curb further transmission. Here we examine behavioural determinants of HIV transmission and how heterogeneity in sexual behaviour influences the outcomes of this strategy. Using a deterministic model, we perform a systematic investigation into the effects of various mixing patterns in a population of men who have sex with men (MSM), stratified by partner change rates, on the elimination threshold and endemic HIV prevalence. We find that both the level of overdispersion in the distribution of the number of sexual partners and mixing between population subgroups have a large influence on endemic prevalence before introduction of ART and on possible long term effectiveness of ART. Increasing heterogeneity in risk behavior may lead to lower endemic prevalence levels, but requires higher coverage levels of ART for elimination. Elimination is only feasible for populations with a rather low degree of assortativeness of mixing and requires treatment coverage of almost 80% if rates of testing and treatment uptake by all population subgroups are equal. In this case, for fully assortative mixing and 80% coverage endemic prevalence is reduced by 57%. In the presence of heterogeneity in ART uptake, elimination is easier to achieve when the subpopulation with highest risk behavior is tested and treated more often than the rest of the population, and vice versa when it is less. The developed framework can be used to extract information on behavioral heterogeneity from existing data which is otherwise hard to determine from population surveys.

Author Summary

HIV is endemic in populations of MSM in Western countries. As ART reduces transmission risk, increased testing and treatment rates are expected to lower HIV incidence. However, concerns are that in MSM populations changing risk behavior may counteract the impact of ART on transmission. Using a mathematical model, we investigated how heterogeneity in sexual behavior influences the possible effects of a test-and-treat strategy on HIV prevalence and in particular the prospects of eliminating HIV from these populations. We demonstrated that behavioral heterogeneity plays an important role in determining the impact of ART on reducing HIV transmission. Knowledge of behavioral heterogeneity is key in setting intervention goals in populations of MSM.

Citation: Rozhnova G, van der Loeff MFS, Heijne JCM, Kretzschmar ME (2016) Impact of Heterogeneity in Sexual Behavior on Effectiveness in Reducing HIV Transmission with Test-and-Treat Strategy. PLoS Comput Biol 12(8): e1005012. doi:10.1371/journal.pcbi.1005012

Editor: Katia Koelle, Duke University, UNITED STATES

Received: November 18, 2015; Accepted: June 2, 2016; Published: August 1, 2016

Copyright: © 2016 Rozhnova et al. This is an open access article distributed under the terms of the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original author and source are credited.

Data Availability: All relevant data are within the paper and its Supporting Information files.

Funding: GR and MEK gratefully acknowledge funding by the Aids Fonds, grant number 2013030. The funders had no role in study design, data collection and analysis, decision to publish, or preparation of the manuscript.

Competing interests: The authors have declared that no competing interests exist.

Introduction

Recently, a large trial conducted at various sites in 34 countries provided evidence that starting ART as soon as possible regardless of CD4 count is advantageous for health prospects of HIV infected persons [1]. The WHO's early-release guideline for ART initiation now reflects these findings recommending to start treatment immediately upon an HIV diagnosis [2]. Therefore, the test-and-treat strategy, where a population is tested for HIV regularly and those found positive are treated immediately, may become widely used in countries with a generalized HIV epidemic. Earlier, it was investigated whether and under which circumstances a test-and-treat strategy and a more general strategy of treatment as prevention would eventually lead to elimination of HIV from a population [3-8]. While much discussion has been devoted to the possible influence of high infectiousness during primary HIV infection on expected effects of large scale ART on HIV incidence [9-12], less attention has been directed to behavioural determinants of HIV transmission dynamics and how heterogeneity in sexual behaviour will influence the impact of a test-and-treat strategy on HIV transmission. Models of HIV treatment as prevention already included heterogeneity (e.g. [7]) but there has been no systematic investigation of how results depended on assumptions about it.

Sexual behaviour influences HIV transmission dynamics in various ways. Changes in sexual risk behaviour over time have been observed, first decreasing risk behaviour as a response to the emerging HIV epidemic in the 1980's, and later increasing risk when ART became available at the end of the 1990's. These changes have been especially apparent in populations of MSM [13-16]. More recently, modelling studies showed the impact of changes in risk behaviour of individuals over time on HIV transmission dynamics [17-27].

Here we developed a modeling approach which allows a systematic investigation into the effects of various mixing patterns in populations stratified by rates of partner change on the basic reproduction number, treatment effects and prospects of elimination. We investigated how endemic levels and elimination threshold depend on the level of overdispersion in the distribution of numbers of partners and on a mixing parameter that reflects assortativeness of mixing. We studied how the infection is distributed across population strata in endemic steady state and how this changes with various levels of diagnosis and treatment. We chose baseline parameter values to reflect the HIV epidemic among MSM in Western countries.

Materials and Methods

Model formulation

We considered an extended version of a previous model used to evaluate prospects of elimination of HIV with test-and-treat strategy [5]. Here, we explicitly incorporated risk heterogeneity in sexual activity and mixing between population groups by sexual activity.

The model represents a population of MSM of size N(t) stratified into m risk groups of size Nl(t) with partner change rates cl, l = 1, ..., m. Here is the time-dependent population size that changes due to additional mortality from HIV infection. Individuals remain in their risk group. The population in group l consists of susceptible, Sl(t), infected, Ilk(t), and treated, Alk(t), individuals in stage of infection k = 1, ..., n, where n is the number of disease stages. The population size in group l can then be expressed as .

In each risk group the model describes HIV infection process, disease progression through n stages of infection, birth, background and HIV related mortalities, the uptake of and dropping out of ART (Fig 1). Individuals enter the risk class l at rate as susceptible, where is the initial number of individuals in group l. Susceptibles can become infected with the first stage of HIV infection at rate Jl(t) (force of infection). In the absence of treatment infected individuals progress through n stages of infection with varying durations and infectivities ending with death from HIV. Infected individuals in any stage can be screened and start ART at rate [tau]. When we will consider heterogeneous ART uptake by risk group, we will denote [tau]l the uptake by risk group l. Treated individuals progress through n stages of infection with varying durations too, albeit at smaller rates. The rates of progression from stage k to stage k + 1 for untreated and treated individuals are denoted as [rho]k and [gamma]k, respectively. Treated individuals in stage k revert to an infection of stage k at rate [varphi]. This transition represents leaving the virally suppressed state; this can be due to treatment failure, dropping out of treatment or other reasons. In the following, we will refer to [varphi] as drop out rate. Finally, all classes of individuals are subject to background mortality at rate M.

Fig 1. Schematic diagram of the HIV model.

The model describes HIV infection process, disease progression through n stages of infection, birth, background mortality and HIV related mortality, the uptake of and dropping out of ART. The model assumes that the population is divided into the classes of susceptibles, Sl(t), infected, Ilk(t), and treated, Alk(t), in n stages of infection, k = 1, ..., n. The population is stratified into m risk groups indicated by the label l = 1, ..., m in the subscript, and the diagram describes the dynamics in one of them. The risk groups differ only in their partner change rates and initial group sizes, whilst the disease progression and dropping out of treatment are the same for all risk groups. The diagram shows the situation when treatment uptake is the same by all risk groups. We denote [tau]l the uptake by risk group l when we consider heterogeneous ART uptake by risk group. The interaction between the groups is encoded in the time-dependent force of infection Jl(t) that takes into account mixing between the risk groups.

http://dx.doi.org/10.1371/journal.pcbi.1005012.g001

The model was formulated as a system of differential equations for the number of individuals in different classes as follows (1) (2) (3) (4) (5) where k = 2, ..., n and l = 1, ..., m. For heterogeneous uptake by risk group [tau] has be substituted by [tau]l in Eqs (1)-(5).

The time-dependent force of infection (per year) in risk group l is given by Here [Lambda] is the transmission probability per partnership, [epsilon] is the infectivity for an individual on ART, whereas hk describe the infectivity in stage k of infection. Infectivity is defined as a dimensionless quantity describing the relative contribution of each disease or treatment stage to overall infectiousness.

The m × m mixing matrix M(t) = [Mll'(t)]l,l' ∈ {1,...,m} with the elements Mll'(t) denoting mixing of susceptible in the risk group l with infected in risk group l' is defined as follows (6) where [delta]ll' = 1 if l = l' and [delta]ll' = 0 otherwise.

The mixing parameter 0 ≤ [omega] ≤ 1 describes the degree of assortative mixing by risk level. When [omega] = 0 mixing between the risk groups is fully assortative (like with like), when [omega] = 1 mixing is fully proportionate. Eq (6) means that a proportion (1 - [omega]) of the partnerships are formed only with the individuals of the same risk group l = l', whereas the remaining proportion [omega] of the partnerships is formed with each risk group (l' = 1, ..., m) proportionally to the number of partnerships offered by those risk groups. Mixing between groups was random meaning that we did not incorporate preferential mixing for adjacent risk groups. This method of incorporating mixing between different population subgroups is commonly used in sexually transmitted infections (STI) models [11, 28-31].

The proportion of new susceptible individuals entering each risk group was chosen such that, in the absence of HIV, this proportion would remain constant. In the model the total population size, N(t), however, as well as the population sizes of m risk groups, Nl(t), l = 1, ..., m, are not constant because of additional mortality from HIV infection. Note that the burden of HIV due to HIV related mortality is different per risk group as the groups with the highest number of infected individuals will have more HIV related deaths. S1 and S2 Figs show the time-dependent behavior of the model variables for the default parameter values without and with ART, respectively.

The mathematical model was implemented in Mathematica 9.0.

Model parametrization

For parameterizing the model we chose values that are plausible for describing populations of MSM in Western countries, but it was not our aim to fit the model to a specific population. We used data to choose the order of magnitude for parameters. Estimates for parameter values relating to heterogeneity in sexual behavior were obtained from Rutgers World Population Foundation (WPF) data on MSM collected in 2005-2006 in the Netherlands and previous studies of STI dynamics among MSM in the UK. Parameters relating to disease progression and infectivity were extracted from published literature. The parameters for the model and their baseline values are summarized in S1 Table.

We assumed that the rate of recruitment to the sexually active population, M, equals the death rate. The average duration of sexual activity is 1/M = 45 years [32-34]. We made use of Rutgers WPF data on the number of MSM in the Netherlands as the baseline value for the total population size in the beginning of HIV epidemic, N0 [35, 36]. The model can accommodate any number of risk groups, m. Here we focused on the case m = 6 considered previously in modeling dynamics of Hepatitis B virus in MSM populations in the UK and the Netherlands [33, 34, 36-38]. From these studies we adopted the initial fractions of the population in the 6 risk groups, ql, where ql ≤ 1 for l = 1, ..., 6 and . We calculated the initial numbers of individuals in each risk group, , from the relation .

Any number of HIV stages, n, can be incorporated into the model. Following Refs. [3, 5, 39] we parametrized the model for the case n = 4, because for this choice of n estimates for the rates of transition between infection stages for untreated, [rho]k, and treated, [gamma]k, individuals, as well as for the infectivities of untreated, hk (all parameters for k = 1, 2, 3, 4), and treated, [epsilon], individuals were available. It should be noted that these infectivities were estimated for heterosexual couples but we used them for MSM in the absence of similar estimates for this population. For n = 4, infection stages are primary infection, asymptomatic chronic stage, the last two stages together define the symptomatic AIDS stage which is subdivided into an infectious and a noninfectious period (due to severe illness leading to cessation of sexual activity). In the model the stages of the population under treatment have no biological interpretation. They were chosen in our previous work (Ref. [5]) such that the survival probability has a distribution function that agrees with CASCADE data from the time period after introduction of ART.

The rate of treatment uptake, [tau], and the rate of dropping out from treatment, [varphi], can be varied in the model. We present the results in terms of annual treatment uptake and dropout percentages, [tau]* and [varphi]*, respectively. These were computed from the expression for the probability that an event (treatment or dropping out) takes place within one year as percentage = (1 - e-rate×1 year)100%. In case of heterogeneous uptake by risk group, the uptake percentage by group l, , is computed using the same expression. Unless stated otherwise, the annual dropout percentage was fixed throughout the analysis at 5%.

The mixing parameter quantifying the degree of assortativeness, [omega], is a variable parameter of the model. Data on the mixing between different risk groups are hard to obtain because information about the characteristics of sexual partners is required. In STI modeling studies, including HIV, the value of [omega] has been either assumed to have a certain value [29, 30] or estimated by fitting a model to incidence data using Bayesian analysis [11, 31]. The direct estimates of mixing were obtained in three studies of sexual behavior based on contract tracing among patients of STI clinics in the USA and in three studies based on a survey of the general population in the USA, France and the UK [40-42]. These estimates indicate weak like-with-like mixing. However, they do not generalize to the USA MSM population, and, to our knowledge, no more data on mixing is available for this population in the USA or other countries. In our analysis, [omega] is a free parameter that takes on the whole range of possible values, [omega] ∈ [0, 1].

We fixed [Lambda] at 5% such that HIV prevalence [43] and R0 [5] in our study represent a plausible range of values that is compatible with HIV dynamics in MSM in Western countries.

Estimation of partner change rates

To estimate partner change rates, cl, l = 1, ..., 6, we used Rutgers WPF sexual behavior data for MSM population in the Netherlands [35]. In the Rutgers WPF survey, respondents were recruited via an existing internet panel. Since there were few MSM in the panel, MSM were additionally recruited via banners on websites that are frequently visited by MSM. This resulted in the final list of respondents. Numbers of respondents with certain demographic characteristics (sex, age, education, residence) were matched with the distribution in the population of the Netherlands, such that the survey population was representative in these demographic variables. Younger age groups (15-45 years) and certain ethnic minorities were oversampled. Therefore weighting factors were then used by Rutgers WPF to come to a representative sample with respect to the above variables. From this survey population, we extracted male respondents who reported sex with men. We used data on the self-reported number of partners in the last 6 months including information about a steady partner in the last 6 months and condom use for 176 respondents aged 15 to 70 years. The question respondents answered was: "How often did you use condoms in the last 6 months with casual partners with whom you had anal sex?" A similar question was posed for steady partners. As we estimated the number of new partners, the number of partners was corrected (minus 1) if a person reported a steady partner, and the duration of relationship with this partner was longer than 6 months. Condom use was encoded in a binary variable (0 = always, 1 = not always), where all respondents who reported condomless sex in the last 6 months were grouped into one category. This variable was used as a multiplication factor for the number of partners, so that individuals who always used condoms effectively had 0 partners.

We fitted the probability density function for a Weibull distribution to the WPF data histogram by maximum likelihood method. The resulting Weibull distribution is a continuous probability distribution with two parameters, a shape parameter [alpha] = 0.5 and a scale parameter [Beta] = 1.26, both of which were obtained from the fitting procedure. Fig 2A shows the probability density function and S3 Fig shows the corresponding cumulative distribution function. From this distribution we obtained the mean rate of partner change, c = 2.54 partners per year, and mean rates of partner change, cl, per intervals defined by the initial fractions of the population in the 6 risk groups, ql, l = 1, ..., 6. Fitting a Gamma distribution with two parameters resulted in similar estimates of the mean partner change rate and of the mean rates in different risk groups, see S4 Fig for details.

Fig 2. Distributions in the partner change rate used in the analysis have the same mean rate estimated from the data and different variances.

(A) Probability density function for the Weibull distribution in the partner change rate fitted to the WPF data histogram by maximum likelihood method. The dashed lines indicate the intervals defined by the initial fractions of the population in the 6 risk groups, ql, l = 1, ..., 6, per which mean rates of partner change were estimated. For a better visualization the ranges of the x and y axes differ among the panel A and the panel B that is why not all dashed lines can be seen, see S3 Fig for more detail. (B) Weibull distributions with the same mean of c = 2.54 partners per year and different variances, [sigma]2, were obtained by varying the shape parameter, [alpha] , and the scale parameter, [Beta]. The Weibull distribution that best fits to the data is shown in black.

http://dx.doi.org/10.1371/journal.pcbi.1005012.g002

Following standard theory, the basic reproduction number for the SIR model with demography and constant population size in a population stratified by partner change rates and with proportionate mixing is proportional to [sigma]2/c, where [sigma]2 and c are the variance and the mean of the distribution in the partner change rate [44]. To study the impact of heterogeneity in partner change rates on the dynamics of the model, we fixed the mean rate and varied [alpha] and [Beta] to obtain Weibull distributions with different variances, [sigma]2. The mean rates of partner change for each of the distributions, cl, were then computed as before (S2 Table). The variance of the distributions in our analysis ranged from 1.8 to 63.5 yr-2, see Fig 2B for probability density functions and S3 Fig for the corresponding cumulative distribution functions. The estimate of the variance obtained from the Weibull distribution best-fitting to the data was [sigma]2 = 32.6 yr-2.

Computation of the effective reproduction number

We analyzed the threshold behavior of the model using the next-generation matrix approach. An extensive discussion of this approach has already been given at length in the literature in the context of compartmental epidemic models, and we refer the reader to [45, 46] for formal details. In S1 Text we described the aspects that are important for this paper. In essence, the method allowed us to compute a parameter, known as the basic reproduction number, R0, from Eqs (1)-(5). Following standard theory [44], if R0 > 1 the infection will reach an endemic equilibrium as t [right arrow] ∞, whereas if R0 < 1 the infection cannot spread in a population. Therefore, R0 is a threshold parameter for the model. We distinguish R0 from Re: the effective reproduction number in a population with treatment. Here we use the term "effective reproduction number" to determine the threshold below which HIV cannot persist in a population under treatment. This is different from the reproduction number in the transient phase of the epidemic. The latter is influenced by density dependent effects, the former is not. Before ART was introduced, the HIV epidemic in populations with persistent HIV transmission was characterized by R0 > 1. As we will see, treatment lowers R0 and thus Re ≤ R0. Elimination by treatment occurs if Re < 1. The computation of both quantities is similar, but in terms of interpretation it is more clear if we distinguish these two, see S1 Text for the details.

Lorenz curves

To describe the distribution of infections across risk groups for populations with different levels of heterogeneity and mixing we used a method based on the so-called Lorenz curve. This method was shown to be useful for calibrating STI models that include risk structure [47]. The Lorenz curve is a graphical representation of a cumulative probability distribution, namely it represents the cumulative proportion of HIV infected individuals as a function of the cumulative proportion of the population in different risk groups ranked in the order of their partner change rate. In the model the population sizes of different risk groups change [over time] as a result of the differential burden of HIV in each risk group, so in principle we could use the cumulative proportion either of the initial or of the final population as the x-axis. We checked that Lorenz curves were not affected much by this choice, therefore we used the initial population fractions in our analysis (S1 Table). The skewness in the distribution of HIV infections across the risk groups is measured by the deviation of the Lorenz curve from the diagonal line. The diagonal denotes the symmetric situation, i. e. the situation where every risk group has the same HIV prevalence.

Results

We present results regarding the elimination threshold, as given by the effective reproduction number Re, the endemic prevalence and the Lorenz curves. First, we analyzed the model without treatment for a range of distributions in the partner change rate with the same mean rate and different variances. Then, we introduced ART (and consequently dropping out from ART too) into the model and focused on the distribution of partner change rate fitted to the Rutgers WPF data. The prevalence was computed in the steady state. In the case when there is treatment, prevalence is computed as the percentage of all infected individuals including those on ART.

Analysis without ART

Fig 3 shows the basic reproduction number R0 as a function of the mixing (assortativeness) parameter [omega] in the model without treatment. We observe that R0 has a strong dependence both on [omega] and on the heterogeneity of the population as quantified by the variance in the rate of partner change, [sigma]2. Our model predicts that R0 is below 1 for populations with a low variance and low levels of assortativeness even in the absence of treatment. For the lowest variance used in the analysis, [sigma]2 of 1.8 yr-2, HIV cannot persist in the population for any level of mixing.

Fig 3. The impact of mixing on the basic reproduction number in the model without treatment.

The results are for populations with different variances in partner change rates, [sigma]2, and the mean rate of partner change kept constant. Mixing is proportionate (assortative) for [omega] = 1 ([omega] = 0). The dashed line indicates the threshold value of R0 = 1 below which HIV cannot spread in the population. R0 decreases as the mixing becomes proportionate and the variance gets lower. For low variances and high levels of mixing (i.e. low levels of assortativeness) R0 can be smaller than 1 even in the absence of treatment. For the lowest variance we used in the analysis HIV cannot persist for any level of mixing (the yellow curve).

http://dx.doi.org/10.1371/journal.pcbi.1005012.g003

Fig 4 demonstrates the impact of mixing on prevalence. For a fixed variance, prevalence does not necessarily change monotonically as mixing ranges from proportionate to intermediate to assortative ([omega] ranges from 1 to 0), see Fig 4A. For example, for [sigma]2 = 63.5 yr-2, the total prevalence is highest for proportionate mixing (blue bars) but for a lower variance, [sigma]2 = 32.6 yr-2, prevalence is highest for intermediate levels of mixing (black bars). This nontrivial effect occurs because population sizes of subgroups are not constant due to HIV related mortality (see S5 Fig). For populations with high variance, prevalence is reduced by more than half (from 2.5% to 1.2%) as the mixing changes from proportionate to assortative. The assortativeness parameter quantifies the extent to which different risk groups of the population are coupled. As [omega] decreases HIV is able to persist in a lower number of risk groups but prevalence per risk group gradually gets higher (Fig 4B). The increase in prevalence with an increasing partner change rate is not unexpected and corroborates the concept of a 'core' group.

Fig 4. The impact of mixing on endemic prevalence in the model without treatment.

(A) Total prevalence in the population. The color of the bars corresponds to the variance in the rate of partner change, [sigma]2. Note that prevalence is zero for the parameters where HIV is not able to spread (yellow curves for all levels of mixing and the green curve for proportionate mixing in Fig 3). (B) Prevalence per risk group for a population with a variance in the rate of partner change [sigma]2 of 32.6 yr-2 corresponding to the data. The color of the bars denotes the number of the sexual activity class.

http://dx.doi.org/10.1371/journal.pcbi.1005012.g004

In Fig 5 we plot the Lorenz curves that represent the cumulative proportion of infected individuals as a function of the cumulative proportion of the initial population when the risk groups are ranked in the order of their average number of partners per year. The diagonal line represents the situation in which every risk group would have the same HIV prevalence. The Lorenz curves deviate significantly from the diagonal, indicating that the infection is concentrated in the groups with the highest numbers of partners. This skewness in the distribution of HIV is more pronounced for higher assortativeness of mixing as seen from the comparison of the solid, dashed and dot-dashed curves in the plot.

Fig 5. Lorenz curves in the model without treatment.

The results are for populations with different variances in partner change rates, [sigma]2, and the mean rate of partner change kept constant. The diagonal line represents the situation in which every risk group would have the same HIV prevalence. Lorenz curves deviate from it which means that the distributions of infection across the risk groups for proportionate, intermediate and assortative mixing are skewed with high prevalence in small high risk groups. This effect gets stronger as the mixing becomes more assortative. Note that we did not plot the results for the parameters used in Fig 3 for which HIV is not able to spread (yellow curves for all levels of mixing and the green curve for proportionate mixing).

http://dx.doi.org/10.1371/journal.pcbi.1005012.g005

Impact of ART on HIV dynamics

Treatment is able to decrease the basic reproduction number and eliminate HIV if Re gets below 1. In Fig 6A Re is shown as a function of annual treatment uptake [tau]* and mixing. The model predicts that for proportionate mixing elimination is feasible in populations with an annual treatment uptake above 30%. The range of [omega] where Re < 1 gets wider with increasing [tau]*. Nonetheless, it is not feasible to eliminate HIV from populations with high degree of assortativeness without additional intervention measures even if treatment uptake is as high as 90% annually if treatment uptake is the same in all risk groups.

Fig 6. Effective reproduction number and treatment coverage for different annual treatment uptakes and dropout percentages.

(A) The impact of annual treatment uptake, [tau]*, on the effective reproduction number for a population with the variance in the rate of partner change corresponding to the data. Mixing is proportionate (assortative) for [omega] = 1 ([omega] = 0). The curve plotted for [tau]* = 0.0 is, by definition, R0. The dashed line indicates the threshold value of Re = 1 below which HIV is eliminated from the population. Re decreases with increasing treatment uptake and mixing. (B) Treatment coverage as a function of annual treatment uptake [tau]* for different dropout percentages [varphi]* and the remaining parameters as in A.

http://dx.doi.org/10.1371/journal.pcbi.1005012.g006

To translate these findings into results on treatment coverage in the population required for HIV elimination, we show in Fig 6B the coverage as a function of annual treatment uptake for different dropout percentages. Here, the coverage is defined as percentage of infected individuals who are on treatment in the steady state, --a measure that can be obtained from HIV data on diagnosis and treatment. Note that in our model the infected population includes those who are unaware of their infection. Our results indicate that annual treatment uptake of more than 30% required for elimination corresponds to a coverage of almost 80% if 5% drop out from ART due to treatment failure or other reasons annually. For mixing with a higher degree of assortativeness, treatment coverage has to be even higher. This 80% coverage is in line with the UNAIDS 90/90/90 treatment target according to which 90% of all people living with HIV will know their HIV status, 90% of all people with diagnosed HIV infection will receive sustained ART and 90% of all people receiving ART will have viral suppression by 2020 [48]. The first two objectives lead to a coverage of 90% × 90% = 81% meaning that it may be possible to reach elimination for a realization of these objectives in relatively homogeneous populations, but not in populations with strong heterogeneity in sexual behavior and mixing.

In Fig 7 we show the percentage reduction in the total HIV prevalence due to treatment. The reduction is 100% for proportionate mixing and annual ART uptake above 60%, meaning that for these parameters we can achieve elimination. For other types of mixing patterns elimination is not feasible but the reduction in prevalence is still significant if treatment uptake is sufficiently high. The reduction is 38% and 57% for [tau]* = 30% and intermediate and assortative mixing, respectively, and even higher for higher values of [tau]*. In some cases, however, treatment can have even a slight adverse effect on prevalence, as for [tau]* = 10% and assortative mixing. This happens because this treatment uptake is not sufficient to decrease HIV transmission substantially when different risk groups do not interact whilst the average lifespan of individuals on ART, and thus the total number of infected individuals, increase. In the model with treatment we again find skewness in the distribution of infections among risk groups that gets more pronounced with decreasing [omega], see Lorenz curves in S6 Fig. ART uptake has only modest impact on this distribution, with its shape almost entirely defined by the type of mixing pattern.

Fig 7. Percentage reduction in the total prevalence due to treatment.

The variance in the rate of partner change equals that of the data. The color of the bars corresponds to the annual treatment uptake, [tau]*.

http://dx.doi.org/10.1371/journal.pcbi.1005012.g007

Heterogeneous uptake of testing and treatment

The national survey data on HIV testing and risk behavior in Britain shows that voluntary confidential HIV testing by men is significantly associated with reporting greater numbers of same-sex partners [50]. We thus investigated how heterogeneity in uptake of testing and treatment possibly affect our conclusions regarding the levels of ART necessary for elimination. In Fig 8 we considered progressively higher uptake rates by groups with higher numbers of partners. Specifically, we assumed that uptakes by the lowest and highest risk groups, and , were 10% and 90%, respectively, and uptakes by the remaining 4 groups were equally spaced and increasing from 26% to 74% (indicated in the figure legends). In this case, elimination was possible for populations with values of mixing parameter above 0.8, i. e. for populations with mixing closer to proportionate. Heterogeneous testing and treatment offers much better prospects for HIV elimination than a constant treatment uptake rate of 23.25% computed as an average of uptakes by different groups weighted by their population size, (Fig 8). For this level of uptake elimination was not feasible at all. However, very high treatment uptakes in the highest risk groups amount to even higher treatment coverages in those groups, which would require intense screening programmes. The results for other combinations of treatment uptakes by different risk groups are shown in S7 Fig. There we show that Re has values above 1 (elimination is unfeasible) in a wider range of mixing parameter when ART uptake by highest risk individuals is smaller than by the rest of the population, and vice versa if they are tested and get treated more frequently.

Fig 8. Effective reproduction number for heterogeneous uptake of testing and treatment.

denotes the uptake by risk group l. We considered progressively higher uptake rates by groups with higher numbers of partners. Specifically, we assumed that (lowest risk group), (highest risk group), and uptakes by the remaining 4 groups were equally spaced and increasing from 26% to 74%. The dashed line is R0 before ART. Also shown is Re for a constant treatment uptake rate of 23.25% computed as an average of uptakes by different groups weighted by their population size, . Heterogeneous test-and-treat with increasing treatment uptakes by higher risk groups has a larger impact on Re than homogeneous test-and-treat with a constant average uptake. See also S7 Fig with results for other combinations of treatment uptakes by different risk groups.

http://dx.doi.org/10.1371/journal.pcbi.1005012.g008

Sensitivity analysis

Ratio of primary to chronic infectivity.

The analysis presented so far implied that the ratio of primary to chronic infectivity, h1/h2, is according to estimates obtained by Hollingsworth et al [39], namely the primary phase is about 26 times as infectious as the chronic stage (h1 = 2.76 and h2 = 0.106, S1 Table). A recent study indicates that this ratio may be as low as 5 (h1 = 0.62 and h2 = 0.12, [49]) and that the duration of primary phase is 1.7 months instead of about 3 months. Using lower values for h1/h2 leads to more optimistic model predictions for the prospects of elimination. For the sensitivity analysis we shifted infectivity from primary to chronic infection while retaining a constant total infectivity. For a ratio of h1/h2 = 4.65 with h1 = 0.828 and h2 = 0.178 elimination can be achieved for smaller values of treatment uptake percentage and in a wider range of mixing parameter ([omega] ∈ [0.85, 1] and [omega] ∈ [0.58, 1] for [tau]* of 30% and 90%, respectively; see Fig 9). If we further shortened the duration of primary phase, the results would get even more favorable for elimination prospects. However, for a given R0 lowering h1/h2 would imply an increase in [Lambda]. Had we taken h1/h2 = 4.65 in our main analysis the value of [Lambda] would have had to be higher to correspond to a plausible range of R0. Taking both effects together (lowering ratio h1/h2 and increasing [Lambda]) would lead us back to similar results as described above, as there is a trade-off between the effects of these two parameters [11].

Fig 9. Effective reproduction number as a function of mixing and treatment uptake for low and high ratios of primary to chronic infectivity.

The solid lines are repeated from Fig 6A and correspond to the baseline value of the ratio h1/h2 = 26.04. The dashed lines were obtained for h1/h2 = 4.65 by shifting infectivity from primary to chronic infection while retaining a constant total infectivity. Elimination can be achieved for smaller values of treatment uptake percentage and in a wider range of mixing parameter if the infectivity of the primary phase is lower.

http://dx.doi.org/10.1371/journal.pcbi.1005012.g009

Treatment uptake during primary infection.

So far we have assumed that treatment uptake is independent of infection stage. In Fig 10 we relax this assumption by setting treatment uptake rate during the primary infection to 0 for all risk groups. As expected, elimination is more difficult to achieve in this case, especially if treatment uptake rate is high (compare the dashed and solid lines in Fig 10).

Fig 10. Effective reproduction number as a function of mixing and treatment uptake with and without uptake during the primary infection.

The solid lines are repeated from Fig 6A and correspond to the annual treatment uptake percentage, [tau]*. The dashed lines were obtained when treatment uptake during the primary infection was set to 0 for all risk groups. As expected, elimination is more difficult to achieve in this case and the impact of not treating primary infection increases with increasing treatment uptake.

http://dx.doi.org/10.1371/journal.pcbi.1005012.g010

Population stratification.

In our analysis the number of risk groups, m, and the respective initial fractions, ql, where l = 1, ..., m, can be chosen arbitrarily. This choice then determines estimates of cl, l = 1, ..., m. To illustrate the effect of population stratification on the effective reproduction number we compared our baseline parameter values to the extreme case when all initial population fractions are equal (ql = 1/6 for l = 1, ..., m = 6), see S8 Fig. The equal stratification leads to the situation closer to the homogeneous population because it smoothens out differences in partner change rates in the 6 risk groups. In this case Re is below 1 for a treatment uptake above 10% independently of the mixing pattern.

Discussion

We investigated how heterogeneity in sexual behaviour impacts on model predictions concerning the effects of ART on endemic HIV prevalence and on the prospects of eliminating HIV from a population. Heterogeneity in the model depended on two parameters, namely the variance in the rate of partner change and the mixing between subpopulations with different risk levels. This allowed us to compare populations that have the same average partner change rates, but differ in the way partnerships are distributed in the population.

We found that both parameters had a large influence on the basic reproduction number and endemic prevalence before ART. HIV would not have been able to spread in populations with proportionate mixing and a low level of overdispersion in the distribution of numbers of partners. For realistic MSM populations, where some degree of assortativeness is always present, R0 is above 1 and is higher if high risk individuals preferably mix with other high risk individuals. The distribution of infection across risk groups is skewed with high prevalence in small high risk subgroups. Moreover, this effect gets more pronounced as assortativeness of mixing increases.

The range of variances and any level of assortativeness in the model can reflect the range of MSM sexual behaviours found in Western societies. We used an average partner change rate estimated from MSM sexual behavior survey in the Netherlands. A different value for this parameter estimated from another data set (e.g. UK NATSAL data [51]) would lead to slightly different quantitative predictions but all qualitative conclusions for the dependence of the basic reproduction number on mixing and overdispersion in the distribution of numbers of partners would remain unchanged.

In the model, ART uptake is able to decrease the effective reproduction number below 1 and lead to HIV elimination. For some optimistic scenarios we found that an annual treatment uptake of at least 30% by all risk groups is necessary to eliminate HIV from populations with proportionate mixing. This uptake translates into a treatment coverage of at least 80% of all HIV infected individuals which is in line with the UNAIDS 90/90/90 treatment target to be reached by 2020. Thus we demonstrate that it may be possible to reach elimination for a realization of these objectives in relatively homogeneous populations regardless of heterogeneity in uptake of test-and-treat by risk group, but not in populations with strong heterogeneity in sexual behavior and mixing. For other types of mixing patterns which are more realistic even higher levels of coverage are necessary. For subpopulations with strongly assortative mixing, the model predicts that elimination with test-and-treat strategy is not feasible and additional interventions reducing the number of sexual partners and/or promoting condom use and PrEP uptake have to be applied. These conclusions agree with those of Dodd et al [52] who showed that test-and-treat in a hyper-endemic African setting generates a smaller impact in a population with heterogeneous risk distribution and assortative mixing than in that with random mixing assuming the intervention is implemented in the same way in both populations. In our model this happens because a high risk core group with a lot of within group mixing will enable persistent transmission within this small group. However, in the presence of heterogeneity in ART uptake, elimination will be easier to achieve when the subpopulation with highest risk behavior is tested and treated more often than the rest of the population. HIV elimination will be easier to achieve as well if the infectivity of primary phase is lower and its duration is shorter as was proposed by Bellan et al [49]. When HIV cannot be fully eliminated, the reduction in endemic prevalence will be still significant. Even if a population is almost fully assortative, we expect the reduction to be of about 57% for an annual treatment uptake of 30% by all risk groups and baseline infectivities.

For many populations we have some knowledge of rates of partner change, or at least numbers of partners reported in a given time period, but usually we have much less information on mixing patterns. Nevertheless, we are interested in using mathematical models based on available data for projecting effects of interventions into the future. We therefore need to be aware of the strong influence of heterogeneity on model outcomes. While it is probably not realistic to gather detailed information on sexual network structure for large populations, our modelling approach offers other ways of extracting information on behavioural heterogeneity from existing data. Linking the impact of behavioural heterogeneity with epidemic outcome distributions in a Lorenz curve allows estimation of the parameters which control heterogeneity by fitting the model to the data based on Lorenz curve. This requires collecting data with individual linkage of sexual behaviour and infection status, as was collected for chlamydia infection in the large UK NATSAL studies [51]. Earlier, comparing such data with outcomes for several individual based models was used to compare the ability of different models to correctly reproduce underlying sexual behaviour networks from population level parameters [47].

As with all models, also our approach has limitations. Our model is deterministic and thus it does not take stochastic effects into account. While stochastic fluctuations can play a role in real populations, where one superspreading individual can have large influence on transmission dynamics, here we were interested in mean effects that can occur over a long time period. The model is simplistic in some aspects of diagnosis, testing and ART uptake. In particular, there is no distinction between HIV-infected individuals who diagnosed and undiagnosed. First, second and third line treatments are not incorporated in the model explicitly. The parameters of the model are based on self-reported sexual behavior which might not be a true reflection. Sexual risk behaviour is stratified into a constant number of levels, and individuals remain in the same strata during their life time. That changes of risk behaviour of individuals in various phases of their lives can be important for HIV dynamics has been highlighted in recent work by Alam et al; Henry et al [19-22, 24, 25]. Also, our model does not take partnership duration into account, and therefore does not allow for long term concurrent relationships, which have been debated as a possible driver of HIV transmission in sub Saharan African heterosexual populations [53-55]. Therefore, our model is more amenable for describing HIV epidemics in MSM populations where concurrent partnerships are less influential for HIV transmission dynamics [56]. The data sets were used as an example to choose plausible parameter values, but we did not attempt to formally fit the model to a comprehensive set of available data. A more data driven approach to analyzing the HIV epidemic among MSM under ART in the Netherlands, for whom the degree of mixing has not been measured directly, is a focus of our ongoing work.

To conclude, we developed a modeling approach to investigate the impact of various mixing patterns in a population stratified by rates of partner change on the basic reproduction number, treatment effects and prospects of elimination. Our analysis revealed that both the variance in the rate of partner change and mixing between subpopulations with different risk levels have a large influence on endemic prevalence before introduction of ART and on possible long term effectiveness of ART. The developed framework offers a way of extracting information on behavioral heterogeneity from existing data, particularly assortativeness of a population, which would be otherwise very hard to measure in a population survey. Such information on behavioural heterogeneity should be taken into account when setting intervention goals and for analysis of cost-effectiveness of test-and-treat programmes in populations of MSM.

Supporting Information

S1 Text. Computation of the effective reproduction number.

(PDF)

S1 Fig. Time-dependent behavior of the model variables for the default parameter values without ART.

(PDF)

S2 Fig. Time-dependent behavior of the model variables for the default parameter values with ART.

(PDF)

S3 Fig.

(A) Cumulative distribution function for the Weibull distribution in the partner change rate and the empirical cumulative frequency distribution obtained from WPF sexual behavior data for MSM in the Netherlands. (B) Distributions in the partner change rate used in the analysis have the same mean rate estimated from the WPF data and different variances. (C) Probability density function for the Weibull distribution in the partner change rate fitted to the WPF data histogram by maximum likelihood method.

(PDF)

S4 Fig. Cumulative distribution function for the Weibull (A) and Gamma (B) distributions in the partner change rate and the empirical distribution obtained from the WPF data.

(PDF)

S5 Fig. Reduction in population size due to AIDS-related mortality.

(PDF)

S6 Fig. Dependence of the Lorenz curves on treatment uptake and mixing.

(PDF)

S7 Fig. Re for heterogeneous uptake of testing and treatment.

(PDF)

S8 Fig. Effect of the population stratification on the effective reproduction number.

(A) Stratification used throughout the main text. (B) Equal stratification into the 6 risk groups. (C) Effective reproduction number as a function of mixing and treatment uptake for the population stratification used in A and in B.

(PDF)

S1 Table. Description of the parameters of the model and their baseline values.

(PDF)

S2 Table. The variance, the shape and the scale parameters of the distributions used in the analysis and the respective partner change rates in the 6 risk groups.

(PDF)

Acknowledgments

We thank Ard van Sighem, Daniela Bezemer, Roel Coutinho, Annemarie Wensing, Amy Matser, Ineke Stolte, Daniela De Angelis, Paul Birrell, and Nicola Low for helpful discussions and their comments on the paper.

Author Contributions

Conceived and designed the experiments: GR MEK. Performed the experiments: GR. Analyzed the data: GR. Contributed reagents/materials/analysis tools: MFSvdL JCMH. Wrote the paper: GR MFSvdL JCMH MEK.

References

The INSIGHT START Study Group. Initiation of antiretroviral therapy in early asymptomatic HIV infection. New England Journal of Medicine. 2015;373:795-807. doi: 10.1056/NEJMoa1506816. pmid:26192873

The World Health Organization. Guideline on when to start antiretroviral therapy and on pre-exposure prophylaxis for HIV. September 2015. http://www.who.int/hiv/pub/guidelines/earlyrelease-arv/en/

Granich RM, Gilks CF, Dye C, De Cock KM, Williams BG. Universal voluntary HIV testing with immediate antiretroviral therapy as a strategy for elimination of HIV transmission: a mathematical model. Lancet. 2009;373:48-57. doi: 10.1016/S0140-6736(08)61697-9. pmid:19038438

Kretzschmar ME, Schim van der Loeff MF, Coutinho RA. Elimination of HIV by test and treat: a phantom of wishful thinking? AIDS. 2012;26:247-248. doi: 10.1097/QAD.0b013e32834e1592. pmid:22179229

Kretzschmar ME, Schim van der Loeff MF, Birrell PJ, De Angelis D, Coutinho RA. Prospects of elimination of HIV with test-and-treat strategy. Proceedings of the National Academy of Sciences USA. 2013;110:15538-15543. doi: 10.1073/pnas.1301801110

Hontelez JAC, Lurie MN, Bärnighausen T, Bakker R, Baltussen R, Tanser F, et al. Elimination of HIV in South Africa through expanded access to antiretroviral therapy: a model comparison study. PLoS Medicine. 2013;10:e1001534. doi: 10.1371/journal.pmed.1001534. pmid:24167449

Eaton JW, Johnson LF, Salomon JA, Bärnighausen T, Bendavid E et al. HIV treatment as prevention: systematic comparison of mathematical models of the potential impact of antiretroviral therapy on HIV incidence in South Africa. PLoS Medicine. 2012;9:e1001245. doi: 10.1371/journal.pmed.1001245. pmid:22802730

HIV Modelling Consortium Treatment as Prevention Editorial Writing Group. HIV treatment as prevention: models, data, and questions-towards evidence-based decision-making. PLoS Medicine. 2012;9:e1001259. doi: 10.1371/journal.pmed.1001259. pmid:22802739

Powers KA, Ghani AC, Miller WC, Hoffman IF, Pettifor AE, Kamanga G, et al. The role of acute and early HIV infection and implication for transmission prevention strategies in Lilongwe, Malawi: a modeling study. Lancet. 2011;378:256-268. doi: 10.1016/S0140-6736(11)60842-8. pmid:21684591

Cohen MS, Chen YQ, McCauley M, Gamble T, Hosseinipour MC, et al. Prevention of HIV-1 infection with early antiretroviral therapy. The New England Journal of Medicine. 2011;365:493-505. doi: 10.1056/NEJMoa1105243. pmid:21767103

Eaton JW, Hallett TB. Why the proportion of transmission during early-stage HIV infection does not predict the long-term impact of treatment on HIV incidence. Proceedings of the National Academy of Sciences USA. 2014;111:16202-16207. doi: 10.1073/pnas.1323007111

Volz EM, Ionides E, Romero-Severson EO, Brandt MG, Mokotoff E, Koopman JS. HIV-1 transmission during early infection in men who have sex with men: a phylodynamic analysis. PLoS Medicine. 2013;10:e1001568. doi: 10.1371/journal.pmed.1001568. pmid:24339751

van Griensven GJ, de Vroome EM, Goudsmit J, Coutinho RA. Changes in sexual behaviour and the fall in incidence of HIV infection among homosexual men. BMJ. 1989;298:218-221. doi: 10.1136/bmj.298.6668.218. pmid:2493866

Stolte IG, Dukers NH, Geskus RB, Coutinho RA, de Wit JB. Homosexual men change to risky sex when perceiving less threat of HIV/AIDS since availability of highly active antiretroviral therapy: a longitudinal study. AIDS. 2004;18:303-309. pmid:15075549 doi: 10.1097/00002030-200401230-00021

Jansen IA, Geskus RB, Davidovich U, Jurriaans S, Coutinho RA, Prins M, et al. Ongoing HIV-1 transmission among men who have sex with men in Amsterdam: a 25-year prospective cohort study. AIDS. 2011;25:493-501. doi: 10.1097/QAD.0b013e328342fbe9. pmid:21192230

Heijman T, Geskus RB, Davidovich U, Coutinho RA, Prins M, Stolte IG. Less decrease in risk behaviour from pre-HIV to post-HIV seroconversion among MSM in the combination antiretroviral therapy era compared with the pre-combination antiretroviral therapy era. AIDS. 2012;26:489-495. doi: 10.1097/QAD.0b013e32834f9d7c. pmid:22156971

Bezemer D, de Wolf F, Boerlijst MC, van Sighem A, Hollingsworth TD, Prins M, et al. A resurgent HIV-1 epidemic among men who have sex with men in the era of potent antiretroviral therapy. AIDS. 2008;22:1071-1077. doi: 10.1097/QAD.0b013e3282fd167c. pmid:18520351

Bezemer D, de Wolf F, Boerlijst MC, van Sighem A, Hollingsworth TD, Fraser C. 27 years of the HIV epidemic amongst men having sex with men in the Netherlands: an in depth mathematical model-based analysis. Epidemics. 2010;2:66-79. doi: 10.1016/j.epidem.2010.04.001. pmid:21352777

Alam SJ, Romero-Severson E, Kim JH, Emond G, Koopman JS. Dynamic sex roles among men who have sex with men and transmissions from primary HIV infection. Epidemiology. 2010;21:669-675. doi: 10.1097/EDE.0b013e3181e9e901. pmid:20585251

Zhang X, Zhong L, Romero-Severson E, Alam SJ, Henry CJ, Volz EM, et al. Episodic HIV risk behavior can greatly amplify HIV prevalence and the fraction of transmissions from acute HIV infection. Statistical Communications in Infectious Diseases. 2012;4:1041. doi: 10.1515/1948-4690.1041. pmid:24058722

Romero-Severson EO, Alam SJ, Volz E, Koopman J. Acute-stage transmission of HIV: effect of volatile contact rates. Epidemiology. 2013;24:516-521. doi: 10.1097/EDE.0b013e318294802e. pmid:23689754

Alam SJ, Zhang X, Romero-Severson EO, Henry C, Zhong L, Volz EM, et al. Detectable signals of episodic risk effects on acute HIV transmission: strategies for analyzing transmission systems using genetic data. Epidemics. 2013;5:44-55. doi: 10.1016/j.epidem.2012.11.003. pmid:23438430

White PJ, Fox J, Weber J, Fidler S, Ward H. How Many HIV infections may be averted by targeting primary infection in men who have sex with men? Quantification of changes in transmission-risk behavior, using an individual-based model. Journal of Infectious Diseases. 2014;210 (Suppl 2):S594-S599. doi: 10.1093/infdis/jiu470. pmid:25381380

Henry CJ, Koopman JS. Strong influence of behavioral dynamics on the ability of testing and treating HIV to stop transmission. Scientific Reports. 2015;5:9467. doi: 10.1038/srep09467. pmid:25902018

Romero-Severson EO, Volz E, Koopman JS, Leitner T, Ionides EL. Dynamic variation in sexual contact rates in a cohort of HIV-negative gay men. American Journal of Epidemiology. 2015;182:255-262. doi: 10.1093/aje/kwv044. pmid:25995288

van Sighem A, Jansen I, Bezemer D, De Wolf F, Prins M, Stolte I, et al. Increasing sexual risk behaviour among Dutch men who have sex with men: mathematical models versus prospective cohort data. AIDS. 2012;26:1840-1843. doi: 10.1097/QAD.0b013e3283574df9. pmid:22781219

van Sighem A, Vidondo B, Glass TR, Bucher HC, Vernazza P, Gebhardt M, et al. Resurgence of HIV infection among men who have sex with men in Switzerland: mathematical modelling study. PLoS One. 2012;7:e44819. doi: 10.1371/journal.pone.0044819. pmid:23024766

Garnett GP, Anderson RM. Sexually transmitted diseases and sexual behavior: insights from mathematical models. Journal of Infectious Diseases. 1996;174 (Suppl 2):S150-S161. doi: 10.1093/infdis/174.Supplement_2.S150. pmid:8843245

Herzog SA, Heijne JCM, Scott P, Althaus CL, Low N. Direct and indirect effects of screening for Chlamydia trachomatis on the prevention of pelvic inflammatory disease: a mathematical modeling study. Epidemiology. 2013;24:854-862. doi: 10.1097/EDE.0b013e31829e110e. pmid:24076992

Barnabas RV, Laukkanen P, Koskela P, Kontula O, Lehtinen M, Garnett GP. Epidemiology of HPV 16 and cervical cancer in Finland and the potential impact of vaccination: mathematical modelling analyses. PLoS Medicine. 2006;3:e138. doi: 10.1371/journal.pmed.0030138. pmid:16573364

Gomez GB, Borquez A, Caceres CF, Segura ER, Grant RM, Garnett GP, et al. The potential impact of pre-exposure prophylaxis for HIV prevention among men who have sex with men and transwomen in Lima, Peru: a mathematical modelling study. PLoS Medicine. 2012;9:e1001323. doi: 10.1371/journal.pmed.1001323. pmid:23055836

Johnson A, Wadsworth J, Wellings K, Field J. Sexual attitudes and lifestyles. Oxford: Blackwell Scientific Publications; 1994

Williams JR, Nokes DJ, Medley GF, Anderson RM. The transmission dynamics of hepatitis B in the UK: a mathematical model for evaluating costs and effectiveness of immunization programmes. Epidemiology & Infection. 1996;116:71-89. doi: 10.1017/S0950268800058970

Anderson RM, Medley GF, Nokes DJ. Preliminary analyses of the predicted impacts of various vaccination strategies on the transmission of the hepatitis B virus. In: Proceedings of the conference on the control of hepatitis B: the role of prevention in adolescence. London; 1991. p. 95-130

Rutgers Nisso Group. Sexual health in the Netherlands. Delft: Eburon; 2006

Xiridou M, van Houdt R, Hahné S, Coutinho R, van Steenbergen J, Kretzschmar M. Hepatitis B vaccination of men who have sex with men in the Netherlands: should we vaccinate more men, younger men or high-risk men? Sexually Transmitted Infections. 2013;89:666-671. doi: 10.1136/sextrans-2012-050900. pmid:23812794

Kretzschmar M, de Wit GA, Smits LJM, van de Laar MJW. Vaccination against hepatitis B in low endemic countries. Epidemiology & Infection. 2002;128:229-244. doi: 10.1017/s0950268801006562

Kretzschmar M, Mangen MJ, van de Laar M, de Wit A. Model based analysis of hepatitis B vaccination strategies in the Netherlands. Vaccine. 2009;27:1254-1260. doi: 10.1016/j.vaccine.2008.12.009. pmid:19118594

Hollingsworth TD, Anderson RM, Fraser C. HIV-1 transmission, by stage of infection. Journal of Infectious Diseases. 2008;198:687-693. doi: 10.1086/590501. pmid:18662132

Garnett GP, Hughes JP, Anderson RM, Stoner BP, Aral SO, Whittington WL et al. Sexual mixing patterns of patients attending sexually transmitted diseases clinics. Sexually Transmitted Diseases. 1996;23:248-257. pmid:8724517 doi: 10.1097/00007435-199605000-00015

Analyse des Comportements Sexuels en France Investigators. AIDS and sexual behaviour in France. Nature. 1992;360:407-409. doi: 10.1038/360407a0. pmid:1448162

Johnson AM, Wadsworth J, Wellings K, Bradshaw S, Field J. Sexual lifestyles and HIV risk. Nature. 1992;360:410-412. doi: 10.1038/360410a0. pmid:1448163

Beyrer C, Baral SD, van Griensven F, Goodreau SM, Chariyalertsak S, Wirtz AL et al. Global epidemiology of HIV infection in men who have sex with men. Lancet. 2012;380:367-377. doi: 10.1016/S0140-6736(12)60821-6. pmid:22819660

Anderson RM, May RM. Infectious diseases of humans. Oxford: Oxford University Press; 1991

van den Driessche P, Watmough J. Reproduction numbers and sub-threshold endemic equilibria for compartmental models of disease transmission. Mathematical Biosciences. 2002;180:29-48. doi: 10.1016/S0025-5564(02)00108-6. pmid:12387915

Diekmann O, Heesterbeek JA, Roberts MG. The construction of next-generation matrices for compartmental epidemic models. Journal of the Royal Society Interface. 2010;7:873-885. doi: 10.1098/rsif.2009.0386

Althaus CL, Turner KME, Schmid BV, Heijne JCM, Kretzschmar M, Low N. Transmission of Chlamydia trachomatis through sexual partnerships: a comparison between three individual-based models and empirical data. Journal of the Royal Society Interface. 2012;9:136-146. doi: 10.1098/rsif.2011.0131

Joint United Nations Programme on HIV/AIDS (UNAIDS). 90-90-90--An ambitious treatment target to help end the AIDS epidemic; 2014. http://www.unaids.org/en/resources/documents/2014/90-90-90

Bellan SE, Dushoff J, Galvani AP, Meyers LA. Reassessment of HIV-1 acute phase infectivity: accounting for heterogeneity and study design with simulated cohorts. PLoS Medicine. 2015;12:e1001801. doi: 10.1371/journal.pmed.1001801. pmid:25781323

McGarrigle CA, Mercer CH, Fenton KA, Copas AJ, Wellings K, Erens B, et al. Investigating the relationship between HIV testing and risk behaviour in Britain: National survey of sexual attitudes and lifestyles 2000. AIDS. 2005;19:77-84. pmid:15627036 doi: 10.1097/00002030-200501030-00009

Fenton KA, Mercer CH, Johnson AM, Byron CL, McManus S, Erens B, et al. Reported sexually transmitted disease clinic attendance and sexually transmitted infections in Britain: prevalence, risk factors, and proportionate population burden. Journal of Infectious Diseases. 2005;191 (Suppl 1):S127-S138. doi: 10.1086/425286. pmid:15627223

Dodd PJ, Garnett GP, Hallett TB. Examining the promise of HIV elimination by 'test and treat' in hyper-endemic settings. AIDS. 2010;24:729-735. doi: 10.1097/QAD.0b013e32833433fe. pmid:20154580

Kim JH, Riolo RL, Koopman JS. HIV transmission by stage of infection and pattern of sexual partnerships. Epidemiology. 2010;21:676-684. doi: 10.1097/EDE.0b013e3181e6639f. pmid:20571409

Eaton JW, Hallett TB, Garnett GP. Concurrent sexual partnerships and primary HIV infection: a critical interaction. AIDS Behavior. 2011;15:687-692. doi: 10.1007/s10461-010-9787-8. pmid:20890654

Kretzschmar M, Caraël M. Is concurrency driving HIV transmission in sub-Saharan African sexual networks? The significance of sexual partnership typology. AIDS Behavior. 2012;16:1746-1752. doi: 10.1007/s10461-012-0254-6. pmid:22790850

Leung KY, Kretzschmar MEE. Concurrency can drive an HIV epidemic by moving R0 across the epidemic threshold. AIDS. 2015;29:1097-1103. doi: 10.1097/QAD.0000000000000676. pmid:25870987