The regulatory landscape of retinoblastoma: a pathway analysis perspective
Computational Genomics Division, National Institute of Genomic Medicine (INMEGEN), Mexico City, Mexico
Medical Research Unit in Infectious Diseases, Hospital de Pediatría, CMN SXXI, Instituto Mexicano del Seguro Social, Mexico City, Mexico
Pharmacology Department, CINVESTAV, Mexico City, Mexico
Center for Complexity Sciences, National Autonomous University of Mexico (UNAM), Mexico City, Mexico
Abstract
Retinoblastoma (Rb) is a rare intraocular tumour in early childhood, with an approximate incidence of 1 in 18 000 live births. Experimental studies for Rb are complex due to the challenges associated with obtaining a normal retina to contrast with diseased tissue. In this work, we reanalyse a dataset that contains normal retina samples. We identified the individual genes whose expression is different in Rb in contrast with normal tissue, determined the pathways whose global expression pattern is more distant from the global expression observed in normal tissue, and finally, we identified which transcription factors regulate the highest number of differentially expressed genes (DEGs) and proposed as transcriptional master regulators (TMRs). The enrichment of DEGs in the phototransduction and retrograde endocannabinoid signalling pathways could be associated with abnormal behaviour of the processes leading to cellular differentiation and cellular proliferation. On the other hand, the TMRs nuclear receptor subfamily 5 group A member 2 and hepatocyte nuclear factor 4 gamma are involved in hepatocyte differentiation. Therefore, the enrichment of aberrant expression in these transcription factors could suggest an abnormal retina development that could be involved in Rb origin and progression.
Untitled section
Keywords: retinoblastoma, functional enrichment analysis, master regulator analysis, expression microarrays, transcriptional regulation
Article notes
Untitled section
Received 2022 Jan 10; Accepted 2022 Apr 13; Collection date 2022 May.
1.Background
Retinoblastoma (Rb) is a rare intraocular malignant tumour from early childhood that originates in the retina with an incidence of 1 in 15 000 to 1 in 18 000 live births [1]. However, experimental studies for Rb are difficult due to the challenges associated with obtaining normal retina, the tissue on which this tumour originates, in addition to the rarity of the tumour and the age of the patients. Rb is considered the most robust clinical model of genetic predisposition to develop cancer. Epidemiological and molecular studies of the tumour led to the discovery of the first tumour suppressor gene, the Rb susceptibility gene RB1 [2–5]. This tumour is a cornerstone of cancer research and has been widely studied with different molecular tools. It has been suggested that alterations in specific pathways across different genes result in cancer, so cancer has been proposed as a disease of biochemical pathways [6,7]. Identifying and quantifying altered pathways are essential for discovering possible therapeutic targets and understanding disease development [8–10]. Pathway analysis (PA) has been used to obtain functional insights with omic studies, especially in cancer-related research [11–14]. PA aims to analyse the data obtained via high-throughput technologies to detect relevant sets of genes altered in specific conditions. These approaches rely on the high-level organization of biological functions stored in structured databases (such as Kyoto Encyclopedia of Genes and Genomes (KEGG)) and the implementation of sophisticated statistical methods. Overrepresentation analysis (ORA), functional class scoring (FCS) and pathway topology-based analysis (PTB) are different approaches for PA [15,16].
ORA methods can detect statistically enriched biological processes from a list of differentially expressed genes (DEGs). Statistical significance will be reached if more DEGs belong to a specific pathway than expected by chance. By contrast, FCS methods overcome the limitation of focusing on a limited number of genes per pathway; these methods account for the expression level of all genes in a pathway or gene set to compute a metric that reflects the state of each pathway or gene set. Pathifier [10] is an FCS method that uses a principal curve approach to calculate the deregulation level for each pathway in each sample by measuring the sample's deviation from normal behaviour. This algorithm does not require detailed knowledge about the tested pathways or their interaction networks. This approach can be used to perform a functional mapping for every selected pathway on every given individual, an essential step towards personalized diagnostics and therapeutics. Additionally, in PTB methods, the pathway topology is integrated into the analysis.
A different approach, master regulators analysis (MRA), has been extensively used for the study of cancer [11–14,17–21]. MRA allows for the determination of key transcriptional regulators that could be responsible for the establishment of cellular and organismal phenotypes. It has been proposed that such master regulators control significant interactions, particularly those in regulatory bottlenecks [22]. MRA assumes that the more DEGs in tumour cells compared to healthy tissue are responsible for maintaining the tumour phenotype. Therefore, the analysis seeks transcription factors (TFs) that target DEGs. A transcriptional regulatory network is generated based on the statistical dependency between each TF expression pattern and the rest of the genes. The regulons per TF (the set of targets of a TF) are assembled, and the DEG content per regulon is assessed through gene set enrichment analysis (GSEA) [23]. The TF with a regulon (genes that are either induced or repressed directly or indirectly by the TF) with the highest normalized enrichment score (NES) will become a candidate transcriptional master regulator (TMR) [12,18]. Identifying the TMRs may facilitate the design of targeted interventions, another desirable instance towards personalized medicine.
In this work, we revisited the dataset generated and used by Kapatai et al. [24]. They were able to obtain samples from three normal retinas, and we think this dataset represents an invaluable opportunity to study Rb from a pathway-level perspective. We have done a pathway-level study in Rb. We obtained and analysed the DEGs from comparing the gene expression between Rb samples and the normal retina. Based on the DEGs, we enriched functional processes with an ORA analysis. We also calculated the deregulation score for each pathway with an FCS analysis by applying the Pathifier algorithm on the same dataset. By either identifying the pathways with a higher proportion of DEGs than expected by chance or identifying the pathways with a higher global expression difference compared to normal samples (higher deregulation score), we simplified the complexity of the model to a handful of relevant pathways. This dimensionality reduction makes the interpretation of the relevant pathways easier to manage, as they are thought to be directly associated with biological functions related to the Rb phenotype.
Furthermore, we also identified TF regulators performing an MRA, which we integrated into a pathway-level analysis by interrogating which TMRs control these pathways. We applied all these analyses to the Kapatai et al. dataset since the number of tumour samples and the presence of control samples offer a rare opportunity to study Rb with pathway-focused tools. In this study, we present a broader view of the transcriptomic context associated with the Rb phenotype by evaluating individual gene contributions, collective gene expression patterns via master regulators, and identifying functional consequences by deregulated pathways, thus presenting a clear landscape of the effects of gene deregulation underlying this pathologic phenotype.
2.Methods
2.1.Datasets
The Kapatai et al. [24] dataset with 20 samples of Rb was generated at Birmingham University and kindly provided by the authors. Normal tissue samples were obtained from the GEO database (https://www.ncbi.nlm.nih.gov/geo/): two normal adult retina (GSM607947, GSM607948) and one 96-day human fetal retina (GSM460264). All CEL files were built up through Affymetrix Human Gene 1.0 ST microarrays. Sample identifications match the original publication [24].
2.2.Data analysis
The analysis pipeline was developed in three levels: differential expression, deregulation pathways analysis and master regulation analysis. (figure 1).
2.3.Pre-processing
The CEL files were pre-processed with the Affymetrix Expression Console software. The robust multichip average (RMA) algorithm was used to normalize and summarize probesets [25]. Pre-processed expression data were annotated with annotation files available in Affymetrix. The processing pipeline was performed in R Bioconductor [26,27].
2.4.Differential expression analysis and overrepresentation analysis
We calculated differential expression via the limma Bioconductor R package [28], contrasting normal retina versus Rb. DEGs (log-fold change (LFC) > 2 and adjusted p-value < 0.05) were grouped by hierarchical clustering (HC) (using Euclidean distance as distance method and complete agglomeration as clustering method) and the results plotted in a heatmap with gplots from R package [29]. Overrepresentation analysis was done using the Webgestalt tool [30]. We chose KEGG pathways like enrichment categories, Affymetrix Human Gene 1.0 ST IDs as the reference list and Benjamini and Hochberg as the multiple-testing correction method. We selected pathways with a false discovery rate (FDR) < 0.05 as statistically significant.
2.5.Functional class scoring analysis: pathway deregulation analysis
To determine the pathway deregulation score (PDS), we used the Pathifier algorithm [10], which uses the normalized gene expression matrix. Briefly, Pathifier applies the Hastie and Stueltzle algorithm [31] in a multi-dimensional space. Each axis represents a gene to calculate a principal curve that passes through and reaches all the points. Each point represents each gene's expression values that belong to one specific pathway in one sample. The centroid of the control group is the origin of the curve, and the distance of each point along the curve is the PDS for that sample–molecular pathway combination [10]. This process is reiterated, and a PDS is calculated for each sample and for each pathway considered. We used KEGG as a reference database of molecular pathways (296 metabolic pathways) [32]. The advantage of the KEGG database is that it is large, carefully curated and has integrated genomic and chemical information databases for biological processes [33]. Results for each sample and each pathway were plotted in a heatmap using Euclidean distance as distance metric, and Ward2 as the clustering method using gplots R Bioconductor package [29].
2.5.1.Pathway deregulation score z-score
The PDSs obtained with Pathifier are based on comparing all samples through each single pathway. As different pathways may have different dynamic ranges, each PDS was z-score transformed to allow a fair comparison between different pathways. The mean and standard deviation of the PDSs per pathway were obtained and used to calculate the z-score per pathway per sample. This transformation allowed us to compare deregulation values across all pathways. Next, we obtained the median PDS z-score per pathway, taking into consideration either all tumoral samples or all samples per subgroup obtained with clustering analysis. From these, we can obtain the top 10 deregulated pathways in Rb or in each tumoral subgroup, respectively.
2.6.Master regulators analysis
We used corto (Correlation Tool) implemented in R [34] to generate a list of correlated TFs and target genes. For this analysis, we obtained the list of human TFs from the TFCheckpoint database [35]. We kept the TFs with experimental evidence that were found in the Human Gene 1.0ST microarray. Corto was run with 100 bootstraps, and interactions with a p-value less than 1 × 10−7 were chosen as true interactions [12]. This algorithm generates a regulatory network based on the correlation of expression, and a regulon is composed by all the target genes of a specific TF.
We used MRA, implemented in the package corto, to find the TMRs. MRA requires a regulatory network, a molecular signature and a null model. The network obtained by corto was used as the regulatory network. A t-statistic per gene was obtained from comparing the healthy versus control expression values, and it was transformed into z-scores. These z-scores were considered the molecular signature. Then, the gene labels were reshuffled, and this randomization was repeated 1000 times to obtain the null model. During the MRA, GSEA was run for each regulon using the z-scores as weights, and a NES was obtained. Each regulon NES was then compared with its null model distribution to calculate a p-value. The regulons with the highest NES value, a statistically significant p-value and with at least 20 genes were selected as TMRs [12,18]. Finally, we investigated the pathways regulated by each TMR.
3.Results
3.1.Differentially expressed genes are related to vision
Differential expression analysis allows us to measure quantitative changes in expression levels between two conditions. We found 305 DEGs (absolute LFC > 2) by comparing Rb samples to the normal retina, including 186 under-expressed genes and 119 over-expressed genes in Rb compared with the normal retina (see electronic supplementary material, table S1). Table 1 shows the first 10 genes most differentially expressed between both tissues, all of them with an absolute LFC ≥ 4.5 of them over-expressed (LFC ≥ 4) and 5 of them under-expressed (LFC ≤ −4).
| gene symbol | Ensembl | logFC | name | description |
|---|---|---|---|---|
| SNORD41 | ENSG00000209702 | 5.47 | small nucleolar RNA, C/D Box41 | SNORD41 (small nucleolar RNA, C/D Box 41) is an RNA gene and is affiliated with the snoRNA class.a |
| SNORD115-23 | ENSG00000201331 | −4.65 | small nucleolar RNA, C/D Box115-23 | SNORD115-23 (small nucleolar RNA, C/D Box115-23) is an RNA Gene and is affiliated with the snoRNA class, a class of small non-coding RNAs associated with nucleotide chemical modification, such as methylation and pseudouridylation, by guiding and tethering partner enzymes to specific sites on RNA targets.a |
| SNORA71D | ENSG00000200354 | 4.59 | small nucleolar RNA, H/ACA Box 71D | SNORA71D is an RNA Gene and is affiliated with the snoRNA class, a class of small non-coding RNAs associated with nucleotide chemical modification, such as methylation and pseudouridylation, by guiding and tethering partner enzymes to specific sites on RNA targets.a |
| SNORA23 | ENSG00000201998 | 4.52 | small nucleolar RNA, H/ACA Box 23 | SNORA23 (small nucleolar RNA, H/ACA Box 23) is an RNA Gene and is affiliated with the snoRNA class.a Has been found upregulated in human pancreatic ductal adenocarcinoma [36]. |
| SNORA3A | ENSG00000200983 | 4.46 | small nucleolar RNA, H/ACA Box 3A | SNORA3A (small nucleolar RNA, H/ACA Box 3A) is an RNA gene and is affiliated with the snoRNA class.a |
| CRABP1 | ENSG00000166426 | −4.03 | cellular retinoic acid binding protein 1 | this gene encodes a specific binding protein for a vitamin A family member and is thought to play an important role in retinoic acid-mediated differentiation and proliferation processes. It is structurally similar to the cellular retinol-binding proteins but binds only retinoic acid at specific sites within the nucleus, contributing to vitamin A-directed differentiation in epithelial tissue.a CRABP1 has found that its under-expression in some other cancers such as breast cancer and ovarian adenocarcinoma is related to poor prognosis [37,38]. |
| FABP7 | ENSG00000164434 | −4.30 | fatty acid binding protein 7 | the gene encodes a small, highly conserved cytoplasmic protein that binds long-chain fatty acids and other hydrophobic ligands. The encoded protein is essential in establishing the radial glial fibre in the developing brain. Alternative splicing and promoter usage result in multiple transcript variants encoding different isoforms.a |
| DPP4 | ENSG00000197635 | 4.15 | dipeptidyl peptidase 4 | the DPP4 gene encodes dipeptidyl peptidase 4, which is identical to adenosine deaminase complexing protein-2, and the T-cell activation antigen CD26. It is an intrinsic type II transmembrane glycoprotein and a serine exopeptidase that cleaves X-proline dipeptides from the N-terminus of polypeptides. Dipeptidyl peptidase 4 is highly involved in glucose and insulin metabolism, as well as in immune regulation.a |
| SNORD115-25 | ENSG00000199489 | −4.03 | small nucleolar RNA, C/D Box115-25 | SNORD115-25 (small nucleolar RNA, C/D Box115-25) is an RNA Gene affiliated with the snoRNA class.a |
| NEFL | ENSG00000277586 | −4.00 | neurofilament light | neurofilaments are type IV intermediate filament heteropolymers composed of light, medium, and heavy chains. Neurofilaments comprise the exoskeleton, and they functionally maintain the neuronal calibre. They may also play a role in intracellular transport to axons and dendrites. This gene encodes the light chain neurofilament protein.a |
Unsupervised HC of DEGs showed four large clusters. As expected, this result exactly reproduces the grouping found by Kapatai et al. (figure 2). Cluster 1 is formed by 13 samples and corresponds to subgroup 1 described by Kapatai et al.; cluster 2 is composed of 5 samples and matches subgroup 2, also described by Kapatai et al. Interestingly, two tumoral samples are closely grouped with the control samples cluster (normal-like).
3.3.Pathway analysis reveals two retinoblastoma subgroups
We used the Pathifier algorithm [10] to determine pathway deregulation in Rb. Pathifier is an algorithm that combines the expression level of all genes within a single pathway or biological process, generating a unique deregulation value for each sample. This value is called the PDS. A PDS is calculated for each sample, and for each pathway, it reflects a deregulation value from 0 to 1. A value close to zero indicates no dysregulation with respect to normal tissue, and a value of 1 is assigned to the most deregulated sample or samples.
For each studied pathway, the tumoral samples with gene expression values similar to the average of control healthy tissue will have a PDS close to 0. By contrast, the samples with the highest PDS (close to 1) will have the most significant differences in gene expression levels compared to the average normal tissue. Previous studies have used this methodology to determine pathway deregulation of metabolic pathways in multiple cancer studies [45–49].
We estimated a PDS for each KEGG pathway (N = 296) for every sample. Unsupervised HC shows that deregulation levels are consistent with Rb subgroups. Figure 3a Illustrates that samples are grouped according to the same subgroups obtained previously via differential gene expression analysis (dendrogram at top subgroup 1, subgroup 2 and normal-like). However, the global clustering presents clear differences, grouping all Rb samples in a higher level cluster. On the other hand, pathways are clustered into two high-level clusters: a larger cluster composed of 246 pathways (cluster (i)) and a small cluster composed of 50 pathways (cluster (ii)). We can observe that pathways in cluster (i) tend to be more deregulated in subgroup 1 compared to subgroup 2 and normal-like samples. On the contrary, pathways in cluster (ii) tend to be more deregulated in subgroup 2 compared to subgroup 1 and normal-like samples. Also, we can see that most pathways tend to be less deregulated in normal-like samples compared to subgroup 1 and subgroup 2 samples; this could be caused by contamination with normal retina cells.
The most deregulated pathways in Rb were selected from the median PDS z-score calculation. The top 10 most deregulated pathways (i.e. with the highest median PDS z-score) were steroid hormone biosynthesis, focal adhesion, 2-oxocarboxylic acid metabolism, fatty acid degradation, Rap1 signalling pathway, peroxisome proliferator-activated receptor (PPAR) signalling pathway, morphine addiction, retrograde endocannabinoid signalling, amyotrophic lateral sclerosis (ALS) pathway and arrhythmogenic right ventricular cardiomyopathy (ARVC). These pathways were selected and grouped by unsupervised HC (figure 3b). This clustering shows that the normal-like samples are grouped with the normal retina and the two Rb subgroups are clustered together. Also, we can observe the same two high-level clusters of pathways as in the global deregulation heatmap: a top cluster that is more deregulated in subgroup 1 compared to subgroup 2 and the bottom one which seems to be equally deregulated in both subgroups (figure 3b).
We also determined which are the most deregulated pathways per Rb subgroup. Interestingly, the most deregulated pathways vary across subgroups (electronic supplementary material, table S2).
3.4.Master regulators analysis finds transcriptional master regulators involved in hepatocyte differentiation
We did an MRA to search possible TMRs in Rb. We used the corto algorithm to infer the regulatory network; this algorithm uses Spearman correlation supplemented with the data processing inequality between a list of previously defined TFs and genes to find statistical correlation from expression data [34,45,46]. This algorithm assumes that the most important regulators will be the ones affecting the expression of the highest number of their target genes in the condition of interest.
Figure 4 shows the first 10 TMRs predicted by the MRA and their complete regulons. It can be noticed how the regulons vary widely in size and how there are few cases of coregulation of target genes by more than one TMR, with the notable exception of OR5AP2 and SLC22A10 that are regulated by three different TMRs (pituitary-specific positive transcription factor 1 (POU1F1), E74-like factor (ELF5) and hepatocyte nuclear factor 4 gamma (HNF4G)). The most statistically significant TMRs are presented in table 2 with a short description of their function.
| name | NES | p-value | LFC | function |
|---|---|---|---|---|
| NR5A2 | 6.43 | 1.25 × 10−10 | 1.06* | also known as liver receptor homolog-1 (LRH1), DNA binding zinc finger transcription factor, member of the fushi-tarazu factor 1 |
| HNF4G | 6.17 | 6.68 × 10−10 | 1.03* | a critical factor for hepatocyte differentiation |
| POU1F1 | 5.97 | 2.34 × 10−9 | 1.21* | a transcription factor for a growth hormone that regulates mammalian development |
| GRHL2 | 5.60 | 2.14 × 10−8 | 0.49* | it plays an important role in primary neurulation and in epithelial development [47,48] |
| MITF | 5.58 | 2.45 × 10−8 | 0.36 | it regulates melanocyte and optic cup derived retinal pigment epithelium development |
| HOXA2 | 5.34 | 9.09 × 10−8 | 1.68* | it may regulate gene expression, morphogenesis, and differentiation during embryonic development |
| KLF4 | 5.18 | 2.20 × 10−7 | 0.09 | member of the SP1 family is required for normal development of skin barrier function. Remarkably this is one of the Yamanaka factors involved in the signalling network of pluripotency [49] |
| AR | 5.00 | 5.80 × 10−7 | 0.16 | it functions as a steroid hormone-activated transcription factor. Upon binding the hormone ligand, the receptor dissociates from accessory proteins, translocates into the nucleus, dimerizes, and then stimulates transcription of androgen responsive genes |
| ELF5 | 4.97 | 6.80 × 10−7 | 0.91* | it is a member of the ETS TF family, unique to animals, implicated in early development and cancer progression. It regulates epithelium specific genes and the later stages of terminal differentiation of keratinocytes |
| NR0B1 | 4.79 | 1.67 × 10−6 | 0.59* | it is a transcription factor that acts as a dominant-negative regulator of transcription which is mediated by the retinoic acid receptor |
4.Discussion
In this work, we performed the functional analysis of Rb on three levels: differential gene expression, pathway deregulation and phenotype transcriptional regulation. Differential expression was carried out to characterize the gene level in contrast with the normal retina. We used Pathifier to analyse the biological functions and metabolic pathways level and MRA for the global transcriptional regulation level. In this approach, as a microscope focuses on different depths, we focus the analysis on different biological layers, and we find evidence of alterations at each level revealing the mechanisms behind Rb.
As previously mentioned, the main cluster features revealed by such analyses (figure 3a) are in agreement with previous findings by the group of Kapatai [24]. (i) Differential gene expression analysis and clustering resulted in two well-defined tumour subtypes plus a normal-like subset of samples that clustered together with the normal tissue. (ii) Overrepresentation analysis of differential genes using the KEGG database pointed out biological processes characteristic of retinal function with several downregulated genes and cancer progression related signalling. (iii) Pathway deregulation analysis uncovered both known and novel highly deregulated metabolic and signalling pathways. Interestingly, clustering analysis of the PDS by sample reveals two distinctive clusters with different sets of globally deregulated pathways. These two clusters correspond to a very high degree with the tumour subtypes first identified by Kapatai and co-workers and later confirmed by our differential expression analysis. This is relevant since this strongly suggests that the differences in gene expression profiles and deregulated pathways in both tumour subtypes are an actual consequence of distinctive regulatory programmes. (iv) In connection to this latter point, TMR analysis allowed us to find gene regulatory switches that may be behind such different transcriptional regulation programmes and, ultimately, perhaps behind phenotypic differences in the tumour subtypes.
Most of the DEGs found at this work (6 out of 10) belong to the snoRNA class, a small non-coding RNA class. One well-studied function of snoRNAs is their role in the modification, maturation and stabilization of rRNA. [50]. Emerging evidence has demonstrated significant roles of snoRNAs in cancer [51–53]. Mutations and aberrant expression of snoRNAs have been reported in cell transformation, tumorigenesis and metastasis, indicating that snoRNAs may serve as biomarkers and/or therapeutic targets of cancer [53]. Increasing evidence suggests that snoRNAs are closely associated with the TP53 regulatory pathway. For instance, a recent study identified that snoRNAs and fibrillarin (an enzymatic small nucleolar ribonucleoprotein) were usually over-expressed in human breast and prostate cancers, and this over-expression promoted tumorigenicity in vitro and in vivo [54].
Most top deregulated pathways have been previously associated with cancer. Steroid hormone biosynthesis is a pathway of recognized importance in the establishment and proliferation of various types of cancer [55–58]. Vertebrate steroid nuclear receptors (NRs) include the androgen receptor (AR), mineralocorticoid receptor, progesterone receptor, glucocorticoid receptor and estrogen receptor alpha [59]. There is very strong evidence of steroid hormone signalling acting as cancer drivers, at least in breast and prostate cancer and lymphoma [60–62]. Steroid hormones formed in a particular cell can activate the corresponding intracellular receptors of the NR superfamily, which act as TFs [56]. There is evidence that NR interactions with coactivators and corepressors are distorted in cancer, which ultimately disrupts the NR function [63–65].
On the other hand, focal adhesion plays an essential role in processes such as cell motility, cell proliferation, cell differentiation, regulation of gene expression and cell survival [66]. Closely related to focal adhesion is Rap1 signalling pathway. Ras-associated protein-1 (Rap1), a small GTPase in the Ras-related protein family, is an important regulator of essential cellular functions like formation and control of cell adhesions and junctions, cellular migration and polarization. Rap1 also regulates MAP kinase (MAPK) activity, playing many roles during cell invasion and metastasis in different cancers [67]. A common characteristic among tumour cells is their ability to reprogramme their metabolism [7,68–71]. The 2-oxocarboxylic acid metabolism pathway plays an important role in metabolic regulation [72–74]. Specifically, in cancer, 2-oxocarboxylic acid metabolism pathway seems to help cell proliferation indirectly through the mitochondrial electron transport chain via the metabolism of NAD+ [72,74]. The fatty acid degradation pathway is fundamental through the role that fatty acids play in cellular metabolism. Fatty acids are essential as structural components of the membrane matrix, secondary messengers and can also serve as fuel sources for energy in normal tissue. Lipid metabolism plays a fundamental role in processes such as metastasis, invasion or even tumorigenesis [75,76]. Moreover, particularly in cancer fatty acid, degradation could be related to the generation of substrates necessary for nucleotide biosynthesis rather than to energy generation [77].
PPARs are nuclear hormone receptors that are activated by fatty acids and their derivatives. PPAR protein family plays a role in processes like the clearance of circulating or cellular lipids via the regulation of gene expression involved in lipid metabolism in liver and skeletal muscle, lipid oxidation and cell proliferation [32]. It has been observed that the PPAR pathway can favour tumour progression through mechanisms such as induction of cell proliferation, inhibition of apoptosis, upregulation of VEGF expression, increased PGE2 production and COX expression [78], the latter responsible for the biosynthesis of prostanoids involved in mitosis and the inflammatory response [79].
Morphine addiction and retrograde endocannabinoid (ECS) signalling pathways are also closely related to carcinogenic processes [44]. Gamma-aminobutyric acid (GABA) is part of the morphine addiction pathway, and it is an important neurotransmitter for retinal neurons. It is involved in retinal maturation and development [80]. It has been shown that most Rb tumours are able to generate endogenous GABA via one or several different biosynthetic pathways [81]. Importantly, GABA is related to cancer cell metastasis and proliferation [82]. Retrograde ECS pathway is the pathway that intervenes in the regulation of endogenous cannabinoids that act as retrograde synaptic messengers in the brain. There is increasing evidence that the endocannabinoid system may be an important regulator of tumour cell malignancy [44,83,84]. However, there is contradictory evidence; it can stimulate tumour cell proliferation, angiogenesis and immunosuppression or inhibit tumour cell growth, angiogenesis, invasion and metastasis [84]. Noteworthy, some of these molecules have been recently proposed as possible therapeutic targets [44,83] since there is evidence that some non-psychoactive cannabinoids can have antimetastatic, anti-invasive effects [85,86] and antiproliferative capabilities [87] by inducing apoptosis [88] and autophagy [89,90].
The ALS pathway, which is a progressive and lethal motor disorder, includes in its gene list TP53, several members of the caspase family pathway, TNF and MAP kinases; all of these genes are widely documented as important agents in many cancer phenotypes [91–94]. Furthermore, a 2016 preclinical study showed how the C9orf72 knockout, the most common cause of ALS, caused an unexpected increase in tumours in mice [95]. On the other hand, more than half of the members of the ARVC pathway are known to be related to cancer, i.e. 20 (out of 77 members, 25%) integrins, 19 calcium voltage-gated channel-related proteins, TCF7 and two TCF7-like proteins [96–98].
Metabolic and signalling pathways previously implicated in Rb development were also found. Fatty acid synthesis is suggested to support an essential functional aspect of Rb. Most tumours present high levels of fatty acid synthase (FAS), and higher FAS expression has been correlated with more advanced choroid and optic nerve invasion, high mitotic index, and less differentiated histology [99].
Throughout the different sections of this work, we have presented the results of three different classes of pathway-level analysis derived from the same initial data: expression data from 29 Rb samples and three normal tumour-free retina samples. In table 3, we show the integration of our results. We present the pathways that are enriched in DEGs, the top 10 most deregulated pathways and the crosstalk between any of them and the regulons of each TMR. In what follows, we will discuss how the integration of these different levels of analysis led us to a unified view of the phenomena, from the molecular to the functional and up to the phenotypic levels.
| source | description | DE | FCS | MRA | |||||||||||
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| ER | FDR | RN | N | AR | ELF5 | GRHL2 | HNF4G | HOXA2 | KLF4 | MITF | NR0B1 | NR5A2 | POU1F1 | ||
| DE | olfactory transduction | 3.27 | 9.45 × 10−6 | 212 | 88 | OR10Q1,OR11H6,OR7C2 | OR10A4,OR10H3,OR10K2 (N = 32) | OR4A16,OR4B1,OR4D10 (N = 9) | OR2D3,OR4A16,OR51G2 (N = 18) | OR10AD1,OR10H2,OR10H3 (N = 13) | NA | OR8D2 | OR10G7,OR1L1,OR7C2 | OR10J1,OR10J5,OR10K1 (N = 10) | OR10G8,OR2F1,OR4F15 (N = 12) |
| DE | phototransduction | 21.11 | 3.30 × 10−10 | 22 | 0 | NA | NA | NA | NA | NA | NA | NA | NA | NA | NA |
| PA/DE | retrograde endocannabinoid signalling | 3.68 | 0.040 | 8 | 1 | NA | MAPK9 | NA | NA | NA | NA | NA | NA | NA | NA |
| PA | amyotrophic lateral sclerosis (ALS) | 5.27 | 0.16 | 4 | 10 | C9orf72,DERL1,DNAH9 | DCTN4,DNAH5,SETX,SPG11 | NA | NA | COX4I1,UQCRC2 | NA | NA | PPP3CA | COX4I1 | NA |
| PA | fatty acid degradation | 1.22 | 1 | 7 | 1 | NA | ADH7 | ADH7 | NA | NA | NA | NA | NA | NA | NA |
| PA | steroid hormone biosynthesis | NA | NA | 2 | 12 | NA | UGT1A1,UGT1A10,UGT1A3 (N = 9) | NA | NA | HSD3B1 | NA | NA | NA | NA | CYP3A7,UGT2B4 |
| PA | 2-oxocarboxylic acid metabolism | NA | NA | 1 | 1 | NA | NA | NA | NA | NA | NA | IDH3B | NA | NA | NA |
| PA | PPAR signalling pathway | 1.47 | 1 | 9 | 0 | NA | NA | NA | NA | NA | NA | NA | NA | NA | NA |
| PA | Rap1 signalling pathway | 0.79 | 1 | 6 | 1 | NA | FGF10 | NA | NA | NA | NA | NA | NA | NA | NA |
| PA | focal adhesion | 1.08 | 1 | 5 | 8 | CHAD,ITGA11 | COL6A6,LAMB4,MAPK9 | NA | NA | NA | THBS2 | NA | MYL2 | LAMC2 | NA |
| PA | morphine addiction | 3.02 | 0.673 | 3 | 0 | NA | NA | NA | NA | NA | NA | NA | NA | NA | NA |
| PA | arrhythmogenic right ventricular cardiomyopathy (ARVC) | 2.24 | 1 | 10 | 2 | ITGA11 | NA | NA | NA | NA | NA | NA | NA | NA | CTNNA3 |
Only one pathway is enriched in DEGs and belongs to the top10 most deregulated pathways: retrograde endocannabinoid signalling. This pathway is also presumably regulated by one TMR: ELF5 regulates MAPK9. ELF5 is a transcription factor of the ETS family. This family is important for its role in cell proliferation, cell differentiation, cell development, apoptosis and tissue remodelling [100]. It has been observed that the deregulation of their members results in malignant cellular transformations [101]. Changes in the expression pattern of the genes in the endocannabinoid signalling pathway have been seen to impact disease progression and patient survival in various types of cancer [44].
The phototransduction pathway contains a lot more DEGs than the number expected by chance (highest enrichment ratio); besides, this pathway is one of the most highly deregulated (ranking 22). This fact is interesting, in particular since genes of the PDE6 family previously associated with the action of photoreceptors as signal transducers of luminal stimuli [102,103] have also recently been argued to be significantly over-expressed in human cancers, in particular in breast carcinomas [36]. PDE6A and PDE6G are present among the DEGs of this pathway. No gene in the phototransduction pathway is regulated by any of the proposed TMRs.
The olfactory transduction pathway has a noticeable behaviour as it is enriched in DEGs. It has a substantial number of genes regulated by the proposed TMRs (88 out of 427 total genes). This pathway is also significantly enriched (FDR < 0.05) for 6 of the 10 TMR regulons (ELF5, grainyhead-like transcription factor 2 (GRHL2), HNF4G, homeobox A2 transcription factor (HOXA2), nuclear receptor subfamily 5 group A member 2 (NR5A2) and POU1F1). This pathway comprises many molecules (448) involved in the transduction of chemical signals into the brain, including receptor proteins, olfactory specific Gs-proteins and type III adenylyl cyclases. Although the structure and function of its components in the nasal cavity are well known [104–106], there are still many unknowns about their role in other tissues [107]. However, they have been described as participating in secretion, migration, apoptosis, differentiation and cell growth [107]. Recent studies suggest the possibility of using some olfactory receptors as biomarkers of some types of cancer, such as OR51E1 and OR51E2 [107–111], which were found to be slightly over-expressed in Rb (LFC = 1.36 and 1.45 respectively, adjusted p-value < 0.001) and present in the HOXA2 and ELF5 regulons, respectively.
Although 2-oxocarboxylic acid metabolism is the most deregulated pathway, it is not enriched in DEGs, and it contains only one gene (IDH3B) regulated by a TMR: melanocyte-inducing transcription factor (MITF). Importantly, this interaction has been experimentally validated [112]. Furthermore, IDH3B is an important factor in the crosstalk between the cell cycle and the Krebs cycle [113]. Additionally, IDH3B has been shown to promote the transition between the G1 and S phase stimulating cell proliferation, and its over-expression has been related to poor survival [113]. Moreover, MITF is known to play an essential role during the differentiation of the retinal pigment epithelium and can regulate melanogenesis [114,115].
ALS is a highly deregulated pathway (ranking 4), and 10 of its genes are regulated by the proposed TMRs. Indeed, five different TMRs regulate ALS genes (AR, ELF5, HOXA2, nuclear receptor subfamily 0 group B member 1 (NR0B1) and NR5A2), several of which have been related to different types of cancer (C9orf72, multiple cancers; DERL1, non-small-cell lung cancer, bladder cancer, colon cancer, breast carcinoma; DNAH9, gastric adenocarcinoma; DCTN4, colon adenocarcinoma; COX4I1, skin cancer; UQCRC2, gastric cancer, glioma, testicular cancer, colorectal cancer, hepatocarcinoma; PPP3CA, ovarian cancer) [95,116–130]. Preclinical studies [95] showed how the C9orf72 knockout, the most common cause of ALS, caused the proliferation of unexpected tumours in different mouse tissues. It is important to highlight that C9orf72 appears as a target of the TMR AR in our regulatory network.
Hepatocyte differentiation was not dysregulated in our analysis but NR5A2 and HNF4G were identified as TMRs. NR5A2 is also known as hepatocytic TF, and HNF4G is a critical factor for hepatocyte differentiation. It has been reported that HNF4G can exert a carcinogenic effect by promoting cell proliferation and inhibiting cell apoptosis in lung cancer [131]. Moreover, NR5A2 transcriptional activation promotes pancreatic cancer progression [132].
There is, of course, a close relationship between Rb (the intrusive intraocular mostly pediatric tumour) and Rb (the tumour suppressor protein). It is widely known that Rb, when fully functional, often acts to prevent excessive cell growth by inhibiting cell cycle progression by binding to cell cycle proteins. Rb malfunction has been also linked to other cancers such as sarcomas, astrocytomas, melanomas, as well as other epithelial tumours [133], most noticeably breast [134], lung [135,136] and ovarian cancer [137]. This is unsurprising, since it has been argued that Rb is indeed a master regulator of the cell cycle, which is mutated or functionally inactivated in the majority of human cancers [138,139].
5.Conclusion
In this work, we studied Rb genetic expression at several organizational levels. For each of these levels, the existence of healthy tissue samples was of utmost importance. We identified the individual genes whose expression is different in Rb in contrast with normal tissue. We determined the pathways whose global expression pattern is more distant from the global expression observed in normal tissue. Finally, we identified which TFs regulated the highest number of DEGs and proposed them as TMRs. The integration of all these layers allows us to pinpoint genes and pathways relevant to the Rb phenotype.
Rb has been described as a highly undifferentiated malignancy. The enrichment of DEGs in the phototransduction and retrograde endocannabinoid signalling pathways could be associated with abnormal behaviour of the processes leading to cellular differentiation and cellular proliferation. On the other hand, the TMRs NR5A2 and HNF4G are involved in hepatocyte differentiation. The enrichment of aberrant expression in the regulons of these TFs could suggest an abnormal retina development which could be involved in Rb origin and progression.
Acknowledgements
The authors would like to extend a sincere thanks to Dr Carmel McConville and all the original paper authors for the data kindly provided and for helping us by publishing it in the GEO.
Contributor Information
M. Verónica Ponce-Castañeda, Email: vponces@gmail.com.
Hugo Tovar, Email: hatovar@inmegen.gob.mx.
Data accessibility
Publicly available datasets were analysed in this study. These data can be found in the Gene Expression Omnibus (GEO) database (National Center of Biotechnology and Information, NCBI), series GSE172170 (https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE172170), and the samples GSM460264 (https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSM460264), GSM607947 (https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSM607947) and GSM607948 (https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSM607948).
The electronic supplementary material, table S1, contains the list of 305 DEGs found in this work, with an adjusted p-value < 0.05 and an LFC > 2. The electronic supplementary material, table S2, lists the top 10 deregulated pathways per Rb subgroup. Median PDS z-score per pathway, taking into consideration all tumoral samples and all samples per subgroup [140].
Competing interests
No potential competing interest was reported by the authors.
Funding
This work was supported by CONACYT grant no. 285544/2016 Ciencia Básica, and grant no. 2115 Fronteras de la Ciencia; the National Institute of Genomic Medicine (Mexico). Additional support has been granted by the National Laboratory of Complexity Sciences grant no. 232647/2014 CONACYT. E.H.L. acknowledges additional support from the 2016 Marcos Moshinsky Fellowship in the Physical Sciences.
References
Untitled section
References
- 1.Pandey AN. 2014. Retinoblastoma: an overview. Saudi J. Ophthalmol. 28, 310-315. ( 10.1016/j.sjopt.2013.11.001)
- 2.Knudson AG. 1971. Mutation and cancer: statistical study of retinoblastoma. Proc. Natl Acad. Sci. USA 68, 820-823. ( 10.1073/pnas.68.4.820)
- 3.Balaban G, Gilbert F, Nichols W, Meadows AT, Shields J. 1982. Abnormalities of chromosome #13 in retinoblastomas from individuals with normal constitutional karyotypes. Cancer Genet. Cytogen. 6, 213-221. ( 10.1016/0165-4608(82)90058-9)
- 4.Sparkes R, Murphree A, Lingua R, Sparkes M, Field L, Funderburk S, Benedict W. 1983. Gene for hereditary retinoblastoma assigned to human chromosome 13 by linkage to esterase D. Science 219, 971-973. ( 10.1126/science.6823558)
- 5.Strong L, Riccardi V, Ferrell R, Sparkes R. 1981. Familial retinoblastoma and chromosome 13 deletion transmitted via an insertional translocation. Science 213, 1501-1503. ( 10.1126/science.7280668)
- 6.Hanahan D, Weinberg RA. 2000. The hallmarks of cancer. Cell 100, 57-70. ( 10.1016/s0092-8674(00)81683-9)
- 7.Hanahan D, Weinberg RA. 2011. Hallmarks of cancer: the next generation. Cell 144, 646-674. ( 10.1016/j.cell.2011.02.013)
- 8.Thomas DC, Baurley JW, Brown EE, Figueiredo JC, Goldstein A, Hazra A, Wilson RT, Rothman N. 2008. Approaches to complex pathways in molecular epidemiology: summary of a special conference of the American Association for Cancer Research. Cancer Res. 68, 10 028-10 030. ( 10.1158/0008-5472.can-08-1690)
- 9.Bild AH, et al. 2006. Oncogenic pathway signatures in human cancers as a guide to targeted therapies. Nature 439, 353-357. ( 10.1038/nature04296)
- 10.Drier Y, Sheffer M, Domany E. 2013. Pathway-based personalized analysis of cancer. Proc. Natl Acad. Sci. USA 110, 6388-6393. ( 10.1073/pnas.1219651110)
- 11.Carro MS, et al. 2010. The transcriptional network for mesenchymal transformation of brain tumours. Nature 463, 318-325. ( 10.1038/nature08712)
- 12.Lefebvre C, et al. 2010. A human B-cell interactome identifies MYB and FOXM1 as master regulators of proliferation in germinal centers. Mol. Syst. Biol. 6, 377. ( 10.1038/msb.2010.31)
- 13.Fletcher MNC, et al. 2013. Master regulators of FGFR2 signalling and breast cancer risk. Nat. Commun. 4, 2464. ( 10.1038/ncomms3464)
- 14.Aytes A, et al. 2014. Cross-species regulatory network analysis identifies a synergistic interaction between FOXM1 and CENPF that drives prostate cancer malignancy. Cancer Cell 25, 638-651. ( 10.1016/j.ccr.2014.03.017)
- 15.García-Campos MA, Espinal-Enríquez J, Hernández-Lemus E. 2015. Pathway analysis: state of the art. Front. Physiol. 6, 47. ( 10.3389/fphys.2015.00383)
- 16.Khatri P, Sirota M, Butte AJ. 2012. Ten years of pathway analysis: current approaches and outstanding challenges. PLoS Comput. Biol. 8, e1002375. ( 10.1371/journal.pcbi.1002375)
- 17.Gatta GD, et al. 2012. Reverse engineering of TLX oncogenic transcriptional networks identifies RUNX1 as tumor suppressor in T-ALL. Nat. Med. 18, 436-440. ( 10.1038/nm.2610)
- 18.Lim WK, Lyashenko E, Califano A. 2009. Master regulators used as breast cancer metastasis classifier. In Biocomputing 2009, pp. 504-515. Singapore: World Scientific. ( 10.1142/9789812836939_0048)
- 19.Remo A, et al. 2015. Systems biology analysis reveals NFAT5 as a novel biomarker and master regulator of inflammatory breast cancer. J. Transl. Med. 13, 138. ( 10.1186/s12967-015-0492-2)
- 20.Tapia-Carrillo D, Tovar H, Velazquez-Caldelas TE, Hernandez-Lemus E. 2019. Master regulators of signaling pathways: an application to the analysis of gene regulation in breast cancer. Front. Genet. 10, 352. ( 10.3389/fgene.2019.01180)
- 21.Tovar H, García-Herrera R, Espinal-Enríquez J, Hernández-Lemus E. 2015. Transcriptional master regulator analysis in breast cancer genetic networks. Comput. Biol. Chem. 59(Pt B), 67-77. ( 10.1016/j.compbiolchem.2015.08.007)
- 22.Chan SSK, Kyba M. 2013. What is a master regulator? J. Stem Cell Res. Ther. 3, 1000e114. ( 10.4172/2157-7633.1000e114)
- 23.Subramanian A, et al. 2005. Gene set enrichment analysis: a knowledge-based approach for interpreting genome-wide expression profiles. Proc. Natl Acad. Sci. USA 102, 15 545-15 550. ( 10.1073/pnas.0506580102)
- 24.Kapatai G, Brundler MA, Jenkinson H, Kearns P, Parulekar M, Peet AC, McConville CM. 2013. Gene expression profiling identifies different sub-types of retinoblastoma. Brit. J. Cancer 109, 512-525. ( 10.1038/bjc.2013.283)
- 25.Affymetrix. 2015. Microarray normalization using signal space transformation with probe guanine cytosine count correction. See http://tools.thermofisher.com/content/sfs/brochures/sst_gccn_whitepaper.pdf.
- 26.Gentleman RC, et al. 2004. Bioconductor: open software development for computational biology and bioinformatics. Genome Biol. 5, R80. ( 10.1186/gb-2004-5-10-r80)
- 27.R Core Team. 2020. R: a language and environment for statistical computing. Vienna, Austria: R Foundation for Statistical Computing.
- 28.Ritchie ME, Phipson B, Wu D, Hu Y, Law CW, Shi W, Smyth GK. 2015. Limma powers differential expression analyses for RNA-sequencing and microarray studies. Nucleic Acids Res. 43, e47. ( 10.1093/nar/gkv007)
- 29.Warnes GR, et al. 2009. gplots: various R programming tools for plotting data. See https://CRAN.R-project.org/package=gplots.
- 30.Liao Y, Wang J, Jaehnig EJ, Shi Z, Zhang B. 2019. WebGestalt 2019: gene set analysis toolkit with revamped UIs and APIs. Nucleic Acids Res. 47, W199-W205. ( 10.1093/nar/gkz401)
- 31.Hastie T, Stuetzle W. 1989. Principal curves. J. Am. Stat. Assoc. 84, 502. ( 10.2307/2289936)
- 32.Kanehisa M, Goto S. 2000. KEGG: Kyoto Encyclopedia of Genes and Genomes. Nucleic Acids Res. 28, 27-30. ( 10.1093/nar/28.1.27)
- 33.Tanabe M, Kanehisa M. 2012. Using the KEGG database resource. Curr. Protoc. Bioinform. 38, 1.12.1-1.12.43. ( 10.1002/0471250953.bi0112s38)
- 34.Mercatelli D, Lopez-Garcia G, Giorgi FM. 2020. corto: a lightweight R package for gene network inference and master regulator analysis. Bioinformatics 36, 3916-3917. ( 10.1093/bioinformatics/btaa223)
- 35.Chawla K, Tripathi S, Thommesen L, Lægreid A, Kuiper M. 2013. TFcheckpoint: a curated compendium of specific DNA-binding RNA polymerase II transcription factors. Bioinformatics 29, 2519-2520. ( 10.1093/bioinformatics/btt432)
- 36.Dong H, Claffey KP, Brocke S, Epstein PM. 2013. Expression of phosphodiesterase 6 (PDE6) in human breast cancer cells. Springerplus 2, 680. ( 10.1186/2193-1801-2-680)
- 37.Maßberg D, Simon A, Häussinger D, Keitel V, Gisselmann G, Conrad H, Hatt H. 2015. Monoterpene (−)-citronellal affects hepatocarcinoma cell signaling via an olfactory receptor. Arch. Biochem. Biophys. 566, 100-109. ( 10.1016/j.abb.2014.12.004)
- 38.Kalbe B, et al. 2017. Helional-induced activation of human olfactory receptor 2J3 promotes apoptosis and inhibits proliferation in a non-small-cell lung cancer cell line. Eur. J. Cell Biol. 96, 34-46. ( 10.1016/j.ejcb.2016.11.004)
- 39.Safran M, Rosen N, Twik M, BarShir R, Stein TI, Dahary D, Fishilevich S, Lancet D. 2022. The GeneCads suite. In Practical guide to life science databases (eds Abugessaisa I, Kasukawa T), pp. 27-56. Singapore: Springer.
- 40.Manteniotis S, et al. 2016. Functional characterization of the ectopically expressed olfactory receptor 2AT4 in human myelogenous leukemia. Cell Death Discov. 2, 15070. ( 10.1038/cddiscovery.2015.70)
- 41.Flegel C, Manteniotis S, Osthold S, Hatt H, Gisselmann G. 2013. Expression profile of ectopic olfactory receptors determined by deep sequencing. PLoS ONE 8, e55368. ( 10.1371/journal.pone.0055368)
- 42.Tsai T, et al. 2017. Two olfactory receptors—OR2A4/7 and OR51B5—differentially affect epidermal proliferation and differentiation. Exp. Dermatol. 26, 58-65. ( 10.1111/exd.13132)
- 43.Weber L, et al. 2018. Olfactory receptors as biomarkers in human breast carcinoma tissues. Front. Oncol. 8, 33. ( 10.3389/fonc.2018.00033)
- 44.Pyszniak M, Tabarkiewicz J, Łuszczki JJ. 2016. Endocannabinoid system as a regulator of tumor cell malignancy: biological pathways and clinical significance. Oncotargets Ther. 9, 4323-4336. ( 10.2147/ott.s106944)
- 45.Margolin AA, Nemenman I, Basso K, Wiggins C, Stolovitzky G, Favera R, Califano A. 2006. ARACNE: an algorithm for the reconstruction of gene regulatory networks in a mammalian cellular context. BMC Bioinf. 7, S7. ( 10.1186/1471-2105-7-s1-s7)
- 46.Cover TM, Thomas JA. 2006. Elements of information theory, 2nd edn. New York, NY: Wiley.
- 47.Ming Q, Roske Y, Schuetz A, Walentin K, Ibraimi I, Schmidt-Ott KM, Heinemann U. 2018. Structural basis of gene regulation by the Grainyhead/CP2 transcription factor family. Nucleic Acids Res. 46, gkx1299. ( 10.1093/nar/gkx1299)
- 48.Petrof G, et al. 2014. Mutations in GRHL2 result in an autosomal-recessive ectodermal dysplasia syndrome. Am. J. Hum. Genet. 95, 308-314. ( 10.1016/j.ajhg.2014.08.001)
- 49.Takahashi K, Yamanaka S. 2006. Induction of pluripotent stem cells from mouse embryonic and adult fibroblast cultures by defined factors. Cell 126, 663-676. ( 10.1016/j.cell.2006.07.024)
- 50.Reichow SL, Hamma T, Ferré-D'Amaré AR, Varani G. 2007. The structure and function of small nucleolar ribonucleoproteins. Nucleic Acids Res. 35, 1452-1464. ( 10.1093/nar/gkl1172)
- 51.Williams GT, Farzaneh F. 2012. Are snoRNAs and snoRNA host genes new players in cancer? Nat. Rev. Cancer 12, 84-88. ( 10.1038/nrc3195)
- 52.Nallar SC, Kalvakolanu DV. 2013. Regulation of snoRNAs in cancer: close encounters with interferon. J. Interf. Cytokine Res. 33, 189-198. ( 10.1089/jir.2012.0106)
- 53.Liang J, Wen J, Huang Z, Chen X, Zhang B, Chu L. 2019. Small nucleolar RNAs: insight into their function in cancer. Front. Oncol. 9, 587. ( 10.3389/fonc.2019.00587)
- 54.Su H, Xu T, Ganapathy S, Shadfan M, Long M, Huang THM, Thompson I, Yuan ZM. 2014. Elevated snoRNA biogenesis is essential in breast cancer. Oncogene 33, 1348-1358. ( 10.1038/onc.2013.89)
- 55.Sushko T, Gilep A, Usanov S. 2014. Genetics, structure, function, mode of actions and role in cancer development of CYP17. Anti-Cancer Agents Med. Chem. 14, 66-76. ( 10.2174/187152061131300330)
- 56.Levin ER, Hammes SR. 2016. Nuclear receptors outside the nucleus: extranuclear signalling by steroid receptors. Nat. Rev. Mol. Cell Biol. 17, 783-797. ( 10.1038/nrm.2016.122)
- 57.Africander D, Storbeck KH. 2018. Steroid metabolism in breast cancer: where are we and what are we missing? Mol. Cell Endocrinol. 466, 86-97. ( 10.1016/j.mce.2017.05.016)
- 58.Cheng M, Michalski S, Kommagani R. 2018. Role for growth regulation by estrogen in breast cancer 1 (GREB1) in hormone-dependent cancers. Int. J. Mol. Sci. 19, 2543. ( 10.3390/ijms19092543)
- 59.Picard D. 1998. The role of heat-shock proteins in the regulation of steroid receptor function. In Molecular biology of steroid and nuclear hormone receptors (ed. Freedman LP), pp. 1-18. Boston, MA: Birkhauser.
- 60.van ‘t Veer LJ, et al. 2002. Gene expression profiling predicts clinical outcome of breast cancer. Nature 415, 530-536. ( 10.1038/415530a)
- 61.Heinlein CA, Chang C. 2004. Androgen receptor in prostate cancer. Endocr. Rev. 25, 276-308. ( 10.1210/er.2002-0032)
- 62.Bloomfield CD, Munck AU, Smith KA. 1984. Glucocorticoid receptor levels predict response to treatment in human lymphoma. Prog. Clin. Biol. Res. 142, 223-233.
- 63.Fereshteh MP, Tilli MT, Kim SE, Xu J, O'Malley BW, Wellstein A, Furth PA, Riegel AT. 2008. The nuclear receptor coactivator amplified in breast cancer-1 is required for neu (erbb2/her2) activation, signaling, and mammary tumorigenesis in mice. Cancer Res. 68, 3697-3706. ( 10.1158/0008-5472.can-07-6702)
- 64.Osborne CK, et al. 2003. Role of the estrogen receptor coactivator AIB1 (SRC-3) and HER-2/neu in tamoxifen resistance in breast cancer. J. Natl Cancer Inst. 95, 353-361. ( 10.1093/jnci/95.5.353)
- 65.Chung ACK, Zhou S, Liao L, Tien JCY, Greenberg NM, Xu J. 2007. Genetic ablation of the amplified-in-breast cancer 1 inhibits spontaneous prostate cancer progression in mice. Cancer Res. 67, 5965-5975. ( 10.1158/0008-5472.can-06-3168)
- 66.Eke I, Cordes N. 2015. Focal adhesion signaling and therapy resistance in cancer. Semin. Cancer Biol. 31, 65-75. ( 10.1016/j.semcancer.2014.07.009)
- 67.Zhang YL, Wang RC, Cheng K, Ring BZ, Su L. 2017. Roles of Rap1 signaling in tumor cell migration and invasion. Cancer Biol. Med. 14, 90-99. ( 10.20892/j.issn.2095-3941.2016.0086)
- 68.Yoshida GJ. 2015. Metabolic reprogramming: the emerging concept and associated therapeutic strategies. J. Exp. Clin. Cancer Res. 34, 111. ( 10.1186/s13046-015-0221-y)
- 69.DeBerardinis RJ, Chandel NS. 2016. Fundamentals of cancer metabolism. Sci. Adv. 2, e1600200. ( 10.1126/sciadv.1600200)
- 70.Dang CV, Kim J. 2018. Convergence of cancer metabolism and immunity: an overview. Biomol. Ther. 26, 4-9. ( 10.4062/biomolther.2017.194)
- 71.Koundouros N, Poulogiannis G. 2019. Reprogramming of fatty acid metabolism in cancer. Brit. J. Cancer 122, 4-22. ( 10.1038/s41416-019-0650-z)
- 72.McKeehan WL, McKeehan KA. 1979. Oxocarboxylic acids, pyridine nucleotide-linked oxidoreductases and serum factors in regulation of cell proliferation. J. Cell. Physiol. 101, 9-16. ( 10.1002/jcp.1041010103)
- 73.Cooper AJL, Ginos JZ, Meister A. 1983. Synthesis and properties of the α-keto acids. Chem. Rev. 83, 321-358. ( 10.1021/cr00055a004)
- 74.Titov DV, Cracan V, Goodman RP, Peng J, Grabarek Z, Mootha VK. 2016. Complementation of mitochondrial electron transport chain by manipulation of the NAD+/NADH ratio. Science 352, 231-235. ( 10.1126/science.aad4017)
- 75.Röhrig F, Schulze A. 2016. The multifaceted roles of fatty acid synthesis in cancer. Nat. Rev. Cancer 16, 732-749. ( 10.1038/nrc.2016.89)
- 76.Munir R, Lisec J, Swinnen JV, Zaidi N. 2019. Lipid metabolism in cancer cells under metabolic stress. Brit. J. Cancer 120, 1090-1098. ( 10.1038/s41416-019-0451-4)
- 77.Schoors S, et al. 2015. Fatty acid carbon is essential for dNTP synthesis in endothelial cells. Nature 520, 192-197. ( 10.1038/nature14362)
- 78.Tachibana K, Yamasaki D, Ishimoto K, Doi T. 2008. The role of PPARs in cancer. PPAR Res. 2008, 102737. ( 10.1155/2008/102737)
- 79.Messerschmidt L, Fischer S, Wiedemann P, Bringmann A, Hollborn M. 2019. Osmotic induction of cyclooxygenase-2 in RPE cells: stimulation of inflammasome activation. Mol. Vis. 25, 329-344.
- 80.Lake N. 1994. Taurine and GABA in the rat retina during postnatal development. Vis. Neurosci. 11, 253-260. ( 10.1017/s0952523800001619)
- 81.Sahoo S, et al. 2019. Metabolite systems profiling identifies exploitable weaknesses in retinoblastoma. FEBS Lett. 593, 23-41. ( 10.1002/1873-3468.13294)
- 82.Young SZ, Bordey A. 2009. GABA's control of stem and cancer cell proliferation in adult neural and peripheral niches. Physiology 24, 171-185. ( 10.1152/physiol.00002.2009)
- 83.Hinz B, Ramer R. 2019. Anti-tumour actions of cannabinoids. Brit. J. Pharmacol. 176, 1384-1394. ( 10.1111/bph.14426)
- 84.Tuerxun H, Cui J. 2018. The dual effect of morphine on tumor development. Clin. Transl. Oncol. 21, 695-701. ( 10.1007/s12094-018-1974-5)
- 85.Ramer R, Merkord J, Rohde H, Hinz B. 2010. Cannabidiol inhibits cancer cell invasion via upregulation of tissue inhibitor of matrix metalloproteinases-1. Biochem. Pharmacol. 79, 955-966. ( 10.1016/j.bcp.2009.11.007)
- 86.Ramer R, et al. 2012. Cannabidiol inhibits lung cancer cell invasion and metastasis via intercellular adhesion molecule-1. FASEB J. 26, 1535-1548. ( 10.1096/fj.11-198184)
- 87.Ligresti A, et al. 2006. Antitumor activity of plant cannabinoids with emphasis on the effect of cannabidiol on human breast carcinoma. J. Pharmacol. Exp. Ther. 318, 1375-1387. ( 10.1124/jpet.106.105247)
- 88.Massi P, Vaccani A, Ceruti S, Colombo A, Abbracchio MP, Parolaro D. 2004. Antitumor effects of cannabidiol, a nonpsychoactive cannabinoid, on human glioma cell lines. J. Pharmacol. Exp. Ther. 308, 838-845. ( 10.1124/jpet.103.061002)
- 89.Shrivastava A, Kuzontkoski PM, Groopman JE, Prasad A. 2011. Cannabidiol induces programmed cell death in breast cancer cells by coordinating the cross-talk between apoptosis and autophagy. Mol. Cancer Ther. 10, 1161-1172. ( 10.1158/1535-7163.mct-10-1100)
- 90.Nabissi M, Morelli MB, Amantini C, Liberati S, Santoni M, Ricci-Vitiani L, Pallini R, Santoni G. 2015. Cannabidiol stimulates Aml-1a-dependent glial differentiation and inhibits glioma stem-like cells proliferation by inducing autophagy in a TRPV2-dependent manner. Int. J. Cancer 137, 1855-1869. ( 10.1002/ijc.29573)
- 91.Riancho J, Delgado-Alvarado M, Andreu MD, Paz-Fajardo L, Arozamena S, Gil-Bea FJ, de Munaín AL. 2020. Amyotrophic lateral sclerosis (ALS), cancer, autoimmunity and metabolic disorders: an unsolved tantalizing challenge. Brit. J. Pharmacol. 178, 1269-1278. ( 10.1111/bph.15151)
- 92.Dhillon AS, Hagan S, Rath O, Kolch W. 2007. MAP kinase signalling pathways in cancer. Oncogene 26, 3279-3290. ( 10.1038/sj.onc.1210421)
- 93.Chu WM. 2013. Tumor necrosis factor. Cancer Lett. 328, 222-225. ( 10.1016/j.canlet.2012.10.014)
- 94.Olsson M, Zhivotovsky B. 2011. Caspases and cancer. Cell Death Differ. 18, 1441-1449. ( 10.1038/cdd.2011.30)
- 95.Sudria-Lopez E, et al. 2016. Full ablation of C9orf72 in mice causes immune system-related pathology and neoplastic events but no motor neuron defects. Acta Neuropathol. 132, 145-147. ( 10.1007/s00401-016-1581-x)
- 96.Shang S, Hua F, Hu ZW. 2017. The regulation of β-catenin activity and function in cancer: therapeutic opportunities. Oncotarget 8, 33 972-33 989. ( 10.18632/oncotarget.15687)
- 97.Desgrosellier JS, Cheresh DA. 2010. Integrins in cancer: biological implications and therapeutic opportunities. Nat. Rev. Cancer 10, 9-22. ( 10.1038/nrc2748)
- 98.Xu X, Liu Z, Tian F, Xu J, Chen Y. 2019. Clinical significance of transcription factor 7 (TCF7) as a prognostic factor in gastric cancer. Med. Sci. Monitor 25, 3957-3963. ( 10.12659/msm.913913)
- 99.Camassei FD, et al. 2003. Expression of the lipogenic enzyme fatty acid synthase (FAS) in retinoblastoma and its correlation with tumor aggressiveness. Invest. Ophth. Vis. Sci. 44, 2399-2403. ( 10.1167/iovs.02-0934)
- 100.Oikawa T, Yamada T. 2003. Molecular biology of the Ets family of transcription factors. Gene 303, 11-34. ( 10.1016/s0378-1119(02)01156-3)
- 101.Oikawa T. 2004. ETS transcription factors: possible targets for cancer therapy. Cancer Sci. 95, 626-633. ( 10.1111/j.1349-7006.2004.tb03320.x)
- 102.Zhang X, Feng Q, Cote RH. 2005. Efficacy and selectivity of phosphodiesterase-targeted drugs in inhibiting photoreceptor phosphodiesterase (PDE6) in retinal photoreceptors. Invest. Ophth. Vis. Sci. 46, 3060-3066. ( 10.1167/iovs.05-0257)
- 103.Cote RH. 2004. Characteristics of photoreceptor PDE (PDE6): similarities and differences to PDE5. Int. J. Impot. Res. 16, S28-S33. ( 10.1038/sj.ijir.3901212)
- 104.Morales B, Bacigalupo J. 1996. Chemical reception in vertebrate olfaction: evidence for multiple transduction pathways. Biol. Res. 29, 333-341.
- 105.Giessel AJ, Datta SR. 2014. Olfactory maps, circuits and computations. Curr. Opin. Neurobiol. 24, 120-132. ( 10.1016/j.conb.2013.09.010)
- 106.Park SW. 2014. Understanding the human sensory conduction of smell. Hanyang Med. Rev. 34, 100-106. ( 10.7599/hmr.2014.34.3.100)
- 107.Maßberg D, Hatt H. 2018. Human olfactory receptors: novel cellular functions outside of the nose. Physiol. Rev. 98, 1739-1763. ( 10.1152/physrev.00013.2017)
- 108.Fuessel S, Weigle B, Schmidt U, Baretton G, Koch R, Bachmann M, Rieber EP, Wirth MP, Meye A. 2006. Transcript quantification of Dresden G protein-coupled receptor (D-GPCR) in primary prostate cancer tissue pairs. Cancer Lett. 236, 95-104. ( 10.1016/j.canlet.2005.05.005)
- 109.Weigle B, et al. 2004. D-GPCR: a novel putative G protein-coupled receptor overexpressed in prostate cancer and prostate. Biochem. Biophys. Res. Commun. 322, 239-249. ( 10.1016/j.bbrc.2004.07.106)
- 110.Weng J, Wang J, Cai Y, Stafford LJ, Mitchell D, Ittmann M, Liu M. 2005. Increased expression of prostate-specific G-protein-coupled receptor in human prostate intraepithelial neoplasia and prostate cancers. Int. J. Cancer 113, 811-818. ( 10.1002/ijc.20635)
- 111.Wang J, Weng J, Cai Y, Penland R, Liu M, Ittmann M. 2006. The prostate-specific G-protein coupled receptors PSGR and PSGR2 are prostate cancer biomarkers that are complementary to α-methylacyl-CoA racemase. Prostate 66, 847-857. ( 10.1002/pros.20389)
- 112.Lachmann A, Xu H, Krishnan J, Berger SI, Mazloom AR, Ma'ayan A. 2010. ChEA: transcription factor regulation inferred from integrating genome-wide ChIP-X experiments. Bioinformatics 26, 2438-2444. ( 10.1093/bioinformatics/btq466)
- 113.Wu Q, Zhang W, Xue L, Wang Y, Fu M, Ma L, Song Y, Zhan QM. 2019. APC/C-CDH1-regulated IDH3β coordinates with the cell cycle to promote cell proliferation. Cancer Res. 79, 3281-3293. ( 10.1158/0008-5472.can-18-2341)
- 114.Amae S, et al. 1998. Identification of a novel isoform of microphthalmia-associated transcription factor that is enriched in retinal pigment epithelium. Biochem. Biophys. Res. Commun. 247, 710-715. ( 10.1006/bbrc.1998.8838)
- 115.Jo IS, Sohn KC, Kim CD, Lee YH. 2016. Effect of MITF-M and MITF-A overexpression on the dendrtic formation in melanocytes. Korean J. Phys. Anthropol. 29, 27-34. ( 10.11637/kjpa.2016.29.1.27)
- 116.Mao M, Zhang J, Jiang J. 2018. Overexpression of Derlin-1 is associated with poor prognosis in patients with non-small cell lung cancer. Ann. Clin. Lab. Sci. 48, 29-34.
- 117.Yang F, et al. 2018. MiR-598 suppresses invasion and migration by negative regulation of Derlin-1 and epithelial-mesenchymal transition in non-small cell lung cancer. Cell. Physiol. Biochem. 47, 245-256. ( 10.1159/000489803)
- 118.Wu Z, Wang C, Zhang Z, Liu W, Xu H, Wang H, Wang Y, Zhang W, Wang SL. 2016. High expression of Derlin-1 is associated with the malignancy of bladder cancer in a Chinese Han population. PLoS ONE 11, e0168351. ( 10.1371/journal.pone.0168351)
- 119.Tan X, He X, Jiang Z, Wang X, Ma L, Liu L, Wang X, Fan Z, Su D. 2015. Derlin-1 is overexpressed in human colon cancer and promotes cancer cell proliferation. Mol. Cell Biochem. 408, 205-213. ( 10.1007/s11010-015-2496-x)
- 120.Wang J, Hua H, Ran Y, Zhang H, Liu W, Yang Z, Jiang Y. 2008. Derlin-1 is overexpressed in human breast carcinoma and protects cancer cells from endoplasmic reticulum stress-induced apoptosis. Breast Cancer Res. 10, R7. ( 10.1186/bcr1849)
- 121.Li D, Shi M, Ji H, Chen G, Jiang H, Wang Z. 2016. MicroRNA-181d is a tumor suppressor in human esophageal squamous cell carcinoma inversely regulating Derlin-1. Oncol. Rep. 36, 2041-2048. ( 10.3892/or.2016.5028)
- 122.Zhu C, Yang Q, Xu J, Zhao W, Zhang Z, Xu D, Zhang Y, Zhao E, Zhao G. 2019. Somatic mutation of DNAH genes implicated higher chemotherapy response rate in gastric adenocarcinoma patients. J. Transl. Med. 17, 109. ( 10.1186/s12967-019-1867-6)
- 123.Wang S, et al. 2018. Distinct prognostic value of dynactin subunit 4 (DCTN4) and diagnostic value of DCTN1, DCTN2, and DCTN4 in colon adenocarcinoma. Cancer Manage. Res. 10, 5807-5824. ( 10.2147/cmar.s183062)
- 124.Dang C, Gottschling M, Manning K, O'Currain E, Schneider S, Sterry W, Stockfleth E, Nindl I. 2006. Identification of dysregulated genes in cutaneous squamous cell carcinoma. Oncol. Rep. 16, 513-519. ( 10.3892/or.16.3.513)
- 125.Wang D, Su F, Zhang T, Yang T, Wang H, Yang L, Zhou F, Feng M. 2020. The miR-370/UQCRC2 axis facilitates tumorigenesis by regulating epithelial-mesenchymal transition in gastric cancer. J. Cancer 11, 5042-5055. ( 10.7150/jca.45553)
- 126.Bai Y, et al. 2018. A novel tumor-suppressor, CDH18, inhibits glioma cell invasiveness via UQCRC2 and correlates with the prognosis of glioma patients. Cell. Physiol. Biochem. 48, 1755-1770. ( 10.1159/000492317)
- 127.Shang Y, Zhang F, Li D, Li C, Li H, Jiang Y, Zhang D. 2018. Overexpression of UQCRC2 is correlated with tumor progression and poor prognosis in colorectal cancer. Pathol. Res. Pract. 214, 1613-1620. ( 10.1016/j.prp.2018.08.012)
- 128.Selvam MKP, Agarwal A, Pushparaj PN. 2019. A quantitative global proteomics approach to understanding the functional pathways dysregulated in the spermatozoa of asthenozoospermic testicular cancer patients. Andrology 7, 454-462. ( 10.1111/andr.12620)
- 129.Jang SI, Lee YW, Cho CK, Yoo HS, Jang JH. 2013. Identification of target genes involved in the antiproliferative effect of enzyme-modified ginseng extract in HepG2 hepatocarcinoma cell. Evid. Based Compl. Alt. Med. 2013, 502568. ( 10.1155/2013/502568)
- 130.Xin B, Ji KQ, Liu YS, Zhao XD. 2019. Higher expression of calcineurin predicts poor prognosis in unique subtype of ovarian cancer. J. Ovarian Res. 12, 75. ( 10.1186/s13048-019-0550-0)
- 131.Wang J, Zhang J, Xu L, Zheng Y, Ling D, Yang Z. 2017. Expression of HNF4G and its potential functions in lung cancer. Oncotarget 9, 18 018-18 028. ( 10.18632/oncotarget.22933)
- 132.Guo F, Zhou Y, Guo H, Ren D, Jin X, Wu H. 2021. NR5A2 transcriptional activation by BRD4 promotes pancreatic cancer progression by upregulating GDF15. Cell Death Discov. 7, 78. ( 10.1038/s41420-021-00462-8)
- 133.Fletcher O, Easton D, Anderson K, Gilham C, Jay M, Peto J. 2004. Lifetime risks of common cancers among retinoblastoma survivors. J. Natl Cancer Inst. 96, 357-363. ( 10.1093/jnci/djh058)
- 134.Knudsen ES, Zacksenhaus E. 2018. The vulnerability of RB loss in breast cancer: targeting a void in cell cycle control. Oncotarget 9, 30 940-30 941. ( 10.18632/oncotarget.25797)
- 135.Rubin SM, Sage J. 2019. Manipulating the tumour-suppressor protein Rb in lung cancer reveals possible drug targets. Nature 569, 343-344. ( 10.1038/d41586-019-01319-y)
- 136.Tanaka H, Fujii Y, Hirabayashi H, Miyoshi S, Sakaguchi M, Yoon H, Matsuda H. 1998. Disruption of the RB pathway and cell-proliferative activity in non-small-cell lung cancers. Int. J. Cancer 79, 111-115. ( 10.1002/(sici)1097-0215(19980417)79:2<111::aid-ijc2>3.0.co;2-w)
- 137.Corney DC, Flesken-Nikitin A, Choi J, Nikitin AY. 2008. Role of p53 and Rb in ovarian cancer. Adv. Exp. Med. Biol. 622, 99-117. ( 10.1007/978-0-387-68969-2_9)
- 138.Vélez-Cruz R, Johnson DG. 2017. The retinoblastoma (RB) tumor suppressor: pushing back against genome instability on multiple fronts. Int. J. Mol. Sci. 18, 1776. ( 10.3390/ijms18081776)
- 139.Giacinti C, Giordano A. 2006. RB and cell cycle progression. Oncogene 25, 5220-5227. ( 10.1038/sj.onc.1209615)
- 140.Gómez-Romero L, Alvarez-Suarez DE, Hernández-Lemus E, Ponce-Castañeda MV, Tovar H. 2022. The regulatory landscape of retinoblastoma: a pathway analysis perspective. Figshare. ( 10.6084/m9.figshare.c.5975228)
Associated Data
Data Citations
References
- Gómez-Romero L, Alvarez-Suarez DE, Hernández-Lemus E, Ponce-Castañeda MV, Tovar H. 2022. The regulatory landscape of retinoblastoma: a pathway analysis perspective. Figshare. ( 10.6084/m9.figshare.c.5975228)
Data Availability Statement
Publicly available datasets were analysed in this study. These data can be found in the Gene Expression Omnibus (GEO) database (National Center of Biotechnology and Information, NCBI), series GSE172170 (https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE172170), and the samples GSM460264 (https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSM460264), GSM607947 (https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSM607947) and GSM607948 (https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSM607948).
The electronic supplementary material, table S1, contains the list of 305 DEGs found in this work, with an adjusted p-value < 0.05 and an LFC > 2. The electronic supplementary material, table S2, lists the top 10 deregulated pathways per Rb subgroup. Median PDS z-score per pathway, taking into consideration all tumoral samples and all samples per subgroup [140].