Host Species and Body Site Explain the Variation in the Microbiota Associated to Wild Sympatric Mediterranean Teleost Fishes

Microorganisms are an important component in shaping the evolution of hosts and as such, the study of bacterial communities with molecular techniques is shedding light on the complexity of symbioses between bacteria and vertebrates. Teleost fish are a heterogeneous group that live in a wide variety of habitats, and thus a good model group to investigate symbiotic interactions and their influence on host biology and ecology. Here we describe the microbiota of thirteen teleostean species sharing the same environment in the Mediterranean Sea and compare bacterial communities among different species and body sites (external mucus, skin, gills, and intestine). Our results show that Proteobacteria is the dominant phylum present in fish and water. However, the prevalence of other bacterial taxa differs between fish and the surrounding water. Significant differences in bacterial diversity are observed among fish species and body sites, with higher diversity found in the external mucus. No effect of sampling time nor species individual was found. The identification of indicator bacterial taxa further supports that each body site harbors its own characteristic bacterial community. These results improve current knowledge and understanding of symbiotic relationships among bacteria and their fish hosts in the wild since the majority of previous studies focused on captive individuals.


INTRODUCTION
Symbiosis among microorganisms and hosts are driving forces of evolution.Among microorganisms involved in symbioses with animals, bacteria are the most widespread and ubiquitous, shaping the evolution of their hosts by contributing both positively (basic functions) and negatively (pathogenicity) to their fitness [1].In-depth knowledge of the taxonomic composition of the microbiota, and their relationships with their hosts, can provide insights into both the function and dysfunction of the host organisms [2], due to the important role that host-microbe interactions play in the physiology and performance of animals [3].
Fish constitute about half of living vertebrate species [4], but the majority of studies on microbiota have been conducted in mammals [5,6].High throughput 16S rRNA gene sequencing is being increasingly employed to investigate the fish microbiota (e.g., [7, 8, 9,]), producing high-resolution descriptions of their community structure and diversity [2,10], and providing further understanding of the consequences of symbioses with bacteria [11].Despite progress, important gaps remain in our current knowledge regarding the factors that shape fish microbiota [12], such as environmental or phylogenetic variables, whose influences could vary among species and / or habitats.Furthermore, most of the information available is focused on the gut of captive species used in aquaculture, i.e., fish of economic significance and commercial interest [8,13], or from model organisms such as the zebrafish [14].
Wild teleost fishes exhibit a broad variety of morphologies, physiologies, ecologies and natural histories [11,15].Therefore, they represent a good model group to investigate the factors that shape host-bacterial assemblages in nature [16].Fishes are in intimate contact with the surrounding water, and bacterial colonisation from the environment is thought to be one of the primary mechanisms of microbiota acquisition for fishes [16,17].Hence, environmental factors influencing water microbial reservoirs, such as water salinity, could ultimately influence fish microbiota [11,18].The water environment also hosts high concentrations of potentially pathogenic microorganisms [19] that may enter the fish body.Therefore, the dynamics and diversity of the microbial communities originating from the surrounding water could influence significantly the individuals' fitness [20].
The external surfaces of fishes, such as the skin and gills, are coated in a mucus secretion that hosts an indigenous microbiota dominated by bacteria [10,21,22] that acts as a protective barrier against pathogens [23][24][25][26].The digestive tract of fish also receive water and food that are populated with microorganisms that will undoubtedly affect the resident microbiota.
Besides the surrounding water, host-related factors may also shape fish-associated microbial communities [18,27].In some teleost species, the skin bacterial diversity is very different from that of the surrounding bacterioplankton, and varies among different parts of the host body [28].In addition, different species reared in the same water environment vary in their gut [18,29] and gills [9] microbiota.
We hypothesize that wild teleost microbiota could be influenced by abiotic environmental factors and host-related factors.If host-related factors are more important than abiotic factors in determining the fish microbiota, we would expect to find specific microbiota from each part of the body, and in different fish species, independently of when the fishes were captured.To test this hypothesis, and evaluate the importance of these different factors, multiple individuals of 13 species of teleost were collected at the same location in the shallow Mediterranean Sea on four different occasions and four body sites were sampled (mucus, skin, gills, and gut).The microbial community composition of these samples was characterized by 16S rRNA gene sequencing and their alpha and beta diversity compared.We also tried to determine indicator bacterial taxa that could be representative for each teleost species and body site, being both unique to a given group (exclusivity) and occurring in all sample units within a group (fidelity) [30,31].
This study is novel in characterizing the bacteria associated with wild fish, while the majority of previous studies addressed bacterial assemblages in captive fish.

Fish sampling
Fishes were captured in the NW Mediterranean Sea, Bay of Banyuls (Gulf of Lion, Banyuls-sur-mer, SE France, 42°29'4.618''N,3°8', 35.39''E) on 4 different days in 2017: June 21, June 26, July 18, and October 4, and several physicochemical variables for each sampling were recorded (Table S1).We captured a total of 59 individuals belonging to 5 families and 13 species (see Table 1).For all fish individuals except the Family Gobiidae, one gill net was placed between 0 and 6 m depth during the entire night period, and at sunrise, it was recovered with the fish (42°29'15.073''N,3°7', 49.688''E).Fish were collected dead, handled with gloves and stored into individual plastic bags.
They were immediately brought from the vessel to the laboratory.They were kept at 4 °C until dissection, within the next 48 h, mostly within 6 h.Individuals from the Family Gobiidae were captured with one net placed for just 2 hours close to the seashore.
All thirteen species share nearly the same environment since they were captured in the same place at the bay, and most of these species have similar feeding habits.In general, they are omnivorous, eating mainly small crustaceans, invertebrates and small fishes, but they can also ingest, to different degrees, other particles such as algae and phytoplankton [33].There is an exception for the species S. salpa, in which juveniles are omnivorous (mainly carnivores), but adults are exclusively herbivores [33].The adult stage is reached at about three years, when the body size is longer than 20 cm [34].The length of individuals belonging to S. salpa that were captured were shorter than 20 cm for four of the five individuals (the lengths were 16.8, 17.2, 15.5 and 14 cm).The fifth individual was of adult size (23.8 cm.)This is the only individual, from all fish investigated here, whose feeding would be exclusively vegetarian.For the other 12 species, no significant differences in feeding habits according to age are described [33].
Once in the laboratory, each sample was taken as follows: (i) the skin mucus was scraped with a sterile scalpel from the entire body surface; (ii) when the mucus was completely removed, a 3 cm 2 piece of skin was cut from the central part of the body, close to the lateral line; (iii) then, one gill arch was cut from each side of the body; finally, (iv) 5 cm of the intestine distal part was cut.We cut the final part of the gut because, due to defecation, it reflects the bacterial community of the entire intestinal tract [35,36].
Because of possible differences in bacterial communities among different parts of the same organ [28], samples were always taken from the same part of the body in all the different species, sterilizing the material by flaming before taking each part of the same individual.Samples were immediately placed into sterile Eppendorf tubes and kept at -80 ºC until DNA extraction.In the case of the species Scorpaena notata, there was no mucus on the skin, so the skin piece was cut directly.
At each sampling site, 2 litres of seawater were taken from the same place where the nets were installed and the water stored in a sterile glass bottle.Briefly, in the lab, the water was vacuum filtered through a 47 mm diameter membrane filter with a pore size of 0.2 µm (Fisherbrand, Thermo Fisher Sientific, France).This filter was kept at -80 ºC prior to DNA extraction.

DNA sequencing
A total DNA extraction was performed using the Kit Quick-DNA Fecal/Soil Microbe MiniPrep Kit (Zymo Research, Orange, California), following the manufacturer's instructions.The V4-V5 regions from the 16S rRNA-encoding gene of isolated DNA were amplified by means of PCR, using the universal primers 515F-Y (5'-GTGYCAGCMGCCGCGGTAA) and 926R (5'-CCGYCAATTYMTTTRAGTTT]) [37].The PCR mix contained 5 µl of KAPA 2G, 0.2 µl of each primer, 3.6 µl of ultrapure water and 1 µl of DNA for a final volume of 10 µl.After 3 min of initial denaturation at 95 °C, the following conditions were applied: 22 cycles of 95 °C for 45 s (denaturation), 50 °C for 45 s (annealing) and 68 °C for 90 s (extension), ending with a final extension at 68 °C for 5 min.For each sample, three PCRs were performed in the same conditions, to increase the DNA quantity, while minimizing cycle numbers to avoid PCR bias [38].
The product of each PCR was verified by agarose gel electrophoresis and triplicate reactions were pooled.
A second PCR was performed to attach the Illumina adapters and 8-bp barcodes for multiplexing.They were added in a reaction mix in which barcode sequences were individually added for each sample.The mix for this second PCR contained 12.5 µl of KAPA 2G, 0.5 µl of each barcode primer (Nextera Index Sequences in http://seq.liai.org/204-2/),10.5 µl ultrapure water and 1 µl of DNA for a final volume of 25 µl.PCR conditions were as follows: initial denaturation at 98 °C for 30 s, 8 cycles of 98 °C for 10 s, 60 °C for 20 s, 72 °C for 30 s, and a final extension at 72 °C for 2 min.PCR products were purified using the USB ExoSAP-IT PCR Product Cleanup (Thermofisher, France) kit by incubating the samples at 37 °C for 30 min, followed by 15 min at 85 °C.The concentration of samples was normalized with the SequalPrep Normalization Plate (96 well) Kit (Thermofisher, France).Thereafter, all the samples with different barcode sequences belonging to each run were pooled, and the pooled products were quantified by using the Quant-iT TM PicoGreen (Thermofisher, France).
The normalized amplicons were concentrated by using the Wizard SV Gel and PCR Clean up Kit (Promega, France) and amplicons were sequenced one run in an Illumina® MiSeq sequencer using the 2 × 250 bp protocol by FASTERIS SA, Switzerland.

Sequence analyses
Sequence analysis was performed by using a combination of tools, including the USEARCH v9.0.2132 [39] program, Qiime V.1.9.1 [40] and bash scripts.Briefly, from demultiplexed fastq files obtained from the sequencing center, forward and reverse reads were merged to create consensus sequences in a single fastq file using usearch9 -fastq_mergepairs.Assembled reads were quality trimmed by the usearch9 -fastq_filter command.Primer sequences were removed from raw amplicon sequences using strip_primers_exclude.py.At this point, sequences were dereplicated and denoised using usearch9 -unoise (minampsize 8).As a certain number of artifactual short sequences remained, sequences less than 360 bp were removed with a combination of bash scripts and the filter_fasta.pycommand from Qiime V.1.9.1.OTU tables were generated by usearch9 -usearch_global against the primer-stripped raw sequences with a 97 % identity to OTU representatives ordered by abundance.OTUs were identified using Qiime V.1.9.1 assign_taxonomy.py-m rdp and a previously described modified version of the Greengenes [41] August 2013 database [9].The taxonomy assignments were corrected to contain full paths from phylum to species.Final OTU tables were generated using biomformat.Sequences matching "Archaea", "Eukaryota", "Unassigned", "Chloroplasts" and "Mitochondria" were discarded using filter_taxa_from_otu_table.py.Details of data analysis, and the analysis environment are described in the supplemental methods.Further details on the bioinformatics methods, as well as the samples characteristics, including the number of reads and of sequences, are in the Electronic Supplementary Material.All the resultant sequences were deposited in the Sequence Read Archive (SRA, http://www.ncbi.nlm.nih.gov/sra)(Accession number: PRJNA531247).
Statistical analyses R 3.5.0[42] was used to calculate both α-diversity (Shannon diversity index) and βdiversity matrices of distances among samples.Matrices of genera were used since they contained the most detailed information on the total community, because in many cases species assignment was not possible using the rdp classifier and the sequenced region.
However, in some of the OTUs, the genus classification was not possible, and in those cases we used the most precise taxonomic level that we could identify (for example the family, or the class).The Shannon diversity index based on OTUs was estimated by the R vegan package [43].Due to the high variation of the index values, it was standardized before subsequent analyses.This technique is very useful and provides simplification, since it modifies the scale and considers the mean as zero (e.g., [44]).This variable was fitted on a normal distribution and met the normality criteria.Differences in Shannon's diversity between fish species and body parts were analysed using one-way ANOVA in Statistica 7.1 [45].In addition, a General Linear Mixed Model (GLMM) was applied to identify which variables (fish species, body part and sampling date, the last being considered as a random factor) influence bacterial richness.
To evaluate the effects that the species, body site, individual fish or sampling date would have in explaining the variability of the bacterial communities' composition in teleosts (i.e., β-diversity), we performed a PERMANOVA analysis.First, a distance matrix on resultant OTUs was constructed with the function vegdist in the R package "vegan", which was considered as the dependent variable.Here, bacterial communities of each sample are compared with all the rest, by means of the Bray-Curtis dissimilarity index for pairs of samples [43].The index values ranging from 0 to 1: a value of 0 (i.e., no distance) indicated that two samples share all the bacterial taxa in the same abundances, and 1 (i.e., maximum distance) indicated that they did not share any of the detected OTUs.This distance matrix was considered as the dependent variable in further analysis.Then, the adonis function of the "vegan" R package was used with all the previously mentioned explanatory factors, considering the sampling date as the stratification variable.Bray-Curtis distances were also represented by means of a non-metric multidimensional scaling (NMDS), implemented in the R package phyloseq [46].Because for some of the groups we have less than 3 replicates (Table 1), we performed the analyses with and without those groups, obtaining virtually the same results.We show here the results with N > 3 for each category, and in the supplementary material, the results obtained with the whole dataset (Tables S2, S3 and S4).
Finally, we performed an analysis to determine indicator species of the community in which the patterns of species distribution were compared among different categories (all body sites within each fish species), by using the Indicspecies R package [47,48].The results display only bacterial taxa that are significantly (p < 0.05) prevalent in a given category while absent in the rest.A target taxon is considered as "indicator" if it is exclusive to a given group (exclusivity) and occurs in all sample units within a group (fidelity).In this way, the analysis is based on the abundance within particular groups (exclusivity), and on the relative frequency of species within groups (fidelity; [31,49]).
However, groups may be composed of a single or multiple categories (such as the case where indicator species are present in all of the samples of both categories but absent from the rest).An indicator value index is assigned between a species and each group, identifying the group with the highest association value.Then, randomization methods (permutation tests) are used to test the statistical significance of values [31,49].

Bacterial taxonomic composition of the fish and water samples
We obtained a total of 2,964,227 sequences assigned to bacteria (i.e., filtering out reads belonging to Archaea, Eukarya and unidentified reads).After filtering, the remaining OTUs belonged to 181 genera, spread over 108 families, 67 orders, 37 classes and 13 phyla.The most abundant bacterial group in the fish samples was Proteobacteria, when considering both body site (Table 2) or host species (Table 3), and this group was also the most abundant in the water column.The second most abundant group in the water samples was Bacteroidetes (34.1 %), while in fishes this group was not highly represented, except in S. salpa (16.1 %).In fishes, the other abundant bacterial phyla were Fusobacteria and Firmicutes, that were rare in the surrounding water (0.7 and 0.3 % respectively).Of the most abundant taxa (> 1% of total sequences) Proteobacteria showed a relative abundance of around 80 % in the mucus, skin and gills, although its abundance was only 37.6 % in the intestine.The other most abundant taxa in the intestine were Tenericutes, Firmicutes (25.9 %), and Fusobacteria (13.7 %), and these groups were either absent from (Tenericutes), or showed lower abundances in the other three body sites.
In detail, the majority of the OTUs retrieved from fishes belonged to the family Vibrionaceae (Proteobacteria, class Gammaproteobacteria), and were mainly assigned to the genera Vibrio and Photobacterium, followed by the genus Cetobacterium (Fusobacteria; Table 2).The relative abundances of these three genera were similar in mucus, skin and gills, although in the latter tissue, Cetobacterium were more abundant than Photobacterium (Table 2).In the intestine, the most abundant group was the order Clostridiales (Firmicutes), and the Genus Cetobacterium (Fusobacteria).Within the family Vibrionaceae, the most abundant genera in the intestine included Aliivibrio and Photobacterium.Vibrio was clearly less abundant in the intestine compared to external surfaces.The proportion of each bacterial group (for groups that represent more than 1% of the total abundance) varied among different species (Table 3), although Photobacterium and Vibrio were predominant in the majority of the species.

Alpha-diversity of bacterial communities
There were significant differences in levels of Shannon α-diversity when comparing all the different fish species and water (F12,174 = 2.16, p = 0,015, Figure 1).However, when removing the water from the analyses, differences among species were no longer significant (F11,171 = 1.54, p = 0.12).This result implied that differences were among the water and the rest of species, with water displaying a higher α-diversity than all the fish species (Fig. 1).When separating by the water and different body sites, we found significant differences among groups, with again water showing the highest α-diversity, and then the mucus (F4,82 = 20.39,p < 0.0001, Fig. 2).Post-hoc tests revealed no significant differences among the water and the mucus samples (p = 0.28), while both categories were different from the other three groups (all p < 0.001).Therefore, in this case, when removing the water from the analyses, significant differences remained due to differences between the mucus and the rest of body sites (F3,179 = 23.13,p < 0.0001).
No differences were found when comparing the skin, gills and intestine (all p > 0.1, Fig. 2).Results with all the dataset are presented in Table S2.
The GLMM results suggest that fish species (F11,165 = 1.82, p = 0.05) and body site (F3,165 = 23.84,p < 0.0001), but not the sampling day (F3,165 = 0.85, p = 0.47), influenced the Shannon diversity index.Results were the same with the entire dataset (Table S3).The α-diversity comparisons among the body sites in each species, and in all the species for each body site, are presented in Fig. S1 and S2.

Beta-diversity
The PERMANOVA results based on Bray-Curtis distances (i.e., β-diversity, Fig. S3) determined that both factors, fish species and body site, explained the variation in bacterial community composition among the samples (Table 4).This means that individual fish within a given species have significantly more similar bacterial communities than with fish from other species, which is also the case with the different body sites.However, individuals and sampling date did not have any significant effect on the structure of the microbiota, which means that individuals from the same species share their microbiota, independently of the sampling date (Table 4).This result clearly discards the possibility of any contamination among fishes at the same sampling date.Tests done with the whole dataset were similar (Table S4).

Indicator bacteria
Several indicator bacterial groups were significant (all p < 0.05, listed in Supplementary material, Table S5).Results displayed a high variation among different fish species regarding indicator microbial species.Sarpa salpa hosts the highest bacterial diversity, and also the highest number of indicator bacteria (95 in total, see Table S5).In addition, only this species hosted the phylum Bacteroidetes, which was relatively abundant in the water.No indicator bacterium was found in Gobius cruentatus, nor in the water samples.
In all fish species except S. salpa, the gut samples did not share any indicator bacteria with the other three body sites.For the other body sites, gills and mucus shared indicator bacteria in five host species, the gills with the skin in two species, and the mucus with the skin in seven species.In six fish species, the mucus, skin and gills shared several indicator bacterial taxa (Table S5).
The body site with the highest number of indicator bacteria was the mucus (Table S5).
However, in the species S. notata that lacks external mucus, indicator bacteria from skin were by far the most numerous in this part of the body (18 taxa against nine in S. salpa, and one in G. bucchichi).In the rest of the species, no indicator bacteria were found exclusively in the skin.
An overview of indicator bacteria from all fish species showed that mucus had the highest number of exclusive bacterial taxa (30 taxa), followed by skin (15 taxa, Fig. 3).Although there were many bacterial taxa that were shared between the different categories, the highest number was shared by mucus, skin and gills, whereas the gut shared a lower number of indicator bacteria with the other three body sites (Fig. 3).

DISCUSSION
Our results confirm that Proteobacteria was the predominant phylum in water and teleost fish in the Mediterranean Sea, as it was indicated by previous studies [50].However, the prevalence of Bacteroidetes, Fusobacteria and Firmicutes differed among the water and fish samples.
Both environmental and/or host-associated factors may shape fish microbiota at different levels [9,51,52,53].Previous studies found that the mucus microbiota may be highly variable and dynamic, depending on the environmental conditions including seasonality [54,55].However, other studies suggested that genetic factors have greater effect than environmental factors in shaping mucus bacterial assemblages [56,57].In our study, we have found that individuals from the same species shared significantly more bacterial communities than they did with individuals of the other species.This was observed despite being captured in the same place and having overall similar feeding habits, which suggests species-specific assemblages.These results are in accordance with previous works claiming that the autochthonous microbes are not a passive reflection of their habitat communities [11], i.e., fish have their own microbiota.We did not detect any effect of the sampling date, although the variation in abiotic factors was moderate among our four samplings (see Table S1).Therefore, other factors may influence teleosts' microbiota, such as seasonal environmental variation throughout the year (e.g.related to temperature or salinity changes), which should be tested in future studies.
Bacterial communities from a given body site were significantly more similar among themselves than with communities from other body sites.Although the majority of bacteria from mucus, skin and gills belong to the Family Vibrionaceae, there are differences in the relative prevalence of genera among the gills with respect to mucus and skin.The highest bacterial diversity in the body was found in the outermost external layer, i.e., the mucus, which is in intimate contact with water, thus a continuous exchanging of bacteria among them could be expected.The mucus acts as a protective barrier for fishes, since it may impede the penetration of potentially pathogenic microorganisms inside the body [24].Our findings that the mucus had a significantly higher diversity than the rest of the body (including the skin), is in accordance with this protection hypothesis, suggesting that environmental bacteria remain trapped here.Moreover, the fact that S. notata, the only fish species without mucus from our sample, had the highest diversity of indicator species in the skin, also supports this hypothesis since in this case, the skin is the direct interface between the fish and the environment.However, we cannot exclude that since the mucus is rich in polysaccharides (e.g.[58,59]), it may also constitute a food source that attracts heterotrophic bacteria, which may in turn lead to a more diverse community.
Pathogenic bacteria can enter the host fish through the skin, gills or gastrointestinal tract, and the integrity of these physical and immunological barriers determine the outcome of host-pathogen interactions (reviewed in [13]).Balanced and complex interplays within the mucus layer are thus key to disease resistance [60], and are essential for supporting host health and fitness [61].
Differences in the relative abundance of bacteria between different sites on the body are more pronounced between the gut and the other three external body parts (mucus, skin and gills, Table 2).For instance, the phylum Proteobacteria was found at lower relative prevalence in the intestine than in the other three parts, and relative abundance of Firmicutes and Fusobacteria also differ.In addition, we found that indicator bacteria can be shared between the gills, mucus and skin (see Table S5); however, almost no indicator bacteria were shared by the gut and those three external body sites.Therefore, our results suggest that, although we found that body site explains a significant fraction of the variance observed in bacterial communities over all the body, it is mostly because the intestine differs markedly from the rest of the body sites.Overall, these results were not unexpected given the very different physico-chemical conditions (oxygen, pH and organic substrate levels) occurring in guts compared to external body parts.Bacterial community composition has been previously shown to vary in the gut depending on the life cycle stage, diet, environment, and region of the gastrointestinal tract, and even varies greatly among individuals within the same species [2,5,10,16,61,62].In our study, we only examined the distal part of the gut, and individuals from the same species were adults, so we did not find any significant within-species effect, while the fish species explained most of the variance in bacterial assemblages.
There are some groups that predominate in the gut of marine fish across different studies, such as Proteobacteria, Firmicutes, Actinobacteria and Bacteroidetes.Also, Fusobacteria, Clostridia, Bacilli and Verrucomicrobia have been identified (reviewed in [2,10]), as well as Mycoplasma [5], a group that we found exclusively in the gut samples of teleosts.There is some congruence in the identity of bacterial taxa we found in the gut (including Vibrio, Clostridium and Mycoplasma) with other marine fish belonging to different species, including some from aquaculture (e.g.Atlantic salmon, [5]).The presence of similar bacterial taxa in the gut-associated microbiota across different fish species, populations or geographic locations, suggests that these microbes are important contributors to host gut functions, such as digestion, nutrient absorption and immune response [63].In general, although the intestine ecosystem is expected to harbour a dense population of microbes, sequence-based analyses have demonstrated that it comprises the lowest phylogenetic diversity compared to external parts of the body [2].This was also observed in the present study when comparing α-diversity (Shannon diversity index) to the external mucus layer.Previous works also found a higher bacterial diversity in external organs such as skin or gills than in the gut [56,60], which could be due to a reflection of environmental diversity in the exterior of the body, whereas the gut may offer more stable habitats leading to a more specialized microbial communities [60].
There was no evidence that sampling day had an effect on the microbiota, which is in accordance with the difference between fish and water microbiota.That is, the identity of the fish species was more important in the structuring of bacterial communities than environmental conditions at the moment of their capture.This observation suggests that our results are not biased by the sampling date, and that fishes display their characteristic microbiota independently of the external conditions, at least between these two time points.Different samplings in the same location and with the same species could be performed throughout the year to verify an effect of seasonality on the bacterial diversity of our fish community as described in other aquatic systems (e.g.[54]).
In conclusion, we found that the thirteen Mediterranean teleost species studied here, living in the same environment and sharing most ecological traits, differed in their bacterial microbiota composition.This suggests that host taxonomic status mainly shapes fish microbiota.In addition, we found a characteristic microbial community in different parts of the body, indicating that microbiota are also influenced by local characteristics of their animal-associated microhabitats.The absence of a significant effect of individuals on the bacterial communities that they host suggests that the level of intra-specific variation is significantly lower than the level of inter-tissue and inter-species variations in shaping fish microbiota.Interestingly, the mucus cover showed the highest bacterial diversity, which supports the hypothesis that it is a barrier between the fish and its environment.That is, bacteria may be retained in this layer but do not reach the skin, and therefore do not penetrate into the body.This physical barrier may thus help to impede pathogenic infections.
Microbial communities associated with fishes are key factors in host physiology, ecology and evolution [1,64].Therefore, the knowledge of the factors shaping microbiota may help predict how changes in abiotic or biotic conditions affect bacterial assemblages and their functions.Our results help to clarify these factors in wild fish, highlighting the importance of fish species and body sites.Further studies should investigate the possible functions of bacterial taxa on their hosts, and also investigate environment-induced variations through time-series samplings.

Flavobacteriales(
Bacteroidetes phylum) were abundant in the seawater samples (30 %) as were the Rhodobacterales (family Rhodobacteraceae) and Rickettsiales (family Pelagibacteriaceae) orders of the Alphaproteobacteria class, showing respectively relative abundances of 18.7% and 14.5%.The class Gammaproteobacteria constituted 7.2 %, with Vibrio representing 1.5 % of total reads of this class and the Synechococcaceae family (Cyanobacteria) accounted for 4.2 of the total sequences.

Figure 2 .
Figure 2. Standardized Shannon diversity index of different body parts and water.Vertical bars denote 0.95 confidence intervals.

Figure 3 .
Figure 3. Venn diagram representing how many indicator bacterial taxa are exclusive or shared among the different body sites across all the fish species.More information on the specific indicator bacteria is given in TableS5.

Table 4 .
Results of PERMANOVA on factors explaining bacterial communities composition (β-diversity).Significant results are in bold.