Characterization of gut microbiota and metabolites in renal transplant recipients during COVID-19 and prediction of one-year allograft function

Abstract

Background

The gut-lung-kidney axis is pivotal in immune-related kidney diseases, with gut dysbiosis potentially exacerbating the severity of Coronavirus disease 2019 (COVID-19) in recipients of kidney transplant. This study aimed to characterize the gut microbiome and metabolome in renal transplant recipients with COVID-19 pneumonia over a one-year follow-up period.

Methods

A total of 30 renal transplant recipients were enrolled, comprising 17 with COVID-19 pneumonia, six with mild COVID-19, and seven without COVID-19. Fecal samples were collected at the onset of infection for gut microbiome and metabolome analysis. Generalized Estimating Equations (GEE) model and Latent Class Growth Mixed Model (LCGMM) were employed to dissect the relationships among clinical characteristics, laboratory tests, and gut microbiota and metabolites.

Results

Four microbial phyla (Deferribacteres, TM7, Fusobacteria, and Gemmatimonadetes) and 13 genera were significantly enriched across three recipients groups, correlating with baseline inflammatory response and allograft function. Additionally, 52 differentially expressed metabolites were identified, with seven significantly correlating with eight altered microbiota genera. LCGMM revealed two distinct classes of recipients, with those suffering from COVID-19 pneumonia exhibiting significantly elevated serum creatinine (Scr) trajectories over the one-year period. GEE further identified 12 genera and 181 metabolites closely associated with these trajectories; a multivariable model incorporating gut metabolites of 1-Caffeoylquinic Acid and PMK was found to effectively predict one-year allograft function.

Conclusions

Our study indicates a possible interaction between the composition of the gut microbiota and metabolites community and COVID-19 in renal transplant recipients, particularly in relation to disease severity and the prediction of one-year allograft function.

Introduction

Kidney transplantation is considered the optimal therapy for end-stage renal diseases [1]. Recipients are obligated to a lifelong regimen of immunosuppressive therapy, which not only diminishes the serological response to Severe Acute Respiratory Syndrome Coronavirus 2 (SARS-CoV-2) vaccination but also significantly augments the susceptibility and severity of COVID-19 [2]. Recent studies have underscored the propensity for SARS-CoV-2 infection to precipitate severe pneumonia, acute respiratory distress syndrome, multi-organ failure, hospitalization, and, in extreme cases, mortality among immunocompromised individuals [3]. Moreover, the pathophysiological perturbations induced by COVID-19 are posited to exert enduring effects on survivors [4, 5], with evidence indicating that individuals discharged from hospital care continue to exhibit compromised pulmonary diffusion capacities and aberrant chest imaging findings six months post-infection [6, 7]. Nonetheless, the long-term prognosis of renal allografts in transplant recipients who contract SARS-CoV-2 infection remains an enigma that warrants further elucidation.

The gut microbiota plays a pivotal role in modulating host immune responses and has been posited as a potential modulator of the host immune system [8, 9]. Substantial evidence has demonstrated that complications and adverse events in kidney transplantation are associated with altered gut microbiota profiles, which exhibit a bidirectional relationship with nephropathy [10, 11]. Specifically, gut dysbiosis may trigger immune responses, thereby accelerating the progression of transplant rejection and infections, including but not limited to urinary tract infections and infectious diarrhea [12]. Previous studies have indicated that severe SARS-CoV-2 infection depletes immunomodulatory gut microorganisms, contributing to persistent multisystem inflammation and symptoms post-viral clearance in some patients [13, 14], as well as exhibiting a reduction in fecal microbial diversity and an enrichment of opportunistic pathogens [15]. Subsequently, altered gut microbiota showed impact metabolites, with lipid metabolomics studies in COVID-19 patients revealing significant changes, suggesting a potential interplay between gut micriobiota, metabolism and COVID-19 [16]. The potential relationship between alternation of gut microbiota and recipients with COVID-19 was rarely determined. Recently, one series of study focused on the impact of dietary fibre supplementation on the gut dysbiosis and promote vaccine responsiveness in renal transplant recipients [17, 18]. It is reported the safety and tolerability of dietary inulin, the diversity and differential abundance of gut microbiota, and COVID-19 vaccine-specific immune cell populations and responses. However, the changes in gut microbiota and metabolites in allograft renal transplant recipients in the context of COVID-19 remain to be elucidated.

In this study, we collected fecal specimens from recipients, stratified into cohorts with and without a history of COVID-19 infection. We employed 16 S rRNA sequencing to delineate the taxonomic composition of the gut microbiota and untargeted metabolomics analysis to profile the metabolic constituents. Our objectives were to elucidate the characteristics of both the gut microbiota and metabolites, and to ascertain their correlations with the severity of COVID-19 manifestations and the functional prognosis of renal allografts over a one-year post-transplant period.

Materials and methods

Subjects recruitment and samples collection

The study was conducted following the declaration of Helsinki. Ethical approval for this study was granted by the Ethics Committee of the First Affiliated Hospital of Nanjing Medical University (2023-SRFA-007). Written informed consent from all transplant recipients was obtained, including the consent of clinical demographic and sample collection, gut microbiota sequencing and Metabolomic analysis, outpatient follow-up, as well as potential risks and benefits. The study was strictly limited to the living-related transplantations of kidney donors to lineal or collateral relatives not beyond the third degree of kinship, or transplantations of kidney donors from cadaveric allograft donors after cardiac death.

All recipients should perform the SARS-CoV-2 test by quantitative reverse transcription PCR (qRT-PCR) based on the nasopharyngeal swabs, as well as the blood test of SARS-CoV-2 IgG between December 30th 2022 and February 4th 2023 in our single center (Kidney Transplant Center of the First Affiliated Hospital of Nanjing Medical University). Patients were classified into three groups based on the severity of SARS-CoV-2 infection as reported previously [19, 20] as COVID-19 pneumonia group, COVID-19 group, and Control group. A total of 30 renal transplant recipients were recruited and completed a one-year outpatient follow-up. The inclusion criteria were described as follows: [1] for COVID-19 group and COVID-19 pneumonia group, recipients were diagnosed as COVID-19 when both qRT-PCR and IgG tests of SARS-CoV-2 were positive; while for Control group, recipients with negative tests were considered as no infection; [2] for COVID-19 pneumonia group, patients were enrolled with moderate to severe if pneumonia with fever and respiratory tract symptoms, with or without respiratory rate ≥ 30 breath per minute, oxygen saturation ≤ 93% when breathing ambient air or PaO₂/FiO₂ ≤ 300 mmHg; [3] for COVID-19 group, recipients were enrolled with mild symptoms if there were no radiographic indications of pneumonia; [4] for Control group, patients were enrolled with stable allograft function which was defined as serum creatinine (Scr) level was less than 120 µmol/L for at least three months post-transplant; [5] patients with willingness to participate in the study; [6] patients with age ranging from 18 years old to 60 years old; [7] all patients underwent renal transplant surgery and follow-up in our center with similar eating habits. Exclusion criteria were: [1] patients who did not meet inclusion criteria; [2] participated in other clinical trials; [3] with active viral infections, such as HIV, chronic hepatitis B and C virus, cytomegalovirus, BK polyomavirus; [4] pregnant women; and [5] who failed to perform long-term follow-up in our center.

For COVID-19 pneumonia group, samples were collected on the first day of hospitalization before the administration of anti-infection protocol; while samples were obtained on the first day of arise of fever and respiratory tract symptoms before any medication in COVID-19 group, and on the subsequent first day of the report of negative results in Control group. All fresh samples were collected in collection tubes containing preservative media and stored at -80℃ instantaneously.

Data collection and follow-up

Clinical characteristics of enrolled recipients were obtained from medical records, including age, sex, Body Mass Index (BMI), transplant duration, history of COVID-19 vaccination, immunosuppressive protocol, past medical history, donor of allograft, as well as the symptoms of SARS-CoV-2 infection. A one-year outpatient follow-up was performed for each subject, and laboratory examination (blood routine, urine routine, serum biochemical indexes, serum concentration of immunosuppressive agents, and immunological indexes) and medical records were collected every three months.

DNA extraction and gut microbiota sequencing

We utilized the CTAB protocol to isolate the complete genomic DNA from the samples. The concentration and purity of the DNA were assessed using a 1% agarose gel electrophoresis. Subsequently, the DNA was diluted to a standard concentration of 1 ng/µL with sterile deionized water. Specific regions of the 16 S rRNA gene, namely the V3-V4 hypervariable regions, were targeted for amplification using tailored primers [341 F (5’-CCTAYGGGRBGCASCAG-3’) and 806R (5’-GGACTACNNGGGTATCTAAT-3’)], each equipped with a unique barcode. The PCR amplification was performed in a 15 µL reaction volume containing Phusion^® High-Fidelity PCR Master Mix (New England Biolabs), with 2 µM of each primer, and approximately 10 ng of template DNA. The thermal cycling profile involved an initial denaturation step at 98 °C for 1 min, followed by 30 cycles of denaturation at 98 °C for 10 s, annealing at 50 °C for 30 s, and extension at 72 °C for 30 s, concluding with a final extension at 72 °C for 5 min. The PCR products were then mixed with an equal volume of 1X TAE buffer and subjected to electrophoresis on a 2% agarose gel to verify the amplification. The PCR products were pooled in equal proportions by density. Subsequently, the mixed PCR products were purified using a Universal DNA Purification Kit (TianGen, China). Sequencing libraries were constructed using the NEB Next^® Ultra DNA Library Prep Kit (Illumina, USA), incorporating unique index codes as per the manufacturer’s guidelines. The quality of the sequencing libraries was evaluated using an Agilent 5400 (Agilent Technologies Co Ltd., USA). Finally, the libraries were sequenced on an Illumina sequencing platform, yielding 250 bp paired-end reads.

Extensive analysis of 16 S rRNA data

Our study’s analytical approach was guided by the “Atacama soil microbiome tutorial” from Qiime2docs, with additional customized script adjustments (https://docs.qiime2.org/2019.1/). Initially, the raw FASTQ data files were converted into a format compatible with the QIIME2 software using the qiime tools import function. Each sample’s demultiplexed sequences underwent a stringent quality filtration process, trimming, denoising, and merging. The QIIME2 dada2 plugin facilitated the identification and removal of chimeric sequences, culminating in the creation of an amplicon sequence variant (ASV) feature Table [21]. For taxonomic classification, the QIIME2 feature-classifier plugin was utilized to map the ASV sequences against the refined GREENGENES 13_8 99% database, focusing on the V3V4 region delineated by the 341 F/806R primer pair. This alignment yielded a comprehensive taxonomy Table [22]. To discern variations in bacterial abundance across samples and groups, the Kruskal Wallis test was applied. Subsequent analysis incorporated LEfSe, employing linear discriminant analysis (LDA) to discern species with significant abundance differences, as indicated by the LDA score [23]. Diversity indices were computed using QIIME2’s core-diversity plugin. Beta diversity was assessed and graphically represented through principal coordinate analysis (PCoA) and nonmetric multidimensional scaling (NMDS) [24]. The “plsda” function from the R package “mixOmics” was employed for partial least squares discriminant analysis (PLS-DA) [25]. Redundancy analysis (RDA), facilitated by the R package “vegan”, was conducted to analyze the relative abundances of microbial species across various taxonomic levels [26]. Spearman’s rank correlations were calculated to identify co-occurrence patterns among dominant taxa, visualized through network plots. The potential functional profiles of microbial communities, based on KEGG Ortholog (KO), were predicted using PICRUSt [27]. Unless otherwise specified, the analysis parameters were set to their default values.

Metabolomic analysis

Ultra-high performance liquid chromatography-tandem mass spectrometry (UHPLC-MS/MS) was conducted utilizing a Vanquish UHPLC system from ThermoFisher, Germany, interfaced with an Orbitrap Q ExactiveTM HF mass spectrometer, also from Thermo Fisher, Germany. The samples were loaded onto a Hypersil Gold column (100 × 2.1 mm, 1.9 μm) and subjected to a linear gradient elution over a 17-minute period at a flow rate of 0.2 mL/min. For the positive ionization mode, the mobile phases consisted of eluent A (water containing 0.1% formic acid) and eluent B (methanol). Conversely, for the negative ionization mode, eluent A was a solution of 5 mM ammonium acetate adjusted to pH 9.0, and eluent B remained methanol. The gradient elution profile was as follows: starting at 2% B for 1.5 min, increasing to 85% B over 3 min, ramping up to 100% B over 10 min, then decreasing back to 2% B by 10.1 min, and finally maintaining 2% B for 12 min. The Orbitrap Q ExactiveTM HF mass spectrometer was configured to operate in both positive and negative polarity modes. The operational parameters included a spray voltage of 3.5 kV, a capillary temperature set at 320 °C, sheath gas flow rate at 35 psi, auxiliary gas flow rate at 10 L/min, an S-lens RF level of 60, and an auxiliary gas heater temperature of 350 °C.

The raw data acquired from the UHPLC-MS/MS system was meticulously processed using Compound Discoverer 3.1 (CD3.1, ThermoFisher). This software facilitated peak alignment, peak picking, and quantification for each detected metabolite. Post-processing, the peak intensities were standardized against the total spectral intensity to ensure consistency across samples. This normalized data was instrumental in deducing the molecular formulae of the metabolites, considering additive ions, molecular ion peaks, and fragment ions. The matching of these peaks was crucial for achieving accurate qualitative and relative quantitative outcomes. For statistical analysis, a combination of software tools was employed, including R (version R-3.4.3), Python (version 2.7.6), and the CentOS operating system (release 6.6). Metabolites with a coefficient of variation (CV) in the relative peak areas of quality control (QC) samples exceeding 30% were excluded to maintain data integrity. This rigorous approach led to the identification and relative quantification of the metabolites. Annotation of the metabolites was carried out utilizing several reputable databases: the KEGG database (https://www.genome.jp/kegg/pathway.html) for pathway information, the HMDB database (https://hmdb.ca/metabolites) for detailed metabolite profiles, and the LIPIDMaps database (http://www.lipidmaps.org/) for lipid-related data. These databases provided a comprehensive resource for the accurate annotation and understanding of the metabolites detected in the study.

Statistical analysis

Continuous variables are presented as mean ± standard deviation (SD) or median (quartile), along with group comparison by Student’s t-test or Wilcoxon rank-sum test. Categorical variables were compared using $\:{\chi\:}^{2}$ test or Fisher’s exact test. For data involving more than two groups, we used analysis of variance (ANOVA), Kruskal Wallis test, and mixture effect model for group comparisons. A SAS macro was developed based on the lifetest procedure to perform log-rank test and Kaplan-Meier estimation for time-to-event endpoint in the follow-up. We used generalized estimating equations (GEE) models to analyze the relationship between dynamic biochemical indicators during one-year follow-up and baseline omics (microbiota and metabolomics) biomarkers of patients. Additionally, we included information such as sex, age, and BMI in the model for covariate adjustment, with p-values adjusted by false discovery rate (FDR) correction.

Because the serum creatinine (Scr) level can effectively reflect patients’ renal function status, we utilized a Latent Class Growth Mixed Model (LCGMM) to fit Scr trajectories of patients based on one-year Scr measurement. LCGMM divides heterogeneous populations by estimating latent classes and fits individual curves through a mixed effects model. Longitudinal Scr measurement was treated as linear or nonlinear functions of time (the number of days between each measurement date and the enrollment date), represented by time, the squared term of time, and explored over 2–4 potential latent groups. The optimal number of groups and best-fitting shapes were determined based on the Bayesian Information Criterion, ensuring an acceptable proportion and posterior probability for each group. Additionally, we developed a predictive model for patients’ one-year creatinine values based on their clinical baseline data and pre-selected omics biomarkers associated with Scr. During the analysis process, we employed a logarithmic transformation method to normalize some of the measured values, ensuring they conform to statistical requirements. These statistical analyses were conducted using R (version 4.3.2) and SPSS (version20.0, Chicago, IL, USA).

Results

Study population and clinical parameter-based evaluation of renal transplant recipients

A total of 30 allograft renal transplant recipients, with an average transplantation duration of 7.77 years, were enrolled in this study (Fig. 1A). Fecal samples were collected within 24 to 72 h after the recipients reported the COVID-19-related symptoms for 16 S rRNA microbiome and untargeted metabolomics analysis at the onset of infection (Fig. 1A). The cohort included 17 recipients in the COVID-19 pneumonia group, six in the COVID-19 group and seven in the Control group; and gut microbiota composition and metabolites identified by gut microbiome and metabolome detection were further explored across the three groups, followed by the association analysis of Scr fluctuations during the one-year follow-up and the establishment of a prediction model for renal allograft prognosis (Fig. 1B).

Fig. 1

Study design and allograft function alteration in renal transplant recipients with COVID-19 pneumonia during one-year follow-up. (A) Timeline of one-year follow-up in 30 recipients. 0 represents the day recipients underwent SARS-CoV-2 infection; positive and negative number represent the day after or before SARS-CoV-2 infection; (B) Flow diagram of this study; (C) Fluctuation of serum creatinine (Scr) in recipients across three groups during one-year follow-up; (D-F) Time-to-event analysis of renal transplant recipients in three groups based on the principle of Scr declined more than 10% (D), more than 15% (E), and more than 20% (F), respectively

The basic demographics of the recipients in the three groups are presented in Table 1. Briefly, no significant differences in most of clinical variables were observed among the groups, except for age (P = 0.002), sex (P = 0.048), and COVID-19 vaccination status (P = 0.025). Recipients with COVID-19 pneumonia reported significantly more symptoms of cough and expectoration (Cough, P = 0.002; Expectoration, P = 0.040). Furthermore, we carried out a one-year outpatient follow-up for enrolled renal transplant recipients to explore the prognosis of allograft function (Fig. 1A). One recipient from the COVID-19 pneumonia group died within one month of SARS-CoV-2 infection, and one recipient from each of the COVID-19 pneumonia and COVID-19 groups was lost to follow-up. Laboratory examination comparisons were presented in Supplemental Table 1, with six parameters (i.e., Scr, lymphocyte counts, neutrophil counts, eosinophil counts, serum albumin, and serum glucose) showing significant differences among the three groups during the follow-up. Notably, Scr level remained stable in recipients from the COVID-19 and Control groups, but dramatically declined in those from the COVID-19 pneumonia group (P = 0.0018, Fig. 1C). Further stratification based on Scr decline showed a remarkable decrease in renal allograft function in the COVID-19 pneumonia group across all stratifications (Fig. 1D-F).

Table 1 Clinical and demographic features of allograft renal transplant recipients at the beginning of infection

Variation of gut microbiota compositions in renal transplant recipients and its association with baseline clinical biomarkers

Initially, we compared gut microbiota composition across the three groups. At the phylum level, the top 20 average relative abundance of microbial phyla detected in patients were summarized in Fig. 2A. Four phyla (i.e., Deferribacteres, TM7, Fusobacteria, and Gemmatimonadetes) were significantly different among the groups (P_FDR<0.05; Supplemental Table 2). Besides, Deferribacteres were characteristic of the Control group (LDA = 2.29, P < 0.001), while Fusobacteria (LDA = 3.65, P < 0.001), Gemmatimonadetes (LDA = 2.34, P < 0.001), Proteobacteria (LDA = 4.33, P = 0.028), and Actinobacteria (LDA = 3.90, P = 0.044) were characterized in the COVID-19 groups, and no characteristic phylum in the COVID-19 pneumonia group (Fig. 2B).

Fig. 2

Alternations in gut microbiota composition from renal transplant recipients among three groups. (A) Average relative abundance of microbial phyla detected in fecal samples from three groups (COVID-19 pneumonia, COVID-19 and Control group); (B) Characteristic microbiota phyla identified by linear discriminant analysis (LDA) analysis; (C) Average relative abundance of microbiota genera detected in fecal samples from three groups; (D) Characteristic microbiota genera identified by LDA analysis; (E) Diversity comparison of gut microbiota in three groups by Alpha diversity analysis; (F) Beta diversity analysis of microbial community composition in samples from three groups using Partial Least Squares Discriminant Analysis (PLS-DA); (G) Association of gut microbiota composition alteration and blood inflammatory markers by Canonical Correspondence Analysis (CCA) analysis; (H) Association of gut microbiota composition alteration and renal allograft function biomarkers by CCA analysis. *P < 0.05; **P < 0.01; ***P < 0.001

At the genus level, the top 20 relative abundance of microbial genus among the three groups were shown in Fig. 2C, with 13 of genera remarkably prevalent (P_FDR < 0.05; Supplemental Table 3). Nine genera were characteristic of the COVID-19 pneumonia group, five in the COVID-19 group, and seven in the Control group (LDA > 4; Fig. 2D). Further at the diversity level, we found that the alpha diversity of gut microbiota was significantly lower in the COVID-19 pneumonia group compared to the COVID-19 (P = 0.036) and Control (P = 0.0068) groups, but no difference between the COVID-19 and Control groups (Fig. 2E). Similar distributions were observed in beta diversity (Fig. 2F).

Subsequently, we analyzed relationships between gut microbiota genera and baseline clinical biomarkers across the three groups. Significant associations were observed with six baseline blood inflammation markers [i.e., blood white cells count (WBC), neutrophil count (NE), lymphocyte count (LY), monocyte count (MN), C-reactive protein (CRP), procalcitonin (PCT); P = 0.021], with the top 30 gut genera showing dominant impacts (Fig. 2G). Similar associations were observed with four markers of allograft function [i.e., Scr, blood urea nitrogen (BUN), estimated glomerular filtration rate (eGFR), urine protein (Upro); P = 0.023; Fig. 2H].

Relationship of gut microbiota composition changes and corresponding derived metabolites in renal transplant recipients

At the gut microbial metabolic level, we identified a total of 1348 metabolites in fecal samples of renal transplant recipients (Supplemental Table 4). Figure 3A shows the top 20 gut metabolites across the three groups. Globally, we observed three distinct clusters in accordance with the three groups (Fig. 3B). Among 1348 gut metabolites, 52 were significantly expressed across the groups (P_FDR < 0.05; Fig. 3C, Supplemental Table 4) and enriched in KEGG pathways (P_FDR < 0.05; Fig. 3D).

Fig. 3

Alternations in gut metabolites from renal transplant recipients among three groups. (A) Top 20 altered gut metabolites among three groups presented by biological function; (B) Multivariate statistical analysis of PLS-DA model was applied to examine the difference of gut metabolites among three groups; (C) 52 differential metabolites of 1348 total gut metabolites across three group; (D) metabolic enrichment analysis based on KEGG Compound Reaction Network; (E) Spearman’s correlation analysis of gut metabolites and microbiota

Correlation analysis of 13 genera and 52 metabolites that differed across three groups revealed significant correlation pairs among eight altered microbiota genera and seven differential metabolites (Fig. 3E). Interestingly, we observed that mycophenolic acid (MPA), a key immunosuppressive agent in renal transplantation [28], was significantly associated with seven gut microbiota genera.

Effect of altered gut microbiota and metabolites on one-year allograft function in renal transplant recipients with COVID-19

To investigate the changes in renal function of recipients after COVID-19 infection, we conducted a trajectory analysis on Scr fluctuations during follow-up (Fig. 4A). According to the selection criteria, a quadratic curve with two potential groups was optimal for Scr (Supplemental Table 5), which stratified the Scr trajectory as low-stable and high-stable (Fig. 4B). The low-stable group (n = 20) consisted of recipients who remained relatively stable allograft function during follow up, while the high-stable group (n = 7) started with elevated Scr and experienced fluctuations (Fig. 4C-D).

Fig. 4

Altered gut microbiota and metabolites in one-year allograft function in renal transplant recipients with COVID-19. (A). Scr fluctuation of each recipient from three groups during one-year follow-up; (B-C). A two-class model (B) with a quadric slope for Scr alternation (C) by LCGMM analysis; (D). Presentation of recipients from three groups and their corresponding class; (E). Top 20 metabolic pathways selected by KEGG metabolic enrichment analysis; (F). Scatter plot of model fitting values and true values by multivariate linear regression modeling

Subsequently, we used GEE to investigate the alteration of the gut microbiota and metabolites alongside the change of renal function. For gut microbiota, 12 differential microbial genera were found to be associated with Scr fluctuation during the one-year follow-up (P_FDR < 0.01; Supplemental Table 6). Notably, the genus of Clostridium contributed to the progression of COVID-19 pneumonia and persistently influenced the one-year changes in allograft function. Additionally, 181 metabolites were statistically associated with Scr changes (P_FDR<0.05; Supplemental Table 7), which were mainly enriched in metabolic pathways of Steroid hormone biosynthesis, Arginine biosynthesis, Valine, leucine and isoleucine biosynthesis, Cysteine and methionine metabolism, Histidine metabolism (Fig. 4E).

Further, we constructed a Scr prediction model for the one-year post-infection. We conducted single-factor regression and LASSO screening on candidate markers within GEE for subsequent multivariate linear regression modeling. After variable selection and model comparison (Supplemental Fig. 1; Supplemental Tables 8–9), we developed a prediction model which performed high predictive value for long-term graft function (Table 2; Fig. 4F).

Table 2 Results of multivariate linear regression modeling using four significant metabolites and clinical demographics

Discussion

To the best of our knowledge, this is the first clinical study to characterize the gut microbiota and metabolites alternations in renal transplant recipients with COVID-19 pneumonia. In this study, we observed significant impairment in allograft function in recipients with COVID-19 pneumonia during a one-year follow-up. We characterized substantial alternations in gut microbiota across three groups of COVID-19 pneumonia, COVID-19, and no-infection, which were significantly associated with inflammatory response and allograft function at baseline. Moreover, changes in gut metabolites were characterized, with lipid metabolism-related metabolites induced by critical microbiota. The LCGMM analysis revealed two distinct trajectories of one-year allograft function, and then GEE models indicated 12 genera and 181 differential gut metabolites between the two trajectories. Finally, a Scr prediction model harboring two gut metabolites, 1-Caffeoylquinic Acid and PMK, were found to have high predictive value for one-year allograft function.

The gut microbiota is crucial for regulating colonic ACE2 expression, which coronaviruses use to enter host cells, and their alterations could be associated with inflammatory cytokines production, thus contributing to the gut-lung axis pathology and the severity of COVID-19 pneumonia [29,30,31]. Additionally, the potential mechanism of intestinal microbiota imbalance affecting the severity of COVID-19 infection included conditional pathogenic bacteria are recognized by innate immune receptors, conditional pathogenic bacteria and toxin translocation, and the reduction of symbiotic bacteria inhibits the recruitment of immune cells [15]. Our findings suggested that gut microbiota may modulate host immune response and then potentially influence disease severity and renal allograft prognosis. Gut dysbiosis, characterized by decrease in Deferribacteres, TM7, and Fusobacteria, with an increase in Gemmatimonadetes, was observed in recipients with COVID-19 pneumonia, aligning with previous studies that Fusobacteria phylum was negatively associated with SARS-CoV-2 infection [32, 33].

The baseline alteration of the gut microbiome is likely to be important for the prognosis of allograft function during the one-year follow-up recovery period from COVID-19 disease in renal transplant recipients. Modifying gut microbiota with antibiotics administered after severe ischemic kidney injury in mice was reported to accelerate the recovery of kidney function [34, 35]. Gut dysbiosis also directly or indirectly contributed to the progression from acute kidney injury to chronic kidney diseases (CKD) [36]. Similarly, we observed that renal allograft function during COVID-19 pneumonia was associated with gut dysbiosis among the three groups; notably, the impact of which persists in the follow-up of recipients with COVID-19 pneumonia, indicating the crucial role of gut dysbiosis in allograft function. Particularly, Peptostreptococcaceae identified by the GEE model may be responsible for the kidney-gut-microbiota interaction as previous studies reported that Peptostreptococcaceae was significantly enriched in patients with CKD receiving hemodialysis or post-transplant diabetes mellitus recipients [37,38,39]. Besides, the kidney-gut-microbiota interaction has a bidirectional influence between kidney disease and gut dysbiosis, with gut dysbiosis potentially potentiating systemic inflammation through inflammatory cell accumulation, cytokines production, and T helper imbalance [40, 41]. Specifically, the gut microbiota affects kidney health by regulating the host’s immune response. Dysbiosis of gut microbiota can exacerbate immune imbalance, promote the production of pro-inflammatory cytokines, trigger systemic inflammatory responses, and accelerate the progression of kidney disease and related cardiovascular complications [42]. Further functional study is required to investigate the biological mechanism of Peptostreptococcaceae in the renal allograft function.

Numerous studies have pointed out the relationship between lipid metabolism and gut microbiota, especially in obese patients [43, 44]. Likewise, we found a significant relationship between gut microbiota composition alterations during COVID-19 pneumonia and lipid metabolism involved in gut metabolites. Furthermore, blood lipid of TG, HDL, and LDL were significantly correlated with gut microbiota alterations in the duration. It is well recognized that short-chain fatty acids, bacterial metabolites derived from fermentation of fibers in the colon, are important for host metabolism and are used as substrates for energy production, lipogenesis, gluconeogenesis and cholesterol synthesis [45,46,47]; as well, perillartine, a compound that may regulate lipid metabolism and lipid transport in hepatocytes [48], was found in high concentration in the intestines of recipients with COVID-19. These evidence suggested that lipid metabolism relevant gut metabolites may accelerate the progression of COVID-19 pneumonia. Notably, MPA, a widely applied immunosuppressive agent, was characterized to be significant correlated with differential gut microbiota of Enterococcus, Streptococcus, Lactobacillus, Faecallibacterium, Bifidobacterium, Roseburia, and Blautia in this study. Gibson et al. [49] has summarized the association between immunosuppressive agents with the changes in gut flora, and combinations of MPA and tacrolimus associated with an increase in colonization of Escherichia coli and Enterococcus sp in solid organ transplantation was observed, which was consistent with our results. Considering the immunosuppressive efficacy, MPA fluctuation in gut could directly result in the changes of peripheral MPA concentrations, contributing to the imbalance of immune status in recipients and following infection complications, such as COVID-19.

The potential role of gut microbiota and metabolites in recipients with COVID-19 could emerge the application of a microbiome-based risk profile to identify the individual with higher risk of impaired immune response [50]. In this study, we developed a gut metabolites-based model to predict disease severity and one-year allograft function, offering novel insights for transplant clinicians to identify recipients at risk of severe infectious diseases and poor prognosis early. Among the two gut metabolites included in this predictive model, 1-Caffeoylquinic Acid and PMK, 1-Caffeoylquinic Acid demonstrated efficacy as an inhibitor against target molecules, acting as the promising chemotherapeutic molecules against the NF-κB precursor protein and prevention of breast cancer [51] and antifungal activity [52].

To be noted, several limitations have been identified in this study. Firstly, the Peptostreptococcaceae was identified to be closely related with the long-term change of renal allograft function, while the biological function still remains to be explored in vitro. Then, the primary shortage of this study is the limited case number across three groups, which restricted the reliability of the conclusion. In fact, we have enrolled all recipients who met inclusion criteria as much as possible without sample calculation. Our further prospective study would perform the sample size calculation to test the alternation of gut microbiota and metabolites, as well as the prediction model. Finally, limited to the case number, potential confounding factors contributing to the strong relationship of four significant phyla with the one-year allograft function should be further determined.

In conclusion, it is the first time to profile the altered characteristics of the gut microbiota and metabolites in allograft renal transplant recipients with COVID-19, including four phyla (i.e., Deferribacteres, TM7, Fusobacteria, and Gemmatimonadetes) and 52 metabolites significantly enriched in recipients with COVID-19 pneumonia. Moreover, altered gut microbiota and metabolites were observed to be closely correlated with one-year allograft function during follow-up, based on which a Scr prediction model was established for one-year allograft prognosis. These findings would enhance our understanding of the molecular mechanisms in immunosuppressed renal transplant patients against virus infection.