Response to Comment on “Circadian rhythms in the absence of the clock gene Bmal1”

Response to Comment on “Circadian rhythms in the absence of the clock gene Bmal1”
Spread the love

Abstract

Ness-Cohn et al. claim that our observations of transcriptional circadian rhythms in the absence of the core clock gene Bmal1 in mouse skin fibroblast cells are supported by inadequate evidence. They claim that they were unable to reproduce some of the original findings with their reanalysis. We disagree with their analyses and outlook.

In their Comment, Ness-Cohn et al. (1) reanalyzed some of the contents of our paper, which reported 24-hour oscillations of the transcriptome, proteome, and phosphoproteome in skin fibroblasts and liver slices lacking Bmal1 (2). The authors focused on transcriptome oscillations found in mouse skin fibroblasts (MSFs). They speculate that cycling genes detected in Bmal1 knockouts are due to “noise” rather than circadian oscillations. However, we disagree with their position and discuss how selected parameters affect the analysis.

The numbers shown in figure 1E of Ray et al. (2) reflect the results of the 37°C rhythmic transcriptome data that were generated as a part of the temperature compensation experiment. During the review process, these additional experiments were performed and provided more consistent results than the initial experiment, likely attributable to the better (2-hour) sampling resolution used in the 37°C temperature compensation experiment. The authors of this Comment will be aware of such sampling resolution effects, given their recent work on the topic (3). Thus, the number of cycling transcripts was updated according to those data, while the heat maps were plotted for the corresponding transcript profiles in the original 72-hour time-course experiment (to keep them comparable with heat maps for the liver tissue time course). Unfortunately, during the revision of the paper, we did not make this clear in the legend. We apologize for this confusion.

The number of rhythmic transcripts (897 in Bmal1–/–, FDR et al. in the 37°C dataset (2-hour resolution) is very similar to the number we detected [see the bar plots in figure 1B of (1)], as also shown in figure 1E of (2). Indeed, we cannot see any notable difference in the number of identified rhythmic transcripts in their and our analysis for any of our MSF datasets (Fig. 1A). Moreover, in their own analysis, Ness-Cohn et al. detected thousands of rhythmic transcripts in the other Bmal1 knockout datasets at FDR 2000 rhythmic transcripts in AM and PM datasets in Bmal1 knockouts [see figure 1B of (1)]. We do not think such a large number of rhythmic transcripts is likely to be attributable to “noise.”

Fig. 1 Analysis of Bmal1 knockout transcriptomics data using various parameters, and the overlap of the different datasets.

(A) Comparative analysis of the number of identified rhythmic transcripts (RAIN independent method, FDR et al. (1) and Ray et al. (2). All the data processing and analysis parameters are the same in both pipelines, except that we used a different peak_border (0.1, 0.9), which provided a slightly different (usually lower) number of rhythmic transcripts relative to Ness-Cohn et al. (B) The number of rhythmic transcripts detected in Bmal1–/– fibroblasts at 37°C (2-hour sampling, log2-transformed FPKM, RAIN, FDR et al. (initial version of their analysis) and various FPKM thresholds after log2 transformation (RAIN longitudinal method). Comparison to Ray et al. parameters [all FPKM> 0, RAIN independent method, peak_border (0.1, 0.9)] is shown. (C) Overlap among the rhythmic transcripts (RAIN, FDR Bmal1+/+ (wild type, WT) and Bmal1–/– MSF experiments performed at 37°C. Bmal1–/– MSF datasets show a substantial number of transcripts (521/3790) rhythmic in at least two of four datasets. (D) Phase difference distribution for the rhythmic transcripts (RAIN, FDR Bmal1+/+ and Bmal1–/– fibroblasts; 83 rhythmic transcripts in AM and PM datasets exhibited antiphase rhythms in Bmal1–/– MSFs. (E) Oppositely phased abundance profiles (log2 transformed FPKM) of a few representative rhythmic genes not shown in Ray et al. (RAIN, FDR F) LimoRhyde analysis indicating differential rhythmic transcriptome in AM and PM datasets.

” data-hide-link-title=”0″ data-icon-position=”” href=”https://science.sciencemag.org/content/sci/372/6539/eabf1930/F1.large.jpg?width=800&height=600&carousel=1″ rel=”gallery-fragment-images-195886060″ title=”Analysis of Bmal1 knockout transcriptomics data using various parameters, and the overlap of the different datasets. (A) Comparative analysis of the number of identified rhythmic transcripts (RAIN independent method, FDR 0, RAIN independent method, peak_border (0.1, 0.9)] is shown. (C) Overlap among the rhythmic transcripts (RAIN, FDR

<figcaption id="F1-caption">
<span>Fig. 1</span> <span>Analysis of&nbsp;<em>Bmal1</em>&nbsp;knockout transcriptomics data using various parameters, and the overlap of the different datasets.</span><p id="p-6">(<strong>A</strong>) Comparative analysis of the number of identified rhythmic transcripts (RAIN independent method, FDR et al. (<em>1</em>) and Ray <em>et al</em>. (<em>2</em>). All the data processing and analysis parameters are the same in both pipelines, except that we used a different peak_border (0.1, 0.9), which provided a slightly different (usually lower) number of rhythmic transcripts relative to Ness-Cohn <em>et al</em>. (<strong>B</strong>) The number of rhythmic transcripts detected in <em>Bmal1</em><sup>&ndash;/&ndash;</sup> fibroblasts at 37&deg;C (2-hour sampling, log<sub>2</sub>-transformed FPKM, RAIN, FDR et al. (initial version of their analysis) and various FPKM thresholds after log<sub>2</sub> transformation (RAIN longitudinal method). Comparison to Ray <em>et al</em>. parameters [all FPKM&gt; 0, RAIN independent method, peak_border (0.1, 0.9)] is shown. (<strong>C</strong>) Overlap among the rhythmic transcripts (RAIN, FDR Bmal1<sup>+/+</sup> (wild type, WT) and <em>Bmal1</em><sup>&ndash;/&ndash;</sup> MSF experiments performed at 37&deg;C. <em>Bmal1</em><sup>&ndash;/&ndash;</sup> MSF datasets show a substantial number of transcripts (521/3790) rhythmic in at least two of four datasets. (<strong>D</strong>) Phase difference distribution for the rhythmic transcripts (RAIN, FDR Bmal1<sup>+/+</sup> and <em>Bmal1</em><sup>&ndash;/&ndash;</sup> fibroblasts; 83 rhythmic transcripts in AM and PM datasets exhibited antiphase rhythms in <em>Bmal1</em><sup>&ndash;/&ndash;</sup> MSFs. (<strong>E</strong>) Oppositely phased abundance profiles (log<sub>2</sub> transformed FPKM) of a few representative rhythmic genes not shown in Ray <em>et al</em>. (RAIN, FDR F) LimoRhyde analysis indicating differential rhythmic transcriptome in AM and PM datasets.</p>  </figcaption>

In the initial version of their analysis (see http://tiny.cc/gntutz), Ness-Cohn et al. analyzed our data without using any FPKM threshold. Thresholding is a standard practice in RNA-seq data analysis and can be executed by a number of methods, including removal of zero-containing transcripts or the use of a cutoff above a mean expression level for each transcript (such as mean FPKM> 0.5) (4). In Ray et al. (2), data were filtered by removing any transcript that had a zero in any of the time points. This was comparable to various mean FPKM thresholds that could have been used (Fig. 1B) (1). The absence of this standard data-processing step led to differences in the numbers of transcripts identified in all of Ness-Cohn et al.’s initial analyses (see their code at http://tiny.cc/dntutz).

Ness-Cohn et al. compared FDR-adjusted P values between different experiments using rank correlation, and they suggest that because there is low correlation among the datasets, this indicates low reproducibility. This methodology is based on their own work investigating the reproducibility of different rhythmicity detection algorithms (e.g., JTK_cycle versus RAIN) in circadian transcriptome datasets using uncorrected P values (3). This correlation metric has not been validated for between-study comparisons, nor for the comparison of FDR-adjusted P values (as opposed to uncorrected ones), which is how they use it here. Also, they did not compare the Ray et al. data with other datasets (i.e., datasets with good concordance by their own criteria). This is important because there is a large body of literature indicating minimal overlap among the genome-scale circadian datasets (5).

We disagree with Ness-Cohn et al.’s concern about the numbers of transcripts that are rhythmic among the different 37°C datasets. To define consistent transcripts across datasets, we used an alternative metric, requiring that a transcript must be rhythmic in at least two datasets. In our analysis, we identified hundreds of transcripts as rhythmic (RAIN FDR

We disagree with the authors’ observation that there are only very few circadian transcripts that exhibit nearly opposite phases in Bmal1–/– MSF “AM” and “PM” datasets. Again, after FPKM-thresholding the data, n=83 transcripts (RAIN, FDR Bmal1–/– MSF AM/PM datasets (Fig. 1, D and E). We also analyzed the AM and PM rhythmic transcriptome datasets using LimoRhyde (6). This framework assesses the extent to which two genome-scale rhythmic datasets differ from one another. LimoRhyde indicates a substantial difference between AM and PM datasets in Bmal1–/– MSFs (438 transcripts at FDR

As far as we are aware, analysis of temperature compensation and/or opposite entrainment using RNA-seq has not been performed before. Therefore, we always analyzed Bmal1+/+ (i.e., wild-type) MSFs in every experiment as a reference and observed levels of consistency among the datasets similar to what we observed in Bmal1–/– cells (2). Consequently, we see no justification in Ness-Cohn et al.’s claim that the rhythmic transcriptome in Bmal1–/– MSFs is “noise” when it is accepted that there are oscillations in Bmal1+/+ MSFs. Furthermore, they do not provide evidence showing superior reproducibility experimentally in similar datasets.

<

p id=”p-12″>Several independent lines of evidence demonstrate molecular and metabolic oscillations in the absence of Bmal1 (7, 8) or an alternative clock knockout model [Cryptochrome‐deficient mice (9)], supporting our observations. To establish that transcriptome-level oscillations observed in Bmal1–/– cells are “noise” rather than true circadian rhythmicity, the authors could perform identical experiments within an environment that they think is completely devoid of external cues to verify that there are absolutely no 24-hour molecular oscillations in cultured Bmal1–/– cells under constant conditions.

Acknowledgments: Funding: A.B.R. acknowledges funding from the Perelman School of Medicine, University of Pennsylvania, and the Institute for Translational Medicine and Therapeutics (ITMAT), Perelman School of Medicine, University of Pennsylvania. Author contributions: Conceptualization and writing, A.B.R., S.R., and U.K.V.; all authors agreed on the interpretation of data and approved the final version of the manuscript; analysis and visualization, A.B.R., S.R., and U.K.V.; funding acquisition, A.B.R. Competing interests: The authors declare no competing interests. Data and materials availability: RNA-seq data are available on the Gene Expression Omnibus (accession number GSE111696 and GSE134333). Analysis scripts are available on GitHub (https://github.com/ReddysLab/Bmal1Paper).