Genome-wide identification and expression profile analysis of nuclear factor Y family genes in Sorghum bicolor L. (Moench)

Summary

Members of the plant Heme Activator Protein (HAP) or NUCLEAR FACTOR Y (NF-Y) are trimeric transcription factor complexes composed of the NF-YA, NF-YB and NF-YC subfamilies. They bind to the CCAAT box in the promoter regions of the target genes and regulate gene expressions. Plant NF-Ys were reported to be involved in adaptation to several abiotic stresses as well as in development. In silico analysis of Sorghum bicolor genome resulted in the identification of a total of 42 NF-Y genes, among which 8 code for the SbNF-YA, 19 for SbNF-YB and 15 for the SbNF-YC subunits. Analysis was also performed to characterize gene structures, chromosomal distribution, duplication status, protein subcellular localizations, conserved motifs, ancestral protein sequences, miRNAs and phylogenetic tree construction. Phylogenetic relationships and ortholog predictions displayed that sorghum has additional NF-YB genes with unknown functions in comparison with Arabidopsis. Analysis of promoters revealed that they harbour many stress-related cis-elements like ABRE and HSE, but surprisingly, DRE and MYB elements were not detected in any of the subfamilies. SbNF-YA1, 2, and 6 were found upregulated under 200 mM salt and 200 mM mannitol stresses. While NF-YA7 appeared associated with high temperature (40°C) stress, NF-YA8 was triggered by both cold (4°C) and high temperature stresses. Among NF-YB genes, 7, 12, 15, and 16 were induced under multiple stress conditions such as salt, mannitol, ABA, cold and high temperatures. Likewise, NF-YC 6, 11, 12, 14, and 15 were enhanced significantly in a tissue specific manner under multiple abiotic stress conditions. Majority of the mannitol (drought)-inducible genes were also induced by salt, high temperature stresses and ABA. Few of the high temperature stress-induced genes are also induced by cold stress (NF-YA2, 4, 6, 8, NF-YB2, 7, 10, 11, 12, 14, 16, 17, NF-YC4, 6, 12, and 13) thus suggesting a cross talk among them. This work paves the way for investigating the roles of diverse sorghum NF-Y proteins during abiotic stress responses and provides an insight into the evolution of diverse NF-Y members.

Keywords:

Biology and life sciences – Genetics – Gene expression – Genomics – Biochemistry – Nucleic acids – Plant science – Organisms – Eukaryota – Plants – Grasses – Computational biology – Physical sciences – Research and analysis methods – Animal studies – Experimental organism systems – Model organisms – Plant and algal models – proteins – DNA-binding proteins – Plant pathology – Database and informatics methods – Bioinformatics – Sequence analysis – Gene regulation – RNA – Non-coding RNA – Sequence motif analysis – Transcription factors – Regulatory proteins – Physics – Classical mechanics – Ecology and environmental sciences – Plant ecology – Plant-environment interactions – Ecology – Plant physiology – Plant defenses – Plant resistance to abiotic stress – Natural antisense transcripts – MicroRNAs – Genome complexity – Brassica – Introns – Mechanical stress – Arabidopsis thaliana – Sorghum – Thermal stresses

Introduction

Nuclear Factor Y (NF-Y), also known as heme activator protein (HAP) or CCAAT-binding factor (CBF) is a ubiquitous, complex, heterotrimeric transcription factor. It is evolutionarily conserved in all plants with three distinct subunits called NF-YA or HAP2, NF-YB or HAP3/CBF-A and NF-YC or HAP5/CBF-C [1]. The assembly of NF-Y is complex and occurs both in cytoplasm and nucleus. While NF-YA and NF-YC family members have a nuclear localization signal (NLS), NF-YB members generally lack the same and hence cannot be transported to the nucleus [2]. The NF-YA subunits are localised to the nucleus and bind with varying affinities to the CCAAT cis-elements in the promoter regions of the target genes [3, 4]. On the other hand, NF-YB and NF-YC subunits contain the conserved Histone Fold Domain (HFD) or Histone Fold Motif (HFM) and help in protein-DNA and protein-protein interactions [5, 4]. The HFD domain of NF-YB/YC is formed by 3 α-helices separated by 2 loops [6]. The α1-helices of both NF-YB and YC contain putative DNA-binding domains [4]. Via the HFDs, NF-YB and NF-YC form a heterodimer [5], which is critical for the translocation of NF-YB from the cytoplasm to the nucleus [7]. Several members of the NF-Y subfamilies play a vital role not only in a wide array of developmental processes but in tolerance to abiotic stresses as well. For example, they are involved in embryogenesis [8], ABA response and seed germination [9], abiotic stress tolerance [10, 11, 12], flowering time [13], primary root elongation [14], photosynthesis [15], endosperm development [16], and photomorphogenesis [17]. On the other hand, in leguminous plants, they are the key regulators of symbiotic root nodule development [18]. Ni et al. [19] reported that GmNF-YA3, a target gene of miR169, is a positive regulator of plant tolerance to drought stress in A. thaliana. Alam et al. [20] demonstrated that overexpression of OsHAP2E gene confers tolerance to drought and salt stresses with enhanced photosynthesis and tiller (stems produced in grass plants) number in comparison with wild-type rice plants. Transgenic rice showed resistance to Magnaporthe oryzae and Xanthomonas oryzae infections. While overexpression of NF-YA5 conferred drought stress tolerance in Arabidopsis [21], NF-YA1 in Arabidopsis resulted in post germinative growth arrest under salt stress [11]. Further, incorporation of several NF-YB and NF-YC genes improved drought stress tolerance in diverse plants like Arabidopsis, maize, poplar and rice [10, 22, 23, 24]. Thus, evidence is accumulating that NF-Y subunits act as key regulators of drought stress tolerance. In plants, NF-Y gene families comprise several paralogs. In Triticum aestivum, 37 paralogs have been described (10 NF-YA, 11 NF-YB, 14 NF-YC and 2 Dr1) [25]; in rice 28 NF-Ys (10 NF-YA, 11 NF-YB and 7 NF-YC) [26]; but later [27] reported 34 members in the same species. In Arabidopsis thaliana, NF-Y members exclude AtNF-YC11, B12, B13 (NC2 subfamily) and AtNF-YC10, C13, and B11 (Dpb3/4 subfamily), but include AtNF-YC12. Accordingly, Petroni et al. [3] pointed out that Arabidopsis contains 30 members of the NF-Y family, 10 from each family (NF-YA, NF-YB and NF-YC) though 36 were originally reported Siefers et al [13]. So, in the updated scheme, 30 members of NF-Y have been considered by Zhao et al. [28] in A. thaliana. In Brachypodium distachyon 36 (7 NF-YA, 17 NF-YB, and 12 NF-YC) [29]; in Brassica napus 33 (14 NF-YA, 14 NF-YB, 5 NF-YC) [30]; in Setaria italica 39 (10 NF-YA, 15 NF-YB and 14 NF-YC) [31]; in Glycine max 68 (21 NF-YA, 32 NF-YB, 15 NF-YC) [32]; in Prunus mume 29 [33]; in Ricinus communis 25 (6 NF-YA, 12 NF-YB and 7 NF-YC) [34]; in Citrus sinensis and C. clemantia 22 (6 NF-YA, 11 NF-YB, 5 NF-YC) [35] were characterised.

Sorghum bicolor is a semi-arid and the second most important staple food grain crop. It provides feed, fodder and fuel and shows genetic diversity [36]. Being a C₄ photosynthetic plant, it is adapted to moderate drought and high temperature. But, salinity and drought coupled with high temperature limit the production and yield stability in sorghum. Enhancing the final yields and productivity of crop plants is especially challenging to the researchers due to the unpredictable nature of drought stress conditions during the growing season and complex drought stress biology [37]. Identification and expression of various transcription factors for abiotic stress tolerance using qRT-PCR and their validation by overexpression or knockouts is therefore critical for developing improved crop varieties with tolerance to water limited conditions. Members of NF-Y subfamilies impart tolerance to a very wide spectrum of both biotic and abiotic stresses as mentioned above. The number of NF-Y genes that exist in sorghum and their detailed biological roles for multiple stress tolerance and ABA-responsiveness remains unexplored. In this study, we identified 42 NF-Y genes using in silico approaches and examined their expression patterns under salt, drought, ABA, cold and high temperature stresses. Our gene expression studies reveal that majority of NF-Y genes (39) exhibited response to high temperature stress. A large number of them (24) were also expressed under multiple stresses like cold, salt (22) and drought (20). Further, 20 SbNF-Ys showed upregulation under ABA stress which indicates their role in ABA-related pathway. Keeping in view of the aforementioned reasons, we aimed to understand how the SbNF-Y members regulate abiotic stresses in an ABA-dependent or independent manner which would further delve into investigating their detailed roles during stress.

Materials and methods

Plant material and abiotic stress treatments

Sorghum bicolor variety BTx623 is an agronomically important inbred line. It is a model variety with known genome sequencing information and moderately tolerant to drought stress. The gene space of the sorghum genome sequence has also been updated by resequencing [38]. Keeping these criteria in mind, seeds of S. bicolor variety BTx623 were obtained from ICRISAT, Patancheru, Hyderabad, and sown in pots filled with 5 kg of black soil and seedlings were grown in glass house conditions at 28/20°C day/night temperatures. Sixty-day-old seedlings were subjected to 200 mM NaCl solution, 200 mM mannitol solution, and 100 μM ABA for 4 h separately. Cold stress was imposed by keeping the plants at 4°C and high temperature stress by exposing the plants to 40°C for 4 h. Corresponding controls (without any treatment) were maintained under identical conditions. After 4 h of exposure, roots, stems and leaves were collected from treated and control plants and snap frozen in liquid nitrogen and stored at -80°C for subsequent use. Three biological and three technical replicates were used for qRT-PCR analysis.

Identification and characterization of NF-Y transcription factors

NFY gene sequences of Oryza, Zea, Setaria were retrieved from plantTFDB (http://planttfdb.cbi.pku.edu.cn/) (S1 Table) database and searched against Sorghum bicolor genome in Gramene database (http://www.gramene.org/) to find out their homologs. Genscan (http://genes.mit.edu/GENSCAN.html) program was used to retrieve the gene and their respective protein sequences. The identified putative Sorghum nuclear factors were scanned using HMMER (https://www.ebi.ac.uk/Tools/hmmer/search/hmmsearch) corresponding to the Pfam database and queried against the Oryza, Setaria and Zea. The identified Sorghum NF-Ys were confirmed by searching against Oryza, Setaria and Zea genomes in Gramene database (http://www.gramene.org/). Based on homology, the identified putative protein sequences were subjected to Motif Search (http://www.genome.jp/tools/motif/) analysis to check the reliability and to identify their conserved domains [39]. The identified NF-Y genes were mapped to their respective chromosomes based on the information provided in the Gramene Genome Database by employing MapInspect software (https://mapinspect.software.informer.com/). Gene Structure Display Server (http://gsds.cbi.pku.edu.cn) software was used for obtaining the NF-Y gene structures—exons, introns, and untranslated sequence regions (UTRs) based on the alignments of their coding sequences [40]. MEME software was employed to analyze new sequence patterns and their significance [41]. The software helps to identify the nature of motifs by setting different default parameters, number of motifs from 1–10 with a motif width of 5–50, and the number of motif sites from 5–10.

NF-Y protein analysis, prediction of potential cis-regulatory elements, identification of miRNA target sites and phylogenetic analysis of NF-Ys

Molecular weight (MW), isoelectric point (pI), and GRAVY (grand average of hydropathy) of NF-Ys were identified for all NF-Y proteins by using ProtParam of Expasy tools (http://web.expasy.org/protparam) [42]. Phosphorylation sites of proteins were predicted using NetPhos3.1 software of Expasy tools [43]. Subcellular localization of NF-Y proteins was carried out by WOLFPSORT program (http://wolfpsort.org/) [44]. To predict the putative cis-acting elements of NF-Y promoter regions, 2000 bp genomic sequences upstream of start codons were analysed using PLANTCARE software [45]. The pSRNATarget software [46] was employed to identify the potential miRNA target sites in identified SbNF-Ys. Finally, the neighbor-joining (NJ) phylogenetic tree was constructed with the NF-Y protein sequences of Sorghum bicolor with the plants as shown in S1 Table using MEGA 6.2 software [47]. The NJ is a recursive algorithm, a fast method which is suited for large datasets and does not require ultra metric data and permits correction for multiple substitutions. The Poisson correction, pairwise deletion, and bootstrap value (1,000 replicates) parameters were used to draw the NJ phylogenetic tree.

Phylogenetic divergence and co-expression analysis

Gene duplication events were found [48, 49] using phylogenetic tree based on 70% similarity and 80% coverage of sequences aligned. PAL2NAL program [50] was followed for finding out synonymous and non-synonymous substitutions rates. Protein-protein interaction (PPI) map of NF-Y proteins was generated from the STRING database [51].

RNA isolation and quantitative real-time PCR analysis

From the stress exposed and control (without any stress) samples, total RNA was isolated using Macherey-Nagel NucleoSpin RNA plant kit by following the instructions given in the manual. To eliminate any genomic DNA contamination in the RNA samples, the purity of RNA was checked using Eppendorf BioPhotometer. Two micrograms of RNA sample was used as template for first strand cDNA synthesis using RevertAid First Strand cDNA Synthesis Kit (#K1622, Thermo Scientific EU, Reinach, Switzerland). To find out the relative gene expression levels of SbNF-Ys, 2X Applied Biosystems (ABI) Master Mix with gene specific primers was used (S2 Table). For qRT-PCR analysis, thermal cycling conditions of 95°C for 5 min followed by 40 cycles of 95°C for 30 s, 57°C for 30 s and 72°C f or 30 s were applied to the ABI 7500 real-time PCR system (Applied Biosystems, Foster City, CA, USA). Expression of SbNF-Y genes in control and treated samples was normalized with EIF4a (Eukaryotic Initiation Factor 4A) and PP2A (protein phosphatase2A subunit A3) reference genes [52]. qRT-PCR was carried out with three biological and three technical replicates for each sample. The PCR reaction specificity was confirmed by melting curve analysis of the amplicons. Comparative 2-DDCT method [53] was used to calculate the relative quantities of each transcript.

Results

Identification and characterization of SbNF-Y transcription factors, motif analysis and subcellular localization

A total of 42 homologous genes comprising 8 NF-YA, 19 NF-YB and 15 NF-YC from the whole genome of sorghum were identified and confirmed. Later, they were crosschecked by using the HMM profile and searching SbNF-Ys against Oryza, Setaria and Zea for further confirm to check their reliability (S3 Table). The predicted 8 NF-YA, 19 NF-YB and 15 NF-YC genes were named as SbNF-YA1 to SbNF-YA8, SbNF-YB1 to SbNF-YB19 and SbNF-YC1 to SbNF-YC15 respectively. Based on the presence of conserved NF-YA, NF-YB and NF-YC domains, the predicted SbNF-Y family of proteins was considered for identification as a member. The exon-intron structures of all the 42 annotated NF-Y genes were analysed (Fig 1). While 17 exons and 16 introns (highest) were detected in SbNF-YA2, 3 exons (least) and 2 introns were found in SbNF-YA4 gene. Among the SbNF-YB family members, SbNF-YB11 showed a maximum of 16 exon and 15 intron regions, and 5 of the members displayed 1 intron. A maximum of 18 exons and 17 introns were noticed on SbNF-YC8. Also, six of the members were intronless and no member exhibited one intron (Fig 1). The sub-cellular localization of SbNF-Y proteins based on consensus sequence showed a majority of them to be localized to nucleolus and chloroplast although a few of them localized to cytoplasm, mitochondria and plastids (Table 1). All the NF-Ys showed nuclear localization signals (NLS); NF-YA holding LRRR sequence (motif 2 in Fig 2A), KRK motif in NF-YB (motif 1 in Fig 2B) and KRR in NF-YC (motif 1 in Fig 2C). Though they contain nuclear localization signals, their subcellular localizations were different. Majority of the SbNF-YAs showed chloroplast as the important target site. Few of them have been found localized in chloroplasts (NF-YA1, NF-YA7 and NF-YA8), mitochondria (NF-YA4) and plastid (NF-YA6). The number of phosphorylation sites in each NF-Y protein is represented in the S3 Table. All the NFYs exhibited higher number of PKC than CK1, CK2, and PKA types. The PKC number is higher in NF-YA subfamily members than in NF-YB (S4 Table). No transmembrane helices were observed except in NF-YA3 protein.

The identified SbNF-Y genes encoded polypeptides with amino acids ranging from 130 to 1430 and pI values varied from 4.26 to 10.83. Characteristically, they showed DNA binding domains. Molecular weights of the proteins ranged from 10.21 to 83.52 kDa (Table 1). Among the SbNF-YA subfamily members, RKPYLHESRHLHAMKRARGSGGRFLNTKQ and EEPIYVNAKQYNAILRRRQARAKLEAZNK large contiguous motifs were found ubiquitous, while rest of the 8 large contiguous motifs showed variability in their distribution (see Figs 2A and S1). Similarly, SbNF-YB proteins revealed the uniform presence of one, highly conserved large contiguous motif, i.e. AKETVQECVSEFISFVTGEASDKCQREKRKTINGDDLLWALATLGFEDYY (Figs 2B and S2). On the other hand, analysis of NF-YC proteins revealed a highly conserved large contiguous motif APVVFAKACEFIQELTLRAWHEENKRRTLQKSSDIAAAIARTEVYDFL (see Figs 2C and S3). Motif analysis of complete SbNF-Y family representing conserved motifs (S4 and S5 Figs) reflect typical diagnostic features for different subunits of NF-Y family proteins in general and hence provide confirmatory identification of SbNF-Y proteins from the sorghum genome.

Phylogeny, divergence, and physical genome mapping of SbNF-Ys

The phylogenetic tree of SbNF-Y proteins was constructed using MEGA 6.2 software. It showed 2 clades, which are subdivided into 4. The phylogenetic analysis displayed a total of 11 paralogous duplication events of which 3 tandem/segmental (SbNF-YB4/SbNF-YB6; SbNF-YB18/SbNF-YB19; SbNF-YB12/SbNF-YC11) and remaining regional duplication events. Interestingly, NF-YBs exhibited paralogous events with SbNF-YCs (SbNF-YB13/SbNF-YC1; SbNF-YB12/SbNF-YC11), indicating their evolutionary relatedness (Figs 3 and S6 and Table 2). To find out the orthologous and the evolutionary relationships of SbNF-Ys, NF-YA, NF-YB and NF-YCs have been compared with other plant genomes (see Fig 4A, 4B and 4C respectively). Not surprisingly, a majority of the identified SbNF-Ys showed orthologous relationship with Zea, few with Setaria and one with Hordeum (Sorbi009G166200 (SbNF-YB15)/MLOC_36879.2). Of the 8 SbNF-YAs, 5 showed orthologous events with Zea and 1 with Setaria (Fig 4 and S5 Table). The major subfamily SbNF-YBs exhibited 11 orthologous events, of which 6 showed with Zea, 4 with Setaria and 1 with Hordeum (Fig 5 and S6 Table). On the other hand, SbNF-YCs showed 10 orthologous events, of which 8 with Zea and 2 with Setaria (Fig 6 and S7 Table).

The identified SbNF-Ys were distributed across all the 10 chromosomes. A maximum number of 7 genes each were located on chromosome 1, and 7, 5 genes each on 2, 3, and 9, 4 each on 8 and 10, 3 genes on 4, and 1 each on 5 and 6 chromosomes (Fig 7 and Table 1). Among the NF-YA subfamily, SbNF-YA1, SbNF-YA2 and SbNF-YA3 genes were located on chromosome 1. Chromosomes 2 and 8 have two genes each located on them, i.e. SbNF-YA4, SbNF-YA5 and SbNF-YA7, SbNF-YA8 respectively. Chromosome 4 is having only SbNF-YA6 localized on it. Among the SbNF-YB genes, a maximum of 4 genes each were located on chromosomes 3 and 9. While SbNF-YB4, SbNF-YB5, SbNF-YB6 and SbNF-YB7 were observed on chromosome 3, SbNF-YB13, SbNF-YB14, SbNF-YB15 and SbNF-YB16 were noticed on chromosome 9. Three genes each SbNF-YB10, SbNF-YB11, SbNF-YB12 and SbNF-YB17, SbNF-YB18, SbNF-YB19 were seen on chromosomes 7 and 10 respectively. Chromosomes 2 and 4 have 2 genes each, i.e., SbNF-YB2 and SbNF-YB3 on 2, SbNF-YB8 and SbNF-YB9 on 4. Chromosome 1 has only SbNF-YB1 located on it. Majority of the SbNF-YC genes were located on chromosome 7, and it accommodates 4 genes (SbNF-YC8, SbNF-YC9, SbNF-YC10, and SbNF-YC11). While chromosome 1 accommodates three genes (SbNF-YC1, SbNF-YC2, and SbNF-YC3), chromosome 8 contains SbNF-YC12 and SbNF-YC13. One gene each SbNF-YC4, SbNF-YC5, SbNF-YC6, SbNF-YC7, SbNF-YC14, and SbNF-YC15 was noticed on chromosomes 2, 3, 5, 6, 9, 10 respectively (Fig 7 and Table 1).

Estimation of non-synonymous and synonymous substitution rates

The non-synonymous (d_N) and synonymous (d_S) substitution (d_N/d_S) rates were calculated for genes which showed duplication events within Sorghum as paralogs and between other genomes as orthologs. The 11 paralogs (S6 Fig) exhibited the d_N/d_S between 0.0010 (SbNF-YB18/SbNF-YB19)-93.9760 (SbNF-YC8/SbNF-YC12). Of the 11 paralogs, only 5 showed the purifying/stabilizing selection (<1), while the remaining exhibited positive/Darwinian selection (>1) (Table 2). The SbNF-YA orthologs exhibited d_N/d_S rates ranging from 0.4000 (Sorbi001G340200/Seita.9G367200) and 28.6781 (Sorbi001G486000/Zm2G000686). This indicates that 2 were following purifying selection and the remaining 4 positive selection (S5 Table). The d_N/d_S rates of SbNF-YB orthologs varied from 0.0031 (Sorbi001G338700/Seita. 9G365700) and 10.2227 (Sorbi010G119200/ZM2G167576). While majority of them (7) evolved through purifying selection, the remaining evolved by positive selection mechanism (S6 Table). The orthologs of SbNF-YCs showed d_N/d_S rates between 0.0797 (Sorbi002G241500/Seita. 2G247700) and 1.9761 (Sorbi007G063200/ZM2G311316) (S7 Table).

Promoter analysis

Analysis of promoter sequences revealed cis-acting elements such as ABA-responsive (ABRE), drought-responsive (DRE, DPBF, MYB and MYC), heat shock-responsive (HSE), and low temperature-responsive (LTR) elements. Aside, methyl jasmonic acid- (MeJA-RE), salicylic acid- (SARE), and defence-responsive elements (TC-rich repeats), associated with biotic stress were detected. Majority of the genes in SbNF-YA family have G-BOX (CACGTG) cis-acting elements (S8 Table). In the SbNF-YB family, G-BOX and Sp1 were observed as major cis-acting elements except in SbNF-YB18 (S9 Table). On the other hand, Skn-1 motif is the most common cis-acting element in SbNF-YC family. The least expressed cis-acting elements in SbNF-YC are AUXRR-CORE and GCN4 motif I (S10 Table).

Protein-protein interaction (PPI) prediction analysis

The PPI mapping of SbNF-Ys showed that a cohort of proteins involved in various cellular, metabolic and molecular pathways are associated with miRNA surveillance pathway, DNA replication, base excision, nucleotide excision repair pathway, purine and pyramidine metabolism. They interacted with core histones, calcineurins, kelch motifs, serine-threonine phosphatases, histone lysine N-methyl transferase, and metal-dependant phosphatase (S7 and S8 Figs).

In silico prediction of miRNA target sites

The SbNF-YAs exhibited different miRNA target sites such as sbi-miR169, sbi-miR5389, sbi-miR6225, sbi-miR5568, sbi-miR6220, sbi-miR5567 and sbi-miR6232, SbNF-YBs showed sbi-miR5565, sbi-miR5568, sbi-miR6232, sbi-miR6220, sbi-miR821, sbi-miR437, sbi-miR528, sbi-miR395, sbi-miR169, sbi-miR171 and sbi-miR172. Likewise, SbNF-YCs showed sbi-miR6232, sbi-miR6235, sbi-miR395, sbi-miR437, sbi-miR395, sbi-miR6220, sbi-miR6225, sbi-miR6227, sbi-miR5569, sbi-miR5568, sbi-miR164, sbi-miR156, sbi-miR160, sbi-miR6218, sbi-miR6230, and sbi-miR821. Interestingly, all of them are known to be associated with translation and cleavage events (S11, S12 and S13 Tables).

Gene expression analysis of SbNF-Ys in different tissues treated with diverse abiotic stresses

Expression of all the 42 NF-Y family of genes was studied at the transcriptional level in different tissues and abiotic stress conditions besides ABA, and the heat map is presented in Fig 6A and 6B. Over all, the gene expressions were higher in leaf tissues in comparison with stem and root (Fig 8A and S8 Table). Among the 8 SbNF-YA genes, NF-YA6 appeared to be associated with salt, drought (imposed by mannitol), cold and high temperature stresses. It was also strongly triggered by ABA. On the other hand, NF-YA1 was induced by multiple stresses like salt, mannitol and high temperature stresses, but NF-YA3 and NF-YA5 were expressed only under high temperature. While salt and high temperature influenced the expression of NF-YA7, NF-YA8 was upregulated by cold and high temperature stresses. Among the 19 SbNF-YB members, SbNF-YB7, B12, B15 and B16 were strongly induced by different stresses like salt, mannitol, ABA, cold and high temperature. In contrast, seven out of fifteen SbNF-YC members, YC1, YC3, YC4, YC7, YC8, YC9, and YC13 were expressed only under high temperature stress. Further, five members (YC6, YC11, YC12, YC14 and YC15) were upregulated by multiple stresses including ABA (Fig 8B and S14 Table).

Discussion

Identification and structural analysis of SbNF-Y genes

It is observed from this study that the number of genes are variable in each of the three distinct subfamilies of NF-Y (NF-YA, NF-YB and NF-YC) across different taxa. NF-Ys are evolutionarily conserved in eukaryotes and each subunit is encoded by a single gene in yeast and animals [4]. But, the same is encoded by a family of genes varying from 8 to 39 in plants. A total of 33 NF-Ys were reported in Brassica napus (14 NF-YA, 14 NF-YB, 5 NF-YC) [30]; while 39 in Setaria italica (10 NF-YA, 15 NF-YB and 14 NF-YC) [31]; whereas 68 in Glycine max (21 NF-YA, 32 NF-YB, 15 NF-YC) [32]. Multiple members of NF-Y subunits in plants reflect the redundancy and differentiated functions of these genes which need to be explored by expression profiling. Pereira et al. [35] identified 22 NF-Y genes in Citrus sinensis and C. clementina (6 NF-YA, 11 NF-YB and 5 NF-YC). Xu et al. [54] pointed out that such a low number of NF-Y genes found in Citrus genome could be due to the whole genome duplication events occurred in A. thaliana when compared to Citrus. It appears that subunit YB has more number of genes in comparison with YA and YC members in most of the species including the present report. Single genes, but with multiple splicing isoforms (that encode NF-Y subunits) is generally noticed in the yeast and mammals. Contrary to this, in higher plants, multigene families are noticed which can encode each subunit. Such multiple genes are vital for plant systems for tissue specific expressions at various stages of growth and development. Further, such a subunit combination can assist the plant systems in performing diverse roles during stress/development. Malaviya et al. [55] searched the plant transcription factor database (Plant TFDB, (http://plantfdb.cbi.pku.edu.cn/) version 3.0 Jin et al. [56]), for identifying NF-Y genes in sorghum. They identified a total of 33 NF-Y transcription factors comprising 8 NF-YA, 11 NF-YB, and 14 NF-YC subunits in sorghum using Plant TFDB. In contrast, in the present study, a total of 42 NF-Y genes, among which 8 code for SbNF-YA, 19 for SbNF-YB and 15 for SbNF-YC subunits were identified. This discrepancy is because, in the present study, sorghum genome sequence available in the public domain has been searched.

Koralewski and Krutovsky [57] pointed that finding out exon-intron organization is crucial since it provides an insight into evolutionary relationships among genes and organisms. Malviya et al. [55] reported no introns in 18 of the genes (out of 33), and 5 of them have only one intron. They reported 2 in NF-YA3, 5 in NF-YA5, 4 in NF-YA6, 3 in NF-YA8, 4 in NF-YB1, 5 in NF-YC4 and 3 in YC7. Interestingly, in the present study, introns were absent in 12 out of 42 NF-Y TFs. While NF-YC18 contained 17 (highest number), YA2 16, YB11 15, YA8 13, YB16 11, YC10 10, YA7, YB15, YC2, YC7 and YC15 8 introns each. Only YB4, YB8, YB14, YB17, and YB18 (in all 5) contained one intron. Like in S. bicolor, single intron genes were also noted in Medicago truncatula which lead to alternative spliced variants of NF-YA1 [58]. Similarly, one intron in the 5′-UTRs of the NF-YA members was observed in A. thaliana, O. sativa, C. sinensis [28, 27, 35]. This suggests that such a post-transcriptional regulatory mechanism is retained among NF-YA genes. Single intron NF-Ys were not observed in SbNF-YC subtype in the present analysis. Chen et al. [59] reported that most of the NF-YB contained only one exon, and the genes from the same clade displayed a similar motif pattern in Gossypium hirsuyum. Chu et al. [60] reported 5 exons and 4 introns (6 genes) or 6 exons and 5 introns (2 genes) in CaNF-YA gene family members in Cicer arietinum. Further, they noticed 1 to 6 exons in CaNF-YB family, and 7 intronless out of 11 members in the CaNF-YC family. They reported 1 intron in NF-YC1, and 3 in NF-YC9. This suggests that a post-transcriptional regulatory mechanism is retained among NF-YA genes. Thus, the presence of multiple exon/intron gene organizations have been found in all the NF-Y family members in other species like B. napus [30], and S. lycopersicum [61] also. This infers that the presence of exon/intron is an attribute and typical of NF-Ys in higher plants. Loss or gain of spliceosomal introns led to the progress in our understanding of the molecular mechanisms associated with intron evolution and variation in gene function [62]. Fusion of exons and intron loss, might play a key role in the evolution of larger families like NF-Ys. Further, several members of the NF-YB and NF-YC have been found without any introns like in S. bicolor [55], Ricinus cummunis [34] and chickpea [60]. Introns are essential parts of all eukaryotic genes. In eukaryotic systems, introns are known to execute several functions like exon shuffling [63], gene expression alterations [64] and also tune the evolutionary rate of genes [62].

Motif identification and chromosomal localization of SbNF-Ys

NF-Y proteins display both conserved and non-conserved regions in Arabidopsis and others. Such conserved sequences may be vital for DNA interactions at CCAAT sites as pointed out by Siefers et al. [13], Romier et al. [2], and Testa et al. [65]. Hahn et al. [66] demonstrated that the yeast CCAAT box factor is a heteromer that contains HAP2 and HAP3 proteins. Xing et al. [67] showed that HAP2 is a 21 residue region with 3 histidines and arginines. Both SbNF-YB and SbNF-YC proteins have histone domains, but not SbNF-YAs. Besides, they also contain centromere kinetochore components and chromatin reorganizing domains. These residues are conserved in all the 8 SbNF-YA proteins (present study) as well as in Oryza sativa, and Triticum aestivum [23, 25]. Romier et al. [2] demonstrated that NF-YC/NF-YB sub-complex interacts through histone fold motifs. The role of the alpha C-helix of NF-YC appears to be vital for trimerization as well as a target for regulatory proteins like that of MYC and p53. It looks that heterotrimeric NF-Y proteins recognize the CCAAT regulatory elements represented in promoter and enhancer regions and modulate the genes. Steidl et al. [68], Liu and Howell [69] pointed out that NF-YB and NF-YC form a dimer in the cytoplasm and then translocated to the nucleus to interact with that of NF-YA to form a heterotrimer complex. Further, it has been demonstrated that bZIP28 and NF-Y transcription factors are activated by endoplasmic reticulum stress and assemble into a transcriptional complex to regulate downstream stress response genes in A. thaliana [69]. Alpha helix transmembrane spans with average hydrophobicity were predicted in 12 of the NF-Y proteins in S. bicolor. Anchoring of protein to glycosylphosphatidyl inositol (GPI) via the C-terminal attachment was predicted in three of the NF-YB proteins namely SbNF-YB1, YB3 and YB6. This appears rational since earlier gene fusion experiments conducted by Caras et al. [70] demonstrated that the C-terminal signal sequence has GPI-anchoring residues. The distribution of NF-Y genes appears to be widespread among different chromosomes. While OsHAP genes were dispersed on 11 out of the 12 rice chromosomes [26], in S. bicolor, they are distributed on 10 chromosomes.

Phylogenetic assessment, divergence and promoter analysis

Among the NF-YA family members, A2 and A5 appeared on the same clade indicating that they are closer to each other compared to others. While Malviya et al. [55] found that SbNF-YB8 was closer to SbNF-YA and SbNF-YC proteins, we could not observe such a correlation. On the contrary, B12 was observed closer to C11 and B13 to C1 in the present study than YA family members. It appears therefore YB and YC members might have close correlations in comparison with other members. Malviya et al. [55] noticed several ortholog and paralog groups through the phylogenetic analysis of SbNF-Y proteins along with 36 Arabidopsis and 28 rice NF-Y proteins. Malviya et al. [55] reported that Sorghum NF-Y family gene expansion is due to segmental duplication events. It appears that SbNF-Y genes retained their function even after duplication. Generally, gene family expansion occurs through segmental, tandem duplications, and transposition events [71]. In the present investigation, 11 paralogs were observed due to 3 regional duplications, and 8 segmental duplications, inferring that segmental duplications are responsible for SbNF-Y gene family expansion. Six duplication events were observed in SbNF-YB family, and this is a large number when compared to other subfamilies. SbNF-YB4, B5, B13 and B14 were phylogenetically distinct from other SbNF-YBs, and might have formed by recent duplications. SbNF-Ys exhibited 20 orthologous events with Zea, 7 with Setaria and 1 with Hordeum, which indicates their monocot ancestors. The synonymous (d_S) and nonsynonymous (d_N) substitutions reveal the selective pressure on duplicated genes. In the present study, phylogenetic relationships and ortholog predictions displayed that sorghum has additional NF-YB genes with unknown functions in comparison with Arabidopsis. The synonymous (d_S) and nonsynonymous (d_N) substitutions reveal the selective pressure on duplicated genes. Nekrutenko et al. [72] pointed out that greater than 1 d_N/d_S value represents positive selection, less than 1 functional constraint, and equal to 1 neutral selection. In the present study, it appears that majority of the duplicated genes evolved through purifying selection. The phytohormone-responsive cis-acting elements make the plants to tolerate various environmental changes. The ABRE play an important role in ABA signalling and abiotic stress tolerance. In the present investigation, large number of ABA-responsive elements were observed in majority of NF-Ys besides Skn elements that participate in endosperm expression [73]. Further, the presence of light-responsive elements like SP1, I-Box, and G-BOX indicate their roles in the regulation of gene responses to light. Interestingly, all the elements are rich with heat shock elements (HSE), which indicates their diverse roles in various stress response mechanisms.

miRNA analysis and protein-protein interactions

It is known that stress-responsive miRNAs target the transcription factors, which regulate the plant growth and development. The rationale behind finding out miRNA target sites is to know if any miRNAs associated in the regulation of SbNF-Ys exist in the genome. Identifying the target sites would subsequently help us in elucidating the regulation of SbNF-Ys during salt, drought and high temperature stress conditions. The miRNAs may also involve in gene networks regulated by transcription factors like NF-Ys. Identifying the interactions between miRNAs and transcription factors like NF-Ys will serve to screen their roles in stress tolerance, signal transduction, different developmental stages and synthesis of secondary metabolites, which will help to develop desired phenotypes with stress tolerance. While Fang et al. [74] reported targeting of the NAC mRNA by miRNA for abiotic stress responses, Stief et al. [75] noticed down regulation of heat stress memory by another miRNA. In the present investigation, miR169 identified was known to participate in post transcriptional regulation [76, 77, 19, 78]. Furthermore, miR169 and NF-YA5 knockout plants showed hypersensitivity to drought indicating their importance in drought tolerance [21]. Overexpression of miR169c in tomato enhanced the drought tolerance by reducing stomatal opening [79]. Therefore, in silico screening for miRNAs and their validation for abiotic stress response is highly important especially in the context of non-coding RNAs playing a gamut of regulatory roles. In addition, the PPI analysis revealed that they interact with calcineurins, the calcium sensors that usually confer spatial specificity in Ca²⁺ signalling, and play important roles in abiotic stress tolerance [80]. NF-Ys also participate in circadian clock and flowering time regulation, serine/threonine phosphatases and metal-dependant phosphatases and control the dephosphorylation of phosphoprotein substrates [81].

Gene expression analysis in different sorghum tissues under abiotic stress conditions

Analysis of NF-Y gene expressions by qRT-PCR indicated tissue-specific and stress-inducible expression profile. NF-YA5, A6, B7, B12, B15, B16, C6, C11, C12, C14 and C15 revealed significant differential expression patterns in response to the abiotic stresses in S. bicolor. Such a tissue-specific expression pattern was earlier noticed in several plants [30, 82]. This may indicate a sub-functionalization of different members in specific tissues under different abiotic stress conditions. Pereira et al. [35] pointed out that CsNF-YA2, CsNF-YB5/11 and CsNF-YC2/3 could form potential complexes in the citrus fruit. Many NF-Y genes were reported to be associated with both biotic and abiotic stresses. Xu et al. [83] reported high expression of BnNF-YA10 and BnNF-YB3, BnNF-YB7, BnNF-YB10 and BnNF-YB14 under NaCl stress. Under polyethylene glycol treatment, expression of BnNF-YA9, 10, 11 and 12 genes increased in B. napus. Malviya et al. [55] performed in silico gene expression analysis under abiotic stress conditions using rice transcriptome data. This revealed several of the sorghum NF-Y genes are associated with salt, drought, cold and temperature stresses. Since such an analysis is based on rice transcriptome database, this cannot give accurate results. But, in the present study, detailed gene expression studies were carried out and the results indicate that SbNF-YA1, 2, and 6 are upregulated under 200 mM salt and 200 mM mannitol stresses. NF-YA7 has been found associated with high temperature (40°C) stress, but NF-YA8 is triggered by both cold (4°C) and high temperature stresses. Among NF-YB genes, 7, 12, 15, and 16 are induced under multiple stress conditions such as salt, mannitol, ABA, cold and high temperatures. Likewise, NF-YC 6, 11, 12, 14, and 15 have been found enhanced significantly in a tissue specific manner under multiple abiotic stress conditions. Thus, the present analysis revealed that several of the NF-Ys are implicated in abiotic stresses and also modulated by ABA. Such a modulation of the NF-Ys by ABA was not shown by Malviya et al. [55]. Zhang et al. [84] found that many Physcomitrella patens NF-Y genes were responsive to abiotic stresses through ABA-dependent or independent pathways. In the present study, several genes were upregulated when treated with ABA, indicating that they are ABA-dependent. It has been observed from the present study that majority of the mannitol (drought)-inducible genes were also induced by salt, high temperature stresses and ABA. Few of the high temperature stress-induced genes are also induced by cold stress (NF-YA2, 4, 6, 8, NF-YB2, 7, 10, 11, 12, 14, 16, 17, NF-YC4, 6, 12, and 13). Seki et al. [85] noticed that drought-inducible genes are also inducible by salt stress and ABA treatments in A. thaliana. Ha et al. [86] observed that diverse transcription factor families modulate plant responses to abiotic stresses independent of ABA or dependent of ABA [87, 88]. Several members of the TFs also function in both ABA-dependent and independent ways [89–91]. Interestingly, such a crosstalk can be achieved via indirect interactions between TFs and cis-elements present in the same promoter regions of the target genes [92].

Quach et al. [32] reported involvement of soybean NF-Y genes in specific developmental stages and also stress responses. In Prunus mume, Yang et al. [33] observed high expression of PmNF-YA1/2/4/5/6, PmNF-YB3/4/8/10/11/13, and PmNF-YC1/2/4/5/6/8 under osmotic stress and ABA. In citrus, CsNF-YA5 and CsNF-YB1/2/4/5/11 were found upregulated by drought stress [35]. Such a finding was proved later by overexpression of AtNF-YB1 in Arabidopsis and its ortholog ZmNF-YB2 in maize which showed enhanced drought tolerance [93]. Similarly, overexpression of osmotic and ABA-inducible NF-YB genes PwNF-YB3 from Picea and PdNF-YB7 from poplar in Arabidopsis exhibited improved drought tolerance activity [94, 95]. Transgenic rice plants harbouring bermudagrass NF-YC gene showed tolerance under drought [95]. NF-Y genes participate in stress tolerance mechanism by interacting with other stress inducible genes like antioxidants. The connection between NF-Ys and antioxidants was observed in previous reports; CsNF-YA5 [35], AtNF-YA5 interacts with glutathione S-transferase, peroxidases and an oxidoreductase [94] and SiNF-YA1 enhance the activity of superoxide dismutase, peroxidase and catalase [94]. Expression profiles exhibited by paralogous SbNF-Y genes in different tissues of sorghum under stress treatments suggest a clear functional redundancy among this gene family members. It is interesting to study how these NF-Ys regulate the expression of downstream genes that perform a wide spectrum of functions. Siefers et al. [13] pointed out that some transcription factors control gene expression by binding to cis-regulatory elements as individual subunits. But, it also appears that others are deployed in a combinatorial fashion both spatially and temporally.

Conclusions

Genome-wide screening revealed the existence of a total of 42 NF-Y genes (8 SbNF-YA, 19 SbNF-YB and 15 SbNF-YC subunit members) in Sorghum bicolor. In silico analysis of promoters revealed that they comprise many stress-related cis-elements such as ABRE and HSE indicating their role in salt, drought and high temperature stress responsiveness. The tissue specific expression of NF-Y transcription factors under salt, drought, ABA, cold and high temperature indicated their role in multiple stress tolerance. In view of this, we firmly believe that our studies have allowed identifying the candidate genes for further validation under an array of abiotic stress conditions in a crop species.

Compliance with ethical requirement

Authors do not have any other interests that influence the results and discussion of this paper. The authors have read the Journal’s policies and the authors of this paper have the following competing interests. RP is the President & CEO of Genomix Molecular Diagnostics Pvt Ltd., Kukatpally, Hyderabad, India. RP is the CEO of Genomix CARL Pvt. Ltd., Andhra Pradesh, India, but does not receive a salary in this capacity. RP is the President and CEO of Genomix Biotech Inc., 2620 Braithwood Road, Atlanta, GA 30345, USA, but does not receive any salary in this capacity. There are no patents, or products in development or marketed products associated with this research to declare. This does not alter our adherence to PLOS ONE policies on sharing data and materials.

Supporting information

S1 Fig [ppt]
1–10 MEME identified motif sequences of .

S2 Fig [ppt]
1–10 MEME identified motif sequences of .

S3 Fig [ppt]
1–10 MEME identified motif sequences of .

S4 Fig [ppt]
Distribution of 1–10 MEME identified SbNFY-A, B, and C conserved motifs.

S5 Fig [ppt]
1–10 MEME identified motif sequences of , , and .

S6 Fig [pptx]
Gene duplications of .

S7 Fig [ppt]
Network analysis of SbNF-Ys.

S8 Fig [ppt]
Functional partners of SbNFYs and their role in interaction network.

S1 Table [doc]
List of plants searched against .

S2 Table [doc]
gene specific primers used in the gene expression analysis.

S3 Table [xlsx]
The reliability of identified .

S4 Table [doc]
Number of phosphorylation sites in NFYs.

S5 Table [doc]
Non-synonymous to synonymous substitution ratios of orthologs.

S6 Table [doc]
Non-synonymous to synonymous substitution ratios of orthologs.

S7 Table [doc]
Non-synonymous to synonymous substitution ratios of orthologs.

S8 Table [doc]
Conserved -acting elements in promoters.

S9 Table [doc]
Conserved -acting elements in promoters.

S10 Table [doc]
Conserved -acting elements in promoters.

S11 Table [doc]
analysis of miRNAs for .

S12 Table [doc]
analysis of miRNAs for .

S13 Table [doc]
analysis of miRNAs for .

S14 Table [doc]
Native and relative expression values of .