Soft tissue sarcomas (STS) are a rare and biologically heterogeneous group of malignancies arising from mesenchymal tissue, accounting for less than 1% of all malignant tumors. Despite low overall incidence, STS spans at least 100 distinct histological subtypes, each with its own biological behavior, prognostic profile, and response to treatment. This diversity creates substantial clinical complexity: patients can range from pediatric cases of rhabdomyosarcoma to elderly patients with undifferentiated pleomorphic sarcoma (UPS), and optimal management differs meaningfully across subtypes. Recurrence rates remain high even after curative-intent resection, and metastatic STS carries a poor prognosis with limited effective systemic options.
Immunotherapy and the response problem: Immune checkpoint inhibitors have transformed outcomes in several solid tumors, and there is interest in extending these benefits to STS. However, individual responses to PD-1 or CTLA-4 blockade in sarcoma vary considerably and are largely unpredictable with current tools. The SARC028 trial of pembrolizumab in STS and bone sarcoma provided modest overall responses, underscoring that identifying the subset of patients who will genuinely benefit from immunotherapy is a critical unmet need.
The stemness hypothesis: Cancer stem cells (CSCs) are a subpopulation of tumor cells with self-renewal and differentiation capacity. In STS specifically, CSCs have been identified across multiple subtypes including rhabdomyosarcoma, synovial sarcoma, and fibrosarcoma. Dysregulated stemness pathways, including Hedgehog, Hippo, and Notch signaling, have been documented in STS. High stemness expression is associated with tumor metastasis, treatment resistance, and poor prognosis, in part through immune evasion mediated by PD-L1 upregulation and suppression of NK cell-mediated killing.
This 2022 study, published in Frontiers in Immunology by Gu et al. from Wuhan University, set out to comprehensively characterize tumor stemness in STS using multi-omic data from 568 patients, build machine learning-derived stemness subtypes, develop a prognostic scoring index, and determine which patients are most likely to respond to immunotherapy.
The study drew on two publicly available STS cohorts. The primary cohort was TCGA-SARC, which provided 259 STS patients with RNA-sequencing data (FPKM values), DNA methylation data from the 450K array, somatic mutation data (MuTect2 Variant Aggregation and Masking), and copy number alteration (CNA) data, all downloaded from the UCSC Xena browser. The validation cohort was GSE21050 with 309 STS patients, whose gene expression data were generated from microarray and previously normalized using the GCRMA (GC-Robust Multi-Array Analysis) algorithm. FPKM values from TCGA-SARC were converted to TPM for comparability, and batch effects between cohorts were removed using the "sva" package in R 4.0.3. Combined, the overall cohort contained 568 STS samples.
Identifying stemness-related signatures with WGCNA: Rather than relying directly on the existing mRNAsi stemness index (calculated as Spearman correlations between a one-class logistic regression model's weight vector and each sample's gene expression), the authors used Weighted Gene Co-expression Network Analysis (WGCNA) to find genes most strongly associated with stemness in STS specifically. A matrix of 3,622 stem cell-associated genes (sourced from 26 human stem cell gene sets via StemChecker) was analyzed. A soft-thresholding power of beta = 8 (R2 = 0.90) was selected to construct a scale-free co-expression network, yielding six gene modules. The blue module, containing 560 genes, showed the highest positive correlation with mRNAsi (cor = 0.56, p = 6e-22), with strong alignment between module membership and gene significance (cor = 0.55, p = 1.4e-45).
Prognostic SRS selection and clustering: Univariate Cox regression analysis (p less than 0.01) was applied to the 560 blue-module genes, identifying 64 prognostic stemness-related signatures (SRSs). These 64 SRSs were then used as input for consensus clustering (K-means algorithm, ConsensusClusterPlus package, repeated 50 times for stability). Consensus matrix and empirical cumulative distribution function (CDF) plots guided selection of k = 3 as the optimal number of stemness subtypes, as this produced the cleanest separation while offering more interpretive granularity than k = 2.
Multi-omic characterization pipeline: To characterize differences between stemness subtypes, the authors deployed a battery of computational tools: immune and stromal scores via ESTIMATE; immune cell quantification via immunedeconv (using 7 algorithms including TIMER, xCell, MCP-counter, CIBERSORT, EPIC, quanTIseq, and IPS); 29 immune gene sets scored with ssGSEA; pathway scoring via GSVA (using KEGG background gene sets from MSigDB v7.4); and differential gene expression via limma (|logFC| greater than 2, FDR less than 0.05). DNA methylation analysis, somatic mutation profiling with maftools, and CNA comparison using chi-square tests were also performed.
Consensus clustering of the 568 STS patients produced three subtypes: cluster A (n=162, highest stemness), cluster B (n=227, intermediate stemness), and cluster C (n=179, lowest stemness). Kaplan-Meier survival analysis confirmed a clear hierarchy: cluster C had the best prognosis, followed by B, with cluster A showing the worst outcomes. The relationship between stemness and metastatic propensity was striking: 41% of cluster A patients developed tumor metastasis, approximately twice the rate seen in cluster C (roughly 20%). These findings were validated separately in both the TCGA-SARC and GSE21050 cohorts.
Immune microenvironment differences: ESTIMATE scores revealed that cluster C, the low-stemness subtype, had significantly higher immune and stromal scores compared to clusters A and B (p less than 0.05). Across all seven immune deconvolution algorithms, cluster C showed enrichment in immune effector populations including B cells, dendritic cells (DC), Th1 cells, CD8+ T cells, activated NK cells, gamma delta T cells, and M1 macrophages. All of these cell types have recognized roles in direct or indirect tumor killing. In contrast, cluster A was characterized by higher proportions of immunosuppressive populations including Th2 cells, M0 macrophages, regulatory T cells (Treg), and cancer stem-like cells. ssGSEA scores from 29 immune gene sets confirmed that both innate and adaptive immunity were broadly activated in cluster C and suppressed in cluster A.
Metabolic pathway differences: GSVA-based pathway analysis showed that low-stemness cluster C was enriched in tumor immune pathways and multiple metabolic programs, including drug metabolism cytochrome P450, histidine metabolism, tryptophan metabolism, and fatty acid metabolism. High-stemness cluster A, by contrast, showed activation of cell cycle, cell division, DNA replication, and mismatch repair pathways, consistent with a more proliferative and genomically unstable phenotype.
Differentially expressed genes: Limma analysis identified 150 DEGs across the three subtypes (|logFC| greater than 2, FDR less than 0.05). GO and KEGG enrichment analyses confirmed that these DEGs were concentrated in cell cycle regulation, cell division, and DNA metabolism, reflecting the mechanistic basis for the observed proliferative advantage and genomic instability of high-stemness STS.
Somatic mutation data from 235 TCGA-SARC samples were analyzed to characterize genomic differences among stemness subtypes. Gene mutation rates decreased progressively from high-stemness to low-stemness groups: 78.79% (52/66) of cluster A samples harbored gene mutations, compared to 73.53% (75/102) in cluster B and 61.19% (41/67) in cluster C. Tumor mutation burden (TMB) followed the same gradient (cluster A greater than B greater than C, p less than 0.05). Specific mutations of interest identified across subtypes included ATRX (a chromatin remodeler implicated in homologous recombination) and MUC16, both of which emerged as potential stemness regulatory targets in STS.
Copy number alterations: CNA analysis comparing the highest-stemness (cluster A) and lowest-stemness (cluster C) groups identified 184 genes with significant copy number differences (p less than 0.0001). Functional annotation of these genes revealed involvement in fat cell proliferation, Notch binding, regulation of nucleotide metabolic processes, and TGF-beta signaling, all pathways with known roles in cancer stemness and tumor biology. High-stemness cluster A patients had the highest overall copy number burden, and gistic scores were elevated across 22 chromosomes, with particularly notable differences in chromosomes 1-4, 7, 9, 13, 15, 17, and 19.
DNA methylation: Differential methylation analysis between the highest and lowest stemness groups identified 2,054 differentially methylated CpG sites (|logFC| greater than 0.25, adjusted p less than 0.05). Counterintuitively, more sites with higher methylation levels were found in the lowest-stemness (cluster C) group. Pathway enrichment of differentially methylated genes pointed to skin and epidermis development, arginine biosynthesis, and steroid hormone biosynthesis. Three stemness-related methylation-driven genes were identified by the MethyMix package: CPXM2 (cor = -0.603, p = 8.039e-09), CYP1B1 (cor = -0.565, p = 1.025e-07), and DES (cor = -0.801, p = 3.909e-18), all showing strong negative correlations between gene expression and methylation level. CYP1B1 has been previously linked to cancer stem cell behavior in head-and-neck carcinoma; DES and CPXM2 have been associated with poor prognosis in colorectal and gastric cancers, respectively.
To translate the stemness subtype concept into a quantitative clinical tool, the authors constructed the Stemness Prognostic Index (SPi). The development process began by identifying the intersection of the 64 prognostic SRSs and the 150 DEGs from the subtype analysis, yielding 57 candidate genes. These were subjected to LASSO (Least Absolute Shrinkage and Selection Operator) regression with 1,000-fold cross-validation using the glmnet package, which selected 16 optimal genes with non-zero coefficients to form the final model.
The SPi formula: SPi = -0.011 * S100A2 + 0.396 * RAD54L + 0.160 * AURKB + 0.104 * TRIP13 + 0.089 * MAD2L1 + 0.030 * KIF15 - 0.048 * CKS2 - 0.278 * CEP55 + 0.027 * PBK + 0.058 * TK1 - 0.050 * PRC1 - 0.277 * OIP5 - 0.079 * UBE2T + 0.134 * SLC2A1 - 0.145 * SERPING1 - 0.061 * IGF1. Several of these genes (RAD54L, AURKB, TRIP13, MAD2L1, KIF15) have known roles in cell cycle control, DNA damage response, and chromosomal segregation, consistent with the proliferative biology of high-stemness sarcoma.
Prognostic performance: SPi was trained in TCGA-SARC and validated in both GSE21050 and the combined overall cohort. Kaplan-Meier analysis confirmed that low-SPi patients had significantly better prognosis than high-SPi patients in all three cohorts. Importantly, SPi dramatically outperformed the raw mRNAsi index: mRNAsi achieved AUC of only 0.560, 0.537, and 0.483 at 1, 3, and 5 years respectively in TCGA-SARC, making it essentially non-informative for prognostic stratification. SPi was confirmed as an independent prognostic factor for STS in both univariate and multivariate Cox regression analyses, adjusting for clinical covariates.
Clinical correlates of SPi: High SPi was significantly associated with patient age greater than 65 years, UPS histological subtype, and presence of metastases (p = 0.014). High-SPi patients were less likely to respond to treatment, while low-SPi patients showed greater treatment responsiveness. SPi correlated positively with mRNAsi (cor = 0.24, p = 9.7e-05), confirming biological consistency with the stemness concept while providing substantially improved clinical discriminatory power.
A central goal of the study was to determine whether tumor stemness could guide immunotherapy decisions in STS. Three independent analytical frameworks were used to assess immunotherapy responsiveness across stemness subtypes and SPi groups. First, T-cell Inflammatory Score (TIS), calculated using GSVA, identifies patients more likely to benefit from immune checkpoint inhibitors based on inflammatory gene expression. Second, Tumor Immune Dysfunction and Exclusion (TIDE) scores, which have been shown to outperform PD-1 expression and TMB as predictors of checkpoint inhibitor response, were used, where lower scores indicate greater predicted responsiveness. Third, SubMap analysis was applied using a 47-patient melanoma cohort treated with CTLA-4 blockade and PD-1 blockade, which provides a cross-cancer estimate of immunotherapy sensitivity.
Stemness subtypes and immunotherapy: Cluster C (low stemness) showed higher TIS scores and lower TIDE scores compared to cluster A, indicating predicted greater benefit from immune checkpoint blockade. SubMap analysis confirmed that cluster C patients were more sensitive to PD-1 blockade. The immune advantage of low-stemness STS was consistent across all three analytical methods and replicated in both the TCGA-SARC and GSE21050 cohorts, though the TIDE score difference between clusters A and C did not reach statistical significance (p = 0.095).
SPi as an immunotherapy predictor: In parallel analyses, low-SPi patients showed higher TIS scores and lower TIDE scores than high-SPi patients, and SubMap analysis confirmed that low-SPi STS patients had significantly greater predicted sensitivity to PD-1 blockade (Benjamini-Hochberg corrected p less than 0.05). High-SPi patients showed enrichment in DNA repair, E2F targets, glycolysis, mTORC1 signaling, and Wnt/beta-catenin pathways, all of which are associated with immune evasion and treatment resistance. Low-SPi patients showed activation of immune pathways, KRAS signaling, and TNFA signaling via NF-kB, consistent with an immune-activated microenvironment.
Adding MSI to SPi: The authors observed a significant correlation between SPi and microsatellite instability (MSI) (cor = 0.25, p = 7e-05). Reasoning that MSI independently informs immunotherapy response, they combined MSI status with SPi to create a four-group classification: MSIhigh-SPihigh, MSIhigh-SPilow, MSIlow-SPihigh, and MSIlow-SPilow. Kaplan-Meier analysis showed markedly different prognoses across these four groups (p less than 0.001). Critically, MSIlow-SPilow patients were most predicted to respond to immunotherapy, while MSIlow-SPilow patients showed the highest TIS scores and greatest sensitivity to both PD-1 and CTLA-4 blockade (adjusted p = 0.04 for MSI contribution to immunotherapy prediction).
Sample size and retrospective design: The authors explicitly acknowledge that the total cohort of 568 patients represents a small sample, particularly given the ambition of a multi-omic analysis across 100+ STS histological subtypes. The low incidence of STS makes assembling large prospective cohorts inherently difficult, and the authors argue that the combination of two well-established public datasets and multiple computational validation methods partially compensates for this limitation. Nevertheless, subtype-specific analyses (e.g., SPi performance within leiomyosarcoma vs. UPS specifically) are constrained by small within-subtype sample sizes.
Immunotherapy validation in a proxy cohort: The SubMap immunotherapy prediction analysis was performed using a 47-patient melanoma cohort rather than a STS-specific immunotherapy dataset. Melanoma is biologically distinct from STS, and while SubMap is designed to allow cross-cancer comparison, extrapolating immunotherapy response predictions from melanoma patients to sarcoma patients introduces uncertainty. No prospective STS immunotherapy cohort with molecular profiling was available for direct validation, which is the most significant gap in the study.
mRNAsi as a starting point: The mRNAsi index used to anchor the WGCNA analysis was originally derived from pluripotent stem cell gene expression in a broad pan-cancer context, not specifically optimized for mesenchymal tumors. The authors note its poor prognostic performance in STS (AUC 0.48-0.56 at 1-5 years), which they use to motivate SPi, but this also means the WGCNA analysis identifying SRSs was guided by an imperfect phenotypic trait. The derived SPi addresses this empirically, but the biological alignment between mRNAsi-based WGCNA modules and true sarcoma stemness biology remains partially assumption-dependent.
Lack of functional experimental validation: The study is entirely computational. While the statistical methods are rigorous and multi-layered, the identified stemness subtypes, methylation-driven genes (CPXM2, CYP1B1, DES), and SPi genes have not been functionally validated in cell line or animal models. In particular, the causal roles of these methylation-driven genes in regulating STS stemness and immune evasion remain to be demonstrated through wet-lab experiments.
Prospective validation in STS immunotherapy trials: The most immediate next step is to validate SPi and the MSI-SPi combined classifier in prospective STS cohorts receiving immune checkpoint inhibitors. Future trials of PD-1 blockade, CTLA-4 blockade, or combination immunotherapy in STS should include baseline molecular profiling to enable SPi scoring and stratified analyses. If low-SPi and MSIlow-SPilow patients show consistently superior outcomes on immunotherapy, SPi could inform patient selection in trial design and, eventually, routine clinical practice.
Functional characterization of stemness-driven immune evasion: The observation that high-stemness STS is simultaneously associated with more somatic mutations, higher CNA burden, and paradoxically lower immune infiltration warrants mechanistic investigation. The authors speculate that high mutation burden in stemness-high tumors may generate "invalid antigens" or impair antigen presentation rather than promoting immunogenicity, but this requires direct testing. In vitro models using CSC-enriched STS lines and co-culture systems with immune effector cells could delineate how stemness drives immune exclusion despite genomic instability.
Targeting stemness to reverse immune resistance: The authors highlight that the reversibility and plasticity of cancer stem cells create potential therapeutic vulnerabilities. Agents targeting the Hedgehog, Notch, or Wnt/beta-catenin pathways, all activated in high-stemness STS, could in principle reduce stemness and thereby restore immune sensitivity. The methylation-driven genes identified in this study (CYP1B1, CPXM2, DES) represent additional candidate targets for epigenetic intervention using demethylating agents or targeted inhibitors. Combining stemness-modulating agents with immune checkpoint blockade is a rational strategy worth exploring in STS preclinical models.
Applicability to other rare sarcoma subtypes: The authors note that their framework, applying multi-omic machine learning to define stemness subtypes and build prognostic indices, provides a replicable template for other rare tumor types beyond the STS subtypes covered here. Extending this approach to bone sarcomas (osteosarcoma, Ewing sarcoma) or to specific STS histotypes with larger available datasets (leiomyosarcoma, gastrointestinal stromal tumor) would enable more subtype-specific tools than what is currently achievable with the 568-patient combined cohort.