Molecular evolutionary analysis of the current Bundibugyo virus disease outbreak in DRC and Uganda

[mahan.ghafari]

Authors: Gina Cuomo-Dannenburg and Mahan Ghafari

Note on data use

This is an informal analysis and is not intended as a formal preprint or manuscript in its current form. The aim of this post is to provide a preliminary molecular evolutionary analysis of the current Bundibugyo virus disease outbreak, based on the 10 available outbreak-associated Bundibugyo virus (BDBV) sequences available as of 28 May 2026.

Specifically, we use these sequences to reconstruct the recent evolutionary history of the outbreak, estimate the time to the most recent common ancestor (tMRCA) of the sampled viruses, and derive preliminary estimates of the epidemic growth rate and corresponding reproduction number.

The sequences analysed here are collated in the Pathoplexus SeqSet: PP_SS_2146.2 | Pathoplexus. The outbreak sequences were first described and shared through this Virological post, and we encourage readers to consult that post for additional details on sequence generation and data release.

We present this post in line with the data-sharing agreements for Restricted-Use Data in Pathoplexus. We are grateful to all laboratories and researchers who generated and rapidly shared the sequence data that made this analysis possible. We also thank the Pathoplexus team for rapidly enabling BDBV sequence sharing.

Summary

We analysed 10 available BDBV sequences from the Democratic Republic of Congo and Uganda associated with the current 2026 Bundibugyo outbreak.

Our main findings are:

1. The inferred tMRCA of the sampled viruses is sensitive to the assumed substitution rate. Across six alternative substitution-rate assumptions, the median posterior tMRCA estimates range from late January to late March 2026, with wide 95% HPD intervals (Figure 2; Table 2). We therefore interpret these estimates as rate-conditioned sensitivity analyses rather than a single definitive estimate of the outbreak origin time.

2. The sampled viruses are already genetically diverse. We identified 23 unique mutations across the 10 sequences. This level of diversity suggests that the sampled viruses do not represent the very beginning of transmission, but instead reflect an outbreak that had already been circulating, at least partly undetected, for several weeks before sequences became available.

3. The estimated epidemic doubling time is relatively slow but highly uncertain. Although the posterior means imply positive growth under all substitution-rate assumptions, the 95% HPD intervals are wide and include scenarios with little or no growth. These estimates should therefore be interpreted as exploratory, rate-conditioned quantities rather than reliable standalone estimates of epidemic growth or transmissibility.

4. The corresponding reproduction number (R_0) estimates are centred above 1, with mean estimates ranging from 1.31 to 1.55 (Table 2), although substantial uncertainty remains because of the limited number of sequences currently available.

5. These estimates should be interpreted cautiously. The analysis is based on only 10 sequences collected over a short sampling window. The estimates of doubling time and R_0 are therefore likely to be affected by sampling bias and should be viewed as preliminary.

Data curation and alignment

We downloaded the outbreak-associated BDBV sequences from Pathoplexus, using the sequence set available as of 28 May 2026: PP_SS_2146.2 | Pathoplexus.

Sequences were aligned to the reference genome FJ217161.1 using MAFFT v7.490 with the -keeplength flag. The alignment was then manually inspected in AliView.

We performed ancestral node reconstruction using IQ-TREE v3.0.1 with the -anr flag and the best-fitting substitution model, HKY+F+I. We used this analysis to reconstruct the consensus sequence of the most recent common ancestor (MRCA) of the sampled outbreak sequences. This reconstructed MRCA sequence was then used as the root for the outbreak mutation reconstruction shown below.

Outbreak mutation reconstruction

The 10 available outbreak sequences are already genetically divergent from one another. The earliest sampled sequence, PP_006XHKB.2|DRC|2026-05-03, is closest to the reconstructed MRCA sequence, with two mutations from the root. All other sequences carry between three and five mutations from the reconstructed MRCA (Figure 1).

This pattern suggests that the sampled sequences do not capture the earliest phase of the outbreak. Instead, they are consistent with a larger, partly unsampled transmission chain that had already been ongoing before the available sequences were generated.

Figure 1 | Manually reconstructed mutation tree of ten available outbreak-associated Bundibugyo virus sequences as of 28 May 2026. The topology was constructed from the observed single nucleotide polymorphisms between the reconstructed MRCA sequence and the outbreak-associated sequences, rather than inferred using a phylogenetic dating program. Mutations are annotated relative to the reconstructed MRCA genome. Bold mutations represent non-synonymous changes. *The mutation T18905G could not be resolved at the amino-acid level because position 18906 is missing following alignment.

Across the 10 sequences, we identified 23 unique mutations – see supplementary CSV file contains the full list of mutations across the genome. Of these, 11 fall within coding regions of the genome. Table 1 lists the coding-region mutations, including their reference and alternative codons, amino-acid changes, mutation type, and coding region.

Five mutations occur in the glycoprotein region: T6657C, C7049T, T7405C, T7447C, and T7975C. Four of these five glycoprotein mutations are non-synonymous.

Table 1 | Coding-region mutations identified among the ten outbreak-associated BDBV sequences. The table lists unique mutations with their reference (ancestor) and alternative (derived) codons and amino acids, mutation type, and coding region. NP = nucleoprotein; GP_full = full glycoprotein; sGP = soluble secreted glycoprotein; ssGP = alternative secreted glycoprotein; L = RNA polymerase.

The reconstructed ancestral sequence should be interpreted as the ancestor of the currently sampled sequences, not necessarily the true ancestor of the outbreak. Because only a small number of sequences are currently available, and because these sequences may be biased spatially and temporally, the reconstructed MRCA may not fully represent the circulating viral diversity.

Estimates of tMRCA, epidemic growth rate, and reproduction number

We estimated the time to the most recent common ancestor (tMRCA) and the number of epidemic doublings per year using the Delphy web interface.¹ For all analyses, we used an exponential-growth coalescent model, reflecting the expectation of a growing epidemic. An earlier Bayesian phylogenetic analysis of these same outbreak sequences was posted on Virological. Our analysis differs in several respects, including the use of Delphy, manual reconstruction of outbreak mutations, and exploration of multiple fixed substitution-rate assumptions. We included this rate-sensitivity analysis because, as seen for Zaire ebolavirus during and after the West African outbreak, early outbreak-based substitution-rate estimates can be higher than later estimates^3-6, which can in turn affect inference of the tMRCA, epidemic growth rate, and corresponding reproduction number. We therefore view the results presented here as a complementary reanalysis rather than a replacement of the earlier work.

We considered six substitution-rate assumptions. In five analyses, we fixed the substitution rate using previously published Ebola virus estimates. These included one estimate from a previous Bundibugyo virus outbreak² and several estimates from the 2013-2016 West African Ebola virus outbreak.^3-6 In the sixth analysis, we estimated the substitution rate directly from the 10 outbreak sequences analysed here.

All Delphy analyses were run until the minimum effective sample size across parameters was greater than 100, and convergence was assessed by visual inspection of the parameter traces. The tMRCA and number of doublings per year are reported as posterior summaries, using the median and 95% highest posterior density (HPD) intervals from Delphy.

We then translated the estimated epidemic growth rate into an estimate of the basic reproduction number, R_0, using:

R_0 = \exp\left(\frac{D \cdot \log_e(2)}{365.25} \cdot T_g\right),

where D is the number of doublings per year and T_g is the generation time.
We calculated the generation time as:

T_g = T_{\text{inc}} + \frac{T_{\text{inf}}}{2},

where T_{\text{inc}} is the incubation period and T_{\text{inf}} is the infectious period.

We used systematically collected Ebola virus disease parameters⁷ available through {epireview}⁸ to identify previous estimates of epidemiological delays for Bundibugyo virus disease.^9,10 This gave a generation time estimate of 9.6 days. The uncertainty intervals for R_0 were obtained by propagating the 95% HPD interval from the estimated number of doublings per year.

The posterior tMRCA estimates across substitution-rate assumptions are shown in Figure 2. Each estimate should be interpreted conditionally on the assumed substitution rate used in that analysis. The resulting tMRCA estimates vary substantially across rate assumptions, illustrating that early phylogenetic dating from this small sequence set is sensitive to the clock rate. Although there is substantial variation among substitution-rate assumptions, all analyses place the median posterior tMRCA no later than 23 March 2026.

The estimated growth rates and R_0 values are likely affected by sampling bias. The sequences included here were collected between 03 May 2026 and 20 May 2026, and therefore represent a short time window around the identification and declaration of the BVD outbreak.¹¹ Situation reports^12,13 and independent epidemiological analyses¹⁴ suggest that, even at the time these sequences were generated, the outbreak was already substantially larger than the available sequence data alone would indicate.

Table 2 | Substitution rates and corresponding estimates of tMRCA, doublings per year, and reproduction number. All substitution rates, except the final estimate, were fixed using previous Ebola virus substitution-rate estimates. These include an estimate from a previous Bundibugyo virus outbreak² and estimates from the 2013-2016 Ebola virus outbreak in West Africa.^3-6 The final row shows an analysis in which the substitution rate was estimated directly from the 10 outbreak sequences analysed here. All analyses were performed in Delphy¹ and run until the minimum effective sample size was greater than 100. The reproduction number was estimated from the number of doublings per year using a generation time estimate derived from published incubation and infectious-period estimates for Bundibugyo virus disease.^7-10

Figure 2 | Posterior estimates of the time to the most recent common ancestor across substitution-rate assumptions. We used Delphy¹ to analyse the ten outbreak-associated sequences under an exponential-growth coalescent model. Five analyses used fixed substitution rates from previous Ebola virus studies, while one analysis estimated the substitution rate directly from the outbreak sequences. Points show the posterior mean tMRCA, and error bars show the 95% highest posterior density (HPD) interval. The uncertainty in the top estimate reflects the additional uncertainty introduced by jointly estimating the substitution rate from the same set of outbreak sequences.

A full list of all mutations and their corresponding amino-acid changes, including the sequences in which they occur, is available here: mutation_bdbv.csv.zip (1.5 KB)

Bibliography

1. Varilly P, Schifferli M, Yang K, et al. Delphy: scalable, near-real-time Bayesian phylogenetics for outbreaks. 2025; 2025.03.25.645253.

2. Hulseberg CE, Kumar R, Di Paola N, et al. Molecular analysis of the 2012 Bundibugyo virus disease outbreak. Cell Reports Medicine 2021; 2: 100351.

3. Gire SK, Goba A, Andersen KG, et al. Genomic surveillance elucidates Ebola virus origin and transmission during the 2014 outbreak. Science 2014; 345: 1369–1372.

4. Erratum for the Report “Mutation rate and genotype variation of Ebola virus from Mali case sequences” by T. Hoenen, D. Safronetz, A. Groseth, K. R. Wollenberg, O. A. Koita, B. Diarra, I. S. Fall, F. C. Haidara, F. Diallo, M. Sanogo, Y. S. Sarro, A. Kone, A. C. G. Togo, A. Traore, M. Kodio, A. Dosseh, K. Rosenke, E. de Wit, F. Feldmann, H. Ebihara, V. J. Munster, K. C. Zoon, H. Feldmann, S. Sow. Science 2015; 348: aac5674.

5. Hoenen T, Safronetz D, Groseth A, et al. Mutation rate and genotype variation of Ebola virus from Mali case sequences. Science 2015; 348: 117–119.

6. Holmes EC, Dudas G, Rambaut A, Andersen KG. The evolution of Ebola virus: insights from the 2013–2016 epidemic. Nature 2016; 538: 193–200.

7. Nash RK, Bhatia S, Morgenstern C, et al. Ebola virus disease mathematical models and epidemiological parameters: a systematic review. Lancet Infectious Diseases 2024; 24: e762–e773.

8. Naidoo T, Nash R, Morgenstern C, et al. epireview: Tools to update and summarise the latest pathogen data from the Pathogen Epidemiology Review Group (PERG). 2025. GitHub - mrc-ide/epireview · GitHub.

9. Choi B, Busch S, Kazadi D, et al. Modeling outbreak data: analysis of a 2012 Ebola virus disease epidemic in DRC. BIOMATH 2019; 8: ID: 1910037.

10. Wamala JF, Lukwago L, Malimbo M, et al. Ebola hemorrhagic fever associated with novel virus strain, Uganda, 2007–2008. Emerging Infectious Diseases 2010; 16: 1087–1092.

11. World Health Organization. Disease Outbreak News: Bundibugyo Virus Disease, Democratic Republic of the Congo and Uganda. 2026. Ebola disease caused by Bundibugyo virus, Democratic Republic of the Congo & Uganda.

12. World Health Organization, African Region. Bundibugyo Virus Disease Outbreak; Weekly External Situational Report 01. 2026. https://iris.who.int/server/api/core/bitstreams/b6e1e783-91c3-43c8-ab90-16ceaa9948f0/content.

13. World Health Organization, African Region. Bundibugyo Virus Disease Outbreak; Weekly External Situational Report 02. 2026. DSpace.

14. McCabe R, Ebbarnezh L, Okware S, et al. Estimation of the size of the Ebola outbreak caused by Bundibugyo virus in the Democratic Republic of the Congo: May 20, 2026 update. Imperial College London, 2026. DOI: https://doi.org/10.25560/13005307.

[arambaut]

Hi Gina & Mahan,

A couple of thoughts on this analysis -

1. Although the Hulseberg et al rate might seem the most relevant, as it is for BDBV, those authors are using the divergence between the previous two outbreaks and thus they are essentially getting a rate for the evolution of the virus in the reservoir species (it is notable this is the slowest rate of them all). It is problematic enough that we are using rates from one virus species to calibrate another, differences in host species is likely to be the greater effect.
2. Some of the rates are not independent - Hoenen et al (once they fixed their erroneously low rate in the erratum) is the same data as Gire et al with the addition of 2 (epidemiologically linked) genomes from Mali.
3. Finding the maximum overlapping interval in the HPDs of the TMRCA under these different rates as the best estimate for the TMRCA is artificially reducing the uncertainty in the estimate by clipping the HPDs - i.e., the more variation there is in the collated rates from the literature, the smaller this overlap will be and the appearance of precision will increase.
4. Furthermore, the inclusion of the estimates for the low rates will bias this considerably back in time – excluding Hulseberg et al and Ghafari et al (an unpublished analysis) would make the HPD overlap region considerably more recent and wider.
5. Or equally, if the Hulseberg et al rate is the closest to reality (because it is the same virus species, despite being a different host species) then including the faster EBOV rates would bias the estimate towards the present day.

Our similar analysis, done with a slightly larger dataset and using two of the same rate estimates from the literature (Gire et al and Holmes et al to represent a range for rates observed at different stages in the same epidemic), takes the approach of considering the union of the HPDs for these rate possibilities. This is, in my view, a more appropriate representation of the bounds of our knowledge of the TMRCA at this stage. These rates from the literature are not independent estimates of the same quantity but rather alternative possibilities as to what the true underlying rate for this epidemic could be.

[ginacuomo-dannenburg]

Dear Andrew,

Thank you so much for your comments. The second part of our analysis where we explored a range of substitution rates and how they affect estimates of R0 and tMRCAwere not intended to be independent estimates but rather an exploration of a plausible range of values of substitution rates and how they influence key estimates for the early phase of the epidemic. We took these estimates to represent rates from varied time frames, datasets, context, and in the case of Ghafari et al a meta-analysis of published zaire rate estimates (using the median of all published rates for this analysis)

We agree that the Hulseberg rate likely represents the reservoir rate, but still feel that this represents a valid rate to explore the potential values of substitution rates and the impact on estimates of R0 and tMRCA for these sequences. The fact that it is estimated across two outbreaks is valid and the inferred rate is similar with our metadata of Zaire rates (Ghafari et al., In preparation). We would expect rates to be lower considering the longer timespan over which the rate was estimated in Hulseberg. At this stage of the outbreak and given the limited previous transmission of BVD, it is hard to determine the likely value of the true substitution rate. The intention was to explore the range of substitution rates and how they may influence estimates of R0 and tMRCA. As you have mentioned, the higher substitution rates are associated with higher estimates of R0, a shorter doubling time and therefore a more recent tMRCA.

The point you make about the addition of more estimates reducing the amount of shared 95% HPDs is correct and indeed may have given a false impression of a tighter or more accurate estimate of tMRCA across the estimates. Table 2 and Figure 2 both show the full uncertainty in these estimates of tMRCA but we agree that this could be communicated more clearly in the results. We will revise Figure 2 to remove the shaded region and revise text to reflect this change.

[arambaut]

I think the most appropriate way of presenting these estimates is the HPD interval for the tMRCA conditioned independently on each of the rate estimates being the actual underlying rate. This is the approach we took in https://virological.org/t/1032/6

I still don’t think the BDBV rate is appropriate to use - the issue is that it is evolving in a different host species. Although the value is lower this is not necessarily because it is over a longer timespan, it could be because the rate in bats is different for all sorts of molecular and cellular reasons (or both).

As we also mention in the caveats to our post that we decline to report our estimated growth rates or derived metrics like R0 at the present because the HPDs include both zero and implausibly high rates. A mean or median of such an interval is potentially misleading. I believe it is better not to report it at this stage - if the HPDs aren’t meaningful then neither are the point estimates.

[mahan.ghafari]

Thanks, Andrew — this is very helpful, and I agree that the cleanest way to frame the analysis is as tMRCA estimates conditioned independently on each assumed substitution rate, which is what we have done in Table 2. As Gina said, we will revise the figure and text accordingly and remove the shared-HPD shading.

On the BDBV rate, I take your point that the Hulseberg et al estimate may reflect evolution in a different host context, and that reservoir-associated rates may not be directly comparable to rates during human-to-human transmission. We included it because, to our knowledge, it is the only available BDBV-specific estimate and because the aim was exploratory: to show how sensitive tMRCA and phylodynamic estimates are to plausible alternative clock rate assumptions. We don’t know whether the substitution rate in the reservoir differs from that in human-to-human transmission in this case, and therefore whether the Hulseberg et al. estimate is inappropriate for this outbreak. We will make this caveat more explicit.

On the growth rate and R_0 estimates, I agree that the HPDs are very wide and include scenarios with shrinking outbreak size, so these estimates should not be interpreted as reliable point estimates on their own. We had intended the table to present the full uncertainty for each assumed rate and associated parameters. We will revise the wording to make clear that these are exploratory, rate-conditioned quantities, and that the current sequence data are insufficient for reliable independent inference of growth rate or R_0.

Our intention was not to present this as a definitive or independent phylodynamic estimate, but to add a complementary exploratory analysis: reconstructing the outbreak-associated mutations from the available genomes, showing how sensitive tMRCA estimates are to alternative substitution-rate assumptions, and making explicit how these choices propagate into inferred growth and R_0. We agree that the latter quantities are highly uncertain at this stage, but we think showing that uncertainty transparently is useful, particularly as a caution against over-interpreting early sequence-only phylodynamic estimates.

More broadly, one of the points we hoped to highlight is precisely that, with only a small number of sequences, phylodynamic inference alone is not yet sufficient to robustly infer the outbreak tMRCA or growth dynamics. We are planning to update the analysis in the coming days to link the existing genomic estimates more explicitly with the epidemiological data.

[erik.volz]

Mahan & Gina - Minor suggestions for improving on the R0 estimation which will come into play when more sequences become available: I think a coalescent model for growth rates combined with renewal equations is an excellent and robust approach, but the exponential generation time assumption is not great for Ebola. Note that a plain BD-sky model would also depend on exponential Tg and would lead you to under-estimate R. There was a recent paper about this and some suggested workarounds: https://doi.org/10.1016/j.epidem.2026.100910 You might use a structured model, but more simple would be to estimate growth and apply the renewal eqn directly with one of the published Gamma Tg’s.

[mahan.ghafari]

Thanks for flagging this, Erik — very useful. And it makes sense, and we’ll update the analysis to take this into account.

[ginacuomo-dannenburg]

Hi Erik,

Thank you so much for a great suggestion of how to improve this analysis. Just to clarify, in the original analysis we used only a central estimate of the generation time for Bundibugyo because there was no generation time estimate available.

Below, we extended this original work by using the method you suggested. We use the daily growth rate as estimated by converting the number of doublings per year which are direct outputs from the Delphy analysis. There are currently no published estimates of generation time for BVD and therefore we used the generation time estimates from the WHO Ebola Response Team, with a mean of 15.3 days and standard deviation of 9.3 days. This estimate was based on the first 9 months of the outbreak in West Africa and is assumed to be representative of a non-isolating cohort.

We convert the growth rate r to an estimate of R_0 using the below renewal equation, as described in Green, Ferguson and Cori, 2022:

R_0 = \frac{(r+b)^a}{b^a}

where r is the epidemic growth rate, and a and b are the shape and rate of the Gamma distribution of the generation time, respectively.

Below is an updated version of Table 1 to reflect this new analytic equation for the reproduction number (only the reproduction numbers have changed). Note that this increases the uncertainty of the reproduction number estimates. In particular the upper interval value increases, and the central estimate also increases.

Substitution rate	tMRCA	Doublings/year	R_0	Source	Species
0.84 \times 10^{-3} [fixed]	2026-01-31 [2025-12-22, 2026-03-21]	14.9 [-2.4, 32.2]	1.49 [0.93, 2.23]	Hulseberg et al., Cell Rep Med, 2021	Bundibugyo
0.91 \times 10^{-3} [fixed]	2026-02-06 [2025-12-15, 2026-03-15]	14.9 [-2.4, 32.2]	1.54 [0.82, 2.56]	Ghafari et al., [In preparation] Median of Zaire substitution rates	Zaire
1.2 \times 10^{-3} [fixed]	2026-02-26 [2026-01-25, 2026-03-31]	19 [-7, 51]	1.65 [0.81, 3.26]	Holmes et al., Nature, 2016 (mean of rates for 2013-2016 EBOV outbreak)	Zaire
1.9 \times 10^{-3} [fixed]	2026-03-23 [2026-03-01, 2026-04-12]	24.0 [-10, 67]	1.85 [0.74, 4.33]	Gire et al., Science, 2014 (early outbreak estimate)	Zaire
1.99 \times 10^{-3} [0.06, 4.15]	2026-02-17 [2025-09-10, 2026-04-27]	21.0 [-17, 64]	1.73 [0.60, 4.11]	Fitted to the 10 sequences from this outbreak	Bundibugyo

The number of sequences available to analyse has increased substantially following an updated analysis available here: Genomic epidemiology of the ongoing 2026 Bundibugyo Virus Disease outbreak in the Democratic Republic of the Congo.

Here we simply extend the previous analysis based on the sequences that were available at the time, although this approach would be applicable to an expanded set of genomic sequences to estimate the reproduction number from the genomic data available and would increase the reliability of these estimates as discussed earlier in this discussion thread.

References:
WHO Ebola Response Team et al., 2014 https://www.nejm.org/doi/full/10.1056/NEJMoa1411100
Green, Ferguson and Cori et al. 2022 https://doi.org/10.1098/rsif.2021.0429