1. Human Biology and Medicine
Download icon

Transcriptomic analyses reveal rhythmic and CLOCK-driven pathways in human skeletal muscle

  1. Laurent Perrin
  2. Ursula Loizides-Mangold
  3. Stéphanie Chanon
  4. Cédric Gobet
  5. Nicolas Hulo
  6. Laura Isenegger
  7. Benjamin D Weger
  8. Eugenia Migliavacca
  9. Aline Charpagne
  10. James A Betts
  11. Jean-Philippe Walhin
  12. Iain Templeman
  13. Keith Stokes
  14. Dylan Thompson
  15. Kostas Tsintzas
  16. Maud Robert
  17. Cedric Howald
  18. Howard Riezman
  19. Jerome N Feige
  20. Leonidas G Karagounis
  21. Jonathan D Johnston
  22. Emmanouil T Dermitzakis
  23. Frédéric Gachon
  24. Etienne Lefai
  25. Charna Dibner  Is a corresponding author
  1. University Hospital of Geneva, Switzerland
  2. University of Geneva, Switzerland
  3. Institute of Genetics and Genomics of Geneva, Switzerland
  4. CarMeN Laboratory, INSERM U1060, France
  5. Nestlé Institute of Health Sciences, Switzerland
  6. Ecole Polytechnique Fédérale de Lausanne, Switzerland
  7. University of Bath, United Kingdom
  8. University of Nottingham, United Kingdom
  9. Edouard Herriot University Hospital, France
  10. NCCR Chemical Biology, University of Geneva, Switzerland
  11. University of St Mark and St John, United Kingdom
  12. Nestlé Research Centre, Switzerland
  13. University of Surrey, United Kingdom
Research Article
  • Cited 3
  • Views 1,645
  • Annotations
Cite as: eLife 2018;7:e34114 doi: 10.7554/eLife.34114

Abstract

Circadian regulation of transcriptional processes has a broad impact on cell metabolism. Here, we compared the diurnal transcriptome of human skeletal muscle conducted on serial muscle biopsies in vivo with profiles of human skeletal myotubes synchronized in vitro. More extensive rhythmic transcription was observed in human skeletal muscle compared to in vitro cell culture as a large part of the in vivo mRNA rhythmicity was lost in vitro. siRNA-mediated clock disruption in primary myotubes significantly affected the expression of ~8% of all genes, with impact on glucose homeostasis and lipid metabolism. Genes involved in GLUT4 expression, translocation and recycling were negatively affected, whereas lipid metabolic genes were altered to promote activation of lipid utilization. Moreover, basal and insulin-stimulated glucose uptake were significantly reduced upon CLOCK depletion. Our findings suggest an essential role for the circadian coordination of skeletal muscle glucose homeostasis and lipid metabolism in humans.

https://doi.org/10.7554/eLife.34114.001

Introduction

Circadian rhythms are daily cycles of most bodily processes driven by a system of intrinsic biological clocks organized in a complex hierarchical manner. This mechanism ensures a temporal coordination of physiology and behavior with a near 24 hr period of rest-activity and feeding-fasting cycles, thus providing the organism with an evolutionary conserved advantage (Albrecht, 2012; Spoelstra et al., 2016). In mammals, the circadian system is driven by a central pacemaker, situated in the paired suprachiasmatic nuclei (SCN) of the hypothalamus, and by secondary oscillators located in peripheral organs. The SCN clock is readjusted on a daily basis, mainly by light inputs coming from the retina. In turn, the central pacemaker orchestrates peripheral clocks through a combination of neuronal, endocrine, and metabolic signaling pathways (Saini et al., 2015). As a result, metabolic processes in the liver, skeletal muscle, and other organs are subject to daily oscillations (Asher and Sassone-Corsi, 2015) with the SCN keeping these rhythms in appropriate synchrony with each other.

Skeletal muscle is responsible for ~70% of glucose uptake resulting from ingested carbohydrates (DeFronzo et al., 1981; Gachon et al., 2017), and perturbations in glucose sensing and metabolism in this organ are strongly associated with insulin resistance in type 2 diabetes (T2D) (Muoio and Newgard, 2008). In rodents, it has been established that the skeletal muscle clock plays an essential role in maintaining proper metabolic homeostasis, with skeletal muscle pathologies stemming from clock disruption via deletion of the core clock component Bmal1 (Andrews et al., 2010; Harfmann et al., 2015). Disruption of muscle insulin sensitivity by modulating glucose uptake, with a reduction in Glut4 mRNA and protein levels, has been reported in two muscle-specific Bmal1 knockout (KO) models (Dyar et al., 2014; Harfmann et al., 2016). In humans, diurnal variations in glucose tolerance have been described (Kalsbeek et al., 2014), although the molecular mechanism responsible for such variations remains largely unexplored. Feeding-fasting cycles accompanying rest-activity rhythms represent major timing cues in the synchronization of peripheral clocks, including skeletal muscle oscillators (Dibner and Schibler, 2015; Wehrens et al., 2017). Although several studies have reported that in murine models 3.4% to 16% of the skeletal muscle transcriptome is expressed in a circadian manner (McCarthy et al., 2007; Miller et al., 2007; Dyar et al., 2014; Zhang et al., 2014; Hodge et al., 2015), it is unclear to what extent the muscle circadian transcriptome is regulated by feeding-fasting rhythms and additional central synchronizers, and how local muscle clocks contribute in generating these transcript oscillations. Cell-autonomous circadian clocks operating in human primary skeletal myotubes (hSKM) established from muscle biopsies have been recently characterized (Perrin et al., 2015; Hansen et al., 2016; Loizides-Mangold et al., 2016). Importantly, proper function of these cellular oscillators is indispensable for the normal secretion of interleukin 6 (IL-6), interleukin-8 (IL-8), the monocyte chemotactic protein 1 (MCP-1) and additional myokines, regulating skeletal muscle insulin sensitivity and inflammation (Perrin et al., 2015).

In order to dissect the impact of the endogenous circadian clock on skeletal muscle gene transcription from external factors and their reciprocal influence, we performed a genome-wide transcriptome analysis by high-throughput RNA sequencing (RNA-seq) in skeletal muscle biopsies collected from human subjects placed under a controlled laboratory routine, as well as in cultured hSKM synchronized in vitro in the presence of a functional or compromised clock. An important overlap between genes exhibiting rhythmic patterns in tissue biopsies and in synchronized hSKM was observed. Expression patterns of genes related to insulin response, myokine secretion, and lipid metabolism were strongly altered in the absence of a fully functional clock in vitro. These transcriptional changes had important functional outputs, with basal as well as insulin-stimulated glucose uptake and lipid metabolism being altered by perturbation of the oscillators operative in hSKM. Altogether, these results strongly suggest that cell-autonomous skeletal muscle clocks drive rhythmic gene expression and are indispensable for proper insulin response, lipid homeostasis, and myokine secretion by the skeletal muscle.

Results

Diurnal rhythms of gene expression in human skeletal muscle under controlled laboratory routine

To assess rhythms of gene expression in human skeletal muscle, RNA samples derived from vastus lateralis biopsies taken every 4 hr across 24 hr from 10 healthy volunteers were analyzed by total RNA-seq (see Supplementary file 1-table S1 for in vivo donor characteristics, n = 10). Sample collection was performed under controlled laboratory routine by implementing a protocol designed to minimize the effect of confounding factors (see Materials and methods and [Loizides-Mangold et al., 2017]). In total 13,377 genes were quantified at the exonic level (Figure 1—source data 1), of which 9211 genes were quantified at the intronic level as well. To identify genes with coordinated rhythmic expression, we used a mixed linear model with harmonic terms across the 10 individuals at the pre-mRNA (intronic signal) and mRNA (exonic signal) levels. This method allowed for the identification of 5748 rhythmic genes that were rhythmic at the pre-mRNA or/and mRNA level with a False Discovery Rate (FDR) of 5% (Figure 1A). When rhythmicity levels were further analyzed, it became apparent that 4792 genes showed rhythmic transcription at the intronic pre-mRNA level (Figure 1A, upper and middle left panels). However, from these rhythmic pre-mRNA transcripts, only 57% were also rhythmic at the mRNA level (R-I.R-E, upper panels Figure 1A), likely because of the longer half-life of their mRNA. Indeed, approximation of mRNA half-life by the exon/intron ratio showed that among these rhythmic pre-mRNA transcripts, those that are only rhythmic at the pre-mRNAs level (R-I) have a longer half-life compared to those that are rhythmic at the mRNA level (R-E and R-I.R-E, Figure 1B). The R-I.R-E group of genes was enriched in circadian clock genes and genes encoding RNA-binding proteins, whereas the R-I group was enriched for genes encoding proteins involved in mRNA translation, mitochondrial activity, TCA cycle, and lipid metabolism (Figure 1—source data 2). In parallel, around 10% of the quantified transcriptome (956 genes) were only rhythmic at the mRNA level (R-E, Figure 1A, lower panels), likely through post-transcriptional regulation and in particular mRNA degradation (LuckLück et al., 2014; Wang et al., 2018). This group was enriched in genes encoding proteins involved in ribosome biogenesis and protein transport (Figure 1—source data 2). Amplitude distribution suggested that among the genes that were rhythmic at the pre-mRNA level, those with higher amplitude of transcription had a greater probability to be rhythmic at the mRNA accumulation level (R-I.R-E, Figure 1C). Taken together, our results highlight the high rhythmicity of gene expression in human skeletal muscle, even under controlled laboratory routine. However, the number of genes being qualified as significantly rhythmic at the pre-mRNA and/or mRNA level was strongly dependent on the threshold level that was applied (Figure 1D).

Figure 1 with 1 supplement see all
Rhythmic gene expression in human skeletal muscle.

(A) Heat map showing genes rhythmic at the pre-mRNA and mRNA level (R-I.R-E: upper panel), at the pre-mRNA level only (R-I: middle panel), and at the mRNA level only (R-E: lower panel). Standardized relative gene expression is indicated in green (low) and red (high) and ordered by their respective phase. (B) mRNA half-life proxy by exon/intron ratio showing lower stability for genes with rhythmic mRNA (R–E) profiles. (C) Amplitude distribution of genes that are rhythmic only at the mRNA level (R-E, blue), the pre-mRNA level (R-I, red), or rhythmic for both (R-I.R-E). Genes with higher amplitude of transcription at the pre-mRNA level have a higher probability to be rhythmic at the mRNA level (R-I.R-E). (D) Number of genes in each group in relation to the -log10 BH corrected p-value; dashed line indicates threshold of 0.05. (E) Phase distribution at the pre-mRNA and mRNA level for the three groups described in (A). (F) Phase distribution for genes activated by acute muscle exercise (red), inflammation (blue), or both (green). (G) Temporal expression of core clock components, and (H) key muscle transcription factors. N = 10 human muscle biopsy donors. (I) Phase distribution of predicted rhythmic DNA motif activity.

https://doi.org/10.7554/eLife.34114.002

As previously reported in mice, rhythmic gene transcription was distributed into two phases of transcript accumulation (04:00 and 16:00, Figure 1E). The afternoon peak (16:00) was enriched in genes related to muscle contraction and mitochondrial activity (Figure 1—source data 2), whereas homologous genes in rodents were shared between the active but also the resting phase (McCarthy et al., 2007; Miller et al., 2007; Hodge et al., 2015). In contrast, among the genes highly activated in the middle of the night (04:00), many were associated with inflammation and immune response (Figure 1F).

Among the rhythmic genes, we observed high amplitude oscillations for the core clock genes ARNTL (BMAL1), NPAS2, CLOCK, PER2, PER3, CRY2, NR1D1 (REV-ERBα) and RORA, which were well synchronized among the 10 individuals (Figure 1G). In addition, several transcription factors associated with muscle metabolism and physiology like TFEB, a key regulator of lysosomal biogenesis and autophagy that also regulates glucose homeostasis and oxidative phosphorylation (Mansueto et al., 2017), KLF13, which plays a role in cardiac muscle cell development (LavalleeLavallée et al., 2006), and KLF15, an important transcriptional regulator of muscle lipid metabolism (Haldar et al., 2012), showed an oscillatory profile. Moreover, also PPARD, the most abundant PPAR isoform in skeletal muscle and master regulator of muscle mitochondrial function (Jordan et al., 2017), and MYOD1, the key regulator of myogenesis and direct target of BMAL/CLOCK (Andrews et al., 2010) were rhythmically expressed in human skeletal muscle, likely orchestrating the temporal muscle transcriptome (Figure 1H). In addition to these transcription factors, genes involved in the secretion of myokines, glucose homeostasis and lipid metabolism displayed rhythmic transcription (Figure 1—figure supplement 1A,B and C, respectively).

To gain more insight into the rhythmic transcriptional regulation observed in this dataset, we performed a DNA-binding motif enrichment analysis to identify those with rhythmic activity. As shown in Figure 1I, the 16:00 peak is enriched in MEF2 and MYOD1 motifs, in phase with MYOD1 expression, and both proteins synergize to activate gene expression (Taylor and Hughes, 2017). In parallel, the 4:00 peak is enriched in MYC and AP1 families of transcription factors, both downstream of the MAP kinase pathway activated by exercise (Aronson et al., 1997) or wound-induced inflammation (Aronson et al., 1998).

Cell-autonomous circadian clocks regulate functional gene expression in hSKM

To assess the impact of cell-autonomous circadian clocks operative in hSKM on skeletal muscle gene transcription and function, we used our previously developed model of efficient siCLOCK-mediated clock disruption (Perrin et al., 2015). RNA-seq was conducted on siCLOCK-transfected hSKM obtained from two male donors matched for age and BMI (Supplementary file 1-table S1, donors 1 and 2, in vitro part). Human primary myoblasts were cultured and differentiated in vitro into myotubes, transfected with siRNA, synchronized in vitro with a forskolin pulse, with subsequent sample collection every 2 hr during 48 hr for RNA-seq analysis (Figure 2—figure supplement 1 and Materials and methods). CLOCK expression was reduced by at least 80% upon siCLOCK depletion, as assessed by RNA-seq and quantitative real-time PCR (Figure 2A and Figure 2—figure supplement 2A). In parallel, bioluminescence profiles derived from hSKM transduced with a lentiviral Bmal1-luciferase vector were monitored as described previously (Perrin et al., 2015). As expected, the circadian amplitude of Bmal1-luciferase bioluminescence reporter profile became dampened upon siCLOCK compared to siControl and non-transfected counterparts (Figure 2—figure supplement 2B).

Figure 2 with 2 supplements see all
Disruption of the circadian oscillator impacts on functional gene expression hSKM.

Differential gene expression analysis between hSKM bearing a disrupted (siCLOCK) or intact (siControl) circadian clock. Comparison of the median gene expression at all analyzed circadian time points between the two groups. A total of 16,776 genes were detected by RNA-seq as expressed above the threshold of counts per million (CPM) >3. (A) Core clock genes; (B) 15,446 genes remained unchanged (dark blue), and 1330 genes exhibited a significantly different expression level upon siCLOCK-mediated clock disruption (light blue), with 588 being up-regulated (orange) and 742 down-regulated (grey) (FDR <0.05 and p-value <0.05). Altered genes comprised those related to glycerophospholipid and triglyceride metabolism, storage and transport (C) inositol phosphate metabolism (D) and sphingolipid metabolism and storage (E). (F) Relative levels of PC, PE, PI, PS, Cer GlcCer, SM, CL and TAG, analyzed by mass spectrometry based lipidomics in hSKM cells transfected with either siControl (orange bar) or siCLOCK (blue bar). Total phosphatidylcholine (PC) and glycosylceramide (GlcCer) levels are significantly decreased or increased, respectively, upon siCLOCK transfection. Data are mean ± SEM, N = 4 (# p-value <0.05). (*) for FDR <0.05, (**) for FDR <0.01, (***) for FDR <0.001.

https://doi.org/10.7554/eLife.34114.006

We first performed a differential analysis of the global gene expression profile across all 25 time points, starting from 0 hr and until 48 hr following synchronization. Out of the 16,776 mapped genes, the median values for all the time points, reflecting the overall expression levels of 1330 genes (8%), were significantly altered in siCLOCK-transfected hSKM compared to their control counterparts, with 742 being downregulated and 588 being upregulated (Figure 2B; Figure 2—source data 1). As expected, core clock gene expression was affected, with NR1D1 (REVERBα), NR1D2 (REVERBβ), PER3, TEF, BHLHE41 (DEC2), and DBP being significantly downregulated, and CRY1 being upregulated upon CLOCK depletion (Figure 2A).

Functional clocks operative in hSKM are required for proper lipid metabolism and response to insulin

Remarkably, genes encoding for proteins essential for vesicle formation including SNAREs, solute transporters, and Rab GTPases exhibited significantly altered expression levels upon CLOCK depletion (Supplementary file 1-table S2). Additional genes involved in secretion pathways and exhibiting altered mRNA expression levels upon CLOCK depletion are listed in Supplementary file 1-table S3. Using the Panther classification system (Mi et al., 2017) for gene ontology (GO) term analysis, overrepresentation of genes associated with the regulation of nucleotide metabolism, transcription and RNA processing, as well as muscle contraction were identified within the significantly down- and/or upregulated genes (Supplementary file 1-table S4, and Figure 2—source data 2). Furthermore, enrichment analysis using the Reactome pathway database was performed on the down- and/or upregulated genes. Of note, overrepresentation of genes related to muscle contraction, regulation of gene transcription, and cellular responses to stress and membrane trafficking were also identified (Supplementary file 1-table S5, Figure 2—source data 3).

In addition, 42 transcripts involved in lipid metabolism were affected by CLOCK disruption. These comprised genes related to glycerophospholipid and triglyceride metabolism as well as lipid storage and transport (Figure 2C), in addition to those regulating inositol phosphate (Figure 2D) and sphingolipid metabolic pathways (Figure 2E). Importantly, the observed modifications in gene expression level were in an accord with significant alterations in absolute lipid metabolite levels, resulting in total phosphatidylcholine levels being downregulated, and glycosylceramide levels being upregulated in the absence of a functional myotube clock (Figure 2F, Supplementary file 1-table S1, donors 7–10 for in vitro part). The first matching a reduction in lysophosphatidylcholine symporter 1 (MFSD2A) and phosphatase LPIN1 levels (Figure 2C and F), and the latter matching the transcriptome outcome for UDP-glucose ceramide glucosyltransferase (UGCG) the key enzyme of de novo glucosylceramide biosynthesis (Figure 2E and F).

Our differential analysis in human muscle cells demonstrates that genes involved in insulin-stimulated and contraction-induced glucose uptake, comprising TBC1D13, TBC1D4 (AS160), 14-3-3θ (YWHAQ), RAB35, STX6, and PDPK1 (PDK1), were significantly downregulated upon siCLOCK (Supplementary file 1-table S2), highlighting the pleiotropic effect of the skeletal muscle CLOCK gene/protein in regulating glucose uptake in response to insulin and/or to contraction.

To determine the protein levels of candidate genes identified by RNA-seq, hSKM cells established from five additional donor biopsies (for donor characteristics see Supplementary file 1-table S1) were transfected by siRNA targeting CLOCK. Matching the changes observed by RNA-seq, CLOCK mRNA was reduced by siCLOCK as determined by RT-qPCR (Figure 3A), leading to a reduction in CLOCK protein expression by 74% (Figure 3B). Moreover, the expression of the 14-3-3θ protein, a key regulator of GLUT4 translocation (Sakamoto and Holman, 2008; Kleppe et al., 2011), was decreased by about 28% under these conditions (Figure 3C), matching the decrease in its gene expression (Supplementary File 1-table S2). In contrast, AS160 protein levels did not change significantly (Figure 3—figure supplement 1) despite a reduction at the mRNA expression level (Supplementary file 1-table S2).

Figure 3 with 1 supplement see all
Basal and insulin-induced glucose uptake by hSKM is downregulated in the absence of a functional circadian clock.

(A) CLOCK mRNA was measured in hSKM cells transfected with siControl or siCLOCK by RT-qPCR and normalized to the mean of 9S and HPRT. CLOCK expression was reduced by 91 ± 2% (mean ± SEM, N = 4; (***) p-value <0.001) in siCLOCK-transfected cells. Protein levels of CLOCK (B) and 14-3-3θ (C) were assessed by western blot. Left panel: representative western blot; right panel protein quantification over all monoplicates (mean ± SEM, N = 5). CLOCK and 14-3-3θ protein levels were reduced by 75 ± 2%, and 28 ± 8%, respectively. (D) Glucose uptake rates (in pmol/mg.min) measured in the absence (basal) or presence (insulin) of insulin (1 hr, 100 nM) in siControl or siCLOCK-transfected cells. Note significant reduction of basal (31 ± 3%) and insulin-stimulated glucose uptake (28 ± 3%). Data are mean ± SEM of four independent experiments using myotubes from four different donors (same as for A-C). (*) p-value <0.05, (**) p-value <0.01, (***) p-value <0.001.

https://doi.org/10.7554/eLife.34114.012

Finally, we analyzed the impact of clock disruption on the ability of hSKM to take up glucose in response to insulin. The assessment of glucose uptake, using a radioactive glucose analogue, demonstrated an increase in glucose uptake upon insulin stimulation in non-synchronized siControl and siCLOCK-transfected myotubes (siControl: basal vs. insulin p-value = 0.019; siCLOCK: basal vs. insulin p-value = 0.017). Importantly, we observed a significant decrease in both basal (30%), and insulin-stimulated (27%) glucose uptake in siCLOCK-transfected myotubes, as compared to their siControl counterparts (Figure 3D left and right panels, respectively).

RNA-seq reveals rhythmically expressed genes in cultured hSKM synchronized in vitro

We next aimed at identifying genes that exhibited rhythmic profiles in hSKM synchronized in vitro. The existing algorithms JTK_CYCLE (Hughes et al., 2010) and CosinorJ (Mannic et al., 2013) do not allow for a satisfactory analysis of datasets containing large differences in amplitude, observed among the two cycles in our datasets. We therefore developed a novel algorithm, adapted for the analysis of our RNA-seq datasets, comprising high-resolution analysis of two full cycles (25 time points over 48 hr) following in vitro synchronization (see Materials and methods for details). Our in vitro algorithm was validated on the temporal expression profiles of key core clock genes from two published large-scale time series (Atger et al., 2015; Petrenko et al., 2016) and on the model dataset cycMouseLiverRNA from the MetaCycle R package. Briefly, after conversion of the raw data to log2 RPKM values and filtering for log2 RPKM >0, temporal patterns of the resulting 12,985 genes were grouped into 16 models (Figure 4A, Figure 4—source datas 1 and 2), with models 1 to 15 comprising 994 rhythmic genes (7.65%), and model 16 comprising 11,991 non-rhythmic genes (92.35%).

Figure 4 with 3 supplements see all
Temporal gene expression analysis in human skeletal myotubes bearing a disrupted or functional circadian clock.

(A) Total of 12,985 genes were identified by RNA-seq as expressed above log2 RPKM >0. Genes were subjected to the circadian analysis, adapted for high-resolution datasets over 48 hr. Genes were categorized as rhythmic or non-rhythmic (NR) (left diagram) and rhythmic genes (994) were grouped into 15 models (right panel). Genes that were non-rhythmic in either one of the 15 models (11,991 genes) are represented in model 16. (B) Heat maps for genes assigned to models 1 to 4. Relative expression is indicated in green (low) and red (high). (C) Individual temporal expression profiles of core clock genes ARNTL, NR1D1, NR1D2, CRY1, CRY2, PER1, PER2 and PER3 in siControl or siCLOCK-transfected cells. (D) Upper panel: Representative examples for genes assigned to models 1–4. Lower panel: Circadian amplitude quantification of siControl and siCLOCK in models 1–4.

https://doi.org/10.7554/eLife.34114.014

Because the number of rhythmic genes exhibited larger variations between the two analyzed cell cultures, established from two different human individuals than between siControl and siCLOCK (Figure 4—figure supplement 1A), we proceeded with a deeper analysis of models 1 to 4, which group together genes that are rhythmic in the siControl situation for both donor cell cultures. According to our analysis, model 1 comprises 73 genes, classified as rhythmic in both donors upon siControl and siCLOCK, models 2 and 3 include 39 and 34 genes, respectively, that are rhythmic in both siControl donors and in one siCLOCK donor respectively, and model 4 comprises 44 genes that are rhythmic in siControl, but not in siCLOCK (Figure 4B). Circadian core clock genes clustered to model 1, as they exhibited a rhythmic mRNA profile in siControl and a flattened, but still rhythmic, profile upon siCLOCK (Figure 4C and Figure 4—figure supplement 2), indicating a presence of partially functional circadian oscillator. Of note, classification of a temporal gene profile as rhythmic by our algorithm did not take into consideration amplitude alterations, like those generated by siCLOCK-treatment, as long as the temporal profile was qualified as circadian. As amplitude values were indeed often blunted upon siCLOCK treatment, to quantify such amplitude changes a log10 transformation was applied, providing approximation to a normal distribution using a paired t-test. In agreement with our previous publication (Perrin et al., 2015), the amplitude of mRNA accumulation was significantly decreased in siCLOCK samples (Supplementary file 1-table S6, Figure 4—figure supplement 1B).

In summary, 190 genes were qualified as rhythmic in the two analyzed cell cultures, and were clustered into models 1–4 (Figure 4—source datas 1 and 2), as exemplified in Figure 4D (upper panels) and in Figure 4—figure supplement 1C. Importantly, similarly to core clock genes, also these functional genes exhibited a blunted circadian amplitude upon clock disruption (Figure 4D, lower panels, Figure 4—figure supplement 1B). For instance, CAMKK1, classified as rhythmic in both siControl and siCLOCK conditions (model 1), exhibited a significant circadian amplitude reduction upon siCLOCK (Figure 4D, Supplementary file 1-table S6). In addition, SERPINE1, a myokine whose secretion by myotube cells was reduced upon clock disruption (Perrin et al., 2015), presented lower amplitude in siCLOCK-transfected cells (Supplementary File 1-table S6). Panther database analysis for genes assigned to models 1–4 suggested enrichment for a number of GO term and Reactome pathways, comprising cell cycle and mitotic regulation (Figure 4—figure supplement 3, Supplementary file 1-tables S4-5, Figure 4—source datas 3 and 4).

Comparative analysis of diurnal rhythms of gene expression in human skeletal muscle tissue and cultured hSKM

Consequently, we compared rhythmic gene expression between muscle biopsy and cultured hSKM cells (Figure 5A). Among the 190 transcripts that were identified as rhythmic in hSKM cells (Figure 4A, models 1–4, Figure 5—source data 1), 14 transcripts were excluded as they were representing non-protein coding sequences or pseudogenes. Additional 26 genes, associated with mitotic cell cycle functions, and further 17 genes related to cell proliferation, adhesion and differentiation, were only found in hSKM and are likely a consequence of incomplete myotube differentiation in cell culture, as opposed to fully differentiated muscle tissue. Cultured muscle cells lack the in vivo environment and the chemical communication that exists within the tissue. Notably, the absence of neuronal connections further limits the final differentiation of cultured myotubes (for review see [Aas et al., 2013]).

Comparison between the circadian transcriptome of human skeletal muscle and human primary myotubes.

(A) Scatter plot, representing the amplitude of expression in relation to the corrected p-value for genes that were rhythmic in vivo (human muscle biopsies). Genes that were also rhythmic in vitro (hSKM, models 1–4) are colored in red. Blue dots represent genes with a p-value <0.01 and log2 amp >0.5. (B) Phase distribution plot of genes rhythmic in muscle biopsies and primary myotubes shows enrichment at the 04:00 time point. (C) Examples of genes, involved in glucose homeostasis and muscle regeneration, that are rhythmic in vivo and in vitro (RNA-seq data, N = 10).

https://doi.org/10.7554/eLife.34114.022

Among the remaining 133 genes, 99 were expressed in human muscle biopsies. Within this group of overlapping genes, 35 were also qualified as rhythmic at the mRNA level (q-val <0.05) in muscle tissue (Figure 5—source data 1). The genes rhythmic in both in vivo and in vitro datasets included the core clock components (Figure 5A), as well as functional genes that were enriched at the 04:00 time point (Figure 5B). Interestingly, genes implicated in glucose homeostasis (KLF11, GFPT2, NAMPT) and in muscle regeneration (PLAU, PLD1, PIM1), were rhythmic both in vivo and in vitro, suggesting an important role for the circadian clock in regulating muscle physiology (Figure 5C).

Discussion

This study presents the first large-scale circadian transcriptome RNA-seq analysis in muscle biopsies from multiple volunteers and in hSKM cells synchronized in vitro with 2 hr resolution over 48 hr, in the presence of a functional or disrupted cell-autonomous clock, with subsequent analysis of its impact on gene expression. Moreover, we demonstrate that CLOCK depletion in cultured primary skeletal myotubes led to significant changes in gene expression (Figure 2), and related physiological outputs, comprising the regulation of basal and insulin-stimulated glucose uptake, lipid homeostasis (Figure 3), and myokine secretion, as summarized in Figure 6. These results provide new insights into the targets of the molecular clock in human skeletal muscle, previously only studied in rodents (McCarthy et al., 2007; Miller et al., 2007; Dyar et al., 2014; Zhang et al., 2014; Dyar et al., 2015). Finally, to dissect the effects of the cell-autonomous endogenous clock from SCN signals and entrainment cues, this dataset was compared to the diurnal transcriptome of human skeletal muscle collected in form of serial muscle biopsies across 24 hr (Figures 1 and 5).

Schema summarizing impact of clock disruption on muscle metabolic function.

Clock disruption leads to impaired insulin sensitivity and decrease in glucose uptake (1), causes a dysregulation of genes involved in vesicle trafficking (2) and impacts lipid metabolism and lipid metabolite oscillations (3) as reported in Loizides-Mangold et al. (2017).

https://doi.org/10.7554/eLife.34114.024

Comparison between the circadian transcriptome of synchronized myotube cells in vitro, and human muscle tissue collected in vivo

Our in vitro myotube system allows us to explore the transcriptional regulation of muscle target genes without confounding effects of the SCN, rest-activity and feeding-fasting cycles (Harfmann et al., 2015). Regarding the rhythmic analysis of in vitro RNA-seq data, larger variations were observed between the two donors than between siControl and siCLOCK conditions (Figure 4—figure supplement 1), likely due to the genetic heterogeneity among human individuals. The low number of subjects therefore represents a limitation of our study, despite high time-resolution of 2 hr for sample collection conducted over 48 hr, that resulted in as many as 25 time points per myotube donor. We therefore concentrated on genes, which were rhythmic in both donors in siControl condition, irrespectively of their rhythmicity disruption by siCLOCK treatment. In total, 994 circadian genes (7.65% of the global transcriptome) were rhythmic in at least one of the four models (models 1–4, Figure 4A–B), exceeding the value found for U2OS cells, exhibiting 1.5% of oscillating gene transcripts (Krishnaiah et al., 2017). When compared with the diurnal transcriptome of human skeletal muscle biopsies, the percentage of rhythmic genes was considerably lower, likely due to the cell culture situation where the circadian amplitude is gradually lost (Figure 2—figure supplement 2B) in the absence of entrainment (Hughes et al., 2009), or due to effects driven by the SCN or behavioral rhythms rather than by the local peripheral clock. Moreover, we cannot exclude that discrepancies between the in vivo and in vitro circadian datasets are in part also influenced by the fiber type composition of vastus lateralis and gluteus maximus, as demonstrated by myosin isoform analysis (Loizides-Mangold et al., 2017). Additional differences with respect to gene rhythmicity between the two datasets may stem from potential differences due to the different algorithms employed for the data analyses.

Among the 35 genes classified as rhythmic in cell culture and in skeletal muscle tissue, were genes involved in glucose metabolism and in muscle regeneration, including PLAU (LluisLluís et al., 2001) and PIM1 kinase (Fischer et al., 2009), along with core clock components such as NR1D1 and NR1D2 (Figure 5A), previously identified as the most rhythmic transcripts across all human and mouse datasets (Laing et al., 2015). Phospholipase PLD1 (Teng et al., 2015), involved in intracellular membrane trafficking and maintenance of glucose homeostasis in human skeletal muscle (Huang et al., 2005) was rhythmically expressed in vivo and in vitro. Moreover, oscillatory genes in both datasets included NAMPT, KLF11, and GFPT2, the latter controlling the flux of glucose into the hexosamine pathway, tightly linked to hyperglycemia and insulin resistance (Coomer and Essop, 2014). The expression of NAMPT, a key regulator of NAD+ synthesis and muscle maintenance (Frederick et al., 2016), was previously shown to be regulated by CLOCK and BMAL1 in complex with SIRT1 (Ramsey et al., 2009; Garten et al., 2015). Importantly, the diurnal rhythm of secreted NAMPT is disturbed by sleep loss, and positively associates with impairment of postprandial glucose metabolism (Benedict et al., 2012). The transcription factor KLF11, a glucose-inducible regulator of insulin transcription and secretion, that is a member of the Krüppel-like family of transcription factors proposed as circadian (Yoshitane et al., 2014), was found to be regulated by the circadian clock in mouse kidney and epididymal fat tissue (Guillaumond et al., 2010), and is possibly involved in postprandial glucose metabolism in skeletal muscle (Neve et al., 2005).

Comparison of our in vitro dataset to the results published on U2OS cells (Hughes et al., 2009; Krishnaiah et al., 2017) revealed that among the 190 genes that were rhythmic in human skeletal myotubes, 30 were also rhythmic in U2OS. Among those were members of the core-clock machinery, and multiple genes involved in cell cycle progression and mitosis. Moreover, the sulfatase ARSJ, involved in glycosphingolipid metabolism, EXOSC8, a regulator of mRNA stability, LRRC16A, involved in actin filament organization, and TUBA1C, encoding for a structural constituent of the cytoskeleton, as well as E2F1 were rhythmic in skeletal myotubes and in U2OS cells. Interestingly, the RING finger domain protein-encoding gene TRIM47 exhibited a rhythmic expression in human skeletal myotubes, in human muscle biopsies, and in U2OS cells.

A comparison of our in vivo dataset (exonic signals) to published temporal gene expression databases of mouse skeletal muscle (Dyar et al., 2014; Zhang et al., 2014) revealed 107 common rhythmic genes between mouse and human skeletal muscle. When comparing the circadian phases of core-clock components in our database to the temporal profiles of the corresponding genes in rodents (Dyar et al., 2014; Zhang et al., 2014), we observed a phase shift of 8–10 hr. This result is in good agreement with a phase shift observed between peripheral clocks in nocturnal versus diurnal species, which is indeed typically smaller than 12 hr (Mure et al., 2018). The question remains at what level such a phase-shift between nocturnal and diurnal species occurs, and why it is not exactly 12 hr in peripheral organs (Mure et al., 2018). Further comparative studies, conducted in the same type of tissue and with the same methodology, will be required to explore this fundamental issue.

An important observation was the strong induction of genes associated with inflammation and immune response in human muscle in the early morning hours (04:00) (Figure 1F), 16 hr after sampling of the first biopsy. We cannot fully rule out that repeated muscle sampling contributed to this signature, as previously reported for repeated biopsy sampling of a single muscle via the same skin incision site over 25 hr (Friedmann-Bette et al., 2012). However, clinical sampling was optimized to minimize this effect, as serial vastus lateralis biopsies were sampled across alternating limbs and from separate skin incision sites, each proximal to the previous (Van Thienen et al., 2014), not excluding the possibility that circulating molecules may diffuse an inflammatory signal between limbs (Catoire et al., 2012). Importantly, this immune signature was restricted to a single time point in the early morning hours, and thus likely does not exclusively result from responses to muscle injury, which would have further increased at the last time point (08:00). Given that inflammatory cytokines have been described as myokines and important regulators of muscle physiology, this could thus represent a true signature with relevant outcomes for muscle physiology.

One limitation of the comparison between RNA-seq datasets obtained for in vivo and in vitro skeletal muscle and hSKM samples in the present study stems from different analyses methods applied for the two datasets. The algorithms applied here for these two datasets were chosen to optimally fit each dataset differing in the number of time points and the time span of samples collection. Indeed, cellular samples were collected every 2 hr over 48 hr resulting in 25 time points, whereas muscle tissue biopsies were collected every 4 hr over 24 hr, resulting in six time points only, due to obvious practical limitations of repetitive muscle tissue biopsy sampling from the same individuals.

Effect of CLOCK depletion on myotube gene transcription and core clock gene regulation

Efficient clock disruption in adult hSKM cells via siRNA-mediated CLOCK knockdown by our previously validated protocol led to significant changes in gene expression (Figure 2A and B) (Perrin et al., 2015; Petrenko et al., 2016; Loizides-Mangold et al., 2017). Most core clock genes were downregulated upon siCLOCK transfection, in addition to a flattening of the Bmal1-luciferase profile (Figure 2—figure supplement 2B), consistent with our previous data for hSKM and human pancreatic islets (Perrin et al., 2015; Saini et al., 2016). However, despite the observed amplitude blunting, core clock components still presented remnant circadian expression profiles that can likely be attributed to the partial downregulation of CLOCK, and to compensatory mechanisms that occur to maintain the circadian machinery (DeBruyne et al., 2007) (Figure 2A, Figure 2—figure supplement 2A), leading to the observation that the effect on absolute gene expression was more pronounced than the effect on mRNA rhythmicity (Figure 4C and D). Although our established experimental system for cellular clock disruption mediated by efficient CLOCK knockdown proved highly useful to study transcriptional and functional outputs in cultured human primary cells, one should keep in mind that core-clock genes also perform clock unrelated functions. The same holds true for genetic mouse models, where different core-clock gene KO strains exhibit distinct phenotypes. To discriminate between clock-related and unrelated effects of the CLOCK gene knockdown, alternative methods for the circadian clock perturbation will be required.

Muscle fiber type parameters are affected in response to siCLOCK

Gluteus maximus is a slow muscle characterized by high levels of MYH7 expression, fatigue resistance, and slow speed contraction as well as an oxidative metabolic type (Talbot and Maves, 2016). Although the levels of key transcription factors regulating fiber-type-specific genes, including MYOD1, NFATC1, SIRT1 and PPARGC1A were not significantly altered upon siCLOCK, we identified an upregulation of multiple genes characteristic for type I slow fibers, as well as a downregulation of genes associated with type II fast fibers (Supplementary file 1-table S7). Whether these changes affect mitochondrial activity stays to be further explored. Importantly, a recent study has demonstrated that Bmal1 KO myotubes display reduced mitochondrial respiration and a reduced expression of hipoxia-inducible factor 1 (HIF1) target genes (Peek et al., 2017). In agreement, we also observed reduced expression of the HIF target gene VEGFA upon clock disruption (Supplementary file 1-table S7), supporting the hypothesis that clock - HIF interactions play an important role in the glycolytic capacity of skeletal muscle. In addition, we also observed an upregulation of myosin light chain kinase MLCK (MYLK) that contributes to force generation by myofilaments. Taken together these observations reinforce the hypothesis that clock disruption induces a global switch in the genetic program towards slow type I muscle fibers, as it was previously suggested in muscle-specific Bmal1 KO mice (Hodge et al., 2015).

Muscle clock alteration impairs glucose uptake in response to insulin

Skeletal muscle is responsible for 70–80% of insulin-stimulated glucose uptake (DeFronzo and Tripathy, 2009). Importantly, we observed a 30% decrease in glucose uptake for both basal and insulin stimulated conditions in siCLOCK-transfected hSKM (Figure 3C). Previous studies have reported similar observations in either Bmal1-specific muscle KO, or in Clock∆19 mutant mice (Kennaway et al., 2007; Dyar et al., 2014; Harfmann et al., 2016). Recently, it was shown that cardiomyocyte-specific Bmal1 KO and ClockΔ19 mutant mice exhibit defects in insulin-regulated processes, including over-activation of AKT and mTOR signaling (McGinnis et al., 2017). Although we did not see significant changes in GLUT1 or GLUT4 gene expression levels, our differential analysis highlighted many genes involved in the regulation of the GLUT4 translocation pathway (Supplementary file 1-table S2 and Figure 2—source data 1).

Upon closer analysis of the GLUT4 translocation and recycling pathways, we observed changes upon siCLOCK treatment at each step, with several genes being differentially expressed. Specifically, the enzyme PI4K2A, catalyzing the phosphorylation of phosphatidylinositol (PI) to phosphatidylinositol 4-phosphate (PI4P), was downregulated at the mRNA level, which may result in decreased PIP2 and PIP3 levels (Pullen et al., 1998; Sakamoto and Holman, 2008). Additionally, siCLOCK-depleted cells overexpressed MAPKAP1 (mSIN1), one component of the mTORC2 complex required for AKT phosphorylation (Frias et al., 2006), and CAV-3, essential for PI3K and AKT activity as well as GLUT4 translocation in muscle (Fecchi et al., 2006; Tan et al., 2012). Moreover, the observed reduction of 14-3-3θ (YWAHQ) upon siCLOCK at both the mRNA and protein level may lead to an attenuated inhibition of TBC1D1 and TBC1D4 (AS160), and thus block GLUT4 translocation to the plasma membrane (Ramm et al., 2006; Roach et al., 2007; An et al., 2010; Kleppe et al., 2011; Szekeres et al., 2012). Consistent with this theory, we observed a modest upregulation of the Rab-GTPase-activator TBC1D1, in addition to a downregulation of RGC2 and an upregulation of TPM3 at the mRNA level. In summary, regulation of GLUT4 translocation and recycling pathways may be affected upon clock disruption with important consequences on glucose uptake and insulin sensitivity as summarized in Figure 6.

Muscle clock disruption influences the expression of genes involved in vesicle trafficking

GLUT4 located at the plasma membrane, is endocytosed in clathrin-coated vesicles and further recycled (Leto and Saltiel, 2012; Jaldin-Fincati et al., 2017). We observed that several factors of the clathrin-mediated endocytosis machinery were altered upon CLOCK depletion (Supplementary file 1-table S2), among them FNBP1 (FBP17), EPN2, HIP1, and SYT1. Furthermore, our results suggest that CLOCK depletion impacts on calcium levels in the cytoplasm as SYT1 acts as a calcium sensor, and in the presence of elevated Ca2+ levels promotes the fusion of close membranes (Martin, 2015).

Once GLUT4 is endocytosed, it is transported to early endosomes using RAB5 (Stenmark et al., 1994; Leto and Saltiel, 2012). As RAB5B was upregulated upon siCLOCK, it is suggesting that this recycling step might be increased. Moreover, we found downregulation of TBC1D16, which was shown to regulate RAB4 activity, suggesting a possible increase in the fast remobilization of GLUT4 at the plasma membrane (Goueli et al., 2012).

We have previously demonstrated that the basal secretion of myokines, such as IL-6, IL-8, and MCP-1, exhibits a circadian pattern, which was strongly disrupted in hSKM after CLOCK silencing in vitro (Perrin et al., 2015). Here, our transcriptional analysis showed that key regulators of the exocytosis machinery were altered upon clock disruption (Supplementary file 1-tables S2 and S3). Both VAMP3 and STX6, which are involved in IL-6 secretion in mouse macrophages (Manderson et al., 2007), were downregulated at the mRNA level (Supplementary file 1-tables S2 and S6), confirming previous results that clock disruption impacts on vesicle trafficking and secretion (Marcheva et al., 2010; Saini et al., 2016). Importantly, when compared with results from clock disrupted human islets (Saini et al., 2016) we found that numerous genes involved in hormone secretion by pancreatic islets were affected in a similar manner in hSKM (Supplementary file 1-table S8).

Further downstream, GLUT4 is sent to the late endosome for degradation by the lysosome or targeted to the endosomal recycling compartment (ERC), through its interaction with VAMP3 (Dugani et al., 2008; Rose et al., 2009). PI4K2A, which was downregulated upon clock depletion (Supplementary file 1-table S2), might be involved here as it regulates VAMP3 trafficking to perinuclear membranes (Volchuk et al., 1995; JovicJović et al., 2014). In addition, CAMSAP2, involved in microtubule stabilization (Hendershott and Vale, 2014), and KIF13A, associated with recycling endosome tubules (Delevoye et al., 2014), were also downregulated upon CLOCK disruption (Supplementary file 1-tables S2 and S3). Taken together these results, as summarized in Figure 6, suggest that the muscle clock may play an important regulatory function in the control of the secretion machinery via transcriptional regulation.

Cell-autonomous clock disruption in hSKM might impact energy substrate utilization

The circadian clock has been associated with the control of muscle development and regeneration, as clock mutant mice exhibit defects in muscle metabolism, architecture and composition (for review see [Chatterjee and Ma, 2016; Schiaffino et al., 2016]). Here, we found alterations in the expression of several genes involved in lipid metabolism, calcium handling, electron transport chain, and glucose metabolism (Figure 2C–E, Supplementary file 1-table S7), suggesting a shift in energy substrate utilization upon clock disruption, as has been proposed previously in rodents upon loss of Bmal1 (Hodge et al., 2015; Harfmann et al., 2016). AMP-activated protein kinase, a potent regulator of skeletal muscle fat metabolism (Thomson and Winder, 2009) might be dysregulated upon clock disruption as we observed upregulation of its regulatory subunit PRKAG2 and downregulation of subunit PRKAG3. Previous work reported downregulation of both subunits in Bmal1-specific muscle KO mice (Hodge et al., 2015), suggesting that this gene could be directly controlled by BMAL1.

Clock disruption causes changes in lipid levels as has been described previously for the liver of Per1/2 KO mice (Adamovich et al., 2014). In hSKM, siCLOCK treatment affected several genes involved in lipid metabolic processes, lipid storage and transport (Figure 2C–E), which resulted in total phosphatidylcholine and glycosylceramide level alterations (Figure 2F). Specifically, we found an increase in the long chain fatty acid transporter CD36 and in FABP3, consistent with previous results obtained in mouse skeletal muscle upon clock disruption (Hodge et al., 2015; Schiaffino et al., 2016). In addition, we observed an upregulation of MSTN upon siCLOCK, which could further promote the increase in CD36 and FABP3, leading to impaired glucose uptake (Figure 3D). Interestingly, muscle-specific myostatin (Mstn) KO mice exhibit a reduction of lipid transporters, including FABP3 and CD36, a diminution of lipid oxidation, higher levels of saturated and unsaturated fatty acids, and a decrease of cardiolipin and triglycerides (Baati et al., 2017). Furthermore, downregulation of Mstn in skeletal muscle from type one diabetic mice leads to an increase of Glut1 mRNA and GLUT4 protein levels, promoting insulin-stimulated glucose uptake (Coleman et al., 2016). Altogether, these results confirm previous rodents studies and indicate a shift in substrate utilization in skeletal muscle from carbohydrates to lipids with impact on muscle metabolism and glucose homeostasis (Dyar et al., 2014; Dyar et al., 2015; Hodge et al., 2015; Harfmann et al., 2016).

Conclusions

In summary, our study provides (1) a comparison between rhythmic transcriptome databases obtained from human muscle tissue and cultured primary cells derived from muscle biopsies, and (2) the identification of pathways regulated by CLOCK in skeletal muscle, involved in glucose uptake, myokine secretion, and lipid metabolism (Figure 6). Human primary cells cultured in vitro have been used as a model to study human disease and metabolism (Aas et al., 2013; Saini et al., 2015). In combination with tissue analysis as presented here, primary cell culture constitutes a powerful model to study human metabolism, and warrants further analyses in additional metabolically active tissues in physiological conditions, and in the context of metabolic diseases.

Materials and methods

Key resources table
Reagent type
(species) or resource
DesignationSource or referenceIdentifiersAdditional information
Antibodyanti-AS160 (C69A7) Rabbit mAbCell SignalingCat. #2670
RRID:AB_2199375
1:1000; for western blot; primary Ab
Antibodyanti-14-3-3 θ Rabbit polyclonalCell SignalingCat. #9638
RRID:AB_2218251
1:200; for western blot; primary Ab
Antibodyanti-CLOCK(H276) Rabit polyclonalSanta Cruz BiotechnologyCat. sc-25361
RRID:AB_2260802
1:200; for western blot; primary Ab
Antibodyanti-actin Rabbit polyclonalSigma-AldrichCat. A2066
RRID:AB_476693
1:1000; for western blot; primary Ab
Antibodygoat anti-rabbit-IgG HRPSigma-AldrichCat. A8275
RRID:AB_258382
1:3000; for western blot; secondary Ab
Recombinant DNA reagentBmal1-luciferase
(luc) reporter
Liu et al., 2008; PMID:18454201
Sequence-based reagentON-TARGETplus Non-targeting PoolDharmaconD-001810-10-20
Sequence-based reagentON-TARGETplus human CLOCK siRNA SMARTpoolDharmaconL-008212-00-0020Target Sequences: CAACUUGCACCUAUAAAUA CGACAGGACUGGAAACCUA GAACAACGGACACGCAUGA CUAGAAAGAUGGACAAAUC
Peptide, recombinant proteinNANANANA
Commercial assay or kitSuperSignal West Pico Chemiluminescent SubstrateThermo Fisher ScientificProd. #34080
Commercial assay or kitQuant-iTª RiboGreenª RNA Assay KitThermo Fisher ScientificR11491
Commercial assay or kitRNeasy Mini kitQiagenRef # 74104
Commercial assay or kitTruSeq Stranded Total RNA Library Prep Kit with Ribo-Zero
Gold Set A (48 samples, 12 indexes)
IlluminaRS-122–2301
Commercial assay or kitTruSeq Stranded Total RNA Library Prep Kit with Ribo-Zero Gold Set B (48 samples, 12 indexes) Indexes onlyIlluminaRS-122–2302
Commercial assay or kitTruSeq RNA Library Prep Kit v2IlluminaRS-122–2001/RS-122–2002
Commercial assay or kitHiSeq PE Cluster Kit V4 - cBotIlluminaPE-401–4001
Commercial assay or kitHiSeq SBS Kit V4 250 cycle kitIlluminaFC-401–4003
Commercial assay or kitKAPA HiFi HotStart
ReadyMixPCR Kit
Kapa BioSystems
(Roche)
KK2602
Commercial assay or kitQuant-iT PicoGreen dsDNA Assay KitThermo Fisher ScientificP7589
Commercial assay or kitLabChip DNA High Sensitivity Reagent KitPerkin ElmerCLS760672
Chemical compound, drugForskolinSigma-AldrichF6886
Chemical compound, drugLuciferinProlume LTD#260150
Chemical compound, drug2-deoxy-[3H]-D-glucosePerkinElmerNET328A001MCSpecific Activity: 5–10 Ci (185-370GBq)/mmol, 1mCi (37MBq)
InsulinSigma-AldrichI9278
Chemical compound, drugPotassium phosphate monobasicSigma-AldrichP5655
Chemical compound, drugHiPerFect transfection reagentQiagenCat No./ID: 301705
Chemical compound, drugTert-butyl methyl etherSigma-Aldrich#20256
Chemical compound, drugMethylamine solutionSigma-Aldrich#534102
Chemical compound, drugMethanol LC-MS CHROMASOLVFluka
(Thermo Fisher Scientific)
#34966
Chemical compound, drugWater LC-MS CHROMASOLVFluka
(Thermo Fisher Scientific)
#39253
Chemical compound, drugChloroform, stabilized with ethanol, for HPLCACROS Organics
(Thermo Fisher Scientific)
#390760010
Chemical compound, drug12:0 PC (DLPC)Avanti Polar Lipids#850335
Chemical compound, drug17:0-14:1 PEAvanti Polar LipidsLM1104
Chemical compound, drug17:0-14:1 PIAvanti Polar LipidsLM1504
Chemical compound, drug17:0-14:1 PSAvanti Polar LipidsLM1304
Chemical compound, drug12:0 SM (d18:1/12:0)Avanti Polar Lipids#860583
Chemical compound, drugC17 Ceramide (d18:1/17:0)Avanti Polar Lipids#860517
Chemical compound, drugC8 Glucosyl(ß) Ceramide (d18:1/8:0)Avanti Polar Lipids#860540
Software, algorithmRstudioRstudioR version 3.3.1
Software, algorithmPrism 5GraphPadNA
Software, algorithmExcel 2016MicrosoftNA
Software, algorithmLumiCycleActimetricsNA
Software, algorithmSTAR: ultrafast universal RNA-seq alignerDobin et al., 2013
PMID:
23104886
Software, algorithmedgeR packageRobinson et al., 2010
PMID:
19910308
edgeR version 3.16.5
Software, algorithmlme4 R packageBates et al., 2015
DOI:10.18637/jss.v067.i01
Software, algorithmlmtest R packageZeileis et al., 2002
DOI:10.18637/jss.v007.i02
Software, algorithmTopGO R packagehttps://bioconductor.riken.jp/packages/3.3/bioc/vignettes/topGO/inst/doc/topGO.pdf
Software, algorithmGEMToolshttp://gemtools.github.io/GEMTools v1.7.1
Lipid data analyzer IIIGB-TUG Graz University; PMID:21169379LDA v.2.5.1
Software, algorithmImage LabBio-Rad
OtherRIPA bufferSigma-AldrichCat# R0278

Human skeletal muscle biopsies

10 healthy volunteers were recruited for the in vivo study (see Supplementary file 1-table S1 for donor characteristics). One week prior to the scheduled laboratory visit, participants had to adhere to a consistent daily feeding and sleeping routine. Participants arrived in the laboratory at 19:00 hr on the day prior to the testing day and ingested one standardized meal that first evening. Participants remained for the duration of their stay in a semi-recumbent position. During the waking hours of the testing day, they were given mixed-macronutrient meal-replacement solutions at hourly intervals (Resource, Nestlé, Switzerland) to ensure energy balance. The laboratory was free from natural light and with artificial lighting standardized to 800 lux in the direction of gaze, ambient temperature maintained between 20 and 25°C and noise levels tightly regulated. Participants were not permitted to sleep during waking hours when lights were on (i.e. 07:00-22:00 hr) and wore eye masks whilst trying to sleep during lights-out (i.e. 22:00-07:00 hr). Anesthetic administration (1% lidocaine w/o epinephrine) and skin/fascia incisions for this procedure (Bergstrom, 1962) were completed within the hour prior to sleep such that night-time samples could be acquired with minimal discomfort. Six 100 mg biopsy samples were acquired from the vastus lateralis muscle at 4 hr intervals (12:00, 16:00, 20:00, 24:00, 04:00 and 08:00 hr) and immediately snap frozen under liquid nitrogen. Samples were taken from each leg in alternating and ascending order with skin incisions separated by 2–3 cm. The study was conducted in accordance with the Declaration of Helsinki and with the approval of the Health Research Authority (NRES Committee South West; 14/SW/0123) and the Swiss Commission cantonal (Canton Vaud) d’éthique de la recherche (Cer-VD). For further details see Loizides-Mangold et al. (2017).

Human muscle RNA-sequencing and data analysis

Total Stranded RNA-Seq (in vivo muscle samples): RNA was quantified with Ribogreen (Life Technologies, Carlsbad, CA) and quality was assessed on a Fragment Analyzer (Advances Analytical). Sequencing libraries were prepared from 250 ng total RNA using the TruSeq Stranded Total LT Sample Prep Kit (Illumina, San Diego, CA) with the Ribo-Zero Gold depletion set. The procedure was automated on a Sciclone NGS Workstation (Perkin Elmer, Waltham, MA). The manufacturer’s protocol was followed, except for the PCR amplification step. The latter was run for 13 cycles with the KAPA HiFi HotStart ReadyMix (Kapa BioSystems, Roche, Switzerland). This optimal PCR cycle number was evaluated using the Cycler Correction Factor method as described previously (Atger et al., 2015). Purified libraries were quantified with Picogreen (Life Technologies) and the size pattern was controlled with the DNA High Sensitivity Reagent kit on a LabChip GX (Perkin Elmer). Libraries were then pooled by 24, and each pool was clustered at a concentration of 8 pmol on 8 lanes of v4 paired-end sequencing flow cells (Illumina). Sequencing was performed for 2 × 125 cycles on a HiSeq 2500 strictly following Illumina’s recommendations.

Paired-end reads were mapped to the Homo sapiens genome (GRCh38/hg38) using STAR (Dobin et al., 2013) with parameters “--alignIntronMin 20 --alignIntronMax 1000000 --GTF (option –sjdbGTFfile). Mapped reads were quantified in intronic and exonic regions. For each gene, we defined a gene body by merging all the respective transcripts using BEDtools (Quinlan and Hall, 2010). A region was defined as exonic if it occurs in a least one of the transcripts while an intronic region had to be shared between all the transcripts. We assigned uniquely mapping paired-reads to exonic regions (exon/exon and complete exon reads) or intronic regions (intron-exon and complete intron reads) considering reads orientation. Genes with less than two intronic reads or 10 exonic reads on average were discarded. Intronic and exonic reads count were normalized using edgeR (Robinson et al., 2010) by the library scaling factor from the exonic regions and the respective intronic and exonic length (rpkm()). Genes with less than −2 RPKM (log2) at the exonic level were discarded. Genes with less than −3 RPKM (log2) at the intronic level were reported as NA for the intronic quantification.

Rhythmicity was assessed with a linear mixed-effects model using lmer() function from the lme4 R package applied on the log2 normalized reads count. A standard harmonic regression with a 24 hr period was fitted with a donor-dependent random effect on the baseline:

yID, t ~cos2π24 t+ sin2π24t+(1|ID)

where yID,t is the log2 normalized reads count for patient ID at time t. This full model was compared to the null model (without the harmonic terms) using the likelihood ratio test function ltest() from the lmtest R package. The p-values were computed from a chi-squared distribution and were adjusted using the Benjamini-Hochberg procedure.

Gene ontology analysis was performed using the TopGO R package . Enrichment analysis for GO terms derived from ‘Biological Process’ was performed for the genes rhythmic in the three groups (R-I.RE, R-I, and R-E) and in the different phase bins. GO terms with p-value <0.05, a minimum number of 3 genes, and less than 200 annotated genes were considered.

Transcription factors enrichment analysis

Predictions of transcription factor binding sites (MotEvo) and promoter regions were downloaded from http://swissregulon.unibas.ch/sr/downloads (database Homo sapiens, hg19:FANTOM5) (Pachkov et al., 2013). Transcription factor binding sites were assigned to their corresponding genes using the promoter regions table. Genes rhythmic at the intronic level, and with amplitude larger than 0.5 (log2), were grouped according to their phase in 4 hr bins. All the genes expressed in the dataset were used as a background and a hypergeometric test was computed for the over-representation of transcription factor binding sites in the different bins. Transcription factor binding sites with -log10(p-value) >104 and belonging to at least five genes were reported and annotated in Figure 1I.

Study participants, skeletal muscle tissue sampling and primary cell culture

Biopsies from the Gluteus maximus muscle were derived from donors with their informed consent (see Figure 2—figure supplement 1 and Supplementary file 1-table S1 for donor characteristics). The experimental protocol (‘DIOMEDE’) was approved by the Ethical Committee SUD EST IV (Agreement 12/111) and performed according to the French legislation (Huriet’s law). All donors had HbA1c levels inferior to 6.0% and fasting glycemia inferior to 7 mmol/L, were not diagnosed for type 2 diabetes (T2D), neoplasia or chronic inflammatory diseases, and not doing shift work. Primary skeletal myoblasts were purified and differentiated into myotubes according to the previously described procedure (Agley et al., 2015; Perrin et al., 2015). After reaching confluence, myoblasts were differentiated into myotubes during 7–10 days in DMEM supplemented with 2% FBS.

siRNA transfection and lentiviral transduction

Human primary myoblasts were differentiated into myotubes as described above. Cells were transfected 24 hr (RNA-seq) or 72 hr (western blot, glucose uptake) prior to experiment with 20 nM siRNA targeting CLOCK (siCLOCK), or with non-targeting siControl (Dharmacon, Lafayette, CO), using HiPerFect transfection reagent (Qiagen, Hilden, Germany) following the manufacturer's protocol. To produce lentiviral particles, Bmal1-luciferase (Liu et al., 2008) lentivectors were transfected into 293 T cells using the polyethylenimine method (for detailed procedure see [Pulimeno et al., 2013]). Myoblasts were transduced with the indicated lentiviral particles with a multiplicity of infection (MOI) = 3 for each, grown to confluence, and subsequently differentiated into myotubes.

In vitro skeletal myotube synchronization and real-time bioluminescence recording

Primary myotubes were synchronized with 10 µM forskolin (Sigma-Aldrich, St. Louis, MO) for 60 min at 37°C in a cell culture incubator, then the medium was changed to a phenol red-free recording medium containing 100 µM luciferin (Prolume LTD, Pinetop, AZ), and cells were transferred to a 37°C light-tight incubator as previously described by us (Pulimeno et al., 2013). Bioluminescence from each dish was continuously monitored using a Hamamatsu photomultiplier tube (PMT) detector assembly. Photon counts were integrated over 1 min intervals. Luminescence traces are shown as detrended data.

Glucose uptake measurement

Human myotubes treated with siControl or siCLOCK as described before were serum-starved for 3 hr then incubated with 2-deoxy-[3H]-D-glucose for 15 min. Incubations were performed with or without insulin stimulation (1 hr, 100 nM). After incubation, the medium was removed prior to cell lysis in 0.05 M NaOH. Cell content radioactivity was determined by liquid scintillation counting (Perkin Elmer, 2900TR) and protein content was quantified by using the Bradford protein assay. Glucose transport is expressed in pmol/mg.min (Chanon et al., 2017).

Protein analysis

Human myotubes transfected with siControl or siCLOCK for 24 to 72 hr, were lysed in RIPA buffer. Protein extracts (8 µg) were analyzed by SDS-PAGE and immunoblotted to 0.45 µm nitrocellulose membrane or 0.2 µm PVDF membrane using a wet system (Bio-Rad, Hercules, CA) according to the manufacturer’s instructions. Membranes were blocked and incubated with primary and secondary antibodies in 5% BSA/TBS-T 0.5% or 5% BSA/TBS-T 0.1%. Primary and secondary antibodies were used at the following dilutions: anti-TBC1D4/AS160 (1/1000, Cell Signaling, Danvers, MA, #2670S), anti-14-3-3θ (1/200, Cell Signaling, #9638S), anti-CLOCK (1/200, Santa Cruz Biotechnology, Santa Cruz, CA, H-276) and anti-ACTIN (1/1000, Sigma-Aldrich, A2066), anti-rabbit-HRP (1:3000, Sigma-Aldrich A8275). Signals were visualized using SuperSignal West Pico Chemiluminescent Substrate (Thermo Fisher Scientific, Waltham, MA). For protein quantification, five donors were analyzed but only the representative western blot result of one donor is shown.

Lipidomics

The lipidomics analysis was performed as described in Loizides-Mangold et al. (2017). Briefly, human primary myotubes were harvested from one confluent 10 cm dish (~1.5×106 cells) and lipid extracts were prepared using the MTBE protocol (Matyash et al., 2008). Total phosphorus was determined as described in (Loizides-Mangold et al., 2017). Phospho- and sphingolipid were analyzed by mass spectrometry using on a TSQ Vantage Triple Stage Quadrupole Mass Spectrometer (Thermo Fisher Scientific) equipped with a robotic nanoflow ion source (Nanomate HD, Advion Biosciences, Ithaca, NY), using multiple reaction monitoring (MRM). Lipid concentrations were calculated relative to the relevant internal standards and then normalized to the total phosphate content of each total lipid extract. Triacylglyceride analysis was performed by mass spectrometry analysis on a hybrid ion trap LTQ-Orbitrap XL mass spectrometer (Thermo Fisher Scientific) equipped with a micro LC binary pump UFLC-XR (Shimadzu, Kyoto, Japan). Lipid identification was carried out with the Lipid Data Analyzer II (LDA v. 2.5.1, IGB-TUG Graz University) (Hartler et al., 2011).

hSKM mRNA extraction and quantitative PCR analysis

Differentiated myotubes were synchronized by 10 µM forskolin, collected every 2 hr during 48 hr (0–48 hr), deep-frozen in liquid nitrogen and kept at −80°C. Total RNA was prepared using RNeasy Mini kit (Qiagen). 0.5 µg of total RNA was reverse-transcribed using Superscript III reverse transcriptase (Invitrogen) and random hexamers, and PCR-amplified on a LightCycler 480 (Roche). Mean values for each sample were calculated from the technical duplicates of each RT-qPCR analysis, and normalized to the mean of two housekeeping genes (HPRT and 9S), which served as internal controls. Primers used for this study are listed in Supplementary file 1-table S9.

hSKM RNA-sequencing and data analysis

Total RNA was prepared from primary human skeletal myotubes from two donors, transfected either with siControl or siCLOCK, synchronized with a forskolin pulse and collected every 2 hr during 48 hr in duplicates, using RNeasy Mini Kit (Qiagen). Total RNA libraries were prepared from 300 ng of RNA following customary Illumina TruSeq v2 protocols for next generation sequencing. PolyA-selected mRNAs were purified, size-fractioned, and subsequently converted to single-stranded cDNA by random hexamer priming. Following second-strand synthesis, double-stranded cDNAs were blunt-end fragmented and indexed using adapter ligation, after which they were amplified and sequenced according to protocol. RNA libraries were sequenced with a HiSeq2000 sequencer producing 49pb paired-end reads. Standard quality checks for material degradation (Bioanalyzer and TapeStation, Agilent Technologies, Santa Clara, CA) and concentration (Qubit, Life Technologies) were done before and after library construction, ensuring that samples are suitable for sequencing.

Paired-end reads were mapped to the human genome (version GRCh37/hg19) with GEMTools (v1.7.1) using GENCODE v19 as gene annotation. RNA-seq reads were subsequently filtered for correct orientation of the two ends, a minimum mapping quality score of 150 and allowing a maximum of 5 mismatches in both ends. GENCODE annotation v19 was used to assign filtered reads to their corresponding exons and genes. For each gene, we processed the exons from all transcripts, which were quantified by considering only filtered reads as above, in which both ends map to exons of the same gene. The gene counts were incremented non-redundantly, meaning reads overlapping two exons were counted once to the total count sum per gene.

The differential gene expression analysis was performed with the R package edgeR (Robinson et al., 2010). First, transcripts expressed lower than three counts per million (CPM) and noninformative (e.g. non-aligned) were filtered. To minimize the log-fold changes between the samples for most genes, a set of scaling factors for the library sizes were estimated with the trimmed mean of M-values (TMM) method (Robinson and Oshlack, 2010). The dispersion was estimated with the quantile-adjusted conditional maximum likelihood (qCML) method. Once the dispersion estimates are obtained, we performed the testing procedures for determining differential expression using the exact test (Robinson and Smyth, 2008).

Regarding the rhythmic analysis, homemade algorithm was developed to analyze these RNA-seq data. In short, raw data were transformed to log2 reads per kilobase per million mapped reads (RPKM) as described previously (Atger et al., 2015), then only transcripts with log2 RPKM >0 for each of the fourth conditions (two subjects, siControl or siCLOCK) were kept avoiding big variability for weakly expressed transcripts. The 48 time points of each condition were used to define a local regression function (LOESS). This step allows smoothing the curve and reducing local variability. The function was then used to calculate 10 different measures (maximum and minimum slopes, first and second extremum, minimum-maximum ratio, autocorrelation, measure of scattering, residues on the loess function, residues on a linear function and period). These features were used to classify gene expression patterns in four different groups: rhythmic genes (category ‘circadian’), genes that show only one peak at the beginning of the time course (category ‘one peak’), linearly (category ‘linear’) and scatteredly expressed genes (category ‘cloud’). The algorithm attributes a probability to each transcript per condition. To be classified in one category, this probability must be the highest value and superior to 0.5 in at least one category. If no probabilities are superior to 0.5 for the four categories, transcripts are grouped into model 16 (non-rhythmic). The 11 major circadian genes, including ARNTL (BMAL1), NR1D1 (REVERBα), NR1D2 (REVERBβ), PER1, PER2, PER3, CRY1, CRY2, NPAS2, TEF and BHLHE41, were selected to train a random forest model. The same number of genes for the three other groups were also integrated in the training dataset. This dataset was then passed to the training algorithm for random forests and gene conditions that were assigned to one of these categories with a high score (0.9) were integrated in the training dataset. This procedure was repeated until 500 curves per group were identified. The last model was kept to classify the whole dataset.1485 curves from 994 transcripts were identified as rhythmic among 12,985 transcripts. Altogether, transcripts were grouped into 16 models.

Bioluminescence and statistical data analysis

Bioluminescence data were analyzed with the Actimetrics LumiCycle analysis software (Actimetrics LTD, Wilmette, IL). RNA-seq data and qPCR data analysis were performed using RStudio, GraphPad Prism five and Excel 2016. Panther analyses were performed using the PANTHER version 12.0 released on 10.07.2017. Statistical analyses were performed using Student’s t-test. Differences were considered significant for (*) p-value <0.05, (**) p-value <0.01, and (***) p-value <0.001. Exact p-values and raw data for Figures 2 and 3 are listed in Supplementary file 2.

Mycoplasma test for primary cultures

Since primary cultures, established from human skeletal muscle tissue biopsies were used in this study, mycoplasma contamination tests were conducted only once for each primary myotube culture. To do so, 100 μl of culture medium were taken 48 hr following the last trypsinization, boiled at 95°C for 5 min, and centrifuged for 10 s at 14,000 rpm. PCR was performed on 5 μl of thus processed samples, using a mix of primers listed in Supplementary file 1-table S9. The PCR program was 5 min at 95°C, followed by 30 cycles with 95°C 30 s, 60°C 30 s, 72°C 30 s, and a final elongation at 72°C for 10 min. PCR products were separated on a 1.5% agarose gel.

References

  1. 1
  2. 2
  3. 3
  4. 4
  5. 5
  6. 6
  7. 7
  8. 8
  9. 9
  10. 10
  11. 11
  12. 12
    Fitting Linear Mixed-Effects Models Usinglme4
    1. D Bates
    2. M Mächler
    3. B Bolker
    4. S Walker
    (2015)
    Journal of Statistical Software, 67, 10.18637/jss.v067.i01.
  13. 13
  14. 14
    Muscle Electrolytes in Man - Determined by Neutron Activation Analysis on Needle Biopsy Specimens - Study on Normal Subjects, Kidney Patients, and Patients with Chronic Diarrhoea
    1. J Bergstrom
    (1962)
    Scandinavian Journal of Clinical and Laboratory Investigation, 14.
  15. 15
  16. 16
    Glucose uptake measurement and response to insulin stimulation in in vitro cultured human primary myotubes
    1. S Chanon
    2. C Durand
    3. A Vieille-Marchiset
    4. M Robert
    5. C Dibner
    6. C Simon
    7. E Lefai
    (2017)
    Journal of Visualized Experiments, 10.3791/55743, 28671646.
  17. 17
  18. 18
  19. 19
  20. 20
  21. 21
  22. 22
  23. 23
  24. 24
  25. 25
  26. 26
  27. 27
  28. 28
  29. 29
  30. 30
  31. 31
  32. 32
  33. 33
  34. 34
  35. 35
  36. 36
  37. 37
  38. 38
  39. 39
  40. 40
  41. 41
  42. 42
  43. 43
  44. 44
  45. 45
  46. 46
  47. 47
  48. 48
  49. 49
  50. 50
  51. 51
  52. 52
  53. 53
  54. 54
  55. 55
  56. 56
  57. 57
  58. 58
  59. 59
  60. 60
  61. 61
  62. 62
  63. 63
  64. 64
  65. 65
  66. 66
  67. 67
    PI(4,5)P₂-binding effector proteins for vesicle exocytosis
    1. TF Martin
    (2015)
    Biochimica Et Biophysica Acta (BBA) - Molecular and Cell Biology of Lipids 1851:785–793.
    https://doi.org/10.1016/j.bbalip.2014.09.017
  68. 68
  69. 69
  70. 70
  71. 71
  72. 72
  73. 73
  74. 74
  75. 75
  76. 76
  77. 77
  78. 78
  79. 79
    Parallel measurement of circadian clock gene expression and hormone secretion in human primary cell cultures
    1. V Petrenko
    2. C Saini
    3. L Perrin
    4. C Dibner
    (2016)
    Journal of Visualized Experiments, 117, 10.3791/54673.
  80. 80
  81. 81
    Phosphorylation and activation of p70s6k by PDK1
    1. N Pullen
    2. PB Dennis
    3. M Andjelkovic
    4. A Dufner
    5. SC Kozma
    6. BA Hemmings
    7. G Thomas
    (1998)
    Science 279:707–710.
  82. 82
  83. 83
  84. 84
  85. 85
  86. 86
  87. 87
  88. 88
  89. 89
  90. 90
  91. 91
  92. 92
  93. 93
  94. 94
  95. 95
    Inhibition of rab5 GTPase activity stimulates membrane fusion in endocytosis
    1. H Stenmark
    2. RG Parton
    3. O Steele-Mortimer
    4. A Lütcke
    5. J Gruenberg
    6. M Zerial
    (1994)
    The EMBO Journal 13:1287–1296.
  96. 96
  97. 97
  98. 98
  99. 99
  100. 100
  101. 101
  102. 102
  103. 103
  104. 104
  105. 105
  106. 106
  107. 107
    strucchange: An R Package for Testing for Structural Change in Linear Regression Models
    1. A Zeileis
    2. F Leisch
    3. K Hornik
    4. C Kleiber
    (2002)
    Journal of Statistical Software, 7, 10.18637/jss.v007.i02.
  108. 108

Decision letter

  1. Joseph S Takahashi
    Reviewing Editor; Howard Hughes Medical Institute, University of Texas Southwestern Medical Center, United States

In the interests of transparency, eLife includes the editorial decision letter and accompanying author responses. A lightly edited version of the letter sent to the authors after peer review is shown, indicating the most substantive concerns; minor comments are not usually included.

Thank you for submitting your article "Transcriptomic analyses reveal rhythmic and CLOCK–driven pathways in human skeletal muscle" for consideration by eLife. Your article has been reviewed by two peer reviewers, and the evaluation has been overseen by a Reviewing Editor and Fiona Watt as the Senior Editor. The following individuals involved in review of your submission has agreed to reveal his identity: Joseph Bass (Reviewer #2).

The reviewers have discussed the reviews with one another and the Reviewing Editor has drafted this decision to help you prepare a revised submission

Summary:

In this manuscript titled 'Transcriptomic analyses reveal rhythmic and CLOCK–driven pathways in human skeletal muscle', the authors performed RNA–seq analysis in muscle tissue biopsies (in vivo) and cultures of synchronized primary differentiated myofibers (in vitro) from human subjects, allowing them to compare cell autonomous and non–cell autonomous rhythms of transcription. Of note, the authors found that rhythmic transcription of genes involved in glucose transport (i.e. GLUT4 regulation), lipid metabolism and immune response were present in both the in vivo and in vitro samples. However, many other genes which were rhythmic in vivo did not appear rhythmic in vitro, suggesting that there are non–cell autonomous factors controlling rhythmic transcription in human muscle. In addition, since total RNA was used, both pre–mRNA and mature mRNAs were quantified. Indeed, the authors identified genes which oscillate at both the transcriptional and post–transcriptional levels. Finally, the dependence of rhythmic transcripts on the core molecular clock was determined using siRNA–mediated knock down of CLOCK in the primary myotube cultures. siCLOCK myotubes displayed reduced expression of many genes which were found to be rhythmic in vitro. The authors also determined that siCLOCK myotubes display reduced insulin–dependent and –independent glucose uptake, suggesting reduced activity of the GLUT4–transport mechanism. Overall, the finding that the circadian clock controls glucose and lipid metabolic gene expression in primary human myofibers is important and provides further evidence of the role of the molecular clock in the regulation of muscle function and glucose homeostasis.

The daily transcriptome analysis of the human muscle biopsies taken at 4 h intervals across 24 h from 10 healthy individuals serve as a valuable and important resource, especially in view of the current attempts to map circadian rhythmicity in humans. In this regard, a more rigorous comparison between available detests of mouse muscle transcriptome (e.g. McCarthy et al., 2007; Miller et al., 2007; Dyar et al., 2014; Zhang et al., 2014; Hodge et al., 2015) with this new data would provide interesting insight regarding the difference in muscle rhythmic gene expression between mice and humans, and likely between nocturnal and diurnal species, respectively.

The major strength of this study lies on the expression data obtained from human healthy individuals in a well–controlled manner, this section by itself is an important resource.

The in vitro experiments primarily corroborate previous studies done with cultured cells (e.g. Krishnaiah et al., 2017) that show dramatic reduction in rhythmic gene expression in culture compared to organs in vivo.

Essential revisions:

1) The in vivo biopsy samples were taken from the vastus lateralis muscle whereas the in vitro samples were from gluteus maximus. Therefore, is it possible that the observed gene expression differences observed between the in vitro and in vivo sample sets are explained by the fact that they are derived from different muscles? One option to address this would be to compare gene expression profiles from biopsies of VL and GM muscles – do they differ in a similar way to the in vitro versus in vivo gene sets? This may beyond the scope here.

2) The authors observed gene expression changes with siCLOCK that support increased mitochondrial oxidative/type I muscle remodeling. This result is consistent with what was observed in RNA–sequencing studies in Bmal1KO mouse muscle (Hodge et al., 2015). However, others have provided evidence for reduced mitochondrial oxidative metabolism in Bmal1KO muscle cells. Therefore, it is possible that the gene signature does not equate with the metabolic phenotype. For instance, gene expression changes may be compensatory rather than causal? In the absence of such analysis, discussion should acknowledge the interest in future functional bioenergetic profiling (e.g., respirometry with mitochondrial fuel substrates).

3) The mathematical model used to detect rhythmicity should be validated on an existing dataset or a demo dataset in order to compare their results with established and widely used algorithms that detect rhythmic profiles, such as JTK_CYCLE (Hughes et al., 2010) and Meta Cycle (Wu et al., 2016) that incorporates ARSER, JTK_CYCLE and Lomb–Scargle.

4) It would be interesting to compare the genes identified as rhythmic in cultured U2OS cells (Krishnaiah et al., 2017) with the ones identified in the current study, even though these are completely different cells. It might shed light on genes, aside from the core clock, that globally maintain rhythmicity in culture.

5) The authors conclusion that the cell–autonomous circadian clock has an essential role in coordinating muscle glucose homeostasis and lipid metabolism in humans should be revised as it is not supported by their results. Their finding that knockdown of clock in vitro affected the overall expression of ~8% of all genes with genes related to glucose and lipid metabolism suggests that the CLOCK protein itself modulate their expression and not necessarily the oscillator. Furthermore, the result that only few genes maintain rhythmicity both in vivo and in vitro, most of them are core clock components, does not support it either.

https://doi.org/10.7554/eLife.34114.036

Author response

Summary:

In this manuscript titled 'Transcriptomic analyses reveal rhythmic and CLOCK–driven pathways in human skeletal muscle', the authors performed RNA–seq analysis in muscle tissue biopsies (in vivo) and cultures of synchronized primary differentiated myofibers (in vitro) from human subjects, allowing them to compare cell autonomous and non–cell autonomous rhythms of transcription. Of note, the authors found that rhythmic transcription of genes involved in glucose transport (i.e. GLUT4 regulation), lipid metabolism and immune response were present in both the in vivo and in vitro samples. However, many other genes which were rhythmic in vivo did not appear rhythmic in vitro, suggesting that there are non–cell autonomous factors controlling rhythmic transcription in human muscle. In addition, since total RNA was used, both pre–mRNA and mature mRNAs were quantified. Indeed, the authors identified genes which oscillate at both the transcriptional and post–transcriptional levels. Finally, the dependence of rhythmic transcripts on the core molecular clock was determined using siRNA–mediated knock down of CLOCK in the primary myotube cultures. siCLOCK myotubes displayed reduced expression of many genes which were found to be rhythmicin vitro. The authors also determined that siCLOCK myotubes display reduced insulin–dependent and –independent glucose uptake, suggesting reduced activity of the GLUT4–transport mechanism. Overall, the finding that the circadian clock controls glucose and lipid metabolic gene expression in primary human myofibers is important and provides further evidence of the role of the molecular clock in the regulation of muscle function and glucose homeostasis.

The daily transcriptome analysis of the human muscle biopsies taken at 4 h intervals across 24 h from 10 healthy individuals serve as a valuable and important resource, especially in view of the current attempts to map circadian rhythmicity in humans. In this regard, a more rigorous comparison between available detests of mouse muscle transcriptome (e.g. McCarthy et al., 2007; Miller et al., 2007; Dyar et al., 2014; Zhang et al., 2014; Hodge et al., 2015) with this new data would provide interesting insight regarding the difference in muscle rhythmic gene expression between mice and humans, and likely between nocturnal and diurnal species, respectively.

We thank the reviewers for this insightful comment as indeed a comparison to the published mouse databases was lacking. Following the reviewers’ suggestion, we first compared genes rhythmic at the exonic level in human skeletal muscle based on our in vivo RNAseq data analysis to previously published transcriptional data performed on mouse skeletal muscle biopsies. We chose for comparison the datasets published by Zhang and colleagues (Zhang et al., 2014) representing temporal RNAseq and microarray analysis conducted in mouse skeletal muscle under a constant darkness regimen, and the work by Dyar and colleagues (Dyar et al., 2014) where the analysis was conducted by microarray under a 12h light–dark regimen. It should be noted that our study was done in human skeletal muscle biopsies obtained from a different muscle type as compared to the published rodent studies (vastus lateralis in our study, versus gastrocnemius, soleus and tibia anterior in the above–mentioned mouse studies). As shown by Venn diagram in Author response image 1, 107 genes were qualified as rhythmic in all three databases.

Zooming in on the specific differences between human and mouse data in skeletal muscle, we found that in human skeletal muscle rhythmic transcription was distributed into two sharp phases of transcript accumulation with the majority of genes peaking at 04:00, and a second smaller peak appearing at 16:00. This bimodal pattern of gene expression has been observed for rodents as well, where the majority of genes peak at CT23 (Zhang et al., 2014) or CT18 (McCarthy et al., 2007; Miller et al., 2007), and minor peaks of enrichment were observed at CT8 or CT8–10, respectively. These data suggest that an overall bimodal expression profile might be kept between rodents and humans, with the major peak time corresponding to late night – early morning hours, prior to the light onset. Keeping in mind that mice are nocturnal animals and humans are diurnal, this major transcriptional peak seems to correspond to the end of the active phase in mice and to the end of the rest phase in humans. The meaning of this difference stays to be explored but we cannot exclude that our experimental procedures also contributed to the observed difference.

With regard to specific genes, we observed in humans a peak of expression for genes related to muscle contraction and mitochondrial activity during the active phase (Supplemental dataset 2). For nocturnal animals, such as rodents, rhythmic genes related to muscle contraction however, do not exclusively peak during the active phase at night, but cluster into two peaks that correspond to the active and the resting phases (McCarthy et al., 2007; Zhang et al., 2014). McCarthy and colleagues report that two genes (Myh3, Actc1) related to muscle contraction peak at CT2 and one gene (Myh10) is enriched at CT14. According to CIRCA (Zhang et al., 2014), Myh2, Myh6, and Myh7 peak during the resting phase in mice, and Myh9 peaks towards the end of the active phase. Further differences were observed for genes involved in lipid metabolism. Elov5, which was enriched at the end of the resting phase in humans (04:00), as shown in our recent manuscript (Loizides–Mangold et al., 2017), reached its peak of expression in rodents at the end of the active period (CT22–24, (Zhang et al., 2014; Hodge et al., 2015)). Taken together these results suggest that in rodent skeletal muscle, rhythmic transcriptional activity might be differently regulated compared to human skeletal muscle with regard gene expression at the active/absorptive phase or resting phase.

Whereas the circadian phase for core–clock genes is kept between diurnal and nocturnal species at the SCN level, a phase–shift is observed in peripheral organs. Curiously, the observed shift in peripheral clocks among diurnal and nocturnal species is typically smaller than 12 hours and corresponds to rather 8–10 hours (Mure et al., 2018). In good agreement with previous studies, when comparing the circadian phases of core–clock genes in our database to the temporal profiles of the corresponding rodent data (Dyar et al., 2014; Zhang et al., 2014), we observed a similar phase shift of 8–10 hours. The question remains at what level such a phaseshift between nocturnal and diurnal species occurs, and why it is not exactly 12 h in peripheral organs (Mure et al., 2018). Further comparative studies, conducted in the same type of tissue and the same methodology, will be required to explore this fundamental issue.

Author response image 1
Comparison of our in vivo vastus lateralis rhythmic dataset with mouse skeletal muscle rhythmic data (Dyar et al., 2014; Zhang et al., 2014).

Blue = our in vivo RNAseq data; red = (Dyar et al., 2014) soleus; green = (Zhang et al., 2014). The numbers presented in the Venn diagram correspond to the total number of rhythmically expressed genes for each segment (exons only). For the human in vivo dataset (our work), this corresponds to R–E + R–I.R–E, with the cutoff described in the manuscript (Figure 1). For (Dyar et al., 2014) (microarray) and (Zhang et al., 2014) (microarray) rhythmic genes were identified using a Benjamini–Hochberg q–value <0.2 in the JTK_Cycle algorithm.

https://doi.org/10.7554/eLife.34114.029

We therefore revised the corresponding sentence in the Results and Discussion sections as follows:

Results section: "As previously reported in mice, rhythmic gene transcription was distributed into two phases of transcript accumulation (04:00 and 16:00, Figure 1E). The afternoon peak (16:00) was enriched in genes related to muscle contraction and mitochondrial activity (Figure 1–source data 2) whereas homologous genes in rodents are shared between the active but also the resting phase (McCarthy et al., 2007; Miller et al., 2007; Hodge et al., 2015).”

Discussion section: “On the other hand, a comparison of our in vivo dataset (exonic signals) to published temporal gene expression databases of mouse skeletal muscle (Dyar et al., 2014; Zhang et al., 2014) revealed 107 common rhythmic genes between mouse and human skeletal muscle. When comparing the circadian phases of core–clock components in our database to the temporal profiles of the corresponding genes in rodents (Dyar et al., 2014; Zhang et al., 2014), we observed a phase shift of 8 – 10 hr. This result is in good agreement with a phase shift observed between peripheral clocks in nocturnal versus diurnal species, which is indeed typically smaller than 12 h (Mure et al., 2018). The question remains at what level such a phase–shift between nocturnal and diurnal species occurs, and why it is not exactly 12 h in peripheral organs (Mure et al., 2018). Further comparative studies, conducted in the same type of tissue and with the same methodology, will be required to explore this fundamental issue.”

The major strength of this study lies on the expression data obtained from human healthy individuals in a well–controlled manner, this section by itself is an important resource.

The in vitro experiments primarily corroborate previous studies done with cultured cells (e.g. Krishnaiah et al., 2017) that show dramatic reduction in rhythmic gene expression in culture compared to organs in vivo.

Essential revisions:

1) The in vivo biopsy samples were taken from the vastus lateralis muscle whereas the in vitro samples were from gluteus maximus. Therefore, is it possible that the observed gene expression differences observed between the in vitro and in vivo sample sets are explained by the fact that they are derived from different muscles? One option to address this would be to compare gene expression profiles from biopsies of VL and GM muscles – do they differ in a similar way to the in vitro versus in vivo gene sets? This may beyond the scope here.

We entirely agree with the reviewers that the muscle origin might have contributed to the observed gene expression differences. To address this important point, we have previously quantified the expression of myosin isoforms MYH1, MYH2, MYH4 and MYH7 in cDNA samples obtained from human gluteus maximus and from vastus lateralis biopsies (Loizides–Mangold et al., 2017). As expected, significant differences in myosin isoform composition between vastus lateralis and gluteus maximus were observed. The dominant myosin isoform in vastus lateralis was MYH2 followed by low levels of MYH1 and MYH7, indicative of mixed fast 2A/2X fibers. The dominant isoforms for gluteus maximus were MYH2 and MYH7, followed by MYH1 indicative of mixed slow and fast 2A fibers. The expression of MYH4, specific for the fastest 2B fiber type, was very low (Harrison et al., 2011).

To further emphasize this valuable point made by the reviewers, we have added the following sentences to the Discussion section:

“Moreover, we cannot exclude that discrepancies between the in vivoand in vitro circadian datasets are in part also influenced by the fiber type composition of vastus lateralis and gluteus maximus, as demonstrated by myosin isoform analysis (Loizides–Mangold et al., 2017).”

To further dissect the impact of the skeletal muscle source on skeletal muscle gene expression, we compared the number of genes commonly expressed in our in vitro and in vivo RNAseq datasets, with published genomic data on human biopsies from quadriceps and gluteus maximus muscle (Eisenberg et al., 2008)., We detected that the number of genes differentially expressed between our in vivo and in vitro datasets is substantially higher than the number of genes differentially expressed between the two biopsies originating from gluteus maximus and quadriceps. Indeed, in the published dataset (Eisenberg et al., 2008) a total of 954 genes (875+79) were differentially expressed. This number was significantly lower compared to the number of genes differentially expressed between skeletal myotubes derived from gluteus maximus tissue biopsies and vastus lateralis tissue. However, 875 of these genes were shared between the published data, and the genes that were differentially expressed in vitro and in vivo in our datasets comparison (gluteus in vitro vs vastus lateralis in vivo). This overlap is significant and cannot be explained by a random overlap between two datasets–. That means that there is an effect of the muscle type origin on the transcriptomic profile of gluteus maximus in vitro vs vastus lateralis in vivo. Given the overwhelming number of genes differentially expressed between gluteus maximus (in vivo) and vastus lateralis (in vitro), it is tempting to speculate that the in vitro transition has an even stronger impact on gene expression than the origin of the muscle tissue biopsy. However, the effect of the muscle type on gene expression is significant, and should be taken into account when interpreting the data.

2) The authors observed gene expression changes with siCLOCK that support increased mitochondrial oxidative/type I muscle remodeling. This result is consistent with what was observed in RNA–sequencing studies in Bmal1KO mouse muscle (Hodge et al., 2015). However, others have provided evidence for reduced mitochondrial oxidative metabolism in Bmal1KO muscle cells. Therefore, it is possible that the gene signature does not equate with the metabolic phenotype. For instance, gene expression changes may be compensatory rather than causal? In the absence of such analysis, discussion should acknowledge the interest in future functional bioenergetic profiling (e.g., respirometry with mitochondrial fuel substrates).

Indeed, a recent publication by the Bass group (Peek et al., 2017) has shown that disruption of Bmal1 in skeletal myotubes and in fibroblasts reduces mitochondrial respiration and anaerobic glycolysis through interaction of the circadian clock with the Hypoxia Inducible Factor (HIF) pathway. This finding indicates that rhythmic oxygen levels can reset circadian clocks, something that has been also suggested by the Asher lab (Adamovich et al., 2017). We agree with the reviewers that it cannot be excluded that the observed increase in slow twitch muscle gene expression (Supplementary file 1–table S7), does not necessarily equate to increased mitochondrial activity. Moreover, in agreement with the data published by Peek et al., we also observe a reduction in the expression of the HIF target gene VEGFA upon clock disruption (Supplementary file 1), supporting the hypothesis that clock – HIF interactions play an important role in the glycolytic capacity of skeletal muscle.

Following this constructive remark, we have now modified the corresponding paragraph in the Discussion section, emphasizing the impact of oxygen sensing on skeletal muscle function in addition to muscle type I remodeling:

"[…]we identified an upregulation of multiple genes characteristic for type I slow fibers, as well as a downregulation of genes associated with type II fast fibers (Supplementary file 1). Whether these changes affect mitochondrial activity stays to be further explored. Importantly, a recent study has demonstrated that Bmal1 KO myotubes display reduced mitochondrial respiration and a reduced expression of Hipoxia Inducible Factor (HIF) target genes (Peek et al., 2017). In agreement, we also observed reduced expression of the HIF target gene VEGFA upon clock disruption (Supplementary file 1), supporting the hypothesis that clock – HIF interactions play an important role in the glycolytic capacity of skeletal muscle."

3) The mathematical model used to detect rhythmicity should be validated on an existing dataset or a demo dataset in order to compare their results with established and widely used algorithms that detect rhythmic profiles, such as JTK_CYCLE (Hughes et al., 2010) and Meta Cycle (Wu et al., 2016) that incorporates ARSER, JTK_CYCLE and Lomb–Scargle.

The well–known harmonic regression (COSINOR, (Nelson et al., 1979) was used to assess rhythmicity of the in vivo RNA–seq data. In order to deal with the fact that time points are repeated measures from the same patient, we used a mixed linear model. This includes the usual fixed effects from the harmonic regression model (or COSINOR) and a random effect (patient–specific) on the intercept that deals with the subject–to–subject variation and dependency of the repeated measures. This allows to, not violating the assumption of independent error terms behind the linear regression theory. A similar method has been used to characterize rhythmic transcript in human blood samples (Arnardottir et al., 2014).

Following the pertinent suggestion by the reviewers regarding the rhythmic analysis of the in vitro samples, we evaluated the temporal expression profiles of key core clock genes, such as ARNTL, NR1D1, NR1D2, PER1, PER2, PER3, CRY1, CRY2, NPAS2, TEF and BHLHE41, extracted from two published large–scale time series (Atger et al., 2015; Petrenko et al., 2017). The expression profiles of all genes listed above, extracted from these two databases, were identified by the algorithm developed and employed by us for the in vitroanalysis as circadian.

We further evaluated our in vitroalgorithm on the model dataset cycMouseLiverRNA from the MetaCycle R package. 20 core–clock and functional genes, listed in Revision table 1, were classified as circadian by our algorithm. The algorithm maps time series gene expression to one of the 4 categories: CIRCA: circadian; LIN: linear; ONE: one peak followed by a flat expression profile; NUA: scattered (from nuage, French word for cloud). The algorithm provides a probability for each category, with the gene expression profile being assigned to the category with the highest probability, presented in the right–most column.

Revision Table 1. Evaluation of core–clock and functional genes from the model dataset cycMouseLiverRNA by our in vitro algorithm. Note that all the genes qualified as circadian by the MetaCycle R package were also confirmed by our analysis method.

NameCIRCALINNUAONEPrediction
Avpr1a_1418603_at0.86200.0860.052CIRCA
Cirbp_1416332_at0.8400.0840.076CIRCA
Clock_1418659_at0.8600.0860.054CIRCA
Elovl5_1437211_x_at0.86200.0860.052CIRCA
Fkbp5_1448231_at0.86600.0860.048CIRCA
Hist1h1c_1416101_a_at0.86800.0860.046CIRCA
Hsp90aa1_1437497_a_at0.7640.0060.1920.038CIRCA
Lipg_1421262_at0.79200.0760.132CIRCA
Nr1d1_1426464_at0.8600.0860.054CIRCA
Nr1d2_1416958_at0.8600.0860.054CIRCA
Nr1h3_1450444_a_at0.86400.10.036CIRCA
Per1_1449851_at0.87600.080.044CIRCA
Per2_1417602_at0.82400.0760.1CIRCA
Rgs16_1426037_a_at0.86600.0860.048CIRCA
Rorc_1425792_a_at0.8600.0860.054CIRCA
Scap_1433520_at0.86200.0860.052CIRCA
Tef_1424175_at0.8600.0860.054CIRCA
Tsc22d3_1420772_a_at0.86200.0860.052CIRCA
Tspan4_1448276_at0.8600.0860.054CIRCA
Tubb5_1416256_a_at0.86200.0860.052CIRCA

Of note, the commonly used JTK_CYCLE algorithm (Hughes et al., 2010) did not prove reliable for our data set, most probably due to the large amplitude differences between the two circadian cycles, which were observed in our experimental settings. A high number of genes qualified as “circadian” by JTK_CYCLE algorithm, were clearly false–positives. The example profiles of the gene expression profiles qualified as circadian by JTK_CYCLE are presented in Author response image 3.

Author response image 3
Examples of the gene expression profiles from our database qualified as circadian by JTK_CYCLE.
https://doi.org/10.7554/eLife.34114.031

4) It would be interesting to compare the genes identified as rhythmic in cultured U2OS cells (Krishnaiah et al., 2017) with the ones identified in the current study, even though these are completely different cells. It might shed light on genes, aside from the core clock, that globally maintain rhythmicity in culture.

We thank the reviewers for this constructive comment and have subsequently compared our dataset to the results published for U2OS cells by (Krishnaiah et al., 2017). Among the 190 genes that were rhythmic in human skeletal myotubes, 30 were also rhythmic in U2OS cells. Among those were members of the core–clock machinery such as ARNTL, CRY1–2, NR2D1–2, TEF, PER1–3, and multiple genes involved in cell cycle progression and mitosis. Beyond those also the sulfatase ARSJ, involved in glycosphingolipid metabolism, EXOSC8, a regulator of mRNA stability, LRRC16A, involved in actin filament organization, and TUBA1C, a structural constituent of the cytoskeleton, as well as E2F1 were rhythmic in both skeletal myotubes and U2OS cells. Interestingly, the RING finger domain protein–encoding gene TRIM47 exhibited a circadian expression in human skeletal myotubes, in human muscle biopsies, and in U2OS cells. The candidate genes that we highlighted in our study for being highly relevant to muscle physiology, and which exhibited circadian profiles in both synchronized myotubes and in skeletal muscle biopsies (KLF11, GFPT2, NAMPT, PLAU, PLD1, PIM1), were all significantly expressed in U2OS cells. However, those genes did not exhibit rhythmic profiles in U2OS cells. KLF11 had a rhythmic tendency, although it did not reach statistical significance to qualify as circadian according to JTK–cycle (adjusted p–value of 0.266437).

To incorporate this information the following part was added to the Discussion section:

“Comparison of our in vitro dataset to the results published on U2OS cells (Hughes et al., 2009; Krishnaiah et al., 2017) revealed that among the 190 genes that were rhythmic in human skeletal myotubes, 30 were also rhythmic in U2OS. Among those were members of the core–clock machinery, and multiple genes involved in cell cycle progression and mitosis. Moreover, the sulfatase ARSJ, involved in glycosphingolipid metabolism, EXOSC8, a regulator of mRNA stability, LRRC16A, involved in actin filament organization, and TUBA1C, encoding for a structural constituent of the cytoskeleton, as well as E2F1 were rhythmic in skeletal myotubes and in U2OS cells. Interestingly, the RING finger domain protein–encoding gene TRIM47 exhibited a rhythmic expression in human skeletal myotubes, in human muscle biopsies, and in U2OS cells.

5) The authors conclusion that the cell–autonomous circadian clock has an essential role in coordinating muscle glucose homeostasis and lipid metabolism in humans should be revised as it is not supported by their results. Their finding that knockdown of clock in vitro affected the overall expression of ~8% of all genes with genes related to glucose and lipid metabolism suggests that the CLOCK protein itself modulate their expression and not necessarily the oscillator. Furthermore, the result that only few genes maintain rhythmicity both in vivo and in vitro, most of them are core clock components, does not support it either.

The authors thank the reviewers for this insightful commentary. Indeed, the current interpretation of the in vitro results might be misleading. We subsequently changed the respective parts in the Discussion section and toned down the message regarding the importance of the cell–autonomous regulation of human skeletal muscle and highlighting the distinction between circadian clock–related and unrelated effects of CLOCK gene knockdown. The revised paragraphs read as following:

Abstract: “Our findings suggest an essential role for the circadian coordination of skeletal muscle glucose homeostasis and lipid metabolism in humans.”

Discussion section: “Moreover, we demonstrate that CLOCK depletion in cultured primary skeletal myotubes led to significant changes in gene expression (Figure 2), and related physiological outputs, comprising the regulation of basal and insulin stimulated glucose uptake, lipid homeostasis (Figure 3), and myokine secretion, as summarized in Figure 6.”

Discussion section: “Though our established experimental system for cellular clock disruption mediated by efficient CLOCK knockdown proved highly useful to study transcriptional and functional outputs in cultured human primary cells, one should keep in mind that core–clock genes also perform clock unrelated functions. The same holds true for genetic mouse models, where different core–clock gene KO strains exhibit distinct phenotypes. To discriminate between clock–related and unrelated effects of the CLOCK gene knockdown, alternative methods for the circadian clock perturbation will be required.”

https://doi.org/10.7554/eLife.34114.037

Article and author information

Author details

  1. Laurent Perrin

    1. Division of Endocrinology, Diabetes, Hypertension and Nutrition, Department of Internal Medicine Specialties, University Hospital of Geneva, Geneva, Switzerland
    2. Department of Cell Physiology and Metabolism, Faculty of Medicine, University of Geneva, Geneva, Switzerland
    3. Diabetes Center, Faculty of Medicine, University of Geneva, Geneva, Switzerland
    4. Institute of Genetics and Genomics of Geneva, Geneva, Switzerland
    Contribution
    Conceptualization, Data curation, Formal analysis, Validation, Investigation, Visualization, Methodology, Writing—original draft, Writing—review and editing
    Competing interests
    No competing interests declared
  2. Ursula Loizides-Mangold

    1. Division of Endocrinology, Diabetes, Hypertension and Nutrition, Department of Internal Medicine Specialties, University Hospital of Geneva, Geneva, Switzerland
    2. Department of Cell Physiology and Metabolism, Faculty of Medicine, University of Geneva, Geneva, Switzerland
    3. Diabetes Center, Faculty of Medicine, University of Geneva, Geneva, Switzerland
    4. Institute of Genetics and Genomics of Geneva, Geneva, Switzerland
    Contribution
    Conceptualization, Data curation, Formal analysis, Supervision, Validation, Investigation, Methodology, Writing—original draft, Project administration, Writing—review and editing
    Contributed equally with
    Stéphanie Chanon and Cédric Gobet
    Competing interests
    No competing interests declared
    ORCID icon 0000-0001-9233-2974
  3. Stéphanie Chanon

    CarMeN Laboratory, INSERM U1060, Oullins, France
    Contribution
    Data curation, Investigation, Visualization, Methodology
    Contributed equally with
    Ursula Loizides-Mangold and Cédric Gobet
    Competing interests
    No competing interests declared
  4. Cédric Gobet

    1. Nestlé Institute of Health Sciences, Lausanne, Switzerland
    2. School of Life Sciences, Ecole Polytechnique Fédérale de Lausanne, Lausanne, Switzerland
    Contribution
    Resources, Data curation, Software, Formal analysis, Validation, Investigation, Methodology
    Contributed equally with
    Ursula Loizides-Mangold and Stéphanie Chanon
    Competing interests
    Is a full-time employee of the Nestlé Institute of Health Sciences SA.
  5. Nicolas Hulo

    1. Institute of Genetics and Genomics of Geneva, Geneva, Switzerland
    2. Service for Biomathematical and Biostatistical Analyses, Institute of Genetics and Genomics in Geneva, University of Geneva, Geneva, Switzerland
    Contribution
    Resources, Data curation, Software, Formal analysis, Supervision, Methodology
    Competing interests
    No competing interests declared
    ORCID icon 0000-0003-2640-636X
  6. Laura Isenegger

    Service for Biomathematical and Biostatistical Analyses, Institute of Genetics and Genomics in Geneva, University of Geneva, Geneva, Switzerland
    Contribution
    Data curation, Software, Formal analysis, Validation
    Competing interests
    No competing interests declared
  7. Benjamin D Weger

    Nestlé Institute of Health Sciences, Lausanne, Switzerland
    Contribution
    Investigation, Methodology
    Competing interests
    Benjamin D Weger: Is a full-time employee of the Nestlé Institute of Health Sciences SA.
  8. Eugenia Migliavacca

    Nestlé Institute of Health Sciences, Lausanne, Switzerland
    Contribution
    Data curation, Investigation, Methodology
    Competing interests
    Eugenia Migliavacca: Is a full-time employee of the Nestlé Institute of Health Sciences SA.
  9. Aline Charpagne

    Nestlé Institute of Health Sciences, Lausanne, Switzerland
    Contribution
    Data curation, Software, Formal analysis, Methodology
    Competing interests
    Aline Charpagne: Is a full-time employee of the Nestlé Institute of Health Sciences SA.
  10. James A Betts

    Department for Health, University of Bath, Bath, United Kingdom
    Contribution
    Data curation, Formal analysis, Investigation, Methodology, Writing—review and editing
    Competing interests
    Has been a consultant for PepsiCo (Quaker) and Kellogg's.
  11. Jean-Philippe Walhin

    Department for Health, University of Bath, Bath, United Kingdom
    Contribution
    Data curation, Formal analysis, Investigation, Methodology, Writing—review and editing
    Competing interests
    No competing interests declared
  12. Iain Templeman

    Department for Health, University of Bath, Bath, United Kingdom
    Contribution
    Data curation, Formal analysis, Investigation, Methodology, Writing—review and editing
    Competing interests
    No competing interests declared
  13. Keith Stokes

    Department for Health, University of Bath, Bath, United Kingdom
    Contribution
    Data curation, Formal analysis, Investigation, Methodology, Writing—review and editing
    Competing interests
    No competing interests declared
  14. Dylan Thompson

    Department for Health, University of Bath, Bath, United Kingdom
    Contribution
    Data curation, Formal analysis, Investigation, Methodology, Writing—review and editing
    Competing interests
    No competing interests declared
  15. Kostas Tsintzas

    MRC/ARUK Centre for Musculoskeletal Ageing, School of Life Sciences, University of Nottingham, Nottingham, United Kingdom
    Contribution
    Data curation, Formal analysis, Investigation, Methodology, Writing—review and editing
    Competing interests
    No competing interests declared
  16. Maud Robert

    Department of Digestive and Bariatric Surgery, Edouard Herriot University Hospital, Lyon, France
    Contribution
    Validation, Investigation, Methodology, Writing—review and editing
    Competing interests
    No competing interests declared
  17. Cedric Howald

    1. Institute of Genetics and Genomics of Geneva, Geneva, Switzerland
    2. Department of Genetic Medicine and Development, Faculty of Medicine, University of Geneva, Geneva, Switzerland
    Contribution
    Resources, Data curation, Software, Formal analysis, Validation
    Competing interests
    No competing interests declared
  18. Howard Riezman

    Department of Biochemistry, NCCR Chemical Biology, University of Geneva, Geneva, Switzerland
    Contribution
    Conceptualization, Funding acquisition, Investigation, Writing—review and editing
    Competing interests
    No competing interests declared
    ORCID icon 0000-0003-4680-9422
  19. Jerome N Feige

    1. Nestlé Institute of Health Sciences, Lausanne, Switzerland
    2. School of Life Sciences, Ecole Polytechnique Fédérale de Lausanne, Lausanne, Switzerland
    Contribution
    Data curation, Formal analysis, Supervision, Validation, Investigation, Methodology, Writing—review and editing
    Competing interests
    Jerome N Feige: Is a full-time employee of the Nestlé Institute of Health Sciences SA.
  20. Leonidas G Karagounis

    1. Experimental Myology and Integrative Biology Research Cluster, Faculty of Sport and Health Sciences, University of St Mark and St John, Plymouth, United Kingdom
    2. Institute of Nutritional Science, Nestlé Research Centre, Lausanne, Switzerland
    Contribution
    Conceptualization, Data curation, Formal analysis, Supervision, Validation, Methodology, Writing—review and editing
    Competing interests
    Is an employee of Nestec Ltd.
  21. Jonathan D Johnston

    Faculty of Health and Medical Sciences, University of Surrey, Guildford, United Kingdom
    Contribution
    Data curation, Formal analysis, Supervision, Validation, Investigation, Methodology, Writing—review and editing
    Competing interests
    No competing interests declared
    ORCID icon 0000-0001-8083-9794
  22. Emmanouil T Dermitzakis

    1. Institute of Genetics and Genomics of Geneva, Geneva, Switzerland
    2. Department of Genetic Medicine and Development, Faculty of Medicine, University of Geneva, Geneva, Switzerland
    Contribution
    Conceptualization, Resources, Data curation, Software, Supervision, Validation
    Competing interests
    No competing interests declared
  23. Frédéric Gachon

    1. Nestlé Institute of Health Sciences, Lausanne, Switzerland
    2. School of Life Sciences, Ecole Polytechnique Fédérale de Lausanne, Lausanne, Switzerland
    Contribution
    Conceptualization, Data curation, Software, Formal analysis, Supervision, Validation, Investigation, Methodology, Writing—review and editing
    Contributed equally with
    Etienne Lefai
    Competing interests
    No competing interests declared
    ORCID icon 0000-0002-9279-9707
  24. Etienne Lefai

    CarMeN Laboratory, INSERM U1060, Oullins, France
    Contribution
    Conceptualization, Data curation, Formal analysis, Supervision, Funding acquisition, Validation, Investigation, Methodology, Writing—original draft, Project administration, Writing—review and editing
    Contributed equally with
    Frédéric Gachon
    Competing interests
    No competing interests declared
  25. Charna Dibner

    1. Division of Endocrinology, Diabetes, Hypertension and Nutrition, Department of Internal Medicine Specialties, University Hospital of Geneva, Geneva, Switzerland
    2. Department of Cell Physiology and Metabolism, Faculty of Medicine, University of Geneva, Geneva, Switzerland
    3. Diabetes Center, Faculty of Medicine, University of Geneva, Geneva, Switzerland
    4. Institute of Genetics and Genomics of Geneva, Geneva, Switzerland
    Contribution
    Conceptualization, Supervision, Funding acquisition, Validation, Investigation, Methodology, Writing—original draft, Project administration, Writing—review and editing
    For correspondence
    Charna.dibner@hcuge.ch
    Competing interests
    No competing interests declared
    ORCID icon 0000-0002-4188-803X

Funding

Schweizerischer Nationalfonds zur Förderung der Wissenschaftlichen Forschung (Sinergia Grant No. CRSII3-154405)

  • Howard Riezman
  • Etienne Lefai
  • Charna Dibner

Schweizerischer Nationalfonds zur Förderung der Wissenschaftlichen Forschung (31003A-166700)

  • Charna Dibner

Fondation Privée HUG

  • Emmanouil T Dermitzakis
  • Charna Dibner

Fondation Ernst et Lucie Schmidheiny

  • Charna Dibner

Société Académique de Genève

  • Charna Dibner

Biotechnology and Biological Sciences Research Council

  • Jonathan D Johnston

The funders had no role in study design, data collection and interpretation, or the decision to submit the work for publication.

Acknowledgements

The authors thank Jacques Philippe and Sylvain Sardy for constructive discussions; Pamela Pulimeno for help with conducting the in vitro around-the-clock sample collection, RNA extractions and RT-qPCR experiments; Svetlana Skarupelova and Camille Saini for help with RNA extractions and Anne-Marie Makhlouf for lentivirus preparation; Luciana Romano and Deborah Beilser for their help with in vitro RNA sequencing (University of Geneva); Marc Moniatte and Jonathan Paz Montoya (PCF, EPFL) for help with the lipid analyses; and Ondine Walter (NIHS) for the ethical and logistical management of human samples. This work was funded by the Sinergia Swiss National Science Foundation (Grant No. CRSII3-154405 to HR, CD, EL), the Swiss National Science Foundation (Grant No. 31003A-166700 (CD), the Fondation Privée de HUG, Fondation Ernst et Lucie Schmidheiny, the Société Académique de Genève (CD) and by the United Kingdom Biotechnology and Biological Sciences Research Council Grant BB/I008470/1 (JDJ).

Ethics

Human subjects: The in vivo study was conducted in accordance with the Declaration of Helsinki and with approval of the Health Research Authority (NRES Committee South West; 14/SW/0123) and the Swiss Commission Cantonal (Canton Vaud) d'éthique de la recherche (Cer-VD). Muscle biopsies were obtained from donors with their informed consent. The experimental protocol ('DIOMEDE') for the in vitro study was approved by the Ethical Committee SUD EST IV (Agreement 12/111) and performed according to the French legislation (Huriet's law).

Reviewing Editor

  1. Joseph S Takahashi, Howard Hughes Medical Institute, University of Texas Southwestern Medical Center, United States

Publication history

  1. Received: December 6, 2017
  2. Accepted: April 4, 2018
  3. Version of Record published: April 16, 2018 (version 1)

Copyright

© 2018, Perrin et al.

This article is distributed under the terms of the Creative Commons Attribution License, which permits unrestricted use and redistribution provided that the original author and source are credited.

Metrics

  • 1,645
    Page views
  • 302
    Downloads
  • 3
    Citations

Article citation count generated by polling the highest count across the following sources: Scopus, Crossref, PubMed Central.

Download links

A two-part list of links to download the article, or parts of the article, in various formats.

Downloads (link to download the article as PDF)

Download citations (links to download the citations from this article in formats compatible with various reference manager tools)

Open citations (links to open the citations from this article in various online reference manager services)

Further reading

    1. Human Biology and Medicine
    2. Neuroscience
    Ming Song et al.
    Research Article Updated
    1. Cancer Biology
    2. Human Biology and Medicine
    Haley Hieronymus et al.
    Research Article Updated