Multiple myeloma (MM) is a clonal plasma cell malignancy that arises within the bone marrow (BM), leading to disruptions in normal hematopoiesis, immune dysfunction, osteolytic bone lesions, anemia, and renal impairment1. Despite the significant advances in treatment, such as proteasome inhibitors, immunomodulatory agents, monoclonal antibodies, and chimeric antigen receptor T-cell (CAR-T) therapies, MM remains largely incurable with most patients eventually relapsing2.
Autologous hematopoietic stem cell transplantation (auto-HSCT) serves as a cornerstone of MM treatment, significantly prolonging progression-free and overall survival in eligible patients3. High-dose chemotherapy following auto-HSCT eliminates both malignant and normal hematopoietic cells, facilitating timely immune reconstitution from infused autologous stem cells4. Although previous studies have evaluated immune reconstitution following the therapy, most have relied on standard clinical assays, which lack the resolution to capture cellular heterogeneity and dynamic transcriptional changes5–7. Sudha et al.8 recently reported the distinct immune microenvironment features of MM. A major challenge in MM management is the persistence of minimal residual disease (MRD)-subclinical populations of malignant plasma cells that can evade therapy and trigger relapses. MM cells exhibit heightened endoplasmic reticulum (ER) stress and an activated unfolded protein response (UPR) in the pre-transplant setting, driven by excessive immunoglobulin production9. However, these stress pathways persist or resolve high-dose melphalan following auto-HSCT remain unclear, which may provide key insights into MRD biology.
In this study single-cell RNA sequencing (scRNA-seq) was used to generate a comprehensive atlas of the BM microenvironment in MM patients before and after treatment with high-dose melphalan following auto-HSCT. Small residual malignant plasma cells in BM identified by scRNA-seq might be a potential risk factor for MM recurrence and treatment with high-dose melphalan following auto-HSCT contributed to reconstitution of pre-B cells, T cells, NK cells, and hematopoietic stem and progenitor cells (HSPCs) in the BM of MM patients. The findings herein provide new insight and a comprehensive understanding for the single-cell atlas of the BM in MM patients with therapy.
scRNA-seq atlas of the BM in MM patients
Bone marrow mononuclear cells (BMMCs) from healthy donors [HD] (n = 3, GSE120221), untreated MM patients (n = 3, GSE189460), and MM patients treated with high-dose melphalan following auto-HSCT MM patients (n = 3) were analyzed to explore plasma cell remodeling, immune reconstitution, and HSC recovery after high-dose melphalan following auto-HSCT. Patient information is listed in Tables S1–4. The overall workflow is displayed in Figure 1A. The uniform manifold approximation and projection (UMAP) identified 22 clusters and 11 major cell types (Figure 1B–F) after stringent quality control (Figure S1). Cell type proportions differed notably among the three groups, including pre-B cells, B cells, monocyte cells, and T cells. HSPCs were significantly reduced in the BM of MM patients but substantially restored after therapy, while plasma and proliferating cells were markedly enriched in the BM of MM patients but decreased in patients treated with high-dose melphalan following auto-HSCT (Figure 1G, H).
Identification of scRNA-seq atlas of BM in MM patients. (A) Schematic overview of the study design and single-cell RNA-seq workflow. (B) UMAP plot shows unsupervised clustering of 43,997 high-quality cells into 22 clusters. (C) Dot plot displaying the average expression and percentage of cells expressing canonical marker genes across 11 major cell types. (D) UMAP visualization colored by cell type annotations. (E, F) UMAP plots with cells colored by sample identity (E) and clinical group (HD for normal, MM for tumor, and auto-HSCT for treatment) (F). (G) Stacked bar chart comparing the relative proportion of cells across the three groups. (H) Box plots show the percentage of each cell type in HD (healthy donors), MM (multiple myeloma), and auto-HSCT (autologous hematopoietic stem cell transplantation) from (G), respectively. The mean ± SD of at least three experiments is shown. DCs, dendritic cells; HSPCs, hematopoietic stem and progenitor cells.
Elimination and transcriptional remodeling of plasma cells by therapy
A hallmark of MM is clonal expansion of malignant plasma cells in the BM, which disrupts hematopoiesis and driving disease progression. Controlling these cells is key to therapeutic success and long-term remission. Traditional assessments reveal residual disease but overlook transcriptional and clonal dynamics in MM patients receiving high-dose melphalan following auto-HSCT. To address this issue, single-cell transcriptomics was used to characterize plasma cell abundance, pseudotime trajectory, and molecular programs before and after transplantation.
Quantitative comparison revealed a dramatic difference in the relative proportion of cells among the three groups (Figures 1G, 1Ha,b, and S2A,B). Plasma cells made up 1.31% of BMMCs in HD but expanded to 48.01% in MM patients, of which the relative proportion of malignant cells was 22.03%, reflecting clonal infiltration. Plasma cells dropped to 0.21% in MM patients receiving high-dose melphalan following auto-HSCT, of which the relative proportion of malignant cells was 0.09%, indicating an effective elimination of malignant cells. In this study the small sample size limited the generalizability and statistical robustness of the findings, which will be validated in more samples.
Pseudotime trajectory analysis was performed using Monocle 2 to assess plasma cell differentiation and transcriptomic dynamics. Plasma cells from HD and MM patients treated with high-dose melphalan following auto-HSCT predominantly clustered in early pseudotime branches, indicating an immature state, while the MM group exhibited a greater distribution, indicating transcriptional remodeling and altered differentiation (Figure S3A). Gene expression and the trajectory were examined next and six distinct clusters of the top 50 differential gene expression (DEG) were identified (Figure S3B). These results indicated that the plasma cells from auto-HSCT samples share an identical developmental trajectory and functions as the plasma cells from HD.
Pairwise DEGs analysis was performed to assess the transcriptomic changes in plasma cells. MM plasma cells exhibited 3247 upregulated and 235 downregulated genes compared to HD controls (Figure S3C). Auto-HSCT samples demonstrated 1993 upregulated and 748 downregulated genes compared to the MM group (Figure S3D). Functional enrichment of upregulated DEGs in MM revealed a strong association with ER stress pathways, including response to ER stress and protein processing in ER (Figure S4A–C), reflecting a high immunoglobulin load. In contrast, these ER stress-related pathways were significantly downregulated after transplantation (Figure S4D–F), demonstrating the reduced proteotoxic stress and clearance of malignant clones. Gene Ontology (GO) and HALLMARK enrichment analyses revealed that the response to ER stress and unfolded protein response were enriched in MM (Figure S3E, F), while diminished in MM patients receiving high-dose melphalan following auto-HSCT (Figure S3G, H), indicating suppression of the UPR following effective disease control.
MM vs. HD upregulated genes and auto-HSCT vs. MM downregulated genes were intersected to pinpoint the potential molecular targets of residual disease clearance, yielding 640 overlapping genes (Figure S5A) with 515 genes normalized to a healthy level in MM patients treated with high-dose melphalan following auto-HSCT (Figure S5B). Among these genes, 25 were ER stress-related (among 258 annotated ER stress-related genes, such as HSPA5, ATF3, CREBZF, and MANF; Figure S5C) and 14 were linked to the UPR among 113 UPR-related genes, such as ATF3, ATF6, EIF4E, and WFS1 (Figure S3I), showing the potential roles in MM pathogenesis and potential markers for MRD. In addition, 77 genes were downregulated in MM but restored in MM patients treated with high-dose melphalan following auto-HSCT (Figure S5D), including 41 normalized to healthy levels (Figure S5E). The overlap with UPR genes highlighted RPS14, a ribosomal biogenesis gene, selectively re-expressed in MM patients treated with high-dose melphalan following auto-HSCT (Figure S3J). Taken together, the findings indicated that these gene sets (HSPA5, ATF3, ATF6, WFS1, and RPS14) reflect stress adaptation in MM plasma cells and the reversal after treatment, offering potential biomarkers for MRD monitoring and therapy. The B cell receptor (BCR) spectrum in patients after auto-HSCT was analyzed. The histogram plot showed the frequency of use of the top 20 groups in IGHV (Figure S3K) and IGHJ (Figure S3L). Sudha et al. recently reported the distinct immune microenvironment features of MM using scRNA-seq8. However, we focused on the identification of the immune reconstitution and hematopoiesis restoration in the system.
In summary, the significance of using scRNA-seq was investigated in a visual abstract. We conclude that the small residual malignant plasma cells were identified by scRNA-seq in the BM of MM patients, which might be a potential risk factor for MM recurrence. The reconstitution of pre-B cells, T cells, NK cells, and HSPCs in BM was demonstrated in the BM of MM patients treated with high-dose melphalan following auto-HSCT. The findings herein provide new insights into the single-cell atlas in the BM of treated MM patients.
Supporting Information
Conflicts of interest statement
No potential conflicts of interest are disclosed.
Author contributions
Conceived and designed the analysis: Guowen Wang, Guang Yang, Xiaodong Zhang, Yafei Wang.
Collected the data: Qian Li, Yaling Wang, Hongfeng Yuan, Qiaomei Cai, Lina Zhao, Li Lin.
Contributed data or analysis tools: Shuang Gao, Lin Chen, Zhiying Zhang, Jing Ma, Su Liu, Zeng Cao, Haifeng Zhao.
Performed the analysis: Qian Li, Yaling Wang, Hongfeng Yuan, Qiaomei Cai, Yufei Wang, Shihui Li, Pan Lv, Huihui Zhang.
Wrote the paper: Qian Li, Yaling Wang, Hongfeng Yuan, Qiaomei Cai.
Data availability statement
The available scRNA-seq data used in this study are available in the Gene Expression Omnibus under accession code GSE304661.
- Received November 10, 2025.
- Accepted January 8, 2026.
- Copyright: © 2026, The Authors
This work is licensed under the Creative Commons Attribution-NonCommercial 4.0 International License.









