Longitudinal multi-omics analysis uncovers the altered landscape of gut microbiota and plasma metabolome in response to high altitude

Abstract

Background

Gut microbiota is significantly influenced by altitude. However, the dynamics of gut microbiota in relation to altitude remains undisclosed.

Methods

In this study, we investigated the microbiome profile of 610 healthy young men from three different places in China, grouped by altitude, duration of residence, and ethnicity. We conducted widely targeted metabolomic profiling and clinical testing to explore metabolic characteristics.

Results

Our findings revealed that as the Han individuals migrated from low altitude to high latitude, the gut microbiota gradually converged towards that of the Tibetan populations but reversed upon returning to lower altitude. Across different cohorts, we identified 51 species specifically enriched during acclimatization and 57 species enriched during deacclimatization to high altitude. Notably, Prevotella copri was found to be the most enriched taxon in both Tibetan and Han populations after ascending to high altitude. Furthermore, significant variations in host plasma metabolome and clinical indices at high altitude could be largely explained by changes in gut microbiota composition. Similar to Tibetans, 41 plasma metabolites, such as lactic acid, sphingosine-1-phosphate, taurine, and inositol, were significantly elevated in Han populations after ascending to high altitude. Germ-free animal experiments demonstrated that certain species, such as Escherichia coli and Klebsiella pneumoniae, which exhibited altitude-dependent variations in human populations, might play crucial roles in host purine metabolism.

Conclusions

This study provides insights into the dynamics of gut microbiota and host plasma metabolome with respect to altitude changes, indicating that their dynamics may have implications for host health at high altitude and contribute to host adaptation.

Video Abstract

Introduction

Living at high altitude presents various challenges for humans, but millions of lowlanders seek out the experience each year as hikers, skiers, or mountaineers. The Qinghai-Tibet Plateau (Tibetan Plateau) in China is a representative high-altitude environment with an average altitude of 4500 meters (m). It is mainly characterized by decreased barometric pressure, hypoxia, low temperature, and intense ultraviolet light, which significantly affect human health. A previous study on the gut microbiota of Chinese Han living in the plain and Tibetan Plateau revealed that high altitude may contribute to shaping human gut microbiota [1]. Furthermore, Tibetans, as the indigenous inhabitants of the Tibetan Plateau, exhibited significant variations in their gut microbial structure with increasing altitude [2], indicating that altitude has an impact on the gut microbiota of both the Han and Tibetan.

Generally, the journey from lowland to high altitude and back to lowland can be divided into three stages: acute response, acclimatization, and deacclimatization. The acute response stage occurs within a week of arriving at high altitude and is often accompanied by acute mountain sickness (AMS), which is the most prevalent ailment experienced by lowlanders at high altitude. Symptoms of AMS typically resolve 3 months after arriving at high altitude during the acclimatization stage. Despite significant advancements in comprehending the molecular mechanisms underlying altitude acclimatization, the determinants of successful acclimatization remain poorly understood [3]. Upon descending to lower altitudes following high-altitude acclimatization, lowlanders gradually lose their hypoxia tolerance over time and encounter alterations in hemoglobin and hormone levels, known as high-altitude deacclimatization [4]. The body undergoes distinct impairments during different stages. However, it remains unclear whether dynamic gut microbial ecosystems are characterized by different stages and participate in altitude sickness.

Beyond the mere identification of gut microbial taxa associated with altitude adaptation and acclimatization, efforts have been made to elucidate the potential functions of signature microbes. Changes in the gut microbiota and metabolites contribute to the development of altitude-induced cardiac hypertrophy in rats during hypobaric hypoxia challenge [5], implying a possible link between altered gut microbes after high altitude exposure and the occurrence of altitude sickness. Therefore, it is imperative to unravel changes in gut microbiota and blood metabolites as well as their associations following ascent to high altitude.

To gain a more profound understanding of the dynamics of gut microbiota after ascending to the high altitude and address the aforementioned knowledge gaps, we conducted a comprehensive longitudinal study integrating multidimensional datasets on the gut microbiome, plasma metabolome, and clinical indices of 138 Tibetans and 472 Hans from three places in China: Kashgar Prefecture (1700 m), Hotan Prefecture (1300 m), and Ngari Prefecture (Tibetan Plateau, 4300 m). All participants were divided into six groups according to altitude, duration of residence, place, and ethnicity. This study aims to provide a comprehensive analysis of gut microbiota, plasma metabolome, and clinical indices in the context of high altitude.

Methods

Stool and plasma samples collection and clinical testing

The fresh stool samples were immediately frozen in liquid nitrogen and stored at − 80 °C. DNA extraction was performed using a magnetic stool DNA kit (TIANGEN, China) according to the manufacturer’s instructions, with the addition of special grinding beads to effectively lyse complex components in fecal samples. Fasting venous blood collected with EDTA-K₂ was centrifuged at 4000 rpm for 15 min to separate plasma, which was then stored at − 80 °C. Plasma clinical indices were analyzed using a hematology analyzer (Cobas 6000; Roche, USA) at Chinese PLA General Hospital. All participants voluntarily participated in this study and provided signed informed consent forms before sample collection. This study was approved by the Ethics Committee of the Chinese PLA General Hospital and conducted in accordance with national and institutional ethical guidelines (S2020-517–01).

Metagenomic sequencing

Quality control

Shotgun paired-end sequencing of all DNA samples was conducted on a HiSeq 2500 instrument (Illumina, USA), generating reads with a fragment length of 350 bp and a read length of 150 bp at Novogene Company (China). High-quality reads were obtained by filtering out reads containing low-quality bases (quality value ≤ 38), ambiguous “N” bases, adapter contamination, and human genome sequences from raw data using Readfq (v8.0) and Bowtie2 (v2.2.4) [6]. The human genome (GRCh37) was used as reference for decontamination. Clean reads were assembled into scaftigs (continuous sequences within scaffolds) using SOAPdenovo software (v2.04) [7] with default parameters. Unused reads were selected through mapping against scaftigs for subsequent mixed assembly using Bowtie2 under identical parameters.

Taxonomy and function annotation and abundance profiling

The assembled scaftigs were utilized for gene prediction using MetaGeneMark (v2.1) [8] with default parameters. CD-HIT (v4.5.8) [9] was employed to obtain a unique initial gene catalogue, which was then mapped to the clean data in order to acquire the gene catalogue (Unigenes). The unigenes were aligned to the NCBI NR (v2018.01) database using DIAMOND (v0.9.9) [10]. The alignment results were used for conducting the lowest common ancestor algorithm, enabling us to obtain taxonomic hierarchy where the abundance of a taxonomic level in one sample equaled the sum of annotated genes assigned to it. Furthermore, we aligned the gene catalogue with Kyoto Encyclopedia of Genes and Genomes (KEGG) database (v2018.01) using DIAMOND in order to acquire functional hierarchy and their relative abundance.

Microbiota diversity comparison

Microbiota diversity analysis was conducted using R package vegan (v2.5–6). The Shannon diversity index was used to assess α-diversity, while principal coordinate analysis (PCoA) based on Bray–Curtis distance and analysis of similarity (Anosim) were conducted to evaluate differences in microbial structure and composition across groups. Hierarchical clustering analysis was performed based on the mean abundance of species in each group.

Genome assembly

We performed genome assembly for Alistipes sp.CAG:435 and Alistipes sp.CAG:514 by mapping reads from each sample obtained from Tibetan4k_NP onto a reference genome. The mapped reads were assembled using SOAP denovo (v2.04), SPAdes (v3.10.0), and ABySS (v1.3.7) respectively. Subsequently, the contigs generated from these assemblies were integrated utilizing Contig Integrator for Sequence Assembly (CISA). The preliminary assembly results underwent optimization through GapCloser (v1.12). Gene prediction was carried out utilizing GeneMarkS (v4.17). The gene function analysis methods employed were consistent with those used in the functional annotation of the sequencing dataset. Based on the assembly outcomes, the results with the best quality were selected for functional enrichment analysis. Sequence alignment and gene enrichment analyses were conducted through online tool BlastKOALA.

Widely targeted mass spectrometry

Plasma samples preparation for ultra-performance liquid chromatography

Plasma samples (50 μL) were transferred into Eppendorf tubes, followed by the addition of pre-chilled ice-cold methanol containing 1 µg/mL of 2-chlorophenylalanine as an internal standard. The mixture was stirred for 3 min and centrifuged at 12,000 r/min for 10 min at a temperature of 4 °C. The resulting supernatant was further centrifuged at the same speed and temperature for another 5 min, and the new supernatant was retained for ultra-performance liquid chromatography-tandem mass spectrometry (UPLC-MS/MS) analysis.

UPLC-MS/MS conditions

Metabolites were separated using an ACQUITY UPLC HSS T3 C18 column (1.8 µm, 2.1 × 100 mm, Waters, MA, USA). The mobile phase consisted of ultrapure water containing 0.04% acetic acid (solvent A), and acetonitrile containing 0.04% acetic acid (solvent B). Metabolites were eluted using the following gradient: at 0 min, A:B ratio was 95:5, (v/v); at 11 min, it changed to a ratio of 5:95; at 12 min, it returned to a ratio of 95:5; and finally at 14 min, it remained constant at a ratio of 95:5. The flow rate of the mobile phase was set to be 0.4 mL/min with a column temperature maintained at 40 °C and an injection volume of 2 μL. The mass spectrometry conditions were as follows: electrospray ionization temperature was set to 500 °C with positive voltage applied as 5500 V and negative voltage as − 4500 V; gas I pressure was 55 psi while gas II pressure was 60 psi; curtain gas pressure was 25 psi; high collision-activated dissociation mode was used. In triple quadrupole, each ion pair underwent scanning based on optimized declustering potential and collision energy.

Plasma metabolites data analysis and quality control

Mass spectrometry data were processed using Analyst 1.6.3. Qualitative and quantitative analyses of the metabolites in the samples were conducted based on the local metabolic database. To ensure repeatability, a quality control sample was prepared by mixing sample extracts and injected after every 10 analyzed samples during instrumental analysis.

Statistical analysis and visualization

The normality test was employed to assess the conformity of the data with a normal distribution, subsequently guiding the selection between parametric and non-parametric tests. Differential abundance of phyla, genera, species, and KEGG functions was assessed using a two-tailed Wilcoxon rank sum test between groups, and the resulting p values were adjusted for multiple testing (q value) using the false discovery rate (FDR) by Benjamin-Hochberg method. Species with an abundance greater than 0.01% in at least half of the samples in any group and a q value less than 0.05 were retained. Principal component analysis (PCA) was conducted to evaluate differences in the composition of the plasma metabolome and clinical indices across groups. The differential analysis of plasma metabolites and clinical indices was performed using a two-tailed Wilcoxon rank sum test, with adjusted p values by FDR less than 0.05 considered significant. The fold change of species, metabolites, and clinical indices was calculated as the ratio of means between each group. Statistical analyses were performed in R (v4.1.0). Differential features were visualized using R package ggplot2 (v2.0.0) and Python package GraPhlAn (v0.9.7) [11].

Machine learning

The eXtreme Gradient Boosting (XGBoost) [12] was used to classify the Han and Tibetan populations, as well as the Han population residing at different altitudes based on the abundance of gut microbes. The area under the receiver operating characteristic curve (AUC) and 95% confidence interval were used to evaluate the performance of the models. All species were first used as raw features for XGBoost for feature selection. To find the most appropriate feature number for each iteration, we fed XGBoost with different numbers of features and evaluated its performance using the mean AUC score of 10 times tenfold cross-validation. The selected features are used to construct the final classifier models. For each group of samples, 70% were randomly selected as the training dataset, and the remaining 30% were selected as the testing dataset using a stratified sampling method. Stratified sampling is a sampling method where a population is divided into similar groups called strata, and then samples are taken from each group. This method is often used to address imbalanced datasets, ensuring that the training and test datasets maintain the same distribution of class labels as the original dataset. We used the train_test_split function from the sklearn.model_selection module (v1.3.0) to conduct stratified sampling. All functions were performed using the Python package scikit-learn (v0.21.3). The significance of the AUC was estimated by a permutation test using the R package sigr (v1.1.4).

Pathway analysis

KEGG orthology (KO) genes that were differentially abundant and involved in metabolism were annotated and displayed in the heatmap, with clustering applied to those appearing within the same module. Differentially abundant genes significantly enriched in the same KEGG module were manually constructed by modifying KEGG pathways.

Effect size analyses

We performed an effect size analysis to determine the potential interplay between omics datasets. To assess the proportion of variance of clinical indices or metabolome explained by gut microbiome or clinical indices, firstly, adonis analysis in the R package vegan was used to estimate the “one-to-all” effect size (R²) between each single variable of the secondary omic (gut microbiome) to the whole original omic dataset (plasma metabolome) [13]. Variables with significant effects were retained (p < 0.05; 999 permutations). Then, the Pearson correlation coefficient between variables was calculated to eliminate redundant variables, and variables with coefficients greater than 0.4 were removed.

Correlation analysis

Correlation analysis was performed by calculating Spearman’s rank correlation coefficient (rho) using R package psych (v2.3.3), as well as Pearson correlation coefficient (r) using the stat_cor function in R package ggpubr (v0.6.0). The p values obtained from Spearman’s rank correlation results were adjusted for multiple testing (q value) using FDR via Benjamin-Hochberg method. The results were visualized using R packages ggplot2, circlize [14], and cytoscape (v3.8.2) [15].

Germ-free animal experiments and plasma uric acid test

For the animal experiments, 10-week-old germ-free C57BL/6J male mice (n = 15) obtained from Gempharmatech Co., Ltd., were randomly assigned to three groups: The E. coli group (n = 5) was gavaged with Escherichia coli, the K. pneu group (n = 5) was gavaged with Klebsiella pneumoniae, and the control group (n = 5) was gavaged with saline. Animal research was approved by the Institutional Animal Care and Use Committee (GPTAP20220429-01) of Gempharmatech Co., Ltd., and was certified by the Association for Assessment and Accreditation of Laboratory Animal Care. Animal experiments were performed strictly in accordance with the Guide for the Care and Use of Laboratory Animals published by the Chinese PLA General Hospital. Mice were gavaged with 200 µL of either Escherichia coli or Klebsiella pneumoniae (10⁹ c.f.u. mL⁻¹) three times every 2 days. On the first and seventh days after the gavage experiment, fecal samples were collected to verify the results of colonization by staining microscopy with crystal violet dye. Two weeks after the gavage experiment, plasma was collected and stored at − 80 °C. Mouse plasma UA levels were detected using a mouse UA ELISA kit (BlueGene Biotech) by measuring the intensity of the color spectrophotometrically at 450 nm in a microplate reader, according to the recommended procedure.

Results

Study cohorts

We collected stool samples from 610 healthy young men of Tibetan and Han ethnicities from three places in China (Table 1), and the basic information of each participant is shown in Supplementary Table 1. The Han individuals residing in the plains (altitude < 1700 m) were recruited from Kashgar Prefecture (Han1k_KP) and Hotan Prefecture (Han1k_HP), while all Tibetans living on the Tibetan Plateau (altitude > 4300 m) were from Ngari Prefecture (Tibetan4k_NP). The remaining Han participants were divided into three groups based on duration of residence and altitude, forming a time-altitude cohort. That included those who ascended to the Tibetan Plateau from Kashgar Prefecture for the first time to live for 1 week (Han4k_1w, the acute response stage) and 6 months (Han4k_6m, the acclimatization stage), and those who had resided on the plateau for more than 6 months before descending back to Kashgar Prefecture for 3 months (Han4k_d3m, the deacclimatization stage).

Table 1 Study cohorts and group information

To eliminate the confounding effects of dietary and lifestyle factors, we followed a fixed cohort of Han Chinese population as the time-altitude cohort. The diet primarily consisted of traditional Chinese cuisine, including rice, wheat, vegetables, and meat. Most of the dishes were served daily while others followed a weekly cycle. All participants ate together in the same canteen, as required by their employer, and at the same time for meals, work, and breaks. Tibetan participants were recruited from the plateau and then lived with the Han cohort together, and the same diets and lifestyles were maintained. In addition, in order to eliminate the influence of other confounding factors, participants for whom specimens were collected had to have not taken antibiotics in the 3 months before sampling. Given that alcohol is a risk factor for altitude sickness, all participants had not consumed alcoholic beverages since recruitment. The baseline characteristics of all participants are presented in Table 2. A detailed description of our study workflow can be found in Fig. 1A and B.

Table 2 The baseline characteristics of each group

Fig. 1

Study design and the dynamics of gut microbiota with altitude and residence time. A The collection of samples and the application of analytical methods. B The study cohort design and groups. C and D PCoA using the Bray–Curtis distance calculated based on the relative abundance of species. The 95% confidence ellipses are shown by circles. E The Bray–Curtis distance between the microbiota of the Han populations and Tibetan4k_NP. Wilcoxon rank sum test was used to evaluate the statistical significance of Bray–Curtis distance to Tibetan4k_NP between Han4k_6m and other groups. “ +  +  + ” above Han4k_6m indicates a significantly smaller Bray–Curtis distance for Han4k_6m compared to the other three groups, with all q values being less than 0.001. F The dendrogram for hierarchical clustering is constructed based on the Bray–Curtis distance between groups. G The relative abundance of top 10 genera in each group

Fluctuation of gut microbiota with altitude

Tibetans, the largest indigenous group residing in the Tibetan Plateau, have developed a range of adaptive characteristics to cope with the challenging environment. The Shannon index revealed no significant difference in terms of α-diversity between Han1k_KP and Han1k_HP, but both were lower than that of Tibetan4k_NP (Supplementary Fig. 1A). Furthermore, a principal coordinate analysis (PCoA) plot demonstrated distinct clusters among Han1k_KP, Han1k_HP, and Tibetan4k_NP over the first two dimensions (Anosim test, R = 0.28, p < 0.001, Fig. 1C), while no distinction was observed between Han1k_KP and Han1k_HP (Anosim test, R = 0.065, p = 0.062). These findings indicate a significant separation in microbial composition and structure between the Tibetan and plain-dwelling Han populations.

In the time-altitude cohort, there was no significant difference in α-diversity among Han1k_KP, Han4k_1w, and Han4k_d3m. However, compared to Han4k_1w, Han4K_6m exhibited a significantly decreased α-diversity (Supplementary Fig. 1B). PCoA results indicated that both Han4k_6m and Tibetan4k_NP showed significant dissimilarities from other groups (Anosim test, R = 0.183, p < 0.001, Fig. 1D). To further analyze similarities in microbial composition between the time-altitude cohort group and Tibetan4k_NP, the Bray–Curtis distance from other groups to Tibetan4k_NP was calculated. A clear concave pattern was observed, and Han4k_6m had the smallest community distance to Tibetan4k_NP (p < 0.001, Fig. 1E). Hierarchical cluster analysis also revealed that Han4k_6m was the closest group to Tibetan4k_NP, followed by Han4k_d3m (Fig. 1F). These results suggest that after ascending to plateau, the microbial composition of Hans gradually converged towards that of Tibetans over time, and it retained some plateau-specific characteristics even after descending back to plains for 3 months.

Among the dominant phyla, Bacteroidetes, Chlamydiae, and Euryarchaeota were significantly enriched in Tibetan4k_NP compared to Han1k_KP and Han1k_HP. Conversely, Firmicutes, Proteobacteria, Actinobacteria, and Verrucomicrobia displayed reversed patterns (Supplementary Table 2). A higher abundance of Bacteroidetes and a lower abundance of Firmicutes in the Tibetan have been previously reported [2]. The analysis of the top 10 most abundant genera (Fig. 1G, Supplementary Table 3) revealed that only Prevotella, a diverse genus of gram-negative anaerobic bacteria from Bacteroidetes and known for its proficiency in producing short-chain fatty acid (SCFA) propionate [16], was significantly enriched in Tibetan4k_NP compared to Han1k_KP and Han1k_HP. On the other hand, Escherichia, Streptococcus, and Faecalibacterium showed contrasting trends. In the time-altitude cohort, it was observed that Prevotella and Bacteroides showed a significant increase in abundance in Han4k_6m compared to Han1k_KP. This suggests that the increased Prevotella is the dominant taxon leading the gut microbiota of Han4k_6m to be similar to that of the Tibetan populations.

Specific microbial signatures of acute response, acclimatization, and deacclimatization stages

To further investigate the microbial signatures, differential analysis of species was performed across groups (Fig. 2A, Supplementary Table 4). Only 39 species exhibited differential abundance between Han1k_KP and Han1k_HP, while consistent variations were observed in 186 species between Tibetan4k_NP and both Han1k_KP and Han1k_HP. Among these, 111 species were enriched in Tibetan4k_NP, whereas 75 species showed enrichment in both Han1k_KP and Han1k_HP. The top 20 significantly enriched species in Tibetan4k_NP are presented (Fig. 2B). In addition to Prevotella spp., Alistipes sp. CAG:435, and Alistipes sp. CAG:514 were significantly enriched in Tibetan4k_NP. Strain genome assembly results revealed that both Alistipes sp. CAG:435 and Alistipes sp. CAG:514 exhibited gene functions associated with thermogenesis for environmental adaptation (Supplementary Fig. 1C and D). Moreover, the butyrate-producer Butyrivibrio crossotus, which is inversely associated with obesity [17], displayed a substantially higher abundance (28-fold) in Tibetan4k_NP compared to Han1K_KP group. The top 20 significantly enriched species in both Han1k_KP and Han1k_HP showed that opportunistic pathogens, such as Kluyvera ascorbate, Klebsiella pneumoniae, Salmonella enterica, Escherichia coli, and Shigella, were significantly increased (Fig. 2C).

Fig. 2

Distinct taxonomic signatures associated with the acute response, acclimatization, and deacclimatization. A In total, 296 differentially abundant species are displayed in a phylogenetic tree, primarily clustered within the phyla Firmicutes, Bacteroidetes, and Proteobacteria. The outer circles indicate species that are significantly enriched (orange) or depleted (green) in the group before “vs.” The innermost circle shows the average relative abundance of species across all samples. B The differentially enriched species in Tibetan4k_NP compared to Han1k_KP with the top 20 significance levels were displayed. The statistical significance of Tibetan4k_NP compared to both Han1k_KP and Han1k_HP is indicated above the box plots using regular plus signs. C The differentially enriched species in both Han1k_KP and Han1k_HP compared to Tibetan4k_NP with top 20 significance. The statistical significance of Han1k_KP and Han1k_HP compared to Tibetan4k_NP is indicated above the box plots using regular plus signs. In B and C, the statistical significance of Han1k_KP compared to Han1k_HP is indicated by italic plus signs. D Significantly differential species between Han1k_KP and Han4k_w. E Heat maps of significantly differential species increased in Han4k_6m and Tibetan4k_NP compared to Han1k_KP, while decreased in Han4k_d3m compared to Han4k_6m. F Heat maps of significantly differential species decreased in Han4k_6m and Tibetan4k_NP compared to Han1k_KP, while increased in Han4k_d3m compared to Han4k_6m. Heat maps are scaled by row. The statistical significance of species between two groups is marked by plus signs. The q value < 0.05 denoted as “ + ”; the q value < 0.01 denoted as “ +  + ”; the q value < 0.001 denoted as “ +  +  + ”

To identify specific species associated with the acute response stage, a comparison was conducted between Han4k_1w and Han1k_KP, resulting in the identification of only 14 differential species (Fig. 2D). Notably, Escherichia coli, Klebsiella pneumoniae, Kluyvera ascorbate, and Shigella were significantly decreased in Han4k_1w. Conversely, the abundances of Ruthenibacterium lactatiformans and Clostridium leptum were higher in Han4k_1w. Furthermore, Hungatella hathewayi, which has been positively correlated with circulating taurine concentration in humans [18], showed significant enrichment in Han4k_1w.

To further identify the species associated with acclimatization and deacclimatization, differential analyses were conducted between Han1k_KP and Han4k_6m, as well as between Han4k_6m and Han4k_d3m. A total of 176 differential species were identified during the acclimatization stage, out of which 81 showed significant increase in abundance in Han4k_6m (Fig. 2A “Han4k_6m vs Han1k_KP”). During the deacclimatization stage, a total of 228 differential species were identified, among which 120 species were significantly enriched in Han4k_d3m (Fig. 2A “Han4k_d3m vs Han4k_6m”). Further attention was given to those species that consistently varied between Han1k_KP or Han4k_d3m and Han4k_6m or Tibetan4k_NP.

In total, there were 51 species that showed significant increase in both Han4k_6m and Tibetan4k_NP compared to that of Han1k_KP, while decreased in Han4k_d3m compared to Han4k_6m (Fig. 2E). Prevotella copri, the most prevalent member within Prevotella genus, exhibited a significant increase in both Han4k_6m and Tibetan4k_NP (Supplementary Fig. 1E). Some species known for their ameliorative effects on diseases, such as Bacteroides xylanisolvens alleviating smoking-related non-alcoholic steatohepatitis [19], Bacteroides faecis improving the epithelial barrier integrity [20], and Parabacteroides distasonis alleviating inflammatory arthritis, obesity, and metabolic dysfunctions [21, 22], showed significant increase during acclimatization stage. Furthermore, immune-related species, such as Bacteroides intestinalis enhancing host immune response [23], Bacteroides ovatus and Prevotella buccalis inducing IgA production [24, 25], and Odoribacter splanchnicus promoting intestinal Th17 cell development [26], were significantly increased. Additionally, the enrichment of crucial players known for their abilities to degrade complex molecules like cellulose, hemicellulose, and xylan, such as Bacteroides cellulosilyticus, Prevotella bryantii, and Paraprevotella xylaniphila, were also observed.

In contrast, 57 species showed significant decreases in both Han4k_6m and Tibetan4k_NP compared to Han1k_KP, but increased in Han4k_d3m (Fig. 2F). Notably, during the deacclimatization stage, Escherichia coli, Klebsiella pneumoniae, Kluyvera ascorbate, and Shigella showed significant increases. Faecalibacterium prausnitzii, one of the top 10 most abundant species (Supplementary Fig. 1E) and a crucial butyrate-producing bacterium in human colon, significantly increased in Han4k_d3m. Short-chain fatty acid producers, such as Subdoligranulum variabile, Coprococcus comes, Holdemania filiformis, Clostridium leptum, Anaerotruncus colihominis, Butyricicoccus pullicaecorum, and Marvinbryantia formatexigens, were also significantly elevated. Species belonging to Bifidobacterium like Bifidobacterium adolescentis, Bifidobacterium pseudocatenulatum, Bifidobacterium longum, and lactic acid producers like Lactococcus lactis, Ruthenibacterium lactatiformans, and Enterococcus faecium were all significantly increased. Blautia, a new functional genus, exhibited an increase in species abundance, such as Blautia wexlerae, Blautia obeum, and Blautia producta. Additionally, Dorea longicatena and Tyzzerella nexilis, as potential biomarkers of obesity and type 2 diabetes [27, 28], as well as Clostridioides difficile and Vibrio parahaemolyticus known to cause diarrhea, were all found to be significantly increased.

Identification of Tibetan and Han populations based on gut microbiota

We subsequently developed an eXtreme Gradient Boosting (XGBoost) classifier using a dataset consisting of 138 Tibetan and 472 Han individuals from all cohorts. To determine the optimal feature combination for classifying Tibetan and Han populations across all species, we conducted ten iterations of tenfold cross-validation on the training data. Employing a selection of 24 species, our model achieved an AUC of 0.96 (p < 1e-05) in the test data, effectively discerning between Tibetan and Han individuals (Supplementary Fig. 2A). Furthermore, a new model established to classify 232 Hans residing on the plain (Han1k_KP, Han1k_HP) and 63 Hans residing on the plateau (Han4k_6m). It attained a test data AUC value of 0.96 (p < 1e-05) based on 17 species features (Supplementary Fig. 2B).

Characteristics of plasma metabolome and clinical indices

To further quantify clinical indices and plasma metabolites, we conducted clinical testing and widely targeted metabolome analyses on 418 and 403 plasma samples, respectively. The demographic and clinical characteristics of each group are presented in Table 2. Based on 65 clinical indices, principal component analysis (PCA) revealed significant differences across groups (Fig. 3A, Anosim test, R = 0.14, p < 0.001). During the acute response stage, Han4k_1w showed significantly elevated levels of ProBNP, FPSA, TNT-HS, Insulin, NH3, APOA, TRIGL, CHO, HDL, T4, and A-HCV, whereas CRPHS, CK, ALT, CORT, LIPC, AST, and MYO levels were decreased. During the acclimatization stage, Han4k_6m had significantly higher levels of CA153, LDIP, TESTO, E2, AST, CPEPTID, APOB, FT4, T4, HDL, TP, and ALB compared to Han1k_KP and Han4k_d3m, whereas GLU, PCT, FERR, PRL, PROG, ProBNP, and CRPHS were decreased. Tibetan4k_NP specifically had higher GGT, HE4, ProGRP, TP1NP, and TSH than Han1k_KP, while no significant differences were observed between other groups. Notably, IL6, BILT, AMY-L, E2, TESTO, FT4, and FT3 were significantly higher in Han4k_1w, Han4k_6m, and Tibetan4k_NP than in Han1k_KP, whereas UA, GLU, and PHO were lower.

Fig. 3

Plasma clinical indices and metabolome analysis. A PCA for plasma clinical indices. B PCA for plasma metabolome. Anosim test based on Bray–Curtis distance was conducted to measure the dissimilarity across groups. C Significantly differential plasma metabolites in both time-altitude cohort and Tibetan4k_NP. The plus signs in column Han1k_KP, Han4k_1w, Han4k_6m, Han4k_d3m, and Tibetan4k_NP represent significance of the comparison between Han1k_KP and Han1k_HP, between Han4k_1w and Han1k_KP, between Han4k_6m and Han1k_KP, between Han4k_d3m and Han4k_6m, and between Tibetan4k_NP and Han1k_KP, respectively. The q values < 0.05 denoted as “ + ”; the q values < 0.01 denoted as “ +  + ”; the q values < 0.001 denoted as “ +  +  + .” D The top 10 representative species with the highest explanation of the plasma clinical indices. E The top 10 representative species with the highest explanation of the plasma metabolome. F The top 10 representative plasma clinical indices with the highest explanation of the plasma metabolome. The total variance in clinical indices or metabolome explained by gut microbiome or clinical indices is displayed by the gray column

A total of 629 plasma metabolites were quantified, and PCA showed significant separation of the plasma metabolome among groups (Fig. 3B, Anosim test, R = 0.17, p < 0.001). The 59 metabolites with significant differences in both time-altitude cohort and Tibetan4k_NP are displayed (Fig. 3C). Among these, 41 metabolites were significantly elevated in Han4k_1w, Han4k_6m, and Tibetan4k_NP compared to Han1k_KP but decreased in Han4k_d3m compared to Han4k_6m. In contrast, 8 metabolites exhibited opposite trends. Notably, lactic acid, D-mannose and D-sucrose, and D-trehalose levels increased significantly after ascending to high altitude. Additionally, sphingosine-1-phosphate and sphinganine-1-phosphate levels were upregulated by high altitude exposure. Sphingosine-1-phosphate has been reported to confer adaptive and protective advantages against hypobaric hypoxia during acute or sub-chronic phases [29]. Plasma taurine (2-aminoethanesulfonic acid) levels also showed a significant increase at high altitude, which also has been shown to have protective effects against hypoxia [30, 31]. Inositol elevated in the frontal cortex of hypoxia-treated rats [32] was increased at high altitude. Conversely, metabolites involved in bile acid metabolism including glycochenodeoxycholic acid, glycodeoxycholic acid, and glycoursodeoxycholic acid exhibited significant decreases in Han4k_1w, Han4k_6m, and Tibetan4k_NP.

Altered plasma clinical indices and metabolome are associated with gut microbiota at high altitude

An effect size analysis was conducted to evaluate the proportion of variance in clinical indices or metabolome explained by gut microbiome or clinical indices. The total effect size of gut microbiota on plasma clinical indices accounted for 48.4% of the variance and the top 10 species with the highest explanation were displayed (Fig. 3D). The total effect size of gut microbiota on plasma metabolome accounted for 29.5% of the variance (Fig. 3E). However, compared to that of the gut microbiota, the effect size of plasma clinical indices on metabolome was smaller, and only 13.2% metabolome variance could be explained by clinical indices (Fig. 3F). The top 10 clinical indices with the highest explanation of the metabolome showed that GLU, TP1NP, ALB, UA, and IL6 exerted significant effects on the variance of plasma metabolome. These results suggest that at high altitude, gut microbiota plays an important role in determining both plasma clinical indices and metabolic landscape.

Correlation analysis was performed in both plain populations (Han1k_HP, Han1k_KP) and plateau populations (Han4k_1w, Han4k_6m, Han4k_d3m, and Tibetan4k_NP). In the plain populations, we identified 82 correlations between 16 clinical indices and 48 species (Supplementary Fig. 3A). IL6 exhibited the highest degree of connectivity and was negatively correlated with opportunistic pathogens such as Shigella spp., Salmonella enterica, and Escherichia coli, while positively correlated with Blautia spp. In contrast, in plateau populations, a total of 134 correlations were found between 18 clinical indices and 100 species (Supplementary Fig. 3B). TP1NP and UA were the two nodes with the highest degree of connectivity that correlated with 50 and 25 species respectively. Interestingly, no correlation was observed between TP1NP levels and gut microbiota in plain populations. Notably, at high altitude however, TP1NP that is considered a specific marker for bone formation showed positive correlation with Prevotella sp. CAG:474 and Alistipes sp. CAG:51. Additionally, it was also positively correlated with microbial ascorbate (M00114), phosphatidylethanolamine (M00093), and tetrahydrofolate (M00841) biosynthesis (Supplementary Fig. 4). Furthermore, plasma GLU levels were positively associated with F-type ATPase (M00158), which generates adenosine triphosphate as the universal major energy source, indicating possible regulation by gut microbiota at high altitude.

In plain populations, 538 correlations between gut microbiota and plasma metabolome showed 118 metabolites that were significantly correlated with 184 species (Supplementary Fig. 3C). Methyl 3-hydroxyphenylacetate and acetaminophen glucuronide were the two nodes with the largest degrees. In plateau populations, 2446 correlations showed that 79 metabolites were significantly correlated with 204 species (Supplementary Fig. 3D). Notably, 760 correlations between 137 specific nodes were exclusively identified in plateau populations, and most species belonged to Prevotella. Guanosine 5′-diphosphate (GDP), terephthalic acid, and thymidine were the three largest degree nodes, and their associations with gut microbiota were only identified in plateau populations. Increased correlations suggest that the associations between plasma metabolites and gut microbiota are significantly enhanced at high altitude. Correlations between creatine and gut microbiota were specifically identified in plateau populations and were positively correlated with Blautia hansenii and Clostridioides difficile. Creatine can enhance corticomotor excitability and cognitive performance during oxygen deprivation [33]. We also found that creatine levels were significantly higher in Han4k_1w than in Han1k_KP through plasma metabolome analysis (Supplementary Fig. 6F). Caffeic acid, which is a novel hydroxylase-2 inhibitor that upregulates hypoxia-inducible factors and provides protection against hypoxia [34], was positively correlated with Firmicutes bacterium CAG:110 and Oscillibacter sp. CAG:241_62_21 in both plain and plateau populations. The caffeic acid levels in Tibetan4k_NP were also significantly higher than Han1k_KP (Supplementary Fig. 6G). Hyodeoxycholic acid and glycoursodeoxycholic acid both were positively correlated with Clostridium sp. AT4 and Ruminococcus gnavus, while negatively correlated with Alistipes spp. and Prevotella spp. The plasma levels of sphingosine-1-phosphate were negatively correlated with Streptococcus thermophilus and Dorea longicatena, while displaying the most associations with microbial module pathways (Supplementary Fig. 5).

Altered microbial functions in acclimatization and deacclimatization stages

To further investigate the changes in microbial functions during acclimatization and deacclimatization stages, differential function analysis was conducted between Han1k_KP and Tibetan4k_NP, Han1k_KP and Han4k_6m, and between Han4k_6m and Han4k_d3m. A total of 524 Kyoto Encyclopedia of Genes and Genomes (KEGG) modules were annotated, among which 418 modules were identified as differentially abundant between Tibetan4k_NP and Han1k_KP (Supplementary Table 5). Notably, we observed a significant enrichment of the cAMP signaling module (M00695), which plays an important role in protecting against hypoxia-induced impairment [35], in Tibetan4k_NP. Additionally, the incomplete reductive citrate cycle of methanogenic archaea (M00620) was significantly increased in Tibetan4k_NP, and it was further validated by the enrichment of Methanobacteriaceae (Supplementary Fig. 1F). Previous studies have showed that methanogens can convert pyruvate to the necessary biosynthetic intermediates in response to anaerobic or microaerophilic growth conditions [36]. We speculate that the enrichment of methanogens in Tibetan populations may help their adaptation to high altitude.

Differential analysis identified 353 differential KEGG modules between Han1k_KP and Han4k_6m, as well as 393 differential KEGG modules between Han4k_6m and Han4k_d3m (Supplementary Table 5). To further investigate the metabolism-related microbial functions, we selected and presented 146 significantly differential KEGG orthology genes (KO genes) between Han4k_6m and Han1k_KP, as well as between Han4k_d3m and Han4k_6m (Fig. 4A). Among these genes, a total of 145 KO genes showed reverse changes in abundance in Han4k_6m and Han4k_d3m, indicating their possible regulation by altitude. Microbial pathway modules enriched with differentially abundant KO genes were manually constructed (Fig. 4B). Notably, the KO genes involved in glycolysis (M00001) showed a significant decrease in Han4k_6m. Specifically, the rate-limiting enzyme like pyruvate kinase (PK) was significantly decreased in Han4k_6m but increased in Han4k_d3m. Furthermore, two rate-limiting enzymes involved in the citrate cycle (M00009), namely 2-oxoglutarate dehydrogenase (OGDH) and isocitrate dehydrogenase (IDH3), were decreased in both Han4k_6m and Tibetan4k_NP. Additionally, downregulation of genes involved in ubiquinone biosynthesis, a crucial component of bacterial electron transfer chains, and heme biosynthesis, a precursor to hemoglobin necessary for oxygen transport in the bloodstream, were observed in Han4k_6m and Tibetan4k_NP. These results suggest that, similar to that observed in Tibetans, the gut microbial aerobic respiration was also inhibited during acclimatization to high altitude in Han populations.

Fig. 4

Altitude-associated changes in microbial genes summarized in KO genes and KEGG pathway modules. A Heat maps of 146 significantly differential KO genes between Han4k_6m and Han1k_KP and Han4k_d3m and Han4k_6m. KO genes annotated in the same modules are clustered and marked. Heat maps are scaled by row. B Representative KO genes appearing in A are shown in pathway modules modified from KEGG pathway maps “Glycolysis,” “Lysine biosynthesis,” “Threonine biosynthesis,” “Shikimate pathway,” “Tryptophan biosynthesis,” “Menaquinone biosynthesis,” “Ubiquinone biosynthesis,” Heme biosynthesis,” “GABA biosynthesis,” “Trans-cinnamate degradation,” and “Ascorbate degradation.” Each box in a pathway represents a KO gene and is marked in red for elevation in Han4k_6m or in blue for depletion in Han1k_KP and Han4k_d3m. Black dots in each pathway represent intermediate metabolites. The statistical significance between two groups is marked by plus signs. The q values < 0.05 denoted as “ + ”; the q values < 0.01 denoted as “ +  + ”; the q values < 0.001 denoted as “ +  +  + ”

During the acclimatization stage, Han4k_6m showed significant decreases in KO genes involved in methionine biosynthesis (M00017), methionine salvage pathway (M00034), threonine biosynthesis (M00018), and valine/isoleucine biosynthesis (M00019), indicating a reduction in the gut microbial capacity to synthesize these amino acids. Moreover, microbial genes involved in the shikimate pathway (M00022) for folate and aromatic amino acid synthesis were notably decreased, as were those associated with beta-oxidation (M00087) and triacylglycerol biosynthesis (M00089). Four genes responsible for γ-aminobutyric acid (GABA) production were also downregulated in Han4k_6m and Tibetan4k_NP. Bacteria have been shown to produce and/or consume GABA, which may affect host physiology [37]. Interestingly, plasma metabolome analysis revealed elevated levels of GABA in Han1k_HP and Han1k_KP but reduced levels in Tibetan4k_NP (Supplementary Fig. 6D). Ascorbate degradation was significantly reduced in Han4k_6m, which is vital for regulating transcription factor hypoxia-inducible factor activity [38]. Specifically, tryptophan biosynthesis, the precursor of the neurotransmitter serotonin, was significantly enriched in Han4k_6m, while the levels of tryptamine, another bacterial-derived tryptophan catabolite that stimulates the release of serotonin by enterochromaffin cells [39], showed a significant increase in Han4k_6m (Supplementary Fig. 6E).

Gut microbiota regulates host purine metabolism at high altitude

We observed a significant upregulation of amidophosphoribosyl transferase (purF), the rate-limiting enzyme in inosine monophosphate (IMP) biosynthesis, and IMP cyclohydrolase (purH) in Han4k_6m and Tibetan4k_NP, while observing downregulation in Han1k_KP and Han4k_d3m (Fig. 5A). Additionally, we found a notable decrease in eight genes involved in purine degradation in Han4k_6m and Tibetan4k_NP, while they increased in Han1k_KP and Han4k_d3m. These alterations suggest an elevation of IMP biosynthesis but inhibition of its degradation at high altitude. We speculate that the aberrant microbial purine metabolism at high altitude may be associated with host blood UA levels.

Fig. 5

Gut microbiota fluctuated with altitude is associated with purine metabolism. A Similar to Fig. 4B, significantly differential KO genes involved in purine metabolism are shown in pathway modules modified from KEGG pathway maps “Inosine monophosphate biosynthesis” and “Purine degradation.” Bar plots show relative gene abundances in different groups. B Plasma metabolite levels analyzed by widely targeted mass spectrometry in different groups. C Plasma UA levels detected by clinical test. The statistical significance between two groups is marked by plus signs. The q values < 0.05 denoted as “ + ”; the q values < 0.01 denoted as “ +  + ”; the q values < 0.001 denoted as “ +  +  + .” D Pearson correlation analysis between relative abundance of genes and purine metabolites levels. E Pearson correlation analysis between relative abundance of genes and plasma UA levels

The metabolome analysis revealed a significant increase in the levels of adenosine, inosine, and hypoxanthine in Han4k_6m and Tibetan4k_NP (Fig. 5B). Clinical testing results demonstrated that plasma UA levels were significantly higher in Han1k_KP compared to other groups (Fig. 5C). Conversely, Tibetan4k_NP exhibited significantly lower UA levels than other groups. Furthermore, the UA levels of Han4k_1w, Han4k_6m, and Han4k_d3m were significantly lower than those observed in Han1k_KP. These results indicate a significant reduction in blood UA levels of both Tibetan and Han populations after ascending to high altitude in our study cohorts.

Pearson’s correlation analysis revealed positive associations between the abundance of purF and plasma adenosine, inosine, and hypoxanthine levels (Fig. 5D), indicating the involvement of KO genes in IMP biosynthesis. Furthermore, the KO genes involved in purine degradation, namely ygeT, ygeU, yagR, yagS, yagT, and allB, were found to be positively correlated with plasma UA levels (Fig. 5E). To further investigate the key species associated with plasma UA levels, a correlation analysis was conducted between species and KO genes (Supplementary Fig. 6H, rho > 0.6). The results demonstrated strong positive correlations between seven species including Shigella dysenteriae, Shigella sonnei, Shigella flexneri, Klebsiella pneumoniae, Kluyvera ascorbate, Salmonella enterica, and Escherichia coli with xanthine dehydrogenase (yagR, yagS, ygeU, yagT). Notably, Escherichia coli in the intestine secretes xanthine dehydrogenase which is an important rate-limiting enzyme responsible for oxidative metabolism of purines [40]. These seven species also exhibited positive correlations with allantoinases (allB, hpxB, EC:3.5.2.5, K01466, K16842). Klebsiella pneumoniae was positively correlated with 2-oxo-4-hydroxy-4-carboxy-5-ureidoimidazoline decarboxylase (hpxQ). Positive correlations between Escherichia coli, Klebsiella pneumoniae, and plasma UA levels were also observed (Supplementary Fig. 6AB).

To establish the causal relationship between Escherichia coli or Klebsiella pneumoniae and plasma UA levels, germ-free mice administered via gavage with these two species and plasma UA levels were analyzed using a mouse UA ELISA kit. The experimental design is illustrated in Fig. 6A. In the control group, no significant difference was observed in plasma UA levels before and after administration (Fig. 6B). As anticipated, gavage with Klebsiella pneumoniae and Escherichia coli resulted in a substantial increase in plasma UA levels (Fig. 6CD). Hence, the elevated abundance of Escherichia coli and Klebsiella pneumoniae in the gut contributes to elevated host plasma UA levels.

Fig. 6

Host plasma UA levels regulated by Escherichia coli (E. coli) and Klebsiella pneumoniae (K. pneu). A Animal experimental design for germ-free mice with administration of E. coli or K. pneu. B Paired t test for plasma UA levels of control group between before and after gavage with saline. C Paired t test for plasma UA levels of experiment group between before and after gavage with K. pneu. D Paired t test for plasma UA levels of experiment group between before and after gavage with E. coli. The p values < 0.05 denoted as “ + ”; the p values > 0.05 denoted as “ns”

Discussion

Our study reveals the dynamics of human gut microbiota during acute response, acclimatization, and deacclimatization to high altitude, providing new evidence for the resilience of microbial communities under altitude-induced perturbations. With prolonged residence time on the plateau, the gut microbial structure of Han individuals gradually approached that of Tibetans, suggesting that there is a representative enterotype in high-altitude populations. Given the adaptive enterotype of Tibetans to high altitude, “assimilation” of gut microbiota towards Tibetan may aid Han individuals in acclimatizing to high altitude. Similar to Tibetan populations, the abundance of Prevotella significantly elevated during the acclimatization stage in Han populations, suggesting its potential as a representative genus in high-altitude populations. A hypoxic environment can lead to an increase in anaerobes in the gut [41]. Prevotella spp., as obligate anaerobes, have a competitive advantage over aerobes and thus being more likely to be abundant in human gut at high altitude. Additionally, in rural populations who follow a more pre-industrial and traditional lifestyle, Prevotella tends to dominate the gut microbiota [42]. We speculate that the abundant Prevotella in high-altitude populations is partly associated with low industrialization levels on the Tibetan Plateau. Our study suggests that after acclimatizing to high altitude, the gut microbiota of the Hans shares similar characteristics with those of Tibetans. However, further mechanistic studies are required to determine whether this feature affects host adaptation at high altitude.

As previously reported, host genetic effects on gut microbiota are universal [43]. Given the prolonged high-altitude habitation of Tibetans, natural selection has led to an enrichment of unique genotypes. Notably, EPAS1 exhibits a highly distinctive haplotype structure exclusively found in Tibetans and Denisovans, with a very low frequency among Han Chinese [44]. Based on large-scale Tibetan whole-genome sequencing data, a high-confidence list of 4320 variants and 192 genes that have undergone selection in Tibetans has been characterized [45]. The variation in gut microbiota is influenced by both environmental and genetic factors. However, multiple studies consistently indicate that environmental effects exert a greater impact than additive genetic effects on gut microbiota [43, 46]. In our study, Han and Tibetan populations cohabitated under identical altitudes, diets, and lifestyles conditions. It is plausible that genetic effects potentially contributed to shaping the distinct characteristics of gut microbiota between these two populations despite their prolonged cohabitation. Therefore, we postulated that despite prolonged cohabitation, the gut microbiota of the two populations would exhibit similarities while still retaining discernible distinctions.

Many SCFA producers were elevated in Tibetan and Han populations after acclimatization to high altitude. Among them, Prevotella spp. are proficient in producing SCFA propionate from arabinoxylans and fructo-oligosaccharides [16]. The increased production of SCFA could provide more efficient energy intake for high-altitude populations compared to plain populations. Alistipes sp. CAG:435 and Alistipes sp. CAG:514, two other highly enriched SCFA producers found in Tibetans with genomic functions favoring thermogenesis for environmental adaptation, may be beneficial for adaptation to the low temperature of Tibetan Plateau as well. Cold-exposed microbiota increased concentrations of SCFA and serum ghrelin levels, thus increasing energy intake and thermogenesis in host adaptation to cold [47]. Therefore, the increased SCFA producers resulting from chronic exposure to high altitude may contribute to host acclimatization to such conditions.

Firmicutes increased significantly in the Han population, while Bacteroidetes increased significantly in the Tibetan population. Previous study has shown that Firmicutes is positively associated with obesity while Bacteroidetes is negatively associated with it [48]. In Tibetans, we also found an increase of microbes known to ameliorate type 2 diabetes or obesity, such as Butyrivibrio crossotus, Prevotella copri, and Parabacteroides distasonis. Moreover, a higher ratio of Prevotella to Bacteroides was observed in subjects who showed improved GLU metabolism [49]. Given that Tibetans have lower incidence rates of diabetes and obesity [50, 51], we speculate that this may be related to their unique microbial composition. Moreover, the decreased plasma GLU levels of the Han after ascending to the high altitude might be partly attributed to their gut microbiota approaching that of the Tibetan (Supplementary Fig. 6C). Consequently, future studies focusing on associations between gut microbiota and blood GLU levels among Tibetans may identify more potential probiotics for treating type 2 diabetes.

We identified 51 species that were elevated and 57 species that decreased at high altitude. These species were validated through comparisons across different cohorts and were found to be independent of ethnicity, suggesting their critical responses to altitude. Most of the elevated species belonged to prevotella and Bacteroides, indicating that they are more adaptive to high altitude. Conversely, we found a decrease in some opportunistic pathogens, such as Escherichia coli, Klebsiella pneumoniae, Shigella, Salmonella enterica, and Clostridioides difficile, in high altitude populations. This suggests that populations residing at higher altitudes may face a lower risk of infection from these pathogens compared to those living in lowland areas.

We identified 41 significantly enriched metabolites in Han population after ascending to high altitude, which may contribute to acclimatization. Notably, sphingosine-1-phosphate and taurine, both of which have protective effects against hypobaric hypoxia exposure, were significantly increased in Tibetan and Han populations after ascending to the plateau. Moreover, sphingosine-1-phosphate was positively correlated with 58 microbial module pathways, such as the PTS system, valine/isoleucine biosynthesis, fatty acid biosynthesis, and leucine biosynthesis; thus, the associations between sphingosine-1-phosphate levels and gut microbiota at high altitude warrant further exploration.

The reduced microbial biosynthesis of acetyl-CoA, pyruvate, citrate, and heme in high altitude populations suggests that aerobic respiration of gut microbiota was inhibited. Oxygenation by the bloodstream is closely related to the gut O₂ environment [52], thus may result in a decrease in aerobic bacteria and an increase in anaerobic bacteria [53]. Additionally, we observed a decrease in GABA biosynthesis-related genes in Han4k_6m and Tibetan4k_NP. Downregulated GABA levels in plasma were also observed in the Tibetans. GABA is the primary inhibitory neurotransmitter in the brain and its levels are associated with anxiety and depression [54]. Although the abundance of fecal GABA producers such as Bacteroides was negatively correlated with depression signatures [55], it remains unclear whether gut-derived GABA can cross the blood–brain barrier and improve mood [56]. Exhilaratingly, two GABA-producing Lactobacillus strains have been shown to reduce depression-like behavior in mice [57]. At high altitude, one’s mood tends to change and the exacerbation of depression was concluded [58], and that may be related to a decrease in gut microbes harboring genes responsible for GABA biosynthesis.

Remarkably, IMP biosynthesis in purine metabolism was enhanced at high altitude, while purine degradation was inhibited. Plasma inosine, a purine nucleoside known for its neuroprotective effects against hypoxia [59], increased at high altitude. Hypoxia and reoxygenation have also been reported to elevate inosine levels [60]. Inosine exerts significant cytoprotective action under hypoxic conditions and surprisingly serves as a more efficient energy source than equivalent amounts of GLU [61]. Therefore, the elevated inosine levels observed at high altitude are crucial for protecting against hypoxia and providing energy, making it a potential therapeutic agent for acclimatization to high altitude. IMP is one of the precursors involved in inosine synthesis. Further verification is required to determine whether increased IMP synthesis by gut microbiota promotes plasma inosine production.

Our study identified several species strongly positively associated with microbial genes involved in purine degradation, such as Shigella, Klebsiella pneumoniae, Kluyvera ascorbate, Salmonella enterica, and Escherichia coli, which were decreased at high altitude. The study has provided evidence that Escherichia coli is capable of secreting xanthine dehydrogenase(XOD) [40], a crucial rate-limiting enzyme responsible for the oxidative metabolism of purine. Although these species are opportunistic pathogens, our study indicates that they may play crucial roles in purine degradation at high altitude. The prevalence of hyperuricemia is higher in Tibetan compared to lowlanders [62, 63]. However, its occurrence is influenced by various factors including dietary and cardiovascular risk factors rather than solely being determined by altitude alone [64]. Our research further demonstrates that plasma UA levels are closely associated with alterations in gut microbial composition at high altitude.

In previous study, stool samples were collected at four time points throughout the year, with 3-month intervals between each collection, and found a distinct seasonal cycle in the gut microbiota of the Hadza [65]. Moreover, it was discovered that significant alterations in gut microbiota can occur within a 3-month period following seasonal changes. This study suggests that changes in the external environment are sufficient to alter the gut microbiota within 3 months. Additionally, it is commonly observed that after residing at high altitudes for approximately 3 months, symptoms associated with altitude sickness subside and individuals reach a state where their respiration and pulse rates resemble those experienced at lower altitudes. Consequently, sampling intervals of 1 week, 3 months, and 6 months provide sufficient data for investigating short-term responses to high altitude. However, a comprehensive investigation into the long-term adaptation of gut microbiota to chronic high altitude exposure may necessitate a duration of 1 year or more.

This study has several limitations. Firstly, our study exclusively focused on investigating the impact of high altitude at 4300 m, without considering multiple altitude gradients. Given the variations in oxygen concentration, ultraviolet radiation, and air temperature across different altitudes, it is plausible that the influence on gut microbiota may exhibit varying magnitudes. Consequently, our study does not provide a comprehensive representation of the changes in gut microbiota across all high-altitude gradients. Secondly, our study primarily focused on the short-term response of gut microbiota to high altitude, excluding individuals who had resided at high altitudes for more than a year or had prolonged exposure to such conditions. It is important to note that chronic mountain sickness may develop with long-term exposure to high altitude, thus necessitating further investigation into the alterations in gut microbiota under prolonged exposure, which were beyond the scope of this study. Thirdly, our study revealed a tendency for the gut microbiota of Han to converge towards that of Tibetan after residing at high altitudes for 6 months. However, it is important to note that while Tibetan serves as an exemplary model of successful adaptation to high altitude, the convergence in gut microbiota composition with Tibetan does not necessarily imply facilitation of high-altitude adaptation in Han. Further investigations are warranted to validate these findings. In our future studies, greater emphasis will be placed on elucidating the mechanisms through which gut microbiota contributes to host adaptation to high altitude.

Conclusions

In conclusion, our study unveiled dynamic fluctuations in gut microbiota of the Han as they migrated from low altitude to high altitude and back again. After residing at high altitude for 6 months, the composition of their gut microbiota gradually converged towards that of the Tibetan. Altitude exerted a significant impact on both plasma metabolome and clinical indices, with a substantial portion of their variance being explained by alterations in gut microbiota. Notably, key players involved in purine degradation in gut microbiota, such as Escherichia coli and Klebsiella pneumoniae, played crucial roles in regulating host plasma UA levels at high altitude.