Comparative transcriptomic analysis reveals the relationship between cold-acclimation and vernalization

The relationship between cold acclimation and vernalization has been a hanging debate. Cold-acclimation is a cold-responding process that enable plants adapt to cold temperature in a rapid time. While vernalization is that plants ower by quantitative accumulation of cold exposure. 4 ℃ is both the acclimation temperature and vernalization temperature. We conducted transcriptomic analyses during 42 days 4 ℃ vernalization process to explore the whole process of transcriptional dynamics during vernalization and responding situation of short-term cold exposure in Arabidopsis Phosphorylation an important role in long-term temperature drop by time-course analysis. Independency cold-acclimation Cold-acclimation suppress owering after vernalization initiation by continuously blocking photoperiod and gibberellin (GA) pathway. FLM inhibited owering in a FLM-β-upregulation and FLM-δ-independent way. Autonomous and aging pathway positively responded but did little contribution to owering.

regulating CBF expression [31]. When temperature is reaching non-freezing cold, ICE1 (INDUCER OF CBF EXPRESSION1) directly binds to prompter of CBF, activating its expression, then results in a further promotion of COR expression downstream, thus enhancing cold tolerance of plants [32,33,34].
One notion revolving around distinction between cold-acclimation and vernalization is time lag, which is approximately 10 days long [35]. In other words, the initiation of vernalization is 10 days later than that of cold acclimation. Thus, the debate that whether the cold responding pathway and vernalization are independent and what role cold acclimation play after vernalization initiation were born. For the purpose of investigating this mystery, we conducted a transcriptome analysis during the whole cold-acclimation and accompanying vernalization processes.

Transcriptional dynamics of vernalization in Arabidopsis thaliana
Cold-acclimation is a cold-responding process that enable plants adapt to cold temperature in a rapid time. While vernalization is that plants ower by quantitative accumulation of cold exposure. At issue is whether overlapping and interactions exist between two pathways. With the aim of pro ling the whole picture of transcriptional dynamics during vernalization process and responding situation of cold shock in Arabidopsis Col-0 with a functional FRI allele. We set up 8 time points of samples, 4 of which (0 d, 14 d, 28 d, 42 d) are designed based on time sequence of whole cold responding process, while 4 of which (0.5 h, 1 d, 29 d, 30 d) are designed for sake of exploring cold stress genes. Notably, among them 29 d samples experienced one day 22℃ recovery, 30 d samples were exposed 1 d cold following the recovery, which were used to examine the validation of short-term responding gene's expression. BGISEQ-500 was applied to detecting differential expression genes (DEG). There are total of 31744 DEGs (including Noncoding RNA) identi ed. The differential genes of short-term response 0 h VS 1 d, 28 d VS 29 d, 29 d VS 30 d and long-term response 0 hVS14 d, 0VS28 d, 0VS42 d were respectively plotted by Wayne diagram.
There are 693 common genes among both groups ( Figure 1B). KEGG enrichment showed that they mainly enriched in plant hormone signal transduction, especially brassinosteroid (BR). This result suggests BR may play a most signi cant role in hormone responding of stress response, which is wellcharacterized in 00905 pathway map ( Figure 1E) (Table S1). There are 5651 common genes in the longterm combination ( Figure 1C). KEGG enrichment mainly focused on photosynthesis, metabolic and circadian rhythm ( Figure 1D). This suggests that long-term responding changed the whole metabolism process and disrupted the photosynthesis, which is characterized by map 00195 ( Figure 1F) ( Table S2). Enrichment of long-term response and short-term response showed the difference between plant stress response and long-term acclimation. At the initiation stage of cold exposure, only 2709 genes had an expression change within 0.5 hour, which is almost 1/4 of 1 d samples ( Figure 1A). Also, the greatest number of DEGs are detected in 1 d samples than any other points, which means one day is the most drastic responding timing of cold-shock. While genes from 30 d samples exhibited no such signi cant difference though they also experienced one day cold exposure after recovery. Acclimation maybe the explanation since plants had already exposed cold for 28 days. Point 14 d, 28 d and 42 d had similar number of differential expression genes, suggesting that plants maintain high responding level to continuing cold exposure ( Figure 1A).
Time-course analysis was conducted by clustering all genes from different time points in purpose of investigating expression dynamics ( Figure 2). Genes in cluster 3 have a quick response to cold within 0.5 h ( Figure 2C). Gene ontology (GO) indicates that these genes are sensitive to stress environment, which can be well characterized by enrichment of wounding response, defense response and ethylene-activated signaling pathway ( Figure 3C). Genes in Cluster 5 and cluster 7 had similar expression pattern ( Figure 2E, 2G). Both of them responded immediately to changing temperature at point 0.5 h and 29 d, and have a relative longer adaptation period than cluster 3 ( Figure 2E). GO analysis showed that up-regulated genes in these two clusters were enriched in DNA repair and RNA modi cation, the reason of which maybe is that such lasting response can cause DNA impair and replication increase ( Figure 2E, 2G). Genes in cluster 2 and cluster 8 were downregulated during cold exposure and are sensitive to raising temperature ( Figure 2B, 2H). GO results reveals that cluster 2 were enriched in reductive pentose-phosphate cycle, oxidation-reduction process, response to cytokinin and cluster 8 were enriched in cell division and cell cycle. This indicates that life cycle, energy consumption and oxidative activity are active in relative high temperature and decrease in low temperature ( Figure 2B). Genes in cluster1 maintained high expression after one-day cold exposure (Figure 2A). These genes are mainly enriched in phosphorylation-related processes including MKK, CIPK9 and intracellular protein transport (Figure 2A). This suggests that phosphorylation have an important role in long-term temperature change. Genes enriched in cold respond and cold-acclimation were from cluster 4 ( Figure 2D). Compared to cluster 6 referring to regulation of ower development, the former expression level reached the peak at 14 d then dropped smoothly afterwards, the latter counterpart reached the peak at 42 d ( Figure 2F). Also, genes in cluster 4 behaved sensitively to temperature change, but genes in cluster 6 had no activities like that. These results can preliminarily show the connections and differences between short-term and long-term response to cold ( Figure 2F, 2G) (Table S3).

Cold-acclimation and vernalization are independent
Short-term response and long-term response as above-mentioned are well characterized by coldacclimation and vernalization, respectively. These two processes are deemed to initiate with a time difference for about 10 days. To explore the overlapping and interactions between cold-acclimation and vernalization, we focused on genes involved in the relative pathways. Pivotal gene of cold-responding pathway CBF was induced to peak within 0.5 h, and had a smooth and steady expression pattern afterwards ( Figure 3A). Other genes well known in respond to cold also showed a typical uptrend of expression. ICE1, COR, LTI, RAB, KIN, OST1 and BES1 all signi cantly ascended during response process within one day ( Figure 3C). While FLC, central gene of vernalization process, was originally suppressed at around 14 days together with the induction of VIN3, which is considered as rst gene in vernalization pathway. ( Figure 3A). NTL8 was recently found to upregulate VIN3 under long-term cold [36], which expression pattern was similar to VIN3. PRC2 genes behaved a normal rising trend. VAL1 showed an obvious uptrend of expression during vernalization, while VAL2 expression level fall to normal level after a slight increase. Notably, almost all genes in vernalization pathway changed their expression level to regulate FLC expression after 14 d ( Figure 3B). There is an overlapping gene HOS1 between coldacclimation and vernalization ( Figure 3C). The instant rise of FLC and decrease of ICE1 within one day are attributed to it, which promotes FLC transcription through chromatin remodeling and negatively regulates the CBF genes by mediates the ubiquitination and proteasomal degradation of ICE1 [37,38]. This may cause antagonistic responding between two pathways. However, this uptrend of FLC and downtrend of ICE1 did not maintain too long. After one-day cold exposure, expression of these two genes went into a conserve expression trend, which means such interactive activity would not affect the initiation timing of two pathways. This indicates that cold-acclimation plays a predominant role in cold responding process. Consisting with the time-lag notion, rapid responding genes did not participate in quantitative responding process. On an objective note, cold-acclimation and vernalization are independent in overlapping-gene regulatory perspective (Table S4).

Cold acclimation suppress owering after vernalization initiation
The signi cance of vernalization is to allow plants to get through the cold winter and ower in the warm spring. Thus, cold is supposed to recognized as winter signal to inhibiting other owering pathways. The impacts of cold on owering can be divided into two types. One is the direct effect on owering that plants perceive the environmental signal thus exhibiting early owering in warmer temperature and late owering in lower temperature [39]. The other is the indirect effect on owering by cold acclimation process in order to obtain cold tolerance. Genes in cold-acclimation are induced within a short period of time and keep stable in a high expression level afterwards. So, what role does cold-acclimation play after plants initiating vernalization response. It has been reported that mutant shi1 (Col background), which is insensitive to cold, can ower after 95 days of 4℃treatment, but Col was delayed extremely longer than 95 days [40]. Interestingly, our analysis showed that many owering genes upregulated during the cold exposure ( Figure 4). Indirect effect involving cold acclimation is associated with photoperiod and gibberellin pathway. The genes related to circadian rhythm had a quick response to cold ( Figure 4A). Almost all genes involved in circadian clock had an expression increase at 1 d and decrease at 29 d. Among them, CCA1/LHY, core genes of feedback loop, had similar expression pattern [41]. They were induced within one day and had a perfect high-regulating response to cold. PRR7/9, TOC1 as the part of the morning loop and evening loop respectively also were up-regulated [42]. While TCP transcription factor CHE (CCA1 HIKING EXPEDITION) participating in an additional module by which CCA1 and CHE reciprocally repress each other showed a severe down-regulation [41]. Phytochromes PHYA and PHYB were upregulated during cold exposure, PHYB's expression rose a little bit slower than PHYA. Another type of photo-receptor Cryptochrome, with ability of sensing blue light and UV-A also showed a coincident upregulation. BBX proteins containing one or two B-box domains in N-terminal regions, provides a key link between the circadian clock and photoperiod [43]. Almost all members in BBX family had an expression increment. However, the expression level of central gene of photoperiod CO was found no distinct difference after cold treatment, which capable to activate FT expression [44]. This may be because CCA1/LHY contributed to cold tolerance by promoting CBF while facilitated expression of CDF at the same time [45,46]. Besides this, GI showed the same expression level at 42 d as 0 d after ascending the peak at 14 d. Also, FKF1, which can interact with CO to stabilize it [47,48], gradually decreased after a quick induction at 0.5 h. Therefore, all increase shown in phoptoperiotic pathway actually had done nothing to owering. These results suggest that when temperature falling to non-freezing cold, circadian rhythm gives preference to initiate cold acclimation, despite upregulation of certain genes may suppress owering (Table S4).
Another pathway associates to both cold acclimation pathway and owering is gibberellin (GA). GA20ox and GA3ox which play an important role in GA activation merely had a slight increase. While key gene for GA deactivation GA2ox6 was distinctly upregulated. Negative regulator of GA signaling SPY also had an expression increment [49]. Other genes participating in GA signaling barely changed. These results indicate that in order to accumulate DELLA proteins to acquire cold tolerance, GA pathway was also inhibited during vernalization by reducing active GA form. Cold acclimation suppressed owering by continuously block photoperiod and gibberellin pathway after vernalization initiation ( gure 4g) ( in Arabidopsis, which all serve as oral repressor in vernalization [50]. However low temperature did not inhibit their expression ( Figure 4C). As a matter of fact, all of these homologous have an uptrend of expression under 4℃ environments. This maybe because MAF1 also participates in thermosensory pathway. FLM-β as a function form of FLM interact with SVP to suppress SOC1 expression [51]. As temperature rising to 27℃, the ratio of two transcripts FLM-β/δ is reduced, thus decreasing the repression effect of SOC1. Recent research found that MAF2 regulates owering time in the same way [52]. Interestingly, under 4℃ FLM-β was the only signi cantly increasing transcript, and FLM-δ had a relative stable expression pattern ( Figure 4C), which means FLM repress owering in a FLM-βupregulation and FLM-δ-independent way. MAF2/MAF3 also regulate owering by changing expression of sole transcript. Additionally, their interacting protein SVP were also upregulated distinctly to cooperate with them to suppress owering. This indicates that MAFs served in thermosensory pathway to inhibit owering rather than vernalization. However, there does existing upregulated genes that can promote owering. SPLs from age pathway were up-regulated in juvenile period, with the result that SOC1 expression was induced ( Figure 4D) (Table 16). At the same time, expression level of genes in autonomous pathway like FCA, FLD, FLK, LD and FY had been increased during the cold treatment, which partially repressed FLC expression ( Figure 4E) (Table 17). However, these pathways failed to activate the downstream owering genes such as AP1, LFY. The reason of this maybe is that current expression increase of SOC1 and decrease of FLC are not su cient to initiate owering transition. These results suggest that cold signi cantly suppress owering by directly upregulating owering repressor FLM in thermosensory pathway and continuously blocking photoperiod and GA pathway. It strikes the right note that, plants choose to acquire adaptation ability of temperature drop rstly then continue reproduction (Table 4).

Alternative splicing mediation during vernalization
Alternative splicing is a ubiquitous post-transcriptional RNA processing by which multiple transcripts can be generated by a single gene. It has been reported that temperature is closely associate with AS process.
Several kinds of mechanisms of alternative splicing have been reported, including skipped exon (SE) (particular exon is excluded from mature mRNA), mutually exclusive exons (MXE) (choice between two constitutive exons), alternative 3'/5'splicing site (A3SS/A5SS) (distinct 3'or 5'splicing sites are generated in resulting isoforms) and retained intron (RI) [53]. IR is a predominant form of alternative splicing in plants, which main function is to generate transcripts with premature termination codons (PTCs), thus leading to nonsense mRNA decay (NMD) [54]. Total of 1540 deferentially alternative splicing genes (DAS gene) were identi ed, which made up 4.85% of all DEG. In the grand scheme of things, the proportion of IR is decreasing during cold exposure, while MXE came to appearance from 1-day treatment, which is the AS form of FLM. MXE plays a role in generating variety kinds of protein without changing structure of them, by which adjust their functions. A3SS form also showed a signi cant increase ( Figure 5A). Another reported alternative splicing regulation in owering is CCA1, which respond to cold by increasing functional CCA1-α [55]. The mechanism of it is A5SS. These results suggest that plants attempt to improve production of proteins to deal with changing environment cues more e ciently ( Figure 5B).
When view the perspective of particular responding pathways, it was found that AS mechanism enormously affect the whole cold exposure process. Activation of calcium channel is incidental to initiation of cold responding. Cyclic nucleotide-gated channels (CNGCs) are Ca 2+ -permeable channels [56]. Full-length transcript of both CNGC5 and CNGC6 barely worked, CNGC6 was even undetectable (Figure 5b). Line chart shows that CNGC5 responded to cold by ratio variation of NM_203224.1 and NM_125179.5, while three transcripts of CNGC6--NM_127960.3, NM_001335885.1 and NM_001335886.1 all had an expression increase, which means they all had a role in responding process. Phosphorylation functioned through whole cold treatment. CIPK3 (CBL-interacting protein kinase3) is crucial for Ca 2+ signal transduction among different abiotic stresses [57]. Full-length transcript NM_179761.2 was found no effect during responding process. However, NM_179763.4 isoform showed a signi cant expression increase, which may be the functional form of CIPK3. With temperature signal transducting in vivo, upstream gene of cold-acclimation ICE1 was induced [58]. ICE1 worked by changing the ratio of NM_001035699.2 and NM_001035697.2, thus leaving them reaching the same expression level. COR15A and COR27 are activated by CBF. Both transcripts of COR15A--NM_180040.3 and NM_129815.5 upregulated to enhance cold tolerance. While COR27A only functioned in full-length form. Hormones are main types of temperature responders. Transcription factor BES1 plays an important role in BR regulation. Full-length of BES1 was also less contributory than splice variants. It was shown that NM_202135.3 is the main responding transcript and NM_202134.1 may react to temperature rise. EIN4 is one of the ethylene receptors [59]. At the early stage of the cold treatment, transcription level of NM_202489.2 and NM_111329.4 are almost the same. However, after long-term treatment when view against NM_202489.2, NM_111329.4 increased drastically. ABA responding protein RAB8 showed a similar regulating pattern as EIN4. Expression level of NM_115221.4 keep going up during the cold exposure, while NM_180365.3 reach the bottom at 42 d.
Circadian rhythm plays a signi cant role in cold-acclimation. PHYA is a putative signal receptor in this process. Transcripts of PHYA all contributed to responding process since they all showed an uptrend of expression. While, splice variant NM_001335612.1 of PHYB showed a distinct induction rather than main transcript. Blue light receptor CRY2 responded to cold in a transcripts ratio change-dependent way. Full length of evening loop genes PRR7, PRR9 had a similar rising expression trend. While splice variants of PRR9 may participate more in responding process.
Vernalization process initiates as FLC suppression. VRN5 showed a curious down-regulation during vernalization ( Figure 3B). However, the reason of this is that spice variant NM_001338682.1 decreased to ensure full transcript NM_113351.3 working effectively. VAL2 was also found having no expression level change ( Figure 3B), the reason of which revolves around the ratio change of NM_001342124.1 and NM_119353.4 to respond cold. ATX1 and ATX2 which both have positive effect on FLC expression exhibited converse expression pattern. Low temperature causes the reduction of ATX1 immediately at 0.5 h and keeps decreasing afterwards. However, after short-term stress response in 0.5 h, 1 d, 29 d, 30 d, ATX2 always maintained a high expression level in long-term cold exposure in 14 d, 28 d and 42 d ( Figure  3B). It was found that ATX1, ATX2 experienced a different way of alternative splicing mediation. ATX1 had a signi cant ratio change of two transcripts, while AS did not affect full length transcript of ATX2.
We classed these analyzed genes in three types: one is ratio change of different transcripts, the next is functional transcript induction, and last one is synergistic work of all transcripts ( Figure 5B). Among them, ICE1, EIN4, CRY2, RAB8, VRN5, VAL2, ATX1 belong to type one, which had a proportion change of two transcripts. This type of genes may have antagonistic functions between two transcripts. CNGC5, CIPK3, COR27, PRR7, PRR9, ATX2 belong to type two, whose certain transcript showed a distinct expression rise. This type of genes which nonfunctional transcripts may involve in nonsense mRNA decay, thus affecting the potency of functional transcript. CNGC6, BES1, PHYA, CPR15A, FCA belong to type three, which transcripts all had response abilities. This type of genes may generate multiple proteins to respond temperature more e ciently. Therefore, alternative splicing responds to temperature change in a delicate regulatory way.
To further verify the genes transcript pro les obtained from RNA-Seq results, we selected six of DEGs: FLC, VIN3, CBF2, CBF2, MAF1, SOC1 to con rm the validation by qRT-PCR assay. The outcome in each case was consistent with the RNA-Seq assay ( Figure 6) (Table S5).

Association between cold-acclimation and vernalization
The relationship between cold acclimation and vernalization has been a hanging mystery. How come two successive pathways have no interactions under same cold treatment. The heart of the contradiction is the initial response time. Cold-acclimation is a rapid responding process in order to adapt the changing temperature as soon as possible, key requirement of which is quickness. Vernalization is a quantitative process in order to overwinter and ower in a warm spring, key requirement of which is accumulation. This led to the impossibility of overlapping regulation between two pathways and the priority of coldacclimation in responding to cold. Time difference for about 10 days of two pathways' initiation veri ed this hypothesis ( Figure 3A). Besides recent study found that NTL8 being upstream regulator of VIN3 does not function to cold stress genes like COR47 [36]. This further con rmed that vernalization is independent with cold acclimation process. Thus, another question is presented. What role does cold-acclimation play after initiation of vernalization. Photoperiod and gibberellin (GA) as crosstalk between cold acclimation and owering characterized this interplay. Circadian rhythm oscillated instantly in response to cold. CCA1/LHY as core oscillators of circadian clock contributed to both CBF and CDF induction, upregulation of which mains this oscillation did not contribute to owering, merely for the sake of promoting cold acclimation. External application of GA can enhance the activity of the LFY promoter and directly initiate the oral switch [60]. DELLA proteins acting as negative regulators in gibberellin (GA) signal transduction [49], was found accumulating under low temperature, thus contributing to CBF-mediated cold acclimation and late owering [61]. Activation genes of GA did not exhibit obvious uctuation while deactivation genes showed a distinct expression increase, which means GA pathway was also occupied by cold acclimation thus inhibiting the owering. These results suggest that cold acclimation suppress owering after vernalization initiation. low temperature can also directly affect owering. Cold did not suppress homologues of FLC. It was found that they all worked in AS-mediated thermosensory pathway. Functional transcript FLM-β was induced to interact with SVP which also had an uptrend expression level, thus repressing owering ( Figure 7). However, there are some expression change in owering pathways actually contributory to oral transition. Endogenous autonomous and aging pathways did have a positive effect on owering during cold exposure. Autonomous genes all had an expression increment, which partially participate in repressing FLC expression. Added to this, almost all SPL members were upregulated, which can promote SOC1 expression. Nevertheless, it is the lonely voice, not su cient to initiate oral switch. It was reported that Col ower after 200 days cold exposure under long days [40]. This suggest that plants would accumulate these owering promoting genes in a relative slower pace under non-freezing cold. After thoroughly adapting to temperature change, they will ower despite the 22℃ recovery after vernalization process (Figure 7).
We focused more on cold exposure process at different stages rather than 22℃ recovery period of vernalization. It was found that WDR5a and ATX1 both facilitate FLC expression by depositing H3K4me3 in its locus, and both of them were downregulated during the vernalization ( Figure 3B). While recently Xi et al (2020) explored and reported the full dynamic of vernalization. Continuous activity of H3K4me3 was shown during vernalization which is not consistent with this [62]. Interestingly, ATX1 experienced alternative splicing mediation during cold exposure. Downregulation of it is due to expression decrease of a splicing variant to ensure a full function of trimethylated. And this indicates that ATX1 may be the main type of methyltransferase during vernalization.

Cold sensors of cold-acclimation and vernalization
What is the cold sensor in plants remains an open question. To explore this, we set a criterion that sensor should respond faster than CBF (0.5 h) and respond to temperature change in 29 d, 30 d quickly and accurately. Some potential cold sensors are proposed recently [63]. First is CNGCs (cyclic nucleotidegated channels), which are Ca2+-permeable channels mostly localizing at plasma membrane [64]. It was found that CNGC5, CNGC6 deferentially expressed during cold exposure and both were mediated by alternative splicing. However, functional transcript of CNGC5 and CNGC6 both start responding after 1d, which is later than CBF responding time (Figure 5b). Secondly, phytochromes served as photoreceptors in plants regulate the photomorphogenesis [65]. PhyB was recently reported as thermosensor in Arabidopsis [65]. However, PHYB also barely responded before 1d and showed no apparent expression uctuation at 29d, 30d. While PHYA upregulated rapidly to cold at 0.5 h, and have a quick response to 29d, 30d temperature change. In tomato FR and R perceived by phyA and phyB have antagonistic effect to cold tolerance, and phyA plays a predominant role [66,67], the mechanism of which may function in the same way in Arabidopsis. But whether phyA is the cold sensor remains lack of evidence.
Whether sensors of cold-acclimation and vernalization are the same is also in debate. Helliwell et al proposed that the initial response to a temperature change is physical in vernalization, by which chromosomes loops are gradually disrupted during cold treatment. Also, transcription activity would fall into lowest level within one week [68]. However, FLC did have a short-term elevation in response to cold at 1d in a HOS1-dependent way (Figure 2A) [37]. Added to this, it was shown that FLC expression maintained high expression at 10-day cold in another transcriptomic analysis recently [64], which means this postulation is questionable. Additionally overexpression of CBFs in Col background increases the expression level of FLC in a vernalization-independent way [69]. Therefore, FLC may share the same cold sensor for downstream of CBF in response to short-term cold and for the sensor of vernalization initiation. As for cold sensor of vernalization process, recent research reported that growth condition can serve as thermo sensor. Temperature falling leads to slower growth and slowing growth leads to reduced NTL8 dilution, thus inducing VIN3 expression [36]. Therefore, there may be no speci c sensor serving as vernalization inducer.
Alternative splicing mediation during cold exposure Temperature has always been considered a key factor in regulating alternative splicing. Recent research found that plants would occur a rapid and dynamic AS mediation in response to short-term cold [70].
With the aim of exploring AS regulation under long-term cold, we explored the overall trend of change and the regulation of individual genes during 42 days cold treatment. It is found that the AS had been responding through the entire cold exposure. From an overall perspective, proportion of ve mechanisms changed during the cold exposure in order to generate more functional transcripts, thus enhancing the cold tolerance more e ciently. From the perspective of particular genes, some of them which total expression has not changed may be greatly up-regulated in a single transcript. In general, genes in coldacclimation and vernalization respond to cold mainly by three ways: ratio change of different transcripts, functional transcript induction, and synergistic work of all transcripts, which enable plants addressing environmental signal more delicately and e ciently.

Plant materials and growth conditions
A functional FRI locus from the Sf2 line (CS6209) was introgressed into Col to construct FRI-Col by the R. Amasino lab [71], and the FRI-Col seeds were provided by Dr. Yuehui He (Shanghai Center for Plant Stress Biology, Chinese Academy of Sciences. Seeds were surface-sterilized by 75% ethyl alcohol for 1 min, 10% sodium hypochlorite for 15mins, washed 6 times by sterile water, and strati ed at 4°C for 2 days before being sown on 1/2MS medium. Seedlings were grown in growth chambers at 22

RNA extraction and RNA-Seq library construction
Total RNA was extracted using RNAiso Plus (TaKaRa). Total RNA is processed by mRNA enrichment method or rRNA removal method. mRNA enrichment: Enrich the mRNA with polyA tail using magnetic beads with OligodT. Add appropriate amount of interruption reagent to the obtained mRNA to fragment it under high temperature conditions, use the interrupted mRNA as a template to synthesize the rst-strand cDNA, and then con gure the two-strand synthesis reaction system to synthesize the second-strand cDNA, and use the kit to purify the recovered End repair, add a base "A" to the 3'end of the cDNA and connect the linker, then select the size of the fragment, and nally perform PCR ampli cation; the quality of the constructed library is checked and sequenced after passing.
Different expression genes and different alternative splicing genes analysis High-throughput sequencing was performed using the BGISEQ-500 platform. After several data processing steps (including removal of adaptor sequences, null reads, and low-quality reads), pure reads were obtained from the original sequence. After getting clean reads, we use HISAT to align clean reads to the reference genome sequence (GCF_000001735.4_TAIR10.1). Bowtie2 [73] was used to align Clean reads to the reference gene sequence, and then use RSEM [74] calculate the gene expression level of each sample.
rMATS was used to detect differential splicing genes between different samples and splicing events of the samples themselves. rMATS is a software for differential alternative splicing analysis of RNA-Seq data. It uses the rMATS statistical model to quantify the expression of alternative splicing events for different samples, and then calculates P value with a likelihood-ratio test to indicate that the two groups of samples are in IncLevel (Inclusion Level) The difference in level, lncLevel and using Benjamini Hochberg algorithm to correct the p value to obtain the FDR value.

Consent for publication
Not applicable.

Availability of data and materials
The datasets generated during the current study are supposed to submit to the NCBI repository.
(We are sorry that we failed to upload the clean data because of the server connection, we will try to acquire the accession number as soon as possible.) Figure 1 Comparison between short-term and long-term of cold responding (A) Bar graph showing total numbers of differentially upregulated (orange) and downregulated (green) genes of T0h VS 0.5h, 1d, 14d, 28d, 29d, 30d, 42d. (B) Venn diagram showing common genes and unique genes in group T0h Vs T0.5h, T0h Vs T1d, T28d Vs T29d, T29 Vs T30d and group T0h Vs T14d, T0h Vs T28d, T0h Vs T42d. (C-D) KEGG terms of short-term cold responding and long-term cold responding (E-F) KEGG pathway maps of short-term cold responding and long-term cold responding Time-course analysis of dynamic expression change during vernalization (A-H) Clustering was performed using the time course analysis software Mfuzz [72]. Shown on the left were heat maps that expression level of genes in cluster1-8 were normalized to Log2, red and blue colors represent up-and downregulated genes in heat map, respectively. Shown on the right were Gene ontology (GO) terms of each cluster.