Ecological filtering and phylogeographic structuring of Psychrilyobacter within two closely related limpet species from the Southern Ocean

Purpose The ecological interdependence between macroorganisms and their microbial communities promotes stable associations over time, potentially leading to their evolutionary co‑diversification. The detection of intricate eco‑evolutionary interactions between animals and their microbiota is challenging, primarily due to complex bacte‑ rial communities related to poorly resolved host population structure. Strikingly, co‑diversification in invertebrates, characterized by generally less complex microbiota, remains largely unexplored. Here, we compared the bacterial communities associated with two distinct lineages of Nacella limpets, a dominant shallow water patellogastropod of the Southern Ocean shores with a well‑described population structure. Our goals were to elucidate the uniqueness of Nacella microbiota, resulting from an ecological filter that selectively favors certain bacterial taxa. Additionally, we aimed to depict the genetic structure of bacterial symbiont seeking evidence of co‑diversification with Nacella . Methods We sequence the V4‑V5 regions of the bacterial 16S rRNA gene in three distinct microenvironments associated with Nacella : rock substrate, radula, and whole intestine. These samples were collected from two popula‑ tions of Nacella deaurata and Nacella concinna , located in the West Antarctic Peninsula and Falkland/Malvinas Islands, respectively. Results We assessed ecological filtering patterns in the limpet microbiota, uncovering unique bacterial communities in both radulas and intestines, with specifically enriched bacterial taxa compared to the surrounding environment. By examining microdiversity patterns of core bacterial taxa, we revealed a deep phylogeographic structure of Psychrilyo-bacter in Nacella intestines. Conclusion We highlight the Southern Ocean limpets of the Nacella genus as a novel and promising model for stud‑ ying co‑diversification between marine mollusks and their resident microbiota.


Findings
The intimate interdependence between macroorganisms and their microbial partners promotes stable associations over evolutionary timescales, which could ultimately result in the simultaneous divergence of both host and microbial symbionts (Moran and Sloan 2015;Moeller et al. 2016).This co-divergence, also known as co-diversification, does not necessarily imply obligatory host-restricted symbionts.For instance, repeated horizontal acquisition of certain microbes, through an ecological filtering of the surrounding environment community shaped by host diet or physiology, can lead to host-restricted microbial lineages (Moran and Sloan 2015;Moran et al. 2019).Nevertheless, the high complexity of microbial communities and the poorly resolved host population structure generally challenge the detection of specific host-microbe associations, particularly in vertebrate models (Alberdi et al. 2022).
In contrast to vertebrates, invertebrates generally possess less complex microbiota, making them intriguing models to comprehend the importance of host evolution in structuring animals-microbes associations (Petersen and Osvatic 2018;O'Brien et al. 2019;Groussin et al. 2020).Among invertebrates, mollusks stand out as ecologically significant and diverse phyla, being Gastropoda the richest class regarding species diversity (Neu et al. 2019).Although descriptions of marine gastropod microbiota are accumulating, they mainly focus on commercial and/or northern hemisphere species, excluding most gastropod diversity (Neu et al. 2019;Maltseva et al. 2021b;Panova et al. 2022).Moreover, the contribution of host diversification in shaping gastropod microbiota assembly remains largely unexplored.This knowledge gap widens in the Southern Ocean, where the true limpet genus Nacella (Nacellidae) thrives as the most conspicuous patellogastropod along intertidal and subtidal shores.Ecologically, species within this genus are characterized as generalist grazers due to their diverse diet, which includes periphyton, macroalgae, and some invertebrates (Rosenfeld et al. 2018).Biogeographically, the genus Nacella contains two major and recently diverged clades, the South American and the Maritime Antarctic/sub-Antarctic Islands, associated with the intensification of the Antarctic Circumpolar Current, which constitutes an oceanographical barrier to species dispersion (González-Wevar et al. 2017).On the premise that host speciation through geographical isolation may trigger its symbionts' diversification (Groussin et al. 2020), we hypothesized that two geographically isolated Nacella species shared specific bacterial symbionts whose phylogeographic structure at fine genetic resolution mirrors the host divergence.In this proof-of-concept study leveraging 16S amplicon sequencing, we aimed to assess (i) the specificity of Nacella species microbiota, often referred to as ecological filtering (Mazel et al. 2018), (ii) the presence of a core microbiome between the South American (Nacella deaurata) and Antarctic (Nacella concinna) species, and (iii) the phylogeographic pattern at fine-genetic resolution (i.e., microdiversity) within core bacterial taxa associated with Nacella.
We sampled two closely related species, N. deaurata (n = 5) and N. concinna (n = 4), along with the rocks they were associated with, from two biogeographic provinces of the Southern Ocean separated by the Antarctic Polar Front: the Falkland/Malvinas Islands (FAL/MAL) and the West Antarctic Peninsula (WAP) (Additional Information S1: Table S1).The radula and the entire intestine were aseptically retrieved from limpet individuals, while rock surfaces were thoroughly scraped to collect the substrate material ("rock" hereafter).The genomic DNA of each microenvironment was extracted with the DNeasy PowerSoil Pro Kit (Qiagen, CA, USA).A touchdown PCR protocol was carried out to amplify the bacterial 16S V4-V5 region using the modified Bakt_341F/Bakt_805R primer pair (Klindworth et al. 2013) (Additional Information S2: Method S1).Subsequently, the PCR products were sequenced on the Illumina Miseq PE300 platform.Sequences were then processed using the Mothur pipeline (Schloss et al. 2009), OTUs were assembled based on 97%-and 99%-identity, and low-abundance OTUs were discarded according to Bokulich et al. (2013).Taxonomic classification of OTUs was performed with Mothur using the SILVA database (v.138) (Quast et al. 2013).Finally, 97%-and 99%-identity OTU tables were rarefied to their corresponding number sequences per sample, specifically 86,714 and 70,390, respectively.All rarefaction curves at both OTU resolution reached saturation, indicating robust coverage of sample bacterial diversity (Additional Information S3: Fig. S1).
By comparing alpha-diversity across microenvironments at 97%-identity OTU resolution, we found a clear reduction in both taxa and phylogenetic diversities from rocks to radulas and intestines (p-values < 0.05 and p-values < 0.001, respectively, Fig. 1a and b).A similar pattern was observed in beta-diversity, with a strong effect of the microenvironment in shaping bacterial community compositions (p-value < 0.001, Fig. 1c, Additional Information S1: Table S2).These findings suggest an ecological filtering of the surrounding microbial community, potentially mediated by host traits such as internal physicochemical properties, and leading to a more constrained and evolutionarily convergent set of bacterial taxa in the intestine (Moran and Sloan 2015).The lack of significant difference in both alpha-and beta-diversity between radulas and intestines (Fig. 1 and Additional Information S1: Table S2) suggests that the most selective filtering probably takes place right after radula scraping, as reported in the freshwater snail Pomacea canaliculate (Li et al. 2019).These findings mark the initial characterization of radula microbiota in marine snails.Moreover, detecting a microbial community in radulas, as diverse as in intestines (Fig. 1a and b), challenges the initial assumption of this microenvironment being inhospitable to bacteria.Additionally, the bacterial communities in radulas and intestines were more heterogeneous than in rocks (p-values < 0.001, Fig. 1d), highlighting the host individual influence in bacterial community composition (i.e.microbiota individualization), as observed in other invertebrate models (Schwob et al. 2020;Gafarova et al. 2023).Alpha-diversity and beta-diversity patterns across microenvironments were globally consistent between 97%-and 99%-identity OTU resolutions, with sole exception being bacterial community dissimilarities among microenvironments (Fig. 1 and Additional Information S3: Fig. S2).Specifically, dissimilarities between radula and intestine were significant at 99%-OTU but not at 97%-OTU (Fig. 1d and Additional Information S3: Fig. S2d).For subsequent comparisons between microenvironments and the phylogeographic analysis of core taxa, we focused on the 97%-identity resolution.
We further explored the community differentiation drivers within microenvironments using Linear discriminant analysis (LDA) effect size (LEfSe) (Yang 2020).We identified 146, 39, and 29 bacterial OTUs significantly enriched in rock, radula, and intestine, respectively (Fig. 2a).A larger number of OTUs sit at the intersection of radula with intestine and rock, suggesting radula as an ecotone.Contrarily, very few OTUs were found at the rock-intestine intersection, further supporting the higher specificity of the intestine microenvironment.The taxonomic assignment of enriched OTUs was aggregated at the genus level, with each microenvironment characterized by its own set of taxonomic groups.Specifically, the rock bacterial community was enriched in Granulosicoccus, Ilumatobacter, Leucothrix, Litoreibacter, Yoonia-Loktanella, the radula in unclassified Bacilli and Mycoplasmataceae genera (Fig. 2b), and the intestine in Psychrilyobacter and Psychromonas (Fig. 2b and Additional Information S1: Table S3).These two latter genera have been consistently found as key components of the gut microbiota of snail taxa like Littorina, Rubyspira, and Haliotis (Aronson et al. 2017;Gobet et al. 2018;Maltseva et al. 2021b), implying their strong adaptation to gut conditions and their significance in the gastropods' ecology.To evaluate the conservativeness level in Nacella-associated bacterial communities and detect bacterial taxa probably important for Nacella ecology, we established a 100% threshold core microbiota for each microenvironment.Of the 130 identified core OTUs, 64 were shared across the three microenvironments.Conversely, 19, 1, and 10 were exclusively observed within rock, radula, and intestine, respectively (Fig. 2c), pointing to lower conservatism of radula and intestine bacterial communities, likely due to host individual-related factors (Pisaniello et al. 2023).Notably, five dominant core taxa (i.e., OTU1-3, OTU8, and OTU10) accounted for 36% of the dataset abundance, including Psychrilyobacter and Psychromonas (γ-Proteobacteria), Sulfitobacter and Ahrensia (α-Proteobacteria), and Bastopirellula (Planctomycetes) (Fig. 2c).Relative abundance patterns of these core OTUs substantially varied according to the host and the microenvironment (Additional Information S1: Table S4), suggesting that setting an arbitrary abundance cut-off in core definition may inadvertently excluded ecologically relevant bacteria in the context of single temporal sampling (Neu et al. 2021).
To explore the microdiversity structuration of these OTUs across microenvironments and Nacella species, we used the Minimum entropy decomposition algorithm (Eren et al. 2015), following the procedure described in a previous work (Schwob et al. 2021).Remarkably, the OTU3 assigned to Psychrilyobacter was the sole OTU displaying unique oligotypes within each Nacella species intestine (Additional Information S1: Table S5).We then reconstructed OTU3 oligotype networks in each microenvironment using the Median Joining method in the PopART software (Leigh and Bryant 2015), estimating genetic and phylogeographic differentiation indices (F ST and Φ ST , respectively) in Arlequin software (Schneider et al. 2000).Intriguingly, unlike rock and radula, OTU3 microdiversity in intestine exhibited strong phylogeographic structure, forming two distinct clusters specific to each Nacella species (Fig. 3).The disparity between differentiation indices in the intestine (F ST < Φ ST , Fig. 3) suggests that the structuration among the oligotypes associated with N. deaurata and N. concinna primarily stemmed from an evolutionary divergence between the oligotypes within each host-species (Maturana et al. 2020).The phylogenetic analysis of these oligotypes using Bayesian inference further supported the hypothesis that these fine-scale genetic clusters represent at least two strongly diverged evolutionary units with different host specificity (Additional Information S2: Method S2 and Additional Information S3: Fig. S2) (Martiny et al. 2023).
Additionally, the observation of identical Nacella species-specific oligotypes in both intestines and rocks indicates that hosts might acquire Psychrilyobacter horizontally by selectively associating with particular strains from their environment.The fine-scale phylogeographic patterns in Nacella intestines suggest host-specific regulatory mechanisms, either passive (e.g., internal condition filtering) or active (e.g., immune system), controlling the association with Psychrilyobacter strains (Cicala et al. 2018;Maltseva et al. 2021b).These host-specific regulatory mechanisms, understudied in Nacella, are likely to vary among species.Previous studies have consistently demonstrated that closely related gastropod species exhibited variations in their diets, proteomics, and metabolomics, thereby triggering the association with divergent microbiomes (Maltseva et al. 2021a(Maltseva et al. , 2021b)).
A recent genomic study of an abalone-associated strain of Psychrilyobacter provided evidence of functions potentially crucial for the host's physiology, such as synthesizing various vitamins or maintaining an acidic gut  (unc., unclassified).c Phylogenetic tree of the 130 OTUs shared across either rock, substrate, or intestine samples.The colors of the tree tips correspond to the core OTUs assignment at the class level.The adjacent heatmaps indicate the membership of the core OTUs to each microenvironment and its relative abundance within the core microbiota dataset pH (Liu et al. 2023).The microdiversity structuration we reported, if indicative of Nacella's metabolism reliance on Psychrilyobacter, may result from stable associations over time and space, potentially culminating in the specialization of Psychrilyobacter strains to their respective hosts and their diversification mirroring the phylogenetic divergence of N. concinna and N. deaurata.Thus, the diversification of Psychrilyobacter strains is likely to be coupled with functional divergence, which remains to be fully explored through genome-resolved metagenomics, transcriptomics, or metabolomics approaches to fully understand their ecological role within Nacella.
In conclusion, our study marks a promising starting point in exploring the intricate evolutionary relationships between Nacella and its microbiome, which raises intriguing questions regarding whether the divergence of metaorganisms in the Southern Ocean may have triggered their microbiome co-diversification.

Fig. 1
Fig.1 Bacterial alpha-and beta-diversity at 97%-identity OTU resolution across Nacella spp.microenvironments.Colors are assigned to the microenvironment.Triangles and circles refer to N. concinna and N. deaurata, respectively.Squares represent interspecies pairwise dissimilarity of beta-diversity.Chao1 (a) and phylogenetic diversity (b) comparisons across microenvironments, independent of Nacella species.c Non-metric multidimensional scaling (NMDS) profiling revealing the bacterial community composition variations among Nacella microenvironments and species.d Comparisons of bacterial community dissimilarities among microenvironments.The effect of the microenvironment factor on sample grouping in the ordination and tested by PERMANOVA, as well as the p-values associated with the Wilcoxon test pairwise comparisons among microenvironments, are provided

Fig. 2
Fig. 2 Shared and specific bacterial taxa across the Nacella microenvironments.a Ternary plot showing the significantly enriched OTUs in each microenvironment.Size is scaled on the LDA scores from the LEfSE analysis.b Compositional barplot revealing the genus taxonomic assignation of the enriched OTUs (unc., unclassified).c Phylogenetic tree of the 130 OTUs shared across either rock, substrate, or intestine samples.The colors of the tree tips correspond to the core OTUs assignment at the class level.The adjacent heatmaps indicate the membership of the core OTUs to each microenvironment and its relative abundance within the core microbiota dataset at 99%-identity OTU resolution across Nacella spp.microenvironments.Colors are assigned to the microenvironment.Triangles and circles refer to N. concinna and N. deaurata, respectively.Squares represent interspecies pairwise dissimilarity of beta-diversity.(a) Chao1 and (b) phylogenetic diversity comparisons across microenvironments, independent of Nacella species.(c) Non-Metric Multidimensional Scaling (NMDS) profiling revealing the bacterial community composition variations among Nacella microenvironments and species.(d) Comparisons of bacterial community dissimilarities among microenvironments.The effect of the microenvironment factor on sample grouping in the ordination and tested by PERMANOVA, as well as the p-values associated with the Wilcoxon test pairwise comparisons among microenvironments, are provided.Fig. S3.Bayesian inference of the phylogenetic relationships among the oligotypes conforming OTU3 assigned to Psychrilyobacter genus.The size of circles in the adjacent bubble plot is scaled on oligotype relative abundance (%).Colors indicate the Nacella species (orange; N. deaurata, green; N. concinna).The support of the tree basal node is provided.