Whole-genome sequencing for the characterization of resistance mechanisms and epidemiology of colistin-resistant Acinetobacter baumannii

Background Multidrug-resistant Acinetobacter baumannii is an important causal pathogen of healthcare-associated infections, and colistin-resistant strains have recently emerged owing to the increased use of colistin. Using next-generation sequencing (NGS), a single whole-genome sequencing (WGS) protocol can identify and type pathogens, analyze genetic relationships among different pathogens, predict pathogenic transmissions, and detect antibiotic resistance genes. However, only a few studies have applied NGS in studying the resistance mechanism and epidemiology of colistin-resistant A. baumannii. This study aimed to elucidate the resistance mechanism of colistin-resistant A. baumannii and analyze its molecular epidemiology through WGS. Materials and methods The subjects in this study were patients who visited a university hospital between 2014 and 2018. Thirty colistin-resistant strains with high minimum inhibitory concentrations were selected from various patient samples, and WGS was performed. Comparative genomic analysis was performed for the 27 colistin-resistant A. baumannii strains using a colistin-susceptible strain as the reference genome. Results The WGS analysis found no mutation for lpxA, lpxC, lpx D, pmrA, pmrB, and mcr1, the genes known to be associated with colistin resistance. Fifty-seven coding sequences (CDS) showed differences; they included 13 CDS with known names and functions that contained 21 genes. From the whole-genome multi-locus sequence typing (wgMLST) and single nucleotide polymorphism (SNP) analyses, two major clusters were found for the colistin-resistant A. baumannii strains. However, no differences were observed by the time of detection for each cluster, the samples, the pattern of antibiotic resistance, or the patient characteristics. In the conventional MLST following the Oxford scheme, the typing result showed ST1809, ST451, ST191, ST1837, and ST369 in the global clone 2 (GC2), without any relation with the results of wgMLST and SNP analyses. Conclusion Based on the findings of the resistance gene analysis through WGS and comparative genomic analysis, the potential genes associated with colistin-resistance or CDS were examined. Furthermore, the analysis of molecular epidemiology through WGS regarding colistin-resistant A. baumannii may prove helpful in preventing infection by multidrug-resistant bacteria and controlling healthcare-associated infections.


Introduction
Acinetobacter baumannii is a widely distributed species in the hospital environment, and as an important causal pathogen of healthcare-associated infections, the species causes infection outbreaks in hospitals. If a patient with a severe underlying disease is affected, the reported mortality rate is substantially higher than that of other patients. In particular, the frequency of multidrug-resistant A. baumannii strains with resistance to three or more types of antibiotics has increased considerably lately [1].
The outbreak of A. baumannii is frequently reported worldwide; the reported incidence of multidrug-resistant (MDR) A. baumannii has increased in regions including Europe, North America and Asia. Treatment with broad-spectrum antibiotics often leads to MDR A. baumannii infection in patients with severe underlying diseases or reduced immunity and those who have undergone invasive surgery. This could lead to ventilator-related pneumonia, urinary tract infection, or sepsis, especially in the intensive care unit (ICU), and cause treatment failure, consequent complications, and increased mortality by preventing the selection of a suitable antibiotic for the treatment.
Colistin shows a relatively high rate of cure against serious infections by multidrug-resistant gram-negative bacteria resulting in a steady increase in its use. As a positively charged peptide, colistin acts by binding to the extracellular lipopolysaccharide (LPS) in gram-negative bacilli, displacing Mg 2+ and Ca 2+ crosslinked to it and destabilizing LPS molecules on the cell membrane. This increases the bacterial cell membrane permeability, leading to increased colistin uptake, resulting in the leakage of intracellular contents of the cell and cell death. Thus, a mutation of LPS can lead to resistance by preventing the action of colistin on the cell membrane [2]. Recently, a multidrug-resistant A. baumannii strain with colistin resistance has been reported [3]. The percentage of colistin-resistant A. baumannii is relatively low (1-3%) worldwide [4]; nonetheless, as colistin resistance indicates pan-drug resistance that shows resistance to all antibiotics in most cases, the lack of a suitable antibiotic for treatment has become a serious problem. It is thus crucial in the control of severe healthcare-associated infections that the resistance mechanism of the multidrug-resistant A. baumannii to colistin be elucidated to prevent the spread of the multidrug-resistant A. baumannii.
Next-generation sequencing (NGS) is a method of analysis where the genome is decoded by dividing it into myriads of fragments and combining each of their sequences. Compared to the conventional Sanger sequencing, NGS offers more rapid analysis, and as a single protocol, it can identify and type pathogens, analyze genetic relationships among different pathogens, predict pathogenic transmissions, and detect antibiotic resistance genes; this method may prove to be valuable for testing and studying clinical microbiology. However, only a few studies have applied NGS in studying the resistance mechanism and epidemiology of colistin-resistant A. baumannii [5].
The purpose of this study is to elucidate the resistance mechanism of colistin-resistant A. baumannii and analyze its molecular epidemiology through NGS-based whole-genome sequencing.

A. Subjects
The subjects in this study were patients admitted to Ewha Womans University Mokdong Hospital between 2014 and 2018, from whom various A. baumannii samples were isolated. The VITEK2 system (bioMerieux, Marcy l'Etoile, France) or VITEK MS (bioMerieux, Marcy l'Etoile, France) was used for species identification, and the strains that showed colistin resistance (minimum inhibitory concentration [MIC] � 4 μg/mL) as the result of antibiotic susceptibility using the VITEK2 system, were collected. The duplicate strains were excluded, and the percentage of colistin resistance in A. baumannii between 2014 and 2018 was 1.9% (84/4467). Among them, whole-genome sequencing was performed for the 30 colistin-resistant A. baumannii strains that showed high MIC (� 8μg/mL) on the broth microdilution test, and the strains were thus named S1-S30. The present study was conducted following the approval of the Institutional Review Board of Ewha Womans University Mokdong Hospital (EUMC 2019-09-013-003). As this study used residual samples, the risk to the subjects was extremely low, and it was a retrospective study, so the Institutional Review Board exempted the participants from filling out the consent form. In addition, all patient medical records and specimen information were anonymized, and individual identification or potential identification information of each patient was not included in this study.

Colistin resistance test.
The antibiotic susceptibility test was performed using the N225 CARD (bioMerieux, Marcy l'Etoile, France) of the VITEK2 system for the following antibiotics: amikacin, ampicillin/sulbactam, cefepime, cefotaxime, ceftazidime, ciprofloxacin, colistin, gentamicin, imipenem, meropenem, minocycline, piperacillin, ticarcillin/clavulanic acid, and trimethoprim/sulfamethoxazole. In addition, the MIC for colistin was determined using the GNX3F plate of the Sensititre system (Thermo Fisher Scientific, Waltham, MA, USA) via the commercial broth microdilution kit and the standard broth microdilution. The breakpoint of the Clinical Laboratory Standard Institute (CLSI) M100-ED29 guideline was used to determine colistin resistance and susceptibility; MIC �4 μg/mL was taken to indicate colistin resistance. Reference MICs for colistin were determined using manual broth microdilution (BMD), performed according to the CLSI M07-ED11 guideline. MIC panels were prepared with pure colistin sulfate powder (Sigma-Aldrich, St. Louis, MO, USA) with two-fold dilutions in the range of 0.125-64 μg/mL in cation-adjusted MH broth (Becton-Dickinson, Sparks, MD, USA) in polystyrene plates (Greiner, Frickenhausen, Germany). Results were read after incubation in an aerobic atmosphere at 35 ± 2˚C for 20-24 h. The results were interpreted based on the MICs according to the M100-ED29 guidelines.
B. Comparative genomic analysis. For comparative genomic analysis, the genome sequences with close association with A. baumannii were identified using the NCBI genome database. To compare the colistin-susceptible strains and the colistin-resistant strains, the genome of the colistin-susceptible A. baumannii ATCC 17978 was downloaded from the NCBI genome database. To understand the phylogenetic relationship among the strains to be analyzed, the OrthoANI was calculated and using the Orthologous ANI Tool (ChunLab, Seoul, Korea), the UPGMA dendrogram was produced. For the genome synteny analysis based on the nucleotide unit, the Blast Ring Image Generator (BRIG) was used with default variables [12,13]. To analyze the pan-genome and core-genome, the Comparative Genomics (CG) pipeline of ezbiocloud Apps (https://www.ezbiocloud.net/apps, ChunLab Inc.) was used. To define the pangenome orthologous groups (POGs), two methods were combined; i) reciprocal best hit using the uBLAST with E-value cut-off at 1×10 −6 [14]; ii) independent open reading frame (ORF) method with at least 70% gene coverage as the threshold of the nucleotide sequence [15]. The magnitude of the pan-genome and core-genome was expressed as the number of orthologous genes on the Venn diagram and using the ezbiocloud CG pipeline, the heat-map was produced to indicate the presence of a given gene. For 2D metabolite biosynthetic gene cluster and antibiotic resistance gene analyses, the antiSMASH [16] and Antibiotic Resistance Genes Database (ARDB) web servers [17] were used. If a given gene was presumed to be a virulence gene, a keyword search using the annotated gene or an amino acid sequence homology search for the predicted virulence protein of A. baumannii was performed.
C. Whole-genome multilocus sequence typing. Roary v3.12.0 pipeline was used to produce the core gene set from the 23 complete A. baumannii genomes [18]. The aforementioned core gene was used to produce the species-specific reference genome named A. baumannii. To conduct wgMLST, the core gene set discovered during the production of the species-specific standard sequence was used. For each core gene, the hidden Markov model (HMMs) was produced using the MAFFT v7.3.10 [19] and nhmmer v3.1b2 tool. The core gene was identified using the HMM-based search and was used to determine the number of unique alleles in each gene. For the minimum spanning tree based on wgMLST, coding was performed with the Kruskal algorithm for the in-house JavaScript (ChunLab, Seoul, Korea) [20].
In addition, the data of whole-genome sequencing were used, and based on the data of seven housekeeping genes (gltA, gyrB, gdhB, recA, cpn60, gpi, and rpoD), the sequence types were identified via conventional MLST using the Oxford scheme [21].
D. Single nucleotide polymorphism analysis. To calculate the genome sequence SNVs of the final dataset, the MUMmer v3.23 program was used, while applying the ATCC 19606 strain as the standard sequence [22]. The multiple sequence alignment was edited based on the detected SNPs, and the maximum likelihood phylogenetic tree was drawn using the RAxML v8.2.11 with 1,000 bootstrap re-samplings and GTRCAT as the nucleotide model [23].

C. Statistical analysis
To analyze the clinical characteristics of the clusters 1 and 2 from whole-genome multi-locus sequence typing (wgMLST), Student's t-test or Mann-Whitney U-test was performed for age and the admission duration prior to A. baumannii isolation; Chi-square test or Fisher's exact test was performed for gender, department of admission, underlying disease (diabetes, cancer, chronic renal disease, chronic liver disease, cardiovascular disease, gastrointestinal disease, and respiratory disease), history of admission within one year, history of antibiotic treatment within 28 days of A. baumannii isolation and the type of antibiotics used, ICU admission within 28 days, use of a ventilator within three months, central vein cannulation, tracheotomy, percutaneous abscess drainage, hemodialysis and other invasive procedures, co-infection with the isolation of a different strain from the same or different sample within three days of A. baumannii identification, and resistance to antibiotics. The significance was set at P<0.05. SPSS 23.0 software (IBM Corp, Armonk, NY, USA) was used for all statistical analyses.

Whole-genome sequencing
A. Microbial species identification and resistance gene analysis. The result of analyzing 30 A. baumannii strains using the TrueBac ID system, following the whole-genome sequencing, showed that 27 strains were A. baumannii. In comparison, three strains (S1, S11, and S12) were Acinetobacter colistiniresistens and thus excluded from subsequent analyses of resistance mechanism and epidemiology. The resistance gene search was performed using the CARD database, and no mutation was detected for the genes lpxA, lpxC, lpxD, pmrA, pmrB, and mcr1 that are known to be related to colistin resistance. The resistance gene search result for other antibiotics is presented in Table 1 and S1  Table 2. The WGS analysis found no mutation for lpxA, lpxC, lpx D, pmrA, pmrB, and mcr1, the genes known to be associated with colistin resistance. C. wgMLST and SNP analyses. The phylogenetic trees obtained through wgMLST and single nucleotide polymorphism (SNP) analyses are presented in Figs 1 and 2. The result of the wgMLST analysis showed two major clusters for the colistin-resistant A. baumannii strain ( Fig  1); likewise, the result of SNP analysis showed two major clusters (Fig 2). Neither wgMLST nor SNP analyses found a difference in the time of detection for each cluster, the hospital ward, the samples, and the antibiotic resistance pattern (Tables 3 and 4). The conventional MLST using the data of seven housekeeping genes (gltA, gyrB, gdhB, recA, cpn60, gpi, and rpoD) led to typing ST1809, ST451, ST191, ST1837, and ST369, with no relation to wgMLST or SNP analyses (Table 5).

Discussion
The most critical mechanism in the acquisition of carbapenem resistance by Acinetobacter is the production of carbapenemase that can hydrolyze carbapenem. Additionally, non-enzymatic mechanisms can contribute to carbapenem resistance by altering the outer membrane protein, modifying the affinity or expression of the penicillin-binding protein, or by overexpressing the efflux pump. The genes associated with carbapenem resistance identified in this study are adeI, adeJ, adeK, adeN, OXA-133, OXA-23, OXA-66, and TEM-1. These genes were

PLOS ONE
acquired by activating the LPS modification operon through the mutations of pmrA and pmrB, the two-component regulatory system, while the mutations of the genes essential in lipid A biosynthesis (lpxA, lpxC or lpxD) may also relate [25,26]. In a study based on wholegenome sequencing, 21 colistin-resistant A. baumannii strains were examined, and pmrAB mutation was detected in 71.4%, while LPS mutation was not observed [27]. In another study based on Sanger sequencing, 17 colistin-resistant A. baumannii strains were examined, and lpxA (29.4%), lpxC (29.4%), and lpxD (52.9%) were detected [28]. In a recent study on gut microflora, a novel colistin resistance mechanism was identified, which involved the expression of plasmid-encoded phosphoethanolamine transferase by the mobilized colistin resistance gene (mcr1). The resulting transfer of a phosphoethanolamine to lipid A on the cell membrane subsequently led to reduced colistin affinity, which led to resistance. The mcr1 gene is mainly

PLOS ONE
detected in Escherichia coli or Klebsiella pneumoniae isolated from a urine or blood sample of patients admitted to the hospital. The mcr1 gene in the gut microflora was identified in over 19 countries, and the reported level per country was within the range of 1.4-2% [29]. For A. baumannii, however, the association between microbial resistance and the mcr1 gene has not been reported. In this study, the resistance gene search using the CARD database showed that no mutation was detected for the following genes: lpxA, lpxC, lpxD, pmrA, pmrB, and mcr1.
In a comparative genomic analysis using NGS, the whole-genome sequences are compared with the reference genome, whereby a causal gene for a distinguished phenotype is detected, or the genomic characteristic shared by a specific group is identified. In this study, a comparative genomic analysis of the colistin-resistant strain was performed using the colistin-susceptible strain as the reference genome. The result showed that 57 CDS exhibited differences, among which 13 had a known name and function. These CDS contained 21 genes (cysH, udg, ureC, DNMT1, dcm, tmk, DTYMK, glxK, garK, CHO1, pssA, mdh, udg, mutY, TST, MPST, sseA, cobS, cobV, ENO, and eno) that may have an association with colistin resistance.
Among these genes, tmk and DTYMK are the genes engaged in dTMP kinase production, while being related to pyrimidine metabolism; glxK and garK govern the production of glycerate kinase and are related to glycerate activity; CHO1 and pssA govern the production of CDP-diacylglycerol-serine O-phosphatidyltransferase and participate in cell membrane and lipid biosynthesis; mdh governs the production of malate dehydrogenase and participates in the oxidation of malate to oxaloacetate; Udg governs the production of uracil-DNA glycosylase (UNG), an enzyme that hydrolyzes the N-glycosylic bond between the uracil base and the deoxyribose sugar in the DNA containing uracil, to create an apyrimidinic site; mutY governs the synthesis of adenine glycosylase and plays a role in correcting G-A mispairs and the errorprone DNA synthesis past GO due to an oxidative damage on guanine; TST, MTST, and sseA are the genes for thiosulfate sulfurtransferase that transfer a sulfur moiety from thiosulfate to thiophilic acceptor; cobS and cobV are the genes engaged in adenosylcobinamide-GDP

PLOS ONE
ribazoletransferase production and are related to the biosynthesis of cell membrane, inner membrane, cobalamin, and magnesium; ENO and eno are the genes engaged in the production of phosphopyruvate hydratase that catalyzes the reversible conversion of 2-phosphoglycerate to phosphoenolpyruvate, which is essential in carbohydrate decomposition during glycolysis. While these genes have not been reported to show an association with colistin resistance, there is a potential association. Notably, the genes for the biosynthesis of the cell membrane or lipids, such as CHO1 and pssA, or the genes for the synthesis of the cell membrane and magnesium, such as cobS and cobV, are likely to be associated with colistin resistance as the target of colistin is the bacterial cell membrane and the mechanism of colistin action is to destabilize LPS by removing magnesium to cause cell death. Further studies should thus be conducted regarding these genes and their association with colistin resistance. In addition, there remains a possibility that, among the CDS that differ based on comparative genomic analysis, a novel gene without a known name or function is associated with colistin resistance. This is meaningful in suggesting the possibility that CDS, which shows a difference between the 27 colistin-resistant A. baumannii group and the colistin-susceptible standard strain, may be involved in the resistance of colistin. Further study is needed to define the association these CDS and colistin resistance.
In this study, molecular epidemiology was analyzed using whole-genome sequencing. The result of wgMLST showed two major clusters for the colistin-resistant A. baumannii strain. However, the two clusters exhibited no difference in the time of detection for each cluster, the samples, the pattern of antibiotic resistance, or the clinical epidemiology. In addition, a phylogenetic analysis was performed using SNP, and as in wgMLST, two major clusters were found,

PLOS ONE
while several strains showed a deviation in classification. In the conventional MLST, the typing result showed ST1809, ST451, ST191, ST1837, ST369, and ST110, without any relation with the results of wgMLST and SNP analyses. The multidrug-resistant or carbapenem-resistant A. baumannii strains belong to the global clone 1 (GC1, formerly known as European clone 1 or international clone 1) and the global clone 2 (GC2, formerly European clone 2 or international clone 2) that are distributed worldwide. GC2 is more commonly distributed in Korea, and A. baumannii ST208, ST357, ST75, ST191, ST137, ST138, ST784, ST451, ST229, ST369, ST357, and ST552, have been reported [30][31][32][33][34][35]. This study examined ST191, ST451, and ST369, all of which belong to GC2. The method of typing using the whole genome is known to ensure excellent results in analyzing an outbreak compared to the conventional methods. Of note is the study by Fitzpatrick et al., where a case of hospital outbreak infection with an Acinetobacter species was investigated. For the 148 strains (including 116 A. baumannii strains) isolated from the blood of patients, the following methods were compared: pulsed-field gel electrophoresis (PFGE), repetitive extragenic palindromic-PCR, the conventional method of typing, and the method using whole-genome sequencing (SNP analysis). Using the conventional MLST as the standard method, the conformity was 39% for PFGE, 65% for Rep-PCR, and 100% for whole-genome sequencing [36].
Regarding the molecular epidemiology of A. baumannii, further studies should compare the accuracy of wgMLST and SNP analysis. In this study, the results of wgMLST and SNP typing showed differences, but these methods may be useful in the strain typing analysis of A. baumannii. The differences between wgMLST and SNP analysis may be due to these factors; wgMLST analyses consider recombination as a single mutational event, whereas the SNP calling analysis considers each SNP in the recombinant segment as an additional mutation. However, as the wgMLST scheme includes accessory genes, the presence/absence of these genes might impact the resulting phylogeny [37]. This study analyzed the molecular epidemiology of colistin-resistant A. baumannii for the first time in Korea.
Through whole-genome sequencing for resistance gene analysis and comparative genomic analysis, colistin resistance was potentially related to certain genes or CDS, but no mutation was detected concerning colistin resistance. This study defined the characteristics of colistinresistant A. baumannii isolated from Korea. In addition, although it is not known to be related to colistin resistance to date, candidate genes that may be associated with colistin resistance were presented, thereby serving as a cornerstone of research to elucidate the mechanism of colistin resistance. To conclude, the molecular epidemiology of colistin-resistant A. baumannii was analyzed using whole-genome sequencing, and the findings may prove useful in devising plans to prevent the spread of resistant bacteria and control the infection.