Integrating whole-genome sequencing within the National Antimicrobial Resistance Surveillance Program in the Philippines

National networks of laboratory-based surveillance of antimicrobial resistance (AMR) monitor resistance trends and disseminate these data to AMR stakeholders. Whole-genome sequencing (WGS) can support surveillance by pinpointing resistance mechanisms and uncovering transmission patterns. However, genomic surveillance is rare in low- and middle-income countries. Here, we implement WGS within the established Antimicrobial Resistance Surveillance Program of the Philippines via a binational collaboration. In parallel, we characterize bacterial populations of key bug-drug combinations via a retrospective sequencing survey. By linking the resistance phenotypes to genomic data, we reveal the interplay of genetic lineages (strains), AMR mechanisms, and AMR vehicles underlying the expansion of specific resistance phenotypes that coincide with the growing carbapenem resistance rates observed since 2010. Our results enhance our understanding of the drivers of carbapenem resistance in the Philippines, while also serving as the genetic background to contextualize ongoing local prospective surveillance.

A ntimicrobial resistance (AMR) is an increasingly serious threat to global public health and the economy that requires concerted action across countries, government sectors and non-government organizations 1 . Without AMR containment, an adverse impact on medical costs, global gross domestic product, livestock production and international trade is expected by 2050, and the sustainable development goals for 2030 are less likely to be attained 2 . The Global Action Plan developed by the World Health Organization (WHO) to tackle AMR highlights the need to strengthen our understanding of how resistance develops and spreads, and the underlying resistance mechanisms 3 . One of the pillars of this objective is the national surveillance system for AMR.
An advanced example of a national surveillance system is the Antimicrobial Resistance Surveillance Program (ARSP), the national laboratory-based surveillance of the Philippine Department of Health (DOH). The ARSP began in 1988 with eight hospital-based laboratories (sentinel sites) and the Antimicrobial Resistance Surveillance Reference Laboratory (ARSRL) in Metro Manila and has since expanded to 24 sentinel sites and two gonorrhoea surveillance sites in 16 regions of the country (Supplementary Fig. 1). Surveillance encompasses common pathogens of public health importance in the Philippines, where infectious diseases represented nine out of the ten leading causes of morbidity 4 , and four out of the ten leading causes of infant mortality 5 in 2013-2014. Results of culture-based identification and antimicrobial susceptibility are stored centrally at the ARSRL within the clinical microbiology management software WHONET (Fig. 1a 6 ). A report summarizing the resistance trends is published annually and disseminated to local, national and international surveillance stakeholders 7 . For policy change, surveillance data generated by ARSP are used by the DOH to develop clinical practice guidelines and determine the panel of antibiotics to be included in the national formulary. Furthermore, the ARSP has been contributing data to international AMR surveillance programmes since the 1980s, including the Global Antimicrobial Surveillance System (GLASS) 8 since 2015. Importantly, the ARSP data has informed the Philippine National Action Plan to combat AMR, one of the key requirements for alignment with the WHO Global Action Plan.
The Philippines has seen a steady increase in resistance rates for several key pathogen-antibiotic combinations in the last 10 years, including carbapenem-resistant organisms (Fig. 1b). The genetic mechanisms underlying carbapenem resistance, one of the biggest therapeutic challenges in the treatment of antimicrobialresistant infections, include increased upregulation of efflux pumps, decreased uptake by altered expression/loss of porin function, and acquisition of hydrolytic enzymes-carbapenemases 9 . Whole-genome sequencing (WGS) of bacterial pathogens can identify distinct clonal lineages (strains) on phylogenetic trees, known AMR mechanisms (genes or mutations) and the vehicles for the acquired AMR genes (mobile genetic elements, MGEs) which, in turn, enables enhanced detection and characterization of high-risk clones (Fig. 1a). WGS is routinely used for infectious disease epidemiology in several high-income countries around the world, where it has improved outbreak investigations and epidemiological surveillance 10,11 and enhanced our knowledge of the spread of antimicrobial-resistant strains and their resistance mechanisms 12,13 . Elucidation of resistance mechanisms and their context can be critical for effective infection control. For example, upregulation of efflux pumps and loss of porin function by mutation are vertically transmitted, while acquired carbapenemases carried in transmissible plasmids or integrative conjugative elements have the potential for horizontal dissemination between strains and species, thus necessitating enhanced infection control measures 14 . Moreover, diverse carbapenemase classes, variants, flanking MGEs and associated plasmid incompatibility groups have been described, requiring multiple biochemical and molecular tests for identification. For example, the class B New Delhi metallo-beta-lactamase (NDM) is found worldwide, represented by over ten different variants that are usually upstream of intact or truncated ISAba125 elements, within plasmid backbones of different incompatibility groups, across multiple lineages of E. coli-as well as of other carbapenem-resistant organisms 15,16 . WGS can determine the interplay between these different components, i.e. the gene, the vehicle and the strain, thus maximizing the epidemiological benefit derived from its cost.
Implementation of WGS within existing national surveillance systems in low-and middle-income countries (LMICs) has the potential to enhance local prevention and control of resistant infections in a sustainable and equitable manner. Integration of WGS into routine surveillance of AMR can be facilitated by international collaboration focused on transfer of expertize and ownership. Our collaboration (See and Sequence, Fig. 1a) aimed to implement WGS within the ARSP via a multi-faceted approach that included a large initial retrospective sequencing survey, technology transfer, utilization of user-friendly web applications and capacity building in laboratory and bioinformatic procedures for local prospective sequencing. We linked what has traditionally been used as the operational unit of laboratory surveillance, the resistance profile (RP), with the operational unit of genomic epidemiology now provided by WGS, genetic relatedness. Here we provide exemplars from the retrospective sequencing survey that highlight how the granular view of strain-gene-vehicle in carbapenem-resistant populations at the local, national and global operational scales can be leveraged for surveillance of AMR and public health.

Results
Integration of laboratory and genomic surveillance. The Philippine ARSP collects bacterial isolates and stores the associated clinical and epidemiological data in the WHONET software 6 . Antibiograms are stored as RPs, which represent the diversity of AMR phenotypes in the country. We implemented WGS within the ARSP via ongoing technology transfer, capacity building in and binational collaboration. In parallel, we conducted a large retrospective sequencing survey to characterize bacterial populations, dissect resistance phenotypes of key bug-drug combinations and provide genomic context for local prospective surveillance. Here, we focus on the carbapenem-resistant organisms, which have been classified as of critical priority for the development of new antibiotics by the WHO.
Operational unit of laboratory surveillance-increasing the burden of specific resistance phenotypes. Carbapenem-resistant P. aeruginosa, A. baumannii, K. pneumoniae and E. coli were all first isolated by the ARSP between 1992 and 1994. Since the year 2000, resistance rates to imipenem and meropenem have remained below 30% for all organisms except for A. baumannii, which has seen a steady rise in resistance rates between 2009 and 2017, reaching 56% (Fig. 2a). This coincides with the expansion of two RPs, first a striking expansion of true extremely drugresistant (XDR) and possible pan drug-resistant (PDR) RP-1, and later of true-XDR RP-2 (full profiles can be found in Fig. 2b). Yearly resistance rates for P. aeruginosa have oscillated between 10 and 25% since the year 2000, and have doubled between 2010 and 2013 (Fig. 2a). Two RPs have seen the largest combined expansion between 2011 and 2013, possible-PDR RP-3 and true-XDR RP-4 (Fig. 2b). The yearly resistance rates to imipenem for Enterobacteriaceae were~3% in 2011 but have since steadily increased to 5% for E. coli and 11% for K. pneumoniae, respectively (Fig. 2a). This coincides with the expansion of three possible-XDR RPs, RP-5, RP-6 and RP-7 (Fig. 2b). The observed expansion of specific RPs could be driven by the rapid dissemination of a discrete genetic lineage (strain) carrying AMR determinants or, alternatively, of transmissible MGEs, which shuttle resistance determinants across different genetic lineages.
Operational unit of genomic surveillance-genetic diversity and mechanisms underpinning carbapenem resistance phenotypes. The sequences of 805 genomes linked to WHONET data were obtained from isolates belonging to the four bacterial pathogens in the WHO critical list (Table 1), which were collected in 2013 and 2014 by between 13 and 18 ARSP sites. The isolates sequenced were non-susceptible to carbapenems and/or to extended-spectrum cephalosporins (Table 1), with the exception of P. aeruginosa, for which susceptible isolates were available and also sequenced.
The distribution of the carbapenem RPs highlighted in Fig. 2 was not concordant with the major clades observed on the phylogenetic trees (Fig. 3). Instead, the same RP was usually observed in multiple genetic lineages across the tree. Yet, the distribution of the multi-locus sequence types (STs), a molecular genotyping method inferred from the sequences of seven housekeeping genes, was largely concordant with the major clades observed on the trees inferred from whole-genome singlenucleotide polymorphisms (SNPs) (Fig. 3), confirming that the homoplasic distribution of RPs was not an artefact of the tree topology. The distribution of the same carbapenem RPs across multiple genetic lineages was observed both at the national and local levels ( Supplementary Fig. 2 20 10 0 1 9 8 8 1 9 8 9 1 9 9 0 1 9 9 1 1 9 9 2 1 9 9 3 1 9 9 4 1 9 9 5 1 9 9 6 1 9 9 7 1 9 9 8  Fig. 1 Implementing whole-genome sequencing (WGS) for AMR surveillance in the Philippines. a ARSP workflow and enhanced detection of high-risk clones by WGS. Isolates collected by sentinel sites are tested for susceptibility to antibiotics (open squares: susceptible, solid squares: resistant). The data are stored as resistance profiles in WHONET and summaries of resistance trends are shared yearly with surveillance stakeholders. WGS of bacterial isolates and interpretation with web applications like Microreact and Pathogenwatch provide information on genetic relatedness (strains), known AMR determinants (mechanisms), and the mobile genetic elements (MGE, vehicles) for their dissemination, thus allowing us to detect high-risk clones. b Detail of trends in antimicrobial resistance in the Philippines. Yearly resistance rates for key bug-drug combinations based on phenotypic data collected by sentinel sites. EBC Enterobacteriaceae (K. pneumoniae, E. coli, Salmonella enterica), PAE P. aeruginosa, ABA A. baumannii (Acinetobacter spp. before the year 2000), SAL Salmonella enterica, SAU S. aureus, ESBL extended-spectrum beta-lactamase production suspected, or non-susceptible to the following antibiotics IPM imipenem, CIP ciprofloxacin, OXA oxacillin, WGS box: period covered by the retrospective sequencing survey. NATURE COMMUNICATIONS | https://doi.org/10.1038/s41467-020-16322-5 ARTICLE NATURE COMMUNICATIONS | (2020) 11:2719 | https://doi.org/10.1038/s41467-020-16322-5 | www.nature.com/naturecommunications assignments increased with the number of isolates representing a RP. Thus, the genomic data suggest that the expansion of specific carbapenem RPs in the Philippines is driven by the horizontal dissemination of resistance determinants across diverse genetic backgrounds, rather than by the sole expansion of single resistant genetic clones, and that resistance to multiple antibiotics is acquired en bloc via MGEs, as previously observed for carbapenem-resistant organisms 16 .
Carbapenemases are a diverse collection of hydrolytic enzymes 16 , and we identified representatives of classes A, B and C carbapenemases in the Philippines (Supplementary Fig. 3). Class B NDM-1 and Verona integron-borne metallo-beta-    lactamase VIM-2 were the most prevalent and geographically disseminated in the Enterobacteriaceae and P. aeruginosa, respectively, while the same was true for class D OXA-23 in A. baumannii. Importantly, different carbapenemase genes/variants (or combinations thereof) were found underlying the same RP in all four organisms (Fig. 3).
Local scale of surveillance-plasmid-driven hospital outbreak of carbapenem-resistant K. pneumoniae. Fifty-seven percent (N = 33) of the K. pneumoniae genomes with the possible-XDR resistance phenotype AMP FOX CAZ CRO FEP IPM AMC TZP GEN AMK CIP SXT (RP-6) were isolated from a single hospital (MMH) but belonged to 12 different STs, with almost half of them (N = 15) placed within ST340 ( Supplementary Fig. 4a). The phylogenetic analysis of these 15 genomes in the wider context of ST340 showed three main lineages, with the 15 possible-XDR isolates from MMH forming a tight cluster (clade III, Fig. 4a). The shorter branch lengths in the imipenem-resistant clade III (6.5 ± 3.7 pairwise SNP differences) compared with the imipenemsusceptible clade II (26.8 ± 6.9 pairwise SNP differences, Mann-Whitney U test z-score −6.265, p value = 3.71 × 10 −10 , Supplementary Fig. 4b), and the isolation dates spanning 12 months, suggest a rapid expansion of this clone, which coincides with the acquisition of the bla NDM-1 gene (Fig. 4a). Epidemiological data showed that the samples were mostly from hospital-acquired infections (HAI, N = 14) and from neonates (N = 12) (Fig. 4a). These observations triggered a retrospective epidemiological investigation that revealed that ten of the isolates originated from patients of the neonatal intensive care unit (NICU), all of which had umbilical catheters and were on mechanical ventilators. The remaining cases were either from paediatric wards (N = 4) or the male ward (N = 1), suggesting wider dissemination of this high-risk clone within the hospital environment.
The genetic diversity underlying the 33 possible-XDR isolates from MMH (12 STs) prompted us to investigate the hypothesis of dissemination of carbapenem resistance in the hospital by a plasmid carrying bla NDM-1 . We identified a 101,540 bp IncFII plasmid with bla NDM-1 , rmtC and sul1 (p13ARS_MMH0112-3, Table 2, Supplementary Note 1 and Supplementary Fig. 5a, b) in a representative isolate from ST340 clade III. We found the entire p13ARS_MMH0112-3 plasmid sequence (i.e., ≥95% of the length) represented in 27 genomes from 9 different STs, including 20 samples from 6 different STs isolated from patients under 1 year old (Fig. 4b). A second plasmid identified in the representative isolate carried 13 AMR genes, including the extended-spectrum beta-lactamase gene bla CTX-M-15 (p13ARS_MMH0112-2, Table 2, Supplementary Note 2 and Supplementary Fig. 5a, c).
Altogether, our findings suggest that the burden of carbapenemresistant K. pneumoniae infections in hospital MMH was largely linked to plasmid p13ARS_MMH0112-3 with bla NDM-1 , circulating within diverse genetic lineages, and leading to outbreaks in high-risk patient populations. Hospital authorities were informed of these findings, and measures for infection control were implemented, including designating a separate multi-drug resistance organism room for cohorting, active surveillance upon identification of any new carbapenem-resistant K. pneumoniae from the NICU, and referral of any new carbapenem-resistant K. pneumoniae from the NICU to ARSRL for sequencing.
National scale of surveillance-regional circulation of a successful K. pneumoniae lineage. The possible-XDR resistance phenotype AMP FOX CAZ CRO FEP IPM AMC TZP CIP SXT (RP-5) was represented by 76 K. pneumoniae isolates from 14 different STs. Seventy-one percent (N = 54) belonged to ST147, an international epidemic clone that was found in 11 sentinel sites in all three island groups in the Philippines. The phylogeny of the ST147 genomes showed that carbapenem non-susceptible isolates (78.8%, N = 63) were found in clades arising from three out of four deep branches of the tree ( Five plasmid sequences obtained from five different isolates representing carbapenem-resistant clusters in clades III and IV showed that the NDM genes were carried within different variants of the insertion sequence ISAba125 ( Supplementary  Fig. 8a) on different plasmid backbones, two of which showed high sequence similarity to international plasmids (Table 3 and Supplementary Figs. 5b and 8b). The distribution of plasmids harbouring bla NDM-1 and bla NDM-7 matched the strong phylogeographic signal in the terminal branches of the tree (Fig. 5). Within clade III-B, a cluster of genomes from sentinel sites MMH and STU was characterized by plasmid p14ARS_MMH0055-5 carrying bla NDM-7 , while another cluster of genomes from VSM and NMC was distinguished by plasmid p13ARS_VSM0593-1 with bla NDM-1 . Plasmids p13ARS_GMH0099 and p14ARS_VSM0843-1, both carrying bla NDM-1 , were found in different subclades of clade IV, one of local circulation in hospital VSM, and another one showing regional expansion across five sentinel sites (Fig. 5). Taken together, the phylogeographic signal (strains) in combination with the distribution of bla NDM variants (mechanism) and plasmids (vehicle) revealed local and regional patterns of circulation of carbapenem non-susceptible ST147 within the Philippines, which could not be identified based solely on the RPs and MLST information. Furthermore, comparison with global ST147 genomes indicated that clade IV may be unique to the Philippines ( Supplementary Fig. 6).
International scale of surveillance-high-risk clone of E. coli ST410 with NDM-1 and OXA-181. Recent reports of E. coli    distribution of carbapenemase genes and variants into three different clones carrying bla NDM-4 (N = 1), bla NDM-1 (N = 8) or both bla NDM-1 and bla OXA-181 (N = 5). The phylogenetic analysis with global E. coli ST410 genomes confirmed that these are independent clones, as they were found interspersed with international genomes within a major high-risk lineage ( Fig. 6b and Supplementary Note 5). K. pneumoniae carrying bla NDM-1 and bla OXA-181 have been previously reported 21,22 , as well as E. coli ST410 harbouring bla NDM-5 and bla OXA-181 from Egypt 23 , Denmark and the UK 18 , but, to our knowledge, this is the first report of E. coli carrying both bla NDM-1 and bla OXA-181 , which is likely to have disseminated between two sentinel sites (NMC and VSM). We identified five plasmids in the representative strain 14ARS_NMC0074 (Table 4 and Supplementary Fig. 10). The IncX3 plasmid carrying bla OXA-181 (p14ARS_NMC0074-5) was identical to plasmid pAMA1167-OXA-181 isolated previously from an E. coli ST410 strain 24 (Supplementary Fig. 10). Mapping short reads to p14ARS_NMC0074-5 showed that this plasmid is the main vehicle of bla OXA-181 in the Philippines, as in the international genomes (Fig. 6b 18 ). However, we did not identify any known plasmids similar to the IncA/C2 plasmid with bla NDM-1 (p14ARS_NMC0074-2), though it shared~90% of its backbone with the IncA/C2 plasmid described above from K. pneumoniae ST147 strain 13ARS_VSM0593 (Supplementary Fig. 10). This plasmid backbone was found in E. coli ST410 genomes from the Philippines, but not in international genomes (Fig. 6b).
Altogether, our results show that that the Philippine E. coli ST410 genomes represent diverse lineages of the global circulating population, with evidence of frequent international introductions. These lineages are characterized by a diverse repertoire of carbapenemase genes and variants, amassed via a combination of conserved international plasmids and locally circulating plasmids.

Discussion
National networks of laboratory-based surveillance are a key pillar within the Global Action Plan to combat AMR, with the reference laboratory playing a central role in the dissemination of surveillance data to local, national and international stakeholders. The Philippines ARSP surveillance data have shown increasing resistance trends for key bug-drug combinations (Figs. 1b and  2a). In this study, the analysis of the dynamics of RPs showed the concomitant expansion of specific resistance phenotypes (Fig. 2b). By complementing laboratory data with WGS and linking the operational units of laboratory and genomic surveillance, we revealed a diversity of genetic lineages, AMR mechanisms and vehicles underlying the expansion of carbapenem resistance phenotypes (Figs. 1a and 3). The combination of these three components vastly improved our understanding of the drivers of carbapenem resistance in the Philippines at different geographical scales, from an plasmid-driven local outbreak of K. pneumoniae ST340 (Fig. 4), to the expansion of a K. pneumoniae ST147 clone carrying NDM-1 in an IncA/C2 plasmid without any close match in international genomes (Fig. 5), to the independent introductions of international epidemic clone E. coli ST410 that can acquire plasmids of local circulation (Fig. 6).
A crucial aspect of AMR surveillance is the detection of the emergence and spread of high-risk clones. In traditional laboratory surveillance, RPs can be analyzed with spatiotemporal algorithms to detect hospital outbreaks 25 or dissemination of phenotypic subpopulations between hospitals 26 . Cluster detection based on RPs depends on consistent and complete antibiograms within and across sentinel sites. During the course of our project it became apparent that there were inconsistencies in the panels of antibiotics tested by the sentinel sites (missing data or different antibiotics tested), which led to the exclusion of samples from the analysis of the dynamics of RPs (Fig. 2). These inconsistencies would also hamper cluster detection based on RPs, or the early detection of novel emerging RPs that could act as indicators of novel high-risk clones. The joint discussions of these observations by the two collaborating teams led to an effort to reinforce the standardization of antibiotic testing across the ARSP sentinel sites, which was coordinated by the ARSRL.
Even with complete and comprehensive susceptibility testing, RPs can only provide a coarse view of the spread of AMR, and sometimes lead to clusters of disparate genetic relatedness ( Fig. 3 and Supplementary Fig. 2). Integration of WGS into laboratorybased AMR surveillance can substantially improve the detection of high-risk clones by providing a high-resolution picture of genetic lineages (strains), AMR mechanisms (genes/mutations) and vehicles (MGEs, Fig. 1a). We identified high-risk clones within known international epidemic clones at the local, national and international scales by linking clonal relatedness, geographic clustering, epidemiological data and gene and plasmid content with interactive web tools 27 .
At the local scale, we identified a plasmid-driven NICU outbreak of carbapenem-resistant K. pneumoniae. The epidemiological data captured in WHONET was key to support the interpretation of the genomic findings, triggering a retrospective investigation. This previously undetected outbreak was traced to an IncFII plasmid carrying bla NDM-1 in the genetic context of ST340 (Fig. 4a). Yet, the endemic IncFII plasmid was found across multiple wards, and genetic backgrounds (STs, Fig. 4b), indicating a major role in the burden of carbapenem resistance K. pneumoniae in this hospital. A second IncFIB-IncFII plasmid found only within the ST340 lineage might have contributed to the persistence and transmission of this clone in the NICU, by carrying several genomic islands with a role in survival in the host and in the environment ( Supplementary  Fig. 5). ST340 is a member of the drug-resistant clonal complex 258 and has been reported to cause outbreaks worldwide 28 . Our findings were disseminated to the local stakeholders (i.e., hospital) via a forum with NICU staff, paediatricians and ARSRL representatives. Infection control strategies were implemented, which ultimately bolstered the infection control team of this hospital. Control of healthcare-associated infections is crucial to containing

LR697124
The replicons and AMR genes were identified with the PlasmidFinder and ResFinder databases, respectively. The match accession indicates the accession code when at least one match of more than 90% query sequence coverage and 99% identity was found in the NCBI nucleotide database.
NA not available.
the spread of AMR 29 . The roadmap from sequence data to actionable data for infection control and prevention within a hospital has been mapped by multiple use cases 30,31 and supported by recent studies on implementation and cost-effectiveness in high-resource settings 32,33 . This study makes a case for extending this roadmap to LMICs that have infection control capacity in place. At the national scale, the combined information on strain, carbapenemase gene/variant and vehicle shows a granular picture of carbapenem-resistant K. pneumoniae ST147 that uncovers the geographic distribution of high-risk clones. ST147 is an international epidemic clone that causes carbapenem-resistant infections mediated by NDM, OXA-48-like and VIM carbapenemases 28,34 . Notably, ST147 was not reported in a recent study of K. pneumoniae from seven healthcare facilities across South and South East Asia 35 , which did not include the Philippines. Within the seemingly well-established genetic background of clade IV (Fig. 5), a high-risk clone characterized by the presence of bla NDM-1 within plasmid p14ARS_VSM0843-1 displayed local clonal expansion in one sentinel site, while another high-risk   Table 1. This interactive view is available at https://microreact.org/project/ARSP_KPN_ST147_2013-14/4c33ace7. The maximum-likelihood tree was inferred from 809 SNP positions identified by mapping the genomes to reference MS6671 (LN824133.1 [https://www.ncbi. nlm.nih.gov/nuccore/LN824133.1]) and masking regions corresponding to mobile genetic elements and recombination. The scale bar shows the number of SNPs per variable site. The distribution of plasmids with NDM genes was inferred by mapping the short reads of the genomes to the complete plasmid sequences, and a match was counted when the reads covered at least 95% of the sequence length with at least 5x depth of coverage. The full data are available at https://microreact.org/project/ARSP_KPN_ST147_2013-14.
NATURE COMMUNICATIONS | https://doi.org/10.1038/s41467-020-16322-5 ARTICLE NATURE COMMUNICATIONS | (2020) 11:2719 | https://doi.org/10.1038/s41467-020-16322-5 | www.nature.com/naturecommunications clone characterized by the presence of bla NDM-1 within plasmid p13ARS_GMH0099 had also disseminated geographically across five different sites (Fig. 5). Clonal expansion followed by geographical dissemination is usually attributed to the acquisition of an AMR determinant 36 . The extent of the geographical dissemination of the different high-risk clones could be attributed to the different plasmid backbones (vehicle), but differences in the genetic background (strain) could also be at play. The presence of a robust AMR surveillance network such as the ARSP is key for the detection of geographical dissemination of high-risk clones at a regional/national level, and for the dissemination of the information to national stakeholders. In this respect, the ARSRL is currently evaluating formats for the inclusion of genomic data into their annual surveillance reports. The reference laboratory directors can also alert hospitals across the country to establish concerted infection control measures, as well as relay this information to the national DOH for the formulation of evidencebased guidelines.
At the global scale, we identified several high-risk clones of E. coli ST410 circulating the in the Philippines, with the evidence of frequent international transmission, in agreement with a previous report 18 . Previous reports of carbapenem-resistant E. coli ST410 from the Western Pacific and South East Asia regions described the presence of bla OXA-181 in China 17,37 , bla NDM-1 in Singapore 38 , bla NDM-5 in South Korea 39 , or the combination of bla NDM-5 and bla OXA-181 in Myanmar 40 and South Korea 39 . The repertoire of carbapenemases in the retrospective collection of ST410 from the Philippines seemed to be in tune with the genes/variants circulating locally, as most isolates carried the prevalent variant bla NDM-1 . Of note was the high-risk clone carrying both bla NDM-1 and bla OXA-181 , a combination hitherto unreported in E. coli ST410, and acquired through the combination of a plasmid of local circulation and a conserved international plasmid, respectively. This suggests that E. coli ST410, a successful pandemic lineage, can not only easily disseminate, but it can also adjust the complement of carbapenemases by acquiring endemic plasmids along the way. Global AMR surveillance networks, such as the WHO GLASS 8 , are paramount to detect the emergence and monitor the spread of resistance at the international level and inform the implementation of targeted prevention and control programmes. The ARSRL has been supplying laboratory surveillance data to GLASS since 2015, and presented their experience with WGS at a WHO technical consultation on the application of WGS for national, regional and global surveillance in 2020, thus becoming an international contributor in both fields.
Collecting information on AMR that can be rapidly transformed into action requires harmonized standards, especially at the national and international levels 8 . For laboratory data, WHONET serves this purpose in the ARSP and in over 2000 clinical, public health, veterinary and food laboratories in over 120 countries worldwide 26 , while GLASS serves this purpose at the global level. Genomic data are amenable to standardization 29 , and national public health agencies 41 and international surveillance networks 42 are adopting diverse schema to identify and name genetic lineages. However, a global standard system for defining a cluster has not been implemented beyond the level of discrimination of MLST. Likewise, different databases of AMR mechanisms [43][44][45] may differ in content and nomenclature. Thus, the standardization of genomic data, as well as the provision of platforms, for the uptake of genomic information is crucially moving forward, and for ongoing AMR surveillance.
Containment of AMR at a global level requires an international concerted effort. High-risk clones have the propensity to disseminate rapidly, and genomics can improve the detection of their emergence and spread. This highlights the increasing need to build equitable partnerships to facilitate ownership transfer of genomic epidemiology capacity (operational, analytical and interpretational) to enhance national AMR surveillance programmes within low-resource settings. The binational partnership of this project led to a large retrospective survey of bacterial pathogens that has provided the first in-depth genomic view of the AMR landscape in the Philippines and established contextual data for ongoing local prospective sequencing. In parallel, through training and transfer of expertize in laboratory procedures and bioinformatics, the open exchange of data, and collective interpretation of results, we have expanded the existing capacity of a national reference laboratory with WGS focused on action for public health. WGS commenced at the ARSRL in 2018 with the Illumina MiSeq equipment available locally. A new dedicated bioinformatics server was installed at ARSRL for sequence data storage and analysis. Collective data interpretation was aided by data sharing via interactive web tools such as Microreact (www.microreact.org) 27 and Pathogenwatch (www. pathogen.watch) to leverage the expertize from both CGPS and ARSRL. The results of the retrospective survey were presented to the representatives of the sentinel sites during the ARSP annual meetings and at international conferences, at first by the members of CGPS, followed by presentations by the members of the ARSRL in subsequent years 46,47 . The ARSRL staff conducted its first outbreak investigation using WGS in July 2019, with

NA LR697132
The replicons and AMR genes were identified with the PlasmidFinder and ResFinder databases, respectively. The match accession indicates the accession code when at least one match of more than 90% query sequence coverage and 99% identity was found in the NCBI nucleotide database.
NA not available.   Table 1. This interactive view is available at https://microreact.org/project/ARSP_ECO_ST410/088ba65b. b Philippine isolates (orange nodes) in global context. This interactive view is available at https://microreact.org/project/ARSP_ECO_ST410_GLOBAL/c701506e. The maximum-likelihood trees were inferred from 703 (A) and 2851 (B) SNP positions, respectively, identified by mapping the genomes to reference AMA1167 (CP024801.1) and masking regions corresponding to mobile genetic elements and recombination. The scale bars show the number of SNPs per variable site. The distribution of plasmids with carbapenemases genes in (b) was inferred by mapping the short reads of the genomes to the complete plasmid sequences, and a match was counted when the reads covered at least 95% of the sequence length with at least 5× depth of coverage. The full data are available at https://microreact.org/project/ARSP_ECO_ST410 (a) and https://microreact.org/project/ARSP_ECO_ST410_GLOBAL (b).
sequencing, bioinformatics, and reporting to the DOH conducted locally, with a turnaround time of 8 days from sample receipt to report preparation. While progress has been made, significant challenges remain to bring genomics into routine surveillance in low-resource settings. These include, but are not limited to, challenges in the supply chain and procurement of WGS reagents and equipment, vastly differing costs between high-and low-income settings, shortage of skilled local bioinformaticians due to both limited access to training and difficulties in staff retention, and the lack of platforms to feedback the actionable genomic data to sentinel sites. In addition, it was expected that the timeline from implementation of a disruptive technology such as WGS to full adoption by all the relevant stakeholders was to exceed the period of one academically funded project. The continued collaboration with the ARSRL (via an NIHR-funded project) is building on the foundations laid by this study to tackle some of the challenges mentioned above, and to improve on the translation of genomic data into public-health action and policy. We laid the foundations for the sustainable implementation of WGS and genomic epidemiology within national surveillance networks in LMICs, which can be extended to other locations to tackle the global challenge of AMR.

Methods
Ethics statement. Ethical approval is not applicable as per the Research Institute for Tropical Medicine-Institutional Review Board. This study uses archived bacterial samples processed by ARSP. No identifiable data were used in this study.
ARSP workflow and data. The ARSP implements AMR surveillance on aerobic bacteria from clinical specimens. The programme collects culture and antimicrobial susceptibility data from its 24 sentinel sites and 2 gonorrhoeae surveillance sites. The sentinel sites participate in an external quality assessment scheme conducted by the reference laboratory to ensure quality of laboratory results. Case finding is based on priority specimens sent routinely to the hospitals' laboratories for clinical purposes, and is thus based on the diagnostic practices of the clinicians, as well as the resources available to the patient to cover the culture and susceptibility testing. Sentinel sites submit monthly results of all culture and susceptibility tests to the ARSRL via WHONET, a free Windows-based database software for the management and analysis of microbiology laboratory data with a special focus on the analysis of antimicrobial susceptibility test results 6 . Multi-drug-resistant phenotypes are classified according to recommended standard definitions 48 . Select isolates are referred to the ARSRL for confirmatory testing if they fall in one of the three following categories: (1) organisms of public health importance regardless of their resistance phenotype, such as Salmonella enterica ser. Typhi, Streptococcus pneumoniae and Vibrio cholerae; (2) organisms with emerging resistance, such as carbapenem resistance; (3) organisms with rare susceptibility patterns, such as vancomycin-resistant Staphylococcus aureus or ceftriaxone-resistant Neisseria gonorrhoeae. At the ARSRL, all referred isolates are re-identified using standardized culture, and susceptibility tested by disc diffusion or minimum inhibitory concentration (MIC) determined using the most current guidance for antibiotic test panels and breakpoints recommended by the Clinical and Laboratory Standards Institute. Serotyping for S. pneumoniae, H. influenzae, Salmonellae and V. cholerae is conducted to further characterize the isolates.
ARSP data and bacterial isolates included in this study. Phenotypic data collected since 1988 from antimicrobial susceptibility tests were summarized with WHONET v5.6 to compute yearly resistance rates for key pathogen-antibiotic combinations. Only the first bacterial isolate per patient per species per year was included in the calculations.
RPs derived from the antibiograms were summarized per organism with WHONET v5.6 6 , and their relative abundance visualized with Tableau Desktop. Data collected between 2009 and 2017 were included in the analysis as the number of sentinel sites collecting data remained relatively stable during this period (22-24 sites). Isolates with missing data for any antibiotic on the panels listed on Table 1 were excluded.
The bacterial pathogens included in the See and Sequence study were those responsible for the majority of antimicrobial-resistant infections in the Philippines. These included carbapenemase-producing Pseudomonas aeruginosa and Acinetobacter baumanii, carbapenemase-producing and ESBL-suspect Escherichia coli and Klebsiella pneumoniae, methicillin-resistant Staphylococcus aureus, Neisseria gonorrhoeae, non-typhoidal Salmonella, and Salmonella enterica ser. Typhi. Linked epidemiological data included location and date of specimen collection, type of specimen, type of patient (in or outpatient) and sex and age of the patient. We utilized a proxy definition for infection origin whereby patient first isolates collected in the community or on either of the first 2 days of hospitalization were categorized as community infection isolates, while isolates collected on hospital day 3 or later were categorized as HAI isolates 8 .
In this article, we focus on carbapenemase-producing organisms. All carbapenemase-producing K. pneumoniae and E. coli isolates referred to, confirmed, and banked by the ARSRL in 2013-2014 were selected for the retrospective sequencing survey. Approximately 100 isolates of carbapenemaseproducing P. aeruginosa and A. baumannii each were selected according to the following criteria: (1) referred to ARSRL in 2013-2014; (2) complete RP (i.e., no missing data); (3) overall prevalence of the RP in the ARSP data (including both referred and non-referred isolates); (4) geographical representation of different sentinel sites. The number of isolates included from each sentinel site was proportional to their relative abundance and estimated from (n/N) × 100 (rounded up), where n is the total number of isolates from one site, and N is grand total of isolates; (5) when both invasive and non-invasive isolates representing a combination of RP, sentinel site and year of collection were available, invasive isolates (i.e., from blood, or cerebrospinal, joint, pleural and pericardial fluids) were given priority. In addition,~100 isolates of ESBL-producing E. coli and K. pneumoniae each were included as per the criteria above. Pan-susceptible isolates were only available for P. aeruginosa isolates and were also included in the retrospective sequencing survey.
Lyophilisates of the selected bacterial isolates were re-suspended in 0.5 mL of tryptic soy broth, plated onto MacConkey agar plates and incubated at 35-37°C for 18-24 h. The species identity of the revived bacterial isolates was confirmed by characterization of the colony morphology and conventional biochemical tests (e.g., oxidase, motility and sugar fermentation using triple sugar iron agar test). Phenotypic re-testing was performed on the resuscitated isolates based on the phenotype stored in WHONET. Carbapenemase production in resuscitated isolates of Enterobacteriaceae was tested with the modified Hodge test, and confirmed using the imipenem MIC method with imipenem Etest strips (BioMerieux). P. aeruginosa and A. baumanii were screened for metallo-beta-lactamase production using the EDTA synergy test between imipenem disk (10 µg) and EDTA (750 µg/ mL), which was confirmed by imipenem/imipenem inhibitor Etest strips. E. coli and K. pneumoniae isolates were screened for ESBL with the double-disk synergy test between aztreonam (30 µg) and amoxicillin/clavulanic acid (20/10 µg) and confirmed by using the gradient diffusion method with cefotaxime/cefotaxime with clavulanic acid or ceftazidime/ceftazidime with clavulanic acid Etest strips.
Whole-genome sequencing, assembly and annotation. A. baumanni, E. coli, K. pneumoniae and P. aeruginosa bacterial strains were shipped to the Wellcome Sanger Institute streaked onto Nutrient Agar butts contained in screw-cap cryovials, where DNA was extracted from a single colony of each isolate with the The replicons and AMR genes were identified with the PlasmidFinder and ResFinder databases, respectively. The match accession indicates the accession code when at least one match of more than 90% query sequence coverage and 99% identity was found in the NCBI nucleotide database.
NA not available.
QIAamp 96 DNA QIAcube HT kit and a QIAcube HT (Qiagen; Hilden, Germany). DNA extracts were multiplexed and sequenced on the Illumina HiSeq platform (Illumina, CA, USA) with 100-bp paired-end reads. The quality of the sequence data was assessed for median base quality >30, mapping coverage >80% of the length of a reference genome sequence of the same species (A. baumannii ATCC 17978, E. coli UPEC_ST131, K. pneumoniae subsp. pneumoniae Ecl8 and P. aeruginosa LESB58), and >80% of the sequence reads assigned to the corresponding species with Kraken v0.10.6-a2d113dc8f 49 and a reference a database of all viruses, archaea and bacteria genomes in RefSeq 50 and the mouse and human reference.
Annotated assemblies were produced using the pipeline described in 51 . For each sample, sequence reads were used to create multiple assemblies using VelvetOptimiser v2.2.5 52 and Velvet v1.2 53 . An assembly improvement step was applied to the assembly with the best N50 and contigs were scaffolded using SSPACE v2.0 54 and sequence gaps were filled using GapFiller v1.11 55 . Automated annotation was performed using PROKKA v1.5 56 and a genus-specific database from RefSeq 50 . The quality of the assemblies was assessed for total size, number of contigs, N50 and GC content appropriate for each species, with all genomes passing quality control characterized by assemblies of <150 contigs and N50 > 60,000.
Phylogenetic analysis. Evolutionary relationships between isolates were inferred from SNPs by mapping the paired-end reads to the reference genomes of A. baumannii A1 (CP010781), E. coli EC958 (HG941718), K. pneumoniae subsp. pneumoniae (AP006725) or P. aeruginosa LESB58 (FM209186), using the Burrows Wheeler aligner v0.7.12 57 to produce a BAM file. PCR duplicate reads were identified using Picard v1.92 58 and flagged as duplicates in the BAM file.
Variation detection was performed using samtools mpileup v0.1.19 59 with parameters -d 1000 -DSugBf and bcftools v0.1.19 60 to produce a BCF file of all variant sites. The option to call genotypes at variant sites was passed to the bcftools call. All bases were filtered to remove those with uncertainty in the base call. The bcftools variant quality score was required to be >50 (quality < 50) and mapping quality >30 (map_quality < 30). If not all reads gave the same base call, the allele frequency, as calculated by bcftools, was required to be either 0 for bases called the same as the reference, or 1 for bases called as an SNP (af1 < 0.95). The majority base call was required to be present in at least 75% of reads mapping at the base (ratio < 0.75), and the minimum mapping depth required was 4 reads, at least two of which had to map to each strand (depth < 4, depth_strand < 2). Finally, strand_bias was required to be <0.001, map_bias < 0.001 and tail_bias < 0.001. If any of these filters were not met, the base was called as uncertain.
A pseudogenome was constructed by substituting the base call at each site (variant and non-variant) in the BCF file into the reference genome, and any site called as uncertain was substituted with an N. Insertions with respect to the reference genome were ignored and deletions with respect to the reference genome were filled with Ns in the pseudogenome to keep it aligned and the same length as the reference genome used for read mapping. For sequence-type-specific phylogenies, MGEs and recombination regions detected with Gubbins v1.4.10 61 were masked from the alignment of pseudogenomes produced by mapping against the genomes of K. pneumoniae strain CAV1217 (ST340, CP018676.1), K. pneumoniae strain MS6671 (ST147, LN824133.1) or E. coli strain AMA1167 (ST410, CP024801.1). Alignments of only variant positions were inferred with snpsites v2.4.1 62 and were also used to compute pairwise SNP differences between isolates (https://github.com/simonrharris/pairwise_difference_count).
Maximum-likelihood phylogenetic trees were generated with RAxML v8.2.8 63 based on the generalized time reversible model with GAMMA method of correction for among-site rate variation and 100 bootstrap replications. Trees were midpoint rooted.
The phylogenetic trees, genomic predictions of AMR and genotyping results were visualized together with the ARSP epidemiological data using Microreact 27 .
In silico predictions of AMR determinants and plasmid replicons. Known antibiotic resistance determinants were identified using Antibiotic Resistance Identification By Assembly v.2.6.1 (ARIBA) 64 . For acquired genes an in-house curated database comprising 2316 known resistance gene variants was used (https://figshare.com/s/94437a301288969109c2) 65 . Full-length assembled matches with identity >90% and coverage >5× were considered as evidence of the presence of the gene in the query genome. Further manual inspection of the results excluded 13 matches that passed the above filters but that were suspected low-level contamination. Mutations in ompK35 and ompK36 genes were detected with ARIBA and the Comprehensive Antimicrobial Resistance Database 44 and the output was inspected for truncations, interruptions and frameshifts.
Plasmid replicons were identified with ARIBA and the PlasmidFinder database 66 . Full-length assembled matches with 100% identity were considered as evidence of the presence of the replicon gene.
MLST calls from P. aeruginosa genome sequences were manually curated as we noted the presence of two genes exhibiting 100% sequence identity to different acsA alleles in the pubmlst database (accessed April 2019), which confounded the in silico predictions. This was not privy to the genomes is this study, since the genome sequence of reference strain NCGM2_S1 (accession AP012280), known to belong to ST235, showed two genes with locus tags NCGM2_1665 (nt position 1808730-1810685) and NCGM2_0806 (nt positions c879796-881733) that are annotated as acsA and acsB, respectively, and are 100% identical to alleles acsA_38 and acsA_225, respectively (https://pubmlst.org/data/alleles/paeruginosa/acsA.tfa, accessed Apr 2019). Furthermore, the same was only observed for the acsA allele both from assemblies (MLST_check) and raw reads (ARIBA-MLST), and therefore unlikely the result of contamination with other strains. Allele acsA_38 defines ST235, but allele acsA_225 defines a novel ST. Thus, when two different acsA alleles were detected in a P. aeruginosa genome we selected the one associated with a known ST.
The presence of AMR genes and plasmid replicons in the plasmid sequences was detected with ResFinder 45 and PlasmidFinder 66 , respectively. Matches with identity larger than 98% were reported. Detailed annotations of plasmids were obtained with MARA 76 . Visual comparisons between plasmid sequences were created with BRIG v1.0 77 with default BLAST v2.7.0 parameters.
The distribution of plasmids across draft genomes was inferred by mapping short reads to the plasmid sequences with smalt v0.7.4 with identity threshold y = 96% and random alignment of reads with multiple mapping positions of equal score (-r 1). The BAM file was sorted and indexed with samtools v0.1.19-44428cd 59 , which was also used to and generate a pileup (mpileup -DSugBf) of the alignment records in BCF format. A pseudo-plasmid sequence in fastq format was constructed with the bcftools view and vcfutils.pl vcf2fq programmes, and the sequence length coverage (%) of the reference plasmid was computed.
Reporting summary. Further information on research design is available in the Nature Research Reporting Summary linked to this article.