RNA-Seq Analysis of Developing Grains of Wheat to Intrigue Into the Complex Molecular Mechanism of the Heat Stress Response

Heat stress is one of the significant constraints affecting wheat production worldwide. To ensure food security for ever-increasing world population, improving wheat for heat stress tolerance is needed in the presently drifting climatic conditions. At the molecular level, heat stress tolerance in wheat is governed by a complex interplay of various heat stress-associated genes. We used a comparative transcriptome sequencing approach to study the effect of heat stress (5°C above ambient threshold temperature of 20°C) during grain filling stages in wheat genotype K7903 (Halna). At 7 DPA (days post-anthesis), heat stress treatment was given at four stages: 0, 24, 48, and 120 h. In total, 115,656 wheat genes were identified, including 309 differentially expressed genes (DEGs) involved in many critical processes, such as signal transduction, starch synthetic pathway, antioxidant pathway, and heat stress-responsive conserved and uncharacterized putative genes that play an essential role in maintaining the grain filling rate at the high temperature. A total of 98,412 Simple Sequences Repeats (SSR) were identified from de novo transcriptome assembly of wheat and validated. The miRNA target prediction from differential expressed genes was performed by psRNATarget server against 119 mature miRNA. Further, 107,107 variants including 80,936 Single nucleotide polymorphism (SNPs) and 26,171 insertion/deletion (Indels) were also identified in de novo transcriptome assembly of wheat and wheat genome Ensembl version 31. The present study enriches our understanding of known heat response mechanisms during the grain filling stage supported by discovery of novel transcripts, microsatellite markers, putative miRNA targets, and genetic variant. This enhances gene functions and regulators, paving the way for improved heat tolerance in wheat varieties, making them more suitable for production in the current climate change scenario.

Introduction

Wheat (Triticum aestivum L.) is the most widely grown cereal crop, accounting for roughly 220 million hectares of agricultural land. Sixty percent of calorie requirement is met by wheat, maize, and rice while wheat alone meets about 25% of protein requirement of the world’s population (Sheoran et al., 2019; Singh et al., 2019). Several abiotic factors affect wheat production, including heat, drought, and salinity. In the present climate scenario, heat stress (HS) is a main constraint affecting the normal growth and productivity of crops, including wheat (Siebert and Ewert, 2014; Fahad et al., 2017), thereby posing a major threat to the global food and nutritional security. Elevated temperature above a threshold level negatively affects all the growth stages. The extent of yield loss depends on the sensitivity of the plant and the phenological/growth stage at which it is exposed to heat stress (Akter and Islam, 2017).

Among all the phenological stages in wheat, the reproductive stage (anthesis and seed development) is considered to be the most sensitive stage for HS (Kumar et al., 2018). It is estimated that during the grain filling stage, per degree temperature rise over and above 18–22°C, decreases grain filling rate and grain yield by 5% and 3–4%, respectively (Gupta et al., 2012). Based on this it was speculated that exposure of 32–38°C temperature during the grain filling stage may cause about 50% yield potential reduction in wheat (Awlachew et al., 2016).

HS tolerance is a very complex process that involves cellular pathway cross-talk. In the early stages of grain filling, HS affects the starch biosynthesis and assimilates translocation leading to poor yield and quality of wheat produce (Farooq et al., 2011; Kumar et al., 2013; Shamloo et al., 2017). Plants modulate their molecular machinery at the transcriptional, post-transcriptional, and epigenetic levels in order to regulate the differential expression of HS-related genes, photosynthetic machinery, cell membrane thermostability, induction of heat shock proteins (HSPs), stress hormones (ABA, ethylene), reactive oxygen species (ROS), compatible solutes (proline, glycine betaine, and sugars) accumulation, and non-coding RNAs synthesis (Farooq et al., 2011; Grover et al., 2013; Goswami et al., 2014; Lal et al., 2021; Zhao et al., 2021).

Many transcriptome studies exist for elucidating molecular mechanisms involving overlapping and distinct regulatory transcriptional mechanism of abiotic stresses response in model plants, as well as in some crops at a specific plant development stage (Coolen et al., 2016; Kang and Yeom, 2018; Cohen and Leach, 2019; Iquebal et al., 2019; Li et al., 2019; Kang et al., 2020; Rangan et al., 2020).

Wheat being a hexaploid with large and complex genome requires modern NGS-based approaches to elucidate tissue and growth stage-specific heat-responsive gene expressions. Wheat transcriptome profiling can reveal the differential gene expression, genome annotations, regulatory factors, molecular markers and expression quantitative trait loci (eQTLs) as well as their sequence variants, controlling the traits of importance (Morozova and Marra, 2008; Berkman et al., 2012; Lal et al., 2021). In wheat, RNA-seq has been adopted mainly to identify new and conserved transcripts associated with abiotic, biotic stress, and nutrient responsive genes (Kumar et al., 2015; Iquebal et al., 2019; Rangan et al., 2020). It is accurate, rapid, and comparatively cheaper and can be applied to non-model plant systems to extract novel genetic information (Unamba et al., 2015). De novo transcriptome assembly may be utilized to study the temporal and spatial gene expression of non-model organisms, which is an otherwise difficult task in the absence of complete genome sequence information (Grabherr et al., 2011).

The expression of genes in developing wheat grain is key in determining the ultimate composition of nutritional properties and yield (Henry et al., 2016). Moreover, response to HS is very complex and governed by numerous genes and their regulated expression (Mwadzingeni et al., 2016). Thus, it is of great importance to study the expression of various crucial genes contributing in grain filling and their regulation under HS. Transcriptomics-based expression analysis utilizing RNA-seq is a rapid, sensitive, and accurate way to identify crucial genes, their transcript level variations along with the regulating factors in wheat (Chu et al., 2021). In the present study, we have utilized NGS-based de novo assembly approach of RNA-seq to unearth the novel and conserved heat-responsive genes transcripts, transcription factors, gene regulatory network, and their complex interactions in response to heat stress during the initial phases of grain filling stage in heat-tolerant wheat variety K7903 (Halna). Further, the study revealed the discovery of putative molecular markers (SSRs, SNPs, and Indel markers), prediction of transcription factors (TFs), and microRNA and their putative targets. The findings of the study will enrich knowledge-based utilization of heat-tolerant varieties in the crossing block for genetic advancement of wheat.

Materials and Methods

Plant Material, Heat Treatment, and Sample Collection

Healthy seeds of a heat stress-tolerant late sown wheat variety, K7903 (Halna), were procured from the germplasm unit of ICAR-Indian Institute of Wheat and Barley Research (IIWBR) Karnal, India. The seeds were surface sterilized with 0.1% HgCl_2, followed by three subsequent washing with sterile distilled water. Four wheat seeds were sown in 10-cm pots in the soil in the net house (15 h light, 9 h dark; 40% RH). A total of 18 such pots were sown and labeled. At seven days post-anthesis (DPA), that is, approximately 3-month-old wheat plants, 9 out of 18 pots were shifted carefully to Temperature Controlled Phenotyping Facility (TCPF), equipped with automatic sensing of the outside environmental temperature and maintained inside temperature up by +5°C (Sharma et al., 2019). Other environmental factors including humidity, light intensity etc. were maintained similar to outside. The remaining nine pots in the net house served as control. The main tiller spike of each wheat plant in both, control (normal temperature) and TCPF (under HS) conditions was tagged, and the spikelets (wheat grains) were harvested from the middle of the spike at 0-, 24-, 48-, and 120-h intervals after starting from 7 DPA from three main tiller spike and pooled. At each sampling point, four spikelets from the middle of spikes from all the plants were harvested and frozen immediately in liquid nitrogen in a 50 ml sterile plastic tube. The samples taken from controlled or normal conditions were labeled as C1, C2, C5, and C7, while samples from heat stress conditions were labeled as R1, R2, R5, and R7 with two biological replicates.

RNA Isolation, Library Construction, and Sequencing

The developing wheat grain samples were carefully taken out of the liquid nitrogen without thawing and grounded to a very fine powder using sterile pestle and mortar in liquid nitrogen. Thereafter, total RNA was extracted using RNeasy RNA isolation kit (Qiagen, Germany) according to the manufacturer’s instructions. DNA was removed by digestion with RNase-free DNase (Qiagen, Germany), and RNA was purified and concentrated using RNeasy column (Qiagen, Germany). RNA quality was evaluated by 1% agarose gel electrophoresis for 28 S/18 S rRNA band intensity (2:1) and Agilent 2,100 Bioanalyser and quantified by using NanoDrop Spectrophotometer (Thermo Fisher Scientific, United States). RNA samples were immediately stored in −80°C and subsequently sent for further processing, library construction, and sequencing to AgriGenome Labs Pvt. Ltd., Hyderabad, Telangana, India.

Pre-processing, de novo Assembly, and Expression Analysis

Paired-end (100 bp × 2) read of control and treated samples were generated using Illumina Hiseq 2,500. The quality assessment of reads was performed using FastQC tool.¹ Cutadapt version 1.8.1 (Martin, 2011) and Sickle version 1.33² were used for removal of adaptors and low-quality reads. Reads having average quality score ≤ 30 were removed. Further, rRNA data from samples were removed with the SILVA database (Quast et al., 2012).³ Transcriptome assembly of wheat was performed using trinity (at default k-mer length of 25), which uses de Bruijn graphs and dynamic programming algorithm (Haas et al., 2013). Finally, the CD-HIT-EST package was employed to remove redundant sequences from trinity assembly (Fu et al., 2012). For the expression analysis, Bowtie 2 (Langmead et al., 2009) was used for alignment of reads, followed by RSEM (RNA-Seq by Expectation–Maximization; Li and Dewey, 2011) for calculation of expression values in the form of FPKM (Fragments per kilobase million reads). RSEM implements improved quantification method having superior accuracy than other methods. Differential gene expression analysis between two biological replicates of each pair was performed using DESeq 1.16.0 program (Anders and Huber, 2010) with parameters, such as Log₂FC ± 2 and adjusted value of p <0.05.

Homology Search, Functional Characterization, and Gene Ontology

Blastx⁴ program was used for sequence similarity search against NCBI non-redundant protein database with expected threshold e-value 10⁻⁵ and identity cutoff 40% (Camacho et al., 2009). Contig Annotator Pipeline (CANoPI) was used for gene ontology, pathways analysis of differential expressed genes, and de novo transcriptome assembly transcripts whose length was ≥200 bp and FPKM ≥1.0. We also identified transcription factors from the sets of differentially expressed genes using PlantTFDB 4.0 (Jin et al., 2016). The miRNA target prediction from differential expressed genes was performed by psRNATarget server (Dai and Zhao, 2011).

Identification and Validation of Simple Sequences Repeat Markers

Simple Sequences Repeats were identified from de novo transcriptome assembly of wheat using perl script of MISA (MIcroSAtellite identification tool) with default parameters, such as ten repeating units for mononucleotides, six repeating units for dinucleotide, and five for tri, tetra, penta, and hexanucleotides (Thiel et al., 2003). Also, the Primer3 tool was used for generating three sets of primers of identified markers (Untergasser et al., 2012).

From the mined SSRs markers, fifteen were randomly selected for validation in 20 different wheat genotypes panel. Genomic DNA was extracted from seedlings by CTAB method. PCR amplification was carried out in BioRad S1000 Thermocycler in 25 μl reaction mixture containing 2.5 μl 10 × buffer,0.5 μl of 10 mM dNTPs, 0.5 μl of 10 μM each reverse and forward primer, 0.5 μl of Taq polymerase (5 U/μl), 60 ng template DNA and nuclease-free water. PCR amplification was performed using the program: 94°C for 5 min for initial denaturation, 30 cycles of 1 min at 94°C, 50 s at 55–63°C, 72°C for 1 min, and a final amplification of 5 min at 72°C. PCR products were analyzed using 3% agarose gel electrophoresis.

Identification of Single Nucleotide Polymorphism and Indels

Single Nucleotide Polymorphism (SNP) and Indels (insertion/deletions) were identified from a sample of prior treatment, that is, C1 and R1 using two references, viz., de novo transcriptome assembly and wheat genome release-31.⁵ Burrows-Wheeler Aligner (BWA; Li and Durbin, 2009) and SAM tools (Li et al., 2009) was used for alignment and calling of variants. Filters, such as read depth coverage >10, quality score > 30, and at least 2 SNPs in 50 bps, were applied to identify significant SNPs and Indels.

Validation of RNA-seq Analysis by q-PCR

RNA-seq results were confirmed by performing q-PCR of 10 randomly selected differentially expressed transcripts (DETs). cDNA preparation was performed using Superscript® III First-Strand Synthesis System (Invitrogen, UK) following manufacturer’s instructions. Prior to q-PCR, the PCR conditions were optimized by performing successful routine PCR. qPCR mixture contained 1 μl of cDNA (10 ng/μl), 6.0 μl of 2 × SYBR Green Master Mix (Thermo Scientific), and 1 μl (10 μM) each of forward and reverse primer of the selected DET and nuclease-free water to a final 12.0 μl total reaction volume. The quantitative reaction was done on Bio-Rad CFX96™ Real-Time PCR System (Bio-Rad, United States). The qPCR program consisted of 95°C for 5 min, then 40 cycles of 94°C for 15 s, 58°C for 30 s, and 72°C for 30 s and a final melt curve step from 65°C to 95°C with a rise of 0.5°C for 5 s. The relative expression level study of target transcripts was done by comparative 2^−ΔΔCt method (Livak and Schmittgen, 2001) using GAPDH as endogenous control.

Results

Pre-processing and de novo Assembly

In this study, a total of 16 samples of control and treated data were generated at four different developmental crop stages viz., 7 days post-anthesis, 0, 24, 48, and 120 h after treatment. Sample names C1, C2, C5, and C7 belonged to the control condition, whereas R1, R2, R5, and R7 represented treated samples. Both control and treated samples had biological replicates. The generated transcriptome data at different developmental crop stages were pre-processed for quality assessment. After trimming and removal of rRNA reads, 253,046,594 (20.5%) poor quality reads were dropped and finally, 979,825,192 (79.47%) high-quality reads were retained (from overall generated data) for de novo transcriptome assembly and further downstream analysis. Trinity generated 2,302,239 transcripts with a minimum read length of 200 and GC content of 50.24%. Finally, 1,696,570 transcripts were obtained after removal of redundant sequences with N50 value ≥473 bp, GC content 50.48%, and transcript lengths ranging between 200 to 18,870 bp (Figure 1).

FIGURE 1

Figure 1. (A) Length distribution, (B) GC content of assembled transcript in control and treated wheat sample at four different developmental stages (at 0, 24. 48 and 120 h of heat stress treatment).

Functional Characterization and Gene Ontology

For homology search and functional annotation, we filtered transcripts having sequence length ≥ 200 bp and FPKM ≥1.0 to finally arrive to a total of 319,879 transcripts for further analyses. The criteria of e-value 10⁻⁵ and % identity cutoff 40% were employed for searching known sequences against NCBI non-redundant protein database. Out of 319,879 transcripts, 152,842 had similarities with other known proteins in the database. Blastx results showed top hits with Aegilops tauschii (49204), Triticum urartu (42504), and Hordeum vulgare (Figure 2). The identified proteins were further annotated using UniProt database, pathways database, PROSITE, InterProscan, TIGRFAMS, Pfam, PANTHER, and OrthoDB. Finally, a total of 105,127 transcripts were found to have significant Blastx hit and UniProt information (Supplementary File 1).

FIGURE 2

Figure 2. Bar graph representing the species distribution of annotated DEGs. Maximum number of annotated DEGs was represented by Aegilops tauschii (49204) followed by T. uratu (42504). Triticum aestivum falling at the fourth place had approximately 10,000 annotated DEGs.

The complete annotations were categorized into Biological Processes (BP), Molecular Functions (MF), and Cellular Components (CC), having 4,345, 3,873, and 2,516 GO terms, respectively. In BP, we found full transcripts represented in DNA integration (5,688 transcripts), while nucleic acid binding (8379) an integral component of membrane (7433) in MF and CC, respectively (Figure 3; Supplementary File 2).

FIGURE 3

Figure 3. Bar graph depicting gene ontology (GO) classification and enrichment of top 20 GO terms represented under biological (DNA integration: 5688 transcripts), molecular (nucleic acid binding: 8379 transcripts), and cellular components (an integral component of membrane: 7433).

Abundance Estimation and Differential Expression Analysis

Reads of all the sixteen samples were mapped separately on the de novo transcriptome assembly, containing 1,696,570 transcripts and ~ 62.48–95.39% (average ~ 83.79%) of reads mapped. Also, expression values of each sample were calculated in the form of FPKM, to identify differentially expressed genes (Figure 4).

FIGURE 4

Figure 4. Graphical representations of DEGs with FPKM value ≥1 of 16 libraries prepared in control and treated sample at four different developmental stages. Samples contained 1,696,570 transcripts and ~ 62.48–95.39% (average ~ 83.79%) of reads mapped.

Differential expression analysis was performed and categorized into 3 sets. Set 1 contained the DEGs sets of control samples, such as C1 vs. C2, C2 vs. C5, and C5 vs. C7, whereas set 2 contained treated samples, such as R1 vs. R2, R2 vs. R5, and R5 vs. R7. In the 3^rd set, DEGs were of control and treated samples, such as C1 vs. R1, C2 vs. R2, C5 vs. R5, and C7 vs. R7. The DEGs were identified from all the paired sets. In total, 109,455 DEGs were obtained, out of which 95,966 were upregulated and 13,489 were downregulated (Table 1).

TABLE 1

Table 1. Detailed list of differentially expressed genes obtained in all the three sets.

In control samples, that is, set 1, no common DEGs were found in all the three comparisons, such as C1 vs. C2, C2 vs. C5, and C5 vs. C7. But, several common DEGs were obtained in the comparison set C1 vs. C2 and C2 vs. C5, that is, 12,332, while 1 and 10 DEGs were common in comparison “C1 vs. C2 and C5 vs. C7” and “C2 vs. C5 and C5 vs. C7,” respectively. A total of 806, 34,670, and 86 DEGs were found unique in C1 vs. C2, C2 vs. C5, and C5 vs. C7, respectively. Similarly, in Set 2 (treated samples), no common DEGs were obtained. In Set 3, which contains the DEGs of control vs. treated samples, 23 DEGs were found common in between C2 vs. R2 and C5 vs. R5 while 48,904, 23, and 158 DEGs were unique in C2 vs. R2, C5 vs. R5, and C7 vs. R7, respectively (Figure 5).

FIGURE 5

Figure 5. Venn diagram representing two maximum common DEGs in set 1(A)-{(C₂ vs. C₅) vs. (C₅ vs. C₇)} and set 3(B)-{(C₂ vs. R₂) vs. (C₅ vs. R₅)}. 10 DEGs were common in set 1(A) and 23 DEGs were common in set 3(B).

Functional Characterization of DEGs

Functional characterization of DEGs was performed using Blastx against NCBI nr database. For identification of transcriptional factors, PlantTFDB v4.0 was used. Transcripts to the tune of 43,861, had hits with other known proteins in the databases (Table 2).

TABLE 2

Table 2. Number of DEGs in each set that showed similarity with other known proteins in the various databases.

MicroRNA Target Prediction

The miRNA target prediction was made by psRNATarget server for DEGs sets against 119 mature miRNAs of T. aestivum. Out of 10 sets, only differential expressed genes of three sets were involved in miRNA target prediction, that is, C1 vs. C2, C2 vs. C5, and C2 vs. R2. In set C1 vs. C2, we found 24 miRNAs that targeted 65 transcripts. Maximum 12 transcripts were targeted by “tae-miR1117,” followed by 5 transcripts each in “tae-miR1130b-3p” and “tae-miR1134.” In set C2 vs. C5, 29 miRNAs were found targeting 125 transcripts with a maximum number of transcripts targeted by “tae-miR1117,” that is, 51, followed by 9 and 6 transcripts in “tae-miR1128” and “tae-miR1130b-3p,” respectively. Lastly, in DEG set C2 vs. R2, 31 miRNAs were found targeting 154 transcripts and maximum hits were targeted by “tae-miR1117,” that is, 55, followed by 9 transcripts each in “tae-miR1130b-3p”and “tae-miR1131” (Supplementary File 5). Twenty miRNAs were found common in C1 vs. C2, C2 vs. C5, and C2 vs. R2, whereas 5 and 3 miRNAs were found unique in C2 vs. C5 and C2 vs. R2, respectively (Figure 6).

FIGURE 6

Figure 6. Venn diagram of miRNAs found in DEG sets C₁ vs. C₂, C₂ vs. C_5, and C₂ vs. R_2. Among all the comparison sets, 20 miRNAs were found to be common.

In silico Identification of SSR Markers

A total of 98,412 simple sequences repeats (SSRs) were identified from 1,696,570 transcripts of wheat de novo transcriptome assembly, while 4,261 repeats were present in compound formation. Out of 98,412 markers, there was abundance of mononucleotides, that is, 37,716, followed by di- (27697), tri- (29659), tetra- (2773), penta- (381), and hexa (186)nucleotide repeats (Supplementary File 6). In dinucleotide repeats, we found the maximum number of repeats was GA, that is, 4,601, followed by AG (3247) and CT (3101). In the case of trinucleotides, the maximum was GGC repeats, that is, 1799, followed by GCG (1487) and GAG (1393; Table 3).

TABLE 3

Table 3. Details of SSR markers obtained from de novo transcriptome assembly as well as unique differential expressed genes of all the stages.

Validation of SSR Markers

Twenty highly diverse randomly selected wheat genotypes were used for validation of SSR. PCR primers were designed for randomly selected 20 SSR loci (Table 4), and PCR was amplified using genomic DNA as a template (Figure 7). Such SSR markers are extremely important in crop improvement programs for heat tolerance and wheat breeding not an exception.

TABLE 4

Table 4. Details of randomly selected SSR markers (thirteen dinucleotide and two trinucleotide repeats) along with their primer sequence and melting temperature ranging from 57.17°C to 60.99°C.

FIGURE 7

Figure 7. Gel pictures depicting the expression of identified genic SSRs in highly diverse wheat genotypes (1-UP2425, 2-WR544, 3-SONARA64, 4-K7903 (HALNA), 5-WH730, 6-DBW14, 7-RAJ4014, 8-DBW71, 9-AKAW1071, 10-RAJ3765, 11-NIAW34, 12-NW1014, 13-K9465, 14-K9644, 15-HD2733, 16-K9107, 17-PBW502, 18-DBW17, 19-RAJ4083, and 20-WH542) for heat stress.

Identification of Single Nucleotide Polymorphism and Indels

Single Nucleotide Polymorphism (SNPs) and insertion/deletion (Indels) were identified from pooled reads of C1 and R1 against de novo transcriptome assembly of wheat and wheat genome Ensembl version 31. We found a total of 107,107 variants against the de novo transcriptome assembly of wheat, containing 80,936 SNPs and 26,171 Indels. Maximum number of variants (38) were found in transcript ID “TRINITY_DN523265_c0_g1_i3 (WD repeat and FYVE domain-containing protein 3)” and interesting all 38 were SNPs (Supplementary File 7). Similarly, we obtained 110,787 variants against the wheat genome, which contains 99,791 SNPs and 10,996 Indels. Top variants identified from chromosome number 3B, that is, 12,977 (11,589 SNPs and 1,388 Indels), followed by 5,788 and 5,771 variants in chromosome 2B and 2A, respectively (Supplementary File 7).

Validation of RNA-seq Analysis by q-PCR

The RNA-seq results were verified by performing q-PCR of ten randomly selected differentially expressed transcripts (DETs; Table 5). The qPCR expression pattern results of these DETs were found to follow the same trend as the RNA-seq expression pattern results, confirming the accuracy of the RNA-seq data obtained in this investigation (Figure 8).

TABLE 5

Table 5. Differentially Expressed Transcripts (DETs) and primers used for qPCR.

FIGURE 8

Figure 8. Comparative expression analysis of selected DEGs using RNA-seq and qPCR. Data are means of two independent biological replicates (p ≤ 0.05, n = 2). Error bars represent the means ± SD (n = 2). Accuracy of the RNA-seq data getting confirmation with the qPCR expression pattern results.

Discussion

Global wheat production is highly vulnerable to climate change and environmental factors, such as heat stress and water deficit, which have severe adverse impact on the wheat crop productivity and quality (Pareek et al., 2020; Rangan et al., 2020). About 57% wheat cultivated area in developing nations is affected by heat stress, thus it can be considered as primary abiotic factor threatening world wheat production and food security (Kosina et al., 2007). Wheat is more sensitive to the heat stress during reproductive (anthesis and post-anthesis) stages and grain filling period which is very crucial for sustaining overall crop yield and quality (Rangan et al., 2020). Heat stress response is very complex mechanism that involves interplay and modulation of numerous biochemical and molecular machinery for sustaining yield and quality (Mwadzingeni et al., 2016). It is of great importance to study the expression of various crucial genes and their regulation in wheat, expressed during grain filling stage under the heat stress. RNA-seq-based next-generation sequencing (NGS) technology provides comparatively economic and efficient way of identifying crucial genes, their transcript level variations along with the regulating factors with high accuracy (Chu et al., 2021).

In the present work, we identified a total of 115,656 genes in heat stress-tolerant wheat genotype, that is, Halna (K7903), exposed to heat stress. Result indicated differential expression of 109,455 DEGs in different sets, out of which 95,966 were upregulated and 13,489 were downregulated (Table 1). The detailed information of DEGs with log2 fold change, adjusted value of p, and their functional characterization is presented in Supplementary File 3. Supplementary File 4 provides the information of transcriptional factors which were obtained in each set. A total of 43,861 transcripts had hits with other known proteins in the databases (Table 2). These DEGs narrowed down to 309 which were involved in many critical processes viz., signal transduction, starch synthetic pathway, antioxidant pathway, and heat stress-responsive conserved and uncharacterized putative genes that play an essential role in maintaining the grain filling rate at high temperature. Categorization of the annotations found their place into Biological Processes (BP), Molecular Functions (MF) and Cellular Components (CC), corresponding to 4,345, 3,873, and 2,516 GO terms, respectively. DNA integration (5,688 transcripts), nucleic acid binding (8379), and integral membrane components (7433) were observed in BP, MF, and CC, respectively (Figure 3; Supplementary File 2). Differential gene expression of transcripts functionally involved in diverse cellular functions is an indicator that the heat stress response in wheat is very complex and it simultaneously regulates different biochemical and molecular pathway sustaining the survival and developmental processes in wheat. Among differentially expressed transcripts, TRINITY_DN521405_c2_g2_i1 (Avenin-like a4) which is cysteine-rich storage proteins expressed during early stages of seed development (3–22 DPA) exclusively in developing endosperm of wheat and its related plant species. Expression of this gene during early grain development plays very significant role in ensuring wheat processing quality (Ma et al., 2013; Chen et al., 2016; Zhang et al., 2018). Jiang et al. (2012) conducted a proteomics-based study and reported that the avenin-like proteins expression is highly upregulated under drought stress and may have an essential role in seed quality. Similarly, Zhang et al. (2020) suggested that avenin-like proteins are multifunctional including role during disease resistance. In our study, we found very high expression of avenin-like proteins which suggests their important role during seed development/grain filling under heat stress. To the best of our knowledge, this is the first report of expression study of avenin-like protein under heat stress in wheat.

TRINITY_DN483169_c0_g1_i6 (Alpha-amylase inhibitor 0.19) are amylase inhibitors, found in cereal grains, including wheat and some legume seeds, which inhibit the target α-amylase enzymes thus protect the starchy seed reserves from degradation into simpler oligosaccharides (Petrucci et al., 1976; Franco et al., 2000; Dombrowski, 2003). Their role in protecting plant amylase inhibitors in biotic stress, particularly insect pests, is well established (Feng et al., 1996; Franco et al., 2000; Rane et al., 2020). Alpha-amylase inhibitors inhibit grain starch degradation and protect the plant against oxidative stress (Hajheidari et al., 2007). Reports suggest that the differential expression of amylase inhibitor at protein level in different cultivars during different grain filling stages play an essential role during grain filling under abiotic stress (Jiang et al., 2012; Zhang et al., 2014). The upregulation of α-amylase inhibitors transcripts under heat stress in our study also confirms their role in heat stress tolerance in wheat during grain filling.

TRINITY_DN518599_c1_g1_i2 (Aspartic proteinase oryzasin-1) are proteolytic enzymes with two highly conserved aspartate residues and play a crucial role in stress responses in plants (Alam et al., 2014). APs has been reported to play an important role in protein processing in response to pathogen attack, programmed cell death, water stress, drought, and other environmental stimuli in many plant species at different developmental stages in various plant organs (Runeberg-Roos et al., 1994; De Carvalho et al., 2001; Guevara et al., 2005; Mendieta et al., 2006; Simões et al., 2007; Muñoz et al., 2010; Yao et al., 2012; Rawlings and Salvesen, 2013). Asakura et al. (1995) reported that a rice aspartic proteinase, oryzasin 1 transcripts level was much more abundant from just after flowering stage up to final seed formation suggesting its active role in seed storage protein processing. Yao et al. (2012) reported overexpression of APs under water-deficient conditions in Zea mays and their role in drought avoidance through ABA-dependent signaling pathway. They also reported that transgenic Arabidopsis overexpressing gene VlAP17, encoding a Group C aspartic protease, enhanced tolerance to salt and drought stress during seed germination, seedling, and maturation. The transgenic Arabidopsis VlAP17 overexpression also increased ABA levels, a reduction in average stomatal aperture size, and elevated expression levels of stress response genes involved in the ABA-dependent pathway well as higher activities of several antioxidases: superoxide dismutase, catalase, and peroxidase. Recently, Sebastián et al. (2020) reported that APA1 encoding AP in known to be highly expressed under water stress in response to salinity, cold, and drought stress. Our study also found that AP expression increases in response to heat stress in plant seed tissue, suggesting their crucial role in response to heat stress.

TRINITY_DN326896_c0_g1_i1 (Elongation factor 1-alpha) Elongation factor 1α- is a multifunctional protein Transcript elongation factors (EFs) play essential role in mediating critical cellular processes related to cellular growth, proliferation, and cell differentiation by interacting with other cellular proteins (Bukovnik et al., 2009; Zheng et al., 2014). Its high expression during heat stress conditions in animals and plants has been reported and thus suggested its essential role in survival under stress conditions (Buckley, 2006; Shamovsky et al., 2006; Bukovnik et al., 2009). They further mentioned its role in heat stress in wheat as accumulation was high in cultivars with better heat tolerance. Zheng et al. (2014) characterized a transcript elongation factor gene in wheat through expression based association analysis, comparing near-isogenic lines and by overexpression in Arabidopsis and reported its role in regulating yield-related traits involved in vegetative growth and reproductive development. In our study, the upregulated expression of the EF also suggests its role during the grain filling stage in wheat. Metabolic proteins of wheat, mainly involved in glycolysis, carbohydrate metabolism, and the stress response, show differential expression profiles in response to elevated temperature conditions (Laino et al., 2010).

TRINITY_DN490415_c0_g1_i8 (α-gliadin) was found to be upregulated in the present study. Yang et al. (2011) has reported that the expression of many storage proteins is under heat or water deficit, viz., α-gliadin, γ-gliadin, low molecular weight glutenin, and globulins in the seeds altered. Further, the abundance of glidins was found to increase under heat stress from anthesis up to 10 DPA. Du Pont et al. (2006) also studied that amount of α-gliadins increases in response to elevated temperature conditions during endosperm development in wheat. Our results are in agreement of the above studies.

TRINITY_DN425597_c0_g2_i2 (heat shock protein 83-like) belongs to the heat shock protein group, which are well known for their protective role in response to elevated temperatures in almost all organisms, including plants. These proteins perform multiple functions, such as protein folding, translocation across membranes, facilitation of protein–protein interaction, and preventing protein aggregation and regulation of synthesis of other stress-related genes during heat stress conditions (Åkerfelt et al., 2010; Roy et al., 2019). Conner et al. (1990) reported the heat-inducible expression of AtHS83 in Arabidopsis thaliana. Felsheim and Das (1992) reported the differential expression of HSP83 member genes in response to heat stress and photoperiod duration in cotyledons. Disruption of genes encoding HSP83 has been shown to have detrimental effects on the growth and survival in yeast at high temperature (Borkovich et al., 1989). Our results also agree with the above study and we can conclude that heat shock proteins also have an essential role in wheat during heat stress.

TRINITY_DN385005_c2_g1_i2 (Zinc transporter), a transmembrane zinc transporter is reported to confer heat tolerance during grain filling stage by maintaining optimum cellular zinc concentration required for the proper working of an enzyme related to stress tolerance carbohydrate metabolism and normal seed filling process under heat stress (Kobayashi et al., 1998). They also reported zinc transporter’s role in normal flowering and fruit development. Moreover, the Zn is crucial for plant growth and development as it is required by various enzymes involved in carbohydrate metabolism, cell membrane stability and integrity, protein synthesis, auxin synthesis regulation, ribosomal stabilization, and cytochrome synthesis (Marschner, 1995). Zn also plays role in the stress tolerance mechanism of plants by regulating gene expression of various stress-related genes (Cakmak, 2000). Homeostasis in cellular zinc concentration improves the plant tolerance to heat and salt stress in various plant species by improving the water uptake and transport (Peck and McDonald, 2010; Tavallali et al., 2010; Disante et al., 2011). Sarwar et al. (2019) reported the role of exogenously applied zinc and other micronutrients in heat stress alleviation in cotton. Specific Zn transporters maintain the optimum cellular concentration through influx, efflux, and compartmentalization. Thus, we can also suggest that upregulation in gene expression of zinc transporter 6 in our study contribute to heat tolerance in wheat.

Therefore, the transcripts validated randomly using q-PCR in the present study are multifunctional. Most of them are reported to be overexpressed in response to multiple stress conditions including heat, drought, pathogen, or insect attack. Some genes are very crucial for maintaining seed quality parameter viz., number of different seed-specific protein content and carbohydrate metabolism specifically under heat stress. Some transcripts also play protective role during heat stress ensuring proper protein folding, protein translocation, and enzyme functioning. Other transcripts play significant role in regulating seed developmental processes under heat stress by regulating molecular processes including gene expression, protein synthesis etc. Further, these transcripts help in better understanding of complex molecular and biochemical mechanism of heat stress response in wheat.

In addition, we also identified SSRs from wheat de novo transcriptome assembly and characterized them. The repeat pattern of identified SSRs indicated an abundance of mononucleotides followed by di, tri, tetra, penta, and hexanucleotide repeats (Supplementary File 6). In dinucleotide repeats, we found the maximum number of repeats as GA followed by AG and CT. In the case of trinucleotides, GGC repeats were most abundant followed by GCG and GAG (Table 3). Additionally, we found a total of about one hundred thousand variants against the de novo transcriptome assembly of wheat, containing SNPs and Indels. Thirty-eight variants were reported in transcript ID “TRINITY_DN523265_c0_g1_i3 representing WD repeat and FYVE domain-containing protein 3” incidentally all of them being SNPs (Supplementary File 7). Similarly, we obtained one hundred and ten thousand variants against the wheat genome, containing SNPs and Indels. Maximum variants were identified from chromosome number 3B followed by chromosomes 2B and 2A (Supplementary File 7). Being high throughput and high density in nature coupled with widespread distribution throughout the genome, SNP markers are popular and important for their utilization in population genetic diversity studies, gene mapping, and high-density linkage map construction (Unterseer et al., 2014; Carrillo-Perdomo et al., 2020). Therefore, the genic-SSR markers and SNP markers reported in this study along with previously developed markers from other laboratories will be a powerful resource for molecular taxonomic studies and construction of a reference molecular map of the wheat genome. Since such genic-SSR and SNP markers represent the gene-rich regions of the genome, many of these can be utilized in marker-assisted breeding (Dutta et al., 2011). Also because of conservation of genic sequences, these markers have a greater chance of transferability across species, compared to genomic SSR markers which show high polymorphism but are less conserved between species (Wang et al., 2021).

The miRNAs are known as key players in regulating plant gene expression via targeting important genes involved in growth and developmental during heat stress (Xin et al., 2010; Kumar et al., 2015). DEGs sets using 119 mature miRNAs were subjected to miRNA target prediction for T. aestivum. Out of total ten sets, only differential expressed genes of three sets were involved in miRNA target prediction. miRNAs were found with a maximum number of transcripts targeted by “tae-miR1117.” Also in DEG set C2 vs. R2, 31 miRNAs were found with maximum hits were targeted by “tae-miR1117” (Supplementary File 5). Common miRNAs in C1 vs.C2, C2 vs. C5, and C2 vs. R2 were 20 in number while there were only 5 and 3 miRNAs that were found unique in C2 vs. C5 and C2 vs. R2, respectively (Figure 6). These tiny miRNAs are very crucial to modulate the physiological, biochemical, and molecular responses by targeting key genes during heat stress.

Taking together, the present study enriches our current understanding of known heat response mechanisms during the grain filling stage in wheat supported by discovery of novel transcripts, microsatellite markers, putative miRNA targets, and genetic variant. Further work on functional validation of these genes, miRNAs, SSRs, and SNPs could pave the way for developing improved heat tolerance in wheat varieties, making them more suitable for production in the current climate change scenario.

Conclusion

Climatic change has affected the quantity and quality of major food crops and thus poses a significant threat to global food security. Initial few days of grain filling are crucial for the net yield of the crop. Heat stress during this period drastically reduces the overall yield. Our de novo transcriptomics-based study of the heat stress-tolerant wheat genotype K7903 (Halna) reveals that the upregulation of key genes during initial grain filling stages might be playing crucial role in heat tolerance mechanism. K7903 (Halna) might have genetic and epigenetic (miRNA expression) control mechanism which might be contributing toward avoidance of the detrimental effect of heat stress and thereby maintaining yield and grain quality. The significant findings could further be functionally validated using reverse genetic tools either in wheat or any other model plants which might assist molecular breeders to design suitable strategies to develop heat-tolerant wheat genotypes. In this study, we identified heat stress-responsive miRNA and their target genes, which have critical role in gene regulation during grain filling under heat stress conditions. This valuable information will enrich our knowledge about involvement of different key genes and their expression pattern regulating the grain filling process when exposed to high temperature conditions. This knowledge can be further utilized in identification, characterization, and breeding strategies to develop heat stress-tolerant wheat varieties.

Data Availability Statement

The datasets presented in this study can be found in online repositories. The names of the repository/repositories and accession number(s) can be found at: https://www.ncbi.nlm.nih.gov/, PRJNA813970, SAMN26531246, SAMN26531247, SAMN26531249, SAMN26531250, SAMN26531251, SAMN26531252, SAMN26531253, SRR18273340, SRR18273339, SRR18273338, SRR18273337, SRR18273336, SRR18273335, SRR18273334, SRR18273333.

Author Contributions

RT and JD conceived the theme of study. SP, SJ, UA, MI, and RS analyzed the data. SP, JD, MI, PS, RS, and RT drafted the manuscript. SP, NR, and SS conducted the field work. OG and DK were instrumental in formatting the tables and graphs. RT, JD, AR, GS, DK, and MI provided overall guidance and edited the manuscript. All authors read and approved the final version of the manuscript.

Funding

This work was supported by CABin grant (F. No. Agril. Edn.4–1/2013-A&P) of ICAR-IASRI, New Delhi in alignment with in-house project of ICAR-IIWBR, Karnal.

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.

Acknowledgments

The authors are thankful to the Director, ICAR-Indian Institute of Wheat and Barley Research, Karnal for providing research facility. Authors are thankful to Indian Council of Agricultural Research, Ministry of Agriculture and Farmers’ Welfare, Govt. of India for providing financial assistance in the form of CABin grant (F. no. Agril. Edn.4-1/2013-AandP) as well as Advanced Super Computing Hub for Omics Knowledge in Agriculture (ASHOKA) facility at ICAR-IASRI, New Delhi, India. The authors are also grateful for the support and guidance received from the department of Biotechnology, Chaudhary Devi Lal University, Sirsa, Haryana, India.

Supplementary Material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fpls.2022.904392/full#supplementary-material

Footnotes

1. ^http://www.bioinformatics.babraham.ac.uk/projects/fastqc/

2. ^https://github.com/najoshi/sickle

3. ^https://www.arb-silva.de/

4. ^https://blast.ncbi.nlm.nih.gov/Blast.cgi

5. ^ftp://ftp.ensemblgenomes.org/pub/plants/release-31/fasta/triticum_aestivum/dna/

References

Åkerfelt, M., Morimoto, R. I., and Sistonen, L. (2010). Heat shock factors: integrators of cell stress, development and lifespan. Nat. Rev. Mol. Cell Biol. 11, 545–555. doi: 10.1038/nrm2938

PubMed Abstract | CrossRef Full Text | Google Scholar

Akter, N., and Islam, M. R. (2017). Heat stress effects and management in wheat. A review. Agron. Sustain. Dev. 37, 1–17. doi: 10.1007/s13593-017-0443-9

CrossRef Full Text | Google Scholar

Alam, M. M., Nakamura, H., Ichikawa, H., Miyao, A., Hirochika, H., Kobayashi, K., et al. (2014). Response of an aspartic protease gene OsAP77 to fungal, bacterial and viral infections in rice Rice. 7:9. doi: 10.1186/s12284-014-0009-2

PubMed Abstract | CrossRef Full Text | Google Scholar

Anders, S., and Huber, W. (2010). Differential expression analysis for sequence count data. Genome Biol. 11:R106. doi: 10.1186/gb-2010-11-10-r106

PubMed Abstract | CrossRef Full Text | Google Scholar

Asakura, T., Watanabe, H., Abe, K., and Arai, S. (1995). Rice aspartic proteinase, oryzasin, expressed during seed ripening and germination, has a gene organization distinct from those of animal and microbial aspartic proteinases. Eur. J. Biochem. 232, 77–83. doi: 10.1111/j.1432-1033.1995.tb20783.x

PubMed Abstract | CrossRef Full Text | Google Scholar

Awlachew, Z. T., Singh, R., Kaur, S., Bains, N. S., and Chhuneja, P. (2016). Transfer and mapping of the heat tolerance component traits of Aegilops speltoides in tetraploid wheat Triticum durum. Mol. Breed. 36, 1–15. doi: 10.1007/s11032-016-0499-2

CrossRef Full Text | Google Scholar

Berkman, P. J., Lai, K., Lorenc, M. T., and Edwards, D. (2012). Next generation sequencing applications for wheat crop improvement. American J. Bot. 99, 365–371. doi: 10.3732/ajb.1100309

PubMed Abstract | CrossRef Full Text | Google Scholar

Borkovich, K. A., Farrelly, F. W., Finkelstein, D. B., Taulien, J., and Lindquist, S. (1989). hsp82 is an essential protein that is required in higher concentrations for growth of cells at higher temperatures. Mol. Cell. Biol. 9, 3919–3930.

Google Scholar

Buckley, B. A. (2006). The cellular response to heat stress in the goby Gillichthys mirabilis: a cDNA microarray and protein-level analysis. J. Exp. Biol. 209, 2660–2677. doi: 10.1242/jeb.02292

PubMed Abstract | CrossRef Full Text | Google Scholar

Bukovnik, U., Fu, J., Bennett, M., Prasad, P. V., and Ristic, Z. (2009). Heat tolerance and expression of protein synthesis elongation factors, EF-Tu and EF-1α, in spring wheat. Funct. Plant Biol. 36, 234–241. doi: 10.1071/FP08266

PubMed Abstract | CrossRef Full Text | Google Scholar

Cakmak, I. (2000). Role of zinc in protecting plant cells from reactive oxygen species. New Phytol. 146, 185–205. doi: 10.1046/j.1469-8137.2000.00630.x

PubMed Abstract | CrossRef Full Text | Google Scholar

Camacho, C., Coulouris, G., Avagyan, V., Ma, N., Papadopoulos, J., Bealer, K., et al. (2009). BLAST+: architecture and applications. BMC Bioinform. 10:421. doi: 10.1186/1471-2105-10-421

PubMed Abstract | CrossRef Full Text | Google Scholar

Carrillo-Perdomo, E., Vidal, A., Kreplak, J., Duborjal, H., Leveugle, M., Duarte, J., et al. (2020). Development of new genetic resources for faba bean (Vicia faba L.) breeding through the discovery of gene-based SNP markers and the construction of a high-density consensus map. Sci. Rep. 10:6790. doi: 10.1038/s41598-020-63664-7

PubMed Abstract | CrossRef Full Text | Google Scholar

Chen, X. Y., Cao, X. Y., Zhang, Y. J., Islam, S., Zhang, J. J., Yang, R. C., et al. (2016). Genetic characterization of cysteine-rich type-b avenin-like protein coding genes in common wheat. Sci. Rep. 6, 1–12. doi: 10.1038/srep30692

CrossRef Full Text | Google Scholar

Chu, C., Wang, S., Paetzold, L., Wang, Z., Hui, K., Rudd, J. C., et al. (2021). RNA-seq analysis reveals different drought tolerance mechanisms in two broadly adapted wheat cultivars 'TAM 111′ and 'TAM 112′. Sci. Rep. 11:4301. doi: 10.1038/s41598-021-83372-0

PubMed Abstract | CrossRef Full Text | Google Scholar

Cohen, S. P., and Leach, J. E. (2019). Abiotic and biotic stresses induce a core transcriptome response in rice. Sci. Rep. 9:6273. doi: 10.1038/s41598-019-42731-8

PubMed Abstract | CrossRef Full Text | Google Scholar

Conner, T. W., La Fayette, P. R., Nagao, R. T., and Key, J. L. (1990). Sequence and expression of a HSP83 from Arabidopsis thaliana. Plant Physiol. 94, 1689–1695. doi: 10.1104/pp.94.4.1689

PubMed Abstract | CrossRef Full Text | Google Scholar

Coolen, S., Proietti, S., Hickman, R., Davila Olivas, N. H., Huang, P. P., Van Verk, M. C., et al. (2016). Transcriptome dynamics of Arabidopsis during sequential biotic and abiotic stresses. Plant J. 86, 249–267. doi: 10.1111/tpj.13167

PubMed Abstract | CrossRef Full Text | Google Scholar

Dai, X., and Zhao, P. X. (2011). psRNATarget: a plant small RNA target analysis server. Nucl. acids res. 39, W155–W159. doi: 10.1093/nar/gkr319

PubMed Abstract | CrossRef Full Text | Google Scholar

De Carvalho, M. H. C., d’Arcy-Lameta, A., Roy-Macauley, H., Gareil, M., El Maarouf, H., Pham-Thi, A. T., et al. (2001). Aspartic protease in leaves of common bean (Phaseolus vulgaris L.) and cowpea (Vigna unguiculata L. Walp): enzymatic activity, gene expression and relation to drought susceptibility. FEBS Lett. 492, 242–246. doi: 10.1016/S0014-5793(01)02259-1

PubMed Abstract | CrossRef Full Text | Google Scholar

Disante, K. B., Fuentes, D., and Cortina, J. (2011). Response to drought of Zn-stressed Quercus suber L. seedlings. Env. Exp. Bot. 70, 96–103. doi: 10.1016/j.envexpbot.2010.08.008

CrossRef Full Text | Google Scholar

Dombrowski, J. E. (2003). Salt stress activation of wound-related genes in tomato plants. Plant Physiol. 132, 2098–2107. doi: 10.1104/pp.102.019927

PubMed Abstract | CrossRef Full Text | Google Scholar

Du Pont, F. M., Hurkman, W. J., Vensel, W. H., Chan, R., Lopez, R., Tanaka, C. K., et al. (2006). Differential accumulation of sulfur-rich and sulfur-poor wheat flour proteins is affected by temperature and mineral nutrition during grain development. J. Cereal Sci. 44, 101–112. doi: 10.1016/j.jcs.2006.04.003

CrossRef Full Text | Google Scholar

Dutta, S., Kumawat, G., Singh, B. P., Gupta, D. K., Singh, S., Dogra, V., et al. (2011). Development of genic-SSR markers by deep transcriptome sequencing in pigeon pea [Cajanus cajan (L.) Millspaugh]. BMC Plant Biol. 11:17. doi: 10.1186/1471-2229-11-17

PubMed Abstract | CrossRef Full Text | Google Scholar

Fahad, S., Bajwa, A. A., Nazir, U., Anjum, S. A., Farooq, A., Zohaib, A., et al. (2017). Crop production under drought and heat stress: plant responses and management options. Front. Plant Sci. 8:1147. doi: 10.3389/fpls.2017.01147

PubMed Abstract | CrossRef Full Text | Google Scholar

Farooq, M., Bramley, H., Palta, J. A., and Siddique, K. H. (2011). Heat stress in wheat during reproductive and grain-filling phases. Crit. Rev. Plant Sci. 30, 491–507. doi: 10.1080/07352689.2011.615687

CrossRef Full Text | Google Scholar

Felsheim, R. F., and Das, A. (1992). Structure and expression of a heat-shock protein 83 gene of Pharbitis nil. Plant Physiol. 100, 1764–1771. doi: 10.1104/pp.100.4.1764

PubMed Abstract | CrossRef Full Text | Google Scholar

Feng, G. H., Richardson, M., Chen, M. S., Kramer, K. J., Morgan, T. D., and Reeck, G. R. (1996). α-Amylase inhibitors from wheat: amino acid sequences and patterns of inhibition of insect and human α-amylases. Insect Biochem. Mol. Biol. 26, 419–426. doi: 10.1016/0965-1748(95)00087-9

PubMed Abstract | CrossRef Full Text | Google Scholar

Franco, O. L., Rigden, D. J., Melo, R., Bloch, C. Jr., Silva, C. P., and Grossi de Sá, M. F. (2000). Activity of wheat αamylase inhibitors towards bruchid αamylases and structural explanation of observed specificities. Eur. J. Biochem. 267, 2166–2173. doi: 10.1046/j.1432-1327.2000.01199.x

PubMed Abstract | CrossRef Full Text | Google Scholar

Fu, L., Niu, B., Zhu, Z., Wu, S., and Li, W. (2012). CD-HIT: accelerated for clustering the next-generation sequencing data. Bioinformatics 28, 3150–3152. doi: 10.1093/bioinformatics/bts565

PubMed Abstract | CrossRef Full Text | Google Scholar

Goswami, S., Kumar, R. R., and Rai, R. D. (2014). Heat-responsive microRNAs regulate the transcription factors and heat shock proteins in modulating thermo stability of starch biosynthesis enzymes in wheat (Triticum aestivum L.) under the heat stress. Aust. J. Crop. Sci. 8:697

Google Scholar

Grabherr, M. G., Haas, B. J., Yassour, M., Levin, J. Z., Thompson, D. A., Amit, I., et al. (2011). Full-length transcriptome assembly from RNA-Seq data without a reference genome. Nat. Biotechnol. 29, 644–652. doi: 10.1038/nbt.1883

PubMed Abstract | CrossRef Full Text | Google Scholar

Grover, A., Mittal, D., Negi, M., and Lavania, D. (2013). Generating high temperature tolerant transgenic plants: achievements and challenges. Plant Sci. 206, 38–47. doi: 10.1016/j.plantsci.2013.01.005

CrossRef Full Text | Google Scholar

Guevara, M. G., Almeida, C., Mendieta, J. R., Faro, C. J., Veríssimo, P., Pires, E. V., et al. (2005). Molecular cloning of a potato leaf cDNA encoding an aspartic protease (StAsp) and its expression after P. infestans infection. Plant Physiol. Biochem. 43, 882–889. doi: 10.1016/j.plaphy.2005.07.004

PubMed Abstract | CrossRef Full Text | Google Scholar

Gupta, P. K., Balyan, H. S., Gahlaut, V., and Kulwal, P. L. (2012). Phenotyping, genetic dissection, and breeding for drought and heat tolerance in common wheat: status and prospectus. Plant Breed. 6, 5–168. doi: 10.1002/9781118358566.ch2

CrossRef Full Text | Google Scholar

Haas, B. J., Papanicolaou, A., Yassour, M., Grabherr, M., Blood, P. D., Bowden, J., et al. (2013). De novo transcript sequence reconstruction from RNA-seq using the trinity platform for reference generation and analysis. Nat. Protoc. 8, 1494–1512. doi: 10.1038/nprot.2013.084

PubMed Abstract | CrossRef Full Text | Google Scholar

Hajheidari, M., Eivazi, A., Buchanan, B. B., Wong, J. H., Majidi, I., and Salekdeh, G. H. (2007). Proteomics uncovers a role for redox in drought tolerance in wheat. J. Proteome Res. 6, 1451–1460. doi: 10.1021/pr060570j

PubMed Abstract | CrossRef Full Text | Google Scholar

Henry, R. J., Rangan, P., and Furtado, A. (2016). Functional cereals for production in new and variable climates. Curr. Opin. Plant Biol. 30, 11–18. doi: 10.1016/j.pbi.2015.12.008

PubMed Abstract | CrossRef Full Text | Google Scholar

Iquebal, M. A., Sharma, P., Jasrotia, R. S., Jaiswal, S., Kaur, A., Saroha, M., et al. (2019). RNAseq analysis reveals drought-responsive molecular pathways with candidate genes and putative molecular markers in root tissue of wheat. Sci. Rep. 9, 13917–13918. doi: 10.1038/s41598-019-49915-2

PubMed Abstract | CrossRef Full Text | Google Scholar

Jiang, S. S., Liang, X. N., Li, X., Wang, S. L., Lv, D. W., Ma, C. Y., et al. (2012). Wheat drought-responsive grain proteome analysis by linear and nonlinear 2-DE and MALDI-TOmass spectrometry. Int. J. Mol. Sci. 13, 16065–16083. doi: 10.3390/ijms131216065

PubMed Abstract | CrossRef Full Text | Google Scholar

Jin, J., Tian, F., Yang, D. C., Meng, Y. Q., Kong, L., Luo, J., et al. (2016). PlantTFDB 4.0: toward a central hub for transcription factors and regulatory interactions in plants. Nucleic Acids Res. 45, D1040–D1045. doi: 10.1093/nar/gkw982

PubMed Abstract | CrossRef Full Text | Google Scholar

Kang, W. H., Sim, Y. M., Koo, N., Nam, J. Y., Lee, J., Kim, N., et al. (2020). Transcriptome profiling of abiotic responses to heat, cold, salt, and osmotic stress of Capsicum annuum L. Sci. Data. 7:17. doi: 10.1038/s41597-020-0352-7

PubMed Abstract | CrossRef Full Text | Google Scholar

Kang, W. H., and Yeom, S. I. (2018). Genome-wide identification, classification, and expression analysis of the receptor-like protein family in tomato. Plant Pathol. J. 34, 435–444. doi: 10.5423/PPJ.OA.02.2018.0032

PubMed Abstract | CrossRef Full Text | Google Scholar

Kobayashi, A., Sakamoto, A., Kubo, K., Rybka, Z., Kanno, Y., and Takatsuji, H. (1998). Seven zincfinger transcription factors are expressed sequentially during the development of anthers in petunia. Plant J. 13, 571–576. doi: 10.1046/j.1365-313X.1998.00043.x

PubMed Abstract | CrossRef Full Text | Google Scholar

Kosina, P., Reynolds, M., Dixon, J., and Joshi, A. (2007). Stakeholder perception of wheat production constraints, capacity building needs, and research partnerships in developing countries. Euphytica 157, 475–483. doi: 10.1007/s10681-007-9529-9

CrossRef Full Text | Google Scholar

Kumar, R. R., Goswami, S., Sharma, S. K., Kala, Y. K., Rai, G. K., Mishra, D. C., et al. (2015). Harnessing next generation sequencing in climate change: RNA-Seq analysis of heat stress-responsive genes in wheat (Triticum aestivum L.). J. Int. Biol. 19, 632–647. doi: 10.1089/omi.2015.0097

PubMed Abstract | CrossRef Full Text | Google Scholar

Kumar, R. R., Goswami, S., Singh, K., Dubey, K., Rai, G. K., Singh, B., et al. (2018). Characterization of novel heat-responsive transcription factor (TaHSFA6e) gene involved in regulation of heat shock proteins (HSPs)—A key member of heat stress-tolerance network of wheat. J. Biotechnol. 279, 1–12. doi: 10.1016/j.jbiotec.2018.05.008

PubMed Abstract | CrossRef Full Text | Google Scholar

Kumar, R. R., Sharma, S. K., and Goswami, S. (2013). Characterization of differentially expressed stress associated proteins in starch granule development under heat stress in wheat (Triticum aestivum L.). Ind. J. Biochem. Biophys. 50, 126–138.

Google Scholar

Laino, P., Shelton, D., Finnie, C., De Leonardis, A. M., Mastrangelo, A. M., Svensson, B., et al. (2010). Comparative proteome analysis of metabolic proteins from seeds of durum wheat (cv. Svevo) subjected to heat stress. Proteomics 10, 2359–2368. doi: 10.1002/pmic.200900803

PubMed Abstract | CrossRef Full Text | Google Scholar

Lal, M. K., Tiwari, R. K., Gahlaut, V., Mangal, V., and Kumar, A. (2021). Physiological and molecular insights on wheat responses to heat stress. Plant Cell Rep. 41, 1–18. doi: 10.1007/s00299-021-02784-4

CrossRef Full Text | Google Scholar

Langmead, B., Trapnell, C., Pop, M., and Salzberg, S. L. (2009). Ultrafast and memory-efficient alignment of short DNA sequences to the human genome. Genome Biol. 10:R25. doi: 10.1186/gb-2009-10-3-r25

PubMed Abstract | CrossRef Full Text | Google Scholar

Li, B., and Dewey, C. N. (2011). RSEM: accurate transcript quantification from RNA-Seq data with or without a reference genome. BMC Bioinform. 12, 1–16. doi: 10.1186/1471-2105-12-323

PubMed Abstract | CrossRef Full Text | Google Scholar

Li, H., and Durbin, R. (2009). Fast and accurate short read alignment with burrows-wheeler transform. Bioinformatics 25, 1754–1760. doi: 10.1093/bioinformatics/btp324

PubMed Abstract | CrossRef Full Text | Google Scholar

Li, H., Handsaker, B., Wysoker, A., Fennell, T., Ruan, J., Homer, N., et al. (2009). The sequence alignment/map format and SAM tools. Bioinformatics 25, 2078–2079. doi: 10.1093/bioinformatics/btp352

PubMed Abstract | CrossRef Full Text | Google Scholar

Li, M., Sui, N., Lin, L., Yang, Z., and Zhang, Y. (2019). Transcriptomic profiling revealed genes involved in response to cold stress in maize. Funct. Plant Biol. 46, 830–844. doi: 10.1071/FP19065

PubMed Abstract | CrossRef Full Text | Google Scholar

Livak, K. J., and Schmittgen, T. D. (2001). Analysis of relative gene expression data using real-time quantitative PCR and the 2^−ΔΔCT method. Methods 25, 402–408. doi: 10.1006/meth.2001.1262

CrossRef Full Text | Google Scholar

Ma, F., Li, M., Li, T., Liu, W., Liu, Y., Li, Y., et al. (2013). Overexpression of avenin-like b proteins in bread wheat (Triticum aestivum L.) improves dough mixing properties by their incorporation into glutenin polymers. PLoS One 8:e66758. doi: 10.1371/journal.pone.0066758

PubMed Abstract | CrossRef Full Text | Google Scholar

Marschner, H. (1995). Mineral Nutrition of Higher Plants. 2nd Edn. London: Academic Press.

Google Scholar

Martin, M. (2011). Cutadapt removes adapter sequences from high-throughput sequencing reads. EMBnet. J. 17, 10–12. doi: 10.14806/ej.17.1.200

CrossRef Full Text | Google Scholar

Mendieta, J. R., Pagano, M. R., Munoz, F. F., Daleo, G. R., and Guevara, M. G. (2006). Antimicrobial activity of potato aspartic proteases (StAPs) involves membrane permeabilization. Microbiol. 152, 2039–2047. doi: 10.1099/mic.0.28816-0

PubMed Abstract | CrossRef Full Text | Google Scholar

Morozova, O., and Marra, M. A. (2008). Applications of next-generation sequencing technologies in functional genomics. Genomics 92, 255–264. doi: 10.1016/j.ygeno.2008.07.001

CrossRef Full Text | Google Scholar

Muñoz, F. F., Mendieta, J. R., Pagano, M. R., Paggi, R. A., Daleo, G. R., and Guevara, M. G. (2010). The swaposin-like domain of potato aspartic protease (StAsp-PSI) exerts antimicrobial activity on plant and human pathogens. Peptides 31, 777–785. doi: 10.1016/j.fitote.2009.10.004

PubMed Abstract | CrossRef Full Text | Google Scholar

Mwadzingeni, L., Shimelis, H., Dube, E., Laing, M. D., and Tsilo, T. J. (2016). Breeding wheat for drought tolerance: progress and technologies. J. Integr. Agric. 15, 935–943. doi: 10.1016/S2095-3119(15)61102-9

CrossRef Full Text | Google Scholar

Pareek, A., Dhankher, O. P., and Foyer, C. H. (2020). Mitigating the impact of climate change on plant productivity and ecosystem sustainability. J. Exp. Bot. 71, 451–456. doi: 10.1093/jxb/erz518

PubMed Abstract | CrossRef Full Text | Google Scholar

Peck, A. W., and McDonald, G. K. (2010). Adequate zinc nutrition alleviates the adverse effects of heat stress in bread wheat. Plant Soil 337, 355–374. doi: 10.1007/s11104-010-0532-x

CrossRef Full Text | Google Scholar

Petrucci, T., Rab, A., Tomasi, M., and Silano, V. (1976). Further characterization studies of the α-amylase protein inhibitor of gel electrophoretic mobility 0.19 from the wheat kernel. Biochim. Biophys. Acta. Proteins Struc. 420, 288–297. doi: 10.1016/0005-2795(76)90320-2

PubMed Abstract | CrossRef Full Text | Google Scholar

Quast, C., Pruesse, E., Yilmaz, P., Gerken, J., Schweer, T., Yarza, P., et al. (2012). The SILVA ribosomal RNA gene database project: improved data processing and web-based tools. Nucleic Acids Res. 41, D590–D596. doi: 10.1093/nar/gks1219

PubMed Abstract | CrossRef Full Text | Google Scholar

Rane, A. S., Joshi, R. S., and Giri, A. P. (2020). Molecular determinant for specificity: differential interaction of α-amylases with their proteinaceous inhibitors. BBA Gen. Sub. 1864:129703. doi: 10.1016/j.bbagen.2020.129703

PubMed Abstract | CrossRef Full Text | Google Scholar

Rangan, P., Furtado, A., and Henry, R. (2020). Transcriptome profiling of wheat genotypes under heat stress during grain-filling. J. Cereal Sci. 91:102895. doi: 10.1016/j.jcs.2019.102895

CrossRef Full Text | Google Scholar

Rawlings, N. D., and Salvesen, G. (2013). Handbook of Proteolytic Enzymes. Vol. 1>. United States: Academic Press.

Google Scholar

Roy, S., Mishra, M., Dhankher, O. P., Singla-Pareek, S. L., and Pareek, A. (2019). “Molecular chaperones: key players of abiotic stress response in plants” in Genetic Enhancement of Crops for Tolerance to Abiotic Stress: Mechanisms and Approaches. Vol. 1. eds. N. Singh, S. Mittal, and N. Thirunavukkarasu (Cham: Springer), 125–165.

Google Scholar

Runeberg-Roos, P., Kervinen, J., Kovaleva, V., Raikhel, N. V., and Gal, S. (1994). The aspartic proteinase of barley is a vacuolar enzyme that processes probarleylectin in vitro. Plant Physiol. 105, 321–329. doi: 10.1104/pp.105.1.321

PubMed Abstract | CrossRef Full Text | Google Scholar

Sarwar, M., Saleem, M. F., Ullah, N., Ali, S., Rizwan, M., Shahid, M. R., et al. (2019). Role of mineral nutrition in alleviation of heat stress in cotton plants grown in glasshouse and field conditions. Sci. Rep. 9, 1–7. doi: 10.1038/s41598-019-49404-6

CrossRef Full Text | Google Scholar

Sebastián, D., Fernando, F. D., Raúl, D. G., and Gabriela, G. M. (2020). Overexpression of Arabidopsis aspartic protease APA1 gene confers drought tolerance. Plant Sci. 292:110406. doi: 10.1016/j.plantsci.2020.110406

PubMed Abstract | CrossRef Full Text | Google Scholar

Shamloo, M., Babawale, E. A., Furtado, A., Henry, R. J., Eck, P. K., and Jones, P. J. (2017). Effects of genotype and temperature on accumulation of plant secondary metabolites in Canadian and Australian wheat grown under controlled environments. Sci. Rep. 7, 1–13. doi: 10.1038/s41598-017-09681-5

CrossRef Full Text | Google Scholar

Shamovsky, I., Ivannikov, M., Kandel, E. S., Gershon, D., and Nudler, E. (2006). RNA-mediated response to heat shock in mammalian cells. Nature 440, 556–560. doi: 10.1038/nature04518

PubMed Abstract | CrossRef Full Text | Google Scholar

Sharma, D., Chandra Pandey, G., Mamrutha, H. M., Singh, R., Singh, N. K., Singh, G. P., et al. (2019). Genotype–phenotype relationships for high‐temperature tolerance: An integrated method for minimizing Phenotyping constraints in wheat. Crop Sci. 59, 1973–1982. doi: 10.2135/cropsci2019.01.0055

CrossRef Full Text | Google Scholar

Sheoran, S., Jaiswal, S., Kumar, D., Raghav, N., Sharma, R., Pawar, S., et al. (2019). Uncovering genomic regions associated with 36 agro-morphological traits in Indian spring wheat using GWAS. Front. Plant Sci. 10:527. doi: 10.3389/fpls.2019.00527

PubMed Abstract | CrossRef Full Text | Google Scholar

Siebert, S., and Ewert, F. (2014). Future crop production threatened by extreme heat. Environ. Res. Lett. 9:041001. doi: 10.1088/1748-9326/9/4/041001

CrossRef Full Text | Google Scholar

Simões, R., Faro, D., and Bur, F. C. (2007). Characterization of recombinant CDR1, an Arabidopsis aspartic proteinase involved in disease resistance. J. Biol. Chem. 282, 31358–31365. doi: 10.1074/jbc.M702477200

PubMed Abstract | CrossRef Full Text | Google Scholar

Singh, R., Iquebal, M. A., Mishra, C. N., Jaiswal, S., Kumar, D., Raghav, N., et al. (2019). Development of model web-server for crop variety identification using throughput SNP genotyping data. Sci. Rep. 9, 5122–5129. doi: 10.1038/s41598-019-41204-2

PubMed Abstract | CrossRef Full Text | Google Scholar

Tavallali, V., Rahemi, M., Eshghi, S., Kholdebarin, B., and Ramezanian, A. (2010). Zinc alleviates salt stress and increases antioxidant enzyme activity in the leaves of pistachio (Pistaciavera L.‘Badami’) seedlings. Turk. J. Agr. Forest. 34, 349–359. doi: 10.3906/tar-0905-10

CrossRef Full Text | Google Scholar

Thiel, T., Michalek, W., Varshney, R., and Graner, A. (2003). Exploiting EST databases for the development and characterization of gene-derived SSR-markers in barley (Hordeum vulgare L.). Theor. Appl. Genet. 106, 411–422. doi: 10.1007/s00122-002-1031-0

PubMed Abstract | CrossRef Full Text | Google Scholar

Unamba, C. I. N., Nag, A., and Sharma, R. K. (2015). Next generation sequencing technologies: The doorway to the unexplored genomics of non-model plants. Front. Plant Sci. 6:1074. doi: 10.3389/fpls.2015.01074

PubMed Abstract | CrossRef Full Text | Google Scholar

Untergasser, I., Cutcutache, T., Koressaar, J., Ye, B. C., Faircloth, R. M., and Rozen, S. G. (2012). Primer3enew capabilities and interfaces. Nucleic Acids Res. 40:e115. doi: 10.1093/nar/gks596

PubMed Abstract | CrossRef Full Text | Google Scholar

Unterseer, S., Bauer, E., Haberer, G., Seidel, M., Knaak, C., Ouzunova, M., et al. (2014). A powerful tool for genome analysis in maize: development and evaluation of the high density 600 k SNP genotyping array. BMC Genomics 15:823. doi: 10.1186/1471-2164-15-823

PubMed Abstract | CrossRef Full Text | Google Scholar

Wang, C., Liu, R., Liu, Y., Hou, W., Wang, X., Miao, Y., et al. (2021). Development and application of the Faba_bean_130K targeted next-generation sequencing SNP genotyping platform based on transcriptome sequencing. Theor. Appl. Genet. 134, 3195–3207. doi: 10.1007/s00122-021-03885-0

PubMed Abstract | CrossRef Full Text | Google Scholar

Xin, M., Wang, Y., Yao, Y., Xie, C., Peng, H., Ni, Z., et al. (2010). Diverse set of microRNAs are responsive to powdery mildew infection and heat stress in wheat (Triticum aestivum L.). BMC Plant Biol. 10:123. doi: 10.1186/1471-2229-10-123

PubMed Abstract | CrossRef Full Text | Google Scholar

Yang, F., Jørgensen, A. D., Li, H., Søndergaard, I., Finnie, C., Svensson, B., et al. (2011). Implications of high-temperature events and water deficits on protein profiles in wheat (Triticum aestivum L. cv. Vinjett) grain. Proteomics 11, 1684–1695. doi: 10.1002/pmic.201000654

PubMed Abstract | CrossRef Full Text | Google Scholar

Yao, X., Xiong, W., Ye, T., and Wu, Y. (2012). Overexpression of the aspartic protease ASPG1 gene confers drought avoidance in Arabidopsis. J. Exp. Bot. 63, 2579–2593. doi: 10.1093/jxb/err433

PubMed Abstract | CrossRef Full Text | Google Scholar

Zhang, Y., Hu, X., Islam, S., She, M., Peng, Y., Yu, Z., et al. (2018). New insights into the evolution of wheat avenin-like proteins in wild emmer wheat (Triticum dicoccoides). Proc. Natl. Acad. Sci. 115, 13312–13317. doi: 10.1073/pnas.1812855115

PubMed Abstract | CrossRef Full Text | Google Scholar

Zhang, Y., Hu, X., Juhasz, A., Islam, S., Yu, Z., Zhao, Y., et al. (2020). Characterising avenin-like proteins (ALPs) from albumin/globulin fraction of wheat grains by RP-HPLC, SDS-PAGE, and MS/MS peptides sequencing. BMC Plant Biol. 20, 19–45. doi: 10.1186/s12870-020-2259-z

PubMed Abstract | CrossRef Full Text | Google Scholar

Zhang, X. M., Zhao, X. Q., Feng, C. X., Liu, N., Feng, H., Wang, X. J., et al. (2014). The cloning and characterization of a DEAD-box RNA helicase from stress-responsive wheat. Physiol. Mol. Plant Pathol. 88, 36–42. doi: 10.1016/j.pmpp.2014.07.004

CrossRef Full Text | Google Scholar

Zhao, J., Lu, Z., Wang, L., and Jin, B. (2021). Plant responses to heat stress: physiology, transcription, noncoding RNAs, and epigenetics. Int. J. Mol. Sci. 22:117. doi: 10.3390/ijms22010117

PubMed Abstract | CrossRef Full Text | Google Scholar

Zheng, J., Liu, H., Wang, Y., Wang, L., Chang, X., Jing, R., et al. (2014). TEF-7A, a transcript elongation factor gene, influences yield-related traits in bread wheat (Triticum aestivum L.). J. Exp. Bot. 65, 5351–5365. doi: 10.1093/jxb/eru306

PubMed Abstract | CrossRef Full Text | Google Scholar