Identification of the Molecular Clockwork of the Oyster Crassostrea gigas

Investigations on circadian behaviors and expression of clock genes were performed on 160 Pacific oysters C. gigas (diploid; 70 ± 1 mm shell length; mean ± SE) from oyster farm source ("Port du Rocher", La Teste de Buch, Arcachon bay, France, Lat. 44°38'N, Long. 1°7'O). Experiments were performed in an isolated room equipped with anti-vibrating benches to minimize external influences on animal behavior at the Marine Station of Arcachon (France) from February to April 2013. Oysters were split into two 150-L tanks and maintained in natural ([Chla] = 0.1 ± 0.07 μg·l^-1) and oxygenated seawater of stable composition (T = 17.6 ± 0.1°C; pH = 7.9 ± 0.1; salinity = 34.9 ± 0.2 ‰; mean ± SD). Physical parameters of seawater were monitored with a R301 pH meter (Consort, Belgium) and a Cond 330 I conductivity probe (WTW, Germany). Following 10 days acclimation to lab conditions, oysters were maintained under L:D 10:14 cycle (light phase from ZT 0 to ZT 10 and dark phase from ZT10 to ZT 24) for 15 days. Oysters were sampled during light phases (ZT 1 in days 14 and 15, ZT 5, and ZT 9) and dark phases (ZT 11, ZT 15 and ZT 23 in days 13 and 14) starting on day 13 of L:D exposure (8 sampling times, S1 Fig). Remaining oysters were thereafter exposed to constant darkness for 15 additional days. Additional samplings were performed at circadian times CT 1 (day 14 and 15), CT 5, CT 9, CT 11, CT 15 and CT 23 (day 13 and 14) starting on day 13 of D:D exposure (8 sampling times). Synthetic diagram of experimental timeframe and sampling times was provided in S1 Fig. At each sampling time, gill tissue was individually dissected from 9 oysters under natural light during light phases or under dim red light during dark phases. Tissues were preserved in RNA later (Qiagen) at 4°C overnight and then transferred at -80°C until RNA extraction.

Total RNA was extracted from individual samples using TRI^® Reagent (Invitrogen, Carlsbad, CA, USA). Total RNA quantity and quality were assessed by spectrophotometry (OD260, OD280) and 5 μg total RNA was individually submitted to reverse transcription using oligo dT17 and Moloney murine leukaemia virus (M-MLV) reverse transcriptase (Promega, Madison, WI, USA).

Putative C. gigas clock sequences were identified through a local combination of tBLASTn and BLASTp searches of CDS, EST and genome [15] databases of C. gigas using Pfam conserved domains as well as homology with clock sequences identified in other organisms [6,16]. CgCry1 was previously described and complete sequence was retrieved from Mat et al. [9]. Similarly, ROR and Rev-Erb homologs in oyster were retrieved from Vogeler et al. [17].

Full cDNA of candidates in C. gigas, including UTR, were obtained by RACE using specific primers (S1 Table) and methods described by Scotto Lavino et al [18,19]. Briefly, total RNA was extracted from gill tissue of oysters. Reverse transcription for the determination of 5' cDNA ends were performed with the SuperScript II reverse transcriptase (Invitrogen, USA) and RT primers (S1 Table). Similarly, 3' cDNA ends were amplified by PCR using reverse transcriptase and Qt primer. Sequences were amplified by PCR using Qo / specific Ro and Qo / specific Fo primer combinations. PCR products were thereafter used as template in nested PCR using Qi / specific Ri and Qi / specific Fi primer combinations (S1 Table) according to Scotto Lavino et al [18,19]. PCR products were separated on agarose gel and purified using Wizard^® SV gel and PCR clean up system (Promega, Madison, USA). Purified products were ligated into pGEM-T vector (Promega, MAdison, USA) and used to transform DH5α bacteria (Invitrogen, USA). Bacteria were cultured in Luria-Bertani broth medium containing 100 μg·ml^-1 ampicillin and plasmids containing inserts were extracted and sequenced by extension from both ends using T7 and SP6 universal primers.

Complete cDNA candidates including UTR were mapped on C. gigas genome and presence of canonical E-box motif (CACGTG) was searched within a sequence of two kilobases upstream of each transcription start site for each gene. In presence of E-box motif, the site for each motif was annotated based on position upstream of the start-site according to Reitzel et al. [20].

The cDNA sequence and deduced amino acid sequence of candidates were analyzed and compared using BLAST algorithm (http://blast.ncbi.nlm.nih.gov/Blast.cgi↗) and the Expert Analysis System (http://www.expasy.ch/↗). Rooted phylogenetic trees of the different families of proteins were generated from sequence alignments by Maximum Likelihood method using Mega 6 software. Statistical confidence of inferred phylogenetic relationships were assessed by bootstraps of 1000 replicates.

Real time PCR reactions were performed on individual samples with Brillant III Ultra Fast SYBR Green QPCR Master Mix kit (Stratagene, USA) and a final concentration of 100 nM for each primer according to manufacturer’s instructions. Primer sets of clock candidates were designed from full cDNA sequence, excepted for CgCry1, CgRev-Erb and CgROR [9,17]; stability of the expression of glyceraldehyde-3-phosphate dehydrogenase (GADPH) and elongation factor (EF1) as reference genes were verified and EF1 was used as reference gene (S1 Table). PCR efficiency (E) was assessed for each primer pair by determining the slope of standard curves obtained from serial dilution analysis of cDNA from different experimental samples. Reactions were initiated with activation of the DNA polymerase at 95°C for of 10 min followed by amplification of the target cDNA (40 cycles: denaturation at 95°C for 30 s, annealing at 58°C for 30s and extension at 72°C for 30 s). Reaction specificity was controlled using melting curve step from 95°C to 60°C (decrease of temperature of 0.5°C every 10 s).

The comparative Ct method (2^-ΔΔCt method, [21]) was used to determine transcript levels of candidates. Expression levels were normalized (ΔCt) in each sample using EF1 (AB122066↗) sequence as housekeeping genes and ΔCt values from all samples were subtracted from the relative transcription level of candidates of each sample (ΔΔCt). Results are given as the mean (2^-ΔΔCt) and standard deviation (n = 9).

Valve activity of oysters was continuously monitored during L:D and D:D exposures (i.e. 30 days). Fifteen oysters were equipped with HFNI (High Frequency—Non Invasive) valvometers to assess rhythms of oysters through variations of hourly valve opening duration according to Tran et al. [12,13]. Double-plotted actograms (each line represents 2 days) were produced with Chronos-Fit 1.06 [22]. Activity levels above the average of the day were represented by a black section, while values below the 24-h average were indicated by a white section. Chronobiological analyses were performed using the software Time Series Analysis Serie Cosinor 6.3 (http://www.euroestech.com↗). Data were processed and analyzed to validate a rhythm in oysters [23,24]. Briefly, the quality of the data set was assessed by controlling the absence of randomness using the autocorrelation diagram and the absence of a stationary character by a Partial Autocorrelation Function (PACF) calculation [25]. The periodicities in the recorded data were tested with the spectral method of the Lomb and Scargle periodogram [26]. These methods provide a threshold of probability (p-value = 0.95) defining the limit below which the signal can be considered as "noise". The confidence interval of the period was determined using the method of Halberg [27]. When a period was statistically validated, rhythmicity was modeled with the Cosinor model, which used a cosine function calculated by regression [28,29]. The model for a given period was written as: Y(t) = A.cos(2πt/τ + φ) + M + ϵ(t) where A was the amplitude, φ the acrophase, τ the period, M the mesor and ϵ the relative error [24]. Two tests provided the validation of the model: the elliptic test [29] had to be rejected and the probabilities (p-value) for the null amplitude hypothesis had to be lower than 0.05. For a set of data, many significant periodicities could exist. Identification of secondary periodicities was performed by the reinjection of residues from the previously calculated Cosinor model in a repeated procedure. Validation of significant rhythms were accomplished when the whole procedure, i.e. checks the quality of data, significant period by spectral analysis and the statistical validation of the Cosinor model, were achieved.

Statistical analyses were performed using SigmaStat (Systat Software Inc., California, USA). ANOVA treatments that generated p-values below 0.05 were followed by a Holm-Sidak post-hoc test comparing different conditions. Non-parametric Kruskal-Wallis one-Way ANOVA on rank was performed whenever normality or homoscedasticity of data were not met. Correlations of transcript expressions between clock candidates were analyzed during LD and DD regimes using Spearman. Differences were considered statistically significant at p < 0.05.

Molecular approaches allowed the identification of homologs for most of the known invertebrates and mammals clock genes in the oyster C. gigas (S2 Fig). Some of which were previously characterized such as CgCry1 [9] or identified such as CgRev-erb and CgROR [17]. For instance, Voleger et al. [17] analyzed nuclear receptor genes from the genome of C. gigas and demonstrated the presence of homologs of NR1F and NR1D that clustered with human ror and Rev-erb genes respectively. Full sequences of clock candidates (CgClock, CgPeriod, CgCry2, CgBmal, CgTim) and other cryptochrome/photolyase homologs (Cg6-4photolyase, CgpCry) were obtained by RACE, deposited into the GenBank database (see S2 Table for accession numbers) and deduced amino acid sequences were used for further analyses of protein structure and phylogenetic relationships (see Methods). Additionally, full mRNA sequences including UTR were mapped in the genome of C. gigas [15] to identify E-box motifs. In mammals and insects, Clock:Cycle/Bmal heterodimers regulate transcription of downstream genes by binding to E-box motifs [30,31]. E-box motif was designed based on consensus sequences [32,33]. Analyses demonstrated the presence of E-box motifs in the genome of C. gigas near the transcription initiation sites of CgCry2 (-36 bp), CgCry1 (-274 bp) and CgPeriod (-86 bp) but no identical match were found for Cg6-4photolyase, CgBmal, CgClock and CgpCry. Absence of complete sequences for CgROR and CgRev-erb (retrieved from another study) did not allow the search of motifs and mapping of CgTim 5’UTR on C. gigas genome was unsuccessful.

Clock and Cycle/Bmal are PAS-bHLH transcription factors. C. gigas ortholog of Cycle/Bmal was named CgBMAL based on its clustering with other Bmal sequences (Fig 1). Both CgCLOCK and CgBMALl contains specific features of the family with highly conserved PAS-A, PAS-B, PAC and bHLH domains as well as NLS and NES signals [34,35]. Moreover, phylogenetic analyses demonstrated a good discrimination between Clock and Cycle/Bmal proteins with, globally, a clustering of CgCLOCK with vertebrate clock proteins and a closer distance of CgBMAL with Cycle/Bmal protein from arthropods than vertebrate ones. More specifically, close phylogenetic relationships of CgCLOCK and CgBMAL were observed with their counterpart from the marine worm Platynereis dumerilii (Fig 1).

Complete sequence of CgPeriod (4216 bp) revealed an open reading frame of 3945 bp encoding 1315 amino acids exhibiting specific characteristics of PERIOD proteins such as the two tandemly organized PAS domains [35–37] (S2 Fig). Phylogenetic analysis of PERIOD proteins from vertebrates and invertebrates demonstrated a clustering of CgPERIOD with PERIOD proteins from other marine mollusks and annelids compared to vertebrate and arthropod counterparts (Fig 2).

Sequence belonging to the Timeless/Timeout family was also identified in C. gigas. CgTim (4750 bp) was composed of an open reading frame of 2802 bp corresponding to 934 amino acids. Surprisingly, CgTIM was not phylogenitically related to the annelid Timeless (P. dumerilii) but clustered with Timeless proteins from insects and crustaceans (Fig 3).

Molecular approaches led to the characterization of 3 sequences belonging to the Cryptochrome/Photolyase family, in addition to the previously characterized CgCry1 [9]. All candidates contained DNA photolyase and FAD binding domains. Phylogenetic analysis of Cryptochrome/photolyase proteins from vertebrates and invertebrates generated a specific dichotomy of the cryptochromes with the presence of several groups previously identified by Oliveri et al. [16]: 6–4 photolyase, vertebrate-type Cry (vcry or cry2), insect-type Cry (dcry or cry1), plant-like cry and cry DASH. Within each sub group, oyster sequences were closely related to P. dumerilii homologs when both sequences were represented (Fig 4). Further analysis of CRYPTOCHROME/PHOTOLYASE proteins organization in C. gigas demonstrated the structural homology of Cg6-4PHOTOLYASE and CgCRY1 with their counterpart in Drosophila. Similarly, CgCRY2 was composed of highly conserved repressor domains (RD) involved in the function of transcriptional repression of CLOCK:CYCLE/BMAL heterodimer [34,35] (Fig 5). More specifically, insert of the Fig 5 showed a high amino acid conservation of CgCRY2 in the sequence alignment corresponding to the nuclear localization signal (NLS) within the RD2b region that is necessary for mammalian-type cry nuclear localization and subsequent repression of CLOCK:CYCLE/BMAL transcription complex. CgPCRY differed from other plant cryptochromes and more generally from CRYPTOCHROME/PHOTOLYASE proteins by the presence of an extended N-terminal fragment before the conserved photolyase/FAD domain (Fig 5).

Investigations of transcriptional variations of clock candidates and Cryptochrome/Photolyase family members in C. gigas under L:D (10:14) and D:D regimes led to the observation of different profiles. Expression for most of candidates peaks during the photophase under L:D regime with a maximum expression at ZT 1 for Cg6-4photolyase, CgClock and CgPeriod whereas CgCry1, Cgpcry, CgTim, CgCry2 and CgBmal tended to peak at ZT 5 (Fig 6). It is interesting to note that CgRev-erb peaks twice, during the photophase (ZT 5) and during the scotophase (ZT 23). Exposure of oyster to constant darkness was associated with significant modulation of gene expressions (Fig 6). For instance, Cg6-4photolyase, CgCry1, CgpCry, CgTim CgClock and CgPeriod tended to increase transcription levels during transitional illumination regime corresponding to light on/off and despite the absence of light shifting under D:D regime (CT 23–1 and CT 9–11). On the contrary to these profiles, CgCry2, CgBmal, CgRev-erb and CgROR exhibited opposite trends under D:D regime with a decrease of expression during the times associated to transitional light regime under L:D (Fig 6).

Correlation analyses at individual level revealed the presence of two groups of genes with a positively correlated transcription level during L:D and D:D regimes (Tables 1 and 2). For instance, transcription levels of Cg6-4photolyase, CgCry1, CgpCry, CgTim, CgClock and CgPeriod followed similar trends of variations with ρ ranging from 0.4 to 0.944 (p < 0.001, Table 1). Similarly, expression of CgCry2, CgBmal CgRev-erb and CgROR followed similar modulations with ρ ranging from 0.273 to 0.891 (p < 0.02, Table 2). However, correlation modifications associated to changes of light regimes were observed. During L:D regime, CgCry2 was positively correlated to the group composed of CgPeriod and CgClock (p < 0.02, ρ 0.264–0.37) whereas D:D regime abolished these relationships. Additionally, D:D regime was associated to negative correlations of Cgclock, CgPeriod, CgCry1 and Cg6-4photolyase with CgBmal and CgRev-erb (p < 0.04, Table 2).

Oyster subjected to a L:D (10:14) regime for 15 days exhibited a daily valve activity, with higher opening activity at the end of the photophase and at the beginning of the scotophase (Fig 7). Spectral analysis of the population revealed an oscillation of oyster behavior at a period of 24.0 h. Further analyses at individual level indicated that 93.3% of oysters exhibited rhythmicity with 92.9% oysters presenting circadian behavior of 24 h. Following exposure to constant darkness for 15 additional days, several behavioral modifications were noticed. For instance, circadian rhythmicity was still observed at a population level (period of 20.6 h) but ultradian period of 12.5 h was also detected. At the individual level, the proportion of oysters exhibiting rhythmicity dropped to 80% and within these animals, only 58.3% had a circadian behavior under constant darkness (Fig 7).

Data presented in the manuscript were submitted to GenBank under the references KX371074 - KX371079. Details were provided in the table S2.