Thermopriming triggers splicing memory in Arabidopsis

Thermopriming induces genome-wide differential gene expression and alternative splicing patterns, and establishes a ‘splicing memory’ that helps plants to survive subsequent and otherwise lethal heat stress.


Introduction
Plants adapt to ever-changing environmental conditions by fine-tuning their molecular responses at the genome, epigenome, transcriptome, epitranscriptome, metabolome, and proteome levels to ensure survival and reproductive success (Hatfield and Prueger, 2015;Burgess et al., 2016). Plants have molecular mechanisms to survive mild stress conditions, which constitute a basal stress-tolerance response. The level of basal stress tolerance varies across plant species and genotypes (Hilker et al., 2016). Intriguingly, plants can acquire tolerance to lethal levels of stress through the establishment of 'molecular stress memory' of previous exposure to mild or severe transient stress (Sani et al., 2013). Brief exposure to stress establishes a new cellular state that differs from the state of unexposed, naïve plants (Bruce et al., 2007). Therefore, the responses (and thus tolerance) to future stress differ between primed and nonprimed plants (Mauch-Mani et al., 2017).
Stress memory allows primed plants to respond robustly and quickly to subsequent exposure to such stresses, which aids in their recovery (Prime-A-Plant Group, 2006;Hilker et al., 2016). Priming could in future be applied in the field to make crops more stress resistant and productive (Liu and Charng, 2012), making the identification of factors controlling priming and stress-induced memory in plants an important focus of study. Priming-induced stress memory functions somatically within a generation; whether such a memory is transgenerational remains largely unknown.
Priming and the establishment of stress memory can help plants to survive a variety of abiotic stress conditions, including drought, heat, and salt, as well as biotic stress (Conrath, 2011;Baurle, 2016). The maintenance of acquired thermotolerance is crucial for successful priming and tolerance to subsequent exposure to otherwise lethal temperatures (Ikeda et al., 2011). The acquisition and maintenance of thermotolerance are two separate processes (Jaskiewicz et al., 2011;Lämke et al., 2016a); the maintenance of thermotolerance is referred to as primingmediated heat-stress memory (Baurle, 2016).
Due to global warming, heat stress poses an increasingly serious threat to plant survival and productivity, and therefore to agriculture and hence food security (Bita and Gerats, 2013). At the molecular level, increased temperature causes protein denaturation or misfolding, thereby compromising cellular viability (Morimoto and Santoro, 1998;Mittler et al., 2012). Heat-stress responses are conserved across eukaryotic species including plants, animals, and fungi, and include the activation of heat-shock factors (HSFs) and the subsequent synthesis of heat-shock proteins (HSPs) (Swindell et al., 2007;Liu et al., 2009Liu et al., , 2011Richter et al., 2010;Nishizawa-Yokoi et al., 2011). HSFs function as transcriptional activators of HSPs, and some HSPs positively regulate HSFs (Scharf et al., 2012). HSPs function as molecular chaperones to protect the proteome against heat stress (Zhang et al., 2010;Wu et al., 2013;Weng et al., 2014;Tang et al., 2016). HSPs are grouped into families based on their molecular size, including HSP70, HSP90, and HSP101, as well as families of small HSPs (Larkindale and Vierling, 2008;Haslbeck and Vierling, 2015). Most plant genomes encode more than 20 HSFs with distinct functions, whereas animals and yeast genomes encode just a few.
The Arabidopsis thaliana genome contains more than 20 HSF genes, including eight whose roles in heat stress responses have been confirmed. One of them, HSFA2, contributes to the establishment of heat stress memory in primed plants (Charng et al., 2007;Banti et al., 2010). HSFA2 activates the chaperonelike protein HEAT STRESS ASSOCIATED 32 (HSA32), which may help maintain cellular homeostasis at high temperature (Charng et al., 2006;Lin et al., 2014). Other genes with a role in the somatic heat stress memory include HSP70 and HSP101, which are up-regulated in response to heat shock (Stief et al., 2014). HSA32 is required for the maintenance of heat stress memory, and therefore it is important to investigate how this protein differentially regulates the targets involved in its acquisition and maintenance (Baurle, 2016).
The molecular mechanisms of stress memory remain unclear, but different categories of stress-inducible genes may play different roles in memory. Heat stress-inducible genes are categorized into three groups. The activation of the first group is not sustained after relief of the stress, and these genes are not thought to play a major role in conferring heat stress memory. For example, expression of the heat stress-inducible genes HSP70 and HSP101 is induced by heat stress but subsequently declines after relief of the stress (Zhang et al., 2010;Wu et al., 2013;Lin et al., 2014;Weng et al., 2014). The activation of the second group is sustained after relief of the stress; these genes play a major role in heat stress-mediated memory (Stief et al., 2014). For example, some HSP genes are induced during the heat stress memory phase, including HSP21, HSP22.0, and HSP18.2 (Sedaghatmehr et al., 2016). The activation of the third group is not sustained after stress relief, but these genes are readily activated by heat stress, and their expression is activated at higher levels after the second exposure to stress conditions (Lämke et al., 2016b). Efforts have focused on identifying the genes that establish and maintain heat stress memory based on their enrichment after the relief of the stress. The levels and duration of such induction vary across plant species and genotypes.
The interplay between transcriptional and post-transcriptional regulation of gene expression is mediated at different levels and epigenetic factors are thought to play a role in establishing heat stress memory (Avramova, 2015;Hu et al., 2015;Brzezinka et al., 2016). Although many studies have focused on transcriptional regulation of stress responses, post-transcriptional responses may also play an important role. Alternative splicing (AS) is a post-transcriptional regulatory mechanism involving the processing of precursor mRNA (pre-mRNA) through the joining of exons in different patterns, resulting in the production of different mRNAs and thus different protein isoforms (Wang et al., 2015). Alternative splicing enriches transcriptome and proteome diversity and also regulates gene expression through multiple mechanisms (Chen and Manley, 2009). Alternative splicing is widespread in photosynthetic eukaryotes, with nearly 61% of intron-containing genes in Arabidopsis exhibiting AS isoforms (Reddy et al., 2013;Laloum et al., 2018). Several types of AS events play key roles in gene regulation, representing an important post-transcriptional regulatory mechanism in response to different environmental conditions (Filichkin et al., 2015). The constitutive pre-mRNA splicing of some transcripts is blocked in human and Drosophila cells under heat stress (Bond, 1988;Biamonti and Caceres, 2009;Shalgi et al., 2014).
Most AS in plants involves intron retention (IR) (Iida et al., 2004;Filichkin et al., 2010); therefore, AS could also be used as a regulatory mechanism to monitor functional versus nonfunctional isoforms and to ensure the production of the correct isoforms at the right levels and times. Specifically, under stress conditions, IR could occur transiently to allow the plant to produce high levels of mRNA upon stress relief AlShareef et al., 2017;Ling et al., 2017). Early reports have indicated that RNA splicing is inhibited during heat-shock responses (Chang et al., 2014), for example the splicing factor SRp38 mediates splicing repression in response to heat shock (Shalgi et al., 2014). How heat stress affects constitutive and alternative splicing is currently unclear, and it remains to be determined whether heat stress can alter the patterns of AS and produce certain pre-mRNA isoforms for proper responses to lethal temperatures. In other eukaryotes, including Drosophila and mice, HSF1 and HSF2 produce different AS isoforms and exhibit tissue-specific and temperature-responsive expression patterns (Lal et al., 2015). Different thermosensitive AS responses have been reported in Neurospora (Diernfellner et al., 2005). Heat stress-inducible AS isoforms of HSFA2 in Arabidopsis may regulate heat-stress responses (Sugio et al., 2009). One of the splice isoforms of HSFA2 has been shown to be important for autoregulation of HSFA2 expression (Liu et al., 2013). In Arabidopsis, heat stress-induced AS of the precursor of miR400 resulted in the accumulation of primary transcripts and reduced levels of mature transcripts (Yan et al., 2012). A genome-wide analysis of heat stress-induced AS in Physcomitrella patens indicated that AS plays an important regulatory role in moss stress responses (Chang et al., 2014).
Previous studies have focused on quantitative changes in gene expression in response to heat stress after the establishment of priming-induced heat stress memory (Charng et al., 2007;Brzezinka et al., 2016;Lämke et al., 2016b). Because AS of pre-mRNA enhances transcriptome and proteome diversity in response to a variety of internal and external cues, including heat stress, we reasoned that AS might contribute to the co-ordination and regulation of gene expression in response to high temperature. Heat-responsive genes exhibit altered AS patterns, prompting us to investigate whether AS co-ordinates the regulation of gene expression during heat priming to establish heat stress memory. Therefore, we proposed that AS may contribute to heat stress-induced memory and we tested this hypothesis using our newly developed heat acclimation protocol, in which heat stress-induced priming leads to memory establishment that allows plants to survive a subsequent and otherwise lethal heat shock. In particular, we also investigated the genome-wide differential gene expression and AS patterns that mediate the establishment of priming and memory, and identified a group of genes with sustained activation levels that are implicated in the establishment of heat stress-induced memory in order to establish a link between AS and heat stress-induced memory.

Plant material and growth conditions
Seeds of wild-type Arabidopsis thaliana Col-0 were surface-sterilized with 10% bleach for 10 min and incubated at 4 °C for 2 d. The seeds were plated on half-Murashige and Skoog (MS) medium agar plates supplemented with 1% sucrose and incubated in a growth chamber (Model CU36-L5, Percival Scientific, Perry, IA, USA) under a 16-h photoperiod (white light, ~75 μmol m −2 s −1 ) and 8-h dark conditions at 22 °C for germination and seedling growth.

Heat priming and heat-shock experiment
The heat-priming platform was designed and developed based on previous studies (Larkindale and Vierling, 2008). Twelve-day-old seedlings were divided into two sets. Set 1 was exposed to a four-step temperature-changing environment as follows. (1) Heat acclimation phase: temperature steadily increased from 22 to 45 °C within 6 h and held constant at 45 °C for 1.5 h.
(4) Recovery/Memory phase: normal growth conditions. Seedlings in Set 2 were grown under normal conditions for an additional 4 d (i.e. to 16-dold), followed by exposure to the conditions described in steps (3) and step (4). Plant samples were collected at different time-points as specified below.

RNA extraction and RNA-seq
Total RNA was isolated from the plants using an Oligotex mRNA Midi Kit (70042, Qiagen Inc., Valencia, CA, USA). The RNA-seq libraries were constructed using an Illumina Whole Transcriptome Analysis Kit following the standard protocol (Illumina, HiSeq system) and sequenced on the HiSeq platform to generate high-quality paired-end reads.

RNA-seq data analysis and gene functional classification
Two replicate samples per time-point were used, except for time-point 3, where only one sample was used for analysis. The annotated Arabidopsis gene models were downloaded from TAIR10 (https://www.arabidopsis.org/). TopHat (Version 2.0.10) was used for alignment and to predict splice junctions (Trapnell et al., 2009). Gene expression levels (fragments per kilobase of transcript per million mapped reads, FPKM) were calculated using Cufflinks (Version 2.0.0). The differentially expressed genes (DEGs) were identified using Cufflinks and the limma package in R (www.r-project.org). Very strict criteria were used to define DEGs: they had to simultaneously show more than 1.8-fold up-/down-regulation in both replicates, and P-values calculated by limma had to be less than 0.05. To filter out false-positive junctions, strict criteria (i.e. an overhang size of more than 20 bp and at least two reads spanning the junctions) were set as cut-off values . JuncBASE was used to annotate all AS events based on the input genome co-ordinates of all annotated exons and all confidently identified splice junctions (Brooks et al., 2011). Fisher's Exact Tests were used to identify differential representation of each type of AS event. For all types of AS events, those with a P-value <0.001 were identified as significantly different. For IR, Fisher's Exact Tests were performed on the intron-read counts and corresponding exon-read counts between different samples. In addition, IR events uniquely identified in any sample were considered significant only if there was at least five-fold coverage of support and the P-values of these events were assigned to zero. For alternative 5´-splice sites (SSs), 3´-SSs, and exon-skipping events, Fisher's Exact Tests were performed on the comparisons of junction-read counts and the corresponding exon-read counts between samples. Gene ontology (GO) classifications were performed with the DAVID software (https://david.ncifcrf.gov/, (Dennis et al., 2003). GO network analysis was performed with Exploratory Gene Association Networks (EGAN, http://akt.ucsf.edu/EGAN/).

RT-PCR and RT-qPCR
Reverse-transcription PCR (RT-PCR) and reverse-transcription quantitative PCR (RT-qPCR) were performed as previously described . Briefly, for RT-qPCR, digestion of contaminating DNA in total RNA samples was performed after RNA extraction using a RNase-Free DNase Set (Invitrogen, cat. no. 18068-015) following the manufacturer's protocol. The total RNA was reverse-transcribed using a SuperScript First-Strand Synthesis System for RT-qPCR (Invitrogen) to generate cDNA. The qPCR was performed using Power SYBR Green PCR Master Mix (Invitrogen, cat. no. 4367659) under the following conditions: 95 °C for 10 min, followed by 40 cycles of 95 °C for 15 s, 60 °C for 1 min.

Establishment of a heat stress-induced priming platform in Arabidopsis
To investigate the transcriptional and post-transcriptional regulation of heat-stress priming and memory, we developed a comprehensive heat stress-induced priming platform in Arabidopsis (see Fig. 1A). Briefly, we grew Arabidopsis seedlings for 12 d and collected some plants at the 4-leaf stage before exposure to heat stress (time-point 1, TP1). Subsequently, we initiated priming by exposing one set of plants (Set 1) to a gradual increase in temperature from 22 to 45 °C over the course of 6 h. We collected one sample after 3 h when the temperature reached 33.5 °C (TP2) and one sample after 6 h when the temperature reached 45 °C (TP3). These samples were used to study early responses to mild and severe heat stress, respectively. The remaining plants were incubated at 45 °C for 1.5 h, at which point another sample was taken (TP4). The remaining plants were transferred to 22 °C and a further two samples were collected, one after 2 d (TP5) and one after 4 d (TP6), prior to second exposure to high temperature (45 °C). The first exposure to heat was intended to induce priming, whereas samples TP5 and TP6 were used to investigate the recovery phase and the maintenance of somatic memory.
Subsequently, Set I (primed) and Set II (non-primed) were used for the second heat shock to test the effects of priming on the acquisition and maintenance of thermotolerance. We collected one sample from set II before exposure to heat shock (TP9) to serve as a control for TP6 plants, in which the expression of memory genes in response to heat stress would be maintained. The two sets were exposed to high temperature for 1.5 h and two samples were collected from Set I (TP7) and Set II (TP10). The plants were then maintained at 22 °C for 2 d, after which we collected samples from Set I (TP8) and Set II (TP11). These two samples were used to investigate the differential responses of primed versus non-primed plants. Plants were photographed and analysed phenotypically 6-10 d after treatment. Set I plants (primed for stress tolerance by an initial and gradual exposure to heat stress) performed much better than Set II plants (non-primed plants; Fig. 1B); this experiment was repeated more than three times with similar results. This priming approach involved different phases, including a priming phase, a recovery from priming phase, a second exposure phase, and a recovery phase after exposure to high temperature.
The second exposure to 45 °C constituted severe heat shock in Arabidopsis. However, primed plants were well adapted to 45 °C after the establishment of stress memory by the priming regime. Thus, our priming platform allowed the acquisition and maintenance of thermotolerance to be investigated at the genome scale, including the transcriptional and post-transcriptional levels of genes that were differentially regulated in their expression and/or splicing patterns to establish and maintain heat stress memory and thermotolerance in primed plants.

Differential gene expression during different phases of heat-stress priming
After establishing the priming method, we investigated the genome-wide transcriptional changes induced by priming and their role in the establishment and maintenance of heat stress memory and priming. As described above, we collected samples at different time-points before, during, and after priming, recovery, and a second exposure to high temperature, and we isolated RNA and performed RNA-sequencing (RNA-seq) to investigate the role of transcriptional changes in inducing priming and maintaining heat stress memory and subsequently conferring thermotolerance.
Our comprehensive analysis of transcriptomes allowed us to identify differential gene expression patterns at all the different time-points. For example, to identify genes that function in the early response to high temperature, we compared the control TP1 plants with TP2 plants. To identify heat stress-response genes that function during severe conditions, we compared TP1 versus TP3 and TP4 plants. To identify all heat-responsive genes, we compared TP1 plants with TP2, TP3, and TP4 plants.
The list of genes that were responsive to heat stress in our study was compared with publicly available microarray data for heat and abiotic stress-responsive genes using Genevestigator (https://genevestigator.com/gv/). The data were consistent with prior results, as shown in the heat map in Supplementary  Fig. S1 at JXB online. Furthermore, GO analysis showed that these differentially expressed genes (DEGs) are involved in responses to heat and other abiotic stresses.

Identification of differentially regulated genes potentially controlling heat-stress memory
Upon exposure to high temperature, thousands of genes were up-or down-regulated, as observed by comparisons between TP1 and TP2, and between TP3 and TP4. However, during the recovery phase after heat stress, the DEGs that maintained the same enrichment patterns were likely to function in maintaining heat stress memory. To identify such genes, we compared TP6 plants (4 d after recovery from heat stress) to TP9 (plants not primed and not exposed to heat stress). DEGs in TP6 plants should include genes whose expression is up-or down-regulated during heat stress and maintained after the recovery phase. Such regulation should be unique to primed plants and not be present in non-primed plants. Therefore, we identified genes whose expression changed by at least 1.8-fold after recovery from heat stress in priming.
Approximately 2600 genes were differentially regulated at TP6 compared to TP9, including 982 up-regulated genes and 1642 down-regulated genes. This list of genes should include those that control the acquisition and/or maintenance of heat stress memory, thereby allowing primed plants to survive subsequent exposure to high temperature. GO analysis indicated that the up-regulated genes are involved in biological processes including responses to heat, high-light stress, water deprivation, and other abiotic/biotic stresses, whereas the down-regulated genes are involved in photosynthesis, chloroplast organization, chlorophyll biosynthetic process, and light stimulus ( Fig. 2A, B).
It is possible that some DEGs between TP6 and TP9 may not be related to the establishment of heat stress memory. To address this, we determined the overlap between heat-responsive genes identified during the priming phase (TP1 versus TP2, TP3, and TP4) and the putative heat stress memory genes identified above by comparing TP6 with TP9 plants. Using the heatmap.2 function in R, we found that TP6, TP1, and TP9 DEGs clustered together, as did TP3 and TP4 (Fig. 2C, D). To determine whether such DEGs in TP6 are associated with the heat stress response, we used Genevestigator to compare our data with the publicly available microarray data. This analysis revealed that DEGs during the recovery phase were associated with the heat stress response (see Supplementary Fig. S2). We also performed GO analysis of these DEGs, which potentially control or help maintain heat stress-induced memory (Fig. 2E, F). Our data showed that 'response to heat', 'response to hydrogen peroxide', and 'response to high light intensity' were the most enriched terms in up-regulated genes, whereas 'chloroplast organization', 'chloroplast RNA processing', and 'chloroplast RNA modification' were the most enriched terms in down-regulated genes. Finally, we used EGAN to construct gene networks of heat-responsive genes ( Supplementary Fig.  S2), and the results showed that a group of memory genes are heat and abiotic stress-responsive. Overall, these results showed that some important heat-responsive genes maintained their expression patterns during the 4-d post-priming phase, which would help the primed plants sustain subsequent heat shock exposure.
Heat stress-primed versus non-primed plants respond differently to a second exposure to high temperature Next, we investigated why primed plants survived the second exposure to heat stress better than non-primed plants. We initially clustered the gene expression profiles using two methods: principal component analysis (PCA) and distance matrix (Fig. 3A, B). The TP1, TP2, TP5, TP6, TP8, and TP9 samples clustered together, which was expected since these plants either were not exposed to heat stress, were exposed to only moderate heat stress, or had recovered from heat stress. Intriguingly, TP11 plants (non-primed plants exposed to heat stress) clustered with TP3, TP4, TP7, and TP10 plants, suggesting that the gene expression profile of TP11 is more like that of plants exposed to heat stress, even after a 2-d of recovery. This result suggests that although heat stress was relieved in the latter samples, the gene expression profiles of these plants mimicked those of plants under heat stress conditions. Equally intriguing, and in contrast to TP11 plants, the gene expression profile of TP8 plants was similar to those of non-primed control plants, suggesting that in primed plants gene expression patterns were reset back to normal once the stress was relieved, but those of non-primed plants were not reset. This would explain why the primed plants survived exposure to high temperature while the non-primed plants did not. There were thousands of DEGs found between TP8 and TP11. Furthermore, we compared these DEGs to known heat-responsive genes and found that they substantially overlapped (~50%), further confirming the high similarity of gene expression profiles between TP11 and plants under severe heat stress conditions. We generated heatmap profiles and they supported this conclusion (see Supplementary Fig. S3).

Transcriptional regulation of a subset of genes may mediate heat stress-induced priming memory
Our RNA-seq data analysis revealed that a subset of transcriptionally regulated genes may mediate the acquisition and maintenance of heat stress memory. We validated the expression patterns of a subset of these genes using qRT-PCR. First, we validated the differential expression of a number of HSPs (HSP101, HSP22, HSP18.5, HSP21, HSP70B, and HSP18.2) and HSFs (HSF32, HSFA2, and HSFB2b) at all time-points. Our qRT-PCR data corroborated the RNA-seq data (Fig. 4A, B). Next, we validated expression of several other heat-responsive genes. The qPCR results demonstrated that many heat-responsive genes, including those that are involved in non-coding RNA-dependent pathways, were differentially regulated. The TRANS-ACTING SIRNA PRECURSORS 1a (TAS1a), TAS1b, and TAS2 were down-regulated by heat stress, whereas their target genes, HEAT-INDUCED TAS1 TARGET1 (HTT1), HTT2, and HTT3 were induced during heat stress. Moreover, targets of miR156 (SPL2 and SPL11) were down-regulated by heat stress (Fig. 4C-E).

Alternative splicing regulates plant responses to heat stress and contributes to priming
Several reports have implicated AS in regulating plant responses to biotic and abiotic stress, prompting us to investigate whether it is involved in the acquisition and/or maintenance of heat stress memory. To analyse changes in pre-mRNA splicing during the priming or recovery phase, or upon exposure to the second heat stress phase, we employed our recently developed pipeline  and performed three steps: prediction of splice junctions, filtering of false-positive junctions, and annotation of AS events. RNA-seq reads were mapped onto the Arabidopsis genome (TAIR 10) to predict the splice junctions using the TopHat software, which identifies the splice junctions of exonic and intronic sequences (Trapnell et al., 2009). Our alignment generated 148 701 splice junctions from the 21 libraries at the 11 time-points for the Set I and Set II samples after removing the false-positives containing short overhangs and low coverage, as previously described . When we compared these junctions with gene annotations from TAIR 10, we found that 70% of were previously annotated and 30% were novel junctions. These data sets were used for comparisons with annotated genes to identify all AS events at all time-points from Set I and Set II, including 127 555 IR events, 6574 alternative 5′-SS events, 19 577 alternative 3′-SS events, 19 mutually exclusive exon events, and 956 co-ordinate cassette exon events (see Supplementary Table S1). Of the IR events, 38 131 had at least five reads supporting the event and >80% of them were found at all time-points. Furthermore, our data indicated that 53% of the intron-containing genes exhibited AS under normal conditions (TP1), 46% under mild stress (TP2), 46% under severe stress (TP4), 55% under lethal stress (TP10), 50% at the recovery phase (TP6), and 52% upon a second exposure to lethal stress (TP7, Supplementary  ). The distance matrix shows that time-points TP1, TP2, TP5, TP6, TP8, and TP9 grouped together, and that TP3, TP4, TP7, TP10, and TP11 grouped together. See Fig. 1 for the temperature regime.
We analysed the different AS events between control plants (TP1) and plants in the priming phase subjected to gradual mild (TP2) or severe heat stress (TP3 and TP4). We observed a significant increase in IR during the priming phase in plants subjected to severe heat stress in TP1 versus TP3, and in TP1 versus TP4 comparisons (Fig. 5A, B). However, during the recovery phase, the number of IR events and genes involved in these events significantly decreased, indicating that heat stress significantly repressed the splicing machinery, leading to the accumulation of pre-mRNAs with retained introns (Fig. 5A, B). GO analysis of the genes that underwent IR in response to heat revealed that this set of genes was enriched in biological process categories involving abiotic stress, including heat stress, RNA splicing, protein folding, and transport (Fig. 5C). These biological processes are modulated during heat episodes to help plants survive the transient stress. In addition, comparisons between primed and non-primed plants revealed that the number of genes with differentially retained introns increased in plants exposed to lethal heat stress. Surprisingly, this number increased more sharply in the recovery phase 2 d after exposure to lethal heat stress (Fig. 5D).

IR contributes to the establishment of heat-stress memory
Next, we investigated whether AS patterns were sustained and remained in plants during the recovery phase (memory establishment phase) and thereby whether they could mediate heat stressinduced memory. We compared the differential intron retention (DIR) events (which significantly increased under heat stress) of TP6 versus TP9 plants. Primed TP6 plants still maintained some DIR events. GO analysis of this group of genes indicated that they were associated with abiotic stresses, response to light stimulus, as well as RNA splicing (Fig. 6A). We then investigated whether the DIR events in TP6 were present in TP2, TP3, and TP4 in order to corroborate our hypothesis that such events are essential for the plant response to heat stress. Interestingly, our data revealed DIR events from a group of 36 genes (Fig. 6B, C, Supplementary Table S3).
To understand the effect of heat shock-and priming-induced changes on IR, we analysed the functional categories and cellular processes regulated by such genes. We identified more than 6000 genes with DIR events during heat shock (TP10) and 6242 genes during the priming phase (TP4), and 6192 genes during the second exposure phase (TP7), when compared to control plants (TP1). A subset of 1180 genes with DIR events was unique to the priming phase; some of them might be responsible for the establishment and maintenance of heat stress memory. Functional analysis of these genes using the DAVID software indicated that they are involved in biological processes including protein transport processes, ribosome biogenesis, rRNA processing, and regulation of translational initiation (Fig. 6E). Interestingly, the functional category 'response to heat stress' was markedly enriched among the DIR genes, and such genes were also enriched in the categories 'heat stress responses' and 'memory establishment', primarily after exposure to lethal temperatures (Fig. 6F).

Differentially expressed and alternatively spliced genes are co-regulated in response to heat-stress priming
Next, we investigated whether the DEGs during different phases of heat-stress priming (acquisition) and maintenance of heat-stress memory were co-regulated. Our RNA-seq analysis revealed that 4483, 11 105, 10 963, 6264, and 4842 genes were differentially regulated at TP2, TP3, TP4, TP5, and TP6, respectively, compared with the control (TP1), and 3233 and 2624 genes were differentially regulated at TP5 and TP6, respectively, compared with TP9 ( Fig. 7A) (fold change >1.8 and P<0.05). Intriguingly, when we compared genes with DIR events with differentially regulated genes, we found significant overlap between DEGs and DIR genes during different phases. The number of overlapping genes increased when plants suffered from increasing temperature (from TP2 to TP3 and TP4), with more than 65% of DIR genes being DEGs. The number of overlapping genes then decreased during the recovery phase. Interestingly, the overlapping genes between primed and non-primed plants during the recovery phase after lethal heat shock sharply increased (TP8 vs TP11), suggesting that these two processes are, to a certain extent, co-regulated to respond to heat stress (Fig. 7A). It should be noted that differential expression of many genes may also be regulated by AS.
Nonetheless, many DIR genes were also DEGs that were co-regulated during the priming phase. GO analysis revealed that these genes, especially those that underwent IR in response to heat, are involved in abiotic stress responses, including heat stress, protein folding and transport, as well as regulation of mRNA processing (Fig. 7B, C). Interestingly, during the recovery and maintenance of the memory phase, we found 1264 genes with DIR events and 246 of these genes were differentially expressed, indicating that the transcriptional and post-transcriptional processes were co-regulated. To investigate whether these subsets of differential AS and DEGs were involved in heat-stress memory acquisition and maintenance, and thereby in facilitating plant survival upon subsequent exposure to heat stress, we examined their expression/splicing patterns during priming. Our data revealed that 112 of 246 genes are either up-or down-regulated at the transcriptional level and also showed DIR during priming (TP4) (see Supplementary Fig. S4).

Heat-stress priming establishes splicing memory
Because our AS analysis revealed that responses to heat stress were regulated post-transcriptionally via effects on pre-mRNA splicing and that some AS patterns were associated with heatstress memory, we asked whether AS contributes to the establishment of heat-stress memory and thereby helps plants survive their next exposure to stress. We compared the AS events in heat-stress primed (TP8) versus non-primed plants (TP11).
Interestingly, the primed plants (TP8) were capable of efficient splicing and produced splicing patterns similar to those of control plants not exposed to high temperature. By contrast, the non-primed plants (TP11) exhibited high levels of IR and were therefore still producing splicing variants mimicking heat-stress conditions (see Supplementary Fig. S5). Thus, the primed plants appeared to maintain splicing memory, and after relief from the second exposure to stress they 'remembered' to undergo efficient splicing, to produce splicing patterns similar to those of control plants under non-stressful conditions in order to support growth and development.
To validate our RNA-seq data, we performed RT-PCR for several key genes to determine their splicing patterns in primed and non-primed plants. Our RT-PCR data corroborated our RNA-seq data, showing that non-primed plants were not capable of efficient splicing to produce patterns similar to those of control plants under non-stressful conditions, as evidenced by the higher levels of IR. By contrast, the primed plants were capable of producing similar splicing patterns to those of the control (see Supplementary Fig. S6).
Heat priming differentially affects pre-mRNA splicing of HSFs. Members of HSF class B, namely HSFB1 and HSFB2b, repress heat-inducible genes, including HSFs under non-heat conditions and in the attenuating period (Ikeda et al., 2011), whereas heat-induced transcript factors, class A HSFs, including HSFA2 and HSFA7a, play an important role in thermotolerance (Larkindale and Vierling, 2008;Meiri and Breiman, 2009). Interestingly, in our RT-PCR experiments HSFB1, HSFB2a, and HSFB2b showed a dramatic increase in IR in both the heat-priming and heat-shock phases, whereas splicing of HSFA2, HSFA7a, and HSFA7b was much less affected in the priming phase (Fig. 8A, Supplementary Fig. S7). Subsequently, although expression of Class A and B HSFs was up-regulated by heat priming/heat shock, plants had increased levels of class A but few functional class B transcripts in the priming phase, whereas the transcript levels of both classes of proteins were clearly affected by heat shock.
Heat priming induces specific isoforms of serine-/arginine-rich (SR) genes, which are the key regulators of AS, and probably modulates pre-mRNA splicing of HSPs. From our previous qPCR results (Fig. 4A, B) and RNA-seq data, we found that most HSPs and HSFs were highly induced by heat priming (TP3 and TP4) and heat shock (TP7 and TP10). We then tested whether the splicing patterns of HSP and HSF pre-mRNAs were altered during the heat-priming/heatshock program. We found that pre-mRNA from a subset of HSP genes, including HSP21, HSP101, HSP70.10, HSP70.6, HSP90.5, and HSP100.3, were alternatively spliced, mostly due to IR, in response to heat stress/priming. Interestingly, the levels of intron-retained isoforms were much higher in the heat-shock compared to the heat-priming phase, except for HSP70.17. Moreover, HSP70.10 expressed mostly the constitutively spliced isoform in the heat-priming phase, but produced two or more isoforms during other phases. By contrast, HSP90.6 produced more isoforms in the priming phase, compared to the other phases (Fig. 8B).
In contrast to the HSPs, our RT-PCR analysis shows that pre-mRNAs encoding SR proteins tended to be alternatively spliced in the heat-priming process. In addition to retained intron splicing of transcripts, SR30, SR45a, SR34, and RS41 underwent other AS events, including alternative first-exon and alternative last-exon events, which may have resulted in addition/deletion of a functional domain of proteins. However, those genes were affected less in splicing during heat shock when compared with the control conditions (Fig. 8C, D). Thus, it remains to be confirmed whether AS of SR pre-mRNAs agrees with what we found in HSPs genes, namely that those splicing processes were less affected during heat priming. Our RT-PCR data on splicing corroborated our RNA-seq data, and these results suggested that heat stress leads to splicing repression and significant levels of IR. Furthermore, heat priming and heat shock induced different levels of IR. During the memory establishment phase (TP5 and TP6), splicing returned to control-like levels. The second exposure to heat stress in primed plants (TP7) and the first exposure to heat stress in non-primed plants led to significant levels of IR, indicating splicing repression. However, after relief from the stress, the primed plants were capable of efficient splicing, unlike nonprimed plants, which exhibited splicing repression mimicking stress conditions. These findings suggested that it might be possible that heat-stress priming could be established at the post-transcriptional level through the maintenance of splicing memory. Efficient splicing after heat-stress relief is crucial for plant survival and the completion of their life cycle.

Discussion
Global climate change is having a substantial impact on agriculture worldwide, resulting in losses in crop productivity, and it poses a serious threat to global food security. Plants live in a dynamic environment with fluctuating temperature regimes, including both predictable and unpredictable seasonal changes. Prior exposure to mild heat stress primes plants, allowing them to acquire and maintain heat-stress memory, leading to increased thermotolerance. Therefore, priming can be viewed as a form of somatic stress memory, which is crucial for plant survival under subsequent stress conditions. Understanding the molecular mechanisms underpinning the acquisition and maintenance of heat stress is of paramount importance to both basic and applied research. Specifically, which factors control the acquisition and maintenance of heat-stress memory? How is this memory established, and how long does such somatic memory persist?
Priming conditions induce transcriptional changes that are sustained during the recovery or the memory-establishment phase (Ding et al., 2012). Several key components function in heat stress-induced priming; however, no previous study has attempted to identify these factors at a genome-wide scale (González-Schain et al., 2016). In this study, we established a comprehensive heat stress-induced priming platform that allowed us to investigate transcriptional and post-transcriptional regulation in response to priming and subsequent exposure to high temperature. Our platform included a priming phase, a stress-recovery and memory-establishment phase, and a stress-challenge phase. During the priming phase, we subjected plants to a gradual increase in temperature and collected samples after the temperature reached 33.5 °C. This sample (TP2) was used to identify the early-response genes and to determine whether the differential regulation of these genes was sustained during severe heat stress (TP3 and TP4). Our results indicated that genes responsive to heat and other abiotic stresses, including drought and cold, and hormones such as abscisic acid, were differentially regulated. Under severe heat stress, the transcriptional patterns involved up-regulation of HSF and HSP genes, as previously reported (Charng et al., 2007;Mittler et al., 2012).
Next, we analysed the priming and recovery phases in order to identify genes that underpin the acquisition and maintenance of heat-stress memory. Interestingly, clustering analysis of gene expression profiles showed that genes expressed during the recovery phase clustered with control samples (TP1 and TP9) compared to plants subjected to severe stress during the priming phase (TP3, TP4). These results indicated that in the recovery phase most gene expression profiles were reset to the expression patterns observed under normal conditions; some DEGs retained their transcriptional level in order to increase plant tolerance to subsequent exposure to heat stress; these DEGs constitute the memory genes.
Plants acquire thermotolerance through the activation of HSPs, the central proteins that confer tolerance to high Fig. 8. Differential splicing pattern of HSF genes, HSP genes and SR genes at different time-points. (A-C) RT-PCR was performed using primers that flank intron(s). (A) Splicing patterns of heat-shock factor (HSF) genes. HSFA2, HSFA7a, and HSFA7b expressed fewer intron-retained isoforms after heat priming (TP3 and TP4) when compared to after heat shock (TP7 and TP10). HSFB1, HSFB2a, and HSFB2b accumulated more intron-retained isoforms after heat priming when compared to after heat shock. (B) Splicing patterns of heat-shock protein (HSP) genes. genes show different ratios of intron retention after heat priming and heat shock. Most genes, except for HSP90.6 and HSP70.17, expressed fewer intron-retained isoforms after heat priming (TP3 and TP4) when compared to after heat shock (TP7 and TP10). (C) Splicing patterns of a group of serine-/arginine-rich (SR) genes. Specific of isoforms of SR genes were induced in heat priming (TP3 and TP4). (D) Gene structure and intron retention of interesting regions from the nine SR genes in (C) as illustrated using the IGV program (http://software.broadinstitute.org/software/igv/). Coverage of RT-PCR-amplified fragments are indicated by the shaded bars at the top of each panel. (This figure is available in colour at JXB online.) temperature (Larkindale and Vierling, 2008). We identified genes encoding HSFs and HSPs among the genes controlling heat-stress memory. For example, HSP21 was induced and its expression maintained during the recovery phase, indicating that it plays a role in heat stress-induced memory. The accumulation of HSP21 was recently shown to increase the thermomemory capacity of Arabidopsis, thereby increasing thermotolerance (Sedaghatmehr et al., 2016). Likewise, HSFA2 is required for the maintenance of heat stress-induced memory, and we found that HSFA2 was still induced at the end of the recovery phase (TP6) (Lämke et al., 2016a(Lämke et al., , 2016b. Furthermore, the small Arabidopsis HSP gene At-HSP17.6A is regulated by heat, osmotic, and drought stress, and plants overexpressing this gene exhibit enhanced stress tolerance (Sun et al., 2001). Moreover, the plant-specific heat-response protein gene HSA32 was up-regulated and maintained in primed plants (TP6), indicating its involvement in the maintenance of heat-stress memory. This result corroborates the finding that HSA32 is required for the maintenance of acquired thermotolerance after a long recovery following acclimation treatment (Charng et al., 2006). Heat-stress memory is established during the first exposure to high temperature.
Epigenetic factors help control priming responses and help maintain priming memory (Berry and Dean, 2015), and epigenetic regulation would ultimately affect the transcriptional patterns of genes. Therefore, the expression patterns identified in this study could have involved epigenetic regulation. Such patterns must be maintained to help confer thermotolerance to heat stress. For example, sustained H3 Lysine 4 tri-methylation and demethylation persist, even after active transcription of these loci has dropped to normal levels (Lämke et al., 2016a). Therefore, it is possible that such loci are epigenetically marked for hyper-induction upon recurring exposure to heat stress. Indeed, when we compared primed and non-primed plants subjected to a second exposure to heat stress, we identified genes that were hyper-induced in the primed plants. However, these genes did not sustain high transcript levels during the recovery phase of memory establishment and maintenance after the first exposure to heat stress. Therefore, such genes could be epigenetically marked for hyper-induction under recurring heat stress.
Alternative splicing is a key mechanism that eukaryotes use to produce diverse transcriptomes and proteomes to support cellular functions, including during abiotic stress responses (Chen and Manley, 2009;Reddy et al., 2013). The role of AS in establishing priming memory, recovery, and resetting is largely unexplored. AS patterns ensure differential regulation of protein abundance to allow plants to overcome heat stress and complete their life cycle. The targeted production of specific transcript isoforms may underpin priming memory and the resetting of such memory to allow an organism to develop thermotolerance. Splicing repression under heat-stress conditions has been observed in multiple eukaryotic systems and might, therefore, represent a conserved mechanism across eukaryotic species (Chang et al., 2014). Nuclear retention of unspliced transcripts under heat stress may prevent the production of aberrant proteins or peptides, thereby reducing the molecular burden on chaperones and the proteasome machinery (Reddy et al., 2013;Staiger and Brown, 2013). Recent studies indicate that IR has an important role in fine-tuning gene expression depending on the physiological or developmental status of cells or organisms (Boothby et al., 2013;Braunschweig et al., 2014). The unspliced transcripts may represent a cellular reservoir that is ready for splicing upon the relief of heat-stress conditions. In support of this, studies using a fern system showed that partially spliced transcripts are stored in nuclear speckles in association with an exon-junction complex protein and splicing of these stored intron-containing mRNAs is triggered by spermatogenesis during specific developmental stages (Boothby and Wolniak, 2011;Boothby et al., 2013). The subset of genes whose splicing patterns were unaffected by heat stress could result from increased stability of pre-made and processed mRNA or from the maintenance of an efficient splicing process. Perhaps the splicing of these unaffected genes occurs in molecular bodies enriched with splicing factors to ensure efficient splicing under heat stress (Ali et al., 2003). These unaffected genes were enriched in the molecular-function categories pertaining protein folding, which are commonly induced under heat-stress conditions.
We investigated whether the repression of splicing was indiscriminate or whether subsets of genes were differentially spliced under high-temperature conditions. Because IR is the predominant form of AS in plants, we investigated the patterns of IR throughout our priming platform. Specifically, we investigated whether primed plants would maintain IR during the recovery and memory-establishment phases and whether these primed plants would respond differently after the second exposure to heat stress at the AS level. We found that plants responded to heat stress by accumulating significant levels of IR events, indicating splicing repression. During the recovery phase, these IR events declined to normal levels. However, after the second exposure to heat stress, the response of primed plants was drastically different from that of non-primed plants, which accumulated higher levels of IR products, indicating significant splicing repression. Intriguingly, the primed plants did not respond to the second exposure to heat stress by accumulating high levels of IR; these primed plants 'remembered' to undergo splicing and produced splicing patterns similar to those of control plants under normal conditions. Therefore, primed plants can maintain a 'splicing memory', whereby after relief of the stress they exhibit efficient splicing and splicing patterns that support plant growth and development. This splicing memory is a key factor conferring thermotolerance to primed plants. We validated the splicing patterns of several genes from primed and non-primed plants after exposure to lethal heat stress. Our data corroborated the RNA-seq data, showing significant IR retention after the exposure to lethal heat-stress in non-primed plants, whereas primed plants exhibited efficient splicing to support plant growth and development and therefore survived the exposure to lethal heat stress. This is the first report demonstrating that heat-stress priming results in splicing memory, which allows the primed plants to execute efficient splicing upon relief of the stress. Transcripts of key genes regulating plant responses to heat stress were alternatively spliced, and these AS patterns were different between primed and non-primed plants. Moreover, the expression of specific subsets of genes controlling key processes of pre-mRNA splicing, including genes encoding splicing-machinery proteins and other regulatory proteins (such as SR and hnRNP), was fine-tuned by heat-stress priming to ensure proper responses and survival under transient, lethal high-temperature conditions (Palusa et al., 2007;Reddy and Shad Ali, 2011;Yeap et al., 2012;Palusa and Reddy, 2015). Co-transcriptional AS might be controlled epigenetically to produce the AS patterns needed for proper responses to ensure plant survival. In mammalian systems, it has been shown that chromatin organization involving epigenetic modifications and the rate of transcription regulate AS (Naftelberg et al., 2015). The fact that the response to heat stress is transcriptionally regulated at the chromatin level (Kumar and Wigge, 2010;Strenkert et al., 2011;Sawarkar and Paro, 2013) suggests that heat stress-induced AS patterns might be epigenetically controlled, a notion that requires further investigation. Heat stress is known to alter the localization of splicing regulators, which may also be responsible for observed changes in splicing patterns (Ali et al., 2003).
Further studies are needed to explore how splicing memory is acquired and maintained in primed plants, and to identify the molecular mechanisms that underpin the interplay between heat stress, transcriptional regulation, epigenetics, and splicing memory. Under natural conditions, plants are faced with a variety of combined abiotic stress conditions. Our findings open up exciting possibilities for future research on the priming of single and combined stress factors, and further elucidation of the mechanisms that contribute to the induction of splicing memory by heat-stress priming might have practical applications for increasing abiotic stress tolerance.

Supplementary data
Supplementary data are available at JXB online. Fig. S1. Heatmaps showing changes in gene expression in response to heat stress. Fig. S2. Functional categories of heat memory with interactions of heat-memory genes. Fig. S3. Comparisons of DEGs between primed and nonprimed plants after heat shock, and heat-responsive genes. Fig. S4. Gene expression pattern of heat-priming memoryrelated DEGs and DIR-related genes. Fig. S5. Differential AS events and genes between primed and non-primed plants. Fig. S6. Splicing pattern of heat-responsive genes. Fig. S7. IVG images of HSF genes. Table S1. Annotated AS events. Table S2. AS events and related genes at different time-points. Table S3. List of genes with differential intron retention between various time-points.