Ras-Mediated Deregulation of the Circadian Clock in Cancer

Summary

Living systems possess an endogenous time-generating system – the circadian clock - accountable for a 24 hours oscillation in the expression of about 10% of all genes. In mammals, disruption of oscillations is associated to several diseases including cancer. In this manuscript, we address the following question: what are the elicitors of a disrupted clock in cancer? We applied a systems biology approach to correlate experimental, bioinformatics and modelling data and could thereby identify key genes which discriminate strong and weak oscillators among cancer cell lines. Most of the discriminative genes play important roles in cell cycle regulation, DNA repair, immune system and metabolism and are involved in oncogenic pathways such as the RAS/MAPK. To investigate the potential impact of the Ras oncogene in the circadian clock we generated experimental models harbouring conditionally active Ras oncogenes. We put forward a direct correlation between the perturbation of Ras oncogene and an effect in the expression of clock genes, found by means of mathematical simulations and validated experimentally. Our study shows that perturbations of a single oncogene are sufficient to deregulate the mammalian circadian clock and opens new ways in which the circadian clock can influence disease and possibly play a role in therapy.

Introduction

All mammalian cells hold an internal circadian clock able to generate daily-endogenous rhythms with a period of approximately 24 hours. Circadian clocks are evolutionary conserved and regulate the expression of about 10% of all genes [1]–[3]. This time-generating mechanism enables the organism to react to external clues, to anticipate environmental changes and to adapt molecular and behavioural processes to specific day-times with the advantage of separating incompatible metabolic processes.

In mammals, the circadian system is hierarchically organized into two major levels of regulation including a main clock, located within the suprachiasmatic nucleus (SCN) and peripheral oscillators [4], [5]. The peripheral clocks can be found in almost all cells in the body. These are able to respond to and synchronize to output signals of the SCN clock, thereby assuring time-precision of molecular processes throughout the organism [6], [7].

Interconnected genetic networks of transcriptional and translational steps drive the oscillator in each individual clock within a cell [8], [9]. The network system can be represented by a core of two main feedback loops: the RORs/Bmal1/REV-ERBs loop and the PERs/CRYs loop which generate oscillations [10], [11]. The core-clock network regulates a series of clock-controlled genes (CCGs) with relevant functions in several cellular and biological processes. CCGs are involved in metabolism, detoxification, cell cycle, cell growth and survival, DNA damage responses and the immune system [12]–[15].

Malfunctions of the circadian clock have been reported in the context of many disorders [16]–[19]. Epidemiological studies have shown that increasing nocturnal light and overnight shift work coincided with a steady increase in the incidence of cancer [20]–[22]. Moreover, the long-term survival rate of patients with metastatic colorectal cancer was about 5-fold higher in those patients with normal clock compared to patients with a severely disrupted clock [23]. Additionally, accelerated malignant growth was observed in mice with an ablated SCN or subjected to experimental chronic jet-lag [24].

Hence, the clock regulation of molecular processes has severe consequences on therapy optimization, and timing of drug intake in cancer [25]–[28]. Cancer patients with altered circadian rhythms have a poorer prognosis [29] and chronotherapy, the administration of anti-cancer drugs at specific times of the day, can improve treatment efficacy as chemotherapeutics may act differently on their targets depending on the time of administration [25], . Furthermore, the rhythmic delivery of cancer therapeutics in colorectal cancer increased the efficacy of oxaliplatin on treatment and patient survival [32], [33].

A molecular basis for the effects of circadian rhythms on cancer patients might be provided by the association of several core-clock genes with cancer promoting mechanisms such as DNA damage [34] [35], metabolism [17] and cell cycle [36]–[38]. Cell cycle check point regulators such as Wee1 (G₂- M transition), Myc (G₀- G₁ transition), and cyclin D1 (G₁- S transition) have been shown to be under the direct regulation of the circadian clock and could represent one way in which the circadian clock regulates cell division [3], [18], [28]. The histone deacetylase sirtuin 1 (SIRT1), a key regulator of metabolism, has recently been identified as a core-clock component [39], [40]. The PER1 and Timeless proteins interact with proteins involved in DNA damage response and Per1 overexpression suppresses growth of human cancer cell lines [26], [41]. Expression of Per1 and Per2 is downregulated in colon, breast and endometrial carcinoma [41], [42]. Per2 is also downregulated in several human lymphoma cell lines and in non-small-cell lung cancer tissues [26], [43]–[45]. In addition, mutations in the clock gene Npas2, a paralog of Clock, have been associated with increased risk of breast cancer and non-Hodgkin's lymphoma [27], [46]. Furthermore, mutations in the gene Clock were found in colon cancer cell lines [47].

These results suggest the existence of a strong cross-regulation between the components of the circadian clock and proto-oncogenes or tumour suppressors. The circadian clock may act as tumour suppressor, whereas a disturbed clock might render the organism more cancer prone. However, a comprehensive view of how cancer genes and clock can influence each other is missing. This motivated us to investigate the mechanism by which the circadian clock might be altered in cancer models.

In the present manuscript, we provide a systems biology approach for the investigation of the circadian clock in several cancer cell lines including colon and skin. Surprisingly, we found strong and weak circadian oscillators within the same type of cancer. We recovered a set of genes which allow the discrimination between the two types of oscillators. Using a theoretical bioinformatics approach, we extended the core-clock network to include a larger set of clock/controlled genes involved in several biological processes and present its interconnection to the discriminative set of genes. Furthermore, we analysed the connections of such discriminative list of genes to cancer pathways among which is the RAS/MAPK pathway. Using experimental data and mathematical modelling, we provide evidence for a putative connection between both systems. Our work provides novel evidence pointing to RAS oncogene being one of the modulators of the mammalian circadian clock.

Results

Cancer cell lines show a rich variety of circadian phenotypes

To investigate a possible link between the circadian core-clock oscillator and tumour-associated pathways, we tested the circadian properties of colorectal cancer cell lines well characterized for their genetic properties and oncogenic pathways. The impact of the circadian systems is apparent for this cancer type, because chronotherapy has shown promising results in colon cancer patients [32], [33].

We studied the oscillation dynamics in the colon carcinoma cell lines HT29, RKO, SW480, LIM1215, CaCo2 and HTC116 using a live-cell imaging approach based on ectopic expression of a luciferase reporter construct driven by the 0.9 kb Bmal1 promoter fragment. As a control we analysed the human osteosarcoma cell line U2OS, a widely used in vitro model to study properties of the mammalian circadian clock. We define a cell line with a clear circadian period and an amplitude variation of at least 20% as a strong oscillator. Surprisingly, the clock properties were very diverse among the colon cancer cell lines (Table 1), which showed strong (Figure 1B, C) and weak to no-oscillation phenotypes (Figure 1D–G). The strongly oscillating cell lines HCT116 and SW480 exhibited shorter doubling times (unpublished observations) indicating an association of perturbations of the circadian clock and effects on the cell cycle and in agreement with previous studies [48]. Moreover, we observed differences in the mRNA levels of the core-clock genes. In Figure 1H, we plotted the relative fold change to Bmal1 for each clock gene and cell line. The highest value for Cry1 is observed in a weak oscillator cell line (CaCo2), while the highest value for Cry2 is observed in the strong oscillator cell lines. For Per1 and Per2, the highest fold change is observed in the strong oscillators. For Clock and Npas2, weak oscillators show the highest relative change to Bmal1. These data clearly show a correlation between expression levels of core-clock genes and the cell oscillator phenotype.

Microarray analysis reveals a set of genes discriminating strong and weak oscillators

To further investigate potential differences between strong and weak oscillators at the gene expression level, we performed transcriptome analysis for each cell line. A leave-one-out cross validation strategy identified a set of differentially expressed genes that discriminate between strong and weak oscillators (Figure 2A). For each cell line we excluded both replica samples once and determined the 100 most significantly expressed probe sets, allowing the identification of two clock groups using a moderated t-test to select the top genes by confidence. The resulting list of 45 best p-value discriminative genes (Table 2 and Figure 2A) allows positioning of U2OS accordingly within the strong oscillator cluster, although this cell line was not previously used for the list generation. Additionally we tested 5 other colon cancer cell lines (Colo205, SW620, SW403, HKe3, HKe-clone8 with/without mifepristone) and two keratinocyte cell lines (HaCaT and A5RT3). Data is shown in Figure S5 and Text S4. The predefined list of 45 discriminative genes could be used to correctly classify seven out of eight cell lines (Text S4). Using the binomial test, we calculated the probability of observing the same or better classification result by random chance. According to this experiment our classifier performs significantly better than randomly expected (p = 0.03516).

Moreover, previous data [49] shows that the siRNA-dependent -knockdown of the majority of these 45 genes confers a circadian phenotype in U2OS cell lines, underlining their potential importance as regulators of the circadian system.

The methodology used to generate the list of 45 discriminative genes, including their p-values, is explained in detail in Text S1 and the resulting p-values for the short list of 45 genes are additionally given in Table 2.

The composition of the set of top genes indicates that phenotypic circadian clock differences are reflected by gene expression differences both in genes of the core network (Figure 1H), but also in additional genes not directly associated with circadian clock functions.

The mammalian circadian clock – an extended regulatory network of a time-generating mechanism

To explore a potential correlation of the discriminative set of genes to the mammalian circadian clock we extend the known core-clock network, which encompasses several genes and proteins interconnected by positive and negative feedback loops [10], to a layer of next neighbours. We scanned a total of 21.4 million abstracts by text mining. After careful curation of the text mining results, we assembled a comprehensive genetic network for the mammalian circadian clock (Figure 2B). First, we gathered a network containing the currently known core-elements of the circadian pathway including the 14 genes: Per1,2,3, Cry1,2, Bmal1,2, Rev-Erbα,β, Rorα,β,γ, Clock and Npas2. In the next step, we added 16 elements reported to be directly interacting to this core set [13]. Subsequently, we used our text mining software GeneView [50] to extract all reported interactions between (i) any two elements of the network and (ii) new elements with direct connections to the 14-element core. As a validation of the performance of GeneView, we analysed the part of the network containing common elements to our previously published network [13]. The software provided evidence for 85% of the interactions described (manuscript in preparation). Additionally, we found 17 new elements in the outer shell and 108 novel interactions, after curation, supported by 132 PubMed references. The enriched circadian-core network contains a total of 47 elements and 229 interactions (Figure 2B, Text S2). This network provides an improved level of specificity regarding interactions of core-clock and clock-controlled genes, in comparison to previously published networks [12], [13]. It includes both protein-protein and DNA-protein interactions and selectively assembles elements which are as well reported to be able to influence the core-clock Figure 2B.

To analyse whether and how the network might convey circadian information to specific output pathways, we performed a detailed analysis of the KEGG pathways associated with the different elements as well as of the data set obtained by text mining analysis (Text S2). We found that many of the genes in the periphery of the network play important roles in pathways frequently deregulated in cancer. Examples are the Wnt pathway (CSNK1ε, CSNK2α, p300, GSK3β, βTRC), the TGF-β signalling pathway (CBP, p300, TNFα), the Jak-STAT signalling pathway (CBP, p300, IFNα) and the MAPK signalling pathway (CBP, p300, GSK3β, PPARγ, TNFα). We also found many genes involved in cell cycle and repair mechanisms (NONO, PPAR1, GSK3β) and in immune defense (CREB, GSK3β, AMPK, TNFα, PGC-1α). Moreover, many of the new network elements were also found to link the circadian system to metabolism and xenobiotics detoxification mechanisms and as such are also potentially relevant in terms of therapy and drug response (PPARα, TEF, ALAS1, AhR, CAR, HLF, E4BP4, DEC2, DEC1, DBP).

Altogether, these results clearly support the view that the circadian clock and oncogenic pathways are strongly connected.

A comprehensive network of circadian regulation in tumourigenesis

To investigate how the identified set of discriminative genes links cancer genes and circadian regulators we again used the text mining software GeneView to create interaction networks between i) the core-clock genes (Figure 2B) and the 45 discriminative genes (Table 2, Text S1); ii) a set of known colon cancer-related genes ([51], Table 3) and the 45 discriminative genes and iii) between all three sets of genes. Only interactions involving elements of two different gene sets were considered (Figure 2B, C). These interactions include both protein-protein as well as DNA-protein interactions. All interactions extracted by the text mining pipeline as well as 391,434 interactions contained in the STRING database version 9.0 [52] were collected. In total 646 Interactions were found to be involved in the assembly of the network (Figure S1). We found 184 connections between our discriminative set of genes and the genes of interest (clock genes and cancer genes, Figure 2C, D). To test whether the number and counts of interactions point to a specific connection between the discriminative genes and the clock/cancer genes, we now tested the significance of the number of connections by comparing to networks build from random gene sets of the same size (45 genes). We created a total of 50 such random sets where each gene of the set was randomly selected from a bin containing genes with the same number of PubMed citations as a gene in the discriminative set. We obtained an average of 165 interactions, standard deviation 28 yielding a p-value of 4.4e-5 when using a non-parametric Wilcoxon-Test. This indicates that a specific connection between the discriminative genes and the clock/cancer genes exists, which points to the relevance of this gene set in the clock-cancer context.

Furthermore, from the network analysis, we found that 20 out of the 45 discriminative genes were associated with clock genes (Figure 2C), 27 were found to be associated with cancer-related genes (Figure 2D) among which 18 intersect with the set of clock genes (Table S1). Discriminative genes associated with cancer pathways (Table 2) include IFNGR2 (involved in the Jak-STAT pathway) [53], PITX2 (TGF-β pathway), RFWD2 (p53 signalling), PPARγ (Wnt pathway). Moreover, several genes are also involved in the MAPK/RAS pathway: LOXL2, PPARD and CTSB are RAS target genes; Rab6 is a RAS family member; SPARC is also targeted by the RAS pathway and it was also shown to be a key modulator of extra-cellular matrix (ECM) remodelling, it affects cell proliferation and differentiation and it was recently reported to downregulate VEGF and thereby suppressing angiogenesis. In addition, we searched among our genes of interest for circadian properties by evaluating their expression profiles in published microarray data. We found 83% of the cancer related genes and 24% of the discriminative genes to show measurable oscillations in gene expression in several tissues and cell lines (Text S3).

Taken together, these discriminative genes are likely relevant for the analysis of possible clock malfunctions in a cancer model as a clear connection of these 45 genes to core-clock genes and to distinct cancer-related genes exists.

Functional coupling of Ras oncogene signalling and the circadian clock

The link of the circadian network to distinct cancer related pathways such as RAS/MAPK, Wnt and Jak/STAT led us to investigate the connection of clock regulation and signalling pathways in an experimental model. As a model pathway, we chose RAS/MAPK signalling, which is one of the most frequently altered signalling pathways in human cancer. KRAS mutations have a high prevalence in colorectal cancers. However, the colon cancer cells analysed exhibited a highly variable genetic background and thus introduced an extra level of complexity when used as models for functional analysis and for studying deregulation of the circadian clock. Therefore, we used a well-established in vitro epithelial model, in which cellular transformation is triggered by the RAS oncogene, that simulates carcinogenesis development, for investigating the influence of RAS transformation on the circadian clock [54]. Human HaCaT skin keratinocytes (and their derivatives represent the different steps of malignant epithelial conversion from the immortal state (HaCaT), to benign (HaCaT I7, class I tumours), advanced (HaCaT II4, class II tumours) and metastatic states (HaCaT A5RT3). The different HaCaT lines were lentivirally transduced with the Bmal1-driven luciferase reporter and the oscillation dynamics was monitored for 5 days. Non-transformed human keratinocytes exhibited an average period of τ = 23.4±0.4 hours (mean ± SEM, n = 5) with a strong oscillation that persisted for several days (Figure 3A). HaCaT I7 and HaCaT II4 periods were not significantly different compared to the immortalized cell line HaCaT (τ = 23.18±0.19 hours, n = 5 and τ = 22.96±0.3 hours, n = 3, respectively). However, the metastatic cell line HaCaT A5RT3 showed a significantly longer period of τ = 24.93±0.2 hours (n = 5), than that observed in normal HaCaT cells (p<0.05, Student's t-test) (Figure 3B). Moreover, HaCaT A5RT3 showed a significantly delayed phase of approximately 1.5 hours (17.9±0.15 hours) compared to HaCaT cells (16.37±0.21 hours; n = 5, p<0.05, Student's t-test). HaCaT II4 showed the opposite effect, having an advanced phase of approximately 0.85 hour significantly earlier than in normal keratinocytes (15.52±0.35 hours; n = 3 p<0.05, Student's t-test) (Figure 3C). A larger amplitude was observed in both HaCaT I7 and II4 (1.49±0.04, and 1.52±0.05, respectively) compared to normal keratinocytes (1.36±0.03) or HaCaT A5RT3 (1.33±0.02) (Figure 3A). To rule out possible secondary effects due to cell growth alterations caused by transfection with the BMAL1-promoter-driven luciferase (BLP) construct, proliferation of cells was measured by determining increased conductance of monolayers over time, using XCelligence technology (Roche) (unpublished observations).

To investigate if the observed HaCaT A5RT3 phenotype is indeed due to a different phase of the oscillator, we carried out a temperature entrainment assay with metastatic HaCaT A5RT3 and immortal HaCaT cells. This technique has the advantage that the cells do not respond to a single pulse of an external signal, but are entrained to follow an environmental cue (temperature) along a period of time. Bioluminescence measurement for at least 3 days following the entrainment revealed that HRAS transformed human keratinocytes A5RT3 indeed have a different circadian phenotype in comparison to normal keratinocytes (Figure 3D). The first peak after release of the cells to a constant temperature was used to determine the phase of entrainment. H-Ras-transformed HaCaT A5RT3 cells showed an advanced phase, of approximately 6 hours (Figure 3D). These data demonstrated that Ras transformation can induce a phase shift in the in vitro model system which is consistent with the observed period and phase perturbations in the circadian phenotype of A5RT3 cells.

To validate these results in an inducible RAS-system, we took the rat fibroblast cell line, 208F, and two derivative clones: IR2 and IR4. 208F cells are preneoplastic rat immortal fibroblasts. IR2 and IR4 cell lines are clones obtained by stable transfection with the H-Ras oncogene (G12V), which is under the lac regulatory control harbouring a Ras IPTG-inducible promoter [55]. All cell lines exhibited oscillations with a period of approximately 24 hours. The effect of the overexpression of H-Ras in rat fibroblasts was a clear phase shift in IPTG-treated cells (Figure 4). The levels of RAS and phosphorylated ERK were analysed via Western Blot and are depicted in Figure S6. The results obtained with the HaCaT system and the rat fibroblast cell lines could be reproduced with colorectal cancer cells HKe3 and are presented in Figure 5. The HKe3 cells are derived from HCT116 colorectal cancer cell lines, but KRAS has been disrupted by genetic recombination (HKe3; [56]). HKe3 cells exhibit a period of τ = 25.3±0.59 hours. Furthermore, we used the Hke3-clone 8, in which we introduced a conditional KRASV12 oncogene. In the absence of KRAS induction Hke3 clone 8 exhibits a period of (τ = 25.1±0.3 hours), similar to HKe3. Most interestingly, we observe a clear period phenotype upon activation of the KRAS oncogene (τ = 37.9±0.96 hours). This resembles our previous observations in HaCaT and IR2 cells.

With this data we could confirm our previous results and clearly showed that upon induction of RAS a larger period phenotype could be measured.

Clock genes are differentially expressed in human keratinocytes

To investigate the clock specific changes at the gene expression level in human keratinocytes, we performed real-time PCR for the clock genes Cry1, Bmal1, Per2, Rev-Erbα and Clock. The measurement was started 24 hours after synchronization to avoid influence of immediate early gene response that may not reflect an accurate effect of the oscillator. We measured the expression of the five clock genes in HaCaT and HaCaT A5RT3 cells over the course of 24 hours and all showed clear oscillations (Figure 6A).

We found marked differences in the levels of Bmal1, Per2, Cry1 and Clock mRNAs, while Per1 and Rev-Erbα mRNAs were largely similar in both cell lines. Consistent with the live-cell oscillation dynamics (Figure 3), Bmal1 transcript levels oscillate with a phase about 6 hours shorter in HaCaT cells than in HaCaT A5RT3 cells. While Per2 mRNA levels are strongly reduced in H-Ras transformed HaCaT cells, Cry1 and Clock gene expression is markedly increased in H-Ras transformed keratinocytes compared to the HaCaT cells, where Cry1 gene expression remained at low levels.

Overall, gene expression in both cell lines is divergent, especially for Cry1, where mRNA levels increased and for Per2 with decreased levels, in H-Ras transformed keratinocytes. This indicates that oncogenic RAS is likely to directly impinge onto the regulation of the circadian clock, however by yet unknown means.

Perturbations in BMAL activity lead to changes in the period

To unravel potential underlying mechanisms of RAS-mediated clock alterations we used our previously developed mathematical model for the mammalian circadian clock and carried out a control coefficient analysis over all model parameters and analysed the effects on period and magnitude [11]. We filtered for parameters for which a perturbation could induce antipodal changes in the magnitude of Per and Cry and at the same time an increase in the period. By perturbing a parameter involved in the transcription of the Cry gene (kt2, [11]) such an effect could be simulated (Figure 5B). Combinations of several parameters could cause a similar effect. However, in our mathematical model-system, the parameter described was the only which could, on its own, cause the observed phenotype. Interestingly, in our model the parameter kt2 (Cry transcription equation [11]) regulates BMAL1/CLOCK-mediated transcription on Cry proposing one possible way of how RAS might interfere with the circadian rhythms: weakening BMAL1's role as a transcriptional activator (perturbation of transcription activation for all clock genes [11]). Figure 7 predicts a gene expression profile in agreement with the experimental data suggesting that RAS activation perturbs the clock possibly via modulating BMAL1/CLOCK transactivation activity.

RAS/MAPK pathway activation is responsible for the phenotype

To test our hypothesis regarding the regulatory effect of RAS on the circadian clock, we first investigated the potential influence of the MAPK/RAS signalling pathway on the circadian phenotype in HaCaT cells. In our model the activation of RAS/MAPK signalling (60% reduction of the parameter which regulates BMAL1-mediated transcription, for each gene) predicts an increase of the period (τ = 24.1 hours), while inhibition of the RAS/MAPK pathway (60% increase in the parameter which regulates BMAL1-mediated transcription, for each gene) led to a shorter period phenotype (τ = 21.4 hours), as shown by the in silico expression profiles of Bmal1 (Figure 8A). To test the predicted effect of reducing RAS/MAPK-mediated signalling on the circadian phenotype experimentally, we treated synchronized HaCaT cells with the MEK inhibitor U0126 (Figure 8B, C). Indeed, U0126-treated cells showed a shorter period compared to vehicle-treated cells (Figure 8D). In HaCaT A5RT3, the MAPK pathway seems to impinge on the length of the period. This can be seen upon the comparison of HaCaT and A5RT3 and now much stronger when we inhibit the MAPK pathway in A5RT3 (Figure S7).

Together, our theoretical and experimental data indicate that the activity of the RAS/MAPK modulates the circadian period (and thereby also the entrained phase), possibly by influencing the transcriptional activity of the CLOCK/BMAL1.

Discussion

Perturbations of the circadian clock have been described in different pathological conditions. In model systems, genetic manipulations and read-out of the circadian clock were found associated with distinct phenotypes such as insulin resistance and obesity [57], premature aging [58], increased tumour development [28] and altered stem cell homeostasis [59]. This diversity of phenotypes reflects the fact that a large number of genes controlling crucial cellular processes such as the cell cycle, DNA repair and metabolic processes are known to be regulated in a circadian manner. Despite these clear functional associations between clock and disease, there are few mechanistic studies investigating possible molecular processes involved in the interconnection between cancer and clock in a systematic way. In this manuscript, we set out to extensively analyse molecular alterations, in colon and skin cancer, which could bring us closer to understand the processes and mechanism involved in how the deregulation of the circadian clock takes place in carcinogenesis.

De-regulation of clock at multiple levels: Gene expression, signalling transduction

In this study, we showed that the clock is perturbed differentially within the same cancer type and observed a rich variety of clock phenotypes. To explain these effects, we correlated the gene expression profiles with the clock phenotype of the cells. This analysis revealed a set of genes discriminating between strong and weak oscillator, thus functioning as a clock phenotype predictor. In fact, this set of genes also allowed the correct classification of the osteosarcoma test cell line U2OS. Moreover, previous data [49] showed that the siRNA dependent-knockdown of the majority of these 45 genes confers a circadian phenotype in U2OS cell lines, underlining their potential importance as regulators of the circadian system. To what extend this set of genes exhibits robustness as a clock phenotype predictor beyond our experimental setup is currently unknown. Of note, Caco2 cells have been described by others as a good oscillators, yet the experimental conditions were different [60]. Thus, nature of the predictive gene set might undergo some alterations when alternative conditions are used.

The microarray data revealed important clock genes such as Nono (a known splicing factor and recently described as connecting the clock to the cell cycle via PER [36]), Rac (involved in proliferation) and CkIIα to be highly expressed in all cell lines. At the same time, the tumour suppressor gene Per3 appears weakly expressed in all cell lines as expected. In addition, we identified a characteristic differential expression of several genes encoding epigenetic regulators, signalling molecules and transcriptional regulators not previously connected to circadian rhythm. CBX7 is an essential component of the polycomb repressive complex 1 (PRC1) involved in the control of histone methylation at tumour suppressor loci such as p16 [61]. CHD, encodes a chromodomain regulator of chromatin remodelling and was lost in approximately 50% of colorectal cancers [62]. The YY1 interacting chromatin remodelling complex factor INO80 [63] and the ubiquitous epigenetic regulator and insulator protein CTCF, are associated with the weak oscillators. CTCF controls long-range chromatin interactions and functions in the establishment and maintenance of epigenetic signatures. This renders CTCF a potentially important factor also for controlling circadian genes [64]–[66]. A link between the circadian clock, energy metabolism and epigenetic (re)programming has been established earlier [67]–[69]. Our new observations indicate that epigenetic events involving DNA methylation, histone modification and chromatin remodelling might also induce differential circadian oscillation in tumour cells.

Differential gene expression between strong and weak oscillators was also seen for the signalling molecules LOXL2, SPARC, CTSB, IFNGR, WASF3, GNG11 and the metabolism-associated PPARD gene [70]. High expression of these genes was strongly associated with the strong oscillator phenotype, but low expression with the weak oscillators. In contrast, high expression of FOXA1 and PDHX were associated with the weak oscillator phenotype. The RAS pathway target genes LOXL2, SPARC and CTSB exert functions in the remodelling of the extracellular matrix, thereby favouring invasion and metastasis of tumour cells [71]–[73]. WASF3 [74] and FOXA1 [75] are downstream targets of the anti-apoptotic PI3K signalling pathway and have been shown to control actin polymerisation and invasion as well as differentiation of secretory gut epithelial cells, respectively. GNG11, a member of the heterotrimeric G-protein family, is involved in senescence induction via environmental stimuli [76] and has been reported recently to be deregulated in TGFIIR knock-out epithelial cells capable of increased metastasis [77]. PPARD and PPARB expression is upregulated in colorectal cancer [78] and the gene was found activated by the K-RAS pathways in rat intestinal epithelial cells [79]. The protein mediates activation of PI3K signalling via PTEN deregulation and can enhance anti-apoptotic signalling and cell survival [70]. Taken together, these observations show that targets of RAS/MAPK-dependent signalling are associated with a certain oscillator phenotype. With the exception of LOXL2, the RAS target genes that play a role in cancer invasion and metastasis are downregulated in the weak oscillators, thus being associated with the bad prognosis phenotype. In addition, PI3K signalling might play a functional role in clock deregulation as indicated by differential expression of WASF3 and FOXA1. Recently it was shown that GSK3β is able to phosphorylate and destabilize CLOCK [80]. Overexpression of oncogenic RAS and subsequent activation of PI3K/Akt resulted in GSK3β inhibition and in an increased stabilization of CLOCK. Comparing clock gene expression between immortal and RAS transformed HaCaT cell, we also detected increased levels of Clock mRNA. To what extent the elevated Clock mRNA levels are due to RAS/MAPK or RAS/PI3K signalling in our cells needs to be investigated.

Furthermore, we correlated the set of discriminative genes found to an extended network of the mammalian circadian clock. From our assembled network it is evident that the core-clock genes directly regulate a set of 33 clock-associated genes involved in the Wnt, the TGF-β, Jak-STAT and MAPK signalling pathway. We also found many of the clock-associated genes to be involved in important biological processes such as cell cycle and repair mechanisms, immune defense, metabolism and the xenobiotics detoxification mechanisms. All these processes and pathways are often found to be deregulated in cancer. The assembled network established in this study provides a novel extension of the circadian system to output clock-associated genes which can potentially be relevant in terms of drug targets or even in the prognosis of cancer. We also suggest that a feedback from the clock-associated genes to the core-clock exists. This might explain the fact that the deregulation of cancer-related pathways can lead to perturbations of the circadian clock leading to the observed phenotype of disturbed sleep patterns as observed in cancer patients. The presented interconnections could as well help to understand why disruption of normal rhythms caused by perturbations of the clock can increase cancer incidence among night shift workers.

RAS/MAPK can directly alter the circadian clock: Novel insights from modelling

We present in the manuscript several pieces of evidence that hint to effects of the MAPK/RAS pathway in perturbation of the clock phenotype. From the bioinformatics analysis of the extended circadian clock network, to the association of genes in the discriminative set found from the microarray analysis, there is support pointing to the MAPK/RAS pathway as an important regulator in the circadian system.

We extended our study of the connection between oncogenic pathways and the circadian clock with a model system which allowed us - in a controlled way - to explore the role of MAPK/RAS for the circadian oscillator. The starting point for the model was the comparison of the expression levels of 6 clock genes, which revealed downregulation of Per2 and increased expression of Cry1 and to a lesser extend of Clock in the RAS-transformed A5RT3 cells. Downregulation of Per2 is in line with previous reports for colorectal and skin cancer [81], [82] indicating that Per2 deregulation is correlated with malignant transformation in different systems. Interestingly, upregulation of Clock and CKIε expression was recently found in colorectal cancer tissues as compared to the adjacent normal mucosa and in these tissues the differential expression of Bmal1 and Per1 was correlated with metastasis [83]. We used a mathematical model for the mammalian circadian clock which we previously had implemented to develop a hypothesis on the mechanism of RAS/MAPK mediated alterations in the circadian clock. We hypothesised that RAS perturbation alone is sufficient to induce the observed altered phenotype, which corresponds to a deregulation in the circadian clock. This assumption could be verified experimentally. Although our mathematical model was optimized for a clock working under a normal scenario of non-pathogenic conditions, we were able to simulate the cancer phenotype observed experimentally, by perturbing a single parameter. The perturbed parameter is involved in the transcriptional regulation of clock genes, which is carried out by BMAL1 as a transcriptional activator. Published data on cells from the chick pineal gland demonstrated that activated ERK is able to phosphorylate BMAL1. This was confirmed in vitro, showing that p-ERK phosphorylates BMAL1 at multiple sites thereby reducing the transactivation capacity of BMAL1:CLOCK [84]. Also in other systems (not related to cancer) indirect correlations between the RAS/MAPK pathway and the circadian clock were reported [85]–[87]. Our observations in the HaCaT and A5RT3 cell lines support this data and suggest that RAS/MAPK activation perturbs the clock by interfering with clock-gene expression. Upon interference with RAS/MAPK signalling in both cells using a MEK inhibitor, we observed the opposite period phenotype as upon RAS overexpression. These results were consistent in vitro and in silico. Furthermore, using either immortal rat 208F cell lines harbouring an inducible H-Ras oncogene, or colon epithelial cells Hke3 cells harbouring an inducible K-Ras oncogene, we obtained similar results upon transient induction of oncogenic Ras. The deregulation of the circadian clock in the HKe3 cells following Ras induction suggests a stronger effect of oncogenic Ras in the colonic epithelial cells as compared to the keratinocytes, although in both cells a period lengthening phenotype can be observed. This indicates that while an influence of RAS/MAPK onto the circadian clock is visible in different cell types, we cannot easily compare different tissue types with each other. While in colorectal cancer cells the circadian clock gets strongly disrupted upon Ras induction (HKe3 clone 8) and gain of metastatic potential (SW620 cells), the metastatic keratinocytes A5RT3 still oscillate. Therefore, the RAS-mediated deregulation of clock rhythm is not necessarily connected to the metastatic potential of the cells, but can rather be considered a specific response of cells towards RAS oncogene-mediated tumourigenesis.

Together, our findings shed light into the mechanism by which the circadian clock is perturbed in carcinogenesis via a single oncogenic pathway, the RAS/MAPK pathway.

Materials and Methods

Text mining approach

Information about protein-protein and protein-DNA interactions were taken from our freely available GeneView repository of facts extracted from PubMed abstracts or full texts [50]. In a nutshell, protein names were recognized using the state-of-the-art tool GNAT [88] and were mapped to Entrez-Gene identifiers. All co-occurring protein citations in a sentence were classified using a custom support vector machine (SVM) kernel [89]. This SVM-kernel achieved very good results in a comprehensive evaluation of nine machine learning kernels for interaction extraction from text [90], with F1-measures ranging between 56.2 and 76.1, depending on target corpus. To account for species specificity we mapped mammalian gene identifiers to HomoloGene cluster [91]. Sentences containing potentially novel interactions were ranked by confidence of the SVM-kernel (distance to the hyperplane). In total, we extracted 15,197,637 protein pairs in 1,372,877 different PubMed articles. Of those we classified 3,921,267 (25%) as interacting. The general quality of this repository was estimated by random evaluation of 181 predicted interactions between circadian elements. From these a domain expert classified the majority of 135 (74.6%) as correct and 46 (25.4%) as incorrect. We also added to the network all relevant interactions described in the STRING database, which provides a regularly updated high quality compendium of protein-protein interactions (PPI) from several PPI-repositories, such as Kegg, Bind, Mint, and IntAct.

Microarray experiments

a) Experimental setup

Cells were all trypsinized and replated on a Friday (8–9:30 am), medium was completely removed and fresh serum was added to make sure that all cells get the mitogenic stimulus at the same time. Cells were plated in 10-cm dishes at a density of 5×10⁵ cells for the 48 h time point to guarantee for constant growth over time and prevent density-dependent growth arrest. On Monday 8 a.m., cells were refed with fresh medium and serum and in cases of K-Ras-induction, the relevant substance was added. Each following day, the cells were collected after 24 hours (microarray time 0) or two days later (microarray time 48).

b) Analysis

Gene expression data (Affymetrix Human Gene 1.0 ST microarrays) were analyzed using the statistical programming environment R. The Affymetrix Human Gene 1.0 ST microarray contains 32,321 probesets. 21,307 probesets are associated to 19,942 distinct human genes. Samples are grouped by type (circadian/non-circadian) and time point (0 h/48 h). Quality of all arrays was determined using array QualityMetrics [92]. The arrays were subsequently normalized using robust multichip average (RMA) [93]. Differential Expression was determined by moderated t-test using R-package limma [94]. The moderated t-test was described as a robust method especially useful for small sample sizes and performed better or at least on-a-par with other commonly used statistically methods [95]–[98].

Detailed information is provided in Text S1. We found no significant difference between the time points 0 h/48 h for the same cell line (lowest Benjamin-Hochberg corrected p-value, 0.99). Consequently, these can be considered as biological replicates.

To find a set of best-discriminating genes, we followed a leave-one-out cross validation strategy. For each cell line we excluded both samples (0 and 48 hour time points) once and determined the 100 most significant probe sets, between the remaining strong/weak oscillator phenotypes. We used moderated t-test to select the most discriminative genes ranked by confidence (p-value). Following this procedure we obtained 6 lists. The quality of each list was evaluated as follows: 1) we generated a heatmap using expression profiles for all six cell lines using the current list of 100 discriminating probes. 2) For heatmaps which correctly cluster the excluded cell line into the expected group, results are retained. We found this procedure to be more robust than global comparison of all arrays, because outliers were filtered out through the intersection. This procedure was repeated for all six different cell lines. 3) Finally, we generated a discriminative list of 45 genes by taking the intersection between all useful discriminative gene-lists (Text S1). The microarray dataset was submitted to GEO with the ref: GSE46549 and will be released upon publication.

c) Validation

For validation purpose we normalized the cell lines using frozen Robust Multi-array Analysis fRMA, which has better characteristics when working incrementally with new batches of arrays [99]. To further reduce batch effects we replaced the cell lines CaCo2 and HCT116 with respective biological replicates. We then generated a heat map for each cell line using the previously defined 45 genes.

Modelling data

All simulations were done using a mathematical model for the mammalian circadian clock which we previously developed [11]. The plots in Figures 6B and 7A were produced using XPPAUT, version 5.85 (http://www.math.pitt.edu/~bard/xpp/xpp.html) and the XPPAUT subsystem AUTO.

Circadian data

Periods and amplitudes were estimated by fitting the cosine wave function via the ChronoStar analysis software [100]. For visualization, data were smoothened by a 4 hours-running average. Different basal luciferase levels from raw data were included by the fold change in luciferase activity relative to GFP controls for de-trended and smoothed data.

Cell culture

HCT116 human colorectal carcinoma cells were maintained in McCoy's 5A culture medium. HaCaT keratinocytes and its derivatives HaCaT I7, II4 and A5RT3, the osteosarcoma U2OS, the human colon carcinoma cell lines HT-29, RKO, SW480, LIM1512, CaCo2 and the rat fibroblast cells 208F and IR2 were cultured in Dulbecco's modified medium (Gibco). HKe3 cells were cultivated in D10 with G418, HKe3 clone 8 cells in D10 with G418, Zeocin and puromycin (all Sigma). For K-Ras induction mifepristone was added at 2.5 µM. All culture media were supplemented with 10% fetal calf serum and 1% penicillin-streptomycin (Biochrom). Stable-transduced cell populations were selected and maintained in medium containing puromycin (Sigma). For live-cell bioluminescence recording (evaluation/analysis/monitoring), cells were maintained in phenol red-free DMEM (Gibco) containing 10% fetal calf serum, 1% penicillin-steptomycin and 0.1 mM D-Luciferin (PJK). Cell morphology and density were controlled by light microscopy. All cells were incubated at 37°C with 5% CO₂ atmosphere.

Crypt cultures

Crypt organoid cultures from Per2 transgenic mice [101] for bioluminescence imaging were established as previously described [102]. Briefly, crypts from the proximal half of the small intestine were embedded in Matrigel (BD) at a density of app. 75 crypts/25 µl. For LumiCycle analysis 3 µl of 10 mM D-Luciferin were added per 25 µl of Matrigel with crypts and 10 drops of 28 µl volume each were plated per dish. Crypt organoids were cultured in Advanced DMEM/F12 without phenol red (Gibco/Life Technologies) supplemented with 0.1 mM D-Luciferin as well as additives and growth factors.

Lentivirus production

Lentiviral elements containing a BMAL1-promoter-driven luciferase (BLP) were generated as previously described [103]. HEK293T cells were seeded in 175 cm² culture flasks and co-transfected with 12.5 µg packaging plasmid psPAX, 7.5 µg envelope plasmid pMD2G and 17.5 µg BMAL1-promoter (BLP) luciferase expression plasmid using the CalPhos mammalian transfection kit (Clontech) according to the manufacturer's instruction. To harvest the lentiviral particles, the supernatant was centrifuged at 4100× g for 15 min to remove cell debris and passed through a 45 µm filter. The lentiviral particles were concentrated to obtain a retentate by passing the supernatant through an Amicon Ultra centrifugal filter device (Millipore). The obtained retentate, containing 100× concentrated recombinant lentivirus, was stored at −80°C.

Transduction with lentiviral vectors

All cell lines were seeded in 12-well plates and incubated at standard cell culture conditions until they reached 80% of confluency. Culture medium was removed and cells were washed once with 1× PBS. Cells were covered with DMEM containing 8 µg/ml protamin sulphate (Sigma-Aldrich) and 300 µl of supernatant or 20× retentate (diluted in DMEM) of the corresponding lentivirus (BLP: BMAL1-promoter luciferase reporter). The next day, the medium was replaced with selection medium (DMEM supplemented 10% fetal calf serum, 1% penicillin-streptomycin and puromycin,) was added to obtain stable transduced cells and incubated at 37°C with 5% CO₂ atmosphere.

Synchronization and measurement of circadian rhythms

To guarantee that confluency is reached only by the end of the measurements, each cell line was specifically tested both prior to the live-cell bioluminescence experiments and then again following transduction with the lentivirus harbouring the Bmal1 reporter.

For bioluminescence measurement, 1–5×10⁵ cells were plated in 35 mm dishes or 1×10⁴ cells in 96-well plates and synchronized by a single pulse of 1 µM dexamethasone (Sigma) for 1 hour. Next, cells were washed once with 1× PBS and cultured with phenol-red-free DMEM supplemented with 0.1 mM D-Luciferin and puromycin or hygromycin. Bmal1-promoter-(BLP)-reporter activity was measured, using a photomultiplier tube (PMT)-based system (Hamamatsu Photonics) in a LumiCycle instrument (Actimetrics) or in a plate luminometer (TopCount equipment, Perkin Elmer). Raw luminescence data were de-trended by dividing luminescence counts by the 24 hours running average. Phase, period and amplitude were analysed using the ChronoStar analysis software [100].

Temperature entrainment

1–5×10⁵ cells were plated in 35 mm dishes and synchronized with 1 µM dexamethasone for 1 hour. Subsequently, cells were washed once with 1× PBS and cultured with phenol red-free DMEM supplemented with 0.1 mM luciferin and puromycin. Cells were subsequently entrained to temperature cycles (5 d: 10 hours 37°C constant; 2 hours 37°C ramp; 33°C constant; 2 hours 33°C ramp; 37°C constant). BLP-reporter activity was measured in the PMT-based system. Detrending of the time series was performed by dividing the luminescence counts by 24 hours running averages. Phase, period and amplitude were analysed using the ChronoStar analysis software [100].

RNA isolation

Cells were seeded in 35 mm dishes with a density of 1–5×10⁵ cells/dish and synchronized for 1 hour with dexamethasone. 24 hours later, cells were collected every three hours and total RNA was extracted. The measurement was started 24 hours after synchronization to avoid an influence of the immediate early gene response that may not reflect an accurate result of the oscillator. Total RNA was isolated using the RNeasy extraction kit (Qiagen), including DNAse digestion, according to manufacturer's instructions. Cells were lysed with 350 µl RLT buffer (with 1% β-mercaptoethanol) and lysate was homogenized with a syringe. Total RNA concentration was determined by OD₂₆₀ using an Ultrospec 3000 UV spectrophotometer and stored at −80°C.

cDNA synthesis and quantitative real-time PCR

Total RNA was isolated using the RNeasy extraction kit (Qiagen), according to manufacturer's instructions. Total RNA concentration was determined by OD₂₆₀ using an UV spectrophotometer Ultrospec 3000. 1 µg of total RNA was reverse-transcribed to cDNA with 2 U/µl RevertAid M-MuLV reverse transcriptase and 500 µM random 15-mer primers (Fermentas).

For quantitative SYBR-green based real-time PCR we used 20 ng of cDNA and the commercial primers from QuantiTect primer assays for Gapdh, Bmal1, Per1, Per2, Rorα, Cry1, Cry2 and Clock (Qiagen). For each primer, a non-template control (NTC) was always included. The reaction was performed in an ABI Prism 7000 SDS thermocycler (Applied Biosystems). The obtained data were analysed using the ABI Prism 7000 System SDS software (Applied Biosystems). Each sample was measured in triplicates. The expression levels were normalized to those of Gapdh (ΔCT) and calibrated to the value at time 0 of the control cell line (ΔΔCT). The relative quantification was performed using the 2_−ΔΔCT method.

Western blot

Total protein lysates for immunoblotting were obtained using SDS lysisbuffer (10 mM Tris-HCl [pH 7.5], 1% SDS, and 2 mM EDTA) supplemented with protease inhibitors (1 mM PMSF, 50 mM NaF, 50 µg/ml leupeptin, 1 mM orthovanadate, 4 µg/ml aprotinin). Protein concentration was determined by the Pierce BCA Protein Assay Kit (Thermo Scientific #23227F). 30 µg of protein were resolved by SDS-PAGE and transferred to a nitrocellulose membrane. Antibodies used were phospho-ERK1/2 (Thr202/Y204; Cell Signaling #4370), Ras (Thermo Scientific #1862335) and GAPDH (Ambion Am4300). Secondary antibodies were IRDye 800CW goat anti-rabbit IgG (LI-COR Bioscience 926-32211), IRDye 800CW goat anti-mouse IgG (LI-COR Bioscience 92632210) and IRDye 680 goat anti-mouse IgG (LI-COR Bioscience 926-322220). Blots were developed using the LI-COR Odyssey system.

MEK inhibition

1–5×10⁵ cells were plated in 35 mm dishes and cultured overnight. The next day the cells were synchronized by a single pulse of 1 µM dexamethasone (Sigma) for 1 hour after which a MEK inhibitor (UO126 final concentrations between 5 and 20 µM, Promega) was added to the culture and incubated for 7 days. Circadian rhythmicity was determined as described above.

Supporting Information

Attachment 1

Attachment 2

Attachment 3

Attachment 4

Attachment 5

Attachment 6

Attachment 7

Attachment 8

Attachment 9

Attachment 10

Attachment 11

Attachment 12

Attachment 13