Single-cell sequencing analysis related to sphingolipid metabolism guides immunotherapy and prognosis of skin cutaneous melanoma

Background We explore sphingolipid-related genes (SRGs) in skin melanoma (SKCM) to develop a prognostic indicator for patient outcomes. Dysregulated lipid metabolism is linked to aggressive behavior in various cancers, including SKCM. However, the exact role and mechanism of sphingolipid metabolism in melanoma remain partially understood. Methods We integrated scRNA-seq data from melanoma patients sourced from the GEO database. Through the utilization of the Seurat R package, we successfully identified distinct gene clusters associated with patient survival in the scRNA-seq data. Key prognostic genes were identified through single-factor Cox analysis and used to develop a prognostic model using LASSO and stepwise regression algorithms. Additionally, we evaluated the predictive potential of these genes within the immune microenvironment and their relevance to immunotherapy. Finally, we validated the functional significance of the high-risk gene IRX3 through in vitro experiments. Results Analysis of scRNA-seq data identified distinct expression patterns of 4 specific genes (SRGs) in diverse cell subpopulations. Re-clustering cells based on increased SRG expression revealed 7 subgroups with significant prognostic implications. Using marker genes, lasso, and Cox regression, we selected 11 genes to construct a risk signature. This signature demonstrated a strong correlation with immune cell infiltration and stromal scores, highlighting its relevance in the tumor microenvironment. Functional studies involving IRX3 knockdown in A375 and WM-115 cells showed significant reductions in cell viability, proliferation, and invasiveness. Conclusion SRG-based risk signature holds promise for precise melanoma prognosis. An in-depth exploration of SRG characteristics offers insights into immunotherapy response. Therapeutic targeting of the IRX3 gene may benefit melanoma patients.


Introduction
Melanoma is a highly aggressive skin tumor with substantial implications for individuals affected by skin cancer, leading to significant mortality rates (1).The median survival period for this condition spans from 6 to 9 months, with a corresponding 5-year survival rate of less than 5% (2).Here is a significant global surge in the prevalence of melanoma.Despite its relatively lower occurrence in Western countries, the Chinese Society of Clinical Oncology has reported a consistent annual rise in melanoma incidence in China, with rates ranging from 3% to 5% (3).In recent decades, notable progress has been made in the treatment of melanoma.The conventional approach of surgery, radiotherapy, and chemotherapy has evolved into a comprehensive multi-modal strategy (4).This modern approach combines surgical interventions with immunotherapy and targeted therapy, resulting in remarkable improvements in patient survival rates (5).However, there remains a subset of melanoma patients who do not experience benefits from immunotherapy or molecular-targeted therapy in the management of their disease.Despite continuous advancements in treatment modalities, the prognosis for melanoma patients remains unfavorable (6).Cancer development is intricately linked to the complex signaling transduction network within cells, where multiple signaling pathways interconnect and collectively regulate the biological phenotypes of tumor cells.Exploring the underlying mechanisms that drive melanoma progression is crucial to enhance therapeutic strategies and improve the prognosis of individuals affected by this condition.
Both genetic instability and the tumor microenvironment play significant roles in the initiation and advancement of tumors (7).Recent studies have uncovered that manipulating the activation of oncogenes or silencing tumor suppressor genes can regulate metabolic homeostasis and impact the promotion or suppression of cancer (8).Furthermore, certain metabolic enzymes possess the capacity to act as oncogenes or tumor suppressor genes, thus playing a role in tumor progression.Cellular energy requirements predominantly rely on glucose and fatty acid metabolism.Although previous investigations have predominantly focused on glucose metabolism in tumor cells, recent research has unveiled the involvement of abnormal lipid metabolism in influencing the aggressive phenotypes observed in various tumor types (9).Sphingolipids, a class of phospholipids that encompass phosphate groups, are abundantly present in cellular membranes and play essential biological functions in organisms (10).Dysregulated sphingolipid metabolism has the potential to influence the composition and functionality of cellular membranes, leading to the facilitation of increased proliferation, invasiveness, and metastasis in cancer cells (11,12).In contrast, sphingolipid molecules play crucial roles as essential regulators in various medically significant biological processes, including cellular signaling and the orchestrated self-destruction process known as apoptosis (13).Activation of acid sphingomyelinase (A-SMase) can be triggered by diverse pro-inflammatory and pro-apoptotic stimuli, thereby playing a role in the induction of apoptosis in tumor cells in response to various anti-tumor therapies.Additionally, A-SMase has been implicated in immune and inflammatory processes.The research conducted by E. Clementi and C. Perrotta highlights the significance of A-SMase as a crucial factor dictating the behavior of melanoma cells (14).
Dysregulated expression of critical enzymes involved in specific sphingolipid synthesis pathways has been linked to the development and progression of various malignancies.For instance, genetic mutations affecting the PI3K catalytic subunit alpha (PIK3CA) gene have been associated with breast cancer, ovarian cancer, cervical cancer, and other tumor types.Another significant factor impacting the susceptibility to breast cancer is the presence of genetic variations in the sphingomyelin transferase 1 (SMT1) gene within the neural sphingolipid synthesis pathway (15,16).Furthermore, sphingolipids possess the capacity to regulate cellular signaling pathways, thereby impacting tumor cell proliferation, advancement, and resistance to therapeutic agents.Perturbations in phosphatidylinositol (PI) metabolism, such as excessive activation, can result in heightened stimulation of the PI3K/AKT/mTOR pathway, a critical regulator of crucial cellular processes, including cellular growth and apoptosis (17).However, the precise role and underlying mechanisms of sphingolipid metabolism in melanoma remain poorly elucidated, urging the need for additional investigation to advance our understanding of these mechanisms.
Single-cell technology is a high-throughput approach extensively employed for the analysis of individual cells in medical research.Its widespread adoption and continuous progress have been noteworthy in recent times.By examining the gene expression profiles of individual cancer cells, this methodology enables the discovery of distinctive epigenetic characteristics inherent to each cancer cell.Single-cell studies have not only revealed the impact of several individual signaling pathways on tumor cell proliferation, metastasis, and drug resistance, but they also hold revolutionary significance in providing comprehensive and profound insights into the complexity of signaling networks within tumor cells and the functional and regulatory mechanisms of key signaling networks involved in tumor cell proliferation, metastasis, and drug resistance (18).These findings will contribute to the development of novel chemotherapy drugs and targeted treatment strategies.
In cancer research, risk profiles are widely used to predict prognostic outcomes.Pei S, Zhang P et al. used genes associated with sphingolipid metabolism to characterize genes strongly associated with survival in patients with breast and lung cancer (11,12).In addition, for the risk profile constructed by SKCM, these prognostic models were shown to be more accurate than traditional methods in predicting clinical outcomes (19).In the field of SKCM research, the molecular regulation of sphingolipid metabolism has not been fully revealed.Therefore, we included sphingolipid metabolization-related genes in the construction of risk profiles to estimate novel strategies for predicting prognosis in patients with SKCM.
In this study, we utilized scRNA-seq and transcriptome data obtained from publicly available databases to identify distinct subsets of melanoma based on SRGs.Subsequently, these SRGs were used to establish risk factors capable of predicting melanoma prognosis.Furthermore, we conducted an in-depth analysis to explore the molecular characteristics derived from SRGs and their clinical relevance.We also investigated the role of signaling pathways in cancer cell proliferation, metastasis, and drug resistance, as well as the effectiveness of immune therapy including immune checkpoint pathways such as PD-1/PD-L1 and CTLA-4, and the activation of immune cells such as NK cells and tumor-associated macrophages in the immune microenvironment.This innovative study provides a groundbreaking perspective on the prognostic stratification of melanoma, facilitating personalized treatment approaches and improving clinical outcomes for patients with melanoma.

Acquisition of original patient data
The scRNA-seq data specific to SKCM were obtained from the Gene Expression Omnibus (GEO) database, with the accession number GSE123139.Subsequently, two cohorts, namely GSE19234 and TCGA data, were selected for subsequent analysis.To ensure data quality, genes expressed in less than three cells or a single cell containing fewer than 250 genes were excluded from the analysis.The Seurat R package's PercentageFeatureSet function was employed to assess the proportion of ribosomal RNA (rRNA) and mitochondria present in the dataset.As a result of this preprocessing, a total of 2725 cells were retained and utilized for further investigation.
The transcriptomic information about SKCM, as well as clinical details, were obtained from the TCGA database.Subsequently, samples lacking outcome status or survival information were excluded from the analysis, resulting in a dataset comprising 300 SKCM samples earmarked for external validation.Additionally, for the training cohort, 44 tumor samples from GSE19234 were selected after eliminating untracked samples sourced from the GEO database.To categorize individual cells into distinct subgroups, the FindNeighbors and FindClusters functions were employed.To reduce the dimensional complexity of the dataset, the RunUMAP function was utilized for UMAP dimension reduction.Within our dataset, we focused our investigation on four genes (PSAP, APOE, ASAH1, DEGS1) associated with sphingolipid metabolism by analyzing their respective gene expression profiles.

Identification of expression SRGs
The Seurat package was employed to re-analyze scRNA-seq data derived from melanoma samples to assess SRGs.Initially, cells expressing fewer than 250 or more than 6000 genes were excluded.The remaining expressed genes were then subjected to log-normalization.To account for batch effects, the FindIntegrationAnchors function was applied.Subsequently, UMAP was utilized for dimensionality reduction with a resolution of 0.1, considering 30 principal components.FindNeighbors and FindClusters functions were utilized to classify cells into distinct subgroups, using a dimensional parameter of 30 and a resolution of 0.1.RuntUMAP was employed for further reduction of UMAP dimensions.Marker genes, namely PSAP, APOE, ASAH1, and DEGS1, were used to annotate the SRGs.Additionally, the SRGs underwent re-clustering using the FindClusters and FindNeighbors methods.The FindAllMarkers tool was employed to compare different clusters and identify marker genes for each cluster within the SRGs data, taking into account input, logFC, and adjusted p-value parameters.

Hub genes identification according to SRGs
We utilized the scale method provided by the "Limma" R package to normalize gene expression profiles in our study.Prognosisassociated key genes were selected based on the criteria of |log2 (FoldChange)|>1 and a false discovery rate (FDR) of<0.05.Cox regression analysis was applied to screen marker genes from seven clusters that were associated with prognosis.To reduce the gene set, we implemented the LASSO technique.To construct an SRG-derived risk profile, we performed multivariate Cox regression analysis using the stepwise regression approach to minimize redundancies.The risk score, computed using a specific formula: Risk score = o n k=1 Coef (k) Â Expr(k), incorporates the regression coefficients (Coef(k)) and the expression levels (Expr(k)) of the prognostic model genes.Zero-mean normalization was employed to categorize patients as either low or high-risk.The predictive capacity of the risk signature was evaluated using timeROC software to analyze the receiver operating characteristic (ROC) curves.

Developing a novel nomogram
We developed a novel nomogram to predict the prognosis of melanoma by considering the risk signature and clinicopathological characteristics.Both univariate and multivariate Cox regression analyses were conducted to analyze the association between various variables and prognosis.Variables with p-values<0.05 were selected and included in the multivariate Cox regression model.To assess the accuracy of the prognostic predictions made by the model, a calibration curve was constructed.

Cluster analysis
Through an iterative process, a partitioning scheme consisting of K clusters was determined by minimizing the loss function associated with the clustering outcomes.The K-means clustering method was applied to group melanoma patients based on 11 modeling genes.

Assessment of immune landscape
To assess the correlation between the risk signature and tumor immune microenvironment (TIME), a combination of algorithms, including CIBERSORT, EPIC, MCPCOUNTER, and TIMER, were employed in the evaluation process.The R package "estimate" was utilized to calculate stromal scores, immune scores, and estimate scores, which represent the combined scores of stromal and immune components.Additionally, the CIBERSORT algorithm was utilized to analyze the distribution of 22 distinct subtypes of immune cells, providing insights into the heterogeneity of the immune response within the tumor microenvironment.Furthermore, a comprehensive study was conducted to explore the relationship between the genes comprising the risk signature and the immunological score, shedding light on the important role of these genes in immune-related functions.

The analysis of the impact of immunotherapy
To evaluate the predictive potential of our risk profile in predicting the response to immune checkpoint blockade therapy, we conducted an assessment of its efficacy using transcriptomic data and corresponding clinical information from patients enrolled in the IMvigor210 dataset.These patients were treated with anti-PD-L1 therapy.In addition, we incorporated transcriptomic data from a separate cohort of melanoma patients in the GSE78220 dataset who had previously been treated with anti-PD-1 checkpoint inhibitors.

Cell culture
The WM-115 and A375 cell lines, obtained from the Cell Resource Center of Shanghai Life Sciences Institute, were cultured in a DMEM medium (Gibco BRL, USA).The cells were maintained at 37°C with 5% CO 2 and supplemented with 10% fetal bovine serum (FBS) sourced from Gibco BRL, USA.

Cell transfection
Two distinct small interfering RNAs (siRNAs) specifically designed to target IRX3 were developed by Ribobio (Guangzhou, China).Transfections were performed using Lipofectamine 3000 (Invitrogen, USA) according to the manufacturer's instructions (20).The siRNA sequences for IRX3 can be found in Supplementary Table 1.

Patients and tissue samples
A cohort consisting of 20 melanoma tissues and paired normal tissues was utilized for qPCR analysis.The tissues included in this study were pathologically verified at the Department of Plastic Surgery, First Affiliated Hospital, Anhui Medical University (AHMU) in China, during the period from 2020 to 2023.Prior approval for conducting this study was obtained from the Ethical Committee of the First Affiliated Hospital of Anhui Medical University.

RT-qPCR analysis
RNA extraction from cell lines was carried out using TRIzol (Thermo, 15596018) following standard protocols.Subsequently, cDNAs were synthesized using the PrimeScriptTM RT kit (Vazyme, R232-01).To quantify gene expression, SYBR qPCR Master Mix (Vazyme, Q111-02) was employed on the Roche LightCycler 480 (Roche, GER), and data analysis was performed using the 2 −DDCt method.The specific primer sequences, provided by Tsingke Biotech (Beijing, China), are available in Supplementary Table 1.For normalization, GAPDH was utilized as the internal reference gene.

The experiment of cell-counting-kit-8 assay
Cells were plated in 96-well plates at a density of 1 × 10 3 cells per well.Following that, the plates were incubated in darkness at 37°D C for 2 hours with CCK-8 labeling reagent (A311-01, Vazyme).The assessment of cell viability was carried out by measuring the absorbance at 450 nm using an enzyme-linked spectrophotometer (A33978, Thermo) at time intervals of 0, 24, 48, 72, and 96 hours.

The experiment of colony formation
A cohort comprising 1000 cells was transfected and cultured in 6-well plates for approximately 14 days.After 2 weeks, the cellular clones were visually examined without magnification.Following that, the cells were washed and fixed using a 4% paraformaldehyde (PFA) solution for 15 minutes.Subsequently, the cells were subjected to crystal violet staining (Solarbio, China) for 20 minutes, and the samples were air-dried at room temperature.Finally, quantification of cells per well was conducted.

Wound healing
The transfected cells were cultivated in 6-well plates and incubated in a cell incubator until reaching a confluency level of 95%.A 200ml pipette tip was employed to create a straight scratch across the cell monolayer.Following the removal of unattached cells and debris using PBS, the cells were transferred to a serum-free culture medium.Subsequently, photographs were captured at the same location both at 0 hours and 48 hours, and the width of the scratch was measured using Image J software.

Transwell
Transwell chambers were employed to perform cell migration and invasion assays.A total of 2×104 cells per well were seeded in the upper compartment using a 200 mL serum-free medium.To assess the migratory and invasive abilities of the cells, the upper region of the chamber was treated with Matrigel solution (BD Biosciences, USA) in some cases, while it remained untreated in others.The lower chamber was filled with 600 mL of complete medium.After incubating for 48 hours, the chambers were retrieved.The cells were fixed with 4% PFA and then stained with 0.1% crystal violet (Solarbio, China).Subsequently, cell counting was performed using a light microscope.The migrated cells were captured in photographs and quantified.

Apoptosis assay
The apoptotic rate was assessed utilizing an Annexin V-APC/PI Apoptosis Detection Kit provided by Multisciences, China, and further analyzed using a flow cytometry system manufactured by BD Biosciences, USA.The proportions of apoptotic cells at early and late stages were evaluated to determine the apoptotic rate.

Statistical analysis
R software version 4.1.3was used for biological analysis, while GraphPad Prism version 8.0 was employed for analyzing experimental data in the field of medicine.The mean values and standard deviations of the outcomes were obtained from three separate studies.Pairwise comparisons between two groups were conducted using Student's t-tests, while comparisons involving more than two groups were analyzed using one-way ANOVAs followed by Tukey's test (*P<0.05,**P<0.01,***P<0.001).

Screen the SRGs
Figure 1 illustrates the flow chart outlining the progression of the study.A total of 2725 cells were obtained after the completion of scRNA-seq data analysis.Following log-normalization and dimensionality reduction, a total of 14 distinct subpopulations were identified in the analysis.Subsequently, based on a literature review, we selected four genes (PSAP, APOE, ASAH1, and DEGS1) that were most closely associated with sphingolipid metabolism and designated them as marker genes for sphingolipid metabolism (Figure 2A).These four key genes were then utilized to identify cells actively involved in sphingolipid gene sets using the AUCell exploration Threshold function.Based on the median AUC scores, the cells were categorized into high-sphingolipid-AUC and lowsphingolipid-AUC groups, which were visualized using the "ggplot2 R" tool.From the 14 cell subsets, the cell groups exhibiting high sphingolipid metabolic activity were selected and labeled as high sphingolipid metabolic cells.Subsequently, cluster analysis was performed once again on the selected cells, resulting in the further division of the high sphingolipid metabolism cell population into 7 subgroups (Figure 2B).Marker genes for each of the seven cell populations were analyzed, and bubble diagrams and volcano plots were used to visually represent the top five marker genes for each cell population (Figures 2C, D).Histograms were employed to demonstrate the distribution of these seven clusters within each cohort (Figure 2E).

Associations between SRG clusters and prognosis
To investigate the prognostic implications of SRG clusters, we initially utilized the GSE123139 dataset to calculate the ssGSEA scores of the marker genes.These marker genes were identified as the DEGs within the seven high sphingolipid metabolic clusters.Intriguingly, our analysis revealed that all seven clusters exhibited significantly higher scores in tumor samples compared to normal samples (Figure 3A).Subsequently, we employed the survminer R package to classify the samples from the TCGA dataset of SKCM into two groups based on their high and low scores of SCRGs.The classification was achieved by determining the optimal cut-off value and minimizing repetition.Importantly, we observed significant differences among all seven clusters between the high-and low-SRGs score groups (Figure 3B).For more comprehensive information regarding the relationship between SRG clusters and clinical characteristics, please refer to Supplementary Figure 1.

Identification of SRGs
To establish a prognostic signature for SKCM, we conducted a comparative analysis between normal and tumor samples to identify DEGs.From these DEGs, we identified marker genes that were significantly associated with gene clusters related to prognosis.To assess the prognostic value of each gene, we performed The flowchart of our study.
univariate Cox regression analysis and identified 10 genes associated with protective factors and 1 gene associated with risk values.To streamline the gene selection process and minimize gene redundancy, we employed Lasso-Cox regression analysis (Figure 4A).Using a stepwise regression method following multivariate Cox regression analysis, we created a risk signature that included eleven genes: IRX3, PLA2G2D, GBP1P1, FCGR2A, GALM, FERMT3, IGKJ5, IL15, IDO1, CMAHP, and HIVEP3 (Figures 4B, C).The risk scores for each sample were calculated based on the expression of these model genes and their corresponding Cox regression coefficients.Z-mean normalization was performed to compute the risk score for each sample, and patients were then classified into high-risk and low-risk clusters.Survival analysis using the Kaplan-Meier method was conducted in both the GSE19234 and TCGA cohorts, demonstrating that patients in the high-risk clusters had a worse prognosis compared to those in the low-risk clusters (Figure 4D).The model exhibited commendable AUC values in both cohorts, indicating its excellent predictive ability.To enhance the precision of our prognostic model, we integrated clinicopathological characteristics and risk scores into univariate and multivariate Cox regression analyses.This was done to reduce redundancy and improve accuracy.Our indepth analysis revealed a strong independent correlation between the risk signature and prognosis in SKCM, with statistical significance indicated by a p-value of less than 0.001 (Figure 4E).In conclusion, our findings highlight the importance of the risk signature as a valuable prognostic tool for SKCM, providing valuable insights into patient outcomes.

Nomogram development and pathway enrichment analysis
In addition, we have developed an innovative nomogram (Figure 5A) that combines the T-stage, N-stage, and risk score to provide a comprehensive prediction of survival outcomes.This nomogram demonstrated a strong predictive capacity for actual survival outcomes (Figure 5B).To further investigate the functional relevance of the eleven genes included in the risk profile, a gene set enrichment analysis was conducted.Interestingly, these genes showed significant associations with nine pathways (Figure 5C).Among these genes, IL15 exhibited a higher immune score in the low-risk group, while the immune scores of the other 10 genes were comparatively lower in the low-risk group compared to the highrisk group (Figure 5D).The innovative nomogram, along with the functional insights provided by the gene set enrichment analysis, contributes to a better understanding of the prognostic implications of the identified genes in SKCM.

Immune infiltrations landscape and risk gene-immunity association
Through our analysis, we have identified several protective genes (GBP1P1, FCGR2A, GALM, FERMT3, IGKJ5, IL15, IDO1, CMAHP, and HIVEP3) that exhibit a positive relationship with various immune infiltration cells.Conversely, the risk gene IRX3 is negatively associated with these immune cells (Figures 6A, B).Furthermore, correlation analysis with immune cells has revealed strong associations between the model genes and neutrophils and fibroblasts (Figure 6C).These findings are further supported by the positive correlation observed between the risk genes and the immune score, as well as the stromal score (Figure 6D).Pathway analysis has highlighted the close relationship between the model genes and key immune signaling pathways, including the JAK-STAT signaling pathway, intestinal immune network, B cell receptor signaling pathway, and toll-like receptor signaling pathway.Interestingly, a negative association has been observed between the risk gene IRX3 and each of these pathways (Figure 6E).

Immunotherapy response prediction of risk signature
To assess the predictive significance of our immune-checkpoint treatment signature, we evaluated its performance in two separate cohorts: GSE78220 and IMvigor210.Specifically, we focused on the outcomes observed in the IMvigor210 cohort, which consisted of 348 patients treated with anti-PD-L1 receptor blockers.These outcomes included partial response (PR), complete response (CR), progressing disease (PD), and stable disease (SD).In Figure 7A, we observed that the high-risk group had a higher proportion of patients with PD/SD compared to the low-risk group.This suggests that the high-risk group had a significantly worse outcome.On the other hand, patients who achieved a complete or partial response (CR/PR) had lower risk scores compared to those with stable disease or progressive disease (SD/PD).This finding indicates that patients with lower risk scores were more likely to exhibit a favorable treatment response.To validate our findings, we also analyzed the GSE78220 cohort.The results obtained were consistent with those from the IMvigor210 cohort.Patients who showed a partial or complete response had decreased risk scores and were less likely to be categorized as high-risk (Figure 7B).Interestingly, these distinct risk groupings were observed primarily in patients with Stage I+II disease, rather than those with Stage III+IV disease (Figure 7C).These findings underscore the potential of our immune checkpoint treatment signature as a predictive tool in T-cell immunotherapy.They suggest that the risk scores derived from this signature can help identify patients who are more likely to have a positive treatment response and may guide personalized treatment decisions in the context of immunecheckpoint blockade therapy.

Analysis of immunological infiltrations and consensus clustering
We employed unsupervised consensus clustering to investigate molecular subtypes based on the expression of SRGs (Signature-Related Genes) that comprise the risk signature.The TCGA cohort was divided into two clusters, with a k-value of 2 determined as the optimal clustering stability.The distribution of these distinct clusters is visualized in Figures 8A, B through a ridge plot.Cluster 1 (C1) exclusively comprised individuals from the lowrisk group, whereas cluster 2 (C2) consisted of both high-risk and low-risk patients, as illustrated in the Sankey diagram (Figure 8C).Survival analysis revealed that patients belonging to the C1 group exhibited significantly better outcomes compared to those in the C2 group (Figure 8D).Further analysis involved calculating TME (Tumor Microenvironment) scores for the different clusters.The C2 cluster demonstrated elevated immune, stromal, and estimate scores, indicating a more prominent presence of immune and stromal components within the tumor microenvironment.In contrast, the C1 cluster exhibited higher tumor purity (Figures 8E, F).Additionally, an examination of immune checkpoint inhibitors revealed a substantial correlation between the C2 cluster and the expression of a majority of immune checkpoints, as depicted in Supplementary Figure 2.This finding highlights a significant association between the C2 cluster and the immune checkpoint expression profile.
Taken together, these results provide insights into the molecular subtypes identified through unsupervised consensus clustering.They demonstrate that the C1 cluster, consisting of low-risk patients, is characterized by superior outcomes and higher tumor purity.On the other hand, the C2 cluster appears to be associated with increased immune and stromal activity, as well as a correlation with immune checkpoint expression.Such findings contribute to our understanding of the diverse tumor microenvironment and its i m p l i c a t i o n s f o r p a t i e n t p r o g n o s i s a n d p o t e n t i a l therapeutic strategies.

Drugs sensitivity
Following an evaluation of the efficacy of various chemotherapeutic drugs in distinct clusters, notable differences in drug response were observed.In cluster 2 (C2), patients displayed higher IC50 values when exposed to Bicalutamide, FH535, and Imatinib, indicating a reduced sensitivity to these drugs (Figure 9A).Conversely, individuals categorized under cluster 1 (C1) exhibited more favorable response rates to ATRA, Gefitinib, and other specific drugs (Figure 9B).These findings highlight that patients in C2 have a diminished response to certain chemotherapeutic agents, suggesting a potential resistance or reduced effectiveness of these drugs in this cluster.Conversely, patients in C1 demonstrate better responses to specific drugs, indicating a potential therapeutic benefit for these individuals.The evaluation of drug response in the different clusters provides valuable insights into the variability of treatment efficacy within distinct molecular subtypes.Understanding these differences can aid in the development of personalized treatment approaches and the selection of appropriate

Upregulated mRNA level of IRX3
To further validate the expression of IRX3 in melanoma, we selected 20 pairs of melanoma and corresponding normal tissues for qPCR verification.In tumor tissues, a significant upregulation of IRX3 was observed (Supplementary Figure 3A).To explore the function of IRX3 in melanoma, two melanoma cell lines, namely A375 and WM-115, were chosen for further experimental validation.Firstly, we silenced IRX3 in A375 and WM-115 cells using two siRNAs and further confirmed its knockdown efficiency through qPCR (Supplementary Figure 3B).

Silencing IRX3 inhibits proliferation, invasion, and metastasis while promoting apoptosis in melanoma cells
To investigate the role of IRX3 in melanoma, we performed a colony formation assay on A375 and WM-115 melanoma cells in the NC and si-IRX3 groups (Figures 10A, B).The results showed that silencing of IRX3 led to smaller colonies in both A375 and WM-115 cells, suggesting that IRX3 silencing inhibits melanoma cell proliferation.This result was further confirmed by the CCK-8 assay (Figure 10C).To explore the effect of IRX3 on melanoma cell migration and invasion, we conducted scratch assay and transwell assay.The results showed that silencing of IRX3 significantly inhibited the invasion ability of A375 and WM-115 cells (Figures 10D-G).Apoptosis plays a crucial role in the malignant behavior of many tumors.We further investigated the effect of IRX3 on tumor cell apoptosis.Apoptosis assay revealed that silencing of IRX3 promoted apoptosis in both A375 and WM-115 cell lines (Figures 10H, I).Inhibition of melanoma cell proliferation, invasion, and migration was observed upon silencing of IRX3, while simultaneously promoting apoptosis.

Discussion
Melanoma, a malignancy known for its aggressive behavior, is influenced by disruptions in lipid metabolism that can significantly influence its onset and advancement (1).Extensive investigation in the field of melanoma has focused on the exploration of multiple genes and metabolites implicated in sphingolipid metabolism.Notably, melanoma showcases aberrant expression or modified activity of specific enzymes that play a pivotal role in the regulation of sphingolipid metabolism (21,22).An illustrative example involves Ceramide synthase (CERS), a crucial component within the sphingolipid metabolism pathway.Perturbations in the functionality of CERS can lead to the accumulation of specific sphingomyelin species in melanoma cells (23).Moreover, the progression and advancement of melanoma have been linked to other enzymes participating in sphingolipid metabolism.Sphingosine kinase (SPHK) and glycosphingolipid transferase (GSLT) are among these enzymes that have been implicated in the pathogenesis of melanoma (24,25).Furthermore, extensive research studies have highlighted a robust association between the advancement and prognosis of melanoma in patients and the levels of sphingolipid metabolism products, such as sphingomyelin and ceramide, found within their biological systems (26).As an example, elevated concentrations of serum sphingomyelin have been linked to decreased overall survival in individuals afflicted with melanoma (19).
Our study focused on exploring potential associations between sphingolipid metabolism and melanoma to explore the molecular characteristics derived from SRGs and their clinical relevance.We also investigated the role of signaling pathways in cancer cell proliferation, metastasis, and drug resistance, as well as the effectiveness of immune therapy including immune checkpoint pathways such as PD-1/PD-L1 and CTLA-4, and the activation of immune cells such as NK cells and tumor-associated macrophages in the immune microenvironment.We leveraged resources such as the Human Protein Atlas to identify four melanoma-specific genes involved in sphingolipid metabolism.Through single-cell RNAsequencing analysis, we identified additional genes associated with sphingolipid metabolism.Differential analysis was conducted, followed by univariate Cox regression, lasso regression, and multivariate Cox regression.This comprehensive approach enabled the identification of a novel risk signature comprising eleven genes (IRX3, PLA2G2D, GBP1P1, FCGR2A, GALM, FERMT3, IGKJ5, IL15, IDO1, CMAHP, and HIVEP3).Notably, the risk score derived from the signature demonstrated independent predictive value for overall survival, as observed in both univariate and multivariate Cox regression models.The accuracy of the prediction was further confirmed by excellent agreement between predicted and actual outcomes using a nomogram for melanoma patients.Moreover, our in vitro experiments revealed that the silencing of IRX3 in melanoma cells resulted in the suppression of cell proliferation, invasion, and migration while promoting apoptosis.
Immunotherapy has emerged as a promising therapeutic strategy in the field of medicine.It aims to counteract the immune evasion mechanism exhibited by tumors, thereby activating the patient's immune cells to target and eradicate cancer cells (27).The products generated through lipid metabolism play a vital role in shaping the immunological microenvironment by impacting diverse aspects of immune cell biology.These products possess the capacity to directly modulate the proliferation, differentiation, function, and activity of immune c e l l s .B y e x e r t i n g c o n t r o l o v e r t h e i m m u n o l o g i c a l microenvironment, they substantially contribute to the behavior of immune cells (28,29).Several lipid metabolism products, such as fatty acids, triglycerides, and cholesterol, have a significant impact on the differentiation and function of various immune cells, including T cells, B cells, and macrophages.As a result, they influence the immune response of the body towards infections and tumors.Moreover, leukotrienes, categorized as lipid metabolism products, actively participate in the chemotaxis and adhesion mechanisms of immune cells, thereby regulating their migration and aggregation (30).In the context of the immune microenvironment, various cells present, such as tumor cells, macrophages, and lymphocytes, are known to produce specific lipid metabolites that play a key role in regulating immune responses (31).In recent times, there has been a remarkable increase in enthusiasm regarding the utilization of immunotherapy as a viable therapeutic approach for patients diagnosed with melanoma.The predominant immunotherapy treatments employed for melanoma include anti-CTLA-4 and anti-PD-1 antibodies, which function by targeting and inhibiting specific receptors found on T cells.Anti-CTLA-4 antibody therapy amplifies T cell functionality by obstructing the CTLA-4 receptor, resulting in escalated T cell assault on malignant cells (32).The application of anti-PD-1 antibody therapy impedes the PD-1 receptor found on T cells, facilitating persistent cancer cell targeting (33).Immunotherapy has exhibited promising efficacy in the management of melanoma; nevertheless, its appropriateness varies among individuals.Hence, it is crucial to employ rigorous screening and evaluation protocols before initiating immunotherapy, aiming to determine its suitability for each patient.The advent of PD-1/PD-L1 inhibitors has introduced novel avenues and challenges in melanoma treatment.We aspire for the utilization of PD-1/PD-L1 inhibitors to assume a pivotal role in the adjuvant therapy of high-risk melanoma (newly diagnosed), thereby furnishing surgeons with precise treatment strategies through pertinent clinical investigations.Ultimately, our objective is to enhance patients' quality of life and improve the disease prognosis (34,35).
The identified signature establishes a link between the IRX3 gene and an adverse prognosis in melanoma.The protein derived from the IRX3 gene assumes a critical role in embryonic development and normal physiological processes in adults (36).Numerous studies have substantiated a notable correlation between the presence, advancement, and outcome of malignant tumors and the expression of the IRX3 gene.Notably, investigations have revealed that heightened expression of IRX3 is strongly associated with diminished survival rates and heightened vulnerability to metastatic disease in individuals diagnosed with melanoma (37).IIRX3 has emerged as a promising diagnostic indicator for several types of cancer, including prostate, colorectal, and gastric cancer.Its expression level can serve as a valuable tool for evaluating tumor prognosis and treatment responsiveness.By detecting the degree of IRX3 expression, clinicians can assess the potential outcomes and responses to therapeutic interventions in patients with these malignancies (38).The role of IRX3 in tumor management, particularly in melanoma, is of utmost significance.Extensive research has elucidated that downregulating IRX3 expression significantly reduces the invasiveness and proliferation of melanoma cells.Moreover, it enhances the sensitivity of these cells towards chemotherapy agents.These findings highlight the therapeutic potential of targeting IRX3 as a means to attenuate tumor progression and improve treatment outcomes in melanoma patients (39, 40).The overexpression of IRX3 in melanoma indicates its plausible involvement in tumorigenesis and disease advancement.These findings align with previous investigations on IRX3 in colorectal cancer and breast cancer.Nevertheless, our study pioneers the examination of IRX3's distinctive function in melanoma.The suppressive impact of IRX3 gene silencing on melanoma cell behavior implies the potential utility of targeting IRX3 as a therapeutic strategy with promising prospects.To sum up, our research successfully developed a robust diagnostic and prognostic model for melanoma while unveiling the upregulation of IRX3 in this particular disease.The potential therapeutic value of targeting IRX3 in melanoma treatment holds great promise.However, additional investigations are warranted to gain a comprehensive comprehension of the underlying mechanisms through which IRX3 influences melanoma development.Furthermore, exploring the clinical applicability of IRX3 targeting necessitates further exploration.

Conclusions
In conclusion, using the sphingolipid-related model can adequately classify patients for prognosis and immunological assessment of patients in melanoma.Our research findings could provide valuable insights into detecting and treating melanoma patients.

2
FIGURE 2 Identification of SRG clusters based on SKCM patient scRNA data.(A) UMAP plots of the expression of the sphingolipid metabolism marker genes and the distribution of 14 clusters.(B) UMAP plots of distributions of 7 high sphingolipid metabolizing cells after clustering.(C) Top-5 marker gene expression of subgroups on a volcano map.(D) Top-5 marker gene expression of subgroups represented in a bubble diagram.(E) Calculating cell numbers as well as neighboring tissue and subgroups in cancer.
FIGURE 3 Based on SRGs clusters, the GSVA analysis.(A) ssGSVA score comparison between tumor samples and normal samples based on each cluster.(B) K-M curves of the high and low SRG score groups in the SRG clusters.****P < 0.0001.

4 A
FIGURE 4 A brand-new risk signature built using several SRGs.(A) Each independent variable's trajectory and distribution for the lambda.(B) Circle plot showing each gene in the risk signature.(C) The multivariate Cox coefficients for each gene in the risk signature.(D) K-M and ROC curves of the risk signature in GSE19234 and TCGA cohort.(E) Results of univariate and multivariate Cox regression analysis based on risk score and clinicopathologic features.

5
FIGURE 5 Creation of a new nomogram that incorporates the risk signature.(A) Construction of the nomogram integrating the T, N-stage, and risk score.(B) Decision curve for nomogram.(C) The gene set enrichment analysis was performed.(D) The immune scores of the 11 genes.***P < 0.001, ****P < 0.0001.

6
FIGURE 6 Analysis of immunological infiltrations.(A) The correlation analysis between risk genes and immunity.(B, C) Correlations between eleven hub genes and 22 immune-related cells.(D) Correlations between the eleven genes and immune score, stromal score, and estimate score.(E) The correlation analysis between eleven hub genes and signaling pathway.*P < 0.05, **P < 0.01, ***P < 0.001.

7
FIGURE 7 Prediction of immunotherapy response (A) Prognoses differences between risk subgroups in the IMvigor210 cohort; Differences between immunotherapy responses; Distribution of immunotherapy response.(B) Prognostic differences between risk subgroups in the GSE78220 cohort; Differences between immunotherapy responses; Distribution of immunotherapy response.(C) Prognostic differences between risk subgroups based on early-stage (stage I-II) and late-stage (stage III-IV) patients in the two cohorts.*P < 0.05, ****P < 0.0001.

8
FIGURE 8 Consensus based on the expression of SRGs.(A) A Stratification of melanoma patients into two clusters according to the consensus clustering matrix (k = 2).(B) PCA depicted the distribution for clusters.(C) The Sankey diagram of the connection between clusters and high-low risk groups.(D) Survival analysis based on the two clusters.(E, F) ESTIMATEScore, SromalScore, ImmuneScore, and TumorPurity difference between two clusters.***P < 0.001.
FIGURE 10 Silencing IRX3 Inhibits Proliferation, Invasion, and Metastasis while Promoting Apoptosis in Melanoma Cells.(A) A colony formation assay was performed on A375 and WM-115 melanoma cells in the NC and si-IRX3 groups.Smaller colonies were observed in the si-IRX3 group, indicating that IRX3 silencing inhibits melanoma cell proliferation.(B) Quantification of colony formation assay results showed a decrease in colony size in the si-IRX3 group compared to the NC group.(C) CCK-8 assay further confirmed the inhibitory effect of IRX3 silencing on melanoma cell proliferation.(D) Quantification of scratch assay results showing a decrease in wound closure percentage in the si-IRX3 group compared to the NC group.(E) Scratch assay revealed a decrease in 17 the migration ability of A375 and WM-115 cell in the si-IRX3 group compared to the NC group.(F) Transwell assay demonstrated a decrease in the invasion ability of A375 and WM-115 cells in the si-IRX3 group compared to the NC group.(G) Quantification of transwell assay results shows a decrease in the number of invading cells in the si-IRX3 group compared to the NC group.(H) Apoptosis assay revealed an increase in apoptosis in both A375 and WM-115 cell lines upon IRX3 silencing.(I) Quantification of apoptosis assay results shows an increase in the percentage of apoptotic cells in the si-IRX3 group compared to the NC group.*P < 0.05, **P < 0.01, ***P < 0.001.