Gene expression datasets for GDM and healthy pregnancies were obtained from the Gene Expression Omnibus (GEO) database (https://www.ncbi.nlm.nih.gov/geo). The bulk RNA-seq dataset (GSE203346) includes data from 21 GDM patients and 22 healthy pregnant women (Supplementary Table 1) (20). Moreover, scRNA-seq data (GSE173193) from placental tissue samples of two GDM patients and two healthy pregnant women were utilized (Supplementary Table 2) (21). A comprehensive list of 1, 136 Mito-RGs was extracted from the MitoCarta3.0 database (22) (https://www.broadinstitute.org/mitocarta).

Differentially expressed Mito-RGs between the GDM and control groups were identified using the R package “limma” (version 3.65.1) on bulk RNA-seq data, with the criteria set at |log2FoldChange (FC)| > 0.585 and adjusted P-value < 0.05. Volcano plots and heatmaps were used to visualize the results by the R packages “ggplot2” (version 3.5.2) and “pheatmap” (version 1.0.12). The biological functions of Mito-RGs were further clarified through Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) enrichment analyses, utilizing the R package “clusterProfiler” (version 4.14.4). GO analysis offered functional annotation for enriched differential genes across biological processes, cellular components, and molecular functions. KEGG pathway analysis utilized a hypergeometric test to identify the genetic functions based on gene enrichment in specific pathways. Gene set enrichment analysis (GSEA) was performed to detect pathways associated with mitochondria. Genes were ranked by log2FC, with pathways deemed significant if their normalized enrichment score (NES) > 1 or < -1, and their P-value was less than 0.05.

Three machine learning algorithms—least absolute shrinkage and selection operator (LASSO) regression, random forest (RF), and extreme gradient boosting (XGBoost)—were employed to identify Hub Mito-RGs. The R package “glmnet” (version 4.1.8) was utilized to conduct LASSO regression, identifying genes with non-zero coefficients by optimizing lambda parameters via 10-fold cross-validation. The RF algorithm was applied using the R package “RandomForest” (version 4.7-1.2) to estimate and rank the importance of Mito-RGs, with the top 10 genes exhibiting the highest MeanDecreaseAccuracy scores selected for further analysis. For XGBoost, the R package “xgboost” (version 1.7.9.1) was utilized to iteratively optimize feature weights, with hyperparameter tuning implemented through grid search combined with 10-fold cross-validation to ensure robust performance. The top 10 genes ranked by Gain importance scores were extracted as candidate features. The genes identified by all three algorithms were designated as hub Mito-RGs, ensuring reliable identification of key Mito-RGs associated with GDM.

The Mito-RGs scores were calculated using multivariate logistic regression, with gene expression levels weighted by their regression coefficients and standardized via z-score normalization. Patients were divided into high-score (z-score > 0) and low-score (z-score ≤ 0) categories. Spearman correlation analysis was conducted to evaluate relationships between Mito-RGs scores and clinical variables, including maternal age, body mass index (BMI), oral glucose tolerance test (OGTT) results, and gestational age. A nomogram was developed using the R package “rms” (version 6.9.0), incorporating Mito-RGs scores alongside clinical parameters. The discriminatory ability of the nomogram was evaluated by calculating the area under the ROC curve (AUC). Calibration curves and decision curve analysis (DCA) were employed to evaluate its calibration accuracy and clinical utility in predicting GDM risk.

The ImmuCellAI algorithm was used to assess the infiltration levels of different immune cell types in placental tissues from both GDM and healthy controls (23). This tool analyzed bulk RNA-seq transcriptional profiles to determine the proportions of immune cell subsets, including T cells, B cells, macrophages, dendritic cells (DCs), and natural killer (NK) cells, in the placental microenvironment. Immune cell infiltration differences between the GDM and control groups were quantified using log2FC for immune cell proportions. A positive log2FC indicated increased infiltration in the GDM group, whereas a negative log2FC reflected reduced infiltration compared to controls. The R package “pheatmap” was utilized to create a heatmap, illustrating the immune cell infiltration landscape and emphasizing variations in cell type abundance.

For scRNA-seq data analysis, initial quality control (QC) was conducted to filter out cells with fewer than 100 expressed genes or those with over 10% mitochondrial gene expression, ensuring high-quality data. Expression matrices were normalized using log-normalization, and highly variable genes were selected for downstream analyses. Dimensionality reduction was performed using principal component analysis (PCA), and graph-based clustering was subsequently conducted with the R package “Seurat” (version 5.1.0). Clusters were visualized in low-dimensional space using the uniform manifold approximation and projection (UMAP). Cell type annotation was performed with SingleR (version 2.8.0) by matching cluster expression profiles against well-established reference datasets. Expression levels of hub genes were visualized across the identified clusters using feature plots and violin plots, showing their specific abundance patterns. The Mito-RGs score at the single-cell level was calculated using the “AddModuleScore” function in Seurat, allowing aggregated expression values to be mapped onto the UMAP projection for spatial distribution analysis. The Mito-RGs scores were mapped onto the UMAP projection through feature plots to reveal their spatial distribution across cell subtypes. Biological functional differences across cell types were examined using gene set variation analysis (GSVA) enrichment analysis and single-cell GSEA (24) with predefined hallmark gene sets. Additionally, the “CellChat” package (version 2.2.0) was applied to infer intercellular communication networks based on ligand-receptor interactions, with a focus on mitochondrial-related pathways. Ligand-receptor pairs mediating key biological functions were identified, and interaction intensity and directionality were visualized to uncover the underlying cellular communication patterns.

Placental tissues were collected from pregnant women with GDM and corresponding healthy controls. All samples were obtained from a teaching hospital affiliated with the Shandong First Medical University. GDM was diagnosed based on the International Association of Diabetes and Pregnancy Study Groups (IADPSG) criteria, which include fasting blood glucose (FBG) levels of at least 5.1 mmol/L, 1-hour postprandial glucose levels of at least 10.0 mmol/L, or 2-hour postprandial glucose levels of at least 8.5 mmol/L during OGTT. Placental tissues from the maternal-fetal interface were collected following vaginal delivery. The study adhered to the Declaration of Helsinki and received approval from the Ethics Committee of the First Affiliated Hospital of Shandong First Medical University (Approval Number: LCYJLL2024-S127).

The in vivo phase of this study, encompassing all animal experimentation and tissue collection, was carried out over a five-month period from March 15, 2024, to July 20, 2024. Female C57BL/6-J mice aged six weeks (16–18 g) were purchased from Jinan Pengyue Experimental Animal Breeding Co., Ltd. (Shandong, China). Mice were provided unrestricted access to food and water and maintained under controlled conditions: a temperature of 22 ± 1°C, relative humidity of 55 ± 1%, and a 12/12 h light/dark cycle. After one week of acclimation, they were randomly divided into a control group (n=8) and a high-fat diet (HFD) group (n=8). A GDM model was established by feeding the HFD group a diet containing 60% of calories from fat (Research Diets) to induce insulin resistance. The HFD regimen was initiated six weeks prior to mating and continued throughout the gestation period. Control mice received a low-fat diet (LFD) with 10% of calories derived from fat (Research Diets). GDM modeling was confirmed by establishing a significantly higher glucose AUC in the HFD group compared to the control group. On gestational day 18 (G18), a glucose tolerance test was performed. After a 12-hour fast, FBG levels were measured using a glucometer. A 20% glucose solution was then administered via oral gavage, and tail vein blood glucose concentrations were subsequently measured at 30, 60, 90, and 120 minutes post-gavage. Maternal body weight was recorded from pre-pregnancy to the end of gestation. The successful induction of GDM was confirmed by sustained hyperglycemia, insulin resistance, and distinct weight gain patterns compared to controls. After establishing the GDM model, we ensured that three mice from each group were selected for final analysis. After a 12-hour fast, all mice were euthanized between June 10 and June 25, 2024, and samples including blood and placental tissue were collected. Euthanasia was carried out under deep isoflurane anesthesia followed by cervical dislocation, consistent with the American Veterinary Medical Association (AVMA) Guidelines. All animal procedures were performed in accordance with the ARRIVE (Animal Research: Reporting In Vivo Experiments) guidelines and approved by the Institutional Animal Care and Use Committee of the First Affiliated Hospital of Shandong First Medical University (Approval Number: LCYJLL2024-S128).

IHC staining was conducted to evaluate the protein expression of dehydrogenase/reductase member 2 (DHRS2), syntaxin 17 (STX17), and translocase of inner mitochondrial membrane 44 (TIMM44). Human and mouse placental tissues were fixed in 4% paraformaldehyde, dehydrated, cleared with xylene, and embedded in paraffin for histological analysis. Sections underwent deparaffinization, rehydration, and antigen retrieval using citrate buffer (pH 6.0). Endogenous peroxidase was blocked with 3% hydrogen peroxide, and nonspecific binding was reduced with 5% bovine serum albumin. Overnight incubation of tissue sections at 4 °C was performed with primary antibodies: DHRS2 (Proteintech, 15735-1-AP), STX17 (Affinity Biosciences, DF12483), and TIMM44 (Affinity Biosciences, DF12332). Sections were washed, incubated with HRP-conjugated secondary antibodies, developed using diaminobenzidine, and counterstained with hematoxylin. Protein expression was semi-quantitatively evaluated using an IHC scoring method incorporating staining intensity and the extent of positive staining. Intensity (I) was graded from 0 to 3 (0, none; 1, weak; 2, moderate; 3, strong). The proportion of positively stained cells (A) was categorized on a 1–4 scale (1, <10%; 2, 10–50%; 3, 50–90%; 4, >90%). The overall IHC score was computed as I×A.

R software (version 4.3.1) was used for statistical analyses. Continuous variables are expressed as the mean ± standard deviation (SD) and categorical data are presented as n (%). For comparisons between two groups, the Student’s t-test was applied for continuous variables and the Chi-square test for categorical variables. Spearman’s rank correlation was used for correlation analyses. Logistic regression models were used for univariate and multivariate analyses, with results expressed as odds ratios (OR) and 95% confidence intervals (95% CI). P-values for multiple comparisons were adjusted using the Benjamini–Hochberg method to control the false discovery rate (FDR), with an adjusted P-value < 0.05 considered statistically significant.

The flow diagram of this study is illustrated in Figure 1. This study examined the gene expression dataset GSE203346, consisting of 21 placental tissue samples from women with GDM and 22 from healthy controls. By intersecting the dataset with the MitoCarta3.0 database, 29 differentially expressed Mito-RGs were identified, all of which exhibited significant upregulation in GDM samples compared to controls. The distribution of these Mito-RGs is visualized in the volcano plot (Figure 2A), and their expression patterns were further illustrated using a hierarchical clustering heatmap (Figure 2B). Functional enrichment analysis of the identified Mito-RGs revealed key molecular processes and pathways central to the pathophysiology of GDM. GO analysis revealed significant enrichment in processes like oxidoreductase activity, fatty acid binding, and small molecule metabolism (Figure 2C). KEGG pathway analysis identified the participation of metabolic pathways, including the citrate cycle (TCA cycle), fatty acid metabolism, and glutathione metabolism (Figure 2D). GSEA highlighted the downregulation of mitochondrial-related pathways in GDM, such as mitochondrial membrane organization, mitochondrial outer membrane permeabilization linked to programmed cell death, and mitochondrial respiratory chain complex assembly (Figure 2E).

In the LASSO regression analysis, the coefficient shrinkage process was visualized (Figure 3A). Through 10-fold cross-validation, the optimal λ value was identified, leading to the selection of 22 genes from the 29 differentially expressed Mito-RGs as potential GDM biomarkers with minimal error (Figures 3B, C). The RF analysis ranked the variables by importance and selected the top 10 genes as potential candidate biomarkers for GDM (Figure 3D). The stability and classification performance of the RF model were verified by plotting the out-of-bag (OOB) error and classification error rates with an increasing number of decision trees (Figure 3E). In addition, the XGBoost analysis ranked the importance of the genes and identified the top 10 important contributing genes to GDM (Figure 3F). A Venn diagram integrating the results from all three machine learning algorithms identified DHRS2, STX17, and TIMM44 as hub Mito-RGs in GDM (Figure 3G). The expression levels of the three hub genes were markedly elevated in GDM placental tissues relative to controls (Figure 3H). ROC curve analysis demonstrated the diagnostic performance of each hub gene for GDM, yielding AUC values of 0.738 for DHRS2, 0.649 for STX17, and 0.688 for TIMM44 (Figure 3I). A composite Mito-RGs score, derived from the three hub genes, was significantly elevated in the GDM group relative to the controls (Figure 3J). Furthermore, the AUC for the composite Mito-RGs score was 0.725 (Figure 3K), indicating its predictive capability. Correlation analysis via a heatmap demonstrated positive association between the Mito-RGs score and BMI (Figure 3L).

Figure 4A illustrates the results of logistic regression analyses assessing potential predictors of GDM. In both univariate and multivariate analyses, there was a statistically significant association between Mito-RGs scores and GDM risk, suggesting that Mito-RGs scores are a strong independent predictor for GDM. To facilitate clinical application, a nomogram was constructed by integrating key variables, including gestational week, maternal age, BMI category, and the Mito-RGs score (Figure 4B). The nomogram allows for individualized GDM risk assessment by assigning points to each predictor variable, summing the total score, and deriving the corresponding risk probability. The ROC analysis demonstrated that the nomogram had a favorable diagnostic performance, achieving an AUC of 0.761 (Figure 4C). The calibration curve demonstrated strong concordance between predicted and observed GDM risk probabilities, validating the reliability of the nomogram (Figure 4D). DCA demonstrated that the nomogram provides substantial net clinical benefits across a wide range of risk thresholds compared to other models, particularly in the ranges of 0.08-0.22 and 0.55-0.69 (Figure 4E).

Immune cell infiltration differences between GDM and control placental tissues were analyzed using the ImmuCellAI algorithm. The heatmap and bar chart highlighted significant alterations in immune cell profiles in GDM samples, characterized by increased infiltration of DCs, macrophages, NK cells, and CD8⁺ T cells (Figure 4F). These findings suggest an abnormal activation of immune function in GDM patients, characterized by changes in placental immune cell composition.

scRNA-seq analysis was conducted to confirm the expression and localization of hub Mito-RGs in GDM. After QC (Supplementary Figure 1), the total number of cells analyzed included 4959 and 4366 cells from control samples, and 6114 and 4366 cells from GDM samples. The GSE173193 dataset was dimensioned into 15 distinct clusters (Figure 5A). Notably, a clear separation between cell populations from GDM and control samples was observed (Figures 5B, C), indicating substantial differences in cellular composition. Cluster-specific marker gene identification (Figure 5D) and SingleR-based annotation (Figure 5E) delineated major cell populations, including epithelial cells, endothelial cells, monocytes, macrophages, neutrophils, myelocytes, NK cells, B cells, common myeloid progenitors (CMPs), and tissue stem cells. Quantitative analysis suggested altered cell-type proportions in the GDM group compared with controls, with a trend toward higher proportions of epithelial cells, macrophages, NK cells, and B cells and a lower proportion of tissue stem cells (Figures 5F, G), although substantial inter-individual variability was observed. Feature plots illustrated DHRS2, STX17, and TIMM44 were predominantly expressed in epithelial cells, macrophages, tissue stem cells, and NK cells (Figure 5H). Moreover, the Mito-RGs score was mapped across cell subtypes, showing pronounced enrichment in epithelial cells, tissue stem cells, NK cells, and macrophages (Figure 5I).

Single-cell GSVA enrichment analysis indicated notable activation of lipid metabolism, glycolysis, oxidative phosphorylation, and proliferation pathways in epithelial cells. Similarly, macrophages exhibited activation in apoptosis, inflammatory responses, and metabolic pathways. Meanwhile, biological processes related to proliferation were significantly activated in NK cells (Figure 6A). The results of single-cell GSEA further confirmed that glycolysis, oxidative phosphorylation, hypoxia, apoptosis, and reactive oxygen species-related pathways were predominantly activated in epithelial cells (Supplementary Figures 2A, B). Meanwhile, multiple inflammatory activation-related pathways, including antigen processing and presentation, IL6-JAK-STAT3 signaling pathway, inflammatory response, and TNF-α signaling via NF-κB were predominantly enriched in macrophages (Supplementary Figure 2C). The cell-cell communication network analysis highlighted strong interactions among epithelial cells, macrophages, and monocytes in GDM (Figure 6B). Key ligand-receptor pathways were identified, including TGFB1-TGFBR1/TGFBR2 signaling between B cells and other cell types (Figure 6C), FN1-ITGA5/ITGB1 signaling between epithelial and other cells (Figure 6D), and LAMA5-CD44 signaling between macrophages and other cell types (Figure 6E).

A total of 12 pregnant women (mean age: 29.73 ± 1.07) were finally included in our study and divided into the GDM group (n = 6) and control group (n = 6). The clinical characteristics of patients are summarized in Supplementary Table 3. To validate the findings, human placental tissues were analyzed via IHC analysis. The study found a notable upregulation in DHRS2 and STX17 expression in placental samples from women with GDM compared to controls, whereas TIMM44 levels showed no significant difference between the groups (Supplementary Figures 3A–C).

A GDM mouse model was developed to further validate these findings. Body weight monitoring showed consistently higher weights in GDM mice compared to controls, with statistical significance reached by week 6 (P-value = 0.0213; Figure 6F). OGTT results revealed significantly elevated blood glucose levels in GDM mice at 30, 60, 90, and 120 minutes post-glucose administration compared to controls, while FBG levels showed no significant difference (Figure 6G). IHC analysis of murine placental tissues further validated the upregulation of hub Mito-RGs. In GDM placentas, DHRS2 and TIMM44 expression levels were significantly elevated compared to controls, while STX17 exhibited an upward trend that did not achieve statistical significance (Figures 6H–J).