In our previous study, we tested the susceptibility of 168 non-duplicate CR-KP clinical strains collected from different patients during June 2014–December 2018 in four tertiary hospitals in Beijing, China (Ni et al., 2021). Among them, 98 CR-KP strains from three tertiary hospitals were used in this study. Most strains (53/98, 54.1%) were collected from the lower respiratory tract specimens, followed by the blood samples (31/98, 31.6%). If multiple strains were isolated from the same patient, only the first strain was chosen. Of the 98 strains, 88 (89.8%) produced the KPC-2 type carbapenemase, and 10 produced the NDM-1 type carbapenemase. A total of eight K. pneumoniae ST types were classified. ST11 was the most prevalent ST (69/98, 70.4%), followed by ST37 (13/98, 12.3%), ST395 (8/98, 8.2%), ST15 (3/98, 3.1%), ST48 (2/98, 2.0%), ST17 (1/98, 1.0%), ST23 (1/98, 1.0%), and ST290 (1/98, 1.0%). All strains were identified by the VITEK 2 Compact System (bioMérieux, Marcy-l’Étoile, France). Escherichia coli ATCC 25922 was used as a quality control strain in the MIC determination. The MICs for colistin and other antibiotics were determined by standard broth microdilution methods according to the recommendations of Clinical and Laboratory Standards Institute [CLSI], 2020. The breakpoint for colistin was based on the criteria established by the CLSI: The MIC ≤ 2 mg/L is judged as susceptible, and >2 mg/L as resistant.

For screening the heteroresistance in CR-KP, population analysis profiling (PAP) assays were performed as previously described (Band et al., 2019). Strains were tested in three times to ensure the accuracy of the entire experiment. In brief, 0.5 McFarland overnight culture bacterial suspensions were serially diluted to bacterial suspensions from 10⁸ to 10⁶ CFU/ml, and 100 μl aliquots with different cell densities were plated onto Mueller–Hinton agar (MHA) plates with increasing concentrations of colistin (from the MIC to 128 mg/L). After incubation at 37°C for 48 h, strains were classified as heteroresistant if the number of colonies growing at 2.0- or 4.0-fold that of the breakpoint concentration was at least 0.0001% (1/10⁶) of those growing on antibiotic-free plates.

In addition, a modified PAP assay was used for examining the efficacy of colistin combined with tetracyclines (tigecycline or minocycline) or aminoglycosides (amikacin or gentamicin) on the suppression of heteroresistance. Bacterial suspensions (100 μl aliquots with cell densities of 10⁸ CFU/ml) were plated onto MHA plates with increasing concentrations of colistin (from the 1 × MIC to 128 mg/L) and a fixed concentration (1 × MIC) of the other drug. Following 48 h incubation at 37°C, the number of colonies was counted.

The time-kill assays for colistin monotherapy and combination therapy were performed as previously described (Ni et al., 2021). The concentration of colistin and other antibiotics used in time-kill assays was 1 × MIC of each tested strain that was obtained from the susceptibility testing. Tubes containing 10 ml of freshly prepared cation-adjusted Mueller–Hinton broth (CAMHB) (Difco, Detroit, MI, United States) supplemented colistin alone or combined with other antibiotics were inoculated with Strains to a final density of approximately 5 × 10⁶ CFU/ml, and incubated with shaking at 37°C. Then, 100 μl aliquots of the culture sample were removed at 0, 3, 6, 9, 12, and 24 h, serially diluted and plated on antibiotic-free MHA plates for determining the total number of cells, or plated on agar plates containing different concentrations of colistin for determining the number of resistant subpopulations.

We selected one clinical strain K65 and its four resistant clones (named C1–C4) obtained from the PAP. After overnight culture, a 1:100 dilution of saturated culture was mixed into fresh Mueller–Hinton broth and shaken at 37°C for 3 h. After that, cultures were diluted in CAMHB to reach turbidity equal to 0.5 McFarland. Next, a 1:100 dilution of each 0.5 McFarland culture was mixed into CAMHB and shaken at 37°C continuously. The absorbance of the bacterial medium at 600 nm was measured every hour. Each strain was cultured in triplicate and the means absorbance was calculated.

Total bacterial RNA was extracted using the RNAprep pure cell/bacteria kit (Tiangen, Beijing, China) according to manufacturer instructions. cDNA was synthesized using the TIANScript cDNA kit (Tiangen, Beijing, China), and the expression of phoP and pmrD was quantified using a 7500 real-time PCR system (Applied Biosystems, Foster City, CA, United States). Housekeeping gene rrsE was used as an internal reference. Relative quantification of target genes (phoP and pmrD) was performed with the 2^–ΔΔCt method, normalized to rrsE and expressed as a relative fold change to rrsE expression.

We performed the WGS and transcriptome analyses for the parent strain K65 and its four heteroresistant clones mentioned above. For next-generation sequencing, the sequencing library was constructed by adding index codes to attribute sequences for each sample, with an average insert size of 474 bp. Paired-end reads of 150 bp were generated using the Illumina HiSeq platform at Shanghai Majorbio Bio-Pharm Technology Co., Ltd. Raw reads were filtered to remove low-quality reads (<Q20). Clean data then were assembled de novo using SOAPdenovo2 (Luo et al., 2012)¹. Coding sequences of the genome were predicted using GeneMarkS, Glimmer², and Prodigal software. The phylogenetic tree of K65 and four heteroresistant clones was constructed by comparsion with 31 house-keeping genes (dnaG, frr, infC, nusA, pgk, pyrG, rplA, rplB, rplC, rplD, rplE, rplF, rplK, rplL, rplM, rplN, rplP, rplS, rplT, rpmA, rpoB, rpsB, rpsC, rpsE, rpsI, rpsJ, rpsK, rpsM, rpsS, smpB, and tsf). Multilocus Sequencing Typing (MLST) analysis with the 7 housekeeping genes gapA, infB, mdh, pgi, phoE, rpoB, and tonB was performed to determine sequence types (STs) using the PubMLST database. The total RNA of each sample was extracted using the RNAprep Pure Bacteria Kit (Tiangen, China) following the manufacturer’s protocol. Ribosomal RNA (rRNA) was removed using the Ribo-Zero kit (Epicentre). RNA-Seq was undertaken using the Illumina HiSeq platform mentioned above. Several databases were used to predict gene functions, such as NR, Swiss-Prot, Pfam, EggNOG, Gene Ontology, and Kyoto Encyclopedia of Genes and Genomes. The RSEM software³ was used for the quantitative analysis of gene expression levels. DEGseq software was used to analyze the differentially expressed genes between groups by screening the threshold value of | log2FC| ≥ 1 and p-adjust < 0.05.

Among the 98 CR-KP strains, 2 (2/98, 2.0%) and 69 (69/96, 71.9%) showed colistin resistance (MIC > 2 mg/L) and colistin heteroresistance, respectively. The PAP results for five representative strains are shown in Figure 1A. For testing the stability of highest degrees of resistance in heteroresistant subpopulations, two heteroresistant colonies of each strain growing on MHA plates containing the highest drug concentrations were randomly selected from the consecutive five-generation cultures on antibiotic-free plates. These heteroresistant subpopulations showed stable high-MIC resistance to colistin. Their MICs were at least 8-fold higher than those of the parent strains, with MIC₅₀ and MIC₉₀ both of ≥16 mg/L (Figure 1B). The distribution of colistin heteroresistance at different MIC values of colistin and nine other antimicrobial agents is shown in Figure 2. The proportion of heteroresistant strains is uniformly distributed, suggesting no distinct associations between colistin heteroresistance and antibiotic resistance.

We then examined the time-kill kinetics of colistin at 1 × MIC against all the strains. For heteroresistant strains, the number of cells decreased initially, but regrowth occurred quickly; for non-heteroresistant strains, the number of cells quickly decreased in the first 6 h without regrowth. The time-kill curves of the four representative strains are shown in Figure 3A. We further monitored resistance in the regrowth population. Four strains (two with colistin heteroresistance and two without heteroresistance) were randomly selected. Then, 100 μl aliquots of the culture sample were extracted at 12 and 24 h in the time-kill assays, serially diluted, and plated on MHA plates containing 2 or 4 mg/L colistin to determine the proportion of resistant subpopulations. As shown in Figure 3B, compared with the heteroresistant strains (K39 and K65), no resistant subpopulations of K28 and K46 recovered at 12 and 24 h (cell density 10⁹–10¹⁰ CFU/ml) in the control group (without colistin exposure). This result further confirmed that these strains were “purely susceptible” without heteroresistance. In the heteroresistant strains (K39 and K65), exposure to colistin significantly increased the frequency of resistance at 12 and 24 h.

The whole genetic profile of four randomly selected resistant clones for K65 (C1–C4) were screened and compared with the parent strain. The parent strain and clones belonged to ST395, and the phylogenetic tree of 31 house-keeping genes showed high homology levels (Supplementary Figure 1A). To determine the difference in fitness costs between the susceptible strains and resistant subpopulations, we determined the growth curves of the parent strain (K65) and its four resistant clones (C1–C4). As shown in Supplementary Figure 1B, the growth rate had no significant differences between the resistant clones and parent strain, indicating that the resistant phenotype did not affect bacterial fitness.

Some specific genes were found in the resistant clones but not in K65. Further analysis showed these clones shared several specific genes involved in gene transfer and recombination: traV encoding lipoprotein in the type IV conjugative transfer system, rdgC encoding recombination-associated protein RdgC, insB encoding IS1 family transposase and umuC encoding Y-family DNA polymerase (Supplementary Table 1). The two-component regulatory system crrAB was not identified, and no mutations in phoP, phoQ, pmrB, pmrD, pmrA, and mgrB were identified in the resistant clones.

We then compared the transcriptomic profiles of K65 and four resistant clones (C1–C4). Figure 4A showed the upregulated and downregulated genes in the resistant clones; Figure 4B showed the number of differentially expressed genes in the resistant clones. GO analysis revealed that the differentially expressed genes were mainly involved in biological processes, such as LPS biosynthesis and metabolism (Figure 4C). Further analysis revealed that resistant clones shared a similar pattern characterized by increased phoPQ and pmrD expressions levels relative to the parent strain (Figure 4D). The qRT-PCR of more strains and their resistant clones (n = 20) verified the findings (Figure 4E). Figure 4F shows the model of resistance mechanisms in the colistin-resistant subpopulations.

We further examined the efficacy of colistin combined with tetracyclines (tigecycline or minocycline) or aminoglycosides (amikacin or gentamicin) in suppressing heteroresistance. The modified PAP analysis showed that after another drug with 1 × MIC was added, no colonies recovered from agar plates containing 2 mg/L colistin in most strains (93/96). The time-kill curves verified the results of modified PAP for combination therapy. Four strains were randomly selected for the assays. As shown in Figure 5, for the heteroresistant strains (K39 and K65), combination therapies could constantly inhibit the growth of CR-KP, and no regrowth occurred within 24 h of treatment. For the non-heteroresistant strains (K28 and K46), colistin monotherapy could completely kill all the cells within 6 h, and the addition of the other drugs did not enhance colistin activity.