Key Points

• CYSLTR1 is downregulated in breast tumors.

• TNBC exhibited less CYSLTR1 than Luminals and HER2 subtypes.

• Low CYSLTR1 expression was associated with worse survival in breast cancer patients.

• Low CYSLTR1 expression was associated with worse survival in TNBC.

Introduction

It is widely known that the severity of breast cancer (BCa) results from a multitude of extrinsic and intrinsic factors, including tumor heterogeneity, which has been identified as the most relevant cause of poor outcome in patients with different subtypes of BCa (1). BCa can be stratified in hormone-dependent tumors, with receptor of human epidermal growth factor 2 (HER2) overexpression or triple negative (TNBC) according to immunohistochemical (IHC) staining for estrogen receptor (ER), progesterone receptor (PR), HER2, and the cell proliferation marker Ki-67. Another form of classification widely used in clinical practice is based on the transcriptomic profiles in Luminal A, Luminal B, HER2+, basal-like, normal-like, and claudin-low (2, 3, 4). This molecular classification has been confirmed by several research groups in different populations of patients with BCa (5, 6, 7). Patients with TNBC do not benefit from hormone therapy or targeted therapies commonly used in luminal and HER2+ cases. This lack of therapeutic options increases the chances of tumor recurrence, leading to a high mortality rate (8, 9).

On the other hand, BCa is strongly associated with inflammation and the release of signaling molecules derived from arachidonic acid, such as leukotrienes, as well as G protein-coupled receptors in the tumor microenvironment (TMev), which results in mediation of allergic, infectious, and inflammatory reactions through a phosphatidylinositol-calcium second messenger cascade (10, 11, 12). Among these messengers, Cysteinyl Leukotriene Receptor 1 (CYSLTR1) is implicated in mediating bronchoconstriction and asthma, and its dysregulation may be of concern in inflammation-related neoplasms (13). For example, in colorectal tumor cells, overexpression of CYSLTR1 is associated with proliferation, survival, and migration, as well as a poor prognosis in patients with colorectal adenocarcinoma (14). Furthermore, in patients with breast tumors, high expression levels of CYSLTR1 and low levels of CYSLTR2 were correlated with high mortality rates in univariate analyses for 144 patients (15). Another study suggests that CYSLTR1 is positively correlated with clinical features, such as tumor size, histologic type, lymph node metastasis, and TNM staging in a BCa population of 90 subjects (16).

Data mining represents a useful aproach to strengthen the knowledge and status of malignant neoplasms. With the era of omics, there is an increasing amount of genomic, transcriptomic, proteomic, and epigenetic data generated by high performance technologies available in public databases. Many of the studies deposited on these platforms have helped to characterize intrinsic cancer subtypes, predict survival, and therapeutic responses, generating a large amount of molecular biomarkers for BCa. There are only a few studies, with limited samples, that have explored the role of CYSLTR1 in women diagnosed with BCa. Therefore, the central objective of this study was to explore the status of CYSLTR1 according to expression levels and its potential prognostic value in BCa using datasets deposited in public repositories.

Materials and Methods

UALCAN and GENT2

UALCAN is a user-friendly online platform that provides easy access to OMICS cancer data. Thus, it allows for easy expression profiling of possible biomarkers, in associations with survival and gene regulation data, rendering a robust profile analysis (17). With this tool, we identified the difference between the CYSLTR1 expression levels of normal and breast tumor tissues. Moreover, in order to confirm our results, we accessed GENT2, a new tool focused on the expression analysis of normal and tumor tissue samples (18).

cBioPortal

The TCGA database (Firehose Legacy) was accessed through the cBioPortal platform to select mRNA expression Z-scores related to 1.108 samples (log RNA Seq V2 RSEM) with a ±2 threshold (19, 20). Clinical pathological data were obtained and cross-linked with CYSLTR1 expression data. Male cases (n = 16) and those who had no information of CYSLTR1 levels (n = 4) were excluded, resulting in 1.088 patients to be assessed.

bc-GenExMiner

The bc-GenExMiner v.5 is a microarray and RNA-seq data-mining tool containing data of BCa patients only. Three analysis modules were explored: Correlation, expression, and prognosis (21, 22). For this study, we considered only RNA-seq data, excluding samples from TCGA.

Kaplan–Meier Plotter

Kaplan–Meier (KM) Plotter is a publicly available platform that hosts data of 21 different types of cancer and contains Affymetrix gene signatures (probes of 20.129 genes) of 3.421 patients (23). For this study, we selected the best probe option corresponding to the CYSLTR1 gene: 230866_at; P2RY10 gene: 236280_at; and XCR1 gene: 221468_at. The overall survival (OS) and relapse-free survival (RFS), adjusted for 120 months’ total follow-up time, were available for all of them. The patients were also stratified by high and low expression of the target gene as the best cut-off between the lower and upper quartile was selected. Analyses were performed according to all deposited cases and only with the basal-like subtype, considering the prognostic value and its impact on poor clinical outcome.

Metascape

Metascape is a user-friendly tool for omics data analysis (24). Here, we accessed CYSLTR1 co-expressed genes previously obtained on bc-GenExMiner for interaction analysis. The protein network data was downloaded and analyzed on Cytoscape v.8.0 (25).

Geo Database

The Geo Database is a microarray and RNA-seq data deposit platform. In order to analyze the expression profile of CYSLTR1 in different subtypes of breast tumors, we accessed the GSE76275 and GSE96058 files (26).

Gene Co-expression Network

The co-expression analysis was conducted using RNA-seq data of TNBC from the bc-GenExMiner v.5 database. A correlation value >0.7 was used as a cut-off, then the data was accessed using the String platform to generate CYSLTR1 co-expressed genes network data and to export it to Cytoscape v8.0 software to select the genes with close interactions with CYLSTR1. In addition, the co-expressed gene list was also accessed using the Metascape software in order to conduct enrichment analyses.

Comparative Toxicogenomics Database

The Comparative Toxicogenomics Database (CTD) is a publicly available tool for manually curated information about chemical interactions with genes, proteins, and chemical relationships with diseases (27). The CTD was accessed to obtain potential drugs capable of interacting positively or negatively with CYSLTR1.

mirTarBase Repository

By using the mirTarBase database we accessed the prediction of experimentally validated miRNAs targeting CYSTLR1 and significant co-expressed genes (28). In addition, we carried out survival analyses of both genes and best-predicted miRNAs in TNBC population.

Statistical Analysis

For platforms with integrated statistical capabilities, the analyzes were performed as described in the topics above. For additional data, analyses were conducted with the Statistical Package for Social Sciences (SPSS) version 25.0 (IBM Inc., Armonk, NY USA) or GraphPad v.7 (California, USA). The chi-square or Fisher’s exact test was applied to compare categorical variables. For univariate and multivariate analysis, the Cox regression method was used. All groups were tested for Gaussian distribution. The Mann–Whitney or t test was used to assess the difference between two groups, and ANOVA (analysis of variance) or Kruskal–Wallis for more than two groups. For survival analysis, survival curves were performed by KM method and compared using log-rank test; additionally, Cox regression univariate and multivariate were performed calculating hazard ratio (HR) with 95% confident interval. A significance level of 5% was adopted.

Results

CYSLTR1 is Downregulated in Breast Tumors and Correlated with Clinical Pathological Parameters

Our study evaluated data from different platforms (Firehose Legacy, cBioPortal and TCGA) which together represented a massive cohort of 1.097 tumor and 114 non-tumor breast samples. Our findings showed that samples from patients with breast tumors had low levels of CYSLTR1 mRNA compared to adjacent healthy tissues (p<0.01) (Figure 1A) and, in a larger cohort, we observed the same profile (p = 0.01) (Supplementary Figure 1A). In addition, significant associations were observed between differential expression of CYSLTR1 with patient age (p = 0.01), histological subtype (p<0.0001), TP53 mutational status (p<0.0001), ER status (p<0.0001), PR status (p<0.0001), and molecular subtype (p<0.0001) (Table 1).

Using the TCGA dataset, we performed analyses to identify the expression profile of CYSLTR1 according to clinicopathological parameters. Initially, we observed significant differences between the histological subtypes, where in invasive ductal carcinoma tumors showed low expression of CYSLTR1 when compared to invasive lobular carcinoma (p<0.0001) (Figure 1B). We also describe expression patterns in accordance with the PAM50 subtype classification, where the basal-like type exhibited a decreased transcriptional distribution of CYSLTR1 compared to the other subtypes (Basal-like vs HER2, p<0.0001; Basal-like vs Luminal A, p<0.0001; Basal-like vs Luminal B, p<0.0001; Basal-like vs Normal-like, p<0.0001) (Figure 1C). In addition, patients whose tumors were negative for hormone receptors (ER and PR) and HER2 had lower levels of CYSLTR1 compared to those with positive expression for these receptors (ER+ vs ER-, p<0.0001; PR+ vs PR-, p<0.0001; HER2- vs HER2+, p = 0.02) (Figures 1D-F).

In order to confirm our findings in a larger cohort, we performed an analysis on bc-GenExMiner database. Consequently, a similar profile was observed in terms of expression levels according to the PAM50 classification (p<0.0001) (Supplementary Figure 1B), as well as estrogen (p<0.0001) (Supplementary Figure 2A), progesterone (p<0.0001) (Supplementary Figure 2B), and HER2+ receptors (p<0.0001) (Supplementary Figure 2C). Moreover, within the TCGA cohort, we assessed the possible differences between Basal-like and non-Basal-like (p<0.0001) (Supplementary Figure 2D), TNBC and non-TNBC (p<0.0001) (Supplementary Figure 2E), Basal-like and TNBC vs non- Basal-like and non-TNBC (p<0.0001) (Supplementary Figure 2F); thus, we confirmed that CYSLTR1 transcription levels were downregulated in samples with negative expression for hormone receptors.

Interestingly, the Basal-like immune-suppressed (BLIS) samples showed lower CYSLTR1 expression compared to Basal-like immune-activated samples (p = 0.004) (Figure 2B), while luminal androgen receptor (LAR) (p = 0.001) (Figure 2B) and mesenchymal (MES) samples (p<0.0001) (Figure 2B) exhibited higher levels when compared to BLIS. Similarly, Basal-like 1 and Basal-like 2 triple-negative tumors levels were lower in the cohort from TCGA/UALCAN (p<0.01 and p<0.001, respectively) (Supplementary Figure 3).

Low CYSLTR1 Expression was Associated with Worse Prognosis

Reduced CYSLTR1 mRNA expression levels were significantly correlated with unfavorable prognosis for both OS (Figure 3A; Supplementary Figure 4A-C) and RFS (Figure 3C; Supplementary Figure 4D-F) for all intrinsic BCa subtypes, but especially in basal subtype (Figure 3B, D).

We employed the GSE96058 dataset to execute univariate and multivariate regression analyses. The low expression of CYSLTR1 was an independent factor associated with lower OS in women with Bca (HR = 1.40, p = 0.002) (Table 2). Tumor size, lymph node status, and age were also related to high risk of the disease.

Gene-interaction and Enrichment Analyses

A list of correlated genes (cut-off ≤ or ≥0.7) within the basal subtype is available in Supplementary Table 1. Among the CYLSTR1 co-expressed genes, it was mainly observed that P2RY10 and XCR1 proteins interact directly with CYSLTR1 (Figure 4A-B). The Gene Ontology enrichment analyses demonstrated that several genes co-expressed with CYLSTR1 are involved in the immune system response and immune cell processing and activation (Figure 4C).

Identification and Prognostic Value of Predicted Genes and MicroRNAs

According to KM plotter repository, P2RY10 and XCR1 were assessed for RFS and OS of transcripts. P2RY10 demonstrated a lower but significant expression associated to poor outcome in all subtypes, as well as in the basal subtype [Supplementary Figure 5A-B (RFS) and 5C-D (OS), respectively]. Moreover, by using the mirTarBase repository, three microRNAs: has-miR-335-5p, has-miR-3130-3p, and has-miR-3607-3p were identified as potential regulatory elements of CYSLTR1, P2RY10 and XCR1, respectively. However, OS in the same TNBC patients according to miRNAs (Figure 5) and transcript expression levels (Supplementary Figure 6) were not significantly associated.

Modulation of CYSLTR1 Expression

Through the CTD, we obtained a list of six drugs capable of interacting with CYSLTR1. Leukotrienes C4, D4, and E4 can bind to CYSLTR1 and increase its activity (Figure 6). Other effects might occur depending on the drug in use; for example, the administration of leukotrienes C4 and E4 results in an abundance of calcium, while D4 increases the expression of widely studied proteins such as interleukin (IL)-6, tumor necrosis factor (TNF), and CXCL8. However, Montelukast, Pobilukast, and Zafirlukast (composts of leukotrienes receptors antagonist, LTRAs) induced a reduction in CYSLTR1 protein activity.

Discussion and Conclusion

CYSLTR1 belongs to the cysteinyl leukotriene synthesis pathway and codes for a transmembrane protein receptor that when coupling with many ligands, triggers inflammation-related signaling which leads to a determined phenotype or disease state (15, 29, 30). Yet, there is no consistent evidence to link CYSLTR1 with underlying Bca pathogenesis or even with prognostic values in individuals with aggressive breast tumors.

According to our findings, breast tumor tissue samples showed reduced levels of CYSLTR1 compared to healthy tissue samples. To date, little is known about the profile of CYSLTR1 transcripts in BCa. A study performed by Wang et al. (16) using the RT-qPCR technique, showed that CYSLTR1 was significantly upregulated in tumor samples (n = 90) vs. paraneoplastic breast tissues (n = 30) (16). However, we have to be careful when comparing our findings to this data due to the sample type and size, and approach utilized. Additionally, we observed that a decrease in CYSLTR1 transcripts leads to an unfavorable survival outcome in patients with TNBC tumors, being the first study that evaluated two different datasets and with a relevant sample size. However, it is necessary to reinforce the idea of working with different tumor stages and treatment cohorts to better understand this possible relationship of CYSLTR1 as a potential biomarker in cancer.

Our univariate prognostic analysis according to the Cox proportional hazards regression model confirmed the results observed in the Kaplan–Meier curves as a function of the differential expression of CYSLTR1. Furthermore, the high CYSLTR1 expression group remained an independent prognostic factor in relation to the risk of cancer-specific death, when adjusted for age, tumor size, and lymph node involvement. Our results differ from the study by Magnusson et al. (15) who did not observe a statistical association between the differential expression of CYSLTR1 and the prognosis of patients with BCa. We have to emphasize that the study evaluated the immunoreactivity of the CYSLTR1 protein and had a small set of samples (n = 139).

Considering that CYLSTR1 gene expression in the basal subtype was significantly decreased when compared to the other subtypes, we were led to carry out an in-depth investigation into this clinically more aggressive molecular subtype of BCa. Therefore, we evaluated the expression patterns of CYSLTR1 in the four stable TNBC subtypes, characterized by the expression of distinct molecular profiles that present different prognoses, proposed through studies by the Burstein and Lehman groups (26, 31). Our results showed that CYSLTR1 is consistently expressed in the MES subgroup. Here, we hypothesize that it is possible that CYSLTR1 is more actively involved in epithelial mesenchymal transition and angiogenesis, than in the processes of tumor differentiation and immune activity. This may be supported by the role of cysteinyl leukotrienes (cys-LTs) since they are pro-inflammatory mediators that modulate vascular leakage, permeability and microvasculature response via other leukotriene molecules (32, 33, 34). Furthermore, CYSLTR1 transcripts were also expressed to a greater extent in tumors of the LAR subgroup, which is highly expressed on the nuclear androgen receptor receptor. Consequently, as CYSLTR1 is linked to kinase activity (35, 36), this could lead us to hypothesize that CYSLTR1 could participate in blocking androgen-dependent signaling and PIK3. To date, no work has focused on studying its possible role in tumorigenesis in TNBC cases. Thus, our work brings unprecedented data about this aggressive type of BCa.

Our analyses of signaling and enrichment pathways for CYSLTR1 have indicated some immunological mechanisms related to inflammation in which toll-like family genes and cytokines may participate. It is noteworthy that two protein-coding genes: P2RY10 and XCR1 exhibited a positive correlation with CYSLTR1 in a TNBC dataset. Eosinophils can be found in the TMev, as they secrete different types of leukotrienes as part of the induction of inflammatory processes (10). Furthermore, they also generate significant amounts of platelet activating factor and promote the production of characteristic cytokines such as TNFα and IL-5 (37, 38). P2RY10 is a G-couple protein receptor that participates in the inflammatory response, stimulated by many molecules such as chemokines, lysophospholipids and prostanoids. Its biological role has not been fully elucidated, but it may participate in eosinophil maturation and eosinophilopoiesis in vitro (38). On the other hand, Yang et al. (39) suggest that XCR1 may act as a progression factor in ER-responsive Bca cell lines through the MAPK/ERK and PI3K/AKT/mTOR pathways that promote migration and invasion by significantly decreasing the protein level of β-catenin (40). Regarding the possible prognostic role, patients with TNBC cases who had high gene expression of XCR1 and P2RY10 exhibited a trend towards greater survival, further confirmed by an independent dataset.

As for epigenetic mechanisms of gene regulation, has-miR-355-5p showed a certain tendency to downregulate CYSLTR1 expression in a TNBC setting. To date, there is no evidence that describes consistent associations between these genes and their Bca-targeted miRNAs. Therefore, we suggest IHC studies to unravel mechanisms underlying survival and immunological processes in TNBC.

Regarding possible drug interactions, we observed that Montelukast, Zafirlukast and Pobilucast played a role in reducing CYSLTR1 expression levels in our in silico experiments. Based on the above, Suknuntha et al.'s (30) group observed that MDA-MB-231 BCa cells, when treated with Montelukast and Zafirlukast molecules, can inhibit cell proliferation and apoptosis, but only Zafirlukast can induce cell cycle arrest. On the other hand, leukotrienes appear as possible positive modulators of CYSLTR1 expression. Both strategies are promising and need to be carefully investigated.

Study Limitations

Some limitations of the analysis performed here must be acknowledged. First, we employed different expression analysis methods compared to other studies examining CYSLTR1 in Bca. Second, as seen in in silico analyses, RNA-seq based expression data were not complemented with protein data to corroborate our findings. Third, many studies available in public databases were deficient in clinicopathological information, and the most important, TNBC studies only accounted for up to 15% of the Bca population, so it is difficult to reach significant conclusions. Nonetheless, based on our findings, we can provide insights into the possible role of CYSLTR1 in BCa disease survival, particularly in TNBC cases.

Our study showed that CYSLTR1 is transcriptionally less expressed in breast tumors compared to adjacent tissue. Additionally, among the tumor subtypes, TNBC had lower levels of CYSLTR1. Low CYSLTR1 expression was associated with worse survival in BCa patients and especially in TNBC. CYSLTR1 is co-expressed with genes that participate in the adaptive immune response and lymphocyte activation. Finally, we suggest that CYSLTR1 may not be working alone, but with linked proteins and miRNAs that could serve as new possible targets for other therapies in Bca, especially in TNBC subtypes.

Ethics Committee Approval: Not necessary.

Informed Consent: Not necessary.

Peer-review: Externall and internally peer-reviewed.

Authorship Contributions

Concept: M.P.F.C., D.R.d.B.; Design: M.P.F.C., D.R.d.B.; Data Collection or Processing: A.G.C., M.P.F.C., D.R.d.B.; Analysis or Interpretation: A.G.C., M.P.F.C., D.R.d.B.; Literature Search: A.G.C., M.P.F.C., D.R.d.B.; Writing: A.G.C., M.P.F.C., D.R.d.B., G.Á.d.G., J.M.R.S.L., R.G.d.N., M.T.F., R.M.L.

Conflict of Interest: No conflict of interest was declared by the authors.

Financial Disclosure: The authors declared that this study received no financial support.