Surprising diversity of new plasmids in bacteria isolated from hemorrhoid patients

Strain acquisition and genome preparation

The samples were taken from various regions of patients with hemorrhoids in the anorectal department of Qilu Hospital of Shandong University (Qingdao). After the patient signs the paper-based informed consent form, seven different sites (healthy buttock skin, anal skin, the outer side of the hemorrhoid, the inner side of the hemorrhoid, healthy anal gland, feces, and the hemorrhoid tissue) on the patient were swabbed with a sterile cotton swab and then immersed in the preservation solution. A total of seven hemorrhoid patients were recruited for the isolation of bacterial strains. After being transported to the laboratory at a low temperature, the solution was cultured on tryptic soy broth (TSB) medium. Dominant bacterial strain was selected for enrichment culture, and the genome was extracted using Tiangen Bacteria DNA Kit (Tiangen Biochemical Technology Co., Ltd., Beijing, China) and stored at −20 °C. Agarose gel electrophoresis is used to check the DNA quality. Ultraviolet spectrophotometry is used to assess the purity of DNA samples, while the fluorescent method is employed to determine the concentration of DNA samples.

Library construction and nanopore sequencing

The DNA library was prepared according to the rapid library construction method. The barcode kit SQK-RBK114.96 was used. The third-generation sequencer used was Nanopore P2solo with R10.4.1 flow cell. Basecalling was performed using Dorado version 0.5.3 (https://github.com/nanoporetech/dorado/) to obtain raw data. NanoQC version 0.9.4 and NanoPlot version 1.42.0 were used for data quality control, and Chopper version 0.8.0 was used to filter the data and generate clean data (De Coster & Rademakers, 2023). The Flye assembler version 2.8.1-b1676 was utilized to assemble the genome (Kolmogorov et al., 2019). This sequence was then corrected twice using Medaka version 1.12.0 (https://github.com/nanoporetech/medaka) to acquire the final whole-genome sequence for further analysis. Quast version 5.0.2, CheckM version 1.0.1, and BUSCO version 5.2.2 were used to evaluate the assembly quality, contamination, and integrity of the genome, respectively (Parks et al., 2015; Mikheenko et al., 2018; Manni et al., 2021). Ineligible samples were resequenced until higher quality assembly results were obtained for all samples. Genome assembly was performed with the Prokaryotic Genome Annotation Pipeline (Li et al., 2021). Genome and plasmid circularity were determined with Flye.

Bioinformatics

Taxonomy of all isolated strains was identified to the species level using GTDB-Tk version 2.1.1 (Chaumeil et al., 2022). Antibiotic resistance genes (ARGs) were annotated with AMRFinder version 3.11.26 (Feldgarden et al., 2021). MLST analysis of the genomes was performed with MLST version 2.23.0 (https://github.com/tseemann/mlst). BURST was used to group sequence types (Jolley, Bray & Maiden, 2018). The PubMLST database (https://pubmlst.org/) was utilized for sequence type analysis (Jolley & Maiden, 2010). Plasmids were typed with plasmidfinder version 2.1.1 (Carattoli & Hasman, 2020). Phylogenetic analysis was performed with SNP-based genome phylogenetic method using Snippy version 4.6.0 (https://github.com/tseemann/snippy). A maximum likelihood tree was constructed using FastTree version 2.1 (Price, Dehal & Arkin, 2010). The phylogenetic trees were visualized with the iTOL tool (Letunic & Bork, 2021). Plasmid maps were generated with SnapGene version 7.2.1. Single Nucleotide Polymorphism (SNP) calculations were performed using Snippy version 4.6.0.

Genomic data

The genomic data involved in this work can be found at GenBank under accession number PRJNA1119654.

Ethics

The study was conducted in accordance with the Declaration of Helsinki and approved by the Medical Ethics Committee of Qilu Hospital of Shandong University (Qingdao) (approval number KYLL-2023045, approval date Aug 11, 2023).

Sample collection, bacteria isolation, and whole genome sequencing

Samples were taken from seven hemorrhoid patients who were subject to surgery in Qilu Hospital (Qingdao). Seven samples were taken from each patient, each from the healthy buttock skin, anal skin, the outer side of the hemorrhoid, the inner side of the hemorrhoid, healthy anal gland, feces, and the hemorrhoid tissue, totaling 49 samples. One bacterial strain was isolated from each sample by plating on TSB plates. Third-generation nanopore sequencing was performed with these strains, generating high quality whole genome sequences. Taxonomic classification was performed with whole genome sequences, leading to species-level classification.

An average of 3.3 contigs were assembled for each genome, including 44.4% (72 contigs) that are circular. Considering 32 of the strains contain plasmids, the low number of contigs for each genome and the high percentage of circular contigs suggest the good level of completeness of the genomes.

The 49 bacterial strains isolated include 34 Escherichia coli strains, six Klebsiella pneumoniae strains, four Enterococcus faecalis strains, two Klebsiella quasipneumoniae strains, two Enterobacter hormaechei strains, and one Morganella morganii strain (Table 1). The majority of the strains (89.8%) are Gram-negative bacteria belonging to a single order Enterobacterales. It is worth noting that all of the isolated strains are opportunistic pathogens commonly found with strong multidrug-resistant phenotypes. A maximum likelihood phylogenetic tree was constructed for all the E. coli strains isolated in this work (Fig. 1). The 34 strains were categorized into 13 sequence types, with ST10, ST1415, and ST4238 forming a single group (according to BURST analysis). In general, it can be observed that strains with the same sequence type are phylogenetically close, suggesting they are the same strain. The same strains are mostly from the same patient, except for two cases: patient Z22 (strains Z22-2, Z22-4, and Z22-5) and patient Z135 (all seven strains) share the same E. coli strain, whereas patient Z81 (Z81-7) and Z226 (Z226-7) share the same E. coli strain.

High resolution identification of plasmids from isolated strains

One particular strength of nanopore sequencing, or any 3^rd generation sequencing, is the ability to generate long reads, thus making recovery and assembly of plasmid sequences more feasible than 2^nd generation short-read high-throughput sequencing methods. The significantly improved single-strand accuracy of Q20+ chemistry, in combination with the R10.4.1 flow cell enables nanopore sequencing to generate long, high-quality (Q value > 18) reads (Zhang et al., 2023), which permits highly accurate discovery of plasmids. In this work, a total of 31 unique plasmids with an average size of 69.2 kb were found, of which 22 (71.0%) were assembled into circular form, indicating the completeness of plasmid sequences (Table 2). The high level of complete plasmid sequences identified suggests the power of Nanopore long-read high-throughput sequencing.

This work identifies 15 new unique plasmids, nearly half of the total of 31 unique plasmids found, suggesting that there are numerous new plasmids in nature that have not been adequately investigated.

The structures of the 15 newly identified plasmids were subjected to more detailed analysis (Fig. 2). One initial strong impression is the high number of genes coding for integrase/transpose/recombinase. Take pZ14232, for example, out of the 56 protein-coding genes this plasmid has, 21 (37.5%) genes are integrase/transposase/recombinase-coding genes. This is true for all the large plasmids found in this work. Considering all new plasmids are from Enterobacteriales, this finding suggests the strong potential for plasmids of Enterobacteriales to load/unload cargo genes for modification of functions. This also suggests the plasticity of genomes of these opportunistic pathogens to varying environments during co-evolution with humans.

The sizes of the newly identified plasmids range from 0.9 to 172.2 kb. All plasmids, except for pZ2221, carry the rep gene essential for plasmid replication, supporting their annotation as plasmids. pZ2221 does not encode a Rep protein homologous to known Rep proteins. Considering the coverage of pZ2221 is only nine, we cannot rule out the possibility that it is part of a larger plasmid, and the circularity of this plasmid is simply an assembly error.

Four of the remaining 14 plasmids (pZ2253, pZ17322, pZ22613, and pZ22663) that carry at least one rep gene are relatively small in size and simple in structure. It is interesting to note that pZ2253, pZ22613, and pZ22663 carry only the rep gene and mobilization-related genes, suggesting they are small and mobile plasmids. In particular, their predicted copy numbers are respectively 45, six, and 13, respectively, indicating they have high copy numbers. This may imply they may play important unknown physiological roles in Enterobacteriales.

Three of the plasmids (pZ8171, pZ14261, and pZ21751) are conjugative plasmids that harbor the conjugative operons belonging to type IV secretion systems. Two of them (pZ8171 and Z21751) are multidrug-resistant plasmids, making them novel dangerous carriers of antimicrobial resistance that can disseminate within bacterial communities. Two additional plasmids, pZ14231 and pZ14232, are multidrug-resistant plasmids that do not encode type IV secretion systems but have mobility-related genes. They may also be involved in antimicrobial resistance dissemination.

These newly identified plasmids also carry other physiological functional gene clusters. Besides the widespread toxin-antitoxin systems that can be found on almost every large plasmid, pZ21751 has a colicin biosynthesis and immunity gene cluster making it a clear virulence factor. Additionally, pZ14211, one of the two plasmids in multidrug-resistant E. coli Z142-1, encodes a plasmid-borne type II secretion system. This protein complex is capable of exporting folded proteins, including toxins (Shaliutina-Loginova, Francetic & Doležal, 2023). The presence of these gene clusters suggests their role in virulence, host-microbe interactions, and microbe-microbe interactions. Metal resistance gene clusters were found on pZ14251 (mercury resistance cluster) and pZ14261 (copper/silver resistance cluster). Interestingly, pZ14231 carries a nitrate respiration gene cluster, whereas pZ14231 and pZ14251 carry cellulose biosynthesis gene clusters, granting the host bacteria additional feats.

Antibiotic resistance determinants

The presence of ARGs in each strain was analyzed with AMRFinder (Table 3). Twenty-five strains were found to contain at least one ARG. When not considering strains carrying ARGs that are part of the core genomes, acquired resistance mechanisms can be found in at least nine strains. Taking into account that the strains were isolated without supplementation of antibiotics in the media, the “lower-than-expected” acquired resistance prevalence (0% for Klebsiella strains, 29.4% for E. coli) suggests that the antibiotic resistance problem has not developed beyond rescue. It needs to be noted that, however, out of the seven multidrug-resistant plasmids identified, four are new plasmids. This reflects how little we know about multidrug-resistance plasmids, in consistence with the finding that only about half of the unique plasmids found in this work are previously identified.

Class 1 integrons can be found in four strains (Fig. 3). Two integrons were found on the chromosome of M. morganii Z22-3, whereas one integron was found in pZ21751 in E. coli Z217-5, pZ21761 in E. coli Z217-6, and pZ21771 in E. coli Z217-7. The structures of these integrons have been previously observed. It is worth noting that the presence of integrons in M. morganii strains has not been widely reported, which could have been a rare event.

Plasmid transmission and modification events

By comparing the strains found in this work, it is evident that plasmid transmission and modification events are quite common. Through SNP-based genomic phylogenetic analysis, we identified that E. coli Z22-2, Z22-4, Z22-5 strains are the same strain. They all carry a 139.1-kb multidrug resistance plasmid (pZ2222 in E. coli E22-2). However, E. coli Z22-2 carries three additional plasmids: pZ2221, pZ2223, and pZ2224, while E. coli Z22-5 carries an additional plasmid, pZ2253. Notably, pZ2224 is an IncFII-type multidrug resistance plasmid that can grant bacteria with additional functions. The variation in plasmid composition within the same strain from the same patient indicates the rapid spread of plasmids, leading to genetic adaptability.

Similarly, in patient Z22, E. coli strains Z22-1, Z22-6, Z22-7 are the same strain. However, strains Z22-6 and Z22-7 carry a plasmid that strain Z22-1 does not possess. In patient Z173, K. quasipneumoniae strains Z173-2 and Z173-4 are the same strain, but strain Z173-2 carries a pZ17322 plasmid that is not found in strain Z173-4; in patient Z226, E. coli strains Z226-5 and Z226-6 are the same strain, but strain Z226-6 carries two additional plasmids, pZ22662 and pZ22663, compared to strain Z226-5.

Plasmid modification events were also observed in this study. In patient Z217, E. coli strains Z217-5, Z217-6, and Z217-7 are the same strain, and the plasmids they harbor, pZ21751, pZ21761, and pZ21771, are the same plasmid. However, pZ21751 lost an ermB gene compared to pZ21761 and pZ21771, which could potentially make E. coli Z217-5 sensitive to erythromycin.

It needs to be noted here that very few SNPs were found between the chromosomes of the same strains. E. coli Z22-1 and Z22-6 have 28 and 0 SNPs in comparison with strain Z22-7. E. coli Z217-5 and Z217-6 have 175 and 141 SNPs in comparison with strain Z217-7. This is negligible (mutation rate < 0.01%) compared to genomic changes incurred with plasmid transmission and modification, which are frequently observed in this work. Therefore, we believe changes in plasmids are the primary source of genetic plasticity in the case of Enterobacteriales that are primarily studied in this work. Taking into account that a large portion of plasmids found in this work are new, we believe genetic plasticity has been underestimated.