The life of nearly all creatures on Earth follows the rhythm of day and night. For example, in fruit flies, darkness and light dictate when the insects feed, rest, move or mate. This is possible thanks to the circadian clock, an internal program which is synchronized with the environment to tell cells in the body when to perform certain roles.

In fruit flies, the structure that keeps the body clock ticking is formed of about 150 ‘circadian neurons’, which are divided into several subgroups. In these cells, a complex genetic programis at work, with networks of genes being ‘switched on’ in a cyclical way. To understand how this program works, scientists need to know which genes are turned on and when, as well as which proteins are created based on the information contained in these genes.

This can be difficult because one gene does not necessarily code for only one protein. Indeed, when a gene is turned on, it gets copied into a pre-messenger RNA (pre-mRNA), which is formed of several modules. The pre-mRNA can then go through a process called alternative splicing that shuffles or removes the different modules. This means that one gene can give rise to different pre-mRNA molecules that will each serve as a template to build a distinct protein. Until now, there have been few studies that examine the different types of pre-mRNAs found in circadian neurons, and how these change with the time of day.

Here, Wang, Abruzzi et al. extract three subgroups of circadian neurons, and one subgroup of non-circadian neurons, from the brain of fruit flies. The pre-mRNAs are isolated, and then a new computational method, called JUM, identifies, counts and categorizes the pre-mRNA molecules in the different groups of neurons.

This analysis reveals hundreds of previously unknown pre-mRNA molecules, many of which differed extensively between the types of brain cells. When comparing circadian and non-circadian neurons, Wang, Abruzzi et al. show that the circadian cells had more pre-mRNAs that code for proteins that help the cell communicate with other neurons. Finally, many genes in the circadian neurons use alternative splicing to turn on different types of pre-mRNA molecules at different times of the day in a cyclical way; this suggests that these pre-mRNAs might be participating in the genetic circadian program.

Many human disorders, such as certain forms of insomnia, emerge when the circadian clock is thrown off balance. The results reported by Wang, Abruzzi et al. show that alternative splicing may be an overlooked mechanism that shapes this complex program.

Organisms ranging from cyanobacteria to mammals contain circadian clocks that synchronize physiology and behavior to time-of-day environmental cues, such as light and temperature. The fruit fly, Drosophila, is no exception and exhibits robust circadian clocks that control a diverse range of behaviors such as activity, sleep, mating and feeding. The core molecular machinery of the circadian clock is a set of transcription factors that act together to drive cycling or time-of-day dependent transcription of at least 20% of genes in the genome (Abruzzi et al., 2011). The heterodimeric transcription factor CLOCK/CYCLE (CLK/CYC) binds to the promoters of many direct target genes and activates their transcription late in the day. Two of these direct target genes, Timeless (Tim) and Period (Per), encode transcriptional repressors that re-enter the nucleus in the late night and act as part of a negative feedback loop to turn off CLK/CYC transcription (Edery et al., 1994; Naidoo et al., 1999; Sehgal et al., 1994). In the early morning, the photoreceptor, Cryptochrome (CRY) and the kinase, Doubletime (Dbt) contribute to the degradation of Tim/Per (Kim et al., 2007; Kloss et al., 1998; Emery et al., 1998; Stanewsky et al., 1998; Ceriani et al., 1999; Hunter-Ensor et al., 1996; Syed et al., 2011). Without these repressors, CLK/CYC-mediated transcription begins again restarting the daily cycle of transcription.

In the Drosophila brain, this molecular clock resides in a neural network of ~150 neurons. These neurons are divided into seven subgroups including three groups of dorsal neurons (DNs; DN1, DN2, and DN3), three groups of lateral neurons (ventral and dorsal lateral neurons; large and small LNvs and LNds, respectively) and the lateral posterior neurons (LPNs) with discrete behavior functions; reviewed in (Peschel and Helfrich-Förster, 2011). These neurons make up a small portion of the Drosophila brain (~0.1%) and head (~0.05%). As a result, it is difficult to profile these specialized neurons as part of studies that focus on Drosophila brain and head tissues (Hughes et al., 2012; Rodriguez et al., 2013).

To learn more about the function of these neurons within the circadian neuronal circuit, three subgroups (DN1s, LNds, and LNvs) of these neurons as well as one non-circadian outgroup (Dopaminergic or tyrosine-hydroxylase (TH) expressing neurons) were labeled with GFP using neuron-specific drivers and manually isolated (Abruzzi et al., 2015; Abruzzi et al., 2017). Illumina RNA-seq datasets from these neurons were analyzed to examine the global mRNA landscape at the transcription level. This analysis revealed many cycling transcripts whose abundance changes with time of day that were not identified in previous studies of brain or head tissues (McDonald and Rosbash, 2001; Claridge-Chang et al., 2001; Wijnen et al., 2006; Rodriguez et al., 2013). In addition, many neuron-specific transcripts including novel circadian neuropeptides were revealed.

Besides transcriptional level gene regulation, post-transcriptional level regulation, such as alternative pre-mRNA splicing, is also known to be essential for the normal functioning of the nervous system. Alternative pre-mRNA splicing (AS) is a major gene regulatory mechanism that enables a single gene locus to produce populations of often functionally distinct mRNAs in a tissue- or cell-type-specific manner, greatly diversifying the metazoan transcriptome. The nervous system is well recognized to exhibit extraordinarily complex and diverse AS patterns in a variety of metazoan organisms (Li et al., 2007; Wang et al., 2008; Irimia et al., 2014; Li et al., 2015). Importantly, these diversified AS patterns have been shown to play important roles in modulating numerous neuronal states and activities (Vuong et al., 2016; Norris and Calarco, 2012). For example, many neuronal receptors such as the NMDA receptors undergo AS that results in functionally distinct proteins that are regulated in different ways to directly affect neuronal membrane depolarization and action-potential firing (Chandrasekar, 2013).

Recent findings reported several cases where AS is involved in the modulation of the circadian clock in plants, mice and Drosophila (Sanchez et al., 2010; Petrillo et al., 2011; McGlincy et al., 2012; Hughes et al., 2012; Preußner et al., 2014). However, a comprehensive understanding of AS profiles and their role in circadian clock regulation remains largely unexplored. As most studies profiling the Drosophila neuronal transcriptome originate from head samples, it has been difficult to detect cycling of alternatively spliced transcripts because heterogeneous tissues can often mask cell-type specific transcript cycling (Kula-Eversole et al., 2010; Nagoshi et al., 2010; Abruzzi et al., 2017). In addition, currently most AS analysis software for short-read Illumina RNA-seq data depend on pre-annotated libraries of known spliced transcripts, making it difficult to detect and quantitate novel neuronal AS events that prevail in neuronal tissues.

Here, we have comprehensively analyzed the global AS profiles in the transcriptomes of manually dissected DN1, LNd, and LNv circadian neurons, as well as the non-circadian, TH control neurons described above, using a new computational tool called JUM (the Junction Usage Model) (Wang and Rio, 2017a). JUM is able to identify, quantitate and categorize tissue-specific AS patterns from RNA-seq datasets without the need for or use of prior genome or transcriptome annotations. Using JUM, we identify hundreds of previously unannotated pre-mRNA isoforms. They differ extensively among the four neuronal subgroups and were totally missed when RNA-seq data from the heterogeneous fly head sample was analyzed using identical methods. Gene ontology analysis of differentially spliced mRNAs in circadian neurons versus non-circadian neurons revealed that they are enriched for transcripts encoding potassium channels, a hallmark of neuronal excitability. In addition, we discovered many new alternatively spliced variants and cycling AS patterns in transcripts encoding proteins of the molecular clock including the photoreceptor Cryptochrome and the kinase Shaggy. Importantly, we also identified a large set of transcripts that exhibit cycling patterns of alternative splicing throughout the day-night cycle, potentially contributing to the normal functioning of the circadian clock. This study highlights and reinforces the idea that post-transcriptional processes like AS impact cell type specification, diversity and function in the brain.

Tissues of the nervous and germline systems, such as brain, testes and ovaries, have more complex transcriptomes than other cell types due to extensive alternative pre-mRNA splicing or AS (Wang et al., 2008; Pan et al., 2008). The nervous system especially exhibits vast numbers of AS isoforms, many of which are novel and are only beginning to be comprehensively identified (Li et al., 2007; Irimia et al., 2014). This increase in transcript isoform complexity likely contributes to the specification and functional diversity of cell types within the nervous system.

Here, we have applied a novel computational algorithm called JUM to characterize the transcript isoform diversity generated by alternative splicing in three circadian neuronal subtypes (LNv, LNd and DN1), as well as a non-circadian dopaminergic neuron population (TH neurons) of the Drosophila central nervous system. JUM can comprehensively analyze, quantitate and compare tissue- or cell-type-specific AS patterns without requiring a priori annotations of known transcripts or transcriptomes (Wang and Rio, 2017a). Our analysis revealed a previously unappreciated diversity and complexity of alternatively spliced transcript isoform patterns in these four neuronal subtypes, suggesting that they contribute to neuronal identity, connectivity, activity and circadian functions. This is because many of these novel, previously undetected and unannotated isoforms were unique to a given neuronal population and occurred in transcripts from genes implicated in neuronal activity or circadian rhythms.

For example, the kinase Shaggy and the blue light photoreceptor cryptochrome play central roles in circadian clock regulation and have novel AS patterns in discrete subsets of the circadian neurons. In addition, nine different transcripts involved in potassium transport undergo differential AS in circadian neurons compared to non-circadian neurons. These transcripts encode six different potassium channels (Figure 5). Many of these genes have a complex organization known to encode populations of functionally distinct proteins isoforms, which change the activation kinetics as well as calcium sensitivity of the channels (Johnson et al., 2011). Neuronal firing is known to play a key role in the circadian circuit with recent studies illustrating that different subgroups of circadian neurons have characteristic time-of-day neuronal firing patterns (Flourakis et al., 2015; Liang et al., 2016; Guo et al., 2017). Although it is not yet fully understood which potassium channels play a critical role in each circadian neuron subgroup, several channel pre-mRNAs that undergo differential splicing in circadian neurons impact circadian behavior and sleep, such as slowpoke (Slo [Jaramillo et al., 2004]), Shaker (Sh; [Cirelli et al., 2005; Pimentel et al., 2016]) and Hyperkinetic (Hk; [Fogle et al., 2015]). It is therefore likely that AS adds diversity and distinct physiological properties to these protein isoforms, which then impacts neuron-specific firing patterns. From a more general perspective, AS augments transcriptional regulation in giving different circadian neurons individual identities and distinct functions.

Approximately 5% of the AS events identified in circadian neurons also undergo time-of-day dependent changes in alternative splicing (cycling splicing). It is important to note that all our experiments were conducted under 12 hr of light and 12 hr of dark conditions, making it impossible to distinguish between light and clock control. Nonetheless, these data indicate that splicing adds a dramatic layer of gene regulation to diurnal changes in gene expression. Moreover, many of the cycling AS transcripts show constant overall mRNA levels, which suggests the existence of neuron-specific splicing factors that are expressed or activated only at specific times of the day. Indeed, we have identified several candidate cycling neuron-enriched transcripts that encode RBPs that may help to drive cycling AS patterns.

A recent trend in biological research is to generate transcriptome profiles from single cells. For example, this strategy is part of the ‘human cell atlas’ project aimed at personalized genomic medicine or the ‘brain initiative’ project to generate profiles of all neurons in the mouse brain (Jorgenson et al., 2015; Regev et al., 2017). One recent study was able to obtain about 20M sequence reads per isolated human iPS cell but only managed to analyze splicing patterns for the most highly expressed genes (Song et al., 2017). Our study in contrast used ~100 isolated Drosophila neurons for each of the four neuron subtypes along with judicious use of both oligo-dT and random hexamer priming of the cDNA libraries. This strategy obtained about 10–30M sequence reads for each sample, including substantial information from the 5’ ends of transcripts, and JUM was able to detect and classify a large number of previously unannotated pre-mRNA isoforms. Many of them are missing from our fly head RNA-seq data assayed and analyzed in parallel, indicating that these new isoforms are cell-type specific. Not surprisingly, the novel isoforms from the three circadian neuron groups fall into many gene ontology (GO) categories associated with specific circadian clock activity and function.

Taken together, the work presented here indicates that the number of alternative splicing events that take place in neuronal tissues is grossly underestimated, even though publically-funded genome projects, such as the NIH modENCODE projects deeply sequenced transcriptomes from a variety of Drosophila tissues and developmental stages. This is despite the appreciation of how much AS occurs in the nervous system, for example recent comprehensive analysis of splicing patterns through deep sequencing of ~50 mouse and human tissues revealed about 2500 neuronally-regulated alternative splicing events (Irimia et al., 2014; Li et al., 2015). We therefore suggest that these events will need to be comprehensively evaluated by much deeper sequencing than is currently afforded by most contemporary single cell RNA-seq studies and by AS analysis software like JUM that is not constrained by a priori knowledge of known splicing events.

All data generated or analysed during this study are included in the manuscript and supporting files.

1. 1. Wang Q
   2. Taliaferro JM
   3. Klibaite U
   4. Hilgers V
   5. Shaevitz JW
   6. Rio DC

   (2016) The PSI-U1 snRNP interaction regulates male mating behavior in Drosophila

   Publicly available at the NCBI Gene Expression Omnibus (accession no: GSE79916).

   https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE79916
2. 1. Abruzzi K
   2. Chen X
   3. Nagoshi E
   4. Zadina A
   5. Rosbash M

   (2015) RNA-seq analysis of Drosophila clock and non-clock neurons reveals neuron-specific cycling and novel candidate neuropeptides

   Publicly available at the NCBI Gene Expression Omnibus (accession no: GSE77451).

   https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE77451

1. 1. Abruzzi K
   2. Chen X
   3. Nagoshi E
   4. Zadina A
   5. Rosbash M

   (2015) RNA-seq profiling of small numbers of Drosophila neurons

   Methods in Enzymology 551:369–386.

   https://doi.org/10.1016/bs.mie.2014.10.025

   • PubMed
   • Google Scholar
2. 1. Abruzzi KC
   2. Rodriguez J
   3. Menet JS
   4. Desrochers J
   5. Zadina A
   6. Luo W
   7. Tkachev S
   8. Rosbash M

   (2011) Drosophila CLOCK target gene characterization: implications for circadian tissue-specific gene expression

   Genes & Development 25:2374–2386.

   https://doi.org/10.1101/gad.178079.111

   • PubMed
   • Google Scholar
3. 1. Abruzzi KC
   2. Zadina A
   3. Luo W
   4. Wiyanto E
   5. Rahman R
   6. Guo F
   7. Shafer O
   8. Rosbash M

   (2017) RNA-seq analysis of Drosophila clock and non-clock neurons reveals neuron-specific cycling and novel candidate neuropeptides

   PLoS Genetics 13:e1006613.

   https://doi.org/10.1371/journal.pgen.1006613

   • PubMed
   • Google Scholar
4. 1. Ceriani MF
   2. Darlington TK
   3. Staknis D
   4. Más P
   5. Petti AA
   6. Weitz CJ
   7. Kay SA

   (1999) Light-dependent sequestration of timeless by cryptochrome

   Science 285:553–556.

   https://doi.org/10.1126/science.285.5427.553

   • PubMed
   • Google Scholar
5. 1. Ceriani MF
   2. Hogenesch JB
   3. Yanovsky M
   4. Panda S
   5. Straume M
   6. Kay SA

   (2002) Genome-wide expression analysis in Drosophila reveals genes controlling circadian behavior

   The Journal of Neuroscience 22:9305–9319.

   https://doi.org/10.1523/JNEUROSCI.22-21-09305.2002

   • PubMed
   • Google Scholar
6. 1. Chandrasekar R

   (2013) Alcohol and NMDA receptor: current research and future direction

   Frontiers in Molecular Neuroscience 6:14.

   https://doi.org/10.3389/fnmol.2013.00014

   • PubMed
   • Google Scholar
7. 1. Cirelli C
   2. Bushey D
   3. Hill S
   4. Huber R
   5. Kreber R
   6. Ganetzky B
   7. Tononi G

   (2005) Reduced sleep in Drosophila Shaker mutants

   Nature 434:1087–1092.

   https://doi.org/10.1038/nature03486

   • PubMed
   • Google Scholar
8. 1. Claridge-Chang A
   2. Wijnen H
   3. Naef F
   4. Boothroyd C
   5. Rajewsky N
   6. Young MW

   (2001) Circadian regulation of gene expression systems in the Drosophila head

   Neuron 32:657–671.

   https://doi.org/10.1016/S0896-6273(01)00515-3

   • PubMed
   • Google Scholar
9. 1. Derr A
   2. Yang C
   3. Zilionis R
   4. Sergushichev A
   5. Blodgett DM
   6. Redick S
   7. Bortell R
   8. Luban J
   9. Harlan DM
   10. Kadener S
   11. Greiner DL
   12. Klein A
   13. Artyomov MN
   14. Garber M

   (2016) End sequence analysis toolkit (ESAT) expands the extractable information from single-cell RNA-seq data

   Genome Research 26:1397–1410.

   https://doi.org/10.1101/gr.207902.116

   • PubMed
   • Google Scholar
10. 1. Dobin A
    2. Davis CA
    3. Schlesinger F
    4. Drenkow J
    5. Zaleski C
    6. Jha S
    7. Batut P
    8. Chaisson M
    9. Gingeras TR

    (2013) STAR: ultrafast universal RNA-seq aligner

    Bioinformatics 29:15–21.

    https://doi.org/10.1093/bioinformatics/bts635

    • PubMed
    • Google Scholar
11. 1. Edery I
    2. Zwiebel LJ
    3. Dembinska ME
    4. Rosbash M

    (1994) Temporal phosphorylation of the Drosophila period protein

    PNAS 91:2260–2264.

    https://doi.org/10.1073/pnas.91.6.2260

    • PubMed
    • Google Scholar
12. 1. Emery P
    2. So WV
    3. Kaneko M
    4. Hall JC
    5. Rosbash M

    (1998) CRY, a Drosophila clock and light-regulated cryptochrome, is a major contributor to circadian rhythm resetting and photosensitivity

    Cell 95:669–679.

    https://doi.org/10.1016/S0092-8674(00)81637-2

    • PubMed
    • Google Scholar
13. 1. Fernández MP
    2. Berni J
    3. Ceriani MF

    (2008) Circadian remodeling of neuronal circuits involved in rhythmic behavior

    PLoS Biology 6:e69.

    https://doi.org/10.1371/journal.pbio.0060069

    • PubMed
    • Google Scholar
14. 1. Flourakis M
    2. Kula-Eversole E
    3. Hutchison AL
    4. Han TH
    5. Aranda K
    6. Moose DL
    7. White KP
    8. Dinner AR
    9. Lear BC
    10. Ren D
    11. Diekman CO
    12. Raman IM
    13. Allada R

    (2015) A conserved bicycle model for circadian clock control of membrane excitability

    Cell 162:836–848.

    https://doi.org/10.1016/j.cell.2015.07.036

    • PubMed
    • Google Scholar
15. 1. Fogle KJ
    2. Baik LS
    3. Houl JH
    4. Tran TT
    5. Roberts L
    6. Dahm NA
    7. Cao Y
    8. Zhou M
    9. Holmes TC

    (2015) CRYPTOCHROME-mediated phototransduction by modulation of the potassium ion channel β-subunit redox sensor

    PNAS 112:2245–2250.

    https://doi.org/10.1073/pnas.1416586112

    • PubMed
    • Google Scholar
16. 1. Gorostiza EA
    2. Depetris-Chauvin A
    3. Frenkel L
    4. Pírez N
    5. Ceriani MF

    (2014) Circadian pacemaker neurons change synaptic contacts across the day

    Current Biology 24:2161–2167.

    https://doi.org/10.1016/j.cub.2014.07.063

    • PubMed
    • Google Scholar
17. 1. Guo F
    2. Chen X
    3. Rosbash M

    (2017) Temporal calcium profiling of specific circadian neurons in freely moving flies

    PNAS 114:E8780–E8787.

    https://doi.org/10.1073/pnas.1706608114

    • PubMed
    • Google Scholar
18. 1. Hansen HT
    2. Rasmussen SH
    3. Adolph SK
    4. Plass M
    5. Krogh A
    6. Sanford J
    7. Nielsen FC
    8. Christiansen J

    (2015) Drosophila imp iCLIP identifies an RNA assemblage coordinating F-actin formation

    Genome Biology 16:123.

    https://doi.org/10.1186/s13059-015-0687-0

    • PubMed
    • Google Scholar
19. 1. Huang daW
    2. Sherman BT
    3. Lempicki RA

    (2009a) Systematic and integrative analysis of large gene lists using DAVID bioinformatics resources

    Nature Protocols 4:44–57.

    https://doi.org/10.1038/nprot.2008.211

    • PubMed
    • Google Scholar
20. 1. Huang daW
    2. Sherman BT
    3. Lempicki RA

    (2009b) Bioinformatics enrichment tools: paths toward the comprehensive functional analysis of large gene lists

    Nucleic Acids Research 37:1–13.

    https://doi.org/10.1093/nar/gkn923

    • PubMed
    • Google Scholar
21. 1. Hughes ME
    2. Grant GR
    3. Paquin C
    4. Qian J
    5. Nitabach MN

    (2012) Deep sequencing the circadian and diurnal transcriptome of Drosophila brain

    Genome Research 22:1266–1281.

    https://doi.org/10.1101/gr.128876.111

    • PubMed
    • Google Scholar
22. 1. Hughes ME
    2. Hogenesch JB
    3. Kornacker K

    (2010) JTK_CYCLE: an efficient nonparametric algorithm for detecting rhythmic components in genome-scale data sets

    Journal of Biological Rhythms 25:372–380.

    https://doi.org/10.1177/0748730410379711

    • PubMed
    • Google Scholar
23. 1. Hunter-Ensor M
    2. Ousley A
    3. Sehgal A

    (1996) Regulation of the Drosophila protein timeless suggests a mechanism for resetting the circadian clock by light

    Cell 84:677–685.

    https://doi.org/10.1016/S0092-8674(00)81046-6

    • PubMed
    • Google Scholar
24. 1. Irimia M
    2. Weatheritt RJ
    3. Ellis JD
    4. Parikshak NN
    5. Gonatopoulos-Pournatzis T
    6. Babor M
    7. Quesnel-Vallières M
    8. Tapial J
    9. Raj B
    10. O'Hanlon D
    11. Barrios-Rodiles M
    12. Sternberg MJ
    13. Cordes SP
    14. Roth FP
    15. Wrana JL
    16. Geschwind DH
    17. Blencowe BJ

    (2014) A highly conserved program of neuronal microexons is misregulated in autistic brains

    Cell 159:1511–1523.

    https://doi.org/10.1016/j.cell.2014.11.035

    • PubMed
    • Google Scholar
25. 1. Jaramillo AM
    2. Zheng X
    3. Zhou Y
    4. Amado DA
    5. Sheldon A
    6. Sehgal A
    7. Levitan IB

    (2004) Pattern of distribution and cycling of SLOB, slowpoke channel binding protein, in Drosophila

    BMC Neuroscience 5:3.

    https://doi.org/10.1186/1471-2202-5-3

    • PubMed
    • Google Scholar
26. 1. Johnson BE
    2. Glauser DA
    3. Dan-Glauser ES
    4. Halling DB
    5. Aldrich RW
    6. Goodman MB

    (2011) Alternatively spliced domains interact to regulate BK potassium channel gating

    PNAS 108:20784–20789.

    https://doi.org/10.1073/pnas.1116795108

    • PubMed
    • Google Scholar
27. 1. Jorgenson LA
    2. Newsome WT
    3. Anderson DJ
    4. Bargmann CI
    5. Brown EN
    6. Deisseroth K
    7. Donoghue JP
    8. Hudson KL
    9. Ling GSF
    10. MacLeish PR
    11. Marder E
    12. Normann RA
    13. Sanes JR
    14. Schnitzer MJ
    15. Sejnowski TJ
    16. Tank DW
    17. Tsien RY
    18. Ugurbil K
    19. Wingfield JC

    (2015) The BRAIN initiative: developing technology to catalyse neuroscience discovery

    Philosophical Transactions of the Royal Society B: Biological Sciences 370:20140164.

    https://doi.org/10.1098/rstb.2014.0164

    • Google Scholar
28. 1. Katz Y
    2. Wang ET
    3. Silterra J
    4. Schwartz S
    5. Wong B
    6. Thorvaldsdóttir H
    7. Robinson JT
    8. Mesirov JP
    9. Airoldi EM
    10. Burge CB

    (2015) Quantitative visualization of alternative exon expression from RNA-seq data

    Bioinformatics 31:2400–2402.

    https://doi.org/10.1093/bioinformatics/btv034

    • PubMed
    • Google Scholar
29. 1. Kim EY
    2. Ko HW
    3. Yu W
    4. Hardin PE
    5. Edery I

    (2007) A Doubletime kinase binding domain on the Drosophila PERIOD protein is essential for its hyperphosphorylation, transcriptional repression, and circadian clock function

    Molecular and Cellular Biology 27:5014–5028.

    https://doi.org/10.1128/MCB.02339-06

    • PubMed
    • Google Scholar
30. 1. Kloss B
    2. Price JL
    3. Saez L
    4. Blau J
    5. Rothenfluh A
    6. Wesley CS
    7. Young MW

    (1998) The Drosophila clock gene double-time encodes a protein closely related to human casein kinase Iepsilon

    Cell 94:97–107.

    https://doi.org/10.1016/S0092-8674(00)81225-8

    • PubMed
    • Google Scholar
31. 1. Kula-Eversole E
    2. Nagoshi E
    3. Shang Y
    4. Rodriguez J
    5. Allada R
    6. Rosbash M

    (2010) Surprising gene expression patterns within and between PDF-containing circadian neurons in Drosophila

    PNAS 107:13497–13502.

    https://doi.org/10.1073/pnas.1002081107

    • PubMed
    • Google Scholar
32. 1. Li H

    (2011) A statistical framework for SNP calling, mutation discovery, association mapping and population genetical parameter estimation from sequencing data

    Bioinformatics 27:2987–2993.

    https://doi.org/10.1093/bioinformatics/btr509

    • PubMed
    • Google Scholar
33. 1. Li Q
    2. Lee JA
    3. Black DL

    (2007) Neuronal regulation of alternative pre-mRNA splicing

    Nature Reviews Neuroscience 8:819–831.

    https://doi.org/10.1038/nrn2237

    • PubMed
    • Google Scholar
34. 1. Li YI
    2. Sanchez-Pulido L
    3. Haerty W
    4. Ponting CP

    (2015) RBFOX and PTBP1 proteins regulate the alternative splicing of micro-exons in human brain transcripts

    Genome Research 25:1–13.

    https://doi.org/10.1101/gr.181990.114

    • PubMed
    • Google Scholar
35. 1. Liang X
    2. Holy TE
    3. Taghert PH

    (2016) Synchronous Drosophila circadian pacemakers display nonsynchronous Ca²⁺ rhythms in vivo

    Science 351:976–981.

    https://doi.org/10.1126/science.aad3997

    • PubMed
    • Google Scholar
36. 1. Love MI
    2. Huber W
    3. Anders S

    (2014) Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2

    Genome Biology 15:550.

    https://doi.org/10.1186/s13059-014-0550-8

    • PubMed
    • Google Scholar
37. 1. McDonald MJ
    2. Rosbash M

    (2001) Microarray analysis and organization of circadian gene expression in Drosophila

    Cell 107:567–578.

    https://doi.org/10.1016/S0092-8674(01)00545-1

    • PubMed
    • Google Scholar
38. 1. McGlincy NJ
    2. Valomon A
    3. Chesham JE
    4. Maywood ES
    5. Hastings MH
    6. Ule J

    (2012) Regulation of alternative splicing by the circadian clock and food related cues

    Genome Biology 13:R54.

    https://doi.org/10.1186/gb-2012-13-6-r54

    • PubMed
    • Google Scholar
39. 1. Meredith AL
    2. Wiler SW
    3. Miller BH
    4. Takahashi JS
    5. Fodor AA
    6. Ruby NF
    7. Aldrich RW

    (2006) BK calcium-activated potassium channels regulate circadian behavioral rhythms and pacemaker output

    Nature Neuroscience 9:1041–1049.

    https://doi.org/10.1038/nn1740

    • PubMed
    • Google Scholar
40. 1. Nagoshi E
    2. Sugino K
    3. Kula E
    4. Okazaki E
    5. Tachibana T
    6. Nelson S
    7. Rosbash M

    (2010) Dissecting differential gene expression within the circadian neuronal circuit of Drosophila

    Nature Neuroscience 13:60–68.

    https://doi.org/10.1038/nn.2451

    • PubMed
    • Google Scholar
41. 1. Naidoo N
    2. Song W
    3. Hunter-Ensor M
    4. Sehgal A

    (1999) A role for the proteasome in the light response of the timeless clock protein

    Science 285:1737–1741.

    https://doi.org/10.1126/science.285.5434.1737

    • PubMed
    • Google Scholar
42. 1. Norris AD
    2. Calarco JA

    (2012) Emerging roles of alternative pre-mRNA splicing regulation in neuronal development and function

    Frontiers in Neuroscience 6:122.

    https://doi.org/10.3389/fnins.2012.00122

    • PubMed
    • Google Scholar
43. 1. Pan Q
    2. Shai O
    3. Lee LJ
    4. Frey BJ
    5. Blencowe BJ

    (2008) Deep surveying of alternative splicing complexity in the human transcriptome by high-throughput sequencing

    Nature Genetics 40:1413–1415.

    https://doi.org/10.1038/ng.259

    • PubMed
    • Google Scholar
44. 1. Peschel N
    2. Helfrich-Förster C

    (2011) Setting the clock--by nature: circadian rhythm in the fruitfly Drosophila melanogaster

    FEBS Letters 585:1435–1442.

    https://doi.org/10.1016/j.febslet.2011.02.028

    • PubMed
    • Google Scholar
45. 1. Petrillo E
    2. Sanchez SE
    3. Kornblihtt AR
    4. Yanovsky MJ

    (2011) Alternative splicing adds a new loop to the circadian clock

    Communicative & Integrative Biology 4:284–286.

    https://doi.org/10.4161/cib.4.3.14777

    • PubMed
    • Google Scholar
46. 1. Petsakou A
    2. Sapsis TP
    3. Blau J

    (2015) Circadian rhythms in Rho1 activity regulate neuronal plasticity and network hierarchy

    Cell 162:823–835.

    https://doi.org/10.1016/j.cell.2015.07.010

    • PubMed
    • Google Scholar
47. 1. Pimentel D
    2. Donlea JM
    3. Talbot CB
    4. Song SM
    5. Thurston AJF
    6. Miesenböck G

    (2016) Operation of a homeostatic sleep switch

    Nature 536:333–337.

    https://doi.org/10.1038/nature19055

    • PubMed
    • Google Scholar
48. 1. Preußner M
    2. Wilhelmi I
    3. Schultz AS
    4. Finkernagel F
    5. Michel M
    6. Möröy T
    7. Heyd F

    (2014) Rhythmic U2af26 alternative splicing controls PERIOD1 stability and the circadian clock in mice

    Molecular Cell 54:651–662.

    https://doi.org/10.1016/j.molcel.2014.04.015

    • PubMed
    • Google Scholar
49. 1. Regev A
    2. Teichmann SA
    3. Lander ES
    4. Amit I
    5. Benoist C
    6. Birney E
    7. Bodenmiller B
    8. Campbell P
    9. Carninci P
    10. Clatworthy M
    11. Clevers H
    12. Deplancke B
    13. Dunham I
    14. Eberwine J
    15. Eils R
    16. Enard W
    17. Farmer A
    18. Fugger L
    19. Göttgens B
    20. Hacohen N
    21. Haniffa M
    22. Hemberg M
    23. Kim S
    24. Klenerman P
    25. Kriegstein A
    26. Lein E
    27. Linnarsson S
    28. Lundberg E
    29. Lundeberg J
    30. Majumder P
    31. Marioni JC
    32. Merad M
    33. Mhlanga M
    34. Nawijn M
    35. Netea M
    36. Nolan G
    37. Pe'er D
    38. Phillipakis A
    39. Ponting CP
    40. Quake S
    41. Reik W
    42. Rozenblatt-Rosen O
    43. Sanes J
    44. Satija R
    45. Schumacher TN
    46. Shalek A
    47. Shapiro E
    48. Sharma P
    49. Shin JW
    50. Stegle O
    51. Stratton M
    52. Stubbington MJT
    53. Theis FJ
    54. Uhlen M
    55. van Oudenaarden A
    56. Wagner A
    57. Watt F
    58. Weissman J
    59. Wold B
    60. Xavier R
    61. Yosef N
    62. Human Cell Atlas Meeting Participants

    (2017) The human cell atlas

    eLife 6:e27041.

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

    • PubMed
    • Google Scholar
50. 1. Robinson JT
    2. Thorvaldsdóttir H
    3. Winckler W
    4. Guttman M
    5. Lander ES
    6. Getz G
    7. Mesirov JP

    (2011) Integrative genomics viewer

    Nature Biotechnology 29:24–26.

    https://doi.org/10.1038/nbt.1754

    • PubMed
    • Google Scholar
51. 1. Rodriguez J
    2. Tang CH
    3. Khodor YL
    4. Vodala S
    5. Menet JS
    6. Rosbash M

    (2013) Nascent-Seq analysis of Drosophila cycling gene expression

    PNAS 110:E275–E284.

    https://doi.org/10.1073/pnas.1219969110

    • PubMed
    • Google Scholar
52. 1. Sanchez SE
    2. Petrillo E
    3. Beckwith EJ
    4. Zhang X
    5. Rugnone ML
    6. Hernando CE
    7. Cuevas JC
    8. Godoy Herz MA
    9. Depetris-Chauvin A
    10. Simpson CG
    11. Brown JW
    12. Cerdán PD
    13. Borevitz JO
    14. Mas P
    15. Ceriani MF
    16. Kornblihtt AR
    17. Yanovsky MJ

    (2010) A methyl transferase links the circadian clock to the regulation of alternative splicing

    Nature 468:112–116.

    https://doi.org/10.1038/nature09470

    • PubMed
    • Google Scholar
53. 1. Schneider CA
    2. Rasband WS
    3. Eliceiri KW

    (2012) NIH image to ImageJ: 25 years of image analysis

    Nature Methods 9:671–675.

    https://doi.org/10.1038/nmeth.2089

    • PubMed
    • Google Scholar
54. 1. Sehgal A
    2. Price JL
    3. Man B
    4. Young MW

    (1994) Loss of circadian behavioral rhythms and per RNA oscillations in the Drosophila mutant timeless

    Science 263:1603–1606.

    https://doi.org/10.1126/science.8128246

    • PubMed
    • Google Scholar
55. 1. Siepel A
    2. Bejerano G
    3. Pedersen JS
    4. Hinrichs AS
    5. Hou M
    6. Rosenbloom K
    7. Clawson H
    8. Spieth J
    9. Hillier LW
    10. Richards S
    11. Weinstock GM
    12. Wilson RK
    13. Gibbs RA
    14. Kent WJ
    15. Miller W
    16. Haussler D

    (2005) Evolutionarily conserved elements in vertebrate, insect, worm, and yeast genomes

    Genome Research 15:1034–1050.

    https://doi.org/10.1101/gr.3715005

    • PubMed
    • Google Scholar
56. 1. Song Y
    2. Botvinnik OB
    3. Lovci MT
    4. Kakaradov B
    5. Liu P
    6. Xu JL
    7. Yeo GW

    (2017) Single-Cell alternative splicing analysis with expedition reveals splicing dynamics during neuron differentiation

    Molecular Cell 67:148–161.

    https://doi.org/10.1016/j.molcel.2017.06.003

    • PubMed
    • Google Scholar
57. 1. Spitzer M
    2. Wildenhain J
    3. Rappsilber J
    4. Tyers M

    (2014) BoxPlotR: a web tool for generation of box plots

    Nature Methods 11:121–122.

    https://doi.org/10.1038/nmeth.2811

    • PubMed
    • Google Scholar
58. 1. Stanewsky R
    2. Kaneko M
    3. Emery P
    4. Beretta B
    5. Wager-Smith K
    6. Kay SA
    7. Rosbash M
    8. Hall JC

    (1998) The cryb mutation identifies cryptochrome as a circadian photoreceptor in Drosophila

    Cell 95:681–692.

    https://doi.org/10.1016/S0092-8674(00)81638-4

    • PubMed
    • Google Scholar
59. 1. Stoiber MH
    2. Olson S
    3. May GE
    4. Duff MO
    5. Manent J
    6. Obar R
    7. Guruharsha KG
    8. Bickel PJ
    9. Artavanis-Tsakonas S
    10. Brown JB
    11. Graveley BR
    12. Celniker SE

    (2015) Extensive cross-regulation of post-transcriptional regulatory networks in Drosophila

    Genome Research 25:1692–1702.

    https://doi.org/10.1101/gr.182675.114

    • PubMed
    • Google Scholar
60. 1. Syed S
    2. Saez L
    3. Young MW

    (2011) Kinetics of doubletime kinase-dependent degradation of the Drosophila period protein

    Journal of Biological Chemistry 286:27654–27662.

    https://doi.org/10.1074/jbc.M111.243618

    • PubMed
    • Google Scholar
61. 1. Thorvaldsdóttir H
    2. Robinson JT
    3. Mesirov JP

    (2013) Integrative genomics viewer (IGV): high-performance genomics data visualization and exploration

    Briefings in Bioinformatics 14:178–192.

    https://doi.org/10.1093/bib/bbs017

    • PubMed
    • Google Scholar
62. 1. Vuong CK
    2. Black DL
    3. Zheng S

    (2016) The neurogenetics of alternative splicing

    Nature Reviews Neuroscience 17:265–281.

    https://doi.org/10.1038/nrn.2016.27

    • PubMed
    • Google Scholar
63. 1. Wang ET
    2. Sandberg R
    3. Luo S
    4. Khrebtukova I
    5. Zhang L
    6. Mayr C
    7. Kingsmore SF
    8. Schroth GP
    9. Burge CB

    (2008) Alternative isoform regulation in human tissue transcriptomes

    Nature 456:470–476.

    https://doi.org/10.1038/nature07509

    • PubMed
    • Google Scholar
64. Preprint

    1. Wang Q
    2. Rio D

    (2017a) The junction usage model (JUM): a method for comprehensive annotation-free differential analysis of tissue-specific global alternative pre-mRNA splicing patterns

    bioRxiv.

    https://doi.org/10.1101/116863

    • Google Scholar
65. Software

    1. Wang Q
    2. Rio D

    (2017b) JUM_Pipeline, version v1.3.12

    Github.

    https://github.com/qqwang-berkeley/JUM
66. 1. Wang Q
    2. Taliaferro JM
    3. Klibaite U
    4. Hilgers V
    5. Shaevitz JW
    6. Rio DC

    (2016) The PSI-U1 snRNP interaction regulates male mating behavior in Drosophila

    PNAS 113:5269–5274.

    https://doi.org/10.1073/pnas.1600936113

    • PubMed
    • Google Scholar
67. 1. Wijnen H
    2. Naef F
    3. Boothroyd C
    4. Claridge-Chang A
    5. Young MW

    (2006) Control of daily transcript oscillations in Drosophila by light and the circadian clock

    PLoS Genetics 2:e39.

    https://doi.org/10.1371/journal.pgen.0020039

    • PubMed
    • Google Scholar
68. 1. Wijnen H
    2. Naef F
    3. Young MW

    (2005) Molecular and statistical tools for circadian transcript profiling

    Methods in Enzymology 393:341–365.

    https://doi.org/10.1016/S0076-6879(05)93015-2

    • PubMed
    • Google Scholar

Thank you for submitting your article "Diversity and circadian cycling of alternative splicing in the transcriptomes of isolated Drosophila neuron populations" for consideration by eLife. Your article has been reviewed by three peer reviewers, and the evaluation has been overseen by a Reviewing Editor and a Senior Editor. The reviewers have opted to remain anonymous.

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 from Wang et al., the authors re-analyzed previously published RNA-seq datasets for alternative splicing (AS) patterns comparing three neuronal populations in Drosophila melanogaster that are critical for circadian rhythms (DN1, LNd, and LNv) with a control group of TH neurons. They have also analyzed these same populations over the course of two days. Overall, this is a coherent paper, applying sophisticated analytical method to high throughput RNA-seq data to uncover the novel findings.

Using the comprehensive computational method JUM (Junction Usage model) to analyze splicing patterns revealed a striking number of AS junctions in the aforementioned neurons absent from the more standard Drosophila head RNA preparation. A significant portion (~15%) of these were exclusive to the individual neuron types and almost entirely were previously unannotated (~95%). Each neuron type was also shown to contain a broad range of unique AS junctions when compared to the TH group, and that these genes could impact neuronal function (ex. N-syb, Shab). GO analysis on the AS genes from the circadian neurons showed globally that these patterns are more prominent in potassium channels, ion exchange pumps (sodium-potassium ATPases), and, more importantly, functions such as vesicle exocytosis that are critical for circadian function. Finally, the authors conducted a 2-day time course dissection of these same neurons every 4 hours and analyzed AS splicing again using JUM. This experiment found unique AS junctions in all circadian neurons that cycled through the day/night; however, the overall level of mRNA from most of these genes did not change, demonstrating that the control of gene function is most likely due to AS rather than transcriptional control.

There are interesting insights from these experiments into the much greater role of alternative splicing in circadian rhythms, however, there is too much emphasis on JUM, which has already been described elsewhere, in Wang and Rio, 2017. Considering this, the extensive description of the tool, nomenclature, and usage is not necessary and detracts from the biological aspect of the manuscript. That is not to say that JUM is not worth a good amount of emphasis – it is impressive, but a revision is needed to achieve a better balance. Substantial portions of the Results section and relevant Figures (Figure 1A in particular) could be removed.

Essential revisions:

1) The authors indicate that they have generated RNA-seq data from individual neuronal populations as described in Abruzzi et al., 2017 but it is not clear if the flies were maintained in 12:12 light:dark or constant darkness during the experiment. This distinction could impact the interpretation of circadian analyses since under 12:12 light:dark many gene expression changes and AS events may be driven by light sensing as opposed to an intrinsic clock.

2) The authors consistently refer only to the number of AS junction sites, without providing the number of affected genes. While the number of events is of interest because it gives a sense of how pervasive alternative splicing is, it how many genes are affected is necessary to get a sense of the scope and selectivity of this phenomenon. For example, in the Results section the authors explain the AS junctions in two genes, Cry and CG10483, but do not provide a perspective – how many of the genes with these AS structures have functions relevant to the different neurons in this study and of neuronal populations over other cell types?

3) There were only 10 AS events found in all three groups of neurons (Figure 3B). Was any further analysis done on these? These might be genes that are particularly important as well.

4) The authors could point out the use of other algorithms to identify AS changes – why rely only on JUM?

5) A critical point of this study is the identification of "novel" AS junctions. Therefore, it will be important to have some quality control metrics to show the identified junctions are real rather than technical artifacts or alignment errors. It will be helpful to clearly describe parameters related to accuracy of exon junction identification (e.g., size of overlap on each side of the junction), and if possible provide an estimate of accuracy (e.g. using simulated synthetic exon junction reads by swapping 5' and 3' parts).

6) Related to the point above, it is also important to clearly state the criteria used to define a "novel" AS junction in the main text. As of now, the authors stated in the method section that "The novelty of each AS junction was compared against the library of annotated junctions in the UCSC genome browser transcriptome annotation (genome 422 version: FB2017_05)." It is unclear what type of annotations was used here (e.g., which gene models, does that include EST sequences?)

7) In subsection “Circadian neurons present specific alternative pre-mRNA splicing patterns”, the authors identified AS junctions differentially spliced in each type of circadian neurons compared to TH neurons. They then examined the overlap of the lists to claim most of these junctions are specifically spliced in one subtype. To be more rigorous, one should perform differential splicing analysis between different types of circadian neurons, as the limited overlap could also be due to limited statistical power.

Similar issue for comparison of cycling exons in Figure 6B.

8) While cell type- and time of day-dependent AS patterns are artfully presented here, a major weakness in this manuscript is that the underlying mechanisms driving AS or their consequences on protein expression are not studied. Regarding upstream mechanism, it would help to know if expression of any RNA binding proteins correlate with the observed AS events across cell types or times of day. Regarding the relative impact of AS on neuronal identity and function, a first pass analysis could be to compare the relative numbers and enriched ontologies between mRNA expression and AS differences between cell types and in the circadian data. A more compelling experiment would be to specifically inhibit an AS event within one of the neuronal populations to test its effect on cell function.

9) One potential concern is the technical/biological variation between replicates, which is expected from the small number of cells used to generate RNA-seq libraries and visible from the examples shown in the figures. It will be helpful to have some data to show the reliability/reproducibility of observed differential splicing (e.g., correlation of splicing pattern of replicates, or even independent RT-PCR validation if possible).

10) The observation of alternative exons that are differentially spliced in specific types of neurons or time-of-the-day during the circadian clock is very interesting. Is there any indication of corresponding changes in splicing factors that potentially regulate these AS events?

11) It will be helpful to provide additional functional analysis of identified AS events with differential splicing, in addition to GO. This reviewer understands that detailed functional dissection of specific events might be beyond the scope of this manuscript, but it should be feasible to perform statistical analysis of cross-species conservation, preservation of reading frame, exon size, etc., which is frequently distinct for exons with cell type specific splicing.

12) In Figure 4B, although the authors claim, "in TH and LNd neurons the Shab transcripts preferentially utilize an alternative terminal exon, […] while in LNv and DN1 neurons, the long Shab isoform that encodes all transmembrane regions is dominant", it is not intuitively showing that by looking at the Sashimi plots. The authors should order the samples similarly to Figure 4A, and also may better prove their point by grouping the LNv and DN1 samples or highlighting the relevant differences in the figure read numbers.

Summary:

[…] There are interesting insights from these experiments into the much greater role of alternative splicing in circadian rhythms, however, there is too much emphasis on JUM, which has already been described elsewhere, in Wang and Rio, 2017. Considering this, the extensive description of the tool, nomenclature, and usage is not necessary and detracts from the biological aspect of the manuscript. That is not to say that JUM is not worth a good amount of emphasis – it is impressive, but a revision is needed to achieve a better balance. Substantial portions of the Results section and relevant Figures (Figure 1A in particular) could be removed.

The reviewers suggested that there is too much emphasis on the computational method (JUM) used to analyze the RNA-seq datasets in this work, since it has already been described in detail in our JUM method paper deposited in bioRxiv. We have modified the manuscript accordingly and only kept a concise description of the nomenclature of the method for easy reference of the readers.

Essential revisions:

1) The authors indicate that they have generated RNA-seq data from individual neuronal populations as described in Abruzzi et al., 2017 but it is not clear if the flies were maintained in 12:12 light:dark or constant darkness during the experiment. This distinction could impact the interpretation of circadian analyses since under 12:12 light:dark many gene expression changes and AS events may be driven by light sensing as opposed to an intrinsic clock.

We (Rosbash lab) have been focusing on LD (light-dark) behavior since beginning these experiments with microarrays a decade ago with the first publications in 2010. More relevant of course are the LD data (Abruzzi et al., 2017) on which this analysis is based. So, we are of course constrained without generating this massive amount of data de novo in constant darkness. We took this strategy at the outset, which we still endorse, because fly behavior is richer in LD than in DD, which allows for a more interesting comparison between transcripts, their presence in specific neurons as well as their cycling, and behavior. Because the reviewer is certainly correct, and light is probably germane (light alone or an interaction between light and the clock), we have also avoided the term "circadian" for these transcripts but refer to them as "cycling" or "oscillating" given their once/day or diurnal pattern. Nonetheless, we now clarify this point explicitly in the Discussion section.

2) The authors consistently refer only to the number of AS junction sites, without providing the number of affected genes. While the number of events is of interest because it gives a sense of how pervasive alternative splicing is, it how many genes are affected is necessary to get a sense of the scope and selectivity of this phenomenon. For example, in the Results section the authors explain the AS junctions in two genes, Cry and CG10483, but do not provide a perspective – how many of the genes with these AS structures have functions relevant to the different neurons in this study and of neuronal populations over other cell types?

We have included additional data throughout the paper that address the number and identity of the transcripts undergoing alternative splicing. The number of genes affected has been added to Figure 1B,D and Figure 3A,B. Gene ontology analysis supporting Figure 1 and Figure 6 has been added as Figure 1—figure supplement 2 and Figure 6—figure supplement 3, as well as Supplementary file 2 and Supplementary file 3. Descriptions of these new data have been added to subsection “There are many novel AS patterns in neuron subpopulations”; subsection “Identification of cycling alternatively spliced isoforms in circadian neurons”.

3) There were only 10 AS events found in all three groups of neurons (Figure 3B). Was any further analysis done on these? These might be genes that are particularly important as well.

We did a detailed analysis of the significantly changed AS events comparing TH neurons to all three groups of circadian neurons. Interestingly, the specific AS changes occur in genes that are functionally enriched in regulating neuronal plasticity, remodeling and synaptic transmission (Pten, PRL-1, Pdp-1, Sap47, Rim, mmd, CG9945 and RpS3A). It is likely that AS serves as an independent layer of gene regulation to affect the diversity of synapses and synaptic plasticity in the circadian neurons. We added the results from this analysis per the reviewer’s request in subsection “Circadian neurons present specific alternative pre-mRNA splicing patterns”. In addition, we made minor revisions to the numbers shown in Figure 3B. The total overlapping differentially spliced AS events from all three groups of circadian neurons versus THs should be 8 instead of 10. Previously some header lines in the gene list files were counted towards the number of genes. These have been corrected, but the conclusions stay the same.

4) The authors could point out the use of other algorithms to identify AS changes – why rely only on JUM?

We chose the JUM algorithm to analyze the neuronal subpopulation datasets because of JUM’s unique feature for analyzing differential AS patterns completely independent of a priori knowledge of the host transcriptome annotation. It is well known that tissues such as neurons and gonads possess exceptionally diverse and dynamic cellular AS profiles associated with various cellular states and functions, greatly exceeding our current knowledge of the transcriptome annotations, especially in non-human model organisms such as Drosophila. Almost all other currently available methods for AS analysis rely on pre-built annotated AS splice junction libraries or transcriptome annotations, which result in vast underestimates of the tissue-specific AS complexity (Wang and Rio, 2017). We have specified the rationale for choosing JUM for AS analysis in subsection “There are many novel AS patterns in neuron subpopulations” in the revised manuscript. We have also explained the rationale and showed the superior performance of JUM on this particular point that the reviewer raised in detail in the preprint of the JUM paper (Wang and Rio, 2017).

5) A critical point of this study is the identification of "novel" AS junctions. Therefore, it will be important to have some quality control metrics to show the identified junctions are real rather than technical artifacts or alignment errors. It will be helpful to clearly describe parameters related to accuracy of exon junction identification (e.g., size of overlap on each side of the junction), and if possible provide an estimate of accuracy (e.g. using simulated synthetic exon junction reads by swapping 5' and 3' parts).

JUM accepts output files from the third-party RNA-seq alignment software tool STAR that performs read mapping and AS junction profiling with the following additional filter settings to control for the quality of AS junction detection and quantification: (1) JUM uses the two-pass alignment approach that has been shown to greatly improve novel splice junction quantification (Veeneman et al., 2016); (2) JUM only considers unambiguously mapped reads over splice junctions for junction profiling and quantification; and (3) JUM applies stringent quality filters for AS junctions. An AS junction will be taken for downstream analysis only if it receives more than 10 unambiguously mapped reads in every one of the biological replicates of the sample. These settings make sure that the AS junction analyzed by JUM has been observed consistently and reproducibly from the RNA-seq experiments. The accuracy of JUM-quantified novel AS junctions has also been tested by multiple rounds of qRT-PCR and RT-PCR validation experiments using experimental data from mouse embryonic cortical neurons and male fruit fly heads and the validation rates stay ~100% so far (Wang et al., 2013, Wang et al., 2016, Wang and Rio, 2017). In addition, JUM uses STAR as the aligner because STAR is among the best short-read alignment software developed so far for both comprehensiveness and accuracy (Baruzzo et al., 2017, Engstrom et al., 2013) and has been tested in numerous both simulated and experimental datasets on multi aspects of alignment performance. Given the extensive testing of STAR referenced above, estimating the performance of STAR using synthetic RNA-seq reads is beyond the scope of this paper.

We have indicated the mapping metrics used in this analysis for the data presented in this manuscript in subsection “RNA-seq data mapping for JUM” in the revised version as follows: “the minimum overhang length for splice junctions on both sides are set to be a default of 30 bp for non-canonical motif junctions, and 12 for the canonical GT/AG (CT/AC), GC/AG (CT/GC), AT/AC (GT/AT) motif junctions, respectively. The minimum uniquely mapped read count per junction is set to be the default 3 for non-canonical motif junctions and 1 for the canonical motif junctions. The minimum allowed distance to other junctions’ donor/acceptor is set to be the default 10 for non-canonical motif junctions and 0, 5, 10 for the canonical GT/AG (CT/AC), GC/AG (CT/GC), AT/AC (GT/AT) motif junctions, respectively. Maximum gap allowed for junctions are set to be follows: junctions supported by 1 read can have gaps <= 50000b; if supported by 2 reads then gaps <= 100000b; if supported by 3 reads then gaps <= 200000; if supported by more than 4 reads then gaps <= alignIntronMax.”

6) Related to the point above, it is also important to clearly state the criteria used to define a "novel" AS junction in the main text. As of now, the authors stated in the method section that "The novelty of each AS junction was compared against the library of annotated junctions in the UCSC genome browser transcriptome annotation (genome 422 version: FB2017_05)." It is unclear what type of annotations was used here (e.g., which gene models, does that include EST sequences?)

The AS junctions identified in this study were compared to the RefSeq gene transcriptome annotation database and it does not include ESTs. We have specified this in the manuscript, subsection “Profiling of AS junctions”.

7) In subsection “Circadian neurons present specific alternative pre-mRNA splicing patterns”, the authors identified AS junctions differentially spliced in each type of circadian neurons compared to TH neurons. They then examined the overlap of the lists to claim most of these junctions are specifically spliced in one subtype. To be more rigorous, one should perform differential splicing analysis between different types of circadian neurons, as the limited overlap could also be due to limited statistical power.

Similar issue for comparison of cycling exons in Figure 6B.

We have performed differential AS analysis between the circadian neurons DN1, LNd and LNv. We found vastly different AS patterns among each circadian neuron subpopulation – a total of 454, 778, and 325 AS events are significantly differentially spliced between DN1 versus LNd neurons, DN1 versus LNv neurons and LNd versus LNv neurons, respectively. This result supports our findings in Figure 3 and in the main text that each circadian neuron group possesses its own distinct AS profiles. We have added the results from this analysis in subsection “Circadian neurons present specific alternative pre-mRNA splicing patterns” and in Figure 3—figure supplement 2. Furthermore, we have conducted more analyses to explore the potential regulatory mechanisms that contribute to the unique AS profiles in circadian neuron subpopulations (detailed more in our response to #8 below). The limited overlap in differentially spliced AS events compared to TH neurons in the circadian neuron groups could be due to the observation that each circadian neuron group expresses its own unique set of RNA-binding protein/splicing factors compared to TH neurons (results added in Figure 3B). For example, we found that two LNv-specific splicing regulators mub and Syp targets ~36% of the AS events that are specifically differentially spliced in LNv neurons versus TH neurons. We have now included all of these additional results in subsection “Circadian neurons present specific alternative pre-mRNA splicing patterns”. It also needs to be noted that the reason for our approach to first perform differential AS analysis between each circadian neuron subgroup and the non-circadian TH neuron and then examining the overlap of each comparison was because we intended to profile for AS events that specifically affect circadian versus non-circadian functions. A direct comparison among the circadian subpopulations does not provide information that directly leads to our original intention and will also reveal AS events that are specialized in other functions in the circadian neurons.

Regarding Figure 6B, it is not trivial to compare two sets of 12 timepoint data (2 days 6 timepoints each) and ask whether the cycling tracks from two circadian neuron groups are statistically different from one another, as cycling determination involves multiple parameters including the shape of the curve, period length, amplitude requirements as well as minimum read requirements. However, there are several possible explanations for the neuron-specific cycling alternative splicing. Some of the alternative splicing events that cycle are neuron-specific (~40% in DN1, 14% in LNd and 25% in LNvs), so they not only cycle but are also only found in one subgroup of neurons. A smaller portion of alternative splicing events are present in all neuron groups yet cycle only in 1 (~20% in all neuron groups). In addition, we also found that each circadian neuron subpopulation possesses its own unique profile of cycling RNA-binding proteins/splicing factors that can potentially contribute to the neuron-specific cycling alternative splicing observed (Figure 7B). We have now included these additional results in subsection “Identification of cycling alternatively spliced isoforms in circadian neurons”.

8) While cell type- and time of day-dependent AS patterns are artfully presented here, a major weakness in this manuscript is that the underlying mechanisms driving AS or their consequences on protein expression are not studied. Regarding upstream mechanism, it would help to know if expression of any RNA binding proteins correlate with the observed AS events across cell types or times of day. Regarding the relative impact of AS on neuronal identity and function, a first pass analysis could be to compare the relative numbers and enriched ontologies between mRNA expression and AS differences between cell types and in the circadian data. A more compelling experiment would be to specifically inhibit an AS event within one of the neuronal populations to test its effect on cell function.

We agree that understanding the mechanisms behind neuron group-specific splicing is important and interesting. Inhibition of a particular AS event in one of the neuronal populations would be technically challenging. However, we have addressed this question as the reviewer suggests using a more extensive analysis of the RNA-seq data from the isolated neuronal populations, as well as publicly available databases of the targets of RNA-binding proteins (through CLIP and RIP experiments).

First, we profiled RNA binding proteins (RBPs) that are differentially expressed in each of the circadian neurons compared to the non-circadian TH neurons, respectively (New Figure 3B, and subsection “Circadian neurons present specific alternative pre-mRNA splicing patterns”). Interestingly, each circadian neuronal subpopulation expresses its own unique set of differentially-expressed RBPs compared to TH neurons, with very limited overlap, which is correlated with the observed unique differentially spliced AS patterns in each of the circadian neuronal populations compared to the AS patterns in TH neurons. We further profiled the targets of several of these RBPs using publicly-available data (Stoiber et al., 2015, Hansen et al., 2015) and found that these neuron-specific RBPs potentially target many AS events observed in each circadian neuronal subpopulation versus TH neurons (subsection “Circadian neurons present specific alternative pre-mRNA splicing patterns”). For example, Syb and mub, two RBPs that are specifically differentially-expressed in LNv neurons compared to TH neurons, target 23 and 32 LNv-specific differentially spliced AS events compared to TH neurons, respectively, with 5 overlapped AS events (~36% of total LNv neuron specific differentially spliced AS events). Although not definitive, this analysis supports the possibility that circadian neuron-specific RBP expression could drive the observed circadian neuron-specific AS patterns compared to TH neurons.

Second, we have performed a similar analysis to identify RBPs involved in splicing that cycles in the circadian neuron subgroups (New Figure 7B and Discussion section). Interestingly, each of the circadian neuron groups have distinct sets of cycling splicing-related transcripts that could contribute to the highly specific cycling splicing patterns observed in each circadian neuron group (Figure 7B). For example, the RNA-binding protein Qkr54B cycles in LNd neurons and one of its reported targets is Sgg (Stoiber et al., 2015), which displays cycling AS also in LNd neurons (Figure 6D). Furthermore, we examined the phase of all cycling AS structures in the three circadian neuron groups. Interestingly, each circadian neuron group exhibits distinct cycling phase patterns (Discussion section). Cycling AS structures in LNd cells peak throughout the day, while both DN1 and LNv cells showed higher levels of cycling AS structures at certain times of the day. DN1 cycling AS structures peak maximal in the late night and LNv cycling AS structure peaks show a bimodal distribution in mid-day and mid-night (Figure 7A and Figure 7—figure supplement 1). Moreover, some (but not all) of the cycling splicing factors found in each neuronal group showed a coordinated phase with their targeted cycling AS events. For example, the cycling RBP Qkr54B mRNA peaks 3-4 hours prior to the time of day dependent cycling inclusion of exon b in its targeted AS transcript Sgg (Figure 6D).

Third, we have compared the enriched gene ontologies between mRNAs and AS events that are differentially expressed or spliced, respectively, in each of the circadian neuron groups versus TH neurons. These comparisons suggest that neuron-specific differences in gene expression and splicing are impacting neuronal identity differently. For example, only one gene ontology category, integral membrane components, is found in common between the differentially spliced AS events and differentially expressed mRNA events in LNv neurons. However, neuropeptide receptors are distinctively enriched in differentially expressed mRNAs while potassium ion channels are distinctively enriched in differentially spliced AS splicing events in LNv neurons.

Lastly, we have used a similar approach to compare GO terms of neuron-specific cycling mRNAs and cycling AS structures. There is very little overlap between the types of genes that cycle in each of the circadian neuron groups at the mRNA level and at the level of splicing. The one exception is that the GO term chemical synaptic transmission is found in both cycling mRNAs and cycling AS structures in LNd neurons. One interesting observation is that the top enriched GO term for mRNAs that undergo cycling alternative splicing in DN1 neurons are mRNA binding proteins themselves (Figure 6—figure supplement 3 and Supplementary file 3). This list includes 8 known factors that impact splicing and the time-of-day dependent alternative splicing of these transcripts could influence AS profiles in DN1s.

These results above suggest that the neuronal-specific RBP expression and cycling patterns can be an important factor that drive the circadian neuron-specificAS profiles and cycling observed. That said, many other factors can contribute to the AS changes and patterns observed in each of the neuronal subpopulations. Possible mechanisms include, but not limited to, differential signaling pathway inductions, post-translational modifications and activation of splicing regulators and chromatin modifications. We believe these can be interesting and crucial topics for future studies that beyond the scope of this manuscript. The results above also suggest that AS and transcriptional regulation are two independent layers of gene regulation in specifying neuronal identify and function.

9) One potential concern is the technical/biological variation between replicates, which is expected from the small number of cells used to generate RNA-seq libraries and visible from the examples shown in the figures. It will be helpful to have some data to show the reliability/reproducibility of observed differential splicing (e.g., correlation of splicing pattern of replicates, or even independent RT-PCR validation if possible).

For software tools like MISO and MATS that can only perform pair-wise comparisons and cannot handle biological replicates, a correlation analysis of splicing patterns from replicates is indeed necessary. However, JUM integrates biological variability from multiple replicates into its core statistical model to make sure that the reported significantly changed AS events are supported by all replicates consistently and reproducibly. JUM achieves this by modeling the reads that map to an AS junction as a negative binomial distribution and fitting a mean-variance function through an inferred dispersion parameter with an empirical Bayesian approach across the observed variability from all replicates. This approach has been described previously, optimized recently and been widely applied in software tools analyzing RNA-seq experiments with replicates, and has been proven to greatly improve the quality of analysis by taking biological replicates into account (Anders et al., 2012, Love et al., 2014, Lu et al., 2005, McCarthy et al., 2012, Robinson et al., 2010, Robinson and Smyth, 2007, Robinson and Smyth, 2008). As a result, the nature of the model design in JUM already takes into account the reliability/reproducibility issues raised here by the reviewer. We also want to point out that we have consistently validated all JUM-reported AS events using multiple rounds of qRT-PCR in mouse embryonic cortical neurons and male fruit fly heads (Wang et al., 2013, Wang et al., 2016, Wang and Rio, 2017). This experimental validation further supports the accuracy of JUM’s predictions.

10) The observation of alternative exons that are differentially spliced in specific types of neurons or time-of-the-day during the circadian clock is very interesting. Is there any indication of corresponding changes in splicing factors that potentially regulate these AS events?

We have addressed this question in detail in our response to #8, above.

11) It will be helpful to provide additional functional analysis of identified AS events with differential splicing, in addition to GO. This reviewer understands that detailed functional dissection of specific events might be beyond the scope of this manuscript, but it should be feasible to perform statistical analysis of cross-species conservation, preservation of reading frame, exon size, etc., which is frequently distinct for exons with cell type specific splicing.

We have performed statistical analyses to compare cross-species conservation and preservation of reading frames in the cassette exons that are alternatively spliced versus the ones that are not alternatively spliced between DN1, LNd and LNv neurons versus the non-circadian TH neurons, respectively. For each comparison of a circadian neuron subpopulation versus TH, we compared the PhastCons conservation score of each identified cassette exon across 27 species of insects and also the size of the cassette exon to see if the alternative inclusion of the cassette exon changes the reading frame of the transcript. We found that neuron-subgroup-specific, alternatively spliced cassette exons are somewhat more conserved (although not statistically significant) and tend to preserve the reading frame of the transcript, indicating these alternatively spliced cassette exons are under more evolutionary selection to preserve their functions, which may be crucial to the identity and function of the neuronal subpopulations. These results are listed in Figure 3—figure supplement 3 and subsection “Circadian neurons present specific alternative pre-mRNA splicing patterns” in the revised manuscript.

12) In Figure 4B, although the authors claim, "in TH and LNd neurons the Shab transcripts preferentially utilize an alternative terminal exon, […] while in LNv and DN1 neurons, the long Shab isoform that encodes all transmembrane regions is dominant", it is not intuitively showing that by looking at the Sashimi plots. The authors should order the samples similarly to Figure 4A, and also may better prove their point by grouping the LNv and DN1 samples or highlighting the relevant differences in the figure read numbers.

We have highlighted the relevant differences in the figure read numbers in Figure 4B as the reviewer suggested.

Author details

1. Qingqing Wang

   1. Department of Molecular and Cell Biology, University of California, Berkeley, Berkeley, United States
   2. Center for RNA Systems Biology (CRSB), University of California, Berkeley, Berkeley, United States
   3. California Institute for Quantitative Biosciences (QB3), University of California, Berkeley, Berkeley, United States

   Contribution

   Software, Formal analysis, Validation, Investigation, Visualization, Methodology, Writing—original draft, Writing—review and editing

   Contributed equally with

   Katharine C Abruzzi

   Competing interests

   No competing interests declared

   "This ORCID iD identifies the author of this article:" 0000-0002-0836-1367

2. Katharine C Abruzzi

   1. Department of Biology, Howard Hughes Medical Institute, Brandeis University, Waltham, United States
   2. National Center for Behavior Genomics, Brandeis University, Waltham, United States

   Contribution

   Data curation, Formal analysis, Validation, Investigation, Writing—original draft, Writing—review and editing

   Contributed equally with

   Qingqing Wang

   For correspondence

   katea@brandeis.edu

   Competing interests

   No competing interests declared

   "This ORCID iD identifies the author of this article:" 0000-0003-3949-3095

3. Michael Rosbash

   1. Department of Biology, Howard Hughes Medical Institute, Brandeis University, Waltham, United States
   2. National Center for Behavior Genomics, Brandeis University, Waltham, United States

   Contribution

   Conceptualization, Supervision, Funding acquisition, Writing—review and editing

   Competing interests

   No competing interests declared

   "This ORCID iD identifies the author of this article:" 0000-0003-3366-1780

4. Donald C Rio

   1. Department of Molecular and Cell Biology, University of California, Berkeley, Berkeley, United States
   2. Center for RNA Systems Biology (CRSB), University of California, Berkeley, Berkeley, United States
   3. California Institute for Quantitative Biosciences (QB3), University of California, Berkeley, Berkeley, United States

   Contribution

   Conceptualization, Supervision, Funding acquisition, Writing—original draft, Writing—review and editing

   For correspondence

   don_rio@berkeley.edu

   Competing interests

   No competing interests declared

   "This ORCID iD identifies the author of this article:" 0000-0002-4775-3515

• 2,850

  Page views

• 438

  Downloads

• 17

  Citations

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