Article Text

## Abstract

**Objectives** To propose a Bayesian approach to uncertainty analysis of sexually transmitted infection (STI) models, which can be used to quantify uncertainty in model assessments of policy options, estimate regional STI prevalence from sentinel surveillance data and make inferences about STI transmission and natural history parameters.

**Methods** Prior distributions are specified to represent uncertainty regarding STI parameters. A likelihood function is defined using a hierarchical approach that takes account of variation between study populations, variation in diagnostic accuracy as well as random binomial variation. The method is illustrated using a model of syphilis, gonorrhoea, chlamydial infection and trichomoniasis in South Africa.

**Results** Model estimates of STI prevalence are in good agreement with observations. Out-of-sample projections and cross-validations also show that the model is reasonably well calibrated. Model predictions of the impact of interventions are subject to significant uncertainty: the predicted reductions in the prevalence of syphilis by 2020, as a result of doubling the rate of health seeking, increasing the proportion of private practitioners using syndromic management protocols and screening all pregnant women for syphilis, are 43% (95% CI 3% to 77%), 9% (95% CI 1% to 19%) and 6% (95% CI 4% to 7%), respectively.

**Conclusions** This study extends uncertainty analysis techniques for fitted HIV/AIDS models to models that are fitted to other STI prevalence data. There is significant uncertainty regarding the relative effectiveness of different STI control strategies. The proposed technique is reasonable for estimating uncertainty in past STI prevalence levels and for projections of future STI prevalence.

- Mathematical model
- sexually transmitted infection
- statistics
- STD surveillance
- uncertainty analysis

## Statistics from Altmetric.com

Mathematical models of sexually transmitted infections (STIs) serve a number of purposes.1 However, the results of mathematical models can only be accepted as an adequate approximation to reality if their key assumptions are supported by empirical data. This requirement is difficult to meet when modelling STIs, as there is much uncertainty regarding the natural history and efficiency of transmission for most STIs, and parameters such as the average duration of untreated infection and the probability of transmission per contact are very difficult to estimate reliably.2 3

Several modellers have attempted to demonstrate the veracity of their models by showing that the model outputs are consistent with observed levels of STI prevalence in the population being modelled.4–7 However, this only partly resolves the problem of parameter uncertainty, as there will usually be many different combinations of parameters that give a similar degree of correspondence to observed prevalence, and these different parameter combinations—although they might be considered equally plausible a priori—do not necessarily produce the same conclusions for users of the model outputs. In addition, it is possible that the true prevalence of the STI in the population may differ from that observed as a result of bias in the sample or imperfect sensitivity and specificity of the test used to detect STIs, or random sampling error. A demonstration of similarity between model estimates and observations does not therefore guarantee that all model outputs are realistic, nor does dissimilarity between model estimates and observations necessarily imply that model outputs are unrealistic.

A more formal approach to dealing with parameter uncertainty is to use uncertainty analysis techniques such as Latin hypercube sampling or Monte Carlo simulation.8–10 This involves specifying probability densities to represent uncertainty regarding key model parameters, randomly sampling from these densities and running the model for each sampled parameter combination in order to estimate the range of uncertainty around the model outputs. When applying such techniques to models of STIs, one might wish to consider uncertainty subject to certain output constraints being met—typically, the constraint is that model estimates of STI prevalence must be roughly consistent with observed STI prevalence. These constraints are usually specified in terms of likelihood functions or ‘sum of squares’ criteria, and several approaches have been proposed in the case of HIV modelling.11–15 However, relatively little work has been done in assessing uncertainty in models of STIs other than HIV. Although it would be possible to apply the techniques developed for HIV to other STI, STIs prevalence data are typically more limited than HIV prevalence data, particularly in developing countries, and this makes it necessary to combine STI data from multiple sources.16 17 In addition, the sensitivity and specificity of the tests used for other STIs are typically much lower than those used for HIV,18 19 and variability in diagnostic accuracy is therefore an important factor to consider when defining the likelihood function.

This paper proposes a Bayesian approach to uncertainty analysis of models of STIs other than HIV. The proposed technique is applied to a model of syphilis, gonorrhoea, chlamydial infection and trichomoniasis in South Africa. It is shown that the method can be used to assess uncertainty with respect to STI model predictions, to estimate STI prevalence at a national level from limited sentinel surveillance data and to make inferences about STI natural history and transmission parameters.

## Methods

The Bayesian approach described below follows three steps: (1) specification of prior distributions to represent uncertainty regarding input parameters; (2) specification of a likelihood function to represent ‘goodness of fit’ to observed levels of STI prevalence, for a given parameter combination and (3) calculation of the posterior distribution. The technique is applied to a deterministic model of STIs in South Africa, which produces estimates of STI prevalence in the South African national population, as well as in various subpopulations (pregnant women, women attending family planning clinics, commercial sex workers and men and women in households). The model assumptions about sexual behaviour and HIV have been have been fixed at the mean values estimated in a previous analysis,20 21 and the STI prevalence data used to define the likelihood function are presented in a review of South African STI prevalence data.22 A summary description of the model and the data is provided in the supplementary appendix (available online only).

### Prior distributions on input parameters

Prior distributions are specified to represent uncertainty around the key model input parameters. These key parameters, represented by parameter vector *ϕ*, include the probability of STI transmission in a single act of unprotected sex with an infected partner, the proportion of cases that become symptomatic, the average length of time spent in different infection states (if the STI is not treated), the average duration of immunity following recovery from infection and the proportion of STI cases that were correctly treated before the adoption of syndromic management protocols. Beta distributions are specified for probabilities and proportions, whereas gamma distributions are specified for average durations, as beta and gamma distributions relate to variables defined on the ranges [0, 1] and [0, ∞), respectively. The means and 95% confidence intervals (CIs) of these prior distributions are shown in table 1 (for gonorrhoea, chlamydial infection and trichomoniasis) and table 2 (for syphilis). The sources on which these prior distributions are based are referenced in the supplementary appendix (available online only).

### Likelihood for observed STI prevalence

To define a likelihood function, it is necessary to allow for differences between observations in terms of sample type (antenatal clinic, family planning clinic, sex worker or household survey), diagnostic type and time. It is also necessary to allow for variation in prevalence between communities, because all observations are specific to particular communities and none of the data are nationally representative. Suppose that *Y _{i}* is the number of individuals testing positive in the

*i*th study. It is assumed that

*Y*is binomially distributed with parameters

_{i}*n*and

_{i}*θ*, where

_{i}*n*is the number of individuals tested in the

_{i}*i*th study, and

*θ*is the prevalence that one would expect to observe in the

_{i}*i*th study. To allow for additional variation between observations due to variation in prevalence between communities and variation in diagnostic accuracy, it is assumed that

*θ*). The number of individuals testing positive in the

_{i}∼beta(α_{i}, β_{i}*i*th study therefore follows a

*beta*-binomial distribution, that is,

Parameters *α _{i}* and

*β*are estimated using the method of moments, by first estimating the mean and variance of

_{i}*θ*. Suppose that

_{i}*M*(

*s*, ϕ) represents the model estimate of the prevalence of the STI in sample type

_{i}, t_{i}*s*(pregnant women, women attending family planning clinics, commercial sex workers, men in households or women in households), in year

_{i}*t*, given a set of input parameters that is denoted by the parameter vector ϕ. If

_{i}*ρ*is the true prevalence of the STI in the

_{i}*i*th study population, then a plausible model for the distribution of

*ρ*might be

_{i}*b*term represents the difference (on the logit scale) between the prevalence in the

_{i}*i*th study population and the prevalence that would be expected in a nationally representative sample. The logit transformation is introduced to normalise the variance of the

*b*terms. As the standard deviation of the study effects,

_{i}*σ*, is unknown, a prior distribution is placed on this parameter. A gamma prior is used for

_{b}*σ*, with its mean and standard deviation based on the standard deviations of the logit-transformed STI prevalence levels before fitting the model to the data. The same mean (0.30) and standard deviation (0.15) are used for all STIs (see tables 1 and 2).

_{b}Equation (2) can be re-expressed as follows:*θ _{i}* is then calculated, taking into account the sensitivity and specificity of the test used in the

*i*

^{th}study,

*Se*and

_{i}*Sp*respectively23:

_{i}The mean and variance of *θ _{i}* are then obtained by noting that

The means and variances of the *Se _{i}* and

*Sp*terms are estimated based on reviews of the sensitivity and specificity of various diagnostic tests,18 24–26 and are shown in the supplementary appendix (available online only). For a given set of input parameters ϕ, a third-order Taylor approximation to

_{i}*ρ*(about the point

_{i}*b*=0) is used to approximate the mean and variance of

_{i}*ρ*. Having obtained the mean and variance of

_{i}*θ*using equations (5) and (6), the

_{i}*α*and

_{i}*β*parameters are calculated using the method of moments, and the beta-binomial likelihood for the

_{i}*i*th observation is then computed.

### Posterior calculation

By Bayes' theorem, the posterior density is proportional to the product of the joint prior density and the likelihood function. As the distribution of this product cannot be evaluated analytically, the posterior distribution is evaluated numerically using sampling importance resampling.27 This was implemented by randomly sampling 20 000 parameter combinations from the joint prior, running the model for each parameter combination and calculating the associated likelihood, and then drawing a resample of 500 parameter combinations from the initial set of 20 000, using the likelihood values as sample weights. This sample of 500 parameter combinations is a sample from the posterior distribution and is used to generate the results presented.

### Validation

To assess the validity of the model predictions and the model estimates, the model is refitted after having omitted the most recent observations and refitted after having randomly omitted half of the observations. The resulting estimates of STI prevalence in 2005 are compared with those obtained in the initial analysis. In addition, the omitted observations are compared with the distributions of estimated prevalence levels, expressing the former as percentiles of the latter and plotting a histogram of the percentiles, similar to the verification rank histogram.28 If the model has good validity, the histogram of percentiles should be approximately uniform. The model fitting procedure is also validated by comparing the model estimates of syphilis prevalence with nationally representative survey results,29 which have been excluded from the likelihood calculation so that they can be used for independent validation.

## Results

The model estimates of trends in STI prevalence in the 15–49 years age group are shown in figure 1, together with observed levels of STI prevalence in ‘low-risk’ groups (antenatal attenders, family planning clinic attenders and men and women in households). Model estimates are generally consistent with observations. However, the 95% CIs around the model estimates are narrower than the ranges of variation in observed prevalence levels, because the former reflect only the uncertainty around the population mean, whereas the latter reflect the variability of the study effects, the variation in diagnostic accuracy, differences in sample type, as well as random binomial variation. Model estimates of STI prevalence in women tend to be similar when comparing pregnant women, family planning clinic attenders and women in households, but are substantially higher in commercial sex workers (see table 12 in the supplementary appendix, available online only).

The observations used in defining the likelihood (figure 1) do not suggest any clear trend in STI prevalence over time, but the model fits a downward trend as a result of the model assumption that condom usage and the use of syndromic management protocols have both increased since the mid-1990s. In the case of syphilis, it is possible to validate this modelled decline in prevalence by comparison with levels of syphilis prevalence measured in nationally representative antenatal surveys.29 Figure 1a shows that the modelled decline in syphilis prevalence in women is consistent with that observed in the antenatal clinic surveys (open diamonds), although the nationally representative antenatal data have not been used in calculating the likelihood function.

Tables 1 and 2 show the posterior means and 95% CIs for the parameters for which prior distributions were specified. In most cases, prior and posterior means are similar, but the posterior distributions tend to have a lower variance, especially in the case of the standard deviation of study effects and the average duration of immunity parameters. For syphilis and gonorrhoea, the posterior mean female-to-male transmission probability is substantially greater than the prior mean, and for trichomoniasis and gonorrhoea, the posterior mean male-to-female transmission probability is greater than the prior mean. The posterior mean of the standard deviation of study effects is substantially greater than the prior mean in the case of trichomoniasis, indicating that there is relatively high inter-study variation in trichomoniasis prevalence.

Figure 2 compares the estimates of STI prevalence in 2005 for three analyses: (1) using all the available data, that is, the same results as presented in figure 1; (2) omitting the most recent data and (3) randomly omitting half of the data. The model estimates in the second and third analyses are not significantly different from those in the first analysis. The distribution of the omitted observations, expressed as percentiles of the model estimates, is also not significantly different from uniform in either analysis (results not shown).

To illustrate the importance of uncertainty analysis to users of STI model outputs, we consider three possible strategies for reducing the prevalence of syphilis: doubling the rate at which individuals with syphilis symptoms seek treatment (from initial rates of 0.23 per week in women and 0.57 per week in men); increasing the proportion of private practitioners using syndromic management protocols (from 45% to 80%); and increasing the proportion of pregnant women who are screened and treated for syphilis (from 60% to 100%). If the syphilis parameters are set at the prior means in table 2 (which produce estimates of syphilis prevalence consistent with the observations shown in figure 1a), the projected reduction in syphilis prevalence by 2020, for the three strategies, is 12%, 4% and 6%, respectively. However, using the 500 parameter combinations in the posterior sample, the estimated average reductions are 43% (95% CI 3% to 77%), 9% (95% CI 1% to 19%) and 6% (95% CI 4% to 7%), respectively. Although improved antenatal screening is the least effective strategy on average, it is the most effective strategy for 7% of the posterior parameter combinations and is more effective than expanding syndromic management in the private health sector for 23% of the posterior parameter combinations. The absolute and relative effectiveness of the different strategies are therefore both subject to significant uncertainty.

## Discussion

This paper extends uncertainty analysis techniques for fitted HIV/AIDS models to models that are fitted to other STI prevalence data. The proposed approach differs from approaches that have been described for HIV models12–14 in two important respects. First, allowance has been made for uncertainty regarding diagnostic accuracy in the definition of the likelihood function; this is necessary in view of the variable sensitivity and specificity of the tests commonly used to detect STIs. Second, data from a number of sources have been used in defining the likelihood function: antenatal clinics, family planning clinics, community household surveys and samples of commercial sex workers. As STI prevalence data may be limited in many developing countries, it is necessary to make use of these multiple data sources when fitting the model. Ideally, the STI model should be capable of estimating different prevalence levels for each of these sample types, in order to control for the biases associated with different sample types, and this is a valuable feature of the model that has been used in this analysis. Nationally representative data could also be incorporated into the analysis, by setting the *b _{i}* term in respect of the national data to zero, by setting the Var(

*ρ*) term in equation (6) to zero, and by making appropriate adjustment to the

_{i}*n*term in equation (1) to take account of the survey design.

_{i}The Bayesian approach that we have proposed has three key uses. First, it can be used to assess uncertainty regarding STI model predictions. For example, it has been shown that the relative effectiveness of different syphilis treatment strategies is subject to significant uncertainty, and considering only a single combination of syphilis parameters that appears to give outputs consistent with survey data may lead to a false sense of confidence in the relative effectiveness of different strategies. Although cost-effectiveness analyses of STI interventions have typically considered uncertainty regarding cost and efficacy parameters,30–33 this analysis shows that uncertainty regarding STI transmission and natural history parameters may also contribute significantly to the overall uncertainty regarding the cost-effectiveness of a particular intervention.

The second use of the proposed Bayesian approach is the estimation of STI prevalence at a national level from limited sentinel surveillance data. The informal validation of the model in the case of syphilis (figure 1a) confirms that model estimates of trends in syphilis prevalence are reasonably consistent with nationally representative data. The more formal validation results suggest that the method is reasonable for estimating past STI prevalence levels and for short-term projections of future STI prevalence (<5 years beyond the most recent observation), but for longer-term projections there are insufficient data to validate the method. Although other studies have estimated STI prevalence at national or regional levels from sentinel surveillance data,16 17 these analyses have not been based on dynamic models of STI transmission and have therefore been forced to assume that STI prevalence levels are stable over time. In addition, these studies have not allowed for variability in diagnostic accuracy. It would be useful to assess whether the method we have proposed produces plausible results when applied in other countries, such as China, where sentinel surveillance data have suggested increases in syphilis prevalence over the past decade.34

The third use of the proposed Bayesian approach is the estimation of STI transmission and natural history parameters, which are notoriously difficult to estimate empirically. Although the posterior distributions are generally not substantially different from the prior distributions, and generally have high variance, there are a few cases in which the STI prevalence data are informative and induce a posterior distribution that differs substantially from the prior. For example, the unusually high observed prevalence of trichomoniasis in South African women (figure 1g) leads to a substantially higher posterior estimate on the male-to-female transmission probability than was assumed a priori, on the basis of a model fitted in other settings.5

A potential criticism of the model that has been used is that it assumes individuals are temporarily immune to re-infection following spontaneous resolution of gonorrhoea and trichomoniasis, contrary to expert opinion that such immunity does not exist.35–37 Attempts were made to fit the model to the gonorrhoea, trichomoniasis and chlamydial infection prevalence data on the assumption of no immunity following recovery, but this produced significantly poorer fits (results not shown). This points to a source of uncertainty not explicitly addressed in this analysis: uncertainty regarding the choice of model structure. Further work is required to evaluate more formally the appropriateness of different model structures, using techniques such as Bayesian model averaging.38

A further limitation of this analysis is that it does not allow for correlation between measurements in the same geographical location. Although it is unrealistic to assume that prevalence measurements in the same location would be independent of one another, regressions of observed STI prevalence levels on various covariates did not reveal residual patterns that differed significantly between locations. This suggests that modifying the method to allow for correlation between measurements in the same location would probably not change the results substantially. A related concern is that sentinel surveillance systems are often not geographically representative. In South Africa, for example, there is a strong surveillance bias towards urban areas and towards the KwaZulu-Natal and Gauteng provinces,22 which could potentially bias model estimates of STI prevalence. There is thus a need to improve the collection of STI prevalence data, with a particular focus on surveillance in regions that have previously been undersampled.

In conclusion, the Bayesian approach to uncertainty analysis of STI models provides a useful means of integrating previous knowledge about STI epidemiology and STI prevalence data, allowing for more reliable estimation of STI prevalence levels and better evaluation of the impact of STI control strategies.

### Key messages

Model estimates of the impact of STI control programmes are subject to significant uncertainty, which needs to be assessed through formal uncertainty analysis.

The uncertainty analysis approach proposed here allows for the estimation of STI prevalence at a national level from sentinel surveillance data obtained using different laboratory methods.

The approach can be used to integrate data from different sample types (antenatal clinics, family planning clinics, households and surveys of sex workers).

The Bayesian approach to uncertainty analysis allows for the integration of previous knowledge about STI epidemiology into model parameter estimation.

## Acknowledgments

The authors would like to thank Debbie Bradshaw for helpful comments on earlier versions of this paper.

## References

## Supplementary materials

## Web Only Data sti.2009.037341

**Files in this Data Supplement:**

## Footnotes

Funding This work was partly funded by a seed grant from the Center for Statistics and the Social Sciences at the University of Washington, and by the National Institute of Child Health and Development through grant no R01 HD054511. The funders did not participate in the study design, analysis and interpretation of data, writing of the paper, or the decision to submit the paper for publication.

Competing interests None.

Provenance and peer review Not commissioned; externally peer reviewed.

## Linked Articles

- Whistlestop tour