RNA-Seq Analysis of Gene Expression Changes Related to Delay of Flowering Time under Drought Stress in Tropical Maize

Few studies have reported on the flowering time mechanism of tropical maize under short-day conditions. Drought, another important factor that affects flowering time, has been reported to delay the silking date in tropical maize. However, due to the lack of genetic information related to flowering in maize, the mechanism by which drought delays flowering is unclear. To further understand this process, we analyzed drought-responsive genes using RNA sequencing and identified genes related to flowering time, including contigs from de novo assembly. The results revealed changes in the expression of flowering-time genes, including INDETERMINATE1 (ID1), Heading date 3a (Hd3a), CONSTANS-like genes, and ZEA MAYS CENTRORADIALES8 (ZCN8), which are known to be crucial factors in flowering. In particular, Hd3a, CONZ1, and ZCN8, which have been reported to accelerate flowering under short-day conditions, were downregulated by drought stress. Changes in gene expression appear to play an important role in changes in flowering time under drought. These expression profiles will help to further understand the flowering-time genes of tropical maize and the delayed flowering time resulting from drought.


Introduction
Maize (Zea mays L.) is one of the most important crops worldwide, with extensive applications in food, feed, and biofuel [1]. Unlike crops such as rice or wheat, maize is an allogamous plant, which involves a more complex genome and a diverse combination of traits that help the species to adapt to various environments. Flowering time is a very important factor for allogamous plants. The pollination of maize occurs when pollen shed from the tassel during anthesis is captured by stigmas (silk) of the ear. Therefore, for maize, pollen shedding and silk emergence must occur at approximately the same time in order to increase the probability of pollination and harvesting [2,3]. Furthermore, drought is an important factor that influences flowering, and drought stress has been reported to increase anthesis-silking interval (ASI) and negatively affect the fertilization rate, number of kernels, seed quality, and weight [4,5]. Due to the relationship between flowering synchrony and drought, the ASI is considered to be a secondary trait associated with stress tolerance and the harvest index under drought [5][6][7][8].
With the release of the maize reference genome [9], the design of nested association mapping (NAM) [10], and the development of technologies such as genome-wide association studies (GWAS), RNA sequencing, and bioinformatics tools, information regarding the diversity of the genome and stress response system has become accessible. Consequently, extensive research efforts have been devoted to understanding how maize responds to drought stress. The reported molecular responses of plants to drought stress involve hormones [11,12], gene networks [13][14][15], and transcription factors [16][17][18][19]. The results of previous studies have also suggested that the drought response system of maize differs between tissues and development stages [19][20][21][22][23]. In addition, several papers on the drought stress response during the flowering stage have identified genes using GWAS [24][25][26][27], gene expression profiles [19,20,[28][29][30], and several transcription factors [31][32][33][34]. Many studies have revealed genes related to flowering time in rice and Arabidopsis, but relatively few flowering-time genes of maize are known. Several studies have attempted to identify the role of maize flowering-time genes using Arabidopsis. As a result, differences in flowering time pathways were identified between dicot and monocot plants and between shortday (SD) and long-day (LD) conditions. For example, in rice under SD conditions, the increased expression of Heading date 3a (Hd3a) induced Early heading date 1 (Ehd1), which promoted flowering time, but this mechanism was not activated in dicot plants [35,36]. Similarly, GIGANTEA (GI) delayed flowering time in rice, but it exerted the opposite effect in Arabidopsis [37,38]. In addition, Heading date 1 (Hd1) promoted flowering in SD but repressed it in LD [39][40][41]. Since the discovery of the maize gene INDETERMINATE 1 (ID1) [42], various genes related to flowering time have been identified. ID1 is a zinc-finger transcription factor and is known to promote flowering time [43]. Maize contains a large number of ZEA CENTRORADIALIS (ZCN) family genes that have similar functions to those of AtFT and AtTFL [44]. ZmZCN8 is highly similar to AtFT1 and interacts with DELAYED FLOWERING1 (ZmDLF1), an FD-like bZIP protein. The expression of ZmZCN8 was increased in the reproductive stage of the tropical line under SD, but it was unchanged in temperate maize [31,45]. ZmCONZ1, an ortholog of AtCO and OsHd1, is also known to regulate flowering by affecting the expression of ZCN8 [31,46]. However, there is no evidence that explains how ZmCONZ1 activates ZmZCN8 expression in SD conditions, and thus, the mechanism remains unclear. The role of Hd3a has been reported to induce flowering time only under SD conditions. Under SD conditions, the expression of Hd3a induces Ehd1, which promotes flowering [36]. However, most studies are based on B73, which is a temperate inbred line grown in long-day conditions. There is relatively little information on the tropical maize line, although the genes involved in flowering time are known to differ depending on the day length [29,34,47,48] and maize genotype [18,49,50]. Furthermore, the complete flowering pathway under SD conditions is unknown.
In this study, we used RNA-seq analysis to investigate the transcriptome changes involved in the flowering time of a tropical maize line under drought conditions. The results of RNA-seq were analyzed using bioinformatics tools to narrow down the potential candidate genes related to changes in flowering time. The results of this study provide new information on the genetic characteristics of tropical maize and can improve the current understanding of the relationship between drought and flowering time.

Plant Growth and Drought Conditions
Maize cv.(cultivar) Ki11 (Zea mays), which was used as the material in this study, is a drought-tolerant inbred line developed in 1982 at Kasetsart University and is one of the NAM parental lines. Ki11 is a tropical maize line that is less affected by drought than its temperate counterparts [51][52][53]. The plants were grown in pots (37 cm × 37 cm × 22 cm) in a greenhouse (26-28 • C) under short-day (12 hours of light and 12 hours of dark) conditions until the top of the tassel was visible. When the tassel was visible, half of the plants were exposed to drought stress conditions. The well-watered pots were maintained at over −0.2 MPa (15~20% soil water content) with irrigation, while drought-treated pots were maintained at less than −1.5 MPa (5~8% soil water content). Days to anthesis (DA), days to silking (DS), and the anthesis-silking interval (ASI) were measured according to Buckler's method [24]. When pollen shedding began, the top collared leaves were harvested from plants in each replicate group and each condition (well-watered/control or drought stress). All harvested leaves were immediately frozen using liquid nitrogen and were transferred to a deep freezer (−80 • C) for storage.

Total RNA Isolation and RNA Sequence
Total RNA was isolated from harvested leaves of cv. Ki11 in each condition with biological replication. The RNA of leaves was extracted using the RNeasy Plant Mini Kit (Qiagen), and the quality of RNA was assessed using a 2100 Bioanalyzer RNA Nanochip (Agilent Technologies). The RNA-seq analysis was performed using the Illumina HiSeq platform (Marcrogen) to obtain the expression profiles of the plants.

Mapping and Differential Gene Expression Analysis
The raw reads were filtered with FASTX-Toolkit and then aligned to the maize reference genome (B73_RefGen_v3 annotation build (5b+)) using default parameters in Tophat2 [54]. The aligned sequences were assembled using Cufflinks and compared with the annotated reference genome using Cuffmerge [55]. Unmapped reads were imported to CLC Genomics Workbench (v7.5.1 CLC Bio). For the de novo assembly, we employed the following parameters: the mismatch count and maximum gap were set to 2, the insertion and deletion costs were set to 3, the minimum contig length was set to 200 bp, and the length fraction and similarity parameters were set to 0.5 and 0.9, respectively. The sequence reads were submitted to NCBI Short Read Archive (SRA) under accession numbers SRP082534 (B73) and SRP135093 (Ki11).
Fragments per kilobase of transcript per million mapped reads (FPKM) values were calculated for each gene. Differentially expressed genes (DEGs) under drought stress were identified by Cuffdiff. DEGs were defined by a false discovery rate (FDR) [56] less than 0.05 and an absolute log2 fold change value greater than 1. In addition, differential usage of exons and alternative splicing events were detected using DEXSeq [57] with an FDR under 0.05 and threshold absolute fold change over 2.

Gene Ontology Analysis and Functional Annotation
The selected DEGs included contigs from de novo assembly results, and genes with differentially expressed isoforms (DEIs) were mapped to Gene Ontology (GO) classifications using Blast2GO [58]. For enrichment analysis, the genes were assessed using the BiNGO plugin for Cytoscape [59] with the Bonferroni correction method (p-value correction with a cut-off of 0.05). To identify transcription factor domains in proteins, all translated sequences were matched against the Pfam database (PlantTFDB v. 3.0).

Quantitative Real-Time Polymerase Chain Reaction (qRT-PCR) Validation
The RNA isolated from each group (well-watered/control or drought stress) was used to construct a cDNA library. The RNA was reverse-transcribed using the PrimeScript First Strand cDNA Synthesis kit (Takara), and qRT-PCR was performed using the CFX Connect™ Real-Time PCR Detection System (Bio-Rad). Gene expression in each condition was calculated using the CT method. The expression of genes was normalized using CT values for the reference genes (18 s rRNA and Actin1) to determine relative fold changes. All genes and contigs were validated using three independent biological replicates in each sample. To determine relative fold changes, gene expression was normalized to CT values for the reference genes. The CT values for each gene under each condition were calculated using the ∆∆ CT method, as previously described [60].

Changes in Anthesis-Silking Interval (ASI) Due to Drought Stress in Ki11
ASI is a key trait that is significantly affected by drought stress. Under normal conditions, days to silking (DS) and days to anthesis (DA) have intervals of 3 days. The ASI is increased by a delay in silking in tropical maize under drought conditions [61]. We found that the ASI was also increased by a delay in silking in Ki11 (Table 1). Under drought, anthesis was delayed by about one day, while silking was delayed by about five days. The delayed silking dates under drought increased the ASI, which can result in reproductive failure.

Analysis of Drought-Responsive Transcriptome
RNA samples were subjected to the RNA-seq pipeline. The number of total reads was much higher in the control, but the rate of concordant pairs was similar between the two conditions. Of the total reads, 71.5% could be mapped to the reference genome (Zea mays AGPv3.31) as unique and multiple matches (Table S1). A total of 38,263 genes and 9900 contigs were expressed under both conditions, and 1002 genes and 771 contigs passed cut-off values and were considered to be DEGs under drought conditions. In addition, principal component analysis (PCA) was performed to compare the responses of B73 and Ki11 to drought stress ( Figure S1). Independent replicates were clustered together for each cultivar. Then, to identify genes that changed at the transcriptional level under drought but not at the expression level, RNA-seq data were analyzed with DEXSeq. We found 234 differently expressed isoforms (DEIs) of 175 genes under drought conditions. Among DEGs, 454 genes were upregulated, while 548 genes were downregulated. A total of 118 genes were not associated with annotations but were sequence matched to the reference genome. These genes are represented by cufflinks ID. In de novo assembly results, upregulated contigs outnumbered downregulated contigs. A total of 465 contigs were upregulated under drought conditions, and 306 contigs were downregulated. Among genes with DEIs, 12 were DEGs, while 163 only changed at the transcriptional level under drought stress ( Figure 1, Table S2). The statistical significance of DEGs is illustrated in a volcano plot in Figure S2, in which DEGs with markedly up-and downregulated expression levels in response to drought stress are shown in red and blue, respectively. increased by a delay in silking in tropical maize under drought conditions [61]. We found that the ASI was also increased by a delay in silking in Ki11 (Table 1).

Analysis of Drought-Responsive Transcriptome
RNA samples were subjected to the RNA-seq pipeline. The number of total reads was much higher in the control, but the rate of concordant pairs was similar between the two conditions. Of the total reads, 71.5% could be mapped to the reference genome (Zea mays AGPv3.31) as unique and multiple matches (Table S1). A total of 38,263 genes and 9900 contigs were expressed under both conditions, and 1002 genes and 771 contigs passed cut-off values and were considered to be DEGs under drought conditions. In addition, principal component analysis (PCA) was performed to compare the responses of B73 and Ki11 to drought stress ( Figure S1). Independent replicates were clustered together for each cultivar. Then, to identify genes that changed at the transcriptional level under drought but not at the expression level, RNA-seq data were analyzed with DEX-Seq. We found 234 differently expressed isoforms (DEIs) of 175 genes under drought conditions. Among DEGs, 454 genes were upregulated, while 548 genes were downregulated. A total of 118 genes were not associated with annotations but were sequence matched to the reference genome. These genes are represented by cufflinks ID. In de novo assembly results, upregulated contigs outnumbered downregulated contigs. A total of 465 contigs were upregulated under drought conditions, and 306 contigs were downregulated. Among genes with DEIs, 12 were DEGs, while 163 only changed at the transcriptional level under drought stress ( Figure 1, Table S2). The statistical significance of DEGs is illustrated in a volcano plot in Figure S2, in which DEGs with markedly up-and downregulated expression levels in response to drought stress are shown in red and blue, respectively.

Gene Ontology Enrichment Analysis of Differentially Expressed Genes and Transcripts
Consensus sequences of drought-responsive genes and transcripts were subjected to Blast2Go to identify gene functions. GO enrichment analysis was performed on the Blast results to investigate biological functions. The results of GO Slim enrichment analyses were divided into two categories: drought-related genes and flowering-time-related genes. The results reveal a strong enrichment of functions involved in the response to stress and response to abiotic stress (GO: 006950 and GO: 009628) ( Figure 2 and Table S3).
Appl. Sci. 2021, 11, x FOR PEER REVIEW 5 of 13 Figure 1. Venn diagram displaying drought-responsive genes. Summary of differentially expressed genes (DEGs), contigs from de novo assembly results, and differentially expressed isoforms (DEIs) from alternative splicing events under drought conditions.

Gene Ontology Enrichment Analysis of Differentially Expressed Genes and Transcripts
Consensus sequences of drought-responsive genes and transcripts were subjected to Blast2Go to identify gene functions. GO enrichment analysis was performed on the Blast results to investigate biological functions. The results of GO Slim enrichment analyses were divided into two categories: drought-related genes and flowering-time-related genes. The results reveal a strong enrichment of functions involved in the response to stress and response to abiotic stress (GO: 006950 and GO: 009628) (Figure 2 and Table S3). In addition, functions such as response to external stimulus (GO: 009605) that might contribute to drought tolerance were observed in the biological process category. GO terms for the response to stress were also observed in the molecular function category.
The gene expression pattern of maize under drought stress has been reported to differ according to the tissue and development stage [19,22,23,62]. We identified genes in Ki11 that are also involved in the drought response in B73. For example, the drought-responsive genes GRMZM2G031033, GRMZM2G366532, GRMZM2G031033, GRMZM2G059799, GRMZM2G126812, GRMZM2G122816, AC208204.3_FG006, and GRMZM2G154628, which have been reported as hub genes in B73 under drought stress, were also observed in Ki11 [19]. Moreover, in contigs, we detected orthologous genes that have been reported as drought-responsive genes in B73 (GRMZM5G833699, GRMZM2G059799, GRMZM2G079440, and GRMZM2G059799). In addition, several contigs (orthologs of GRMZM2G063340, GRMZM2G133802, and GRMZM2G132093) and genes (GRMZM5G858784, GRMZM5G884600, and GRMZM2G159724) reported to be In addition, functions such as response to external stimulus (GO: 009605) that might contribute to drought tolerance were observed in the biological process category. GO terms for the response to stress were also observed in the molecular function category.
The gene expression pattern of maize under drought stress has been reported to differ according to the tissue and development stage [19,22,23,62]. We identified genes in Ki11 that are also involved in the drought response in B73. For example, the drought-responsive genes GRMZM2G031033, GRMZM2G366532, GRMZM2G031033, GRMZM2G059799, GR-MZM2G126812, GRMZM2G122816, AC208204.3_FG006, and GRMZM2G154628, which have been reported as hub genes in B73 under drought stress, were also observed in Ki11 [19]. Moreover, in contigs, we detected orthologous genes that have been reported as drought-responsive genes in B73 (GRMZM5G833699, GRMZM2G059799, GRMZM2G079440, and GRMZM2G059799). In addition, several contigs (orthologs of GRMZM2G063340, GR-MZM2G133802, and GRMZM2G132093) and genes (GRMZM5G858784, GRMZM5G884600, and GRMZM2G159724) reported to be directly or indirectly involved in drought stress from genome-wide association results were identified as drought-responsive genes. Blast analysis characterized these genes as involved in the response to stress (GO: 0006950) and response to stimulus (GO: 0009628). Furthermore, the expression profiles of 66 genes involved in flowering time were annotated with GO terms for reproduction (GO: 0000003). The flowering-time-related genes include ZCN, FLOWERING LOCUS T, MADS-box, and CONSTANS. We also found genes that were not associated with GO terms for reproduction but could affect flowering time. Of the drought-responsive genes, four genes and five contigs were also DEGs with potential effects on flowering time. These genes were obtained by comparing candidate genes involved in flowering time in maize with the genes of other plants, including cloned and homologous flowering-time genes [27] (Table S2). Three genes (GRMZM2G018627, GRMZM2G126018, and GRMZM2G348863) involved in the response to drought stress may have roles analogous to those of chlorophyll a-b binding protein (CAB), SQUAMOSA PROMOTER BINDING LIKE PROTEIN (SPL) 17, and flowering-promoting factor 1 (FPF1), respectively. One gene (GRMZM2G546097, homolog of GRMZM2G126018) and five contigs (contig854, homolog of GRMZM2G018627; contig4549, homolog of GRMZM2G431157; contig4777 and contig312, homolog of GRMZM2G142718; contig6971, homolog of GR-MZM2G373928; and contig7875, homolog of GRMZM2G013478) that were identified as responsive to drought are expected to play roles analogous to those of CAB, DNA-binding one zinc finger (Dof ), and Heading date 3a (Hd3a). In addition, we detected DEGs that were previously reported as drought-responsive genes in B73 under short-day and drought conditions [30]. By comparing the two datasets, we identified genes (GRMZM2G094990, GRMZM2G13501, GRMZM2G012717, and GRMZM2G120619) and contigs (two homologs of GRMZM2G142718) that are responsive to drought and potentially change flowering time. These genes potentially have roles analogous to those of Dof, CAB, and CONSTANS (CO)-like 16. Several genes that are known to regulate the flowering time under SD (shortday) conditions were identified as DEGs. INDETERMINATE 1 (ID1) plays a role in the activation of Early heading date 1 (Ehd1) and is an important factor in the initiation of flowering [19,63]. In our study, GRMZM2G011357, corresponding to the ID1 gene, was upregulated by drought. The MIKC-type MADS genes, including ZmMADS1 (ZmMADS56 and GRMZM2G070034), are known to be closely related to identified flowering time activators [34]. ZmMADS56 is a homolog of the rice gene OsMADS50, which is reported to be a flowering activator [64,65]. Genes corresponding to the MIKC_MADS family (GR-MZM2G026223 and GRMZM2G070034) were not among DEGs but were found to undergo alternative splicing in response to drought. In rice, Ehd1 promotes Hd3a and accelerates flowering under SD. In our study, Ehd1 (GRMZM2G096171) and Hd3a (contig6971 and GRMZM2G103666) were downregulated under drought conditions. Similarly, Zea mays1 (CONZ1, GRMZM2g026024) [46] is also known to affect ZCN8 (Zea CENTRORADIALIS, collinear orthologs of Hd3a) [66] in maize, and its expression decreases under drought.

Transcription Factors That Mediate Drought Tolerance and Flowering Time
Several genes (342 contigs and 50 sequences without annotations) could not be mapped to GO terms, even after Blast analysis. All DEGs were further analyzed using a plant transcription factor database (PlantTFDB v. 4.0) to predict their function. Transcription factor motifs were identified in a total of 689 genes, 45 of which were identified as transcription factors. Of these transcription factors, bHLH, C2H2, C3H, NAC, MYB, WRKY, ERF, and bZIP are the most frequently reported to have roles in development and the drought stress Several genes (342 contigs and 50 sequences without annotations) could not be mapped to GO terms, even after Blast analysis. All DEGs were further analyzed using a plant transcription factor database (PlantTFDB v. 4.0) to predict their function. Transcription factor motifs were identified in a total of 689 genes, 45 of which were identified as transcription factors. Of these transcription factors, bHLH, C2H2, C3H, NAC, MYB, WRKY, ERF, and bZIP are the most frequently reported to have roles in development and the drought stress response of various plants. These TFs are known to be affected by drought and may also affect flowering time (Figure 3) [15,17,32,[68][69][70][71]. Our results included many transcription factors, such as GATA, Dof, CO, and MADS families, that are known to be associated with flowering time (Table S5). Recently, a study on the expression pattern of each transcription factor in B73 was carried out, and a total of five hub genes (GRMZM2G113779, SPL; GRMZM2G387528, bHLH; GRMZM2G017349, bHLH; GRMZM2G011357, C2H2; GRMZM2G025642, NAC) were identified. Among the five hub genes, three genes (GRMZM2G113779, GRMZM2G387528 and GRMZM2G011357) were expressed in the leaves of B73. Our results confirmed the expression of three hub genes (GRMZM2G161673, homolog of GRMZM2G113779; contig1816, homolog of GRMZM2G387528; GRMZM2G011357) under drought conditions. GRMZM2G113779 is a hub gene that encodes SPL. The SPL family has been reported to contain major plant-specific TFs involved in flower development [72][73][74]. Interestingly, our results showed that the identified gene was down- Our results included many transcription factors, such as GATA, Dof, CO, and MADS families, that are known to be associated with flowering time (Table S5). Recently, a study on the expression pattern of each transcription factor in B73 was carried out, and a total of five hub genes (GRMZM2G113779, SPL; GRMZM2G387528, bHLH; GRMZM2G017349, bHLH; GRMZM2G011357, C2H2; GRMZM2G025642, NAC) were identified. Among the five hub genes, three genes (GRMZM2G113779, GRMZM2G387528 and GRMZM2G011357) were expressed in the leaves of B73. Our results confirmed the expression of three hub genes (GRMZM2G161673, homolog of GRMZM2G113779; contig1816, homolog of GR-MZM2G387528; GRMZM2G011357) under drought conditions. GRMZM2G113779 is a hub gene that encodes SPL. The SPL family has been reported to contain major plant-specific TFs involved in flower development [72][73][74]. Interestingly, our results showed that the identified gene was downregulated under drought conditions. Other hub genes that are associated with drought tolerance were upregulated under drought conditions.

qRT-PCR Validation of Differentially Expressed Transcripts
To confirm our RNA-seq results, the expression of drought-responsive genes was validated by using quantitative real time PCR (qRT-PCR). In total, 11 genes were selected from the identified DEGs that were predicted to delay flowering time. The patterns of expression from RNA-seq were replicated by the qRT-PCR results for all analyzed genes ( Figure 4). regulated under drought conditions. Other hub genes that are associated with drought tolerance were upregulated under drought conditions.

qRT-PCR Validation of Differentially Expressed Transcripts
To confirm our RNA-seq results, the expression of drought-responsive genes was validated by using quantitative real time PCR (qRT-PCR). In total, 11 genes were selected from the identified DEGs that were predicted to delay flowering time. The patterns of expression from RNA-seq were replicated by the qRT-PCR results for all analyzed genes ( Figure 4).

Discussion
The flowering time of tropical maize is known to be more dependent on photoperiods than temperate maize lines [47]. The reference genome B73 is a temperate line that was developed under long-day conditions, whereas Ki11 is a tropical line of maize and has adapted to short-day conditions. Ki11 is a drought-tolerant cultivar; its ASI is about 2.22 days under normal conditions [24] and experiences small changes under drought [51,53]. The ASI of tropical maize is increased by drought stress, resulting in a delay in silking [61]. As in previous studies, the ASI in this study increased by about 4 days under drought conditions, manifesting as a delay in silking (Table 1). To determine whether changes in the ASI were due to changes in the expression of drought-responsive genes, we analyzed genetic changes through RNA-seq and bioinformatics.
To reveal candidate genes that may cause changes in the ASI under drought, we analyzed previous reports on drought and flowering-time genes, including homologs in other plants [27], expression profiles of drought stress-responsive genes at each stage and in each tissue [19], and changes in the ASI under drought and short-day conditions [30]. We compared the DEGs of Ki11 to those of B73 under the same conditions. A total of 130

Discussion
The flowering time of tropical maize is known to be more dependent on photoperiods than temperate maize lines [47]. The reference genome B73 is a temperate line that was developed under long-day conditions, whereas Ki11 is a tropical line of maize and has adapted to short-day conditions. Ki11 is a drought-tolerant cultivar; its ASI is about 2.22 days under normal conditions [24] and experiences small changes under drought [51,53]. The ASI of tropical maize is increased by drought stress, resulting in a delay in silking [61]. As in previous studies, the ASI in this study increased by about 4 days under drought conditions, manifesting as a delay in silking (Table 1). To determine whether changes in the ASI were due to changes in the expression of drought-responsive genes, we analyzed genetic changes through RNA-seq and bioinformatics.
To reveal candidate genes that may cause changes in the ASI under drought, we analyzed previous reports on drought and flowering-time genes, including homologs in other plants [27], expression profiles of drought stress-responsive genes at each stage and in each tissue [19], and changes in the ASI under drought and short-day conditions [30]. We compared the DEGs of Ki11 to those of B73 under the same conditions. A total of 130 genes with changes in expression in response to drought were found to be common to both cultivars. Among DEGs, five genes (GRMZM2G012717, GRAMZM2G142718, GR-MZM2G135013, GRMZM2G094990, and GRMZM2G104549) appear to have the potential to change the ASI. Four of these genes (excluding GRMZM2G142718) were expressed in Ki11, and GRMZM2G142718 was found to be expressed in two homologous forms (contig312 and contig4777). GRMZM2G104549 was upregulated under drought in B73 but downregulated in Ki11, and other genes showed the same expression pattern. Thus, few flowering-time genes showed similar patterns between B73 and Ki11 under drought and SD conditions, suggesting that the flowering mechanisms of tropical and temperate lines are independent.
Furthermore, we compared the flowering-time genes of rice and maize under SD. Our results suggest that changes in the expression of florigens may cause a delay in flowering times. ID1 (GRMZM2G011357), which promotes the expression of ZCN8, is reported to be a regulator of flowering time; the overexpression of ZCN8 causes early flowering, while its downregulation results in late flowering [31,44,75,76]. Our data show that the ID1 gene was upregulated under drought. Other major genes related to flowering time were downregulated by drought. Two ZCN8 genes (GRMZM2G179264 and GRMZM2G103666), which are crucial factors for flowering time, were downregulated under drought conditions. The expression of ZCN8 has been reported to increase at the reproductive stage and interact with DELAYED FLOWERING 1 (ZmDLF1) under SD conditions in tropical maize [31,45], but ZmDLF family genes were not detected among DEGs. Another flowering-time gene, Hd3a (contig6971, homolog of GRMZM2G373928), was also downregulated under drought.
Hd3a is known to be activated by HD1 and to accelerate flowering time under SD [77]. Similarly, Dof genes were mostly downregulated under drought. Dof is a repressor of CONSTANS and CO-like transcription factor, a regulator of photoperiodic flowering [78,79]. These transcriptome changes due to drought stress may delay silking and increase the ASI of Ki11. However, the detailed flowering time mechanism in SD conditions remains unclear, so we could not perform a more detailed analysis. Based on our results and the known flowering pathway under SD conditions, a decrease in the expression level of Ehd1 and CONZ1 may lead to a decrease in the expression level of ZNC8, which may cause a delay in the flowering time. Although we cannot confirm this possibility because of the insufficient understanding of the genetic pathway of anthesis and silk emergence [63], a decrease in ZCN8 expression may be indirect evidence of the mechanism underlying the delay of silk production.
With the development of more advanced technologies and bioinformatics, flowering time mechanisms are being rapidly characterized. We aimed to use recent data to analyze flowering-time genes and elucidate the mechanism by which drought increases the ASI. However, a difference in the flowering time mechanisms between temperate and tropical maize was evident, and research on the flowering time of tropical maize is not yet sufficient to fully clarify these pathways. We investigated changes in genes related to the flowering time in this study and identified factors involved in the delayed flowering time of tropical maize under drought. Although we cannot confirm the exact mechanism because of the current insufficient understanding of the genetic pathway of anthesis and silk emergence, our results provide candidate genes for further study, which will help to improve the current understanding of the flowering time of tropical maize and its delay under drought.
Supplementary Materials: The following are available online at https://www.mdpi.com/article/10 .3390/app11094273/s1, Figure S1: Principal component analysis under control and drought stress conditions; Figure S2: Volcano plot of the distribution of differentially expressed genes (DEGs); Table S1: Summary of RNA-seq data from control and drought samples; Table S2: Summary of Blast results of drought-responsive genes; Table S3: Enrichment analysis of drought-responsive genes and transcripts; Table S4: Results of alternative splicing; Table S5: Results of Blast in a plant transcription factor database.