Pathogenic strains of Shewanella putrefaciens contain plasmids that are absent in the probiotic strain Pdp11

Bacterial strains and growth conditions

Shewanella putrefaciens strains SH4, SH6, SH9, SH12 and SH16 isolated from diseased eels were kindly provided by Dra. Esteve C. from the University of Valencia (Valencia, Spain) (Esteve, Merchán & Alcaide, 2017). Two saprophytic isolates (SdM1 and SdM2) of S. putrefaciens were obtained from environmental sources (Seoane et al., 2019) and one probiotic strain, S. putrefaciens Pdp11 CECT 7627, was isolated from the skin of healthy Sparus aurata L. (Díaz-Rosales et al. 2009).

All S. putrefaciens strains were grown on tryptic soy agar plates (Oxoid Ltd., Basingstoke, UK) supplemented with NaCl (1.5%) (TSAs) for 24 h at 23 °C. Then, one or two colonies of each strain were picked and grown in 10 mL tubes of tryptic soy broth (Oxoid Ltd., Basingstoke, UK) added with NaCl (1.5%) (TSBs) for 24 h at 23 °C on shaking at 80 rpm (ELMI DOS-20M Digital Orbital Shaker, USA). Given the possibility that, if present, plasmids could be integrated into the bacterial chromosome, strains were cultured under different growth conditions (temperature, incubation time, growth medium and freeze-thaw) to favour their excision. For this, one or two colonies of each strain were picked and cultured in 10 mL tubes of TSBs and minimal (M9) media, and incubated at 23 °C or 4 °C for 24 and 48 h on shaking at 80 rpm. Cultures inoculated in parallel with TSBs medium containing glycerol (20%) were subjected to a freeze-thaw cycle (Pesaro et al., 2003) for 24 h at −80 °C prior to incubation.

Plasmid DNA isolation

One or two colonies of the S. putrefaciens strains grown under different growth conditions were picked from the pure culture and grown in 10 mL of TSBs (Oxoid Ltd., Basingstoke, UK) for 24 h at 23 °C under agitaton at 80 rpm (ELMI DOS-20M Digital Orbital Shaker). As a positive control, Escherichia coli V157, a strain harboring seven plasmids (Macrina et al., 1978), was grown on Luria-Bertani agar plates (LB) (Becton, Dickson and Company, Le Pont de Claix, France) for 24 h at 37 °C. The cultures were centrifuged at 8,000×g for 5 min and the pellet was used for plasmid DNA (pDNA) isolation using the GeneJet Plasmid Miniprep kit (Thermo Fisher, Waltham, MA, USA) following the manufacturer’s protocol. Plasmid DNA integrity was checked by agarose gel (0.8%, w/v) electrophoresis in the presence of RedSafe nucleic acid staining solution (InTRON Biotechnology, Seongnam-Si, South Korea). The 0.2–10 kb Hyperladder molecular weight marker (Bioline, Taunton, MA, USA) was used to check plasmid size. The pDNA was stored at −20 °C for further processing.

Plasmid DNA amplification and sequencing

Rolling circle amplification (RCA) was performed using the TempliPhi 100 amplification kit (GE Healthcare, Chicago, IL, USA) for each isolated plasmid following the manufacturer’s instructions. As a positive control for the reaction, pUC19 was used. DNA from each RCA-amplified plasmid was digested with EcoRI and EcoRV (Takara, San Jose, CA, USA), separately. They were then ligated employing T4 DNA ligase (Thermo Fisher Scientific, Waltham, MA, USA) to pBluescript SK (+) (pBSK S+) (Addgene, London, UK) previously treated with shrimp phosphatase alkaline (New England Biolabs, Ipswich, MA, USA). Plasmid DNA was transformed into CaCl₂-treated E. coli DH5α competent cells, as described by Sambrook, Fritsch & Maniatis (1989). Recombinant bacteria were selected using plates containing ampicillin (100 μg/ml) and X-gal (20 μg/ml), and colonies resuspended in 50 µl of DEPC-Water (Sigma, St. Louis, MI, USA) and tested by PCR amplification. Reaction mixture contained 2 µl of bacterial suspension, 62,5 U of Taq Accustart II Trough Mix (Boimerieux, Marcy-l’Étoile, France), 20 pmol of M13R primer 5′-CAGGAAACAGCTATGAC-3′ and 20 pmol of M13F primer 5′-TGTAAAACGACGGCCAGT-3′ in a final volume of 10 µL. The PCR program was: (1) 94 °C 3 min, (2) 94 °C 30 s, (3) 50 °C 30 s, (4) 72 °C 1 min/kb. Steps (2) to (4) were repeated during 35 cycles. Amplicons were checked by agarose electrophoresis gel as above. PCR products were sent to Macrogen (Madrid, Spain) for Sanger sequencing.

Plasmid DNA sequence assembly

Plasmids were assembled de novo using a workflow integrating Sanger and Illumina reads (Chevreux, 2005; Bankevich et al., 2012). Plasmid-specific Illumina sequences used here were obtained by Seoane et al. (2019) using a sequencing library with the Nextera protocol for 2 × 300 bp PE (raw reads available in BioProject PRJNA510237).

First, the Sanger plasmid sequences were pre-processed with BBDuk (Bushnell, 2015; Kechin et al., 2017) a tool developed for quality filtering and adapter trimming using k-mers matching, where k-mers are substrings of length “k” contained in the original nucleotide sequence. To remove the sequence of the pBSK S+ vector, different values of the two clean-up parameters k and hdist were applied, where “k” is the length of the k-mers and “hdist” is the “hamming distance” between two k-mers. We used a k = 20 and hdist = 1 for pSH12 and a k = 22 and hdist = 2 for pSH4 (Fig. 1A). Two sets of vector-free Sanger reads were generated when quality trimming with BBDuk was applied using “trimq” parameter. This parameter is the minimum average quality required in a sequence window, if the average is below the threshold, the nucleotides are trimmed from the read. For each plasmid, a low (10) and a high (17) stringency quality values were applied. The low stringency threshold generate longer reads that are used to capture Illumina reads from the BioProject PRJNA510237, denoted as Capture read set (Fig. 1B). With the high strigency threshold, a high quality read set for the final assembly is generated, denoted as Assembly read set (Fig. 1C).

Next, Illumina plasmid reads were captured with the sanger Capture read set using bowtie 2 (Langmead & Salzberg, 2012; Langmead, 2017). These reads were identified and subsequently extracted with Samtools (Li et al., 2009) by selecting pairs of Illumina reads in which at least one member had aligned to the Capture read set (Fig. 1D). Then, Sanger-Illumina hybrid assembly was performed with two assemblers, Mira (Chevreux, 2005) and Spades (Bankevich et al., 2012) using default parameter values. The aforementioned Illumina reads and the Assembly read set, obtaining two possible and complete pDNA assemblies per strain (Fig. 1E). The two sequences per strain were recircularized by the MARS method (Ayad & Pissis, 2017) (Fig. 1F) to make possible the next step. The plasmid sequences of each strain, they were aligned by the Clustal method (Sievers & Higgins, 2014) of the Seaview program (Galtier, Gouy & Gautier, 1996) (Fig. 1G), leading to a unique plasmid consensus sequence per strain.

Assembly validation

Alignment of the sequence set of SH4 and SH12 plasmids using Mira and Spades programs revealed the presence of unreliable regions in both plasmids. Validation was carried out by PCR and primers were designed for each non-validated plasmid sequence using Primer3 (Untergasser et al., 2012). The composition of the reaction mixture was done as above, except that 10 pmol of the primers pSH4_R 5′<CGGATTGAATGGCTGGCTGGACTG>3 (reverse) and pSH4_F 5′<ATACCAAACGCCCACAGA>3′ (direct) were used to validate the assembly of pSH4, and of the primers pSH12_R 5′<GGCTCCACCCTTACCCAAAAAA>3′ (reverse) and pSH12_F <5′GCGAGCCCCTCCATGATTTT>3′ (direct) to validate the assembly of pSH12. The PCR program was: (1) 94 °C 3 min, (2) 94 °C 30 s, (3) 66 °C 30 s, (4) 72 °C 1 min/kb, with 35 repeat cycles of steps (2) through (4). PCR products were checked by agarose gel electrophoresis using the 100 bp molecular weight marker Hyperladder (Bioline, Taunton, MA, USA). Finally, the plasmid assemblies of strains SH4 and SH12 were mapped back to the Sanger assembly read set sequence, the Illumina sequences captured from Sanger capture read set sequence and the PCR product obtained from validation.

Plasmid annotation and functional characterization

Once complete consensus sequences were obtained for each plasmid, sequences were searched for on the PLSDB web server (Galata et al., 2019), a resource containing numerous plasmid records collected from the NCBI nucleotide database. Next, gene prediction was performed for annotation of bacterial operons and open reading frames (ORFs) using Softberry Software (http://www.softberry.com/berry.phtml) including the fgenesB tool (Solovyev et al., 2011) with “Bacterial generic” as the closest organism; the Prokaryotic GeneMark. Hmm version 2 program (Besemer, Lomsadze & Borodovsky, 2001) with “Shewanella putrefaciens CN_32” as the selected species. Consensus sequences were also analysed with ORF finder at NCBI website (https://www.ncbi.nlm.nih.gov/orffinder/) with the search parameters: 150 nucleotides of minimal ORF length, genetic code “Bacterial, Archaeal and Plant Plastid”, and “ATG and alternative initiation codons” as ORF start codon to be used. For each of the identified ORFs, their amino acid sequence was obtained and queried in blastp (NCBI, Bethesda, Maryland, USA), to obtain clues about the ORF functions (Altschul et al., 1990; Gish & States, 1993). PHYRE 2 version 2.0, was also used for subsequent protein prediction and modelling (Kelley et al., 2015). Finally, Dfast (Tanizawa, Fujisawa & Nakamura, 2018) was used for ORF protein identification and functional annotation. Plasmid maps were drawn with pDRAW32 version 1.1.146 (AcaClone software, Greeley, CO, USA) (https://www.acaclone.com/).

Plasmid to genome comparison

To detect the degree of similarity between the genomes of the SH4 and SH12 strains and their respective plasmids, a comparison between them was carried out using Sibelia (Minkin et al., 2013). The complete putative pDNA sequences and their corresponding ORFs were also searched by blastn against the probiotic genome and all non-probiotic strains (NPS) genomes previously described by Seoane et al. (2019).

Phylogenetic study

The repB protein sequences of pSH12 and pSH4 were searched in Blast (NCBI, Bethesda, Maryland, USA). The repB protein sequences of plasmids from different bacterial species were downloaded in FASTA format in a single export file, and alignment was performed by the Clustal method (Sievers & Higgins, 2014). A phylogenetic tree was constructed by PhyML employing the Maximum Likelihood method (Edwards & Cavalli-Sforza, 1964) using Seaview software (Galtier, Gouy & Gautier, 1996).

Ethical statement

To conduct the research, we used bacterial strain Shewanella putrefaciens (strains SH4, SH6, SH9, SH12, SH16, SdM1, SdM2 and Pdp11), which is not considered a human or animal sample. Therefore, no ethics approval was needed, and no informed consent was required.

Identification of plasmids

S. putrefaciens strains SH4, SH6, SH12, SH16, SH19, Pdp11, SdM1 and SdM2 were analysed to confirm the presence or absence of plasmids under different growth conditions (Table S1). After growth under optimal conditions, only SH4 and SH12 strains were found to contain one plasmid around 3 kb in size (Fig. 2). The same results were obtained when other growth conditions were used (Table S1). E. coli V157 strain, used as a positive control, showed all seven plasmids previously described in the literature (Macrina et al., 1978).

Plasmid assembly and sequence validation

The consensus sequences of pSH4 and pSH12 obtained using Mira-Spade workflow from Sanger and Illumina reads showed a size of 3,003 bp for pSH4 and 2,990 bp for pSH12. However, unvalidated regions were detected in both plasmids (Fig. 3) because the Capture read set sequence did not sufficiently capture the Illumina sequences in pSH4 (Fig. 3A) and pSH12 (Fig. 3B), rendering these regions unreliable. Moreover, the Assembly read set sequence did not support this region in the alignment of pSH4 (Fig. 3A), unlike in pSH12 (Fig. 3B). From the PCR product obtained, a new alignment using Sanger, Illumina and PCR product sequences was carried out for each plasmid and mapped. After this, the unvalidated region of pSH4 was supported by a total of 22 Illumina and one Sanger sequences as was the PCR amplicon, while that of pSH12 was invalid. Thus, fully assembled circular sequences showed 3,003 bp in pSH4 but 2,960 bp in pSH12, 43 pb less due to the absence of that unreliable region (Fig. 3B). The GC content of both plasmids was similar, 42.69% (pSH4) and 43.04% (pSH12). The raw sequence information from Seaview program and agarose gel of PCR amplicons are showed in Fig. S1.

ORF analysis of pSH4 and pSH12

No hits sequences of the plasmids were available on the PLSDB web server (Galata et al., 2019), so a comprehensive ORF analysis was carried out. Analysis using the ORF finder revealed 10 possible ORFs in pSH4 and nine in pSH12. However, after removing nested predictions (i.e., overlapping ORFs) a total of seven ORFs remained in pSH4 and six in pSH12 (Table S2). The gene prediction tools, fgenesB, GeneMark and ORF finder, confirmed that only four and three ORFs in pSH4 and pSH12, respectively, were common and highly similar in both plasmids (Table 1).

Both plasmid sequences resulted in the presence of the same ORFs, except for the presence of an extra 43 bp region in pSH4 that introduced a stop codon that split in two the equivalent of ORF5 in pSH12 giving rise to ORF5 and ORF6 in pSH4. The localization and orientation of all ORFs of pSH12 and pSH4 are shown in the circular map in Fig. 4. Since unique EcoRI and EcoRV cleavage sites were found in both plasmids, they were used for cloning and sequencing.

Prediction of functional proteins

Based on blastp and Phyre 2, ORF5 and ORF6 of pSH4, and ORF5 of pSH12 presented a high identity to a plasmid replication initiator protein of the repB family of Shewanella algae (Table 2). The proteins encoded by these ORFs were predicted to show structural analogy with a crystal structure of the plasmid initiator protein of the repB family, with a nucleic acid as ligand and functions such as DNA replication and plasmid copy control (Table 2). Similarly, the proteins encoded by ORF7 of pSH4 and ORF6 of pSH12 showed analogy with the Type II toxin-antitoxin system PemK/MazF family toxin from Avibacterium paragallinarum, presenting DNA as ligand and functions of cell growth inhibitor by cleaving mRNA with high sequence-specificity as part of the bacterial stress response, and as toxic component for plasmid maintenance (Table 2). No reliable predictions were obtained for the remaining ORFs that showed similarities to hypothetical or unknown proteins.

Genomic and plasmid comparison between strains

Comparison of the genomes of SH4 and SH12 and their respective plasmids yielded high similarity at the plasmid level with 94.24% coverage (Fig. S2), differing only by an additional 43 pb region in pSH4 as explained above. At the same time, genome-wide comparison of SH4 and SH12 strains also showed high similarity composed of 46 and 58 scaffolds, respectively with 99.42% coverage between each other (Fig. S2). These percentage values support the idea that SH4 and SH12 are different strains belonging to the same species.

In addition, a blastn alignment was carried out to check for a possible insertion of these plasmids or any of their ORFs in the complete genome of the other strains under study. A total absence of plasmids was found in the rest of NPS genomes, as wells as in the genome of the probiotic S. putrefaciens Pdp11.

Phylogenetic analysis of the repB protein

The repB protein sequences found in pSH12 (ORF5) and pSH4 (ORF5 and ORF6) were compared extensively and a phylogenetic analysis was performed (Fig. 5) to look for similar disruptions in the repB sequences and their possible effects on functionality. The disruption in the repB protein in pSH4 gave rise to ORF5 and ORF6, two short repB sequences encoding 127 and 186 amino acids, respectively, indicating that the ORF6 was the longest conserved region of the repB protein in pSH4, close to the only 228 amino acid protein encoded by ORF5 in pSH12. Phylogenetic analysis showed that the pSH12 protein was conserved and structurally similar to other repB proteins described in Shewanella species (Fig. 5). In addition, the repB proteins of both plasmids also were closely related to the repB protein of Pseudoalteromonas sp, Thiomicrospira sp. and other Shewanella species (Fig. 5).