Abstract Compositional variation of the gut microbiota across host allopatric populations can reflect both adaptation and stochasticity since the time of separation. Major factors shaping this variation include the host phylogeographic and demographic history, the microbiota inheritance, environmental inputs and dispersal of bacteria. Here we explored the impact of these factors in driving gut community diversity in seven allopatric populations of the omnivorous lizard Podarcis lilfordi from the Menorcan coastal islets, all descending from an ancestral mainland population. Using 16S rRNA Illumina sequencing, we showed that ‘islet’ and ‘age’ (time since islet separation from mainland) were the only significant variables in microbial community clustering, suggesting a partial islet-restricted diversification following these lizards phylogeography. Despite a significant variation, islets/populations were characterized by a remarkably low bacterial uniqueness (2.4% of total OTUs) and a minor differential enrichment of taxa, indicating a negligible impact of local inputs and important host common constraints. Overall, the extant pattern of similarity/dissimilarity among islets is compatible with partial retention of the ancestral mainland microbial pool, with differences among islets potentially explained by a differential loss of bacteria following population fragmentation and bottlenecks (i.e. ecological drift). While more quantitative data are needed to validate this hypothesis, this study unveils the importance of considering both neutral and niche-driven processes in driving contemporary patterns of gut metacommunity diversity. island, drift, dispersal, inheritance, population, allopatry INTRODUCTION A multitude of factors shape the gut microbiota composition and challenge its stability over ecological and evolutionary times. Most critical among them are the host constraints (e.g. gut morphology, physiology and immune system) (Rawls et al.2006; Ley et al.2008; Benson et al.2010; Moeller et al.2016) and environmental inputs (e.g. diet and parasites) (Muegge et al.2011; David et al.2014; Delsuc et al.2014; Aivelo and Norberg 2016; Kohl et al.2017). Although diet and host genotype have been consistently found to be the most influential determinants, large part of animal intraspecific microbial variation remains unexplained by adaptive variables. Hence, the role of neutral processes (drift and dispersal) is gaining more attention (Costello et al.2012; Lankau, Hong and Mackie 2012; Nemergut et al.2013; Zeng et al.2015; Burns et al.2016; Martinson, Douglas and Jaenike 2017; Shoemaker, Locey and Lennon 2017). Population-level studies of gut microbiotas dynamics might represent the right scale to explore the influence of both neutral and selective forces on metacommunity structure. While substantial interindividual microbial variation might exist, populations do typically carry discrete gut metacommunities (Yatsunenko et al.2012; Song et al.2013), resulting from a restricted, although not exclusive circulation of gut bacteria favored by parental transmission (Funkhouser and Bordenstein 2013), preferential dispersal by shared ecological niche (e.g. common diet), reproductive strategy (including parental care) and social behavior (Nemergut et al.2013; Song et al.2013; Ren et al.2016). How gut metacommunities respond to restrictions in ‘free’ circulation among host members, for instance in the event of population split, remains to date poorly understood/investigated in natural systems (Lankau, Hong and Mackie 2012), owing to the concomitant effect of adaptive variables (ecological traits). Allopatric but ecologically similar populations found in islands might provide excellent natural systems to explore the impact of geographic barriers and host demographic history on gut microbiota diversity. In the event of fragmentation of an original insular population following isolation of some specimens into novel islets, distinct factors can act over time to shape the microbiota compositional structure (for simplicity, hereafter use of microbial/microbiota will refer to bacteria only). These include the conservation of taxa present in the original pool by stable transmission through host generations (i.e. microbiota inheritance) (Funkhouser and Bordenstein 2013; Zeng et al.2015), selective and stochastic changes in taxa relative abundances (the latter known as ecological drift) (Chase and Myers 2011; Lankau, Hong and Mackie 2012), acquisition of novel taxa from local pools (i.e. local exposure) and dispersal of taxa among populations (Nemergut et al.2013). Of these, ecological drift during and following population isolation and local exposure should drive microbiota divergence among populations (Lankau, Hong and Mackie 2012), whereas microbiota inheritance and dispersal of common taxa across islands should maintain/increase microbial similarity among populations (Zeng et al.2015). Disentangling the contribution of each of these processes on extant patterns of similarity/dissimilarity of gut microbiota would provide a better understanding of the compositional stability and dynamics of these microbial communities over spatial and temporal changes. We explored this scenario using allopatric populations of the lizard species Podarcis lilfordi inhabiting coastal islets around Menorca, one of the Balearic Islands (Fig. 1). The species originated during the re-flooding of the Mediterranean at the end of the Messinian salinity crisis, about 5.33 Mya (Krijgsman et al.1999), with a first cladogenic event (2.88–2.66 Ma) leading to the separation between P. lilfordi (currently present in coastal islets of Mallorca and Menorca, and Cabrera) and P. pityusensis (Ibiza, Formentera and adjacent coastal islets) (Brown et al.2008; Terrasa et al.2009). To date, endemic populations of P. lilfordi inhabit 16 of the 19 Menorca islets (Pérez-Mellado 1989; Terrasa et al.2009) (Fig. 1). According to mtDNA data, these Menorca populations were largely panmictic during the Pleistocene (Terrasa et al.2009); the following sea level rise during the past 10 000 years (Holocene) (Vecchi et al.2016) led to coastal fragmentation and progressive isolation of islet populations up to date (Pretus, Marques and Pérez-Mellado 2004; Brown et al.2008; Vecchi et al.2016). Archaeological data indicate that the mainland stock became extinct ∼2000 years ago subsequent to the introduction of competitive species and foreign predators (Alcover, Moyá-Solà and Pons-Moyá 1981; Sanders 1984; Pretus, Marques and Pérez-Mellado 2004), while islets provided valuable refugia. Extant populations have been geographically and reproductively isolated for the past 2000 years, although time since separation most likely ranges between 4000 and 9000 years, according to bathymetry with respect to mainland and Holocene sea level curves (Vecchi et al.2016). While morphologically diverse, suggesting incipient speciation (Salvador 2009), mtDNA data show low genetic differentiation among P. lilfordi populations, in accordance with their recent segregation (Terrasa et al.2009). The presence of additional vertebrates (birds, rabbits and humans) is restricted to few islets (Pérez-Mellado 1989; Salvador 2009). Several studies indicate that a dietary shift occurred with variable intensity in the different populations, from an ancestral insectivorous niche to a current omnivorous niche, including the uptake of a large proportion of arthropods and vegetable matter, with seasonal variations (Saez and Traveset 1995; Traveset and Riera 2005; Perez-Cembranos, Leon and Perez-Mellado 2016). Figure 1. View largeDownload slide Map of Menorca and sampled islets (circled). Population size according to Pérez-Mellado et al. (2008) is provided below each islet (no estimates available for Ravl, although expected population size is likely below 1000 individuals, personal observations). Bars below each islet correspond to a 50-mt scale. Islet abbreviations: AdGr (Addaia Gran), Aire (Aire), Colm (Colom), PrCv (Porros de Cavalleria), Ravl (Ravells), Rei (Rei) and Tosq (Tosqueta). Figure 1. View largeDownload slide Map of Menorca and sampled islets (circled). Population size according to Pérez-Mellado et al. (2008) is provided below each islet (no estimates available for Ravl, although expected population size is likely below 1000 individuals, personal observations). Bars below each islet correspond to a 50-mt scale. Islet abbreviations: AdGr (Addaia Gran), Aire (Aire), Colm (Colom), PrCv (Porros de Cavalleria), Ravl (Ravells), Rei (Rei) and Tosq (Tosqueta). Overall, two aspects make these populations a rare system for studying natural processes shaping the vertebrate gut microbiota: their genetic and trophic similarity, which allowed us to detract the differential effect of these two major factors affecting the microbiota, and their distribution into several allopatric populations descending from a single ancestral population. Here we aimed to address some key questions: Was the geographic isolation and subsequent exposure to distinct local pools of bacteria sufficient to drive significant divergence among enteric communities in these populations? Do their gut microbial communities still carry an identifiable signature of their past history, i.e. a core component that has been putatively inherited from the ancestral mainland population? What is the actual impact of ecological drift and dispersal in explaining the current microbiota variation? To explore these questions, we sampled whole intestines from P. lilfordi specimens from seven Menorca islets and performed comparative analyses of their gut microbial community structure. We specifically partitioned their gut microbiota into shared and variable components and investigated forces responsible for shaping contemporary patterns of similarity/dissimilarity. MATERIALS AND METHODS Sample collection and morphometric body measures A total of 72 individuals from 7 of the 16 populated islets surrounding the major Menorca island (3 to 18 individuals per population/islet) were sampled in summer 2015, between 27 July and 1 August (Fig. 1 and Table S1, Supporting Information). Summer represents the season when trophic resources are most limited for these lizard populations. Sampled islets were chosen to vary along a gradient of bathymetry (from deep to shallow channels between islets and Menorca) and distance to the coastline (together providing a proxy for time since population isolation), and for their different exposure to the open sea (including both islets found in open sea and within protected bay areas). Specimens were collected using food traps containing fruits or fruit juices. For each specimen, the following morphometric measurements were taken: number of ventral-row scales (VS), weight (W), snout-vent length (SVL), head width (HW) and head length (HL). Of the 72 specimens, 33 adult males were chosen for microbiota sequencing (five for population, except for Rei, for which three specimens were sampled), while the remaining specimens were immediately released (Table S1, Supporting Information). Adulthood stage was assigned when SVL ≥ 51 mm for females and SVL ≥ 55 mm for males, according to the mean values for the smallest subspecies as described in Salvador (2009). Sex was assigned primarily according to number of VS (<27 for males) and validated by body size and femoral pores (males are larger than females and show pores with visible lipophilic compounds (Salvador 2009) (Table S1, Supporting Information). Statistical correlations were assessed through the cor.test function in the R Stats package version 3.2.2. Specimen sampling and manipulation were carried out in accordance with the ethics guidelines and recommendations of the Species Protection Service (Department of Agriculture, Environment and Territory, Government of the Balearic Islands), under permit No. CEP 31/2015 given to LB and JLP. Microbiota DNA extractions and 16S rRNA Illumina sequencing Sampled specimens for the microbiota study were kept in individual plastic containers for up to 6 h, and then transferred to –20°C until death. The full digestive tracts were dissected with sterile instruments and preserved in RNAlater. Samples were immediately stored at 4°C for 24 h and then kept at –20°C until further processing. Sex was then carefully reassessed, revealing the inclusion of three adult females in the final sample set (Table S1, Supporting Information). Under a sterile hood, the entire digestive tracts were wiped from RNAlater using a Kimwipes filter paper (Sigma-Aldrich) and lengths measured with a Vernier Caliper (Table S1, Supporting Information). Stomachs and intestines were separated according to their different tissue characteristics (stomachs had visibly ticker walls). Stomach contents were removed for examination under a stereoscope. The intact intestines were transferred to tubes with a 800-μL lysis buffer (see Baldo et al.2015), chopped using sterile scissors for facilitating the lysis process, vigorously shacked and incubated at 70°C for 30 min with gentle agitation. After addition of 0.3 g of 0.1-mm zirconia-silica beads, tube contents were homogenized with a Precellys Evolution bead beating instrument (Bertin Technologies) at 5500 rpm for 2 × 45 s, and centrifuged at 4°C for 5 min at 16 000 × g. Supernatant was stored in a separate tube. Additional 300 μL of lysis buffer was added to the pellet, and samples were homogenized again in the Precellys instrument for 2 × 30 s. Stored supernatant was placed back into the tube with beads, 20 μL of proteinase K (to a final of 10 mg/mL) was added and samples incubated at 55°C for 1 h at 450 rpm. Extracted DNA was purified according to Baldo et al. (2015). Briefly, residues were eliminated with ammonium acetate and DNA precipitated with isopropanol and cleaned with standard ethanol washes followed by Qiagen columns (DNeasy Blood and Tissue kit). The region V3-V4 (∼460 bp) of 16S rRNA was amplified using the following non-barcoded primers including TruthSeq adapters: S-D-Bact-0341-b-S-17 (5΄-CTACACGACGCTCTTCCGATCT CCTACGGGNGGCWGCAG-3΄) and S-D-Bact-0785-a-A-21 (5΄-CAGACGTGT GCTCTTCCGATCT GACTACHVGGGTATCTAATCC -3΄). Amplicons from individual specimens were barcoded and pooled according to Baldo et al. (2017). The final library was sequenced on Illumina MiSeq v3 (600-cycle cartridge, 300 paired-end reads) at the Centre for Genomic Regulation in Barcelona (Spain). Sequence analyses After demultiplexing, paired-end reads were merged and quality-filtered using the scripts make.contigs and screen.seqs (maxambig = 0, maxhomop = 8, maxlength = 490) implemented in mothur v.1.25.0 (Schloss et al.2009). Chimeras were identified using uchime_ref in vsearch v1.4.4 (VSEARCH GitHub repository) against the Chimera Slayer reference database (‘gold’) in the Broad Microbiome Utilities and discarded. Sequences were input into the QIIME pipeline (Caporaso et al.2010), using Macqiime 1.9.1., and clustered into operational taxonomic units (OTUs) at 97% similarity through the command pick_open_reference_otus.py. Taxonomy was assigned with the UCLUST method against the Greengenes gg_13_8_otus with RDP classifier (80% confidence). Mitochondrial and chloroplast-derived OTUs were discarded, retaining only bacterial OTUs with ≥10 total sequences. Sequences were aligned with PyNAST default parameters (Caporaso et al.2010) using the reference database ‘core_set_aligned.fasta’ from Greengenes website and default lanemask. Alignments were further trimmed to ensure overlapping ends, and a phylogenetic tree was built with FastTree (Price, Dehal and Arkin 2009). Alpha diversity measures (phylogenetic diversity (PD) whole tree and Shannon Index) were estimated on multiple rarefied biom tables to an even depth of 38 160 reads (10 iterations). Statistical differences in alpha diversity across islets were calculated with a non-parametric two-sample t-test (999 permutations). Beta diversity and PERMANOVA analyses Unweighted Unifrac and binary Jaccard distance matrices were built on the rarefied OTU table (38 160 sequences). Statistical correlations between microbiota distances and categorical variables (‘islet’, ‘age’, presence of ‘rats’ and stable colonies of ‘seagulls’, and ‘exposure’ to open sea, see Table S1, Supporting Information) were assessed by permutational multivariate analysis of variance (PERMANOVA, Anderson 2001; Anderson, Gorley and Clarke 2008) for PRIMER v7 (Clarke and Gorley 2015). All factors were considered fixed save for ‘islet’, which was defined as a random factor and included in all designs as nested within the levels of fixed factors. Categorical ‘age’ provides an approximation of each islet isolation time and was largely assigned according to bathymetry (as depth of the shallowest channel between islets and main island) choosing a threshold of ≥9 m to discriminate between old and more recently separated islets, in accordance with the speed of sea level rising, which slowed down after 7500 yr BP (Vecchi et al.2016) (Table S1, Supporting Information). Bathymetry and distance from mainland were also combined to define different age categories, grouping Aire and PrCv (old) against all remaining islets (recent) (see the Results section). To summarize microbial metacommunity diversity by islet, OTU counts, with and without previous rarefaction, were binned by islet and the data were used to calculate ‘uniqueness’ by islet (i.e. taxa and OTUs occurring only in one islet to the exclusion of all other islets). To test whether microbiota distances correlated with the geographic distances among islets, the original OTU table was filtered for OTUs above 50 total counts, reads were binned by islet and rarefied to the minimum depth per islet. Jaccard and unweighted Unifrac distance matrices were built on this OTU table and tested against geographic (i.e. linear) distances across islets with a Mantel test. Islet distances from mainland were instead approximated by the ‘age’ category and tested with PERMANOVA. Taxa and OTUs enrichment among islets To statistically compare differences in taxa and OTU relative abundances among islets, we used the script group_significance.py in QIIME with the non-parametric Kruskal–Wallis test. In all cases, a false discovery rate (FDR) correction was applied to assess significant differences. Differentially enriched OTUs (i.e. indicator OTUs) were estimated on the original OTU table rarefied to 38 160 reads. For differentially enriched taxa (i.e. indicator taxa), OTU counts were binned to each taxon level. Results were filtered for OTU/taxa with total sum of means across islets greater than 10. Core analyses We used the non-rarefied OTU table to estimate core OTUs and taxa (from phylum to species, estimated by binning OTU counts at each taxon level, as above) at different depths (shared by 100%, 90% and 80% of total specimens), retaining those taxa and OTUs with representation in all islets. Core OTUs were searched with batch BLASTN against the local nt database (cut-off e-value = 30), retrieving the first two best hits. Outputs were parsed using in-house perl scripts. Metadata associated with each BLASTN hit was extracted from Genbank records with Biopython scripts. The 134 core OTUs present in at least 90% of the specimens were used to explore a potential signature of within-population diversification. For this purpose, all reads associated with core OTUs were pooled and reclustered at 98% and 99% of similarity (now referred as to 98% and 99% OTUs or more generally as to sub-OTUs), picking the most abundant as representative and retaining OTUs with at least 50 reads total and presence in at least two specimens. A binary Jaccard distance matrix was built on each sub-OTU abundance table, rarefied to the minimum depth coverage across samples (13 600 and 11 000 for 98% and 99% OTUs, respectively), and significance of clustering was tested with PERMANOVA, as for 97% OTUs (see above). The 23 core OTUs that occurred in all specimens (representing the 100% core) were analyzed individually for their sub-OTU contribution to the islet-based structuring. Specifically, reads from each OTU were reclustered at 99% similarity and significance tested with permutational analysis of variance using adonis (999 permutations) with the script compare_categories.py in QIIME. RESULTS Morphometrics and stomach content The large majority of the captured specimens were adult males, putatively because they were less shy in approaching food sources, although a skewed sex ratio cannot be ruled out (Table S1, Supporting Information). Males were on average significantly larger than females (weight, W was 8.2 g ± 2.17 for males, 5.4 g ± 1.55 for females) and had longer SVL (65.2 mm ± 4.20 for males, 58.6 mm ± 3.68 for females) and larger HW/SVL and HL/SVL values (Table S1, Supporting Information). Considering males only, Addaia Gran (AdGr) hosted the smallest lizards, which showed among the largest HL and HW/SVL ratios, while Aire hosted the biggest and longest specimens. These results are compatible with previous findings (Salvador 2009). Overall, SVL positively correlated with W (Pearson's r = 0.677, P < 0.0001), HW (r = 0.6735581, P < 0.0001) and HL (r = 0.863, P < 0.0001), but not with intestine length (r = 0.132, P = 0.333). Stomach content analyses indicated widespread and abundant presence of arthropods, particularly insects (especially ants) in most specimens and all islets, compatible with a dominant role of insectivory during the dry season. Sequence summary We characterized the gut microbiota of 33 total specimens (30 males and 3 females), with three to five representatives per population (Table S1, Supporting Information). After quality filtering, removal of the three OTUs found in the blank control (water only) and OTUs with less than 10 reads total, we obtained a total of 2 014 276 high-quality reads and an average of 60 081 reads per specimen (ranging between 38 163 and 106 206), corresponding to 6380 OTUs. These were assigned to 17 phyla, 32 classes (6378 annotated OTUs), 58 orders (6369), 99 families (4975) and 147 genera (2321). Rarefaction curves based on the number of observed OTUs indicated that we captured most of the sample diversity, although we did not reach saturation at the chosen rarefied depth of 38 160 sequences (Fig. S1, Supporting Information). The taxonomic composition of the P. lilfordi gut microbiota is highly conserved Podarcis lilfordi showed the typical vertebrate phyla composition, with three phyla representing altogether 97% of the microbiota: Firmicutes (average 65% across islets), Bacteroidetes (18%) and Proteobacteria (14%) (Fig. 2). The proportion of Firmicutes was highly comparable across islets, while the representation of both Bacteroidetes and Proteobacteria varied to a greater extent. Firmicutes and Bacteroidetes were represented by six and five orders, respectively, while being dominated by a single order in all islets: Clostridiales (Firmicutes) and Bacteroidales (Bacteroidetes) (Fig. 2). Two families within both phyla showed the highest relative abundances: Lachnospiraceae and Ruminococcaceae within Firmicutes, and Porphyromonadaceae and Bacteroidaceae within Bacteroidetes (Fig. 2). Patterns of family relative abundance within islet were remarkably conserved, especially for Firmicutes. On the other hand, Proteobacteria were more diverse, including a larger number of represented orders (24) and families (33). Among the most represented taxa within this phylum were four orders (Desulfobrionales, Enterobacteriales, Campylobacteriales and Pseudomonadales) and three families (Helicobacteraceae, Desulfovibrionaceae and Enterobacteriaceae), present in all islets (Fig. 2). Proteobacteria varied widely in their quantitative profiles across islets at all taxonomic levels. Nonetheless, alpha diversity estimates for this phylum were not significantly different among islets (two-sample t-test P > 0.05, both Shannon and PD whole tree for all islet pairwise). Figure 2. View largeDownload slide Order and family-level taxonomic composition of the three major phyla, Firmicutes, Bacteroidetes and Proteobacteria, as proportions of OTUs. Values were estimated on the rarefied dataset and averaged by islet. Legends list only the 10 most abundant taxa. Figure 2. View largeDownload slide Order and family-level taxonomic composition of the three major phyla, Firmicutes, Bacteroidetes and Proteobacteria, as proportions of OTUs. Values were estimated on the rarefied dataset and averaged by islet. Legends list only the 10 most abundant taxa. As a whole, gut communities did not show significantly distinct taxonomic diversity and richness across islets, according to both Shannon and PD-whole tree (Fig. S2, Supplementary; P > 0.05 at both indexes for any islet pairwise), although a visual trend indicates Porros de Cavalleria (PrCv) as the islet hosting the largest diversity. Diversity and richness did not correlate with intestinal length (Pearson's r = –0.092, P = 0.49 for PD whole tree and –0.194, P = 0.15 for Shannon). Overall, the major taxonomic structure of the P. lilfordi gut microbiota was remarkably conserved across islets, with Proteobacteria content representing the most variable quantitative trait. Islet as the major determinant in the P. lilfordi microbiota clustering According to both unweighted Unifrac and Jaccard distances on 97% OTUs, ‘islet’ was the main significant factor accounting for 13.5% of the community variation, as measured by partition into components of variation (PERMANOVA, P = 0.0001, Table 1). There were no significant effects of age (but see below), exposure or presence of rats (P > 0.05). The presence of seagull colonies was barely significant, accounting for only 3.5% of microbiota variation. The islet Rei was grouped with the oldest islets according to bathymetry, but it presents several unique features, including its location within a natural harbor at short distance from the mainland and in a context of a largely anthropized environment (the islet hosted a paleochristian basilica in the 6th century and a hospital in the 19th century). When both bathymetry and distance from the mainland were used as grouping factors, thus excluding Rei as an ‘ancient’ islet and retaining only Aire and PrCv, ‘age’ became a significant factor (P = 0.0001, Table 1). The only interaction that could be tested (age by exposure) was not significant. Microbiota distances across islets did not correlate with geographic distances across islets (Mantel test, r = 0.00588, P = 0.980 for Jaccard and r = –0.02849, P = 0.906 for Unifrac), but perhaps indirectly with their distances from the mainland, here approximated by the ‘age’ category (see the Methods section). The ‘islet’ and ‘age’ effects are also illustrated by a principal coordinate analysis based on unweighted Unifrac distances, depicting the specimens clustering according to the islet origin (color-coded, Fig. 3). The two oldest islets, Aire and PrCv, clearly diverged from the remaining islets along PC2. Figure 3. View largeDownload slide Principal coordinate analysis according to unweighted Unifrac distances. Distances were calculated for OTUs at 97% similarity using the rarefied dataset (38 160 sequences). Dots correspond to single specimens. Specimens within each islet are connected to the median coordinates for each islet. Figure 3. View largeDownload slide Principal coordinate analysis according to unweighted Unifrac distances. Distances were calculated for OTUs at 97% similarity using the rarefied dataset (38 160 sequences). Dots correspond to single specimens. Specimens within each islet are connected to the median coordinates for each islet. Table 1. Results of PERMANOVA analyses, separately for the whole dataset and core OTUs after reclustering at 98 and 99% similarity. Whole dataset Core (90% of specimens) Unweighed Unifrac Jaccard 97% Jaccard sub-OTU 98% Jaccard sub-OTU 99% Model Term P-value Variation P-value Variation P-value Variation P-value Variation estimate estimate estimate estimate Islet Islet 0.0001 0.028 0.0001 0.035 0.0001 0.029 0.0001 0.031 Age + islet(A) Age N.S. [0.0001] 0.004 [0.011] N.S. [0.0001] 0.006 [0.012] N.S. [0.0001] 0.005 [0.012] N.S. 0.003 [0.008] Islet(A) 0.0001 0.024 [0.022] 0.0001 0.029 [0.028] 0.0001 0.024 [0.023] 0.0001 0.028 [0.027] Exposure + islet(E) Exposure N.S. 0.000 N.S. 0.000 N.S. 0.000 N.S. 0.000 Islet(E) 0.0001 0.028 0.0001 0.034 0.0001 0.029 0.0001 0.031 Age + exposure + islet(A, E) Age N.S. 0.002 [0.010] N.S. 0.005 [0.012] N.S. [0.0286] 0.002 [0.012] N.S. 0.001 [0.008] Exposure N.S. 0.003 [0.000] N.S. 0.003 [0.000] N.S. 0.003 [0.001] N.S. 0.002 [0.001] A*E N.S. 0.003 N.S. 0.002 N.S. 0.003 N.S. 0.003 Islet(A, E) 0.0001 0.022 0.0001 0.027 0.0001 [N.T.] 0.022 0.0001 0.026 Rats + islet(R) Rats N.S. 0.001 N.S. 0.001 N.S. 0.002 N.S. 0.002 Islet(R) 0.0001 0.027 0.0001 0.034 0.0001 0.028 0.0001 0.030 Seagulls + islet(S) Seagulls 0.0465 0.007 0.0388 0.008 0.0299 0.008 0.0340 0.007 Islet(S) 0.0001 0.024 0.0001 0.030 0.0001 0.025 0.0001 0.027 Residual 0.178 0.201 0.161 0.175 Whole dataset Core (90% of specimens) Unweighed Unifrac Jaccard 97% Jaccard sub-OTU 98% Jaccard sub-OTU 99% Model Term P-value Variation P-value Variation P-value Variation P-value Variation estimate estimate estimate estimate Islet Islet 0.0001 0.028 0.0001 0.035 0.0001 0.029 0.0001 0.031 Age + islet(A) Age N.S. [0.0001] 0.004 [0.011] N.S. [0.0001] 0.006 [0.012] N.S. [0.0001] 0.005 [0.012] N.S. 0.003 [0.008] Islet(A) 0.0001 0.024 [0.022] 0.0001 0.029 [0.028] 0.0001 0.024 [0.023] 0.0001 0.028 [0.027] Exposure + islet(E) Exposure N.S. 0.000 N.S. 0.000 N.S. 0.000 N.S. 0.000 Islet(E) 0.0001 0.028 0.0001 0.034 0.0001 0.029 0.0001 0.031 Age + exposure + islet(A, E) Age N.S. 0.002 [0.010] N.S. 0.005 [0.012] N.S. [0.0286] 0.002 [0.012] N.S. 0.001 [0.008] Exposure N.S. 0.003 [0.000] N.S. 0.003 [0.000] N.S. 0.003 [0.001] N.S. 0.002 [0.001] A*E N.S. 0.003 N.S. 0.002 N.S. 0.003 N.S. 0.003 Islet(A, E) 0.0001 0.022 0.0001 0.027 0.0001 [N.T.] 0.022 0.0001 0.026 Rats + islet(R) Rats N.S. 0.001 N.S. 0.001 N.S. 0.002 N.S. 0.002 Islet(R) 0.0001 0.027 0.0001 0.034 0.0001 0.028 0.0001 0.030 Seagulls + islet(S) Seagulls 0.0465 0.007 0.0388 0.008 0.0299 0.008 0.0340 0.007 Islet(S) 0.0001 0.024 0.0001 0.030 0.0001 0.025 0.0001 0.027 Residual 0.178 0.201 0.161 0.175 N.S.: not significant (P > 0.05). N.T.: no test performed. Rarefaction depth used: whole dataset (38 160 sequences); core sub-OTU 98% (13 600), sub-OTU 99% (11 000). In brackets are estimates for ‘age’ when excluding Rei from the ‘old’ category (see Table S1, Supporting Information). Residual are common to all models. View Large Despite a main significant islet effect, a large fraction of the P. lilfordi gut microbiota was common to all populations: according to Jaccard pairwise distances across islets, mean distance between any two islets was only 0.38, with the Aire and PrCv representing the closest pair (mean Jaccard = 0.28). Hereafter, we investigated both shared and variable microbiota components driving the above pattern. Low microbiota uniqueness and differential enrichment among islets To explore the variable component responsible for the observed clustering and distances, we considered both taxa/OTU that are islet-specific in terms of presence/absence (i.e. uniqueness) and those that are differentially enriched. Taxa and OTUs unique to one islet can either represent novel gut acquisitions from the local pool, contamination by allochthonous bacteria, or simply loss of these taxa in all other islets. Total number of OTUs per islet, without rarefaction, ranged between 3013 and 4816 for a comparable sampling effort (five individuals per islet, except for Rei, with three individuals sampled) (Table 2). A total of 1093 OTUs (17% of total OTUs) were shared among all islets (although not all specimens). Uniqueness by islet was remarkably low, accounting for only between 1.1% and 5.2% of an islet's OTU content (average 2.4%) and showing limited distribution (1–2 individuals). This suggests that we either failed to sample them in most specimens or, more likely, that they represent part of the allochthonous component. Of the 666 total unique OTUs detected (Table 2), none corresponded to unique taxa when checked against the unfiltered OTU table including all reads. In all cases, most unique OTUs occurred in Aire, with Tosqueta (Tosq) hosting the lowest uniqueness. Overall, at least 95% of each islet OTU content was shared between two or more islets. The same analyses performed on the rarefied OTU table returned highly comparable results. Table 2. OTU uniqueness and shared component. Number of OTUs Islet Totala Unique (%) Shared %b Tosq 3287 37 (1.13) 98.87 Rei 3013 61 (2.02) 97.98 Ravl 3370 71 (2.11) 97.89 Colm 4334 82 (1.89) 98.11 PrCv 4816 87 (1.81) 98.19 AdGr 4133 110 (2.66) 97.34 Aire 4219 218 (5.17) 94.83 Average 3881.7 95.1 (2.40) 97.60 Total Unique Sharedb Commonc 90% of specimens All islets 6380 666 5714 1093 134 Number of OTUs Islet Totala Unique (%) Shared %b Tosq 3287 37 (1.13) 98.87 Rei 3013 61 (2.02) 97.98 Ravl 3370 71 (2.11) 97.89 Colm 4334 82 (1.89) 98.11 PrCv 4816 87 (1.81) 98.19 AdGr 4133 110 (2.66) 97.34 Aire 4219 218 (5.17) 94.83 Average 3881.7 95.1 (2.40) 97.60 Total Unique Sharedb Commonc 90% of specimens All islets 6380 666 5714 1093 134 aof OTUs was estimated from a non-rarefied dataset, considering OTUs > 10 total sequences. bby two or more islets. cin all islets, but not all specimens. View Large As for presence/absence, differences in relative abundances of taxa/OTUs among islets were also minor. According to enrichment analyses, only a few taxa (Table S2, Supporting Information) and OTUs (total of 30, Table 3) showed a significant differential representation among islets (kw, FDR P < 0.05). Among indicator taxa (Table S2, Supporting Information), the phylum Synergistetes and its only recognized family Synergistaceae were a signature of Aire and PrCv, with only a few reads detected in AdGr and Colom (Colm). Aire had a particularly high representation of this family in terms of reads (376) and number of specimens (all). The islet of Rei showed a particular enrichment in Oceanospirillales and Halomonadaceae, while Ravells (Ravl) was enriched in the genus Turicibacter. Table 3. Differentially enriched OTUs among islets (Kruskal Wallis, FDR P < 0.05). Mean of reads Finest taxonomy BLASTN Hits OTU_IDa Aire AdGr Rei Tosq Ravl Colm PrCv Specificity Phylum Order Family Genus/species Accession number Identity (%) Hostb N.Ref.OTU895 516.8 0 0 0 0 0 0 Aire Firmicutes Clostridiales Peptococcaceae GQ502520 91.8 Termite 721569 26.8 0 0 0 0 0 0 Aire Firmicutes Clostridiales Ruminococcaceae CU926622 98.9 175037 22.2 0 0 0 0 0 0 Aire Firmicutes Clostridiales Christensenellaceae FJ362815 97.7 Human NC.Ref.OTU185722 20.8 0 0 0 0 0 0 Aire Firmicutes Clostridiales FJ535004 96.8 NC.Ref.OTU27689 20.6 0 0 0 0 0 0 Aire Bacteroidetes Bacteroidales Bacteroidaceae Bacteroides HQ745181 96.7 Human NC.Ref.OTU224538 18.2 0 0 0 0 0 0 Aire Firmicutes Clostridiales Peptococcaceae GQ502520 90.8 Termite NC.Ref.OTU27099 17.2 0 0 0 0 0 0 Aire Firmicutes Clostridiales LC028698 97.3 Buffalo NC.Ref.OTU212111 15 0 0 0 0 0 0 Aire Firmicutes Clostridiales Ruminococcaceae GQ294881 96.4 Mouse 566434 14 0 0 0 0 0 0 Aire Firmicutes Clostridiales FJ535004 97.1 NC.Ref.OTU43839 10.4 0 0 0 0 0 0 Aire Firmicutes Clostridiales Ruminococcaceae AB218341 92.3 Dugong NC.Ref.OTU59527 0 0 0 0 0 0 10.6 PrCv Firmicutes Clostridiales AF371819 96.8 Swine 688934 0 0 0 0 25.2 0 0 Ravl Proteobacteria Enterobacteriales Enterobacteriaceae Serratia KY124195 98.3 NC.Ref.OTU106014 0 0 0 0 22.8 0 0 Ravl Proteobacteria Enterobacteriales Enterobacteriaceae FR848151 98.1 Millipedes 182970 128.6 0 0 0 0 0 664.2 Aire, PrCv Verrucomicrobia Verrucomicrobiales Verrucomicrobiaceae Akkermansia LT629973 99.1 Python N.Ref.OTU217 46.4 0 0 0 0 0 66.2 Aire, PrCv Firmicutes Clostridiales KU506346 94.4 Pig N.Ref.OTU651 28 0 0 0 0 0 73.4 Aire, PrCv Bacteroidetes Bacteroidales [Barnesiellaceae] HQ792288 91.5 Human NC.Ref.OTU124139 108 0 0 0 0 0 6 Aire, PrCv Firmicutes Clostridiales AB606338 95.7 Mouse NC.Ref.OTU110653 67.6 0 0 0 0 0 0.4 Aire, PrCv Firmicutes Clostridiales JQ248100 97.3 Chicken NC.Ref.OTU229167 2.2 0 0 0 0 0 16.8 Aire, PrCv Firmicutes Clostridiales Ruminococcaceae HM123952 95.5 Human 3970651 1041.6 0 0 0 0 23 0 Aire, Colm Bacteroidetes Bacteroidales Porphyromonadaceae Parabacteroides JX457157 99.3 Cockroach NC.Ref.OTU12308 15.6 0 0 0 0 0.2 0 Aire, Colm Bacteroidetes Bacteroidales Porphyromonadaceae Parabacteroides JX457157 96.7 Cockroach NC.Ref.OTU231957 7.4 0 12 0 0 0 0 Aire, Rei Firmicutes Clostridiales Ruminococcaceae CU919465 94.4 NC.Ref.OTU25986 0 0 0 0 100.4 0 3.2 Ravl, PrCv Firmicutes Clostridiales AB606272 94.3 Mouse 347529 25 6.2 0 0 329.2 5.2 26.6 Several Firmicutes Turicibacterales Turicibacteraceae Turicibacter LT624967 99.1 727945 8.4 2.6 0 0 11 1.4 2.2 Several Firmicutes Clostridiales Eubacteriaceae Anaerofustis LT223622 99.3 Human NC.Ref.OTU117227 1.4 0 0 0 0.2 0 16.4 Several Firmicutes Clostridiales Lachnospiraceae Coprococcus EU779145 96.6 Sheep NC.Ref.OTU63406 7 42.2 0 0 0 0 54.2 Several Bacteroidetes Bacteroidales Bacteroidaceae Bacteroides GU105959 92.8 Human N.Ref.OTU255 10 0 0 0.2 0 14.2 7.2 Several Firmicutes Clostridiales [Mogibacteriaceae] Anaerovorax JQ084279 94.8 Mouse N.Ref.OTU41 274.4 134 3.7 0 0 0 86.2 Several Bacteroidetes Bacteroidales Porphyromonadaceae Parabacteroides HQ767608 94.3 Human N.Ref.OTU544 251.6 0.2 46.7 0 0 0 51.6 Several Bacteroidetes Bacteroidales Porphyromonadaceae P. distasonis HQ797049 97.4 Human Mean of reads Finest taxonomy BLASTN Hits OTU_IDa Aire AdGr Rei Tosq Ravl Colm PrCv Specificity Phylum Order Family Genus/species Accession number Identity (%) Hostb N.Ref.OTU895 516.8 0 0 0 0 0 0 Aire Firmicutes Clostridiales Peptococcaceae GQ502520 91.8 Termite 721569 26.8 0 0 0 0 0 0 Aire Firmicutes Clostridiales Ruminococcaceae CU926622 98.9 175037 22.2 0 0 0 0 0 0 Aire Firmicutes Clostridiales Christensenellaceae FJ362815 97.7 Human NC.Ref.OTU185722 20.8 0 0 0 0 0 0 Aire Firmicutes Clostridiales FJ535004 96.8 NC.Ref.OTU27689 20.6 0 0 0 0 0 0 Aire Bacteroidetes Bacteroidales Bacteroidaceae Bacteroides HQ745181 96.7 Human NC.Ref.OTU224538 18.2 0 0 0 0 0 0 Aire Firmicutes Clostridiales Peptococcaceae GQ502520 90.8 Termite NC.Ref.OTU27099 17.2 0 0 0 0 0 0 Aire Firmicutes Clostridiales LC028698 97.3 Buffalo NC.Ref.OTU212111 15 0 0 0 0 0 0 Aire Firmicutes Clostridiales Ruminococcaceae GQ294881 96.4 Mouse 566434 14 0 0 0 0 0 0 Aire Firmicutes Clostridiales FJ535004 97.1 NC.Ref.OTU43839 10.4 0 0 0 0 0 0 Aire Firmicutes Clostridiales Ruminococcaceae AB218341 92.3 Dugong NC.Ref.OTU59527 0 0 0 0 0 0 10.6 PrCv Firmicutes Clostridiales AF371819 96.8 Swine 688934 0 0 0 0 25.2 0 0 Ravl Proteobacteria Enterobacteriales Enterobacteriaceae Serratia KY124195 98.3 NC.Ref.OTU106014 0 0 0 0 22.8 0 0 Ravl Proteobacteria Enterobacteriales Enterobacteriaceae FR848151 98.1 Millipedes 182970 128.6 0 0 0 0 0 664.2 Aire, PrCv Verrucomicrobia Verrucomicrobiales Verrucomicrobiaceae Akkermansia LT629973 99.1 Python N.Ref.OTU217 46.4 0 0 0 0 0 66.2 Aire, PrCv Firmicutes Clostridiales KU506346 94.4 Pig N.Ref.OTU651 28 0 0 0 0 0 73.4 Aire, PrCv Bacteroidetes Bacteroidales [Barnesiellaceae] HQ792288 91.5 Human NC.Ref.OTU124139 108 0 0 0 0 0 6 Aire, PrCv Firmicutes Clostridiales AB606338 95.7 Mouse NC.Ref.OTU110653 67.6 0 0 0 0 0 0.4 Aire, PrCv Firmicutes Clostridiales JQ248100 97.3 Chicken NC.Ref.OTU229167 2.2 0 0 0 0 0 16.8 Aire, PrCv Firmicutes Clostridiales Ruminococcaceae HM123952 95.5 Human 3970651 1041.6 0 0 0 0 23 0 Aire, Colm Bacteroidetes Bacteroidales Porphyromonadaceae Parabacteroides JX457157 99.3 Cockroach NC.Ref.OTU12308 15.6 0 0 0 0 0.2 0 Aire, Colm Bacteroidetes Bacteroidales Porphyromonadaceae Parabacteroides JX457157 96.7 Cockroach NC.Ref.OTU231957 7.4 0 12 0 0 0 0 Aire, Rei Firmicutes Clostridiales Ruminococcaceae CU919465 94.4 NC.Ref.OTU25986 0 0 0 0 100.4 0 3.2 Ravl, PrCv Firmicutes Clostridiales AB606272 94.3 Mouse 347529 25 6.2 0 0 329.2 5.2 26.6 Several Firmicutes Turicibacterales Turicibacteraceae Turicibacter LT624967 99.1 727945 8.4 2.6 0 0 11 1.4 2.2 Several Firmicutes Clostridiales Eubacteriaceae Anaerofustis LT223622 99.3 Human NC.Ref.OTU117227 1.4 0 0 0 0.2 0 16.4 Several Firmicutes Clostridiales Lachnospiraceae Coprococcus EU779145 96.6 Sheep NC.Ref.OTU63406 7 42.2 0 0 0 0 54.2 Several Bacteroidetes Bacteroidales Bacteroidaceae Bacteroides GU105959 92.8 Human N.Ref.OTU255 10 0 0 0.2 0 14.2 7.2 Several Firmicutes Clostridiales [Mogibacteriaceae] Anaerovorax JQ084279 94.8 Mouse N.Ref.OTU41 274.4 134 3.7 0 0 0 86.2 Several Bacteroidetes Bacteroidales Porphyromonadaceae Parabacteroides HQ767608 94.3 Human N.Ref.OTU544 251.6 0.2 46.7 0 0 0 51.6 Several Bacteroidetes Bacteroidales Porphyromonadaceae P. distasonis HQ797049 97.4 Human Taxonomy in brackets is not officially recognized (under revision). aOTUs were filtered for total sum of means > 10. bPresence in gut or feces samples. Blanks correspond to putative environmental hits. View Large Among indicator OTUs (Table 3), the majority belonged to the order Clostridiales (Firmicutes). Except for six putative environmental OTUs, all others had best hits to gut/feces of animals, typically vertebrates. Of them, 13 were single islet-specific, 10 being uniquely found in Aire. One Aire-specific OTU (N.Ref.OTU895) was present in high abundance (516 mean reads) in most specimens and had best BLASTN hit with low identity (91.8%) to the termite gut. Six OTUs were specific to Aire and PrCv, three of them with large mean read abundance in both islets (Table 3). OTU-182970, corresponding to a member of the genus Akkermansia, was particularly abundant in PrCv and had the best BLASTN hit to the reticulated python gut. This is overall the variable component responsible of the islet-based compositional structure of the P. lilfordi microbiota. Hereafter, we focused on the large shared component, the core microbiota. Core composition of the P. lilfordi gut microbiota We considered core taxa those with presence in at least 90% of specimens and all islets (calculated on the non-rarefied dataset). The core comprised 7 phyla, 12 classes, 13 orders, 25 families, 27 genera and 7 species (for species we considered presence in at least 80% of specimens) (Table 4). Notably, the only represented order within Cyanobacteria was YS2, a recognized member of the novel non-photosynthetic phylum Melainabacteria (Di Rienzi et al.2013) that occurred at very low abundances. Among core genera, most abundant were Bacteroides (average 10% of total reads per specimen) and Parabacteroides (4.6%), followed by Oscillospira, Dorea, Ruminococcus, Desulfovibrio and Coprococcus (all above 1%) (Fig. S3, Supporting Information). Table 4. List of core taxa (90%) by their finest taxonomic resolution, with number of core OTUs per taxon shown in parentheses (total = 134). Phylum Class Order Family Genus Species1 Firmicutes (118) Erysipelotrichi (5) Erysipelotrichales (5) Erysipelotrichaceae (5) Coprobacillus (1) PSB-M-3 Holdemania Clostridia (113) Clostridiales (113) Lachnospiraceae (37) Ruminococcus (8) [Ruminococcus] gnavus Blautia (2) Blautia producta Dorea (10) Coprococcus Roseburia Ruminococcaceae (51) Oscillospira (24) Anaerotruncus Eubacteriaceae (1) Pseudoramibacter (1) Anaerofustis Clostridiaceae Clostridium Dehalobacteriaceae Dehalobacterium Christensenellaceae (4) Christensenella (1) [Mogibacteriaceae] (3) Anaerovorax Peptostreptococcaceae (1) Peptococcaceae (3) Veillonellaceae [Tissierellaceae] Bacilli Lactobacillales Lactobacillaceae Lactobacillus Leuconostocaceae Bacillales Bacillaceae Bacillus Actinobacteria (6) Coriobacteriia (6) Coriobacteriales (6) Coriobacteriaceae (6) Adlercreutzia (1) Proteobacteria (6) Deltaproteobacteria (6) Desulfovibrionales (6) Desulfovibrionaceae (6) Desulfovibrio (5) Bilophila Gammaproteobacteria Enterobacteriales Enterobacteriaceae Alphaproteobacteria Campylobacterales Helicobacteraceae Epsilonproteobacteria Bacteroidetes (4) Bacteroidia (4) Bacteroidales (4) Bacteroidaceae (3) Bacteroides (3) Bacteroides fragilis Bacteroides ovatus Porphyromonadaceae (1) Parabacteroides (1) Parabacteroides distasonis Parabacteroides gordonii Rikenellaceae Rikenella [Odoribacteraceae] Odoribacter Verrucomicrobia Verrucomicrobiae Verrucomicrobiales Verrucomicrobiaceae Akkermansia Akkermansia muciniphila Tenericutes Mollicutes Anaeroplasmatales Anaeroplasmataceae Anaeroplasma RF39 Cyanobacteria (Melainabacteria) 4C0d-2 YS2 Phylum Class Order Family Genus Species1 Firmicutes (118) Erysipelotrichi (5) Erysipelotrichales (5) Erysipelotrichaceae (5) Coprobacillus (1) PSB-M-3 Holdemania Clostridia (113) Clostridiales (113) Lachnospiraceae (37) Ruminococcus (8) [Ruminococcus] gnavus Blautia (2) Blautia producta Dorea (10) Coprococcus Roseburia Ruminococcaceae (51) Oscillospira (24) Anaerotruncus Eubacteriaceae (1) Pseudoramibacter (1) Anaerofustis Clostridiaceae Clostridium Dehalobacteriaceae Dehalobacterium Christensenellaceae (4) Christensenella (1) [Mogibacteriaceae] (3) Anaerovorax Peptostreptococcaceae (1) Peptococcaceae (3) Veillonellaceae [Tissierellaceae] Bacilli Lactobacillales Lactobacillaceae Lactobacillus Leuconostocaceae Bacillales Bacillaceae Bacillus Actinobacteria (6) Coriobacteriia (6) Coriobacteriales (6) Coriobacteriaceae (6) Adlercreutzia (1) Proteobacteria (6) Deltaproteobacteria (6) Desulfovibrionales (6) Desulfovibrionaceae (6) Desulfovibrio (5) Bilophila Gammaproteobacteria Enterobacteriales Enterobacteriaceae Alphaproteobacteria Campylobacterales Helicobacteraceae Epsilonproteobacteria Bacteroidetes (4) Bacteroidia (4) Bacteroidales (4) Bacteroidaceae (3) Bacteroides (3) Bacteroides fragilis Bacteroides ovatus Porphyromonadaceae (1) Parabacteroides (1) Parabacteroides distasonis Parabacteroides gordonii Rikenellaceae Rikenella [Odoribacteraceae] Odoribacter Verrucomicrobia Verrucomicrobiae Verrucomicrobiales Verrucomicrobiaceae Akkermansia Akkermansia muciniphila Tenericutes Mollicutes Anaeroplasmatales Anaeroplasmataceae Anaeroplasma RF39 Cyanobacteria (Melainabacteria) 4C0d-2 YS2 Taxonomy in brackets is not officially recognized (under revision). Not all core OTUs were identified below order level. aPresence in 80% of the specimens and all islets. View Large At OTU level, 23 OTUs were present in 100% of the specimens, 134 at 90% and 303 at 80%. Four phyla hosted the 134 core OTUs at 90% prevalence, Firmicutes, Proteobacteria, Actinobacteria and Bacteroidetes (Table 4), with the majority of core OTUs (118) belonging to the phylum Firmicutes and, within this, to the order Clostridiales (113), families Ruminococcaceae (51), Lachnospiraceae (37), and genera Oscillosphira (24), Ruminococcus (8) (both Ruminococcaceae) and Dorea (10, Lachnospiraceae) (Table 4). Despite the fact that Bacteroidetes and Proteobacteria were among the most abundant phyla, they hosted only four and six core OTUs, respectively. All Proteobacteria OTUs belonged to the family Desulfovibronaceae, five of which were classified to the genus Desulfovibrio. All six Actinobacteria core OTUs were Coriobacteriaceae. Overall, these 134 core OTUs accounted for only 2.1% of the total number of OTUs detected, with individual core OTUs typically occurring at low abundances, ranging on average across all specimens between 0.01% and 2.7% of total reads per specimen (Table S3, Supporting Information). The most abundant core OTU was N.Ref.OTU202 (2.7%, classified as Bacteroides), followed by OTU-113 542 (1.6%, Desulfovibrio), N.Ref.OTU980 (1.6%, Erysipelotrichaceae) and OTU-194 286 (1.5%, Parabacteroides gordonii). Altogether, the 134 core OTUs contributed 34.8% of the total reads across all samples and 35.1% average reads per specimen. According to best BLASTN hits and taxonomic classification (Table S3, Supporting Information), core OTUs (i) typically belonged to bacteria known as anaerobic or strictly anaerobic, suggesting a limitation in their habitat range and potential for dispersal; (ii) were not classified as pathogens at present abundance; and (iii) in most cases matched to stool/gut from humans and other vertebrates, rarely invertebrates (representing the main food source). These features suggest large specificity to the vertebrate gut, making these OTUs good candidates to investigate microbiota heritability across host generations and the effect of long-term geographic isolation. Sub-OTU resolution recovers islet signal To address the above scenario, we exploited the sub-OTU diversity of the core (sequence variants within each OTU as inferred by re-clustering of sequences at 98 and 99% similarity). If core OTUs represent a conserved component inherited from the ancestral microbial pool and then ‘confined’ in each islet after population splits, we would expect the sub-OTU variants to carry a signature of their islet-restricted diversification, including local drift and mutation/recombination events. On the contrary, if acquisition of this core occurred by more recent dispersal across islets, the sub-OTU variants should present a stochastic distribution or correlate with geographic distances, but not with islet structure. For this purpose, we used individual and pooled core OTUs, which carry no information for discrimination among specimens/islets, and reclustered all reads at 98 and 99% similarity (the latter approximating the haplotype level). PERMANOVA analyses on all previously tested categorical variables indicated ‘islet’ as the main significant factor (Table 1, P = 0.0001), explaining 15% of the core microbiota sub-OTU variance, a proportion comparable to the one explained by clustering of the whole dataset at 97% similarity. According to Jaccard distance matrices built on de novo 98% and 99% OTUs, distances were significantly shorter within islets than between islets (two-sided Student's two-sample t-test, Bonferroni-corrected P < 0.001). Moreover, mean Jaccard distances across islets mirrored the 97% OTU-based distances: in all cases, the closest islet pair was Aire-PrCv (mean Jaccard = 0.25 and 0.30 at 98% and 99% clustering, respectively). Again, microbiota distances did not correlate with geographic distances (Mantel r = –0.08838, P = 0.752 for 99% clustering and r = –0.10154, P = 0.674 for 98% clustering). To understand the individual contribution of core OTUs to this pattern, the 23 OTUs with presence in 100% of the specimens were individually analyzed (99% OTU clustering, Table S4, Supporting Information). Intra-OTU variation analysis revealed that 18 of the 23 OTUs carried a significant signal by islet (adonis, P < 0.05), explaining between 14% and 23% of the total sub-OTU variation. The five OTUs with no significance had their best BLASTN hits to putative environmental locations, human skin and a novel bacterium from the human gut (i.e. Emergencia timonensis). DISCUSSION Here we provide the first characterization of the P. lilfordi gut microbiota, exploring its compositional variation in adult males from seven isolated populations inhabiting coastal islets around Menorca. Taking advantage of the peculiar phylogeographic history of these populations and the expected minor diversifying role of two crucial variables in microbiota studies (host genotype and diet), we performed a preliminary investigation of the forces driving current patterns of microbiota similarity/dissimilarity among populations, including inheritance, ecological drift, dispersal and local environmental inputs. Geographic separation and age drive the diversity and structuring of the Podarcis microbiota Compatible with a significant spatial and temporal effect, islet was the most significant factor in structuring the Podarcis microbiota, accounting for 15% of the observed microbiota variation (Table 1). Aire in particular, with its geographic remoteness and morphological diversity (it hosts a unique melanic form) (Salvador 2009), presented the most distinct microbial community as indicated by Unifrac distances (Fig. 3) as well as by taxa/OTU uniqueness and enrichment (Tables 2 and 3). Among the fixed factors, within which islet was nested, only ‘age’ was significant when considering Aire and PrCv as the oldest islets. Most of the microbiota variation remained unexplained by all other factors tested, including geographic distances among islets. Consistent with a minor although significant islet effect, microbiota beta-diversity distances across populations were relatively short (average pairwise Jaccard distance = 0.38) and partly associated with time since islet/population isolation. Although the phylogeographic history of the Menorcan populations is far from being fully resolved (and out of the scope of this study), Aire and PrCv represent the two first islets/populations to separate from mainland about 8000–9000 years (Pérez-Mellado 1989), in agreement with both bathymetry and curves of sea level rise in the Mediterranean during Holocene (Vecchi et al.2016). Despite being the most geographically distant islets within our sample set (Fig. 1), their gut microbial composition showed the shortest pairwise Jaccard distances and presented several unique commonalities (Table 3 and Table S2, Supporting Information). Retention of ancestral taxa in these two islets or convergence driven by similar environmental opportunities (open sea and scarce vegetation) could both be invoked. However, none of the ecological variables tested with PERMANOVA were strictly associated with these two islets (Table S1, Supporting Information), largely excluding convergence by niche selection. In terms of gut community similarities across islets, OTU uniqueness per islet/population was remarkably low, with the large majority of OTUs (98% on average, Table 2) occurring in more than one islet. Quantitative differences in common taxa and OTUs representation across islets were also minor (Table 3 and Table S2, Supporting Information) and associated with bacteria belonging to core taxa (i.e. order Clostridiales and family Ruminococcaceae, Table 4), suggesting minor functional community differences among islets. While, to some extent, this large microbial membership and structural similarity might result from common host genotypic constraints and lack of major dietary differences (all populations were largely omnivorous), important differences in ecological niches exist among islets: some show rich vegetation and large anthropic impact (Colm and Rei), while others are characterized by more natural environments, with scarce vegetation and large exposure to the open sea (Aire and PrCv). Overall, this suggests that the striking compositional similarity observed is, at least in part, independent from the environmental inputs/pressures. Ancestral retention versus dispersal effect To explore the origin of the shared taxonomic content among populations, we looked at the nature of core OTUs (defined as presence in 90% of the specimens and all islets). Core OTUs could either represent a conserved inherited component since divergence from the mainland population or an independent acquisition after isolation through taxa dispersal across islets followed by gut selective retention (common niche selection). Excluding recent host immigration events (at least during the past 2000 years), as supported by bathymetric data and the degree of morphological and pigmentation diversity among these populations (Salvador 2009), dispersal should be largely restricted to other biotic or abiotic vectors able to cross the sea barrier. Invertebrates consumed with the diet represent the most likely biotic vectors, particularly flying insects (i.e. all populations are predominantly insectivorous, especially during the summer season) and occasional mollusks and crustaceans captured along the shores. According to best BLASTN hits and taxonomic classification, however, the majority of core OTUs typically matched gut/stool samples of vertebrates, with rare invertebrate hits (Table S3, Supporting Information). This suggests that despite invertebrates (mostly insects) represent the largest component of these lizards’ diet, they are not a major source of gut bacteria. This is in accordance with previous findings showing that invertebrate communities did not significantly contribute to the lizard gut microbiota (Ren et al.2016). Likewise, vertebrates (including rats, rabbits, seagull colonies and humans) are improbable vectors of pervasive dispersal of shared OTUs, as their presence is largely restricted to a few islets (Table S1, Supporting Information; except for Colm and Rei, none of the islets hosts stable human settlements). Aside from animal vectors, recent studies have shown that the phyllosphere microbiome might play a substantial role in shaping the gut microbiota in both lizards (Kohl et al.2017) and desert woodrats (Kohl and Dearing 2012); in this respect, P. lilfordi populations are known to consume vegetable matter (Salvador 2009; Perez-Cembranos, Leon and Perez-Mellado 2016), a finding that is also confirmed by stomach content analysis from this study, revealing plant material and flowers in the majority of samples, and by the presence of several fiber-degrader taxa in the gut (see below). While our current lack of microbial samples from the surrounding environment does not allow us to test the impact of the phyllosphere in these lizards’ microbiota, both best BLASTN matches and large differences in vegetation cover among islets make this scenario improbable. A second important observation is that the 98% and 99% OTU-based structure of the core presented a significant clustering by islet and age (Table 1), which is compatible with an islet-restricted diversification of core OTUs (putatively driven by loss of variants, mutation and recombination events). This hypothesis clearly needs to be corroborated by strain-level analyses. Moreover, these data might suffer from undersampling in some populations/individuals, despite the high sequence depth coverage per specimen obtained (38 160 sequences after rarefaction). Finally, whether any of these sub-OTU variants are under selection or exist in a transitory neutral state remains to be investigated. Indirect evidence of gut microbiome heritability Clearly, the presence of shared OTUs, that are unlikely the result of dispersal, indirectly stands for the inheritance of a microbial component over relatively large temporal and spatial scales. Up to date, levels and mechanisms of gut microbiota inheritance in lizards, and reptiles in general, remain unclear (Colston 2017). Gut microbiotas from Anolis lizards showed a remarkably high intraspecific variation, with only 10% on average of shared OTUs between any two conspecifics in nature (Ren et al.2016), suggesting low microbiota heritability in these lizards. Although their large generalist diet was proposed as a major diversifying factor, this is unlike our findings for P. lilfordi, also a generalist species, where any two conspecifics within an islet share on average 31% of their OTUs and populations are allopatric. Methodological differences, such as the use of full intestines in our study versus fecal samples for Anolis, could be partly responsible for the observed pattern. In contrast, for viviparous lizards of genera Liolaemus and Phymaturus, recent evidence indicated around 34% of microbiome heritability in captivity (Kohl et al.2017), an estimate that nicely matches the percentage of core microbial sequences observed in Podarcis (34.8%). While the Podarcis core summed up to more than one-third of the microbiota, individual core OTUs typically occurred at low abundances (mean = 0.26 ± 0.4% over total reads per specimen, Table S3, Supporting Information), suggesting that, while being more vulnerable to stochastic loss (Chase and Myers 2011), their persistence over time is sustained by important selective forces and putative high inheritability. This is compatible with the general understanding that core gut microbiotas in vertebrates are rather stable in the absence of differential adaptive pressures (such as type of metabolism) (Ley et al.2008; Lozupone et al.2012). In this respect, the overall taxonomic composition of the P. lilfordi gut microbiota mirrored that found in other vertebrates (Fig. 2), with the majority of core families and genera detected (Table 4) representing stable residents of human and mouse guts (Xiao et al.2015), and being consistently found in other lizards (Ren et al.2016; Kohl et al.2017). Core taxa were particularly enriched in bacteria known as active fiber degraders, including the families Ruminococcaceae and Lachnospiraceae (Fig. 2, both Clostridiales) (Biddle et al.2013) and the genera Bacteroides, Blautia, Akkermansia, Oscillospira and Desulfovibrio (Fig. S3, Supporting Informtion). Bacteroides, for instance, which was the most conspicuous genus in our samples, also makes up a very large proportion of the mammalian gut where it is an active degrader of plant material (Ley et al.2008; Xiao et al.2015). Oscillospira, also highly abundant in humans, where it shows high heritability and is associated with leanness, is largely enriched in plant-fed lizards (Kohl et al.2017). Overall, the P. lilfordi core composition is compatible with major conserved functions in vertebrates and the consumption of plant material, suggesting common selective pressures for its maintenance across host generations. Nonetheless, considering that parental care in reptiles is rare and not documented in the genus Podarcis (oviparous), means of bacteria recycling at each new host generation remain unclear. Coprophagy, cannibalism or other types of social behavior have been proposed as mechanisms for bacteria transmission (Colston and Jackson 2016; Perez-Cembranos, Leon and Perez-Mellado 2016; Colston 2017); however, they might suffer from much stochasticity, compromising their ability to preserve a metacommunity structure over time. Covering the eggs with a bacterial layer by means of fecal contamination during oviposition through the single urogenital opening present in lizards (the cloaca) might represent a more reliable mechanism for a stable microbiota transmission. Overall, although our results provided important indirect evidence for gut microbiome heritability in P. lilfordi, confirming previous data for lizards (Kohl et al.2017), this topic definitely demands more empirical investigation. Islet microbiotas as nested subsets of the ancestral mainland microbial pool: a hypothesis To summarize, our main findings indicated that (i) there is a large taxonomic microbiota conservation among all seven allopatric populations studied, down to OTU level, with a reduced impact of local exposure (i.e. uniqueness is negligible); (ii) the minor variation detected is explained by ‘islet’ and ‘age’, supporting a role of the host phylogeographic history in shaping microbiota composition; (iii) dispersal of OTUs among islets appears to be low; and (iv) core variation at 98%–99%-OTU level carries a significant ‘islet’ and ‘age’ signature, compatible with an islet-restricted evolution. Altogether, these findings suggest that persistence of taxa from the mainland stock and ecological drift might have played a substantial role in shaping current pattern of microbiota similarity/dissimilarity among P. lilfordi populations. This intriguing hypothesis clearly requires more extensive quantitative data, targeting additional specimens, islets and microbial markers, as well as ad hoc statistical tools to test drift in host-associated microbiota across populations, a field still in its infancy. Although overlooked in current literature, host population demographic history is likely to represent a crucial factor in explaining current gut metacommunity diversity, as host drift could, to some extent, cause microbiota ecological drift (Orrock and Watling 2010; Chase and Myers 2011). Up to date, only few studies have attempted to address this topic. In the Galápagos archipelago, processes of ecological drift were proposed to explain intraspecific differences in marine iguanas (Lankau, Hong and Mackie 2012). Likewise, in two Anolis lizard species, the complex pattern of similarity/dissimilarity observed, only partially attributed to geography and with no differences among ‘ecomorphs’ (Ren et al.2016), could be potentially explained by chance variation driven by drift and dispersal. In our study system, population densities varied greatly across islets (Fig. 1), with Aire presenting one of the highest population densities of lizards/vertebrates worldwide (Pérez-Mellado et al.2008). Although past population demography for these islets is poorly known (Pérez-Mellado et al.2008), the largest populations (Aire, PrCv and Colm) carried slightly higher diversities (PD whole tree) compared to smaller islets (Tosq, Ravl and Rei), though differences were not significant (Fig. S2, Supporting Information). The small population found in Tosqueta, in particular, which showed the strongest genetic drift among all 16 islets according to unpublished microsatellite data (personal observations), carried the least diverse microbial community and lowest uniqueness, despite the extensive vegetation cover and ample opportunities for local inputs in this islet. These observations, although preliminary, point to an interesting correlation between population demographic history and drift effect in gut metacommunities, which could potentially explain a large part of the extant variation, a topic that should definitely warrant further investigation. CONCLUSIONS Overall, the minor dietary and genotypic differences among our sampled populations (suggesting weak diversifying selection) and the geographically restricted areas to which they are confined, reducing the exposure to microbial pools from other vertebrates, have most likely contributed to the large compositional stability of P. lilfordi gut microbiota, ultimately revealing the strength of drift. We envision that, owing to a progressive within-islet microbiota evolution, in a context of no host immigration events, the core OTU signature of these populations’ past host legacy will gradually disappear, while resemblance at higher taxonomic levels will persist for a longer time. Integrating our data with a fine-scale resolution of P. lilfordi phylogeographic history and a bacterial strain-level characterization that targets core OTUs (under development) might offer a way to disentangle the dynamics of this microbiota component over such short evolutionary time scales. Finally, extending the sampling to all remaining populations from Balearic Islands, including the closely related species P. pityusensis, will provide important comparative data for a most integrative understanding of the gut microbiota diversification in these lizards. Data accessibility Raw sequences are available in the Sequence Read Archive (SRA) database at NCBI under BioProject ID PRJNA429123PRJNA429123. SUPPLEMENTARY DATA Supplementary data are available at FEMSEC online. Acknowledgements We thank Sam Pons, Jordi García-Darocas and Montserrat Farré for aid in sampling, and the Animal Biology department (University of Barcelona) for use of equipment. FUNDING This work was supported by funding from the Spanish Ministerio de Economía y Competitividad to LB (RYC-2012-11443). Conflicts of Interest. None declared. REFERENCES Aivelo T, Norberg A. Parasite-microbiota interactions potentially affect intestinal communities inwildmammals. J Anim Ecol 2017, doi: 10.1111/1365-2656.12708. Alcover JA, Moyá-Solà S, Pons-Moyá J. Les quimeres del passat. Els vertebrats fòssils del Plio-Quaternari de les Balears i Pitiüses . Palma, Espanya: Editorial Moll, 1981. Anderson MJ. A new method for non-parametric multivariate analysis of variance. Austral Ecol 2001; 26: 32– 46. Anderson MJ, Gorley RN, Clarke KR. PERMANOVA+ for PRIMER: Guide to Software and Statistical Methods . Plymouth, UK: PRIMER-E, 2008. Baldo L, Pretus JL, Riera JL et al. Convergence of gut microbiotas in the adaptive radiations of African cichlid fishes. ISME J 2017; 11: 1975– 87. Google Scholar CrossRef Search ADS PubMed Baldo L, Riera JL, Tooming-Klunderud A et al. Gut microbiota dynamics during dietary shift in eastern african cichlid fishes. PLoS One 2015; 10: e0127462. Google Scholar CrossRef Search ADS PubMed Benson AK, Kelly SA, Legge R et al. Individuality in gut microbiota composition is a complex polygenic trait shaped by multiple environmental and host genetic factors. P Natl Acad Sci USA 2010; 107: 18933– 8. Google Scholar CrossRef Search ADS Biddle A, Stewart L, Blanchard J et al. Untangling the genetic basis of fibrolytic specialization by lachnospiraceae and ruminococcaceae in diverse gut communities. Divers Distrib 2013; 5: 627– 40. Google Scholar CrossRef Search ADS Brown RP, Terrasa B, Perez-Mellado V et al. Bayesian estimation of post-Messinian divergence times in Balearic Island lizards. Mol Phyl Evol 2008; 48: 350– 8. Google Scholar CrossRef Search ADS Burns AR, Stephens WZ, Stagaman K et al. Contribution of neutral processes to the assembly of gut microbial communities in the zebrafish over host development. ISME J 2016; 10: 655– 64. Google Scholar CrossRef Search ADS PubMed Caporaso JG, Kuczynski J, Stombaugh J et al. QIIME allows analysis of high-throughput community sequencing data. Nat Methods 2010; 7: 335– 6. Google Scholar CrossRef Search ADS PubMed Chase JM, Myers JA. Disentangling the importance of ecological niches from stochastic processes across scales. Philos T Roy Soc B 2011; 366: 2351– 63. Google Scholar CrossRef Search ADS Clarke KR, Gorley RN. PRIMER v7: User Manual/Tutorial . Plymouth, UK: PRIMER-E, 2015. Colston TJ. Gut microbiome transmission in lizards. Mol Ecol 2017; 26: 972– 4. Google Scholar CrossRef Search ADS PubMed Colston TJ, Jackson CR. Microbiome evolution along divergent branches of the vertebrate tree of life: what is known and unknown. Mol Ecol 2016; 25: 3776– 800. Google Scholar CrossRef Search ADS PubMed Costello EK, Stagaman K, Dethlefsen L et al. The application of ecological theory toward an understanding of the human microbiome. Science 2012; 336: 1255– 62. Google Scholar CrossRef Search ADS PubMed David LA, Maurice CF, Carmody RN et al. Diet rapidly and reproducibly alters the human gut microbiome. Nature 2014; 505: 559– 63. Google Scholar CrossRef Search ADS PubMed Delsuc F, Metcalf JL, Wegener Parfrey L et al. Convergence of gut microbiomes in myrmecophagous mammals. Mol Ecol 2014; 23: 1301– 17. Google Scholar CrossRef Search ADS PubMed Di Rienzi SC, Sharon I, Wrighton KC et al. The human gut and groundwater harbor non-photosynthetic bacteria belonging to a new candidate phylum sibling to Cyanobacteria. eLife 2013; 2: e01102. Google Scholar CrossRef Search ADS PubMed Funkhouser LJ, Bordenstein SR. Mom knows best: the universality of maternal microbial transmission. PLoS Biol 2013; 11: e1001631. Google Scholar CrossRef Search ADS PubMed Kohl KD, Brun A, Magallanes M et al. Gut microbial ecology of lizards: insights into diversity in the wild, effects of captivity, variation across gut regions and transmission. Mol Ecol 2017; 26: 1175– 89. Google Scholar CrossRef Search ADS PubMed Kohl KD, Dearing MD. Experience matters: prior exposure to plant toxins enhances diversity of gut microbes in herbivores. Ecol Lett 2012; 15: 1008– 15. Google Scholar CrossRef Search ADS PubMed Krijgsman W, Hilgen FJ, Raffi Y et al. Chronology, causes and progression of the Messinian salinity crisis. Nature 1999; 400: 652– 5. Google Scholar CrossRef Search ADS Lankau EW, Hong PY, Mackie RI. Ecological drift and local exposures drive enteric bacterial community differences within species of Galapagos iguanas. Mol Ecol 2012; 21: 1779– 88. Google Scholar CrossRef Search ADS PubMed Ley RE, Lozupone CA, Hamady M et al. Worlds within worlds: evolution of the vertebrate gut microbiota. Nat Rev Microbiol 2008; 6: 776– 88. Google Scholar CrossRef Search ADS PubMed Lozupone CA, Stombaugh JI, Gordon JI et al. Diversity, stability and resilience of the human gut microbiota. Nature 2012; 489: 220– 30. Google Scholar CrossRef Search ADS PubMed Martinson VG, Douglas AE, Jaenike J. Community structure of the gut microbiota in sympatric species of wild Drosophila. Ecol Lett 2017; 20: 629– 39. Google Scholar CrossRef Search ADS PubMed Moeller AH, Caro-Quintero A, Mjungu D et al. Cospeciation of gut microbiota with hominids. Science 2016; 353: 380– 2. Google Scholar CrossRef Search ADS PubMed Muegge BD, Kuczynski J, Knights D et al. Diet drives convergence in gut microbiome functions across mammalian phylogeny and within humans. Science 2011; 332: 970– 4. Google Scholar CrossRef Search ADS PubMed Nemergut DR, Schmidt SK, Fukami T et al. Patterns and processes of microbial community assembly. Microbiol Mol Biol R 2013; 77: 342– 56. Google Scholar CrossRef Search ADS Orrock JL, Watling JI. Local community size mediates ecological drift and competition in metacommunities. Proc Biol Sci 2010; 277: 2185– 91. Google Scholar CrossRef Search ADS PubMed Perez-Cembranos A, Leon A, Perez-Mellado V. Omnivory of an insular lizard: sources of variation in the diet of Podarcis lilfordi (Squamata, Lacertidae). PLoS One 2016; 11: e0148947. Google Scholar CrossRef Search ADS PubMed Pérez-Mellado V. Estudio ecológico de la lagartija balear Podarcis lilfordi (Günther, 1874) en Menorca. Rev de Menorca 1989; 53: 455– 511. Pérez-Mellado V, Hernández-Estévez JA, García-Díez T et al. Population density in Podarcis lilfordi (Squamata, Lacertidae), a lizard species endemic to small islets in the Balearic Islands (Spain). Amphibia-Reptilia 2008; 29: 49– 60. Google Scholar CrossRef Search ADS Pretus JL, Marqués R, Pérez-Mellado V. Holocene sea level rise 975 and range of fragmentation of Podarcis lilfordi on Minorcan islets: a vicariance scenario reviewed through a mtDNA tree. In Pérez-Mellado V, Riera N, Perera A, (eds). The Biology of lacertid lizards: Evolutionary and Ecological Perspectives . Institut Menorquí d'Estudis . 2004; 8: 279– 291. Price MN, Dehal PS, Arkin AP. FastTree: computing large minimum evolution trees with profiles instead of a distance matrix. Mol Biol Evol 2009; 26: 1641– 50. Google Scholar CrossRef Search ADS PubMed Rawls JF, Mahowald MA, Ley RE et al. Reciprocal gut microbiota transplants from zebrafish and mice to germ-free recipients reveal host habitat selection. Cell 2006; 127: 423– 33. Google Scholar CrossRef Search ADS PubMed Ren T, Kahrl AF, Wu M et al. Does adaptive radiation of a host lineage promote ecological diversity of its bacterial communities? A test using gut microbiota of Anolis lizards. Mol Ecol 2016; 25: 4793– 804. Google Scholar CrossRef Search ADS PubMed Saez ET, Traveset A. Fruit and nectar feeding by Podarcis lilfordi (Lacertidae) on Cabrera archipelago (Balearic islands). Herpetol Rev 1995; 26: 121– 3. Salvador A. Lagartija balear – Podarcis lilfordi (Günther, 1874). In Salvador A, Marco A., (eds). Museo Nacional de Ciencias Naturales, Madrid . Enciclopedia Virtual de los Vertebrados Españoles, 2009. http://www.vertebradosibericos.org. Sanders ED. Evidence Concerning Late Survival and Extinction of Endemic Amphibia and Reptilia from the Bronze and Iron Age settlement of Torralba d’En Salort (Alayor, Menorca) . Palma de Mallorca: Monografies Cientìfiques, Palma de Mallorca. 1984. Schloss PD, Westcott SL, Ryabin T et al. Introducing mothur: open-source, platform-independent, community-supported software for describing and comparing microbial communities. Appl Environ Microb 2009; 75: 7537– 41. Google Scholar CrossRef Search ADS Shoemaker W, Locey KJ, Lennon JT. A macroecological theory of microbial biodiversity. Nat Ecol Evol 2017; 1: 107. Google Scholar CrossRef Search ADS PubMed Song SJ, Lauber C, Costello EK et al. Cohabiting family members share microbiota with one another and with their dogs. eLife 2013; 2: e00458. Google Scholar PubMed Terrasa B, Perez-Mellado V, Brown RP et al. Foundations for conservation of intraspecific genetic diversity revealed by analysis of phylogeographical structure in the endangered endemic lizard Podarcis lilfordi. Divers Distrib 2009; 15: 207– 21. Google Scholar CrossRef Search ADS Traveset A, Riera N. Disruption of a Plant-Lizard seed dispersal system and its ecological effects on a threatened endemic plant in the Balearic islands. Conserv Biol 2005; 19: 421– 31. Google Scholar CrossRef Search ADS Vecchi M, Marriner N, Morhange C et al. Multiproxy assessment of Holocene relative sea-level changes in the western Mediterranean: Sea-level variability and improvements in the definition of the isostatic signal. Earth-Sci Rev 2016; 155: 172– 97. Google Scholar CrossRef Search ADS Xiao L, Feng Q, Liang S et al. A catalog of the mouse gut metagenome. Nat Biotechnol 2015; 33: 1103– 8. Google Scholar CrossRef Search ADS PubMed Yatsunenko T, Rey FE, Manary MJ et al. Human gut microbiome viewed across age and geography. Nature 2012; 486: 222– 7. Google Scholar PubMed Zeng Q, Sukumaran J, Wu S et al. Neutral models of microbiome evolution. PLoS Comput Biol 2015; 11: e1004365. Google Scholar CrossRef Search ADS PubMed © FEMS 2017. All rights reserved. For permissions, please e-mail: email@example.com
FEMS Microbiology Ecology – Oxford University Press
Published: Feb 1, 2018
It’s your single place to instantly
discover and read the research
that matters to you.
Enjoy affordable access to
over 18 million articles from more than
15,000 peer-reviewed journals.
All for just $49/month
Query the DeepDyve database, plus search all of PubMed and Google Scholar seamlessly
Save any article or search result from DeepDyve, PubMed, and Google Scholar... all in one place.
Get unlimited, online access to over 18 million full-text articles from more than 15,000 scientific journals.
Read from thousands of the leading scholarly journals from SpringerNature, Elsevier, Wiley-Blackwell, Oxford University Press and more.
All the latest content is available, no embargo periods.
“Hi guys, I cannot tell you how much I love this resource. Incredible. I really believe you've hit the nail on the head with this site in regards to solving the research-purchase issue.”Daniel C.
“Whoa! It’s like Spotify but for academic articles.”@Phil_Robichaud
“I must say, @deepdyve is a fabulous solution to the independent researcher's problem of #access to #information.”@deepthiw
“My last article couldn't be possible without the platform @deepdyve that makes journal papers cheaper.”@JoseServera