Background {#Sec1}
==========

Vibrios represent a large subgroup of the class *Gammaproteobacteria* and are Gram-negative, usually motile rods. Generally, members of the *Vibrio* species are mesophilic and chemoorganotrophic, and have a facultatively fermentative metabolism \[[@CR1]\]. Vibrios are ubiquitously distributed in water columns of marine and estuarine water, marine sediments, and aquaculture settings worldwide \[[@CR2]\]. It has been shown that vibrios are characterized by remarkable biodiversity, having evolved to develop complex lifestyles \[[@CR3]\]. Many studies have indicated that vibrios have enormous biomass in the ocean, and comprises more than 10% of the aquatic culturable bacterial population \[[@CR4]--[@CR7]\]. To date, more than 110 species are recognized in the genus *Vibrio*, many of which have been identified in recent years.

As representative obligate heterotrophs, vibrios utilize a wide range of carbon sources for energy \[[@CR8]\], an interesting feature of which is the utilization of chitin. Chitin is the most abundant polysaccharide in the ocean, and the second most abundant polysaccharide in nature after cellulose \[[@CR9]\]. Chitin is composed of repeating units of the monomer *N*-acetylglucosamine (GlcNAc), and this material is an abundant carbon-nitrogen-energy source for marine microorganisms. Cycling of chitin is mainly by chitinolytic bacteria \[[@CR10]\]. Many *Vibrio* species, such as *V. cholerae* \[[@CR10]\], *V. harveyi* \[[@CR11]\] and *V. furnissii* \[[@CR12]\], are able to use chitin as their sole carbon source. A group of enzymes is required for vibrios to successfully convert chitin to monomer \[[@CR10]\], in which process chitinases function to hydrolyze chitin into small oligomers. Chitinases have been received much attention and were believed to play a central role in the chitin degrading processes by bacteria \[[@CR13]\]. Many physiological characteristics of vibrios were associated with the capacity of chitin-degradation, including chemotaxis, cell division, biofilm formation and natural competence \[[@CR14]\].

The genomic features of vibrios in the past decades have been well studied, and with the increasing number of sequenced species, there are growing number of reports which focus on the genus rather than on specific members of the genus. The diversity and phylogenetics of vibrios has long been a matter of interest, and multilocus sequence analysis (MLSA) and constructing supertrees using housekeeping genes from whole genome sequencings has been used in the analysis \[[@CR15]\], which was believed to be the most powerful tool until alternative genome-based methods were developed. Sawabe et al. \[[@CR16]\] have studied the evolutionary history of vibrios using 8 housekeeping genes, and suggested the method to be very efficient and reliable. While comprising a considerable pan-genome \[[@CR15], [@CR17]\], the genome plasticity of vibrios aided by horizontal gene transfer (HGT) was considered to be an important method of acquiring new genes and accelerating evolution \[[@CR18]\], and this rapid change in the composition of the genomes mainly related to cell envelope modification and secondary metabolites, helping vibrios to succeed in competition with other organisms \[[@CR19]\]. The wide range of metabolisms in vibrios, which facilitates them to adapt to a wide array of environmental conditions has drawn much attention, especially the capability of degrading chitin \[[@CR8], [@CR12]\]. Nevertheless, although many studies mentioned above have been carried out to investigate the diversity, evolution and ecology of vibrios, many fundamental works are still needed in order to discover the underlying genetic basis, and this is now possible as the novel genomes and productive tools have recently become available.

For a long time, draft genomes have been employed to perform genomic studies. However, a complete sequence is better than partial genomes to examine all the genes or study the complete hereditary information of an organism \[[@CR20]\]. The advent of the improved sensitivity and precision of complete genome sequencing technology provides a primary important path toward understanding the genome features of organisms. Lukjancenko et al. \[[@CR17]\] noticed that the core gene numbers of *Vibrio* chromosome in complete genomes dropped sharply to only half of that when considering draft genomes, which was partly due to the incomplete information stored in the draft genomes. Besides, the rRNA operon number of draft *Vibrio* genomes is only one in most cases, whereas the complete genomes have seven or eight copies generally. Given the increased availability of complete genome sequences in the public databases and the development of comparative genomic analysis approaches, a more integrated and comprehensive view of the genus *Vibrio* becomes possible. There were more than 1500 genomes of vibrios available in the GenBank database at the time of this study, 18 of which were complete genomes sequences of the *Vibrio* species.

In this study, we sequenced two new complete *Vibrio* genomes, and presented a broad comparative analysis of these genomes across the genus with 18 other complete genomes of *Vibrio* spp. available in the GenBank database (<http://www.ncbi.nlm.nih.gov>). Firstly, the phylogeny of the genus *Vibrio* was inferred using core-genome. Secondly, the ancestral reconstruction and the estimation of emerge time of vibrios were performed to provide new knowledge about the evolutionary history of the genus. Thirdly, various environmental adaptation strategies, especially the ability of the chitin metabolism, were investigated to help understand how vibrios adapt to a variety of habitat niches. Finally, all the *Vibrio* genomes, including those draft sequences in the public database, along with a set of other representative marine bacteria, were employed to provide a comprehensive assessment of genome features. Our analysis helps in the understanding of the genetic characteristics of *Vibrio* and serves as a guide for further research.

Methods {#Sec2}
=======

Bacterial strains {#Sec3}
-----------------

The two strains *V. rotiferianus* B64D1 and *V. mediterranei* QT6D1 were isolated from bottom seawater of the Bohai Sea (119.04°E, 38.23°N, 17.5 m water depth) in August 2015 and the East China Sea (122.0051°E, 26.8957°N, 92 m water depth) in October 2015, respectively. No animal ethics approval was required, and the field sampling procedures met local guidelines. Both strains were isolated using thiosulfate citrate bile salts sucrose (TCBS) agar (Hopebio, China) and demonstrated strong capacity to degrade chitin when growing on chitin agar plates. The strains were routinely grown aerobically on marine 2216E agar (MA; Becton Dickinson) at 28 °C. Genomic DNA of the two strains were extracted using the phenol-chloroform-isoamylic alcohol extraction protocol described by Marmur \[[@CR21]\], and the 16S rRNA genes were sequenced to validate the strains.

Genome sequencing and annotation {#Sec4}
--------------------------------

Whole-genome shotgun sequencing of *V. rotiferianus* B64D1 and *V. mediterranei* QT6D1 was carried out with paired-end sequencing by HiSeq 4000 (Illumina) and long sequencing by PacBio RS II Sequencer (Pacific Biosciences, Menlo Park, CA). The long sequence was assembled using Canu v1.1 \[[@CR22]\], and Pilon v1.16 \[[@CR23]\] was used to polish PacBio assemblies according to Illumina data. Misassemblies were identified and corrected using REAPR \[[@CR24]\] and manual inspection. Open reading frame (ORF) prediction and annotation were carried out using the RASTtk pipeline \[[@CR25]\]. tRNA-encoding genes and rRNA operons were found by using tRNAscan \[[@CR26]\] and RNAmmer \[[@CR27]\], respectively. Restriction-modification (R-M) systems were identified by searching the RAST annotation, and CRISPRfinder \[[@CR28]\] was used to search for CRISPR (clustered regularly interspaced short palindromic repeats) regions. When the CRISPR regions was recognized as questionable CRISPRs, the presence of Cas (CRISPR-associated protein) genes was considered evidence for a functional CRISPR system.

The two complete genomes sequenced in this study and a further 18 complete genomes of representative strains in distinct species of the genus *Vibrio* were collected together as the "small dataset" (Additional file [1](#MOESM1){ref-type="media"}: Table S1), while all the 1582 *Vibrio* genomes available in GenBank except those with low qualities (Contig N50 \< 10,000 bp) were retrieved to make up the "large dataset". *Vibrio* species included in the small dataset were isolated from various sources, such as oceans, estuarine environments, humans, marine fishes and marine molluscs, which enabled the provision of a holistic profile of the genus. In particular, *V. fischeri* (*Allivibrio fischeri*) and *V. damselae* (*Photobacterium damselae* subsp. *damselae*) had been considered involved in the genus *Vibrio* for a long time. Although they have been reclassified as other closely related genus \[[@CR29], [@CR30]\], we decided to involve them in this study to acquire richer information on the evolution and forces for genome changes of the group due to their relatively unique living strategies.

Pan/core genomes and phylogenetic analysis {#Sec5}
------------------------------------------

The program GET_HOMOLOGUES v3.0.3 \[[@CR31]\], with three popular clustering algorithms, i.e. bidirectional best hit (BDBH), COGtriangles \[[@CR32]\] and OrthoMCL \[[@CR33]\], was used for clustering orthologous genes and identifying core- and pan-genomes within the small dataset mentioned above under default parameter values. Sequences of all the single-copy orthologues were aligned individually using the MAFFT software with the "linsi" algorithm \[[@CR34]\], and the alignments were concatenated into one extended alignment. Subsequently, TrimAl v1.2 \[[@CR35]\] was applied to remove ambiguous and unreliable fragments from the concatenated alignment with the parameter "-gt 1". The phylogenetic relationship of the concatenated data was inferred using RAxMLv8.2.4 \[[@CR36]\] based on the Maximum-likelihood (ML) algorithm with the PROTGAMMALGF evolutionary model and 1000 bootstrap replicates, after Prottest 3.4 \[[@CR37]\] was performed to estimate the best-fit substitution model. The genome of *Shewanella denitrificans* OS217 was chosen to serve as an outgroup. The phylogenetic network of the concatenated data was constructed using SplitsTree v4.14.6, with a neighbor net drawing and Jukes-Cantor correction \[[@CR38]\]. A pan-genome tree was constructed according to the presence/absence of all homologous genes using RAxML. Besides, we identified the specific genome (unique gene pool of a species) of the two new strains we sequenced in this study among the noncore region, in order to investigate the unique properties of them.

A chronogram for 23 taxa RAxML phylogeny was constructed using the method of penalized likelihood \[[@CR39]\], implemented in r8s v1.81 \[[@CR40]\]. The estimated origin time of cyanobacteria (2700\~ 3500 million years ago, MYA) \[[@CR41]\] and the divergence time between *Escherichia coli* and *Salmonella enterica* (140 MYA) \[[@CR42]\] were used as calibration points. The reference ML tree was built based on 332 single-copy core orthologues shared among all the 23 bacterial species. The optimal smoothing parameter was determined through a cross-validation procedure that involved pruning terminal branches.

Reconstruction of gene gain and loss events {#Sec6}
-------------------------------------------

Reconstruction of ancestral gene content and the gene gain and loss during the evolution of the small dataset was conducted using the phylogenetic gain--loss-duplication model employed in the program COUNT \[[@CR43]\]. The number of paralogs in each genome for all the identified gene families was taken into account in the reconstruction, and the topology of the core-genome tree was used as the reference tree. All the genes acquired or lost during the evolutionary process of the genus *Vibrio* were categorized based on searching against the COG database ([http://www.ncbi.nih.gov/COG/](http://www.ncbi.nih.gov/COG)) using RPS-BLAST program. SIGI-HMM \[[@CR44]\] was used to infer whether genes in the genomes were acquired through HGT.

Average nucleotide identity (ANI) calculation {#Sec7}
---------------------------------------------

ANI between any two genomes of the small dataset was calculated using a custom Perl script adapted from ANI.pl (<https://github.com/jhbadger/scripts/blob/master/ANI.pl>), following the algorithm described by Goris et al. \[[@CR45]\], and the resulted matrix was clustered and visualized using R packages pheatmap software \[[@CR46]\].

Analysis of genes related to chitin degradation {#Sec8}
-----------------------------------------------

Functional proteins related to chitin metabolism in all annotated genomes were determined according to the chitin pathway defined for *V. cholerae* \[[@CR8]\] and by being searched against all the annotations. Among these predicted chitinolytic enzymes, the evolutionary history of chitinase A in the small dataset was inferred using the Neighbor-Joining method \[[@CR47]\] implemented in MEGA7 \[[@CR48]\], and the phylogenetic tree was further tested for the selection pressure among amino acid sites using several codon-based ML procedures as implemented in the codeml tool of PAML v4.9 \[[@CR49]\]. The conserved domains of chitinase A were defined by the NCBI batch CD-search online tool (<https://www.ncbi.nlm.nih.gov/Structure/bwrpsb/bwrpsb.cgi>).

Genomic comparison with other marine bacteria {#Sec9}
---------------------------------------------

Thirteen complete genomes of common marine bacteria (Additional file [2](#MOESM2){ref-type="media"}: Table S2) were retrieved from GenBank to make comparisons with *Vibrio* genomes. These genomes are from four different classes including *Cyanobacteria*, *Alphaproteobacteria*, *Gammaproteobacteria*, and *Flavobacteria*, which are considered to account for up to 80% of the total bacterioplankton in the ocean \[[@CR50]\]. Basic features of these 13 genomes, including genome sizes, 16S rRNA gene copy numbers, the number of genes derived from HGT process and genes encoding chitinase, were inferred according to the same methods as for the *Vibrio* species described above.

Database accession numbers {#Sec10}
--------------------------

The complete genome sequences of *V. rotiferianus* B64D1 and *V. mediterranei* QT6D1 have been deposited in NCBI GenBank server under the accession number CP018311 to CP018312 and CP018308 to CP018310, respectively.

Results {#Sec11}
=======

General features of *Vibrio* species {#Sec12}
------------------------------------

The general information of the 20 complete *Vibrio* genomes is summarized in Table S1. Each genome contains two chromosomes, and some species have one or more plasmids. The genome size of vibrios ranges from 4.09 to 6.32 Mb, with a mean size of 5.14 Mb. Most the other marine bacteria compared in this study have a lower genome size (3.62 Mb on average) compared to vibrios (Fig. [1a](#Fig1){ref-type="fig"}). Vibrios have a wide range of G + C content, ranging from 38.3% in *V. fischeri* to 50.7% in *V. furnissii*. The ORF number and orthologue number range from 3807 to 5920 and 3669 to 5593, respectively. Strikingly, vibrios contain a much higher number of 16S rRNA gene copies (7--15) than most of the other marine bacteria (3 on average) used for comparison in our study (Fig. [1b](#Fig1){ref-type="fig"}), which may reflect the capability of rapid growth of *Vibrio* species. In addition, the genome size of the large dataset ranges from 2.9 to 6.7 Mb, and the G + C content ranges from 38.0% to 57.2%, which suggests a higher diversity of *Vibrio* genomes. However, due to the missing information of those draft genomes, there is some bias definitely existing in the large dataset.Fig. 1Comparison of genomic characteristics between the 20 vibrios with complete genomes and other representative marine bacteria. Box-and-whisker plot showing data distribution and interquartile ranges for each group with the median value displayed as a line in each box plot. Outliers are shown as circles, and mean values are shown as saltires in the figure. **a** Comparison of the genome sizes. **b** Comparison of the copy numbers of 16S rRNA gene. **c** Comparison of the percentage of HGT genes. **d** Comparison of the number of chitinase coding genes

R-M and CRISPR systems {#Sec13}
----------------------

Microbes rely on diverse defense mechanisms that protect them from mobile elements or foreign DNA implicated in the reduction of HGT. The common bacterial defense systems include R-M systems and CRISPR systems \[[@CR51]\]. We identified these two defense systems in the small dataset (Additional file [1](#MOESM1){ref-type="media"}: Table S1), and found that almost all the *Vibrio* genomes contained the R-M system except for *V. breoganii* FF50 and *V. tritonius* JCM 16456. Specifically, all the vibrios containing R-M systems have type I R-M system, while type II R-M system was additionally harbored by *V. cholerae* KW3, *V. rotiferianus* B64D1, *V. tubiashii* ATCC 19109 and *V. vulnificus* FORC 017, and type III R-M system was additionally harbored by *V. coralliilyticus* RE98. Eight out of the 20 genomes had CRISPR systems (Additional file [1](#MOESM1){ref-type="media"}: Table S1), while the cas gene numbers (1 to 8) and the lengths of direct repeats (DRs) of each CRISPR (23 to 51) were distinct from each genome. The identification of DRs in the CRISPR would be quite difficult if performed in draft genomes, since the short and repetitive sequencing reads would be difficult to assembly correctly.

Core and pan-genome {#Sec14}
-------------------

To determine the global gene repertoire and the genes conserved across all the vibrios, the power law model of Tettelin et al. \[[@CR52]\] was used to describe the pan- and core-genomes based on the 20 complete *Vibrio* genomes. All the genes were grouped into 21,844 families, which represented the *Vibrio* pan-genome. The progression of the pan- and core-genomes after ten random samplings was calculated (Additional file [3](#MOESM3){ref-type="media"}: Figure S1). The number of novel gene families in the pan-genome gradually increased with no limitation when more genomes were considered, leading to an open pan-genome. Nevertheless, the number of conserved gene families constituting the core genome decreased slightly, and the extrapolation of the curve indicated that the core genome reached a minimum of 1630, and would remain relatively constant, even as many more genomes were added.

Specific genomes of V. rotiferianus B64D1 and V. mediterranei QT6D1 {#Sec15}
-------------------------------------------------------------------

To investigate the unique characters of the two complete genomes sequenced by us, we identified the specific genome of each genome, i.e., orthologous families that contain genes exist in each of these two genomes but absent in other 18 previously sequenced genomes. In the total 21,844 pan-genome families, 397 specific genes were possessed by *V. rotiferianus* B64D1, whereas 871 were possessed by *V. mediterranei* QT6D1. The COG functional annotation of the specific genomes (Additional file [4](#MOESM4){ref-type="media"}: Figure S2) showed that genes related to carbohydrate transport and metabolism (G) and transcription (K) are extremely abundant (53 and 40, respectively) in *V. mediterranei* QT6D1, while the most abundant COG category in the specific genome of *V. rotiferianus* B64D1 is general function prediction only (R) and cell wall/membrane/envelope biogenesis (M). Furthermore, it is notable that there is a chitinase gene (BSZ05_18380) annotated by the RAST server involved in the specific genome of *V. mediterranei* QT6D1, which may contribute to its strong degradation capability for chitin.

Phylogenetic and divergence times analysis of vibrio species {#Sec16}
------------------------------------------------------------

To investigate the phylogenetic relationship of *Vibrio* species, a robust core-genome tree was inferred based on the concatenated 1003 single-copy orthologues proteins shared by all the 20 complete genomes as shown in Fig. [2a](#Fig2){ref-type="fig"}. Two main clades could be determined in the tree, and these two clades could be further proved by the split network analysis (Fig. [2b](#Fig2){ref-type="fig"}). Clade I consisted of five species: *V. harveyi*, *V. campbellii*, *V. rotiferianus*, *V. alginolyticus* and *V. parahaemolyticus*, which was the so-called *Vibrio* core group \[[@CR53]\] or Harveyi clade, and Clade II consisted of *V. cholerae*, *V. mimicus*, *V. fluvialis* and *V. furnissii*, which formed the Cholerae clade proposed by Sawabe et al. \[[@CR54]\]. The pan-genome tree (Additional file [5](#MOESM5){ref-type="media"}: Figure S3) also had a relatively high bootstrap values at most of the nodes similar to the core-genome tree. However, the topologies of the pan- and core-genome trees were different, which might be attributed to the dispensable genes. The phylogenetic tree based on 16S rRNA genes was also inferred (Additional file [6](#MOESM6){ref-type="media"}: Figure S4), whereas it was not as reliable as the core-genome tree insofar as the bootstrap values were extremely low.Fig. 2Phylogenetic relationships of the 20 vibrios with complete genomes. The Harveyi clade and the Cholerae clade are highlighted in red and yellow respectively. **a** Core-genome tree reconstructed by RAxML software. The tree is rooted using *S. denitrificans* OS217. The number at each node denotes the bootstrap value based on 1000 replicates, and the scale bar indicates the number of substitutions per site. **b** Phylogenetic network reconstructed by SplitsTree4 software

To explore the divergence times of vibrios, several calibration points were included in the analysis. The most recent common ancestor (MRCA) of vibrios (i.e., the emergence time of *V. damsela*) was estimated to have emerged around 536 MYA, whereas *V. fischeri* was estimated to have emerged around 385 MYA (Fig. [3](#Fig3){ref-type="fig"}).Fig. 3Chronogram of the 20 vibrios with complete genomes. The value near each internal branch is the estimated emerge time for that branch. Nodes with fossil record corrections are indicated with an asterisk

ANI values {#Sec17}
----------

ANI values calculated through pairwise comparison of the 20 *Vibrio* genomes further elucidated their genetic relatedness. ANI has been considered the best approach to measure genetic distance between genomes, as it represents a mean of identity/similarity values between homologous genomic regions shared by two genomes, which enables minimization of bias brought about by variable evolutionary rates and HGT events \[[@CR55]\]. It was demonstrated that ANI values within a genus were generally higher than 62%, and 95--96% ANI has been recommended as the delimitation criterion for different species \[[@CR56]\]. The ANI values between *Vibrio* species ranged from 70% to 88%, which was within the range of the classical boundary for a genus. Species groups could be derived after clustering according to the rows (Additional file [7](#MOESM7){ref-type="media"}: Figure S5), among which the Harveyi clade and Cholerae clade were consistent with the core-genome tree (Fig. [2a](#Fig2){ref-type="fig"}). Members of the Harveyi clade shared the highest ANI across the genus ranging from 80% to 88%, indicating extremely close genetic relatedness within the Harveyi clade. Furthermore, the clades possessing high ANI values could serve as an explanation for the relatively short terminal branch lengths of their members in the phylogenetic tree, as the ANI represented the genetic distances between genomes.

Gene gain and loss {#Sec18}
------------------

To gain insight into the origin and the variations of *Vibrio* genomes, we investigated the gene gain/loss events that might have occurred during the evolution of the genus *Vibrio*. The gene family numbers showed a noticeable upward trend during the evolutionary history, and the first time that gene loss events exceeded gene gain events was present at the MRCA of *V. anguillarum* 90--11-286 and *V. mimicus* ATCC 33654 (Fig. [4](#Fig4){ref-type="fig"}). This internal node led to a cluster that was mainly comprised by the Cholerae clade, and all the members at terminal nodes in this cluster had relatively small orthologues families (\< 4600), which may be explained by the fact that gene loss events posed more significant impact on this cluster in the evolutionary process. On the contrary, other members in this phylogenetic tree experienced more and continuous gene gain events. However, no matter whether the species have suffered more gene loss events in a few internal nodes, the gene gain events in the evolutionary process were more numerous than the gene loss events in most instances, suggesting that vibrios have experienced remarkable gene family expansion in the evolutionary history. This generally expansive trend made the orthologue number of vibrios increase by approximately 50%, i.e., from 2828 gene families in MRCA of vibrios to 4547 gene families on average as currently predicted (Additional file [1](#MOESM1){ref-type="media"}: Table S1). According to the COG functional annotation (Fig. [5](#Fig5){ref-type="fig"}) of the gain/loss genes, other than the poorly characterized genes (general function prediction only and function unknown, COG categories R and S), genes related to transcription (K) formed the most abundant category of gene gains, followed by genes involved in carbohydrate transport and metabolism (G), whereas the most abundant category of gene losses was associated with signal transduction mechanisms (T). Gene gains were present in every COG category except for extracellular structures (W) and nuclear structure (Y, not applicable to prokaryotic COGs), whereas gene losses in the history were not present in RNA processing and modification (A), extracellular structures (W), nuclear structure (Y) and cytoskeleton (Z, not applicable to prokaryotic COGs).Fig. 4Ancestral genome content reconstruction of the 20 vibrios with complete genomes. The numbers of gene gain and loss are shown above each branch, and the numbers of gene present at each node are shown below the corresponding branch. The area of the circle at each node indicate the genome size of the node. Nodes with more gene gain events (net gains \> net losses) are indicated in yellow; Nodes with more gene loss events (net losses \> net gains) are indicated in green; Nodes with significant genome expansions (\> 10% net gains) are indicated in red; Nodes with significant genome expansions and reductions at the same time (\> 10% net gains and \> 10% net losses) are indicated in blue. The bar chat on the right indicates the orthologues number of each terminal nodeFig. 5COGs classification of putative gene gains and losses occurred in the evolutionary history of vibrios**.** Designations of functional categories: \[A\] RNA processing and modification, \[B\] chromatin structure and dynamics, \[C\] energy production and conversion, \[D\] cell cycle control and mitosis, \[E\] amino acid metabolism and transport, \[F\] nucleotide metabolism and transport, \[G\] carbohydrate metabolism and transport, \[H\] coenzyme metabolism, \[I\] lipid metabolism, \[J\] translation, \[K\] transcription, \[L\] replication and repair, \[M\] cell wall/membrane/envelope biogenesis, \[N\] Cell motility, \[O\] post-translational modification, protein turnover, chaperone functions, \[P\] Inorganic ion transport and metabolism, \[Q\] secondary metabolites biosynthesis, transport and catabolism, \[R\] general functional prediction only, \[S\] function unknown, \[T\] signal transduction, \[U\] intracellular trafficking and secretion, \[V\] Defense mechanisms, \[W\] Extracellular structures, \[Y\] nuclear structure, \[Z\] cytoskeleton

HGT is a significant driver of gene innovation in the life of bacteria \[[@CR57]\]. A search for genes with aberrant codon usage statistics was carried out to predict possible genomic islands, referred to as alien genes. Surprisingly, the average percentage of alien genes among these vibrios was predicted to be 7.3%, with some species containing up to 10% alien genes. The class *Gammaproteobacteria* was the major donor of HGTs, contributing as much as 72.1% alien genes, followed by *Sphingobacteria*, which contributed 6.5% of the alien genes. This percentage of alien genes was significantly higher than that of the other compared marine bacteria, with nearly all of the genomes only involved having no more than 5% alien genes, except for species from *Marinobacter* and *Rhodobacter*, which had 9.96% and 9.25% alien genes, respectively (Fig. [1c](#Fig1){ref-type="fig"}).

Chitin degradation related genes {#Sec19}
--------------------------------

Many *Vibrio* species have been reported as capable of degrading chitin. Chitin degradation is achieved by a complex pathway including multiple chitinases \[[@CR58]\]. Seventeen different protein families related to chitin metabolism were retrieved from *Vibrio* genomes (Additional file [8](#MOESM8){ref-type="media"}: Figure S6). Eighteen out of the 20 vibrios (except for *V. breoganii* FF50 and *V. tritonius* JCM 16456) contained chitinase (EC 3.2.1.14), and all the 20 vibrios contain β-N-acetylhexosaminidase (EC 3.2.1.52). Moreover, \~ 99% of the genomes in the large dataset carried both these two genes, which have the ability to completely hydrolyze chitin to monomer GlcNAc \[[@CR59]\]. In contrast to this abundance, these genes seldom exist in the compared marine bacteria, except for strains from genus *Pseudoalteromonas* and *Shewanella*, which carried eight and two chitinase encoding genes in their genomes, respectively (Fig. [1d](#Fig1){ref-type="fig"}).

Chitinases distributed in these 18 *Vibrio* species come from 15 protein families, which can be classified into four categories (Fig. [6](#Fig6){ref-type="fig"}). Class I encompasses a family GH18 chitinase domain, while Class II possesses a similar chitinase domain as Class I but with an insertion domain in it. Class III encompasses a family GH19 chitinase domain, which was considered to be rare in bacteria \[[@CR60]\]. However, all the 18 *Vibrio* genomes carried GH19 chitinase according to our analysis, indicating that GH19 chitinase is abundant in vibrios. Class IV does not contain any typical chitinase domain, but could be annotated as chitinase using blastp against the nr database (<ftp://ftp.ncbi.nlm.nih.gov/blast/db/>), which may result from the relatively low conservation of these chitinases. Most studies on diversity of chitin-degrading enzymes in aquatic environments focused on chitinase A, which belonged to Class I chitinase mentioned above, since this enzyme was thought to be conserved in Proteobacteria and may be a potential indicator of chitinoclastic ability \[[@CR8]\]. Gene *chiA* was present in 17 species of the genus *Vibrio*, almost all of which were single-copy, except *V. coralliilyticus*, which has two different copies of *chiA*. The phylogenetic tree of 18 *chiA* genes was constructed as shown in Fig. [7](#Fig7){ref-type="fig"}. The topology of the *chiA* tree was very similar to that of the core-genome tree (Fig. [2a](#Fig2){ref-type="fig"}), indicating the *chiA* gene evolved along with the *Vibrio* genus. Intriguingly, one of the two copies of *V. coralliilyticus* was shown to have a relatively distant relationship with other orthologues, suggesting its different evolutionary origin. Furthermore, there were more synonymous (Ks) changes than non-synonymous (Ka) changes in the evolutionary process of *chiA*, which indicated the *chiA* gene had been subjected to significant purifying (negative) selection.Fig. 6Distribution of Subgroups of enzyme Chitinase involved in the 20 *Vibrio* species with complete genomes. Each of the columns corresponds to a chitinase cluster. The number within the box indicates the copy numbers of that gene cluster. Class I contains a family 18 chitinase domain; Class II contains a family 18 chitinase domain with a chitinase insertion domain; Class III contains a family 19 chitinase domain; Class IV contains no significant chitinase domain based on NCBI conserved domain database. Columns belonged to the same class indicate they share similar component patterns of structure domains but failed to be clustered into the same protein familyFig. 7Phylogenetic tree built based on *chiA* carried by 17 complete genomes of *Vibrio*. The phylogenetic tree was obtained using the neighbor-joining method. The number at each node denotes the bootstrap value based on 1000 replicates. The '*V. coralliilyticus* RE98 (2)' branch represents the second copy of *chiA* involved in strain *V. colralliilyticus* RE98

Discussion {#Sec20}
==========

Genomic characters of V. rotiferianus B64D1 and V. mediterranei QT6D1 {#Sec21}
---------------------------------------------------------------------

The specific genome of an organism among a related group represent its unique character. We identified the specific genome of each of the two complete genomes we sequenced and found they occupied many unique orthologues families which absent from the previously genes pool of the complete genomes of vibrios. *V. rotiferianus* B64D1 comprised many unique genes related to cell wall/membrane/envelope biogenesis (COG category M), e.g., six orthologues families were annotated as glycosyltransferase, which may play an important role as the first level of protection against environmental stress. Up to 53 orthologues families in the specific genome of *V. mediterranei* QT6D1 were assigned into the COG category of carbohydrate transport and metabolism (G), among which 6 orthologues families were annotated as the component of ABC-type sugar transport system, which may indicate the potential specific ability for *V. mediterranei* QT6D1 to utilize different kinds of sugar in ambient environments. Besides, 46 orthologues families in the specific genome of *V. mediterranei* QT6D1 were assigned into the COG category of transcription (K), and 31 out of these 46 families were transcriptional regulator. The abundant unique transcriptional regulators present in the genome of *V. mediterranei* QT6D1 could enable it to regulate metabolic processes to adapt complex habitats. Furthermore, it is interesting that *V. mediterranei* QT6D1 carried a putative chitinase gene in its specific genome, which might explain the high activity to degrade chitin showed by the strain in our chitin degradation study.

Genomics and the evolutionary dynamics of *Vibrio* species {#Sec22}
----------------------------------------------------------

The genus *Vibrio* encompasses one of the most diverse group of heterotrophic marine bacteria \[[@CR61], [@CR62]\]. We investigated 20 *Vibrio* species with complete genomes (18 species were from the GenBank and two were from this study) and all the 1582 available draft *Vibrio* genomes in the GenBank. A set of common marine bacterial species which are widely distributed in the ocean was employed to make comparisons with vibrios. The numbers of predicted ORFs in *Vibrio* species differed from \~ 3800 to \~ 6000. Given the close phylogenetic relationship of these bacteria according to the ANI analysis, such a difference suggested substantial gene loss and/or gain in their evolution. Variable gene composition was the primary mode of adaptive evolution in microorganisms, leading to the enhancement of their physiological and ecological capabilities \[[@CR63]\]. We therefore performed ancestral state reconstruction of vibrios and found generally expansion of the orthologues number throughout the whole evolutionary history of this group. Genome expansion evidence has been found in a few other groups of marine bacterioplankton, such as the genus *Glaciecola* \[[@CR64]\], whereas the genome reduction (also known as the genome streamlining) has been reported in groups such as SAR11 \[[@CR65]\] and the genus *Prochlorococcus* \[[@CR66]\]*.* Many putative gene gains belonged to the function of transcription and carbohydrate transport and metabolism, which enabled vibrios to quickly response to environmental changes and use many different kinds of substrates to survive under various conditions. The function of signal transduction mechanisms has experienced the most gene loss, despite it has also experienced more than three times of gene gain in the evolution. The constantly changing of the components of this function may indicate vibrios absorb and discard genes as needed, in response to their environments. The massive genome content changes helped vibrios succeed within a niche, and the process of HGT can promote their evolution. We found that the percentage of alien genes in *Vibrio* genomes was significantly more than most of the other compared marine bacteria, suggesting that alien genes provided vibrios with a genetic basis of diversity and helped them to survive under their respective niche conditions. In addition, the HGT events appear to occur among bacteria in close relationship, as more than 70% alien genes of vibrios were estimated to come from the same class *Gammaproteobacteria.*

Bacteria have to change their genomes to adapt to variable conditions, and greater niche diversity will require larger pan-genomes. Given the varying marine environments in which the species distributed, it was understandable that vibrios had a large amount of pan-genome contents, and this feature was in accordance with previous reports of some other marine bacteria such as the genus *Shewanella* \[[@CR67]\] and *Glaciecola* \[[@CR64]\]. In addition, the pan-genome tree was inferred based on homologous genes distributed in the whole-genome scales. The taxa relationships in the pan-genome tree differed from that in the core-genome tree, because the pan-genome tree focused on those genes that were not present in every genome, while on the contrary the core-genome tree focused on those genes that were present in all the genomes. Such a substantial difference between the two trees indicated that the dispensable genes, which may contain genes obtain from the HGT processes, played a crucial role in vibrios' genomic composition and evolution.

Phylogeny and divergence times of *Vibrio* species {#Sec23}
--------------------------------------------------

*Vibrio* species have considerable biodiversity and extremely close phylogenic relationships within the genus \[[@CR16]\]. 16S rRNA gene sequences could only provide a relatively low discriminating power to elucidate the evolutionary relationship of the group. The availability of complete genomes of *Vibrio* species enables a more definitive analysis of their relationships. Though highly plastic genomes may suffer lots of recombination and HGT events, it has been proved that the phylogenetic signal was still very strong within the core-genome \[[@CR68]\]. We constructed a core-genome tree of the genus *Vibrio* from concatenated proteins of all the homologues shared by every species, and the result showed improved resolution and enormously increased robustness of phylogenetic analyses with nearly all the nodes having the bootstrap value of 100. The availability of a robust phylogeny of the genus *Vibrio* enabled a reliable basis to conduct the following and future analysis.

The species *V. fischeri* might have occurred at 385 MYA, corresponding to the Devonian period (\~ 419.2--358.9 MYA). This was a period that the cephalopods (the order which squid belong to) has been estimated to divergent \[[@CR69]\]. Besides, fish reached substantial diversity during this time, and many ray-finned fishes (Actinopterygii) occurred during the Devonian period \[[@CR70]\]. Based on the phylogenetic tree, *V. fischeri* was a early member of vibrios to emerge, and it is a marine luminous bacterium that lives symbiotically with squids and fish \[[@CR71]\]. Therefore, we could deduce that the rise of vibrios might have coincided with the occurrence of cephalopods and the flourishing of fish.

Environmental adaptation strategies of vibrios {#Sec24}
----------------------------------------------

Free-living bacteria encompass large genomes in most cases \[[@CR72]\], and vibrios occupy relatively large genomes compared to most other common marine bacteria in this study, which enable them to search for resources (chemotaxis), change strategies to exploit patchy nutrients and target metabolites such as complex polysaccharides, which require complex enzyme repertoires.

The capacities of microorganisms to adapt to various environments can be laterally indicated by their growth rates. To the best of our knowledge, vibrios have the fastest growth rate of marine bacteria, and *V. natriegens* has been reported to have a maximal generation time of less than 10 min \[[@CR73]\]. It has been considered that the number of rRNA operons may dictate the rate of ribosomes synthesis and the response to favorable changes in growth environments \[[@CR74]\]. Therefore, the extremely rapid growth rate of vibrios may be explained by the much higher number of 16S rRNA gene copies, which was, on average, more than three times of other compared marine bacteria in this study.

Bacteria have developed protective systems including R-M systems and CRISPR systems to defend external DNA invasion either by plasmids or virus \[[@CR51]\]. Most vibrios in our analysis possessed at least one type of R-M system, and few of them also had CRISPR systems. These protective systems found in *Vibrio* genomes indicate their ability to defend themselves in marine environment and also suggest multiple encounters with foreign DNA.

Chitin metabolism is essential for carbon and nitrogen recycling in marine ecosystems. The capability of vibrios to use chitin provided possibilities for them to form biofilms on the surface of chitin, and the formation of these biofilms could affect many other biological functions \[[@CR75]\]. Therefore, the chitinase family played an essential role in the growth and development of vibrios, which may explain the expansion of this family during the evolutionary process of the genus. We found that the numbers of chitinase genes in the genomes of *Vibrio* species (up to 9 and 16 in the small and large dataset, respectively) are much higher than that in the other compared marine bacteria (only two strains from *Gammaproteobacteria* had eight and two, respectively, while others did not carry chitinase genes), yet only two of which exist in the predicted MRCA of the genus *Vibrio*. As an expansion in size of a particular gene family suggested its importance for the adaptation of the species during evolution \[[@CR76]\], the remarkable expansion of chitinase gene demonstrated the importance of the chitin catabolism process in vibrios.

The chitinolytic cascade in *Vibrio* species was attributed to a bunch of chitin degrading related genes. We found that \~ 99% of *Vibrio* strains contained at least one chitinase coding gene in their genomes in the large dataset. However, a previous chitinase enzyme activity survey found that only \~ 80% (37 out of 47) of *Vibrio* strains had the ability to grow on chitin substrates \[[@CR8]\], which may be because some components required in the complex pathway of chitin-degrading were absent in the genomes. Studies on the nonpathogenic marine bacterium *V. furnissii* showed that chitin utilization was a complex process involving at least three steps: chitin sensing, attachment, and degradation \[[@CR77]\], and it has been demonstrated that the lack of essential chitin degradation elements led to the loss of chitin-degrading activity \[[@CR78]\]. Our analysis showed that the chitinase coding genes suffered HGT and duplication during the evolutionary history of vibrios, which may result in mutations in the sequences and followed by the loss of the chitin-degrading activity. In addition, the selective pressure calculation and phylogenic analysis of the gene *chiA* showed that the gene has undergone significant purifying selective force and may evolve along with most other protein families in vibrios. It could be deduced that there has been evolutionary pressure to make *chiA* conserve the ancestral state, and members of the *chiA* family may evolve independently by silent nucleotide substitution and may be subject to evolution by a birth-and-death process at the DNA level, which implies new genes are generated by gene duplication, and some duplicate genes may maintain in the genome, while others may be deleted or dysfunction through deleterious mutations \[[@CR79]\]. A previous study showed that vibrios tend to acquire new genetic material more easily when they attach to chitin, and this process could help vibrios to grow in a previously hostile niche by slightly changing their genome content \[[@CR80]\]. We noticed that the chitinase genes in *Vibrio* genomes presented random distribution patterns, and does not seem to correlate with phylogenetic relatedness, which may because the HGT or duplication events occurred in the *Vibrio* evolutionary history were randomly. Intriguingly, in the 19 *Vibrio* complete genomes, the two genomes of which does not possess genes encoding chitinase coincide with the two does not have R-M systems, suggesting that the R-M systems might be deleted from the genomes as the systems would be less useful for vibrios which lose chances to intake external DNA. This confirmed the important role of chitin metabolism in vibrios, and the fact that chitinase genes were abundant in *Vibrio* genomes potentially explained the reason why we found an extremely high rate of HGT occurring in vibrios.

Conclusions {#Sec25}
===========

This extensive comparative genomics study based on complete genomes revealed the evolutionary history and ecological adaptation strategies of the widely distributed marine group vibrios. The phylogenetic history of vibrios was inferred according to the core-genome tree, and the extremely gene gain events were found during the evolutionary process from 536 MYA, which may be largely driven by HGT. The steady core-genome and tremendous pan-genome made vibrios comprise necessary genetic materials and flexible accessory genes to be cosmopolitan. Most vibrios had defensive mechanisms including R-M systems and CRISPR systems to resist the invasion of foreign DNA. Furthermore, chitin-degrading related genes were carried in nearly all the *Vibrio* genomes, and the number of chitinase genes (up to 16) have extremely expanded from only two exist in the most recent ancestor of vibrios. Our study provided a comprehensive view of vibrios at genomic level, shedding light on the process that the evolution drives the abundant heterotrophic bacteria to be diversity and adapt to various marine environments.

Additional files
================

 {#Sec26}

Additional file 1:**Table S1.** General features of the 20 complete genomes of *Vibrio* species analyzed. (DOCX 48 kb) Additional file 2:**Table S2.** General features of the 13 marine bacteria complete genomes used to compare with *Vibrio*. (DOCX 34 kb) Additional file 3:**Figure S1.** Pan and core plot of the 20 *Vibrio* species with complete genomes. Accumulation curves for (A) the number of genes in common or (B) total number of genes are plotted. The grey dots represent 10 different random input sequences of genomes, and the extrapolated limiting value for core genes is shown as a dashed line. Equations which can be used to fit the curves are shown above the plot respectively. (PDF 1023 kb) Additional file 4:**Figure S2.** COGs classification of the specific genomes of *V. rotiferianus* B64D1 and *V. mediterranei* QT6D1. Designations of functional categories: \[A\] RNA processing and modification, \[B\] chromatin structure and dynamics, \[C\] energy production and conversion, \[D\] cell cycle control and mitosis, \[E\] amino acid metabolism and transport, \[F\] nucleotide metabolism and transport, \[G\] carbohydrate metabolism and transport, \[H\] coenzyme metabolism, \[I\] lipid metabolism, \[J\] translation, \[K\] transcription, \[L\] replication and repair, \[M\] cell wall/membrane/envelope biogenesis, \[N\] Cell motility, \[O\] post-translational modification, protein turnover, chaperone functions, \[P\] Inorganic ion transport and metabolism, \[Q\] secondary metabolites biosynthesis, transport and catabolism, \[R\] general functional prediction only, \[S\] function unknown, \[T\] signal transduction, \[U\] intracellular trafficking and secretion, \[V\] Defense mechanisms, \[W\] Extracellular structures, \[Y\] nuclear structure, \[Z\] cytoskeleton. (PDF 822 kb) Additional file 5:**Figure S3.** Pan genome tree. The tree was created based on the presence or absence of gene clusters in the 20 complete *Vibrio* genomes. The number at each node denotes the bootstrap value based on 1000 replicates. The color red and yellow are corresponding to the core genome tree, suggesting the discrepancy between core and pan genome trees. (PDF 194 kb) Additional file 6:**Figure S4.** Neighbor-joining phylogenetic tree of the 20 vibrios genomes based on 16S rRNA gene. The number at each node denotes the bootstrap value based on 1000 replicates. *S. denitrificans* OS217 was used as the outgroup. Bar, 0.01 substitutions per site. (PDF 186 kb) Additional file 7:**Figure S5.** Heatmap presentation of pairwise average nucleotide identity of the 19 *Vibrio* species with complete genomes. The genomes are hierarchical clustered according to the values of rows. The clusters present in blue and orange are corresponding to the core genome tree. (PDF 212 kb) Additional file 8:**Figure S6.** Genes related to the chitin-degrading process in 20 vibrios with complete genomes. Each column indicates a chitin metabolism-related gene family, with the family name indicating the predicted function. The number in the box indicates the copy number of that gene family in the corresponding genome. (PDF 282 kb)

ANI

:   Average nucleotide identity

HGT

:   Horizontal gene transfer

MA

:   Marine 2216E agar

ML

:   Maximum-likelihood

MLSA

:   Multilocus sequence analysis

MYA

:   Million years ago

ORF

:   Open reading frame

TCBS

:   Thiosulfate citrate bile salts sucrose

**Electronic supplementary material**

The online version of this article (10.1186/s12864-018-4531-2) contains supplementary material, which is available to authorized users.

The authors thank professor Brian Austin from University of Stirling for critical reading of the manuscript.

Funding {#FPar1}
=======

This work was supported by the National Natural Science Foundation of China through grants 41730530, 41476112, 41521064 and 41506154.

Availability of data and materials {#FPar2}
==================================

The data sets supporting the results of this article are available at the GenBank database under the accession numbers CP018311 to CP018312 and CP018308 to CP018310. All data used for the phylogenetic tree of this article are available in the TreeBase (<http://treebase.org/treebase-web/home.html>) repository, under the accession URL <http://purl.org/phylo/treebase/phylows/study/TB2:S22249>. The materials generated during the current study are available from the corresponding author on reasonable request.

X-HZ designed the study, interpreted the data and thoroughly revised the manuscript. HL and MY performed the bioinformatics analysis and wrote the manuscript. XW isolated the strains and sequenced the genomic DNA and analyzed the genomic data. All authors edited and approved the final version of the manuscript.

Ethics approval and consent to participate {#FPar3}
==========================================

Not applicable, as this study did not involve human or animal subjects, and the field sampling procedures met local guidelines.

Consent for publication {#FPar4}
=======================

Not applicable.

Competing interests {#FPar5}
===================

The authors declare that they have no competing interests.

Publisher's Note {#FPar6}
================

Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.
