Full-length transcriptome analysis of maize root tips reveals the molecular mechanism of cold stress during the seedling stage
BMC Plant Biology volume 22, Article number: 398 (2022)
As maize originated in tropical or subtropical zones, most maize germplasm is extremely sensitive to low temperatures during the seedling stage. Clarifying the molecular mechanism of cold acclimation would facilitate the breeding of cold tolerant maize varieties, which is one of the major sustainability factors for crop production. To meet this goal, we investigated two maize inbred lines with contrasting levels of cold tolerance at the seedling stage (IL85, a cold tolerant line; B73, a cold sensitive line), and performed full-length transcriptome sequencing on the root tips of seedlings before and after 24 h of cold treatment.
We identified 152,263 transcripts, including 20,993 novel transcripts, and determined per-transcript expression levels. A total of 1,475 transcripts were specifically up-regulated in the cold tolerant line IL85 under cold stress. GO enrichment analysis revealed that 25 transcripts were involved in reactive oxygen species (ROS) metabolic processes and 15 transcripts were related to the response to heat. Eight genes showed specific differential alternative splicing (DAS) in IL85 under cold stress, and were mainly involved in amine metabolism. A total of 1,111 lncRNAs were further identified, 62 of which were up-regulated in IL85 or B73 under cold stress, and their corresponding target genes were enriched in protein phosphorylation.
These results provide new insights into the molecular mechanism of cold acclimation during the seedling stage in maize, and will facilitate the development of cultivars with improved cold stress tolerance.
Maize (Zea mays L.) is one of the most important cereal crops. However, most maize germplasm is sensitive to low temperature, as it originated in tropical or subtropical zones [1, 2]. At high latitudes and altitudes, there is a high risk of exposure to cold stress for maize, which reduces plant growth, and grain yield and can even cause death [3, 4]. Especially at the seedling stage, the root, as the fundamental organ for water and nutrient absorption, would affect the growth and development of maize seedlings at suboptimal temperatures . Therefore, understanding the root response to cold stress at the molecular level is of great importance for developing cold-resistant crops, which is one of the major sustainability factors for crop production.
Assessing gene expression patterns based on transcriptome sequencing data is an important way to elucidate the molecular mechanism of cold acclimation in plants. Many studies have shown the effectiveness of transcriptome analysis in understanding how altered expression contributes to complex traits, such as cold tolerance [6, 7]. Using transcriptome analysis, more than 100 regulons of CBF transcription factors were found to be rapidly induced under cold stress . In addition to genes related to the CBF cold response pathway, genes associated with signal transduction, membrane part, response stimulus, plant hormone signal transduction, abiotic stress response and zeatin biosynthesis, were significantly differentially expressed in plants under cold stress according to a comparison of transcriptomes between genotypes with contrasting levels of cold tolerance [9,10,11]. In maize, a number of cold response genes have been identified, and the functions of these genes were found to be largely similar to those of other plants [12,13,14,15,16,17,18,19,20]. Although thousands of cold-responsive genes have been identified in plants, most of these studies mainly used the Illumina sequencing platform, which generates short reads (< = 150 bp), making it difficult to quantify expression at the transcript level.
Full-length RNA sequencing has recently become commercially available recently, and can be performed by Oxford Nanopore Technologies (ONT), Pacific Biosciences (PacBio) and others [21, 22]. These technologies enable single-molecule sequencing of complete individual RNA molecules or cDNA, and enabled the generation of full-length cDNA/RNA reads for transcripts up to 15 kb . Thus, this technique permits efficient transcript assembly, the quantification of transcript levels and alternative splicing (AS) events, and a novel transcript identification [24, 25]. Using the ONT sequencing platform, Yu et al. (2021) identified 1390 transcripts in P. damicornis in response to different growing environments . Utilizing the PacBio sequencing platform, a number of novel transcripts and genes were identified and 56 differentially expressed genes were found under drought/heat stress in rice . Most maize genes contain multiple exons and introns, which can produce several transcripts . Identifying the changes in transcript abundance under cold stress is crucially important for researchers. To date, few studies have focused on the full-length transcriptome in response to cold stress in maize. Therefore, to advance the knowledge of molecular mechanisms under cold stress, it is necessary to assess global full-length transcriptome changes under cold stress during the seedling stage in maize.
In this study, we screened two maize inbred lines with contrasting levels of cold tolerance at the seedling stage. To characterize different transcriptional responses in the root tips of the two lines under cold stress, the ONT sequencing platform was used on the 2 cm samples of root tips from seedlings before and after 24 h of cold treatment (10 °C). Due to the long-read length of ONT Nanopore-based RNA sequencing, the accuracy of transcript assembly was greatly improved, and per-transcript expression abundance was quantified and used for further differentially expressed transcript (DET) and alternative splicing (AS) event analysis. The results of this study provide important information regarding the molecular basis of cold acclimation and the future of cold tolerance breeding in maize.
Effects of cold stress on the electrolyte leakage of roots
Based on germplasm screening, IL85 and B73 were regarded as cold tolerant and sensitive lines, respectively. Electrolyte leakage frequently serves as an indicator of cold stress, as the membrane is the first organelle to sense cold stress. To evaluate the electrolyte leakage caused by cold stress, the root tips of the seedlings that were grown at the control temperature for 7 days, and then transferred to cold conditions for 24 h were used to measure the relative electrolyte conductivity (REC) (Fig. 1A). Compared with the control, the REC of root tips of both IL85 and B73 increased after cold stress, as expected. The increase in REC in B73 was more pronounced (approximately 3 times) than that in IL85 (Fig. 1B). These results revealed that IL85 was more resistant to cold than B73.
Full-length sequencing and transcript expression profiles of 12 samples
To characterize the transcriptomic difference between IL85 and B73 under cold stress, seedlings grown under control or cold conditions were used for root tips collection, library preparation and full-length transcriptome sequencing. A total of 55 Gb of clean data was produced with 3,890,121 to 6,541,618 reads. After removing rRNA sequences, 3,759,176 to 6,309,821 reads were obtained. These reads had an N50 of 695 – 1,265 bp, a mean length of 661 – 1,089 bp, and a maximum length of 14,970–233,119 bp. The number of full-length reads, with primers at both ends of the read, ranged from 3,297,036 to 5,744,254, with an N50 of 552 – 1,122 bp, a mean length of 460–881 bp, and a maximum length of 6,950–11,677 bp, among which 93.7% to 96.9% were mapped to the reference genome (Table 1). Thanks to the long reads produced by the ONT sequencing platform, we identified 25,309 novel transcripts and 4,316 novel genes (Table S1, Dataset 1 and Dataset 2).
Principal component analysis (PCA) and hierarchical clustering analysis based on the counts per million (CPM) of all transcripts in 12 samples revealed a major genotypic effect as well as an effect of treatment. Replicate samples clustered together and showed a reasonably high correlation as expected (Fig. 2). Furthermore, both control and cold-treated samples of B73 trended to cluster together, as did those of IL85 (Fig. 2B).
DET identification and characterization
To assess the transcriptome changes under cold stress, DETs with at least a twofold difference in expression and a p value less than 0.01 were identified for two comparisons, BCK v.s BCD and ICK v.s ICD. A total of 3,655 and 2,598 up- or down-regulated DETs were identified in B73. For IL85, 2,712 and 1,693 up- or down-regulated DETs were identified. Further analysis using a Venn diagram showed that both unique and overlapping DETs were detected (Fig. 3A-C). For up-regulated DETs, a total of 2,418 and 1,475 DETs were unique to B73 and IL85, and 1,237 DETs were commonly differentially expressed in both lines. Regarding down-regulated DETs, a total of 1,958 and 1,053 DETs were unique to B73 and IL85, respectively, and 640 DETs were commonly differentially expressed in both lines. Among these DETs, more than 50% were novel transcripts (Fig. 3A-C). Among these DETs, a total of 982 could be recognized as transcription factors (TFs), 591 or 391 of which were up- or down-regulated in both lines after 24 h cold stress (Table S2).
Functional analysis of DETs
The common or specific DETs might participate in the basic or genotype-specific cold response process. To analyse the function of these DETs, GO enrichment analyses were performed. GO enrichment analysis of DETs revealed that these transcripts were related to various functions in the biological process (BP) category (Fig. 3D). For up-regulated transcripts under cold stress, DETs that were common to both B73 and IL85 were mainly related to reactive oxygen scavenging, sulfur compound or amino acid metabolic processes. Transcripts that were specifically differentially expressed in IL85 were also related to reactive oxygen scavenging and other GO terms, such as response to heat and inositol. Transcripts specific to B73 were mainly involved in sulfur compound metabolic process, alpha-amino acid metabolic process and so on. For down-regulated transcripts under cold stress, the DETs that were common to both B73 and IL85 were mainly related to reactive oxygen species scavenging. Furthermore, the same reactive oxygen species scavenging related GO terms were identified for transcripts that were specifically down-regulated in B73. Transcripts that were specifically expressed in IL85 were mainly related to glutathione metabolic process, glutamine metabolic process and response to jasmonic acid.
To identify the molecular pathways underlying cold stress, KEGG pathway enrichment analysis was performed on the specific and common DETs mentioned above. These transcripts were enriched in a total of ten pathways (Fig. 3E). The common up-regulated DETs in both IL85 and B73 were related to six pathways. Transcripts related to the pentose phosphate pathway were specifically up-regulated in IL85 under cold stress. B73 specific up-regulated transcripts were enriched in three pathways, such as tryptophan metabolism. For down-regulated transcripts, the diterpenoid biosynthesis pathway was identified and was unique to B73, but no pathway was detected for IL85.
The expression patterns of key transcripts in IL85 under cold stress were investigated. There were 67 transcripts related to three kinds of functions, including 34 for reactive oxygen species scavenging, 21 for response to heat and 12 for inositol (Table S3). The heatmap revealed that the expression of these transcripts was relatively lower in B73 than in IL85 under both conditions. After 24 h of cold stress, the expression of these transcripts rapidly increased in IL85, while they slightly decreased in B73 (Fig. S1).
Characterization of alternative splicing events
Alternative splicing (AS) is a key step in plant abiotic stress acclimation. Identifying and characterizing AS events enhances our understanding of the biological role of transcript isoform diversity. In this study, five categories of AS events were identified in 12 samples, the most abundant event was intron retention (30.3%-53.3%), and the least abundant event was mutually exclusive exons (0.3%-1.4%). The number of AS events was higher in samples under cold stress than under the control temperature, especially intron retention events (Fig. 4A). The AStalavista tool was also used to detect differential AS (DAS) events between samples with or without cold stress in each inbred line. A total of 662 and 509 genes with DAS (DAS genes) were identified for IL85 and B73, respectively. There were 201 genes that showed DAS in both lines, and a total of 461 and 308 DAS genes were specific to IL85 and B73, respectively (Fig. 4B). Interestingly, the common DAS genes were significantly enriched in seven GO terms in the BP category, and all these terms were related to RNA processing, as expected (Fig. 4C). The DAS genes uniquely identified in IL85 were mainly enriched in three amine metabolic process related GO terms (Fig. 4D), and the specific DAS genes for B73 were not significantly enriched in any GO terms in the BP category (Fig. 4E). KEGG enrichment analysis of these DAS genes was used to explore their molecular pathways. A total of four pathways were detected, and all three DAS gene sets were enriched in the spliceosome. The TCA cycle was found to be a unique pathway in IL85, while cysteine and methionine metabolism was unique to B73 (Fig. 4F). The results presented here suggested that the AS of genes played an important role in the maize cold stress response and that there were obvious genotypic differences.
Eight DAS genes that were unique to IL85 were related to amine metabolic process according to GO enrichment analysis. These genes contained 126 transcripts including 40 novel transcripts (Table S4). The expression of each transcript is shown in Figure S2A. More than one-third of the transcripts were predominantly expressed in IL85 root tips under cold stress, and seven out of eight genes had at least one transcript that was predominantly expressed in ICD samples. To further investigate the AS event, sashimiplots were generated for Zm00001d019698 and Zm00001d053404 using igv tools. Zm00001d019698 encodes a cyclase-like protein 4 with 7 exons, and many reads mapped to the beginning of the third intron in B73 (Fig. S2B), which suggested that different intron retention events occurred between the two genotypes. The gene Zm00001d053404 has 13 transcripts, and encodes a spermidine synthase. Additionally, because of intron retention events, six novel transcripts were identified in B73 under cold stress (Fig. S2C). These results suggested that the genotype- and stress-specific splicing isoforms might play an important role in maize cold tolerance.
Characterization of lncRNAs
Long noncodling RNAs (lncRNAs) were identified by performing CNCI, CPC, CPAT and Pfam analysis, and a total of 1,111 lncRNAs were detected by all four methods (Fig. 5A and Dataset 3). These lncRNAs grouped into four categories with 43 anti-sense lncRNAs, 20 intronic lncRNAs, 542 lincRNAs, and 506 sense lncRNAs (Fig. 5B). The length of most of the lncRNAs was less than 3000 bp (Fig. 5C). The expression profiles of lncRNAs are shown in Fig. 5D. All lncRNAs can be divided into two major groups, those mainly expressed in B73 or IL85, and there were differences in lncRNA expression patterns between treatments within each genotype (Fig. 5D). The expression differences of each lncRNA were calculated by the DEseq2 package. With a cut-off of an absolute value of log2FC > 1 and p-value < = 0.01, a total of 106 lncRNAs showed a significant difference between treatments within each genotype, with 62 up-regulated (Fig. 5E), and 44 down-regulated lncRNAs (Fig. 5F). For up-regulated lncRNAs, 33 and 23 were specifically expressed in IL85 and B73, and 6 lncRNAs were commonly expressed in both IL85 and B73. Among the down-regulated lncRNAs, 20 and 18 were specific to IL85 and B73, respectively, and 6 lncRNAs were common to IL85 and B73.
To explore the potential function of lncRNAs, their target genes were predicted according to their position and complementary sequences. A total of the 9069 genes were predicted to be targets of these 1,111 lncRNAs, including 3,580 cis- targeted and 5,489 trans-targeted genes (Dataset 4). For the 62 up-regulated lncRNAs in IL85 or B73, there were 495 corresponding target genes, which were mainly involved in the phosphorylation process according to GO enrichment analysis (Fig. 5G). There were 572 target genes for the 44 down-regulated lncRNAs, and no significant GO terms or KEGG pathways were identified by enrichment analysis, indicating that down-regulated lncRNAs might not participate in particular biological processes or be involved in more extensive processes.
Validation of the expression of selected genes by Quantitative Real-Time PCR (qRT-PCR)
The expression of ten genes, including five novel genes, four known genes and an endogenous control gene (ZmTub), was quantified with qRT-PCR to validate the results of ONT RNA sequencing (Fig. 6). As expected, there was a similar tendency and a high correlation between the results of qRT-PCR and the ONT RNA sequencing global profiling (Fig. S3), demonstrating the robustness of the RNA sequencing results.
Abiotic stress leads to substantial changes in the transcriptome in plants. To precisely characterize gene expression patterns and molecular features under cold stress, several studies have assessed transcriptome changes at the germination and seedling stages using the Illumina sequencing platform, and hundreds of genes in response to cold stress have been identified [12,13,14,15,16,17,18,19,20]. In maize, approximately 50% of genes contain more than 4 exons, with an average of 3.3 transcripts per gene . Given that mRNA transcripts are straightforward template for translation, identifying the changes in transcript abundance under cold stress is the main task. In this study, we performed full-length RNA sequencing on the root tips of seedlings with or without cold stress in maize. A total of 152,263 transcripts were detected, including 20,993 novel transcripts (Table S1). We identified 3,655 up-regulated DETs in IL85 under cold stress through comparative transcriptome analysis (Fig. 3A), these DETs were mainly involved in the ROS scavenging process (Fig. 3D). Alternative splicing (AS) events are markedly induced in response to environmental stress, allowing rapid adjustment of the abundance and function of key stress-response components . Long read lengths contributed to the identification of AS events. A total of 970 DAS events were detected in response to cold stress, and the overlapping set of DAS involved genes was mainly related to RNA splicing or spliceosome, as expected (Fig. 4D and F). The DAS genes unique to IL85 were related to the amine metabolic process (Fig. 4C). LncRNAs have been shown to be involved in the response to cold stress [30,31,32]. We also identified 62 and 44 lncRNAs that were up- and down-regulated in both lines under cold stress, and the target genes of up-regulated lncRNAs were mainly involved in protein phosphorylation (Fig. 5G), which is consistent with the findings of others . These results suggested that three analysis methods, the identification of DETs, AS events and differential lncRNAs, reflected different aspects of transcriptome response to cold stress, and no single analysis can capture the complete view.
ROS scavenging enzymes played an important role in maize cold acclimation
ROS can be generated as a consequence of electron leakage during photosynthesis and respiration in cells. Superoxide radicals, hydrogen peroxide, hydroxyl radicals, and singlet oxygen are the main ROS formed in response to the reduction of oxygen molecules in planta . The production and removal of ROS must be strictly controlled, because ROS cause lipid peroxidation, which leads to the loss of plasma membrane integrity and electrolyte leakage [34,35,36,37,38]. However, the equilibrium between the production and scavenging of ROS may be perturbed by numerous adverse abiotic stress factors, such as cold stress . To moderate intracellular ROS homeostasis, plant have evolved an enzymatic antioxidative mechanism, including superoxide dismutase (SOD), catalase (CAT), and peroxidase (POD). In this study, a total of 34 transcripts that related to ROS scavenging enzymes were predominantly expressed in IL85 under cold stress, 25 of which encode peroxidase, while transcripts with similar function were significantly down-regulated in B73. This is perhaps the principal reason that the REC of IL85 is much lower than that of B73 under cold stress (Fig. 1B). To illustrate this, we have measured the activities POD, CAT and SOD, the results show that the activities of all the three antioxidant enzymes in IL85 were higher than that in B73, although the activity of POD was also increased after 24 h of cold stress in B73 (Fig. 7A to 7C). As a result, changes in antioxidant enzymes would alter the cellular redox homeostasis. Malondialdehyde (MDA) is the best investigated product of lipid peroxidation induced by oxidative stress , it could be served as an indirect indicator of ROS level. Then the content of MDA was measured, and it was lower in IL85 root tips than that in B73 (Fig. 7D) as expected, which may explain part of the different REC between IL85 and B73 under cold stress (Fig. 1B). It has been revealed that the activities of SOD, POD and CAT increase in X. sorbifolia under cold stress . The increase in antioxidant enzyme activity was the main reason for the enhanced cold tolerance in transgenic rice (GmFAD3A) , as several POD genes were dramatically induced after cold treatment. The expression of these genes increased in the overexpression plants (PtrbHLH) but decreased in the RNAi plants of trifoliate orange . As a result, the accumulation of ROS scavenging enzymes was appreciably alleviated in the transgenic lines of trifoliate orange under cold stress .
Amine metabolism and maize cold tolerance
Polyamines (PAs) are the main forms of amine-containing compounds in plants, and are elevated under cold stress. They have dual effects, such as being protectors and perpetrators of stress damage to the cells. It has been demonstrated that PAs play a role in providing cold tolerance [45, 46]. In this work, we found that amine metabolism related genes showed specific DAS events (Fig. 4C), and some transcripts were predominantly expressed in IL85 under cold stress (Fig. S2A). Compared with the other treatments, putrescine (Put) and spermidine (Spd) levels were much greater under cold stress . Putrescine (Put) has been shown to accumulate in the leaves of tomato but to be higher in chilling-tolerant cultivar than in chilling-sensitive cultivars . The rate of polyamine (spermine and spermidine) synthesis and the levels of transcripts encoding arginine decarboxylase, SAMDC and spermine/spermidine synthase also increased after cold exposure [49,50,51]. A naturally cold-tolerant variety of rice showed higher expression of SAMDC genes . Loss of function of ADC1 and ADC2, two key genes for polyamine metabolism, reduced cold tolerance in Arabidopsis . Increasing the content of PAs contributed to alleviating chilling stress in cucumber and tomato seedlings when melatonin was applied [48, 53]. These results suggested that understanding the role of the PA metabolism at the molecular level will be pivotal for improving maize performance under cold stress conditions.
Polyamine oxidases involved in ROS and PA metabolism
The accumulation of ROS was frequently accompanied by an increase in polyamines during cold acclimation . ROS, such as H2O2, can be produced by polyamine catabolism processes. In plants, polyamine oxidases (PAOs) catalyse the terminal catabolism (TC) of Spm and Spd and produce N-(3-aminopropyl)-4-aminobutanal from Spm and 4-aminobutanal from Spd, along with 1,3-diaminopropane (DAP) and H2O2 [54,55,56], PAO is the most important enzyme in polyamine homeostasis and plays an essential role in growth and developmental processes and responses to abiotic stresses . In this study, we found that at least one of the polyamine oxidase encoding transcripts was more highly expressed in IL85 under cold stress (Fig. 4C). In rice, at least three PAOs (OsPAO4/6/7) were induced by cold treatment, while some PAOs were invariable or even decreased under cold stress [58, 59]. In turn, OsPAO2 and OsPAO6 were strongly induced upon oxidative stress (H2O2) .
In summary, we characterized the transcript changes related to early cold stress and the initiation of cold acclimation using full-length RNA sequencing. By comparing the full-length transcriptome between seedlings subjected to cold stress across two maize inbred lines (IL85 and B73), we suggested that transcripts involved in ROS scavenging enzymes and amine metabolism played an important role in cold acclimation. These results provide new insights into the molecular mechanism of cold acclimation at the seedling stage in maize, and facilitate to development of cultivars with improved cold stress tolerance.
According to our previous experiment, the maize inbred line IL85 was selected as a cold tolerant line, and B73 was used as a cold susceptible line. The seeds of both lines were harvested from the experimental base of the Institute of Nanfan & Seed Industry after self-pollination. Then, the seeds, sterilized with 1% sodium hypochlorite for 5 min, were placed on moist germination paper and another sheet of moist paper was used as a cover. The two pieces of germination paper were rolled together and placed vertically in a sealed plastic bag . After germinating at 25 °C/22 °C with an 8 h light/16 h dark cycle for 7 days, the root tips (bottom ~ 2 cm) of the primary roots were sampled. Then, the remaining paper rolls were transferred to an incubator with a constant temperature of 10 °C with an 8 h light/16 h dark cycle for 24 h, and the root tips were harvested. Each sample was divided into two equal parts. One part was immediately stored in liquid nitrogen for total RNA extraction, and the remaining part was used for electrical conductivity measurement.
Determination of electrolyte leakage
A total of 15 freshly harvested root tips were incubated at room temperature in tubes with 15 mL of distilled water for 20 min. Then, the electrical conductivity (EC1) of the solution was measured using a conductivity meter (Mettler Toledo). The tubes were heated at 100 °C for 15 min. After cooling to room temperature, the electrical conductivity (EC2) was determined. The final relative electrical conductivity (REC) was estimated according to the equation :
RNA extraction, library preparation, and Oxford Nanopore Full-Length Sequencing
Root tips stored in liquid nitrogen were pooled and ground in liquid nitrogen for total RNA extraction using the RNAprep Pure Plant Kit (Tiangen, China). Then, 1 ug of total RNA was prepared for cDNA libraries using a cDNA-PCR Sequencing Kit (SQK-PCS109) according to the manufacturer’s protocol (Oxford Nanopore Technologies, US). The final cDNA libraries were added to FLO-MIN109 flowcells and run on the PromethION platform at Biomarker Technology Company.
Read quality control, alignment, transcript expression quantification and differentially expressed transcript (DET) identification
For each sample, raw reads were first filtered with the cutoff of minimum average read quality score = 7 and minimum read length = 500 bp, and the ribosomal RNAs were removed after mapping to the rRNA database. Second, full-length non-chemiric (FLNC) transcripts were identified by searching for primers at both ends of the reads. Clusters of FLNC transcripts were obtained after mapping to the B73 reference genome (GCA_000005005.6_B73_RefGen_v4) with mimimap2 , and consensus isoforms were obtained after polishing within each cluster by the pinfish tool. The mapped consensus reads were further collapsed by the cDNA_Cupcake package with min-coverage = 85% and min-identity = 90%. A 5’ difference was not considered when collapsing redundant transcripts. Per-transcript expression levels were determined based on reads with a match quality above 5. Expression levels were estimated by read counts per gene/transcript per 10,000 reads mapped (CPM) in each sample, allowing direct comparison across samples.
The DESeq2 package (1.6.3) was employed to identify differentially expressed transcripts (DETs) between two samples . The resulting P-values were adjusted using Benjamini and Hochberg’s approach for controlling the false discovery rate. Genes/transcripts with an FDR < 0.01 and fold_change ≥ 2 were considered significantly differentially expressed.
Transcription factors identification
All maize transcription factors (TFs) were identified in Plant Transcription Factor Database (Plant TFDB, http://planttfdb.gao-lab.org/index.php?sp=Zma)  and GRASSIUS (https://grassius.org/tfomecollection.php) . A gene was regarded as TF if it was in any of databases, then the corresponding transcripts of TF encoding genes were extracted.
Alternative splicing event analysis
Transcripts were validated against known reference transcript annotations with gffcompare . AS events, including intron retention (IR), exon skipping (ES), alternative 5’ splice site (A5SS), alternative 3’ splice site (A3SS) and mutually exclusive exon (MEE), and the differential AS events between two samples were detected by the AStalavista tool .
Four computational approaches, CPC/CNCI/CPAT/Pfam, were combined to sort long noncoding RNA candidates from putative protein-coding RNAs in the transcripts [68,69,70,71]. Putative protein-coding RNAs were filtered out using a minimum length and exon number threshold. Transcripts with lengths greater than 200 nt and with more than two exons were selected as lncRNA candidates and further screened using CPC/CNCI/CPAT/Pfam, which has the power to distinguish protein-coding genes from the noncoding genes. The differential lncRNAs were identified using DESeq2 (1.6.3) as in the identification of DETs.
Functional annotation and enrichment analysis
Gene/transcript functions were annotated based on the Gene Ontology (GO) [72, 73] and Kyoto Encyclopedia of Genes and Genomes (KEGG) databases [74,75,76]. We performed GO and KEGG enrichment analyses of DEGs/DETs based on hypergeometric tests using the phyper function in R software , and the resulting p-values were adjusted using Benjamini and Hochberg’s approach for controlling the false discovery rate. The results of enrichment analysis were visualized by the enrichplot and ggplot2 packages.
Quantitative Real-Time PCR (qRT-PCR) analysis
The expression of ten randomly selected genes, including five novel genes, four known genes and an endogenous control gene, was validated by qRT-PCR. To normalize the relative expression levels for each gene, maize ZmTub served as an endogenous control. The primers used for qRT-PCR were designed by the web-based software Primer3 (https://bioinfo.ut.ee/primer3/) and their details are shown in Table S5. cDNA was synthesized from the total RNA prepared for transcriptome sequencing. qRT-PCR MasterMix (Applied Biological Materials Inc.) was used for expression level analysis according to the manufacturer’s protocol. Three biological replicates were conducted and each biological replicate was technically repeated three times. The 2−ΔCT method was used to calculate the relative abundance of genes.
The activities of antioxidant enzymes and the content of MDA measurement
The root tips were sampled as described above, the activities of antioxidant enzymes (CAT, POD, SOD) were measured according to Elavarthi and Zeng with minor modification [78, 79]. Briefly, approximately 0.1 g root tips tissues were homogenized with 5 μl of Tris–HCl buffer (pH 7.0) containing 20% glycerinum, 1 mM EDTA,1 mmol/L ASA, 1 mmol/L DTT, 1 mmol/L GSH and 5 mmol/L MgCl2 on ice. Then the homogenate was centrifuged at 10000 g for 20 min at 4 °C. The supernatant was collected and used for CAT, POD and SOD activity analysis. CAT activity was defined as the absorbance changes (OD240) per gram of fresh weight per minute, POD activity was expressed the absorbance change (OD470) per gram of fresh weight per minute, and the SOD activity was represented as unit that the amont of SOD causing 50% inhibition of NBT reduction.
The content of MDA was measured using a Solarbio kit (BC0025), the main principle behind this kit was that MDA and TBA could be synthetized trimethine under acidic and high temperature conditions, and the trimethine has a maximum absorbance at 532 nm.
Availability of data and materials
The datasets supporting the conclusions of this article is(are) available at NCBI under SRA accession number PRJNA822071.
Oxford Nanopore technologies
Differentially expressed transcript
Differentially alternative splicing
Long noncoding RNA
Reactive oxygen species
Relative electrolyte conductivity
Principal component analysis
Counts per million
Kyoto encyclopedia of genes and genomes
Quantitative real-time PCR
Alternative 5’ splice site
Alternative 3’ splice site
Mutually exclusive exon
Ethylene diamine tetraacetic acid
Bhosale SU, Rymen B, Beemster GTS, Melchinger AE, Reif JC. Chilling tolerance of Central European maize lines and their factorial crosses. Ann Bot. 2007;100:1315–21.
He F, Shen H, Lin C, Fu H, Sheteiwy MS, Guan Y, et al. Transcriptome Analysis of Chilling-Imbibed Embryo Revealed Membrane Recovery Related Genes in Maize. Front Plant Sci. 2016;7:1978.
Stirling CM, Nie GY, Aguilera C, Nugawela A, Long SP, Baker NR. Photosynthetic productivity of an immature maize crop: changes in quantum yield of CO2 assimilation, conversion efficiency and thylakoid proteins. Plant Cell Environ. 1991;14:947–54.
Leipner J, Fracheboud Y, Stamp P. Effect of growing season on the photosynthetic apparatus and leaf antioxidative defenses in two maize genotypes of different chilling tolerance. Environ Exp Bot. 1999;42:129–39.
Li Y, Wang X, Li Y, Zhang Y, Gou Z, Qi X, et al. Transcriptomic Analysis Revealed the Common and Divergent Responses of Maize Seedling Leaves to Cold and Heat Stresses. Genes. 2020;11:881.
Fowler S, Thomashow MF. Arabidopsis transcriptome profiling indicates that multiple regulatory pathways are activated during cold acclimation in addition to the CBF cold response pathway. Plant Cell. 2002;14:1675–90.
Dasgupta P, Das A, Datta S, Banerjee I, Tripathy S, Chaudhuri S. Understanding the early cold response mechanism in IR64 indica rice variety through comparative transcriptome analysis. BMC Genomics. 2020;21:425.
Park S, Lee C-M, Doherty CJ, Gilmour SJ, Kim Y, Thomashow MF. Regulation of the Arabidopsis CBF regulon by a complex low-temperature regulatory network. Plant J. 2015;82:193–207.
Jia Y, Liu H, Qu Z, Wang J, Wang X, Wang Z, et al. Transcriptome Sequencing and iTRAQ of Different Rice Cultivars Provide Insight into Molecular Mechanisms of Cold-Tolerance Response in Japonica Rice. Rice (N Y). 2020;13:43.
da Maia LC, Cadore PRB, Benitez LC, Danielowski R, Braga EJB, Fagundes PRR, et al. Transcriptome profiling of rice seedlings under cold stress. Funct Plant Biol. 2017;44:419–29.
Pan Y, Liang H, Gao L, Dai G, Chen W, Yang X, et al. Transcriptomic profiling of germinating seeds under cold stress and characterization of the cold-tolerant gene LTG5 in rice. BMC Plant Biol. 2020;20:371.
Frascaroli E, Revilla P. Genomics of Cold Tolerance in Maize. In: Bennetzen J, Flint-Garcia S, Hirsch C, Tuberosa R, editors. The Maize Genome. Cham: Springer International Publishing; 2018. p. 287–303.
Jończyk M, Sobkowiak A, Trzcinska-Danielewicz J, Skoneczny M, Solecka D, Fronk J, et al. Global analysis of gene expression in maize leaves treated with low temperature. II. Combined effect of severe cold (8 °C) and circadian rhythm. Plant Mol Biol. 2017;95:279–302.
Li Z, Hu G, Liu X, Zhou Y, Li Y, Zhang X, et al. Transcriptome Sequencing Identified Genes and Gene Ontologies Associated with Early Freezing Tolerance in Maize. Front Plant Sci. 2016;7:1477.
Lu X, Zhou X, Cao Y, Zhou M, McNeil D, Liang S, et al. RNA-seq Analysis of Cold and Drought Responsive Transcriptomes of Zea mays ssp. mexicana L. Front Plant Sci. 2017;8:136.
Mao J, Yu Y, Yang J, Li G, Li C, Qi X, et al. Comparative transcriptome analysis of sweet corn seedlings under low-temperature stress. The Crop Journal. 2017;5:396–406.
Rymen B, Fiorani F, Kartal F, Vandepoele K, Inze D, Beemster GTS. Cold nights impair leaf growth and cell cycle progression in maize through transcriptional changes of cell cycle genes. Plant Physiol. 2007;143:1429–38.
Sobkowiak A, Jończyk M, Jarochowska E, Biecek P, Trzcinska-Danielewicz J, Leipner J, et al. Genome-wide transcriptomic analysis of response to low temperature reveals candidate genes determining divergent cold-sensitivity of maize inbred lines. Plant Mol Biol. 2014;85:317–31.
Sowiński P, Fronk J, Jończyk M, Grzybowski M, Kowalec P, Sobkowiak A. Maize Response to Low Temperatures at the Gene Expression Level: A Critical Survey of Transcriptomic Studies. Front Plant Sci. 2020;11:576941.
Li X, Hu H, Hu X, Wang G, Du X, Li L, et al. Transcriptome Analysis of Near-Isogenic Lines Provides Novel Insights into Genes Associated with Seed Low-Temperature Germination Ability in Maize (Zea mays L.). Plants. 2022;11:887.
Gupta PK. Single-molecule DNA sequencing technologies for future genomics research. Trends Biotechnol. 2008;26:602–11.
Ameur A, Kloosterman WP, Hestand MS. Single-Molecule Sequencing: Towards Clinical Applications. Trends Biotechnol. 2019;37:72–85.
Stark R, Grzelak M, Hadfield J. RNA sequencing: the teenage years. Nat Rev Genet. 2019;20:631–56.
Tan C, Liu H, Ren J, Ye X, Feng H, Liu Z. Single-molecule real-time sequencing facilitates the analysis of transcripts and splice isoforms of anthers in Chinese cabbage (Brassica rapa L. ssp. pekinensis). BMC Plant Biol. 2019;19:517.
Chao Y, Yuan J, Li S, Jia S, Han L, Xu L. Analysis of transcripts and splice isoforms in red clover (Trifolium pratense L.) by single-molecule long-read sequencing. BMC Plant Biol. 2018;18:300.
Yu X, Yu K, Chen B, Liao Z, Qin Z, Yao Q, et al. Nanopore long-read RNAseq reveals regulatory mechanisms of thermally variable reef environments promoting heat tolerance of scleractinian coral Pocillopora damicornis. Environ Res. 2021;195:110782.
Schaarschmidt S, Fischer A, Lawas LMF, Alam R, Septiningsih EM, Bailey-Serres J, et al. Utilizing PacBio Iso-Seq for Novel Transcript and Gene Discovery of Abiotic Stress Responses in Oryza sativa L. Int J Mol Sci. 2020;21:E8148.
Jiao Y, Peluso P, Shi J, Liang T, Stitzer MC, Wang B, et al. Improved maize reference genome with single-molecule technologies. Nature. 2017;546:524–7.
Laloum T, Martín G, Duque P. Alternative Splicing Control of Abiotic Stress Responses. Trends Plant Sci. 2018;23:140–50.
Zhao X, Li J, Lian B, Gu H, Li Y, Qi Y. Global identification of Arabidopsis lncRNAs reveals the regulation of MAF4 by a natural antisense RNA. Nat Commun. 2018;9:5056.
Cao Z, Zhao T, Wang L, Han J, Chen J, Hao Y, et al. The lincRNA XH123 is involved in cotton cold-stress regulation. Plant Mol Biol. 2021;106:521–31.
Waititu JK, Cai Q, Sun Y, Sun Y, Li C, Zhang C, et al. Transcriptome Profiling of Maize (Zea mays L.) Leaves Reveals Key Cold-Responsive Genes, Transcription Factors, and Metabolic Pathways Regulating Cold Stress Tolerance at the Seedling Stage. Genes (Basel). 2021;12:1638.
Thorpe GW, Reodica M, Davies MJ, Heeren G, Jarolim S, Pillay B, et al. Superoxide radicals have a protective role during H2O2 stress. Mol Biol Cell. 2013;24:2876–84.
Dawood M, Cao F, Jahangir MM, Zhang G, Wu F. Alleviation of aluminum toxicity by hydrogen sulfide is related to elevated ATPase, and suppressed aluminum uptake and oxidative stress in barley. J Hazard Mater. 2012;209–210:121–8.
Shen X, Xiao X, Dong Z, Chen Y. Silicon effects on antioxidative enzymes and lipid peroxidation in leaves and roots of peanut under aluminum stress. Acta Physiol Plant. 2014;36:3063–9.
Yin L, Mano J, Wang S, Tsuji W, Tanaka K. The involvement of lipid peroxide-derived aldehydes in aluminum toxicity of tobacco roots. Plant Physiol. 2010;152:1406–17.
Hussain HA, Men S, Hussain S, Zhang Q, Ashraf U, Anjum SA, et al. Maize Tolerance against Drought and Chilling Stresses Varied with Root Morphology and Antioxidative Defense System. Plants (Basel). 2020;9:E720.
Waqas MA, Wang X, Zafar SA, Noor MA, Hussain HA, Azher Nawaz M, et al. Thermal Stresses in Maize: Effects and Management Strategies. Plants (Basel). 2021;10:293.
Prasad TK, Anderson MD, Martin BA, Stewart CR. Evidence for Chilling-Induced Oxidative Stress in Maize Seedlings and a Regulatory Role for Hydrogen Peroxide. Plant Cell. 1994;6:65–74.
Del Rio D, Stewart AJ, Pellegrini N. A review of recent studies on malondialdehyde as toxic molecule and biological marker of oxidative stress. Nutr Metab Cardiovasc Dis. 2005;15:316–28.
Wang J, Guo J, Zhang Y, Yan X. Integrated transcriptomic and metabolomic analyses of yellow horn (Xanthoceras sorbifolia) in response to cold stress. PLoS One. 2020;15:e0236588.
Wang X, Yu C, Liu Y, Yang L, Li Y, Yao W, et al. GmFAD3A, A ω-3 Fatty Acid Desaturase Gene, Enhances Cold Tolerance and Seed Germination Rate under Low Temperature in Rice. Int J Mol Sci. 2019;20:E3796.
Huang X-S, Wang W, Zhang Q, Liu J-H. A basic helix-loop-helix transcription factor, PtrbHLH, of Poncirus trifoliata confers cold tolerance and modulates peroxidase-mediated scavenging of hydrogen peroxide. Plant Physiol. 2013;162:1178–94.
Huang X-S, Zhang Q, Zhu D, Fu X, Wang M, Zhang Q, et al. ICE1 of Poncirus trifoliata functions in cold tolerance by modulating polyamine levels through interacting with arginine decarboxylase. J Exp Bot. 2015;66:3259–74.
Schaberg PG, Minocha R, Long S, Halman JM, Hawley GJ, Eagar C. Calcium addition at the Hubbard Brook Experimental Forest increases the capacity for stress tolerance and carbon capture in red spruce (Picea rubens) trees during the cold season. Trees. 2011;25:1053–61.
Sun X, Li X, Zhu J, Huang N, Bian X, Li H, et al. Polyamines and ethylene metabolism during cold acclimation in zoysiagrass (Zoysia Japonica Steud.). Acta Physiol Plant. 2020;42:138.
Kovács Z, Simon-Sarkadi L, Szucs A, Kocsy G. Differential effects of cold, osmotic stress and abscisic acid on polyamine accumulation in wheat. Amino Acids. 2010;38:623–31.
Song Y, Diao Q, Qi H. Polyamine metabolism and biosynthetic genes expression in tomato (Lycopersicon esculentum Mill.) seedlings during cold acclimation. Plant Growth Regul. 2015;75:21–32.
Beike AK, Lang D, Zimmer AD, Wüst F, Trautmann D, Wiedemann G, et al. Insights from the cold transcriptome of Physcomitrella patens: global specialization pattern of conserved transcriptional regulators and identification of orphan genes involved in cold acclimation. New Phytol. 2015;205:869–81.
Pillai MA, Akiyama T. Differential expression of an S-adenosyl-L-methionine decarboxylase gene involved in polyamine biosynthesis under low temperature stress in japonica and indica rice genotypes. Mol Genet Genomics. 2004;271:141–9.
Tian A-G, Zhao J-Y, Zhang J-S, Gai J-Y, Chen S-Y. Genomic characterization of the S-adenosylmethionine decarboxylase genes from soybean. Theor Appl Genet. 2004;108:842–50.
Cuevas JC, López-Cobollo R, Alcázar R, Zarza X, Koncz C, Altabella T, et al. Putrescine is involved in Arabidopsis freezing tolerance and cold acclimation by regulating abscisic acid levels in response to low temperature. Plant Physiol. 2008;148:1094–105.
Zhao H, Zhang K, Zhou X, Xi L, Wang Y, Xu H, et al. Melatonin alleviates chilling stress in cucumber seedlings by up-regulation of CsZat12 and modulation of polyamine and abscisic acid metabolism. Sci Rep. 2017;7:4998.
Liu T, Kim DW, Niitsu M, Maeda S, Watanabe M, Kamio Y, et al. Polyamine oxidase 7 is a terminal catabolism-type enzyme in Oryza sativa and is specifically expressed in anthers. Plant Cell Physiol. 2014;55:1110–22.
Cona A, Rea G, Angelini R, Federico R, Tavladoraki P. Functions of amine oxidases in plant development and defence. Trends Plant Sci. 2006;11:80–8.
Angelini R, Cona A, Federico R, Fincato P, Tavladoraki P, Tisi A. Plant amine oxidases “on the move”: an update. Plant Physiol Biochem. 2010;48:560–4.
Yu Z, Jia D, Liu T. Polyamine Oxidases Play Various Roles in Plant Development and Abiotic Stress Tolerance. Plants (Basel). 2019;8:E184.
Sagor GHM, Inoue M, Kusano T, Berberich T. Expression profile of seven polyamine oxidase genes in rice (Oryza sativa) in response to abiotic stresses, phytohormones and polyamines. Physiol Mol Biol Plants. 2021;27:1353–9.
Gholizadeh F, Mirzaghaderi G. Genome-wide analysis of the polyamine oxidase gene family in wheat (Triticum aestivum L.) reveals involvement in temperature stress response. PLoS One. 2020;15:e0236226.
Li X, Wang G, Fu J, Li L, Jia G, Ren L, et al. QTL Mapping in Three Connected Populations Reveals a Set of Consensus Genomic Regions for Low Temperature Germination Ability in Zea mays L. Front Plant Sci. 2018;9:65.
Dionisio-Sese ML, Tobita S. Antioxidant responses of rice seedlings to salinity stress. Plant Sci. 1998;135:1–9.
Li H. Minimap2: pairwise alignment for nucleotide sequences. Bioinformatics. 2018;34:3094–100.
Love MI, Huber W, Anders S. Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2. Genome Biol. 2014;15:550.
Jin J, Tian F, Yang D-C, Meng Y-Q, Kong L, Luo J, et al. PlantTFDB 4.0: toward a central hub for transcription factors and regulatory interactions in plants. Nucleic Acids Res. 2017;45:D1040-5.
Yilmaz A, Nishiyama MY, Fuentes BG, Souza GM, Janies D, Gray J, et al. GRASSIUS: a platform for comparative regulatory genomics across the grasses. Plant Physiol. 2009;149:171–80.
Pertea G, Pertea M. GFF Utilities: GffRead and GffCompare. F1000Res. 2020;9:304.
Foissac S, Sammeth M. ASTALAVISTA: dynamic and flexible analysis of alternative splicing events in custom gene datasets. Nucleic Acids Res. 2007;35(Web Server issue):W297-299.
Kong L, Zhang Y, Ye Z-Q, Liu X-Q, Zhao S-Q, Wei L, et al. CPC: assess the protein-coding potential of transcripts using sequence features and support vector machine. Nucleic Acids Res. 2007;35(Web Server issue):W345-349.
Sun L, Luo H, Bu D, Zhao G, Yu K, Zhang C, et al. Utilizing sequence intrinsic composition to classify protein-coding and long non-coding transcripts. Nucleic Acids Res. 2013;41:e166.
Wang L, Park HJ, Dasari S, Wang S, Kocher J-P, Li W. CPAT: Coding-Potential Assessment Tool using an alignment-free logistic regression model. Nucleic Acids Res. 2013;41:e74.
Mistry J, Chuguransky S, Williams L, Qureshi M, Salazar GA, Sonnhammer ELL, et al. Pfam: The protein families database in 2021. Nucleic Acids Res. 2021;49:D412–9.
Ashburner M, Ball CA, Blake JA, Botstein D, Butler H, Cherry JM, et al. Gene ontology: tool for the unification of biology. The Gene Ontology Consortium Nat Genet. 2000;25:25–9.
Gene Ontology Consortium. The Gene Ontology resource: enriching a GOld mine. Nucleic Acids Res. 2021;49:D325–34.
Ogata H, Goto S, Sato K, Fujibuchi W, Bono H, Kanehisa M. KEGG: Kyoto Encyclopedia of Genes and Genomes. Nucleic Acids Res. 1999;27:29–34.
Kanehisa M. Toward understanding the origin and evolution of cellular organisms. Protein Sci. 2019;28:1947–51.
Kanehisa M, Furumichi M, Sato Y, Ishiguro-Watanabe M, Tanabe M. KEGG: integrating viruses and cellular organisms. Nucleic Acids Res. 2021;49:D545–51.
Yu G, Wang L-G, Han Y, He Q-Y. clusterProfiler: an R package for comparing biological themes among gene clusters. OMICS. 2012;16:284–7.
Elavarthi S, Martin B. Spectrophotometric Assays for Antioxidant Enzymes in Plants. In: Sunkar R, editor. Plant Stress Tolerance: Methods and Protocols. Totowa, NJ: Humana Press; 2010. p. 273–80.
Zeng Q, Ling Q, Hu F, Wu J, Yang Z, Qi Y, et al. Genotypic Differences in Growth and Antioxidant Enzyme Activities Under Cadmium Stress in Sugarcane. Bull Environ Contam Toxicol. 2017;99:607–13.
We thank all the laboratory members for their continuous technical advice and helpful discussions. We thank the BioMarker Technologies Corporation (Beijing, China) for whole-transcriptome libraries preparing and Oxford Nanopore Full-length Sequencing.
This work was co-funded by the special project of Guangdong Academy of Sciences, China (2020GDASYL-20200103073), Basic and Applied Basic Research Fund of Guangdong Province (2021A1515110745), National Natural Science Foundation of China (32072027), China Agricultural Research System (CARS201707) and the Laboratory of Lingnan Modern Agriculture Project (NZ2021014).
Ethics approval and consent to participate
All experimental studies on plants were complied with relevant institutional, national, and international guidelines and legislation.
Consent for publication
The authors declare that they have no competing interest.
Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.
Figure S1. Heatmap of transcripts corresponding to three GO function (ROS scavenging, response to heat and inositol)
Figure S2. Differential AS events corresponding to amine metabolism
Figure S3. Scatter plot revealed the correlation of gene expression characterized by full-length RNAseq and qRT-PCR
Table S1. Statistics of genes and transcripts
Table S2. Transcription factor encoding DETs that related to cold stress
Table S3. ROS scavening, responst to heat and inositol related transcripts
Table S4. Important DETs related to cold stress
Table S5. Primer used for qRT-PCR
DataSet1. Transcript count matrix and annotation
DataSet4. LncRNAs and target gene
About this article
Cite this article
Xuhui, L., Weiwei, C., Siqi, L. et al. Full-length transcriptome analysis of maize root tips reveals the molecular mechanism of cold stress during the seedling stage. BMC Plant Biol 22, 398 (2022). https://doi.org/10.1186/s12870-022-03787-3
- Cold stress
- Full-length transcriptome
- Seedling stage