- 1Department of Biology, Faculty of Science, Bülent Ecevit University, Zonguldak, Türkiye
- 2Department of Zoology, Charles University, Prague, Czechia
- 3Institute of Vertebrate Biology, Academy of Sciences of the Czech Republic, Brno, Czechia
- 4Department of Food Engineering, Izmir Institute of Technology, Izmir, Türkiye
- 5Institute of Ecology, Ilia State University, Tbilisi, Georgia
- 6Institute of Physiology, University Hospital Essen, University of Duisburg-Essen, Essen, Germany
- 7Department of Biology, Dokuz Eylül University, İzmir, Türkiye
The animal gut microbiome acts as a crucial link between the host and its environment, playing a vital role in digestion, metabolism, physiology, and fitness. Using 16S rRNA metabarcoding, we investigated the effect of altitude on the microbiome composition of Anatolian Blind Mole Rats (Nannospalax xanthodon) across six locations and three altitudinal groups. We also factored in the host diet, as well as host microsatellite genotypes and thyroid hormone levels. The altitude had a major effect on microbiome composition, with notable differences in the relative abundance of several bacterial taxa across elevations. Contrary to prior research, we found no significant difference in strictly anaerobic bacteria abundance among altitudinal groups, though facultatively anaerobic bacteria were more prevalent at higher altitudes. Microbiome alpha diversity peaked at mid-altitude, comprising elements from both low and high elevations. The beta diversity showed significant association with the altitude. Altitude had a significant effect on the diet composition but not on its alpha diversity. No distinct altitude-related genetic structure was evident among the host populations, and no correlation was revealed between the host genetic relatedness and microbiome composition nor between the host microbiome and the diet. Free thyroxine (FT4) levels increased almost linearly with the altitude but none of the bacterial ASVs were found to be specifically associated with hormone levels. Total thyroxine (TT4) levels correlated positively with microbiome diversity. Although we detected correlation between certain components of the thyroid hormone levels and the microbiome beta diversity, the pattern of their relationship remains inconclusive.
1 Introduction
The animal microbiome, which constitutes the microbial community within the gastrointestinal system, profoundly influences diverse aspects of the host, encompassing digestion, development, immunity, and energetics (Lindsay et al., 2020; McFall-Ngai et al., 2013; Suzuki, 2017). It also has the potential to enhance the host’s evolutionary capacity and fitness by shaping its phenotype (Alberdi et al., 2021; Henry et al., 2019; Suzuki, 2017; Henry et al., 2021). Humans, laboratory animals, and domestic mammals often serve as valuable models to study the microbiome’s importance in host ecology, evolution, health, and disease (Lin and Zhang, 2017; Lynch and Pedersen, 2016; Thursby and Juge, 2017; Wang et al., 2019; Bäckhed et al., 2007). While studying captive animals and humans is convenient, this approach has limitations, including the effect of controlled captivity conditions on the microbiome (Rogers et al., 2014; Roeselers et al., 2011; Van Leeuwen et al., 2020; Belheouane et al., 2020), lack or loss of heterogeneity in model (or held in captivity) hosts and their microbiome (Wang et al., 2014; Kohl et al., 2014), lack of environmental pressure (Chong-Neto et al., 2022), and variation in lifestyles and diets in humans (De Filippo et al., 2010) make it challenging to extrapolate the significance of findings to natural environments. Thus, employing natural models in research becomes not just a choice but a necessity whenever feasible to reveal the implication of host-microbiota interactions in ecological and evolutionary processes (Hird, 2017; Greyson-Gaito et al., 2020).
The microbiome is influenced by various factors, such as habitat variation (Amato et al., 2013), lifestyle (McKenzie et al., 2017), host diet (Pellizzon and Ricci, 2018; De Filippo et al., 2010; Kuang et al., 2022), host phylogeny (Ley et al., 2008), reproductive status of host (Nuriel-Ohayon et al., 2016), social interactions (Li et al., 2016c; Grieneisen et al., 2017) and environment (Rothschild et al., 2018; Chong-Neto et al., 2022; Ahn and Hayes, 2021). Among these factors, the environment holds significant importance as it can impact the physiology and immune responses of the host, while also playing a pivotal role in shaping dietary preferences and the selection process of the diverse pool of microbes in the surrounding environment. The microbiome can also help the host to adapt to new environments by modulating gene expression associated with nutrient metabolism, enhancing immune function, and providing essential nutrients (Bisschop et al., 2022; Petersen et al., 2023; Henry et al., 2021). In addition, microbiome diversity correlated positively with metabolic rate in a few studies, e.g., in wild and captive giant pandas (Zhu et al., 2011) and in hibernating brown bears (Sommer et al., 2016).
Highland environments are characterized by decreased atmospheric oxygen levels and colder temperatures, which poses various challenges to animals. Specifically, hypoxia tolerance and energy-efficient metabolism among species living at high altitude is the most common adaptation mechanisms (Storz and Moriyama, 2008; Storz et al., 2019; Cheviron et al., 2014; Schippers et al., 2012). Interestingly, the high-altitude adaptations do not seem costly, for example, humans living at high altitude have longer lifespans (Midha et al., 2023; Rogers et al., 2023). Recently, gut microbiome (hereafter GM) with specific diversity, composition, and function have been proposed as an additional factor that could contribute to high-altitude adaptation. For example, strong positive correlation has been demonstrated between the elevation and the proportion of anaerobic bacteria in the GM of wild house mice (Mus musculus) (Suzuki et al., 2018), which suggests that microbiota might contribute to the adaptation to hypoxia by helping to cope with low oxygen levels. Studies on wild pika (Ochotonidae) have demonstrated that an increase in altitude positively affects the diversity of gut microbial community, both at the host individual (alpha-diversity) and the population (beta-diversity) level (Li et al., 2019). Several species of wild ungulates on the Tibetan plateau have more diverse microbiomes, suggesting specific microbiome changes associated with high altitudes that facilitate extracting more energy from the plant diet (Ma et al., 2019).
Despite the growing body of research on the variation of microbiota with the altitude, it is not easy to disentangle which components of this variation represent genuine adaptations to high altitude environments. Apart from the environmental factors that systematically change with the altitude (e.g., oxygen level or temperature), there are multiple other modulators of gut microbiota that can co-vary with the altitude. One of them is host genetics, (Tabrett and Horton, 2020; Dąbrowska and Witkiewicz, 2016), which is a basis of evolutionary adapation (Suzuki et al., 2019). Because host species colonization history often follows the altitudinal gradient, the genome-wide population differentiation is expected to vary with altitude. The change in altitude also affects the composition of the ecological communities (Gale, 2004; Di Musciano et al., 2021; Lee et al., 2021), which in turn shape the diet of animals (De Filippo et al., 2010; Kuang et al., 2022). The diet stands as one of the primary factors affecting the microbiome composition. Finally, seasonal variaion might also affect the microbiome diversity and composition (Guo et al., 2021; Hu et al., 2018; Jiang et al., 2021; Fan et al., 2022).
Some study systems might be more suitable for disentangling the direct drivers of gut microbiota composition and diversity from those that simply co-vary with the altitude. Subterranean rodents, particularly Nannospalax xanthodon (Anatolian Blind Mole Rat, hereafter ABMR) possesses a unique combination of several traits that could offer an opportunity to disentangle such effects. The AMBR is an obligate subterranean rodent species (Arslan et al., 2016), setting it apart from many of its terrestrial and social rodent counterparts. This distinct ecological niche is likely to have profound implications for the composition and dynamics of its gut microbiota. Unlike social rodents that live in close-knit communities, ABMRs are solitary, each constructing an intricate network of underground tunnels and aggressively defending them from conspecifics, except for mating and raising young (Sözen, 2005; Nevo, 1979, 2007). In contrast to social rodents, this solitary lifestyle limits direct individual interactions, and thus slows down the exchange of microbial communities within the population. The subterranean environment, characterized by more stable temperature and humidity, can also be severely hypoxic when compared to the ambient atmospheric O2 levels. This effect could be exacerbated by the high altitude hypoxia, putting pressure on the aerobic microbes and favoring the anaerobic ones (Suzuki et al., 2018). Another potential factor relevant for gut bacteria is the increased concentration of carbon dioxide inside the AMBR tunnel networks, since some gut bacteria may compete for CO2 as a substrate (Nollet and Verstraete, 1996).
AMBR inhabits various regions of Turkey, from the warm Mediterranean coast of the Aegean Sea to the cold alpine climates of the Taurus Mountains and Eastern Anatolia. At the same time, its dispersing ability is expected to be quite limited compared to above-ground rodents. Slow dispersal, combined with small population size, is expected to allow ample time for the various types of environmental adaptation, which may include changes in microbiome composition, to manifest into notable differences between neighboring populations. Finally, the AMBR feeds on a wide range of plant species, predominantly those with nutrient-rich root parts, yet with a high fibre content (Sözen, 2005; Heth et al., 1989; Ülgen and Tavşanoğlu, 2024). This factor alone is expected to contribute to a high diversity of GM, offering a rich material for exploring the roles of the various factors listed above. To date, only three studies addressed the microbiome composition of the representatives of the Nannospalax genus (Sibai et al., 2020; Solak et al., 2023; Kuang et al., 2022), and just two of them were set in a natural environment (Kuang et al., 2022; Solak et al., 2023).
Notably, the ABMR is in fact a taxonomic complex of multiple, cytogenetically distinct, and potentially genetically isolated geographic populations (Arslan et al., 2016). At the same time, many (cyto) genetically uniform populations have continuous distributions occupying a range of diverse habitats. In this study, we focused on one the most striking example from the Central Taurus (Bolkar) mountains, where a single chromosomal race of ABMR (“cilicicus”) with a diploid chromosome number 2n = 58 occurs at the elevations from ca 1,000 m to 3,000 m above sea level (Sözen et al., 2006). Employing a culture-independent, advanced high-throughput metabarcode amplicon sequencing approach, we investigated: (i) the gut microbiota of ABMR using the 16S rRNA, (ii) its diet composition using 18S rRNA, (iii) thyroid hormone levels as an indirect proxy of the metabolic rate and (iv) the host genetic structure using microsatellite markers. The primary aim was to determine whether and how the altitude in conjunction with other factors influences the composition of the GM of ABMR.
2 Methods
2.1 Study setup and sampling
Our sampling site is located in the Central Taurus Mountains in Türkiye, encompassing an altitudinal gradient ranging between 1,000 m asl and 3,000 m asl (Figure 1). The annual mean temperature and mean precipitation at the highest sampling location are 0.16°C and 744.6 mm, respectively, while they are 11°C and 344.7 mm at the lowest sampling location. Partial snow cover is present until June at high altitudes and may limit digging activities and access to water for these animals.
Figure 1. The sampling area in the Central Taurus mountains, Türkiye. The insert on the right shows the geographic location within the Anatolian peninsula.
The sampling was conducted at three different altitude categories, each represented by two locations. The study area is about 460 km2, the distance between farther most localities (Eregli and Kiziltepe) is about 68 kilometers and the average distance between localities (except Eregli) is 8.5 kilometers. To minimize the influence of seasonal variation, the sampling was carried out during the same 2 months (late June to mid-July) in 2019 and 2021. From each locality, at least eight samples were collected, resulting in a total of 65 individuals. However, some samples failed during the PCR step of 16S and 18S rRNA library preparation. Sampling details are shown in Figure 1 and Table 1.
Animals were live-captured by opening the burrow passageways at specific sections and blocking retreat when the animal mended the tunnel (Wertheim and Nevo, 1971). After recording body mass and sex, animals were euthanized and dissected. Approximately 1 ml of whole blood was taken from the heart, kept overnight at ~4\u00B0C, and centrifuged at 7 RPM to separate the serum. The serum samples and entire digestive tracts were placed on dry ice and kept frozen until transferred to a − 80°C freezer within several days. Additionally, a ~ 5 g caecum sample (both the content and the caecum wall) was stored in ethanol immediately upon dissection. The procedure was approved by the Animal Ethics Committee of Bülent Ecevit University (#91330202).
2.2 Genotyping of gut microbiome and diet content
A small (~5 g) piece from the terminal end of the caecum was collected and ~ 0.25 g of the tissue with its content was used for DNA extraction using the DNEasy PowerSoil Kit (Qiagen, Cat No: 47014) following the manufacturer’s protocol.
For the microbiome genotyping we used 17 individuals from low altitude, 15 individuals from middle altitude, and 19 individuals from high altitude, with a minimum of 7 samples per locality, on a total of 51 samples (Table 1). For the diet genotyping we used 12 individuals from low altitude, 14 individuals from middle altitude, and 12 individuals from high altitude, with a minimum of 2 samples per locality, on a total of 38 samples. On diet genotyping, 13 samples failed during the PCR amplification stage.
To amplify the V3-V4 variable region of the bacterial 16S rRNA gene, the universal primers S-D-Bact-0341-b-S-17(5’-CCTACGGGNGGCWGCAG-3′) and S-D-Bact-0785-a-A-21 (5’-GACTACHVGGGTATCTAATCC-3′) (Klindworth et al., 2013) were utilized. The 5′ ends of the primers were extended with inline barcodes to increase the multiplexing capacity (Supplementary Table S1). PCR reactions done in 10 μL volume, with 1x KAPA HIFI Hot Start Ready Mix (Kapa Biosystems, United States) and each primer at 0.2 μM and 4.6 μL of DNA template with the following cycling conditions: initial denaturation at 98°C for 5 min, followed by 30 cycles of 98°C (15 s), 55°C (20 s), and 72°C (40 s), and a final extension at 72°C for 5 min. The dual-indexed Nextera sequencing adapters were ligated during the second PCR, which is the same cycling conditions except it was performed in 20 μL volume, the concentration of each primer was 1 μM, 1.5 μL of the first PCR product diluted × 12.5 was used as a template and the number of PCR cycles was 12. The quantity and the expected length of the PCR products were evaluated by running them on a 1.5% agarose gel and successfully amplified bands were pooled equimolarly, and purified with SpriSelect beads (Beckman Coulter, United States). The resulting pool was subjected to size selection using the Pippin Prep automatic size selection system (Sage Science), targeting an amplicon size window of 520–750 bp. The pool of libraries was sequenced using MiSeq (Illumina, United States) and v3 chemistry (i.e., 2 × 300-bp paired-end reads) at the CEITEC Genomics Core Facility (Brno, Czech Republic).
To amplify 110–150 bp of eukaryotic 18S rRNA, a similar two-step PCR protocol was employed. We modified the primers from (Guardiola et al., 2015). F primer appended with 4–5 bp inline barcode to increase multiplexing possibilities and both primers appended with the “tail sequence” for priming the second PCR (Supplementary Table S2). Two elongation blockers were added to the first PCR to prevent amplification of (1) the mole rat DNA (Rodent_blckr: +G + T + C + C + C + C + C + A + A + C + T/3SpC3/; plus sign stands for LNA modification; 3SpC3 is C3 spacer at the 3′ end) and (2) the intestinal nematode common in rodents, Syphacia spp. (Syphacia_blckr: +T + G + T + C + T + G + A + A + A + T + A + C + T/SpcC3/). The final sequence of primers are: Guard_18S_F_longer: GATYTGTCTGGTTVATTCCGand Guard_18S_R_longer: CATCACAGACCTGTTATYGC. First PCR was performed in 10 ul reaction (1x KAPA HIFI Hot Start Ready Mix, 0.3 μM of each primer, 1uM of Rodent_blckr, 2 uM of Syphacia_blckr and 1.4 ul of caecal DNA), with 3 min at 95°C and 34 cycles of 98°C (20 s), 56.5°C (15 s) and 70°C (15 s), and a final extension at 70°C (30 s). The second PCR was performed in 15 uL (1x KAPA HIFI Hot Start Ready Mix, 0.5 uM of each indexed primer, and 1 uL of the first-PCR product as template), with 3 min at 95°C and 12 cycles of 98°C (20 s), 55°C (15 s) and 72°C (15 s), and a final extension at 72°C (3 min). Each PCR was performed in a technical duplicate. PCRs were quantified using 2% agarose gel electrophoresis and pooled equimolarly. Fragments between 240 and 360 bp were extracted using PippinPrep and sequenced with Illumina Nextseq, 150 bpPE reads at the CEITEC Genomics Core Facility (Brno, Czech Republic).
To account for possible amplification stochasticity, each sample was amplified and genotyped twice. In the subsequent analyses, the data from the duplicates were treated as individual samples.
2.3 Individual genotyping of the host
For the host genotyping, we used 23 individuals from low altitude, 15 individuals from middle altitude, and 27 individuals from high altitude, with a minimum of 5 samples per locality, on a total of 65 samples (Table 1).
The genotyping protocol utilized a combination of seven microsatellite markers sourced from (Popa et al., 2014) and six additional microsatellite markers from (Karanth et al., 2004). Consequently, a total of 13 variable microsatellite markers were employed for the genotyping analysis (see Supplementary Table S3). Multiplex PCR was performed with FAM and HEX fluorescent-labeled primers, using a QIAGEN Multiplex PCR Kit (QIAGEN, Cat. No: 206143). This enabled the amplification of multiple targets within a single reaction. Fragments were differentiated by the respective fluorescent labels as well as by their expected sizes.
PCR reactions were performed in PCR plates in 7 μL reaction volumes (3 μL of Multiplex PCR Buffer, 0.6 μL of Q-solution, 1.5 μl of DNA, 0.5 to 1 μl of primer, and distilled water). Notably, the HEX-labeled primers required higher concentrations (10 picomoles) compared to FAM markers (5 picomoles) for optimal performance.
The PCR temperature profile consisted of an initial denaturation step at 95°C for 15 min, followed by 10 cycles of denaturation at 93°C for 40 s, annealing at 60°C for 40 s (with a decrease of 0.4°C per cycle), and extension at 72°C for 80 s. This was followed by 20 cycles of denaturation at 93°C for 30 s, annealing at 56°C for 40 s, extension at 72°C for 80 s, and concluded with a final extension step at 60°C for 30 min.
For fragment analysis, we employed the GeneScan™ 500 LIZ™ dye Size Standard (Applied Biosystems™, Cat. No: 4322682). Initially, the size standard was diluted with Hi-Di™ Formamide (Cat. no. 4311320) at a ratio of 40 size standard to 1,000 formamide. Subsequently, 10 μL of the diluted standard was added to each well, followed by the addition of 2–3 μL of PCR product. The plate was subjected to a thermal cycler at 95°C for 3 min, followed by immediate chilling. Fragment analysis was performed using the 3130XL Genetic Analyzer (Applied Biosystems™), with allele identification carried out using Gene Mapper Software (version: 5.0).
2.4 Thyroxine levels quantification assays
Using the serum from four localities among 3 altitude groups (Table 1), we measured free fractions of thyroxine and triiodothyronine (fT4 and fT3) as well as total thyroxine and triiodothyronine (TT4 and TT3) using commercial ELISA kits according to manufacturer’s instructions (DRG Diagnostics, Marburg, Germany: FT4 - EIA 2386; FT3 - EIA 2385; TT4 - EIA 4568; TT3 - EIA 4569).
2.5 Bioinformatic analyses
Following demultiplexing and trimming of the raw 16 s rRNA and 18 s rRNA sequencing data using Skewer (Jiang et al., 2014) reads with low quality were eliminated by setting the expected error rate per paired-end read >1 (Jiang et al., 2014).
The quality-filtered reads were denoised with DADA2 software (Callahan et al., 2016), resulting in an abundance matrix containing the number of reads for each amplicon sequence variant (ASV) in each sample. The UCHIME software (Edgar et al., 2011) was employed to identify and remove sequence chimeras, with gold.fna database1 serving as a reference for chimera filtering. For 16S rRNA bacterial ASV annotation in the DADA2 software, the Silva database version 138.1 (updated in March 2021; Quast et al., 2013) was used as a reference. For the 18S rRNA dataset, for each ASV, the top 200 Blastn hits were downloaded from the NCBI nucleotide database (Camacho et al., 2023) and used to construct the custom reference database. The ASV taxonomy was then assigned using the RDP classifier as previously described (Quast et al., 2013). Finally, the phyloseq (McMurdie and Holmes, 2013) package was utilized to construct a comprehensive database containing the Amplicon Sequence Variants (ASVs) table, ASV sequences, taxonomic annotations at phylum and family and genus level (when possible), and phylogeny for both datasets.
2.6 Microbiome and diet
The microbiome database comprised 1,103,094 high-quality sequences grouped in 4841 non-chimeric ASVs. The number of sequences in the microbiome database per sample varied between 8,936 and 30,237. The diet database comprised 2,392,247 high-quality sequences grouped in 140 non-chimeric ASVs. The number of sequences in the diet database per sample varied between 515 and 205,099. For the rarefaction of the ASV table, “phyloseq_mult_raref_avg” function from the metagMisc package was used with 100 iterations, which provides robustness by applying repeated subsampling. Using the minimal sequencing depth as the rarefaction threshold, we ensured even sequencing depth per sample and utilized the down-sampled dataset for further analysis unless otherwise stated. All the statistical analyses were done in R version 4.2.2.2
The procrustes test (Procrustes Rotation of Two Configurations in vegan pack (Oksanen et al., 2024)) was used to compare duplicates and for both datasets procrustes test showed the composition was consistent across all duplicates (number of permutations = 999, procrustes sum of squares (m12 squared) = 0.00414, Correlation in a symmetric Procrustes rotation = 0.9979, p-value = 0.001).
We employed the exact observed number of ASVs, Shannon, and Simpson indices to estimate alpha diversity using the “estimate_richness” function in the phyloseq package (McMurdie and Holmes, 2013). To compare alpha diversities between altitudinal groups, we used the “wilcox.test” function (hereafter WT, stats package) (Bauer, 1972). Per-sample Shannon and Simpson diversity indices were used as response variables in the Generalized Linear Mixed Models with Gaussian distribution (hereafter GLMM, glmmTMB package; Brooks et al., 2017), with the altitude category (low, middle, and high) included as a fixed variable and sampling location included as a random variable.
As a measure of beta diversity, we employed the Bray–Curtis dissimilarity index, which focuses on relative abundances using the “distance” function with a specified “bray” method (phyloseq package). We visualized the between-sample divergence pattern using Principal Coordinate Analysis (PCoA). Furthermore, we applied PERMANOVA (Permutational Multivariate Analysis of Variance Using Distance Matrices, “adonis2” function from the vegan package) to test for differences in the gut microbiota composition between altitudes. Additionally, to take account of the possible effect of the sampling localities, we used the first two principal components of the Bray-Curtis PCoA as response variables in the GLMM analysis, with altitude category (low, middle, and high) included as a fixed variable and sampling location included as a random variable. Furthermore, the effect of fixed variables was tested using likelihood ratio tests (“anova” function, stats package). Finally, we employed the MDMR (Multivariate Distance Matrix Regression, mdmr package (McArtor et al., 2017)) with the Bray-Curtis distance matrix as a response variable, altitude category as a fixed variable, and sampling location as a random variable. Then, for both alpha and beta diversity measures, we tested the effects of sex, sampling year, and body mass using the same statistical approach. We checked for the possible correlation between microbiome and diet composition. First, the samples not present in both datasets were excluded. Then the Mantel Test (ape package; (Mantel, 1967; Paradis and Schliep, 2019)) was employed to check the correlation between diet and microbiome composition.
We used the “DA.kru” function in the DAtest package for calculating Kruskal-Wallis analysis (Russel88/DAtest on GitHub (Russel et al., 2018)) to compare the differential abundances of bacterial phyla and families among different altitudes. The p-values from “DA.kru” function were adjusted by “FDR” method (Benjamini and Hochberg, 1995) as default. To explore the relative abundances of ASVs that significantly differ among altitudinal groups or sampling localities, we used the DESeq2 package (Anders and Huber, 2010). Additionally, we specifically focused on ASVs that exhibit higher abundance in one location or a select few, while being nearly or entirely absent in others. To do that, we employed “estimateSizeFactors” function in the DESeq2 package with LTR test parameter (Likelihood ratio test (chi-square test) for GLMs).
Finally, we predicted high-level bacterial characteristics using BugBase4 (Ward et al., 2017) based on 16S rRNA data, following the online instructions. These phenotypes included Gram Positive, Gram Negative, Biofilm Forming, Pathogenic Potential, Mobile Element Containing, Oxygen Utilizing, Oxidative Stress Tolerant, and Facultatively Anaerobic bacteria. The relative frequencies of predicted phenotypes were then tested for altitudinal difference using Kruskal-Wallis Test (“kruskal.test” in stats package (Hollander et al., 2013)) and the p-values were corrected using “p.adjust” function in stats package with the FDR method.
2.7 Host genetics
The individual genotypes were handled in adegenet package (Jombart, 2008) as a “genind” object. We calculated Nei’s pairwise Fst (Nei, 1973) between all pairs of populations by “pairwise.neifst” function hierfstat package, (Goudet, 2005). Then, we calculated the mean proportions of shared alleles by using the “propShared” function (adegenet package), averaged per respective sampling localities. Then, allelic richness for each altitudinal group was calculated using “allelic.richness” function (hierfstat package) and compared with “t-test” function (stats package). Finally, we run AMOVA (Analysis of Molecular Variance) for altitudinal groups and sampling localities “using poppr.amova” function in poppr package (Kamvar et al., 2014), and generated the p-values using “randtest” function in ade4 package.
STRUCTURE software (V2.3.4, (Pritchard et al., 2000)) was used to reveal the population genetic structure using the Admixture Model with 5,000 burnin, 50,000 MCMC replicates after burnin, 5 iterations, and correlated allele frequencies. The analyses run with and without sampling locality as prior, and a number of inferred genetic clusters (K) is set to a range between 1 to 6. The resulting output is reviewed on the Structure Harvester website (v.0.6.94, (Earl and vonHoldt, 2012)).
Using the “coancestry” function (related package; Pew et al., 2015), we calculated the empirical relatedness between same-population individuals according to Li et al. (1993) to seek for possible correlation between kinship and microbiome composition. Only those samples present in both datasets were included. The pairwise matrix of genetic distances (measured as 1-relatedness) was computed using the “mat_gen_dist” function, utilizing the “PCA” method (graph4lg package; Savary et al., 2021). Subsequently, Euclidean distance matrices were generated for each population. The Bray–Curtis dissimilarity distance matrices were constructed on the individual microbiome data, and the Mantel Test (ape package; Mantel, 1967; Paradis and Schliep, 2019) was employed to check for the correlation between host kinship and microbiome composition. Additionally, “grouprel” function (related package) is employed to calculate average within-group relatedness to compare with expected relatedness in simulated populations with 100 iterations (ie. 100 parent-offspring, 100 full siblings, 100 half siblings, and 100 unrelated).
2.8 Thyroid hormone levels
Concentrations of fT4, TT4, fT3, and TT3 were analyzed with multiple linear regression models with altitude, sex, and weight as explanatory variables. All data were tested for normal distribution using three normality tests (Anderson-Darling, D’Agostino-Pearson omnibus, and Shapiro Wilk). Since fT4 concentrations were not normally distributed, fT4 concentrations were log-transformed before analyses. We first calculated one model with main effects (model 1) and one model with two-way interactions (model 2) to identify the best-fitting model. The Akaike information criterion (AIC) was used to estimate model fit. Based on the best-fitting model, explanatory variables with significant effects on hormone concentrations were identified. These analyses were conducted using GraphPad Prism (vers. 9.3.1, San Diego, CA, United States).
Effect of altitude on hormone levels was tested using ANOVA. Since this approach did not take into account the random effect of sampling altitude, we used GLMMs and MDMRs with alpha diversity metrics as response variables, hormone measurements as fixed, and sampling locality as a random variable. For the beta diversity, we employed GLMMs and MDMRs using two principal components of the Bray-Curtis PCoA as response variables, and hormone levels as fixed variables, and the sampling locality and altitude as a random effect. Finally, we use the DA.test package to see if there is correlation between hormone levels and differential abundance of a bacterial taxa.
3 Results
3.1 Gut microbiome
We successfully genotyped the 16S rRNA amplicons of a total of 52 ABMR caecum samples in duplicates. For subsequent analyses, the data from duplicates were merged based on the consistent coverage of ASVs between duplicates, meaning that only the ASVs found in both duplicates were retained.
The microbiome database was dominated by Firmicutes (53% of all reads), Bacteroidota (36%), and Desulfobacterota (0,7%). Another 9 phyla detected were: Patescibacteria, Cyanobacteria, Elusimicrobiota, Campylobacterota, Proteobacteria, Actinobacteriota, Synergistota, Halobacterota, and Thermoplasmatota had low abundances (<1% reads, Supplementary Figure S1). At the family level, the data was dominated by Muribaculaceae (35%), Lacnospiraceae (28%), Christensenellaceae (~5%), Ruminococcaceae (~5%), Desulfovibrionaceae (~5%), Oscillospiraceae (~5%) and 29 other bacterial families with less than 1% abundance (Supplementary Figure S1). The relative abundances of bacterial phyla were visualized for each sample, sampling location, and altitudinal group (Figures 2A,B) and individual microbiome variability shown in Supplementary Figure S2.
Figure 2. (A) The relative abundance of dominant bacterial phyla among the sampling localities and (B) altitudinal groups (X-axis represents % of the abundance of all reads). (C) Variation in the GM diversity (number of observed ASVs, Shannon, and Simpson diversity indices) across altitudes (p-values are given from WT analysis). (D) PCoA plot on Bray–Curtis dissimilarity metric showing the divergence between gut microbiota from altitudinal groups.
3.2 The Alpha and Beta diversity of microbiome
The alpha diversity of the gut microbiota in ABMR exhibited variations across different sampling locations and altitudinal groups. Specifically, the animals inhabiting the middle altitude demonstrated a more diverse composition of the microbiome compared to those at low and high altitudes (Figure 2C).
In the case of all three alpha diversity metrics (number of observed ASVs, Shannon and Simpson diversity indices), the WT analysis revealed that samples from the middle altitude group exhibited significantly higher levels of among-sample diversity compared to the other altitudinal groups. Moreover, the WT analysis indicated a significant difference between the middle and high altitudes, except for the Simpson diversity index (p-values for each test are shown in Figure 2C). GLMM analyses revealed that the likelihood ratio test comparing the two models was not significant for the Shannon index (ΔD.F = 2, χ2 = 4.6772, p-value = 0.09), and Simpson index (ΔD.F = 2, χ2 = 2.4964, p-value = 0.28), but significant for the number of observed ASVs (ΔD.F = 2, χ2 = 10.571, p-value = 0.005). While WT analyses revealed significantly higher bacterial diversity at middle altitudes using all three alpha diversity metrics, GLMM showed significance only for the number of observed ASVs. It is important to note that, unlike GLMM, WT analyses do not account for the effect of sampling locality.
The PERMANOVA results demonstrated significant differences between the altitudinal groups in terms of Bray-Curtis dissimilarities (F-value = 4.6839, df = 2, p-value = 0.001, and R2 = 0.16). Furthermore, we employed the first two axes of the Bray-Curtis-based PCoA as response variables in a GLMM analysis to assess the altitude’s effect. The first PCoA axis mostly corresponded to the elevation gradient, separating the three altitude categories. A notable exception to this trend was the high-altitude site of Kiziltepe, as indicated by the data points in the top right corner of Figure 2D, which clustered far apart from all other data points. At the same time, the second PCoA axis represented the elevation gradient more clearly, with data points at the bottom of Figure 2D corresponding to low-altitude samples and those at the top corresponding to high-altitude samples. In the GLMM analyses, the likelihood ratio test comparing the two models for the first PCoA axis was not significant (χ2 = 1.6172, df = 2, p = 0.44) which is most probably due to the Kiziltepe samples which clustered far from the rest, especially along the first axis. However, the comparison for the second PCoA axis indicated that the categorical model significantly enhanced the explanation of variation in the ABMR gut microbiota (χ2 = 9.9051, df = 2, p = 0.007). Finally, the results from MDMR analyses also showed a significant effect of altitude (df = 2, p-value =0.0006). All of the above patterns remained significant when the data collected in each year were analysed separately, but at the same time, the sampling year was itself a significant factor in models considering altitude and sampling localities as random variables (GLMM: χ2 = 10.423, df = 1, p-value = 0.001; MDMR: df = 1, p-value =0.00001). Neither GLMM for sex and Shannon diversity considering sampling locality and altitude as random variables (χ2 = 1.6125, df = 1, p-value = 0.20) nor GLMM for sex and beta diversity (i.e., the first PCoA axis; χ2 = 0.0132, df = 1, p-value = 0.90) were significant. Similarly, GLMM for host body mass and Shannon diversity considering altitude and sampling localities as random variables (Z-value = 0.542, SE = 0.00091, p-value = 0.55) and GLMM for host body mass and beta diversity (i.e., the first PCoA axis; χ2 = 1.581, df = 1, p-value = 0.21) found to be insignificant.
3.3 Differential ASV abundance testing
Anayses within the DA.test package revealed a significant effect of the altitude on the relative abundance of bacterial phyla and families. Specifically, the relative abundance of Bacteroidota increased with altitude, whereas Firmicutes and Desulfobacterota were significantly more abundant at low altitudes, and Verrucomicrobiota was only found at low altitudes (Table 2; Figure 2A).
At the bacterial family level, we observed that the relative abundances of Muribaculaceae, and Unclassified Bacteroidales, increased with altitude, while the relative abundances of Lacnospiraceae, and Ruminococcaceae decreased. Additionally, Christensenellaceae was significantly more abundant at high altitude but less abundant at middle altitude, and Akkermansiaceae was only found at low altitude in both localities (Table 3).
The DESeq2 analysis revealed a significant effect of the sampling locality on the relative abundances of 51 ASVs (Table 4; Supplementary Figure S3). Twenty four of them were overabundant in one or in a few localities but significantly underrepresented in all the other localities (Table 4). Among these, Kiziltepe exhibited the highest number of ASVs with significant differential abundance (n = 12), while both high-altitude localities collectively harbored 19 significant ASVs. Intriguingly, no such ASVs were detected in the middle altitude localities. These DESeq2 results correspond to alpha diversity analyses, which also suggest that the middle altitude samples comprise a blend of two other altitudinal groups.
3.4 Bacterial characteristics prediction
We used BugBase, a tool for measuring high-level phenotypes in the microbiome, to assess changes in bacterial phenotypes across different altitudinal groups (Figure 3). The ‘Facultatively Anaerobic’ phenotype was the only one that displayed a significant difference (Kruskal-Wallis fdr adjusted p-value = 0.049), attributed to its prevalence at the high altitude. Note that almost all ASVs that contributed to this phenotype belonged to the bacterial phylum Firmicutes, which itself was more abundant at the high altitude (Figure 2A).
Figure 3. Variation of predicted bacterial phenotypes among the altitudinal groups. Values in the y-axis represent the BugBase bacterial phenotype prediction.
3.5 Diet composition
We successfully genotyped the 18S rRNA amplicons from 38 ABMR caecum samples (12 low, 14 middle, and 12 high altitude) in duplicates. Similar to the 16S results, the data from duplicates were merged based on the same coverage of ASVs between duplicates, and the data from duplicates were treated as individual samples.
The raw dataset was clustered in four phyla: Arthropoda, Basidiomycota, Chordata, and Streptophyta. However, considering the ABMRs are herbivorous rodents (Sözen, 2005), the database was filtered by excluding the host genome and other possible contamination-caused sequences from non-plant phyla, retaining only those belonging to plants (Streptophyta). This resulted in approximately two-thirds of the reads being retained for further analysis. At the order level, the data was dominated by Asterales (42.7%) and Brassicales (21%), Apiales (15.21%), Fabales (10.37%), Asparagales (7.92%), Rosales (6.61%), and as well as 12 other plant orders with less than 5% relative abundance. Due to the relatively short sequence lengths of our 18S rRNA libraries which become increasingly unreliable when assigning low-level taxa and lack of complete reference database focused on plant 18S rRNA, we only resolved the diet taxonomy down to the order level.
The relative abundances of plant orders across altitudinal groups and sampling locations are illustrated in Figures 4A,B. Asterales and Fabales were prevalent across all groups, whereas Apiales were exclusively observed at low and middle altitudes, with significantly higher abundance at the middle altitude (Kruskal-Wallis fdr adjusted p-value = 0.05). In contrast, Asparagales were observed exclusively at the high altitude, consistent with the DA.test results (Kruskal-Wallis fdr adjusted p-value <0.0001).
Figure 4. (A) The relative abundance of dominant plant orders among the sampling localities and (B) altitudinal groups (X-axis represents % of the abundance of all reads). (C) Box-plot showing variation in diet diversity (number of observed ASVs, Shannon, and Simpson diversity indices) between altitudes. Both WT and GLMM showed no significant difference between altitudes. (D) PCoA of Bray–Curtis dissimilarity in diet composition among the individuals, altitudes, and localities.
3.6 The Alpha and Beta diversity of diet
The observed number of ASVs, Shannon, and Simpson diversity measures did not significantly differ by altitude (Figure 4C). This result was supported by both WT and GLMMs (both p-values >0.05). Nevertheless, all three indices displayed the widest observable range at the middle altitude (Figure 4). After examining the variation in Bray–Curtis dissimilarity index of the diet composition with PCoA, PERMANOVA (p-value = 0.001, R2 = 0.14, and F-value = 2.8517), the GLMM for the first PCoA axis (χ2 = 7.6772, df = 2, p = 0.021), and MDMR for PCoA axes (df = 2, p-value = 0.00001) showed a significant effect of the altitude (Figure 4D). In the PCoA, the high-altitude animals clustered separately, while the low and the middle altitudes overlapped to some extent. There was no significant effect of the sex, body mass, and sampling year on the diet of ABMR. The Mantel Test indicated no significant correlation between the diet and the microbiome composition (p-value = 0.20, r = 0.03).
3.7 Host genetics
A total of 65 samples were genotyped at the six loci reported in (Karanth et al., 2004), with 6.81% missing data, but only the samples collected in 2019 were genotyped for the seven loci reported in (Popa et al., 2014), with 25% missing data. The mean allelic richness values calculated in low, middle, and high altitude groups were 1.636348, 1.735554, and 1.603827, respectively. There was no significant difference in mean allelic richness among the altitudinal groups (t-test p-values >0.05). The mean allelic richness at the sampling locality level ranged between 1.352573 (Kiziltepe) and 1.721504 (Madenkoy), and the only significant difference was observed between Eregli and Kiziltepe (t-test p-value = 0.01). The pairwise Fst and proportion of shared alleles between populations are shown in Table 5.
Table 4 presents the mean proportions of shared alleles and the pairwise Fst among populations. The shared alleles proportions ranged from 38% (Eregli and Kiziltepe) to 52% (Madenkoy and Darbogaz), while the Fst varied from 0.025 (Karagol and Madenkoy) to 0.25 (Karagol and Darbogaz). Notably, the Fst between Uluskisla and Madenkoy was negative (−0.025), likely due to the high proportion of missing data in the Madenkoy samples. Altitudinal genetic differences in Fst were 0.04 (low to middle), 0.06 (middle to high), and 0.11 (low to high). Notably, the confidence intervals for the Fst estimates were very wide, spanning zero in all cases except for four values: Darboğaz-Eregli = 0.056 (0.014–0.095), Karagol-Eregli = 0.112 (0.047–0.199), Karagol-Darbogaz = 0.233 (0.020–0.445), Kiziltepe-Karagol = 0.118 (0.042–0.188).
The MANOVA revealed the only significant variation between the sampling localities at the level of 11.3% (p-value = 0.01), but the variation among altitudes (6%, p-value = 0.49), and among samples within the same altitude (1%, p-value = 0.76) were non-significant. Finally, the two runs of STRUCTURE (with and without the prior) both indicated the best model at K = 2 with the Delta K values of 44.840244 and 121.992245, respectively. These results suggest that while our dataset could potentially represent an admixture of two genetically distinct clusters (see Supplementary Table S4; Supplementary Figure S5), this weak genetic structure does not reflect the altitudinal variation.
Our analysis of genetic relatedness at the individual level using related (Pew et al., 2015) revealed that the average observed within-group relatedness across all localities was significantly higher compared to the expected values (p-value = 0.01). At the level of locality, Eregli, Darbogaz, and Kiziltepe showed higher than expected relatedness, while Ulukisla, Madenkoy, and Karagol showed lower than expected relatedness (Supplementary Figure S6). The Mantel tests between the microbiome composition and the individual genetic relatedness in each sampling locality showed no significant correlation (r < 0.09, p-values >0.05).
3.8 Thyroxine hormone levels
To assess the effects of altitude on thyroid hormone levels and to account for additional factors known to influence thyroid hormone concentrations, we constructed multiple linear regression models with altitude, sex, and the body mass as explanatory variables. Altitude as a main effect exerted a significant effect on fT4 concentrations (F = 6.38, p-value = 0.021), while the sex had no effect (F = 0.0001, p-value = 0.98) and the body mass only had a marginal effect (F = 4.06, p-value = 0.059) on fT4. In contrast, none of the given factors had a significant effect on fT3. The main effects did not affect TT4 and TT3 concentrations, but the interaction of sex and altitude exerted a statistically significant effect on total hormone concentrations (Supplementary Table S5). Furthermore, the fT4/TT4 ratio was significantly lower in animals living at low altitudes compared to middle altitudes and showed a trend toward a significantly lower ratio compared to high altitudes (p-value = 0.09, Figure 5).
ANOVA and GLMMs revealed that there is a marginally significant effect of TT4 hormone levels on the Shannon and Simpson alpha diversity metrics. Specifically, TT4 levels showed negative correlation with alpha diversity metrics. The PERMANOVA performed on the Bray-Curtis distance matrix and GLMMs applied to the first axis of PCoA revealed that there is a significant effect of the FT4/TT4 levels on the microbiome beta diversity. However, while the correlation of the Bray-Curtis distance matrix correlated with the FT4/TT4 ratio positively, the PCoA axis showed negative correlation (Table 6; see Supplementary Table S6). Finally, the DA.test results revealed no significant correlation with the relative abundance of bacterial taxa and the levels of thyroid hormones.
4 Discussion
Our study is the first to address the multiple aspects of altitude adaptation in the subterranean Blind Mole Rat (Nannospalax xanthodon) - a biomedically and evolutionarily important model organism. We found that altitude is indeed a key factor shaping the variation observed in (i) microbiome, (ii) diet, and (iii) thyroid hormone levels, whereas (iv) its population genetic structure reflects a fine-scale geographic pattern rather than the elevation gradient. We could not determine statistically whether the biological aspects (i-iv) also interact in a meaningful way among themselves, or if (i-iii) simply co-vary with the altitude. Below, we take advantage of the existing knowledge of the complex nature of the microbiome functioning and mammalian adaptation to high altitudes, in order to hypothesize on the possible mechanisms behind our findings.
4.1 Changes in bacterial taxa and phenotypes among the altitudinal groups
The ratio of Firmicutes/Bacteroidota (F/B) is linked to various host health, disease traits, and also linked with altitude according to most studies (Indiani et al., 2018; Jasirwan et al., 2021; Koliada et al., 2017; Mariat et al., 2009; An et al., 2023). A study on Rhesus macaques indicated that animals inhabiting high altitudes exhibit a higher abundance of Bacteroidota and a lower abundance of Firmicutes compared to those at lower altitudes (Wu et al., 2020). Another study examining the composition of the oral microbiome in Tibetans living at high (2800–3,650 m) and ultra-high (3650–4,500 m) altitudes reported a low F/B ratio in the populations at the highest altitude (Liu et al., 2021). This observation was speculated to be an adaptation to maintain normal blood pressure in hosts, considering the decrease in ambient oxygen and maximal oxygen consumption (VO2max) of the human body with altitude (Squires and Buskirk, 1982), and lower F/B ratio has been associated with a lower VO2max in another study (Durk et al., 2019). We observed the same pattern in AMBR (see Table 2; Figure 2A), which could suggest a role of the low F/B ratio in blood pressure maintenance. At the same time, the extreme specialization of the AMBR to subterranean lifestyle may create a very different context compared to the above listed examples. The ABMR responds to the hypoxia caused by low levels of O2 inside the tunnels rather than to its ambient atmospheric concentration, with the additional factors such as tunnel depth, soil type and permeability, animal own activity, etc. (see the respiratory stress hypothesis in (Arieli, 1979; Darden, 1972)) playing a role. Importantly, some studies have presented conflicting results regarding the link between the F/B ratio and altitude (ie., Li and Zhao, 2015; Zeng et al., 2017).
An increase in the abundance of strictly anaerobic bacteria at higher altitudes is a well-known phenomenon, for example in wild house mice (Suzuki et al., 2018) and humans (Lan et al., 2017), although such effects were mostly observed at much greater elevations (i.e., up to 4,500 m a.s.l.) compared to our study (the highest location was at 2900 m a.s.l.). The previously mentioned factors specific to the subterranean environment and AMBR lifestyle may also explain the lack of significant difference in the abundance of strictly anaerobic bacteria among our altitudinal groups. Even more striking is the complete absence from our 16S rRNA dataset of the anaerobic bacterial genus Prevotella, which was previously found to be strongly associated with increasing altitude in wild house mice (Suzuki et al., 2019), wild pikas (Li et al., 2016b), high altitude yaks, in Tibetan sheep (Zhang et al., 2016), and in humans (Li et al., 2016d; Lan et al., 2017), but (see Li and Zhao, 2015). The Prevotella is also present at ~3.5% in subterranean plateau zokors - another member of the family Spalacidae living at high altitudes (Hu et al., 2023). Our results as well as the previous studies suggest that the absence of Prevotella in Blind Mole Rats GM may actually be a common feature of the wild populations of these rodents: i.e. it was detected neither in the wild-caught N. xanthodon [(Kuang et al., 2022; Solak et al., 2023), nor in N. ehrenbergi (Kuang et al., 2022; Solak et al., 2023)], but Sibai et al. (2020) reported the first appearance of Prevotella in the fecal microbiomes on the AMBR only after they spent a month in captivity. Therefore, the importance of this bacterial taxon in the adaptation to high altitude does not seem to apply to our research system.
We showed that unlike the strictly anaerobic bacteria, the facultatively anaerobic ones were still significantly more abundant at high altitudes (Figure 3, Supplementary Figure S4). This functional group has a unique ability to grow in the presence or in the absence of oxygen and is thus well-adapted to changing environments (André et al., 2021). It is possible that extreme fluctuation of temperature and precipitation regime (including dense snow cover in winter) at high altitude also translates into a wide range of the ambient O2 concentrations experienced by the AMBR, and thus favors the facultatively anaerobic microbial elements.
The Ruminococcaceae bacterial family was found to be more prevalent during periods of limited energy availability in black howler monkeys, potentially compensating for reduced energy intake (Amato et al., 2015). This finding was further supported by (Wu et al., 2020), who reported an increased abundance of Ruminococcaceae in rhesus macaques inhabiting cold, high-altitude environments, suggesting a role in energy-saving. Interstingly, our results revealed not an increase but a decline in the relative abundance of Ruminococcaceae with increasing altitude (Table 3; Figure 2). This result could be attributed to the specifics of biology of the AMBR, in particular to its reduced tolerance of the heat stress (Šumbera et al., 2023). During the summer season when all our samples were collected, the animals at the low altitudes have to deal with much hotter and drier environment compared to the milder and wetter conditions at the higher elevation. This would likely cause a greater cumulative stress, causing the animals to reduce their activity and so exploit any available opportunity for greater energy efficiency, including the adjustment in the microbiome.
Muribaculaceae was reported to be linked with immune function (Sha et al., 2022; Huang et al., 2022) and to be more abundant in hypobaric hypoxic rats (Ma et al., 2023). A study on the plateau pikas found that Muribaculaceae is more abundant in the warm season and at high altitudes (Tang et al., 2023). Similarly, we observed a significant increase in the abundance of Muribaculaceae with increasing altitude, while all our samples were collected during the warmest time of the year. Previously, we showed that the high-altitude ABMRs from the same study area possess a greater innate immune response (Solak et al., 2020). We hypothesize that the higher abundance of Muribaculaceae could have contributed to a stronger immune function in high-altitude ABMRs. We note that apart from Muribaculaceae, the Akkermansiaceae (any many other bacteria) family has also been linked with immunoregulation and we only found it at low altitudes (Portincasa et al., 2024; Cani et al., 2022; Grenda et al., 2022). However, different bacterial species can stimulate different immune responses by activating different immune cells and pathways (Acosta and Alonzo, 2023; Zheng et al., 2020).
4.2 The effect of altitude on microbiome diversity and composition
In the context of laboratory rodent studies simulating high-altitude conditions, the relationship between GM diversity and altitude has yielded inconsistent results. A study on rats by (Tian et al., 2018) reported significantly higher GM alpha diversity at high altitudes, whereas (Zhang et al., 2018) observed slightly lower diversity (albeit not statistically significant). Higher GM alpha diversity at high altitudes was found in wild pikas (Li et al., 2019) and rhesus macaques (Wu et al., 2020; Zhao et al., 2018, 2023). Another study on wild mice ranging from sea level to 4,000 meters asl. Found a slight increase in diversity at high altitudes, although not statistically significant (Suzuki et al., 2018). It should be noted that none of these studies (except Suzuki et al., 2018) were conducted in continuous transects, i.e., the effect of latitude over altitude is unknown. The effect of altitude on the microbiome in the wild animal system thus remains understudied. Consequently, each study was only able to draw conclusions based on the specific biology of the host species and the characteristics of the study system.
In contrast to the above examples, our study design aimed at reducing the possible effect of the local geography and seasonal variation in the AMBR system within its natural habitat. We detected a very strong but non-linear effect of the altitude: the highest alpha-diversity was observed in both populations from the middle altitude group, while the low and high altitudes did not differ significantly (Figure 2C). The reason behind this is that the middle altitude in fact combines all the specific (i.e., significantly over-abundant) ASVs from both low and high altitudes (Table 4). A straightforward interpretation for such a pattern would be the existence of two polarized assemblies of the microbial taxa, at low and high elevations respectively, with the middle altitude serving as a bridge between them. Recall that the obligate subterranean lifestyle of the AMBR implies very slow and gradual dispersal (Rado et al., 1992). If we then assume that the transfer of the GM elements co-occurs with the host, the direct exchange between the low and the high altitudes would seem difficult: i.e. there is no direct pathway for the lowland bacteria to reach highland habitats (and vice versa) without first passing through the middle altitude. In addition, the middle altitude belt is expected to have relatively mild environmental conditions (i.e., less extreme fluctuation in temperature and humidity), which could also favor the higher microbiome diversity at the individual level. Verification of this intriguing hypothesis will require collecting additional data, accounting for multiple potential ecological factors that could act to ‘polarize’ the high and low GM communities.
The altitude is known to exert significant effect on microbial composition at the population level (beta diversity), as reported in previous studies on rodents (Li et al., 2019; Li et al., 2016a; Li et al., 2016b; Suzuki et al., 2018) and humans (Zeng et al., 2020; Lan et al., 2017; Li and Zhao, 2015). Most studies explained this effect as an indirect result of other altitude-related processes, such as variation in the host diet, and in host genotypes. In most cases it was difficult to eliminate the effect of geography / spatial variation while testing for the effect of the altitude. In our study, we attempted to minimize the effect of geography by collecting samples within a small area featuring a steep elevation gradient. A very strong altitude effect on the GM composition revealed was mostly caused by a single high-altitude locality (Kiziltepe) clustering far from the rest of the data (i.e., PCoA plot on Figure 2D). Further analysis using DESeq2 revealed that this particular population harbored the highest number (n = 15) of significantly overrepresented ASVs - at least three times that of the second high-altitude locality Karagöl (n = 4), which is only a few hundred meters lower in terms of absolute elevation (Table 1). The population residing at Kiziltepe had the lowest genetic diversity, but the host genetics is an unlikely factor to explain the microbiome difference there, given the overall lack of correlation between the GM and population genetic structure. Neither could we find any significant difference in diet composition between Kiziltepe and Karagöl, which undermines the role of diet on the microbiome composition, at least in the case of AMBR. In addition, we note that Kiziltepe is situated next to an open mining area active since 1825, primarily extracting minerals such as iron, lead, silver, and gold (Kahya, 2018; Anatolia, 2023). Previous studies have indicated a significant impact of heavy metals on the gut microbiota of mice (Breton et al., 2013). In the absence of other clues, we can suggest that the distinct composition of the gut microbiota observed in Kiziltepe could potentially be attributed to the influence of heavy metal exposure as well as the adaptation to high altitude.
4.3 The effect of diet and host genetics on microbiome composition
The relationship between microbiome and host diet is discussed in several studies (Scott et al., 2013; Conlon and Bird, 2014; Beam et al., 2021; Leeming et al., 2021) with majority of authors finding significant effect of diet in GM. Even though we found significant changes in microbiome (Figure 2) and diet composition (Figure 4) among the altitudes, we did not find any correlation between the GM at the level of ABMR sampling localities and their respective diet composition. Despite the fact that we found no significant link between the altitude and the plant taxonomic diversity in the ABMR diet, the composition of the latter was still affected by the altitude. This is generally expected, given a major effect of the altitude on vegetation. Seasonal, geographical or even cultural variation in the diet has been shown to have a strong influence on microbiome (Amato et al., 2015; Guo et al., 2021); (De Filippo et al., 2010), and we attempted to minimize such effects by performing sampling during the same season and within the same mountain range. Notably, the plant order Asparagales, which includes species with characteristic bulbous high-nutrient content roots, was only present in the AMBR diet at high-altitude localities (Figures 4A,B). The ABMRs are known to prefer the bulbous rhizomes and collect them for storage in their tunnels (Sözen, 2005). It is possible that the ABMRs, being food generalists, also possess a ‘generalist’ gut microbiota equally suitable to help with the digestion of different plant species.
A review on effect of host genetic control over microbiome reported that over 110 different genetic loci were found to be associated with the abundance of specific gut microbes (Bubier et al., 2021) and other studies found significant link between certain genes and the abundance of bacterial taxa in the gut (Org et al., 2015; Zhernakova et al., 2024). In addition, the microbiomes usually remain species-specific despite similar diets and shared habitat (Kaneko et al., 2023). We found no systematic effect of altitude on genetic differentiation, i.e., two distant populations at low altitude differed much more than any populations across the elevation gradient. This result may suggest a relatively recent colonization of the higher elevations in our study area. It was reported that within the same population, the family members or the ones with closer social network often share similar microbiota compared to unrelated individuals (Lee et al., 2011; Yatsunenko et al., 2012; Tims et al., 2013; Raulo et al., 2021). Another study found differences in specific microbial taxa between the social and solitary hyena species (Münger et al., 2018). Unlike in the highly social animals, the direct horizontal transmission of bacteria in the BMR is expected to be more difficult due to their solitary lifestyle. In line with this expectation, no correlation was observed between the host relatedness and microbiome in our data.
Several animal and human studies have found that sex significantly affects GM diversity and composition, though this pattern is not always consistent (Kim et al., 2020). We did not observe any significant differences in GM diversity and composition between female and male animals. This inconsistency extends to the diversity and composition of the diet. Previous studies have discussed the relative strength of sex’s impact compared to other factors, such as host genetic background, age, and diet, with sex generally found to have lower influence on shaping the GM (Kovacs et al., 2011; Elderman et al., 2018; Org et al., 2016). Our findings support the notion that the effect of sex on GM is not universal and environmental factors, host diet and genetics have a stronger effect on GM.
4.4 Thyroid hormone levels
Previously it has been reported that thyroid hormone levels in Middle East blind mole-rats are influenced by climatic conditions (Avivi et al., 2014). Other thyroid hormone parameters were found to be correlated with microbiome diversity in humans (Zheng et al., 2020; Virili et al., 2024). We found that the free thyroxine (fT4) levels increased significantly with altitude after correcting for sex and weight (Figure 5; Supplementary Table S5), while the fT4/TT4 ratio decreased in animals at low altitude. These findings suggest that in the animals living at lower elevations, a lower proportion of total T4 resources is recruited into the biologically-active free fraction. As thyroid hormones are positive regulators of thermogenesis, lower recruitment of T4 into the free form in animals living at lower altitudes might be linked to warmer ambient temperatures (Yau and Yen, 2020; Gerhardt et al., 2023), which reduces thermoregulatory constraints. Overall, these results point toward higher metabolic rate at high altitudes (Supplementary Table S5). As an indicator of metabolic rate, the fT4 levels might also explain the higher abundance of facultatively anaerobic bacteria at high altitudes (see above).
The microbiota can affect thyroid hormone levels by regulating iodine uptake, influencing the availability of essential micronutrients like selenium and zinc, and impacting the absorption of thyroid medications and overall thyroid function (Fröhlich and Wahl, 2019); (Knezevic et al., 2020). We did not find any correlation between the thyroid hormone levels and the relative abundance of bacterial taxa. At the same time, all of the alpha diversity metrics were found to correlate with TT4 levels (Table 6). Specifically, there was negative correlation with TT4 and a positive correlation with fT4/TT4, suggesting that the recruitment of T4 into its biologically-active form might contribute to an increase in the alpha diversity metrics. This is also in line with the negative correlation between TT4 levels and alpha diversity metrics, because circulating fT4 suppresses thyroid hormone synthesis and secretion via a negative feedback loop (Ortiga-Carvalho et al., 2016). Furthermore, PERMANOVA and GLMM on the first principle component of PCoA of Bray-Curtis similarities were found to be positively associated with fT4/TT4 (Table 6), similar to alpha diversity metrics. At this stage, it is not possible to establish whether there is a causal relationship between thyroid hormones and microbiome composition, or if both respond to the altitude in an independent manner. The relationship between metabolism and the GM is well-established (Martin et al., 2019), in particular between the metabolic rate and thyroid hormones (Mullur et al., 2014). In the case of ABMR, facilitating the nutrient breakdown could be necessary to fuel the thyroid hormone-driven thermogenesis - which in turn would put direct selective pressure on the GM landscape at higher altitudes.
5 Conclusion
This study explored the effect of altitude on the gut microbiome composition of ABMRs, across six localities and three altitudinal categories, considering the factors such as diet, thyroid hormone levels, and host genetics. Notable differences in the relative abundance of a number of bacterial taxa at different altitudes may reflect the specific roles these bacteria play in the complex adaptation of the host to the challenging mountain environment. The fact that the abundance of strictly anaerobic bacteria was unaffected by the altitude may reflect the relatively narrow range of absolute elevations in our study (1000–3,000 m a.s.l.), or it may be a specific feature of the host: the latter possibility is partly supported by the fact that the previously identified high altitude-linked genus Prevotella was absent from our data, in line with the previous studies done on the BMR. At the same time, the facultatively anaerobic bacteria were more prevalent in the high-altitude host specimens. We showed that the microbiome alpha-diversity reached its peak at the middle altitude, and that it incorporated elements from both lower and higher elevations. The beta-diversity correlated positively with the altitude. The altitude also affected the host diet composition, but not its alpha-diversity. These observations are unlikely to be caused by the host genetics, since we detected no clear association between the population genetic structure and the altitude, nor there was any correlation between the host relatedness and the microbiome composition nor diet. Thyroid hormone levels, specifically free thyroxine (FT4), increased almost linearly with altitude; however, no specific associations were found between bacterial ASVs and hormone levels. At the same time, the total thyroxine (TT4) levels did show a positive correlation with microbiome diversity. While some correlations between thyroid hormone components and microbiome beta diversity were identified, the nature of these relationships remains unclear.
5.1 Perspectives
Several properties of our study design could have unavoidably acted as caveats of the results presented. In particular, we tried to reduce the effect of seasonal variation on microbiome by sampling in late June to mid-July during both years when the fieldwork was performed. At this time of year, the high altitudes experience humid and colder conditions, but also feature a lusher vegetation compared to the low altitudes, which are both drier and hotter. The effect of season on the microbiome composition and diversity is well-known (Davenport et al., 2014; Hu et al., 2018; Jiang et al., 2021; Fan et al., 2022). While the results presented here offer a comprehensive snapshot of the microbiome, diet and thyroid hormones levels variation during the summer months, we cannot exclude a very different pattern if the samples were collected during other seasons. For example, the conditions at the lower elevation in the Taurus mountains during early spring resemble those found in mid-summer at high altitudes, which may profoundly influence the AMBR activity pattern as well as diet - which in turn may cause a major shift in the data collected across the entire elevation gradient. In the future, it would be preferable to have samples from other seasons in order to gain better insight into the microbiome and physiological dynamics in this system.
Rather than exclusively concentrating on microbial diversity and composition, future studies could also explore the functional profiling of bacterial ASVs specific to each altitude category using high-resolution metabarcoding data (i.e., full 16S rRNA sequencing), or full microbial genomes. Field experiments to measure the in situ metabolic rate of the animals might also help to reveal the interplay between microbiome and host physiology. Additionally, considering the distinct microbiome composition we found at Kiziltepe, further investigation could involve comparing the soil microbiome with the host microbiome to understand potential correlations. For dietary composition analysis, more specific genetic markers that focus on plants (e.g., rbcL, trnL, ITS) could offer higher taxonomic assignment resolution. High-resolution genomic data (e.g., ddRAD-seq or WGS) might unveil more subtle genetic divergence and relatedness patterms among the sampling localities and altitudes. Lastly, it would be extremely interesting to check if the patterns we found in AMBR will also be observed in other terrestrial animal species residing at different elevations in the Taurus mountains of Anatolia.
Data availability statement
16S, 18S rRNA phyloseq databases, and thyroid hormone datasets are available at FigShare (DOI: 10.6084/m9.figshare.25982311), and the microsatellite dataset is available in the Supplementary Tables. The 16S and 18S rRNA raw reads available at the European Nucleotide Archive under the project number of the study PRJEB76873 and the accession numbers for the samples are ERS20304461 - ERS20305352. The R code used to perform the analyses is available at GitHub (https://github.com/mrsolak/nannospalax_altitude).
Ethics statement
The animal study was approved by Animal Ethics Committee of Bülent Ecevit University. The study was conducted in accordance with the local legislation and institutional requirements.
Author contributions
HMS: Resources, Data curation, Formal analysis, Investigation, Methodology, Visualization, Writing – original draft, Writing – review & editing. JK: Formal analysis, Funding acquisition, Methodology, Project administration, Resources, Supervision, Validation, Writing – review & editing. DC: Formal analysis, Funding acquisition, Methodology, Project administration, Resources, Supervision, Validation, Writing – review & editing. ES: Funding acquisition, Project administration, Supervision, Writing – review & editing. LS: Investigation, Data Curation, Formal analysis, Writing – review & editing. MM: Investigation, Data curation, Writing – review & editing. YH: Formal analysis, Funding acquisition, Methodology, Writing – review & editing. FÇ: Resources, Writing – review & editing. FM: Resources, Writing – review & editing. AY: Conceptualization, Project administration, Funding acquisition, Resources, Investigation, Data curation, Methodology, Formal analysis, Resources, Supervision, Validation, Writing – original draft, Writing – review & editing.
Funding
The author(s) declare that financial support was received for the research, authorship, and/or publication of this article. The fieldwork and the laboratory analysis were supported by TÜBITAK grants 117Z596 and 220Z032 (to AY), respectively, while DC and JK were supported by the Czech Science Foundation (19-19307S). This research has been supported by the Ministry of Education, Youth and Sports of the Czech Republic grant talking microbes - understanding microbial interactions within One Health framework (CZ.02.01.01/00/22_008/0004597). The microsatellite genotyping was supported by Bülent Ecevit University (Grant Number: 2019-YKD-84906727-01 to AY). YH was supported by the German Research Foundation (Deutsche Forschungsgemeinschaft, grant number HE 7661/1–1) and the IFORES program of the Faculty of Medicine, University of Duisburg-Essen. The visits of HMS to the Institute of Vertebrate Biology Czech Academy of Sciences and to Charles University were supported by the EMBO Scientific Exchange Grant (#9427) and WAME Research Exchange Scheme (funded by the European Society for Evolutionary Biology), respectively. The numerical calculations reported in this paper were partially performed at TUBITAK ULAKBIM, High Performance and Grid Computing Center (TRUBA).
Acknowledgments
We are grateful to two reviewers for their valuable insights and constructive feedback on the manuscript. We are grateful to Gerald Parolly for his insights on the interpretation of the 18S rRNA diet dataset and to İshak Özel Tekin for helping with the use of ELISA instruments. HMS, AY, FM, and FÇ thank the staff of Pala Remzi restaurant in Çiftehan for helping us to recover after daily fieldwork labor in the Taurus mountains. This study is a part of the Ph.D. dissertation by HMS, supervised by AY and JK.
Conflict of interest
The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.
Publisher’s note
All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.
Supplementary material
The Supplementary material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fmicb.2024.1476845/full#supplementary-material
Footnotes
1. ^Available at: https://drive5.com/uchime/gold.fa
References
Acosta, I. C., and Alonzo, F. (2023). The intersection between bacterial metabolism and innate immunity. J. Innate Immun. 15, 782–803. doi: 10.1159/000534872
Ahn, J., and Hayes, R. B. (2021). Environmental influences on the human microbiome and implications for noncommunicable disease. Annu. Rev. Public Health 42, 277–292. doi: 10.1146/annurev-publhealth-012420-105020
Alberdi, A., Martin Bideguren, G., and Aizpurua, O. (2021). Diversity and compositional changes in the gut microbiota of wild and captive vertebrates: a meta-analysis. Sci. Rep. 11:22660. doi: 10.1038/s41598-021-02015-6
Amato, K. R., Leigh, S. R., Kent, A., Mackie, R. I., Yeoman, C. J., Stumpf, R. M., et al. (2015). The gut microbiota appears to compensate for seasonal diet variation in the wild black howler monkey (Alouatta pigra). Microb. Ecol. 69, 434–443. doi: 10.1007/s00248-014-0554-7
Amato, K. R., Yeoman, C. J., Kent, A., Righini, N., Carbonero, F., Estrada, A., et al. (2013). Habitat degradation impacts black howler monkey (Alouatta pigra) gastrointestinal microbiomes. ISME J. 7, 1344–1353. doi: 10.1038/ismej.2013.16
An, J., Kwon, H., and Kim, Y. J. (2023). The firmicutes/bacteroidetes ratio as a risk factor of breast cancer. J. Clin. Med. 12:216. doi: 10.3390/jcm12062216
Anatolia. (2023). Library of Congress. Available at: https://www.loc.gov/item/a22000917 (Accessed June 6, 2023).
Anders, S., and Huber, W. (2010). Differential expression analysis for sequence count data. Genome Biol. 11:R106. doi: 10.1186/gb-2010-11-10-r106
André, A. C., Debande, L., and Marteyn, B. S. (2021). The selective advantage of facultative anaerobes relies on their unique ability to cope with changing oxygen levels during infection. Cell. Microbiol. 23:e13338. doi: 10.1111/cmi.13338
Arieli, R. (1979). The atmospheric environment of the fossorial mole rat (Spalax ehrenbergi): effects of season, soil texture, rain, temperature and activity. Comp. Biochem. Physiol. A Physiol. 63, 569–575. doi: 10.1016/0300-9629(79)90197-X
Arslan, A., Kryštufek, B., Matur, F., and Zima, J. (2016). Review of chromosome races in blind mole rats (Spalax and Nannospalax). Folia Zool. 65, 249–301. doi: 10.25225/fozo.v65.i4.a1.2016
Avivi, A., Nevo, E., Cohen, K., Sotnichenko, N., Hercbergs, A., Band, M., et al. (2014). They live in the land down under: thyroid function and basal metabolic rate in the blind mole rat, Spalax. Endocr. Res. 39, 79–84. doi: 10.3109/07435800.2013.833216
Bäckhed, F., Manchester, J. K., Semenkovich, C. F., and Gordon, J. I. (2007). Mechanisms underlying the resistance to diet-induced obesity in germ-free mice. Proc. Natl. Acad. Sci. USA 104, 979–984. doi: 10.1073/pnas.0605374104
Bauer, D. F. (1972). Constructing confidence sets using rank statistics. J. Am. Stat. Assoc. 67, 687–690. doi: 10.1080/01621459.1972.10481279
Beam, A., Clinger, E., and Hao, L. (2021). Effect of diet and dietary components on the composition of the gut microbiota. Nutrients 13:795. doi: 10.3390/nu13082795
Belheouane, M., Vallier, M., Čepić, A., Chung, C. J., Ibrahim, S., and Baines, J. F. (2020). Assessing similarities and disparities in the skin microbiota between wild and laboratory populations of house mice. ISME J. 14, 2367–2380. doi: 10.1038/s41396-020-0690-7
Benjamini, Y., and Hochberg, Y. (1995). Controlling the false discovery rate: a practical and powerful approach to multiple testing. J. R. Stat. Soc. Ser. B Stat. Methodol. 57, 289–300. doi: 10.1111/j.2517-6161.1995.tb02031.x
Bisschop, K., Kortenbosch, H. H., van Eldijk, T. J. B., Mallon, C. A., Salles, J. F., Bonte, D., et al. (2022). Microbiome heritability and its role in adaptation of hosts to novel resources. Front. Microbiol. 13:703183. doi: 10.3389/fmicb.2022.703183
Breton, J., Massart, S., Vandamme, P., De Brandt, E., Pot, B., and Foligné, B. (2013). Ecotoxicology inside the gut: impact of heavy metals on the mouse microbiome. BMC Pharmacol. Toxicol. 14:62. doi: 10.1186/2050-6511-14-62
Brooks, M. E., Kristensen, K., Benthem, K. J., Van, M. A., Berg, C. W., Nielsen, A., et al. (2017). glmmTMB balances speed and flexibility among packages for zero-inflated generalized linear mixed modeling. R J. 9:378. doi: 10.32614/RJ-2017-066
Bubier, J. A., Chesler, E. J., and Weinstock, G. M. (2021). Host genetic control of gut microbiome composition. Mamm. Genome 32, 263–281. doi: 10.1007/s00335-021-09884-2
Callahan, B. J., McMurdie, P. J., Rosen, M. J., Han, A. W., Johnson, A. J. A., and Holmes, S. P. (2016). DADA2: high-resolution sample inference from Illumina amplicon data. Nat. Methods 13, 581–583. doi: 10.1038/nmeth.3869
Camacho, C., Boratyn, G. M., Joukov, V., Vera Alvarez, R., and Madden, T. L. (2023). ElasticBLAST: accelerating sequence search via cloud computing. BMC Bioinformatics 24:117. doi: 10.1186/s12859-023-05245-9
Cani, P. D., Depommier, C., Derrien, M., Everard, A., and de Vos, W. M. (2022). Akkermansia muciniphila: paradigm for next-generation beneficial microorganisms. Nat. Rev. Gastroenterol. Hepatol. 19, 625–637. doi: 10.1038/s41575-022-00631-9
Cheviron, Z. A., Connaty, A. D., McClelland, G. B., and Storz, J. F. (2014). Functional genomics of adaptation to hypoxic cold-stress in high-altitude deer mice: transcriptomic plasticity and thermogenic performance. Evolution 68, 48–62. doi: 10.1111/evo.12257
Chong-Neto, H. J., D’amato, G., and Rosário Filho, N. A. (2022). Impact of the environment on the microbiome. J. Pediatr. 98, S32–S37. doi: 10.1016/j.jped.2021.10.001
Conlon, M. A., and Bird, A. R. (2014). The impact of diet and lifestyle on gut microbiota and human health. Nutrients 7, 17–44. doi: 10.3390/nu7010017
Dąbrowska, K., and Witkiewicz, W. (2016). Correlations of host genetics and gut microbiome composition. Front. Microbiol. 7:1357. doi: 10.3389/fmicb.2016.01357
Darden, T. R. (1972). Respiratory adaptations of a fossorial mammal, the pocket gopher (Thomomys bottae). J. Comp. Physiol. 78, 121–137. doi: 10.1007/BF00693609
Davenport, E. R., Mizrahi-Man, O., Michelini, K., Barreiro, L. B., Ober, C., and Gilad, Y. (2014). Seasonal variation in human gut microbiome composition. PLoS One 9:e90731. doi: 10.1371/journal.pone.0090731
De Filippo, C., Cavalieri, D., Di Paola, M., Ramazzotti, M., Poullet, J. B., Massart, S., et al. (2010). Impact of diet in shaping gut microbiota revealed by a comparative study in children from Europe and rural Africa. Proc. Natl. Acad. Sci. U. S. A. 107, 14691–14696. doi: 10.1073/pnas.1005963107
Di Musciano, M., Zannini, P., Ferrara, C., Spina, L., Nascimbene, J., Vetaas, O. R., et al. (2021). Investigating elevational gradients of species richness in a Mediterranean plant hotspot using a published flora. Front. Biogeogr. 13:7. doi: 10.21425/F5FBG50007
Durk, R. P., Castillo, E., Márquez-Magaña, L., Grosicki, G. J., Bolter, N. D., Lee, C. M., et al. (2019). Gut microbiota composition is related to cardiorespiratory fitness in healthy young adults. Int. J. Sport Nutr. Exerc. Metab. 29, 249–253. doi: 10.1123/ijsnem.2018-0024
Earl, D. A., and vonHoldt, B. M. (2012). STRUCTURE HARVESTER: a website and program for visualizing STRUCTURE output and implementing the Evanno method. Conserv. Genet. Resour. 4, 359–361. doi: 10.1007/s12686-011-9548-7
Edgar, R. C., Haas, B. J., Clemente, J. C., Quince, C., and Knight, R. (2011). UCHIME improves sensitivity and speed of chimera detection. Bioinformatics 27, 2194–2200. doi: 10.1093/bioinformatics/btr381
Elderman, M., Hugenholtz, F., Belzer, C., Boekschoten, M., van Beek, A., de Haan, B., et al. (2018). Sex and strain dependent differences in mucosal immunology and microbiota composition in mice. Biol. Sex Differ. 9:26. doi: 10.1186/s13293-018-0186-6
Fan, C., Zhang, L., Jia, S., Tang, X., Fu, H., Li, W., et al. (2022). Seasonal variations in the composition and functional profiles of gut microbiota reflect dietary changes in plateau pikas. Integr. Zool. 17, 379–395. doi: 10.1111/1749-4877.12630
Fröhlich, E., and Wahl, R. (2019). Microbiota and thyroid interaction in health and disease. Trends Endocrinol. Metab. 30, 479–490. doi: 10.1016/j.tem.2019.05.008
Gerhardt, P., Begall, S., Frädrich, C., Renko, K., Hildebrandt, T. B., Holtze, S., et al. (2023). Comparative analysis of thyroid hormone systems in rodents with subterranean lifestyle. Sci. Rep. 13:3122. doi: 10.1038/s41598-023-30179-w
Goudet, J. (2005). Hierfstat, a package for R to compute and test hierarchical F-statistics. Mol. Ecol. Notes 5, 184–186. doi: 10.1111/j.1471-8286.2004.00828.x
Grenda, A., Iwan, E., Chmielewska, I., Krawczyk, P., Giza, A., Bomba, A., et al. (2022). Presence of Akkermansiaceae in gut microbiome and immunotherapy effectiveness in patients with advanced non-small cell lung cancer. AMB Express 12:86. doi: 10.1186/s13568-022-01428-4
Greyson-Gaito, C. J., Bartley, T. J., Cottenie, K., Jarvis, W. M. C., Newman, A. E. M., and Stothart, M. R. (2020). Into the wild: microbiome transplant studies need broader ecological reality. Proc. Biol. Sci. 287:20192834. doi: 10.1098/rspb.2019.2834
Grieneisen, L. E., Livermore, J., Alberts, S., Tung, J., and Archie, E. A. (2017). Group living and male dispersal predict the core gut microbiome in wild baboons. Integr. Comp. Biol. 57, 770–785. doi: 10.1093/icb/icx046
Guardiola, M., Uriz, M. J., Taberlet, P., Coissac, E., Wangensteen, O. S., and Turon, X. (2015). Deep-Sea, deep-sequencing: Metabarcoding extracellular DNA from sediments of marine canyons. PLoS One 10:e0139633. doi: 10.1371/journal.pone.0139633
Guo, N., Wu, Q., Shi, F., Niu, J., Zhang, T., Degen, A. A., et al. (2021). Seasonal dynamics of diet-gut microbiota interaction in adaptation of yaks to life at high altitude. NPJ Biofilms Microbiomes 7:38. doi: 10.1038/s41522-021-00207-6
Henry, L. P., Bruijning, M., Forsberg, S. K. G., and Ayroles, J. F. (2019). Can the microbiome influence host evolutionary trajectories? BioRxiv 2019:237. doi: 10.1101/700237
Henry, L. P., Bruijning, M., Forsberg, S. K. G., and Ayroles, J. F. (2021). The microbiome extends host evolutionary potential. Nat. Commun. 12:5141. doi: 10.1038/s41467-021-25315-x
Heth, G., Golenberg, E. M., and Nevo, E. (1989). Foraging strategy in a subterranean rodent, Spalax ehrenbergi: a test case for optimal foraging theory. Oecologia 79, 496–505. doi: 10.1007/BF00378667
Hird, S. M. (2017). Evolutionary biology needs wild microbiomes. Front. Microbiol. 8:725. doi: 10.3389/fmicb.2017.00725
Hollander, M., Wolfe, D. A., and Chicken, E. (2013). Nonparametric statistical methods. 3rd Edn. Hoboken, NJ: Wiley.
Hu, X., Liu, G., Li, Y., Wei, Y., Lin, S., Liu, S., et al. (2018). High-throughput analysis reveals seasonal variation of the gut microbiota composition within Forest musk deer (Moschus berezovskii). Front. Microbiol. 9:1674. doi: 10.3389/fmicb.2018.01674
Hu, B., Wang, J., Li, Y., Ge, J., Pan, J., Li, G., et al. (2023). Gut microbiota facilitates adaptation of the plateau zokor (Myospalax baileyi) to the plateau living environment. Front. Microbiol. 14:1136845. doi: 10.3389/fmicb.2023.1136845
Huang, J., Liu, D., Wang, Y., Liu, L., Li, J., Yuan, J., et al. (2022). Ginseng polysaccharides alter the gut microbiota and kynurenine/tryptophan ratio, potentiating the antitumour effect of antiprogrammed cell death 1/programmed cell death ligand 1 (anti-PD-1/PD-L1) immunotherapy. Gut 71, 734–745. doi: 10.1136/gutjnl-2020-321031
Indiani, C. M. D. S. P., Rizzardi, K. F., Castelo, P. M., Ferraz, L. F. C., Darrieux, M., and Parisotto, T. M. (2018). Childhood obesity and firmicutes/bacteroidetes ratio in the gut microbiota: a systematic review. Child. Obes. 14, 501–509. doi: 10.1089/chi.2018.0040
Jasirwan, C. O. M., Muradi, A., Hasan, I., Simadibrata, M., and Rinaldi, I. (2021). Correlation of gut Firmicutes/Bacteroidetes ratio with fibrosis and steatosis stratified by body mass index in patients with non-alcoholic fatty liver disease. Biosci. Microbiota Food Health 40, 50–58. doi: 10.12938/bmfh.2020-046
Jiang, F., Gao, H., Qin, W., Song, P., Wang, H., Zhang, J., et al. (2021). Marked seasonal variation in structure and function of gut microbiota in forest and alpine musk deer. Front. Microbiol. 12:699797. doi: 10.3389/fmicb.2021.699797
Jiang, H., Lei, R., Ding, S.-W., and Zhu, S. (2014). Skewer: a fast and accurate adapter trimmer for next-generation sequencing paired-end reads. BMC Bioinformatics 15:182. doi: 10.1186/1471-2105-15-182
Jombart, T. (2008). Adegenet: a R package for the multivariate analysis of genetic markers. Bioinformatics 24, 1403–1405. doi: 10.1093/bioinformatics/btn129
Kahya, A. (2018). Geology and geochemistry of Madenköy (Ulukışla/Niğde) area carbonate-hosted au-ag-Zn±Pb deposits. J. Sci. Eng. 18, 648–663. doi: 10.5578/fmbd.67032
Kamvar, Z. N., Tabima, J. F., and Grünwald, N. J. (2014). Poppr: an R package for genetic analysis of populations with clonal, partially clonal, and/or sexual reproduction. PeerJ 2:e281. doi: 10.7717/peerj.281
Kaneko, C., Shinohara, A., Kikuchi, T., Tokuda, A., Irie, T., Yamada, K., et al. (2023). Distinctly different gut microbiota in Japanese badgers and Japanese raccoon dogs despite sharing similar food habits and environments. Mamm. Biol. 103, 363–373. doi: 10.1007/s42991-023-00362-7
Karanth, K. P., Avivi, A., Beharav, A., and Nevo, E. (2004). Microsatellite diversity in populations of blind subterranean mole rats (Spalax ehrenbergi superspecies) in Israel: speciation and adaptation. Biol. J. Linn. Soc. 83, 229–241. doi: 10.1111/j.1095-8312.2004.00384.x
Kim, Y. S., Unno, T., Kim, B. Y., and Park, M. S. (2020). Sex differences in gut microbiota. World J Mens Health 38, 48–60. doi: 10.5534/wjmh.190009
Klindworth, A., Pruesse, E., Schweer, T., Peplies, J., Quast, C., Horn, M., et al. (2013). Evaluation of general 16S ribosomal RNA gene PCR primers for classical and next-generation sequencing-based diversity studies. Nucleic Acids Res. 41:e1. doi: 10.1093/nar/gks808
Knezevic, J., Starchl, C., Tmava Berisha, A., and Amrein, K. (2020). Thyroid-gut-Axis: how does the microbiota influence thyroid function? Nutrients 12:1769. doi: 10.3390/nu12061769
Kohl, K. D., Skopec, M. M., and Dearing, M. D. (2014). Captivity results in disparate loss of gut microbial diversity in closely related hosts. Conserv. Physiol. 2:cou009. doi: 10.1093/conphys/cou009
Koliada, A., Syzenko, G., Moseiko, V., Budovska, L., Puchkov, K., Perederiy, V., et al. (2017). Association between body mass index and Firmicutes/Bacteroidetes ratio in an adult Ukrainian population. BMC Microbiol. 17:120. doi: 10.1186/s12866-017-1027-1
Kovacs, A., Ben-Jacob, N., Tayem, H., Halperin, E., Iraqi, F. A., and Gophna, U. (2011). Genotype is a stronger determinant than sex of the mouse gut microbiota. Microb. Ecol. 61, 423–428. doi: 10.1007/s00248-010-9787-2
Kuang, Z., Li, F., Duan, Q., Tian, C., Nevo, E., and Li, K. (2022). Host diet shapes functionally differentiated gut microbiomes in sympatric speciation of blind mole rats in upper Galilee, Israel. Front. Microbiol. 13:1062763. doi: 10.3389/fmicb.2022.1062763
Lan, D., Ji, W., Lin, B., Chen, Y., Huang, C., Xiong, X., et al. (2017). Correlations between gut microbiota community structures of Tibetans and geography. Sci. Rep. 7:16982. doi: 10.1038/s41598-017-17194-4
Lee, M. A., Burger, G., Green, E. R., and Kooij, P. W. (2021). Relationships between resource availability and elevation vary between metrics creating gradients of nutritional complexity. Oecologia 195, 213–223. doi: 10.1007/s00442-020-04824-4
Lee, S., Sung, J., Lee, J., and Ko, G. (2011). Comparison of the gut microbiotas of healthy adult twins living in South Korea and the United States. Appl. Environ. Microbiol. 77, 7433–7437. doi: 10.1128/AEM.05490-11
Leeming, E. R., Louca, P., Gibson, R., Menni, C., Spector, T. D., and Le Roy, C. I. (2021). The complexities of the diet-microbiome relationship: advances and perspectives. Genome Med. 13:10. doi: 10.1186/s13073-020-00813-7
Ley, R. E., Hamady, M., Lozupone, C., Turnbaugh, P. J., Ramey, R. R., Bircher, J. S., et al. (2008). Evolution of mammals and their gut microbes. Science 320, 1647–1651. doi: 10.1126/science.1155725
Li, K., Dan, Z., Gesang, L., Wang, H., Zhou, Y., Du, Y., et al. (2016d). Comparative analysis of gut microbiota of native tibetan and han populations living at different altitudes. PLoS One 11:e0155863. doi: 10.1371/journal.pone.0155863
Li, H., Li, T., Beasley, D. E., Heděnec, P., Xiao, Z., Zhang, S., et al. (2016a). Diet diversity is associated with Beta but not alpha diversity of Pika gut microbiota. Front. Microbiol. 7:1169. doi: 10.3389/fmicb.2016.01169
Li, H., Li, T., Yao, M., Li, J., Zhang, S., Wirth, S., et al. (2016b). Pika gut may select for rare but diverse environmental bacteria. Front. Microbiol. 7:1269. doi: 10.3389/fmicb.2016.01269
Li, H., Qu, J., Li, T., Li, J., Lin, Q., and Li, X. (2016c). Pika population density is associated with the composition and diversity of gut microbiota. Front. Microbiol. 7:758. doi: 10.3389/fmicb.2016.00758
Li, C. C., Weeks, D. E., and Chakravarti, A. (1993). Similarity of DNA fingerprints due to chance and relatedness. Hum. Hered. 43, 45–52. doi: 10.1159/000154113
Li, L., and Zhao, X. (2015). Comparative analyses of fecal microbiota in Tibetan and Chinese Han living at low or high altitude by barcoded 454 pyrosequencing. Sci. Rep. 5:14682. doi: 10.1038/srep14682
Li, H., Zhou, R., Zhu, J., Huang, X., and Qu, J. (2019). Environmental filtering increases with elevation for the assembly of gut microbiota in wild pikas. Microb. Biotechnol. 12, 976–992. doi: 10.1111/1751-7915.13450
Lin, L., and Zhang, J. (2017). Role of intestinal microbiota and metabolites on gut homeostasis and human diseases. BMC Immunol. 18:2. doi: 10.1186/s12865-016-0187-3
Lindsay, E. C., Metcalfe, N. B., and Llewellyn, M. S. (2020). The potential role of the gut microbiota in shaping host energetics and metabolic rate. J. Anim. Ecol. 89, 2415–2426. doi: 10.1111/1365-2656.13327
Liu, F., Liang, T., Zhang, Z., Liu, L., Li, J., Dong, W., et al. (2021). Effects of altitude on human oral microbes. AMB Express 11:41. doi: 10.1186/s13568-021-01200-0
Lynch, S. V., and Pedersen, O. (2016). The human intestinal microbiome in health and disease. N. Engl. J. Med. 375, 2369–2379. doi: 10.1056/NEJMra1600266
Ma, Y., Ma, S., Chang, L., Wang, H., Ga, Q., Ma, L., et al. (2019). Gut microbiota adaptation to high altitude in indigenous animals. Biochem. Biophys. Res. Commun. 516, 120–126. doi: 10.1016/j.bbrc.2019.05.085
Ma, Q., Ma, J., Cui, J., Zhang, C., Li, Y., Liu, J., et al. (2023). Oxygen enrichment protects against intestinal damage and gut microbiota disturbance in rats exposed to acute high-altitude hypoxia. Front. Microbiol. 14:1268701. doi: 10.3389/fmicb.2023.1268701
Mantel, N. (1967). The detection of disease clustering and a generalized regression approach. Cancer Res. 27, 209–220
Mariat, D., Firmesse, O., Levenez, F., Guimarăes, V., Sokol, H., Doré, J., et al. (2009). The Firmicutes/Bacteroidetes ratio of the human microbiota changes with age. BMC Microbiol. 9:123. doi: 10.1186/1471-2180-9-123
Martin, A. M., Sun, E. W., Rogers, G. B., and Keating, D. J. (2019). The influence of the gut microbiome on host metabolism through the regulation of gut hormone release. Front. Physiol. 10:428. doi: 10.3389/fphys.2019.00428
McArtor, D. B., Lubke, G. H., and Bergeman, C. S. (2017). Extending multivariate distance matrix regression with an effect size measure and the asymptotic null distribution of the test statistic. Psychometrika 82, 1052–1077. doi: 10.1007/s11336-016-9527-8
McFall-Ngai, M., Hadfield, M. G., Bosch, T. C. G., Carey, H. V., Domazet-Lošo, T., Douglas, A. E., et al. (2013). Animals in a bacterial world, a new imperative for the life sciences. Proc. Natl. Acad. Sci. USA 110, 3229–3236. doi: 10.1073/pnas.1218525110
McKenzie, V. J., Song, S. J., Delsuc, F., Prest, T. L., Oliverio, A. M., Korpita, T. M., et al. (2017). The effects of captivity on the mammalian gut microbiome. Integr. Comp. Biol. 57, 690–704. doi: 10.1093/icb/icx090
McMurdie, P. J., and Holmes, S. (2013). Phyloseq: an R package for reproducible interactive analysis and graphics of microbiome census data. PLoS One 8:e61217. doi: 10.1371/journal.pone.0061217
Midha, A. D., Zhou, Y., Queliconi, B. B., Barrios, A. M., Haribowo, A. G., Chew, B. T. L., et al. (2023). Organ-specific fuel rewiring in acute and chronic hypoxia redistributes glucose and fatty acid metabolism. Cell Metab. 35, 504–516.e5. doi: 10.1016/j.cmet.2023.02.007
Mullur, R., Liu, Y.-Y., and Brent, G. A. (2014). Thyroid hormone regulation of metabolism. Physiol. Rev. 94, 355–382. doi: 10.1152/physrev.00030.2013
Münger, E., Montiel-Castro, A. J., Langhans, W., and Pacheco-López, G. (2018). Reciprocal interactions between gut microbiota and host social behavior. Front. Integr. Neurosci. 12:21. doi: 10.3389/fnint.2018.00021
Nei, M. (1973). Analysis of gene diversity in subdivided populations. Proc. Natl. Acad. Sci. USA 70, 3321–3323. doi: 10.1073/pnas.70.12.3321
Nevo, E. (1979). Adaptive convergence and divergence of subterranean mammals. Annu. Rev. Ecol. Syst. 10, 269–308. doi: 10.1146/annurev.es.10.110179.001413
Nevo, E. (2007). “Mosaic evolution of subterranean mammals: tinkering, regression, progression, and global convergence” in Subterranean rodents (Berlin, Heidelberg: Springer Berlin Heidelberg), 375–388.
Nollet, L., and Verstraete, W. (1996). Gastro-enteric methane versus sulphate and volatile fatty acid production. Environ. Monit. Assess. 42, 113–131. doi: 10.1007/BF00394045
Nuriel-Ohayon, M., Neuman, H., and Koren, O. (2016). Microbial changes during pregnancy, birth, and infancy. Front. Microbiol. 7:1031. doi: 10.3389/fmicb.2016.01031
Oksanen, J., Simpson, G. L., Blanchet, F. G., Kindt, R., Legendre, P., Minchin, P. R., et al. (2024). Vegan: community ecology package. The R Foundation.
Org, E., Mehrabian, M., Parks, B. W., Shipkova, P., Liu, X., Drake, T. A., et al. (2016). Sex differences and hormonal effects on gut microbiota composition in mice. Gut Microbes 7, 313–322. doi: 10.1080/19490976.2016.1203502
Org, E., Parks, B. W., Joo, J. W. J., Emert, B., Schwartzman, W., Kang, E. Y., et al. (2015). Genetic and environmental control of host-gut microbiota interactions. Genome Res. 25, 1558–1569. doi: 10.1101/gr.194118.115
Ortiga-Carvalho, T. M., Chiamolera, M. I., Pazos-Moura, C. C., and Wondisford, F. E. (2016). Hypothalamus-pituitary-thyroid Axis. Compr. Physiol. 6, 1387–1428. doi: 10.1002/cphy.c150027
Paradis, E., and Schliep, K. (2019). Ape 5.0: an environment for modern phylogenetics and evolutionary analyses in R. Bioinformatics 35, 526–528. doi: 10.1093/bioinformatics/bty633
Pellizzon, M. A., and Ricci, M. R. (2018). Effects of rodent diet choice and fiber type on data interpretation of gut microbiome and metabolic disease research. Curr. Protoc. Toxicol. 77:e55. doi: 10.1002/cptx.55
Petersen, C., Hamerich, I. K., Adair, K. L., Griem-Krey, H., Torres Oliva, M., Hoeppner, M. P., et al. (2023). Host and microbiome jointly contribute to environmental adaptation. ISME J. 17, 1953–1965. doi: 10.1038/s41396-023-01507-9
Pew, J., Muir, P. H., Wang, J., and Frasier, T. R. (2015). Related: an R package for analysing pairwise relatedness from codominant molecular markers. Mol. Ecol. Resour. 15, 557–561. doi: 10.1111/1755-0998.12323
Popa, O. P., Chişamera, G. B., Murariu, D., and Popa, L. O. (2014). Development of nuclear microsatellite markers for the lesser blind mole rat Nannospalax leucodon (Rodentia: Spalacidae). Conserv. Genet. Resour. 6, 787–789. doi: 10.1007/s12686-014-0220-x
Portincasa, P., Khalil, M., Graziani, A., Frühbeck, G., Baffy, G., Garruti, G., et al. (2024). Gut microbes in metabolic disturbances. Promising role for therapeutic manipulations? Eur. J. Intern. Med. 119, 13–30. doi: 10.1016/j.ejim.2023.10.002
Pritchard, J. K., Stephens, M., and Donnelly, P. (2000). Inference of population structure using multilocus genotype data. Genetics 155, 945–959. doi: 10.1093/genetics/155.2.945
Quast, C., Pruesse, E., Yilmaz, P., Gerken, J., Schweer, T., Yarza, P., et al. (2013). The SILVA ribosomal RNA gene database project: improved data processing and web-based tools. Nucleic Acids Res. 41, D590–D596. doi: 10.1093/nar/gks1219
Rado, R., Wollberg, Z., and Terkel, J. (1992). Dispersal of young mole rats (Spalax ehrenbergi) from the Natal burrow. J. Mammal. 73, 885–890. doi: 10.2307/1382211
Raulo, A., Allen, B. E., Troitsky, T., Husby, A., Firth, J. A., Coulson, T., et al. (2021). Social networks strongly predict the gut microbiota of wild mice. ISME J. 15, 2601–2613. doi: 10.1038/s41396-021-00949-3
Roeselers, G., Mittge, E. K., Stephens, W. Z., Parichy, D. M., Cavanaugh, C. M., Guillemin, K., et al. (2011). Evidence for a core gut microbiota in the zebrafish. ISME J. 5, 1595–1608. doi: 10.1038/ismej.2011.38
Rogers, G. B., Kozlowska, J., Keeble, J., Metcalfe, K., Fao, M., Dowd, S. E., et al. (2014). Functional divergence in gastrointestinal microbiota in physically-separated genetically identical mice. Sci. Rep. 4:5437. doi: 10.1038/srep05437
Rogers, R. S., Wang, H., Durham, T. J., Stefely, J. A., Owiti, N. A., Markhard, A. L., et al. (2023). Hypoxia extends lifespan and neurological function in a mouse model of aging. PLoS Biol. 21:e3002117. doi: 10.1371/journal.pbio.3002117
Rothschild, D., Weissbrod, O., Barkan, E., Kurilshikov, A., Korem, T., Zeevi, D., et al. (2018). Environment dominates over host genetics in shaping human gut microbiota. Nature 555, 210–215. doi: 10.1038/nature25973
Russel, J., Thorsen, J., Brejnrod, A. D., Bisgaard, H., Sorensen, S. J., and Burmolle, M. (2018). DAtest: a framework for choosing differential abundance or expression method. BioRxiv 2018:802. doi: 10.1101/241802
Savary, P., Foltête, J., Moal, H., Vuidel, G., and Garnier, S. (2021). graph4lg: a package for constructing and analysing graphs for landscape genetics in R. Methods Ecol. Evol. 12, 539–547. doi: 10.1111/2041-210X.13530
Schippers, M.-P., Ramirez, O., Arana, M., Pinedo-Bernal, P., and McClelland, G. B. (2012). Increase in carbohydrate utilization in high-altitude Andean mice. Curr. Biol. 22, 2350–2354. doi: 10.1016/j.cub.2012.10.043
Scott, K. P., Gratz, S. W., Sheridan, P. O., Flint, H. J., and Duncan, S. H. (2013). The influence of diet on the gut microbiota. Pharmacol. Res. 69, 52–60. doi: 10.1016/j.phrs.2012.10.020
Sha, Y., Ren, Y., Zhao, S., He, Y., Guo, X., Pu, X., et al. (2022). Response of ruminal microbiota-host gene interaction to high-altitude environments in Tibetan sheep. Int. J. Mol. Sci. 23:2430. doi: 10.3390/ijms232012430
Sibai, M., Altuntaş, E., Yıldırım, B., Öztürk, G., Yıldırım, S., and Demircan, T. (2020). Microbiome and longevity: high abundance of longevity-linked Muribaculaceae in the gut of the Long-living rodent Spalax leucodon. OMICS 24, 592–601. doi: 10.1089/omi.2020.0116
Solak, H. M., Sezgin, E., Cizkova, D., Kreisinger, J., Çolak, F., Çetintaş, O., et al. (2023). The microbiota of long-living and cancer-free blind mole rat (Nannospalax xanthodon) from the edge of its distribution in northern Anatolia. Commun. Fac. Sci. Univ. Ank. Ser. C 32, 105–118. doi: 10.53447/communc.1281221
Solak, H. M., Yanchukov, A., Çolak, F., Matur, F., Sözen, M., Ayanoğlu, İ. C., et al. (2020). Altitudinal effects on innate immune response of a subterranean rodent. Zool. Sci. 37, 31–41. doi: 10.2108/zs190067
Sommer, F., Ståhlman, M., Ilkayeva, O., Arnemo, J. M., Kindberg, J., Josefsson, J., et al. (2016). The gut microbiota modulates energy metabolism in the hibernating Brown bear Ursus arctos. Cell Rep. 14, 1655–1661. doi: 10.1016/j.celrep.2016.01.026
Sözen, M. (2005). A biological investigation on Turkish Spalax Guldentaedt, 1770 (Mammalia: Rodentia). G. U. J. Sci. 18, 167–181.
Sözen, M., Matur, F., Çolak, E., Özkurt, Ş., and Karataş, A. (2006). Some karyological records and a new chromosomal form for Spalax (Mammalia: Rodentia) in Turkey. Folia Zool. 55, 247–256.
Squires, R. W., and Buskirk, E. R. (1982). Aerobic capacity during acute exposure to simulated altitude, 914 to 2286 meters. Med. Sci. Sports Exerc. 14, 36–40. doi: 10.1249/00005768-198201000-00007
Storz, J. F., Cheviron, Z. A., McClelland, G. B., and Scott, G. R. (2019). Evolution of physiological performance capacities and environmental adaptation: insights from high-elevation deer mice (Peromyscus maniculatus). J. Mammal. 100, 910–922. doi: 10.1093/jmammal/gyy173
Storz, J. F., and Moriyama, H. (2008). Mechanisms of hemoglobin adaptation to high altitude hypoxia. High Alt. Med. Biol. 9, 148–157. doi: 10.1089/ham.2007.1079
Šumbera, R., Lövy, M., Nevo, E., and Okrouhlík, J. (2023). Thermal biology in the upper Galili Mountain blind mole rat (Nannospalax galili) and an overview of spalacine energetics. J. Therm. Biol. 115:103618. doi: 10.1016/j.jtherbio.2023.103618
Suzuki, T. A. (2017). Links between natural variation in the microbiome and host fitness in wild mammals. Integr. Comp. Biol. 57, 756–769. doi: 10.1093/icb/icx104
Suzuki, T. A., Martins, F. M., and Nachman, M. W. (2018). Altitudinal variation of the gut microbiota in wild house mice. Mol. Ecol. 28, 2378–2390. doi: 10.1111/mec.14905
Suzuki, T. A., Phifer-Rixey, M., Mack, K. L., Sheehan, M. J., Lin, D., Bi, K., et al. (2019). Host genetic determinants of the gut microbiota of wild mice. Mol. Ecol. 28, 3197–3207. doi: 10.1111/mec.15139
Tabrett, A., and Horton, M. W. (2020). The influence of host genetics on the microbiome. [version 1; peer review: 2 approved]. F1000Res 9:835. doi: 10.12688/f1000research.20835.1
Tang, X., Zhang, L., Ren, S., Zhao, Y., and Zhang, Y. (2023). Temporal and geographic distribution of gut microbial enterotypes associated with host thermogenesis characteristics in plateau pikas. Microbiol. Spectr. 11:e0002023. doi: 10.1128/spectrum.00020-23
Thursby, E., and Juge, N. (2017). Introduction to the human gut microbiota. Biochem. J. 474, 1823–1836. doi: 10.1042/BCJ20160510
Tian, Y. M., Guan, Y., Tian, S. Y., Yuan, F., Zhang, L., and Zhang, Y. (2018). Short-term chronic intermittent hypobaric hypoxia alters gut microbiota composition in rats. Biomed. Environ. Sci. 31, 898–901. doi: 10.3967/bes2018.122
Tims, S., Derom, C., Jonkers, D. M., Vlietinck, R., Saris, W. H., Kleerebezem, M., et al. (2013). Microbiota conservation and BMI signatures in adult monozygotic twins. ISME J. 7, 707–717. doi: 10.1038/ismej.2012.146
Ülgen, C., and Tavşanoğlu, Ç. (2024). A taxonomic snapshot of belowground organs in plants of Anatolian steppes. Folia Geobot. 58, 231–243. doi: 10.1007/s12224-024-09442-z
Van Leeuwen, P., Mykytczuk, N., Mastromonaco, G. F., and Schulte-Hostedde, A. I. (2020). Effects of captivity, diet, and relocation on the gut bacterial communities of white-footed mice. Ecol. Evol. 10, 4677–4690. doi: 10.1002/ece3.6221
Virili, C., Stramazzo, I., Bagaglini, M. F., Carretti, A. L., Capriello, S., Romanelli, F., et al. (2024). The relationship between thyroid and human-associated microbiota: A systematic review of reviews. Rev. Endocr. Metab. Disord. 25, 215–237. doi: 10.1007/s11154-023-09839-9
Wang, J., Lang, T., Shen, J., Dai, J., Tian, L., and Wang, X. (2019). Core gut bacteria analysis of healthy mice. Front. Microbiol. 10:887. doi: 10.3389/fmicb.2019.00887
Wang, J., Linnenbrink, M., Künzel, S., Fernandes, R., Nadeau, M.-J., Rosenstiel, P., et al. (2014). Dietary history contributes to enterotype-like clustering and functional metagenomic content in the intestinal microbiome of wild mice. Proc. Natl. Acad. Sci. USA 111, E2703–E2710. doi: 10.1073/pnas.1402342111
Ward, T., Larson, J., Meulemans, J., Hillmann, B., Lynch, J., Sidiropoulos, D., et al. (2017). Bugbase predicts organism level microbiome phenotypes. BioRxiv 2017:462. doi: 10.1101/133462
Wertheim, G., and Nevo, E. (1971). Helminths of birds and mammals from Israel: III. Helminths from chromosomal forms of the mole-rat, Spalax ehrenbergi. J. Helminthol. 45, 161–169. doi: 10.1017/S0022149X00007045
Wu, Y., Yao, Y., Dong, M., Xia, T., Li, D., Xie, M., et al. (2020). Characterisation of the gut microbial community of rhesus macaques in high-altitude environments. BMC Microbiol. 20:68. doi: 10.1186/s12866-020-01747-1
Yatsunenko, T., Rey, F. E., Manary, M. J., Trehan, I., Dominguez-Bello, M. G., Contreras, M., et al. (2012). Human gut microbiome viewed across age and geography. Nature 486, 222–227. doi: 10.1038/nature11053
Yau, W. W., and Yen, P. M. (2020). Thermogenesis in adipose tissue activated by thyroid hormone. Int. J. Mol. Sci. 21:3020. doi: 10.3390/ijms21083020
Zeng, B., Zhang, S., Xu, H., Kong, F., Yu, X., Wang, P., et al. (2020). Gut microbiota of Tibetans and Tibetan pigs varies between high and low altitude environments. Microbiol. Res. 235:126447. doi: 10.1016/j.micres.2020.126447
Zeng, B., Zhao, J., Guo, W., Zhang, S., Hua, Y., Tang, J., et al. (2017). High-altitude living shapes the skin microbiome in humans and pigs. Front. Microbiol. 8:1929. doi: 10.3389/fmicb.2017.01929
Zhang, W., Jiao, L., Liu, R., Zhang, Y., Ji, Q., Zhang, H., et al. (2018). The effect of exposure to high altitude and low oxygen on intestinal microbial communities in mice. PLoS One 13:e0203701. doi: 10.1371/journal.pone.0203701
Zhang, Z., Xu, D., Wang, L., Hao, J., Wang, J., Zhou, X., et al. (2016). Convergent evolution of rumen microbiomes in high-altitude mammals. Curr. Biol. 26, 1873–1879. doi: 10.1016/j.cub.2016.05.012
Zhao, J., Yao, Y., Dong, M., Xiao, H., Xiong, Y., Yang, S., et al. (2023). Diet and high altitude strongly drive convergent adaptation of gut microbiota in wild macaques, humans, and dogs to high altitude environments. Front. Microbiol. 14:1067240. doi: 10.3389/fmicb.2023.1067240
Zhao, J., Yao, Y., Li, D., Xu, H., Wu, J., Wen, A., et al. (2018). Characterization of the gut microbiota in six geographical populations of Chinese Rhesus macaques (Macaca mulatta), implying an adaptation to high-altitude environment. Microb. Ecol. 76, 565–577. doi: 10.1007/s00248-018-1146-8
Zheng, D., Liwinski, T., and Elinav, E. (2020). Interaction between microbiota and immunity in health and disease. Cell Res. 30, 492–506. doi: 10.1038/s41422-020-0332-7
Zhernakova, D. V., Wang, D., Liu, L., Andreu-Sánchez, S., Zhang, Y., Ruiz-Moreno, A. J., et al. (2024). Host genetic regulation of human gut microbial structural variation. Nature 625, 813–821. doi: 10.1038/s41586-023-06893-w
Keywords: gut microbiome, diet, thyroid, altitude adaptation, high altitude, blind mole rats, 16S, 18S
Citation: Solak HM, Kreisinger J, Čížková D, Sezgin E, Schmiedová L, Murtskhvaladze M, Henning Y, Çolak F, Matur F and Yanchukov A (2024) Altitude shapes gut microbiome composition accounting for diet, thyroid hormone levels, and host genetics in a subterranean blind mole rat. Front. Microbiol. 15:1476845. doi: 10.3389/fmicb.2024.1476845
Edited by:
Lin Zhang, Hubei University of Chinese Medicine, ChinaReviewed by:
Adrien Assié, Baylor College of Medicine, United StatesXin Dai, Yangzhou University, China
Copyright © 2024 Solak, Kreisinger, Čížková, Sezgin, Schmiedová, Murtskhvaladze, Henning, Çolak, Matur and Yanchukov. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.
*Correspondence: Halil Mert Solak, c29sYWtoQG5hdHVyLmN1bmkuY3o=
†ORCID: Halil Mert Solak, orcid.org/0000-0001-6877-9274
Jakub Kreisinger, orcid.org/0000-0001-9375-9814
Dagmar Cizkova, orcid.org/0000-0001-5031-7792
Efe Sezgin, orcid.org/0000-0002-8000-7485
Lucie Schmiedová, orcid.org/0000-0003-2180-5281
Marine Murtskhvaladze, orcid.org/0000-0002-4647-2287
Yoshiyuki Henning, orcid.org/0000-0002-0166-2204
Faruk Çolak, orcid.org/0000-0003-3985-7864
Ferhat Matur, orcid.org/0000-0001-9488-1408
Alexey Yanchukov, orcid.org/0000-0002-9613-8770