Flatfish intestinal microbiota depend on various host traits, and vary with sediment type and bottom trawling effort

Comprising a diverse and complex ecosystem of bacteria, archaea, fungi, protozoa, and viruses, the intestinal microbiota play a critical role in animal physiology in general, and in fish particularly^1,2. For instance, the intestinal microbiota influence immune functions, homeostasis, food digestion, growth and reproduction^3,4,5,6,7,8. In humans, where intestinal microbiota have been studied most extensively, healthy intestinal microbiota have even been shown to enhance brain health and stress response⁹. Due to their shorter generation times and higher genetic variability, microbiota can respond more rapidly to environmental changes, thereby influencing the adaptive responses of their marine hosts^10,11.

The factors that shape the composition and diversity of these microbial communities through space and time are less well understood, but appear to be habitat, species and tissue specific^12,13,14,15, and correlate with variables such as diet, age, pollutants, salinity or temperature^{3,16,17,18,19}.

As the largest and most diverse group of vertebrates comprising more than 36,000 species^20,21, fish play an essential role in a variety of aquatic ecosystems. Notably, they contribute to maintaining the biodiversity, stability and resilience of those ecosystems. Further, as a fundamental part of aquatic food webs, they ensure the recycling and transport of diverse nutrients, even between aquatic and terrestrial ecosystems²². As consumers in marine food webs, fish contribute to the biological carbon pump by transporting carbon to the seafloor, thereby helping to mitigate the impacts of climate change²³. Beyond providing fundamental ecosystem services, fish account for important economic value and contribute to global nutrition^24,25. Despite the high species diversity, only a small fraction is targeted by commercial fisheries^24,25. For instance, in the North Sea, home to approximately 197 fish species, only a small number of species accounts for the majority of total landings^26,27. Demersal fish, such as flatfish, constitute a significant proportion of these landings, highlighting their importance in the region’s fisheries^27,28. With the rapid decline of fish populations and the collapse of entire commercially exploited stocks, such as the North Sea herring (Clupea harengus) in the 1970’s, gaining a deeper understanding of population dynamics and fish health is increasingly important^24,29,30.

Fishing practices are known to impact the marine environment, with bottom trawling being the most pervasive anthropogenic disturbance of seabed habitats^27,31. Shallow shelf sea regions are especially subject to trawling, with the central and southeastern North Sea being among the most heavily trawled areas in the world³². Besides the removal of flatfish directly, trawling may also indirectly impact the remaining flatfish populations, by altering the communities of seabed organisms that partially form their diet, or by resuspending the upper sediment layers, stirring up fine particulate matter and contaminants³¹. Furthermore, recent evidence suggests that trawling may also alter the microbiota of the sediments^33,34. Both diet and sediment serve as crucial microbial sources for the microbial communities associated with flatfishes and other seafloor-inhabiting holobionts, playing a primary role in shaping the composition and diversity of their intestinal microbiota. However, whether trawling impacts on flatfish diet and sediment microbiota may be indirectly affecting the microbiota of seafloor dwelling holobionts, including flatfishes, remains unknown.

To gain better insight into the flatfish intestinal microbiota and improve our understanding of wild flatfish populations, we analyzed the intestinal microbiota of three flatfish species (Buglossidium luteum, Risso 1810; Limanda limanda, Linnaeus 1758; and Pleuronectes platessa, Linnaeus 1758) that were collected in the southeastern North Sea. Using 16S rDNA metabarcoding, we characterized the prokaryotic communities of 162 fishes. We refer to this community of prokaryotes as intestinal microbiota from here on, based on the definition by Berg et al.³⁵, that identifies the microbiota as microorganisms associated with the host.

Subsequently, we used uni- and multivariate generalized linear models to test the predictive potential of different variables on the structure and diversity of the intestinal microbiota, each representing individual hypotheses. This included a combination of host variables (i.e., species, sex, condition factor, age and weight), and environmental variables (i.e., sediment type and trawling intensity). In addition, we analyzed differential abundances associated with these variables, to identify microbial markers that potentially fulfill important roles in the flatfish intestinal microbiota.

Sampling

Flatfishes were sampled during an expedition with the German Fishery Research Vessel ‘Solea’ from 07.05.2021 to 18.05.2021 (Cruise number: SOL791). The expedition took place in the Natura 2000 area Sylt Outer Reef in the southeastern North Sea (Fig. 1) and was carried out within the frame of a systematic survey on demersal fish and epibenthic fauna. A 2-m beam trawl with a 20 mm mesh size outer net and 4 mm mesh size codend, and a 7-m beam trawl with a 70 mm mesh size forenet and a 20 mm mesh size codend were used. At each station, the two-meter beam trawl had a towing time of five minutes on ground, while the seven-meter beam trawl remained on the bottom for 15 min. Flatfishes were caught exclusively during the day and for the present study the species Pleuronectes platessa (European plaice), Limanda limanda (Common dab) and Buglossidium luteum (Solenette) were selected from the catch. When possible, one male and one female were randomly selected at each station. The fish were placed in plastic bags and frozen at −20 °C on board. After the cruise, the fishes were transferred to the laboratory in cooling boxes, where they were again stored at −20 °C. In total 183 individuals of the species P. platessa, L. limanda and B. luteum were analyzed. Prior to dissection of fish individuals, the sex of the thawed fish was re-determined. Then, the total length of the fish was measured from the snout to the tip of the caudal fin and the weight was recorded. Scissors, forceps, and scalpels were cleaned and sterilized with a burner before the dissection of each fish and all procedures were carried out in a sterile hood.

To determine the age of the fish, the otoliths were removed, cleaned and stored in individually in tubes, and the age was determined at the Thünen Institute of Sea Fisheries (Bremerhaven, Germany). Subsequently, the stomach, intestines, liver and in females also the gonads were removed, and their weight was measured individually. Intestines were then collected in separate tubes and stored at −20 °C for DNA extraction.

DNA extraction and library preparation

To allow optimal DNA extraction from intestinal prokaryotes, the frozen intestines were first cut into small fragments and then crushed inside the tube with sterile scissors, resulting in a well-mixed suspension. From this suspension, subsamples of 100 to 125 mg were used for DNA extraction with the ZYMO Quick-DNA™ Fecal/Soil MicroPrep Kit (D6012; ZYMO Research) according to the manufacturer’s protocol. To verify potential contamination during the extraction process, DNA was also extracted from several blanks. Amplicon library preparation followed the two-step PCR protocol from Gohl et al.³⁶. The 16S-V4 region was amplified using the primers 515 F (S-*-Univ-0515-a-S-19) and 806R (S-D-Arch-0786-a-A-20³⁷), using the Phusion Green Hot Start II High-Fidelity PCR Master Mix (ThermoFisher Scientific), and 1 µl of DNA template per reaction. PCRs were also conducted on Mock communities (ZYMO-D311) and negative controls. The program of the first PCR included an initial denaturation step of 98 °C for 3 min, followed by 30 cycles of 98 °C for 30 s, 55 °C for 30 s and 72 °C for 30 s and a final elongation step of 72 °C for 5 min. PCR products were run through a 1% agarose gel electrophoresis and samples without a visible product underwent five additional PCR cycles. Then, PCR products were diluted 1:100 and used as templates for the second PCR, which was conducted with the same reagents and the indexing primers from Gohl et al.³⁶. This PCR included 10 cycles but was otherwise identical to the PCR 1. The final amplicons were then normalized, purified and pooled using sequelPrep plates (ThermoFisher Scientific). The sequencing, on the Illumina MiSeq platform with the MiSeq Reagent Kit v3 (600-cycle), was done in the German Centre for Marine Biodiversity Research (DZMB; Wilhelmshaven, Germany).

Data preparation and processing

The raw sequencing files were demultiplexed, quality filtered, denoised, purged from chimeras and clustered into OTUs based on the 97% similarity criterion, using the OptiClust algorithm³⁸with the software MOTHUR (v.1.48.0³⁹) and classified with the SILVA database (v138.1⁴⁰). For details on the bioinformatic processing see the dedicated script on https://github.com/gbonthond/flatfish_microbiota. Finally, samples were purged of reads classified as chloroplasts, mitochondria or eukaryotes, as well as those lacking a domain level classification. Samples with fewer than 1,000 remaining reads were excluded from downstream analyses. The final dataset was rarefied by averaging 100 replicated count tables that were subsampled to 3500 read counts.

Statistical analysis

All analyses were performed in R (v4.4.2, https://www.r-project.org/; see for rendered scripts https://github.com/gbonthond/flatfish_microbiota). Figures were assembled and refined using Affinity Photo 2 (v2.6.0, https://www.affinity.studio/photo-editing-software). We considered a total of eight host and environmental variables. The host variables included species, sex, age (in years), length (in cm), total weight (in grams) and Fulton’s condition factor (total weight/length^3 41). The weight and age were transformed with a natural log and the median grain size with a log-base 2. The environmental variables that were considered for the analysis were the median grain size (in µm) as a proxy for habitat type and the trawling intensity in swept area ratio (SAR) per year. Median grain size data was obtained from Bockelmann et al.⁴². As a measure for the trawling intensity, we used fishing intensity of the subsurface (≥ 2 cm) penetrating gears from the OSPAR data & information management system⁴³. Based on spearman correlation coefficients, which revealed a strong correlation between length and weight (Spearman’s rank coefficient = 0.99, Figure S1, and see Figures S2-5 for more details on the host variables), it was decided to use the variable weight as an indicator of size and exclude length from downstream analyses.

To evaluate the impact of the remaining predictors on the diversity of the intestinal microbial community, we calculated the OTU richness and the effective number of OTUs⁴⁴. On both responses, we fitted linear mixed models (LMMs) using the R package lme4 (v1.1–35.5⁴⁵). First a global model was fitted with the same structure as used for the PERMANOVA (i.e., the main effects and all possible interactions with the factor species identity) and with the station identity as a random effect. Model assumptions were assessed using diagnostic plots. Model selection was then performed for both response variables by comparing the global model to all possible simpler models, keeping the random effect fixed. The best model was selected based on the corrected Akaike Information Criterion (AICc).

The community composition of the intestinal microbiota was analyzed with PERMANOVA⁴⁶ with the adonis2 function in the R package vegan (v2.6.8⁴⁷) based on Bray–Curtis distances and 9999 permutations. To test whether predictors had general effects applicable to all species, as well as species-specific effects, all seven predictors were included along with all possible second order interactions with the factor species identity. To visualize community similarity patterns, we conducted non-metric multidimensional scaling (nMDS) based on Bray–Curtis distances, using the R package vegan (v2.6.8).

Subsequently, a differential abundance analysis was conducted to identify microbial markers for each of the predictors. First, the OTU dataset was reduced to OTUs with at least 1% occurrence and at least 0.1% relative abundance. Second, a multivariate generalized linear model (mGLM) was fitted on the reduced community matrix with the package mvabund (v4.2.1⁴⁸). This model included all seven predictors and assumed a negative binomial distribution. Third, the estimated model coefficients were used to identify marker OTUs for each predictor, with coefficients considered significant if the respective 95% confidence intervals did not overlap with zero. To identify species-specific OTUs, the mGLM was fitted three times, each time using a different species as the reference level, allowing us to obtain pairwise estimates of species differences. The two coefficients and standard errors of the differences with both other species were then pooled to obtain estimates and confidence intervals for species-specific OTUs.

Ethics statement

The fishes collected in this study (Limanda limanda, Pleuronectes platessa, and Buglossidium luteum) are not protected under European or German nature conservation legislation and are not considered threatened or endangered. All specimens were obtained from a routine beam-trawl survey forming part of a legally mandated and authorized long-term fisheries monitoring program conducted by the Thünen Institute of Sea Fisheries. This monitoring is carried out under the authority of the Federal Office for Agriculture and Food (BLE), the competent federal authority for the oversight and authorization of federal fisheries research and monitoring under German sea fisheries legislation (Seefischereigesetz) and the framework of the EU Common Fisheries Policy (CFP), using standard commercial fishing gear in areas open to commercial trawling. All individuals used in this study were already dead when they were handed over for tissue sampling and were therefore handled in accordance with relevant guidelines and regulations. No experimental procedures or manipulations beyond standard fisheries monitoring were performed. In accordance with EU Directive 2010/63/EU and German animal welfare legislation (Tierschutzgesetz; Tierschutz-Versuchstierverordnung), such routine fisheries monitoring does not constitute an animal experiment within the meaning of the legislation and therefore does not require ethics approval or informed consent.

After quality filtering a total of 162 samples remained. These belonged to 65 individuals of the species B. luteum, 51 of L. limanda, and 46 of P. platessa, and included 96 females and 66 males. The rarefied community matrix counted 12,676 OTUs. Pseudomonadota is the most abundant phylum followed by Actinomycetota, Planctomycetota, Bacillota, and Verrucomicrobiota. At the family level, Pirellulaceae, Ilumatobacteraceae, Paracoccaceae, Hyphomicrobiaceae, and Actinomarinales were the five most abundant groups, accounting for 36.7% of the total OTU abundance (Fig. 2).

While the intestinal microbiota were generally dominated by a combination of taxa from the phyla Pseudomonadota, Planctomycetota and Actinomycetota, some samples deviated from this pattern and only contained few taxa (e.g., Mycoplasmataceae or Hyphomicrobiaceae), or taxa that were not among the 30 most abundant families.

Diversity

From the comparison of models fitted on the effective number of OTUs (ENOs) the best model yielded only trawling intensity as informative predictor (Table S1). However, the effect of trawling intensity on ENO was not significant (χ²₁ = 3.06, p = 0.080, Table S2) and the model explained < 2% of the variation (R²m = R²c = 0.019). In contrast, for OTU richness, the best model (R²m = 0.151, R²c = 0.225) included the variables species (χ²₂ = 7.33, p = 0.026), age (χ²₁ = 5.05, p = 0.025), weight (χ²₁ = 3.89, p = 0.049), the median grain size (χ²₁ = 7.77, p = 0.005) and the interaction between the median grain size and species (χ²₂ = 6.74, p = 0.034, (Table S1-S2). Post-hoc pairwise comparisons (adjusting p-values to control for familywise error rates with the Holm method) among species, revealed significant differences in OTU richness between B. luteum and P. platessa (t₁₆₇ = −2.525, p = 0.030) and between L. limanda and P. platessa (t₁₅₈ = −2.605, p = 0.030), and identified that of species-specific changes in OTU richness with the median grain size, only L. limanda was significant (t_116.9 = −3.209, p = 0.005, Fig. 3, Table S3).

Community composition

Non-metric multidimensional scaling (nMDS) did not reveal strong clustering patterns, but indicated that community composition varied with age, condition factor, median grain size and trawling intensity (Fig. 4A), and subtly differed among species (Fig. 4B).

These patterns were confirmed by the PERMANOVA analysis, which could explain only 16.9% of the overall variation in community composition, but resolved the variables species (F_2,141 = 1.678, p < 0.001, R² = 0.020), condition factor (F_1,141 = 1.561, p = 0.015, R² = 0.015), the interaction between condition factor and species (F_2,141 = 1.311, p = 0.036, R² = 0.036), age (F_2,141 = 1.834, p = 0.003, R² = 0.011), median grain size (F_1,141 = 3.621, p < 0.001, R² = 0.021) and trawling intensity (F_1,141 = 3.488, p < 0.001, R² = 0.021) as significant predictors of community composition.

Differential abundances

The differential abundance analysis (Fig. 5A) identified various OTUs that were negatively or positively associated with trawling intensity (negative: 29, positive: 24, Fig. 5B), and median grain size (negative: 12, positive 35, Fig. 5C). The host variable with the highest amount of differentially abundant OTUs was the condition factor (negative: 13, positive: 11, Fig. 5D), followed by age (negative: 3, positive: 10, Fig. 5E). For B. luteum (Fig. 5F), L. limanda (Fig. 5G), and P. platessa (Fig. 5H), the relative abundance of six, 10, and 10 OTUs, respectively, increased exclusively with each species. Additionally, 16 OTUs showed significantly lower relative abundance in B. luteum, 12 in L. limanda, and two in P. platessa.

This study provides new insights into the intestinal microbial community of three flatfish species from the Sylt Outer Reef (southeastern North Sea), highlighting the influence of environmental and host predictors on bacterial diversity and community composition. Overall, the predictors evaluated in this study were able to account for a portion of the total variation, but our models also indicated substantial variation to remain unexplained, indicating that other variables, not examined in the present study, importantly contribute as well to microbial composition and diversity in the fish intestines. The significance of sediment properties (i.e., measured by the median grain size) and trawling intensity emphasizes the important role of environmental factors, which is consistent with previous studies^{10,49,50,51,52}. Additionally, we identified several host variables (species, age, condition factor and weight) that contribute to diversity as well as to community composition. However, while we expected that intestinal microbiota would vary between sexes^53,54,55, our study found no significant differences between females and males on either diversity or community composition, suggesting that other factors, such as environmental influences or host genetics, play a more important role in shaping the microbial community.

It is important to emphasize that an interplay of different factors shapes the intestinal microbial community^14,56,57. Some environmental variables such as temperature^58,59, salinity^18,60,61, water microbiota^62,63, pollutants^16,64,65,66 and diet^15,17,54,67, were not considered in this study, although they are known to affect the intestinal microbiota of fish. While at the scale of our sampling area, bottom temperature and salinity vary minimally, they are likely to play an important role at larger spatial scales, across seasons or in coastal zones with substantial freshwater inflow. Diet is known to have a primary influence on the diversity and composition of intestinal microbiota^15,17,54,67. While variation in diet may be partially captured by both host and environmental variables, we did not have direct information of dietary intake of individual flatfishes, which makes other dietary factors likely candidate sources for the variation that this study could not explain. The migratory behavior of flatfish^68,69 further complicates efforts to link environmental factors to variation in their intestinal microbiota, since the location where a fish is captured does not necessarily represent the habitats it has occupied throughout its life. Consequently, life-history influences on the microbiota, such as exposure to different environments, anthropogenic impacts such as contaminated sediments and trawling history, or differences in diet, likely contribute to the variation that our models could not explain^3,64,70.

The presence of several similar bacterial taxa in the intestinal microbiota of one or more fish species from different populations suggests that these taxa play an important role in host intestinal functions⁷¹. As reviewed in the study by Rombout et al.², the phylum Pseudomonadota is the most common in the intestinal microbial community in fish, as it is in B. luteum, L. limanda, and P. platessa in this study comprising about 29% of all reads. Another review by Ghanbari et al.⁷² found that the phyla Pseudomonadota, Bacteroidetes, and Bacillota make up to 90% of the fish intestinal microbiota. For the flatfish species examined here, these phyla account for only ~ 40% of all amplicon reads. Instead, the phyla Actinomycetota and Planctomycetota make up a significantly larger share, with the most abundant three phyla together compromising for approximately 70% of the reads.

Intestinal microbiota vary among flatfish species

It is well known that the intestinal microbiota are species dependent and often vary even among closely related species^52,57,73,74. For instance, Miyake et al.⁷³ showed that intestinal microbial community composition differed even among closely related species within the same genus (Acanthurus). The three species analyzed in this study are more distantly related, all belonging to the same order (Pleuronectiformes): P. platessa and L. limanda are members of the family Pleuronectidae, while B. luteum belongs to the family Soleidae⁷⁵. Given that evolutionary distance is positively correlated with differences in intestinal microbiota, it is expected that these species have distinct intestinal microbiota¹⁵. Moreover, Smith et al.⁵² examined the intestinal microbiota of different populations of three-spined stickleback (Gasterosteus aculeatus) and found that its composition and diversity varied, with more genetically divergent populations having more distinct intestinal microbiota. Our results are in line with this, showing that both community composition and OTU richness varied among the three flatfish species. In addition, we found that changes in composition with the condition factor were species-specific, and changes in OTU richness with sediment type were different among species.

Our study identified six OTUs that were significantly more abundant in B. luteum, and 10 OTUs in both P. platessa and L. limanda. Interestingly, one of the L. limanda specific OTUs was classified to the genus Endozoicomonas, which is recognized as a diverse symbiont, found across a wide range of marine animals, including sponges, corals, mollusks, and fish⁷⁶. While the flatfish species studied here are benthic predators, there are subtle dietary differences. A study by Schückel et al.⁷⁷, examined the dietary overlap among four flatfish species, including the flatfish studied here, throughout the German Bight and found differences in prey selection. This prey resource partitioning among species likely contributes to the differences in prokaryotic community composition and diversity found here among species. Beyond species-specific dietary differences, other unique traits may further contribute to the distinct composition of the intestinal microbiota. Even though all three species live in marine demersal habitats, only L. limanda and P. platessa tolerate brackish water. Since salinity is an important factor influencing the composition of the microbiota this may contribute to differences found in intestinal community composition within species^18,60,61. Moreover, the flatfish species studied here have different migratory behavior^68,69, which influences both physiological changes and environmental exposure. As a result, migration history may therefore also contribute to differences among species^78,79,80.

Microbiota change with age in composition and become more diverse

Microbial community assembly of the intestine in fish can be divided into colonization and persistence⁵³. The interplay between colonization and persistence begins immediately after hatching when microbes start to colonize the intestine^81,82,83. The order in which different microbiota colonize the intestines can have a decisive influence on the developing community, a phenomenon referred to as the priority effect^83,84. In addition to the skin or gills, the intestine acts as an entry point, allowing prokaryotes in the surrounding environment and the first food to colonize it⁸⁵. With each successive food or water intake, the intestinal microbiota of the fish becomes more different from the environment⁸⁶ and diversifies further⁸⁷. Accordingly, the diversity of intestinal microbiota tends to increase with age^15,80,88. However, Yan et al.⁸⁹ found that diversity decreased with age in freshwater fish species. Especially in early life phases, environmental factors can have a particularly strong and long-lasting impact on the composition of the microbiota⁸⁶. Considering that fish are exposed to various factors throughout their life that can influence their bacterial community, it is natural to expect a gradual and consistent change with age. Our results are in line with many other studies, who found that microbial diversity in the intestine increases with age in fishes^19,88, as well as in humans³. In the studied flatfish species, we found a similar increasing trend in OTU richness (but not on the effective species number). Furthermore, community composition shifted with age: three OTUs were significantly less abundant in older individuals, while ten OTUs increased with age.

Intestinal microbiota do not differ between females and males

Our data showed no evidence for sex-related differences in diversity. These findings are somewhat contrasting to other studies on fish intestinal microbiota where more pronounced differences have been detected between females and males^53,54,55. In a study by Martyniuk et al.⁹⁰ a comparison between females and males in zebrafish (Danio rerio) found also no differences in diversity between sexes, but several individual taxa were more common in either males or females. For example, male zebrafish showed higher abundance in the families Erythrobacteraceae and Lamiaceae⁹⁰. In another study by Li et al.⁹² on largemouth bronze gudgeon (Coreius guichenoti), certain taxa were found to differ in abundance between sexes. Pseudomonadota was the most abundant phylum in males, whereas females exhibited dominance of five different phyla. Additionally, Li et al.⁹¹ found significant differences in diversity. Such sex-associated differences may be explained by dietary variations between females and males¹⁷, but could also be driven by physiological factors, such as differences in hormone production^92,93.

Richness decreases with size and composition changes

Few studies have examined the relationship between size and the intestinal microbiota in fish. Since length and weight were highly correlated, we used weight as a proxy for size. Size is linked to other physiological variables such as length or age. Even though a fish’s weight increases with age, especially in early life stages, growth rate reduces with age. This is one of the reasons that makes size and age different from each other and potentially impact the intestinal microbiota differently. Our study was able to disentangle these effects to some extent, as we found that OTU richness increased with age, whereas a decreasing trend in richness was found for weight, after correcting for age. This result, however, contrasts with the general observation made across vertebrates⁵⁷ that larger organisms, which consequently have a longer intestine, exhibit higher intestinal microbial diversity. Likewise, a study by Zhao et al.⁹⁴ on Gymnocypris chilianensis found that diversity was higher in larger individuals (> 300 g) than in smaller ones. These results align with the island biogeography theory⁹⁵, which has been linked to diversity of intestinal microbiota, predicting that as an isolated ecosystem, the intestine, can support greater diversity as it increases in size⁹⁶. Considering the different habitats in the intestines that are associated with distinct microbial communities, this correlation gets even more evident⁹⁷. First of all, the intestine is composed of the foregut, midgut, and hindgut, each hosting distinct microbial communities^{56,98,99,100,101}. Furthermore, allochthonous microbiota are transient and associated with digesta, while autochthonous microbiota colonize the mucosal surface, forming the core community. Interestingly, our study suggests a size-diversity relationship for flatfish intestinal microbiota that is opposite from what island theory predicts. Exceptions, however, exist as for instance a study on largemouth bronze gudgeon found no effect of weight on diversity⁹¹.

Fulton’s condition factor is associated with compositional changes

In this study, we demonstrated a significant effect of the variable ´condition factor´, as well as an interaction between species and condition factor on community composition. However, no significant effect was found on diversity. Among all studied host variables, condition factor was found to have most differently abundant OTUs consisting of 13 negatively differently abundant OTUs and 11 positively differently abundant OTUs. While weight provides a measure of size, the condition factor, as a ratio of weight to length, offers a relative measure that reflects a fish’s health and food availability¹⁰². Also in this case, diet, which is known as an important factor influencing the composition of the intestinal microbiota across vertebrates, and specifically in fish, may be associated with the condition factor. In particular, food availability, or more precisely periods of starvation, are a prime example that affects the condition factor and the intestinal microbiota simultaneously. Consequently, Xia et al.⁶⁷ observed shifts in bacterial communities in Asian seabass (Lates calcarifer), with a significant enrichment of Bacteroidetes and a significant depletion of Betaproteobacteria as a result of starvation. Although comparable studies that analyze the relationship between the condition factor and the intestinal microbiota in fish are rare, community composition of the human intestinal microbiota has been reported to vary with the body mass index, which can be considered equivalent to Fulton’s condition factor^103,104. The relationship between the intestinal microbiota and the condition factor is reciprocal, as some prokaryotes of the intestinal microbiota are able to break down complex sugars and provide essential short-chain fatty acids and energy as well as other nutrients, directly influencing the nutritional and health condition of the fish¹⁰⁵. Due to this fundamental importance of the microbiota for the condition of a fish, the intestinal microbial diversity has been used a biomarker for fish health and metabolic capacity⁶³.

Microbial composition and diversity vary with sediment type

With their bottom-dwelling life-style and a diet, which largely consists of benthic invertebrates⁷⁷, flatfish are in close contact with both the sediment and its microbiota. Sediment properties (i.e., median grain size, mud content and organic matter content) are primary drivers of microbial community composition and diversity of the top sediment layer³³. Moreover, also the organisms that make up the flatfish diet and pass through their intestines vary across habitat types and are thus strongly dependent on sediment properties^106,107,108. Using the median grain size as a simple measure for sediment type (e.g., mud, sand or gravel), we found that microbial diversity and composition within the flatfish intestine is indeed linked with local sediment properties.

That the sediment is an important source of microbiota found in the intestine is also suggested by others. Although grass carp primarily inhabit the mid to upper water column⁸, found that their intestinal microbiota composition mainly originates from the surrounding water and sediment. Species dwelling closer to the bottom sediment, such as flatfish, may be even more influenced by it. Besides the highly significant association with OTU richness and community composition, the median grain size also yielded a large number of differentially abundant OTUs, with 12 OTUs decreasing with median grain size, and 35 OTUs increasing.

Intestinal microbial community composition varies with trawling intensity

We found that trawling intensity explained small but significant changes in intestinal community composition. In total, we detected 24 OTUs to increase with trawling intensity, while 29 OTUs decreased. While trawling activity does not impact the isolated microbial ecosystems in the intestines of individual flatfishes directly, we propose three potential indirect pathways that may explain this observed trend.

First, mobile bottom-contacting fishing gears impact the environment physically. Beam trawls penetrate several centimeters deep into the seafloor, resuspending large amounts of sediment and organic matter, and modifying the seabed morphology¹⁰⁹. This affects both the water column and the seafloor, with which flatfish live in close contact.

Second, trawling alters benthic faunal communities. In the North Sea, up to 70% of benthic invertebrates, including bivalves, polychaetes, echinoderms or ophiuroids, die when dragged by a trawl (reviewed in Eigaard et al.¹¹⁰; Sciberras et al.¹¹¹). Many of these animals naturally influence the remineralization of organic matter and the regeneration of nutrients by microorganisms^112,113, but also make up the flatfish diet, which importantly impacts intestinal microbiota composition and diversity^57,101. Link et al.²⁸ analyzed dietary data of flatfish in the Northwest Atlantic over 25 years and found that the average weight of stomach contents of flatfish decreased in heavily fished areas. This supports that trawling can influence the diet of flatfishes. Moreover, besides changes in the identity of organisms in the flatfish diet, periods of reduced food availability and starvation are known to alter the intestinal microbial community⁶⁷. Therefore, trawling driven changes in flatfish diet may impact the intestinal microbiota in different ways and present a possible cause of the here observed association between intestinal microbial community composition and trawling intensity.

Third, the demersal flatfish species studied here live in close contact with the sediment, which also acts as a source of microbes that colonize their intestines⁸. Sediment microbiota have been shown to vary with bottom trawling effort as well, showing a decrease in alpha diversity and change in overall community composition^33,34. Therefore, trawling related shifts in benthic microbiota may offer another potential pathway through which trawling intensity could indirectly impact the intestinal microbial community of flatfish.

While physical alterations to the seabed morphology, changes in faunal communities and therefore diet, and shifts in sediment microbiota, present interesting hypotheses for how trawling may indirectly affect flatfish intestinal microbiota, they currently remain speculative. Nonetheless, they currently offer the best explanation for the observed association between bottom trawling intensity and intestinal microbiota composition. Moreover, they highlight the need for further investigation and serve as a hypothetical basis for future research on the effects of environmental disturbance caused by fishing activities on fish health and fish intestinal microbiota in particular.

While our study detected a trawling signal in flatfish intestinal microbiota, future work is needed to resolve the causal pathways underlying this pattern. Because trawling operates at large spatial scales and interacts with environmental conditions and fish life history, we recommend sampling designs including gradients of trawling intensity (as in the present study), while also covering nursery and migratory grounds and dietary analyses to evaluate how migration history and diet may mediate trawling effects on intestinal microbiota.

Here, we disentangled how species identity, age, size, condition factor, and environmental factors (i.e., sediment type and trawling intensity), contribute to shaping the intestinal microbiota of the demersal flatfish species B. luteum, L. limanda and P. platessa in the southeastern North Sea. The strong effects of sediment type and minor, but significant, effects of trawling intensity indicate that the environment plays a key role in shaping the intestinal microbiota of these flatfishes. To the best of our knowledge, the present study is the first to detect an association between fish intestinal microbiota and bottom trawling intensity. This adds to previous work, that found benthic microbiota to vary along a trawling gradient³³, and may hint that such effects could extent to higher trophic levels. However, we note that substantial variation could not be explained by our models, indicating that other processes, not captured by the variables examined here, importantly influence fish intestinal microbiota as well. Dietary factors likely account for a substantial portion of the unexplained variation in both composition and diversity of the flatfish intestinal microbiota. Another important factor that we were not able to account for is the migratory history of the flatfish studied^69,79.

While many studies have focused on captive fish due to their economic importance in aquaculture, wild populations remain understudied¹⁴. Given that Ramírez & Romero ⁵¹ found differences in intestinal microbial communities between wild and captive fine flounder (Paralichthys adspersus), and Xie et al.⁵⁷ documented such differences across vertebrates in a meta-analysis, the focus of this study on wild populations contributes to narrowing this knowledge gap. Furthermore, this study helps to compensate for a geographic bias, as most studies on this topic to date have come from North America or East Asia, while Central Europe and Africa tend to be underrepresented¹⁴.

These findings contribute to our understanding of how host variables as well as environmental and anthropogenic processes may directly or indirectly affect host-associated microbial communities, with potential ecological and evolutionary implications. As this study is the first to observe an association of trawling effort and microbial community composition of flatfish intestines, it merits for more research to identify the mechanisms that underly this trend, gain insight into the long-term impacts, and potential consequences for fish health and population dynamics.