Direct size estimates were generated using empiric data collection methods between 2012 and 2014 through Integrated Bio-Behavioral Surveillance Studies (IBBSS) in Zambezi, Ohangwena, Erongo, and Khomas for FSW [19] and in Karas, Oshana, Erongo, and Khomas for MSM (Table S1 in Multimedia Appendix 1) [20]. In the IBBSS, direct estimates were derived using multiple methods of size estimation, including mapping, key informant interviews, unique object identifiers, wisdom of the crowds, literature review, and stakeholder consensus [16,19]. These methods have been described previously [16,19]. Respondent-driven sampling was used to recruit all participants. Briefly, a set of 6-9 seeds from 4 study sites were selected, and each seed used coupons to refer others into the study. Recruitment continued until the desired sample size was reached [19,20]. Eligible FSW were aged 18 years or older and had received monetary payment for sex in the previous 6 months. Eligible MSM were aged 18 years or older and reported engaging in sex with other men in the previous 6 months. Both FSW and MSM were required to have resided in the study region for the last 6 months.

Auxiliary data to inform extrapolated estimates were identified in consultation with implementing partners and key stakeholders of key population programs in Namibia and were accessed through the Data.FI Consortium, a 5-year collaboration funded by the President’s Emergency Fund for AIDS Relief (PEPFAR) through the United States Agency for International Development (USAID) [21]. Specific data types included census data and geographically disaggregated population-level data such as regional literacy rates and HIV prevalence estimates. Auxiliary data sources comprised the Namibia 2011 Census [22], Namibia Population Projections calculated by the Namibia Statistics Agency [23], and the Namibia Population-Based HIV Impact Assessment (NAMPHIA) [24]. Census data and population projections were publicly available, while NAMPHIA data were obtained by request from the Centers for Disease Control and Prevention. Program data sources included indicators collected through PEPFAR’s key populations Data for Accountability, Transparency, and Impact Monitoring (DATIM) efforts [25] as well as indicators reported by the partner-led Key Populations—Strengthening Technical Assistance and Response for sustainable HIV prevention and treatment (KP-STAR) program [26].

Population denominators for direct estimates were extracted from Namibia’s Census Population Projections for the years in which data were collected in order to transform size estimates into standardized proportions. For FSW, this included regional estimates of cisgender women aged between 18 and 49 years; for MSM, this comprised adult cisgender men aged between 18 and 49 years. Age bands for each population were selected to align denominators with the age composition of the direct size estimates (numerators) collected through the initial IBBSS respondent-driven sampling survey.

Auxiliary variables identified for potential inclusion in SAE approaches included population density, employment rates, literacy proportions, and HIV prevalence. Population density, literacy, and employment rates were abstracted from the publicly available Namibia 2011 Census, while HIV prevalence estimates came from NAMPHIA (Tables S2 and S3 in Multimedia Appendix 1 provide comparisons of these characteristics across regions, stratified by sex). Variables were selected for inclusion in analyses based on their availability across on- and off-sample regions and through previously identified associations with key population size in the literature [13,14,27]. Population density and HIV prevalence were ultimately selected, given the variation in these data across regions where direct estimates were available. Program data used to inform priors included counts of the number of individuals who accessed services, reported through both DATIM and KP-STAR.

Each method of direct size estimation produces a single estimate. Given the variation between available direct estimates within regions, we used a consensus estimation approach and a publicly available tool (ie, the consensus estimator tool) to generate a single estimate for each region where direct estimates were collected [28]. This allowed for the weighting of the IBBSS-collected direct estimates both by precision and by quality of study implementation. The methods of direct size estimation were ranked by level of confidence, with successive sampling population size estimation (SS-PSE) having the highest confidence, followed by unique object multiplier, mapping, stakeholder consensus, literature review, key informant interview, and lastly, wisdom of the crowds. This initial ranking scheme was based on the rigor and quality of methods, guided by the World Health Organization’s synthesis of these direct size estimation methods and investigator experience implementing these approaches [16,29-31].

Prior beliefs for the distribution of the derived size estimates were defined based on program data. Quarterly data between fiscal year (FY) 2020 quarter 3 and FY 2021 quarter 2 from KP-STAR and between FY 2016 quarter 4 and FY 2021 quarter 1 from DATIM were used. Program data included the number of FSW or MSM that were reached by programs during each quarter. A 2-tailed, 2-sample t test was performed to identify any differences in the underlying distribution between the KP-STAR and PEPFAR data. The data from PEPFAR and KP-STAR were pooled due to similarity in the underlying distribution, and the minimum, maximum, and median numbers were calculated, in addition to the SD for the pooled numbers from 2016 quarter 4 and 2021 quarter 2. These 4 values, along with the assumed log-normal distribution for population sizes, served as the prior beliefs for the consensus estimator tool in order to derive a single direct estimate for each region . Of note, the data availability limited the accessibility of consistently reported direct estimates across regions, such as missing direct estimates collected through SS-PSE for MSM in Karas and Oshana. The USAID team implemented routine data quality assessments of data from DATM and KP-STAR. This involved a random sample review of data captured in client files in comparison to the implementing partner’s client-level electronic database. A similar comparison was routinely collected at the aggregate level, comparing site-level monthly reports to the quarterly data.

Sensitivity analyses were completed to determine the extent to which adjusting the numerical level of confidence but retaining the rank order affected the weighting scheme and the resulting combined direct estimates. Multiple permutations and values for the confidence levels were assigned to assess the sensitivity of the rankings across the resulting triangulated estimates. The resulting size estimates generated from the consensus estimator were similar across the varying confidence levels as long as rank order was preserved.

We conducted simple imputation, stratified imputation, and regression modeling to determine the regional proportion of female individuals aged between 18 and 49 years and the proportion of male individuals aged between 18 and 49 years that were FSW and cisgender MSM, respectively. The derived direct estimates from the consensus estimator tool were used in the simple imputation and stratified imputation approaches, whereas all the direct estimates were used in the regression approach to maximize the number of input data points.

For the imputation approaches, the consensus-generated direct estimates were divided by the total cisgender female population aged between 18 and 49 years for FSW and the total cisgender male population aged between 18 and 49 years for MSM to derive proportions. Using simple imputation, the mean of these proportions was applied to each of the regions where direct size estimates were not available. In the stratified imputation, HIV prevalence estimates reported by NAMPHIA were used as the stratification variable. Regions were categorized by whether their HIV prevalence fell above or below the median HIV prevalence of all the regions, as reported by NAMPHIA (17.85% among women and 7.45% among men of reproductive age). The mean of the proportions for each respective category was applied to each of the regions where direct size estimates were not available. In sensitivity analyses, simple and stratified imputations were also conducted using all the available direct size estimates to evaluate how the incorporation of program data as prior beliefs in the consensus estimator tool changed the results.

In a supplemental analysis, multivariable Poisson regression models were fit using the full set of available direct estimates to generate size estimate proportions for FSW and MSM. To ensure a maximum number of available inputs were available to fit the regression models, the consensus estimator tool was not used. Generalized estimating equations accounted for clustering by region, and candidate predictors included HIV prevalence, population density, literacy rate, employment, and annual projected population growth, all overall and stratified by sex. Model selection was determined using Akaike information criterion (AIC) and Bayesian information criterion (BIC) values, with the lowest values indicating the best fit.

The resulting estimated proportions from all approaches were applied to the projected 2021 population of female individuals for FSW and male individuals for MSM aged between 15 and 49 years, resulting in the final PSEs.

Imputations were conducted using Microsoft Excel [32], multivariable Poisson regression was performed using Stata (version 16) [33], and maps were constructed using ArcGIS software [34,35].

The process of generating the extrapolated size estimates was an iterative process, requiring collaboration and discussion between researchers and government, programs, and community stakeholders in order to refine and select a final set of estimates. Meetings with implementing partners occurred on multiple occasions, with discussion around how these size estimates could be used programmatically to inform service delivery targets. A brief synopsis of these meetings has been provided in Multimedia Appendix 1.

This study comprised a secondary analysis of deidentified surveillance data and routinely collected programmatic data; patient-level data were not shared with the study team nor included in analyses. The Johns Hopkins Bloomberg School of Public Health Institutional Review Board designated this work as not being a human subjects research study (IRB00007442).

Direct PSEs collected through the IBBSS varied across direct estimation methods and by region for FSW (Table 1). The largest number of FSW was generally seen in Khomas, ranging from 100 (95% CI 50-1700) using key informant interviews to 5240 (95% CI 3373-11,706) using the unique object multiplier method. Khomas also produced the largest PSEs for MSM, ranging from 300 (95% CI 100-1600) using key informant interviews to 2416 (95% CI 850-4000) using stakeholder consensus. After incorporating program data as priors in the consensus estimator tool, the resulting generated PSEs in Khomas were 1480 (95% CI 1099-1799) for FSW and 511 (95% CI 143-958) for MSM. Direct PSEs for each region, including the resulting weighted estimates from the consensus estimator tool, are presented in Table 1.

Extrapolated national estimates for FSW ranged from 4777 (stratified imputation by HIV prevalence) to 13,148 (stratified imputation by population density, preprogrammatic data). For MSM, estimates ranged from 4611 (simple imputation and stratified imputation by HIV prevalence) to 10,171 (stratified imputation by population density). Regional estimates from the imputation approaches are reported in Tables 2 and 3 for FSW and Tables 4 and 5 for MSM. Regional estimates from multivariable Poisson regression models are included in Table S4 in Multimedia Appendix 1.

When program data were incorporated as priors in the imputation models, national extrapolated estimates ranged from 4777 (stratified imputation by HIV prevalence) to 9288 (simple imputation) for FSW. For MSM, these estimates ranged from 2236 (stratified imputation by population density) to 2372 (simple imputation). Regional estimates resulting from the inclusion of programmatic priors are presented by population for each approach in Tables 2-5 and are presented visually in Figures 1 and 2.

After including program data as priors in the simple imputation model for FSW, the estimated proportion of female individuals aged between 15 and 49 years who were FSW increased from 1.9% to 2.8% (Table 3). Among male individuals aged between 15 and 49 years, the proportion of MSM decreased from 1.5% to 0.75% (Table 5).

In models where HIV prevalence informed stratified imputations, the estimated proportion of women who were FSW increased from 2.6% to 4.0% in regions where HIV prevalence was greater than the median; in regions with less than median HIV prevalence, the estimated proportion decreased from 1.4% to 1.2%. Among MSM, the proportion decreased from 2.1% to 0.8% in regions where the HIV prevalence was greater than the median and decreased from 1.3% to 0.7% in regions where the HIV prevalence was less than the median.

In models where population density informed strata, the estimated proportion of women who were FSW increased from 1.1% to 3.4% in regions with a greater than median population density and decreased from 2.8% to 1.0% in regions with less than median population density. Among MSM, the proportion decreased from 1.4% to 0.45% in regions where the population density was greater than the median and decreased from 1.5% to 1.0% in regions where the population density was less than the median.