Muhammad Ikram
Jingfang Chen1†
Yanshi Xia
Kadambot H. M. Siddique
Peiguo Guo- 1Guangdong Provincial Key Laboratory of Plant Adaptation and Molecular Design, International Crop Research Center for Stress Resistance, School of Life Sciences, Guangzhou University, Guangzhou, China
- 2The UWA Institute of Agriculture, UWA School of Agriculture & Environment, The University of Western Australia, Perth, WA, Australia
Flowering Chinese cabbage (Brassica campestris L. ssp. chinensis var. utilis Tsen et Lee, 2n=20, AA) is a vegetable species in southern parts of China that faces high temperatures in the summer and winter seasons. While heat stress adversely impacts plant productivity and survival, the underlying molecular and biochemical causes are poorly understood. This study investigated the gene expression profiles of heat-sensitive (HS) ‘3T-6’ and heat-tolerant (HT) ‘Youlu-501’ varieties of flowering Chinese cabbage in response to heat stress using RNA sequencing. Among the 37,958 genes expressed in leaves, 20,680 were differentially expressed genes (DEGs) at 1, 6, and 12 h, with 1,078 simultaneously expressed at all time points in both varieties. Hierarchical clustering analysis identified three clusters comprising 1,958, 556, and 591 down-regulated, up-regulated, and up- and/or down-regulated DEGs (3205 DEGs; 8.44%), which were significantly enriched in MAPK signaling, plant–pathogen interactions, plant hormone signal transduction, and brassinosteroid biosynthesis pathways and involved in stimulus, stress, growth, reproductive, and defense responses. Transcription factors, including MYB (12), NAC (13), WRKY (11), ERF (31), HSF (17), bHLH (16), and regulatory proteins such as PAL, CYP450, and photosystem II, played an essential role as effectors of homeostasis, kinases/phosphatases, and photosynthesis. Among 3205 DEGs, many previously reported genes underlying heat stress were also identified, e.g., BraWRKY25, BraHSP70, BraHSPB27, BraCYP71A23, BraPYL9, and BraA05g032350.3C. The genome-wide comparison of HS and HT provides a solid foundation for understanding the molecular mechanisms of heat tolerance in flowering Chinese cabbage.
Introduction
Heat stress is an environmental factor affecting crop vegetative growth, performance, and productivity worldwide, even causing plant death (Caers et al., 1985; Song et al., 2014a). According to a crop-based algorithm study, a 1°C rise in seasonal temperature causes 2.5–16% direct yield losses for key crops in tropical and subtropical regions (Lobell et al., 2008). Furthermore, greenhouse gas emissions have increased atmospheric temperatures by 0.3°C per decade (Jones et al., 1999), indicating that extreme temperature and climate change events will threaten food security by decreasing crop production (Battisti and Naylor, 2009; Rosenzweig et al., 2014). Flowering Chinese cabbage (Brassica campestris L. ssp. chinensis var. utilis Tsen et Lee) is a commercial vegetable crop mainly cultivated in southern parts of China (Chen et al., 2017). The demand for flowering Chinese cabbage is increasing due to its rich source of vitamin C, favourable taste, soluble fiber, phenolics, glucosinolates, and other nutritional value (Huang et al., 2017), but southern China’s high summer and autumn temperatures seriously impair crop quality and production (Fan et al., 2017). Heat stress adversely affects plant growth, phenological stages, grain filling, pollen viability, even the structural changes in tissues and cell organelles, loss of leaf water, cell membrane damage, and photosynthetic membranes (Yoshinaga et al., 2005). At a molecular level, heat stress also affects the synthesis of primary and secondary metabolites, antioxidant enzymes, and lipid peroxidation via the production of reactive oxygen species (ROS). However, plants have evolved various molecular and physiological mechanisms during domestication to resist heat stress (Wahid et al., 2007; Hasanuzzaman et al., 2013), such as accumulating different metabolites (antioxidants, osmoprotectants, etc.) and activation of signaling and metabolic pathways. Therefore, understanding the mechanism of signaling cascades and specific genes expressed in response to HT will be beneficial for developing stress-tolerant varieties.
Transcription factors (TFs), signal transduction components, and proteins related to metabolism are responsive to heat stress (Grover et al., 2013). Recently, heat-responsive siRNAs and miRNAs in flowering Chinese cabbage have been reported (Yu et al., 2012; Song et al., 2021), with most up-regulated under high temperature, and their target genes with differential gene expression involved in responses to temperature stimulus, cell membrane, signal transduction, and mitogen-activated protein kinase (MAPK) signaling pathways (Song et al., 2021). In Chinese cabbage, Song et al. (2016) used RNA-seq datasets of heat treatments to identify 14,329 DEGs and 9,687 novel LncRNAs, of which LncRNAs controlled 192 DEGs under heat stress. Transcriptome sequencing analysis of Chinese cabbage (B. rapa ssp. pekinensis) inbred lines revealed that heat stress affects many genes, including those linked with membrane leakage, heat-shock proteins (HSPs), and enzymes that regulate ROS homeostasis (Dong et al., 2015). Heat-shock transcriptional factors (HSFs), which stimulate the expression of HSPs and other stress proteins, play a critical role in modulating stress responses (Swindell et al., 2007; Dong et al., 2015; Liu et al., 2018; Yao et al., 2020). The HSF gene family has been identified in various species, including Arabidopsis (Nover et al., 2001), tomato (Yang et al., 2016), soybean (Li et al., 2014), potato (Dossa et al., 2016), and Chinese cabbage (Song et al., 2014b). HSF and HSP family members are involved in several stress response pathways, including heat, cold, osmotic, and salt, in Arabidopsis (Swindell et al., 2007), with 21 HSFs cloned and analyzed (Nover et al., 2001). In global transcription profiles, many HSFs/HSPs were up-regulated in Brassica napus siliques, with many other heat-responsive marker genes, such as ROF2, MBF1c, DREB2a, and Hsa32, involved in heat resistance in many plants (Yu et al., 2014). Another study found that genes involved in protein protection, biotic stress responses, oxidative stress, programmed cell death, and metabolism were differentially expressed during heat stress (Kotak et al., 2007). Only limited transcriptome data are available for flowering Chinese cabbage under heat stress relative to Chinese cabbage, maize, rice, and Arabidopsis (Song et al., 2014a; Chen et al., 2017; Fan et al., 2017; Huang et al., 2017; Ali et al., 2022).
In recent years, next-generation sequencing technologies have been developed and significantly reduced the cost and increased the efficiency of genome sequencing. Therefore, the entire genome of B. rapa (Chiifu-401-42) was sequenced using Illumina GA II sequencing and annotated (Wang et al., 2011), which provides the foundation to identify the heat-responsive genes in flowering Chinese cabbage. Using RNA-seq, researchers can identify the expressed genes, particularly those with low abundance. To date, microarrays and RNA-Seq projects have been undertaken in many species, including rice, Chinese cabbage, maize, cotton, Arabidopsis, and vegetables, to detect genes responsive to heat and cold stress (Grover et al., 2013; Nakashima et al., 2014; Ikram et al., 2020; Song et al., 2021; Ali et al., 2022). For example, in Arabidopsis, nearly 30% of DEGs associated with abiotic stresses were identified, and 2,409 genes were involved in salt, cold, and drought stresses (Kreps et al., 2002). In wheat, 2% of genes were found to be associated with cold stress (Winfield et al., 2010). In A. mongolicus, 9,309 and 23,419 up- and down-regulated genes were detected under cold stress (Pang et al., 2013). In previous studies, various genes have been functionally characterized for heat tolerance in different crop plants, such as ZmWRKY106 (Wang et al., 2018), PpEXP1 (Xu et al., 2014), hsp26 (Xue et al., 2010), WsSGTL1 (Mishra et al., 2013), and OsOr-R115H (Jung et al., 2021).
Our previous research identified expressed sequence tag-simple sequence repeat (EST-SSR) markers in flowering Chinese cabbage (Chen et al., 2017) and also detected significant SNPs/genes for plant height (Ikram et al., 2022a; Ikram et al., 2022b; Li et al., 2022) and wilt resistance (Lai et al., 2021) in tobacco. The present study used RNA-seq based on the Illumina HiSeq2000 platform to undertake a genome-wide analysis of gene expression in leaves of HT (Youlu-501) and HS (3T-6) flowering Chinese cabbage varieties under heat stress (0, 1, 6, and 12 h) to (1) identify DEGs in HS and HT genotypes after 1, 6, and 12 h of heat stress: (2) annotate these genes based on gene ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) enrichment analysis; (3) identify the expression pattern of genes using Hierarchical clustering analysis; (4) identify TFs associated with the DEGs and their role in heat stress. The findings of this study will provide new insights into the heat tolerance of flowering Chinese cabbage.
Materials and methods
Plant materials and physiological parameters
Two natural flowering Chinese cabbage (Brassica campestris L. ssp. chinensis var. utilis Tsen et Lee) varieties [3T-6 (heat-sensitive) and Youlu-501 (heat-tolerant)] were obtained from the Guangdong Academy of Agricultural Sciences, China. The seeds of both accessions were sterilized and then germinated in a growth chamber of our lab at the International Crop Research Center for Stress Resistance (113.37_E, 23.05_N), Guangzhou University, under the following condition: 28/22°C for 14/10 h (day/night) with 80% relative humidity. Fifteen-day-old seedlings at five-leaf stage were transferred into another chamber at 38/29°C (14/10 h day/night) for heat stress. After the heat-stress treatment, samples were collected at three time points from fully expanded upper leaves of heat-sensitive (HS-1h, HS-6h, and HS-12h) and heat-tolerant (HT-1h, HT-6h, and HT-12h) plants. For the control (CK), leaf samples were collected from HS (HS-CK) and HT (HT-CK) plants under normal conditions at 25°C. The harvested leaf samples were immediately frozen in liquid nitrogen and preserved at –80°C until RNA extraction. The experiment had five replicates. Peroxidase (POD), superoxide dismutase (SOD), and catalase (CAT) were assessed as described in previous studies (Kono, 1978; Aebi, 1984). The fresh plants were taken at different times to calculate the fresh weight in both varieties.
RNA extraction, cDNA library construction, and sequencing
For transcriptome assembly, total RNA was extracted from five biological replicates at CK, HS/HT-1h, HS/HT-6h, and HS/HT-12h using TRIzol reagent (Takara Bio, Ostu, Japan), following the manufacturer’s recommendations. The RNA was pooled into a single sample in equal amounts to find the broad gene library associated with heat resistance. DNase I (Takara Bio, Ostu, Japan) was used to digest genomic DNA, and RNA integrity and quantity were assessed using 1% agarose gel electrophoresis and a microplate spectrophotometer (BioTek Company, USA). cDNA libraries were prepared using RNA with high purity. In brief, mRNAs were purified using poly-T oligo-attached magnetic beads from total RNA. Following fragmentation, first-strand cDNA was synthesized using a random hexamer primer, followed by second-strand cDNA synthesized using DNA Polymerase I and RNase H. The paired-end Illumina sequencing required the purification and ligation of the double-stranded cDNAs to adaptors. After PCR amplification, library quality was verified on the Agilent 2100 Bioanalyzer system before sequencing the cDNA libraries using the Illumina HiSeq2000 system by the Beijing Genomics Institute (BGI) to generate 125 bp paired-end reads.
RNA-seq data analysis and annotation
Before library assembly, the raw reads were processed using SOAPnuke v1.5.2 (https://github.com/BGI-flexlab/SOAPnuke; Cock et al., 2009) to remove low-quality reads with >50% of bases with Q ≤ 20 or >10% unknown (N) bases and adapter sequences. The raw reads with quality check parameters were transformed to clean reads for further analysis. Finally, clean reads were mapped and aligned to B. rapa cv. Chiifu v3.0 reference genome (http://brassicadb.org/brad/; Wang et al., 2011) using HISAT2 v2.1.0 (Kim et al., 2015) and Bowtie2 v2.2.5 (Langmead and Salzberg, 2012). For functional annotation, all assembled gene sequences were aligned to databases, including reference genome (BRAD database), NT, NR, UniProt, Cluster of Orthologous Groups of proteins (COG), KEGG, and GO using Trinotate software (http://trinotate.github.io/) with E value ≤1e-5. The mapped reads were used to measure expression levels in fragments per kilobase of exon per million mapped fragments (FPKM) using RSEM software (Li and Dewey, 2011). Pearson’s correlation coefficients for the eight samples were estimated using the ‘cor’ function in R4.1.0 (http://www.r-project.org/) software based on FPKM values, and a correlation graph was drawn using the ‘ggplot2’ R package (https://cran.r-project.org/web/packages/ggplot2/index.html). The ‘prcomp’ function in R4.1.0 (http://www.r-project.org/) was used to perform principal component analysis (PCA).
Differential expression and cluster analyses
For HS and HT varieties, DEGs between control and heat-stress conditions were calculated using the DESeq2 R package (Love et al., 2014), a reliable and efficient package for detecting DEGs between diverse samples. The significant DEGs between heat stress and control conditions were identified at |log2(fold-change) | > 2 and false discovery rate (FDR) adjusted P < 0.001. Further, significant DEGs in HS or HT were exposed to hierarchical or K-mean clustering based on the hclust function in R4.1.0 (distance: euclidean; method: ward.d), with log2FC subjected as the input. Finally, DEGs were analyzed for soft clustering based on fuzzy c-mean using the R4.1.0 package ‘Mfuzz’ (Kumar and Futschik, 2007).
Gene ontology and KEGG pathway enrichment analysis of differentially expressed genes
GO and KEGG enrichment analysis was undertaken to reveal the functional annotation and pathways associated with up- and down-regulated DEGs in each cluster or group. An online tool (http://www.geneontology.org/) was used to conduct GO enrichment analysis at FDR adjusted P ≤ 0.05, classifying the DEG functions into three categories: cellular component, biological process, and molecular function. Similarly, KEGG pathway enrichment analysis was performed using the KEGG database (https://www.genome.jp/kegg/pathway.html) with significant criteria at an adjusted P ≤ 0.05 (Kanehisa and Goto, 2000).
Validation of RNA-seq result by real-time quantitative PCR
Ten genes were chosen randomly to confirm the reliability of RNA-seq results using RT-qPCR. The RNA was the same as that used for transcriptomic sequencing, with the relative expression calculated using ABIPrism7000 RT-qPCR platform (Applied Biosystems, USA). The PCR reaction volume was 10 mL with 1 μg cDNA, 5 mL SYBR Premix Ex Taq II (Takara, China), and 200 nM primers. Three technical and biological replicates for each sample were measured with the following protocol: 95°C for 5 min, followed by 40 cycles of 95°C for 15 s, 60°C for 60 s, and until 65°C to 97°C with a ramp rate of 0.02°C s–1 for dissociation curve analysis. BraActin was used as housekeeping gene to normalize the data. The relative expression level of all selected genes at each time point was calculated using the 2−ΔΔCT method. The in-house R script was used to conduct a student t-test at P ≤ 0.05 to find the significant differences between time points.
Results
Plant phenotypic and enzymatic activities under heat stress
Heat stress significantly altered plant growth and the reproduction system. At the initial stage of stress, plant leaves exhibited moderate to high wilting symptoms. Twelve hours of heat stress significantly decreased whole plant weight (g) in both species relative to the control (Figure 1A). The heat treatment significantly increased catalase (CAT), superoxide dismutase (SOD), and peroxidase (POD) activities in plant leaves at all time points, relative to the control, and more so in Youlu-501 than 3T-6 (Figure 1A).
Figure 1 Phenotypic effect and overview of RNA sequencing data based on heat-sensitive and heat-tolerant flowering Chinese cabbage varieties. (A) Plant weight and enzymatic activities between control and stress conditions; (B) Number of transcripts and their sequence size; (C) Correlation analysis between heat-sensitive and heat-tolerant varieties at different time points; (D) log10 transformation of FPKM values at 0, 1, 6, and 12 h; (E) Distribution of genes based on FPKM values. The different small letters indicate the significant differences between the time points.
Overview of transcriptomic sequencing data
Libraries were constructed and sequenced on the IlluminaHiseq2000 deep sequencing platform using the PE125 protocol (Table S1). Each library generated 21.94–24.90 million raw reads. After removing joint contamination, low-quality reads, and reads with unknown base (N) content, 21.86–24.81 million clean reads (99.42–99.963%) were collected (Table S1). The clean Q30 base rates of the eight samples ranged from 89.71–91.04%. Of the clean reads, 88.59–92.06% total and 73.06–76.72% unique reads were mapped to the B. rapa reference genome (Table S1), with 46,250 predicted B. rapa genes (Figure 1B). The PCA revealed an overall 86.80% (PC1 = 73.40% and PC2 = 13.40%) variation in gene expression datasets, with the heat-stress conditions (HS-1h/HR-1h, HS-6h/HR-6h, and HS-12/HR-12) clustered or grouped nearby but separated from their corresponding control HS-CK/HR-CK (Figure S1). Pearson’s correlation coefficient analysis based on FPKM values of each sample revealed a significant positive correlation between different time points; for example, HS-6h and HT-6h highly correlated with other time points compared to HS-12h and HR-12h (Figure 1C). The HS and HT varieties significantly differed based on the log10 transformation of FPKM (Figure 1D), with the FPKM values of most genes >10 (Figure 1E). Sixteen genes had FPKM values >2,000 in both varieties at all time points.
Differential expression profiling of two varieties under heat stress
DESeq2 software was used to investigate DEGs with FPKM values >1. During heat stress, 30,842, 32,260, 32,926, 31,540, 32,635, and 33,113 genes were regulated for HS-1h, HS-6h, HS-12h, HT-1h, HT-6h, and HT-12h, respectively (Figures 2A, B). The DEGs were selected using the following criteria of log2 fold-change |log2(foldchange)|≥2 and FDR adjusted P ≤ 0.001. Volcano plots identified the significant up- or down-regulated genes during heat stress at each time point (Figure 2A), with 10,153 (1,115 up- and 9,038 down-regulated), 7,879 (4,438 up- and 3,441 down-regulated), and 7,662 (5,349 up- and 2,313 down-regulated) DEGs in 3T-6 at HS-1h, HS-6h, and HS-12h compared to the control (Figures 2A–C) and 13,704 (882 up- and 12,822 down-regulated), 5,161 (2,206 up- and 2,955 down-regulated), and 7,438 (3,989 up- and 3,449 down-regulated) DEGs in Youlu-501 at HT-1h, HT-6h, and HT-12h compared to the control (Figures 2A–C).
Figure 2 Analysis of differentially expressed genes in flowering Chinese cabbage varieties in response to heat stress. (A) Volcano plots of all expressed genes in heat-sensitive (HS) and heat-tolerant (HT) varieties at different time points relative to CK, with log2FC values drawn against –log10 (FDR) adjusted P-values. Red and blue dots represent up- and down-regulated DEGs based on |log2(foldchange)|≥2 and FDR ≤ 0.001; gray dots represent non-significant DEGs; (B) The Sankey diagram represents the number of common and unique DEGs between each stress point; (C) Number of up and down-regulated DEGs identified under different heat-stress treatments relative to control HS and HT varieties; (D) Common and unique DEGs under different heat-stress treatments in HS variety; (E) Common and unique DEGs in HT under different heat-stress treatments; (F) Common DEGs in all heat-stress treatments in HS and HT.
Subsequently, Sankey and Venn diagrams demonstrated the number of common and uniquely expressed genes at each time point or between two consecutive time points (Figures 2B, D–F), with 5,443 DEGs unique in all samples, and 7,025, 3,626, 2,382, 1,126, and 1,078 DEGs shared between two, three, four, five, and six samples (Figure 2B). Additionally, 2,073 and 2,210 DEGs were common at all time points of HS and HT genotypes, respectively (Figures 2D, E), indicating higher number of genes were regulated in heat-tolerant variety than sensitive. Of these, 1,078 DEGs were co-expressed differently in the HS and HT genotypes (Figure 2F), with 237 and 747 up- and down-regulated, respectively, and 94 differentially regulated. After 6 and 12 h of heat treatment, 328/325 and 330/323 more DEGs were up-regulated in HS/HT. Finally, 3,205 genes were differentially expressed at 1, 6, and 12 h in 3T-6 or Youlu-501 (Figures 2B, D–F).
Functional annotation of DEGs responsive to heat stress
GO enrichment analyses of up- and down-regulated genes in HS and HT at 1, 6, and 12 h categorized the heat stress-related genes associated with key biological processes. Figures 3, S2–S5 show the top 25 GO terms in three categories. In HS, for HS-CK vs. HS-1h, 462 up-regulated genes were mainly involved in ‘response to stimulus (GO:0009628),’ ‘cellular process (GO:0009987),’ ‘primary metabolic process (GO:0044238),’ ‘response to abiotic stimulus (GO:0009628),’ ‘developmental process involved in reproduction (GO:0003006),’ and ‘response to heat (GO:0009408)’ (Figure 3A). In comparison, most of the 4,055 down-regulated genes were involved in ‘tissue development (GO:0009888),’ ‘auxin-activated signaling pathway (GO:0009734),’ ‘cell part (GO:0044464),’ ‘intracellular part (GO:0044424),’ ‘cellular response to abiotic stimulus (GO:0071214),’ ‘cell (GO:0005623),’ ‘intracellular (GO:0005622),’ and ‘cellular process (GO:0009987)’ (Figure 3A). For HS-CK vs. HS-6h and HS-CK vs. HS-12h, 2,046 and 2,275 up-regulated genes were significantly enriched in ‘anion binding (GO:0043168),’ ‘nucleotide binding (GO:0000166),’ ‘response to temperature stimulus (GO:0009266),’ ‘purine nucleotide binding (GO:0017076),’ ‘carbohydrate derivative binding (GO:0097367),’ ‘response to abiotic stimulus (GO:0009628),’ and ‘photosynthesis (GO:0015979)’ (Figure S2). In contrast, 1,568 and 902 down-regulated genes were mainly enriched in ‘metabolic process (GO:0008152),’ ‘detoxification (GO:0071722),’ ‘intracellular organelle (GO:0043229),’ ‘chloroplast (GO:0009507),’ and ‘RNA binding (GO:0003723)’ (Figure S2). In HT, for HT-CK vs. HT-1h, 418 up-regulated genes were involved in ‘binding (GO:0005488)’ ‘cellular process (GO:0009987),’ and ‘primary metabolic process (GO:0044238)’, and 6,626 down-regulated genes were enriched in ‘leaf development (GO:0048366),’ ‘response to stimulus (GO:0009628),’ ‘chloroplast outer membrane (GO:0009707),’ and ‘biosynthetic process (GO:0009058)’ (Figure 3B). Similarly, the highest number of DEGs for HT-6h and HT-12 were involved in ‘response to stimulus (GO:0009628),’ ‘chloroplast (GO:0009507),’ ‘cell periphery (GO:0071944),’ ‘auxin homeostasis (GO:0010252),’ ‘membrane-bounded organelle (GO:0043227),’ and ‘photosynthesis (GO:0015979)’ (Figure S3). The GO enrichment analysis indicated that HT had more enriched DEGs for each term than HS; for example, 1,410 HT and 415 HS genes were involved in ‘response to stimulus’ (Figures 3A, B, S2, S3).
Figure 3 GO and KEGG enrichment analysis of up- and down-regulated genes in heat-sensitive (3T-6) and heat-tolerant (Youlu-501) flowering Chinese cabbage. (A, B) The top 25 GO terms in three categories are listed at HS-1h vs. HS-CK and HT-1h vs. HT-CK at P ≤ 0.05; (C, D) The top 25 KEGG pathways and their corresponding number of genes HS-1h vs. HS-CK and HT-1h vs. HT-CK at P ≤ 0.05. The x-axis represents the enrichment factor, the y-axis represents GO/KEGG terms with the number of genes, and the blue to red color indicates the significance and number of genes.
The KEGG database was used to investigate the key genes involved in KEGG pathways, with 2,389 and 1,453 up-regulated and 3,390 and 4,286 down-regulated genes in HS and HT, respectively, enriched in 39 KEGG pathways (adjusted P ≤ 0.05), including ‘plant-pathogen interaction (ko04626),’ ‘oxidative phosphorylation (ko00190),’ ‘sulfur metabolism (ko00920),’ ‘steroid biosynthesis (ko00100),’ ‘nitrogen metabolism (ko00910),’ ‘plant hormone signal transduction (ko04075),’ ‘brassinosteroid biosynthesis (ko00905),’ ‘MAPK signaling pathway - plant (ko04016),’ ‘circadian rhythm - plant (ko04712),’ ‘carbon metabolism (ko01200),’ and ‘photosynthesis (ko00195)’ (Figures 3C, D, S4, S5). The HT genes were more involved in signal transduction and MAPK pathways than the HS genes (Figures 3C, D).
Cluster analysis of differentially expressed genes
We used gene expression changes (log2FC value) to perform hierarchical clustering for 3205 DEGs (Figures 2D–F) in response to heat stress among flowering Chinese cabbage varieties. The clusters with similar patterns of expression change were identified using the hclust method. As a result, we found three clusters of genes response to heat stress at three time points relative to control and a black line representing the mean changes in expression (Figures 4A–C and Table S2). Cluster 1 contained 1,958 DEGs that were down-regulated in both varieties at 1, 6, and 12 h, and the expression variations of cluster 1 genes in Youlu-501 were comparatively more repressed than in 3T-6 (Figures 4A, B; Table S2). Interestingly, cluster 2 included 691 DEGs, the expression changes were down-regulated at 1h in both varieties compared to the control and up-regulated at 6 and 12 h, but expression patterns of HT were relatively higher than HS at 6 and 12 h stress condition (Figures 4A, B; Table S2). A total of 556 DEGs were found in Cluster 3, and most of them were up-regulated by heat stress in both varieties. The difference in Cluster 3 was that the genes had moderately greater induced expression levels in HS than in HT during heat stress except for 12h (Figures 4A, B; Table S2).
Figure 4 Expression patterns of differentially expressed genes (DEGs) and their functional annotation in heat-sensitive (3T-6) and heat-tolerant (Youlu-501) flowering Chinese cabbage in response to heat stress. (A) Hierarchical clustering analysis indicates DEG expression patterns in 3T-6 and Youlu-501 after 1, 6, and 12 h of heat stress relative to their controls; (B) Cluster analysis shows that the lines reflect expression patterns for each DEG at heat-stress time points. Black lines in each cluster indicate mean changes in DEG expression; (C) Gene ontology enrichment analysis of DEGs in each cluster.
Further, we investigated the biological functions of these clusters using GO enrichment analysis (Figure 4C; Table S2). The cluster 1 DEGs were mainly involved in cellular component organization or biogenesis (91), biological regulation (143), cellular process (394), developmental process (46), metabolic process (379), and response to stimulus (165). The DEGs belonging to developmental process (24), growth (10), response to stimulus (54), localization (22), metabolic process (155), and multi-organism process (17) were also significantly enriched in cluster 2. In addition, DEGs in cluster 3, which were up-regulated in both varieties, were involved in biological regulation (41), carbon utilization (8), biogenesis (28), reproduction (7), and response to stimulus (59). The fold-change expression values were higher in HT Youlu-501 than in HS 3T-6 (Figure 4A; Table S2). In brief, the top 15 potential DEGs expressed in all stages were identified, and their homologous have been reported to play an essential role in heat stress (Table 1).
Table 1 List of top potential differentially expressed genes involved in heat tolerance in flowering Chinese cabbage.
Gene families associated with heat stress
Gene families, including HSPs, cytochrome P450 (CYP), photosystem II (PSII), MAPK, and phenylalanine ammonia-lyase (PAL), play a significant role in plant abiotic and biotic stresses. HSP genes were clustered into two groups: (1) two genes up-regulated at high-stress conditions; (2) initially up-regulated and then down-regulated under high stress (Figures 5A–F; Table S3). Forty-two cytochrome P450 (CYP450) genes were differentially expressed under heat stress, with 16 were up-regulated in both varieties at 1, 6, and 12 h, including BraCYP98A3, BraCYP82G1 BraCYP83A1, BraCYP83A1, BraCYP72A13, and BraCYP94B3 (Figure 5A; Table S3). Of 16 DEGs, three genes (BraCYPBC1, BraCYP97B3, and BraCYP77A4) were up-regulated at 1, 6, and 12 h, while the others were up-regulated after 6 and 12 h stress (Figure 5A). Furthermore, BraA01g002420.3C (BraCYP84A1-X1) was down-regulated under heat stress (Figure 5A), with similar results obtained using RT-qPCR (Figure 5B). Photosynthesis is mainly sensitive to heat stress, and photosystem II (PSII) is the stress-sensitive site with its oxygen-evolving complex. Three PSII genes [PSBQ3 (BraA02g017000.3C), PSBQ2 (BraA03g024820.3C), and Psb27 (BraA08g035330.3C)] were negatively expressed in HS (Figure 5D; Table S3). Another gene, PSBQ2 (BraA05g013730.3C), encoding PSII oxygen-evolving enhancer protein 2, was up-regulated in HS and HT varieties. Moreover, 5 MAPK and 2 PAL gene expressions significantly differed in HS and HT after 1, 6, and 12 h. PAL1 and PAL3 have a redundant role in flavonoid biosynthesis (Table S3). BraA09g032810.3C (MPK-10) gene was down-regulated, while other genes were up- and down-regulated at different time points (Figure 5C).
Figure 5 Heatmaps of up- and down-regulated genes in heat-sensitive (3T-6) and heat-tolerant (Youlu-501) Chinese flowering cabbage related to stress at different times and confirmation of two DEGs through RT-qPCR analysis. (A) Number of DEGs of cytochrome P450 gene family at four heat-stress time points; (B) BraA01g002420.3C (BraCYP84A1-X1 cytochrome P450 84A1 isoform X1) showing similar expression pattern in RNA-seq data; Heatmaps of MAPK (C), photosystem II (D), and HSPs (E) related genes in HS and HT varieties under heat stress; (F) Confirmation of BraHSP70-5 gene identified through RNA-seq using RT-qPCR analysis. Asterisks indicate significance between treatment and control using Student’s t-test at **P < 0.01 and ***P < 0.001.
Differential expression of transcription factors during heat stress
TF families of flowering Chinese cabbage were retrieved from the PlantTFDB database for enrichment analysis. We analyzed the DEGs related to heat-responsive TFs, such as MYBs, ARFs, HSFs, NACs, ERFs, and WRKYs families, which could be involved in regulating genes during heat stress (Figures 6A–I; Table S4), identifying 12 MYBs, 13 NACs, 11 WRKYs, 31 ERFs, 17 HSFs, six ARFs, and 16 bHLHs TFs (Table S4). Apart from partial bHLH TFs family members, heat stress significantly up-regulated most gene family members of MYBs, NACs, WRKYs, and ERFs. For example, BraA03g051480.3C (WRKY53) in HS at 1, 6, and 12 h, with similar results obtained using RT-qPCR (Figures 6H, I). WRKY18 was down-regulated at the initial stress but up-regulated under high-stress conditions in HT (Figure 6I). HSFs encoding HSPs, as direct transcriptional activators of genes regulated by heat stress, showed higher up-regulated expression levels in HT than HS: BraA02g043370.3C encodes heat stress transcription factor B-2a-like (HsfB-2a) and had higher expression levels at 1, 6, and 12 h in HS and HT (Figures 6F, G; Table S4). Similarly, three bZIP TFs (BZIPHY5, BZIPHY5-X1, and BZIP61) were differentially expressed after heat stress. Of 12 MYBs, BraA07g030110.3C (MYB1R1) and BraA03g044990.3C (MYB44) were up-regulated at 1, 6, and 12 h, with the remaining ten down-regulated (Figure 6F; Table S4). Among 13 NACs, nine genes, including five NACs (BraA02g002870.3C, BraA03g053670.3C, BraA07g004020.3C, BraA07g034350.3C, and BraA10g031880.3C) were up-regulated in HS and HT at 1, 6, and 12 h, while BraA06g027220.3C, BraA07g014500.3C, and BraA10g000880.3C were up-regulated under high stress. In addition, most of the auxin (ARFs) and ethylene-responsive transcription factor (ERFs) up-regulated at 12h; for example, e.g., BraA01g030800.3C, BraA06g041230.3C, BraA09g022610.3C, BraA09g022620.3C, and BraA02g031120.3C, were down-regulated at initial stress (1h) and up-regulated at 6 and 12 h in HS and HT (Figures 6A, C–I; Table S4). These transcription factors possibly involve different thermotolerance in the two genotypes (Table S4) and might play essential roles in flowering Chinese cabbage heat resistance during the reproductive stage.
Figure 6 Heatmaps of transcription factor families involved in abiotic stress and validation of DEGs in heat-sensitive (3T-6) and heat-tolerant (Youlu-501) flowering Chinese cabbage using RT-qPCR analysis. The up- and down-regulated genes of TF families: (A) ERFs, (B) NAC, (C) ARFs, (D) bHLHs, (E) HSFs, (F) MYBs, and (H) WRKYs for different high-temperature treatments; (G–I) Confirmation of HsfB-2a and WRKY53 using RT-qPCR. Asterisks indicate significance between treatment and control using Student’s t-test at *P < 0.01, **P < 0.01, and ***P < 0.001.
Validation of RNA-seq data by quantitative real-time PCR
To verify the results of RNA-seq, ten genes with diverse expression profiles, including up-regulated or down-regulated at all stages and/or regulated in HT or HS, were randomly selected for real-time qPCR to measure expression levels. Table S5 lists the primer sequences of selected genes. The relative expression level was calculated with 2−ΔΔCT using BraActin as a reference gene (Table S5). As a result, two TFs genes, WRKY53 and HsfB-2a, had similar expression patterns as those in RNA-seq (Figure 6G, I). In brief, all ten genes exhibited the same expression profile, with high correlation coefficients observed between RNA-seq data and RT-qPCR at HS-1h (r = 0.83, P = 0.0032), HS-12h (r = 0.86, P = 0.0013), HT-1h (r = 0.88, P = 0.0017), and HT-12 (r = 0.86, P = 0.0028), indicating the consistency of RNA sequencing data (Figure 7).
Figure 7 Validation of ten differentially expressed genes in heat-sensitive (3T-6) and heat-tolerant (Youlu-501) flowering Chinese cabbage based on RT-qPCR. The RNA-seq and RT-qPCR values in the form of log2 fold-change were plotted using linear regression analysis at (A) HS-1, (B) HS-12h, (C) HT-1h, and (D) HT-12h. R and P indicate the correlation coefficient and corresponding P-value.
Discussion
Flowering Chinese cabbage, domesticated from Chinese cabbage with flowering stalk, is a major food due to its high nutrient value (Wang et al., 2016). However, heat stress significantly affects the yield of flowering Chinese cabbage. Thus, developing new varieties with high heat resistance is critical for solving this dilemma. In previous studies, stress-associated genes have been reported and functionally annotated in different vegetable crops under biotic and abiotic stress conditions, including high temperature, salt stress, heavy metal stress, drought, nutrient deficiency, bacterial wilt, fungal stress, and UV-B radiation (Gill et al., 2016; Gill et al., 2017; Cai et al., 2020; Yao et al., 2020; Ikram et al., 2022a; Li et al., 2022). In this study, the comparative transcriptome analysis between the heat treatment and control in HS and HT accessions significantly identified 20,680 DEGs (|log2(foldchange)|≥2, FDR ≤ 0.001) (Figures 2, 3), and 3,205 putative candidate genes involved in response to heat stress in flowering Chinese cabbage (Figure 4; Table S2).
Heat stress affects flowering Chinese cabbage growth and productivity and produces ROS in plant cells, causing oxidative damage. In response to stress, a plant defense mechanism against oxidative damage activated antioxidant enzymes, including SOD, POD, and CAT, that efficiently scavenge ROS (Wahid et al., 2007; Hasanuzzaman et al., 2013). In this study, SOD, POD, and CAT activities significantly increased under stress conditions compared to the control, more so in the HT genotype than the HS genotype (Figure 1A), as reported elsewhere (Wang et al., 2016; Song et al., 2019). Wahid et al. (2007) also measured increasing SOD, POD, and CAT activities in response to stress. Using RNA sequencing, we obtained 21.94–24.90 million raw reads from HS and HT accessions at 0, 1, 6, and 12 h, with 73.06–76.72% of the clean reads mapped to the B. rapa reference genome (Table S1) and higher from Chinese cabbage leaf RNA sequencing (Zhao et al., 2013) and flowering Chinese cabbage stalk RNA sequencing (Huang et al., 2017). We compared the expression data of leaves at 1, 6, and 12 h with the control, identifying 10,153, 7,879, 7,662, 13,704, 5,161, and 7,438 DEGs in HS-1h, HS-6h, HS-12h, HT-1h, HT-6h, and HT-12h, respectively (Figures 2A–C). Further, these genes were involved in cell division, cellular process, response to stimulus, primary metabolic process, response to heat, the developmental process involved in reproduction, photosynthesis, detoxification, RNA binding, and auxin homeostasis (Figures 3, S2, S3) in line with previous studies (Wang et al., 2016; Cai et al., 2020; Helal et al., 2021; Lai et al., 2021). In addition, most DEGs were enriched in the MAPK signaling pathway – plant, brassinosteroid biosynthesis, plant–pathogen interaction, plant hormone signal transduction, and circadian rhythm (Figures 3, S4, S5); these results are similar to previous reports (Song et al., 2014a; Song et al., 2021).
Hierarchical clustering analysis identified three clusters containing 1,958, 591, and 556 down-regulated, up- and down-regulated, and up-regulated, respectively (Figure 4; Table S2), involved in response to stimulus, biological regulation, oxidation, reproduction, cellular process, developmental process, and carbon utilization. These findings indicate that heat stress caused oxidative damage, leading to the up-regulation of genes to increase tolerance; for example, BraA07g024750.3C and BraA10g001200.3C encoded abscisic acid receptor PYL9-like and gibberellin 2-beta-dioxygenase 6, respectively (Table 1), and their homologous genes have been functionally characterized to mediate salt stress tolerance in Arabidopsis (Hichri et al., 2016) and cassava (Chang et al., 2020). BraA05g032350.3C encodes bidirectional sugar transporter SWEET2-like and plays a role in the organic substance biosynthetic process (Table 1). Moreover, cluster 3 DEGs were involved in sugar metabolism and proline/polyamine metabolism (Table S2), with the accumulation of sugar and proline detected as osmolytes, enhancing heat stress resistance (Wahid et al., 2007; Kumazaki and Suzuki, 2019; Song et al., 2019). The Sid2-1 mutant was reported in Arabidopsis, significantly increasing soluble sugars to improve drought and heat stress (Kumazaki and Suzuki, 2019). Similarly, heat stress increased proline content in leaves of soybean (Amirjani, 2010), tomato (Rivero et al., 2004), chickpea (Chakraborty and Tongden, 2005), and Arabidopsis (Zhang and Huang, 2013; Gurrieri et al., 2020; Hoermiller et al., 2022). The cluster 2 DEGs (Figure 4; Table S2) are associated with late embryogenesis abundant protein (Magwanga et al., 2018) and ubiquitin (Qian et al., 2020a), which protect sub-cellular and cellular structures from dehydration and oxidative forces. Therefore, the transcriptomic results were similar to enzymatic activities, and the DEGs related to different proteins in flowering Chinese cabbage might be necessary for heat stress tolerance.
In this study, we identified 42 CYPs DEGs involved in primary and secondary metabolism, including BraCYP83A1, BraCYP83B1, BraCYP71A1, BraCYP71B13, BraCYP71B2, BraCYP71B37, and BraCYP71B9 regulated in HS and HT (Figure 5; Table S3). BraCYP83A1 was up-regulated at all time points (Figure 5), with a similar result reported in Arabidopsis that the CYP83A1 was involved in flavonoid metabolism (Bilodeau et al., 1999). Similarly, CYP71A1 genes were up-regulated in Panicum virgatum under heat stress and responsible for the synthesis of indole alkaloid (Li et al., 2013) and down-regulated in Rhazya stricta at high temperatures (Obaid et al., 2016). One gene (CYP71A23) in this study was previously identified using genome-wide association study, which caused pollen sterility in Brassica napus under heat stress (Rahaman et al., 2018). Similarly, five MAPK genes were identified in this study, which are critical regulators of plant growth in stress responses (Agrawal et al., 2002; Rohila and Yang, 2007). TFs are essential regulators of genes, involving many growth, developmental, and stress processes. The published studies have reported that heat stress increased the ethylene response, auxin response, zing finger proteins, MADS-box proteins, AP2 domain, WRKY, and leucine zipper factors (Chen et al., 2012; Zander et al., 2012; Kim et al., 2016; Wang et al., 2016; Gu et al., 2018; Rehman et al., 2021). RNA-seq analysis of HS and HT lines in rice and maize determined ethylene response factors that increased stress tolerance (Thirunavukkarasu et al., 2013; Minami et al., 2018). In this study, 31 ERFs were enriched, with >12 up-regulated at 6 and 12 h, possibly playing a significant role in heat tolerance (Figure 6A). Likewise, five WRKYs (BraWRKY18, BraWRKY53, BraWRKY23, BraWRKY25, and BraWRKY40) were up-regulated in HS and HT at 6 and 12 h (Figure 5; Table S4). Of these, WRKY25 was also reported for heat tolerance in Arabidopsis (Chen et al., 2012), and WRKY26 and WRKY33 were also reported for stress resistance (Chen et al., 2012; Kim et al., 2016; Gu et al., 2018). Moreover, many other TFs, such as MYB, NAC, bHLH, and ARF, were detected (Figure 6; Table S4), indicating that TF families are directly related and enhance heat stress tolerance in flowering Chinese cabbage.
HSPs are expressed more during the initial stage of stress than long-term stress (Singh et al., 2019; Qian et al., 2020b). In this study, only ten DEGs were annotated as HSPs, including BraHSP70-5, BraHSP70-10, BraHSP70-8, and BraHSP70-6 up-regulated at 1h and BraHSP70-6-90-1 and BraHSP70-17 up-regulated at 12h (Figure 5; Table S3). Likewise, Wang et al. (2016) reported that HSP70 and HSPB27 were expressed in B. rapa after heat stress, with HSP70 negatively controlling the heat stress response. HSP70s have been identified in many vegetables in response to heat stress, such as potato, tomato, cabbage, and pepper (Guo et al., 2016; Lee et al., 2017; Liu et al., 2018). Studies have demonstrated that HSFs are key regulators of HSPs under heat stress (Wang et al., 2016; Liu et al., 2018), with 35 HSF genes in A, B, and C groups reported in Chinese cabbage (Song et al., 2014b). However, 17 BraHSFs genes identified in this study, including HsfB1, HsfB-2a, HsfA-7a, HsfB-1, HsfA-4a, HsfA-7b), and HsfA-3, were up-regulated at 1, 6, and 12 h in HS and HT (Figure 6; Table S4). Similarly, Scharf et al. (2012) identified HSFA2 as a major heat stress factor that induces HSPs expression in stressed plants. Moreover, overexpression of HSFA2 significantly increased salt, heat, light, and drought stress tolerance (Guo et al., 2016). A tomato HSFA1-a was identified as a master regulator for heat response in rice (Mishra et al., 2018). Thus, the present study increased our understanding of DEGs, their roles in stress, and BraHSF genes under heat stress.
Conclusion
We compared the transcriptomes of HS and HT varieties of flowering Chinese cabbage under heat stress. Approximately 3,205 genes were differentially expressed, with 1,078 DEGs identified at 1, 6, and 12 h after heat stress in both varieties. Cluster analysis divided these genes into three clusters containing 1,958, 591, and 556 genes, which participated in response to stimulus, cell division, cellular process, heat, programmed cell death, ribosome biogenesis, etc. Finally, 15 potential heat-tolerant genes were identified based on functional annotation and literature search. These results provide useful genetic resources for understanding the heat-tolerance mechanism in flowering Chinese cabbage, and these candidate genes require further functional validation and cloning to determine their actual role in heat tolerance.
Data availability statement
The original contributions presented in the study are publicly available. This data can be found here: NCBI, PRJNA885253.
Author contributions
PG and YX conceived and designed the experiments. JC, YX, MI, and RL performed the experiments and analyzed data. YX, MI, and RL contributed reagents/materials/analysis tools. MI, JC, YX, RL, KS, and PG wrote the paper. All authors contributed to the article and approved the submitted version.
Funding
This work was funded by the Natural Science Foundation of Guangdong Province, China (2019A1515011587).
Acknowledgments
We thank the Guangdong Academy of Agricultural Sciences, China, for providing flowering Chinese cabbage accessions.
Conflict of interest
The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.
Publisher’s note
All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.
Supplementary material
The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fpls.2022.1077920/full#supplementary-material
References
Aebi, H. (1984). Catalase in vitro. Methods Enzymol. 105, 121–126. doi: 10.1016/s0076-6879(84)05016-3
Agrawal, G. K., Rakwal, R., Iwahashi, H. (2002). Isolation of novel rice (Oryza sativa l.) multiple stress responsive MAP kinase gene OsMSRMK2 whose mRNA accumulates rapidly in response to environmental cues. Biochem. Biophys. Res. Commun. 294, 1009–1016. doi: 10.1016/S0006-291X(02)00571-5
Ali, A., Wu, T., Xu, Z., Riaz, A., Alqudah, A. M., Iqbal, M. Z., et al. (2022). Phytohormones and transcriptome analyses revealed the dynamics involved in spikelet abortion and inflorescence development in rice. Int. J. Mol. Sci. 23, 7887. doi: 10.3390/ijms23147887
Amirjani, M. R. (2010). Effect of salinity stress on growth, mineral composition, proline content, antioxidant enzymes of soybean. Am. J. Plant Physiol. 5, 350–360. doi: 10.3923/ajpp.2010.350.360
Battisti, D. S., Naylor, R. L. (2009). Historical warnings of future food insecurity with unprecedented seasonal heat. Science 323, 240–244. doi: 10.1126/science.1164363
Bilodeau, P., Udvardi, M. K., Peacock, W. J., Dennis, E. S. (1999). A prolonged cold treatment-induced cytochrome P450 gene from Arabidopsis thaliana. Plant Cell Environ. 22, 791–800. doi: 10.1046/j.1365-3040.1999.00444.x
Caers, M., Rudelsheim, P., Van Onckelen, H., Horemans, S. (1985). Effect of heat stress on photosynthetic activity and chloroplast ultrastructure in correlation with endogenous cytokinin concentration in maize seedlings. Plant Cell Physiol. 26, 47–52. doi: 10.1093/oxfordjournals.pcp.a076894
Cai, Z., He, F., Feng, X., Liang, T., Wang, H., Ding, S., et al. (2020). Transcriptomic analysis reveals important roles of lignin and flavonoid biosynthetic pathways in rice thermotolerance during reproductive stage. Front. Genet. 11. doi: 10.3389/fgene.2020.562937
Chakraborty, U., Tongden, C. (2005). Evaluation of heat acclimation and salicylic acid treatments as potent inducers of thermotolerance in Cicer arietinum l. Curr. Sci. 89, 384–389.
Chang, Y., Bai, Y., Wei, Y., Shi, H. (2020). CAMTA3 negatively regulates disease resistance through modulating immune response and extensive transcriptional reprogramming in cassava. Tree Physiol. 40, 1520–1533. doi: 10.1093/treephys/tpaa093
Chen, J., Li, R., Xia, Y., Bai, G., Guo, P., Wang, Z., et al. (2017). Development of EST-SSR markers in flowering Chinese cabbage (Brassica campestris l. ssp. chinensis var. utilis tsen et Lee) based on de novo transcriptomic assemblies. PloS One 12, e0184736. doi: 10.1371/journal.pone.0184736
Chen, L., Song, Y., Li, S., Zhang, L., Zou, C., Yu, D. (2012). The role of WRKY transcription factors in plant abiotic stresses. Biochim. Biophys. Acta - Gene Regul. Mech. 1819, 120–128. doi: 10.1016/j.bbagrm.2011.09.002
Cock, P. J. A., Fields, C. J., Goto, N., Heuer, M. L., Rice, P. M. (2009). The Sanger FASTQ file format for sequences with quality scores, and the Solexa/Illumina FASTQ variants. Nucleic Acids Res. 38, 1767–1771. doi: 10.1093/nar/gkp1137
Dong, X., Yi, H., Lee, J., Nou, I.-S., Han, C.-T., Hur, Y. (2015). Global gene-expression analysis to identify differentially expressed genes critical for the heat stress response in Brassica rapa. PloS One 10, e0130451. doi: 10.1371/journal.pone.0130451
Dossa, K., Diouf, D., Cissé, N. (2016). Genome-wide investigation of hsf genes in sesame reveals their segmental duplication expansion and their active role in drought stress response. Front. Plant Sci. 7. doi: 10.3389/fpls.2016.01522
Fan, Z.-Q., Tan, X.-L., Shan, W., Kuang, J.-F., Lu, W.-J., Chen, J.-Y. (2017). BrWRKY65, a WRKY transcription factor, is involved in regulating three leaf senescence-associated genes in Chinese flowering cabbage. Int. J. Mol. Sci. 18, 1228. doi: 10.3390/ijms18061228
Gill, R. A., Ali, B., Cui, P., Shen, E., Farooq, M. A., Islam, F., et al. (2016). Comparative transcriptome profiling of two Brassica napus cultivars under chromium toxicity and its alleviation by reduced glutathione. BMC Genomics 17, 885. doi: 10.1186/s12864-016-3200-6
Gill, R. A., Ali, B., Yang, S., Tong, C., Islam, F., Gill, M. B., et al. (2017). Reduced glutathione mediates pheno-ultrastructure, kinome and transportome in chromium-induced Brassica napus l. Front. Plant Sci. 8. doi: 10.3389/fpls.2017.02037
Grover, A., Mittal, D., Negi, M., Lavania, D. (2013). Generating high temperature tolerant transgenic plants: Achievements and challenges. Plant Sci. 205–206, 38–47. doi: 10.1016/j.plantsci.2013.01.005
Guo, M., Liu, J.-H., Ma, X., Zhai, Y.-F., Gong, Z.-H., Lu, M.-H. (2016). Genome-wide analysis of the Hsp70 family genes in pepper (Capsicum annuum l.) and functional identification of CaHsp70-2 involvement in heat stress. Plant Sci. 252, 246–256. doi: 10.1016/j.plantsci.2016.07.001
Gurrieri, L., Merico, M., Trost, P., Forlani, G., Sparla, F. (2020). Impact of drought on soluble sugars and free proline content in selected arabidopsis mutants. Biol. (Basel). 9, 1–14. doi: 10.3390/biology9110367
Gu, L., Wei, H., Wang, H., Su, J., Yu, S. (2018). Characterization and functional analysis of GhWRKY42, a group IId WRKY gene, in upland cotton (Gossypium hirsutum l.). BMC Genet. 19, 48. doi: 10.1186/s12863-018-0653-4
Hasanuzzaman, M., Nahar, K., Alam, M. M., Roychowdhury, R., Fujita, M. (2013). Physiological, biochemical, and molecular mechanisms of heat stress tolerance in plants. Int. J. Mol. Sci. 14, 9643–9684. doi: 10.3390/ijms14059643
Helal, M., Gill, R. A., Tang, M., Yang, L., Hu, M., Yang, L., et al. (2021). SNP- and haplotype-based GWAS of flowering-related traits in Brassica napus. Plants 10, 2475. doi: 10.3390/plants10112475
Hichri, I., Muhovski, Y., Clippe, A., Žižková, E., Dobrev, P. I., Motyka, V., et al. (2016). SlDREB2, a tomato dehydration-responsive element-binding 2 transcription factor, mediates salt stress tolerance in tomato and arabidopsis. Plant Cell Environ. 39, 62–79. doi: 10.1111/pce.12591
Hoermiller, I. I., Funck, D., Schönewolf, L., May, H., Heyer, A. G. (2022). Cytosolic proline is required for basal freezing tolerance in arabidopsis. Plant Cell Environ. 45, 147–155. doi: 10.1111/pce.14196
Huang, X., Lei, Y., Guan, H., Hao, Y., Liu, H., Sun, G., et al. (2017). Transcriptomic analysis of the regulation of stalk development in flowering Chinese cabbage (Brassica campestris) by RNA sequencing. Sci. Rep. 7, 15517. doi: 10.1038/s41598-017-15699-6
Ikram, M., Han, X., Zuo, J. F., Song, J., Han, C. Y., Zhang, Y. W., et al. (2020). Identification of QTNs and their candidate genes for 100-seed weight in soybean (Glycine max l.) using multi-locus genome-wide association studies. Genes (Basel). 11, 1–22. doi: 10.3390/genes11070714
Ikram, M., Lai, R., Xia, Y., Li, R., Zhao, W., Siddique, K. H. M., et al. (2022a). Genetic dissection of tobacco (Nicotiana tabacum l.) plant height using single-locus and multi-locus genome-wide association studies. Agronomy 12, 1047. doi: 10.3390/agronomy12051047
Ikram, M., Xiao, J., Li, R., Xia, Y., Zhao, W., Yuan, Q., et al. (2022b). Identification of superior haplotypes and candidate genes for yield-related traits in tobacco (Nicotiana tabacum l.) using association mapping. Ind. Crops Prod. 189, 115886. doi: 10.1016/j.indcrop.2022.115886
Jones, P. D., New, M., Parker, D. E., Martin, S., Rigor, I. G. (1999). Surface air temperature and its changes over the past 150 years. Rev. Geophys. 37, 173–199. doi: 10.1029/1999RG900002
Jung, Y. J., Go, J. Y., Lee, H. J., Park, J. S., Kim, J. Y., Lee, Y. J., et al. (2021). Overexpression of orange gene (OsOr-R115H) enhances heat tolerance and defense-related gene expression in rice (Oryza sativa l.). Genes (Basel). 12, 1891. doi: 10.3390/genes12121891
Kanehisa, M., Goto, S. (2000). KEGG: Kyoto encyclopedia of genes and genomes. Nucleic Acids Res. 28, 27–30. doi: 10.1093/nar/28.1.27
Kim, D., Langmead, B., Salzberg, S. L. (2015). HISAT: A fast spliced aligner with low memory requirements. Nat. Methods 12, 357–360. doi: 10.1038/nmeth.3317
Kim, C. Y., Vo, K. T. X., Nguyen, C. D., Jeong, D. H., Lee, S. K., Kumar, M., et al. (2016). Functional analysis of a cold-responsive rice WRKY gene, OsWRKY71. Plant Biotechnol. Rep. 10, 13–23. doi: 10.1007/s11816-015-0383-2
Kono, Y. (1978). Generation of superoxide radical during autoxidation of hydroxylamine and an assay for superoxide dismutase. Arch. Biochem. Biophys. 186, 189–195. doi: 10.1016/0003-9861(78)90479-4
Kotak, S., Larkindale, J., Lee, U., von Koskull-Döring, P., Vierling, E., Scharf, K. D. (2007). Complexity of the heat stress response in plants. Curr. Opin. Plant Biol. 10, 310–316. doi: 10.1016/j.pbi.2007.04.011
Kreps, J. A., Wu, Y., Chang, H. S., Zhu, T., Wang, X., Harper, J. F. (2002). Transcriptome changes for arabidopsis in response to salt, osmotic, and cold stress. Plant Physiol. 130, 2129–2141. doi: 10.1104/pp.008532
Kumar, L., Futschik, M. E. (2007). Mfuzz: A software package for soft clustering of microarray data. Bioinformation 2, 5–7. doi: 10.6026/97320630002005
Kumazaki, A., Suzuki, N. (2019). Enhanced tolerance to a combination of heat stress and drought in arabidopsis plants deficient in ICS1 is associated with modulation of photosynthetic reaction center proteins. Physiol. Plant 165, 232–246. doi: 10.1111/ppl.12809
Lai, R., Ikram, M., Li, R., Xia, Y., Yuan, Q., Zhao, W., et al. (2021). Identification of novel quantitative trait nucleotides and candidate genes for bacterial wilt resistance in tobacco (Nicotiana tabacum l.) using genotyping-by-sequencing and multi-locus genome-wide association studies. Front. Plant Sci. 12. doi: 10.3389/fpls.2021.744175
Langmead, B., Salzberg, S. L. (2012). Fast gapped-read alignment with bowtie 2. Nat. Methods 9, 357–359. doi: 10.1038/nmeth.1923
Lee, S. S., Jung, W. Y., Park, H. J., Lee, A., Kwon, S.-Y., Kim, H.-S., et al. (2017). Genome-wide analysis of alternative splicing in an inbred cabbage (Brassica oleracea l.) line ‘HO’ in response to heat stress. Curr. Genomics 19, 12–20. doi: 10.2174/1389202918666170705151901
Li, B., Dewey, C. N. (2011). RSEM: accurate transcript quantification from RNA-seq data with or without a reference genome. BMC Bioinf. 12, 323. doi: 10.1186/1471-2105-12-323
Li, H., Ikram, M., Xia, Y., Li, R., Yuan, Q., Zhao, W., et al. (2022). Genome-wide identification and development of InDel markers in tobacco (Nicotiana tabacum l.) using RAD-seq. Physiol. Mol. Biol. Plants 28, 1077–1089. doi: 10.1007/s12298-022-01187-3
Liu, J., Pang, X., Cheng, Y., Yin, Y., Zhang, Q., Su, W., et al. (2018). The Hsp70 gene family in Solanum tuberosum: Genome-wide identification, phylogeny, and expression patterns. Sci. Rep. 8, 16628. doi: 10.1038/s41598-018-34878-7
Li, Y.-F., Wang, Y., Tang, Y., Kakani, V. G., Mahalingam, R. (2013). Transcriptome analysis of heat stress response in switchgrass (Panicum virgatum l.). BMC Plant Biol. 13, 153. doi: 10.1186/1471-2229-13-153
Li, P.-S., Yu, T.-F., He, G.-H., Chen, M., Zhou, Y.-B., Chai, S.-C., et al. (2014). Genome-wide analysis of the hsf family in soybean and functional identification of GmHsf-34 involvement in drought and heat stresses. BMC Genomics 15, 1009. doi: 10.1186/1471-2164-15-1009
Lobell, D. B., Burke, M. B., Tebaldi, C., Mastrandrea, M. D., Falcon, W. P., Naylor, R. L. (2008). Prioritizing climate change adaptation needs for food security in 2030. Science 319, 607–610. doi: 10.1126/science.1152339
Love, M. I., Huber, W., Anders, S. (2014). Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2. Genome Biol. 15, 550. doi: 10.1186/s13059-014-0550-8
Magwanga, R. O., Lu, P., Kirungu, J. N., Lu, H., Wang, X., Cai, X., et al. (2018). Characterization of the late embryogenesis abundant (LEA) proteins family and their role in drought stress tolerance in upland cotton. BMC Genet. 19, 6. doi: 10.1186/s12863-017-0596-1
Minami, A., Yano, K., Gamuyao, R., Nagai, K., Kuroha, T., Ayano, M., et al. (2018). Time-course transcriptomics analysis reveals key responses of submerged deepwater rice to flooding. Plant Physiol. 176, 3081–3102. doi: 10.1104/pp.17.00858
Mishra, M. K., Chaturvedi, P., Singh, R., Singh, G., Sharma, L. K., Pandey, V., et al. (2013). Overexpression of WsSGTL1 gene of withania somnifera enhances salt tolerance, heat tolerance and cold acclimation ability in transgenic arabidopsis plants. PloS One 8, e63064. doi: 10.1371/journal.pone.0063064
Mishra, N., Srivastava, A. P., Esmaeili, N., Hu, W., Shen, G. (2018). Overexpression of the rice gene OsSIZ1 in arabidopsis improves drought-, heat-, and salt-tolerance simultaneously. PloS One 13, e0201716. doi: 10.1371/journal.pone.0201716
Nakashima, K., Yamaguchi-Shinozaki, K., Shinozaki, K. (2014). The transcriptional regulatory network in the drought response and its crosstalk in abiotic stress responses including drought, cold, and heat. Front. Plant Sci. 5. doi: 10.3389/fpls.2014.00170
Nover, L., Bharti, K., Döring, P., Mishra, S. K., Ganguli, A., Scharf, K. D. (2001). Arabidopsis and the heat stress transcription factor world: How many heat stress transcription factors do we need? Cell Stress Chaperones 6, 177–189. doi: 10.1379/1466-1268(2001)006<0177:aathst>2.0.co;2
Obaid, A. Y., Sabir, J. S. M., Atef, A., Liu, X., Edris, S., El-Domyati, F. M., et al. (2016). Analysis of transcriptional response to heat stress in Rhazya stricta. BMC Plant Biol. 16, 252. doi: 10.1186/s12870-016-0938-6
Pang, T., Ye, C.-Y., Xia, X., Yin, W. (2013). De novo sequencing and transcriptome analysis of the desert shrub, Ammopiptanthus mongolicus, during cold acclimation using Illumina/Solexa. BMC Genomics 14, 488. doi: 10.1186/1471-2164-14-488
Qian, Y., Cao, L., Zhang, Q., Amee, M., Chen, K., Chen, L. (2020b). SMRT and illumina RNA sequencing reveal novel insights into the heat stress response and crosstalk with leaf senescence in tall fescue. BMC Plant Biol. 20, 366. doi: 10.1186/s12870-020-02572-4
Qian, H., Zhang, Y., Wu, B., Wu, S., You, S., Zhang, N., et al. (2020a). Structure and function of HECT E3 ubiquitin ligases and their role in oxidative stress. J. Transl. Intern. Med. 8, 71–79. doi: 10.2478/jtim-2020-0012
Rahaman, M., Mamidi, S., Rahman, M. (2018). Genome-wide association study of heat stress-tolerance traits in spring-type Brassica napus l. under controlled conditions. Crop J. 6, 115–125. doi: 10.1016/j.cj.2017.08.003
Rehman, A., Azhar, M. T., Hinze, L., Qayyum, A., Li, H., Peng, Z., et al. (2021). Insight into abscisic acid perception and signaling to increase plant tolerance to abiotic stress. J. Plant Interact. 16, 222–237. doi: 10.1080/17429145.2021.1925759
Rivero, R. M., Ruiz, J. M., Romero, L. M. (2004). Importance of n source on heat stress tolerance due to the accumulation of proline and quaternary ammonium compounds in tomato plants. Plant Biol. 6, 702–707. doi: 10.1055/s-2004-821293
Rohila, J. S., Yang, Y. (2007). Rice mitogen-activated protein kinase gene family and its role in biotic and abiotic stress response. J. Integr. Plant Biol. 49, 751–759. doi: 10.1111/j.1744-7909.2007.00501.x
Rosenzweig, C., Elliott, J., Deryng, D., Ruane, A. C., Müller, C., Arneth, A., et al. (2014). Assessing agricultural risks of climate change in the 21st century in a global gridded crop model intercomparison. Proc. Natl. Acad. Sci. U. S. A. 111, 3268–3273. doi: 10.1073/pnas.1222463110
Scharf, K. D., Berberich, T., Ebersberger, I., Nover, L. (2012). The plant heat stress transcription factor (Hsf) family: Structure, function and evolution. Biochim. Biophys. Acta - Gene Regul. Mech. 1819, 104–119. doi: 10.1016/j.bbagrm.2011.10.002
Singh, D., Singh, C. K., Taunk, J., Jadon, V., Pal, M., Gaikwad, K. (2019). Genome wide transcriptome analysis reveals vital role of heat responsive genes in regulatory mechanisms of lentil (Lens culinaris medikus). Sci. Rep. 9, 12976. doi: 10.1038/s41598-019-49496-0
Song, X., Hu, J., Wu, T., Yang, Q., Feng, X., Lin, H., et al. (2021). Comparative analysis of long noncoding RNAs in angiosperms and characterization of long noncoding RNAs in response to heat stress in Chinese cabbage. Hortic. Res. 8, 48. doi: 10.1038/s41438-021-00484-4
Song, X., Li, Y., Liu, T., Duan, W., Huang, Z., Wang, L., et al. (2014a). Genes associated with agronomic traits in non-heading Chinese cabbage identified by expression profiling. BMC Plant Biol. 14, 71. doi: 10.1186/1471-2229-14-71
Song, X., Liu, G., Duan, W., Liu, T., Huang, Z., Ren, J., et al. (2014b). Genome-wide identification, classification and expression analysis of the heat shock transcription factor family in Chinese cabbage. Mol. Genet. Genomics 289, 541–551. doi: 10.1007/s00438-014-0833-5
Song, X., Liu, G., Huang, Z., Duan, W., Tan, H., Li, Y., et al. (2016). Temperature expression patterns of genes and their coexpression with LncRNAs revealed by RNA-seq in non-heading Chinese cabbage. BMC Genomics 17, 297. doi: 10.1186/s12864-016-2625-2
Song, Q., Yang, F., Cui, B., Li, J., Zhang, Y., Li, H., et al. (2019). Physiological and molecular responses of two Chinese cabbage genotypes to heat stress. Biol. Plant 63, 548–555. doi: 10.32615/bp.2019.097
Swindell, W. R., Huebner, M., Weber, A. P. (2007). Transcriptional profiling of arabidopsis heat shock proteins and transcription factors reveals extensive overlap between heat and non-heat stress response pathways. BMC Genomics 8, 125. doi: 10.1186/1471-2164-8-125
Thirunavukkarasu, N., Hossain, F., Mohan, S., Shiriga, K., Mittal, S., Sharma, R., et al. (2013). Genome-wide expression of transcriptomes and their co-expression pattern in subtropical maize (Zea mays l.) under waterlogging stress. PloS One 8, 125, e70433. doi: 10.1371/journal.pone.0070433
Wahid, A., Gelani, S., Ashraf, M., Foolad, M. R. (2007). Heat tolerance in plants: An overview. Environ. Exp. Bot. 61, 199–223. doi: 10.1016/j.envexpbot.2007.05.011
Wang, A., Hu, J., Huang, X., Li, X., Zhou, G., Yan, Z. (2016). Comparative transcriptome analysis reveals heat-responsive genes in Chinese cabbage (Brassica rapa ssp. chinensis). Front. Plant Sci. 7. doi: 10.3389/fpls.2016.00939
Wang, C.-T., Ru, J.-N., Liu, Y.-W., Li, M., Zhao, D., Yang, J.-F., et al. (2018). Maize WRKY transcription factor ZmWRKY106 confers drought and heat tolerance in transgenic plants. Int. J. Mol. Sci. 19, 3046. doi: 10.3390/ijms19103046
Wang, X., Wang, H., Wang, J., Sun, R., Wu, J., Liu, S., et al. (2011). The genome of the mesopolyploid crop species Brassica rapa. Nat. Genet. 43, 1035–1040. doi: 10.1038/ng.919
Winfield, M. O., Lu, C., Wilson, I. D., Coghill, J. A., Edwards, K. J. (2010). Plant responses to cold: transcriptome analysis of wheat. Plant Biotechnol. J. 8, 749–771. doi: 10.1111/j.1467-7652.2010.00536.x
Xue, Y., Peng, R., Xiong, A., Li, X., Zha, D., Yao, Q. (2010). Over-expression of heat shock protein gene hsp26 in Arabidopsis thaliana enhances heat tolerance. Biol. Plant 54, 105–111. doi: 10.1007/s10535-010-0015-1
Xu, Q., Xu, X., Shi, Y., Xu, J., Huang, B. (2014). Transgenic tobacco plants overexpressing a grass PpEXP1 gene exhibit enhanced tolerance to heat stress. PloS One 9, e100792. doi: 10.1371/journal.pone.0100792
Yang, X., Zhu, W., Zhang, H., Liu, N., Tian, S. (2016). Heat shock factors in tomatoes: genome-wide identification, phylogenetic analysis and expression profiling under development and heat stress. PeerJ 4, e1961. doi: 10.7717/peerj.1961
Yao, S., Liang, F., Gill, R. A., Huang, J., Cheng, X., Liu, Y., et al. (2020). A global survey of the transcriptome of allopolyploid brassica napus based on single-molecule long-read isoform sequencing and illumina-based RNA sequencing data. Plant J. 103, 843–857. doi: 10.1111/tpj.14754
Yoshinaga, K., Arimura, S. I., Niwa, Y., Tsutsumi, N., Uchimiya, H., Kawai-Yamada, M. (2005). Mitochondrial behaviour in the early stages of ROS stress leading to cell death in Arabidopsis thaliana. Ann. Bot. 96, 337–342. doi: 10.1093/aob/mci181
Yu, E., Fan, C., Yang, Q., Li, X., Wan, B., Dong, Y., et al. (2014). Identification of heat responsive genes in Brassica napus siliques at the seed-filling stage through transcriptional profiling. PloS One 9, e101914. doi: 10.1371/journal.pone.0101914
Yu, X., Wang, H., Lu, Y., De Ruiter, M., Cariaso, M., Prins, M., et al. (2012). Identification of conserved and novel microRNAs that are responsive to heat stress in Brassica rapa. J. Exp. Bot. 63, 1025–1038. doi: 10.1093/jxb/err337
Zander, M., Chen, S., Imkampe, J., Thurow, C., Gatz, C. (2012). Repression of the arabidopsis thaliana jasmonic acid/ethylene-induced defense pathway by tga-interacting glutaredoxins depends on their c-terminal alwl motif. Mol. Plant 5, 831–840. doi: 10.1093/mp/ssr113
Zhang, Z., Huang, R. (2013). Analysis of malondialdehyde, chlorophyll proline, soluble sugar, and glutathione content in arabidopsis seedling. BIO-PROTOCOL 3, 817. doi: 10.21769/BioProtoc.817
Keywords: heat stress, transcriptome, differentially expressed genes, cluster analysis, flowering chinese cabbage, RT-qPCR
Citation: Ikram M, Chen J, Xia Y, Li R, Siddique KHM and Guo P (2022) Comprehensive transcriptome analysis reveals heat-responsive genes in flowering Chinese cabbage (Brassica campestris L. ssp. chinensis) using RNA sequencing. Front. Plant Sci. 13:1077920. doi: 10.3389/fpls.2022.1077920
Received: 23 October 2022; Accepted: 21 November 2022;
Published: 02 December 2022.
Edited by:
Adnan Mustafa, Brno University of Technology, CzechiaReviewed by:
Muhammad Waheed Riaz, Zhejiang Agriculture and Forestry University, ChinaRafaqat Ali Gill, Oil Crops Research Institute (CAAS), China
Fazal Rehman, South China Botanical Garden (CAS), China
Copyright © 2022 Ikram, Chen, Xia, Li, Siddique and Guo. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.
*Correspondence: Peiguo Guo, guopg@gzhu.edu.cn; Yanshi Xia, xiayanshi922@163.com
†These authors have contributed equally to this work