Plant growth and treatments

Arabidopsis thaliana (Col-0) plants were grown in potting soil in growth rooms at 22°C with a 12-h light photoperiod and light intensity of 180 μmol m^-2 s^-1. Repeated dehydration stresses were performed by air-drying for 2 h followed by a 22 h period of full re-hydration recovery as described [9]. RNA-Seq analyses were performed on rosette leaves from pre-stressed (W) plants, from plants exposed to the first dehydration treatment (S1), and during a third stress (S3) following two stress/recovery cycles. Plants from two independent biological samples were used. After removal from soil and before initiating the stress/recovery cycle, experimental plants were conditioned overnight in humid chambers and used to establish the basal transcript levels (W) for the transcriptome analyses. The transcriptional behavior of the marker memory genes RD29B and RAB18 were monitored in leaf samples taken from fresh plants in soil, from overnight conditioned plants after their removal from soil, and from plants kept under watered conditions throughout the stress/recovery as indicators of the normal dehydration stress behavior of plants in S1 and S2 and internal controls for the whole genome transcriptome analyses.

RNA extraction and RNA-Seq library construction

Leaf tissues were collected and immediately frozen in liquid nitrogen. Total RNA was extracted with Trizol (Invitrogen Inc. Carlsbad, CA, USA), treated with DNase I (Qiagen, Valencia, CA), and purified using Qiagen RNeasy Mini Kit. RNA integrity was confirmed on a Bioanalyzer 2100 using Nano 6000 LabChip (Agilent Technologies, Santa Clara, CA). Complementary DNA sequencing library was prepared from the total RNA using the mRNA-Seq Sample Preparation Kit (Illumina, San Diego, CA). Briefly, poly-adenylated RNA was isolated from 10 μg total RNA by Sera-Mag Magnetic Oligo-dT beads (ThermoFisher Scientific, Waltham, MA). RNA Purified mRNA was fragmented, annealed to high concentrations of random hexamers, and reverse transcribed. Following second strand cDNA synthesis, end repair, and A-tailing, Oligo adapters complementary to sequencing primers were ligated to cDNA fragment ends. Resultant cDNA libraries were size fractionated on an agarose gel, 200 bp fragments excised, and amplified by 15 cycles of polymerase chain reaction. Clusters were generated from the cDNA sequencing library on the surface of a flowcell in the Cluster Station (Illumina) by so-called bridge amplification. Replicates for the watered, S1 and S3 sample libraries were each run on a single lane in a flowcell on an Illumina GAIIx at the Genomics Core Facility at the University of Nebraska-Lincoln.

Reverse transcription and real-time PCR

Total RNA isolation and reverse transcription with oligo(dT) (18418–012, Invitrogen) were performed as described previously [10]. The amounts of individual genes were measured with gene-specific primers by real-time PCR analysis with a CyclerIQ real-time PCR Instrument (Bio-Rad) and SYBR Green mixture (Bio-Rad). The relative expression or amount of specific genes was quantitated with the 2^-ΔΔCt calculation [11], according to the manufacturer’s software (Bio-Rad), where the reference gene was ubiquitin. Primers used in real-time RT-PCR are in Additional file 1: Table S1.

Bioinformatics analysis

Transcriptome sequencing of the watered, S1, and S3 samples yielded a total of 53.4, 76.8, and 78.4 million reads, respectively, summed over the two biological replicates per sample (see Additional file 2: Table S2). The read length for S1 and S3 is 75 bases, while for watered it is 101 bases. To determine the quality of the replicates we performed a least-square simple linear regression for each of the three samples. We calculated the R ² statistic (0.96 ≤ R ²  ≤ 0.99) and slope (1.07 ≤ b ≤ 1.09), which provide measures of goodness-of-fit and correlation, respectively, using the regress function in MATLAB® (version 8.0.0.783 [R2012b]; The MathWorks™) (Additional file 3: Table S3). For use in all further analyses, the Arabidopsis thaliana genome and gene models were downloaded from Illumina’s iGenomes (Ensembl, TAIR10), and gene functional descriptions and Gene Ontology assignments from the Arabidopsis Information Resource (release TAIR10, [12]).

The bowtie (version 2.1.0; [13]) and tophat (version 2.0.8; [14]) packages were used with default parameters to map the RNA sequence reads from watered, S1, and S3 to the genome (see Additional file 3: Table S3), and to determine the expression quantity of known transcripts in each sample. The cuffdiff tool from the cufflinks package (version 2.0.2; [15]) was run with default parameters to calculate expression changes and associated q-values (False Discovery Rate adjusted p-values) for each gene, between the samples S1 and water, and S3 and S1. We further classify genes as being significantly differentially expressed when all three of the following conditions are met: q ≤ 0.05; | log₂(fold change) | ≥ 1; and the FPKM-normalized expression value of at least one sample out of the two needs to be larger than the 25^th percentile. The output files of cuffdiff are further annotated (in-house Perl script) by adding gene functional descriptions and GO classifications, and merged into a master file containing all data for S1 versus water and S3 versus S1 (Additional file 3: Table S3, Additional file 4: Table S4, and Additional file 5: Table S5).

The raw transcriptome sequence files for watered, S1, and S3 have been uploaded, together with gene expression result files, to NCBI’s Gene Expression Omnibus under sequence number GSE48235.

Dehydration stress memory response genes of Arabidopsis thaliana

Watered non-stressed plants and plants subjected to one or three dehydration stresses by exposure to dry air for 2 hours followed by 22 h periods of watered recovery intervals [9] were analyzed. The levels of mRNAs from the leaves of these plants were determined by RNA-Seq. Transcripts from 33,555 A. thaliana genes were identified (Table 1; Additional file 3: Table S3). Of these, 6579 genes (~20% of the genes in the genome) significantly increased or decreased their transcript levels during the first dehydration stress (S1), compared to their transcript levels in watered plants (W), representing the dehydration stress-responsive genes of A. thaliana. Of the genes responding in S1, 3396 genes were up-regulated and 3183 were down-regulated (Table 1; Additional file 4: Table S4). A comparison of the transcript levels during the first (S1) and third dehydration stress (S3) identified 1,963 dehydration responsive genes displaying significantly different amounts of transcripts in S3 versus S1. Our operational criterion for transcriptional memory is that the transcriptional responses to similar stress conditions must be different. Accordingly, these 1,963 genes display dehydration stress transcriptional memory (Table 1; Additional file 4: Table S4).

Four distinct transcriptional memory response-patterns were recognized within the S1-responding fraction. We designate these transcription patterns as [+/+], [-/-], [+/-] and [-/+] memory responses, where the first sign indicates higher (+) or lower (-) transcript levels in S1 relative to the levels in pre-stressed watered plants (W). The second sign indicates higher (+) or lower (-) transcripts in S3 relative to the levels in S1. Non-memory genes are denoted as [+/=] or [-/=] indicating that up-regulated or down-regulated transcript levels occurring in S1, respectively, have similar levels in S3. In addition, 798 genes were up-regulated in S2 ([=/+] W = S1 < S3) and 573 were down-regulated ([=/-] W = S1 > S3). These delayed response patterns are different from the memory responses represented by the four S1-responding memory types. Here, we focus mainly on the transcriptional patterns of the S1-responding gene fraction of Arabidopsis.

Of the four different memory classes, 362 genes denoted as [+/+] memory genes were induced in S1 and induced to higher levels in S3 (Table 1; Figure 1A). All 7 memory genes described in our earlier study [9] were present in the [+/+] memory category in the genome-wide RNA-Seq data (Additional file 4: Table S4). The [-/-] class contained 310 genes with decreased transcript levels in S1 and further reduced levels in S3 (Table 1; Figure 1B). The [-/+] class contained 434 memory genes that were down-regulated in S1, but produced significantly higher transcript levels in S3 (Figure 1C); conversely, 857 [+/-] memory genes, were induced in S1 but had lower levels in S3 (Figure 1D). The latter two categories ‘revise’ their transcriptional behavior in a subsequent stress: after robustly responding in S1, these genes show weaker/no responses in S3, producing transcripts at levels close to their initial watered (pre-stressed) levels. These genes are referred to as ‘revised response’ memory genes. The distribution of transcript levels for the dehydration stress-responsive memory genes in S1 and S3 illustrating the four memory categories is shown (Figure 2A, B). The two non-memory gene classes ([+/=] and [-/=]) consistently provide transcriptional responses in S3 similar to those in S1 (Figure 1E, F).

Functional distribution of Arabidopsis memory genes

The large number of the dehydration stress-memory responsive genes and the diversity of their transcriptional responses raised questions of their functions and of whether there was a preferential association of particular cellular functions with any of the four memory types. To address these questions, genes displaying transcriptional memory were analyzed for biological function according to their Gene Ontology (GO) classification (Table 2; Additional file 5: Table S5).

Genes involved in ABA/abiotic stresses responses are signature for the [+/+] memory genes

Genes implicated in responses to salt, salinity, cold/heat acclimation, and abscisic acid (ABA) constitute about a quarter of the [+/+] memory genes. LEA genes represent 3% among them (Table 2; Additional file 5: Table S5; Figure 3).

Ribosome/protein synthesis, chloroplast, and thylakoid membrane associated functions are signature for [-/-] memory genes

Memory genes encoding ribosomal, as well as chloroplast and photosynthesis proteins are the major constituents of the [-/-] memory subgroup (Additional file 4: Table S4; Additional file 5: Table S5; Figure 3). Encoded proteins are involved in ribosome structure and amino acid biosynthesis, as well as in photosynthesis, in the light harvesting complex (photosystem II), and in responses to light (red/far red light, in particular). The largest number of thylakoid membrane-associated proteins is encoded by [-/-] memory genes.

Functions encoded by [-/+] revised response memory genes

There are no particular types of cellular functions specifically enriched within the [-/+] group of memory genes, although the largest percentages of genes encode chloroplast and membrane (plasma, organellar, and thylakoid membrane) related functions. Many of these are for electron transport, photosystem II assembly, chloroplast organization and re-location, functions similar to those encoded by [-/-] memory genes (Table 2; Additional file 5: Table S5; Figure 3).

Dehydration stress [+/-] memory genes are shared with multiple signaling pathways

The signature [+/-] dehydration memory genes are regulated by multiple signaling pathways, including the ABA, ethylene, auxin (IAA), gibberrellic acid (GA), jasmonic acid (JA), and salicylic acid (SA) pathways (Table 2; Additional file 5: Table S5; Figure 3). A prominent functional group represented by [+/-] memory genes encodes proteins associated with membranes. Although membrane (plasma and organellar)-related genes are highly represented within the entire drought-responding fraction (including the non-memory and all four memory types), genes for tonoplast intrinsic proteins regulating water transport (GAMMA-TIP, PIP2B, TMP-C, PIP2A, TIP2, PIP1A, RD28) and the inward K⁺ channel proteins (KAT1 and KAT2) display, exclusively, [+/-] memory responses (Table 2; Additional file 5: Table S5).

Dehydration stress memory of transcription factor genes

Among the ~1500 Arabidopsis genes encoding transcription factors (TFs) [16], members of about two dozen families have been implicated in responses to drought [17–19]. Among these families, transcriptional memory behavior was found for members of ten of the families (Additional file 6: Table S6). Five families (AP2/ERF, bHLH, Homeodomain-Zip, MYB/Myb-like, and ZF) are represented in all memory categories; GRAS and HSF members display only [+/+] and [+/-] memory, and the CAAT family is represented only by [+/+] memory genes. TF genes are least represented in the [-/-] memory subgroup (Table 2; Additional file 6: Table S6).

[+/+] memory genes for MYB2, MYB112, MYB13, MYB47, ATHB-7, ATHB-14, and ABF2 specifically cluster with ABA/abiotic stress responding genes but are rarely found among other hormonally responsive pathways (Figure 4A); six putative TF genes display [-/-] memory type responses and only three are shared with other response pathways (Figure 4B). TFs with [+/-] or [-/+] revised response memory behavior are broadly shared with multiple abiotic and hormone response networks, outlining distinct clusters (Figure 4C, D).

Biological relevance of the transcription memory behavior

Plants experience more dehydration stress during the day when transpiration exceeds the rate of water uptake and fully, or partially, recover during the night when the rate of water uptake exceeds the greatly diminished loss of water from leaves [20]. This creates a diurnal cycle of oscillating leaf water potential that is probably most pronounced during periods of moderate drought conditions [21, 22]. The repetitive dehydration stress system employed in this study has similarities to this natural diurnal stress. We hypothesize that dehydration stress memory helps plants to prepare for the next day’s stress if they were stressed the day before, despite alleviated stress signals during the night [9].

Genes displaying transcriptional memory are among the genes responding to the first dehydration stress. Altering their expression levels in subsequent stresses, presumably, allows the plant to finely tune its responses to ongoing/recurring dehydration stress. The possible biological relevance of the genes displaying transcriptional memory is considered in the context of four overlapping strategies generally employed by a plant to improve its stress tolerance and/or survival: 1) increased synthesis of membrane protecting, damage-repairing, and detoxifying functions; 2) coordinating photosynthesis and growth under repetitive stress; 3) re-adjusting osmotic and ionic equilibrium to maintain homeostasis [23]; and 4) re-adjusting interactions between dehydration and other stress/hormone regulated pathways [24, 25].

Increased synthesis of protective functions

In general, the proteins encoded by [+/+] memory genes are involved in cell protective roles. These include heat shock proteins and chaperons, proteins involved in the repair of stress-induced damages to membranes, in modifying membrane lipid composition, in regulating membrane fluidity and permeability to toxic ions, dehydrins (RD29B, RAB18, LTi30), and lipid transfer proteins (LTP2, LTP3, LTP4) [26, 27]. LEA (late embryogenesis abundant) proteins functioning as molecular chaperones to maintain membrane structures, ion balance and homeostasis [28, 29] and enzymes for the synthesis of isoleucine (toxin degradation), serine (redox responses), and proline (an osmolyte) are encoded by [+/+] memory genes. Thereby, [+/+] memory genes ensure elevated synthesis of factors critical for cell survival under multiple dehydration stresses (Table 2; Figure 3).

Abscisic acid (ABA) is the key mediator of dehydration signaling and is involved in responses to other abiotic stress response pathways [30, 31]. ABA responsive genes that respond to dehydration are often points of overlap in various abiotic stress pathways [32, 33] (Table 2). The [+/+] memory and [+/=] non-memory genes, generally, encode proteins with similar protective functions. However, the memory component, presumably, allows the plant to optimize its stress tolerance within the context of protection. Only a few [-/-] memory genes and a smaller number of down-regulated [-/=] genes are shared between dehydration/ABA and other abiotic stress responding systems (Table 2; Additional file 5: Table S5; Figure 3).

Photosynthesis and growth during repeated stresses

More than 40% of the [-/-] memory genes encode ribosome/protein synthesis, chloroplast, and thylakoid membrane related proteins indicating repetitive stress results in a stronger attenuation of the pathways for photosystem II electron transport, small RuBisCo subunit and sugar biosynthesis, as well as for ribosome assembly, protein translation, DNA replication and histones (Additional file 4: Table S4; Additional file 5: Table S5; Figure 3). None of the [+/+] memory genes encodes a thylakoid membrane function. Of note, 12 [+/+] memory genes associated with chloroplasts encode functions that are different, even opposite to, the functions encoded by the [-/-] memory genes. Chlorophyll and toxin catabolic enzymes, enzymes for the synthesis of isoleucine (toxin degradation), serine (redox responses), and proline (an osmolyte) are encoded by [+/+] memory genes. Photosynthetic and starch synthesizing activities, in addition to proteins involved in chloroplast organization and relocation, are encoded by [-/-] memory genes. Down-regulated [-/=] non-memory genes encode ribosome- and chloroplast-associated proteins as well (Table 2; Additional file 5: Table S3; Figure 3). The functional implications of the [-/-] transcriptional response patterns (decreased protein synthesis, cell growth and photosynthesis) and the protective functions by the [+/+] memory genes are in agreement with known responses of plants enduring dehydration stress [34–37]. Of note, on the background of a rather even distribution of GO functions encoded by the late responding [=/+] and [=/-] genes, the fraction of genes encoding chloroplast and thylakoid membrane proteins constitute ~ 20% among the [=/-] genes (Table 2; Additional file 5: Table S5).

Re-adjusting cellular homeostasis

The [-/+] and [+/-] revised memory genes respond robustly to the initial stress, but in a subsequent exposure provide responses at levels closer to their non-stressed levels. Presumably, these genes help in restoring homeostasis as the plant adjusts to dehydration stress. Thus, [-/+] memory genes for chloroplast organization, relocation, photosynthetic, light harvesting, and metabolic functions are similar to the functions encoded by the [-/-] memory genes (Additional file 5: Table S5; Figure 3). However, by producing significantly more transcripts in S3 than in S1([-/+]), or returning to pre-stressed (W) levels in S3 ([+/-]), these genes support the ability of plants to restore photosynthesis after temporary water withdrawal [38]. Presumably, the [-/+] genes contribute to an initial decrease in chloroplast/photosynthetic activity in S1, but by reverting to pre-stressed transcription levels in S3, contribute to a subsequent re-adjustment. The [-/+] memory genes GLK1 and GLK2, implicated in regulating the photosynthetic apparatus and chloroplast development in a cell-autonomous manner [39, 40] are candidate genes for a role in such an adjustment.

Within the same GO category, [-/+] and [+/-] memory genes may encode opposing biochemical activities. For example, a number of the [-/+] memory genes encode chlorophyll and carbohydrate (starch) synthesizing enzymes, while chlorophyll and starch degrading activities are encoded by [+/-] memory genes (Additional file 5: Table S5). Juxtaposing these functions would suggest complementary roles during both the initial and repeated stresses. Acting together, these genes are likely to decrease overall chlorophyll and starch biosynthesis during the first stress but to partially restore it during subsequent stresses. These complementary transcriptional patterns and functional roles support the fine coordination occurring between metabolic and energy adjustments during adaptation to drought [41].

Membrane-associated memory genes present another paradigm of biological relevance for the dehydration memory genes. As dehydration may damage membrane integrity, it is not surprising that a large number of genes (more than 750 of memory and non-memory) encode functions related to membranes; more than 260 of the late responding genes encode membrane-associated proteins (Table 2; Additional file 5: Table S5; Figure 3). Proteins regulating osmotic pressure, water balance, and wall modifications have been implicated in plants’ stress responses and environmental adaptation [42]. These functions are encoded by genes of all four-memory response types, suggesting that after providing an up/down transcriptional response in S1, a large number of Arabidopsis genes alter their transcript outputs (further up- or down-regulating or reversing towards pre-stressed levels) in a subsequent exposure to more finely modify cell membrane structure, ion balance and homeostasis.

Dehydration stress memory genes in crosstalk with other response pathways

Individual genes that are regulated by multiple different stimuli represent overlapping points of these plant’s signaling networks [43–48]. The subset of the dehydration stress memory genes that are at these convergent points revealed new aspects of plants’ stress responses. For example the dehydration/JA pathways share the largest number of both memory and non-memory genes co-regulated by these pathways, consistent with reported crosstalk between jasmonic acid (JA) and ABA regulated signaling networks during dehydration [45, 47, 49], salt [50], and cold stresses [51]. More than 200 JA-responsive genes are up-regulated in S1 (including memory and non-memory genes) and about 80 are down-regulated (Table 2; Additional file 5: Table S5). In S3 stress, however, 89 [+/=] and 46 [-/=] genes continued to provide similar transcript levels, seven produced higher ([+/+]), and 16 lower ([-/-]), transcript levels, but 140 ([+/-]) genes did not respond in the subsequent stress. Thereby, in addition to a finer tuning, the altered transcript levels from the memory genes suggest an altered crosstalk between these pathways in S3 compared to S1.

Non-memory ([+/=] and [-/=]) and unidirectional memory ([+/+] and [-/-]) dehydration stress responding genes respond every time the plant experiences dehydration stress. Concurrent with amplified responses from the [+/+] and [-/-] memory genes, however, engagement of the [+/-] and [-/+] genes lessens in S3, suggesting the nature and the dynamics of the interactions between dehydration and other signaling pathways are different during repeated exposures than the interactions occurring in a single exposure.

Lastly, we note that TF encoded by memory genes can be critical for the expression of dependent genes in S1, but that they do not necessarily determine the memory behavior of regulated genes in S3. For example, the [+/-] memory of the MYC2 gene, identified as a master regulator of the crosstalk between the ABA, SA, GA, JA, and auxin signaling pathways [52–54], correlates with the [+/-] memory of a large number of MYC2-dependent genes common for these signaling pathways (Figure 4C; Additional file 6: Table S6). However, the signature gene, RD22, directly regulated by MYC2 [55, 56] is a non-memory gene that is up-regulated in both S1 and in S3 (Additional file 4: Table S4; Additional file 6: Table S6) despite the fairly low levels of MYC2 transcripts in S3. Evidently, another TF activates RD22 in S3. The results suggest that diverse gene-specific mechanisms are involved in regulating the behavior of dehydration stress responding genes and that the memory behavior of individual TFs alone does not necessarily determine or predict the memory or non-memory behavior of their targets.