Skip to main content

Oyster hemolymph is a complex and dynamic ecosystem hosting bacteria, protists and viruses



The impact of the microbiota on host fitness has so far mainly been demonstrated for the bacterial microbiome. We know much less about host-associated protist and viral communities, largely due to technical issues. However, all microorganisms within a microbiome potentially interact with each other as well as with the host and the environment, therefore likely affecting the host health.


We set out to explore how environmental and host factors shape the composition and diversity of bacterial, protist and viral microbial communities in the Pacific oyster hemolymph, both in health and disease. To do so, five oyster families differing in susceptibility to the Pacific oyster mortality syndrome were reared in hatchery and transplanted into a natural environment either before or during a disease outbreak. Using metabarcoding and shotgun metagenomics, we demonstrate that hemolymph can be considered as an ecological niche hosting bacterial, protist and viral communities, each of them shaped by different factors and distinct from the corresponding communities in the surrounding seawater. Overall, we found that hemolymph microbiota is more strongly shaped by the environment than by host genetic background. Co-occurrence network analyses suggest a disruption of the microbial network after transplantation into natural environment during both non-infectious and infectious periods. Whereas we could not identify a common microbial community signature for healthy animals, OsHV-1 μVar virus dominated the hemolymph virome during the disease outbreak, without significant modifications of other microbiota components.


Our study shows that oyster hemolymph is a complex ecosystem containing diverse bacteria, protists and viruses, whose composition and dynamics are primarily determined by the environment. However, all of these are also shaped by oyster genetic backgrounds, indicating they indeed interact with the oyster host and are therefore not only of transient character. Although it seems that the three microbiome components respond independently to environmental conditions, better characterization of hemolymph-associated viruses could change this picture.


Since the proposition of the hologenome theory of evolution [1], the majority of microbiome-related research has focused on its prokaryotic component (for review see [2]). However, recent findings indicate that we should also consider protists and viruses to understand the holobiont ecosystem (the host with its associated microbiota) and the role of microbiota in metaorganism functioning. For instance, the presence of protozoa correlates with higher diversity of the bacterial gut microbiome in humans, suggesting their potential beneficial role in the maintenance of homeostasis [3]. Similarly, the human body (including the circulatory system [4, 5]) is inhabited by highly diverse viral communities that can directly or indirectly influence the balance between health and disease [6], notably through their ability to modulate microbiota by infecting both bacteria and eukaryotes [7, 8]. On the other hand, the bacterial microbiota can shape the outcome of viral infections in mammals [9]. Nevertheless, the role and dynamics of host-associated non-bacterial microbiotas in invertebrates is understudied, with a few exceptions like corals [10, 11], hydra [12] or some insects [13,14,15].

One of the services that microbiota may provide to its host is the protection against pathogens. Numerous examples of microbiota-mediated protection span the tree of life and have recently inspired a concept of co-immunity, reflecting the increasing evidence that a host is protected not only by its own immune system, but also by its microbiota [16, 17]. On the other hand, disease is often associated with dysbiosis [18, 19] and a growing body of evidence supports a shift from a primarily “a pathogen - a disease” paradigm towards a more ecological view of disease (“a dysbiosis - a disease”), as reflected in the recently proposed “ecological Koch’s postulates” [20]. Overall, it is becoming clear that the understanding of a powerful selection pressure such as disease requires an understanding of the holobiont system.

In this respect, Pacific oysters (Crassostrea gigas) represent an interesting model for studying host-microbiome interactions, as for or more than a decade, oyster production has been plagued with a recurring mortality syndrome, called Pacific oyster mortality syndrome (POMS) [21]. POMS affects all coastal regions of France and numerous other producing countries [22, 23]. It coincides with the recurrent detection of Ostreid herpesvirus variants (OsHV-1 μVar) in moribund oysters both in France [24,25,26] and worldwide [22, 27,28,29,30,31]. In addition OsHV-1 μVar, POMS has been associated with various strains of the bacterial genus Vibrio [32]. Among these, populations of Vibrio crassostreae have been repeatedly identified in diseased oysters [33, 34]. We recently demonstrated that the disease is caused by polymicrobial infection, with the invasion of oyster hemocytes by OsHV-1 μVar as the initial and necessary step [35]. Altogether, these data highlight the need to simultaneously investigate the role of all hemolymph-associated microorganisms (bacteria, protists and viruses) in Pacific oyster health and disease.

So far, all studies regarding bivalve microbiota have focused on bacteria (reviewed in [36, 37]). The bacterial communities are tissue-specific: the microbiotas in oyster solid tissue differ from the hemolymph microbiota and both are distinct from the surrounding bacterial plankton communities [38]. The presence of living bacteria in the hemolymph of healthy bivalves [39,40,41] has been considered paradoxical, as the hemolymph contains hemocytes - circulating immune cells with phagocytic activity - and thus plays a key role in the oyster defense mechanisms [42]. However, it has been suggested that some of the resident hemolymph bacteria may contribute to oyster protection by producing antimicrobial peptides [43]. Although the number of descriptive studies of oyster bacterial microbiotas has been growing, few have addressed the factors shaping these communities [37]. Temperature seems to be a key driver of bacterial community structure and dynamics in oyster hemolymph [44], but many other host-intrinsic and environmental factors are potentially involved. Transplantation experiments conducted across different natural environments indeed suggested that the bacterial community structure and dynamics are influenced by a complex interplay of host-related factors, biotic interactions within the microbiota and local environmental conditions that require further exploration [45]. Finally, gill-associated bacterial communities are influenced by oyster genetics [46], but we know nothing about the influence of genetic factors on the hemolymph microbiota.

In this study, we aimed to explore how the environment, oyster genetic background and the hemolymph microbiota (bacteria, protists and viruses) interact with each other in the natural environment before and during a disease outbreak. Specifically, we aimed to: (i) describe the diversity and composition of oyster bacterial, protist and viral hemolymph microbiome, (ii) estimate the relative weight of host genetic background and environmental factors affecting the oyster blood microbiota composition in controlled as well as natural conditions (with low and high POMS risk), and (iii) detect co-occurrence patterns of microbes in oyster hemolymph and to examine their potential associations with the host health. By raising genetically distinct oyster families in controlled hatchery conditions, we were able to estimate the effect of host genetic backgrounds on the hemolymph microbiome, whereas transplantation experiments allowed us to examine the interactions between the environment, host genetic backgrounds and hemolymph microbiome in both health and disease.


Production of biparental oyster families

Five biparental oyster families were produced by in vitro fertilization from wild genitors sampled in farming and non-farming areas as previously described [35]. Briefly, wild genitors used to produce Atlantic families F9 and F15 were collected at two sites approximately 20 km apart, in Logonna Daoulas (lat.: 48.335263 long.: − 4.317922, farming areas) and Dellec (lat.: 48.353970 long.: − 4.566123, non-farming areas), respectively. Wild genitors used to produce Mediterranean families F32 and F44 were collected at two sites approximately 40 km apart, in Vidourle (lat.: 43.553906 long.: 4.095175, non-farming areas) and Thau lagoon (lat.: 43.418736 long.: 3.622620, farming areas) respectively. Genitors coming from aquaculture areas were assumed to be exposed to stronger selection pressure due to mass mortality outbreaks occurring annually at these sites. The last family, F21, was generated from a pair of broodstocks derived from mass selection conducted in the field during four generations in an aquaculture area at the Atlantic coast (La Tremblade, lat 45.781741 long − 1.12191, 450 km from the Mediterranean and around 370 km from the Atlantic sampling sites) [47]. In order to minimize environmental effects on hemolymph microbiota establishment, the oyster families were maintained in hatchery under controlled biosecured conditions before transplantation experiments.

Transplantation experiments and sampling

In order to examine how oyster genetic background and environmental conditions affect the oyster microbiome in health and disease, we performed transplantation experiments before and during a mortality outbreak. Each of the five oyster families was divided into three batches (~ 300 oysters per batch); one was kept in the hatchery (at 16 °C) and the other two were transplanted for 5 days into an oyster farming area, the Thau lagoon (Lat. 43.378888 log. 3.571111), before and during a mortality outbreak. For transportation, oysters were packed in polystyrene containers and kept moist by covering with a damp cloth (duration of transport was less than 40 h). We previously showed that of 5 days of transplantation in the natural environment was sufficient to allow the oysters to be infected and precede massive mortality [48]. The non-infectious period (March 2016), was characterized by low average seawater temperature (10.6 °C), whereas the temperature during the infectious period (May 2016) was 18.7 °C. As microbiota may be affected by temperature stress [44], the oysters were progressively acclimated (2 °C/day) to Thau lagoon temperature before transplantation. Five days after the transplantation, the oysters were transported out of water to IHPE laboratory facilities located at 40 km, and hemolymph was drawn from the adductor muscle using a 1 mL syringe and 23 x ¼ needle. We sampled three replicates (of at least 30 oysters) per oyster family and per condition (hatchery, Thau non-infectious and Thau infectious). Hemolymph aliquots were kept on ice during the sampling. Hemocytes were removed from the hemolymph by filtration (5 μm membranes, Sartorius Minisart). Hemocyte-free hemolymph samples were filtered on 0.2 μm membranes, then the membranes and the filtrates were flash-frozen in liquid nitrogen and stored at − 80 °C. The membranes were subsequently used to extract DNA for bacterial and protist microbiota analysis. The filtrates were used to quantify virus-like particle and to extract DNA for viral microbiota analysis.

Seawater samples (10 L), were collected from each environment, filtered on 20 μm membranes to remove large particles in suspension and microorganisms collected by filtration as for hemolymph samples except that viruses contained in 0.2 μm filtrates were concentrated by iron chloride flocculation [49].

DNA extraction and sequence processing for bacterial and protist metabarcoding

Total genomic DNA was extracted from 0.2 μm filters used to collect bacteria and protists from hemolymph and seawater using the NucleoSpin® Tissue Genomic extraction kit (Macherey-Nagel). Amplification and sequencing were performed by Genome Quebec Company (Genome Quebec Innovation Center, McGill University, Montreal, Canada). PCR amplifications were performed using primers 341F (5′-CCTACGGGNGGCWGCAG-3′) and 805R (5′-GACTACHVGGGTATCTAATCC-3′) targeting the V3-V4 region of the bacterial 16S rRNA genes [50]. For the V1-V2 region amplification of the eukaryotic 18S rRNA genes, we used the eukaryotic forward primer NEW EUK F (5′-ACCTGGTTGATCCTGCCA-3′) adapted from the EukA primer described by Medlin et al. [51] that we shortened by two nucleotides in 3′-end to catch more diversity. Based on a selection of complete sequences covering major eukaryotic lineages extracted from the PR2 database [52], we observed that this primer will target most eukaryotic lineages (93%), including oysters, with the notable exception of some metazoan (human for instance) and Microsporidia (Additional file 1: Table S1). The reverse primer NEW EUK R (5′-GTARKCCWMTAYMYTACC-3′) was newly designed with seven degenerated nucleotides for targeting major eukaryotic lineages excluding metazoans (especially oysters). Metazoans have ≥3 mismatches with this primer, with the notable exception of Cnidaria, Porifera and Ctenophora which have no mismatch (Additional file 1: Table S1). This primer has no mismatch with most of other eukaryotic lineages (about 77% of sequences based on our selection), with few exceptions, notably Parabasalia and Microsporidia. Paired-end sequencing with 250 bp read length was performed on a MiSeq system (Illumina) using the v2 chemistry according to the manufacturer’s protocol. Pre-process analyses and clustering was performed on FROGS pipeline (Find Rapidly OUT with Galaxy Solution) [53]. In brief, after denoising and primer/adapter removal with cutadapt, clustering was performed using SWARM, which uses a two-step clustering algorithm with a threshold corresponding to the maximum number of differences between two operational taxonomic units (OTU) (denoising step d = 1; aggregation distance = 3) [54]. To produce OTU and affiliation tables in the standard BIOM format, we performed an affiliation using Blast + against the Silva 16S rRNA database (release 128, Sept 2016) and the pr2 database (release 13, Sep 2017) for bacteria and protists respectively. OTU relative abundances and taxonomic affiliations of the bacteria and protists found in the seawater and the hemolymph are available in Additional file 2: Table S2 and Additional file 3: Table S3, respectively.

In order to verify that we did not introduce contamination during metabarcoding processes (from DNA extraction to sequencing), we made six independent DNA extractions blank controls using the NucleoSpin® Tissue Genomic extraction kit. As we hypothesized that the most probable source of contamination was of bacterial origin, amplification and sequencing were performed on these 6 samples using the primers targeting V3-V4 region of the bacterial 16S rRNA genes. At the end of the bio informatics analyses we obtained 15 OTUs and none of them was common to the six blank controls. As only two of these OTUs are found in some of the samples and represent on average 0.09% of the reads (45 reads/sample out of 50,000) and as we used abundance-weighted analysis methods, we considered that our DNA extraction and sequencing process do not introduce significant contamination.

DNA extraction and sequence processing for viral metagenomics

Viral particles were collected from the 48 samples (45 hemolymph samples and 3 seawater) by ultracentrifugation (257,000 g for 2 h at 4 °C, Beckman Optima). The pellets were re-suspended in sterile 1× PBS and treated with DNase I (10 U/μL RQ1 DNAse, Promega). Viral DNA was extracted using NucleoSpin® Virus kit (Macherey Nagel), and DNA amplification was run with a GenomiPhi™ V3 Kit (GE Healthcare). DNA concentration was measured using Qubit dsDNA HS assay kit (Invitrogen). Genomes were fragmented (average size 1225 bp for the 48 samples, ranking from 814 bp to 2227 bp), amplified and sequenced on a MiSeq system at the Genomer platform at Station Biologique in Roscoff. The 48 viral libraries corresponding to 92,891,278 paired-end 300 bp sequences (2,109,133 paired-end reads per library +/− 904,066) were cleaned by adapter removal and trimmed of low quality bases using Trimmomatic [55]. To remove host sequences, cleaned reads were mapped (Bowtie 2) to the C. gigas genome (GCA_000297895.1). rRNA contamination was removed using SortMeRNA software [56]. Full metavirome was assembled using MetaSPAdes [57] and coverage virome analyses were performed using the Bowtiebatch and Read2RefMapper scripts developed within the iVirus protocol [58]. Raw demultiplexed sequence data are available at NCBI SRA under the accession number PRJNA381401.

Viral contigs were annotated using VirSorter [59], Phaster [60] and BLAST against the Global Ocean Virome database (E-value <1e− 10). We further employed a procedure that improves short scaffold annotation [61]. Briefly, open reading frames (ORFs) were predicted with Prodigal [62] and aligned against nr database using BLASTP (E-value <1e− 5). We considered an ORF to be of viral origin based on the best blast hit annotation (BBH), or - in the case of bacteriophages - if one of the top 10 hits was a phage protein with one of the following phage-related functions: tail, coat, head, capsid, portal protein. Viral contig abundance and taxonomic affiliation in seawater and hemolymph are available in Additional file 4: Table S4.

VLP quantification and cell sorting

VLP (Virus-Like-Particle) abundances were determined for three 1.5 mL aliquots per sample by flow cytometry (adapted from [63]), using a FACSCanto II cytometer (Becton Dickinson, Franklin Lakes, NJ, U.S.A.) equipped with an air-cooled laser providing 15 mW at 488 nm with the standard filter set-up. Abundances of VLPs were calculated using specific calibrated Becton Dickinson Trucount™ beads. Samples were stained with SYBR Green I (S7563, Invitrogen; 2% final concentration) and incubated at 80 °C for 20 min. After cooling down to room temperature, 1 μL of mixed fluorescent 0.5 μm-diameter beads (Molecular Probes Inc., Eugene, OR, U.S.A.) were added as an internal standard to variable sample volume mixed with 0.02 μm filtered TE buffer [63] to reach a final dilution between 400 and 2000 times. The flow rate of the cytometer was set at low level (acquisition time: 1 min) and artificially induced fluorescence (FL1: 530 nm) was used to detect VLPs. The flow cytometric data were analyzed using the Diva software (Becton Dickinson).

Rarefaction analysis and subsampling of bacterial and protist datasets

We performed a rarefaction analysis of Good’s coverage and alpha diversity indices (Shannon’s H, evenness, Chao 1 and Observed OTUs) in QIIME [64] to determine an adequate subsampling depth for the bacterial and protist datasets. We chose the cutoff of 50,000 reads per sample for bacteria and 8000 for the protists (two samples with 7000+ reads were included as they were). We found high Good’s coverage for the above cutoff values as well as leveling-off of the rarefaction curves for Shannon’s H and evenness (Additional file 5: Figure S1).

Statistical analyses

All statistical analyses were performed in R [65] at the OTU level (or contigs for viruses). Alpha diversity was represented by Shannon’s H and Observed OTUs indices calculated in mothur [66] and analyzed by generalized least squares by maximum likelihood linear models implemented in the package nlme [67]. The variance structure was included in the model in case of the variance heterogeneity as determined by the Levene’s test from the package lawstat [68]. Beta diversity based on Bray-Curtis index was visualized by nonmetric multidimensional scaling - NMDS and statistically analyzed by Permanova, [69] using the package vegan [70].

For both alpha and beta diversity, we first compared the seawater and hemolymph microbiota in each environment. We subsequently focused on the hemolymph microbiota, fitting three different models. The first two models tested the differences between the families grouped either by genitors’ origin or selection pressure experienced by genitors in the hatchery and infectious environment, respectively. The third model tested for differences between the hatchery and each of the natural environments (i.e., hatchery vs. non-infectious and hatchery vs. infectious environment). In this third model, the families were included as random effects in case of alpha diversity or used as blocks to constrain permutations in case of beta diversity. For each model, we first checked if any of the tested factors were significant and, if this was the case, we subsequently inspected the differences between two sets of a priori defined orthogonal contrasts, genitors’ origin and selection pressure experienced by genitors (Additional file 6 Table S5). We also performed indicator species analysis (using R package indicspecies, [71]) on all three microbiotas in order to find taxa significantly changing with environmental conditions and depending on the family origin and phenotype (Additional file 7: Table S6). We included only taxa with mean relative abundance > 1% in the analysis, thus testing 2505 bacterial, 357 protist and 6788 viral taxa. We defined indicators as OTU/taxa with both specificity (component A of the index, probability that oyster belongs to the group if the indicator is found) and sensitivity (component B, probability that the oyster in the group hosts the indicator) ≥ 0.9.

Associations between bacterial and eukaryotic OTUs were inferred from an undirected co-occurrence network. Pairwise scores between OTUs were computed using Spearman’s rank correlations. Only co-occurrences corresponding to correlations with a coefficient (rho) > 0.7 and a statistical significance (P-value) < 0.001 were considered for further analysis. Non-random co-occurrence patterns were tested with the checkerboard score (C-score) under a null model preserving site frequencies [72, 73]. Network characterization for comparisons between hatchery, infectious and non-infectious environments (15 samples in each network) was performed using a set of overall network topological indices (i.e., average node connectivity, average path length, average clustering coefficient and modularity) [74, 75]. All analyses were run using the R packages vegan [70], igraph [76] and WGCNA [77].

Results and discussion

Transplantation experiments of oyster families with different resistance phenotypes to Pacific oyster mortality syndrome (POMS)

To investigate the role of the oyster genetic background and environment in the composition and dynamics of the hemolymph microbiota, we examined five genetically differentiated biparental oyster families with various POMS resistance phenotypes [35] (summarized in Fig. 1a) in controlled (hatchery) and natural (before and during disease outbreak) conditions (Fig. 1b). We observed mortality events only during the infectious period (Fig. 1c). The families produced from genitors collected in non-aquaculture areas had markedly lower survival rates (5.5 and 12.8% for F15 and F32, respectively) than the families from genitors collected in aquaculture areas (80.5% for F9 and 100% for F21 and F44). This difference justifies grouping the families into resistant (S+) and susceptible (S-), and supports our hypothesis that the selection pressure experienced by genitors determines the resistance to POMS.

Fig. 1

Production of oyster families with different level of resistance to Pacific oyster mortality syndrome (POMS). a Origin of genitors and production of the five biparental families. Four pairs of wild genitors were collected either in aquaculture (S+) or non-aquaculture (S-) areas on Atlantic (Atl.) or Mediterranean (Med.) coasts (blue and orange lines, respectively). The fifth pair of genitors (*) comes from a mass selection program (more details in method section). Genitors coming from aquaculture areas were assumed to be exposed to stronger selection pressure due to mass mortality outbreaks occurring annually at these sites. b Transplantation experiments. Oysters from the five families were raised in the hatchery under the same controlled conditions, and subsequently transplanted into the Thau lagoon either before (non-infectious period) or during a mortality outbreak (infectious period). c Survival curves for the five oyster families after transplantation into the Thau lagoon during the infectious period (a minimum of 125 oysters was used for survival tests). Day 5 corresponds to the end of the transplantation, when oysters were brought back to the laboratory. Resistant and susceptible oyster families are indicated by S+ and S-, respectively. **** p < 0.0001 (Mantel-Cox Log-rank test)

Hemolymph hosts a complex microbiota different from the environment and shaped by the oyster genetic background

Our results show that oyster hemolymph hosts a complex community of bacteria, protists and viruses, which differs from the microbial community in the seawater column not only at the OUT level (Additional file 8: Table S7), but also at higher taxonomic levels (Fig. 2a-c). Compared to the seawater from the hatchery tank, oyster hemolymph was characterized by a higher proportion of Proteobacteria (hemolymph: 78.1% ± 1.5 SEM, seawater: 53.4%, Fig. 2a) and Alveolata (hemolymph: 68.4% ± 3.5 SEM, seawater: 17.0%, Fig. 2b), and a lower proportion of Actinobacteria (hemolymph: 2.3% ± 0.2 SEM, seawater: 15%) and viruses of the “other phages” category (hemolymph: 3.1% ± 1.1 SEM, seawater: 23.3%, Fig. 2c). Differences in bacterial, protist and viral seawater and hemolymph communities were also observed after transplantation into a natural environment (Fig. 2d-i). The fact that the hemolymph and seawater microbiota differ already at a low taxonomic resolution indicates that the hemolymph is indeed a distinct environment. Whereas the distinctness of hemolymph bacterial microbiota has been previously reported [45, 78] here we also show that the hemolymph hosts rich and distinct protist and viral microbial communities as well.

Fig. 2

Relative abundance of bacterial classes, protistan phyla and viral taxa in oyster hemolymph and seawater in controlled and natural conditions. a-c) Hatchery, d-f) natural environment during non-infectious and g-i) infectious period. Bacterial classes and protistan phyla with relative abundance < 1% are grouped as “other”. The viral taxonomic affiliation was done on the 200 most abundant scaffolds, which represented 42.6% of the total reads (c, f, i). SW: seawater

To determine if oyster genetic background can influence hemolymph microbiota composition, we compared oysters with different genetic (from the five biparental families) that were maintained in the hatchery from fertilization until the age of one year. The response of hemolymph microbiota to environmental changes is mostly restricted to fine taxonomic scales [44] and we therefore examined the diversity and composition of the hemolymph microbiome at the OTU level (or contig level for viruses). We were particularly interested if the oyster families of the same origin (Atlantic or Mediterranean) and phenotype (resistant S+ and susceptible S-) are more similar to each other despite different genotypes. We therefore defined two sets of a priori contrasts, in order to explore the observed differences from both the aspect of origin and phenotype (Additional file 6 Table S5). Whereas we found no difference in viral alpha diversity between the families, we observed some differences in Shannon’s H indexes for bacterial (F4,10 = 6.00, p = 0.010) and protist (F4,10 = 4.01, p = 0.034) communities (Fig. 3a-c; Additional file 9: Table S8). Specifically, differences in bacterial alpha diversity were restricted to resistant families, with the Atlantic resistant families (F9 and F21) having on average lower Shannon’s H than the Mediterranean family (F44) (0.56, 95%CI: (0.137, 0.985)). However, this difference was driven by very low diversity of the F9 family, which was also significantly lower than that of F21, the other Atlantic resistant family (0.66, 95%CI: (0.290, 1.036)). Conversely, the differences in protist Shannon’s H were restricted to Mediterranean families, where the resistant family (F44) had higher diversity than the susceptible one (F32) (0.27, 95%CI: (0.117, 0.419)). Regarding beta-diversity, Permanova analyses revealed that 16% of bacterial and 10% of viral community variability was explained by selection pressure experienced by genitors (Fig. 4a, d and c, f, Additional file 10: Table S9). On the other hand, we found no significant general effect of selection pressure experienced by genitors on the protist communities (S+/S-), while 12% of variability was explained by genitors’ origin (Atl/Med) (Fig. 4b, e, Additional file 10: Table S9). Our findings support the previously reported association between the host genotype and oyster gill and hemolymph bacterial communities [45, 46], which we now extend to protist and viral microbiota.

Fig. 3

Alpha diversity expressed as Shannon’s H index for bacteria, protists and viruses in the hemolymph and sweater in controlled and natural conditions. a-c) Hatchery, d-f) natural non-infectious and g-i) infectious environment. Bacterial and protistan diversity estimates are based on OTUs, viral on the recovered contigs. The Shannon’s H indices of bacterial (a, d, g), protist (b, e, h) and viral (c, f, i) communities are presented on the same row. Error bars represent standard error of the mean

Fig. 4

Effect of oyster genotype on hemolymph microbiota. NMDS plots represent beta diversity based on Bray-Curtis dissimilarities for the bacterial (a), protistan (b) and viral (c) communities of oysters born and kept in the hatchery under controlled conditions. Bar plots are graphical representations of Permanova results showing the variation explained by oyster genetic background (genitors’ origin and selection pressure experienced by genitors) for the bacterial (d), protistan (e) and viral (f) communities. Both bars represent the same analysis (model: ~ Family), but with two different sets of orthogonal contrasts (see Additional file 6 Table S5). The terms which are not significant are shown as shaded rectangles (for completeness), but should not be interpreted. Table on which the bar plots are based can be found in Additional File 10: Table S9. Bacterial and protistan diversity estimates are based on OTUs, viral on contigs

Dynamics of microbial consortia during transplantation into a natural environment is mainly shaped by environmental factors

To investigate how the host and environmental factors influence oyster hemolymph microbiota, we compared the microbiota of oyster biparental families after transplantation from the hatchery into the natural environment during non-infectious and infectious periods. Transplantation into natural environments differentially affected the alpha diversity of the three microbial components. According to Shannon’s H index, transplantation did not affect the alpha diversity of bacterial microbiota at all (Fig. 3a, d and g, Additional file 11: Table S10). Although the Shannon’s H of F9 family seems to dramatically increase after the transplantation in the non-infectious environment, neither the interaction terms (environment x genitors’ origin or environment x selection-pressure) nor family random effects were significant, likely due to overall high within-group variability (Additional file 11: Table S10). Average diversity of viral communities increased after transplantation into the non-infectious environment (1.49, 95%CI: (0.977, 2.006)), whereas it remained the same after transplantation into the infectious environment (Fig. 3c, f and i, Additional file 11: Table S10). However, we found significant environment x origin (− 0.59, 95%CI: (− 1.167, − 0.016)) and environment x origin x selection-pressure (0.58, 95%CI: (0.007, 1.158)) interactions for viral microbiota, probably reflecting a large diversity drop in the F15 family during the infectious period. Finally, the average alpha diversity of protist microbiota was lower during the non-infectious period (− 0.45, 95%CI: (− 0.704, − 0.204)) and higher during the infectious period (0.72, 95%CI: (0.467, 0.966)) than in the hatchery (Fig. 3b, e and h, Additional file 11: Table S10). Although we previously observed lower diversity of hemolymph bacterial microbiota coupled with proliferation of opportunistic bacterial pathogens of the genus Arcobacter following a Vibrio sp. infection [44], we see no link between OsHV-1 μVar infection and bacterial diversity, suggesting that the oyster bacterial microbiota responds differently to bacterial and viral infection.

For all three community types, NMDS plots based on Bray-Curtis dissimilarities revealed grouping according to the environment type (Fig. 5a-c). Permanova confirmed a huge effect of the environment, with 37, 48 and 35% of the variability in bacterial, protist and viral communities, respectively, explained by the environment (Fig. 5d-f; Additional file 12: Table S11), indicating that the composition and dynamics of hemolymph microbial consortia are mainly shaped by environmental factors. This was also supported by the indicator species analysis, as we identified a number of taxa specific to each environment, but found virtually no evidence of origin or selection-pressure signature taxa (Additional File 7:Table S6). Despite the dominant effect of the environment, genitor origin and selection pressure experienced by genitors each explained a low (2–3%) but significant amount of variability in the structure of all three microbial components, suggesting that oyster-related factors do play a role in the microbiota assembly (Fig. 5d-f; Additional file 12: Table S11).

Fig. 5

Effect of oysters genetic background and environmental conditions on oyster hemolymph microbiota structure in controlled and natural conditions. NMDS plots of Bray-Curtis dissimilarities for the bacterial (a), protistan (b) and viral (c) communities in the hatchery, and natural non-infectious and infectious environments. Bar plots are graphical representations of Permanova results for the bacterial (d), protistan (e) and viral (f) communities. Permanova model fit here is ~ environment*genitors’ origin *selection pressure (table on which the barplots are based can be found in Additional File 12: Table S11). The terms which are not significant are shown as shaded rectangles (for completeness), but should not be interpreted. Each term is described by the factor or interaction (interactions are marked by *) it represents (i.e. environment, environment * genitors’ origin …), followed by a specific contrast after the colon. Each term is more similar to a linear model coefficient (although it describes the variability explained), and is different from the shared amount of variability explained in redundancy analysis. Bacterial and protistan diversity estimates are based on OTUs, viral on contigs

Our results corroborate previous findings showing that the hemolymph bacterial microbiota is mainly influenced by environmental factors [45], and extend these to protists and viruses.

Trans-kingdom co-occurrence analysis shows disruption of microbial network after transplantation

Co-occurrence network properties can be used to study community stability across different conditions [74, 79, 80]. To assess the influence of the environment on network properties, we constructed separate microbial networks for hatchery, non-infectious and infectious periods. Only bacteria and protists were considered in the network, due to incomplete assembly of viral genomes that can lead to generation of “artefactual networks” (as several contigs may correspond to the same virus). Network comparisons are commonly achieved by calculating the average of four network indices for each node (i.e., taxon): (i) connectivity or degree distribution, representing the number of edges of a node toward other nodes; (ii) path length (the shortest path between two nodes); (iii) clustering coefficient, describing how well a node is connected to its neighbors; and (iv) modularity, a measure of the structuration of the network into modules [81].

Overall, we detected significant differences in terms of average connectivity, average shortest path and average clustering coefficient between the environments (Fig. 6a-c, note that x-axes are not on the same scale). Particularly, these results point to a disruption of the microbial network after transplantation into a natural environment during both non-infectious and infectious periods (Fig. 6b and c). Indeed, the higher values for average shortest path in non-infectious and infectious networks may imply lower efficiency and speed in the transmission of information, energy or material in the system, potentially altering the response of the microbial community to environmental perturbations [81, 82]. In addition, the lower connectivity and lower average clustering coefficient indicate a potential perturbation of the network and the loss of important community members [83]. A lower clustering coefficient is indeed an indication of the loss of highly connected nodes (i.e., hubs), which have been related to the concept of keystone species [84, 85].

Fig. 6

Relationship between closeness and betweenness centrality values for microbial co-occurrence networks in controlled and natural conditions. The plots show closeness vs. betweeness centrality values for nodes in the hatchery (a), non-infectious (b) and infectious (c) networks based on bacterial and protist OTUs. Only the 20 nodes with the highest degree value are represented for each network. Color indicates phylum or class-level taxonomic assignment of the nodes. General network properties are specified in the lower left corner of each plot. A superscript ‘a’ indicates a significant difference (p < 0.01) between the hatchery and infectious or non-infectious general network properties. A superscript ‘b’ indicates a significant difference (p < 0.01) between infectious and non-infectious general network properties

Calculation of network indices for individual nodes, particularly those related to the concept of keystone species [86], demonstrated this loss. Specifically, nine OTUs from the Bacteroidetes and Planctomycetes phyla scored among the 20 highest values of betweenness centrality, degree and closeness centrality in the hatchery environment, but disappeared in the non-infectious and infectious networks. The loss of such keystone species, the ‘backbone’ of the community [85], could lead to the fragility of the affected network, ultimately leading to its fragmentation and to secondary extinctions of species relying on the hubs [87]. OTUs playing a pivotal role in the structuration of the three networks were mainly dominated by members of the Proteobacteria phylum. OTU-6 (Gamaproteobacteria, Alteromonadales, Colwelliaceae, Colwellia) and OTU-7 (Alphaproteobacteria, Rhodobacterales, Rhodobacteraceae, Planktomarina) were among the keystone species in the hatchery network, but only OTU-7 was identified in all networks (Fig. 6a-c). Bacteria of the Roseobacter clade are known for their probiotic properties and protective effects against Vibrio infections in marine vertebrates [88] and invertebrates [89].

Hemolymph microbiota dynamics in the context of Pacific oyster mortality syndrome

To evaluate the role of hemolymph microbiota in health and disease, we further investigated the effect of oyster genetic background on β-diversity during the infectious period. Whereas we did not observe any differences in the structure of bacterial communities between the families, we found that selection pressure experienced by genitors explained 10% (p < 0.05) of the variability in protist microbiota (Additional file 13: Table S12). Closer inspection revealed that this effect was mainly due to the difference between two Mediterranean families (F32 and F44, Additional file 13: Table S12, compare S+/S- selection pressure contrast and Med: F44/F32 genitors’ origin contrast). For viruses, oyster genetic background explained 69% of the variation in the community structure (Fig. 7a; Additional file 13: Table S12). Interestingly, the differences between Atlantic resistant families accounted for the major part of the variability explained (F9/F21, 33%), which could be explained by the fact that F9 was produced from wild genitors collected in farming area whereas F21 was obtained by a breeding program using mass selection [47, 90]. Still, a considerable portion of variation was explained by the differences between the resistant and susceptible families (S+/S-: 12%; F32/F44: 12%; F15/F9&F21: 16%), as well as between the susceptible families (F32/F15: 16%). Indicator species analysis showed that susceptible families in the infectious environment were characterized by high abundance of OsHV-1 μVar (Fig. 2i; Additional File 7:Table S6) demonstrating the key role of host genetic background in resistance to OsHV-1 μVar as observed in a previous genome-wide association study [91].

Fig. 7

Viral beta diversity and VLP load quantification in oyster hemolymphs during the infectious period. a Bar plots are graphical representations of Permanova results showing the amount of variation explained by families grouped by genitors’ origin or selection pressure (for additional explanations see the legend to Fig. 4; table on which the barplots are based can be found in Additional File 13: Table S12). b Representative flow cytometry dot plots obtained for resistant (F9, F21 and F44) and susceptible (F15 and F32) families. For families F15 and F32, during the infectious period, we observed an additional VLP population characterized by a higher DNA content (red arrows). c) OsHV-1 μVar loads in the oyster hemolymph of families F15 and F32 (mean SEM)

To accurately quantify OsHV-1 μVar load within the oyster hemolymph, we used a flow cytometry technique to quantify virus-like particles (VLPs). Flow cytometry dot plots show an additional VPL population in susceptible families (red arrows in Fig. 7b), representing 62.4% (± 1.95 SEM) and 20.3% (± 2.56 SEM) of the overall VLP population in families F15 and F32, respectively. This VLP population was sorted from the hemolymph of the two susceptible families and sequenced; in both cases it corresponded to OsHV-1 μVar (Fig. 7c); corroborating the metagenomic analysis results and thus the role of the genetic background in the control of the OsHV-1 μVar dynamics within oyster hemolymph. In addition, high OsHV-1 μVar concentrations, reaching up to 2.8 × 108 viruses per mL hemolymph, were found only in susceptible families, supporting the previous findings that viral replication occurred in hemocytes [92, 93]. Moreover, we recently showed that the hemocyte infection by OsHV-1 μVar followed by intense replication was an initial and necessary step for disease development in susceptible families, while resistant families successfully controlled viral replication [35]. We hypothesize that the intense viral replication within hemocytes leads to the liberation of OsHV-1 μVar into the circulatory system and its spread to all oyster tissues.

Interestingly, we found no association between oyster susceptibility and the occurrence of particular bacterial taxa. However, oyster colonization by opportunistic bacteria was shown, on the base of whole animal analyses, to be a necessary second step of the infectious process leading to the death of oysters [35], supporting previous findings that identified bacteria as important etiological agent of the disease [48]. The bacteria most frequently associated with oyster mortality belong to the genera Vibrio, Arcobacter, Marinobacterium, Fusibacter, Psychromonas and Psychrobium [32, 35, 44, 94]. Although the bacteria belonging to these genera characterized oyster hemolymph microbiota during the infectious period, they were not specifically linked to susceptible animals that experienced mortality (Additional File 7: Table S6). However, previous findings were based on solid tissues, whose composition and dynamics differ from those of the hemolymph microbiota [45], which may explain the difference in response to infection. Alternatively, it is possible that these taxa proliferate only during a particular period of disease that we did not cover with our sampling.


Deciphering the complex interplay between the host, its microbiota and the environment is a vital part of understanding dynamics and outcome of polymicrobial infections, such as the Pacific Oyster Mortality Syndrome. Our study shows that oyster hemolymph is a complex ecosystem strongly influenced by environmental conditions and hosting diverse bacteria, protists and viruses. However, each of these microbiota components is weakly but significantly shaped by oyster genetics, indicating that they indeed interact with the oyster host and that holobiont framework is a helpful concept describing the oyster-microbes community [93]. More studies would be needed to determine if certain microbiota components could be part of the holobiont, potentially affecting fitness of their host, or if they only constitute inconstant environmental component [95].

Interestingly, despite the shared environment, each microbiota components seems to respond differently and independently of the others to environmental conditions. Contrary to previous findings [35, 44, 48, 94], we observed no proliferation of opportunistic bacterial pathogens following the invasion by OsHV-1 μVar, which highlights the need to more closely investigate the mechanisms and dynamics of transkingdom interactions within the hemolymph ecosystem in order to explain this discrepancy.

Finally, the bacterial microbiome displayed stable high-level taxonomic composition and bacteria were the only organisms identified as important in the microbial networks, suggesting a higher stability and more fine-tuned interactions with the oyster host compared to protists. Therefore, it seems that at least some bacteria may be really adapted to the hemolymph conditions, whereas viral and protist communities are primarily transient. However, most of the viral diversity we observe here is dark matter [96] and it is possible that better characterization of oyster-associated viruses may reveal a different picture.

Availability of data and materials

The datasets generated during the current study are available in the Sequence Read Archive repository under BioProject ID PRJNA381401 and SRP105754.



Biological Observation Matrix


Find Rapidly OUT with Galaxy Solution

OsHV-1 μVar:

Ostreid herpesvirus variants


Operational Taxonomic Unit


Pacific Oyster Mortality Syndrome


Nonmetric Multidimensional Scaling


Quantitative Insights into Microbial Ecology


Standard Error of the Mean


Virus Like Particle


Weighted Gene Co-expression Network Analysis


  1. 1.

    Rosenberg E, Koren O, Reshef L, Efrony R, Zilber-Rosenberg I. The role of microorganisms in coral health, disease and evolution. Nat Rev Microbiol. 2007;5:355–62.

    CAS  PubMed  Article  Google Scholar 

  2. 2.

    Rosenberg E, Zilber-Rosenberg I. The hologenome concept of evolution after 10 years. Microbiome. 2018;6:78.

    PubMed  PubMed Central  Article  Google Scholar 

  3. 3.

    Chabe M, Lokmer A, Segurel L. Gut Protozoa: friends or foes of the human gut microbiota? Trends Parasitol. 2017;33:925–34.

    PubMed  Article  Google Scholar 

  4. 4.

    Breitbart M, Rohwer F. Method for discovering novel DNA viruses in blood using viral particle selection and shotgun sequencing. Biotechniques. 2005;39:729–36.

    CAS  PubMed  Article  Google Scholar 

  5. 5.

    Moustafa A, Xie C, Kirkness E, Biggs W, Wong E, Turpaz Y, Bloom K, Delwart E, Nelson KE, Venter JC, Telenti A. The blood DNA virome in 8,000 humans. PLoS Pathog. 2017;13:e1006292.

    PubMed  PubMed Central  Article  CAS  Google Scholar 

  6. 6.

    Rascovan N, Duraisamy R, Desnues C. Metagenomics and the human Virome in asymptomatic individuals. Annu Rev Microbiol. 2016;70:125–41.

    CAS  PubMed  Article  Google Scholar 

  7. 7.

    Keen EC, Dantas G. Close encounters of three kinds: bacteriophages, commensal Bacteria, and host immunity. Trends Microbiol. 2018;26:943–54.

    CAS  PubMed  PubMed Central  Article  Google Scholar 

  8. 8.

    Mirzaei MK, Maurice CF. Menage a trois in the human gut: interactions between host, bacteria and phages. Nat Rev Microbiol. 2017;7:397–408.

    Article  CAS  Google Scholar 

  9. 9.

    Pfeiffer JK, Virgin HW. Transkingdom control of viral infection and immunity in the mammalian intestine. Science. 2016;351:aad5872.

    PubMed  Article  CAS  Google Scholar 

  10. 10.

    Hester ER, Barott KL, Nulton J, Vermeij MJ, Rohwer FL. Stable and sporadic symbiotic communities of coral and algal holobionts. ISME J. 2016;10:1157–69.

    CAS  PubMed  Article  Google Scholar 

  11. 11.

    Peixoto RS, Rosado PM, Leite DCD, Rosado AS, Bourne DG. Beneficial microorganisms for corals (BMC): proposed mechanisms for coral health and resilience. Front Microbiol. 2017;8:341.

    PubMed  PubMed Central  Article  Google Scholar 

  12. 12.

    Bosch TC, Grasis JA, Lachnit T. Microbial ecology in Hydra: why viruses matter. J Microbiol. 2015;53:193–200.

    PubMed  Article  Google Scholar 

  13. 13.

    Gurung K, Wertheim B, Salles JF. The microbiome of pest insects: it is not just bacteria. Entomologia Experimentalis Et Applicata. 2019;167:156–70.

    Article  Google Scholar 

  14. 14.

    Harvey E, Rose K, Eden JS, Lawrence A, Doggett SL, Holmes EC. Identification of diverse arthropod associated viruses in native Australian fleas. Virology. 2019;535:189–99.

    CAS  PubMed  Article  Google Scholar 

  15. 15.

    Pettersson JH, Shi M, Eden JS, Holmes EC, Hesson JC. Meta-Transcriptomic comparison of the RNA Viromes of the mosquito vectors Culex pipiens and Culex torrentium in northern Europe. Viruses. 2019;11:1033.

  16. 16.

    Chiu L, Bazin T, Truchetet ME, Schaeverbeke T, Delhaes L, Pradeu T. Protective microbiota: from localized to long-reaching co-immunity. Front Immunol. 2017;8:1678.

    PubMed  PubMed Central  Article  CAS  Google Scholar 

  17. 17.

    Masson F, Lemaitre B. Symbiosis: protection from within. Elife. 2017;6:e24111.

    PubMed  PubMed Central  Article  Google Scholar 

  18. 18.

    Selber-Hnatiw S, Rukundo B, Ahmadi M, Akoubi H, Al-Bizri H, Aliu AF, Ambeaghen TU, Avetisyan L, Bahar I, Baird A, et al. Human gut microbiota: toward an ecology of disease. Front Microbiol. 2017;8:1265.

    PubMed  PubMed Central  Article  Google Scholar 

  19. 19.

    Sweet MJ, Bulling MT. On the importance of the microbiome and Pathobiome in coral health and disease. Front Mar Sci. 2017;4(9).

  20. 20.

    Vonaesch P, Anderson M, Sansonetti PJ. Pathogens, microbiome and the host: emergence of the ecological Koch's postulates. FEMS Microbiol Rev. 2018;42:273–92.

    CAS  PubMed  Article  Google Scholar 

  21. 21.

    Paul-Pont I, Dhand NK, Whittington RJ. Spatial distribution of mortality in Pacific oysters Crassostrea gigas: reflection on mechanisms of OsHV-1 transmission. Dis Aquat Org. 2013;105:127–38.

    PubMed  Article  Google Scholar 

  22. 22.

    EFSA PoAHW: Oyster mortality. EFSA J 2015, 13:4122.

  23. 23.

    Barbosa Solomieu V, Renault T, Travers MA. Mass mortality in bivalves and the intricate case of the Pacific oyster, Crassostrea gigas. J Invertebr Pathol. 2015;131:2–10.

    PubMed  Article  Google Scholar 

  24. 24.

    Martenot C, Oden E, Travaille E, Malas JP, Houssin M. Detection of different variants of Ostreid Herpesvirus 1 in the Pacific oyster, Crassostrea gigas between 2008 and 2010. Virus Res. 2011;160:25–31.

    CAS  PubMed  Article  Google Scholar 

  25. 25.

    Renault T, Moreau P, Faury N, Pepin JF, Segarra A, Webb S. Analysis of clinical Ostreid Herpesvirus 1 (Malacoherpesviridae) specimens by sequencing amplified fragments from three virus genome areas. J Virol. 2012;86:5942–7.

    CAS  PubMed  PubMed Central  Article  Google Scholar 

  26. 26.

    Segarra A, Pepin JF, Arzul I, Morga B, Faury N, Renault T. Detection and description of a particular Ostreid herpesvirus 1 genotype associated with massive mortality outbreaks of Pacific oysters, Crassostrea gigas, in France in 2008. Virus Res. 2010;153:92–9.

    CAS  PubMed  Article  Google Scholar 

  27. 27.

    Abbadi M, Zamperin G, Gastaldelli M, Pascoli F, Rosani U, Milani A, Schivo A, Rossetti E, Turolla E, Gennari L, et al. Identification of a newly described OsHV-1 microvar from the North Adriatic Sea (Italy). J Gen Virol. 2018;99:693–703.

    CAS  PubMed  PubMed Central  Article  Google Scholar 

  28. 28.

    Burioli EA, Prearo M, Riina MV, Bona MC, Fioravanti ML, Arcangeli G, Houssin M. Ostreid herpesvirus type 1 genomic diversity in wild populations of Pacific oyster Crassostrea gigas from Italian coasts. J Invertebr Pathol. 2016;137:71–83.

    CAS  PubMed  Article  Google Scholar 

  29. 29.

    Burioli EAV, Prearo M, Houssin M. Complete genome sequence of Ostreid herpesvirus type 1 microVar isolated during mortality events in the Pacific oyster Crassostrea gigas in France and Ireland. Virology. 2017;509:239–51.

    CAS  PubMed  Article  Google Scholar 

  30. 30.

    Lynch SA, Carlsson J, Reilly AO, Cotter E, Culloty SC. A previously undescribed ostreid herpes virus 1 (OsHV-1) genotype detected in the pacific oyster, Crassostrea gigas, in Ireland. Parasitology. 2012;139:1526–32.

    CAS  PubMed  Article  Google Scholar 

  31. 31.

    Peeler EJ, Reese RA, Cheslett DL, Geoghegan F, Power A, Thrush MA. Investigation of mortality in Pacific oysters associated with Ostreid herpesvirus-1 mu Var in the Republic of Ireland in 2009. Preventive Veterinary Med. 2012;105:136–43.

    Article  Google Scholar 

  32. 32.

    Lemire A, Goudenege D, Versigny T, Petton B, Calteau A, Labreuche Y, Le Roux F. Populations, not clones, are the unit of vibrio pathogenesis in naturally infected oysters. ISME J. 2014;9:1523–31.

    PubMed  PubMed Central  Article  Google Scholar 

  33. 33.

    Bruto M, James A, Petton B, Labreuche Y, Chenivesse S, Alunno-Bruscia M, Polz MF, Le Roux F. Vibrio crassostreae, a benign oyster colonizer turned into a pathogen after plasmid acquisition. ISME J. 2017;11:1043–52.

    PubMed  Article  Google Scholar 

  34. 34.

    Faury N, Saulnier D, Thompson FL, Gay M, Swings J, Le Roux F. Vibrio crassostreae sp. nov., isolated from the haemolymph of oysters (Crassostrea gigas). Int J Syst Evol Microbiol. 2004;54:2137–40.

    CAS  PubMed  Article  Google Scholar 

  35. 35.

    de Lorgeril J, Lucasson A, Petton B, Toulza E, Montagnani C, Clerissi C, Vidal-Dupiol J, Chaparro C, Galinier R, Escoubas J-M, et al. Immune-suppression by OsHV-1 viral infection causes fatal bacteraemia in Pacific oysters. Nat Commun. 2018;9:4215.

    PubMed  PubMed Central  Article  CAS  Google Scholar 

  36. 36.

    Li Z, Nicolae V, Akileh R, Liu T. A brief review of oyster-associated microbiota. Microbiol Res J Int. 2017;20:1–14.

    Article  Google Scholar 

  37. 37.

    Pierce ML, Ward JE. Microbial ecology of the Bivalvia, with an emphasis on the family Ostreidae. J Shellfish Res. 2018;37:793–806.

    Article  Google Scholar 

  38. 38.

    Lokmer A, Kuenzel S, Baines J, Wegner KM. The role of tissue-specific microbiota in initial establishment success of Pacific oysters. Environ Microbiol. 2016;18:970–87.

    CAS  PubMed  Article  Google Scholar 

  39. 39.

    Garnier M, Labreuche Y, Garcia C, Robert M, Nicolas JL. Evidence for the involvement of pathogenic bacteria in summer mortalities of the Pacific oyster Crassostrea gigas. Microb Ecol. 2007;53:187–96.

    CAS  PubMed  Article  Google Scholar 

  40. 40.

    Olafsen JA, Fletcher TC, Grant PT. Agglutinin activity in Pacific oyster (Crassostrea gigas) hemolymph following in vivo Vibrio anguillarum challenge. Dev Comp Immunol. 1992;16:123–38.

    CAS  PubMed  Article  Google Scholar 

  41. 41.

    Wendling CC, Batista FM, Wegner KM. Persistence, seasonal dynamics and pathogenic potential of Vibrio communities from pacific oyster hemolymph. PLoS One. 2014;9:e94256.

    PubMed  PubMed Central  Article  Google Scholar 

  42. 42.

    Bachere E, Rosa RD, Schmitt P, Poirier AC, Merou N, Charriere GM, Destoumieux-Garzon D. The new insights into the oyster antimicrobial defense: cellular, molecular and genetic view. Fish Shellfish Immunol. 2015;46(1):50–64.

  43. 43.

    Desriac F, Le Chevalier P, Brillet B, Leguerinel I, Thuillier B, Paillard C, Fleury Y. Exploring the hologenome concept in marine bivalvia: haemolymph microbiota as a pertinent source of probiotics for aquaculture. FEMS Microbiol Lett. 2014;350:107–16.

    CAS  PubMed  Article  Google Scholar 

  44. 44.

    Lokmer A, Wegner KM. Hemolymph microbiome of Pacific oysters in response to temperature, temperature stress and infection. ISME J. 2015;9:670–82.

    CAS  PubMed  Article  Google Scholar 

  45. 45.

    Lokmer A, Goedknegt MA, Thieltges DW, Fiorentino D, Kuenzel S, Baines JF, Wegner KM. Spatial and temporal dynamics of Pacific oyster Hemolymph microbiota across multiple scales. Front Microbiol. 2016;7:1367.

    PubMed  PubMed Central  Article  Google Scholar 

  46. 46.

    Wegner KM, Volkenborn N, Peter H, Eiler A. Disturbance induced decoupling between host genetics and composition of the associated microbiome. BMC Microbiol. 2013;13:252.

    PubMed  PubMed Central  Article  Google Scholar 

  47. 47.

    Dégremont L, Nourry M, Maurouard E. Mass selection for survival and resistance to OsHV-1 infection in Crassostrea gigas spat in field conditions: response to selection after four generations. Aquaculture. 2015;446:111–21.

    Article  Google Scholar 

  48. 48.

    Petton B, Bruto M, James A, Labreuche Y, Alunno-Bruscia M, Le Roux F. Crassostrea gigas mortality in France: the usual suspect, a herpes virus, may not be the killer in this polymicrobial opportunistic disease. Front Microbiol. 2015;6:686.

    PubMed  PubMed Central  Article  Google Scholar 

  49. 49.

    John SG, Mendez CB, Deng L, Poulos B, Kauffman AK, Kern S, Brum J, Polz MF, Boyle EA, Sullivan MB. A simple and efficient method for concentration of ocean viruses by chemical flocculation. Environ Microbiol Rep. 2011;3:195–202.

    CAS  PubMed  PubMed Central  Article  Google Scholar 

  50. 50.

    Herlemann DPR, Labrenz M, Jurgens K, Bertilsson S, Waniek JJ, Andersson AF. Transitions in bacterial communities along the 2000 km salinity gradient of the Baltic Sea. ISME J. 2011;5:1571–9.

    CAS  PubMed  PubMed Central  Article  Google Scholar 

  51. 51.

    Medlin L, Elwood HJ, Stickel S, Sogin ML. The characterization of enzymatically amplified eukaryotic 16s-like Rrna-coding regions. Gene. 1988;71:491–9.

    CAS  PubMed  Article  Google Scholar 

  52. 52.

    Guillou L, Bachar D, Audic S, Bass D, Berney C, Bittner L, Boutte C, Burgaud G, de Vargas C, Decelle J, et al. The Protist ribosomal reference database (PR2): a catalog of unicellular eukaryote small sub-unit rRNA sequences with curated taxonomy. Nucleic Acids Res. 2013;41:D597–604.

    CAS  PubMed  Article  Google Scholar 

  53. 53.

    Escudie F, Auer L, Bernard M, Mariadassou M, Cauquil L, Vidal K, Maman S, Hernandez-Raquet G, Combes S, Pascal G. FROGS: find, rapidly, OTUs with galaxy solution. Bioinformatics. 2018;34:1287–94.

    CAS  PubMed  Article  Google Scholar 

  54. 54.

    Mahe F, Rognes T, Quince C, de Vargas C, Dunthorn M. Swarm: robust and fast clustering method for amplicon-based studies. PeerJ. 2014;2:e593.

    PubMed  PubMed Central  Article  Google Scholar 

  55. 55.

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

    CAS  PubMed  PubMed Central  Article  Google Scholar 

  56. 56.

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

    CAS  PubMed  Article  Google Scholar 

  57. 57.

    Nurk S, Meleshko D, Korobeynikov A. Pevzner PA: metaSPAdes: a new versatile metagenomic assembler. Genome Res. 2017;5:824–34.

    Article  CAS  Google Scholar 

  58. 58.

    Bolduc B, Youens-Clark K, Roux S, Hurwitz BL, Sullivan MB. iVirus: facilitating new insights in viral ecology with software and community data sets imbedded in a cyberinfrastructure. ISME J. 2017;11:7–14.

    PubMed  Article  Google Scholar 

  59. 59.

    Roux S, Enault F, Hurwitz BL, Sullivan MB. VirSorter: mining viral signal from microbial genomic data. PeerJ. 2015;3:e985.

    PubMed  PubMed Central  Article  CAS  Google Scholar 

  60. 60.

    Arndt D, Grant JR, Marcu A, Sajed T, Pon A, Liang YJ, Wishart DS. PHASTER: a better, faster version of the PHAST phage search tool. Nucleic Acids Res. 2016;44:W16–21.

    CAS  PubMed  PubMed Central  Article  Google Scholar 

  61. 61.

    Bettarel Y, Halary S, Auguet JC, Mai TC, Van Bui N, Bouvier T, Got P, Bouvier C, Monteil-Bouchard S, Christelle D. Corallivory and the microbial debacle in two branching scleractinians. ISME J. 2018;4:1109–26.

    Article  CAS  Google Scholar 

  62. 62.

    Hyatt D, et al. Prodigal: prokaryotic gene recognition and translation initiation site identification. BMC Bioinformatics. 2010;11:119.

  63. 63.

    Mojica KDA, Evans C, Brussaard CPD. Flow cytometric enumeration of marine viral populations at low abundances. Aquat Microb Ecol. 2014;71:203–9.

    Article  Google Scholar 

  64. 64.

    Caporaso JG, Kuczynski J, Stombaugh J, Bittinger K, Bushman FD, Costello EK, Fierer N, Pena AG, Goodrich JK, Gordon JI, et al. QIIME allows analysis of high-throughput community sequencing data. Nat Methods. 2010;7:335–6.

    CAS  PubMed  PubMed Central  Article  Google Scholar 

  65. 65.

    R Core Team. R: a language and environment for statistical computing. Vienna: R Foundation for statistical Computing; 2017.

    Google Scholar 

  66. 66.

    Schloss PD, Westcott SL, Ryabin T, Hall JR, Hartmann M, Hollister EB, Lesniewski RA, Oakley BB, Parks DH, Robinson CJ, et al. Introducing mothur: open-source, platform-independent, community-supported software for describing and comparing microbial communities. Appl Environ Microbiol. 2009;75:7537–41.

    CAS  PubMed  PubMed Central  Article  Google Scholar 

  67. 67.

    Pinheiro J, Bates D, DebRoy S, Sarkar D, R Development Core Team: nlme: Linear and Nonlinear Mixed Effects Models. R package version 3.1–131. Availlabe at: edition; 2017.

  68. 68.

    Gastwirth JL, Gel YR, Wallace Hui WL, Lyubchich V, Miao W, Noguchi K. lawstat: Tools for Biostatistics, Public Policy, and Law. R package version 3.1 edition; 2017.

    Google Scholar 

  69. 69.

    Anderson MJ. A new method for non-parametric multivariate analysis of variance. Austral Ecol. 2001;26:32–46.

    Google Scholar 

  70. 70.

    Oksanen J, Blanchet FG, Kindt R, Legendre P, Minchin PR, O'Hara RB, Simpson GL, Solymos P, Henry M, Stevens H, Wagner H: vegan: Community Ecology Package. R package version 2.4–3. Available at: edition; 2017.

  71. 71.

    De Caceres M, Legendre P. Associations between species and groups of sites: indices and statistical inference. Ecology. 2009;90:3566–74.

    PubMed  Article  Google Scholar 

  72. 72.

    Gotelli NJ, McCabe DJ. Species co-occurrence: a meta-analysis of J. M Diamond’s assembly rules model. Ecology. 2002;83:2091–6.

    Article  Google Scholar 

  73. 73.

    Stone L, Roberts A. The checkerboard score and species distributions. Oecologia. 1990;85:74–9.

    PubMed  Article  Google Scholar 

  74. 74.

    Jeanbille M, Gury J, Duran R, Tronczynski J, Agogue H, Ben Said O, Ghiglione JF, Auguet JC. Response of Core microbial consortia to chronic hydrocarbon contaminations in coastal sediment habitats. Front Microbiol. 2016;7:1637.

    PubMed  PubMed Central  Google Scholar 

  75. 75.

    Newman MEJ. The structure and function of complex networks. SIAM Rev. 2003;45:167–256.

    Article  Google Scholar 

  76. 76.

    Csardi G, Nepusz T. The igraph software package for complex network research. Inter J Complex Syst. 2005. p. 1695.

  77. 77.

    Langfelder P, Horvath S. WGCNA: an R package for weighted correlation network analysis. Bmc Bioinformatics. 2008;9:559.

    PubMed  PubMed Central  Article  CAS  Google Scholar 

  78. 78.

    Vezzulli L, Stagnaro L, Grande C, Tassistro G, Canesi L, Pruzzo C. Comparative 16SrDNA gene-based microbiota profiles of the Pacific oyster (Crassostrea gigas) and the Mediterranean mussel (Mytilus galloprovincialis) from a shellfish farm (Ligurian Sea, Italy). Microb Ecol. 2018;75:495–504.

    CAS  PubMed  Article  Google Scholar 

  79. 79.

    Faust K, Raes J. Microbial interactions: from networks to models. Nat Rev Microbiol. 2012;10:538–50.

    CAS  PubMed  Article  Google Scholar 

  80. 80.

    Sun MY, Dafforn KA, Johnston EL, Brown MV. Core sediment bacteria drive community response to anthropogenic contamination over multiple environmental gradients. Environ Microbiol. 2013;15:2517–31.

    PubMed  Article  Google Scholar 

  81. 81.

    Deng Y, Jiang YH, Yang Y, He Z, Luo F, Zhou J. Molecular ecological network analyses. BMC Bioinformatics. 2012;13:113.

    PubMed  PubMed Central  Article  Google Scholar 

  82. 82.

    Zhou JZ, Deng Y, Luo F, He ZL, Tu QC, Zhi XY. Functional molecular ecological networks. Mbio. 2010;1:e00169–10.

    PubMed  PubMed Central  Google Scholar 

  83. 83.

    Montoya JM, Pimm SL, Sole RV. Ecological networks and their fragility. Nature. 2006;442:259–64.

    CAS  PubMed  Article  Google Scholar 

  84. 84.

    Steele JA, Countway PD, Xia L, Vigil PD, Beman JM, Kim DY, Chow CET, Sachdeva R, Jones AC, Schwalbach MS, et al. Marine bacterial, archaeal and protistan association networks reveal ecological linkages. ISME J. 2011;5:1414–25.

    PubMed  PubMed Central  Article  Google Scholar 

  85. 85.

    Paine RT. A note on trophic complexity and community stability. Am Nat. 1969;103:91–3.

    Article  Google Scholar 

  86. 86.

    Berry D, Widder S. Deciphering microbial interactions and detecting keystone species with co-occurrence networks. Front Microbiol. 2014;5:219.

    PubMed  PubMed Central  Article  Google Scholar 

  87. 87.

    Sole RV, Montoya JM. Complexity and fragility in ecological networks. Proc Royal Soc B-Biological Sci. 2001;268:2039–45.

    CAS  Article  Google Scholar 

  88. 88.

    Prol-Garcia M, Pintado J. Effectiveness of probiotic Phaeobacter Bacteria grown in biofilters against Vibrio anguillarum infections in the rearing of turbot (Psetta maxima) larvae. Mar Biotechnol. 2013;15:726–38.

    CAS  PubMed  Article  Google Scholar 

  89. 89.

    Karim M, Zhao WJ, Rowley D, Nelson D, Gomez-Chiarri M. Probiotic strains for shellfish aquaculture: protection of eastern oyster, Crassostrea Virginica, larvae and juveniles against bacterial challenge. J Shellfish Res. 2013;32:401–8.

    Article  Google Scholar 

  90. 90.

    Azema P, Lamy JB, Boudry P, Renault T, Travers MA, Degremont L. Genetic parameters of resistance to Vibrio aestuarianus, and OsHV-1 infections in the Pacific oyster, Crassostrea gigas, at three different life stages. Genet Sel Evol. 2017;49:23.

    PubMed  PubMed Central  Article  Google Scholar 

  91. 91.

    Gutierrez AP, Bearf TP, Hooper C, Stenton CA, Sanders MB, Paley RK, Rastas P, Bryrom M, Matika O, Houston RD. A genome-wide association study for host resistance to Ostreid Herpesvirus in Pacific oysters (Crassostrea gigas). bioRxiv. 2017;8(4):1273–80.

  92. 92.

    Martenot C, Gervais O, Chollet B, Houssin M, Renault T. Haemocytes collected from experimentally infected Pacific oysters, Crassostrea gigas: detection of ostreid herpesvirus 1 DNA, RNA, and proteins in relation with inhibition of apoptosis. PLoS One. 2017;12:e0177448.

    PubMed  PubMed Central  Article  CAS  Google Scholar 

  93. 93.

    Morga B, Faury N, Guesdon S, Chollet B, Renault T. Haemocytes from Crassostrea gigas and OsHV-1: a promising in vitro system to study host/virus interactions. J Invertebr Pathol. 2017;150:45–53.

    CAS  PubMed  Article  Google Scholar 

  94. 94.

    King WL, Jenkins C, Go J, Siboni N, Seymour JR, Labbate M. Characterisation of the Pacific oyster microbiome during a summer mortality event. Microb Ecol. 2018;77:502–12.

    PubMed  Article  CAS  Google Scholar 

  95. 95.

    Theis KR, Dheilly NM, Klassen JL, Brucker RM, Baines JF, Bosch TCG, Cryan JF, Gilbert SF, Goodnight CJ, Lloyd EA, et al. Getting the Hologenome Concept Right: an Eco-Evolutionary Framework for Hosts and Their Microbiomes. mSystems. 2016;1(2):e00028–16.

  96. 96.

    Krishnamurthy SR, Wang D. Origins and challenges of viral dark matter. Virus Res. 2017;239:136–42.

    CAS  PubMed  Article  Google Scholar 

Download references


We thank the staff of the Ifremer station of Argenton (LPI, PFOM), and the Comité Régional Conchylicole de Méditerranée (CRCM) for technical support in the collection of oyster genitors, reproduction and transplantation experiments. We also thank Lionel Dégremont (SGMM) for providing us the genitors used to generate the F21 family. We also thank Marc Leroy for technical assistance and Delphine Destoumieux-Garzon, Guillaume Charrière, Marlène Chiarello, Frédérique Le Roux, Camille Cleriss and Yves Desdevises for fruitful discussions.


The present study was supported by: the DHOF program of the UMR5244/IHPE (; the Ifremer project HémoMicrobiote; EMBRC-France, whose French state funds are managed by the ANR within the Investments of the Future (ANR-10-INBS-02); the ANR projects DECIPHER (ANR-14-CE19–0023) and ENVICOPAS (ANR-15-CE35–0004); the EU funded project VIVALDI (H2020 program, n°678589); by Ifremer; CNRS, Université de Montpellier; Université de Perpignan via Domitia; by the IHU Méditerranée Infection, Marseille, France, by the French Government under the «Investissements d’avenir» program (reference: Méditerranée Infection 10-IAHU-03) and by the Région Provence-Alpes-Côte d’Azur and by the European funding FEDER PRIMI. This study is set within the framework of the “Laboratoires d’Excellences (LABEX)” TULIP (ANR-10-LABX-41).

Author information




JME, MC, DLJ, MG, and GY were involved in the study concept and design. DS, ET, PB, GY, CM, DLJ and JME were involved in the collection of samples. DS, CE, TG, SC, PD, GL, DC, LS, BKJ and JME were involved in the experimental work. AL did the statistical analysis. AJC did the network analysis. AL, DS, and JME drafted the manuscript, and all authors contributed to critical revisions and approved the final manuscript.

Corresponding author

Correspondence to J.-M. Escoubas.

Ethics declarations

Ethics approval and consent to participate

Not applicable.

Consent for publication

Not applicable.

Competing interests

The authors declare no conflict of interest.

Additional information

Publisher’s Note

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

Supplementary information

Additional file 1: Table S1. Overview of specificity of primers used to amplify V1-V2 region of 18S rRNA genes.

Additional file 2: Table S2. Bacteria OUT table with sequence tag count per sample and taxonomic affiliations.

Additional file 3: Table S3. Protist OUT table with sequence tag count per sample and taxonomic affiliations.

Additional file 4: Table S 4. Viral contig table with sequence tag count per sample and taxonomic affiliations.

Additional file 5: Figure S1. Sample rarefaction curves of alpha diversity indices for bacterial (a) and protists (b) datasets. (PPTX 995 kb)

Additional file 6: Table S5. Orthogonal contrasts used to examine the differences between the oyster families.

Additional file 7: Table S6. Indicator taxa for the oyster environment, family origin and phenotype.

Additional file 8: Table S7. Generalized least squares by maximum likelihood linear models testing for differences in alpha diversity between the seawater and hemolymph microbiota within each environment.

Additional file 9: Table S8. Generalized least squares by maximum likelihood linear models testing for differences in alpha diversity between the families in the hatchery.

Additional file 10: Table S9. Permanova based on Bray-Curtis dissimilarities testing for differences in beta diversity between the families in the hatchery.

Additional file 11: Table S10. Linear mixed-effects models testing for effects of the environment, genitors’ origin, selection pressure and their interactions on alpha diversity.

Additional file 12: Table S11. Permanova based on Bray-Curtis dissimilarities testing for effects of the environment, genitors’ origin and selection pressure and their interactions on beta diversity.

Additional file 13: Table 12. Permanova based on Bray-Curtis dissimilarities testing for differences in beta diversity between the families in the infectious environment.

Rights and permissions

Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, 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 changes were made. 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

Reprints and Permissions

About this article

Verify currency and authenticity via CrossMark

Cite this article

Dupont, S., Lokmer, A., Corre, E. et al. Oyster hemolymph is a complex and dynamic ecosystem hosting bacteria, protists and viruses. anim microbiome 2, 12 (2020).

Download citation


  • Oyster genetic background
  • Hemolymph microbiota dynamics
  • Early-life microbiota
  • Trans-kingdom interactions
  • Crassostrea gigas
  • Within-host ecosystem