Quantifying aggregated uncertainty in Plasmodium falciparum malaria prevalence and populations at risk via efficient space-time geostatistical joint simulation

Abstract

Risk maps estimating the spatial distribution of infectious diseases are required to guide public health policy from local to global scales. The advent of model-based geostatistics (MBG) has allowed these maps to be generated in a formal statistical framework, providing robust metrics of map uncertainty that enhances their utility for decision-makers. In many settings, decision-makers require spatially aggregated measures over large regions such as the mean prevalence within a country or administrative region, or national populations living under different levels of risk. Existing MBG mapping approaches provide suitable metrics of local uncertainty--the fidelity of predictions at each mapped pixel--but have not been adapted for measuring uncertainty over large areas, due largely to a series of fundamental computational constraints. Here the authors present a new efficient approximating algorithm that can generate for the first time the necessary joint simulation of prevalence values across the very large prediction spaces needed for global scale mapping. This new approach is implemented in conjunction with an established model for P. falciparum allowing robust estimates of mean prevalence at any specified level of spatial aggregation. The model is used to provide estimates of national populations at risk under three policy-relevant prevalence thresholds, along with accompanying model-based measures of uncertainty. By overcoming previously unchallenged computational barriers, this study illustrates how MBG approaches, already at the forefront of infectious disease mapping, can be extended to provide large-scale aggregate measures appropriate for decision-makers.

Figure 1. Progression of footprint algorithm for efficient joint simulation of the space-time Gaussian random field.

The sequence of schematic diagrams shows the algorithm at six different stages. In this schematic, the prediction space is 25 columns by 25 rows by two months. In each diagram the target column to be predicted is marked in red, pixels already predicted to the left or below the target column are shaded, those yet to be predicted are left white. The ‘footprint’ of conditioning data used in each prediction is shaded blue. In this example the full footprint extent is specified to include, in the target month, seven columns to the left of the target column and, in the preceding month, seven columns to the left, seven to the right and the column directly below the target column. This full extent is thinned to include only every second column and row. Diagram 1 (lower left) shows the algorithm at an early stage: having already simulated values in the first three columns of the first month, the target column being simulated is the fourth from left. The full footprint is truncated and consists of only two columns to the left of the target column. As the algorithm scans across this first month, more columns become available to the left and the footprint grows (diagrams 2 and 3). In diagram 4 the algorithm has moved to the second month, and the footprint can now begin to include simulated pixels from the preceding month. In diagram 5 the full footprint is shown, truncated neither to the left nor right. As the algorithm scans further to the right to complete the second month, the footprint becomes truncated once more, this time by the right-hand margin (diagram 6).

Figure 2. Validation results.

A validation procedure generated many sets of aggregated pixels for which a posterior predictive distribution and point estimate of the set-mean PfPR_2–10 could be compared to the true value. In (A) the error between the point estimate and true value is plotted against the size (number of pixels) of each aggregated set (black dots). Also shown are smoothed moving averages of the mean error (green line) and mean absolute error (red line) in relation to set size. (B) is a coverage plot comparing, for aggregated sets of different sizes, the correspondence between predicted probability thresholds (as provided by the modelled posterior predictive distributions of mean PfPR_2–10 in the validation sets) and actual probability thresholds (defined as the observed proportions of true set means exceeding the predicted threshold values).

Figure 3. Simulated global surfaces of PfPR_2–10.

Examples of five of the 500 realisations of PfPR_2–10 generated via the joint simulation algorithm. Each of these maps represents an equally possible ‘reality’ and the full set of 500 provides a model of the probable prevalence at all locations. Because each map is jointly simulated, pixels within any spatial region can be aggregated together to define a regional mean, and the 500 different versions of that mean across the set of maps provides a model of the uncertainty for that mean value. Simulation is constrained to the global limits of stable transmission.

Figure 4. Simulated posterior predictive distributions of PfPR_2–10 at different spatial scales.

The crucial feature of the jointly simulated realisations of PfPR_2–10 was that they could be aggregated over arbitrary spatial and/or temporal regions to generate posterior predictive distributions of mean PfPR_2–10, constituting appropriate models of regional uncertainty. In this example, such predicted distributions are provided at three different scales, predicting mean PfPR_2–10 for the entire African continent (i); across the nations of Ghana (ii), Democratic Republic of Congo (iv) and Kenya (vi); and across a first administrative level unit of those countries (Ashanti Region, Ghana(iii); Kinshasa Province, DRC (v); Nyanza Province, Kenya (vii)).

Figure 5. Predicted national populations at risk of high stable transmission and associated ranked uncertainty metrics.

Map A shows the estimated population at risk of high stable transmission (PfPR_2–10>40%) for each of the 80 P. falciparum malaria endemic countries considered, based on the mean of each posterior predictive distribution. Map B shows the ranked uncertainty associated with each of these national estimates, quantified using the width of each posterior inter-quartile range, and ordered into quintiles such that countries in quintile one have the largest uncertainty and quintile five the smallest. Map C presents the same uncertainty metric but based on national population proportions in this risk class rather than absolute numbers.