Bioinformatic Exploration for Prognostic Significance of Sphingolipid Metabolism-Related Genes in Invasive Ductal Carcinoma Using the Cancer Genome Atlas Cohort

Introduction

Breast cancer is the most common and life-threatening malignancy in females worldwide.¹ Breast carcinoma is the most prevalent malignant type and is classified as carcinoma in-situ and invasive breast cancer.² Invasive ductal carcinoma (IDC) is the most common type of invasive breast cancer, accounting for up to 80% of diagnosed breast cancer cases.³ There are clinical prognostic biomarkers for breast cancer, including size, histological grade, and estrogen receptor, progesterone receptor, and human epidermal growth factor receptor 2 status.⁴ In particular, molecular sub-classification systems such as receptors play an important role in clinical therapeutic strategy.⁵ However, despite the availability of therapeutic strategies for breast cancer based on molecular subtype, breast cancer is still not overcome.¹ Therefore, further studies are required to identify promising molecular biomarkers that can provide new treatment avenues.

Recent advances in genomic profiling using next generation sequencing have made it possible to identify the genetic characteristics of disease, particularly in cancer. Several large-scale cancer genome studies have been conducted, and The Cancer Genome Atlas (TCGA) is a research consortium that may be used to investigate genes in different cancer types.⁶ Moreover, TCGA may be used to investigate specific histological types of cancer, such as female IDC, by utilizing histological data and clinical parameters.

Sphingolipid metabolism is a highly-regulated intracellular process that controls the synthesis and degradation of bioactive lipids, including ceramide and sphingosine-1-phosphate,⁷ which plays an important role in biological processes such as angiogenesis, ageing, cancer biology, degenerative diseases, diabetes, immune responses, and inflammation.⁸ Accumulating evidences revealed that dysregulated sphingolipid metabolism-related genes are implicated in human breast cancer. For example, LAG1 longevity assurance homolog 2 (LASS2) and LASS6 mRNA levels are increased in breast cancer tissues compared with matched normal tissues,⁹ and overexpressed LASS4 and LASS6 reduced cell proliferation.¹⁰ Furthermore, upregulation of sphingosine kinase 1 (SPHK1) was significantly associated with poor prognosis¹¹ and metastasis.¹² Moreover, downregulation of 15-hydroxyprostaglandin dehydrogenase (HPGD) has an unfavorable effect on the overall survival of patients with triple negative breast cancer.¹³ While the altered expression of various sphingolipid metabolism-related genes in breast cancer and their potentials as prognostic factors have been reported in the aforementioned studies, few studies have investigated sphingolipid metabolism-related genes in female IDC.

Therefore, the aim of the present study was to investigate the expression of sphingolipid metabolism-related genes in female IDC, as well as to evaluate their prognostic significance, using gene expression RNAseq data obtained from TGCA breast invasive carcinoma (TCGA BRCA) datasets.

Materials and Methods

Gene Expression Datasets and Cluster Analysis

TCGA BRCA gene expression RNAseq datasets (level 3, dataset ID: TCGA.BRCA.sampleMap/HiSeqV2) and clinical parameters (dataset IDs: TCGA.BRCA.sampleMap/BRCA_clinicalMatrix and survival/BRCA_survival.txt) were downloaded from the UCSC Xena public database (https://xena.ucsc.edu). TCGA BRCA dataset consisted of 1218 samples, including 1097 primary tumor tissues, 7 metastatic tumor tissues, and 114 normal solid tissues (NST). NST were taken from normal tissues adjacent to the tumor. To analyze the RNAseq data of female IDC, female IDC datasets were sorted from TCGA BRCA using clinical parameters. The mRNA expression of the sphingolipid metabolism-related genes¹⁴ was identified from the female IDC dataset. In order to identify genes with ≥2-fold changes (2FC) in mRNA expression levels between IDC and NST, the difference between average values of the two groups was calculated, and genes with a value greater than 1 were selected. This study met the publication guidelines for using TCGA datasets (http://www.cancer.gov/about-nci/organization/ccg/research/structrual-genomics/tcga/ using-tcga/citing-tcga). Cluster 3.0¹⁵ was used for cluster analysis, and samples with statistically similar gene expression were classified into groups. TreeView 1.6 (www.eisenlab.org/eisen) was used to visualize the resulting heat map. The mRNA expression levels of the heat maps were scaled (quantile normalization¹⁶ and median-centered) within columns for visualization.

Survival Analysis

Survival data (death event and survival time) were available for 655 female IDC patients. For survival analysis, the mean gene expression value of the selected sphingolipid metabolism-related genes was used as a cutoff to divide the patients into high- and low-expression groups. Survival analysis was performed using the Kaplan–Meier method, and the log rank test was used to identify statistically significant differences between the two groups.

Statistical Analysis

Statistical analysis was performed using SPSS software (version 25.0; IBM SPSS, Armonk, NY, USA). The Kolmogorov–Smirnov test was performed to assess normality. Differences in mRNA expression levels between groups were analyzed using the Student’s t-test. The associations between clinicopathological parameters and the dysregulated sphingolipid metabolism-related genes were analyzed using Chi-square test or Fisher’s exact test for categorical variables. Correlation analysis between inter-individual mRNA expression levels of the sphingolipid metabolism-related genes was performed using the Spearman correlation coefficient analysis for continuous variables. P<0.05 was considered to indicate a statistically significant difference.

Results

Sphingolipid Metabolism-Related Genes are Dysregulated in Female IDC Compared with NST in TCGA BRCA

The heat map revealed relative mRNA expression levels of various sphingolipid metabolism-related genes in female IDC tissues and NST (Figure 1). The Student’s t-test revealed that a total of 36 sphingolipid metabolism-related genes were significantly altered in female IDC tissues compared with NST. The results are presented as a heat map (P<0.05; Figure 2A). When a ≥2FC in mRNA expression was used a cutoff, sphingomyelin phosphodiesterase acid like 3B (SMPDL3B), beta-1,4-N-acetyl-galactosaminyltransferase 1 (B4GALNT1), lysophosphatidic acid receptor 2 (LPAR2), and LAG1 longevity assurance homolog 2 (LASS2) were significantly upregulated in female IDC, whereas LAG1 longevity assurance homolog 3 (LASS3), lysophosphatidic acid receptor 1 (LPAR1), beta-1,4-galactosyltransferase 6 (B4GALT6), galactose-3-O-sulfotransferase 1 (GAL3ST1), HPGD, ST8 alpha-N-acetyl-neuraminide alpha-2,8-sialyltransferase 1 (ST8SIA1), UDP glycosyltransferase 8 (UGT8), sphingosine-1-phosphate receptor 1 (S1PR1) were significantly downregulated (P<0.001, |2FC| ≥1.0; Figure 2B), compared with NST. These results suggested that the aforementioned 12 dysregulated sphingolipid metabolism-related genes may play a crucial role in pathophysiology of female IDC.

Figure 1 Heat map showing the relative mRNA expression levels of the sphingolipid metabolism-related genes in female IDC tissues and NST obtained from TCGA BRCA cohort. In the data shown in matrix format, each row represents an individual gene and each column represents a single tissue. Each cell in the matrix represents the relative mRNA expression level of a gene feature in an individual tissue. The red and green in the cells reflects relatively high and low expression levels, respectively, as indicated by the scale bar. The samples are sorted into the NST group on the left and the IDC group on the right. Each cells is arranged in descending order of the mean difference between the scaled mRNA expression levels of each gene in the NST and IDC groups. Abbreviations: TCGA BRCA, The Cancer Genome Atlas breast invasive carcinoma; IDC, invasive ductal carcinoma; NST, normal solid tissues.

Figure 2 (A) Heat map showing significantly altered mRNA expression of the sphingolipid metabolism‑related genes in female IDC tissues compared with female NST obtained from TCGA BRCA cohort. In the data shown in matrix format, each row represents an individual gene and each column represents a single tissue. Each cell in the matrix represents the relative mRNA expression level of a gene feature in an individual tissue. The red and green in the cells represent relatively high and low expression levels, respectively, as indicated by the scale bar. The samples are sorted into the NST group on the left and the IDC group on the right. Each cells is arranged in descending order of the mean difference between the scaled mRNA expression levels of each gene in the NST and IDC groups. *Student’s t-test, P<0.05 (NST versus IDC). (B) Relatively altered mRNA expression levels of various sphingolipid metabolism-related genes in IDC samples obtained from TCGA BRCA. **Student’s t-test, P<0.001, |2FC| ≥ 1.0. Abbreviations: TCGA BRCA, The Cancer Genome Atlas breast invasive carcinoma; IDC, invasive ductal carcinoma; NST, normal solid tissues.

KM Survival Analysis Identified SMPDL3B as a Prognostic Biomarker in Female IDC

Among the altered sphingolipid metabolism-related genes in female IDC, the KM survival analysis and Log rank test demonstrated that higher SMPDL3B mRNA expression levels were found to be associated with favorable overall survival in female IDC (Figure 3, P=0.003).

Figure 3 Survival analysis of the dysregulated sphingolipid metabolism-related genes in IDC samples obtained from TCGA BRCA cohort. Kaplan–Meier estimates of female patients with IDC according to the relative mRNA expression values of GAL3ST1, SMPDL3B, LASS3, LASS2, B4GALT6, HPGD, UGT8, ST8SIA1, S1PR1, LPAR1, B4GALNT1, and LPAR2. Abbreviations: TCGA BRCA, The Cancer Genome Atlas breast invasive carcinoma; IDC, invasive ductal carcinoma.

Altered mRNA Expression Levels of Sphingolipid Metabolism-Related Genes are Associated with Clinicopathological Parameters in Female IDC

To investigate the clinicopathological implications of the dysregulated sphingolipid metabolism-related genes, Chi-square test or Fisher’s exact test were performed. The clinicopathological characteristics of female IDC patients of TCGA BRCA are presented in Tables 1 and 2. The SMPDL3B mRNA expression level was significantly associated with age, estrogen receptor (ER) status, progesterone receptor (PR) status, and histological subtypes according to immunohistochemistry (IHC) (P=0.002, P<0.001, P<0.001, and P<0.001, respectively, Table 1). The B4GALNT1 mRNA expression level was significantly associated with T stage, ER status, PR status, and histological subtype according to IHC (P=0.004, P<0.001, P<0.001, and P<0.001, respectively, Table 1). The LPAR2 mRNA expression level was significantly associated with T stage, ER status, PR status, and histological subtype according to IHC (P<0.001, P<0.001, P<0.001, and P<0.001, respectively, Table 1). The LASS2 mRNA expression level was significantly associated with age, ER status, PR status, and histological subtype according to IHC (P=0.003, P<0.001, P<0.001, and P<0.001, respectively, Table 1). Moreover, the GAL3ST1 mRNA expression level was significantly associated with ER status, PR status, and histological subtype according to IHC (P=0.003, P=0.011, and P<0.001, respectively, Table 2). The LASS3 mRNA expression level was significantly associated with age and histological subtype according to IHC (P=0.036 and P=0.041, respectively, Table 2). The B4GALT6 mRNA expression level was significantly associated with ER status, PR status, epidermal growth factor receptor type 2 (HER2), and histological subtype according to IHC (P=0.006, P<0.001, P=0.003, and P<0.001, respectively, Table 2). The HPGD mRNA expression level was significantly associated with age, ER status, PR status, and histological subtype according to IHC (P=0.024, P<0.001, P=0.007, and P<0.001, respectively, Table 2). The UGT8 mRNA expression level was significantly associated with age, N stage, ER status, PR status, HER2 and histological subtype according to IHC (P=0.003, P=0.029, P<0.001, P<0.001, P=0.003, and P<0.001, respectively, Table 2). The ST8SIA1 mRNA expression level was significantly associated with age, T stage, ER status, PR status, and histological subtype according to IHC (P=0.019, P=0.008, P<0.001, P<0.001, and P<0.001, respectively, Table 2). The S1PR1 mRNA expression level was significantly associated with M stage and histological subtype according to IHC (P=0.005 and P<0.001, respectively, Table 2). The LPAR1 mRNA expression level was significantly associated with T stage, ER status, PR status, and histological subtype according to IHC (P<0.001, P<0.001, P<0.001, and P<0.001, respectively, Table 1). Interestingly, we found that lower expression levels of SMPDL3B were 2.676 times more frequent in ER-positive IDC than in ER-negative (odds ratio (OR), 2.676; 95% confidence interval (CI), 1.821 to 3.932; P<0.001). Moreover, we found that lower expression levels of SMPDL3B were 1.898 times more frequent in PR-positive IDC than in PR-negative (OR, 1.898; 95% CI, 1.357 to 2.655; P<0.001).

Table 1 Upregulated mRNA Expression Levels of Sphingolipid Metabolism-Related Genes in Relation to Clinicopathological Parameters of IDC

Table 2 Downregulated mRNA Expression Levels of Sphingolipid Metabolism-Related Genes in Relation to Clinicopathological Parameters of IDC

Altered mRNA Expression Levels of Sphingolipid Metabolism-Related Genes Have Inter-Individual Correlations in Female IDC

Spearman correlation coefficient analysis was performed to explore the correlation among the significantly altered sphingolipid metabolisms-related genes in female IDC tissues. Prior to the correlation analysis, the mRNA expression levels of the 12 significantly dysregulated sphingolipid metabolism-related genes in female IDC samples were retrieved from TCGA BRCA cohort. The resulting analysis uncovered 43 significant correlations among the mRNA expression levels of specific sphingolipid metabolism-related genes in female IDC samples (Table 3).

Table 3 Spearman Correlation Analysis Between Inter-Individual Components of Sphingolipid Metabolism-Related Genes

Discussion

Sphingolipid metabolism regulates diverse biological processes, including breast cancer pathophysiology.¹⁷ IDC of the breast cancer, arises from epithelial cells in the inner lining of the milk ducts, and is the most common type of invasive breast cancer that microscopically penetrates through the epithelium of the ducts into the stroma in breast tissue.^2,3 As the importance of precision medicine becomes apparent, the number of studies using genomic profiles to identify novel prognostic and therapeutic targets is steadily increasing.^18,19 Although methods for precisely targeting or imaging such as peptide-based SPECT radiopharmaceuticals have been studied using the discovered target for diagnosis or treatment of IDC,^20–22 its prognosis is still unfavorable. In the present study, gene expression RNAseq data from TCGA BRCA dataset were used to evaluate whether sphingolipid metabolism-related genes are significantly dysregulated, and whether they may serve as potential prognostic indicators in female IDC.

The present study investigated altered genes with ≥2FC in mRNA expression levels in IDC compared with NST. SMPDL3B, B4GALNT1, LPAR2, and LASS2 were significantly upregulated and the other genes (LASS3, LPAR1, B4GALT6, GAL3ST1, HPGD, ST8SIA1, UGT8, and S1PR1) were significantly downregulated (Figure 2). Among the dysregulated sphingolipid metabolism-related genes, only SMPDL3B had prognostic significance in female IDC (Figure 3). SMPDL3B is one of the lipid raft enzymes that regulates lipid composition and fluidity in the plasma membrane of macrophages.²³ Interestingly, SMPDL3B modulates podocyte migration and apoptosis, and activates integrin in podocytes.²⁴ Moreover SMPDL3B inversely regulates ceramide-1-phosphate (C1P) levels by interacting with ceramide kinase (CERK) in human podocytes.²⁵ In addition, it has been reported that C1P is involved in enhancement of cancer cell growth, migration, and survival.²⁶ However, CERK was downregulated in IDC (Figure 1), suggesting SMPDL3B and C1P are regulated differently in IDC and podocytes. Notably, a recent report demonstrated that increased SMPDL3B was associated improved prognosis in prostate cancer.²⁷ The present study demonstrated that increased mRNA expression levels of SMPDL3B is significantly associated with improved prognosis (Figure 3). These results suggested that SMPDL3B may function as a tumor suppressor gene.

In the present study, we found that SMPDL3B and LPAR2 were upregulated (Figure 2B), and a significant positive correlation was identified between the genes (Table 3). Additionally, higher LPAR2 mRNA expression levels tended to show favorable overall survival in female IDC (Figure 3). However, LPAR2 showed aggressive characteristics of gastric cancer cells, such as inducing cell migration²⁸ and LPAR2-KO had significantly suppressed the invasion of ovarian cancer cells.²⁹ Although more studies on the mechanism of each gene in IDC are needed, the aforementioned results suggested that SMPDL3B and LPAR2 have different regulatory mechanisms and functions. On the other hand, the downregulated sphingolipid metabolism-related genes such as HPGD, S1PR1, and LPAR1, which are inversely correlated with the SMPDL3B, have been identified (Table 3). Additionally, higher HPGD mRNA expression levels tended to indicate a favorable prognosis (Figure 3). However, it was previously shown that HPGD was significantly upregulated and higher HPGD is significantly associated with a poor overall survival and the increased risk of disease relapse in breast cancer.³⁰ In order to determine the dysregulation and function of HPGD in breast cancer, it is necessary to analyze it by breast cancer subtype.

Recent studies reported significant correlations among sphingolipid metabolism-related genes, such as ceramide synthases. For example, significant correlations were determined between LASS2 and LASS4/LASS6, and between LASS4 and LASS6 in breast cancer.⁹ Moreover, significant correlations between LASS2 and LASS4, and between LASS5 and LASS4/LASS6 were identified in various colorectal cancer cohorts.³¹ Interestingly, in the present study, only the correlation between LASS2 and LASS3 was significant in female IDC (Table 3).

Although further studies are required to elucidate the underlying mechanisms involved in the inter-individual correlation analysis in female IDC, our results suggested that networks of sphingolipid metabolism-related genes may be implicated in the pathophysiology of the disease.

Conclusions

To the best of our knowledge, the present study was first to identify the significantly dysregulated sphingolipid metabolism-related gene, SMPDL3B as a novel prognostic biomarker for female IDC using TCGA BRCA datasets. Our results suggested that SMPDL3B plays an important role in the pathophysiology of female IDC, and may serve as a potential prognostic biomarker as well as promising therapeutic target for the disease.