Introduction

Throughout the history of life on Earth, organisms have had to adapt to a constantly changing environment, including the ∼24-hour daily rhythm of light/dark, driving the development of endogenous biological clocks. The circadian clock, which is entrained by external stimuli such as light, enables the organism to anticipate the onset of the light and dark phases and synchronize its physiology and behavior in harmony with the environment. This, in turn, enhances the organism’s fitness and survival13. From single-celled organisms to metazoans, circadian clocks have evolved multiple times, highlighting their importance to living organisms1,2. Despite the fundamental role circadian clocks play in regulating the rhythmicity of living organisms, their evolutionary origin and intricate molecular mechanisms remain ambiguous in early diverging animal lineages, such as cnidaria.

Rhythmic phenomena, including calcification, reproduction, and diel behavior patterns, have been examined in cnidarian species38. While some of these patterns were found to be directly triggered by environmental stimuli, such as light, others persist in the absence of external cues, suggesting the presence of an internally generated and self-sustaining circadian clock6,9,10. At the molecular level, cnidarians possess homologs of putative core pacemaker genes found in bilaterians1113. Several studies have shown that most of these genes display diel expression patterns under light/dark cycles. However, in contrast to most animals, their oscillation generally ceases without the presence of light cues1416. Thus, how the core pacemaker genes orchestrate rhythmic gene expression and circadian behaviors in cnidarians remains unclear.

One of the most studied cnidarian species in the field of chronobiology is the estuarine sea anemone, Nematostella vectensis. Few studies have shown that in diel lighting, the locomotor behavior of Nematostella has a ∼24-h rhythm that is maintained under constant conditions suggesting it is regulated by an endogenous circadian clock6,1719. In accordance, the Nematostella genome codes for conserved core pacemaker genes such as NvClk, NvCycle, and the cryptochromes NvCry1a and NvCry1b12,13. The proposed circadian clock model in Nematostella is composed of the positive transcription factors (bHLH-PAS family), NvCLK and NvCYCLE, that heterodimerize and upregulate light-dependent cryptochrome genes in the feedback loop, and NvPAR-bZIPs in the feed-forward loop, which repress the transcription of the positive elements11,12,19. More recently, the NvCLK-interacting pacemaker, NvCIPC was predicted to act as an additional repressor of the NvCLK:NvCYCLE dimer19. However, in contrast to the free-running oscillation demonstrated for Nematostella behavior6,17,19, transcriptional expression profiles of most candidate genes implicated in the pacemaker do not retain their oscillation period in the absence of light14,15,19.

Here, we employed the CRISPR/Cas9-mediated genome editing system to establish a NvClk mutant (NvClk1-/-) Nematostella. By combining behavioral monitoring and transcriptomic analysis, we aimed to elucidate the role of NvClk in regulating rhythmic locomotor activity and gene expression under varying light conditions. Our study revealed a robust light response pathway capable of compensation and a conserved function of CLOCK as a timekeeper in the absence of a light cue.

Results

Phylogenetic analysis and spatial expression pattern of NvClk

Phylogenetic analysis of NvCLK protein sequences revealed that NvClk is positioned within the cnidarian branch (Fig. 1a). It contains a basic helix–loop–helix (bHLH) DNA binding domain and two Per-Arnt-Sim (PAS) domains, similar to the protein structure found in other animals. PAS domains are crucial structural motifs involved in protein-protein interactions that drive the self-sustaining molecular mechanism underlying the circadian clock20,21.

NvClk1-/- cannot maintain circadian behavior in non-diel light conditions.

(a) Phylogenetic tree showing the evolutionary relationship of CLK orthologs across different animal species. (b) In situ hybridization of NvClk in the WT juvenile with scale bars representing 0.1mm. (c) Schematic representation of the NvClk gene in grey, with the open reading frame (ORF) in dark grey and the conserved protein domains bHLH (yellow) and PAS1 and PAS2 (dark red). The CRISPR generated NvClk1 allele has a +20nt insertion after the PAS1 domain, represented by a black arrowhead. NvCLK dimerizes via its three functional domain with NvCYCLE binding the CACGTG ebox to drive rhythmic transcription. (d-g-j-m). Normalized Movement (a.u), hourly binned over 72h, under different light conditions: 12h light:12h dark, continuous dark (Dark - Dark), continuous light (Light - Light), and 6h light:6h dark. The black line represents the WT and the red line represents the NvClk1-/- mutant. (e-h-k-n) Lomb-Scargle Periodograms for each corresponding light condition. The significant period value (p<0.01) is indicated for each genotype in the top left corner of each graph. (f-i-l-o) Phase detection (Cosinor) and genotype comparison of 24h rhythmic individuals. See number nrhythmic/ntotal on the x-axis indicating for each genotype the number of 24h-rhythmic animals over the total population.

To localize NvClk expression at the polyp stage, in situ hybridization chain reaction (HCRv.3) was performed. Polyps were sampled at ZT10, i.e., peak expression of NvClk6,12,22. NvClk specific expression was observed throughout the animal tissue, and enriched expression was visible in the tentacle endodermis and mesenteries, while no signal was observed in the negative control (Fig.1b, Extended Data Fig.1a). This expression pattern resembled the expression observed at the larvae stage22. To date, functional manipulation of the NvClk gene has not been performed in basal animal lineages including cnidarians, and its function is unknown in cnidarians20,21 (Fig.1a).

Generation of NvClk1-/- Nematostella

To investigate the function of NvClk in Nematostella, we employed the CRISPR-Cas9 system to generate mutants. Based on existing knowledge from mouse and Drosophila models, we hypothesized that NvCLK:NvCYCLE dimer binds to the DNA motif CACGTG within the promoter of rhythmic target genes (Fig. 1c). Guide RNA (gRNA) was synthesized to target a region between the two PAS domains of the NvClk coding sequence (CDS). This gRNA, along with the Cas9 endonuclease, was microinjected into zygotes (Methods). Subsequently, F0 animals were outcrossed with wild type (WT), and the F1 progeny were raised to adulthood. Genotyping of F1 polyps identified 10 different mutated alleles, with six displaying a frame-shift mutation, including one with a 20 bp insertion (NvClk1), resulting in a premature stop codon (Extended Data Fig. 1b). The predicted 203 amino acid truncated protein lacked 459 amino acids, including one co-factor dimerization PAS domain (Fig. 1c, Extended Data Fig. 1b). To obtain homozygous NvClk1-/- polyps, we crossed heterozygous NvClk1 F1 animals. Genotyping of F2 polyps confirmed the expected 25% frequency of NvClk1-/- mutants. Subsequently, we intercrossed NvClk1-/- animals to obtain F3 NvClk1-/- polyps for use in subsequent experiments aimed at assessing the role of NvClk in regulating behavioral and genetic rhythms.

NvClk is necessary to maintain circadian behavior in constant conditions

To assess the impact of the NvClk1 mutation on circadian rhythm, we monitored locomotor behavior of WT and NvClk1-/- polyps under different light conditions (Supp. Table 1). Both the WT and NvClk1-/- populations exhibited a 24-hour periodicity in 12:12 hour Light-Dark (LD) cycles (Fig. 1d,e), with 15 out of 21 WT animals displaying 24-hour rhythmicity compared to only 9 out of 20 NvClk1-/- animals (Fig. 1f, Table 1, Supp. Table 2). The average acrophase for WT polyps (13.3 hours) was significantly lower than for NvClk1-/- polyps (17.3 hours) (Fig. 1f, Table 1).

Summary of rhythmic analysis of individual behavior.

We then investigated locomotor behavior under continuous conditions, namely continuous dark (DD) or continuous light (LL). WT polyps exhibited a 22-hour rhythmic behavior under both constant light conditions, with 17 out of 25 WT polyps displaying a 24-hour rhythm under DD and 7 out of 25 under LL (Fig. 1g-l). In contrast, few NvClk1-/- polyps displayed rhythmic behavior under constant conditions (1 out of 24 in DD and 1 out of 26 in LL) (Table 1). Additionally, we observed an intermediate phenotype in the locomotor behavior of heterozygous polyps for the NvClk1 allele in DD (Extended Data Fig. 1c-f).

To further explore the 24-hour rhythm of NvClk1-/- polyps under LD, we tracked locomotor activity under a 6-hour light: 6-hour dark (LD 6:6) cycle after a regular diel 72-hour entrainment under 12:12 LD. While WT polyps maintained a marginally significant periodicity of 22 hours, NvClk1-/- polyps displayed a 12-hour rhythm at the population level (Fig. 1m-o). Specifically, we identified an important difference of 12- hour rhythmic individual polyps between WT and NvClk1-/- groups (1 out of 25 WT polyps vs. 13 out of 26 NvClk1-/- polyps) (Table 1). Notably, entrainment with LD 6:6 did not lead to a 12-hour rhythm in DD for both WT and NvClk1-/- polyps (Extended Data Fig. 1g-i).

NvClk regulates rhythmic gene expression differentially in response to light conditions

To investigate the underlying molecular correlates of the behavioral phenotype found in NvClk1-/- polyps, we conducted transcriptional profiling. WT and NvClk1-/- polyps were sampled seven times every 4 hours over 24 hours under LD and DD conditions (Fig. 2a). We defined true rhythmic genes using a combination of both RAIN (rhythmicity analysis incorporating non-parametric methods) and JTK_CYCLE algorithms (Fig. 2b, Supp. Table 3). Within our dataset, we identified 119 genes rhythmic under LD and 107 rhythmic genes under DD in WT polyps (Fig. 2b). In NvClk1-/- polyps, we detected 147 rhythmic genes under LD and only 37 under DD (Fig. 2b).

NvClk1-/- shows rhythmic genes reduction in constant darkness with altered rhythmic features.

(a) Overview of the experimental design used to generate RNA-seq data. Polyps were entrained for 72 hours before sampling at 4-hour intervals over a 24-hour period (dark arrows) in both LD and DD cycles. (b) Venn diagram comparing the total number of 24h rhythmic genes identified in WT and NvClk1-/- in LD and DD cycles with a p <0.01 with RAIN and JTK. (c) Average acrophase comparison between rhythmic genes in LD and DD in WT polyps. Mann-Whitney test, p<0.001. (d) Average acrophase comparison between rhythmic genes in LD and DD in NvClk1- polyps. Mann-Whitney test, p>0.05. (3) Average relative amplitude comparison between rhythmic genes in LD and DD in WT polyps. Mann-Whitney test, p<0.05. (d) Average relative amplitude comparison between rhythmic genes in LD and DD in NvClk1-/- polyps. Mann-Whitney test, p>0.05. (c-f) sample size (n) indicated below each boxplots.

The rhythmic genes in WT polyps displayed a higher acrophase under DD compared to LD (17.20h vs 12.93h, Fig. 2c). However, no differences were detected between LD and DD rhythmic genes in NvClk1-/- polyps (Fig. 2d). Similarly, the relative amplitude (rAMP) of DD rhythmic genes was higher in WT polyps compared to LD (0.61 vs 0.43, Fig. 2e), but no rAMP difference was observed between LD and DD rhythmic genes in NvClk1-/- polyps (Fig. 2f).

We performed a clustering analysis on the rhythmic genes using the DPGP model (Dirichlet process Gaussian process mixture model). The average number of genes per cluster between LD and DD conditions in WT polyps did not differ significantly (7.3 vs 7.6, Fig. 3a, Supp. Table 4). However, in NvClk1-/- polyps, the average number of genes per cluster was significantly lower in DD compared to LD (4.1 vs 8.6, Fig. 3b). Specifically, in WT DD polyps, we observed clusters with high numbers of genes peaking at subjective night (Fig. 3c, Supp. Table 4). We did not identify any GO term enrichment to any of the clusters.

NvClk1-/- shows temporal gene clusters recruitment alteration in constant darkness.

(a) Gene per cluster count comparison between LD and DD in WT polyps. Number of cluster indicated below each boxplot. Mann-Whitney, p > 0.05. (b) Gene per cluster comparison between LD and DD in NvClk1-/-. Sample size (n) indicated below each boxplot. Mann-Whitney, p < 0.05. (c) Representative temporal gene cluster for each condition. For both light conditions, from top to bottom cluster peaking at day, then twilight, then at night. Number of gene for each specific cluster indicated on the x-axis.

NvClk regulates temporal expression pattern of pacemaker genes

In line with previous findings in Nematostella12,14, candidate pacemaker genes showed arrhythmic expression under DD conditions (Fig. 4a, Supp. Table 3). However, the altered expression patterns observed in LD NvClk1-/- polyps, compared to LD WT polyps (Fig. 4a, Supp. Table 3), suggest a complex regulatory effect.

NvClk1-/- alters temporal pacemaker genes expression.

(a) Four pacemaker gene are individually plotted showing the read count over 24h in LD and DD in WT (black) and NvClk1-/- (red). The continuous line represents significant rhythmicity (RAIN&JTK p<0.01) while the dashed line indicates no rhythmicity. (b) correlation matrix of candidate pacemaker genes expression in LD for WT on the left and NvClk1-/- on the right. (c) schematic representation of the promoter sequences analyses 5kb upstream the putative ATG. Black boxes represent canonical E-boxes while circadian E-boxes are green. Below the logo motif we used to identify canonical and circadian Ebox. (d) Circadian / Canonical ratio (in %) per condition. Kruskal-Wallis, multiple comparison, a vs b : p<0.05.

To systematically analyze the mutation’s impact on these genes, we constructed a correlation matrix based on their expression levels. In LD WT polyps, genes clustered into two groups: one with NvClk peaking during the day, anticorrelated to the second cluster (NvPar-bzipc and NvCipc) peaking at night. In LD NvClk1-/- polyps, the second cluster contained two additional genes and showed weaker anticorrelation with the NvClk cluster (Fig. 4b).

To go further into the regulatory mechanisms downstream of the pacemaker, we examined the presence of circadian E-box motifs (CACGTG) within 5kb upstream of the predicted ATG of rhythmic genes. We calculated circadian/canonical E-box enrichment to account for the total variation in the number of canonical E-boxes (Fig. 4c). Notably, only the candidate pacemaker genes exhibited a significant enrichment in circadian E-boxes in their promoters (15.9%) compared to the WT (5.6%), NvClk1-/- (4.8%) rhythmic genes, and non-rhythmic genes (6.8%) (Fig. 4d).

NvClk coordinates cell division and neuronal pathways in constant darkness

In addition to the transcriptomic rhythmic analysis, our aim was to identify processes regulated by NvClk that may not necessarily exhibit rhythmicity. To achieve this, we conducted a differential gene expression analysis on the total transcriptome between genotypes under each light condition. Under LD conditions, NvClk1-/- polyps exhibited 457 down-regulated genes and 646 up-regulated genes, with no significant enrichment in GO terms observed (Fig. 5a, Supp. Table 4 and 5). However, in DD conditions, NvClk1-/- displayed 2450 down-regulated genes and 1770 up-regulated genes (Fig. 5b, Supp. Table 4). Notably, the down-regulated genes were predominantly enriched in processes related to mitosis, microtubules, and ciliary/flagellar motility. Conversely, the up-regulated genes showed significant enrichment in processes such as the modulation of another organism’s processes, axonal guidance, and sensory perception (Fig. 5b, Supp. Table 5).

NvClk1-/- disrupts cell-cycle and neuronal pathways in constant darkness.

(a) Volcano plot showing the differential expression of genes (DEG) between WT and NvClk1-/- in LD (Left) and DD (right). Dashed line indicate the threshold used to detect DEG (p.adj<0.01). Red dots indicate down regulated genes and black dots up-regulated genes in NvClk1-/- compare to WT polyps (b) Gene Ontology (GO) terms with with significant fold-enrichment (Bonferroni corrected p-value or p.adjusted <0.01) for the DEG analysis in DD. Down regulated genes in Red while Up regulated genes in Black.

Discussion

Conserved behavioral CLOCK function through animal evolution

Our study provides valuable insights into the evolution of circadian clocks by characterizing the effects of the first Clock mutation in a cnidarian, the sea-anemone Nematostella vectensis. Through our behavioral assays, we have shown that NvClk is essential for maintaining rhythmic locomotor activity in the absence of an entraining light cue. Although the rhythmicity of the NvClk1+/- heterozygote polyps was affected in DD, our results could not discriminate a dominant-negative from a full loss of function to identify the nature of this mutation (Extended Data Fig. 1g-i). Studies in various model organisms further support the importance of CLOCK in regulating circadian locomotion. For instance, homozygous DmClk-/- flies exhibit a loss of circadian locomotion in constant darkness23,24, while the DnClk1a zebrafish mutant displayed a shortened period under the same conditions25. Interestingly, MmClock-/- mice maintain robust circadian locomotion in constant darkness, suggesting potential functional redundancy among mammalian clock genes2629. Overall, these findings support a conserved role of CLOCK in preserving circadian behavioral rhythms in absence of light cues across the distant Nematostella, flies, zebrafish, and mice.

Moreover, the conservation of a 24-hour locomotion rhythm in LD of the NvClk1-/- polyps with a different mean acrophase, revealed a light-response pathway independent of the circadian circuit, consistent with observations in other animal models23,26,30 (Fig. 1f). NvClk1-/- polyps exposed to a 12:12h LD cycle exhibited a 24-hour period, whereas those exposed to a 6:6h LD cycle displayed a 12-hour period. Notably, nearly no WT polyps exhibited a 12-hour rhythm under this condition, suggesting that the circadian clock overrides the light-response pathway (Fig. 6a). While some of the circadian factors can directly sense the light such as CRY proteins31, 29 putative NvOpsin have been identified in the genome which could be involved in the light-response pathway32. Behavioral tracking of NvClk1-/- polyps exposed to different wavelength could help to identify candidates for further functional studies of the light-response pathway.

Summary of NvClk function in the regulation of Nematostella circadian rhythmicity.

Transcriptional rhythmicity plasticity downstream NvClk

At the transcriptomic level, previous studies in Nematostella have shown large changes in the transcriptional profile of many genes after a single day of constant darkness, including the candidate pacemaker genes that were found arrhythmic, despite sustain circadian locomotion12,14,22. Consistent with previous transcriptomic analysis in cnidarian7,14,16,33, most of rhythmic genes identified in LD were different from those identified in DD in the WT polyps. Notably, they displayed higher mean acrophase and larger mean amplitude in DD suggesting a differential regulation in response to light conditions, which has not been investigated in previous cnidarian studies. Additionally, an overlap between our LD rhythmic genes and Leach et al. dataset14 suggests the robustness of pacemaker rhythmic transcription in LD conditions compared to downstream rhythmic genes (Extended data Fig. 2). Overall, these observations suggest a potential plasticity of the pacemaker in selecting specific genes depending on environmental conditions.

Similar to the rhythmic behavior of NvClk1-/- polyps under LD conditions, the identification of 24-hour rhythmic genes suggests potential compensation by the light-response pathway. However, most of rhythmic genes in the NvClk1-/- compared to WT polyps in LD conditions are different. This suggests the light-response pathway may not fully reproduce the normal circadian pattern observed in WT polyps and would require further investigation to understand their recruitment and function. Furthermore, the identification of a reduced number of rhythmic genes in DD NvClk1-/- polyps revealed the importance of NvClk in maintaining molecular rhythm in absence of light. However, the persistence of these rhythmic genes raises questions regarding their origin. They could be the result of the mutation itself, a compensatory mechanism or false positives resulting from our analysis. The study by Ray et al., which report controversial findings regarding RNA rhythmicity in the Bmal1-/- mouse cells (NvCycle orthologous in Nematostella) is an interesting case of rhythmicity analysis to consider3438. Our attempt to identify rhythmic genes faced challenges, particularly in the choice of statistical parameters. Initially, using stringent criteria like BHQ for the JTK method or p.adj with RAIN resulted in a very limited number of genes for analysis (Supp. Table 3). Consequently, we opted for a balanced approach combining both algorithms with a p<0.01 threshold. We made a compromise to identify rhythmic genes with greater confidence, but we are aware that the methodological choices are critical for the results.

The clustering analysis revealed that rhythmic genes can be categorized by expression pattern similarity. This approach defined “temporal gene clusters” as group of seven/eight genes in average in the WT (Fig. 6b). Their discovery poses a fundamental question: How a group of genes are co-regulated in times and space (cell types) by the pacemaker? Their recruitment is disrupted in the DD NvClk1-/- polyps suggesting an essential function of NvClk in absence of light. The combination of published scAtlas39 and multiplexed FISH techniques40 will be essential to further investigate the biological regulation and function of these temporal gene clusters.

NvClk temporally organize pacemaker genes expression to drive rhythmic genes recruitment

Our study reveals that NvClk plays a key role in regulating the temporal transcription of pacemaker candidate genes (Fig.4a). Our analysis identified two clusters of pacemaker genes: One containing NvClk, and a second one containing a potential NvClk inhibitor (NvCipc)41,42. The alteration of clusters composition with a weaker anticorrelation in LD NvClk1-/- polyps might generate a desynchronization of the pacemaker factors availability. Indeed, regulation of rhythmic transcription involved a complex protein-protein-DNA timing interaction. Furthermore, we did not identify any circadian E-boxes enrichment in rhythmic genes between conditions, except for the candidate pacemaker genes. Altogether, this supports the function of NvClk in orchestrating the timing interaction of pacemaker factors to select downstream rhythmic genes, indicating a more complex regulatory landscape at play.

However, one significant unanswered question in our study is the reason for the arrhythmic transcription of putative pacemaker genes in DD. We hypothesize that using whole animals for sampling material might mask oscillating gene expression signals, especially if signals are present in a small number of cells or if tissues exhibit rhythmic gene expression in different phases. Furthermore, it is important to acknowledge a limitation in our interpretation, common in chronobiology: the use of RNA oscillation as a proxy for protein oscillation and function. The development of tools to study the pacemaker factors at the protein level in Nematostella will leverage this limitation in the field.

NvClk regulates processes involved in cell proliferation and the neural system in absence of light

Our study of NvClk suggests coordination of cellular processes, especially in the absence of light. The results of our rhythmic transcriptomic analysis (Fig. 2-3) raised questions regarding indirect effects, and the non-rhythmic function of NvClk. We performed a differential gene expression analysis on the total transcriptome for each light condition. Under LD conditions, NvClk1-/-polyps exhibited significant changes in gene expression but no GO term enrichment was found (Fig. 5a, Supp. Table 5) revealing multiple altered processes we cannot yet identified. In contrast, under DD conditions, NvClk1-/- polyps displayed more pronounced alterations, with more DEGs and enriched GO-terms for down-regulated genes related to mitosis, microtubule organization, and ciliary/flagellar motility, while the up-regulated genes showed enrichment in processes such as the modulation of other organism’s processes, axonal guidance, and sensory perception (Fig. 5b, Supp. Table 5). These results imply that NvClk has non-circadian functions that are dependent on light availability. This is particularly noteworthy considering the expression of core pacemaker genes, which are known to be arrhythmic during larvae stages, potentially involved in developmental processes22.

Taken as a whole, this study provides novel insights into circadian regulation in Nematostella vectensis and sheds light on the evolutionary origin of circadian time maintenance. Our findings indicate that CLOCK function is conserved from Cnidaria to mammals to maintain rhythmicity in absence of light cues. Furthermore, it revealed a light-response pathway able to compensate at both behavioral and molecular level using light cues. This circadian clock mutant opens new avenues for investigating cell-type-specific mechanisms of the circadian clock that drive the molecular and phenotypical oscillations of cnidarians. By further exploring the circadian clock mechanisms in cnidarians, we can gain deeper insights into the evolutionary origins of this critical aspect of biology, enhancing our understanding of how organisms have evolved to keep track of time and adapt to their environment.

Methods

Animal husbandry

Nematostella were grown in 12g.L sea salt water at 17°C, and fed with Artemia salina nauplii three times a week. Spawning of gametes and fertilization was performed according to a published protocol43. In brief, temperature was raised to 25 °C for 9 h and the animals were exposed to strong white light. Three hours after the induction, oocytes were mixed with sperm to allow fertilization.

CRISPR/Cas9 mediated mutagenesis

Genome editing in Nematostella was carried out following established CRISPR/Cas9 protocols, with slight modifications44,45. ZiFiT targeting software (http://zifit.partners.org/)46 was used to select a guide RNA (gRNA) target site within the beginning of the NvClk exon and to design complementary oligonucleotides. To ensure specificity of the gRNA, the selected target site sequence (GGTCCTCTCGTGGACTCTAC) was BLASTed against Nematostella vectensis genome (using JGI expected E value threshold of 0.1 to adjust for short sequences: http://genome.jgi.doe.gov/Nemve1/Nemve1.home.html). To generate the gRNA template, the following oligonucleotides were used:

Oligo 1: 5’- TAGGTCCTCTCGTGGACTCTAC -3’

Oligo 2: 5’- AAACGTAGAGTCCACGAGAGGA -3’

To construct the gRNA expression vector, pDR274 (plasmid # 42250; Addgene) was digested with the BsaI restriction enzyme. Subsequently, the gRNA oligonucleotides were annealed and cloned into the BsaI-digested pDR274 vector. Next, DraI-digested gRNA expression vectors, purified via ethanol precipitation followed by PureLink PCR purification kit (Invitrogen), were transcribed and purified using HiScribe T7 High Yield Transcription Kit (New England BioLabs) and illustra Microspin G-50 Columns (GE Healthcare Life Sciences), respectively. Cas9 recombinant protein with nuclear localization signal (260 ng/μl; PNA Bio, USA) was co-injected with the gRNA (140ng/µl) into Nematostella zygotes. Injected embryos were raised in petri dishes at 22°C under constant darkness with daily water changes.

CRISPR/Cas9 mediated mutagenesis screening

To evaluate genome editing efficiency and mosaicism in F0 animals, genomic DNA flanking the guide sequence was amplified and Sanger sequenced. PCR was performed using two strategies. For the first, PCR reactions were performed using individual injected Nematostella (7-days post-fertilization), directly pipetted into a 25- μl PCR reaction containing a 2x concentration of PCR MasterMix (Tiangen) and 10 pmol of each PCR primer. For the second, genomic DNA was extracted from tissue sampling from live adult animals, after relaxation in 7% MgClphy2 (Sigma-Aldrich), using NucleoSpin Tissue DNA purification kit (MACHERY-NAGEL). Subsequent PCR reactions were performed as above using 50-ng of genomic DNA. The primers used for these reactions (listed below) were designed to amplify a ∼750bp region around the targeted NvClk genomic locus. Mosaicism was determined if sequenced PCR products showed overlapping peaks in their chromatograms. The second strategy, which takes advantage of the ability of Nematostella to fully regenerate within a few days47,48, is the one we refer to in the text hereafter. The injected individuals determined mosaic mutants were raised as F0 founders to sexual maturity and outcrossed with wild-type animals. The progeny of these crosses was raised and individually genotyped as described above. To determine inheritable mutations, sequences were further analyzed using the Tracking of Indels by DEcomposition web-tool (TIDE). TIDE quantifies editing efficiency and identifies the predominant types of DNA insertions and deletions (indel) mutation composition from a heterogeneous PCR product compared to a wild-type sequence49. Different heterozygous mutants were raised to sexual maturity and outcrossed with wild-type animals. The resulting F2 progenies were then raised to sexual maturity and genotyped before spawning for F3. Heterozygous mutants from each F2 progeny were intercrossed to obtain 25% homozygous F3 mutants. All animals used in this study are derived from heterozygous F2 mutants intercrosses, harboring the mutant allele NvClk1. PCR genotyping was performed using the following primer:

Forward 5’- GATAAACACGGGCCGAAGATA -3’

Reverse 5’- CAGTCCACGCTGGTCTAAAT -3’

Determination of NvClk1 F3 mutant genotypes

Genomic DNA was extracted as described above and used for following PCR and electrophoresis-based genotyping. PCR primers (listed below) encompassing the NvClk targeted site were used to produce PCR products of approximately 100 bp. The PCR products were then loaded and migrated by electrophoresis on a 3% Tris-borate- EDTA (TBE) agarose gel supplemented with GelStar Nucleic Acid Gel Stain (Lonza) for approximately one hour. The genotype of each F3 animal was determined by visualizing differences in migration speed of the PCR products caused by the CRISPR/Cas9 genome editing. The homozygous mutant animal (NvClk1−/−) produces only the larger ∼120 bp amplicon while the wild-type animal (NvClk1+/+) produces only the lower ∼100 bp amplicon. Animals heterozygous for the deletion (NvClk1+/−) produce both the larger mutant and the smaller wild-type amplicons. PCR genotyping was confirmed by subsequent DNA sequencing of selected F3 animals.

PCR was performed using the following primer:

Forward 5’- ACCCCACTGAGTGACCTCTT -3’

Reverse 5’- ATACGCCTGCGCTATACACC -3’

Behavioral assays

Locomotor activity of individual Nematostella were monitored using a lab-made setup equipped with an IP Infra-Red camera (Dahua Technology, Hangzhou, China), a white neon illumination (Aquastar t8, Sylvania Lightning Solution,) and constantly illuminated with low-intensity infrared (850nm) LED light. The camera output 1h mp4 movie files which were AVI converted then stitched. The data collection and analysis were carried out by EthoVision XT8 video tracking software (Noldus information technology, Wageningen, Netherlands). Animals were isolated in wells of six-well plates, each of which was manually defined as a tracking ‘arena’ in the EthoVision software. Center- point detection with gray scaling (detection range of 25–77, contour erosion of 1-pixel, high pixel smoothing) was used to monitor movements, which were calculated according to the change in position of the average center pixel each second.

Illumination was provided with an intensity of 17 PPFD (+/−2) and did not significantly affect the experimental temperature (20 °C). The illumination cycles were 12: 12 h Light-Dark, 6: 6h Light-Dark or continuous light (LL). Parameters were optimized to ensure that organisms were detected throughout the entire observation period.

Behavior analysis

The total distance moved was summed in hourly bins and individually normalized min/max by the software GraphPad Prism 9.4. The average and standard errors were calculated for all tested animals based on the normalized values of each hour. The oscillation frequencies of the average population were evaluated based on the average values of each experiment using Fourier analysis-based software LSP with a p<0.01. For individual analysis we used the online platform Discorythm, combining different algorithms including Cosinor, JTK and LSP50. We chose Cosinor it is the one designed to detect efficiently the acrophase.

RNA-seq experimental design

All polyps were isolated in wells of six-well plates. Then, they were subjected to the 12: 12h LD cycle with 17 PPFD (+/−2) light intensity during 72 hours for entrainment in an incubator with a stable temperature at 18°C. Subsequently, the polyps were divided into 2 experimental subgroups: 12: 12h LD and Continuous Dark or DD. Sampling began at 7 am (ZT0) and was performed at 4-hour intervals over 24 hours. At each time point, three or four individual polyps were sampled from each experimental group, immediately snap-frozen in liquid nitrogen, and transferred to − 80°C for storage.

RNA extraction, library preparation, and sequencing

Total RNA was extracted from all sampled polyps (n = 96) using TRIzol reagent (Invitrogen). Purified RNA samples were analyzed using a NanoDrop 1000 spectrophotometer (Thermo Fisher Scientific) to assess RNA quantity and 2200 TapeStation (Agilent) to assess RNA quality (RNA integrity number, >8.5). From each of the 96 samples (4 biological replicates in LD and 3 in DD, with the highest-quality extracts × 4 experimental subgroups × 7 time points, 1.5 μg of RNA was sent for library preparation (INCPM mRNAseq protocol) and sequencing at the Weizmann Institute Sequencing Unit, Israel. The libraries were sequenced using the bulk MARS-seq protocol (Jaitin et al. Science, 2014, Keren-Shaul et al., Nature Protocols,2019) on an Illumina NovaSeq 6000, resulting in an average of 17 million single-end reads of 113 bases per sample.

Bioinformatic analysis

First, the unique molecular identifier (UMI) sequence of each read was extracted and placed within the read1 header file using UMI-tools extract (umi_tools v1.1). Next, the reads were mapped onto the Nematostella genome (NCBI genome GCA_000209225.1) using STAR (v2.6.0a) (Dobin et al., 2013) with default parameters. Mapped reads were then deduplicated based on UMIs using the umi_tools dedup. The mapped reads were sorted by SAMtools (version 1.9). The number of reads per gene were quantified using HTSeq-Count (v0.12.4) (Anders et al., 2015).

Rhythmicity analysis

Rhythmicity in transcript expression was assessed using the RAIN (ref-23) and metac ycle (Wu et al., 2016) packages in R. The RAIN and JTK algorithms from metacycle w ere run separately for each Nematostella genotype in both light conditions (LD and D D), treating them as individual datasets. All replicates (n=3) for each time point within a dataset were analyzed as regular time series to identify transcripts exhibiting daily o scillations. Specifically, we focused on transcripts with a precise 24-hour period, exclu ding those with a range (e.g., 10 to 14 or 20 to 28 hours). To improve the accuracy of identifying true rhythmic genes, only transcripts with a P value < 0.01 in both RAIN an d JTK analyses were deemed confidently cycling transcripts. Genes identified as sign ificant cycling genes were subsequently utilized as input for the DPGP_cluster progra m (McDowell et al., 2018), which clusters genes based on their expression trajectories. Gene clusters comprising 10 or more genes underwent testing for GO term enrichm ent. Heatmaps were generated using the heatmap package (v4.5.5) in R. Venn diagr ams were generated using the web tool Venn diagram (http://bioinformatics.psb.ugent.be/webtools/Venn/) And redraw with Inkscape. Expression plots were generated usi ng GraphPad Prism (V.9.1).

GO term enrichment analysis

After obtaining the differential gene expression results, Gene Ontology (GO) analysis was performed using the R TopGO package (v2.50.0).This analysis aimed to identify significantly enriched biological processes, cellular components, and molecular functions among the differentially expressed genes.The file "nveGenes.vienna130208.GO_annotation_141017.txt" was utilized for GO analysis, and it was obtained from the following source: https://figshare.com/articles/dataset/Nematostella_vectensis_transcriptome_and_gene_models_v2_0/807696. This file contains the set of GO-transcript annotations that served as input for TopGO. The algorithm assigns a significance score to each GO term based on the enrichment p-value and the specificity of the term. In this study, the GO analysis was performed separately for the up-regulated and down-regulated genes in each condition (LD and DD) to identify the specific biological processes and molecular functions that are affected by the NvClk1 mutation.

E-box motif enrichment analysis

Sequences for promoter regions (1000kb upstream ATG) of differentially expressed genes were extracted. We manually identified in the list of motif enrichment all the E- box motifs and Circadian E-box motifs. Boxplot were generated using GraphPad Prism version 9.5.1.

Differential Expression Analysis

Differential expression analysis was performed using R (v4.2.2) Bioconductor package, DESeq2 (v1.38.3)51 . Raw read counts were obtained using I used HTSeq- Count (v0.12.4)52 and then imported into DESeq2 for normalization and statistical analysis. Differentially expressed genes were identified using the Wald test with an adjusted p-value cutoff of 0.05. The analysis was performed on all the time-points pooled of each genotype per light condition. The output of the analysis includes a list of genes with their log2 fold change, p-value, and adjusted p-value. Volcano-plots were generated using GraphPad Prism version 9.5.1.

HCR v.3 in situ hybridization

A custom NvClk (NVE2080, amplifier: B3 and B5) and NvMyhc-st probe set (NVE14552, amplifier: B5) were generated. We used zfHcrt probe set (ZDB-GENE- 040324-1, amplifier: B1 and B3) as a negative control. For HCR on Nematostella juvenile, several alterations were made to a previously described protocol40. Briefly, polyps were plucked and fixed in 4% PFA overnight at 4°C. Polyps were washed 3× in 1× PBS and then dehydrated and permeabilized with 2×5 min washes in 100% methanol. The samples were stored at –20°C overnight. To rehydrate the samples, a series of graded MeOH/PBST washes were used for 5 min each: 75% MeOH:25% 1× PBST, 50% MeOH:50% 1× PBST, 25% MeOH:75% 1× PBST, and finally 2× washes in 100% 1× PBST. To further permeabilize the polyps, samples were incubated in 10 μg/ml Proteinase K diluted in 1× PBST for 10 min. Samples were quickly washed 3× in 1× PBST, and then post-fixed with 4% PFA for 10 min. After post-fixation, samples underwent 3×5 min washes with 1× 2xSSC+0.1% Triton. From now, the following solutions (Pre-hybridization, hybridization and probe wash buffers) were lab-made from the cnidarian-adapted hybridization buffer53. Samples were then pre-hybridized with pre-hybridization buffer at 37°C for 30 min. After pre-hybridization, samples were incubated with 2 pmol of the probe set diluted in hybridization buffer for 16 hr at 37°C. To remove the probe mixture solution, samples were washed 2× for 30 min each with probe wash buffer at 37°C. Samples were washed 2× for 5 min with 5× SSC +0.1% Triton and then treated with probe amplification buffer for 30 min at room temperature. Samples were washed into hairpin amplification buffer containing snap cooled amplifier hairpins and were incubated at room temperature, protected from light, overnight. Samples were then washed with successive 3× 5xSSC +0.1% Triton washes: 2× washes for 15min. Nuclear staining was performed using DAPI 1:1000 in PBST during 1h. Samples were then washed with successive 2× 5xSSC +0.1% Triton washes: 2× washes for 5min. Eventually were slide-mounted in glycerol and stored at 4°C.

Microscopy and Image processing

Samples were imaged using a Zeiss LSM 710 with a x63 oil objective. They were slide- mounted in glycerol. Image manipulation was performed with Fiji54. For the double probes NvClk imaging (Fig1.B), ROIs were generated from each NvClk probes signal and only the ROIs positive for the 2 fluorophores were kept. These ROIs were then used to extract from the original picture the signal considered as true mRNA signal. Figures were then assembled in Inkscape (http://www.inkscape.org/)

Data availability

The RNA-seq data reported in this study have been deposited to the Sequence Read Archive (SRA), under accession PRJNA935092. All data supporting the findings of this study are included in the manuscript and its supplementary files or available upon request from the corresponding authors.

Acknowledgements

We would like to thank the lab of Yehu Moran (The Hebrew University) for their help with advice on applying the CRISPR-Cas9 system in Nematostella vectensis. We thank Ms. Roni Turgeman for her assistance with Nematostella cultures. We thank Dre Julie Batut for hosted RA during the manuscript revision of the manuscript. The research was funded by the Moore Foundation "Unwinding the Circadian Clock in a Sea Anemone" (OL, grant no. 4598), and the Israel Science Foundation (LA, ISF, grant no. 961/19), we also acknowledge German Israeli Foundation GIF Nexus (OL, LA, No. G-1566-413.13/2023). Raphael Aguillon was funded by the Azrieli Foundation. This study represents partial fulfillment of the requirements for a Ph.D. thesis for M. Rinsky at the Faculty of Life Sciences Bar-Ilan University, Israel.

Author Contributions

M.R., R.A., N.S.B., and O.L. designed the research; M.R., R.A., and N.S.B. carried out the experiments; M.R. and R.A. analyzed the data; T.D. performed the bioinformatics;

O.L. and L.A. conceived and supervised the project; M.R., R.A., L.A and O.L. wrote the manuscript. All authors read and approved the final manuscript.

No competing Interests

(a) In situ hybridization HCRv.3 of control probe (zebrafish transcripts) in the WT juvenile with scale bars representing 0.1mm. (b) Nucleotide and translated amino acid sequences of wild-type and NVClk1 alleles. sgRNA target site is boxed in green, PAM site is in bold, a black arrow indicates predicted site of editing. Predicted translation termination site is boxed in black. Insertion is labeled in red. Predicted immature peptide sequences in red. (c) Normalized Movement (a.u), hourly binned over 72h, under continuous dark (Dark - Dark), Blue line represents the NvClk1+/- heterozygote. (d) Individual locomotor amplitude comparison between the three genotype. Kruskal-wallis, multitest comparison *** is p<0.001; * is p<0.05. (e) Lomb-Scargle Periodograms for each corresponding genotypes in constant darkness. The significant period value (p<0.01) is indicated for each genotype in the top left corner. (f) Phase detection (Cosinor) and genotype comparison of 24h rhythmic individuals. See number nrhythmic/ntotal on the x-axis indicating for each genotype the number of 24h rhythmic animals over the total number analyzed. (g) Normalized Movement (a.u), hourly binned over 72h, under continuous dark (Dark - Dark) after 72h of LD 6h: 6h entrainment, Black line represent the WT and red line represents the NvClk1-/-. (h) Lomb-Scargle Periodograms for each corresponding genotypes in constant darkness. The significant period value (p<0.01) is indicated for each genotype in the top left corner. (i). Phase detection (Cosinor) and genotype comparison of 24h rhythmic individuals. See number nrhythmic/ntotal on the x-axis indicating for each genotype the number of 24h rhythmic animals over the total number analyzed.

Venn diagram showing overlaps of rhythmic genes in LD (WT on the left and NvClk1-/- on the right) condition with rhythmic genes from Leach et al., 2019 and the candidate pacemaker genes (based on protein conservation).