- Research article
Comparative transcriptome analyses revealed different heat stress responses in high- and low-GS Brassica alboglabra sprouts
BMC Genomicsvolume 20, Article number: 269 (2019)
Chinese kale (Brassica alboglabra) contains high nutritional elements and functional molecules, especially anticarcinogenic and antioxidant glucosinolates (GS), which was highly affected by environment temperature. To investigate the link of GS biosynthesis with heat stress response in Chinese kale, global transcription profiles of high-GS line (HG), low-GS line (LG), high-GS line under heat stress (HGT) and low-GS line under heat stress (LGT) were analyzed.
Based on three biological replicates of each RNA sequencing data, 3901, 4062 and 2396 differentially expressed genes in HG vs HGT, LG vs LGT and HGT vs LGT were obtained, respectively. GO annotation, KEGG pathway analysis and a comprehensive analysis of DEGs showed a strong correlation between the GS biosynthesis and heat stress response. It was noticed that 11 differentially expressed genes tied to the GS biosynthesis were down-regulated, 23 heat shock transcription factors and 61 heat shock proteins were up-regulated upon the heat treatment. Another two Chinese kale varieties Cuibao and Shunbao with high- and low- GS content respectively, were used to validate the relationship of GS content and heat-response, and the results showed that high-GS content variety were more thermotolerant than the low-GS content one although GS significantly decreased in both varieties under heat stress. In addition, HSP100/ClpB, HSP90, HSP70 and sHSPs were differentially expressed in high- and low-GS varieties. Notably, HSP90 and sHSPs showed an obviously early response to heat stress than other related genes.
The higher heat resistance of high-GS Chinese kale and the sharp decrease of glucosinolate content under heat stress indicated a strong relationship of GS accumulation and heat stress response. Combined with the previous report on the low expression of HSP90 at elevated temperatures in GS-deficient mutant TU8 of Arabidopsis, the differential expression pattern of HSP90 in high- and low- GS varieties and its early heat response implied it might be a key regulator in GS metabolism and heat-resistance in Chinese kale.
Brassica species, included many important crops in the Brassicaceae family, play an important role in agriculture and horticulture fields and contribute greatly to economies and population health worldwide. Brassica alboglabra, also known as Chinese kale or Kai Lan, is widely distributed in Southern China and Southeast Asia . B. alboglabra is preferentially consumed in China owing to its tender bolting stem and leaves, good flavor and high nutritional content in anticarcinogenic compounds and antioxidants, such as total phenolics, carotenoids, and glucosinolates (GS) [2,3,4] .
GS also known as mustard oil glucosides, are a group of nitrogen- and sulfur-rich secondary metabolites, mainly found in the order Capparales. They can be grouped into indole, aliphatic, and aromatic GS according to the characteristic of the amino acids they are derived from, such as Trp (indole GS), Ala, Ile, Leu, Met, or Val (aliphatic GS), and Phe or Tyr (aromatic GS) [5,6,7]. Generally, the biosynthesis of GS include three separate phases: the elongation of aliphatic GSs chain, the formation of the core structure, and the modification of the side chain of GS . In recent years, a wide variety of crucial transcriptional factors and biosynthetic genes involved in the GS metabolism and regulation of B. oleracea , B. alboglabra [2, 3, 10], B. rapa [11, 12] were identified and studied.
GSs, also affected by some abiotic stresses, have been implicated to play an important role in plant abiotic stress response. A few studies have demonstrated that the GS biosynthesis was linked to heat stress [10, 13,14,15]. Ludwig-Müller et al. (2000) found that the deficiency in GS metabolism of Arabidopsis thaliana was remarkably less tolerant upon exposure to elevated temperatures than wild-type plants. Charron and Sams (2004) discovered that high temperatures could enhance higher aliphatic GS contents of Brassicaceae vegetables. Hara et al. (2013) showed that exogenous application of isothiocyanates in Arabidopsis enhanced thermo-tolerance and induced the expression of HSP genes. Guo et al. (2016a) found that heat and hypoxia stresses enhanced the accumulation of aliphatic GSs and sulforaphane in broccoli sprouts.
Higher temperature, one of the most vital abiotic stresses, severely restricted plant growth and development . To survive high temperatures, plants have evolved sophisticated response mechanisms to in the long-term evolution. The response mechanism was orchestrated by modifying biochemical and physiological factors at the cellular and molecular levels. The thermal tolerance involved the expression of a set of heat shock protein (HSP) genes, mediated by effector genes and various heat shock transcription factors (HSFs). HSFs, widely existed in both plants and animals, played an important role in regulation of this heat-induced transcriptional reprogramming . HSPs can act as molecular chaperones, protect cells against heat damage and function in protein folding, as well as the intracellular distribution and degradation of other proteins .
In recent years, RNA sequencing (RNA-Seq) has becoming a popular tool for profiling gene expression and for identifying genes or molecular pathways that are differentially expressed between two or more biological conditions . Recent examples [2, 20, 21] of using RNA-Seq to investigate GS metabolism and heat stress found that genes related to GS metabolism were down-regulated under heat stress in B. napus silique wall , the aliphatic GS and sulforaphane were accumulated in broccoli sprouts under the stresses of heat, hypoxia and heat plus hypoxia , and that the exogenous allyl-isothiocyanate treatment triggered the expression of genes related to GS metabolism, sulphate uptake and assimilation, heat stress response, oxidative stress response, elicitor perception, plant defence and cell death mechanisms .
In this study, the RNA-Seq was adopted to examine the transcriptome of high-GS line (HG), low-GS line (LG), high-GS line under heat stress (HGT) and low-GS line (LG) under heat stress of B. alboglabra. From the transcriptome, differentially expressed genes (DEGs), especially participated in GS biosynthesis and heat stress response, were identified and characterized between HG and LG sprouts. Combining the profiles and contents of GS in another two Chinese kale varieties with high- and low- GS content respectively, and quantitative real-time PCR assay of HSF and HSP genes, we revealed a strong correlation between GS biosynthesis and heat stress response.
Analysis of RNA-Seq data in HG and LG sprouts under heat stress
RNA-Seq was used to examine the transcript profiles in both the high-GS line and low-GS line with three biological replicates at 25 °C and 42 °C, respectively. After the filtration, the reads from RNA-seq were aligned to B. oleracea reference genome, and the mapped ratio ranges from 76.95 to 83.21% (Table 1). Most of the total reads (the uniquely mapped ratio ranges from 71.95 to 80.86%) were uniquely mapped to the B. oleracea genome, whereas a small proportion was mapped to multiple locations in the B. oleracea genome. The mean of GC percentage and Q30 percentage are 47.70 and 88.64%, respectively. Data for untreated high-GS and low-GS line were listed in Guo et al. (2016) .
The sample homogeneity of the four groups with three replications were analyzed. RNA collected from the three replications in every sample had a low biological coefficient of variation (BCV) value, indicating that gene expression in three replications of every sample was homogeneous (Fig. 1a). This suggests that the replicated samples produced data that are acceptable for further analyses. The expression levels of 12 cDNA libraries were delineated by boxplot profiles (Fig. 1b), and a pairwise comparison of gene expression levels of the four groups was shown as Fig. 1c, showing that the samples were highly reproducible and the experimental results were reliable and credible. Gene differential expression analysis identified 8234 DEGs (log2|Fold change| ≥1 and FDR ≤ 0.01) with 2716 DEGs in HG vs LG, 3901 DEGs in HG vs HGT, 4062 DEGs in LG vs LGT and 2396 DEGs in HGT vs LGT respectively (Fig. 1d).
Alternative splicing (AS) analysis after heat stress
AS events had occurred and there was an average of 13,000 AS events in every sample. AS events were classified into 6 types using soft Cufflinks, and the largest and smallest AS event groups were intron retention event and alternative last exon event, respectively (Fig. 2). The total AS events occurred in the samples germinated under heat stress were more than that in the samples germinated at normal temperature, and the statistically significant difference was mainly reflected in intron retention, suggesting that introns could play an important role in environmental responses. These results showed that the probability of AS events occurred more when circumstances changed, and that the organism could orchestrate sophisticated and refined regulatory mechanisms to respond to the altered environment and develop effective mitigation strategies.
GO analysis of DEGs in response to heat stress
GO was used to analyze the functional classifications of the DEGs in HG and LG grown at 25 °C and 42 °C, respectively (Fig. 3). The DEGs were divided based on cellular component (CC), molecular function (MF) and biological process (BP). The largest CC, MF and BP subcategory for the DEGs in HG vs HGT were plasma membrane, ATP binding and response to salt stress, which comprised 33.44, 14.82 and 11.90% of the DEGs in the subcategory, respectively. This data suggested that the high-GS line underwent more exquisite metabolic activities in the plasma membrane upon response to heat stress. The largest CC, MF and BP subcategory for the DEGs in LG vs LGT were nucleus, protein binding and response to chitin, which comprised 34.69, 23.13 and 14.23% of the DEGs in the subcategory, respectively. It suggested that the low-GS line conducted more sophisticated metabolic activities in nucleus response to heat stress. Furthermore, the largest CC, MF and BP subcategory for the DEGs in HGT vs LGT were nucleus, protein binding and response to salt stress, which comprised 40.96, 21.37 and 12.27% of the DEGs in the subcategory, respectively, suggesting the different response to heat stress in HG and LG were enriched in protein functions of nucleus. In summary, these results indicated the different response mechanism between HG and LG when confronted with heat stress.
KEGG analysis of DEGs in response to heat stress
The KEGG Pathway Database was initiated to identify the biological pathways of DEGs, and the DEGs in HG vs HGT, LG vs LGT and HGT vs LGT were found to be mainly enriched in the 104, 102 and 106 KEGG pathways, respectively. The most enriched KEGG pathway in HG vs HGT, LG vs LGT and HGT vs LGT were protein processing in endoplasmic reticulum, plant hormone signal transduction and plant hormone signal transduction, respectively (Fig. 4 and Additional file 1: Figure S1), which coincided with the results of GO analysis. Those DEGs involved in protein processing in endoplasmic reticulum related to translocon, ribosome anchor, protein recognition by luminal chaperones, N-glycan biosynthesis, re-glucosylation, de-glucosylation, lectin associated, accumulation of misfolded proteins, protein targeting, ER-associated degradation and ubiquitin ligase complex. However, compared to the DEGs in HG vs HGT, there was more DEGs in LG vs LGT and HGT vs LGT for plant hormone signal transduction related to tryptophan metabolism, zeatin biosynthesis, diterpenoid biosynthesis, carotenoid biosynthesis, cysteine and methionine metabolism, brassinosteroid biosynthesis, a-linolenic acid metabolism and phenylalanine metabolism. These changes indicated significant changes in protein modification in response to heat stress in HG, and plant hormone signal transduction in response to heat stress in LG.
Furthermore, the significantly enriched 20 KEGG pathways in HG vs HGT, LG vs LGT and HGT vs LGT were general metabolic pathways and biosynthesis of secondary metabolites (Fig. 5 and Additional file 1: Figure S2). Of these, brassinosteroid biosynthesis, glycosaminoglycan degradation, oxidative phosphylation and phagosome uniquely represented in the HG vs HGT (Additional file 1: Figure S2A), endocylosis, ether lipid metabolism, plant hormone signal transduction and plant-pathogen interaction uniquely represented in the LG vs LGT (Additional file 1: Figure S2B). However, the significant enriched 20 KEGG pathways in HGT vs LGT were clearly different from those in HG vs HGT and LG vs LGT (Fig. 5). It should be noted that the enrichment factor of GS biosynthesis was the second smallest and its Q-value was the highest of all the KEGG pathways in HGT vs LGT, indicating that high- and low-GS Chinese kale sprouts did have different heat stress responses ascribed to the GS biosynthesis.
Analysis of DEGs involved in GS biosynthesis
First of all, we found that 25 genes and 10 transcription factors related to GS biosynthesis altogether were differentially expressed (Tables 2 and 3), and almost all DEGs in HG vs HGT and LG vs LGT were down-regulated and accounted for the decline of GS content under heat stress. Notably, the DEGs in HG vs LG showed similar expression pattern with those in HGT vs LGT, suggesting that GS biosynthesis pathway is contributed to the differential responses between HG and LG to heat stress. This finding correlate with the results shown in Fig. 5. Together, these results demonstrated that DEGs of GS biosynthesis were down-regulated during heat stress, and that the difference of GS biosynthesis in both high- and low-GS sprouts was conducive to their different responses to heat stress during growth stage.
Analysis of DEGs involved in heat stress response
To identify genes responsive to heat stress in high-GS and low-GS line, we analyzed the HSFs and HSPs. In total, 23 HSFs and 61 HSPs were differentially expressed in HG vs LG, HG vs HGT, LG vs LGT, and HGT vs LGT (Table 4). Overall, it also showed that FC values of most of HSF and HSP transcripts in HG were slightly lower than that in LG exposed to heat stress, which implied that low-GS B. alboglabra was less thermotolerant to heat stress.
In plants, there were 3 classes in the HSF protein families, classes A, B and C, which were discriminated by peculiarities of their flexible linkers and HR-A/B regions . In this study, 6 class A, 4 class B and 1 class C were differentially expressed under heat shock stress, and the detailed fold change (FC) values of these HSF members were provided (Table 4). Four interesting features were worth mentioning. First, the 3 members (Bo4g190180, Bo2g165560 and Bo4g190160) of class A4a group and one class B4 were found to be constitutively expressed and were not up-regulated by heat stress, therefore, it is speculated that they were generally present in an inactive form. Secondly, in contrast to the expressions of classes B and C, the expressions of class A were different and diverse in response to heat stress, suggested that the class A was paraphyletic and HSFA genes were functionally characterized extensively. Thirdly, the A7, B1 and B2a Hsf genes with highly expressed under heat stress were likely to be important HSF members for keeping the expression of Hsp genes at a high level. Lastly, HsfC1 was down-regulated in response to heat stress, which was different from HSFA and B, indicating distinct regulational modes of HSFs during heat stress response.
HSPs were divided into five conserved classes: HSP100/Clp, HSP90, HSP70, HSP60 and the small heat-shock proteins (sHSPs) . In this work, one ClpB1, 8 HSP90s, 4 HSP70s and 17 sHSPs were identified (Table 4). All sHSPs (except Bo8g007140) were up-regulated expressed in response to heat stress to endow cells with thermotolerance. Notably, the expressions of two HSP17.6 (Bo4g169420 and Bo3g130430) in line JL09 were higher than in JL08 under same condition, suggesting that some HSP17.6 in JL09 may be constitutively expressed at relatively high levels and more susceptible to high temperature induction. However, the expressions of HSP90 and HSP70 were intricate compared with sHSPs: three HSP70s in the same line were up- or down-regulated expressed induced by heat stress, the expression of HSP70 (Bo1g078880) and ClpB1 were not influence by heat shock, which were similar with HSF4a (Bo4g190180, Bo2g165560, and Bo4g190160), while all HSP90 were up-regulated but the expression of HSP90 in HG and LG were differently induced by heat stress. Furthermore, we found that the co-chaperone of HSP70, DnaJ, was also differentially expressed under heat shock conditions.
Different heat response in HG and LG Chinese kale sprouts under heat stress
Another two Chinese kale varieties Cuibao and Shunbao were used to validate the differential heat-stress response in high- and low- GS content sprouts. To examine the GS diversity of both lines, we applied HPLC method to compare the GS profile and the content. As expected, Cuibao and Shunbao sprouts mainly contained aliphatic GS including gluconapin (GNP) and glucoiberin (GIB), and a small quantity of indole GS (IGS). As well, the initial GS contents (0 h) of Cuibao were significantly higher than those of Shunbao (Figs. 6a, b and c). The GS contents (GNP, GIB and IGS) of Cuibao and Shunbao were sharply reduced when grown at 42 °C for 48 h compared to the control (Figs. 6a, b and c), suggesting that GS contents in both lines o were decreased by heat stress. It was also found that after treated with heat stress for 72 h Cuibao showed more thermotolerance than Shunbao (Fig. 6d). These results indicated that high-GS content variety was more resistant to heat stress than low-GS variety.
Expression pattern of heat-resistant related genes in high- and low- GS sprouts
To investigate the expression pattern of HSFs and HSPs in high- and low- GS sprouts, a time course of expression of these genes was performed (Fig. 7). Most of the genes showed relatively high expression in Shunbao at 42 °C relative to Cuibao at 42 °C with little to no major changes over time, except for HSF4a2, DnaJ2, ClpB1, and HSP90. Especially, the expression of HSFA4a1 was delayed compared with HSFA1e, HSFA4a2, HSFA4a3, and HSFB3 in Shunbao at 42 °C. The expressions of HSP70, DnaJ1, DnaJ2 and sHSP2 in Shunbao and Cuibao were later than HSFs, suggesting the regulation of HSF on some of HSPs. However, the expressions of ClpB1 and all HSP90s in Shunbao and Cuibao were rapidly activated in heat treatment, except for HSP90–2 in Cuibao, indicating that there were some HSF-dependent pathways for activated expression of HSPs in response to heat stress.
With global warming growing more severe, heat stress is becoming an agricultural problem in many areas in the world. Heat stress has a significant adverse impact on plant growth and development, which leads to agricultural productivity worldwide . Thus, the molecular mechanism of the heat stress response and breeding of heat-tolerant plants is necessary to protect food production and ensure crop safety. To survive high temperatures, many heat-responsive transcription factors and genes are induced to protect plants from heat shock stress . Since the GS contents of B. alboglabra sprouts greatly changed, this study was initiated with the aim to explore possible crosstalk between GS metabolism and heat shock stress response.
In this work, we performed comparative studies to assess the effects of exposure heat stress of high- and low-GS content sprouts at the transcriptomic levels. Our particular interest was the identification of DEGs involved in GS metabolism in response to heat stress, which could help explain the regulatory mechanism of GS metabolism under heat stress.
GS accumulation under heat stress condition
The mechanism of GS accumulation under high temperatures has not been conclusive, because most of the earlier studies focused on higher temperatures as a way to improve GS contents in Brassica sprouts or cruciferous products. In sprouts, the change of GS contents under high temperatures is increased upon heat treatments. The 7-day old to 11-day old broccoli sprouts grown under three day/night temperature regimes (30 °C /15 °C, 22 °C /15 °C, and 18 °C /12 °C) showed decreased contents of glucoraphanin and total GS contents under higher temperature condition (30 °C /15 °C) . When higher temperature (60 °C) were used to treat broccoli sprouts every 24 h, the aliphatic GS including GRA, GNA, and GER in 3-day, 5-day, and 7-day old sprouts were significantly increased . In plants, it is proven that the concentration of GS in B. oleracea leaves and stems were higher at 32 °C while at 22 °C the GS in roots decreased . The application of GS degradation products ITCs including phenethyl ITC, methyl ITC, and allyl ITC can mitigate the growth inhibition induced by heat stress at 55 °C for 1 h in Arabidopsis. The pre-administration of phenethyl ITC increased the expression level of heat tolerance related heat shock protein (HSP) 70 .
The thermal degradation of GS was first demonstrated in synthetic glucobrassicin . After that, cooking method including boiling, microwaving and canning were used to study the effect of temperature on GS reduction. Severe heat treatments like canning were found to result in 73% degradation of GS in red cabbage [28,29,30]. The kinetic modeling of different kinds of GS degradation in broccoli and Chinese kale after treated by a heating block at 100 °C were established in which it was found that the GS degradation rates were related to their structure and the kinetic parameters that was affected by planting season . The heat treatment for broccoli sprouts powder was carried out at 130 °C and the results showed that GS with different structures has different degradation rates. The methylsulfinyl GS (glucoiberin, glucoraphanin and glucoalyssin) are more stable than the methylsulfanyl GS (glucoiberverin and glucoerucin) and among the methylsulfinyl GS a longer chain (glucoalyssin) results in increased thermal breakdown . It is ambiguous how high temperature affects the GS in Brassica plants. Mutant experiment showed that the expression of HSP 90 were reduced at elevated temperatures in GS-deficient mutant TU8 of Arabidopsis , indicating a potential relation between high temperatures and GS biosynthesis. In our study, the GS content was decreased dramatically and the expression of genes related to GS biosynthesis were down-regulated in Cuibao and Shunbao after treated with the high temperature (42 °C) for 48 h. This is consistent with the reduction of GS related gene expressions in heat-treated Brassica napus’ silique wall at the seed-filling stage . However, in low-GS variety JL-09, the expression of GS biosynthetic genes was up-regulated and sulfur metabolism genes were down-regulated compared with the high-GS variety JL-08 under 42 °C condition, indicating that under heat stress condition, the metabolic flux in low-GS content variety were orienting GS biosynthesis. Maybe the amount of GS is related to the heat response.
Comparison of the transcriptional response in different GS-content Chinese kale sprouts under heat stress
HSFs are the terminal components of signal transduction and mediate the expression of HSPs responsive to both heat stress and a large number of chemical stressors . In plants, there are 3 classes in the HSF protein family (classes A, B, and C). Among these, HSF A genes generally serve as bona fide transcription factors binding to heat shock element (HSE) that are found in Arabidopsis . HSF B members have been demonstrated to act as repressors or co-activator of class A HSFs and were also capable of binding to HSE [34,35,36]. Yet knowledge about the biological role of the HSF B and HSF C members was scarce. However, individual HSF members had species-specific features in regulating genes involved in the heat stress response . In Lycopersicon esculentum, there are at least 25 HSF protein families. HSF A1a acted as a master regulator for triggering the heat response and acquired thermo-tolerance , and HSF B1 functions as a transcription co-activator, cooperating with class A HSFs during heat stress . Whereas, there are 21 HSF family members in Arabidopsis thaliana and no master regulator has been found . Furthermore, three subclass A1 HSFs (A1a, A1b, and A1d) in A. thaliana were functionally redundant for triggering the heat response [41, 42], and HSF B1 was found to act as a repressor [34,35,36, 43].
In this work, 23 HSFs were expressed differentially in B. alboglabra sprouts. It was shown that the expressions of HSF genes in Chinese kale were significantly up-regulated by heat shock stress, suggesting commonalities of HSFs in the heat stress response in these two lines. Unexpectedly, we found that all HSFs except HSFB4 interacted with HSP70 (Bo1g037840), HSP90 (Bo03352s010, Bo2g037660, Bo3g020890, Bo3g021230, Bo5g048850, Bo32505s010, Bo13090s010) and chaperone protein dnaJ (Bo8g102330, Bo1g017440, Bo7g113150, Bo9g149730, Bo1g134560, Bo00581s030, Bo5g131620, Bo6g112630, Bo1g005990, Bo6g115880, Bo5g139610), but no interactions between HSFs. These analyses seemingly suggest that 1) no master regulator, transcription co-activator or repressor in B. alboglabra, which was different from L. esculentum and A. thaliana; 2) HSFs in B. alboglabra functioned redundantly, as in the case of the Arabidopsis Hsf A1 subfamily [41, 44,45,46]; 3) HSF C member was a transcriptional repressor of HSP genes, which has not been shown in any other plant species. Importantly, it showed that the expression patterns of HSF A1e, HSF A4a (Bo7g106320), HSF A8, HSF B3, HSF B4 and HSFC1 were different in high- and low- GS varieties, indicating the different heat shock response between high- and low-GS B. alboglabra sprouts under heat stress.
HSP families, including HSP100/ClpB, HSP90, HSP70, HSP60 and sHSPs, regulated by HSFs, can be induced to respond to an array of physiological and environmental stresses. HSP protect cells from degradation, oxidative stress, hypoxia and thermal stress by catalyzing the refolding of damaged or denatured proteins to prevent protein aggregation and promote protein disaggregation [47, 48]. The action mode of HSFs and HSPs in tomato has been widely researched. The HSF and HSP networks in plants were controlled at the transcriptional level by cooperation of distinct HSF members and by interaction with chaperones. In response to heat shock stress, HSFs not only regulated HSP genes expression, but also were concomitantly subjected to cell-type-dependent feedback regulation through factor-specific physical and functional interactions with chaperones belonging to HSP90, HSP70 and small HSP families . Hahn et al. (2011) found that HSP70 and HSP90 regulated HSF function by direct interactions in tomato, HSP70 repressed the activity of HSF A1, HSP90 affected the abundance of HSF A2 and HSF B1 by modulating HSF A2 transcript degradation involved in regulation of the timing of HSF A2 synthesis .
In this study, HSP100/ClpB, HSP90, HSP70 and sHSPs were differentially expressed and the expression levels of most HSPs in high-GS variety were different from those in low-GS variety, indicating the different actions of HSPs when the two lines exposure to heat stress. It is worth to note that HSP90 was expressed prior to some HSFs, indicating there were some HSF-independent pathways for activated expression of HSP in response to heat stress.
Taken together, this study confirmed the presence of the following transcription features: 1) down-regulation of pathways may either serve as a passive adaptation or an active resistance/protection measure to save energy or reduce consumption, such as, the down-regulated GS pathway under heat shock stress; and 2) the up-regulation of a subset of specific heat-response genes could prompt or rapidly induce a protective process upon heat shock.
Four kinds of Chinese kale (B. alboglabra) including JL08, JL09, Cuibao and Shunbao were used in this study. Seeds were germinated in sterile petri dishes with sterile filter paper and 15 mL of sterile water after disinfected in 0.7% sodium hypochlorite for 30 min, then drained and washed seven times with distilled water. The germinated seedlings were transferred to another sterile petri dish with sterile filter paper. One hundred seedlings and 10 mL sterile water was added per dish every three days after day one. All plants were grown in 16 h day/8 h night at 25 °C for 5 days and then 42 °C for 2 days, respectively. The lines JL08 were of high-GS contents, named as HG, and JL09 were of low-GS contents, named as LG. After treated with high temperature (42 °C), HG and LG were abbreviated as HGT and LGT respectively. JL08 and JL09 were cultivated in the greenhouse of College of Horticulture in Fujian Agricultural and Forestry University and used for transcriptome sequencing as they have been self-crossed for five generations as described in Guo et al. (2016) . To validate the relationship of heat-resistance and GS content in Chinese kale sprouts, another two varieties with different GS content and heat resistance ability, Cuibao and Shunbao (Produced by the Musashino Seed Co., LTD. Tokyo, Japan), were used for GS profiling, content assays and heat resistant related gene assays. Fujian Agriculture and Forestry University has approved the study. The land accessed is not privately owned and protected. No protected species were sampled. No specific permissions were required for these locations/activities.
Measurement of GS content
GS were extracted and analyzed as previously described [51, 52]. Samples (500 mg) were boiled in 3 mL water for 10 min. After transferring the supernatant to a new tube, the residues were washed with 3 mL distilled deionized water (ddH2O), and the combined aqueous extract was applied to a DEAE-Sephadex A-25 (30 mg) column (pyridine acetate form) (Sigma, St. Louis, MO, USA). The column was washed three times with 20 mM pyridine acetate and twice with ddH2O. The GS were converted into their desulpho analogues by overnight treatment with 100 μL of 0.1% (1.4 units) aryl sulphatase (Sigma, St. Louis, MO, USA) added into the column, and the desulpho GSs were eluted with 2 × 0.5 mL water. HPLC analysis was performed using an HPLC system consisting of an Agilent HPLC series chromatograph (Agilent Technologies). A Hypersil C18 column (5 μm particle size, 4.6 mm × 250 mm, Elite Analytical Instruments Co., Ltd., Dalian, China) was used with a mobile phase of acetonitrile and water at a flow rate of 1.0 mL/min. The procedure employed isocratic elution with 1.5% acetonitrile for the first 5 min; a linear gradient to 20% acetonitrile over the next 15 min followed by isocratic elution with 20% acetonitrile for the final 10 min. A 40 μL sample was injected into the column by an autosampler. Absorbance was detected at 226 nm. Ortho-nitrophenyl-β-d-galactopyranoside (Sigma, St. Louis, MO, USA) was used as an internal standard for HPLC analysis. Desulphoglucosinolates were identified by comparison of retention time and quantified by peak area. The GS concentration was expressed as μmol/g fresh weight (fw) of sprouts.
RNA extraction, RNA-Seq library construction and sequencing
Total RNA was extracted from 7-day old Chinese kale sprouts of lines JL-08 and JL-09 using the TRIzol reagent (Invitrogen, Carlsbad, CA) according to the manufacturer’s instructions. The purity and concentration of the total RNA was detected by NanoDrop 1000 spectrophotometer (Thermo Fisher Scientific, Wilmington, DE, USA) and Qubit® 2.0 Flurometer (Life Technologies, CA, USA), respectively. The integrity of the RNA samples was assessed using the RNA Nano 6000 Assay Kit of the Agilent Bioanalyzer 2100 system (Agilent Technologies, CA, USA). The purified RNA was dissolved in RNase-free water and stored at − 80 °C.
The construction of RNA-seq library was performed following procedures similar to previous study [2, 53]. Briefly, after enrichment and purification with oligo (dT)-rich magnetic beads, the mRNA was cleaved into short fragments. Then first- and second-strand cDNA were synthesized using the mRNA fragments as templates. The cDNA was purified by AMPure XP beads and resolved with elution buffer for end reparation and single nucleotide adenine addition, and the short fragments were connected with adapters. Then the suitable fragments were selected as templates for PCR amplification. Finally, RNA-seq libraries were sequenced using an Illumina HiSeq™ 2500 at Biomarker Technologies Corporation (Beijing, China) and submitted to the SRA database (PRJNA427496).
RNA-Seq reads mapping and transcript assembly
After sequencing, raw reads containing sequencing adapters and low-quality reads were removed. Brassica oleracea genome and gene model annotation files were directly downloaded from the genome website (http://plants.ensembl.org/Brassica_oleracea/Info/Index). Index of the reference genome was constructed using Bowtie , and paired-end clean reads were aligned to the reference genome using Top Hat 2 (http://ccb.jhu.edu/software/tophat/index.shtml)  with all parameters set to their default values. Then the Sequence Alignment Map (http://samtools.sourceforge.net) (Li et al., 2009) files were generated by Top Hat2 and subsequently transcripts were assembled by Cufflinks (http://cufflinks.cbcb.umd.edu) .
Transcriptome annotation, novel transcript prediction, and alternative splicing analysis
To annotate the Chinese kale genes, the sequences were compared to Gene Ontology (GO) , Kyoto Encyclopedia of Genes and Genomes (KEGG) , Swiss-Prot  and the NCBI non-redundant protein sequence databases by using the BLASTx algorithm with an e-value threshold of 10− 5 . Prior to differential gene expression analysis for each sequenced library, the read counts were adjusted using the Edge R program package through one scaling normalized factor. Sequence splicing and redundancy was determined by DESeq . The DEGs were identified from the Fragments Per Kilobase of transcript per Million fragments mapped (FPKM) values using the DEGseq R package . The resulting P values were adjusted using the Benjamini and Hochberg’s approach for controlling the false discovery rate (FDR). Genes with an adjusted P-value < 0.01 found by DESeq were assigned as differentially expressed, and the FDR < 0.01 and |log2 (fold change)| ≥1 was set as the threshold for significantly differential expression.
The Cufflinks v2.1.1 Reference Annotation Based Transcript assembly method was used to construct and identify both known and novel transcripts from Top Hat alignment results. Alternative splicing (AS) events were classified into 6 basic types by Asprofile v1.0. The number of AS events in each sample was estimated.
Hierarchical clustering, GO, and KEGG analysis
The overall gene expression patterns were separately analyzed based on the expression levels of the DEGs in lines JL-08 and JL-09, respectively. GO enrichment analysis (http://www.geneontology.org) of the DEGs was implemented by the GOseq R packages , and the corrected P-value < 0.05 was set as the threshold for GO analysis. Enrichment analysis of DEGs in the KEGG (http://www.genome.jp/kegg/) biological pathway was conducted, and we used KOBAS software  to test the statistical enrichment of DEGs in KEGG pathways.
Gene expression validation
Real-time quantitative RT-PCR (qRT-PCR) analyses were performed with independent samples from lines JL-08 and JL-09 grown at 25 °C and 42 °C with the same conditions as those used for RNA-Seq analysis. The cDNAs were reversely transcribed using PrimeScript™ RT reagent Kit with gDNA Eraser (Takara, Otsu, JP), and qRT-PCR analyses were performed on an CFX96 Touch™ Real-Time PCR Detection System (Bio-Rad Laboratories, Inc.) using SYBR® Premix Ex Taq™ II (Tli RNaseH Plus) (Takara, Otsu, JP). The Actin gene from B. alboglabra was used as an internal control to normalize the expression data. Gene expression was evaluated by the 2-ΔΔCt method . The gene-specific primers are listed in Additional file 1: Table S1.
Biological coefficient of variation
distilled deionized water
Differentially expressed genes
False discovery rate
Fragments per kilobase of transcript per million fragments mapped
High-GS line treated with heat stress
Heat shock transcription factors
Heat shock proteins
Kyoto encyclopedia of genes and genomes
Low-GS line treated with heat stress
Real-time quantitative RT-PCR
Qian H, Liu T, Deng M, Miao H, Cai C, Shen W, Wang Q. Effects of light quality on main health-promoting compounds and antioxidant capacity of Chinese kale sprouts. Food Chem. 2016;196:1232–8.
Guo R, Huang Z, Deng Y, Chen X, Xu X, Lai Z. Comparative transcriptome analyses reveal a special Glucosinolate metabolism mechanism in Brassica alboglabra sprouts. Front Plant Sci. 2016;7:1497.
Guo R, Deng Y, Huang Z, Chen X, Xu X, Lai Z. Identification of miRNAs affecting the establishment of Brassica Alboglabra seedling. Front Plant Sci. 2016;7:1760.
Srikanth S, Choong TW, Yan A, He J, Chen Z. An efficient method for adventitious root induction from stem segments of Brassica species. Front Plant Sci. 2016;7:943.
Halkier BA, Gershenzon J. Biology and biochemistry of glucosinolates. Annu Rev Plant Biol. 2006;57:303–33.
Kim JI, Dolan WL, Anderson NA, Chapple C. Indole Glucosinolate biosynthesis limits Phenylpropanoid accumulation in Arabidopsis thaliana. Plant Cell. 2015;27(5):1529–46.
Mach J. Metabolic crosstalk: interactions between the Phenylpropanoid and Glucosinolate pathways in Arabidopsis. Plant Cell. 2015;27(5):1367.
Sonderby IE, Geu-Flores F, Halkier BA. Biosynthesis of glucosinolates - gene discovery and beyond. Trends Plant Sci. 2010;15(5):283–90.
Gao J, Yu X, Ma F, Li J. RNA-Seq analysis of transcriptome and Glucosinolate metabolism in seeds and sprouts of broccoli (Brassica oleracea var. italic). PLoS One. 2014;9(2):e88804.
Guo L, Yang R, Zhou Y, Gu Z. Heat and hypoxia stresses enhance the accumulation of aliphatic glucosinolates and sulforaphane in broccoli sprouts. Eur Food Res Technol. 2016;242(1):107–16.
Seo M-S, Jin M, Chun J-H, Kim S-J, Park B-S, Shon S-H, Kim JS. Functional analysis of three BrMYB28 transcription factors controlling the biosynthesis of glucosinolates in Brassica rapa. Plant Mol Biol. 2016;90(4–5):503–16.
Zhang J, Liu Z, Liang J, Wu J, Cheng F, Wang X. Three genes encoding AOP2, a protein involved in aliphatic glucosinolate biosynthesis, are differentially expressed in Brassica rapa. J Exp Bot. 2015;66(20):6205–18.
Charron CS, Sams CE. Glucosinolate content and myrosinase activity in rapid-cycling Brassica oleracea grown in a controlled environment. J Am Soc Hortic Sci. 2004;129(3):321–30.
Ludwig-Muller J, Krishna P, Forreiter C. A glucosinolate mutant of Arabidopsis is thermosensitive and defective in cytosolic Hsp90 expression after heat stress. Plant Physiol. 2000;123(3):949–58.
Hara M, Harazaki A, Tabata K. Administration of isothiocyanates enhances heat tolerance in Arabidopsis thaliana. Plant Growth Regul. 2013;69(1):71–7.
Jin G-H, Gho H-J, Jung K-H. A systematic view of rice heat shock transcription factor family using phylogenomic analysis. J Plant Physiol. 2013;170(3):321–9.
Kotak S, Larkindale J, Lee U, von Koskull-Doering P, Vierling E, Scharf K-D. Complexity of the heat stress response in plants. Curr Opin Plant Biol. 2007;10(3):310–6.
Hartl FU, Hayer-Hartl M. Protein folding - molecular chaperones in the cytosol: from nascent chain to folded protein. Science. 2002;295(5561):1852–8.
Young MD, Wakefield MJ, Smyth GK, Oshlack A. Gene ontology analysis for RNA-seq: accounting for selection bias. Genome Biol. 2010;11(2):R14.
Kissen R, Overby A, Winge P, Bones AM. Allyl-isothiocyanate treatment induces a complex transcriptional reprogramming including heat stress, oxidative stress and plant defence responses in Arabidopsis thaliana. BMC Genomics. 2016;17:740.
Yu E, Fan C, Yang Q, Li X, Wan B, Dong Y, Wang X, Zhou Y. Identification of heat responsive genes in Brassica napus siliques at the seed-filling stage through transcriptional profiling. PLoS One. 2014;9(7):e101914.
Nover L, Bharti K, Doring P, Mishra SK, Ganguli A, Scharf KD. Arabidopsis and the heat stress transcription factor world: how many heat stress transcription factors do we need? Cell Stress Chaperones. 2001;6(3):177–89.
Kim KK, Kim R, Kim SH. Crystal structure of a small heat-shock protein. Nature. 1998;394(6693):595–9.
Wahid A, Gelani S, Ashraf M, Foolad MR. Heat tolerance in plants: an overview. Environ Exp Bot. 2007;61(3):199–223.
Chauhan H, Khurana N, Agarwal P, Khurana P. Heat shock factors in rice (Oryza sativa L.): genome-wide expression analysis during reproductive development and abiotic stress. Mol Gen Genomics. 2011;286(2):171–87.
Pereira FMV, Rosa E, Fahey JW, Stephenson KK, Carvalho R, Aires A. Influence of temperature and ontogeny on the levels of glucosinolates in broccoli (Brassica oleracea var. italica) sprouts and their effect on the induction of mammalian phase 2 enzymes. J Agric Food Chem. 2002;50(21):6239–44.
Chevolleau S, Gasc N, Rollin P, Tulliez J. Enzymatic, chemical, and thermal breakdown of 3H-labeled glucobrassicin, the parent indole glucosinolate. J Agric Food Chem. 1997;45(11):4290–6.
Oerlemans K, Barrett DM, Suades CB, Verkerk R, Dekker M. Thermal degradation of glucosinolates in red cabbage. Food Chem. 2006;95(1):19–29.
Verkerk R. Evaluation of glucosinolate levels throughout the production chain of Brassica vegetables: Towards a novel predictive modelling approach. PhD. Dissertation. Wageningen, The Netherlands: Wageningen University; 2002.
Verkerk R, Dekker M. Glucosinolates and myrosinase activity in red cabbage (Brassica oleracea L. var. Capitata f. Rubra DC.) after various microwave treatments. J Agric Food Chem. 2004;52(24):7318–23.
Hennig K, Verkerk R, Bonnema G, Dekker M. Rapid estimation of Glucosinolate thermal degradation rate constants in leaves of Chinese kale and broccoli (Brassica oleracea) in two seasons. J Agric Food Chem. 2012;60(32):7859–65.
Hanschen FS, Rohn S, Mewis I, Schreiner M, Kroh LW. Influence of the chemical structure on the thermal degradation of the glucosinolates in broccoli sprouts. Food Chem. 2012;130(1):1–8.
Kotak S, Port M, Ganguli A, Bicker F, von Koskull-Doring P. Characterization of C-terminal domains of Arabidopsis heat stress transcription factors (Hsfs) and identification of a new signature combination of plant class a Hsfs with AHA and NES motifs essential for activator function and intracellular localization. Plant J. 2004;39(1):98–112.
Ikeda M, Ohme-Takagi M. A novel Group of Transcriptional Repressors in Arabidopsis. Plant Cell Physiol. 2009;50(5):970–5.
Zhu X, Thalor SK, Takahashi Y, Berberich T, Kusano T. An inhibitory effect of the sequence-conserved upstream open-reading frame on the translation of the main open-reading frame of HsfB1 transcripts in Arabidopsis. Plant Cell & Environment. 2012;35(11):2014–30.
Ikeda M, Mitsuda N, Ohme-Takagi M. Arabidopsis HsfB1 and HsfB2b act as repressors of the expression of heat-inducible Hsfs but positively regulate the acquired thermotolerance. Plant Physiol. 2011;157(3):1243–54.
Xue G-P, Sadat S, Drenth J, McIntyre CL. The heat shock factor family from Triticum aestivum in response to heat and other major abiotic stresses and their role in regulation of heat shock protein genes. J Exp Bot. 2014;65(2):539–57.
Mishra SK, Tripp J, Winkelhaus S, Tschiersch B, Theres K, Nover L, Scharf KD. In the complex family of heat stress transcription factors, HSfA1 has a unique role as master regulator of thermotolerance in tomato. Genes Dev. 2002;16(12):1555–67.
Bharti K, von Koskull-Doring P, Bharti S, Kumar P, Tintschl-Korbitzer A, Treuter E, Nover L. Tomato heat stress transcription factor HsfB1 represents a novel type of general transcription coactivator with a histone-like motif interacting with the plant CREB binding protein ortholog HAC1. Plant Cell. 2004;16(6):1521–35.
Scharf K-D, Berberich T, Ebersberger I, Nover L. The plant heat stress transcription factor (Hsf) family: structure, function and evolution. Biochimica Et Biophysica Acta-Gene Regulatory Mechanisms. 2012;1819(2):104–19.
Liu H-C, Liao H-T, Charng Y-Y. The role of class A1 heat shock factors (HSFA1s) in response to heat and other stresses in Arabidopsis. Plant Cell & Environment. 2011;34(5):738–51.
H-c L, Y-y C. Acquired thermotolerance independent of heat shock factor A1 (HsfA1), the master regulator of the heat stress response. Plant Signal Behav. 2012;7(5):547–50.
Czarnecka-Verner E, Pan SQ, Salem T, Gurley WB. Plant class BHSFs inhibit transcription and exhibit affinity for TFIIB and TBP. Plant Mol Biol. 2004;56(1):57–75.
Lohmann C, Eggers-Schumacher G, Wunderlich M, Schoffl F. Two different heat shock transcription factors regulate immediate early expression of stress genes in Arabidopsis (vol 271, pg 11, 2004). Mol Gen Genomics. 2004;271(3):376.
Nishizawa-Yokoi A, Nosaka R, Hayashi H, Tainaka H, Maruta T, Tamoi M, Ikeda M, Ohme-Takagi M, Yoshimura K, Yabuta Y, et al. HsfA1d and HsfA1e involved in the transcriptional regulation of HsfA2 function as key regulators for the Hsf signaling network in response to environmental stress. Plant Cell Physiol. 2011;52(5):933–45.
Yoshida T, Ohama N, Nakajima J, Kidokoro S, Mizoi J, Nakashima K, Maruyama K, Kim J-M, Seki M, Todaka D, et al. Arabidopsis HsfA1 transcription factors function as the main positive regulators in heat shock-responsive gene expression. Mol Gen Genomics. 2011;286(5–6):321–32.
Parsell DA, Lindquist S. The function of heat-shock proteins in stress tolerance: degradation and reactivation of damaged proteins. Annu Rev Genet. 1993;27:437–96.
Lindquist S, Craig EA. The heat-shock proteins. Annu Rev Genet. 1988;22:631–77.
Fragkostefanakis S, Roeth S, Schleiff E, Scharf K-D. Prospects of engineering thermotolerance in crops through modulation of heat stress transcription factor and heat shock protein networks. Plant Cell & Environment. 2015;38(9):1881–95.
Hahn A, Bublak D, Schleiff E, Scharf K-D. Crosstalk between Hsp90 and Hsp70 chaperones and heat stress transcription factors in tomato. Plant Cell. 2011;23(2):741–55.
Guo R, Yuan G, Wang Q. Sucrose enhances the accumulation of anthocyanins and glucosinolates in broccoli sprouts. Food Chem. 2011;129(3):1080–7.
Guo R, Yuan G, Wang Q. Effect of sucrose and mannitol on the accumulation of health-promoting compounds and the activity of metabolic enzymes in broccoli sprouts. Sci Hortic. 2011;128(3):159–65.
Lou Q, Liu Y, Qi Y, Jiao S, Tian F, Jiang L, Wang Y. Transcriptome sequencing and metabolite analysis reveals the role of delphinidin metabolism in flower colour in grape hyacinth. J Exp Bot. 2014;65(12):3157–64.
Langmead B, Salzberg SL. Fast gapped-read alignment with bowtie 2. Nat Methods. 2012;9(4):357–U354.
Kim D, Pertea G, Trapnell C, Pimentel H, Kelley R, Salzberg SL. TopHat2: accurate alignment of transcriptomes in the presence of insertions, deletions and gene fusions. Genome Biol. 2013;14(4):R36.
Trapnell C, Williams BA, Pertea G, Mortazavi A, Kwan G, van Baren MJ, Salzberg SL, Wold BJ, Pachter L. Transcript assembly and quantification by RNA-Seq reveals unannotated transcripts and isoform switching during cell differentiation. Nat Biotechnol. 2010;28(5):511–U174.
Ashburner M, Ball CA, Blake JA, Botstein D, Butler H, Cherry JM, Davis AP, Dolinski K, Dwight SS, Eppig JT, et al. Gene ontology: tool for the unification of biology. The gene ontology consortium. Nat Genet. 2000;25(1):25–9.
Kanehisa M, Goto S, Kawashima S, Okuno Y, Hattori M. The KEGG resource for deciphering the genome. Nucleic Acids Res. 2004;32:D277–80.
Apweiler R, Bairoch A, Wu CH, Barker WC, Boeckmann B, Ferro S, Gasteiger E, Huang HZ, Lopez R, Magrane M, et al. UniProt: the universal protein knowledgebase. Nucleic Acids Res. 2004;32:D115–9.
Altschul SF, Madden TL, Schaffer AA, Zhang J, Zhang Z, Miller W, Lipman DJ, Gapped BLAST. PSI-BLAST: a new generation of protein database search programs. Nucleic Acids Res. 1997;25(17):3389–402.
Wang L, Feng Z, Wang X, Wang X, Zhang X. DEGseq: an R package for identifying differentially expressed genes from RNA-seq data. Bioinformatics. 2010;26(1):136–8.
Mao XZ, Cai T, Olyarchuk JG, Wei LP. Automated genome annotation and pathway identification using the KEGG Orthology (KO) as a controlled vocabulary. Bioinformatics. 2005;21(19):3787–93.
Livak KJ, Schmittgen TD. Analysis of relative gene expression data using real-time quantitative PCR and the 2(−Delta Delta C(T)) method. Methods. 2001;25(4):402–8.
This work was supported by the Outstanding Youth Training Program in Colleges and Universities of Fujian Province (K8015060A) for the study and collection of data, Special Financial Grant (No. 2017 T100464) for analysis of data, National Natural Science Foundation of China (No. 31401859, 31772310) for interpretation of data, and Horticulture Postdoctoral Funding (No. 132300155) in writing the manuscript.
Availability of data and materials
The materials of this study were provided by Joint FAFU-Dalhousie Lab at Fujian Agriculture and Forestry University. Correspondence and requests for materials should be addressed to Gefu Wang-Pruski (Gefu.Wang-Pruski@Dal.Ca). The raw sequencing data have been submitted to the NCBI SRA database (PRJNA427496).
Ethics approval and consent to participate
Consent for publication
We declare that the research was conducted in the absence of any commercial or financial relationships that could be a potential conflict of interest.
Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.
Table S1. Primers used in the verification of differential expressed genes in high- and low-GS Chinese kale sprouts under heat stress. Figure S1. (A) KEGG annotation of DEGs in HG vs HGT. (B) KEGG annotation of DEGs in LG vs LGT. Yellow column shows cellular processes, blue column shows organismal systems, purple column shows environmental information processing, green column shows metabolism, and pink column shows genetic information processing. Figure S2. (A) KEGG enrichment of DEGs in HG vs HGT. (B) KEGG enrichment of DEGs in LG vs LGT. The y-axis corresponds to different Q-value of different KEGG pathway, and the x-axis shows the enrichment factor of different KEGG pathway. The different symbols represented for different KEGG pathways. (DOCX 767 kb)