Flavobacterium spartansii (T16^T) was grown aerobically in CYE broth. Escherichia coli DH5α was used for cloning and plasmid propagation. E. coli S17 (λ pir) was used for conjugation. E. coli EC100D pir⁺ was used to recover transposons from F. spartansii. E. coli strains were routinely grown in Luria-Bertani (LB) broth. Liquid cultures were incubated with shaking (ca. 200 rpm) at either 28°C (F. spartansii) or 37°C (E. coli). For solid LB medium, Bacto agar (Difco, Detroit, MI) was added at a final concentration of 20 g/liter with kanamycin (50 μg/ml) or ampicillin (100 μg/ml) added for plasmid selection in E. coli or erythromycin (Em) added (100 μg/ml) for transposon or plasmid selection in F. spartansii. Demonstration of motility was done by culturing F. spartansii cells on PY2 agar medium at 28°C.

Isolation and purification of bacterial genomic DNA was performed with the Wizard Genomic DNA Purification Kit (Promega, CA, USA) according to the manufacturer's instruction. The purity, quality and quantity of purified genomic DNA were assessed using a NanoDrop 2000 UV-Vis spectrophotometer (Thermo Scientific, MA, USA) and a Qubit 2.0 fluorometer (Life Technologies, MA, USA), respectively.

Random transposon mutagenesis in F. spartansii T16^T was performed using pHimarEm1 (Braun et al., 2005). Donor E. coli S17-1 λ pir cells with pHimarEm1 and recipient cells T16^T were grown to the log phase. Cells were collected by centrifugation at 3,800 rpm for 15 min. F. spartansii T16^T and E. coli cells were mixed at 1:1 ratios in CYE medium supplemented with Ca²⁺ and spotted onto CYE agar. After 24-h incubation at 28°C, the mixed cells were scraped from agar plates, washed once with CYE broth, and plated on the CYE agar containing erythromycin (100 μg/ml). Next, erythromycin transconjugants were transferred to PY2 agar (8%) for screening the spreading-deficient mutants. The gliding deficiency mutants were further observed for their motility on glass slides under microscopy. The pHimarEm1 insertion sites in the non-spreading mutants were identified using a published method (Braun et al., 2005). Briefly, genomic DNA was digested by XbaI overnight, purified with Qiagen PCR purification kit, and self-ligated. The ligation mixture was electroporated into E. coli Ec100D⁺ (Epicentre Technologies). Purified plasmid DNA from the Km-resistant colonies was sequenced using the primer Walker85 (TGGGAATCATTTGAAGGTTGG) and Walker86 (TCGGGTATCGCTCTTGAAGGG). The sequences were blasted against F. spartansii T16^T genome to characterize the insertion site. Growth between the selected mutants and WT was compared by determination of the OD600_nm after overnight culturing in PY2 broth (Braun et al., 2005).

Hemolysin production in F. spartansii T16^T was examined by inoculating bacteria on Remel Blood Agar (Thermo Scientific, KS). Hemolytic activity was evaluated following incubation at 28°C for 48 h. The controls for α- and β-hemolysin producers were Elizabethkingia meningoseptica ATCC 13253 and Staphylococcus aureus MSU001, respectively (Chen et al., 2017).

Flavobacterium spartansii isolates T16^T and S12 were tested for antibiotic susceptibility using the Kirby-Bauer disk diffusion method as previously described (Loch and Faisal, 2014b). In brief, cultures grown on Hsu-Shotts medium (48 h) were re-suspended in sterile 0.85% saline and adjusted to OD600_nm of 0.5 in a Biowave CO8000 Cell Density Meter (WPA Inc., Cambridge, UK). The bacterial suspension in duplicate was inoculated onto dilute Mueller-Hinton agar (Hawke and Thune, 1992) without 5% calf serum. Antibiotic-impregnated disks were stamped onto the medium and plates incubated at 22°C for 48 h, at which time the zones of inhibition were measured and means calculated. Antibiotics included polymyxin-B (PB; 300 iu), oxytetracycline (T; 30 μg), trimethoprim-sulfamethoxazole (SXT; 25 μg), erythromycin (E; 15 μg), ampicillin (AMP; 10 μg), florfenicol (FFC; 30 μg), penicillin G (P; 10 iu), and the vibriostatic agent O/129 (O129; 2,4-diamino,6,7-di-isopropyl pteridine; 10 μg). The cutoffs for sensitivity were as follows: PB < 8 mm resistant; 8–12 mm intermediate sensitivity; >12 mm sensitive; T < 15 mm resistant; 15–18 mm intermediate sensitivity; >18 mm sensitive; SXT < 11 mm resistant; 11–15 mm intermediate sensitivity; >15 mm sensitive; E < 14 mm resistant; 14–18 mm intermediate sensitivity; >18 mm sensitive; AMP < 12 mm resistant; 12–13 intermediate sensitivity; >13 mm sensitive; FFC < 16 mm resistant; 16–21 mm, intermediate sensitivity; >21 mm, sensitive; O/129, < 7 mm resistant, >7 mm sensitive (Whitman, 2004).

Next generation sequencing (NGS) libraries were prepared using the Illumina TruSeq Nano DNA Library Preparation Kit. Completed libraries were evaluated using a combination of Qubit dsDNA HS, Caliper LabChipGX HS DNA, and Kapa Illumina Library Quantification qPCR assays. Libraries were combined in a single pool for multiplexed sequencing and this pool was loaded on one standard MiSeq flow cell (v2) and sequencing was performed in a 2 × 250 bp paired end format using a v2, 500 cycle reagent cartridge. Base calling was done by Illumina Real Time Analysis v1.18.54 and the output of RTA was demultiplexed and converted to a FastQ format with Illumina Bcl2fastq v1.8.4 (Klein et al., 2013). NGS libraries were sequenced by Illumina Miseq paired-end sequencing technology at the Research Technology Support Facility (RTSF) at Michigan State University.

The reads were assembled using SPAdes (version 3.9.0. Gene annotation was carried out by NCBI Prokaryotic Genome Automatic Annotation Pipeline (PGAAP 3.3) (http://www.ncbi.nlm.nih.gov/genome/annotation_prok/). Initial prediction and annotation of open reading frames (ORF) and tRNA/rRNA gene prediction were carried out with Glimmer 3 through the Rapid Annotation using Subsystem Technology server (RAST).

The functional categorization and classification for predicted ORFs were performed by RAST server-based SEED viewer (Aziz et al., 2008). The multi-drug resistance genes were predicted by Resistance Gene Identifier (RGI) tool implemented in CARD website (McArthur et al., 2013). Determination of phage elements was conducted in PHASTER (Arndt et al., 2016). Identification of genomic islands was performed in IslandViewer3 using F. johnsoniae genome as the reference (http://www.pathogenomics.sfu.ca/islandviewer). Prophage and Clustered Regularly Interspaced Short Palindromic Repeats (CRISPR) were predicted using the CRISPRfinder (Grissa et al., 2007). Detection and identification of virulence factors were carried out using MvirDB database (Zhou et al., 2007). Some virulence factors found in Bacteroidetes were manually compared among those fish pathogens and environmental isolates (Chen et al., 2005). For genome similarity assessment, digital DNA-DNA hybridization (dDDH) values were computed using web tool GGDC 2.1 (formula 2, identities/HSP length). Further, the average nucleotide identity (ANI) values were computed by ANI calculator (http://enve-omics.ce.gatech.edu/ani/).

The pan genome, core genome, and specific genes of F. spartansii were analyzed by comparison with 13 representative Flavobacterium using EDGAR 2.0 (Blom et al., 2016). The sizes of pan genome and core genomes were approximated using the core/pan development feature.

This Whole Genome Shotgun project has been deposited at DDBJ/ENA/GenBank under the accession MIKE00000000. The version described in this paper is version MIKE00000000. The BioProject designation for this project is PRJNA341797. BioSample accession number is SAMN05728473.

The general genome features of F. spartansii compared to those representative Flavobacterium, and assembly statistics, are presented in Table 1. The assembly of F. spartansii produced 5.35 Mb of sequence across 29 contigs with an N₅₀ of 515,812 bp, the longest sequence of 1,308,521 bp and a G+C content of 36% (Table 1). The genome was predicted to have at least 4,422 protein-coding genes and 130 RNA genes. Among the protein-encoding genes, 2,864 of them could be assigned a putative function, whereas 1,772 genes were predicted to encode hypothetical proteins. At least 2,119 proteins were assigned to 25 different functional categories with 368 subsystems using SEED subsystems by RAST analysis (Figure S1). The annotated genome had 112 genes responsible for bacteriocins (2 genes), resistance to antibiotics and toxic compounds (84 genes), and invasion and intracellular resistance (26 genes; Figure S1).

The ANI and dDDH values of representative Flavobacterium species are summarized in Table 2. F. spartansii and Flavobacterium sp. WG21 have an ANI value of 96.9% and a dDDH value of 72.7% (Table 2), suggesting that they are the same species according to the microbial taxonomy for species delineation (95 and 70% cut-off for ANI and dDDH, respectively) (Goris et al., 2007).

The genome size for ubiquitous genes (core genome) or entire set of genome (pan-genome) amongst the selected Flavobacterium genomes was plotted against the number of genomes (Figure 1). The pan-genome plot showed that the power trend line had not reached the plateau (Figure 1), indicating that Flavobacterium possess an open pan-genome. Core genome analysis showed that the number of shared genes decreased with the addition of the input genomes and was predicted to converge against 1,182 (Figure 1). Actually, the core genome for the selected Flavobacterium in comparison was calculated to be 1,125 CDS per genome; thus, given the assumption that the prediction slightly overpredicts the core genome size, the current core genome represents the Flavobacterium genus quite well. F. spartansii shared at least 3,649, 3,201, 2,996, 2,774, 1,971, and 1,868 CDS with those environmental or opportunistic pathogenic Flavobacterium including Flavobacterium sp. WG211667 (87.2% of its total encoding genes), F. hibernum DSM 12611 (70.1%), F. denitrificans DSM 15936 (71.7%), F. beibuense F44-8 (69.6%) F. daejeonense DSM 17708 (68.2%), and F. johnsoniae UW101 (61.6%) (Figure 2). On the other hand, it shared up to 1,627, 1,651, and 1,712 CDS with the selected pathogenic Flavobacterium species such as F. columnare strains (~81%) and three F. psychrophilum strains (~75%) and F. branchiophlium (~71.3%).

Antibiotic susceptibility assays were conducted using the disk diffusion method (Table 3). Of the 8 tested antibiotics, F. spartansii T16^T was intermediately sensitive to trimethoprim-sulfamethoxazole, polymyxin-B, florfenicol, erythromycin as well as oxytetracycline, while it was resistant to ampicillin, penicillin G, 2,4-diamino and 6,7-di-isopropyl pteridine (Table 3). The resistome in strain S12 was similar to that in strain T16^T. However, compared to strain T16^T, strain S12 was more sensitive to low concentrations of erythromycin and immediate sensitive to ampicillin. Neither of the two F. spartansii isolates (Table 3) were completely sensitive to the two antibiotics (i.e., florfenicol and oxytetracycline) that are approved to treat certain flavobacterial diseases (e.g., bacterial coldwater disease caused by F. psychrophilum, and columnaris disease caused by F. columnare) in food fish (Bowker et al., 2016).

Both the RAST SEED subsystem and CARD were used to identify antibiotic resistance genes within the F. spartansii genome using the default settings. The predicted genes conferring resistance to aminocoumarin, mupirocin, chloramphenicol, fluoroquinolone, tetracycline, rifampin, and β-lactam are summarized in Table S1. Interestingly, up to 10 β-lactam resistance-related genes were found in the F. spartansii T16^T genome including those encoding the putative β-lactamases (8 copies), metallo-β-lactamases (1 copy) and penicillin-binding protein (1 copy). This observation is consistent with resistance to both penicillin and ampicillin in our antimicrobial tests (Table 3). Most β-lactamases in F. spartansii showed very low identity (<50%) to other fish pathogens (F. psychrophilum PG2, F. columnare and F. branchiophilum FL-15) (Table S1). However, a few lactamases have relatively high identity to those in F. johnsoniae and Flavobacterium sp. WG21. No erythromycin resistance related genes were detected in F. spartansii. More drug resistance genes and non-specific efflux pumps (5 hydrophobic/amphiphilic exporters, HAE1 family) possibly conferring antibiotic resistance are also conserved in Flavobacterium sp. WG21 and F. hibernum DSM 12611.

At least 5 prophages were predicted in F. spartansii by PHAST (Table S2). Among them, the prophage 1 (29.4 kb) consisted of 20 CDSs that encodes 12 proteins with known functions and 8 hypothetical proteins. Of interest, prophage 1 may be complete, as it consists of a phage tail, head, portal, integrase, lysin, terminase, and other component proteins involved in phage structure and assembly. The GC content (35.9%) in prophage 1 is very close to the average of GC in the whole genome (35.5%), indicating that prophage 1 may have been integrated into F. spartansii genome long ago. Unlike prophage 1, the remaining 5 prophages seem to be incomplete, as evidenced by the absence of a portion of phage structure proteins (Table S2). However, these regions carried genes encoding some important enzymes which may contribute to flavobacterial fitness under certain conditions. For example, prophage 2 had a DNA polymerase III beta subunit, methylase, NinG recombination, ATPase and Rec T protein, indicating this region may participate in DNA replication and modification. Prophage 3 carried many genes encoding enzymes related to polysaccharide (capsule precursor) metabolism; there were at least two 1,3-N-acetylglucosaminyltransferases, four glycosyl transferases, one GDP-D-mannose 4,6-dehydratase, one alpha-1,2-fucosyltransferase, one fucose synthetase, one dTDP-d-glucose 4,6-dehydratase, one L-fucosamine transferase and one UDP-N-acetylglucosamine 2-epimerase (Table S2). Prophage 4 may be linked to RNA metabolism because ribokinase-like domain-containing protein, RNA polymerase, ribonuclease H and putative ribonuclease 3 were organized as a cluster. Prophage 5 had genes encoding arsenate reductase and metallophosphoesterase, which possibly contribute to heavy metal resistance.

At least 8 genomic islands (GIs) were identified in F. spartansii (Figures 3A,B and Table S3) by both IslandPick and IslandPath-DIMOB methods (Nagamatsu et al., 2015). The GI size ranged from 11.3 to 63.8 Kb (Figure 3A). Genes encoding virulence factors, toxins, DNA metabolism, transposases, regulators, modification and restriction systems and resistance to antibiotics occurred in these GIs (Figure 3 and Table S3), indicating that F. spartansii T16^T possibly acquired these genes, thereby forming GIs favoring adaptation to diverse environments (see below). Among the 8 GIs, GI-8 had the largest size (63.7 Kb). Notably, GI-8 has many genes encoding components of a conjugative transposon in Bacteroidetes (here, named CTnFs), e.g., traN, traM, traK, traJ, traI, and traG. Moreover, the remaining components of this large conjugative transposon CTnFs can be found in GI-7 (33.3 kb) including traD, traC, traE, traF, traG, traI, traJ, traK, traM, traI, and traQ genes (Table S3). We further examined the surrounding regions around GI-7 and found the traA and traC genes which are located upstream of traD (Table 4). The distribution spectrum of the CTnFs-like conjugative transposons were investigated by searching the conserved transfer protein gene traG (a signature gene for conjugative transposons) among the selected Flavobacterium species. No traG gene sequences (cutoff 50%) were found in F. psychrophilum, F. columnare, and F. branchiophilum FL-15. The same results were obtained when using other tra genes to search against their genomes, indicating that the conjugative transposon is absent in these fish pathogens. By contrast, the similar conjugative transposons were observed in F. johnsoniae, Flavobacterium sp. W21, F. hibernum DSM 12611, F. daejeonense DSM 17708, F. denitrificans DSM 15936, and F. beibuense F44-8 (Table 4).

GI-2 had many genes encoding respiratory chain components such as cytochrome C subunits, cytochrome P450, quinol:cytochrome C oxidoreductase, and several other membrane proteins (Table S3). GI-2 had the predicted virulence factor VirJ, possibly participating in the Type IV secretory pathway (Table S3). Immediately downstream of GI-2, there are genes involved in the gliding motility (RemB and SprA). Another large genomic island, GI-6, contained a gene cluster encoding DNA repair proteins and enzymes such as MsrA (preventing oxidative damage), NTP pyrophosphohydrolases, chaperone protein DnaK, deoxyribodipyrimidine photolyase (DNA repair enzyme), DNA alkylation repair enzyme, peroxiredoxin, and DNA-damage-inducible protein D (Table S3). On the same GI, there were transcriptional regulators from GntR, HxlR, and AraC families, showing that the GI-1 may be involved in the defense against oxidative stress and DNA damage repair. Compared to the above 4 GIs, the remaining ones (GI-1, GI-3, GI-4, and GI-5) were smaller (Table S3).

When the proteome of F. spartansii T16^T was used in a BLAST search against MvirDB database, 1,001 virulence factors were predicted (cutoff E-10) (Table S4). Further, 133 putative virulence proteins were found to be highly conserved in F. spartansii T16^T (cutoff 50%) (Table S4). Virulence factor matches in the MvirDB database include well-known proteases, catalase/peroxidase HPI, conjugal transfer proteins (4 Tra components), molecular chaperones (GroEL, DnaJ, and DnaK), translation initiation factor IF-1, transcription regulation factors, antimicrobial resistance proteins, LPS biosynthesis protein WbpP, biopolymer transporter (ExbD and SusD/RagB), ribosomal proteins, carbon and nitrogen metabolism related proteins as well as iron-binding proteins (Table S4).

At least three hemolysin genes were present in F. spartansii (Table 5), which complements observations that strain T16^T has in vitro hemolytic activity (Figure 4). When grown on sheep's blood agar plate for 48 h, colonies of F. spartansii T16^T appeared greenish, suggesting α-hemolytic activity (Figure 4). Moreover, many genes encoding proteolytic enzymes were discovered in the genome of F. spartansii T16^T (Table 5). Further, the amino sequence of a collagenase (BHE19_RS20565, peptidases U32) was conserved in most of the selected Flavobacterium (Table 5).

In F. spartansii, at least 21 genes encoding proteins associated with gliding motility were identified, including gldA, gldB, gldC, gldD, gldE, gldF, gldG, gldH, gldI, gldJ, gldK, gldL, gldM, gldN, SprA, SprB, SprC, SprD, SprE, SprF, and SprT (Table 6). The motility protein sequences were conserved among these Flavobacterium genomes except GldE, SprB, and SprC (Table 6). Genes encoding PorV and RemA also existed in F. spartansii T16^T (Table 6). Therefore, F. spartansii possesses the minimal requirement of the full function of T9SSs (also important for gliding) consisting of GldK, GldL, GldM, GldN, SprA, SprE, and SprT (McBride and Nakane, 2015).

A transposon used for random mutagenesis in flavobacteria, pHimarEm1, was functional in F. spartansii T16^T. Out of 2,000 conjugants, at least 4 mutants were recovered with impaired ability to spread on PY2 agar or glide on the surface of glass slides (Figure 5). For mutant 2F2-2, the transposon inserted 6-bp upstream of the ATG start codon of a gene encoding a hypothetical protein (BHE19_RS13930). For mutants 7B5-5, the insertion of pHimarEm1 transposon was 14 bp downstream of the ATG start codon of an aconitase gene (BHE19_RS00150). In mutant 10E3-2 the transposon inserted 1,380 bp downstream of the ATG start codon of the gldM gene (BHE19_RS21040). For the mutant M69, the transposon inserted at 67 bp downstream of the ATG start codon of EpsM gene (BHE19_RS15305). The difference in biomass determination was negligible between the WT and mutants, indicating that these gliding genes were not critical for cell growth in PY2 broth (Figure 5).

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.