Emergence of mcr-9.1 in Extended-Spectrum-β-Lactamase-Producing Clinical Enterobacteriaceae in Pretoria, South Africa: Global Evolutionary Phylogenomics, Resistome, and Mobilome

Colistin is currently the last-resort antibiotic for difficult-to-treat bacterial infections. However, colistin resistance genes that can move from bacteria to bacteria have emerged, threatening the safe treatment of many bacterial infections. One of these genes, mcr-9.1, has emerged in South Africa in bacteria that are multidrug resistant, further limiting treatment options for clinicians. In this work, we show that this new gene is disseminating worldwide through Enterobacter and Salmonella species through multiple plasmids. This worrying observation requires urgent action to prevent further escalation of this gene in South Africa and Africa.

Subsequently, colistin was reintroduced into clinical medicine as a last-resort agent to manage infections resistant to carbapenems, other ␤-lactams, and other antibiotics (5,6).
Although the WHO has classified ESBL-and carbapenemase-positive Enterobacteriaceae as critical and high-priority pathogens, relatively less attention has been given to rarely isolated pathogenic species such as Citrobacter freundii, Enterobacter hormaechei, Klebsiella variicola, Providencia spp., and Proteus spp. However, these relatively less isolated species can harbor important resistance genes that can be transferred to commonly isolated species such as Klebsiella pneumoniae, Escherichia coli, and Salmonella enterica (4,11,(20)(21)(22). Hence, we investigated the resistome, mobilome, and genomic epidemiology of such less-isolated species to increase our understanding of their resistance mechanisms and epidemiology.
In this report, we describe for the first time the emergence of the mcr-9.1 colistin resistance gene in South Africa and Africa as well as the cooccurrence of mcr-9.1 with ESBLs in Enterobacter hormaechei. We further undertake a global phylogenomic analysis of C. freundii, Enterobacter hormaechei, Klebsiella variicola, and Providencia/Proteus spp., along with a global evolutionary analysis of mcr-9.1. These findings call for advanced One Health genomic surveillance to contain such multidrug-resistant strains from further dissemination.

RESULTS
Demographics, genome characteristics, and resistance phenome. The specimens were mainly sourced from six males and four females, aged between 18 and 61 years and hospitalized at a single tertiary academic hospital in Pretoria, South Africa (Table 1). Except for the carbapenems, cephalothin, and combinations of clavulanate and tazobactam with cefotaxime/ceftazidime and piperacillin, respectively, most of the isolates were resistant to all the ␤-lactams. Moreover, reduced resistance to piperacillintazobactam (in isolates EC009, EC010, and K001) and sensitivity to aztreonam and cefepime (CF004) and cefoxitin (EC015, K001, K006, K130, K063, and PM005) were observed. All the isolates were phenotypically positive for ESBL production, except for CF004 (see Table S1 in the supplemental material; Table 1). None of the isolates was  (Table S1). The draft genome sizes of the genomic sequences ranged from 4 Mb to 6 Mb, with PM005 having the lowest GC content (39.2%) and number of coding sequences (3,542). The genome coverage ranged from 90ϫ to 102ϫ, while CRISPR arrays were identified in only two isolates (K001 and K063) ( Table 1; Data Set S1).
The isolates were all multidrug resistant, and agreements as well as discrepancies between the resistance phenotypes (phenomes) and resistance genes (resistomes) were observed (Table S1). In particular, CF003 and CF004 had only bla CMY and aac(6')-If or fosA7 genes but expressed multidrug resistance to several antibiotics for which no resistance genes were identified. Such observations were also made with the presence of fosA genes in the isolates and the absence of fosfomycin resistance in almost all the strains as well as the presence of chromosomal AmpCs, such as ACT, LEN, and CMY in EC015, K001, K006, K130, K063, and PM005, which were all susceptible to cefoxitin (Table S1). Notable was the absence of colistin resistance in K006 and K130, which had the mcr-9.1 gene; comparative genomic analyses showed that the colistin-susceptible and -resistant strains shared the same chromosomal mutations in mgrB and pmrB ( Table 2). The absence of carbapenem resistance was corroborated by the absence of any carbapenemase gene, while the ESBL-positive phenotypes agreed with the ESBL genes identified in the genomes (Table S1; see Fig. 2).
Resistome, mobilome, and evolutionary epidemiology of mcr-9. CF003 and CF004 had the least number of resistance genes, with no ESBL genes. Dominant ESBLs in the other nine isolates included CTX-M-15, TEM-1, OXA-1, OXA-9, and ACT (Table S1 and Data Set S1). The genetic environment of the contigs on which the AmpC genes such as bla CMY , bla ACT , bla LEN , and bla LAP were located strongly suggested their presence on chromosomes (Fig. 1A); this was confirmed by subjecting the contigs to a BLAST search on GenBank, and the sequence percent identities of the AmpC contigs aligned closely with only chromosomes. Notably, the bla CMY in PM005 was in close synteny with ISEc9 (Fig. 1A, panel III). bla CTX-M-15 and bla TEM-1 ESBLs were almost always found on the same contig, bracketed by ISEc9 and a Tn3 transposon and IS19 and an integrase/recombinase, respectively. In particular, aph(6')-Id::aph(3=')::sul2 resistance genes were commonly found joined to IS19, just after bla TEM-1 in mcr-9-negative E. hormaechei (Fig. 1A). In PM005, K006, and K130, the genetic environment of bla TEM-1 was different (Fig. 1A). The bla OXA gene was commonly found between catB3 and aac(6')-Ib on a class 1 integron as gene cassettes or directly on the chromosome, similar to what was observed for bla SHV, bla LEN , and bla ACT. (Fig. 1A).
A recent allele of the mobile colistin resistance gene, mcr-9.1 (Fig. 1B), was identified in three E. hormaechei strains, the first to be discovered from Africa, dating back to 2013. mcr-9.1 genes were all found in the immediate environment of a cupin fold metalloprotein, WbuC, similar to the bla CTX-M-15 gene. mcr9.1 only occurred once with The reference genome used was Enterobacter hormaechei strain C15117 (PRJNA494598). b -, intrinsic resistance.
bla CTX-M-15 in K063 (Fig. 1A). A nucleotide BLAST search of the contigs, namely, K006 contig 131, K063 contig 53, and K130 contig 60, and subsequent phylogenetics analysis using the fast minimum evolution method showed that the contigs harboring the mcr9.1 gene were most closely aligned to plasmid genomes from different Enterobacteriaceae species, with a few chromosomes such as that of Enterobacter kobei strain DSM 13645 being also closely aligned (Fig. 2). Sequence alignment showed close similarity between the three mcr-9-bearing contigs (Data Set S2), particularly between the contigs of K006 and K130 ( Fig. 2A and B). Notably, both K006 and K130 aligned with very similar nucleotide identities with the same genomes, although minor differences existed in their nucleotide sequences (Data Set S2). Both K006 and K130 were phylogenetically distant from K063 ( Fig. 2A and B), which is further evinced by the different genomes to which they both most closely aligned and the distance of those genomes from each other in the distance trees ( Fig. 2A to C). The evolutionary epidemiology of the mcr-9.1-bearing plasmids or chromosomes is further shown in Fig. 2A to C, with different plasmid types and chromosomes mediating the spread of this gene within and between species throughout the world. The various mcr-9.1-bearing genomes (plasmids or chromosomes) clustered phylogenetically into six clades, herein labeled A to F, with Enterobacter sp. strain 18A13 plasmid pECC18A13-1 DNA from Japan seeming to be the earliest ancestor and Enterobacter hormaechei strain S13 plasmid pSHV12-1301491 from the United States seeming to be the most recent. All the mcr-9.1 contigs in this study fell within phylogenetic cluster B, although the genomes that aligned most closely with these contigs changed their phylogenetic clustering between K006/K130 and K063 ( Fig. 2A to C).
Several integrons bearing gene cassettes of antibiotic resistance genes were identified in various contigs in the isolates (Table 3), with none being detected in the C. freundii strains. PM005 (Providencia alcalifaciens) had a unique integron, In27, while K001 (K. variicola) and four E. hormaechei strains harbored In191. K006, K063, and K130, originating from different patients and having the mcr-9.1 gene, harbored the same integrons and very similar gene cassettes. EC009 and EC010 had the same class 1 integron and gene cassettes and belonged to the same clone ( Table 3). Most of the isolates had multiple plasmid replicons, with ColRNAI being the most common; ColRNAI was also found in the mcr-9.1-positive strains. The isolates were further grouped using their plasmid multilocus sequence types (pMLST). Except for CF004 and PM005, which Pathogenomics of mcr-9.1-and ESBL-Positive Enterobacteriaceae had no pMLST but contained an A/C 2 replicon, IncHI2[ST-1] and IncF subtypes were the main plasmid types. Notably, all the mcr-9.1 strains had an IncHI2[ST-1] plasmid, while two, K006 and K063, coharbored IncF subtype plasmids. Although the resistance genes within strains having the same plasmid types were different, core resistance genes with same genetic environment were present ( Fig. 1A; Table 4).
Evolution of mcr genes. The phylogenetic and evolutionary relationships between the 10 mcr variants are shown in Fig. 1B (see Fig. S1 in the supplemental material), with the mcr-9 genes being clustered into two clusters, I and II. Cluster I of mcr-9 was closely related to that of mcr-10, mcr-3, mcr-7, and mcr-8, while cluster II was closely related to that of mcr-5. Notably, strains containing mcr-1, one of which was used as the reference strain, were closely related to those containing mcr-2 and mcr-6, with mcr-2 and mcr-6 clustering closely together. Particularly, there were more subclades within cluster II than cluster I of the mcr-9 genes, which is due to the sequence differences between the genes in cluster II.
Phylogenomics. The C. freundii strains were not related to any South African or African strains but were related to American and European strains. CF004 and CF003 were not related to each other, but CF004 was of the same clone as C. freundii , and CRE20 (bla CMY ) from the United States (Fig. 3).
The closest evolutionary relative of K001 was a single strain, WUSM_KV_44 (bla LEN-13 , fosA, oqxA, oqxB) from the United States (Fig. 5), while PM005, which was initially identified as Proteus mirabilis by MicroScan but later revised by NCBI's ANI to P. alcalifaciens, was not related to any Providencia sp. genome ( Fig. 6A and B); however, it was closely related to Proteus mirabilis strain Pr2921 from Uruguay [catA, tet(J)] (Fig. 6C).

DISCUSSION
We here report on the first emergence of mcr-9.1 in both South Africa and Africa in three E. hormaechei strains with a rich repertoire of resistomes and mobilomes. The mcr-9.1 strains were collected among clinical Enterobacteriaceae strains during a molecular surveillance procedure to identify ESBL producers in a referral laboratory. It is a worrying observation that all the strains, sourced from different patients and wards within the same hospital, were multidrug resistant. Although the resistomes largely reflected the resistance phenotype, resistance to antibiotics was observed without an underlying resistance gene being identified. Specifically, the C. freundii strains expressed resistance to several antibiotics, although not more than three resistance genes were found in them. This could suggest the presence of an unknown resistance determinant or the use of active efflux. Furthermore, sensitivity to certain antibiotics such as colistin was not corroborated by the presence of mcr-9.1. Sensitivity to colistin in strains having the mcr gene has been reported previously, and mutations in chromosomal colistin resistance genes have been shown to exert a higher resistance MIC than mcr genes (3,(23)(24)(25).
Only a single strain with the mcr gene was colistin resistant, and the identified mutations in mgrB and pmrAB were also found in susceptible strains. Hence, the identified mutations cannot be responsible for colistin resistance. Consequently, other factors identified for colistin resistance or other unknown determinants could be responsible for the observed colistin resistance (3,7,11,22). The absence of colistin resistance in mcr-9.1 strains is not new or surprising, as it has been observed in other studies (23,24,26). However, it should be known that in the presence of certain promoters in some bacterial species, these silent genes can become actively expressed, leading to phenotypic resistance (10, 23). Furthermore, it has been suggested that a reduction in the colistin breakpoints might increase the detection of otherwise colistin-

Sample code (MLST) Integron(s) Cassette arrays (location on contig) a CF003
None None CF004 None None EC001 (ST-459 resistant strains that are classified as susceptible, although this reduction must be done after additional investigations (26). The three strains bearing the mcr-9.1 gene were closely related, but the contigs bearing the mcr-9.1 genes did not have 100% nucleotide sequence identity among themselves (see Data Set S2 in the supplemental material), with the mcr-9.1 contigs' distance trees showing that K006 and K130 were of a more recent common ancestor (Fig. 2). It is worth noting that K006 and K130 were phylogenetically distant but that their mcr-9.1 contigs were of closer nucleotide identity than those of K006 and K063, which were of the same clone (Fig. 4). This observation, plus the close alignment of the mcr-9.1 contigs with plasmid genomes in GenBank (Fig. 2), strongly suggests that the mcr-9.1 gene might have been horizontally, instead of vertically, acquired. It also suggests that either the patients (one was female, and the other was male) from which K006 and K063 were obtained got the mcr-9.1 gene from the same source or that one of the two patients transferred it to the other directly, through a health care worker, or via the hospital environment.
The clustering dynamics of the mcr genes confirm their nomenclature as mcr variants of the same allele clustered together on the same branch. For instance, all mcr-1 strains clustered together, showing their close sequence identity (Data Set S3). The evolutionary distance between the various mcr variants also confirms their sequence differences stemming from the nucleotide polymorphisms between these alleles. These differences in sequence percent identity between the various mcr alleles begs the question of whether these variants evolved from a single ancestor or emerged independently from each other. It can, however, be seen that some variants are from a more recent common ancestor than others, owing to their sequence percent identity ( Fig. 1B; Fig. S1). The phenotypic effect of these variants, in terms of the level of colistin resistance, is yet to be determined to know which of these variants confers the most resistance in terms of MIC levels. Nevertheless, most of these mcr variants were plasmid borne, with a few being chromosomal, explaining their wide distribution in diverse species (11,27).
The global evolutionary trajectory of the mcr-9 gene and genomes shows six major phylogenetic clades, herein labeled A to F. The strains and plasmids in clade A seem to be the earliest ancestors of this gene, with clade F being the last and most recent. Moreover, Enterobacter spp., particularly E. hormaechei, and Salmonella enterica plasmids are the most common hosts of the mcr-9 gene (Fig. 2). This is not surprising, as mcr-9.1 was first identified in Salmonella enterica serovar Typhimurium (23). The diversity of species and plasmids involved in the dissemination of this gene explains its promiscuity and rapid spread around the globe, further corroborating the need to restrict the use of colistin in both veterinary and human medicine (5, 28). The genetic support of the various resistance genes, particularly bla CTX-M-15 , bla TEM-1 , bla OXA , and mcr-9.1, identified in these isolates is not new. In particular, the cooccurrence of bla CTX-M-15 and bla TEM-1 within Tn3 composite transposons, ISEc9, and IS19 on IncF type plasmids is widely reported in both South Africa and worldwide (4,20,21,(29)(30)(31). Further, the presence of the mcr-9.1 gene on IncHI2 plasmids as well as the presence of a cupin fold metalloprotein, as observed here, is common around mcr-9.1 genes (23). The promiscuity of IncF plasmids and the abundance of integrons, trans-  A and B) and Proteus (C) genomes deposited in GenBank and PATRIC were analyzed and drawn into three trees, as shown in panels A, B, and C. Genomes belonging to the same subclade with the closest evolutionary distance are indicated in red-colored text labels. The various phylogenetic clades and subclades are highlighted together and uniquely to show their evolution and epidemiology. The trees were drawn using the maximum-likelihood method in RAxML, using Streptomyces sp. MUSC164 as a reference. Bootstrap values are shown on the branches. posons, and insertion sequences (ISs) in the genomes of these strains are worrying, as they might mobilize and facilitate the more rapid dissemination of these multiple resistance genes to other species and clones (4,31,32).
It is interesting to note that not all the closely related strains around the mcr-9.1 strains harbored the mcr gene. Some of the global strains with very close phyletic clustering with the non-mcr-positive strains contained mcr and, in some cases, carbapenemase (bla NDM , bla OXA-48 , bla VIM , and bla IMP ) genes. Such resistome differences between strains belonging to the same clade suggest the gain and loss of resistance plasmids during the evolutionary and epidemiological trajectory. This could be the case for the other E. hormaechei strains that had no mcr gene. The 11 strains included in this study were phylogenetically distant from any strain from South Africa or Africa but were of very close evolutionary distance to international strains, some of which harbored very rich resistomes, including the cooccurrence of carbapenemase and mcr genes.
Thus, the possibility of these strains having been imported cannot be ruled out, a situation that warrants constant surveillance and screening of medical tourists, particularly those from areas of carbapenemase and mcr gene endemicity (6,33).
Although carbapenemases were not found in these strains, as has been reported elsewhere in mcr-positive strains (34)(35)(36)(37)(38), the cooccurrence of mcr and ESBL genes is worrying, as they restrict therapeutic options. To date, only mcr-1 genes have been  reported in South Africa (6,25,39,40), making the identification of mcr-9.1 in strains as long ago as 2013 a concerning situation. This suggests that other mcr variants may be present in South Africa and Africa, and intensive clinical surveillance is necessary to unearth these and preempt further escalation.

MATERIALS AND METHODS
Clinical bacterial specimens and antibiotic sensitivity testing. An old collection of clinical specimens (n ϭ 72) that was obtained in 2013 at a referral laboratory from patients having drug-resistant infections at an academic teaching hospital were analyzed (20,21,29). These patients suffered from several different diseases, and the specimens were obtained from different body sites. The specimens were cultured on blood agar overnight, and pure colonies were transferred to freshly prepared Mueller-Hinton plates for another 24 h of incubation. The fresh cultures were then used for antimicrobial sensitivity testing and species identification on a MiscroScan Walkaway system. The isolates were further tested for ESBL production using double-disk synergy testing (41). Enterobacteriaceae species having reduced sensitivity to cephems and colistin were further selected for genomic analyses.
Whole-genome sequencing. Eleven isolates belonging to Enterobacteriaceae species that are rarely isolated at the referral laboratory were cultured for 24 h to determine the reservoirs of resistance genes in less-isolated pathogens or commensals. Their genomic DNA was subsequently isolated using the ZR Fungal/Bacterial DNA MiniPrep kit (Zymo Research, USA). Briefly, 200-bp libraries were prepared from the genomic DNA (gDNA) and size-selected using 2% agarose gel and Pippen prep (Sage Science, Beverly, MA, USA). The libraries were barcoded and pooled for sequencing on an Ion Proton system (Thermo Fisher, Waltham, MA, USA).
Bioinformatic and phylogenomic analyses. Cutadapt 2.8 was used to remove sequence adapters, artifacts, and poor sequences using default parameters. The SPAdes assembler was used to assemble the trimmed raw reads de novo using default parameters. Reads below 200 nucleotides were removed, deposited at GenBank, and annotated with PGAP (42). The resistomes of the isolates were determined from the Isolates Browser database of NCBI (https://www.ncbi.nlm.nih.gov/pathogens/isolates#/search/), while the mobilomes were manually curated from the gff files from GenBank. The isolates' multilocus sequence typing (MLST) sequence types were identified with the MLST database at CGE (http://cge.cbs .dtu.dk/services/MLST/). INTEGRALL (http://integrall.bio.ua.pt/) was used to identify the integrons and gene cassettes within the genomes, while PlasmidFinder (http://cge.cbs.dtu.dk/services/PlasmidFinder/) was used to identify the plasmid replicons.
Contigs bearing the mcr-9.1 gene were subjected to a BLAST search (using nucleotide BLAST) to identify genomes with the closest nucleotide identity; these genomes were used to draw distance trees using the fast minimum evolution method at a maximum sequence difference of 0.75. The genomes of Citrobacter freundii, Enterobacter hormaechei, Klebsiella variicola, Providencia spp., and Proteus mirabilis were curated from PATRIC (https://www.patricbrc.org) and filtered to remove poor genomes; i.e., genomes with fewer than 1,000 consensus genes for alignment with the total genome collection were discarded. The filtered genomes were then aligned using CLUSTAL and run through RAxML to draw maximum-likelihood trees (using a bootstrap reassessment of 1,000ϫ), which were annotated using Figtree (http://tree.bio.ed.ac.uk/software/figtree/). All mcr-9 genes and selected RefSeq mcr-1 to mcr-10 genes were downloaded from GenBank and aligned with MUSCLE. The aligned file was then used to draw a maximum likelihood tree using PhyML. The tree was thereafter annotated using Figtree.

SUPPLEMENTAL MATERIAL
Supplemental material is available online only. FIG S1, PDF file, 0.03 MB.