Spatio-temporal estimates of HIV risk group proportions for adolescent girls and young women across 13 priority countries in sub-Saharan Africa

Data

We analysed nationally-representative household survey data from 13 countries: Botswana, Cameroon, Kenya, Lesotho, Malawi, Mozambique, Namibia, South Africa, Eswatini, Tanzania, Uganda, Zambia and Zimbabwe. We included surveys conducted in these countries between 1999 and 2018 in which women were interviewed about their sexual behaviour and sufficient geographic information was available to locate survey clusters to health districts. Demographic and Health Surveys (DHS) [12], AIDS Indicator Surveys (AIS) [13], Population-based HIV Impact Assessment (PHIA) [14] surveys, and the Botswana AIDS Impact Survey 2013 (BAIS) [15] were included.

For each survey, we classified female respondents aged 15-29 years into one of four behavioural risk groups according to reported sexual risk behaviour in the past 12 months. These risk groups were: not sexually active, one cohabiting sexual partner, non-regular or multiple sexual partner(s), and AGYW who report transactional sex (Table 1). In the case of inconsistent responses, women were categorised according to the highest risk group they fell into, ensuring that the categories were mutually exclusive. Exact survey questions varied slightly across survey types and between survey phases (S2 Text). Questions captured information about whether the respondent had been sexually active in the past twelve months, and if so how with many partners. For their three most recent partners, respondents were also asked about the type of partnership (spouse, cohabiting partner, partner not cohabiting with respondent, friend, sex worker, sex work client, and other).

Some surveys included a specific question asking if the respondent had received or given money or gifts for sex in the past twelve months. In these surveys, 2.64% of women reported transactional sex. In surveys without such a question, women almost never (0.01%) answered that one of their three most recent partners was a sex work client. Due to this incomparability across surveys, we did not include surveys without a specific transactional sex question when estimating the proportion of the population who engaged in transactional sex. We focused on estimating the proportion of women who reported transactional sex at a district level, and subsequently adjusted these proportions to align to national estimates for the number of female sex workers.

We used estimates of population, people living with HIV (PLHIV) and new HIV infections stratified by district and age group from HIV estimates published by UNAIDS that were developed using the Naomi model [11]. The model synthesises data from multiple sources to produce subnational estimates of indicators of interest, and has been used by countries as a part of the HIV estimates process supported by UNAIDS. The administrative area hierarchy and geographic boundaries we used correspond to those used for health service planning by countries, exceptions being Cameroon and Kenya where we conducted analysis one level higher at the department and county levels, respectively (Table E in S2 Text). We used the most recent 2022 estimates for all countries, apart from Mozambique where, due to data accuracy concerns, we used the 2021 estimates (in which the Cabo Delgado province is excluded due to disruption by conflict).

Two-stage model for sexual risk group proportions

To estimate the proportion of AGYW in each risk group, we took a two-stage modelling approach. First, we fit a spatio-temporal multinomial model stratified by district, year (1999-2018) and five-year age group (15-19, 20-24, and 25-29) for the proportion of AGYW in three categories: (1) not sexually active, (2) one cohabiting partner, and (3) either non-regular or multiple partner(s), or transactional sex. Combining the two highest risk groups (high and very high) in this way allowed data from all surveys to be included in this first stage model. Second, we fit a spatial logistic regression model separating those who have non-regular or multiple partner(s) from those who reported transactional sex, stratified by district, and five-year age group, and using only data from the surveys with a specific transactional sex question. As surveys were only available in the years 2013-2018, we assumed the proportion in the very high risk group among those in the two highest risk groups was constant over time. We combined the two models using 1000 samples from each posterior distribution to produce samples for all four risk groups. Finally, we adjusted the samples from the transactional sex category to match age- and country-specific FSW population size estimates. We modified the samples from the non-regular or multiple partner(s) risk group to ensure that after adjustment the risk group proportions still summed to one. FSW population size estimates by age were obtained by disaggregating national 15-49 FSW population size estimates [16] using the FSW age distribution in South Africa from the Thembisa model [17] in combination with country-specific age at first sex distributions [18]. Further technical details are in Section 1 of S1 Text.

We considered four model specifications for the space-age-time multinomial model for the three risk groups. All models included intercepts for each risk group, as well as age, country, and age-country random effects. To account for district-level variation we used spatial random effects consisting of a parameter for each district. We considered alternative model specifications in which the spatial random effects were either independent or spatially correlated such that more information was shared between neighboring districts than those far apart. Similarly, we used temporal random effects to allow variation in risk group proportions over time, and considered alternative model specifications as independent versus first-order auto-regressive where a smooth temporal trend is assumed. To understand the importance of each part of the model we analysed the relative sizes of the variance parameters for each effect.

For the logistic regression model of the proportion engaging in transactional sex among those with non-regular or multiple partner(s), we considered six specifications. Each included an intercept, age and country random effects, and a spatial random effect allowing district-level variation. Both independent and spatially correlated spatial random effects were considered. To improve estimation with sparse data, during model selection we considered alternatives with national-level covariates for either the proportion of men who reported ever having paid for sex or having paid for sex in the last twelve months [19].

We performed inference using the integrated nested Laplace approximation (INLA) [20] algorithm via the R-INLA package [21]. For models with a Gaussian latent field, INLA has comparable accuracy to Markov chain Monte Carlo with realistic, finite samples [22], and is substantially more computationally tractable for high dimensional models like ours, which has 940 districts, 20 years, 3 age groups, and 4 risk groups. It is not possible to directly fit multinomial logistic regression models in R-INLA, so we used the multinomial-Poisson transformation [23]. Details of this approach, including how we used Kronecker products over Gaussian Markov random fields to appropriately define random effects, are provided in Section 1.2 of S1 Text. The best performing model was selected according to the conditional predictive ordinate (CPO) criterion [24], a measure of leave-one-out model performance that can be calculated directly in R-INLA without explicit model refitting [25]. The R [26] code used to implement the models and produce results is available from github.com/athowes/multi-agyw. We used sf [27] for handling of spatial data, orderly [28] for reproducible research, ggplot2 for data visualisation [29] and rticles [30] for reporting via rmarkdown [31].

HIV indicators and prevention prioritisation

Using risk group proportion estimates, and subnational estimates developed using the Naomi model, we calculated HIV prevalence, PLHIV, HIV incidence and number of new HIV infections stratified by district, age group and risk group. Further details are provided in S1 Text.

Data

We included 46 surveys in our analysis (Fig 1, Table C in S2 Text), with a total sample size of 274,970 women aged 15-29 years (103,063 aged 15-19 years, 92,173 aged 20-24 years, and 79,734 aged 25-29 years). Of these, 12 surveys included a specific transactional sex question, with a total sample size of 62,853 (28,753 aged 15-19 years, 26,324 aged 20-24 years, and 7,776 aged 25-29 years¹). The median number of surveys per country was four, ranging from one in Botswana and South Africa to six in Uganda.

• Download figure
• Open in new tab

Fig 1.

Surveys used in our analysis by year, survey type, sample size, and whether the survey included a question about transactional sex. Details of included surveys are in Table C in S2 Text.

Model selection and model fit

The best fitting multinomial regression model included correlated spatial random effects and independent and identically distributed temporal random effects. The best logistic regression model for transactional sex included correlated spatial random effects and the proportion of men who reported ever paying for sex covariate. Model performance according to the CPO criterion is provided in S1 Text for all models we considered. Direct estimates of the four risk group proportions from these surveys were highly correlated with our modelled estimates at a national-level (Fig E - Q in S2 Text).

Risk group estimates

Fig 2 and Fig 3 show posterior mean estimates for the proportion in each risk group from our final model (provided in S1 Data). In subsequent results, all estimates refer to 2018, the most recent year included in our analysis, unless otherwise indicated.

• Download figure
• Open in new tab

Fig 2.

The spatial distribution (posterior mean) of the AGYW risk group proportions in 2018. Estimates are stratified by risk group (columns) and five-year age group (rows). Countries in grey were not included in our analysis.

• Download figure
• Open in new tab

Fig 3.

National (in white) and subnational (in color) posterior means of the risk group proportions. Estimates are stratified by risk group (columns) and five-year age group (rows).

The median national FSW proportion was 1.1% (95% CI 0.4–1.9) for the 15-19 age group, 1.6% (95% CI 0.6–2.8) for the 20-24 age group and 1.9% (95% CI 0.5–3.5) for the 25-29 age group.

In the 20-24 and 25-29 year age groups, the majority of women were either cohabiting or had non-regular or multiple partner(s). Countries in eastern and central Africa (Cameroon, Kenya, Malawi, Mozambique, Tanzania, Uganda, Zambia and Zimbabwe) had a higher proportion of women in these age groups cohabiting (63.1% [95% CI 35–78.7%] compared with 21.3% [95% CI 10.1–48.8%] with non-regular partner[s]). In contrast, countries in southern Africa (Botswana, Eswatini, Lesotho, Namibia and South Africa) had a higher proportion with non-regular or multiple partner(s) (58.9% [95% CI 43.2–70.5%], compared with 23.4% [95% CI 9.7–39.1%] cohabiting). This clear geographic delineation passes along the border of Mozambique, through the interior of Zimbabwe and along the border of Zambia (Fig 2).

In most districts (57.9%; 95% credible interval [CI] 27.7–79.7) adolescent girls aged 15-19 were not sexually active. The exception was Mozambique, where the majority (64.23%) were sexually active in the past year and close to a third (34.17%) were cohabiting with a partner.

Age group was the most important factor explaining variation in risk group proportions, accounting for 65.9% (95% CI 54.1–74.9%) of total variation. The primary change in risk group proportions by age group occurs between the 15-19 age group and 20-29 age group (Fig 3). The next most important factor was location. Country-level differences explained 20.9% (95% CI 11.9–34.5%) of variation, while district-level variation within countries explained 11.3% (95% CI 8.2–15.3%). Temporal changes only explained 0.9% (95% CI 0.6–1.4%) of variation, indicating very little change in risk group proportions over time. We observed similar variance decomposition results fitting each country individually, and using other model specifications.

Expected infections reached

For any given fraction of AGYW prioritised, substantially more new infections were reached by strategies that included behavioural risk stratification. Reaching half of all expected new infections required reaching 19.4% of the population when stratifying by subnational area and age, but only 10.6% when behavioural stratification was included (Fig 4). The majority of this benefit came from reaching FSW, who were 1.3% of the population but 10.6% of all new infections.

• Download figure
• Open in new tab

Fig 4.

Percentage of new infections reached across all 13 countries, taking a variety of risk stratification approaches, against the percentage of at risk population required to be reached.

Considering each country separately, on average, reaching half of new infections in each country required reaching 14.6% (range 8.7-21.8%) of the population when stratifying by area and age, reducing to 5.1% (range 2.1-13.2%) when behaviour was included. The relative importance of stratifying by age, location and behaviour varied between countries, analogous to the varying contribution of each to the total variance (Fig C in S2 Text). For example, FSW in Kenya were estimated to be 1.1% of the population and close to a third (16.1%) of all new infections, whereas FSW in Tanzania were just 1.2% of the population and 16% of all new infections.

Limitations

Our analysis was subject to challenges shared by most approaches to monitoring sexual behaviour in the general population [38]. In particular, under-reporting of higher risk sexual behaviours among AGYW could affect the validity of our risk group proportion estimates. Due to social stigma or disapproval, respondents may be reluctant to report non-marital partners [40] or may bias their reporting of sexual debut [18,41,42]. For guidance of resource allocation, differing rates of under-reporting by country, district, year or age group are particularly concerning to the applicability of our results; and, while it may be reasonable to assume a constant rate over space-time, the same cannot be said for age, where aspects of under-reporting have been shown to decline as respondents age [43], suggesting that the elevated risks we found faced by younger women are likely a conservative estimate. If present, these reporting biases will also have distorted the estimates of infection risk ratios and prevalence ratios we used in our analysis, likely over-attributing risk to higher risk groups.

We have the least confidence in our estimates for the FSW risk group. As well as having the smallest sample sizes, our transactional sex estimates do not overcome the difficulties of sampling hard to reach groups. We inherent any limitations of the national FSW estimates [16] which we adjust our estimates of transactional sex to match. Furthermore, we do not consider seasonal migration patterns, which may particularly affect FSW size. More generally, we did not consider covariates potentially predictive of risk group proportions (such as sociodemographic characteristics, education, local economic activity, cultural and religious norms and attitudes), which are typically difficult to measure spatially. Identifying measurable correlates of risk, or particular settings in which time-concentrated HIV risk occurs, is an important area for further research to improve risk prioritisation and precision HIV programme delivery.

The efficiency of each stratified prevention strategy depends on the ability of programmes to identify and effectively reach those in each strata. Our analysis of new infections potentially averted assumed a “best-case” scenario where AGYW of every strata can be reached perfectly, and should therefore be interpreted as illustrating the potentially obtainable benefits rather than benefits which would be obtained from any specific intervention strategy. In practice, stratified prevention strategies are likely to be substantially less efficient than this best-case scenario. Factors we did not consider include the greater administrative burden of more complex strategies, variation in difficulty or feasibility of reaching individuals in each strata, variation in the range or effectiveness of interventions by strata, and changes in strata membership that may occur during the course of a year. Identifying and reaching behavioural strata may be particularly challenging. Empirical evaluations of behavioural risk screening tools have found only moderate discriminatory ability [9], and risk behaviour may change rapidly among young populations, increasing the challenge to effectively deliver appropriately timed prevention packages. This consideration may motivate selecting risk groups based on easily observable attributes, such as attendance of a particular service or facility, rather than sexual behaviour.