Integrated transcriptomics, proteomics, and metabolomics to uncover inflammatory-metabolic crosstalk and key molecules in a mouse model of hyperoxaluria-induced kidney injury

Introduction

Hyperoxaluria is widely recognized as a major contributor to renal injury and the progression to end-stage kidney disease (ESKD) (1-3). This disorder is defined by elevated urinary oxalate excretion and encompasses two primary subtypes: primary hyperoxaluria and secondary hyperoxaluria (4). Primary hyperoxaluria arises from rare inherited defects in glyoxylate metabolism, which trigger excessive hepatic oxalate production. In contrast, secondary hyperoxaluria is predominantly driven by increased dietary intake of oxalate or its precursors, alongside alterations in the intestinal microflora that modulate oxalate handling.

The pathological consequences of hyperoxaluria extend beyond the formation of recurrent calcium oxalate (CaOx) kidney stones and nephrocalcinosis. Persistent oxalate overload further induces renal tubular damage, interstitial inflammation, and progressive fibrosis—pathological processes that collectively culminate in irreversible renal dysfunction and eventual ESKD. For primary hyperoxaluria management, FDA-approved small interfering RNA (siRNA) therapies that suppress glycolate oxidase expression, along with pyridoxine supplementation, are employed to curb hepatic oxalate synthesis (5,6). In enteric forms of secondary hyperoxaluria, clinical approaches include adherence to an oxalate-restricted diet and calcium supplementation to limit intestinal oxalate absorption (7). However, a significant therapeutic gap persists: no established treatments exist to mitigate the hyperoxaluria-induced renal inflammation and fibrosis that are tightly linked to the development of renal failure.

Multi-omics techniques help uncover disease mechanisms. A combined transcriptomic and proteomic study has identified altered genes and proteins in rats with kidney stones (8). However, integrated transcriptomic, proteomic, and metabolomic analysis is lacking for mouse models. The molecular basis of hyperoxaluria induced-kidney injury also remains unclear.

In this study, transcriptomic, proteomic, and untargeted metabolomic datasets were integrated and analyzed in a mouse model of hyperoxaluria-induced kidney injury. We aimed to: characterize gene expression and metabolic changes in hyperoxaluria-induced kidney injury, identify enriched biological pathways and key molecular players, and explore crosstalk between transcriptional and metabolic networks. This multi-omics approach offers new insights into the mechanisms of hyperoxaluria-induced kidney injury. We present this article in accordance with the ARRIVE reporting checklist (available at https://tau.amegroups.com/article/view/10.21037/tau-2025-534/rc).

Methods

Experiment design

Male C57BL/6J mice (8 weeks old, 23–25 g) were used. We randomly divided them into two groups: a control group (Control, n=6) and a glyoxylate group (Gly, n=6). All mice lived in controlled conditions. These included a 12-hour light-dark cycle, 22±1 ℃ temperature, and 40–50% humidity. Standard food and water were always available. CaOx crystals were induced via daily intraperitoneal injection of glyoxylate (100 mg/kg) for 7 consecutive days, as modified from previous protocols (9). The Control group mice were administered equal volumes of saline. Seven days later, the mice were fully anesthetized using 4% isoflurane. We then euthanized them via cardiac puncture. Both kidneys were collected from each mouse. One kidney was placed in 4% paraformaldehyde. This fixed kidney was for histology. The other kidney was immediately frozen in liquid nitrogen. Frozen kidneys were used for transcriptomic and metabolomic analysis. A protocol was prepared before the study with registration in Changhai Hospital. All animal experiments were performed under a project license (No. CHEC2025-021) granted by the Ethics Committee of Changhai Hospital, in compliance with the institutional guidelines for the care and use of animals.

Histopathology examination

Fixation was performed on kidney tissue using paraformaldehyde. Subsequent processing included graded ethanol dehydration, embedding in paraffin, and preparation of 4 µm sections. CaOx crystal deposition was detected by Pizzolato staining, as previously described (10). Hematoxylin and eosin (HE) staining was used to examine morphological changes.

Histopathological evaluation included cell lysis, brush border loss, and tubular dilation. For the key outcome of tubular damage scoring, C.L., a trained examiner, was blinded to the group allocation. C.L. independently scored the images without knowledge of which experimental group each sample belonged to, thereby minimizing subjective bias: 0 (none), 1 (<25%), 2 (≥25 and ≤50%), 3 (>50 and ≤75%), and 4 (>75%) (11). Injury scores were averaged from three randomly selected images and analyzed using GraphPad Prism (v10.1.2, GraphPad Software, CA, USA).

Biochemical assay

Renal function was assessed by measuring blood urea nitrogen (BUN) and creatinine (CRE) levels. These measurements were performed with the Beckman AU5800 automated biochemical analyzer.

Transcriptome analysis

Total RNA isolation from kidney tissues was carried out using TRIzol reagent (Invitrogen, Carlsbad, USA). RNA integrity was evaluated using agarose gel electrophoresis, while concentration and purity were quantified with a NanoDrop 2000 spectrophotometer (Thermo Fisher, Waltham, USA). Only samples displaying an A260/A280 ratio exceeding 1.8 and an A260/A230 ratio above 2.0 were selected for subsequent use. Preparation of RNA-sequencing (RNA-seq) libraries was executed using the NEBNext^® Ultra^TM RNA Library Prep Kit for Illumina^® (NEB, Ipswich, USA). Library concentration was determined with a Qubit 4.0 Fluorometer (Thermo Fisher, USA), and quality evaluation was performed via an Agilent Bioanalyzer 2100 system. Sequencing procedures were conducted on the Illumina NovaSeq 6000 platform (2×150 bp, paired-end) in strict accordance with the manufacturer’s protocols.

The initial sequencing datasets underwent processing via Fastp (v0.22.0) to execute adapter trimming and elimination of low-quality reads (12). To verify data validity, quality indicators such as Q20, Q30, and GC content were computed. Gene expression levels were gauged using the metric of fragments per kilobase of transcript per million mapped reads (FPKM), with a comprehensive compilation of FPKM values documented in online table 1 (available at https://cdn.amegroups.cn/static/public/tau-2025-534-Onlinetable1.xlsx).

For the identification of differentially expressed genes (DEGs), the DESeq2 tool was employed. Genes were categorized as significant when exhibiting an adjusted P value below 0.05 and an absolute log₂ fold change exceeding 1. Functional enrichment analyses including Gene Ontology (GO) terms and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathways were performed utilizing the clusterProfiler package.

Proteome analysis

Kidney tissues were harvested from mice and immediately frozen in liquid nitrogen. For protein extraction, tissues were ground into fine powder using a mortar and pestle under liquid nitrogen, followed by lysis in SDT buffer [4% sodium dodecyl sulfate (SDS), 100 mM Tris-HCl, pH 7.6] supplemented with 1% protease inhibitor cocktail. The lysates were sonicated on ice (3 cycles of 30 s on/30 s off) and boiled at 95 ℃ for 15 min to denature proteins. After centrifugation at 14,000 ×g for 40 min at 4 ℃, the supernatant was collected, and protein concentration was determined using the BCA Protein Assay Kit (Bio-Rad, Hercules, USA) according to the manufacturer’s protocol.

Equal amounts of protein (20 µg per sample) were subjected to reduction with 10 mM dithiothreitol at 37 ℃ for 1.5 h, followed by alkylation with 20 mM iodoacetamide in the dark at room temperature for 30 min. Samples were then transferred to 10 kDa cutoff Microcon filter units and washed sequentially with UA buffer (8 M urea in 100 mM Tris-HCl, pH 8.5) and 25 mM ammonium bicarbonate buffer. Trypsin digestion was performed overnight at 37 ℃ with a trypsin-to-protein ratio of 1:50 (w/w), and peptides were collected as filtrates. Digested peptides were desalted using C18 SPE cartridges (Empore^TM, Sigma-Aldrich, St. Louis, USA), concentrated by vacuum centrifugation, and reconstituted in 0.1% formic acid. Peptide concentration was estimated by UV absorbance at 280 nm.

A pooled peptide sample (equal aliquots from all samples) was fractionated into 10 fractions using a High pH Reversed-Phase Peptide Fractionation Kit (Thermo Scientific^TM Pierce^TM) to enhance proteome coverage. Each fraction was desalted using C18 cartridges, reconstituted in 0.1% formic acid, and spiked with iRT calibration peptides (Biognosys, Schlieren, Switzerland) for retention time normalization.

Fractionated peptides were analyzed using a timsTOF Pro mass spectrometer (Bruker Daltonics, Bremen, Germany) coupled with an Evosep One ultra-high performance liquid chromatography (UHPLC) system (Evosep, Odense, Denmark). The liquid chromatography (LC) separation was performed on a reversed-phase column (25 cm × 75 µm i.d.) with a gradient of 6–35% solvent B (0.1% formic acid in acetonitrile) over 87 min at a flow rate of 450 nL/min. For data-dependent acquisition (DDA), the mass spectrometer was operated in parallel accumulation-serial fragmentation (PASEF) mode, with a mass range of m/z 100–1,700. Precursor ions with charge states 2–5 were selected for fragmentation, and dynamic exclusion was set to 30 s. Tandem mass spectrometry (MS/MS) spectra were acquired with 10 PASEF scans per cycle, and collision energy was adjusted based on ion mobility.

Individual peptide samples were analyzed in data-independent acquisition (DIA) mode using the same timsTOF Pro and Evosep One system. The LC gradient was optimized to 60 min, with solvent B ranging from 6% to 35%. DIA scans covered the m/z range 100–1,700, with 4 ion mobility windows per trapped ion mobility spectrometry (TIMS) cycle (100 ms per window). Collision energy was linearly ramped from 20 to 59 eV based on ion mobility (1/K0 =0.85–1.30 Vs/cm²). iRT peptides were spiked into each sample to ensure retention time consistency.

DDA data were processed using Spectronaut^TM 14 (Biognosys) to generate a spectral library, searching against the UniProt database (Mus musculus, accessed July 2023) with added iRT sequences. Parameters included trypsin as the cleavage enzyme (max 1 missed cleavage), carbamidomethylation (C) as a fixed modification, and oxidation (M) and N-terminal acetylation as variable modifications. False discovery rate (FDR) was controlled at <1% for both peptides and proteins. DIA data were quantified using the same software with the generated library, enabling dynamic iRT correction and cross-run normalization. Final results were filtered to FDR <1% (Q value <0.01).

Subcellular localization of differentially expressed proteins (DEPs) was predicted using CELLO (http://cello.life.nctu.edu.tw/). Functional annotation was performed via Blast2GO, incorporating GO terms for biological processes (BP), cellular components (CC), and molecular functions (MF). KEGG pathway enrichment was analyzed using KOBAS software, with significance determined by Fisher’s exact test and Benjamini-Hochberg correction (adjusted P<0.05).

Metabolome analysis

Frozen kidney tissues (100 mg) were finely minced on dry ice and homogenized using 200 µL of water in conjunction with ceramic beads. Metabolite extraction was performed with 800 µL of a methanol/acetonitrile mixture. The dried samples were reconstituted in 100 µL of an acetonitrile/water solution, subjected to a subsequent centrifugation step, and prepared for LC-MS analysis.

Metabolite separation was achieved via UHPLC using a BEH Amide column. The mobile phases consisted of 25 mM ammonium acetate/ammonium hydroxide (A) and acetonitrile (B) under gradient elution conditions. Mass spectrometry (MS) analysis was conducted on a TripleTOF 6600 (AB Sciex) instrument in both positive and negative electrospray ionization modes. MS/MS data (m/z 25–1,000) were acquired using information-dependent acquisition with a collision energy of 35±15 V.

MzXML format conversion for raw MS files was executed via ProteoWizard MSConvert (v3.0), followed by processing in XCMS. CentWave was employed for peak detection. CAMERA was utilized to annotate isotopes and adducts. Features were retained only when detected in over 50% of samples within at least one group. Metabolite identification relied on matching accurate mass (within 10 ppm) and MS/MS spectra against an in-house authentic standards database.

Data normalization was conducted using sum intensity, with analyses performed in R via the ropls package. Multivariate approaches involved Pareto-scaled principal component analysis (PCA) and orthogonal partial least-squares discriminant analysis (OPLS-DA). Model quality was validated through 7-fold cross-validation and response permutation testing. Variable importance in projection (VIP) scores from OPLS-DA ranked feature contributions, while Student’s t-test identified group disparities. Metabolites with VIP >1 and P<0.05 were deemed significantly altered. Pearson correlation analysis explored variable relationships.

Integrative transcriptome and metabolome analysis

To investigate the association between gene expression and metabolite profiles, integrative analysis was conducted using the OmicsPLS R package. Bidirectional orthogonal partial least-squares (O2PLS) analysis was employed to detect shared variation between transcriptomic and metabolomic datasets (13). Pearson’s correlation analysis was performed with the igraph R package to construct an omics interaction network. Subsequently, Spearman correlation was applied to identify key nodes, emphasizing gene–metabolite pairs with absolute correlation coefficients of at least 0.8. The resulting correlation network was visualized and analyzed using igraph in R version 4.2.3.

RT-qPCR validation

To corroborate RNA-seq results, six DEGs (Cxcl2, Havcr1, Slc22a7, Areg, Retnla, and Lcn2) were randomly chosen for quantitative real-time polymerase chain reaction (qRT-PCR) analysis. Total RNA extraction from kidney tissues was performed using TRIzol reagent (Invitrogen), followed by complementary DNA (cDNA) synthesis with a reverse transcription kit (Vazyme, China). qRT-PCR was performed with a LightCycler 480 System (Roche, Basel, Switzerland). Gene expression levels were normalized to Gapdh and quantified via the 2^−ΔΔCt method. Primer sequences are detailed in Table S1.

Statistical analysis

All statistical analyses were performed using SPSS 26.0 (SPSS Inc., Chicago, IL, USA) and GraphPad Prism (v10.1.2). Data are presented as mean ± standard deviation. First, the Kolmogorov-Smirnov test was used to assess the normality of data distribution. For normally distributed variables, differences were analyzed using either a two-tailed Student’s t-test (for two groups). For non-normally distributed variables, the Mann-Whitney U test was applied to compare differences between two groups. A P value <0.05 was considered statistically significant.

Results

Establishment of a mouse model of kidney injury induced by hyperoxaluria

A mouse model of kidney injury induced by hyperoxaluria was established through intraperitoneal glyoxylate injection. Kidney samples were then subjected to transcriptomic and metabolomic analyses (Figure 1A). Histological assessment confirmed kidney injury in the model group. HE staining revealed typical pathological changes, including tubular dilation, and cell lysis (Figure 1B). Pizzolato staining further demonstrated marked CaOx crystal deposition in renal tissues (Figure 1C). Compared with the control mice, model mice exhibited significantly elevated serum CRE and BUN levels (Figure 1D). These findings confirmed the successful induction of hyperoxaluria-induced kidney injury in mice.

Figure 1 Establishment and validation of the glyoxylate-induced kidney stone mouse model. (A) Schematic diagram of the experimental design for transcriptomic and metabolomic analysis in kidney tissues. (B) HE staining images from the control group and the Gly group. Original magnification, × 4 and ×20. (C) Pizzolato staining for CaOx crystal deposition in kidney tissues. Red arrow: intratubular crystals; Orange arrow: interstitial crystals. Original magnification, × 4 and ×20. (D) Serum levels of BUN and CRE in the Control group and the Gly group. N=5. BUN, blood urea nitrogen; CaOx, calcium oxalate; CRE, creatinine; Gly, glyoxylate; HE, hematoxylin and eosin; qPCR, quantitative polymerase chain reaction.

Transcriptomic profiling of the kidney in hyperoxaluria mice

We first examined global gene expression patterns across all samples. Overall expression profiles were consistent within each group (Figure S1A). Pearson correlation analysis confirmed strong reproducibility among biological replicates (Figure S1B). PCA showed tight clustering within groups and clear separation between control and model mice (Figure S1C). These results confirmed robust transcriptional changes following glyoxylate treatment.

In total, 5,099 DEGs were identified. Among them, 2,546 genes were upregulated and 2,553 were downregulated (Figure 2A,2B). Il36a emerged as the most significantly upregulated gene, although its role in nephrolithiasis remains unknown [online table 2 (available at https://cdn.amegroups.cn/static/public/tau-2025-534-Onlinetable2.xlsx)]. Several known markers of renal injury and inflammation—including Lcn2, Havcr1, and Ccl2—also showed marked upregulation, indicating substantial kidney damage at the molecular level. These findings suggest that glyoxylate-induced hyperoxaluria trigger widespread transcriptional reprogramming, especially in inflammatory and injury-related pathways.

Figure 2 Overview of transcriptomic analysis. (A) Bar chart, and (B) volcano plot of DEGs in the kidney of calculi mouse. (C) GO terms enrichment of all DEGs. (D) GO enrichment analysis for up-regulated and down-regulated DEGs. (E) KEGG enrichment analysis of all DEGs. (F) KEGG enrichment analysis for up-regulated DEGs. (G) KEGG enrichment analysis for down-regulated DEGs. BP, biological process; CC, cellular component; Gly, glyoxylate; GO, Gene Ontology; KEGG, Kyoto Encyclopedia of Genes and Genomes; MF, molecular function.

GO enrichment analysis of all DEGs highlighted significant BP (Figure 2C). Enriched terms in the BP category included “response to external stimulus”, “immune response”, and “defense response” reflecting activation of systemic and local immune defenses. In the MF category, enriched terms including “oxidoreductase activity”, “transmembrane transporter activity”, and “signaling receptor binding” suggested modified enzymatic and receptor functionalities in the kidney. Within the CC category, enrichment for “extracellular space”, “cell periphery”, and “extracellular region” highlighted alterations in extracellular structural organization and membrane architecture. Further analysis separated upregulated and downregulated genes (Figure 2D). Upregulated genes were primarily associated with “cell periphery” and “extracellular region”, suggesting roles in membrane remodeling and extracellular matrix changes during crystal formation. In contrast, downregulated genes were enriched in metabolic processes, particularly “small-molecule metabolic process”, indicating suppressed renal metabolic activity in response to crystal-induced injury.

KEGG pathway enrichment analysis for all DEGs uncovered notable activation of immune and metabolic pathways (Figure 2E). Pathways showing significant enrichment encompassed “metabolic pathways”, “cytokine-cytokine receptor interaction”, and “complement and coagulation cascades”. These results suggest that kidney crystal formation triggers a strong inflammatory response and disrupts normal metabolic processes. Activation of the complement and coagulation systems likely reflects tissue injury and immune activation. Meanwhile, dysregulation of cytokine signaling indicates altered immune communication, and changes in metabolic pathways suggest a disturbed renal metabolic state. Analysis of upregulated DEGs (Figure 2F) showed strong enrichment in immune-related pathways, including “cytokine-cytokine receptor interaction”, “TNF signaling pathway”, “IL-17 signaling pathway”, and “PI3K-Akt signaling pathway”. These pathways are involved in inflammation, immune cell recruitment, and tissue repair, reinforcing the idea that glyoxylate-induced injury drives a robust immune response. In contrast, downregulated DEGs (Figure 2G) were predominantly enriched in metabolism-related pathways. These included “metabolic pathways”, “valine, leucine, and isoleucine degradation”, and “oxidative phosphorylation”. The suppression of these pathways suggests impaired energy production and amino acid metabolism, reflecting compromised renal function.

Overall, the GO and KEGG enrichment analyses demonstrate that kidney stone formation induces broad transcriptional changes. These changes affect immune responses, membrane and extracellular processes, and core metabolic functions in the kidney.

Proteomic profiling of the kidney in hyperoxaluria mice

To further dissect the molecular alterations underlying glyoxylate-induced hyperoxaluria, we conducted proteomic analysis on kidney tissues from control and glyoxylate-treated mice. For protein identification, we got a total of 28,456 unique peptides, mapping to 5,687 proteins (Figure S1D). PCA of the proteomic data (Figure S1E) revealed clear separation between the control and glyoxylate groups. Tight clustering within each group indicated good reproducibility of the proteomic profiles, validating the reliability of subsequent analyses.

A total of 2,966 DEPs were identified, with 2,095 upregulated and 871 downregulated in the glyoxylate-treated mice compared to controls [Figure 3A; online table 3 (available at https://cdn.amegroups.cn/static/public/tau-2025-534-Onlinetable3.xlsx)]. Subcellular localization analysis of DEPs (Figure 3B) showed that the majority of altered proteins were localized in the nucleus (1,306 proteins) and cytoplasm (895 proteins), followed by extracellular space (388 proteins), mitochondria (513 proteins), and plasma membrane (431 proteins). This distribution suggests widespread dysregulation across multiple cellular compartments, with prominent changes in nuclear and cytoplasmic processes, potentially reflecting transcriptional and translational reprogramming, as well as cytoskeletal or metabolic perturbations during crystal formation. The volcano plot (Figure 3C) visualized these changes, with red dots representing upregulated proteins and dark blue dots representing downregulated proteins.

Figure 3 Overview of proteomic analysis. (A) Bar chart, (B) pie chart and (C) volcano plot of DEPs in the kidney of calculi mouse. (D) GO enrichment analysis for up-regulated and down-regulated DEPs. (E) KEGG enrichment analysis for up-regulated DEPs. (F) KEGG enrichment analysis for down-regulated DEPs. DEPs, differentially expressed proteins; GO, Gene Ontology; KEGG, Kyoto Encyclopedia of Genes and Genomes.

GO enrichment analysis of DEPs further categorized their biological functions (Figure 3D). For the upregulated proteins, they were primarily associated with terms related to nuclear functions, such as “nucleus”, “nucleic acid metabolic process”, “nuclear part”, “RNA metabolic process”, and “RNA processing”. These terms suggest that processes related to nucleic acid metabolism, RNA processing, and gene expression in the nucleus are significantly activated in the kidney of hyperoxaluria mice. For the downregulated proteins, they were enriched in functions like “monocarboxylic acid metabolic process”, “oxidoreductase activity”, “oxoacid metabolic process”, “carboxylic acid metabolic process”, “oxidation-reduction process”, and “mitochondrion”. These findings indicate that metabolic processes, especially those related to carboxylic acid and oxidation-reduction, as well as mitochondrial functions, are suppressed in the kidney during hyperoxaluria. This aligns with transcriptomic data, reinforcing the activation of nuclear-related and stress responses alongside suppression of core metabolic and mitochondrial functions in the kidney.

For downregulated proteins (Figure 3E), pathways such as “spliceosome”, “cell cycle”, “ATP-dependent chromatin remodeling”, “DNA replication”, “RNA degradation”, “NF-kappa B signaling pathway”, “osteoclast differentiation”, and “apoptosis” were also involved. These pathway analyses indicate that multiple cellular processes, including RNA processing, cell cycle, and immune-related pathways, are perturbed in the kidney during hyperoxaluria. For the downregulated proteins (Figure 3F), “oxidative phosphorylation” shows the most significant changes. Other pathways like “peroxisome”, “fatty acid degradation”, “citrate cycle (TCA cycle)”, “PPAR signaling pathway”, and “lipoic acid metabolism” were also enriched.

Collectively, the proteomic data reveals a coordinated shift in protein expression, with upregulation of proteins related to nuclear nucleic acid metabolism and gene expression, and downregulation of metabolic and mitochondrial-associated proteins. The KEGG pathway analyses further highlight the involvement of diverse cellular processes in hyperoxaluria. These changes complement the transcriptomic findings, reinforcing the interplay between inflammation, metabolic dysfunction, and cellular damage in hyperoxaluria.

Metabolomic characterization of the kidney in hyperoxaluria mice

PCA in both negative and positive ion modes clearly distinguished between control and glyoxylate-treated mice (Figure S2A,S2B). Samples from each group clustered closely, confirming data consistency and highlighting clear metabolic differences between normal and hyperoxaluria-affected kidneys. OPLS-DA further improved group separation (Figure S2C,S2D), indicating substantial metabolic shifts. Permutation tests (Figure S2E,S2F) showed high R²Y and satisfactory Q² values, validating the robustness and predictive performance of the OPLS-DA models. A pie chart of classified metabolites (Figure 4A) showed that most identified compounds fell into several major superclasses. Lipids and lipid-like molecules accounted for ~18.7%, organic acids and derivatives for ~18.6%, and organoheterocyclic compounds for ~15%. These results suggest broad metabolic reprogramming involving diverse biochemical categories.

Figure 4 Overview of metabolomic analysis. (A) Pie chart of DEMs categorized by superclasses. (B) Bar charts showing the number of up-regulated and down-regulated DEMs in negative (left) and positive (right) ion modes. (C) Top 15 up-regulated and down-regulated DEMs in negative ion modes. (D) Top 15 up-regulated and down-regulated DEMs in positive ion modes. (E) KEGG pathway enrichment of all DEMs. (F) Differential abundance scores of enriched pathways. DEMs, differentially expressed metabolites; FC, fold change; KEGG, Kyoto Encyclopedia of Genes and Genomes; NEG, negative; POS, positive.

In total, 255 differentially expressed metabolites (DEMs) were detected under negative ion mode—131 were upregulated and 124 downregulated. Under positive ion mode, 117 DEMs were identified—27 were upregulated and 90 downregulated [Figure 4B; online tables 4,5 (available at https://cdn.amegroups.cn/static/public/tau-2025-534-Onlinetables4-5.xlsx)]. Bar charts show detailed profiles of DEMs in both ion modes (Figure 4C,4D). Notably, in the negative ion mode, syringic acid, Pro-Thr, and N-propyl gallate were among the most upregulated metabolites, while rosiglitazone, D-glutamine, and adipic acid showed marked downregulation. In the positive ion mode, 2-pivaloyl-1,3-indandione was highly upregulated, whereas propranolol and caffeine were significantly reduced. These shifts point to disrupted amino acid metabolism, lipid regulation, and xenobiotic processing during kidney crystal formation.

KEGG pathway enrichment of all DEMs revealed several significantly affected metabolic routes (Figure 4E). Notable pathways included “cAMP signaling pathway”, “biosynthesis of amino acids”, and “beta-alanine metabolism”, all showing significant enrichment as indicated by high −log₁₀(P value) scores. These pathways suggest that kidney crystal formation disrupts essential signaling and biosynthetic functions.

Differential abundance scores (Figure 4F) were calculated to quantify the degree of alteration for each enriched pathway. Among them, the cyclic adenosine monophosphate (cAMP) signaling pathway and vitamin digestion and absorption pathways displayed prominent shifts, highlighting their relevance to the metabolic reprogramming observed in hyperoxaluria-affected kidneys. Together, these results emphasize the complexity of metabolic dysregulation in hyperoxaluria-induced kidney injury and point to specific pathways that may mediate disease progression.

Integrative transcriptome and proteome analysis

A total of 32,563 genes were identified in the transcriptomic profiling, while 6,225 proteins were detected in the proteomic analysis. Among these, 6,225 genes were quantified at both the mRNA and protein levels, representing the overlapping entities between the two omics layers (Figure 5A). This intersection provided a foundation for exploring the consistency and discrepancy in gene expression patterns across different regulatory levels. The Spearman correlation coefficient between transcript and protein expression ratios was calculated to assess the overall consistency. A moderate positive correlation (R=0.6778) was observed, indicating a partial alignment between mRNA and protein abundance changes (Figure 5B). The scatter plot further visualized this relationship. Based on the differential expression status, the overlapping genes were categorized into distinct groups (Figure 5C). A total of 375 genes exhibited concurrent upregulation at both mRNA and protein levels, while 418 genes were consistently downregulated. In contrast, 1,698 genes showed upregulation only at the protein level without significant mRNA changes, and 445 genes were downregulated specifically at the protein level. Additionally, 8 genes were upregulated at the mRNA level but downregulated at the protein level, and 22 genes displayed the opposite trend, highlighting the complexity of post-transcriptional regulation.

Figure 5 Integrated transcriptomic and proteomic analysis. (A) Venn diagram of DEGs and DEPs. (B) Quadrant plot of transcript-protein expression trend consistency. (C) Venn diagram of overlapping classes of DEGs and DEPs. (D) GO enrichment analysis of up-regulated DEGs and DEPs. (E) KEGG enrichment analysis of up-regulated DEGs and DEPs. (F) GO enrichment analysis of down-regulated DEGs and DEPs. (G) KEGG enrichment analysis of down-regulated DEGs and DEPs. BP, biological process; CC, cellular component; DEGs, differentially expressed genes; DEPs, differentially expressed proteins; GO, Gene Ontology; KEGG, Kyoto Encyclopedia of Genes and Genomes; MF, molecular function; NDEGs, non-differentially expressed genes; NDEPs, non-differentially expressed proteins.

GO enrichment analysis of the 375 consistently upregulated genes revealed significant enrichment in BP terms related to cell-substrate adhesion, actin filament organization, wound healing, and chemotaxis (Figure 5D). CC terms were dominated by collagen-containing extracellular matrix, cell leading edge, and endocytic vesicle, while MF terms highlighted actin binding, integrin binding, and GTPase activity. KEGG pathway analysis further indicated that these genes were primarily involved in regulation of “actin cytoskeleton”, “focal adhesion”, and “phagosome” (Figure 5E). Notably, pathways associated with infectious diseases (e.g., “salmonella infection”, “tuberculosis”) were also enriched, reflecting potential inflammatory responses.

For the 418 genes with consistent downregulation, BP terms were enriched in metabolic pathways, including fatty acid metabolism, purine nucleotide metabolism, and organic acid catabolism (Figure 5F). CC terms such as mitochondrial matrix and peroxisome were prominently represented, while MF terms focused on oxidoreductase activity, transmembrane transporter activity, and cofactor binding. KEGG analysis demonstrated that these downregulated genes were significantly involved in energy metabolism pathways (“oxidative phosphorylation”, “carbon metabolism”), and lipid metabolism (“fatty acid degradation”, “PPAR signaling pathway”) (Figure 5G). These findings suggest a broad suppression of metabolic processes and mitochondrial function in the model of hyperoxaluria-induced kidney injury.

In summary, the integrated transcriptomic and proteomic analysis uncovered a set of genes with consistent expression changes at both mRNA and protein levels, which are primarily involved in immune activation and metabolic regulation.

Integrative transcriptome and metabolome analysis

To examine the relationship between gene expression and metabolic changes, an integrated analysis was performed using Spearman correlation (|r|>0.8). Correlation networks (Figure 6) revealed strong interactions between key DEGs and DEMs. Notably, Areg, Cxcl2, Havcr1, Lcn2, Retnla, Il36a, Mmp7, and Slc22a7 emerged as central transcriptional nodes connected to metabolic alterations. On the metabolite side, compounds such as 1-(3,4-dichlorophenyl)-3,3-dimethylurea and 1-heptanoyl-sn-glycero-3-phosphorylcholine stood out as key nodes. These findings suggest a coordinated molecular response involving inflammation, injury, and metabolic stress in glyoxylate-induced kidney injury. Full network details are available in Tables S2-S4.

Figure 6 Gene expression and metabolite abundance correlation network diagram. Green squares: genes. Teal squares: metabolites. Red lines: positive correlations. Blue lines: negative correlations.

We performed O2PLS analysis to integrate transcriptomic and metabolomic datasets. The joint loading plot (Figure 7) illustrates the top 50 variables contributing to the shared variance between the two omics layers, providing insights into the core molecular networks associated with hyperoxaluria-induced kidney injury (Tables S5,S6).

Figure 7 O2PLS loadings plot of DEGs and DEMs. DEGs, differentially expressed genes; DEMs, differentially expressed metabolites; O2PLS, orthogonal partial least-squares.

Notably, 2,4-dichlorophenol, a xenobiotic compound linked to oxidative stress, and 1-(3,4-dichlorophenyl)-3,3-dimethylurea, a known metabolic disruptor, displayed high positive loading values, suggesting their potential involvement in renal pathophysiological processes. These metabolites clustered closely with transcripts such as Lcn2, a marker of renal injury, indicating a coordinated response to cellular stress. On the other end of the loading spectrum, lipids including 1-palmitoyl-2-docosahexaenoyl-sn-glycero-3-phosphocholine and cis-4,7,10,13,16,19-docosahexaenoic acid (DHA) showed negative loadings, implying dysregulated lipid metabolism in the kidney. These lipid species were spatially associated with transcripts like Havcr1, which is involved in tubular damage, highlighting a potential link between lipid peroxidation and renal tubular dysfunction.

Overall, the O2PLS model captures the interconnectedness of transcriptional regulation and metabolic shifts, with the joint loading patterns emphasizing the roles of oxidative stress, lipid metabolism dysfunction, and renal injury responses as central hubs in the molecular landscape of kidney crystal formation. These findings provide a multidimensional framework for understanding the pathophysiology of hyperoxaluria-induced kidney injury and identifying potential therapeutic targets.

Verification of the DEGs by qRT-PCR

To verify the accuracy of RNA-seq data, six DEGs (Areg, Cxcl2, Havcr1, Lcn2, Retnla, and Slc22a7) were selected from the correlation network for qRT-PCR validation. The expression patterns obtained by qRT-PCR (Figure 8) were consistent with the transcriptomic results, confirming the accuracy of the sequencing data and supporting the identified gene expression changes in the hyperoxaluria model.

Figure 8 qRT-PCR validation of DEGs. N=3. DEGs, differentially expressed genes; qRT-PCR, quantitative real-time polymerase chain reaction.

Discussion

In this study, we integrated transcriptomic, proteomic, and untargeted metabolomic analysis to systematically characterize the molecular perturbations in a mouse model of hyperoxaluria-induced kidney injury. Our multi-omics approach revealed coordinated changes in gene expression, protein abundance, and metabolite profiles, shedding new light on the pathophysiological processes governing crystal formation, renal injury, and inflammatory responses.

Transcriptomic profiling identified 5,099 DEGs, with significant upregulation of inflammatory pathways such as cytokine-cytokine receptor interactions, tumor necrosis factor (TNF) signaling, and IL-17 signaling. These findings corroborate earlier studies highlighting the role of immune activation in hyperoxaluria (14,15). Notably, markers of renal injury, including Lcn2 and Havcr1, were highly upregulated. Downregulated genes were predominantly involved in metabolic processes, particularly oxidative phosphorylation and amino acid degradation, suggesting a shift in renal energy metabolism during crystal formation.

Of particular interest, both transcriptomic and proteomic analysis revealed significant downregulation of fatty acid metabolism pathways, including PPAR signaling, fatty acid degradation, and mitochondrial β-oxidation. This suppression of lipid metabolic processes aligns with recent studies demonstrating that impaired fatty acid oxidation contributes to renal tubular injury and crystal adhesion (16-18). Furthermore, peroxisomal dysfunction, reflected in reduced activity of fatty acid-processing enzymes, could disrupt lipid homeostasis, promoting lipid peroxidation and oxidative stress, which are known to accelerate CaOx crystal formation (19-21). The proteomic data further supported this observation, showing decreased expression of key enzymes involved in lipid catabolism, such as acyl-CoA dehydrogenases and carnitine palmitoyltransferases. Given that fatty acids serve as a major energy source for renal tubular cells, their dysregulated metabolism may exacerbate cellular stress and promote crystal retention, further aggravating crystal formation (18,22).

Metabolomic analysis revealed 372 DEMs, with significant disruptions in lipid metabolism, amino acid biosynthesis, and energy pathways. The downregulation of metabolites involved in the TCA cycle and oxidative phosphorylation aligns with the transcriptomic and proteomic data, reinforcing the notion of impaired mitochondrial function in crystal-forming kidneys. Notably, the cAMP signaling pathway emerged as a key player, consistent with its established role in renal tubular function and inflammation (23-25). Our metabolomic profiling uncovered profound disruptions in amino acid, lipid, and energy metabolism, which we propose are both a consequence and a driver of CaOx crystal deposition. Amino acid metabolism, in particular, showed widespread alterations: downregulation of branched-chain amino acid (BCAA) degradation pathways (e.g., valine, leucine, and isoleucine) and perturbations in glutamate and glutamine metabolism. BCAAs are critical for energy production and protein synthesis in renal tubules; their dysregulation may impair tubular repair capacity, rendering the kidney more susceptible to crystal-induced damage (26,27). Lipid metabolism was another key target of disruption, with significant changes in phospholipids and fatty acids. For instance, decreased levels of DHA, an anti-inflammatory omega-3 fatty acid, and increased pro-inflammatory lipid mediators (e.g., lysophosphatidylcholines) suggest a shift toward a pro-inflammatory lipid environment. This is consistent with transcriptomic data showing upregulation of genes involved in lipid peroxidation, further supporting a role for lipid-derived oxidative stress in tubular injury. Importantly, our integrated OPLS-DA analysis identified strong correlations between these lipid perturbations and inflammatory genes (e.g., Lcn2), highlighting a bidirectional relationship where lipid dysregulation amplifies inflammation.

Although this research provides significant insights into the molecular mechanisms of glyoxylate-induced hyperoxaluria, several limitations require acknowledgment. First, the study was conducted at a single time point, limiting insights into the temporal dynamics of molecular changes during crystal formation. Future time-course studies could clarify whether metabolic perturbations precede inflammatory activation. Second, the metabolomic analysis focuses solely on kidney tissue. Serum or urine metabolomics might reveal systemic metabolic shifts contributing to hyperoxaluria pathogenesis. Third, all findings are based on a mouse model. Species-specific differences in metabolism and immune responses may limit the direct applicability to human disease. Fourth, the functional roles of key molecules (e.g., Il36a) remain to be explored experimentally. Targeted interventions, such as inhibition of TNF signaling or modulation of cAMP pathways, could further elucidate their contributions to crystal formation. Lastly, although our study uses a glyoxylate-induced mouse model of CaOx crystal deposition, which differs from human kidney stone formation as mice do not naturally form CaOx stones, the molecular mechanisms we identified may partially resemble the early pathological processes of human kidney stone formation. For example, in human kidney stone formers, interactions between CaOx crystals and renal tubular cells can also trigger local inflammatory responses and metabolic changes, which may contribute to stone formation and progression. However, it should be noted that human kidney stone formation is a more complex process involving multiple factors such as genetic background, dietary habits and urinary environment. These factors are not fully reflected in our mouse model. Thus, the findings of this study provide a reference for understanding the molecular basis of CaOx crystal-related renal injury, but further studies are needed to verify the relevance of these mechanisms in human kidney stone patients.

Conclusions

In summary, our multi-omics approach reveals that hyperoxaluria-induced kidney injury is characterized by a dynamic interplay between inflammatory responses and metabolic dysfunction, driven by transcriptional and post-transcriptional reprogramming. The suppression of fatty acid metabolism pathways, along with mitochondrial dysfunction, may play a critical role in promoting renal injury and crystal retention. These findings underscore the potential of integrated omics strategies to unravel complex disease mechanisms and inform precision medicine approaches for hyperoxaluria prevention and treatment.

Acknowledgments

None.

Footnote

Reporting Checklist: The authors have completed the ARRIVE reporting checklist. Available at https://tau.amegroups.com/article/view/10.21037/tau-2025-534/rc

Data Sharing Statement: Available at https://tau.amegroups.com/article/view/10.21037/tau-2025-534/dss

Peer Review File: Available at https://tau.amegroups.com/article/view/10.21037/tau-2025-534/prf

Funding: This work was supported by the National Natural Science Foundation of China (No. 82270805).

Conflicts of Interest: All authors have completed the ICMJE uniform disclosure form (available at https://tau.amegroups.com/article/view/10.21037/tau-2025-534/coif). The authors have no conflicts of interest to declare.

Ethical Statement: The authors are accountable for all aspects of the work in ensuring that questions related to the accuracy or integrity of any part of the work are appropriately investigated and resolved. All animal experiments were performed under a project license (No. CHEC2025-021) granted by the Ethics Committee of Changhai Hospital, in compliance with the institutional guidelines for the care and use of animals.

Open Access Statement: This is an Open Access article distributed in accordance with the Creative Commons Attribution-NonCommercial-NoDerivs 4.0 International License (CC BY-NC-ND 4.0), which permits the non-commercial replication and distribution of the article with the strict proviso that no changes or edits are made and the original work is properly cited (including links to both the formal publication through the relevant DOI and the license). See: https://creativecommons.org/licenses/by-nc-nd/4.0/.

References

1. Leumann E, Hoppe B. The primary hyperoxalurias. J Am Soc Nephrol 2001;12:1986-93. [Crossref] [PubMed]
2. Nazzal L, Puri S, Goldfarb DS. Enteric hyperoxaluria: an important cause of end-stage kidney disease. Nephrol Dial Transplant 2016;31:375-82. [Crossref] [PubMed]
3. Robijn S, Hoppe B, Vervaet BA, et al. Hyperoxaluria: a gut-kidney axis? Kidney Int 2011;80:1146-58. [Crossref] [PubMed]
4. Demoulin N, Aydin S, Gillion V, et al. Pathophysiology and Management of Hyperoxaluria and Oxalate Nephropathy: A Review. Am J Kidney Dis 2022;79:717-27. [Crossref] [PubMed]
5. Dutta C, Avitahl-Curtis N, Pursell N, et al. Inhibition of Glycolate Oxidase With Dicer-substrate siRNA Reduces Calcium Oxalate Deposition in a Mouse Model of Primary Hyperoxaluria Type 1. Mol Ther 2016;24:770-8. [Crossref] [PubMed]
6. Garrelfs SF, Frishberg Y, Hulton SA, et al. Lumasiran, an RNAi Therapeutic for Primary Hyperoxaluria Type 1. N Engl J Med 2021;384:1216-26. [Crossref] [PubMed]
7. Knauf F, Asplin JR, Granja I, et al. NALP3-mediated inflammation is a principal cause of progressive renal failure in oxalate nephropathy. Kidney Int 2013;84:895-901. [Crossref] [PubMed]
8. Zhu W, Qiong D, Yanli G, et al. Proteomics and transcriptomics profiling reveals distinct aspects of kidney stone related genes in calculi rats. BMC Genomics 2023;24:127. [Crossref] [PubMed]
9. Okada A, Nomura S, Higashibata Y, et al. Successful formation of calcium oxalate crystal deposition in mouse kidney by intraabdominal glyoxylate injection. Urol Res 2007;35:89-99. [Crossref] [PubMed]
10. Pizzolato P. Histochemical recognition of calcium oxalate. J Histochem Cytochem 1964;12:333-6. [Crossref] [PubMed]
11. Brooks C, Wei Q, Cho SG, et al. Regulation of mitochondrial dynamics in acute kidney injury in cell culture and rodent models. J Clin Invest 2009;119:1275-85. [Crossref] [PubMed]
12. Chen S, Zhou Y, Chen Y, et al. fastp: an ultra-fast all-in-one FASTQ preprocessor. Bioinformatics 2018;34:i884-90. [Crossref] [PubMed]
13. Bouhaddani SE, Houwing-Duistermaat J, Salo P, et al. Evaluation of O2PLS in Omics data integration. BMC Bioinformatics 2016;17:11. [Crossref] [PubMed]
14. Anders HJ, Suarez-Alvarez B, Grigorescu M, et al. The macrophage phenotype and inflammasome component NLRP3 contributes to nephrocalcinosis-related chronic kidney disease independent from IL-1-mediated tissue injury. Kidney Int 2018;93:656-69. [Crossref] [PubMed]
15. Canela VH, Bowen WS, Ferreira RM, et al. A spatially anchored transcriptomic atlas of the human kidney papilla identifies significant immune injury in patients with stone disease. Nat Commun 2023;14:4140. [Crossref] [PubMed]
16. Kang HM, Ahn SH, Choi P, et al. Defective fatty acid oxidation in renal tubular epithelial cells has a key role in kidney fibrosis development. Nat Med 2015;21:37-46. [Crossref] [PubMed]
17. Gewin LS. Sugar or Fat? Renal Tubular Metabolism Reviewed in Health and Disease. Nutrients 2021;13:1580. [Crossref] [PubMed]
18. Hwang S, Chung KW. Targeting fatty acid metabolism for fibrotic disorders. Arch Pharm Res 2021;44:839-56. [Crossref] [PubMed]
19. Zhang D, Zhang YH, Liu B, et al. Role of peroxisomes in the pathogenesis and therapy of renal fibrosis. Metabolism 2025;166:156173. [Crossref] [PubMed]
20. Peroxisomes Vasko R, Injury Kidney. Antioxid Redox Signal 2016;25:217-31. [Crossref] [PubMed]
21. Vasko R, Ratliff BB, Bohr S, et al. Endothelial peroxisomal dysfunction and impaired pexophagy promotes oxidative damage in lipopolysaccharide-induced acute kidney injury. Antioxid Redox Signal 2013;19:211-30. [Crossref] [PubMed]
22. Xu Y, Liang H, Xia K, et al. Umbelliferone attenuates calcium oxalate crystal-induced renal injury and inflammation by attenuating autophagy through the PI3K/AKT pathway. Int Immunopharmacol 2025;150:114250. [Crossref] [PubMed]
23. Beavo JA, Brunton LL. Cyclic nucleotide research -- still expanding after half a century. Nat Rev Mol Cell Biol 2002;3:710-8. [Crossref] [PubMed]
24. Teixeira MM, Gristwood RW, Cooper N, et al. Phosphodiesterase (PDE)4 inhibitors: anti-inflammatory drugs of the future? Trends Pharmacol Sci 1997;18:164-71. [Crossref] [PubMed]
25. Shabb JB. Physiological substrates of cAMP-dependent protein kinase. Chem Rev 2001;101:2381-411. [Crossref] [PubMed]
26. Tavares LP, Negreiros-Lima GL, Lima KM, et al. Blame the signaling: Role of cAMP for the resolution of inflammation. Pharmacol Res 2020;159:105030. [Crossref] [PubMed]
27. Brosnan JT, Brosnan ME. The sulfur-containing amino acids: an overview. J Nutr 2006;136:1636S-1640S. [Crossref] [PubMed]