We first assessed the dynamics of the bacterioplankton community, a factor known to covariate with fish skin mucus communities [12-14]. PERMDISP and ANOSIM tests performed on unweighted UniFrac distances (Table 1), and Kruskal-Wallis tests performed on Simpson index (Supplementary data, Figure S1, Tables S2, S3) showed that phylogenetic structure and alpha diversity of bacterioplankton did not exhibit any time or treatment-specific pattern. This result suggests that bacterioplankton was not significantly associated to the microbial community restructuring observed in clownfish and anemones from physical (PI) and remote (RI) interaction groups from T1 to T3, as detailed below.
Thetayc dissimilarity analyses. Analysis of the anemones’ epithelial microbiota (Figure 2a) shows that prior to contact with clownfish, dissimilarity between test and control groups was minimum (0.12 ± 0.01) (T0: after three weeks of acclimation). At T1, 15 minutes after the clownfish test individuals were transferred from the fish control tank system into their respective two-tanks systems for remote interaction (RI) (i.e. six biological replicates of one anemone tank connected with one fish tank), and after the first 15 min of physical contact between physical interaction (PI) clownfish individuals with their respective anemone (i.e. six biological replicates of physical interaction), dissimilarity between test (physical and remote interaction) and control anemones was significantly higher (0.39 ± 0.05) (Student T test, p < 0.001 ***) relatively to that of T0. Then, the dissimilarity between test and control anemones remained high during the interaction period (0.43 ± 0.03) (T1 to T3). From T4 (one week after PI / RI clownfish individuals were retrieved), to T5 (two weeks after PI / RI clownfish individuals were retrieved), dissimilarity between test and control anemones was significantly lower (0.36 ± 0.03) compared to that of the contact period (Student T test, p < 0.009 ***).
Regarding clownfish skin microbiota (Fig. 2b), the same pattern as that of the anemones occurred from T0 to T3: dissimilarity at T0 between PI / RI clownfish test and control groups was minimum (0.029 ± 0.005) prior to fish contact with their respective anemone. At T1, as soon as PI and RI clownfish test individuals were placed into their respective tank systems, dissimilarity between PI / RI test and control clownfish was significantly higher (0.61 ± 0.08, Student T test, p < 0.001 ***) relatively to that of T0. Then, the dissimilarity between PI / RI test and control clownfish remained high (0.46 ± 0.04) during the interaction period (T1 to T3). From T4 (one week after PI / RI clownfish individuals were retrieved and moved back to the control clownfish water system), to T5 (two weeks after), dissimilarity between PI / RI test and control clownfish groups remained stable (0.47 ± 0.03) and not significantly lower compared to that of the contact period (T1-T2-T3) (Student T test, p = 0.6).
Finally, regarding dissimilarity between fish and anemone microbiota (Fig. 2c), it was similar at T0 (0.41 ± 0.02) in all groups. Then, at T1, T2 and T3, the dissimilarity dropped to 0.16 ± 0.04 in PI and 0.11 ± 0.02 in RI test groups, significantly below the stable dissimilarity (0.56 ± 0.04) observed between fish and anemone control groups (Student T test, p < 0.001***). At T4 and T5, after fish-anemone pairs' separation, the dissimilarity values increased in PI (0.44 ± 0.03) and RI (0.30 ± 0.04) contact groups, but without reaching that of their respective controls (0.50 ± 0.04). This partial recovery is most likely explained by the stability of the dissimilarity between PI / RI test and control clownfish groups from T3 to T5 (Fig. 2b). In addition, the dissimilarity between the host microbiota and the bacterioplankton (Fig. S2) was never significantly different between the test and control groups.
Differential abundance analysis in clownfish and anemone skin microbiota between contact and control groups. De novo ASV abundances in clownfish and anemone were monitored during the whole experiment (from T0 to T5) using differential abundance analysis (DESeq2) to identify bacterial taxa that were mostly associated to fish-anemone epithelial microbiota convergence. To reach this goal, de novo ASV abundances of clownfish and anemone were combined, PI and RI groups were combined as an interaction group, as well as clownfish and anemone controls were combined as a control group. ASVs with log2-normalized fold-change over 1 and Bonferroni corrected p-value < 0.05 were kept (Table S1). At T0, there were only 5 differentially abundant taxa between interaction and control fish-anemone pairs. At T1, after the first 15 minutes of clownfish-anemone interaction, differentially abundant taxa increased to 10 ASVs. At T2, after one week of interaction, differentially abundant taxa doubled to reach 21 ASVs. At T3, after two weeks of interaction, differentially abundant taxa peaked at 30 ASVs. At T4, one week after separation of interaction fish-anemone pairs, the number of differentially abundant ASVs remained at 30. At T5, two weeks after separation of interaction fish-anemone pairs, the number of differentially abundant taxa dropped to 17 ASVs. From T2 to the end of the experiment (T5), three ASVs (2, 49, 177) matching to Cellulophaga tyrosinoxydans strain EM41 (95% identity, 100% coverage, 1.31E-119 to 2.81E-121 e-values) exhibited an interesting dynamic: they peaked at T3, with the three highest Bonferroni corrected p-values, with a fold change ranging from 9 to 12, then decreased gradually after separation of fish-anemone pairs: fold change ranging from 6 to 11 at T4, and from 3 to 8 at T5, relatively to fish-anemone control group. Therefore, these three ASVs related to Cellulophaga tyrosinoxydans were further analyzed in terms of abundance dynamics in water, sea anemone and clownfish mucus.
Differential abundance analysis (DESeq2) on Cellulophaga sp. in water, sea anemone and clownfish. The monitoring of the three ASVs (2, 49, 177) related to C. tyrosinoxydans, which were differentially abundant from T2 to T5 (DSeq2, Fig. 3, Table 3) was decomposed in terms of host community (sea anemone, clownfish, water), experimental groups and time. At T0, Cellulophaga sp. counts were both low and variable across experimental groups and host communities. ASVs with log2-normalized fold-change over 1 and FDR corrected p-values < 0.0001 were kept (Table 2).
Tank system water: From T0 to T1, Cellulophaga sp. dropped in all experimental groups to become undetectable except for FC, where only ASV 2 was significantly higher to PI (8.9 fold change). At T2, Cellulophaga sp. was still undetectable in anemone control group, and dropped under the detection threshold in clownfish control group. On the contrary, Cellulophaga sp. counts increased for the three ASVs both in PI and RI tank water to become statistically higher than in AC and FC groups (8.6 to 13.6 fold changes). At T3, the three ASVs were still undetectable in both control group water, whereas peaking in both PI and RI groups (9.7 to 13.7 fold changes). From T4 to T5, one and two weeks after clownfish retrieving from PI and RI tank systems, the three ASVs counts decreased gradually (5.2 to 13.7 fold changes at T5) in both PI and RI tank system water.
Sea anemone epithelium. From T0 to T1, Cellulophaga sp. counts dropped under the detection threshold in the three experimental groups hosting anemones (AC, PI, RI). At T2, Cellulophaga sp. counts increased for ASVs 49 and 177 in PI and RI to become significantly higher than in control (7.7 and 9.5 fold changes). At T3, the counts of the three ASVs peaked in PI and RI groups, and were still significantly higher than in control (8 and 13 fold changes). From T4 to T5, one and two weeks after clownfish retrieving from PI and RI tank systems, Cellulophaga sp. counts decreased quickly: ASV 2 was no more differentially abundant, and ASVs 49 and 177 dropped from 5.8 to 10.2 fold changes at T4, and were no more significantly different from their control at T5.
Clownfish skin mucus.
At T0, the three Cellulophaga sp. related ASVs counts were low and comparable between the three clownfish groups, which had shared the same tank system water for the three weeks acclimation period. From T0 to T1, Cellulophaga sp. counts remained low and comparable between the three clownfish groups, despite the transfer of PI and RI individuals to their respective PI and RI tank systems hosting anemones. At T2, Cellulophaga sp. counts of ASVs 49 and 177 increased in PI and RI to become significantly higher than in control (7.7 to 9.4 fold changes). At T3, the counts of the three ASVs peaked in PI and RI groups, and were still significantly higher than in control (10.7 to 13.3 fold changes). At T4, one week after PI and RI clownfish were reintroduced into the fish control water system, the counts of the three ASVs remained high and significantly higher than in control fish (7.1 to 12.9 fold change), despite sharing the same tank system water. At T5, two week after PI and RI clownfish reintroduction into the fish control water system, the three ASVs counts remained high only in PI clownfish (5 to 10.6 fold changes), whereas ASVs 2 and 49 dropped drastically in RI clownfish, both of them being no more significantly higher than in control fish.