Executive Summary
In this report, we develop and apply spatio-temporal risk models for predicting the risk of invasion of BDBV to new health zones previously not reporting cases. This work was motivated by the need to prioritise packages of interventions in the DRC for health zones at risk.
We used epidemiological data (confirmed cases by symptom onset) up until and including July 13th together with demographic (population counts, socioeconomic deprivation, GDP per capita, health site count and density) and mobility data (Flowminder relocations from Ituri epicenter and travel times between health zones). Aside from confidential line list data, all other sources are available via our open-access data repository and visualised on our real-time situational awareness dashboard. We used this data to build and test spatiotemporal invasion models to anticipate the short-term spatial transmission risk in health zones that have not yet reported a confirmed case of BDBV.
According to 1-week horizon predictions from the best-performing model on 13th July 2026, the following health zones are estimated to be the 10 most at risk of invasion: Rethy, Mahagi, Jiba, Biringi, Linga, Nyarambe, and Angumu (Ituri), Lubero (Nord Kivu), Watsa and Gombari (Haut Uele). These predictions were derived from a model using a generation time distribution with mean 15.3 days and composite Flowminder mobility data with gravity mobility, demographic covariates, and geographic distances. The inclusion of Flowminder mobility data was consistently found to improve model performance.
Summary Figure: A) Relative risk of invasion, 1-week forecast, on a scale of 0-1, plotted on the log scale (grey = already invaded). The most ‘at risk’ (highest probability of invasion) health zone (Rethy) has a score of 1, and the least ‘at risk’ health zone has a score of 0. B) Ranking of the 20 most ‘at risk’ health zones, coloured by province. The mean posterior probability of invasion is marked by a circle, with confidence intervals marked by the horizontal bars.
Our invasion probability outputs for the next two weeks will be updated daily using the latest data to generate operational information for public health decision-making, and displayed for interactive examination on our dashboard.
Results
Briefly, the BDBV outbreak has been geographically concentrated in north-east DRC up until 13th July 2026, particularly in the health zones of Mongbwalu, Rwampara, Bunia, Nyankunde and Nizi, together with short-term outflows from a subset of these health zones (Fig. 1A). The outbreak spread initially among health zones in Ituri before expansion to other provinces (to date: Nord Kivu, Sud Kivu, Haut Uele and Tshopo) (Fig. 1B).
Figure 1. Spatiotemporal spread of BDBV at the health zone level: A) Cumulative confirmed BDBV cases per health zone (binned choropleth; grey = no cases), overlaid with the Flowminder orange arcs representing short-trip outflows from the principal epicentre (Bunia, Rwampara, Mongbwalu). Approximately 88% of the outflow share is to zones within Ituri, with other outflows reaching Kinshasa and the southernmost provinces. B) Epidemic trajectory by symptom-onset date: cumulative confirmed cases (top) and the number of health zones with ≥1 confirmed case (bottom).
Predicted compared to observed invasions
To evaluate the real-time performance of the model, we performed leave-future-out forecasting for one- and two-week forecast horizons. The Bayesian invasion models generally clustered in terms of ability to predict invasion events (Fig. 2A). Based on performance across metrics and horizons, the ‘best’ performing model included the “median” generation time distribution (mean 15.3 days) and composite Flowminder with covariates (population counts, socioeconomic deprivation, distance to nearest affected health zone, and health site density), gravity mobility and geographic distances, producing an AUC-PR skill of 70.2x (90% CI: 45.7, 109.3). 70.2x AUC-PR skill indicates that the model is 70.2 better than chance for anticipating spatial invasions across the five folds. Models including Flowminder mobility data are mostly comparable in performance but including information on previously reported suspected cases in the location of interest and in other locations (weighted by mobility and generation time) did not lead to improved performance (Fig. 2A and Fig. 2C).
Figure 2. Leave-future-out cross-validation. Leave-future-out cross-validation on the at-risk zones only, pooled over 5 folds (17 invasion events). A) Discrimination across the Bayesian model grid: AUC-PR skill (average precision ÷ base rate; 1 = no skill, dashed line; higher = better), at both horizons. Points are the pooled leave-future-out estimate; bars are the 90% zone-cluster bootstrap confidence interval. B) Real-time prioritisation skill: the share of true invasions caught when the top-K highest-risk zones are actively monitored each week, versus a random watch-list of the same size (dashed) and zones monitored just using inflows to locations from epicentre as per the Flowminder data. C) Mean ranking of invaded zones across models for one- and two-week-ahead forecasts of invasion probability. D) Invasion-ranking evolution: the rank trajectory of the current top-25 watch-list across the successive forecast rounds (fold cutoffs) with one line per health zone coloured by province, and ranks worse than 30 are drawn at the “30+” baseline Predicted vs observed invasion per forecast round (fold), faceted by horizon: mean predicted probability (blue), observed invasion fraction (black), and mean predicted probability at the zones that actually invaded (orange). E) Ranking accuracy: the 1-week-ahead predicted probability of invasion among all at-risk zones that actually invaded (orange) or not invaded (grey).
In a scenario where health officials had targeted the top-10 predicted at-risk health zones (sorted by invasion probability) one week ahead of time over the five folds, ~73 % of the next week’s invasions would have been caught versus ~2 % for a random list of the same size (Fig. 2B). The average rank of the newly invaded zones prior to invasion is 9.7 out of a possible ~ 500 health zones. Highly ranked health zones as of 13th July all experienced increases in their rankings as more zones observed their first confirmed case, with zones in Ituri and Nord-Kivu still dominating the most-at-risk locations (Fig. 2D). Zones with the highest predicted probability of invasion were generally associated with subsequent invasions, either in the fold/round in question or the subsequent round/fold (Fig. 2E). While the best-fitting Bayesian model for discrimination performs well at anticipating future spatial invasions, the same model over-predicts the invasion probability out-of-sample by ~2.3×. This is likely in part because the model is trained on the earlier stages of the outbreak where there are fewer control measures and the invasion parameter estimated from such training weeks could be overestimated, as a unit of import pressure converts to a first case more often in the earlier stages around the epicentre. As a result, the ordering/ranking is a better predictor of future invasion than the estimated probability magnitudes. This means that the models are able to detect which zones are most likely to be invaded better than the models can inform us about the probability of the invasion event. We therefore interpret the models’ absolute probabilities as upper bounds and suggest using the ranking and relative-risk scores for informing decisions instead of raw predicted probabilities.
Communicating forecasted invasion probabilities across space and time
We then used the best-fitting Bayesian model from leave-future-out cross-validation to anticipate future spatial invasion risks. We report the relative rankings of health zones by probability of invasion in Fig. 3A, together with the uncertainty in these predictions in Fig. 3B and Fig. 3C.
For the top-10 rankings of health zones, we visualise the posterior median for the invasion probabilities and their corresponding 90% credible intervals (Crls). We report the relative risk of invasion per health zone alongside their relative vulnerability (see Methods) to communicate a ‘priority score’, i.e. a measure of recommended prioritisation of health zones based on their likelihood of, and vulnerability to, invasion (Fig. 3D and Fig. 3E). The relative risk of invasion remained greatest in health zones of Ituri and Nord Kivu, with higher vulnerability in Haut-Uele.
Figure 3. Forecasted invasion probabilities. A) Invasion-risk rank map: each at-risk zone’s relative risk (RR) on a logarithmic scale among at-risk zones (higher =more likely to record its first case). B) Forecast uncertainty: the width of the 90% posterior credible interval on the invasion probability. C) Highest-risk zones (the top by 1-week risk): posterior mean invasion probability with the 90% credible interval, shown for both horizons in one panel (blue = 1 week, red = 2 weeks). D) Invasion-vulnerability scatter (1 week ahead): each zone’s relative invasion risk (y) against its vulnerability-capacity gap (x), with point size proportional to the composite priority score and colour by province; the upper-right corner represents health zones most likely to be invaded and vulnerable/ under-resourced. E) Bivariate invasion-probability × vulnerability choropleth across DRC for 1-week-ahead and 2-week-ahead risks from 13 July 2026.
Discussion & Caveats
We developed and evaluated simple spatio-temporal invasion models based on a renewal framework to predict future risk of BDBV across the DRC. This work aims to inform current and future risk assessment and surveillance prioritisation during a time when the epidemic is still growing (July 14th, 2026). Future work will focus on estimating where cases are growing, declining, and have stabilised. We are working on implementations of a federated framework where others can upload their models into an INRB-hosted environment to then generate spatial probabilistic risk forecasts. Our model could be improved by incorporating testing effort in each location and healthcare-seeking. For now, we have added contextual information by plotting invasion risk against vulnerability scores (Fig. 3E).
Materials and Methods
Data
Epidemiological data: Daily number of cases at the health zone level in the DRC by symptom onset from the DHIS2 linelist. We imputed onset dates when only sample dates were known based on estimates of the delay distribution for time from symptom onset to sample date. Up to June 13th 2026, there were 2161 confirmed cases across 43 affected zones.
There are some discrepancies in confirmed between the linelist and INSP Situational Reports (SitReps); for example, the health zone of Oicha does not have any confirmed cases according to the linelist. However, according to the SitRep, this health zone has already had confirmed cases. For the purposes of this report, we treated the linelist as the ground truth, with three exceptions: we manually included the presence of cases in Oicha, Makiso-Kisangani, and Lubunga. These health zones are densely populated and well connected, thus the presence of cases in these locations (which has been confirmed by the SitReps) will likely have a strong effect on the predicted risk profiles.
Demographic data: We used the following data at the health zone spatiotemporal resolution: population counts (and log-transformed population counts) from WorldPop1, socioeconomic deprivation from the Climate-Conflict Vulnerability Index2, mean GDP per capita, and health site count and density from GRID33. These data were resampled and aggregated using the DART pipeline4.
Mobility data: To build mobility matrices, we used Flowminder data that estimated proportions of subscribers present in the epicentre (Bunia, Mongbwalu or Rwampara over a 21-day reference period)5,6 later detected in other health zones in the DRC (over a 31-day period until May 24th) and paired these data with mobility models (gravity and/or radiation models) that are used for connections between locations that are not captured by mobility flows from the epicentre. We also incorporated zone-to-zone road travel times calculated using the OpenStreetMap Open Source Routing Machine API7.
Table 1: Data streams used in analysis.
| Source | Description | Treatment | Resulting Metrics |
|---|---|---|---|
| WorldPop | Demographic data on population counts and density (2025 projections) | Aggregated from WorldPop raw data by summing raster cells within each health zone in the shapefile to yield population counts. Density calculated by dividing population counts by polygon area. | Population counts per health zone (number persons) |
| Population density per health zone (number persons per sq km) | |||
| Climate-Conflict Vulnerability Index | Socioeconomic data (as of 2025 Q4) | Resampled and aggregated to health zone level via the darts pipeline (mean; sremapbil) | Socioeconomic deprivation per health zone (0-1 index) |
| Kummu et. al | Mean GDP per capita (as of 2022) | Resampled and aggregated to health zone level via the darts pipeline (mean; sremapbil) | Mean GDP per capita per health zone |
| GRID3 | Healthcare site information (as of 19 December 2025) | Point data aggregated (counts) by health zone using HDX shapefile, yielding healthsite count per health zone. Divided by polygon area to yield healthsite density per health zone. | Number of healthcare sites per health zone |
| Healthcare density per health zone (sites per sq km) | |||
| Flowminder population mobility estimates (April 24 - May 24th) | Proportion of an anonymised mobile subscriber cohort - defined by presence in Bunia, Mongbwalu or Rwampara health zones during 3–23 April 2026 - who were later detected at least once in each of the other DRC health zones during 24 April – 24 May 2026. Note that 79 out of 519 DRC health zones are excluded for data quality reasons and carry no detection value (flag = "No data"). | None | Proportion of an anonymised mobile subscriber cohort. |
Spatial invasion model
Bayesian renewal modelling frameworks
Our primary goal was to estimate short-term spatial risks to health zones that had not reported a confirmed case as of the date of forecast generation. To do so, we developed a suite of simple, mobility-informed Bayesian models based on the renewal equation and specifically calibrated for the task of predicting spatial invasion at one- and two-week horizons.
In this framework, we are estimating for each health zone the likelihood of a first confirmed-case onset, cumulative for one and two-week horizons, updated each day for the following week/two weeks. More formally, for each source zone j that has experienced a confirmed case, we formed a generation-time-weighted g(k) sum of nowcasted incidence Yⁿᶜ at time t, giving Ỹ(t):
Ỹⱼ(t) = Σₖ≥₁ g(k) · Yⁿᶜ(j, t − k),
which means that the (confirmed) infection importation pressure into location i at time t can be written as:
Λᵢ(t) = Σⱼ≠ᵢ W[j, i] · Ỹⱼ(t) = (Wᵀ Ỹ)ᵢ.
Then, to model the invasion probability pᵢ for the binary invasion outcome for health zone i (invadedᵢ), we used a Bernoulli distribution with a complementary-log-log link consisting of i) an intercept β₀, ii) optional standardised covariates xₘᵢ (e.g. population, socioeconomic deprivation, health-site density, OSRM travel time to the nearest already-affected zone), and iii) log Λ as a fixed offset:
invadedᵢ ~ Bernoulli(pᵢ),
cloglog(pᵢ) = log(−log(1 − pᵢ)) = β₀ + Σₘ γₘ xₘᵢ + log Λᵢ.
which makes βᵢΛᵢ the cumulative hazard of a first arrival:
pᵢ = 1 − exp(−exp(β₀ + Σₘ γₘ xₘᵢ) · Λᵢ) = 1 − exp(−βᵢ Λᵢ), with βᵢ = exp(β₀ + Σₘ γₘ xₘᵢ).
Where the parameter βᵢ governs how imported infections generate new infections in a new location i. Note that a standard approach would use the reproduction number R(t) (or an inward reproduction number). However, we found that by calibrating a single free parameter β = exp(β₀) directly to the observed, rare invasion frequency, we better estimated the introduction, establishment, and ascertainment mechanisms, compared to using R(t), which resulted in overconfident, poorly calibrated invasion probabilities (results not shown).
For prediction over a one-/two-week horizon, we accumulate the cumulative hazard μ across horizons (so a first invasion is counted once) up to max horizon h, to give invasion probabilities p:
μᵢ,ₕ = βᵢ · Λᵢ(t+h) ,
pᵢ(≤h) = 1 − exp(−Σ_{h′≤h} μᵢ,ₕ′).
For this Bayesian renewal equation framework, we varied assumptions/data across axes of generation time distributions, mobility matrices, covariates, and observation link function (complementary-log-log link or logit link). We fitted these Bayesian renewal models using cmdstanr in R, with 2,000 iterations (1,000 warmup iterations) across 2 chains.
To characterise the vulnerability, we also combined ranks of health zones by health-site density, health-site count, travel time to the nearest zone with a facility, and Climate Conflict Vulnerability Index (CCVI). socioeconomic deprivation. We then produced priority scores by multiplying the invasion probability by the vulnerability and rescaling to get a value between 0 and 1.
As of 13 July 2026, there were no BDBV-specific estimates of the generation time distribution8. We therefore considered three different generation time distributions (short, medium, and long) from EVD outbreaks, which, measured in units of days, were Gamma probability distributions with mean-standard deviation values of (12.0, 6.5), (15.3, 9.3), and (18.0, 10.5), respectively. We discretised the Gamma distribution and binned it to weekly values to align with our weekly modelling horizons.
Key assumptions that could be addressed in future work include the two-stage nowcasting-invasion modelling approach, which could be done jointly. Likewise, we do not correct for within-location growth from undetected transmission chains and do not incorporate information on time-varying human movement patterns. We did see relatively stable movement patterns throughout the observation period (March 2026 compared to May 2026).
Model evaluation
To evaluate and select from the different tested models, we performed rolling origin cross-validation, specifically focusing on the real-time performance of models for anticipating spatial invasion. This meant that we split the data into folds where we trained up to a historical cut-off date, generated predictions for the subsequent two weeks and evaluated whether our models’ predicted invasion probability (of a health zone observing its first case) aligned with the observed invasion frequency.
Across the dataset, at a one-week forecast horizon, there were 2,457 at-risk zone-weeks with 17 invasion events (base rate 0.69%). As spatial invasions are rare events, we have a highly imbalanced dataset for which standard discrimination metrics such as Area Under the Receiver Operating Characteristic curve (AUC-ROC) can be highly misleading for assessing our models’ abilities for anticipating the spatial invasion events. Therefore, we instead used the area under the Precision-Recall curve (AUC-PR), which plots precision against recall for different probability thresholds. AUC-PR was our primary discrimination metric to ask the question: when the model raises an alarm about a spatial invasion, how often is it correct, and how many true spatial invasions does it detect? Specifically, precision measures the proportion of the model’s predicted invasions that are true invasions, while recall measures the proportion of the true invasions that are detected by the model. To compute AUC-PR skill, we divided AUC-PR by the base rate to measure how much better our models were than random/chance (i.e. a no-skill model). For uncertainty on metrics, we performed a cluster bootstrap resampling of zones (not zone-weeks), to respect within-zone correlation.
For operationally focused evaluation quantities, we also computed the mean invasion probability rank of the newly invaded health zones, recall with a fixed weekly top-K health zones budget, reliability diagrams, log scores, and further precision and calibration metrics.
To select models, we balanced discrimination (AUC-PR and average rank of newly invaded zones) against calibration (logarithmic score) across both one- and two-week horizons.
Authors
Key contributors to data collection, molecular testing, data interpretation, analysis, and writing:
Institut National de Recherche Biomédicale (INRB) and Institut National de Santé Publique (INSP), Kinshasa, Democratic Republic of the Congo, and partners
Dav M. Ebengo (INRB)*
Cathal Mills (Oxford)
Ciara Judge (Oxford)
Pierre Akilimali (INSP; Kinshasa School of Public Health, Faculty of Medicine, University of Kinshasa, Kinshasa, Democratic Republic of the Congo)*
Bernardo Gutierrez (Oxford)
Ellie Bourgikos (Oxford)
Joseph L.-H. Tsui (Oxford)
Tania Bishola Tshitenge (INRB)
Adelard Lofungola (INSP)
Joel Kosianza (WHO)
Abdi Mahamud (WHO)
Benjamin Kanku (INSP)
Etien Koua (WHO)
Olga Ntumba (WHO)
Lorenzo Subissi (WHO)
Nicksy Gumede (WHO)
Marie Roseline Darnycka Belizaire (WHO)
Otim Patrick Cossy Ramadan (WHO)
John Otokoye Otshudiema (WHO)
Thierno Balde (WHO)
Brian Ajong (WHO)
Mamadou Saliou Kalifa Diallo (WHO)
Gianni Donkor (WHO)
Tamayi Mlanda (WHO)
Richy Ngombo (WHO)
Jean-Marie Tshilumbu (WHO)
Martina McMenamin (WHO)
Linus Bengtsson (Flowminder)
Thomas Smallwood (Flowminder)
Romain Goldenberg (Flowminder)
Benjamin Reddy (Oxford)
Simon Cauchemez (Pasteur)
Justus Nsio (AfricaCDC)
Yap Boum (AfricaCDC)
Jeanine Nkakulu (AfricaCDC)
Yenew Kebede (AfricaCDC)
Collins Tanui (AfricaCDC)
Mosoka Fallah (AfricaCDC)
Pascal Adroba Tandele (Laboratoire Provincial de Bunia, Ituri)
Micheline Magapa Chuma (Division Provincial de la Santé, Ituri)
Didier Bompangue Nkoko (INOHA)
Francine Ntoumi (Fondation Congolaise pour la Recherche Médicale)
Jean-Jacques Muyembe Tamfum (INRB)
Olivier le Polain de Waroux (WHO)
Samuel V. Scarpino (Northeastern U)
Moritz U.G. Kraemer (Oxford)
Christian Ngandu (INSP)
Dieudonne Mwamba Kazadi (INSP)
Placide Mbala-Kingebeni (INRB)*
Acknowledgements
We are grateful to all institutions and partners for their support to surveillance efforts in the DRC.
Statement on continuing work and analyses before publication
Please note that this data is based on work in progress and should be considered preliminary. Our analyses are ongoing, and a publication communicating our findings is in preparation. If you intend to use data for similar analyses and before our publication, please contact Prof. Placide Mbala-Kingebeni (INRB, DRC), Prof. Dav Ebengo (INRB, INOHA), or Prof Pierre Akilimali (INSP).
Collaborating institutions and agencies
-
L’Institut National de Recherche Biomédicale (INRB), Democratic Republic of the Congo
-
Institut National de Santé Publique (INSP), Democratic Republic of the Congo
-
Kinshasa School of Public Health, Faculty of Medicine, University of Kinshasa, Kinshasa, Democratic Republic of the Congo
-
Africa Centres for Disease Control and Prevention, Addis Ababa, Ethiopia
-
World Health Organization, Geneva, Switzerland
-
World Health Organization Country Office, Kinshasa, Democratic Republic of the Congo
-
World Health Organization Regional Office for Africa, Brazzaville, Republic of Congo
-
Pandemic Sciences Institute, University of Oxford, Oxford, UK
-
Northeastern University, Boston, MA, USA
References
1. Open Spatial Demographic Data and Research. WorldPop https://www.worldpop.org/ (2020).
2. German Federal Foreign Office, Truth, Beauty, PIK & UniBW. Climate—Conflict—Vulnerability Index (CCVI). CCVI — Climate Conflict Vulnerability Index.
3. GRID. GRID3 COD - Health Facilities v8.0. (2025).
4. Dasgupta, A. et al. Scalable, open-access and multidisciplinary data integration pipeline for climate-sensitive diseases. Wellcome Open Res. 10, 467 (2025).
5. Flowminder Foundation. Population movements from Bunia, Mongbwalu and Rwampara based on privacy-secure analysis of mobileoperator data from Vodacom Congo. (2026).
6. https://www.flowminder.org/resources/publications-reports/drc-reports-publications
7. Project OSRM. https://project-osrm.org/
8. Epidemiological Parameters: Ebola Bundibugyo Virus (BVD). Epidemiological Parameters: Ebola Bundibugyo Virus (BVD) • epireview