A natural history of AMR in Klebsiella pneumoniae: Global diversity, predictors, and predictions of evolutionary pathways

15 Antimicrobial resistance (AMR) is a substantial and growing global health burden. Understanding, 16 and predicting, its evolution in specific pathogens will help responses across scales from individual 17 patient cases to large-scale policy. Data-driven approaches to this question often focus more on 18 genomes and less on the evolutionary dynamics generating these genomes. Here, we use global 19 data on AMR features in Klebsiella pneumoniae with hypercubic transition path sampling 20 (HyperTraPS), a machine learning approach, for Bayesian inference of the evolutionary pathways of 21 AMR in K. pneumoniae (KpAMR) in 102 diUerent countries, territories and areas. We identify 22 “globally consistent” evolutionary behaviours that hold across countries, and “globally divergent” 23 behaviours including carbapenem and fluoroquinolone resistance that vary across countries. We 24 show how these divergent dynamics covary both with public health superregion and drug use policy, 25 and reveal competing evolutionary pathways within and between countries. Using newly-26 sequenced data across several decades from sub-Saharan Africa, we show that this inferred global 27 roadmap of KpAMR evolution successfully predicts prospective evolutionary dynamics. Together, 28 we hope that the ability to characterize and predict evolutionary dynamics of AMR acquisition, 29 connected to socio-economic and drug policy predictors, will help strengthen our understanding of 30 AMR evolution worldwide. 31

32 Antimicrobial resistance (AMR), the resistance of microbial pathogens to the medicines we use to 33 treat them, is a major threat to global health. Antibiotic resistance, the resistance of bacteria to 34 antibacterials, is an important subset of AMR. In 2019 and 2021, 1.27m and 1.14m deaths 35 respectively were attributable to bacterial AMR (many more were associated with AMR), with highest 36 rates in sub-Saharan Africa (Murray et al., 2022; Naghavi et al., 2024). The evolution of genetic and 37 phenotypic features that cause AMR in bacterial pathogens is the focus of tremendous research 38 interest, from basic biology through epidemiology, clinical studies and health policy. 39 Klebsiella pneumoniae (Kp) is a gram-negative bacterium and a major cause of nosocomial 40 infections across settings. (Related species and subspecies of Kp form the KpSC (species complex); 41 in this report we will focus on K. pneumoniae sensu stricto). Kp is found in multiple environmental 42 niches and can both cause opportunistic infections in hospitalized patients and community 43 acquired infections, including severe infections caused by hypervirulent strains (Wyres et al., 2020). 44 Kp has been identified as an important traUicker for antibiotic resistance genes from diUerent 45 ecological niches into several other clinically important pathogenic bacteria (Wyres & Holt, 2018). 46 Kp contributes significantly to the global burden of AMR. In 2021, 45 600 and 36 000 deaths were 47 attributed to carbapenem and fluoroquinolone resistant Kp respectively. In 2021, Kp accounted for 48 .CC-BY 4.0 International licenseperpetuity. It is made available under a preprint (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in The copyright holder for thisthis version posted September 21, 2025. ; https://doi.org/10.1101/2025.09.20.677523doi: bioRxiv preprint 12.7% of all deaths attributable to AMR in the population over five years old and 19.4% of AMR-49 attributable deaths in children under five years (Naghavi et al., 2024). 50 51 Given the substantial health burden associated with the evolution of AMR in Kp, large-scale data-52 driven research has explored the genomic structure of, and AMR features in, Kp populations 53 (Argimón et al., 2021; Hetland et al., 2025; Holt et al., 2015; Munk et al., 2022; Wyres et al., 2020). 54 Here, we consider a complementary perspective: the evolutionary pathways by which these 55 features are acquired (Renz et al., 2024). Many related subquestions here can broadly be unified as: 56 from a global perspective, in what orders are AMR features acquired over diUerent instances of 57 pathogen evolution? For example, given a newly observed isolate in the clinic with a given set of 58 AMR features, can we predict which drug resistance(s) will evolve next (and thereby inform 59 treatment guidelines accordingly)? And can we identify which extrinsic features – from physical 60 environment, demography, drug policy, or others – influence the evolutionary acquisition of AMR? 61 These questions require approaches that can use the expanding volumes of AMR data to make and 62 train descriptive and predictive evolutionary models. However, the large sets of features involved in 63 typical AMR systems are often not addressable using traditional methods from evolutionary biology. 64 For example, an often-studied dataset of Mycobacterium tuberculosis isolates contains 65 susceptibility-resistance information for a panel of 10 diUerent drugs (Casali et al., 2014). 66 Comparative methods like the Mk (Markov k-state) model often struggle for more than around 7 67 potentially coupled characters under study (Johnston & Diaz-Uriarte, 2024). An alternative class of 68

for studying the evolution of multiple characters – evolutionary accumulation modelling or 69 EvAM – has emerged jointly from the cancer and evolution literatures (Diaz-Uriarte & Herrera-Nieto, 70 2022; Renz et al., 2024; Schill et al., 2024), but has traditionally considered independent samples, 71 without the phylogenetic relationships found in bacterial populations (Dewar et al., 2025). While 72 many EvAM approaches have the ability to infer dependencies between multiple interacting 73 characters, neglecting this phylogenetic relatedness can lead to substantial artefacts in analysis 74 (Maddison & FitzJohn, 2015). 75 Hypercubic transition path sampling (HyperTraPS) is an EvAM approach designed to incorporate 76 phylogenetic information into the (Bayesian) inference of evolutionary pathways for multiple 77 interacting characters (Aga et al., 2024; Greenbury et al., 2020; Johnston & Williams, 2016). Other 78 EvAM approaches including TreeMHN (Luo et al., 2023), HyperHMM (Moen & Johnston, 2023), and 79 HyperMk (Johnston & Diaz-Uriarte, 2024), also support phylogenetically coupled data, but we focus 80 on HyperTraPS here because of its flexibility (Aga et al., 2024; Renz et al., 2024). In this article, we 81 harness large-scale global datasets on AMR in Kp (KpAMR) with HyperTraPS to learn the structure 82 and variability of evolutionary pathways of drug resistance in Kp across the globe, as well as the 83 extrinsic factors shaping these pathways. 84 85 .CC-BY 4.0 International licenseperpetuity. It is made available under a preprint (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in The copyright holder for thisthis version posted September 21, 2025. ; https://doi.org/10.1101/2025.09.20.677523doi: bioRxiv preprint Figure 1. Global genomic data on KpAMR. (A) We use 47 721 Klebsiella pneumoniae (Kp) isolates from 86 around the world, with (i) country-specific counts of isolates given by colour and years of isolation given in (ii). 87 (B) UpSet plot giving counts of di\erent KpAMR characters across these isolates. Character names: AGly 88 (resistance to aminoglycosides); Col (colistin); Fcyn (Fosfomycin); Flq (Fluoroquinolones); Gly 89 (Glycopeptides); MLS (Macrolides); Phe (Phenicols); Rif (Rifampin); Sul (Sulfonamides); Tet (Tetracyclines); Tgc 90 (Tigecycline); Tmt (Trimethoprim); Bla_a (β-lactam resistance without extended spectrum or inhibitor-91 resistance); Bla_inhR (resistance to β-lactam/inhibitor combinations); Bla_ESBL (resistance to extended-92 spectrum β-lactams); Bla_ESBL_inhR (resistance to extended-spectrum β-lactam/inhibitor combinations); 93 Bla_Carb (resistance to carbapenems); SHV (SHV β-lactamase with expanded enzyme activity), Bla_chr (SHV 94 alleles conferring resistance to ampicillin), Omp (outer membrane protein). Further details available via 95 Kleborate documentation at https://github.com/klebgenomics/Kleborate. We will use -a: acquired (via 96 horizontal gene transfer); -m: mutation; -chr: chromosomal (intrinsic). 97

98 Local and global inferred pathways of AMR evolution in Klebsiella 99 To learn and compare the likely pathways by which KpAMR is acquired in diUerent countries, we 100 obtained a dataset of 47 721 genomes from 102 diUerent countries, territories and areas, and 101 extracted details of the presence or absence of genes corresponding to each of 22 AMR classes, 102 using PathogenWatch and Kleborate (Argimón et al., 2021) (Lam et al., 2021) (Methods; Fig. 1A-B). 103 We refer to presence/absence of resistance genes for each of these classes as binary “resistance 104 characters” (“character” referring to a particular property of a species in evolutionary biology). The 105 set of characters we consider, following Kleborate, is given in the Fig. 1 caption; in this report we will 106 use the shorthand defined there, grouping genes by drug resistance class (and Lahey class for β-107 lactamases) (Gupta et al., 2014; Lam et al., 2021; Tsang et al., 2024). 108 For each country in our dataset, we next used hypercubic transition path sampling (HyperTraPS), a 109 machine learning approach, to infer the “pathways” of KpAMR character acquisition: that is, the 110 ordered sequences with which characters are acquired over time. Figs. 2A-C illustrate this pipeline 111 for a single country example (Romania, chosen for its moderate number of samples). HyperTraPS 112 takes presence-absence patterns of characters on a phylogeny as input (Fig. 2A) and outputs an 113 evolutionary model describing the pathways of character acquisition most supported by the data 114 (Aga et al., 2024; Greenbury et al., 2020; Johnston & Williams, 2016; Renz et al., 2024). This model 115 can be summarized, for example, as the set of posterior probabilities Pij with which character i is 116 acquired at ordering j in a putative evolutionary process (Fig. 2B). We refer to a plot where Pij is 117 visualised with point size, as in Fig. 2B, a “bubble plot” . The full output of the inference process is a 118 parameterised transition network describing the probability of diUerent sequences of transitions 119 through the state space of possible resistance characters (Fig. 2C). Respecting the fact that many 120 KpAMR characters may be acquired reversibly (for example, via gain and loss of plasmids) and 121 HyperTraPS assumes irreversible dynamics, we confirmed with synthetic control studies that 122 reversible dynamics marginally increase uncertainty in our evolutionary inference but do not lead to 123 systematic errors (Supp. Fig. 1). 124 .CC-BY 4.0 International licenseperpetuity. It is made available under a preprint (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in The copyright holder for thisthis version posted September 21, 2025. ; https://doi.org/10.1101/2025.09.20.677523doi: bioRxiv preprint 125 Figure 2. Inferring evolutionary dynamics from genome data across countries. (A) Example source data 126 from Romania, with 118 Kp isolates connected via a putative phylogeny derived from LIN codes (see Methods) 127 and assigned presence (dark pixel) or absence (white pixel) markers for resistance to each of 22 antibiotic 128 groups. (B) Summary “bubble plot” of evolutionary pathways inferred for Romanian data with HyperTraPS. The 129 area of a bubble gives the posterior probability that a given character (row) is acquired at a given ordering in 130 the accumulative evolutionary process. For example, Bla_chr is inferred to be a likely early acquisition, while 131 Bla_ESBL-a likely occurs later and Tgc_a is inferred likely later still (or never). (C) A subset of the hypercubic 132 transition graph inferred by HyperTraPS for Romanian data. Each edge’s width gives the probability of a 133 transition involving a single character acquisition from the upper state to the lower state. The uppermost state 134 is the ancestral state with no AMR characters. Only early transitions above a threshold probability 0.025 are 135 plotted for clarity. (D) Global plot, taking mean acquisition probability over all countries in the dataset. 136 Probabilities higher than the “null model” of uniform acquisitions are plotted in blue. 137 138 The large majority of genomes with timing information were isolated in the period 2015-2020 (Fig. 139 1A). Our approach here does not directly consider the specific, real-world timings of the 140 evolutionary transitions involved – although, given reliable data on the divergence times of all pairs 141 of isolates, HyperTraPS can readily analyse continuous-time data (Aga et al., 2024; Renz et al., 142 2025). Rather, our approach uses the relative evolutionary timings from the source data phylogeny 143 (for example, Fig. 2A) to report the orderings of events: which acquisitions precede or follow which 144 others, given that characters may influence each other, with the timescale of these events (see 145 Discussion). 146 While Figs. 2A-C demonstrate the approach for a single country, Fig. 2D gives the average 147 acquisition probabilities across all countries in our dataset, and Supp. Fig. 2 gives the 148 corresponding summary plots for every country in the dataset. Resistance to some antibiotic groups 149 .CC-BY 4.0 International licenseperpetuity. It is made available under a preprint (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in The copyright holder for thisthis version posted September 21, 2025. ; https://doi.org/10.1101/2025.09.20.677523doi: bioRxiv preprint are consistently observed across all countries: following very common presence of Bla-chr and 150 SHV-m, resistance characters against betalactams and extended-spectrum betalactams (ESBLs) 151 (Bla-a, Bla_ESBL-a) are often acquired earlier, while resistance to colistin (drug of last resort for 152 multi-drug resistant Kp; Col-a and Col-m) is acquired later. Country-specific diUerences exist in the 153 evolutionary acquisition of other characters, including resistance to carbapenems (Bla_Carb-a). 154 In several instances, bimodality in the inferred evolutionary dynamics is observed. That is, a given 155 character may be acquired at an earlier stage or a later stage in the evolutionary process, but not at 156 intermediate stages. Such bimodality is characteristic of several diUerent evolutionary pathways 157 existing in a system. Examples can readily be seen in Fig. 2B for Romania and Fig. 2D for the global 158 mean behaviour, where SHV-m is likely acquired either very early or at more intermediate stages, 159 with little probability in between. These competing dynamics are directly visible in Fig. 2C, where 160 distinct, canalised pathways exist: one involving SHV-m as the first acquisition, one involving 161 Bla_chr as the first acquisition and SHV-m only following rather later. Across countries (Supp. Fig. 2), 162 bimodality is observed in several characters: some examples are in Bla_chr (Tanzania, Zambia); Bla-163 a (Ireland, S Korea; Bla_ESBL-a (Laos); and SHV-m (Romania and at least 15 other countries across 164 GDB regions). The bimodality in the mutational and chromosomal characters reflects the 165 diUerences in evolutionary dynamics between horizontal-gene-transfer-mediated acquisitions and 166 the other KpAMR characters we consider. 167 Across subsets of countries, interactions between our KpAMR resistance characters are observed 168 (Supp. Fig. 3). For example, in around a quarter of cases, the acquisition of tigecycline resistance 169 (Tgc_a) is inferred to promote the acquisition of glycopeptide resistance (Gly_a), and the acquisition 170 of sulfonamide resistance (Sul_a) is inferred to promote the acquisition of trimethoprim resistance 171 (Tmt_a). In some cases, the acquisition of MLS resistance (MLS_a) is inferred to repress the 172 acquisition of carbapenem resistance (Bla_carb_a), suggesting possible interactions that could be 173 clinically exploited in combination therapies (see Discussion). 174 However, although Supp. Fig. 2 gives a detailed “roadmap” of KpAMR evolution across countries, it 175 is hard to immediately extract global-scale insights from this representation. We therefore next 176 considered how to compare these results in a reduced-dimensionality picture. 177 Global consistency in evolutionary acquisition of resistance to a subset of drug families 178 To proceed, we used principal components analysis (PCA) to embed the high-dimensional set of 179 posterior probabilities in a 3D space (Fig. 3). Here, the axes of highest variance in the inferred 180 evolutionary pathways are identified and used to naturally lay out the individual country outputs. 181 Comparison of country-by-country behaviour along these axes can then reveal major sources of 182 variability (and similarity) in inferred KpAMR evolution across countries. 183 Inferred evolutionary pathways for Kp across countries consistently involve early acquisition of 184 characters conferring resistance to aminoglycosides, sulfonamides, trimethoprim, and betalactams 185 (Agly, Sul, Tmt and Bla characters). Aminoglycosides and sulfonamides were both discovered and 186 entered clinical use early in the mid 19th century. Resistance properties inferred to be acquired at 187 later orderings include extended spectrum beta-lactam and phenicol resistance, acquired via 188 mobile genetic elements (Bla_ESBL-a, Phe-a). Acquisitions inferred to be acquired late across 189 countries include inhR-driven resistance (Bla_inhR-a, Bla_ESBL_inhR-a), and resistance to later- 190 and final-line drugs including colistin, fosfomycin, tigecyclin, and glycopeptides (Col-a, Fcyn-a, Tgc-191 a, Gly-a). These form a set of “globally consistent” evolutionary observations. 192 .CC-BY 4.0 International licenseperpetuity. It is made available under a preprint (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in The copyright holder for thisthis version posted September 21, 2025. ; https://doi.org/10.1101/2025.09.20.677523doi: bioRxiv preprint 193 Figure 3. Global structure of inferred pathways of AMR evolution. (A-B) Principal components analysis 194 (PCA) of country-specific “bubble plots” as in Fig. 2B, reflecting inferred evolutionary dynamics of KpAMR in 195 di\erent countries. Individual countries are plotted on PCAs 1-2 (A) and 2-3 (B), individually and with ellipses 196 demarcating their corresponding GBD region (acronyms below). (C-D) Example “bubble plots” from countries 197 reflecting the extremes of PCA1 (C) and PCA2 (D). In (C), Venezuela reflects poorly characterised evolutionary 198 dynamics, with relatively uniform probabilities everywhere; Nigeria reflects precisely characterised 199 evolutionary dynamics, with individual characters having highly likely, precisely specified orderings. In (D), 200 country-specific di\erences in ordering emerge, including in those characters highlighted with boxes. For 201 example, KpAMR evolution in Argentina is likely to involve Flq-m acquisition earlier than other characters; 202 evolution in Senegal is likely to involve Flq-m acquisition later than other characters. GBD regions: CEEECA, 203 Central Europe, Eastern Europe, Central Asia; HI, High-income; LAC, Latin American and Caribbean; NAME, 204 North Africa and Middle East; SA, South Asia; SEEAO, South-east and East Asia and Oceania; SSAf, Sub-205 Saharan Africa. Zanzibar (dark outline) is from data newly sequenced in this study. Character names given in 206 Fig. 1 caption. 207 208 Country-to-country di@erences in evolutionary acquisition of resistance to carbapenems, 209 fluoroquinolones, tetracycline, and others 210 The general patterns above appear across evolutionary dynamics across countries and regions, with 211 the first principal component of inferred Kp variability (accounting for 23% of the overall variance) 212 simply reporting the precision with which these patterns can be characterised (Fig. 3C, Fig. 4). In 213 countries supporting less precise inference (for example, with limited datasets), the expected 214 ordering of these characters is more uniform and closer to the average null hypothesis of all 215 characters behaving equally; in countries with more precise output the expected orderings diverge 216 to early and late values (Fig. 3C, Fig. 4). 217 218 .CC-BY 4.0 International licenseperpetuity. It is made available under a preprint (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in The copyright holder for thisthis version posted September 21, 2025. ; https://doi.org/10.1101/2025.09.20.677523doi: bioRxiv preprint 219 Figure 4. KpAMR characters determining aspects of global variability in evolutionary dynamics. Each 220 country has an expected acquisition ordering for each KpAMR character. The plots show how these patterns of 221 expected ordering are linked to the principal components of global KpAMR evolution variability in Fig. 3: PCA1 222 (top, grey); PCA2 (centre, red); and PCA3 (bottom, blue). (A) “Globally consistent” characters: those that 223 covary tightly with PCA1. All characters here have either early or late expected orderings across all countries 224 (either left or right of the central vertical line). Their position on PCA1 reflects the precision of the inference 225 output. Less precise outputs (low PCA1 values, like Venezuela in Fig. 3C) will have orderings close to the 226 central “null”; more precise outputs (high PCA1 values, like Nigeria in Fig. 3C) will have more divergent “early” 227 and “late” orderings towards the edges of the plot. (B-C) “Globally divergent” characters: those that covary 228 more tightly with PCA2 (B) or PCA3 (C). Expected orderings in these character sets span a wider range, 229 showing sources of variability in the evolutionary dynamics that are independent of PCA1 (and hence less 230 influenced by observational uncertainty). An alternative visualisation is given in Supp. Fig. 4. Character names 231 given in Fig. 1 caption. 232 233 Once variability in this precision is accounted for (by considering the first principal component), 234 variability corresponding to other diUerences in the inferred dynamics can be extracted. Several 235 characters covary comparatively little with the (precision-aligned) PCA1 and instead show stronger 236 covariance with PCA2 or PCA3, suggesting that there is genuine country-to-country variability in 237 these characters after controlling for diUerences in observational noise (Fig. 3D, Fig. 4, Supp. Fig. 4). 238 PCA2 (accounting for 15% of overall variance) displays strong covariance with carbapenem, 239 fluoroquinolone, MLS, and tetracycline resistance and OMP mutation (Bla_Carb-a, Flq-m, MLS-a, 240 Tet-a, Omp-m,). PCA3 (accounting for 8% of overall variance) displays covariance with Flq-a, Rif-a, 241 SHV-m. None of these “globally divergent” characters display the “null-to-extreme” divergence seen 242 in those characters most closely related to PCA1, suggesting genuine country-to-country 243 diUerences in evolutionary dynamics, rather than observational diUerences, are responsible for this 244 variability. 245 .CC-BY 4.0 International licenseperpetuity. It is made available under a preprint (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in The copyright holder for thisthis version posted September 21, 2025. ; https://doi.org/10.1101/2025.09.20.677523doi: bioRxiv preprint 246 Figure 5. Geographical variation in inferred AMR evolutionary pathways. (A-B) Distributions of country-247 specific PCA projection values across di\erent GBD regions for (A) PCA2 and (B) PCA3. In (A), sub-Saharan 248 African countries show systematically higher positions in PCA2 than other regions. (C) Inferred acquisition 249 orderings of KpAMR characters by GBD region. The left five characters covary with PCA2; the right three covary 250 with PCA3. Carbapenem resistance and OMP mutations are inferred to be acquired substantially later in sub-251 Saharan Africa; rifampicin and fluoroquinolone resistance are inferred to be acquired substantially earlier in S, 252 SE, E Asia and Oceania. New Zanzibar data from this study is shown by large circles. Character names given in 253 Fig. 1 caption. 254 Geographical correlates with Kp evolutionary dynamics 255 To explore the sources of this evolutionary variability in more detail, we asked whether geographical 256 region was linked to diUerent variants of the inferred evolutionary dynamics (after accounting for 257 observational diUerences). We first asked whether countries from diUerent regions showed 258 systematic diUerences in their PCA2 and PCA3 projections – corresponding to diUerent orderings of 259 the “globally divergent” country-specific characters above (Fig. 5). We found that sub-Saharan 260 Africa was a dramatic outlier in PCA2, involving notably later acquisition of carbapenem and 261 fluoroquinolone resistance and OMP mutation than the other regions (Fig. 5A, C). This ordering 262

agrees with the continuous-time history of treatment in the region, where carbapenem and 263 fluoroquinolone use was established later (especially for children) than in other regions. Other 264 regions were broadly consistent in PCA2 but showed substantial diUerences in PCA3. These 265 diUerences correspond mainly to acquired resistance to fluoroquinolone and rifampicin, with 266 substantially earlier inferred acquisition in South-Southeast-East Asia and Oceania, and 267 substantially later in Central and Eastern Europe, Central Asia, and sub-Saharan Africa (Fig. 5B, C). 268 .CC-BY 4.0 International licenseperpetuity. It is made available under a preprint (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in The copyright holder for thisthis version posted September 21, 2025. ; https://doi.org/10.1101/2025.09.20.677523doi: bioRxiv preprint 269 Figure 6. Connecting AMR character acquisition propensity with drug treatment regimens across 270 countries. Median (across observed years 2016-2021) per-capita drug use in di\erent drug classes (J codes), 271 with expected acquisition ordering of associated resistance characters. “Globally consistent” characters Tet-272 a, Bla_inhR-a, Sul-a, AGly-a were associated with PCA1 in the global analysis and display little correlation 273 between evolutionary ordering and drug use levels. “Globally divergent” characters Bla_Carb-a and MLS-a 274 were associated with PCA2 in the global analysis and display stronger correlations, with earlier evolution of 275 resistance linked to higher levels of drug use. Character names given in Fig. 1 caption. 276 Policy predictors of variance in AMR evolutionary dynamics 277 We next asked whether known policy of antibiotic use influence the inferred evolutionary dynamics 278 from our analysis. We gathered statistics on drug use from 2016-2021 in countries available from 279 the WHO GLASS surveillance programme (Global Antimicrobial Resistance and Use Surveillance 280 System (GLASS) Report 2022, 2022) and looked at correlations between per-capita drug use 281 statistics and the inferred evolutionary dynamics of resistance to those drugs (Fig. 6). We found that 282 for those drug types associated with PCA1 (where variability in inferred dynamics is largely 283 determined by observation uncertainty) there was little correlation between drug use and the 284 ordering of resistance acquisition. However, for carbapenems, associated with PCA2 (more non-285 observational country-to-country variability), there was a link between drug usage statistics and 286 inferred resistance evolution, with resistance evolving earlier in countries with high drug use rates. 287 The direction of this correlation was also observed in MLS, also associated with PCA2, though this 288 relationship did not show statistical support at the p < 0.05 level. This pattern of correlations is 289 consistent with the picture above: acquisition of some (PCA1) KpAMR characters is relatively 290 constant across countries and circumstances, while other (PCA2) characters are more plastic and 291 depend more tightly on country- and region-specific drivers. 292 .CC-BY 4.0 International licenseperpetuity. It is made available under a preprint (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in The copyright holder for thisthis version posted September 21, 2025. ; https://doi.org/10.1101/2025.09.20.677523doi: bioRxiv preprint Testing predictions of AMR evolutionary dynamics with newly sequenced, cross-decadal 293 data 294 One motivating application of evolutionary inference to AMR problems is the ability to predict the 295 next steps in evolutionary pathways (Renz et al., 2024). A question of potentially direct clinical 296 relevance is: given a clinically observed isolate with a given set of AMR characters, which 297 character(s) will likely evolve next – and can treatment guidelines be adapted accordingly? 298 Previous studies of AMR with HyperTraPS have demonstrated its ability to predict properties of 299 evolutionary transitions in a withheld subset of the overall dataset, following the common training-300 test machine learning paradigm (Aga et al., 2024). Here, we are in the unique position of using 301 independently obtained and Kp genome data to test the predictions of the evolutionary models 302 trained above. We sequenced Kp isolates from patients with bloodstream infections in 2001-2 and 303 2017-8 in Tanzania and 2015-6 in Zanzibar, collected in previous studies (Blomberg et al., 2007; 304 Moyo et al., 2020; Onken et al., 2015, 2024) (Fig. 7; see Methods). Bioinformatic analysis placed 305 these new genomes in a phylogenetic tree with a subset of existing genomes from the training data, 306 including their AMR profiles (Fig. 7A, Supp. Fig. 5). Using only transitions independent of those in the 307 training data (see Methods), we tested the predictions of the HyperTraPS model trained on data from 308 Tanzania. Specifically, for each new transition, we queried the trained model about likely next steps 309 from the ancestral state, and recorded the predicted ranks of the next steps actually present in the 310 transition (Fig. 7B). By far the most common outcome was that the true next step was ranked first 311 (predicted to be most likely); compared to predictions from an untrained model, the trained model 312 was systematically much better at predicting the true next steps. 313 This analysis demonstrates that within-country predictions of future evolutionary dynamics are 314 possible given a trained HyperTraPS model. We also asked whether the model could predict 315 evolutionary dynamics at the regional scale. We use Zanzibar as a test case here – a region of 316 Tanzania that was not explicitly represented in the original Kleborate training dataset. Our findings 317 would predict KpAMR in Zanzibar to fall into the sub-Saharan Africa set of behaviours, with high 318 values on PCA2 and associated signatures of late carbapenem resistance and others (Fig. 5). To this 319 end, we applied HyperTraPS inference to newly sequenced isolates from Zanzibar. The resulting 320 inferred dynamics clearly place Zanzibar within the set of sub-Saharan Africa behaviours (Fig. 3A-B), 321 with high values on PCA2 arising from later acquisition of carbapenem resistance, and high values 322 on PCA3 arising from later acquisition of rifampicin resistance (Fig. 5). 323 We also explored the representation of diUerent KpAMR characters in the snapshot genome 324 sequences alone, without connecting with evolutionary dynamics. Fig. 7C shows the relative 325 prevalence of acquired KpAMR characters in existing and newly-sequenced data, not accounting for 326 phylogenetic relatedness (and without representation of an explicit non-AMR control group). With 327 this sampling, most characters showed an increase from 2001-2 to 2015-6, particularly in mobile 328 genetic acquisition of fluoroquinolone resistance, while the relative prevalence of acquired 329 phenicol resistance and fluoroquinolone resistance mutations decreased. Proportions of these 330 characters were in broad agreement with known coarser-grained observations from Tanzania 331 (Mapunjo et al., 2025; Sangeda et al., 2021) and existing Tanzanian genomes from Pathogenwatch, 332 but we observed more acquired rifampicin resistance and lower acquired phenicol resistance than 333 in existing samples. Our samples from Zanzibar showed lower relative proportions of 334 aminoglycoside and β-lactam resistance, and higher prevalence of outer membrane protein 335 mutations. 336 .CC-BY 4.0 International licenseperpetuity. It is made available under a preprint (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in The copyright holder for thisthis version posted September 21, 2025. ; https://doi.org/10.1101/2025.09.20.677523doi: bioRxiv preprint 337 Figure 7. Prediction of new evolutionary transitions in KpAMR. (A) Pathogenwatch (grey) and newly 338 sequenced 2001-2 (blue) and 2017-8 (red) genome properties from Tanzania, connected by a phylogeny 339 estimated from average nucleotide identity (ANI). (B) For each independent transition in the set of newly 340 sequenced observations, we used HyperTraPS trained on existing data to predict which characters would 341 likely evolve next given the ancestral state. The bars show the predicted rankings of the characters that were 342 truly observed to evolve in the data. The trained model (blue-green) overwhelmingly predicted top rankings for 343 the characters that were truly acquired; an untrained model (red) predicts intermediate rankings. Character 344 names given in Fig. 1 caption. (C) Proportion of Tanzanian isolates in Kleborate and new datasets (uncorrected 345 for relatedness) where di\erent KpAMR characters are present. 346 347

348 Powerful data-driven approaches to study the evolution of AMR are emerging in parallel with large 349 genomic datasets (Argimón et al., 2021; Hetland et al., 2025; Holt et al., 2015; Munk et al., 2022; 350 Wyres et al., 2020). Most of these studies focus on the current genomic structure of pathogen 351 populations. In this research we have attempted to take a complementary perspective, using large-352 scale genome data while focusing on the (inferred) evolutionary dynamics that produce these 353 genomic patterns. We believe that this new way of working has substantial power to support other 354 data-driven approaches, particularly in the prediction of unseen and future evolutionary behaviours 355 of pathogens (Renz et al., 2025). 356 The core technology we use, HyperTraPS, is an instance of evolutionary accumulation modelling or 357 EvAM (Diaz-Uriarte & Herrera-Nieto, 2022; Diaz-Uriarte & Johnston, 2024). EvAM has previously 358 been used to learn AMR pathways, including in HIV (Beerenwinkel et al., 2005) and multidrug 359 resistance in tuberculosis (Aga et al., 2024; Greenbury et al., 2020; Moen & Johnston, 2023; Renz et 360 .CC-BY 4.0 International licenseperpetuity. It is made available under a preprint (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in The copyright holder for thisthis version posted September 21, 2025. ; https://doi.org/10.1101/2025.09.20.677523doi: bioRxiv preprint al., 2024). One question that has remained is to what extent such results reflect universal, general 361 behaviors, versus a response to the specific selective pressures arising from a given country’s drug 362 use (for example). Our global comparison suggests a form of answer to this question. A subset of 363 the KpAMR characters we consider appear to behave similarly across countries regardless of 364 specific selective pressure, drug regimes, and other diUerences (those aligned with PCA1). Another 365 set of characters (those aligned with PCA2-3) evolves in a way that is more dependent on country- 366 and region-specific influences, which we have shown include public health superregions and drug 367 use levels – although certainly other factors, from population change to the prevalence of other 368 diseases, will also influence these dynamics. While these results are specific to Kp, and the 369 characters we consider, this combination of global and country-specific behaviors may reasonably 370 hold across other pathogens. 371 Inferring evolutionary dynamics of multiple, coupled traits across groups is a challenging problem. 372 There will be countless, multiscale selective influences that impact the evolution of specific 373 lineages of Kp around the globe. Our perspective attempts to “coarse-grain” away from these 374 (unobservable) specific pressures and consider their aggregated eUect on large-scale genomic 375 properties in populations. This perspective is not without challenges. Some of these are: (a) The 376 acquisition of many of the KpAMR characters we consider is reversible (for example, mediated by 377 mobile genetic elements), which challenges an assumption of our modelling. (b) We assume that 378 genomes from a country can be reasonably placed on their own phylogeny. (c) We use datasets 379 which likely exhibit bias towards AMR presence: in most circumstances, due to timing and/or 380 human research focus, genome databases will be enriched for cases with AMR characters acquired. 381 And (d), in most of this analysis, we use relative evolutionary timings rather than the absolute real-382 world (day-month-year) timings associated with observed genomes. For (a), we have shown with a 383 synthetic control study that violating the irreversibility assumption does not dramatically challenge 384 our findings. (We note, in parallel to this point, that simultaneous acquisition of multiple features is 385 very compatible with our approach; the resulting output simply places equal probabilities over all 386 compatible ordering pathways.) For (b), coupling between countries (through, for example, host or 387 pathogen migration) will not aUect the interpretation of our results, which make statements 388 conditioned on the set of pathogens reported within a country. And likewise for (c), the enrichment 389 of AMR characters in our source data does not introduce systematic bias in our results, which report 390 on accumulation of AMR characters, controlling for the existing AMR profile of an ancestral state in 391 each case. Our phylogenetic reconstruction bridges the gap between observed, AMR-enriched 392 genomes and their unobserved, less-enriched ancestors, ensuring that the whole spectrum of 393 character profiles is captured in the underlying model. 394 For (d), the precise timescales of the evolutionary process are in principle addressable using a 395 version of our approach called HyperTraPS-CT (continuous time) (Aga et al., 2024), but this requires 396 a consistent estimation of divergence times between all isolates in the dataset. Importantly, our 397 approach here does not completely neglect a timescale – it is implicit in the phylogenetic 398 relationships included in the source data, and this inclusion both accounts for the temporal 399 ordering of events and “pseudoreplication” due to the relatedness of sampled genomes (Boyko & 400 Beaulieu, 2023; Dewar et al., 2025; Renz et al., 2025). However, this approach based on relative 401 evolutionary timings means that we cannot directly tie diUerences in inferred evolutionary 402 behaviour to known real-world historical events shaping AMR profiles: shifts in drug policy, for 403 example, or outbreaks or other events changing the eUective size of the evolving population. While 404 our results, focusing on the ordering of acquisition events, already demonstrate descriptive and 405 predictive power in KpAMR evolution across countries, a continuous time picture on a smaller scale 406 in the future (perhaps a single country or small set) would enable finer-grained linking between 407 model prediction and the real-world details of policy and disease history. 408 In learning estimates for evolutionary dynamics, HyperTraPS also infers potential interactions 409 between characters, including whether the acquisition of one character makes acquiring another 410 character less likely. Such “repressive interactions” , if truly present, could form the basis of 411 combination therapies: resistance to drug X makes resistance to drug Y less likely. These 412 .CC-BY 4.0 International licenseperpetuity. It is made available under a preprint (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in The copyright holder for thisthis version posted September 21, 2025. ; https://doi.org/10.1101/2025.09.20.677523doi: bioRxiv preprint interactions correspond to the phenomenon of collateral sensitivity in AMR, where resistance to one 413 drug induces sensitivity to another (Pál et al., 2015; Szybalski & Bryson, 1952): for example, in 414 Escherichia coli, beta-lactamase expression is linked to heightened sensitivity to colistin (among 415 other) (Herencias et al., 2024). While not observed consistently throughout all countries in our 416 dataset, a subset of countries showed statistical support for some such “repressive interactions” 417 (Supp. Fig. 3). These include, for example, MLS resistance repressing the acquisition of carbapenem 418 resistance. Such indications could suggest more detailed followup investigation to verify and 419 explore mechanisms behind these possibly exploitable collateral interactions. 420 In conclusion, we hope to show that evolutionary accumulation modelling – inferring the historic 421 evolutionary pathways of AMR acquisition, in addition to considering contemporary properties – can 422 complement and enhance established ways of working in AMR genomic analysis. By considering 423 the dynamics by which global KpAMR patterns have become established, a collection of research 424 directions are opened or expanded, including the similarities and diUerences between countries, 425 links with social and environmental covariates, and predictions of future dynamics. We anticipate 426 that these methods may readily be applied to provide insight into AMR evolution in other pathogens 427 in future. 428

429 Existing Klebsiella data acquisition 430 We obtained all 47 721 Klebsiella pneumoniae sensu stricto records from pathogenwatch as of 431 March 2024 (Argimón et al., 2021). The most common isolation years were 2015-2020, with a peak 432 at 2018 (Fig. 1A). We used AMR features (genes and mutations) reported by drug class as reported 433 by Kleborate version 2.3 (Lam et al., 2021). Kleborate groups genes or mutations known to confer 434 resistance to clinically relevant antibiotic groups and related phenotypes. Betalactamases are 435 further grouped by enzyme activity (Lam et al., 2021). We consider a set of L = 22 drug classes for 436 each genome, dichotomizing Kleborate output by presence of any gene or mutation related to each 437 group. That is, if there is any gene or mutation present for a resistance phenotype, the isolate is 438 considered resistant, otherwise it is considered susceptible. This assumption holds for most 439 resistance phenotypes considered in this study. However, fluoroquinolone resistance (Flq-m) often 440 require multiple mutations to cause resistance. Thus, the presence of these characters cannot be 441 interpreted as a resistance phenotype, but rather an increased MIC (minimum inhibitory 442 concentration). The isolates were divided by country based on metadata, and countries were 443 grouped into superregions as defined by the Global Burden of Disease regional classification 444 system (Rudd et al., 2020). A coarse-grained phylogeny was generated via LIN-codes (Life 445 identification numbers) using lincoding.py supplied by (Hennart et al., 2022). 2 102 genomes were 446 excluded due to missing metadata. Any country with less than two genomes available were 447 excluded due to the inability to construct a tree. The resistance profiles and phylogenetic tree were 448 combined into a phylogenetic tree annotated with resistance to drug classes. 449 Evolutionary pathway inference with hypercubic transition path sampling (HyperTraPS) 450 HyperTraPS (Aga et al., 2024) models an “evolutionary space” containing every possible state of a 451 system involving L binary characters, then infers the probabilities of (Markovian) transitions between 452 states in this space (here, involving the ordered accumulation of diUerent AMR characters) that is 453 most compatible with observed data (Greenbury et al., 2020; Johnston & Williams, 2016). A recent 454 overview of the method in an AMR context is given in (Renz et al., 2024). The source data is a 455 collection of length-L “barcodes” labelling presence or absence of our characters -- resistance to 456 drug classes as defined in Kleborate -- on the tips of a phylogeny. Ancestral state reconstruction 457 (here, assuming that AMR resistance character accumulation is rare and irreversible, so that an 458 ancestor possesses a character if and only if all descendants possess it) is used to infer ancestral 459 states, and thereby construct a set of transitions from the data. These transitions are the input data 460 for the HyperTraPS algorithm. The target of inference is an L x L matrix 𝜃, where 𝜃!! is the base rate 461 for acquiring character i, and 𝜃!" is the influence that acquisition of character j has on the base rate 462 .CC-BY 4.0 International licenseperpetuity. It is made available under a preprint (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in The copyright holder for thisthis version posted September 21, 2025. ; https://doi.org/10.1101/2025.09.20.677523doi: bioRxiv preprint of acquisition of character i (as in mutual hazard networks (Schill et al., 2020)). In this way, positive 463 and negative interactions between the acquisition of diUerent AMR characters is supported, so that 464 (for example) acquisition of resistance to drug X can make acquisition of resistance to drug Y more 465 (or less) likely. The output of the inference process is, most broadly, a transition matrix 𝜆, with 466 elements 𝜆#→% describing the probability with which an isolate in state a will next transition to state 467 b. For each of the 102 countries, territories and areas we ran HyperTraPS on three diUerent random 468 seeds to ensure that the results were consistent across runs (Supp. Fig. 2). Simulations were run for 469 106 steps with 103 random walkers for sampling, employing a penalised likelihood. The US genomes 470 were sub-sampled to a subset of 1000 genomes; three diUerent subsamples clustered closer 471 together than with any other countries in the PCA plot. 472 Synthetic control study for inference of reversible characters 473 To determine the eUect of our assumption of irreversible accumulation of AMR characters (given 474 that the loss of plasmids is not uncommon in Kp evolution (Holt et al., 2015)), we constructed a 475 synthetic control study. A random phylogeny over 128 tips was constructed using a birth-death 476 process with birth rate 1 and death rate 0.5 (chosen empirically to roughly match the phylogenies 477 from our case studies) using phangorn (Schliep, 2011). A synthetic accumulation process reflecting 478 a single evolutionary pathway was simulated on this phylogeny. Reversibility was captured by 479 modelling random loss of the first L/2 characters as a Poisson process with a given rate parameter. 480 The same HyperTraPS pipeline (assuming irreversibility), including ancestral state reconstruction 481 and pathway inference, was run on both the pre-loss and post-loss datasets, and the outputs 482 compared (Supp. Fig. 1). 483 Global structure in inferred pathways 484 The outputs of HyperTraPS can be summarized as a matrix P, where element Pij gives the posterior 485 probability that character i is acquired at ordering j in a putative evolutionary process from an 486 ancestor with no AMR characters towards a final state with all AMR characters. For example, P32 487 gives the probability that the third character is the second to be acquired in such an evolutionary 488 process. We obtained this matrix for each country in our dataset, then performed principal 489 components analysis (PCA) on the set of matrices (Williams et al., 2013). We observed the structure 490 of matrices with “bubble plots” , where an array of points is plotted with areas proportional to Pij. 491 Drug usage data 492 We collected antibiotic consumption data from 57 countries, territories and areas enrolled in the 493 WHO Glass surveillance program in the period 2016 to 2021 (Global Antimicrobial Resistance and 494 Use Surveillance System (GLASS) Report 2022, 2022). Overlaying the consumption data with the 495 countries that had 3 or more available genomes reduced the number to 53. We use the median 496 consumption statistics over the available time points for each country. 497 New Klebsiella data acquisition 498 The 79 newly sequenced Kp isolates in this article are derived from a collection of studies 499 performed in Tanzania and Zanzibar over the last three decades. These comprise bacterial isolates 500 from (a) blood samples from pediatric patients with bloodstream infections in Dar-es-Salaam, 501 Tanzania in 2001-2002; (b) blood samples from adult patients with bloodstream infections in Dar-502 es-Salaam, Tanzania in 2017-2018; (c) fecal samples from the patients in (b); and (d) blood samples 503 from adult and pediatric patients with bloodstream infections in Mnazi Mmoja hospital, Zanzibar in 504 2015-2016. These cultures were sequenced by MicrobesNG (MicrobesNG, Birmingham, UK) using 505 Illumina HiSeq technology. Assembled contigs were provided by the sequencing service, which 506 formed the initial step for the bioinformatics pipeline below. 507 New Klebsiella data analysis 508 The analysis pipeline begins with the contigs from the initial sequencing and assembly process. We 509 first used dnadiQ, which wraps nucmer from the MUMmer 3.0 package (Kurtz et al., 2004), to 510 .CC-BY 4.0 International licenseperpetuity. It is made available under a preprint (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in The copyright holder for thisthis version posted September 21, 2025. ; https://doi.org/10.1101/2025.09.20.677523doi: bioRxiv preprint scaUold the contigs from a given isolate on the Klebsiella reference genome 511 GCA_000240185.2_ASM24018v2 (Liu et al., 2012) and report statistics of the alignment. We 512 retained only isolates with <20% unaligned bases as instances of Klebsiella. Combining the newly 513 sequenced isolates with the existing Tanzanian genomes from Pathogenwatch, we then ran pairwise 514 dnadiQ across the combined set. We extracted the average identity across the aligned regions for 515 each pair and recorded this as the ANI (average nucleotide identity). We then used d = 1-ANI as a 516 distance measure for phylogenetic tree estimation using the unweighted pair group method with 517 arithmetic mean (UPGMA) (Sokal & Michener, 1958) implemented in the phangorn (Schliep, 2011) 518 package in R (R Core Team, 2022). We confirmed that an alternative method, neighbourhood joining 519 (Saitou & Nei, 1987), did not give qualitatively diUerent results for the predictions of AMR evolution. 520 We used AMR profiles of existing and new isolates to reconstruct ancestral states across this 521 complete tree, then identified the set of subtrees that had only newly-sequenced isolates at their 522 tips. As transitions within such subtrees are independent of transitions in the (Pathogenwatch) 523 training data set, we used this set of transitions as our independent test set. 524 Testing predicted dynamics from the trained model on new data 525 We queried the model fitted from the training data to predict the likely next states from each 526 precursor state in this independent set of transitions (Aga et al., 2024; Renz et al., 2024). We ranked 527 the possible next steps from each state from most likely to least likely under the fitted model. We 528 then recorded the ranks corresponding to the actual step(s) corresponding to the observed 529 transition, and compared these to a null model where each possible transition could occur with 530 equal probability. 531 Data and code availability 532 The code base for this study is publically available at https://github.com/StochasticBiology/kp-533 evolution-inference. In addition to the software referenced above, we use R (R Core Team, 2022) 534 with libraries hypertrapsct for evolutionary inference (Aga et al., 2024), dplyr (Wickham et al., 2023), 535 tidyverse (Wickham et al., 2019), countrycode (Arel-Bundock et al., 2018) for data manipulation, 536 lme4 (Bates et al., 2015) for statistics, phytools (Revell, 2012) for phylogenetic analysis, and ggplot2 537 (Wickham, 2016), ggpubr (Kassambara, 2020), ggbeeswarm (Clarke et al., 2016), ggrepel 538 (Slowikowski, 2021), ggupset (Ahlmann-Eltze, 2025) for visualisation. We also use a Python script for 539 LINcoding by Melanie Hennart, which uses numpy (Harris et al., 2020). Authors’ note: we are in the 540 process of uploading the new Kp sequences to a repository; a link will be provided when this 541 process is complete. In the meantime the raw FASTA files (which are the input for the pipeline) can 542 be downloaded from https://osf.io/36r45 . 543

544 This project has received funding from the European Research Council (ERC) under the European 545 Union’s Horizon 2020 Research and Innovation Programme [Grant agreement No. 805046 546 (EvoConBiO) to I.G.J.]. This project was supported by the Trond Mohn Foundation [project HyperEvol 547 under grant agreement No. TMS2021TMT09], through the Centre for Antimicrobial Resistance in 548 Western Norway (CAMRIA) [TMS2020TMT11], This project has received support from ERA-Net: JPI-549 AMR STRESST project (NFR333432). The authors are grateful to the CAMRIA collaboration and the 550 Stochastic Biology Group at UiB for useful discussions. 551 552

553 Aga, O. N., Brun, M., Giannakis, K., Dauda, K. A., Diaz-Uriarte, R., & Johnston, I. G. (2024). HyperTraPS-CT: 554 Inference and prediction for accumulation pathways with flexible data and model structures. PLOS 555 Computational Biology, 20(9), e1012393. 556 Ahlmann-Eltze, C. (2025). ggupset: Combination Matrix Axis for ‘ggplot2’ to Create ‘UpSet’ Plots [R package]. 557 Arel-Bundock, V ., Enevoldsen, N., & Yetman, C. (2018). countrycode: An R package to convert country names 558 and country codes. Journal of Open Source Software, 3(28), 848. https://doi.org/10.21105/joss.00848 559 .CC-BY 4.0 International licenseperpetuity. It is made available under a preprint (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in The copyright holder for thisthis version posted September 21, 2025. ; https://doi.org/10.1101/2025.09.20.677523doi: bioRxiv preprint Argimón, S., David, S., Underwood, A., Abrudan, M., Wheeler, N. E., Kekre, M., Abudahab, K., Yeats, C. A., 560 Goater, R., Taylor, B., Harste, H., Muddyman, D., Feil, E. J., Brisse, S., Holt, K., Donado-Godoy, P ., 561 Ravikumar, K. L., Okeke, I. N., Carlos, C., … NIHR Global Health Research Unit on Genomic 562 Surveillance of Antimicrobial Resistance. (2021). Rapid Genomic Characterization and Global 563 Surveillance of Klebsiella Using Pathogenwatch. Clinical Infectious Diseases, 73(Supplement_4), 564 S325–S335. https://doi.org/10.1093/cid/ciab784 565 Bates, D., Mächler, M., Bolker, B., & Walker, S. (2015). Fitting Linear Mixed-E\ects Models Using lme4. Journal 566 of Statistical Software, 67, 1–48. https://doi.org/10.18637/jss.v067.i01 567 Beerenwinkel, N., Däumer, M., Sing, T., Rahnenführer, J., Lengauer, T., Selbig, J., Ho\mann, D., & Kaiser, R. 568 (2005). Estimating HIV evolutionary pathways and the genetic barrier to drug resistance. The Journal 569 of Infectious Diseases, 191(11), 1953–1960. 570 Blomberg, B., Manji, K. P ., Urassa, W. K., Tamim, B. S., Mwakagile, D. S. M., Jureen, R., Msangi, V ., Tellevik, M. 571 G., Holberg-Petersen, M., Harthug, S., Maselle, S. Y., & Langeland, N. (2007). Antimicrobial resistance 572 predicts death in Tanzanian children with bloodstream infections: A prospective cohort study. BMC 573 Infectious Diseases, 7, 43. https://doi.org/10.1186/1471-2334-7-43 574 Boyko, J. D., & Beaulieu, J. M. (2023). Reducing the biases in false correlations between discrete characters. 575 Systematic Biology, 72(2), 476–488. 576 Casali, N., Nikolayevskyy, V ., Balabanova, Y ., Harris, S. R., Ignatyeva, O., Kontsevaya, I., Corander, J., Bryant, J., 577 Parkhill, J., Nejentsev, S., Horstmann, R. D., Brown, T., & Drobniewski, F. (2014). Evolution and 578 transmission of drug-resistant tuberculosis in a Russian population. Nature Genetics, 46(3), Article 3. 579 https://doi.org/10.1038/ng.2878 580 Clarke, E., Sherrill-Mix, S., & Dawson, C. (2016). ggbeeswarm: Categorical scatter (violin point) plots. CRAN: 581 Contributed Packages. https://cir.nii.ac.jp/crid/1360019690351052672 582 Dewar, A. E., Belcher, L. J., & West, S. A. (2025). A phylogenetic approach to comparative genomics. Nature 583 Reviews Genetics, 26(6), 395–405. https://doi.org/10.1038/s41576-024-00803-0 584 Diaz-Uriarte, R., & Herrera-Nieto, P . (2022). EvAM-Tools: Tools for evolutionar y accumulation and cancer 585 progression models. Bioinformatics, 38(24), 5457–5459. 586 Diaz-Uriarte, R., & Johnston, I. G. (2024). A picture guide to cancer progression and monotonic accumulation 587 models: Evolutionary assumptions, plausible interpretations, and alternative uses 588 (arXiv:2312.06824). arXiv. https://doi.org/10.48550/arXiv.2312.06824 589 Global Antimicrobial Resistance and Use Surveillance System (GLASS) Report 2022 (1st ed). (2022). World 590 Health Organization. 591 Greenbury, S. F ., Barahona, M., & Johnston, I. G. (2020). HyperTraPS: Inferring probabilistic patterns of trait 592 acquisition in evolutionary and disease progression pathways. Cell Systems, 10(1), 39–51. 593 Gupta, S. K., Padmanabhan, B. R., Diene, S. M., Lopez-Rojas, R., Kempf, M., Landraud, L., & Rolain, J.-M. 594 (2014). ARG-ANNOT, a new bioinformatic tool to discover antibiotic resistance genes in bacterial 595 genomes. Antimicrobial Agents and Chemotherapy, 58(1), 212–220. 596 https://doi.org/10.1128/AAC.01310-13 597 Harris, C. R., Millman, K. J., van der Walt, S. J., Gommers, R., Virtanen, P ., Cournapeau, D., Wieser, E., Taylor, J., 598 Berg, S., Smith, N. J., Kern, R., Picus, M., Hoyer, S., van Kerkwijk, M. H., Brett, M., Haldane, A., del Río, 599 J. F., Wiebe, M., Peterson, P ., … Oliphant, T. E. (2020). Array programming with NumPy. Nature, 600 585(7825), 357–362. https://doi.org/10.1038/s41586-020-2649-2 601 Hennart, M., Guglielmini, J., Bridel, S., Maiden, M. C. J., Jolley, K. A., Criscuolo, A., & Brisse, S. (2022). A Dual 602 Barcoding Approach to Bacterial Strain Nomenclature: Genomic Taxonomy of Klebsiella pneumoniae 603 Strains. Molecular Biology and Evolution, 39(7), msac135. https://doi.org/10.1093/molbev/msac135 604 Herencias, C., Álvaro-Llorente, L., Ramiro-Martínez, P ., Fernández-Calvet, A., Muñoz-Cazalla, A., DelaFuente, 605 J., Graf, F. E., Jaraba-Soto, L., Castillo-Polo, J. A., Cantón, R., San Millán, Á., & Rodríguez-Beltrán, J. 606 (2024). β-lactamase expression induces collateral sensitivity in Escherichia coli. Nature 607 Communications, 15(1), 4731. https://doi.org/10.1038/s41467-024-49122-2 608 Hetland, M. A. K., Winkler, M. A., Kaspersen, H. P ., Håkonsholm, F ., Bakksjø, R.-J., Bernho\, E., Delgado-Blas, 609 J. F., Brisse, S., Correia, A., Fostervold, A., Lam, M. M. C., Lunestad, B.-T., Marathe, N. P., 610 Ra\elsberger, N., Samuelsen, Ø., Sunde, M., Sundsfjord, A., Urdahl, A. M., Wick, R. R., … Holt, K. E. 611 (2025). A genome-wide One Health study of Klebsiella pneumoniae in Norway reveals overlapping 612 populations but few recent transmission events across reservoirs. Genome Medicine, 17(1), 42. 613 https://doi.org/10.1186/s13073-025-01466-0 614 Holt, K. E., Wertheim, H., Zadoks, R. N., Baker, S., Whitehouse, C. A., Dance, D., Jenney, A., Connor, T. R., Hsu, 615 L. Y ., Severin, J., Brisse, S., Cao, H., Wilksch, J., Gorrie, C., Schultz, M. B., Edwards, D. J., Nguyen, K. V ., 616 Nguyen, T. V ., Dao, T. T., … Thomson, N. R. (2015). Genomic analysis of diversity, population structure, 617 virulence, and antimicrobial resistance in Klebsiella pneumoniae, an urgent threat to public health. 618 .CC-BY 4.0 International licenseperpetuity. It is made available under a preprint (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in The copyright holder for thisthis version posted September 21, 2025. ; https://doi.org/10.1101/2025.09.20.677523doi: bioRxiv preprint Proceedings of the National Academy of Sciences, 112(27), E3574–E3581. 619 https://doi.org/10.1073/pnas.1501049112 620 Johnston, I. G., & Diaz-Uriarte, R. (2024). A hypercubic Mk model framework for capturing reversibility in 621 disease, cancer, and evolutionary accumulation modelling (p. 2024.06.27.600959). bioRxiv. 622 https://doi.org/10.1101/2024.06.27.600959 623 Johnston, I. G., & Williams, B. P . (2016). Evolutionary Inference across Eukaryotes Identifies Specific Pressures 624 Favoring Mitochondrial Gene Retention. Cell Systems, 2(2), 101–111. 625 https://doi.org/10.1016/j.cels.2016.01.013 626 Kassambara, A. (2020). Ggpubr:“ggplot2” based publication ready plots. R Package Version 0.4. 0, 438. 627 Kurtz, S., Phillippy, A., Delcher, A. L., Smoot, M., Shumway, M., Antonescu, C., & Salzberg, S. L. (2004). 628 Versatile and open software for comparing large genomes. Genome Biology, 5(2), R12. 629 https://doi.org/10.1186/gb-2004-5-2-r12 630 Lam, M. M. C., Wick, R. R., Watts, S. C., Cerdeira, L. T., Wyres, K. L., & Holt, K. E. (2021). A genomic 631 surveillance framework and genotyping tool for Klebsiella pneumoniae and its related species 632 complex. Nature Communications, 12(1), 4188. https://doi.org/10.1038/s41467-021-24448-3 633 Liu, P ., Li, P ., Jiang, X., Bi, D., Xie, Y ., Tai, C., Deng, Z., Rajakumar, K., & Ou, H.-Y. ( 2 0 1 2 ) . C o m p l e t e g e n o m e 634 sequence of Klebsiella pneumoniae subsp. Pneumoniae HS11286, a multidrug-resistant strain 635 isolated from human sputum. Journal of Bacteriology, 194(7), 1841–1842. 636 https://doi.org/10.1128/JB.00043-12 637 Luo, X. G., Kuipers, J., & Beerenwinkel, N. (2023). Joint inference of exclusivity patterns and recurrent 638 trajectories from tumor mutation trees. Nature Communications, 14(1), Article 1. 639 https://doi.org/10.1038/s41467-023-39400-w 640 Maddison, W. P ., & FitzJohn, R. G. (2015). The Unsolved Challenge to Phylogenetic Correlation Tests for 641 Categorical Characters. Systematic Biology, 64(1), 127–136. https://doi.org/10.1093/sysbio/syu070 642 Mapunjo, S., Mbwasi, R., Nkiligi, E. A., Wilbroad, A., Francis, E. N., Msovela, K., Yahya, T., Mpembeni, R., 643 Masunga, E., Nkungu, K., Saitoti, S., Lusaya, E., & Konduri, N. (2025). National consumption of 644 antimicrobials in Tanzania: 2020–2022. JAC-Antimicrobial Resistance, 7(2), dlaf026. 645 https://doi.org/10.1093/jacamr/dlaf026 646 Moen, M. T., & Johnston, I. G. (2023). HyperHMM: E\icient inference of evolutionary and progressive dynamics 647 on hypercubic transition graphs. Bioinformatics, 39(1), btac803. 648 Moyo, S. J., Manyahi, J., Blomberg, B., Tellevik, M. G., Masoud, N. S., Aboud, S., Manji, K., Roberts, A. P ., 649 Hanevik, K., Mørch, K., & Langeland, N. (2020). Bacteraemia, Malaria, and Case Fatality Among 650 Children Hospitalized With Fever in Dar es Salaam, Tanzania. Frontiers in Microbiology, 11, 2118. 651 https://doi.org/10.3389/fmicb.2020.02118 652 Munk, P ., Brinch, C., Møller, F . D., Petersen, T. N., Hendriksen, R. S., Seyfarth, A. M., Kjeldgaard, J. S., 653 Svendsen, C. A., van Bunnik, B., Berglund, F ., Larsson, D. G. J., Koopmans, M., Woolhouse, M., & 654 Aarestrup, F . M. (2022). Genomic analysis of sewage from 101 countries reveals global landscape of 655 antimicrobial resistance. Nature Communications, 13(1), 7251. https://doi.org/10.1038/s41467-022-656 34312-7 657 Murray, C. J. L., Ikuta, K. S., Sharara, F ., Swetschinski, L., Aguilar, G. R., Gray, A., Han, C., Bisignano, C., Rao, P ., 658 Wool, E., Johnson, S. C., Browne, A. J., Chipeta, M. G., Fell, F., Hackett, S., Haines-Woodhouse, G., 659 Hamadani, B. H. K., Kumaran, E. A. P ., McManigal, B., … Naghavi, M. (2022). Global burden of 660 bacterial antimicrobial resistance in 2019: A systematic analysis. The Lancet, 399(10325), 629–655. 661 https://doi.org/10.1016/S0140-6736(21)02724-0 662 Naghavi, M., Vollset, S. E., Ikuta, K. S., Swetschinski, L. R., Gray, A. P ., Wool, E. E., Aguilar, G. R., Mestrovic, T., 663 Smith, G., Han, C., Hsu, R. L., Chalek, J., Araki, D. T., Chung, E., Raggi, C., Hayoon, A. G., Weaver, N. 664 D., Lindstedt, P . A., Smith, A. E., … Murray, C. J. L. (2024). Global burden of bacterial antimicrobial 665 resistance 1990–2021: A systematic analysis with forecasts to 2050. The Lancet, 404(10459), 1199–666 1226. https://doi.org/10.1016/S0140-6736(24)01867-1 667 Onken, A., Moyo, S., Miraji, M. K., Bohlin, J., Marijani, M., Manyahi, J., Kibwana, K. O., Müller, F ., Jenum, P . A., 668 Abeid, K. A., Reimers, M., Langeland, N., Mørch, K., & Blomberg, B. (2024). Predominance of 669 multidrug-resistant Salmonella Typhi genotype 4.3.1 with low-level ciprofloxacin resistance in 670 Zanzibar. PLOS Neglected Tropical Diseases, 18(4), e0012132. 671 https://doi.org/10.1371/journal.pntd.0012132 672 Onken, A., Said, A. K., Jørstad, M., Jenum, P . A., & Blomberg, B. (2015). Prevalence and Antimicrobial 673 Resistance of Microbes Causing Bloodstream Infections in Unguja, Zanzibar. PLoS ONE, 10(12), 674 e0145632. https://doi.org/10.1371/journal.pone.0145632 675 Pál, C., Papp, B., & Lázár, V. (2015). Collateral sensitivity of antibiotic-resistant microbes. Trends in 676 Microbiology, 23(7), 401–407. https://doi.org/10.1016/j.tim.2015.02.009 677 .CC-BY 4.0 International licenseperpetuity. It is made available under a preprint (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in The copyright holder for thisthis version posted September 21, 2025. ; https://doi.org/10.1101/2025.09.20.677523doi: bioRxiv preprint R Core Team. (2022). R: A language and environment for statistical computing. R Foundation for Statistical 678 Computing, Vienna, Austria. 2012. 679 Renz, J., Dauda, K. A., Aga, O. N. L., Diaz-Uriarte, R., Löhr, I. H., Blomberg, B., & Johnston, I. G. (2024). 680 Evolutionary accumulation modelling in AMR: Machine learning to infer and predict evolutionary 681 dynamics of multi-drug resistance (arXiv:2411.00219). arXiv. 682 https://doi.org/10.48550/arXiv.2411.00219 683 Renz, J., Dauda, K. A., Aga, O. N. L., Diaz-Uriarte, R., Löhr, I. H., Blomberg, B., & Johnston, I. G. (2025). 684 Evolutionary accumulation modeling in AMR: Machine learning to infer and predict evolutionary 685 dynamics of multi-drug resistance. mBio, 16(6), e00488-25. https://doi.org/10.1128/mbio.00488-25 686 Revell, L. J. (2012). phytools: An R package for phylogenetic comparative biology (and other things). Methods 687 in Ecology and Evolution, 2, 217–223. 688 Rudd, K. E., Johnson, S. C., Agesa, K. M., Shackelford, K. A., Tsoi, D., Kievlan, D. R., Colombara, D. V ., Ikuta, K. 689 S., Kissoon, N., Finfer, S., Fleischmann-Struzek, C., Machado, F . R., Reinhart, K. K., Rowan, K., 690 Seymour, C. W., Watson, R. S., West, T. E., Marinho, F ., Hay, S. I., … Naghavi, M. (2020). Global, 691 regional, and national sepsis incidence and mortality, 1990–2017: Analysis for the Global Burden of 692 Disease Study. The Lancet, 395(10219), 200–211. https://doi.org/10.1016/S0140-6736(19)32989-7 693 Saitou, N., & Nei, M. (1987). The neighbor-joining method: A new method for reconstructing phylogenetic 694 trees. Molecular Biology and Evolution, 4(4), 406–425. 695 Sangeda, R. Z., Saburi, H. A., Masatu, F . C., Aiko, B. G., Mboya, E. A., Mkumbwa, S., Bitegeko, A., Mwalwisi, Y . 696 H., Nkiligi, E. A., Chambuso, M., Sillo, H. B., Fimbo, A. M., & Horumpende, P . G. (2021). National 697 Antibiotics Utilization Trends for Human Use in Tanzania from 2010 to 2016 Inferred from Tanzania 698 Medicines and Medical Devices Authority Importation Data. Antibiotics, 10(10), 1249. 699 https://doi.org/10.3390/antibiotics10101249 700 Schill, R., Klever, M., Rupp, K., Hu, Y . L., Lösch, A., Georg, P ., Pfahler, S., Vocht, S., Hansch, S., Wettig, T., 701 Grasedyck, L., & Spang, R. (2024). Reconstructing Disease Histories in Huge Discrete State Spaces. KI 702 - Künstliche Intelligenz. https://doi.org/10.1007/s13218-023-00822-9 703 Schill, R., Solbrig, S., Wettig, T., & Spang, R. (2020). Modelling cancer progression using mutual hazard 704 networks. Bioinformatics, 36(1), 241–249. 705 Schliep, K. P . (2011). phangorn: Phylogenetic analysis in R. Bioinformatics, 27(4), 592–593. 706 https://doi.org/10.1093/bioinformatics/btq706 707 Slowikowski, K. (2021). ggrepel: Automatically Position Non-Overlapping Text Labels with’ggplot2’ . R package 708 version 0.9. 1, 2021. 709 Sokal, R. R., & Michener, C. D. (1958). A statistical method for evaluating systematic relationships. University 710 of Kansas Science Bulletin, 38, 1409. 711 Szybalski, W., & Bryson, V . (1952). Genetic studies on microbial cross resistance to toxic agents i. Journal of 712 Bacteriology, 64(4), 489–499. https://doi.org/10.1128/jb.64.4.489-499.1952 713 Tsang, K. K., Lam, M. M. C., Wick, R. R., Wyres, K. L., Bachman, M., Baker, S., Barry, K., Brisse, S., Campino, S., 714 Chiaverini, A., Cirillo, D. M., Clark, T., Corander, J., Corbella, M., Cornacchia, A., Cuénod, A., D’Alterio, 715 N., Di Marco, F ., Donado-Godoy, P ., … KlebNET-GSP AMR Genotype-Phenotype Group. (2024). 716 Diversity, functional classification and genotyping of SHV β-lactamases in Klebsiella pneumoniae. 717 Microbial Genomics, 10(10), 001294. https://doi.org/10.1099/mgen.0.001294 718 Wickham, H. (2016). ggplot2: Elegant Graphics for Data Analysis. Springer-Verlag. 719 https://link.springer.com/book/10.1007/978-0-387-98141-3 720 Wickham, H., Averick, M., Bryan, J., Chang, W., McGowan, L., François, R., Grolemund, G., Hayes, A., Henry, 721 L., Hester, J., Kuhn, M., Pedersen, T., Miller, E., Bache, S., Müller, K., Ooms, J., Robinson, D., Seidel, D., 722 Spinu, V ., … Yutani, H. (2019). Welcome to the Tidyverse. Journal of Open Source Software, 4(43), 723 1686. https://doi.org/10.21105/joss.01686 724 Wickham, H., François, R., Henry, L., Müller, K., & Vaughan, D. (2023). Dplyr: A grammar of data manipulation. 725 R package version 1.1. 2. Computer Software. 726 Williams, B. P ., Johnston, I. G., Covsho\, S., & Hibberd, J. M. (2013). Phenotypic landscape inference reveals 727 multiple evolutionary paths to C4 photosynthesis. eLife, 2, e00961. 728 https://doi.org/10.7554/eLife.00961 729 Wyres, K. L., & Holt, K. E. (2018). Klebsiella pneumoniae as a key tra\icker of drug resistance genes from 730 environmental to clinically important bacteria. Current Opinion in Microbiology, 45, 131–139. 731 https://doi.org/10.1016/j.mib.2018.04.004 732 Wyres, K. L., Lam, M. M. C., & Holt, K. E. (2020). Population genomics of Klebsiella pneumoniae. Nature 733 Reviews. Microbiology, 18(6), 344–359. https://doi.org/10.1038/s41579-019-0315-1 734 735 .CC-BY 4.0 International licenseperpetuity. It is made available under a preprint (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in The copyright holder for thisthis version posted September 21, 2025. ; https://doi.org/10.1101/2025.09.20.677523doi: bioRxiv preprint 736 .CC-BY 4.0 International licenseperpetuity. It is made available under a preprint (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in The copyright holder for thisthis version posted September 21, 2025. ; https://doi.org/10.1101/2025.09.20.677523doi: bioRxiv preprint Supplementary Information 737 738 739 Supplementary Figure 1. Synthetic control study on reversibility. (top) A synthetic dataset involving the 740 progressive accumulation of characters X1 to X10. The output inferred from HyperTraPS accurately reports the 741 accumulation ordering of characters, with ambiguity about the first three characters which are observed 742 across all samples. (bottom) A synthetic dataset involving reversibility. The same accumulation process is 743 simulated, but with the potential for acquired characters to be lost reversibly. The output from HyperTraPS 744 (which assumes irreversibility) is not dramatically di\erent, except for a moderate amplification of the 745 uncertainty in the early characters. The e\ect of reversibility here is to increase uncertainty rather than to 746 introduce any systematic bias in the outputs. 747 748 .CC-BY 4.0 International licenseperpetuity. It is made available under a preprint (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in The copyright holder for thisthis version posted September 21, 2025. ; https://doi.org/10.1101/2025.09.20.677523doi: bioRxiv preprint 749 Supplementary Figure 2. All posterior orderings of AMR character acquisitions. This “global roadmap” 750 summarises “bubble plots” describing inferred evolutionary dynamics of KpAMR characters in di\erent 751 countries. Each row is a country, ordered vertically by descending number of Kp samples. Each column is a 752 KpAMR character; the horizontal axis within columns gives evolutionary orderings (from early to late). The size 753 of a point gives the probability that that character is acquired at that ordering in that country. For example, 754 across the vast majority of datasets (including the top row, USA, for example), Bla_chr has a high probability of 755 early acquisition and Gly_acquired has a high probability of late acquisition. In Cameroon, all characters have 756 similar inferred evolutionary orderings, reflecting the limited data available to make more precise statements. 757 Di\erent colour points give the output of runs with di\erent random seeds, to demonstrate consistency of the 758 estimates. Character names given in Fig. 1 caption. 759 760 .CC-BY 4.0 International licenseperpetuity. It is made available under a preprint (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in The copyright holder for thisthis version posted September 21, 2025. ; https://doi.org/10.1101/2025.09.20.677523doi: bioRxiv preprint 761 Supplementary Figure 3. Inferred interactions between KpAMR characters. An arrow from character X to 762 character Y reports an inferred interaction between the evolutionary acquisition of those characters. Blue 763 arrow: acquiring X makes acquiring Y more likely (“promoting”). Red arrow: acquiring X makes acquiring Y less 764 likely (“repressing”). The number on each arrow gives the proportion of countries in which that interaction was 765 inferred. Character names given in Fig. 1 caption. 766 .CC-BY 4.0 International licenseperpetuity. It is made available under a preprint (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in The copyright holder for thisthis version posted September 21, 2025. ; https://doi.org/10.1101/2025.09.20.677523doi: bioRxiv preprint 767 Supplementary Figure 4. KpAMR characters determining aspects of global variability in evolutionary 768 dynamics (alternative visualisation). Connected to Fig. 4. KpAMR characters (horizontal axis) are plotted by 769 their expected acquisition ordering (vertical axis) for each country (points). The country’s projection on PCA1, 770 PCA2, and PCA3, for those characters that most strongly covary with each PCA axis, are given by the colour of 771 a point. As in Fig. 4, those characters that covary with PCA1 form a consistent spectrum linked to precision: 772 low values of PCA1 correspond to less precise, more uniform timing, while higher values lead to more precise 773 earlier or later timings. Characters covarying with the other PCAs display wider ranges connected to PCA1-774 independent variability. Character names given in Fig. 1 caption. 775 776 .CC-BY 4.0 International licenseperpetuity. It is made available under a preprint (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in The copyright holder for thisthis version posted September 21, 2025. ; https://doi.org/10.1101/2025.09.20.677523doi: bioRxiv preprint 777 Supplementary Figure 5. New KpAMR data from Zanzibar and Dar es Salaam, Tanzania. Phylogenies and 778 KpAMR profiles from new datasets: (A) Zanzibar 2015-6. (B) Tanzania 2001-2. (C) Tanzania 2017-8. (D) 779 Pathogenwatch Tanzania. (E) Circular projection of the phylogeny of Tanzanian isolates in Fig. 7A. Colours: 780 green, 2001-2; blue, 2017-8; red, existing Pathogenwatch data. 781 .CC-BY 4.0 International licenseperpetuity. It is made available under a preprint (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in The copyright holder for thisthis version posted September 21, 2025. ; https://doi.org/10.1101/2025.09.20.677523doi: bioRxiv preprint