Skip to main content
Download PDF

Interactions between gilthead seabream intestinal transcriptome and microbiota upon Enteromyxum leei infection: a multi–omic approach

Animal Microbiome volume 7, Article number: 22 (2025) Cite this article

Abstract

Background

The enteric myxozoan parasite Enteromyxum leei is an important problem in gilthead seabream aquaculture invading the intestinal epithelium and leading to chronic intestinal inflammation, poor food conversion rates, cachexia, and mortalities, with no treatments available, resulting in significant economic losses. It is known that myxozoan infections are affected by factors such as temperature, duration of exposure, stocking densities, and seasonality. Gut microbiota has key effects on host health, including disease resistance and immune system training and development, tightly interacting with the host, affecting systemic and local physiological functions. This study aimed to gain insights into the host–microbiota–parasite interactions integrating metataxonomics, host transcriptomics, and metatranscriptomics within this disease model.

Results

Exposure to E. leei together with temperature and age differences led to alterations in gilthead seabream intestinal microbiota. Samples from 240 g fish kept at 18ºC during a winter trial at 10 weeks post–parasite exposure showed the highest significant changes in their microbial composition with Proteobacteria increasing in abundance from 32.3% in the control group up to 89.8% in the infected group, while Firmicutes and Actinobacteria significantly decreased in relative abundance from 23% and 37.8–2.4% and 1.1%, respectively. After LEfSe analysis, Acinetobacter was identified as the best biomarker for the parasite–exposed group. Parasite exposure also altered the expression of 935 host genes, highlighting genes involved in immune responses such as pathways related to Interleukins, MHCI and Interferons. Microbial transcripts, also showed significant changes upon parasite infection. Integration of the results revealed differential effects on the host induced directly by the parasite or indirectly by parasite–induced microbial shift.

Conclusions

Intestinal microbiota and local host gene expression showed significant changes upon en enteromyxosis. The detected activation of the host immune response was not exclusively linked to the parasite infection but also to changes in microbiota, demonstrating the key role of the different components of the mucosal system during disease. These results provided different datasets of bacterial taxa and microbial and host transcripts that will allow a better understanding of host–microbiota–parasite interactions and can serve as starting points for studying and evaluating mucosal health in aquaculture during parasitosis or other diseases.

Background

Gilthead seabream (Sparus aurata) is one of the most important species for Mediterranean aquaculture constituting together with European seabass (Dicentrarchus labrax), the majority of finfish aquaculture production (79% by 2021) [1,2,3]. However, sub–optimal aquaculture practices pose several concerns, such as disease outbreaks that fish farms often experience resulting in economic losses [4].

In particular, the enteric myxozoan parasite Enteromyxum leei represents a threat to gilthead seabream farming, and it has already led to the abandonment of the production of aquacultured valuable fish such as red porgy (Pagrus pagrus) and sharp snout seabream (Diplodus puntazzo) at some sites [5]. This parasite is responsible for enteromyxosis, a slow–progressing disease starting in the posterior intestine, spreading to other gut segments [6], affecting the intestinal epithelium and producing inflammation, anorexia, cachexia, and eventually death [7]. Myxozoans’ life cycle typically presents a complex cycle involving two hosts, an invertebrate and a vertebrate, as definitive and intermediate hosts. The invertebrate host of E. leei is still unknown, but direct fish–to–fish transmission is also possible increasing the risk of propagation in a farming environment [7, 8]. To date, E. leei cannot be cultured in vitro and no effective preventive or curative treatments are available against enteromyxosis [9]. Nevertheless, in vivo experimental infection models have contributed to filling the knowledge gap on how biotic and abiotic factors interact in this fish–parasite model. Thus, it is now known that factors such as high temperature, increased duration of exposure, and elevated stocking density enhance enteromyxosis [7]. In addition, previous myxozoan–related studies have revealed that seasonality also affects parasitic abundance in water and infection frequencies in fish [10,11,12,13].

Nowadays research on fish microbiomes is exponentially increasing as the link between microbiota and host health has been strengthened by emerging evidence that microbiota plays a key role in the fish health status, including immune response, parasite resistance, food digestion and generally maintaining host homeostasis [4, 14, 15]. Hence, shifts in microbiome composition and activity can be used to detect changes in fish health [16, 17]. Consequently, applying an integrative approach between the sequencing of the 16 S ribosomal RNA (16 S rRNA) gene and different multi–omic methods (i.e. RNAseq, metatranscriptomics, metabolomics, metagenomics, etc.) offers an improved resolution of these communities and unveil host–microbiome–pathogen interactions that can be linked to host health [16]. Thus, this study aimed to gain insights into these interactions by integrating metataxonomics, host transcriptomics, and metatranscriptomics approaches within an E. leei infection model in gilthead seabream.

Materials and methods

Ethics statement

Animals were handled and trials were carried out according to the European Union Council legislation on the handling of experimental fish (Directive 2010/63/EU) and to the current Spanish Royal Decree RD 53/2013. All procedures were approved by the Ethics and Animal Welfare Committee of the Institute of Aquaculture Torre de la Sal (IATS–CSIC, Castellón, Spain), CSIC, and “Generalitat Valenciana” (permit number 2021/VSC/PEA/0194).

Infection trials and samplings

Clinically healthy juvenile gilthead seabream (Sparus aurata) were obtained from a commercial fish farm and transported to the IATS, CSIC facilities (40°5’N, 0°10’E). Upon arrival, a fish subsample (n = 20) was checked by PCR (18 S ribosomal RNA gene) and histological analyses [6, 18] to confirm it was free of intestinal parasites. Fish were kept in UV–treated 5 μm–filtered sea water (open flow), with natural photoperiod and temperature (unless otherwise stated), and fed daily ad libitum with a commercial diet (Biomar). Two experimental infection trials with the intestinal parasite E. leei were conducted using the same batch of fish to avoid differences due to the genetic background, one in summer (S) and the other in winter (W). Infections were performed by effluent exposure as previously described [8]. Briefly, in each experimental trial, 40 fish were distributed in two 500 L tanks to constitute the control (C, n = 20) and recipient (R, n = 20) groups. The R tank received the water from a donor tank containing 20 heavily parasite–infected fish of similar size. After 2 (t1), 5 (t2) and 10 (t3) weeks post–exposure (wpe) 6 C and 6 R fish were sampled. The summer trial was conducted between July and September 2021. Experimental fish had a mean weight of 45 g (± 1.95 SD) at the beginning of the summer trial and the temperature was 23.5 °C, 25.2 °C, 26.8 °C and 26.4 °C at 0, 2, 5 and 10 wpe, respectively. The winter trial was conducted between February and March 2022 and the mean fish weight was 240.4 g (± 35.3 SD). In this trial, the water temperature of all tanks was maintained at constant 18 ± 1 °C to achieve infection, as lower temperature arrests parasite establishment and growth [5].

Sampled fish were fasted for 24 h, sacrificed by overexposure to tricaine methanesulfonate (MS–222, 0.1 g/L; Merck), and spleens and whole intestines were dissected. Spleens were preserved in RNAlater (Invitrogen) for transcriptomic analyses. The posterior part of the intestine (~ 3 cm from the anal ampoule) was separated and a portion of ~ 0.5 cm was preserved in RNAlater. The rest of the posterior intestine was cut open, washed by flushing with sterile PBS to remove intestinal contents and non-adherent material, and mucus was scrapped out from the PBS-soaked tissue using the blunt end of a sterile scalpel. Approximately, 200 and 500 µl of mucus were obtained from summer and winter fish, respectively. Mucus samples were collected in 100 µl aliquots in Eppendorf tubes and immediately frozen at − 80 °C.

Nucleic acid extraction

DNA from 100 µl intestinal mucus was extracted using the High Pure PCR Template Preparation Kit (Roche) with a routinely optimized protocol including a previous lysozyme lysis step [19]. DNA concentration was measured with Nanodrop 2000c (Thermo Scientific) and samples were stored at − 20 °C. Intestinal mucus DNA was used for parasite diagnosis and for Illumina 16 S rRNA gene sequencing for microbiota analysis.

RNA was extracted from 100 µl of intestinal mucus or 60 mg of spleen or posterior intestine tissue. Mucus or tissue samples were homogenized in 0.5 ml of TRI Reagent solution (Invitrogen) using microbial lysis tubes (Qiagen, for the mucus samples) or Lysing Matrix D tubes (MP Biomedicals, for the tissues). Homogenization was performed in a FastPrep homogenizer (MP Biomedicals) using 1 cycle of 30 s at 6 m/s. Then, total RNA was isolated using the MagMAX™-96 for Microarrays Total RNA Isolation Kit (Applied Biosystems) according to the manufacturer’s instructions. RNA concentration and purity were measured using a Nanodrop 2000c (Thermo Scientific). Finally, RNA quality and integrity were measured on an Agilent Bioanalyzer 2100 total RNA Nano series II chip (Agilent). RNA integrity number (RIN) values were between 8 and 10. Samples were stored at − 80 °C until sequencing. Tissue RNA was used for host transcriptome Illumina sequencing, and mucus RNA was used for microbial metatranscriptome Illumina sequencing.

Parasite diagnosis

Parasite 18 S rRNA (DQ448298.1) gene copies in intestinal mucus DNA were estimated by qPCR as previously described [9, 18, 20]. Copy numbers were interpolated from the cycle thresholds (Ct values) using a standard curve generated with known concentrations of a plasmid containing the target gene providing the information of parasite infection intensity. The standard dilutions covered six orders of magnitude. All samples from the same trial were run in the same plate with two DNA sample dilutions. Samples with Ct < 38 were considered positive.

Illumina 16 S rRNA gene sequencing for microbiota

For characterizing the intestinal adherent microbiota, the V3-V4 regions of the 16 S rRNA gene were sequenced by the Illumina MiSeq platform (2 × 300 paired-end run) at the Genomics Unit from the Madrid Science Park Foundation (FPCM, Spain). A total of 60 samples were sequenced, five fish per group (C and R), sampling (2, 5, 10 wpe) and trial (summer and winter). From the six sampled fish, the five with higher DNA quality and concentration were chosen for microbiota sequencing. Details of PCR and amplicon sequencing can be found elsewhere [19]. After sequencing, quality–filtering and preprocessing steps were performed using Prinseq [21]. Then, clean forward and reverse reads were merged with fastq-join [22], and finally, reads were assigned as distinct amplicon sequence variants (ASVs), and then taxonomically assigned with Minimap2 v2.17-r941 [23] using SILVA v138.1 [24] as a reference database.

Illumina RNA sequencing of host transcriptome

The results from the microbiota study showed that the group with the largest differences was from the winter trial at 10 wpe. Therefore, samples from this sampling point (n = 6 C and n = 6 R fish) were chosen for the subsequent analyses. Host RNA from spleen and posterior intestine were sequenced by GeneWiz (Azenta, Germany) on the Illumina NovaSeq 6000 platform with a 2 × 150 nucleotides paired-end (PE) read format. Raw reads were filtered by quality, length and percentage of Ns using Prinseq [21] (-min_len 31, -min_qual_mean 25, -ns_max_p 10, -trim_qual_right 20). Clean reads were mapped against the CSIC gilthead seabream annotated reference genome [25] using hisat2 v2.0.5 [26]. Mapped read counts were quantified using featureCounts v1.5.0-p3 [27].

Illumina metatranscriptome sequencing

Using the mucus RNA samples from the winter trial at 10 wpe, rRNA was removed using the Illumina Ribo-Zero Plus rRNA Depletion Kit (Illumina), which targets both eukaryotic and prokaryotic rRNA. Then, Illumina RNA-seq libraries were prepared by NovoGene (UK) from 500 ng of total ribo–depleted RNA using the Illumina TruSeq™ Stranded Total RNA Library Prep Kit (Illumina) according to the manufacturer’s instructions. RNA-seq libraries were sequenced on an Illumina NovaSeq 6000 platform as a 2 × 250 nucleotides paired–end (PE) read format. Quality filtering was performed with Trimmomatic v0.40 [28] to remove adaptor contamination and reads with > 10% Ns or mean sequence quality < 20. Obtained reads were pre–processed using SortMeRNA [29] to filter out left over rRNAs. Then, samples were mapped against the Ensembl gilthead seabream annotated reference genome (fSpaAur1.1) using STAR 2.7.11b [30]. Unmapped reads (putative microbial reads) were retrieved with SAMtools [31] and assembled from each metatranscriptomic sequence using megahit [32] with default parameters. Finally, retained reads were mapped against the UniRef 90 database using HUMAnN 3.0 [33]. Unique descriptions of bacterial genes with at least 10 counts in at least one sample were used to perform a differential expression analysis using DESeq2 [34]. Differentially expressed transcripts were considered at padj < 0.1.

Data analysis, statistics, and visualizations

For microbiota analysis, R packages phyloseq, vegan, and mixomics [35,36,37] were used. Differences in richness and diversity indexes were determined by the Kruskal–Wallis test using Dunn’s post–test (P < 0.05) and beta diversity was analyzed using permutational multivariate analysis of variance (PERMANOVA) with the non–parametric method adonis and 10,000 random permutations. Nonmetric multidimensional scaling (NMDS) using Bray-Curtis dissimilarity was used to show differences among groups. Linear discriminant analysis (LDA) effect size (LEfSe) was used to identify significant microbial biomarkers among groups at genus level [38]. For normalization, sample depths were normalized by total sum scaling and then made proportional to the total sequencing depth.

For the RNA sequencing data, differentially expressed (DE) transcripts between C and R groups were detected using DESeq2 [34] with a significance level of padj < 0.05. Unique annotated DE transcripts were used to perform pathway enrichment analysis using the Reactome database [39] after converting the gilthead seabream identifiers into their human equivalents, when possible. Although this database is mainly constructed for human genes, a large proportion of our input genes were recognized for the analysis (71%). Therefore, despite the limitations, the use of this tool is considered adequate to obtain an overview of the enriched biological pathways, in lack of a fish-specific implement. Pathways were considered enriched at FDR < 0.05. For visualization, the relationships between enriched categories or pathways according to their shared transcripts were performed using the runGSA function of the piano R package [40] and the resulting networks were constructed using Cytoscape v3.8.2 [41].

Correlation analyses between parasite intensity, microbial biomarkers and DE immune gene counts (C and R, n = 10), or between metatranscriptome DE genes and DE immune gene counts (C and R, n = 12) from each individual for which all correlated data were available were calculated by pairwise Spearman correlation coefficients. Data were normalized using centered log-ratio (clr) transformation and the Spearman correlation coefficient, P values and Benjamini–Hochberg adjusted P values (FDR) were calculated using the cor-test function of the corrplot R package [42]. Significant correlations were considered at FDR < 0.1. Correlation networks were visualized using Cytoscape v3.8.2 [41].

Data Availability

Raw sequence data were uploaded to the Sequence Read Archive (SRA) under Bioproject accession number PRJNA1134388, BioSample accession numbers SAMN42418970-SAMN42419029 for the 16 S rRNA amplicon sequencing, SAMN42419030-SAMN42419053 for the host RNA sequencing, and SAMN42419054-SAMN42419065 for the metatranscriptome sequencing.

Results

Infection results

To determine the infection status of the experimental fish, we performed a qPCR of the E. leei 18 S rRNA gene in the extracted DNA from the posterior intestine of all fish (n = 6 C and n = 6 R in each sampling point). In the summer trial, exposed animals were already infected at 2 wpe, with a prevalence of infection of 83.3%, 100% and 66.7% at 2, 5 and 10 wpe, respectively. The intensity of infection in 18 S rRNA copy numbers ranged from 0 to 105, from 103 to 107, and from 0 to 106 at 2, 5 and 10 wpe, respectively. In the winter trial, no positives were detected at 2 wpe, and prevalence of infection was 100% at 5 and 10 wpe. Intensity of infection in copy numbers ranged from from 102 to 105. Intensity of infection was not significantly different among groups (Kruskal–Wallis + Dunn’s post-test, P > 0.05). All C fish were negative for the parasite, as no amplification of E. leei 18 S rRNA gene was detected. No mortalities were registered along the experimental trials. The individual information for the parasite diagnosis results can be found in Additional file 1.

Enteromyxum Leei exposure and trial impact on intestinal microbiota diversity and composition

After Illumina sequencing of 16 S rRNA V3-V4 amplicons, 5,910,205 clean reads were taxonomically assigned at a mean number of 98,503 reads per sample (Additional file 1). ASVs were assigned to a total of 1546 taxonomies up to the genus level. Alpha diversity indexes (Shannon and Simpson) showed significant differences among groups. However, richness estimates (Chao1 and ACE) did not show important differences regardless of the variable used to perform the comparison (Fig. 1A, Additional file 2). When considering all the fish groups and analyzing only the trial variable (all Summer fish vs. all Winter fish), alpha diversity was significantly higher in the winter trial than in summer (Additional file 2 A). Additionally, considering only the parasite exposure variable (all R fish vs. all C fish) R fish showed a significantly decreased diversity (Additional file 2B).

Considering the parasite exposure variable and trial (Summer vs. Winter trials and R vs. C), C fish from the winter trial had a significantly higher diversity than all the other groups, and no differences appeared between C and R in summer (Additional file 2 C). Finally, taking into account all variables (Summer vs. Winter, R vs. C and t1, t2, t3), C fish from the winter trial were more diverse than all R groups in summer, and no differences were found between C and R in the summer trial at any sampling point. Only at 10 wpe (t3), R fish from the winter trial (WR-t3) showed significantly lower diversity than all C fish (Fig. 1A, Additional file 2D).

Fig. 1
figure 1

Intestinal microbiota alpha and beta diversity. A Boxplot representing the Shannon diversity index of all groups (n = 5). Boxes represent the interquartile range (IQR) between the first and third quartiles and the horizontal line inside the box defines the median. Whiskers represent the lowest and highest values within 1.5 times the IQR from the first and third quartiles. Samples with a relative abundance exceeding those values are represented as points outside the boxes. Different letters indicate significant differences among the groups (Kruskal–Wallis + Dunn’s post-test, P < 0.05). B NMDS ordination plot based on Bray-Curtis distance matrixes of normalized sample counts showing community differences among samples. Groups showing a clearer separation are control (C) and parasite exposed (R) fish in the winter (W) trial at the 10 weeks post-exposure (wpe) sampling (t3) and are highlighted with black borders in the plot. t1, t2, t3 stand for the different samplings performed at 2, 5 and 10 wpe, respectively. W and S stand for the different trials performed in winter and summer, respectively. C and R stand for control and parasite-exposed fish, respectively

Full size image

Beta diversity analyses showed that there were significant differences among groups regardless of the variable (Table 1) except for the parasite exposure variable that yielded no significant differences. The group separation can be visualized in the NMDS plot shown in Fig. 1B. In the NDMS plot, the groups that showed the highest separation were C and R from the winter trial at 10 wpe (WC-t3 and WR-t3).

Table 1 Intestinal microbiota beta diversity (PERMANOVA) in the different comparisons
Full size table

In all groups, three phyla, Proteobacteria, Firmicutes, and Actinobacteriota, accounted for 85.9% up to 98.8% of the total abundance and showed significant changes among groups (Additional file 3 A). The most abundant phylum was Proteobacteria, ranging between 94.8% in the control fish of the summer trial at the last sampling point (SC-t3) and 32.2% in the control fish of the winter trial at the last sampling point (WC-t3). This phylum, while very prevalent in all groups of abundance in the winter trial with a tendency to increase upon parasite exposure (Fig. 2A). Firmicutes and Actinobacteriota were generally more abundant in control fish of the winter trial than in fish from the summer trial. These two phyla showed a significant decrease after long-term parasite exposure concomitant with the increase in Proteobacteria. At the family level, 15 families were detected at a mean abundance of > 10% in at least one of the groups (Additional file 3B). The most abundant families were Vibrionaceae, Moraxellaceae, and Pseudomonadaceae, with a global mean abundance of > 5%. The only family showing statistically significant differences among groups (Fig. 2B) was the most abundant, Vibrionaceae, with abundance values that ranged between 53% in summer and 9% in the winter trial, but with a significant increase in parasite-exposed fish in winter at 10 wpe (WR-t3), reaching summer abundance values. Other families showed high abundance variations but were not significant due to the high individual variability. For example, Moraxellaceae ranged between 22.7% and 0.3%, and Burkholderiaceae ranged from 23.1% to less than 1% in several sampling groups.

Fig. 2
figure 2

Intestinal microbiota taxonomic structure. Bar chart representing the mean relative abundance of bacterial phyla (A) and families (B) in each one of the sampling groups (n = 5). Only the 15 most abundant phyla or families are represented. Different letters indicate significant differences among the groups (Kruskal–Wallis + Dunn’s post-test, P < 0.05). Taxa presenting significant differences are marked with an asterisk in the legend. U. stands for unchacterized. t1, t2, t3 stand for the different samplings performed at 2, 5 and 10 weeks post-exposure, respectively. W and S stand for the different trials performed in winter and summer, respectively. C and R stand for control and parasite-exposed fish, respectively

Full size image

LEfSe analysis of all groups revealed four microbial markers at the genus level (Fig. 3A), Staphylococcus, Photobacterium, Cutibacterium, and an uncharacterized Burkholderiaceae genus. Staphylococcus and Cutibacterium were markers for both winter groups (WC, WR). In particular, Staphylococcus was more abundant in the WR-t2, WC-t1, and WC-t2 groups accounting for the 9.8%, 19.7%, and 19.6% of the total abundance (Additional file 3 C). Cutibacterium was more abundant in the WC-t1 accounting for 12.3% of the total abundance. On the other hand, Photobacterium and the uncharacterized Burkholderiaceae genus were identified as markers for both summer groups (SC, SR). Photobacterium was the most abundant genus in all summer groups, ranging between 68.8% and 19.3% of total sample abundance, while the uncharacterized Burkholderiaceae genus accounted for 23.14% and 6.33% of the total sample abundance from the SC-t1 and SR-t1 groups, respectively.

Fig. 3
figure 3

Microbial biomarkers. Linear discriminant effect size (LEfSe) analysis performed at the level of genus representing the significant bacterial biomarkers for each group (n = 5) and their abundance in normalized counts. The analysis was performed using all experimental groups (A) or only considering samples from the winter trial at 10 weeks post-exposure (wpe) (B) which were used for further analyses. Boxes represent the interquartile range (IQR) between the first and third quartiles and the vertical line inside the box defines the median. Whiskers represent the lowest and highest values within 1.5 times the IQR from the first and third quartiles. Samples with a relative abundance exceeding those values are represented as points outside the boxes. t1, t2, t3 stand for the different samplings performed at 2, 5 and 10 wpe, respectively. W and S stand for the different trials performed in winter and summer, respectively. C and R stand for control and parasite-exposed fish, respectively

Full size image

All previous analyses showed that samples from the winter trial at 10 wpe had the most striking and significant changes in the intestinal microbial community composition between control and treatment. Therefore, those samples (WC-t3 and WR-t3) were chosen for the subsequent analyses. At the phylum level, Proteobacteria increased in relative abundance from 32.3% in the WC-t3 group up to 89.8% in the WR-t3 group, while Firmicutes and Actinobacteria significantly decreased in abundance from 23% and 37.8–2.4% and 1.1%, respectively. After LEfSe analysis, seven microbial biomarkers were identified (Fig. 3B). Acinetobacter (Proteobacteria) was the only genus increasing in relative abundance in the WR-t3 group from 0.5 to 9.3%, while Streptococcus, Staphylococcus, Bacillus (Firmicutes), Pseudarthrobacter, Cutibaterium and Corynebacterium (Actinobacteria) decreased in relative abundance in the WR-t3 group.

Parasite exposure induced a significant change in intestinal transcriptome with only mild changes in spleen

The following results focus exclusively on the groups WC-t3 and WR-t3, which were chosen because they presented the largest differences in microbial composition.

Approximately 682 million paired-end reads were obtained from the RNA sequencing of posterior intestine and spleen whole tissue of control and recipient fish from the winter trial at 10 wpe (~ 28 million reads per sample) (Additional file 1). After bioinformatic analysis, ~ 446 million reads were mapped against the gilthead seabream genome [25], where unique hit counts were associated with 78,238 transcripts. Differential expression analysis revealed 14 differentially expressed (DE) transcripts in spleen (ten up-regulated, four down-regulated) and 935 DE transcripts in posterior intestine (788 up-regulated, 147 down-regulated) in recipient fish when compared to control ones (Additional file 4 A and B). Considering the low number of differences found in spleen, subsequent analyses were performed only with the intestinal data.

The 935 DE intestinal transcripts corresponded to 594 transcripts that were annotated to unique known-protein coding sequences, out of which Reactome pathway analysis recognized 422 transcripts. A total of 88 pathways were significantly enriched in the Reactome analysis (Additional file 4 C) with the highest number of genes belonging to pathways related to Immune system (131 transcripts), followed by Metabolism of RNA (77 transcripts), Disease (48 transcripts) and Cell cycle (35 transcripts) (Fig. 4A). Within the Immune system pathways, enriched pathways could be grouped in four specific categories: MHCI related, Interleukin related, Colony stimulated factor 3 related, and Interferon related, being the later the one with more enriched pathways (Fig. 4B). Remarkably, of the 131 transcripts belonging to Immune system pathways only 10 were downregulated in R fish (Additional file 4B), pointing to a general activation of the immune response.

Fig. 4
figure 4

Cytoscape networks showing the Reactome pathway analysis results. Pathway analysis of all differentially expressed genes (935) resulted in 88 overrepresented pathways that were grouped into 11 parent categories (A). The immune system overrepresented pathways are shown in detail in B. Numbers within nodes represent number of genes in each pathway or category. Edge thickness represent number of shared genes. The node color gradient in A represents the number of enriched pathways within each category. The node color in B represents different groups of pathways: Blue: Interleukin related pathways, Dark green: MHCI related pathways, Light green: interferon related pathways, Purple: colony-stimulating factor 3 signaling related pathways

Full size image

Changes in immune gene expression are not only linked to parasite intensity

A total of 218 significant correlations were found between microbial biomarkers or E. leei and the expression of the 131 transcripts related to the immune system. Microbial markers with significant correlations were Corynebacterium, Cutibacterium, Streptococcus and Staphylococcus, which significantly decreased in abundance in parasite-exposed fish. Considering that most immune genes were up-regulated upon infection, significant correlations between these bacteria relative abundance and immune genes were negative, whereas, the correlations between immune genes and the parasite intensity were mainly positive (Fig. 5A). The correlation network shows that the expression of 47 transcripts was significantly associated with both, bacteria and parasite. However, some transcripts were exclusively associated with the parasite (12 transcripts) or with bacteria (42 transcripts) (Fig. 5B, Additional file 4D).

Fig. 5
figure 5

Cytoscape network showing significant correlations between the expression of immune system related genes (blue squares), and parasite (green oval) or bacterial biomarkers (pink ovals) abundance (A). Red and blue lines represent positive and negative correlations, respectively. Thicker lines represent lower FDR values (FDR threshold < 0.1). B Venn diagram showing the number of immune related transcripts exclusively correlated to the parasite (green), to the bacterial biomarkers (pink), or shared (intersection)

Full size image

Enteromyxum Leei infection induced significant changes in intestinal metatranscriptome

Illumina sequencing of the 12 RNA intestinal mucus samples from the winter trial at 10 wpe (n = 6 C, n = 6 R) yielded a total of 534 M reads (∼44 M reads per sample) (Additional file 1). After quality filtering and an in silico rRNA removal step, 521 M reads remained.

Mapping against the UniRef90 database revealed a total of 29,581 transcripts, which corresponded to 757 unique descriptions with at least 10 reads in one sample. After DESeq analysis, a total of 24 bacterial transcripts were found to be differentially expressed, with only two transcripts down-regulated in the R group (Fig. 6A, Additional file 5). Thirteen of these DE transcripts showed significant associations with immune gene expression (Fig. 6B).

Fig. 6
figure 6

Bacterial metatranscriptome results. A Bar chart representing the differentially expressed bacterial transcripts in Log2 fold change (± standard error of the mean, lfcSE) when comparing parasite-exposed to control fish. B Cytoscape network showing significant correlations between the expression of immune system related genes (blue squares) and bacterial transcripts (pink squares). Red and blue lines represent positive and negative correlations, respectively (P < 0.01, correlation coefficient > 0.76)

Full size image

Discussion

Parasitic infections in fish have been linked to changes in the microbiome, leading to dysbiosis and secondary infections [17, 43]. Dysbiosis or microbial imbalance can lead to alterations in the microbial community function and modify the host-microbiota metabolic crosstalk, resulting in increased infection risk and eventually compromising host health [17]. In this study, our results showed that enteromyxosis and other variables, such as temperature and age are associated with shifts in the intestinal mucosal microbial populations, where significant differences were observed between trials (S-t1,t2,t3 vs. W-t1,t2,t3) and status of infected fish (C vs. R). Regarding the most abundant phyla, the results share similarities with previous research that had characterized gilthead seabream intestinal microbiota [44,45,46].

So far, there are no published studies concerning the crosstalk between mucosal microbiota and fish infected with myxozoan parasites. However, research involving the analysis of gut microbiota of diseased fish has been published, where overall, an increase in pathogen abundance and changes in the microbiome structure were the common responses to the infections [17]. Concerning fish parasitic diseases, some studies have been focused on gut microbiota [47,48,49], where little signs of dysbiosis caused by parasitism were found. However, Fu et al. [47] reported significant changes in the microbial composition at the genus level where the parasite Khawia japonensis (cestoda) showed a significant correlation with several bacterial genera. In other disease models, changes in the Proteobacteria/ Firmicutes ratios were observed in the intestine of diseased fish (suffering from bacterial infections) which had a higher abundance of Proteobacteria [50, 51]. Interestingly, in our study, even when Proteobacteria was prevalent in summer regardless of the parasite exposure, an increase of this phylum was observed during the winter trial (18 °C) along with parasite infection (WR-t3). Among other variables, seasonal changes and water temperature have previously been shown to influence gut microbiota composition and have highlighted some prejudicial impacts [52,53,54]. In particular, a study on Atlantic salmon (Salmo salar) found that higher temperatures were associated with lower abundance of lactic acid bacteria (Firmicutes) and the genus Acinetobacter, and higher abundance of the genus Vibrio (Proteobacteria) [55]. These findings were similar to the results found in gilthead seabream, where the family Vibrionaceae was the most abundant family in all summer groups and in the winter group with the longer parasite exposure (WR-t3). In addition, temperature also affects the immune responses in fish, since it can act as an environmental stressor weakening the host’s immune system [56]. However, we need to consider that, in order to have the same genetic background, fish from the summer and winter group were different size and age, thus these variables can also have an effect (discussed in more detail below).

In this study, several microbial markers were identified by LEfSe analysis, where the genus Photobacterium (Vibrionaceae family) and an unidentified genus of the Burkholderiaceae family were more abundant during the summer trial regardless of the infection status (SC, SR). On the other hand, the genera Cutibacterium and Staphylococcus were more abundant in the winter trial (WC, WR). Some Photobacterium species can act as mutualistic bacteria in the gut aiding with chitin digestion, nevertheless, many Photobacterium species can act as pathogens as well [53]. The family Burkholderiaceae has been reported to be present in the human gut, especially if there are immunosuppression signs [57] since it contains several opportunistic pathogens, as well as primary pathogens for humans and animals [58]. In addition, this family is associated with the degradation of complex organic materials [58, 59]. Concerning the genus Cutibacterium [60], it contains well-known opportunistic human pathogens usually present on the skin. However, it has been commonly identified as highly prevalent in the gut of pre-smolt Atlantic salmon and present in the intestinal microbiomes of farmed chinook salmon (Oncorhynchus tshawytscha), European seabass and gilthead seabream [19, 61,62,63]. This genus is often described as beneficial and it is found in healthy fish [19]. In our case, it was significantly more abundant in all control winter samplings (WC-t1 to WC-t3) and at the first winter sampling point (WR-t1), where R fish were still diagnosed as parasite-negative. Cutibacterium is known for producing propionate, a short-chain fatty acid, vitamins, and linolenic acid [61]. These microbial metabolites induce the strengthening of the epithelial barrier, reduce inflammation, and increase the production of mucus and antimicrobial peptides [61, 64]. Short-chain fatty acids have shown their benefits on gilthead seabream intestinal health [19] and in particular in their resilience against E. leei infections [9, 20]. The genus Staphylococcus has already been reported as dominant in the gut of gilthead seabream [4, 19, 65] and it is known for acting as a pathogenic or mutualistic genus, the latter being linked to host benefits such as pathogen inhibition and immune training [66]. Interestingly, Staphylococcus has been reported as a taxon increasing along host development [65]. Similarly, our results showed that this genus is mainly found in the winter trial, where fish were older than in the summer trial. Related to this, in this study we used fish from the same batch in an attempt to avoid microbiota variations due to genetic background, which have been described as highly relevant in this species [65, 67]. However, this implies that the age of the fish was different in the summer and the winter trials. Therefore, we cannot exclude the age effect on the detected variations among trials. A previous study conducted with gilthead seabream of different age groups (1, 2, and 4 years old) at an intermediate temperature of those used in this study (22 °C) detected Cutibacterium and Staphylococcus as markers for 1-year-old fish [19], coinciding with the results obtained from our fish group of similar age (winter sampling). A decrease in Proteobacteria, linked to an increase in Actinobacteria along the production cycle of gilthead seabream has been previously reported [65]. However, in this previous study, the peak of Actinobacteria abundance was detected in the summer trial. Future experiments should be conducted to specifically determine the effect of each variable.

Overall, fish from the summer trial showed an altered microbiota profile masking the effect of the parasite. E. leei proliferation can be suppressed at sustained high temperatures [45, 68] which may have affected the results of the final summer sampling, where the prevalence of infection is starting to decrease (from 100% in t2 to 66.7% in t3). As mentioned before, WR-t3 samples showed most changes in their intestinal microbial composition upon parasite infection, where Acinetobacter (Proteobacteria) was identified as the best marker for the parasite-exposed group. Acinetobacter spp. are a group of aerobic, non-fermentative, Gram-negative bacteria that have been distributed in diverse environments, including fish gut [69], and they are considered opportunistic pathogens that cause diseases in immunocompromised fish [70]. Interestingly, Piazzon et al. reported the presence of Acinetobacter in healthy 1-year-old fish [19], indicating that this bacteria is normally found in gilthead seabream gut as a potential opportunist. On the other hand, the genera Staphylococcus, Streptococcus, Rhodococcus, Pseudarthrobacter, Lawsonella and Corynebacterium were identified as markers for the C group in the final winter sampling. Briefly, Streptococcus, a genus of of the lactic acid bacteria, has been previously described as part of the intestinal microbiome of several fish species, and is considered a beneficial organism related to a healthy epithelium [71]. The genus Rhodococcus is present in the gills and intestines of some marine fish, showing bacteriocinogenic properties against fish pathogens [72]. In addition, it has been tested as a probiotic for some fish species [73]. Pseudarthrobacter has been described as a conditional fish pathogen [74]. However, in our results, its mean abundance in the control group was not very high, only reaching up to 1.4%. Finally, Lawsonella and Corynebacterium, have been reported in the gill and gut of healthy fish microbiomes, including gilthead seabream [4, 19].

On the host side, enteromyxosis induced larger significant changes in the intestinal transcriptomic profile than in spleen (935 DE genes vs. 14 DE genes, respectively). These results are in line with previous research in intestinal myxozoan infections [75, 76] reporting the highest effects on the intestine, the parasite’s main target tissue. In particular, out of the 422 annotated DE genes recognized in the pathway analysis, the intestinal response against E. leei was characterized by the activation of a local immune response, represented in four specific categories being the “interleukin related pathways” the one containing more DE genes. Previous results on gilthead seabream have also shown significant changes in the expression of interleukins (ILs) in the intestine (and not in the spleen) upon E. leei infection [76, 77]. In our previous studies, we have characterized that E. leei induced a type 1 immune response in the intestine [76] with up-regulation of interferon-stimulated genes such as ifi44, ifi35, and irf1 [78], also found in the current study. Type 1 immune responses are classically known to control intracellular infections activating cytotoxic functions [79]. However, type 1 immunity can be activated by large eukaryotic pathogens and is considered protective [80]. Interferon-related pathways (prototypical of a type 1 response) are linked to the signaling between innate and adaptive immune responses, involved in anti-parasitic response and resistance [81, 82]. The current results, with a significant up-regulation of interferon-related pathways, agree with the previous findings in this host-parasite model. On the contrary, transcriptome analysis of turbot (Scophthalmus maximus) intestines infected with Enteromyxum scophthalmi showed a marked down-regulation of genes related to interferon and MHCI pathways [75]. The different host response against these two Enteromyxum species has already been described [75] and seems to be related to the higher susceptibility of turbot and the chronic nature of infection in gilthead seabream during enteromyxosis. The protective role of type 1 responses has been described in other fish-parasite models. For example, irf1 and ifi35 up-regulation were related to resistance to Myxobolus cerebralis (Myzozoa) in rainbow trout (Oncorhynchus mykiss) and to amoebic gill disease in Atlantic salmon [83, 84].

Changes in the microbiota have already been described to influence host gene expression at local and systemic levels [85]. Therefore, we attempted to decipher whether the elicited immune responses were induced directly by the parasite or indirectly by changes in the microbiota. Correlation analyses showed that the differential expression of immune-related genes cannot be exclusively attributed to the parasite infection or changes in the microbiota, with most of the genes changing being correlated with both variables. Remarkably, the correlation network indicates that the detected activation of the immune response is linked to the decrease of several bacterial taxa, already described above to be common commensals and some with anti-inflammatory activities (i.e. Cutibacterium with the production of propionate). The decrease in these taxa could therefore explain part of the detected inflammatory response. This highlights again the relevance of including microbiota analyses when evaluating host responses to parasite infections, as it has been already described in other animal models including fish [86, 87]. Having pinpointed these marker genes and bacterial taxa, future functional studies should be directed towards elucidating the detailed molecular mechanisms in this interaction.

The current study also included a metatranscriptomic approach that allowed to identify the bacterial genes most differentially expressed during parasite infection. Within these bacterial-specific differentially expressed transcripts, several relevant genes have been identified most having significant correlations with host immune gene expression. For example, two flagellin genes have been linked to the expression of ifnγ. Flagellin, a major structural protein of the flagellum of motile bacteria, is a potent stimulator of inflammatory responses including interferon production and signaling [88], and is being used as vaccine adjuvant in fish aquaculture [89]. CysM, also linked to the expression of ifnγ, is a crucial gene for the synthesis of cysteine, key in bacterial defense mechanisms against host defenses. This gene has been tagged as a target for the development of antimicrobials against pathogens [90]. RihB expression has been positively linked with the expression of a large number of immune genes. It is a ribonucleoside hydrolase, but its exact function remains elusive [91]. Nonetheless, the current results show that RihB might be implicated in the regulation or activation of inflammatory responses in fish intestines. UvrY is part of the bacterial two-component system and has been shown to have a role in gut colonization of pathogens [92]. Interestingly, it can be expressed by Acinetobacter species and has been labeled as an interesting target for novel antibiotic development [93].

Conclusions

Enteromyxum leei infections induced a clear shift in the microbial populations of the intestine. However, the microbiota is highly susceptible to many different variables and, in the current study, we demonstrated that high temperatures and age can have a very strong effect capable of masking the impact of the parasite. The multi-omic approach undertaken allowed us to establish that the host response, classically attributed to the parasite infection, can be induced by the microbial shifts caused by the infection, not only shifts in the taxonomies found in these populations but also changes in their gene expression. Still, knowledge about which component of this complex system is responsible for which specific effect remains elusive. However, the present work presents a series of datasets, bacterial markers and target genes to perform future functional studies to understand this three-way relationship (host-parasite-microbiota) in more detail.

Data availability

All relevant data are within the paper and its additional files. All the sequencing data presented in this study is available at the Sequence Read Archive (SRA), under Bioproject accession number PRJNA1134388, BioSample accession numbers SAMN42418970-SAMN42419029 for the 16 S rRNA amplicon sequencing, SAMN42419030-SAMN42419053 for the host RNA sequencing, and SAMN42419054-SAMN42419065 for the metatranscriptome sequencing.

References

  1. FAO. The state of world fisheries and aquaculture 2024 – Blue transformation in action. Rome: FAO; 2024.

    Google Scholar 

  2. Sicuro B. The evolution of aquaculture in the mediterranean region: an anthropogenic climax stage? PLoS ONE. 2024;19:e0290870.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  3. Mhalhel K, Levanti M, Abbate F, Laurà R, Guerrera MC, Aragona M, et al. Review on Gilthead seabream (Sparus aurata) aquaculture: life cycle, growth, aquaculture practices and challenges. J Mar Sci Eng. 2023;11:11102008.

    Article  Google Scholar 

  4. Quero GM, Piredda R, Basili M, Maricchiolo G, Mirto S, Manini E, et al. Host-associated and environmental microbiomes in an open-sea mediterranean Gilthead sea Bream fish farm. Microb Ecol. 2023;86:1319–30.

    Article  CAS  PubMed  Google Scholar 

  5. Estensoro I, Redondo MJ, Álvarez-Pellitero P, Sitjà-Bobadilla A. Novel horizontal transmission route for Enteromyxum Leei (Myxozoa) by anal intubation of Gilthead sea Bream Sparus aurata. Dis Aquat Organ. 2010;92:51–8.

    Article  PubMed  Google Scholar 

  6. Sitjà-Bobadilla A, Palenzuela O. Enteromyxum species. In: Woo P, Buchmann K, editors. Fish parasites pathobiol prot. CABI; 2012. pp. 163–76.

  7. Picard-Sánchez A, Estensoro I, Del Pozo R, Palenzuela OR, Piazzon MC, Sitjà-Bobadilla A. Water temperature, time of exposure and population density are key parameters in Enteromyxum Leei fish-to-fish experimental transmission. J Fish Dis. 2020;43:491–502.

    Article  PubMed  Google Scholar 

  8. Sitjà-Bobadilla A, Diamant A, Palenzuela O, Álvarez-Pellitero P. Effect of host factors and experimental conditions on the horizontal transmission of Enteromyxum Leei (Myxozoa) to Gilthead sea Bream, Sparus aurata L., and European sea bass, Dicentrarchus labrax (L). J Fish Dis. 2007;30:243–50.

    Article  PubMed  Google Scholar 

  9. Palenzuela O, Del Pozo R, Piazzon MC, Isern-Subich MM, Ceulemans S, Coutteau P, et al. Effect of a functional feed additive on mitigation of experimentally induced Gilthead sea Bream Sparus aurata enteromyxosis. Dis Aquat Organ. 2020;138:111–20.

    Article  CAS  PubMed  Google Scholar 

  10. Alama-Bermejo G, Šíma R, Raga JA, Holzer AS. Understanding Myxozoan infection dynamics in the Sea: seasonality and transmission of Ceratomyxa Puntazzi. Int J Parasitol. 2013;43:771–80.

    Article  PubMed  Google Scholar 

  11. Ishimaru K, Matsuura T, Tsunemoto K, Shirakashi S. Seasonal monitoring of Kudoa Yasunagai from sea water and aquaculture water using quantitative PCR. Dis Aquat Organ. 2014;108:45–52.

    Article  CAS  PubMed  Google Scholar 

  12. Okamura B, Hartikainen H, Schimidt-Posthaus H, Wahli T. Life cycle complexity, environmental change and the emerging status of salmonid proliferative kidney disease. Freshw Biol. 2011;56:735–53.

    Article  Google Scholar 

  13. Jones SRM, Bartholomew JL, Zhang JY. Mitigating Myxozoan disease impacts on wild fish populations. In: Okamura B, Gruhl A, Bartholomew JL, editors. Myxozoan evol ecol Dev. Cham: Springer International Publishing; 2015. pp. 397–413.

    Chapter  Google Scholar 

  14. Diwan AD, Harke SN, Panche AN. Host-microbiome interaction in fish and shellfish: an overview. Fish Shellfish Immunol Rep. 2023;4:100091.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  15. Medina-Félix D, Garibay-Valdez E, Vargas-Albores F, Martínez-Porchas M. Fish disease and intestinal microbiota: A close and indivisible relationship. Rev Aquac. 2023;15:820–39.

    Article  Google Scholar 

  16. Legrand TPRA, Wynne JW, Weyrich LS, Oxley APA. A microbial sea of possibilities: current knowledge and prospects for an improved Understanding of the fish Microbiome. Rev Aquac. 2020;12:1101–34.

    Article  Google Scholar 

  17. Xavier R, Severino R, Silva SM. Signatures of dysbiosis in fish microbiomes in the context of aquaculture. Rev Aquac. 2024;16:706–31.

    Article  Google Scholar 

  18. Palenzuela O, Bartholomew JL. Molecular tools for the diagnosis of Ceratomyxa Shasta (Myxozoa). In: Cunningham CO, editor. Mol diagnosis salmon dis. Dordrecht: Springer Netherlands; 2002. pp. 285–98.

    Chapter  Google Scholar 

  19. Piazzon MC, Naya-Català F, Simó-Mirabet P, Picard-Sánchez A, Roig FJ, Calduch-Giner JA, et al. Sex, age, and bacteria: how the intestinal microbiota is modulated in a protandrous hermaphrodite fish. Front Microbiol. 2019;10:2512.

    Article  PubMed  PubMed Central  Google Scholar 

  20. Piazzon MC, Calduch-Giner JA, Fouz B, Estensoro I, Simó-Mirabet P, Puyalto M, et al. Under control: how a dietary additive can restore the gut Microbiome and proteomic profile, and improve disease resilience in a marine teleostean fish fed vegetable diets. Microbiome. 2017;5:164.

    Article  PubMed  PubMed Central  Google Scholar 

  21. Schmieder R, Edwards R. Quality control and preprocessing of metagenomic datasets. Bioinformatics. 2011;27:863–4.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  22. Aronesty E. Comparison of sequencing utility programs. Open Bioinforma J. 2013;7:1–8.

    Article  Google Scholar 

  23. Li H. New strategies to improve minimap2 alignment accuracy. Bioinformatics. 2021;37:4572–4.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  24. Yilmaz P, Parfrey LW, Yarza P, Gerken J, Pruesse E, Quast C, et al. The SILVA and All-species living tree project (LTP) taxonomic frameworks. Nucleic Acids Res. 2014;42:D643–8.

    Article  CAS  PubMed  Google Scholar 

  25. Pérez-Sánchez J, Naya-Català F, Soriano B, Piazzon MC, Hafez H, Gabaldón T, et al. Genome sequencing and transcriptome analysis reveal recent species-specific gene duplications in the plastic Gilthead sea Bream (Sparus aurata). Front Mar Sci. 2019;6:760.

    Article  Google Scholar 

  26. Kim D, Paggi JM, Park C, Bennett C, Salzberg SL. Graph-based genome alignment and genotyping with HISAT2 and HISAT-genotype. Nat Biotechnol. 2019;37:907–15.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  27. Liao Y, Smyth GK, Shi W. FeatureCounts: an efficient general purpose program for assigning sequence reads to genomic features. Bioinformatics. 2014;30:923–30.

    Article  CAS  PubMed  Google Scholar 

  28. Bolger AM, Lohse M, Usadel B. Trimmomatic: A flexible trimmer for illumina sequence data. Bioinformatics. 2014;30:2114–20.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  29. Kopylova E, Noé L, Touzet H. SortMeRNA: fast and accurate filtering of ribosomal RNAs in metatranscriptomic data. Bioinformatics. 2012;28:3211–7.

    Article  CAS  PubMed  Google Scholar 

  30. Dobin A, Davis CA, Schlesinger F, Drenkow J, Zaleski C, Jha S, et al. STAR: ultrafast universal RNA-seq aligner. Bioinformatics. 2013;29:15–21.

    Article  CAS  PubMed  Google Scholar 

  31. Danecek P, Bonfield JK, Liddle J, Marshall J, Ohan V, Pollard MO et al. Twelve years of samtools and BCFtools. Gigascience. 2021;10.

  32. Li D, Liu C-M, Luo R, Sadakane K, Lam T-W. MEGAHIT: an ultra-fast single-node solution for large and complex metagenomics assembly via succinct de Bruijn graph. Bioinformatics. 2015;31:1674–6.

    Article  CAS  PubMed  Google Scholar 

  33. Beghini F, McIver LJ, Blanco-Míguez A, Dubois L, Asnicar F, Maharjan S, et al. Integrating taxonomic, functional, and strain-level profiling of diverse microbial communities with biobakery 3. Elife. 2021;10:e65088.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  34. Love MI, Huber W, Anders S. Moderated Estimation of fold change and dispersion for RNA-seq data with DESeq2. Genome Biol. 2014;15:550.

    Article  PubMed  PubMed Central  Google Scholar 

  35. McMurdie PJ, Holmes S. Phyloseq: an R package for reproducible interactive analysis and graphics of Microbiome census data. PLoS ONE. 2013;8:e61217.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  36. Oksanen J, Simpson G, Blanchet F, Kindt R, Legendre P, Minchin P et al. vegan: Community ecology package [Internet]. CRAN.R-project; 2022. Available from: https://cran.r-project.org/package=vegan

  37. Rohart F, Gautier B, Singh A, Lê Cao K-A, mixOmics. An R package for ’omics feature selection and multiple data integration. PLoS Comput Biol. 2017;13:e1005752.

  38. Segata N, Izard J, Waldron L, Gevers D, Miropolsky L, Garrett WS, et al. Metagenomic biomarker discovery and explanation. Genome Biol. 2011;12:R60.

    Article  PubMed  PubMed Central  Google Scholar 

  39. Gillespie M, Jassal B, Stephan R, Milacic M, Rothfels K, Senff-Ribeiro A, et al. The reactome pathway knowledgebase 2022. Nucleic Acids Res. 2022;50:D687–92.

    Article  CAS  PubMed  Google Scholar 

  40. Väremo L, Nielsen J, Nookaew I. Enriching the gene set analysis of genome-wide data by incorporating directionality of gene expression and combining statistical hypotheses and methods. Nucleic Acids Res. 2013;41:4378–91.

    Article  PubMed  PubMed Central  Google Scholar 

  41. Shannon P, Markiel A, Ozier O, Baliga NS, Wang JT, Ramage D, et al. Cytoscape: A software environment for integrated models of biomolecular interaction networks. Genome Res. 2003;13:2498–504.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  42. Wei T, Simko V. R package corrplot: Visualization of a correlation matrix [Internet]. 2021. Available from: https://github.com/taiyun/corrplot

  43. Toxqui-Rodríguez S, Riera-Ferrer E, Del Pozo R, Palenzuela O, Sitjà-Bobadilla A, Estensoro I, et al. Molecular interactions in an holobiont-pathogen model: integromics in Gilthead seabream infected with Sparicotyle chrysophrii. Aquaculture. 2024;581:740365.

    Article  Google Scholar 

  44. Ruiz A, Andree KB, Furones D, Holhorea PG, Calduch-Giner JÀ, Viñas M, et al. Modulation of gut microbiota and intestinal immune response in Gilthead seabream (Sparus aurata) by dietary bile salt supplementation. Front Microbiol. 2023;14:1123716.

    Article  PubMed  PubMed Central  Google Scholar 

  45. Piazzon MC, Naya-Català F, Pereira GV, Estensoro I, Del Pozo R, Calduch-Giner JA, et al. A novel fish meal-free diet formulation supports proper growth in Gilthead sea Bream (Sparus aurata) with a reshape of tissue-specific gene expression patterns and gut microbiota, and without worsening intestinal parasite susceptibility. Aquaculture. 2022;558:738362.

    Article  CAS  Google Scholar 

  46. Firmino JP, Vallejos-Vidal E, Balebona MC, Ramayo-Caldas Y, Cerezo IM, Salomón R, et al. Diet, immunity, and microbiota interactions: an integrative analysis of the intestine transcriptional response and microbiota modulation in Gilthead seabream (Sparus aurata) fed an essential oils-based functional diet. Front Immunol. 2021;12:625297.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  47. Fu PP, Xiong F, Feng WW, Zou H, Wu SG, Li M, et al. Effect of intestinal tapeworms on the gut microbiota of the common carp, Cyprinus carpio. Parasit Vectors. 2019;12:252.

    Article  PubMed  PubMed Central  Google Scholar 

  48. Vasemägi A, Visse M, Kisand V. Effect of environmental factors and an emerging parasitic disease on gut Microbiome of wild salmonid fish. mSphere. 2017;2:e00418–17.

    Article  PubMed  PubMed Central  Google Scholar 

  49. Minich JJ, Power C, Melanson M, Knight R, Webber C, Rough K, et al. The Southern Bluefin tuna mucosal Microbiome is influenced by husbandry method, net pen location, and anti-parasite treatment. Front Microbiol. 2020;11:2015.

    Article  PubMed  PubMed Central  Google Scholar 

  50. Tran NT, Zhang J, Xiong F, Wang G-T, Li W-X, Wu S-G. Altered gut microbiota associated with intestinal disease in grass carp (Ctenopharyngodon idellus). World J Microbiol Biotechnol. 2018;34:71.

    Article  PubMed  Google Scholar 

  51. Li T, Long M, Ji C, Shen Z, Gatesoupe F-J, Zhang X, et al. Alterations of the gut Microbiome of largemouth bronze gudgeon (Coreius guichenoti) suffering from furunculosis. Sci Rep. 2016;6:30606.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  52. Hovda MB, Fontanillas R, McGurk C, Obach A, Rosnes JT. Seasonal variations in the intestinal microbiota of farmed Atlantic salmon (Salmo Salar L). Aquac Res. 2012;43:154–9.

    Article  Google Scholar 

  53. Egerton S, Culloty S, Whooley J, Stanton C, Ross RP. The gut microbiota of marine fish. Front Microbiol. 2018;9:873.

    Article  PubMed  PubMed Central  Google Scholar 

  54. Ringø E, Zhou Z, Vecino JLG, Wadsworth S, Romero J, Krogdahl Å, et al. Effect of dietary components on the gut microbiota of aquatic animals. A never-ending story? Aquac Nutr. 2016;22:219–82.

    Article  Google Scholar 

  55. Neuman C, Hatje E, Zarkasi KZ, Smullen R, Bowman JP, Katouli M. The effect of diet and environmental temperature on the faecal microbiota of farmed Tasmanian Atlantic salmon (Salmo Salar L). Aquac Res. 2016;47:660–72.

    Article  CAS  Google Scholar 

  56. Hansen GH, Olafsen JA. Bacterial interactions in early life stages of marine cold water fish. Microb Ecol. 1999;38:1–26.

    Article  CAS  PubMed  Google Scholar 

  57. Zhang Y, Zhao R, Shi D, Sun S, Ren H, Zhao H, et al. Characterization of the Circulating Microbiome in acute-on-chronic liver failure associated with hepatitis B. Liver Int. 2019;39:1207–16.

    Article  PubMed  Google Scholar 

  58. Coenye T. The family Burkholderiaceae. In: Rosenberg E, DeLong EF, Lory S, Stackebrandt E, Thompson F, editors. The prokaryotes. Berlin, Heidelberg: Springer Berlin Heidelberg; 2014. pp. 759–76.

    Chapter  Google Scholar 

  59. Miranda RM, Aguila-Torres P, Aranda CP, Maldonado J, Casado A. Taxonomy and diversity of bacterial communities associated with marine sediments from Chilean salmonid farms. Aquac Res. 2021;52:1605–20.

    Article  CAS  Google Scholar 

  60. Dréno B, Pécastaings S, Corvec S, Veraldi S, Khammari A, Roques C. Cutibacterium acnes (Propionibacterium acnes) and acne vulgaris: a brief look at the latest updates. J Eur Acad Dermatology Venereol. 2018;32:5–14.

    Article  Google Scholar 

  61. Ferreira M, Abdelhafiz Y, Abreu H, Silva J, Valente LMP, Kiron V. Gracilaria gracilis and Nannochloropsis oceanica, singly or in combination, in diets alter the intestinal microbiota of European Seabass (Dicentrarchus labrax). Front Mar Sci. 2022;9:1001942.

    Article  Google Scholar 

  62. Zhao R, Symonds JE, Walker SP, Steiner K, Carter CG, Bowman JP, et al. Salinity and fish age affect the gut microbiota of farmed Chinook salmon (Oncorhynchus tshawytscha). Aquaculture. 2020;528:735539.

    Article  CAS  Google Scholar 

  63. Dehler CE, Secombes CJ, Martin SAM. Seawater transfer alters the intestinal microbiota profiles of Atlantic salmon (Salmo Salar L). Sci Rep. 2017;7:13877.

    Article  PubMed  PubMed Central  Google Scholar 

  64. Pérez Chaia A, Zarate G, Oliver G. The probiotic properties of propionibacteria. Lait. 1999;79:175–85.

    Article  Google Scholar 

  65. Naya-Català F, Piazzon MC, Torrecillas S, Toxqui-Rodríguez S, Calduch-Giner JÀ, Fontanillas R, et al. Genetics and nutrition drive the gut microbiota succession and host-transcriptome interactions through the Gilthead sea Bream (Sparus aurata) production cycle. Biology. 2022;11:1744.

    Article  PubMed  PubMed Central  Google Scholar 

  66. Cogen AL, Nizet V, Gallo RL. Skin microbiota: a source of disease or defence? Br J Dermatol. 2008;158:442–55.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  67. Piazzon MC, Naya-Català F, Perera E, Palenzuela O, Sitjà-Bobadilla A, Pérez-Sánchez J. Genetic selection for growth drives differences in intestinal microbiota composition and parasite disease resistance in Gilthead sea Bream. Microbiome. 2020;8:168.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  68. China M, Nakamura H, Hamakawa K, Tamaki E, Yokoyama H, Masuoka S, et al. Efficacy of high water temperature treatment of myxosporean emaciation disease caused by Enteromyxum Leei (Myxozoa). Fish Pathol. 2014;49:137–40.

    Article  Google Scholar 

  69. Sun YZ, Yang HL, Ling ZC, Chang JB, Ye JD. Gut microbiota of fast and slow growing grouper Epinephelus coioides. Afr J Microbiol Res. 2009;3:713–20.

    Google Scholar 

  70. Bunnoy A, Na-Nakorn U, Srisapoome P. Probiotic effects of a novel strain, Acinetobacter KU011TH, on the growth performance, immune responses, and resistance against Aeromonas hydrophila of Bighead catfish (Clarias macrocephalus Günther, 1864). Microorganisms. 2019;7:613.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  71. Estruch G, Collado MC, Penaranda DS, Tomas Vidal A, Jover Cerda M, Perez Martinez G, et al. Impact of fishmeal replacement in diets for Gilthead sea Bream (Sparus aurata) on the Gastrointestinal microbiota determined by pyrosequencing the 16S rRNA gene. PLoS ONE. 2015;10:e0136389.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  72. Ghotbi M, Kelting O, Blümel M, Tasdemir D. Gut and gill-associated microbiota of the flatfish European plaice (Pleuronectes platessa): Diversity, metabolome and bioactivity against human and aquaculture pathogens. Mar Drugs. 2022;20:573.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  73. Sharifuzzaman SM, Rahman H, Austin DA, Austin B. Properties of probiotics Kocuria SM1 and Rhodococcus SM2 isolated from fish guts. Probiotics Antimicrob Proteins. 2018;10:534–42.

    Article  CAS  PubMed  Google Scholar 

  74. Charo FJ, Mbuthia PG, Bebora LC, Nguta JM. Influence of aquaculture management practices and water quality on bacterial occurrence in fish culture units in Kenya. Int J Fish Aquat Stud. 2023;11:01–7.

    Article  Google Scholar 

  75. Robledo D, Ronza P, Harrison PW, Losada AP, Bermúdez R, Pardo BG, et al. RNA-seq analysis reveals significant transcriptome changes in turbot (Scophthalmus maximus) suffering severe enteromyxosis. BMC Genomics. 2014;15:1149.

    Article  PubMed  PubMed Central  Google Scholar 

  76. Piazzon MC, Estensoro I, Calduch-Giner JA, Del Pozo R, Picard-Sánchez A, Pérez-Sánchez J, et al. Hints on T cell responses in a fish-parasite model: Enteromyxum Leei induces differential expression of T cell signature molecules depending on the organ and the infection status. Parasit Vectors. 2018;11:443.

    Article  PubMed  PubMed Central  Google Scholar 

  77. Pérez-Cordón G, Estensoro I, Benedito-Palos L, Calduch-Giner JA, Sitjà-Bobadilla A, Pérez-Sánchez J. Interleukin gene expression is strongly modulated at the local level in a fish-parasite model. Fish Shellfish Immunol. 2014;37:201–8.

    Article  PubMed  Google Scholar 

  78. Davey GC, Calduch-Giner JA, Houeix B, Talbot A, Sitjà-Bobadilla A, Prunet P, et al. Molecular profiling of the Gilthead sea Bream (Sparus aurata L.) response to chronic exposure to the myxosporean parasite Enteromyxum Leei. Mol Immunol. 2011;48:2102–12.

    Article  CAS  PubMed  Google Scholar 

  79. Yamaguchi T, Takizawa F, Fischer U, Dijkstra JM. Along the axis between type 1 and type 2 immunity; principles conserved in evolution from fish to mammals. Biology. 2015;4:814–59.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  80. Spellberg B, Edwards JE Jr. Type 1/type 2 immunity in infectious diseases. Clin Infect Dis. 2001;32:76–102.

    Article  CAS  PubMed  Google Scholar 

  81. Le Bon A, Tough DF. Links between innate and adaptive immunity via type I interferon. Curr Opin Immunol. 2002;14:432–6.

    Article  PubMed  Google Scholar 

  82. Zou J, Secombes CJ. Teleost fish interferons and their role in immunity. Dev Comp Immunol. 2011;35:1376–87.

    Article  CAS  PubMed  Google Scholar 

  83. Baerwald MR, Welsh AB, Hedrick RP, May B. Discovery of genes implicated in whirling disease infection and resistance in rainbow trout using genome-wide expression profiling. BMC Genomics. 2008;9:37.

    Article  PubMed  PubMed Central  Google Scholar 

  84. Wynne JW, O’Sullivan MG, Stone G, Cook MT, Nowak BF, Lovell DR, et al. Resistance to amoebic gill disease (AGD) is characterised by the transcriptional dysregulation of immune and cell cycle pathways. Dev Comp Immunol. 2008;32:1539–60.

    Article  CAS  PubMed  Google Scholar 

  85. Naya-Català G, Piazzon MC, Fernandes AM, Calduch-Giner JA, Sitjà-Bobadilla A et al. The cross-talk between intestinal microbiota and host gene expression in juveniles of gilthead sea bream (Sparus aurata). Insights in fish feeds for increased circularity and resource utilization. Front Physiol. 2021;12:784265.

  86. Brealey JC, Kodama M, Rasmussen JA, Hansen SB, Santos-Bay L, Lecaudey LA, et al. Host-gut microbiota interactions shape parasite infections in farmed Atlantic salmon. mSystems. 2024;9:e0104323.

    Article  PubMed  Google Scholar 

  87. Niciura SCM, Cardoso TF, Ibelli AMG, Okino CH, Andrade BG, Benavides MV, et al. Multi-omics data elucidate parasite-host-microbiota interactions and resistance to Haemonchus contortus in sheep. Parasit Vectors. 2024;17:102.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  88. Kang W, Park A, Huh J-W, You G, Jung D-J, Song M, et al. Flagellin-stimulated production of interferon-β promotes anti-flagellin IgG2c and IgA responses. Mol Cells. 2020;43:251–63.

    CAS  PubMed  PubMed Central  Google Scholar 

  89. Wangkahart E, Secombes CJ, Wang T. Studies on the use of Flagellin as an immunostimulant and vaccine adjuvant in fish aquaculture. Front Immunol. 2018;9:3054.

    Article  CAS  PubMed  Google Scholar 

  90. Hicks JL, Oldham KEA, McGarvie J, Walker EJ. Combatting antimicrobial resistance via the cysteine biosynthesis pathway in bacterial pathogens. Biosci Rep. 2022;42:BSR20220368.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  91. Petersen C, Møller LB. The RihA, RihB, and RihC ribonucleoside hydrolases of Escherichia coli. Substrate specificity, gene expression, and regulation. J Biol Chem. 2001;276:884–94.

    Article  CAS  PubMed  Google Scholar 

  92. Miki T, Hoshino Y, Sudo N, Ito M, Haneda T, Okada N. UvrY deletion and acetate reduce gut colonization of Crohn’s disease-associated adherent-invasive Escherichia coli by decreasing expression of type 1 fimbriae. Infect Immun. 2022;90:e0066221.

    Article  PubMed  Google Scholar 

  93. Ahmad S, Raza S, Uddin R, Azam SS. Comparative subtractive proteomics based ranking for antibiotic targets against the dirtiest Superbug: Acinetobacter baumannii. J Mol Graph Model. 2018;82:74–92.

    Article  CAS  PubMed  Google Scholar 

Download references

Acknowledgements

The authors thank the staff of IATS, CSIC fish facilities for technical assistance.

Funding

This work has been financed by the project MucosalFrontier (PID2020-115070RA-I00) funded by MCIN/AEI/https://doiorg.publicaciones.saludcastillayleon.es/10.13039/501100011033. Additional funding was obtained from the Generatitat Valenciana AICO2023 funding (CIAICO/2022/144), and the European Union’s Horizon 2020 MSCA-ITN-2020 under grant agreement no. 956697; EATFISH. This publication reflects only the authors’ view and the European Union cannot be held responsible for any use that may be made of the information contained herein. This study is part of the ThinkInAzul programme and was supported by MICIU with funding from the European Union Next Generation EU (PRTR-C17.I1) and Generalitat Valenciana (GVA-THINKINAZUL/2021/022; Principal investigator: A. Sitjá-Bobadilla, IATS, CSIC).

Author information

Authors and Affiliations

  1. Fish Pathology Group, Instituto de Acuicultura Torre de la Sal (IATS, CSIC), Castellón, Spain

    Socorro Toxqui-Rodríguez, Itziar Estensoro, Raquel Del Pozo, Ariadna Sitjà-Bobadilla & M. Carla Piazzon

  2. Nutrigenomics and Fish Endocrinology Group, Instituto de Acuicultura Torre de la Sal (IATS, CSIC), Castellón, Spain

    Socorro Toxqui-Rodríguez, Ricardo Domingo-Bretón & Jaume Pérez-Sánchez

  3. Laboratory of Microbiology, Wageningen University and Research, Wageningen, The Netherlands

    Detmer Sipkema

Authors
  1. Socorro Toxqui-Rodríguez
    View author publications

    You can also search for this author inPubMed Google Scholar

  2. Itziar Estensoro
    View author publications

    You can also search for this author inPubMed Google Scholar

  3. Ricardo Domingo-Bretón
    View author publications

    You can also search for this author inPubMed Google Scholar

  4. Raquel Del Pozo
    View author publications

    You can also search for this author inPubMed Google Scholar

  5. Jaume Pérez-Sánchez
    View author publications

    You can also search for this author inPubMed Google Scholar

  6. Detmer Sipkema
    View author publications

    You can also search for this author inPubMed Google Scholar

  7. Ariadna Sitjà-Bobadilla
    View author publications

    You can also search for this author inPubMed Google Scholar

  8. M. Carla Piazzon
    View author publications

    You can also search for this author inPubMed Google Scholar

Contributions

Conceptualization: MCP, ASB, IE. Data curation and formal analysis: STR, MCP, IE, RDB. Funding acquisition: MCP, ASB, JPS, DS. Investigation: STR, IE, RDB, RDP, MCP. Project administration: MCP, ASB, JPS, DS. Resources: MCP, ASB, JPS, DS. Supervision: DS, ASB, MCP, JPS. Visualization: MCP, STR. Writing - original draft: STR, MCP. Writing - review & editing: all.

Corresponding author

Correspondence to M. Carla Piazzon.

Ethics declarations

Ethics approval and consent to participate

This study does not involve the use of human subjects. The ethic statements regarding animal experimentation are included in the Methods section.

Consent for publication

Not applicable.

Competing interests

The authors declare no competing interests.

Additional information

Publisher’s note

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

Electronic supplementary material

Rights and permissions

Open Access This article is licensed under a Creative Commons Attribution-NonCommercial-NoDerivatives 4.0 International License, which permits any non-commercial use, sharing, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if you modified the licensed material. You do not have permission under this licence to share adapted material derived from this article or parts of it. The images or other third party material in this article are included in the article’s Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article’s Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by-nc-nd/4.0/.

Reprints and permissions

About this article

Check for updates. Verify currency and authenticity via CrossMark

Cite this article

Toxqui-Rodríguez, S., Estensoro, I., Domingo-Bretón, R. et al. Interactions between gilthead seabream intestinal transcriptome and microbiota upon Enteromyxum leei infection: a multi–omic approach. anim microbiome 7, 22 (2025). https://doiorg.publicaciones.saludcastillayleon.es/10.1186/s42523-025-00388-x

Download citation

  • Received:

  • Accepted:

  • Published:

  • DOI: https://doiorg.publicaciones.saludcastillayleon.es/10.1186/s42523-025-00388-x

Keywords

Animal Microbiome

ISSN: 2524-4671