Population genetic and evolution analysis of Vibrio isolated from Turkish fish farms

Aquaculture 562 (2023) 738728 Available online 19 August 2022 0044-8486/© 2022 The Authors. Published by Elsevier B.V. This is an open access article under the CC BY-NC license (http://creativecommons.org/licenses/bync/4.0/). Population genetic and evolution analysis of Vibrio isolated from Turkish fish farms Muhammed Duman a , Noemí Buj´ an b , c , * , Soner Altun a , Jesús L. Romalde c , d , Izzet Burcin Saticioglu e , 1 a Department of Aquatic Animal Diseases, Faculty of Veterinary Medicine, Bursa Uludag University, 16059 Bursa, Turkey b Biology Department, Boston College, Chestnut Hill, MA 02467, USA c Department of Microbiology and Parasitology, CIBUS-Faculty of Biology, Universidade de Santiago de Compostela, Santiago de Compostela 15782, Spain d CRETUS, Universidade de Santiago de Compostela, Santiago de Compostela 15782, Spain e Department of Aquatic Animal Diseases, Faculty of Veterinary Medicine, Erciyes University, 38280 Kayseri, Turkey ARTICLE INFO Keywords: Vibrio spp. MLST Population genetics Aquaculture ABSTRACT The genus Vibrio includes important pathogenic species for human and aquatic organisms such as Vibrio cholera, V. parahaemolyticus, V. vulnificus, V. anguillarum or V. harveyi. At present, Vibrionaceae family consists of >190 described species, classified into nine genera. Vibrio are widespread in shellfish, finfish and marine ecosystems and show resourceful ecologies, which recognized as one of the most diverse bacterial groups for illuminating the genome evolution. In the present study, to clarify the relationship among aquatic species in the genus, a multilocus sequence analysis (MLSA) and typing (MLST) approach was applied to characterize 51 Vibrio isolates from Turkish fish farms, 146 strains deposited in the PubMLST database and 59 type strains from GenBank. For all studied isolates (n =256), diversity analysis, population structure, determination of recombination, demographic history and gene flow were performed based on the MLST scheme. Vibrio isolates, subjected to the study, showed a high diversity within the Vibrio population and also genetic interactions into the genus. We found 17 new described sequence types by MLST analysis that were isolated from rainbow trout, sea bream and sea bass in Turkish fish farms, which clearly indicate that the genes underwent recombination frequently. While predominant sequence types were found in the presented study, differences of genotypes need to be evaluated in a disease situation or preventing measurements. The findings about genetic recombination possibly helps to understand differences of Vibrio infections in fish. Furthermore, elucidating of genetic variability within species shed light on providing effective measurements in aquaculture by vaccine production and drug applications. 1. Introduction Vibrio genus comprises ubiquitous bacteria in the aquatic environment, including essential zoonotic pathogens such as Vibrio cholera, V. parahaemolyticus, and V. vulnificus (Farmer et al., 2015; Johnson et al., 2012). In addition to described zoonotic species, the genus Vibrio includes about 145 species and some of them cause devastating mortalities in numerous cultured fish species, primarily V. anguillarum, V. harveyi, V. ordalii, V. alginolyticus (Farmer et al., 2015; Parte et al., 2020; Woo and Bruno, 2011). Besides the known pathogenic species of Vibrio, there are also numerous newly characterized species even though their zoonotic and pathogenic role are still unknown. In many countries, fish are vaccinated each year for protection against pathogenic agents, including V. anguillarum, V. salmonicida, V. alginolyticus and V. vulnificus; some of them are licensed vaccines (Colquhoun and Lillehaug, 2014); still, there are numerous mortality reports in cultured fish species because of Vibrio whether fish were vaccinated or not (Mechri et al., 2017; Sharma and Dube, 2017; Woo and Bruno, 2011). The success of vaccine development deals with obstacles such as highly variable species, highly close but different pathogenic species or misidentified bacteria (Castiblanco and Anaya, 2015). Thus, the accurate identification and exact characterization of the Vibrio spp. play a * Corresponding author at: Department of Microbiology and Parasitology, CIBUS-Faculty of Biology, Universidade de Santiago de Compostela, Santiago de Compostela 15782, Spain. E-mail address: [email protected] (N. Buj´ an). 1 Present adress: Bursa Uludag University, Vocational School of Keles, 16740, Bursa, Turkey. Contents lists available at ScienceDirect Aquaculture journal homepage: www.elsevier.com/locate/aquaculture https://doi.org/10.1016/j.aquaculture.2022.738728 Received 7 February 2022; Received in revised form 30 May 2022; Accepted 10 August 2022 Aquaculture 562 (2023) 738728 2 key role in protecting fish. DNA-based molecular methods have become more attractive and widely acceptable because of high discriminatory power (Mulet et al., 2010) and one of the most used molecular methods is 16S rRNA gene sequencing for identification for an agent. Together with molecular identification, genetic typing and determination of genetic variabilities lead up effective monitoring pathogen transmission especially Repetitive Extragenic Palindromic sequence (REP-PCR), Random Amplified Polymorphic DNA (RAPD), DNA-DNA hybridization technique (dDDH) and Multilocus sequence typing (MLST) (Mulet et al., 2010; Saticioglu et al., 2018). MLST method allows discrimination of the species of the genus and determines low intraspecific genetic variations by using different housekeeping genes at the same time. To better understand genetic variations in bacteria, a public database (PubMLST), was launched and about 739,319 isolates and >20 million alleles were added recently (PubMLST, 2014). Previously, Vibrio spp. MLST database was constructed with different Vibrio spp., which were isolated from live marine animals collected from Italy in 2007 and reference strains, including 160 Sequence type (ST) and >440 alleles (https://pubmlst. org/organisms/Vibrio-spp) but the evolutionary relationship among the strains was not studied. One of the MLST advantages is the possibility of progressively enriching the database by adding new isolates; therefore, we added 51 new Vibrio isolates from different hosts isolated between 2014 and 2018 in Turkey. The aims of this study are to clarify the identification of 51 Turkish Vibrio isolates by phylogenetic analysis and to elucidate the evolution of genus Vibrio focusing on the clades Anguillarum, Cholera, Harveyi, Mediterranei, Orientalis, Ponticus, Scophtalmi, Splendidus and Vulnificus. 2. Materials and methods 2.1. Bacterial isolates and phenotypic characterization A total of 51 Vibrio isolates (Table 1) were isolated from mainly seabream (Sparus aurata), sea bass (Dicentrarchus labrax) cultures in offshore farms. Additionally, isolates were collected from rainbow trout (Oncorhynchus mykiss) farms using conventional flow-through systems with non-disinfected surface water and spring water located in five different regions of Turkey. The farms were selected due to their high production capacity (the farms that provided samples supplied a total of at least 80% of the production capacity of Turkey). The bacterial isolates were collected during a fish health surveillance program from internal organs of healthy (no external lesions) and diseased fish (moribund and with external lesions) of different species, weight, age and regions. Samples from external lesions, internal organs such as kidney, liver, spleen, and ascites fluid were collected monthly over four years from at least 30 farms (2014–2018). Samples were cultured on blood agar (BA; 5% sheep blood added), brain heart infusion agar (BHIA, Merck, USA), tryptic soy agar (TSA, Merck, USA) and marine agar (MA, Merck, USA). Sampling was carried out following the guidelines for the diagnosis of fish diseases and in consideration of international and national guidelines for animal welfare (OIE, 2015). The biochemical characteristics of all isolates were determined using conventional microbial tests such as colony morphology, color upon Gram staining, bacterial and swarming motility, oxidase and catalase activities, oxidative fermentation (O/F), growth on Thiosulfate-citrate- bile salts-sucrose agar (TCBS, Merck, USA), susceptibility to Vibriostatic agent (0/129, 10 μ g and 150 μ g) different temperature and salinities, hemolysis on blood agar (5% sheep blood added), lysis on gelatin, and Tween 20 and 80 (Farmer et al., 2015). Moreover, the sequences of 146 strains deposited in PubMLST (https://pubmlst.org/organisms/vibrio-spp), as well as the sequences of type strains (n =59) housed in the GenBank database, representative of the clades identified, were used (Supplementary Table S1). It is interesting to note that some strains deposited in the database were excluded for different reasons: i) It is not Vibrio (ES114, ATCC 43305, V63 and V78), ii) It has an extra codon and some software does not work with spaces (Vi_20, Vi_51, Vi_54, Vi_60, Vi_62, Vi_73, Vi_9a and Vi_16a), iii) we used only type strains to known species (40B, 775, ATCC 43996, LMG 23867, NCTC 11218, ATCC 39315, RIMD 2210633 and VM603) and iv) It is duplicated (BAA-1116 and LMG 19703 T). 2.2. DNA extraction, PCR and sequence analysis DNA was extracted with spin column filtration kits according to the manufacturer’s instructions (QIAamp DNA mini kit, 51,306, Hilden, Germany). The amount and purity of DNA in each sample were measured at wavelengths of 260 nm and 260/280 nm with a spectrophotometer (Multiskan Go, Thermo, Vantaa, Finland). The genomic DNA of the bacterial isolates was extracted using a QIAamp DNA mini kit (Qiagen, Hilden, Germany) according to the manufacturer’s instructions. All DNA samples were optimized according to concentration and purity before molecular analysis. All Turkish isolates included in this study were previously identified as Vibrio sp. by 16S rRNA gene sequencing (data not shown). Following the previous MLST database, the housekeeping genes used were atpA (ATP synthase α subunit), gyrB (DNA gyrase β subunit), pryH (Uridylate kinase) and recA (Recombinase A) (https://pubmlst.org/organisms/vibrio-spp; J¨ ackel et al., 2020; Preheim et al., 2011; Rahman et al., 2014; Schwartz et al., 2017). Polymerase chain reaction (PCR) was performed in a total volume of 25 μ l containing 2.5 μ l 10×PCR Buffer, 5 μ l 5×Q-Solution, 0.5 μ l 0.2 mM dNTP mix (10 mM of each), 37.5 pmol of each primer (Supplementary Table S2), 1 μ l 25 mM MgCl 2 , 0.125 μ l HotStarTaq DNA Polymerase (Qiagen, Venlo, Netherlands), 10 ng DNA template, and DNase/RNase-Free Distilled Water. All PCR products were sent for sequence analysis to Macrogen Korea (Republic of Korea). Sequences obtained for atpA, gyrB, pyrH and recA gene regions of 51 Vibrio isolates were deposited in the pubmlst database (https://pubmlst.org/o rganisms/vibrio-spp), allele numbers and sequence types (STs) were generated by pubmlst database. 2.3. Diversity analyses and population structure Both forward and reverse sequences for each strain were trimmed, aligned and analyzed using BioNumerics software (version 7.6, Applied Maths, Sint-Martens-Latem, Belgium). For each MLST locus, the sequences obtained for all isolates were compared, and each unique sequence (differing from the others at one or more nucleotides) was assigned with a new allele number. Each isolate was unambiguously defined by an allele profile or sequence type (ST) derived from the combination of numbers corresponding to the alleles at each of the loci analyzed. The same ST was assigned to multiple strains when they shared the same allelic profile. Allelic profiles were used for subsequent analysis. A descriptive diversity analysis was performed for the individual genes and concatenated sequences (atpA, gyrB, pyrH and recA). The mean G +C content of the DNA along with the d N /d S ratios (where d N is the number of nonsynonymous substitutions per nonsynonymous site, and d S is the number of synonymous substitutions per synonymous site), number of polymorphic sites, indexes of association (I A and I S A ) (Maynard et al., 1993) were calculated using START2 (http://pubmlst. org/software/analysis/start2/). Nucleotide diversity per site ( π ) was estimated using DnaSP v5 (Librado and Rozas, 2009). The clonal structure of the population was established using Phyloviz software by goeBURST 1.2.1 algorithm (Francisco et al., 2009). ClonalFrame v1.1 (Didelot and Falush, 2007) was used to construct a dot graph using a 50% majority rule consensus tree with ten runs of MCMC (Markov Chain Monte Carlo). For phylogenetic analysis, concatenated sequences were aligned and trees were performed for the whole collection with the program MEGA version X (Kumar et al., 2016) by Neighbor-joining (NJ) M. Duman et al. Aquaculture 562 (2023) 738728 3 Table 1 Origin and MLST data of Vibrio isolates employed in this study. Bacteria code Fish species Isolation date Isolated city Isolated region Fish lesion MLSA result Sequence type (ST) V1 live marine animals ND ND ND ND V. alginolyticus 130 V10 Dicentrarchus labrax Jul-2016 Mugla Aegean region ND V. alginolyticus 134 V13 Dicentrarchus labrax Jan-2016 Mugla Aegean region ND V. alginolyticus 136 V16 Dicentrarchus labrax Apr-2017 Mugla Aegean region ND V. alginolyticus 137 V19 Dicentrarchus labrax Dec-2016 Kayseri Central Anatolia ND V. alginolyticus 134 V20 Sparus aurata Feb-2016 Kayseri Central Anatolia ND V. alginolyticus 137 V22 Dicentrarchus labrax Jun-2017 Elazıg Eastern anatolia ND V. alginolyticus 136 V23 Sparus aurata May-2017 Elazıg Eastern anatolia ND V. alginolyticus 134 V26 Dicentrarchus labrax Feb-2017 Mugla Aegean region ND V. alginolyticus 134 V69 Oncorhynchus mykiss May-2011 Rize Black sea region darkening in color, haemorhagia, fin erosion V. alginolyticus 134 V72 Dicentrarchus labrax Jan-2016 Mugla Aegean region No lesion V. alginolyticus 134 V73 Dicentrarchus labrax Jul-2017 ˙ Izmir Aegean region Anemia in liver V. alginolyticus 140 V74 Dicentrarchus labrax Jul-2017 ˙ Izmir Aegean region Anemia in liver V. alginolyticus 140 V75 Oncorhynchus mykiss Feb-2015 Mugla Aegean region No lesion V. alginolyticus 141 V103 Dicentrarchus labrax Jun-2015 Canakkale Maramara region High mortality V. alginolyticus 142 V5 live marine animals Jun-2008 Norway ND ND V. anguillarum 131 V9 Dicentrarchus labrax Jun-2016 Mugla Aegean region ND V. anguillarum 133 V12 Dicentrarchus labrax Jan-2016 Mugla Aegean region ND V. anguillarum 135 V15 Oncorhynchus mykiss May-2017 Mugla Aegean region ND V. anguillarum 133 V18 Oncorhynchus mykiss Dec-2016 Tunceli Eastern anatolia ND V. anguillarum 133 V21 Sparus aurata Jun-2017 Elazıg Eastern anatolia ND V. anguillarum 133 V25 Sparus aurata May-2017 Mugla Aegean region ND V. anguillarum 135 V28 Sparus aurata Jan-2017 Mugla Aegean region ND V. anguillarum 133 V32 Sparus aurata Jan-2017 Mugla Aegean region ND V. anguillarum 133 V76 Oncorhynchus mykiss Aug-2013 Mugla Aegean region Moribund fish V. anguillarum 133 V80 Oncorhynchus mykiss Feb-2015 Mugla Aegean region No lesion V. anguillarum 133 V82 Oncorhynchus mykiss Jan-2015 Mugla Aegean region No lesion V. anguillarum 133 V85 Oncorhynchus mykiss Feb-2011 Rize Black sea region darkening in color, haemorhagia, fin erosion V. anguillarum 133 V88 Oncorhynchus mykiss Feb-2014 Mugla Aegean region Darkening in color V. anguillarum 133 V91 Oncorhynchus mykiss Feb-2014 Mugla Aegean region Darkening in color V. anguillarum 133 V92 Oncorhynchus mykiss Feb-2014 Mugla Aegean region Darkening in color V. anguillarum 133 V95 Oncorhynchus mykiss Nov-2013 Mugla Aegean region Haemorhage in fin base V. anguillarum 133 V96 Oncorhynchus mykiss Nov-2013 Mugla Aegean region Haemorhage in fin base V. anguillarum 133 V98 Dicentrarchus labrax Jun-2015 Canakkale Maramara region High mortality V. anguillarum 133 V99 Oncorhynchus mykiss Apr-2015 Mugla Aegean region No lesion V. anguillarum 133 V100 Oncorhynchus mykiss Apr-2015 Mugla Aegean region No lesion V. anguillarum 133 V101 Dicentrarchus labrax Jun-2015 Canakkale Maramara region High mortality V. anguillarum 133 V102 Dicentrarchus labrax Jun-2015 Canakkale Maramara region High mortality V. anguillarum 133 V104 Dicentrarchus labrax Jun-2015 Canakkale Maramara region High mortality V. anguillarum 133 V109 Feb-2017 Mugla Aegean region No lesion V. anguillarum 133 (continued on next page) M. Duman et al. Aquaculture 562 (2023) 738728 4 method using Jukes-Cantor method and by Maximum Likelihood (ML) algorithm using General Time Reverse model (GTR +G +I), with 1000 bootstrap pseudoreplicates in both cases. For the ML reconstruction, optimal evolution models were estimated from nucleotide data using MEGA X software considering 24 substitution types. The best model was selected using the Bayesian Information Criterion (BIC). A network constructed by the Median Joining method using the program NETWORK 4.1 (Fluxus-engineering) (Bandelt et al., 1999) was displayed to establish the relationship among haplotypes. 2.4. Recombination analyses The phi test for recombination, based on individual loci from the whole strain collection, was performed using SplitsTree v4 (Huson and Bryant, 2006). The P-value (a value of P <0.05 was considered significant herein) indicated the DNA regions exhibiting the strongest evidence of mosaicism. By DnaSP v5 (Librado and Rozas, 2009), the minimal number of recombination events (R min ) were calculated. Structure software (Falush et al., 2007; Falush et al., 2003; Hubisz et al., 2009; Pritchard et al., 2000) was used to identify groups with distinct allele frequencies and detect strains carrying foreign DNA among isolates. The number of clusters (K) was set from 2 to 39 and all runs were replicated five times. The best model of probability, estimated using a burn-in period of 100,000 interactions and an MCMC =50,000, was at twenty-eight groups (K =28). The best model of probability (K = 28) was estimated by webtool Structure Harvester (Earl and VonHoldt, 2012). 2.5. Demographic history and gene flow DnaSP5 (Librado and Rozas, 2009) was employed to calculate the average number of nucleotides per site ( π ) and the number of segregating sites (θ). In addition, the Tajim´ as test (D) (Tajima, 1989) and Fu and Li statistics (F* and D*) (Fu and Li, 1993) were used to confirm if the mutant alleles are selectively neutral (Kimura, 1983). The haplotype structure based on haplotype frequency distribution was evaluated by Fús (F’S) (Fu, 1997) and Strobeck’s (S) (Strobeck, 1987) statistics. Ramos-Onsis R2 test (Ramos-Onsins and Rozas, 2002) assessed the differences between the number of singleton mutations and the average number of nucleotide differences in the population and was used to detect a possible demographic expansion. 3. Results 3.1. Isolation, phenotypic characterization and identification The diseased fish affected by the strains isolated in this study showed darkening in color, haemorrhagia in internal organs, skin and fin base, fin erosion, anemia in the liver and exophthalmia. Isolates were Gramnegative, motile, oxidase and catalase-positive, glucose fermentative, hydrolyzed to gelatin, tween 20 and tween 80 and sensitive to vibriostatic agent (Supplementary Table S3). All Turkish isolates included in this study were primarily identified as Vibrio sp. by 16S rRNA gene sequencing (data not shown), being presumptively assigned to V. anguillarum, V. alginolyticus, V. harveyi, V. diabolicus, V. fluvialis and V. parahaemolyticus (Supplementary Fig. S1). 3.2. Genetic diversity, clustering and phylogenetic analysis We aligned the sequences for each strain at each MLST locus and confirmed that neither indels nor premature stop codons were present in any of the sequences. The DNA fragments used for the analysis are between 570 bp (gyrB) to 462 (recA) (Supplementary Table S2), and the concatenated sequences produced a 2022-bp fragment. The locus with more polymorphic sites was gyrB (131 biallelic, 44 triallelic and 71 tetrallelic sites), in contrast with the atpA (108 biallelic, 45 triallelic and 5 tetrallelic sites). The most frequent allele detected in the four genes were allele 77 in atpA (59%), allele 107 in gyrB (66%), allele 84 in pyrH (63%) and allele 105 in recA (65.6%). The G +C content means of these genes varied from 45.68% to 52.19% (Table 2). The nucleotide diversity per site ( π ) and the number of segregation sites (θ) for the concatenated sequences were 0.12637 and 0.12509, respectively. Values of these estimators for each individual gene in the whole population are presented in Table 3. The d N /d S ratio is a common indicator of selection pressure, with d N /d S >1 indicating positive selection and d N /d S <1 indicating purifying selection for the gene sequence tested. For all loci, the ratio d N /d S was less than one (Table 3). This indicated that all loci were under strong purifying selection, which is optimal for MLST loci. The concatenated sequences were constituted by 201 different alleles or STs (Table 3), 76 of which were sinapomorfic (shared and derived). Although we avoided using some strains from the original MLST, all old STs were checked with the aim to recognize new STs. Among the 51 Turkish Vibrio strains which were isolated in the present study, we determined 7 new STs (103, 105, 106, 107, 108, 109, 110) (Table 2). Table 1 (continued) Bacteria code Fish species Isolation date Isolated city Isolated region Fish lesion MLSA result Sequence type (ST) Dicentrarchus labrax P143 Dicentrarchus labrax Jan-2017 Antalya Mediterranean region Darkening in color, exophthalmos V. anguillarum 133 V145 Dicentrarchus labrax 2014 Mu˘ gla Aegean region No lesion V. anguillarum 133 V146 Dicentrarchus labrax 2014 Mugla Aegean region No lesion V. anguillarum 133 V24 Dicentrarchus labrax May-2016 Elazıg Eastern anatolia ND V. diabolicus 139 V7 live marine animals Apr-2007 Bangladesh Decca ND V. fluvialis 132 V17 Dicentrarchus labrax Aug-2017 Mugla Aegean region ND V. harveyi 138 V110 Sparus aurata Feb-2018 Mugla Aegean region No lesion V. harveyi 143 V111 Sparus aurata Feb-2017 Mugla Aegean region No lesion V. harveyi 144 V147 Sparus aurata 2017 Mugla Aegean region No lesion V. harveyi 145 V148 Sparus aurata 2017 Mugla Aegean region No lesion V. harveyi 146 V6 live marine animals ND ND ND ND V. parahaemolyticus 117 ND: No data, table order was construed according to the alphabetical order of species. M. Duman et al. Aquaculture 562 (2023) 738728 5 Table 2 Sequence type and G +C content of Vibrio isolates. Isolate ST atpA gyrB pyrH recA Mean %GC Vi_2 1 48,47 50,35 49,3 48,05 49,04 Vi_9 2 48,67 47,89 47,9 48,05 48,13 Vi_10 3 47,24 49,3 48,9 45,67 47,78 Vi_11 4 47,03 49,3 49,5 45,89 47,93 Vi_12 3 47,24 49,3 48,9 45,67 47,78 Vi_13 5 47,03 49,12 49,7 46,32 48,04 Vi_14 6 47,03 48,77 49,5 45,89 47,8 Vi_16 7 47,03 48,6 46,51 50 48,03 Vi_18 8 47,85 45,09 46,71 48,27 46,98 Vi_21 9 48,26 51,93 52,69 50,65 50,88 Vi_22 10 46,63 48,77 49,5 45,45 47,59 Vi_23 11 47,44 49,3 49,3 45,89 47,98 Vi_24 12 47,44 49,12 49,5 45,89 47,99 Vi_25 13 47,03 49,12 49,5 46,1 47,94 Vi_26 14 47,44 46,32 47,7 46,54 47 Vi_27d 15 47,65 47,02 49,5 49,35 48,38 Vi_28 16 48,67 49,12 49,5 48,7 49 Vi_29 17 47,44 47,89 48,1 48,7 48,04 Vi_31 18 47,65 47,02 49,5 48,92 48,27 Vi_32 19 47,44 47,72 48,1 48,92 48,05 Vi_33 20 47,85 47,89 49,5 48,7 48,49 Vi_34 21 47,85 47,37 49,3 48,48 48,25 Vi_35 22 47,85 47,72 49,3 49,13 48,5 Vi_36 23 47,44 47,72 47,9 48,7 47,94 Vi_37 24 47,65 46,84 49,5 48,7 48,17 Vi_38 25 47,85 48,42 48,5 46,54 47,83 Vi_39 26 47,85 46,84 49,3 48,92 48,23 Vi_40d 27 47,85 48,77 49,5 48,48 48,65 Vi_41 28 47,85 47,19 48,3 47,62 47,74 Vi_42 29 49,69 46,32 48,7 48,27 48,25 Vi_43 30 48,06 47,37 49,3 48,7 48,36 Vi_44 31 48,06 47,37 49,3 48,92 48,41 Vi_45 32 47,85 46,32 48,7 48,05 47,73 Vi_46 33 47,85 46,32 48,7 48,05 47,73 Vi_47 5 47,03 49,12 49,7 46,32 48,04 Vi_48 32 47,85 46,32 48,7 48,05 47,73 Vi_49 32 47,85 46,32 48,7 48,05 47,73 Vi_50 34 47,85 46,67 49,1 48,27 47,97 Vi_52 35 46,42 48,07 48,3 47,4 47,55 Vi_53 36 46,42 47,02 47,5 47,4 47,09 Vi_55 37 47,24 48,6 49,3 46,32 47,86 Vi_56 38 48,47 47,19 48,5 47,84 48 Vi_57 39 47,24 49,3 49,5 46,1 48,04 Vi_58 40 47,65 48,95 49,3 46,32 48,05 Vi_59 3 47,24 49,3 48,9 45,67 47,78 Vi_61 41 46,83 47,19 47,9 47,4 47,33 Vi_63 42 47,65 47,72 47,9 48,92 48,05 Vi_64 43 47,65 48,95 49,5 48,7 48,7 Vi_66 43 47,65 48,95 49,5 48,7 48,7 Vi_67 44 47,44 48,42 47,5 48,48 47,96 Vi_68 45 47,65 49,12 48,1 47,84 48,18 Vi_69 46 47,85 46,32 48,7 48,27 47,78 Vi_70 47 47,44 47,89 49,3 48,27 48,23 Vi_71 48 47,85 49,12 48,5 46,54 48 Vi_72 49 47,03 49,12 49,5 45,89 47,89 Vi_74 32 47,85 46,32 48,7 48,05 47,73 Vi_79 50 47,03 46,84 47,31 47,19 47,09 Vi_80 51 48,06 47,37 49,3 48,7 48,36 Vi_81 52 47,65 49,3 48,3 46,32 47,89 Vi_1a 53 47,44 47,37 49,5 49,35 48,42 Vi_2a 54 45,6 45,44 47,9 47,4 46,59 Vi_3a 55 47,65 49,47 48,5 47,62 48,31 Vi_4a 56 47,85 49,12 48,1 47,4 48,12 Vi_5a 57 47,65 47,54 49,3 48,7 48,3 Vi_6a 58 48,26 48,77 48,1 44,59 47,43 Vi_7a 59 47,65 48,42 48,3 48,27 48,16 Vi_8a 60 48,06 47,72 49,3 48,92 48,5 Vi_10a 29 49,69 46,32 48,7 48,27 48,25 Vi_11a 61 48,26 47,19 48,9 47,84 48,05 Vi_12a 62 48,06 46,67 48,9 48,05 47,92 Vi_13a 63 47,44 47,02 49,3 48,48 48,06 Vi_14a 64 47,44 47,19 47,9 48,05 47,65 Vi_15a 65 47,03 48,6 46,71 50 48,08 Table 2 (continued) Isolate ST atpA gyrB pyrH recA Mean %GC Vi_17a 66 47,65 47,02 49,7 49,35 48,43 Vi_18a 67 47,65 48,95 48,7 46,54 47,96 Vi_19a 68 47,44 47,54 48,7 48,05 47,94 Vi_20a 69 47,65 46,84 49,5 48,7 48,17 Vi_21a 70 47,65 48,77 48,5 48,05 48,24 Vi_22a 71 47,85 48,42 48,1 46,54 47,73 Vi_23a 65 47,03 48,6 46,71 50 48,08 Vi_24a 72 47,03 50,53 46,71 49,57 48,46 Vi_25a 70 47,65 48,77 48,5 48,05 48,24 Vi_26a 67 47,65 48,95 48,7 46,54 47,96 Vi_27a 73 48,06 48,42 48,7 46,32 47,88 Vi_28a 74 47,85 46,84 49,5 48,7 48,22 Vi_29a 75 48,26 48,77 48,1 46,32 47,86 Vi_30a 76 47,85 47,19 49,5 48,7 48,31 Vi_31a 77 47,03 48,6 46,51 49,78 47,98 Vi_32a 78 48,26 47,37 48,9 48,27 48,2 Vi_33a 79 47,44 47,89 48,1 48,7 48,04 Vi_34a 19 47,44 47,72 48,1 48,92 48,05 Vi_35a 80 47,44 47,54 48,3 48,48 47,94 Vi_36a 81 47,44 47,54 48,1 48,48 47,89 Vi_37a 82 47,85 48,95 48,1 46,97 47,97 Vi_38a 83 46,63 47,37 47,5 47,62 47,28 Vi_39a 84 48,26 47,19 48,7 47,84 48 Vi_40a 85 47,85 47,72 47,9 48,7 48,04 Vi_41a 86 48,26 47,54 48,3 48,7 48,2 Vi_42a 87 47,65 47,19 47,7 47,4 47,49 Vi_43a 88 47,44 47,72 48,1 48,48 47,94 Vi_44a 89 47,85 49,12 48,1 47,4 48,12 Vi_45a 90 47,85 48,77 49,5 48,05 48,54 Vi_46a 91 47,65 49,3 48,5 46,75 48,05 Vi_47a 92 47,44 48,95 48,1 46,97 47,87 Vi_48a 93 47,44 47,37 47,9 48,48 47,8 Vi_49a 94 47,65 49,3 46,51 46,97 47,61 Vi_50a 95 47,44 48,07 47,5 46,1 47,28 Vi_51a 96 46,83 47,72 47,5 45,02 46,77 Vi_52a 97 47,85 48,6 48,5 46,54 47,87 Vi_53a 98 47,85 50,35 47,11 50,87 49,04 Vi_54a 99 48,06 47,19 48,1 48,27 47,91 Vi_55a 100 48,06 48,77 49,3 48,27 48,6 Vi_56a 101 47,65 47,19 47,9 48,48 47,81 Vi_57a 102 47,03 48,6 46,71 50 48,08 Vi_58a 103 47,85 49,3 46,31 47,84 47,82 Vi_59a 104 47,24 48,95 46,11 47,19 47,37 Vi_60a 105 47,44 47,54 47,31 46,1 47,1 Vi_62a 106 48,47 47,19 49,5 47,62 48,19 Vi_63a 107 47,65 49,3 46,91 46,97 47,71 Vi_64a 108 47,03 48,77 45,91 46,97 47,17 Vi_65a 109 47,24 49,12 46,11 47,19 47,41 Vi_67a 110 47,24 49,12 46,11 50 48,12 Vi_68a 111 51,12 52,46 51,5 53,68 52,19 Vi_69a 112 47,85 48,42 48,3 48,27 48,21 Vi_70a 113 46,63 47,02 48,7 45,89 47,06 Vi_71a 114 47,24 49,12 45,91 47,19 47,36 Vi_72a 115 47,24 47,72 47,7 47,4 47,52 Vi_73a 115 47,24 47,72 47,7 47,4 47,52 Vi_74a 116 47,85 48,77 48,3 47,62 48,14 Vi_75a 117 47,85 48,95 48,3 46,97 48,02 Vi_76a 117 47,85 48,95 48,3 46,97 48,02 Vi_77a 117 47,85 48,95 48,3 46,97 48,02 Vi_78a 118 48,26 50,35 49,5 48,92 49,26 Vi_79a 119 46,63 47,19 48,5 45,89 47,05 Vi_80a 120 47,85 47,02 48,1 48,27 47,81 Vi_81a 121 47,24 49,3 46,31 47,19 47,51 Vi_1bd 117 47,85 48,95 48,3 46,97 48,02 Vi_2b 122 47,44 47,37 47,9 45,45 47,04 Vi_3b 123 47,44 47,37 47,11 45,89 46,95 Vi_4b 124 47,65 48,77 49,5 48,27 48,55 Vi_5b 125 47,85 49,47 47,5 47,84 48,17 Vi_6b 126 47,44 47,54 47,11 46,1 47,05 Vi_7b 127 47,44 47,54 47,31 45,89 47,05 Vi_8b 87 47,65 47,19 47,7 47,4 47,49 Vi_9b 128 47,24 48,95 45,71 47,4 47,32 Vi_10b 129 47,65 47,89 49,1 48,48 48,28 V1 130 47,44 47,54 47,9 48,7 47,9 V5 131 47,24 48,77 49,5 45,67 47,8 (continued on next page) M. Duman et al. Aquaculture 562 (2023) 738728 6 Moreover, other new 43 STs were described in this work, composed of some strains employed in the analysis. The most heterogeneous species cluster in V. alginolyticus/diabolicus, including 32 different ST and commonly collected from live marine animals. goeBURST analysis (Fig. 1) revealed 22 clonal complexes (CC, STs related by SLV or/and DLV) and 1202 singletons (STs). Some of these CC Table 2 (continued) Isolate ST atpA gyrB pyrH recA Mean %GC V6 117 47,85 48,95 48,3 46,97 48,02 V7 132 48,88 51,4 51,3 52,81 51,1 V9 133 47,24 49,47 49,1 46,1 47,98 V10 134 47,65 47,54 48,1 48,7 48 V12 135 47,24 49,12 49,5 46,1 47,99 V13 136 47,44 47,37 48,3 48,48 47,9 V15 133 47,24 49,47 49,1 46,1 47,98 V16 137 47,44 47,54 48,1 48,48 47,89 V17 138 47,85 48,6 49,5 48,27 48,55 V18 133 47,24 49,47 49,1 46,1 47,98 V19 134 47,65 47,54 48,1 48,7 48 V20 137 47,44 47,54 48,1 48,48 47,89 V21 133 47,24 49,47 49,1 46,1 47,98 V22 136 47,44 47,37 48,3 48,48 47,9 V23 134 47,65 47,54 48,1 48,7 48 V24 139 47,85 47,02 49,5 49,13 48,38 V25 135 47,24 49,12 49,5 46,1 47,99 V26 134 47,65 47,54 48,1 48,7 48 V28 133 47,24 49,47 49,1 46,1 47,98 V32 133 47,24 49,47 49,1 46,1 47,98 V69 134 47,65 47,54 48,1 48,7 48 V72 134 47,65 47,54 48,1 48,7 48 V73 140 47,44 47,54 47,9 48,05 47,74 V74 140 47,44 47,54 47,9 48,05 47,74 V75 141 58,49 47,54 47,9 48,05 50,5 V76 133 47,24 49,47 49,1 46,1 47,98 V80 133 47,24 49,47 49,1 46,1 47,98 V82 133 47,24 49,47 49,1 46,1 47,98 V85 133 47,24 49,47 49,1 46,1 47,98 V88 133 47,24 49,47 49,1 46,1 47,98 V91 133 47,24 49,47 49,1 46,1 47,98 V92 133 47,24 49,47 49,1 46,1 47,98 V95 133 47,24 49,47 49,1 46,1 47,98 V96 133 47,24 49,47 49,1 46,1 47,98 V98 133 47,24 49,47 49,1 46,1 47,98 V99 133 47,24 49,47 49,1 46,1 47,98 V100 133 47,24 49,47 49,1 46,1 47,98 V101 133 47,24 49,47 49,1 46,1 47,98 V102 133 47,24 49,47 49,1 46,1 47,98 V103 142 47,24 47,19 48,1 46,1 47,16 V104 133 47,24 49,47 49,1 46,1 47,98 V109 133 47,24 49,47 49,1 46,1 47,98 V110 143 47,85 48,6 49,5 48,27 48,55 V111 144 47,65 48,6 49,5 48,27 48,5 P143 133 47,24 49,47 49,1 46,1 47,98 V145 133 47,24 49,47 49,1 46,1 47,98 V146 133 47,24 49,47 49,1 46,1 47,98 V147 145 47,24 48,6 49,5 48,27 48,4 V148 146 47,24 48,6 49,5 48,05 48,35 V alfacsensis CAIM 1831T 147 48,47 47,37 47,7 48,48 48,01 V alginolyticus ATCC 17749T 130 47,44 47,54 47,9 48,7 47,9 V anguillarum NCIMB 6T 131 47,24 48,77 49,5 45,67 47,8 V artabrorum Vb 118T 148 47,24 47,37 49,3 48,05 47,99 V atlanticus Vb 1111T 149 47,44 46,84 46,11 46,97 46,84 V azureus LC2 ¡005T 150 47,44 47,37 47,31 46,75 47,22 V barjaei 3062T 151 47,44 48,42 46,71 50 48,14 V bivalvicida 605T 152 47,03 49,47 48,7 45,67 47,72 V brasiliensis LMG 20546T 153 46,63 48,07 48,3 46,75 47,44 V campbellii CECT 519T 154 48,06 49,12 48,5 48,92 48,65 V celticus CECT 7224T 155 46,83 46,84 47,9 47,62 47,3 V chagasii LMG 21353T 156 50,31 47,54 47,5 46,1 47,86 V cholerae CECT 514T 157 50,92 51,58 50,1 52,16 51,19 V cidicii 2756-81T 158 48,26 51,58 52,89 51,52 51,06 V cincinnatiensis NCTC 12012T 159 46,22 48,77 50,9 49,35 48,81 V crassostreae LMG 22240T 160 48,67 46,67 48,3 48,27 47,98 V cyclitrophicus LMG 21359T 161 47,24 46,49 47,9 46,32 46,99 162 47,85 47,02 49,3 49,13 48,33 Table 2 (continued) Isolate ST atpA gyrB pyrH recA Mean %GC V diabolicus CNCM I- 1629T V fluvialis NCTC 11327T 163 48,67 51,4 51,3 52,81 51,05 V fortis LMG 21557T 164 47,65 49,47 47,5 47,62 48,06 V furnissii NCTC 13120T 165 51,12 52,11 51,5 53,46 52,05 V gallaecicus CECT 7244T 166 46,01 46,32 47,9 43,94 46,04 V gigantis LGP 13T 167 48,67 46,67 48,5 47,19 47,76 V harveyi ATCC 14126T 168 47,65 46,84 48,7 48,05 47,81 V harveyi 1DA3 169 48,26 48,77 49,1 48,05 48,55 V hepatarius LMG 20362T 170 47,44 50,88 47,31 48,27 48,47 V hyugaensis 090810aT 171 48,06 47,02 49,3 48,27 48,16 V ichthyoenteri ATCC 700023T 172 48,06 47,19 46,91 46,32 47,12 V jasicida CECT 7692T 173 48,06 46,67 49,1 48,05 47,97 V kanaloae LMG 20539T 174 47,03 46,84 46,71 47,4 47 V lentus CECT 5110T 175 47,44 46,49 47,5 47,62 47,26 V maritimus CAIM 1455T 176 47,65 50,88 44,71 50,87 48,53 V mediterranei CECT 621T 177 47,24 49,3 45,91 47,19 47,41 V metoecus OP3HT 178 50,51 49,47 49,9 52,81 50,67 V metschnikovii NCTC 8443T 179 46,42 48,6 49,3 49,13 48,36 V mimicus LMG7896T 180 49,08 51,05 50,1 50,87 50,27 V mytili CAIM 528T 181 47,44 47,72 49,3 48,7 48,29 V natriegens LMG 10935T 182 47,44 47,72 47,9 50,22 48,32 V navarrensis CIP 103381T 183 48,26 52,46 54,49 52,81 52,01 V ordalii ATCC 33509T 184 48,06 48,42 49,7 45,45 47,91 V orientalis CECT 629T 185 49,49 48,42 48,3 49,13 48,84 V owensii LMG 25443T 186 49,28 47,02 48,5 47,84 48,16 V panuliri CAIM 703T 187 46,83 47,02 48,5 46,1 47,11 V parahaemolyticus ATCC 17802T 117 47,85 48,95 48,3 46,97 48,02 V pelagius LMG 3897T 188 47,85 49,3 45,91 47,62 47,67 V pomeroyi LMG 20537T 189 47,03 47,02 47,9 46,54 47,12 V ponticus CAIM 1751T 190 48,26 49,12 47,9 47,4 48,17 V rotiferianus LMG 21460T 191 47,85 48,77 47,9 47,84 48,09 V sagamiensis NBRC 104589T 192 47,24 44,39 43,91 47,19 45,68 V salilacus DSG-S6T 193 45,6 48,6 49,9 48,48 48,15 V scophthalmi LMG 19158T 194 48,47 48,25 48,1 46,97 47,95 V sinaloensis LMG 25238T 195 47,44 50 49,3 46,75 48,37 V splendidus LMG 19031T 196 47,24 45,61 47,9 46,75 46,88 V tasmaniensis LMG 20012T 197 47,03 46,32 46,31 47,62 46,82 V thalassae CECT 8203T 198 47,03 48,6 46,31 48,27 47,55 V toranzoniae CECT 7225T 199 47,65 45,26 46,91 46,1 46,48 V tubiashii ATCC 19109T 200 47,24 48,25 49,1 47,84 48,11 V variabilis CAIM 1454T 201 48,47 50,18 46,11 50,22 48,74 V vulnificus ATCC 27562T 118 48,26 50,35 49,5 48,92 49,26 Mean – 47,68 48,3 48,49 47,75 48,05 M. Duman et al. Aquaculture 562 (2023) 738728 7 are composed by the isolates belonging to the same species which a low genetic variability as CC1 V. alginolyticus, V. orientalis CC4 CC5 V. chagasii,CC9 iV. mediterranei,V. fluvialis CC10; V. barjaei CC11and the CC15 composed by Vi_42a, Vi_8b, Vi_72a and Vi_73a probably a new species belonging to Harveyi clade. On the other hand, the species with more genetic diversity are divided in many CC as V. diabolicus CC2, CC6, CC18 and CC21, V. parahaemolyticus CC7, CC8 and CC15, V. anguillarum CC12, CC13, CC19, and CC16 and V. harveyi CC18 and CC20. Other species more close genetically are mixed in the CC3 and CC22, V. owensii, V. jasicida and V. hyugaensis isolates. Most of the Turkish isolates belong to ST133 identified as V. anguillarum. The phylogenetic relationship of all aquatic species within the genus Vibrio was constructed by performing a NJ and ML tree (Supplementary Fig. S1). Before starting with the MLST studies, phylogenetic trees were performed with all Vibrio spp. (data not shown) with the aim of a correct identification of the isolates then, all Vibrio spp. not related with the strains used were deleted. The phylogenomic tree showed a clear identification of the different Vibrio spp. except in the case of V. jasicida/ V. hyugaensis and V. ordalii/V. anguillarum. The isolates were classified in 28 Vibrio spp. but most of them were identified as V. alginolyticus, V. diabolicus, V. parahaemolyticus and V. anguillarum/V. ordalii. Many strains deposited in the pubmlst database were reclassified (Supplementary Table S1). The isolates V103, Vi_2, Vi_28, Vi_41, Vi_67, Vi_2a, Vi_24a, Vi_27a, Vi_38a, Vi_42a, Vi_54a, Vi_56a, Vi_67a, Vi_72a, Vi_73a, Vi_80a and Vi_8b are not clearly identified. Turkish strains (51 of 256) were grouped into six different species V. anguillarum (28 strains), V. alginolyticus (15 strains), V. harveyi (5 strains), V. diabolicus, V. fluvialis and V. parahaemolyticus (one of each). The structure of the Median Joining (MJ) algorithm network (Fig. 2) is similar to the phylogenetic tree, with some exceptions. V. ordalli (h_184) and Vi_22 (h_10) were together in the same branch in V. anguillarum/V. ordalli group, but in MJ network, these strains were separated but in the same group. V. mytili (h_181) and V. natriegens (h_182) were shared branch in the ML tree; however, in the network, V. mytili was close to V. diabolicus (h_162) and V. natriegens to V. parahaemolyticus (h_117). The isolate Vi_2a (h_54) was included in the Table 3 Genetic characteristics, evolutionary variation and demographic history statistics of Vibrio population. Locus Number of alleles Number of polymorphic site (PIM) R min π θ d N /d S D D* F* F S S R 2 atpA 130 208 (186) 57 0,090 0,121 0,033 −0,801 −0,632 -0,851 −22,35 1 0,105 gyrB 160 246 (228) 77 0,137 0,124 0,0329 0,333 0,650 0,577 −21,83 1 0,158 pyrH 133 212 (198) 59 0,133 0,123 0,027 0,263 1051 0,767 −9,14 1 0,156 recA 160 207 (184) 59 0,144 0,133 0,0242 0,262 −0,009 0,154 −29,73 1 0,160 Concatenated 201 873 (796) 253 0,126 0,125 NA 0,033 0,301 0,189 −7,79 1 0,145 PIM, parsimony informative sites; R min , minimal number of recombination events; π , average number of nucleotides differences per site; θ, number of segregating sites; d N /d S ,ratio of mean non-synonymous substitutions per non-synonymous site/mean synonymous substitution per synonymous site; D, Tajima’s D statistic; D* and F*, Fu and Li statistics; F S , Fu’s statistic; S, Strobecks’s statistic; R 2 , Onsis R2 statistic; NA, not available. Fig. 1. Population snapshot obtained by goeBURST analysis. Each number corresponds to the ST. The size of each circle reflects the number of strains belonging to the ST. Single-locus variants (SLVs) are in black (link drawn without recourse to tiebreak rules) and blue (Link drawn using tiebreak rule 1). Double-locus variants (DLVs) are in grey and green (link drawn using tiebreak rule 2). CC, Clonal Complex. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.) M. Duman et al. Aquaculture 562 (2023) 738728 8 Mediterranei clade in the phylogenetic tree while in the network belonged to Scophtalmi clade. Vi_2 (h_1) and Vi_28 (h_16) were identified as Vibrio sp. closely to V. parahaemolyticus and V. alfacsensis but in the network were affiliated to V. harveyi 1DA3 (h_169). The isolates grouped into the two branches of V. harveyi (h_168 and h_169) swiped the positions in the MJ network. 3.3. Recombination and mutation The phi test (Pairwise Homoplasy Test) indicates non-significant statistical evidence for recombination in any gene. However, putative recombination events were detected with the minimal number of recombination events (R min ), 57, 77, 59 and 59 in atpA, gyrB, pyrH, recA, respectively. Moreover, the hypothesis that some recombination could also occur is supporting by the linkage disequilibrium estimated from allelic profiles for the whole isolates collection, I A parameter (2.0039 (P =0.000)) and standardized I S A =0.0721 (P =0.000)), significantly different from zero. Moreover, the presence of homoplasies was shown by the interconnection network structures in the split graphs, based on each gene (Supplementary Fig. S2) and on the concatenated sequences (Fig. 3). To investigate the population structure, we performed population analysis using Structure software. This software classifies individuals for a given number of clusters (K) represented by colors using one color for each ancestor lineage (Fig. 4). We performed Structure analysis by setting K from 2 to 39, executing the simulations five times for each K. The best model estimated by Structure Harvester (Earl and VonHoldt, 2012) was K =28, therefore, 28 genetic ancestors. The admixture of these lineages resulted in 37 profiles (Fig. 4). The Turkish strains showed a predominance of profiles with ≥80% of ancestor lineages (profiles 5, 16, 21, 25, 30 and 36). The isolates V7 (profile 7, V. fluvialis), V24 (profile 8, V. diabolicus), V103 (profile 37, V. alginolyticus) and V147 (profile 22, V. harveyi) were the isolates with more admixture. On the contrary, Italian isolates showed a high degree of admixture. It’s interesting to note that the Italian isolates Vi_42a, Vi_8b Vi_72a and Vi_73a (identified as a possible new Vibrio spp.) belonged to a different ancestry lineage (profile 2, yellow bars). All species encompass many different profiles and/or share the profile with other species except V. vulnificus (profile 4), V. barjaei (profiles 18 and 19, green bars) and V. harveyi ATCC 14126 T (profiles 21, 22 and 23, red bars). The patterns did not show any relation with host, geographical point or isolation year. 3.4. Gene flow and demographic trends The neutral theory of molecular evolution holds that most variation at the molecular level evolutionary changes occurs at the molecular level does not affect fitness; therefore, the genetic variation is explained by random processes (stochastic processes) as genetic drift. The statistics Tajima’s D, Fu and Li´s (F* and D*), and the Ramos-Onsis R 2 were different to cero rejecting the neutral theory. Tajima’s D and F* and D* were positive for all genes (except atpA and recA) and concatenated sequences reflecting a minimal presence of rare polymorphisms, typical of a population in contraction and under balancing selection. On the contrary, the negative values of Fu’s F S , high values of Strobeck’s statistic (S =1) and low values of Ramos-Oasis R 2 test are consistent with deviation from neutrality due to either selection or population expansion (Table 3). 4. Discussion Molecular identification of Vibrio genus by PCR and sequence analysis of one single gene has some incongruences in the species identification, a reliable and conceivable alternative is is identify and typing species by MultiLocus Sequence Analysis (MLSA) and by MultiLocus Sequence Typing (MLST) (Choopun et al., 2002). Similar to the reports, we found that V. alginolyticus and V. anguillarum are the most prominent species isolated from marine fish. The population structure and evolution of bacteria are primarily based on the MLST analysis which is one of the most useful methods. Despite other authors employed several genes to MLSA studies, the MLST scheme for Vibrio spp. (Rahman et al., 2014) was established in four genes (atpA, gyrB, pyrH and recA) on the pubmlst database. The population presented only 77 nonsynonymous changes of the 873 polymorphic sites, indicating a strong selection against amino acid changes, a typical feature of housekeeping genes (P´ erez-Losada et al., 2006). Determination of high number singletons in analyzed strains, about 89.5%, which resulted by the goeBurst scheme, also reflects some genetic discontinuity among the strains. More elaborated analysis by creating the median-joining network showed a tight relationship between some species of the Vibrio genus, despite the presence of missing intermediates in the population structure in this study. On the other hand, it is important to highlight the limited number of isolates available in the GenBank for the recently identified fish pathogenic Vibrio spp.(Sawabe et al., 2013) such as V. anguillarum, V. Fig. 2. Median joining network of Vibrio sp. isolates employed in the study based on the concatenated sequence (atpA, gyrB, pyrH and recA). The size of the circles is related to the number of strains within each ST. Red circles represent putative intermediate sequences. Number in grey circles indicate the nucleotide substitutions. Lines without numbers correspond with one nucleotide substitution. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.) M. Duman et al. Aquaculture 562 (2023) 738728 9 alginolyticus, V. harveyi, and this could cause a perturbation in the genetic structure (Balboa and Romalde, 2013; Bisharat et al., 2007; Bruto et al., 2017; Kim and Cho, 2015; Sawabe et al., 2013; Steinum et al., 2016; Tao et al., 2012; Warner and Oliver, 2008; Wollenberg and Ruby, 2009). V. parahaemolyticus and V. cholera are the most common strains into the Vibrio genus in the world because they cause pandemic disease in humans, and they have been found to have a high degree of genetic diversity within the species (Makino et al., 2003). On the contrast of these species, we found that V. vulnificus, V. crassostreae, V. fischeri, V. anguillarum, V. ordalii, V. tapetis and other pathogenic species reported that have a light linkage disequilibrium and commonly two lineages revealed by MLST analysis. The phylogenetic analysis, on the basis of the concatenated sequences, showed some species in a separated monophyletic branch, which was consistent with a recent taxonomic study on the genus Vibrio (Sawabe et al., 2013). From all studied aquatic Vibrio isolates here, phylogenetic tree showed four big groups differentiating the species V. alginolyticus, V. diabolicus, V. parahaemolyticus and V. anguillarum, each of them grouping isolates from live marine animals including rainbow trout, sea bream and sea bass isolates, and the rest of isolates suggesting that the changes in the sequences can be caused by geographical isolation. Abayneh et al. (2012) reported that the MLSA approach is capable of resolving isolates according to their geographical origin or host. In a recent report, Vibrionaceae has been reported as a very complex group with determining 51 distinct clades including 21 newly defined ones by Jiang et al. (2022) and MLSA was reported that it is still an effective and reliable tool for delineating new species and monophyletic groups/clades. However, our results by MLST analysis did not establish any association between geographical location, isolation source or date, and further research needs to be clarified if any relation between the isolates. The bacteria belonging to the family Vibrionaceae present a high level of genotypic diversity supported by a high recombination/mutation rate, especially because vibrios are natural constituents of marine environments and spend much of their life cycle outside the host and have adapted to fluctuating conditions, including temperature shifts, nutrient limitation, osmotic stress, and predator (Ceccarelli et al., 2019). In this study, the high recombination (R min ), 253 events for concatenated sequence, supports genotypic diversity of aquatic Vibrio spp.. The construction of the Median-Joining network revealed a high number of mutations to separate some species (V. anguillarum, V. alginolyticus, V. parahaemolyticus, V. mediterranei/V. shilonii); therefore, the high number of mutations together with the I S A value showed that, although the impact of recombination is an important evolution factor, Fig. 3. Neighbor-Net graph based on the concatenated sequences of the four housekeeping genes. The formation of an interconnecting network structure was suggestive of recombination. Scale represents the genetic distance (Clades were determined according to Sawabe et al., 2013 and Jiang et al., 2022). M. Duman et al.