Genome-wide identification of the MADS-box gene family in Avena sativa and its role in photoperiod-insensitive oat

Identification of MADS genes in Avena sativa L

To identify the potential MADS-box genes in oat, the genome sequences were obtained from the GrainGenes database (https://wheat.pw.usda.gov/GG3/graingenes_downloads/oat-ot3098-pepsico). The MADS-box protein sequences of Arabidopsis (https://www.arabidopsis.org) and rice (http://rapdb.dna.affrc.go.jp/download/irgsp1.html) were downloaded. Hidden Markov Model-based searches (http://hmmer.janelia.org/), built from these known MADS-box protein sequences, were used to search the MADS-box proteins in oat and as queries to search against the oat protein databases with the BLASTP program with an e-value of 1 × e⁻¹⁰ as the threshold. The software pfamscan and the PFAM A database were used to annotate the domains of the candidate AsMADS sequences, and the sequence containing the PF00319 domain was determined to be the final MADS-box sequence in oat. The information of sequence length, molecular weight, isoelectric points, and the instability index was obtained from the Expsay website (http://web.expasy.org/protparam/) (Panu et al., 2012). The subcellular localization was predicted using the Softberry website (http://www.softberry.com/).

Sequence alignment and phylogenetic analysis

The MADS-box protein sequences of Arabidopsis (https://www.arabidopsis.org) and rice (http://rapdb.dna.affrc.go.jp/download/irgsp1.html) were downloaded and contained 108 and 75 sequences, respectively. Sequence alignment and evolution analysis were carried out on the MADS protein sequences identified in oat, Arabidopsis and rice. Mafft (https://mafft.cbrc.jp/alignment/software/) was used for protein multiple sequence alignment. DNAMAN software was used to align sequences and recheck data result. Both NCBI batch CDD searches and pfamscan were used to double check the conserved domains. TBtools software (https://bio.tools/tbtools) was used for visual display (Chen et al., 2020). MEGA version 7 was used for neighbor-joining (NJ) tree construction based on the alignment of MADS proteins (Kumar et al., 2018). AsMADS genes were classified according to their phylogenetic relationships with the corresponding Arabidopsis (Parenicová et al., 2003) and rice (Arora et al., 2007) MADS genes.

Gene structure, conserved motif analysis, Cis-element and protein structure analysis

Gene structure prediction was performed in GSDS (http://gsds.cbi.pku.edu.cn/). The conserved motifs were analyzed by MEME software (http://meme-suite.org/) (Bailey et al., 2009).The upstream regulatory regions (2 kb from the translation start site) of AsMADS genes were obtained from the whole genome sequence of Avena sativa by TBtools. The cis-elements on the promoter were predicted by the website http://plantregmap.gao-lab.org/. The positions of the binding sites in the physical map of the gene promoters are marked and displayed. Protein secondary structure was investigated by the Prabi website (https://npsa-prabi.ibcp.fr/cgi-bin/npsa_automat.pl?page=npsa_sopma.html), and the 3D structure was predicted on the SWISS-MODEL website (https://swissmodel.expasy.org/interactive).

Chromosome localization and collinearity analysis

According to the position of AsMADS genes on the chromosome, MegaGene2Chrom (http://mg2c.iask.in/mg2c_v2.0/) was used to draw the chromosome physical location map. MCScanX software was used to perform family gene collinearity analysis based on chromosome length information and gene structural annotation information. An intraspecific collinearity map of oat and an interspecific collinearity map of oat, rice and Arabidopsis were constructed.

Transcriptional profile analysis

For AsMADS gene expression analysis, RNA-seq data of Avena sativa L. materials with different photoperiod sensitivities under short-day conditions were obtained from previous laboratory data (An et al., 2020; Becker & Günter, 2003). Short-day treatment-related transcriptome analysis was conducted in the initial stage of panicle development in that study, in which 10 MADS genes were differentially expressed in the panicles and leaves of HQ2 and MSY4. Heatmaps were generated using these data. qPCR data results of 10 AsMADS genes was normalized by z-score. Heatmap is performed in Metware cloud platform tool (https://cloud.metware.cn/#/tools/tool-form?toolId=168). Each bar represents a gene. Red indicates that the gene is up-regulated, and blue indicates that the gene is down-regulated. GO function annotation of MADS genes (https://www.geneontology.org/).

Plant materials and treatments

Similar to the materials used in the previous transcriptome study, MSY4 and HQ2 were used for specific expression analysis at different panicle developmental stages under short-days condition (12 h). The seeding and shading methods were described in previous studies (Yang et al., 2018; An et al., 2018, 2020). The young panicles, including those at the initial stage, elongation stage, spikelet differentiation stage, floret differentiation stage, pistil differentiation stage, and tetrad stage, were collected for expression analysis. Collected tissues were immediately frozen in liquid nitrogen and stored at −80 °C until RNA isolation was performed.

Total RNA extraction and qRT–PCR expression analysis

The total RNA of oat panicles was extracted using a TransZol Up Plus RNA Kit (Transgene, Beijing, China) and then reverse transcribed into cDNA by a TransScript^® One-Step gDNA Removal and cDNA Synthesis Kit. The integrity of the RNA was analyzed by agarose gel electrophoresis, and the concentration and purity of the RNA were detected by a NanoDrop 2000. The cDNA was stored at −20 °C for subsequent qRT‒PCR. Specific primers for candidate genes in the MADS family were designed for qRT‒PCR analysis using the Actin gene of oat as an internal reference gene (Yang et al., 2013). qRT‒PCR was performed with SYBR-Green on a Jena Qtower 2.2 analyzer. According to the manufacturer’s instructions, the 20 μL reaction contained 1 μL cDNA, 400 nM of each primer and 10 μL SYBR Green Mix. The amplification conditions were 95 °C for 1 min and 40 cycles of 95 °C for 10 s, 60 °C for 10 s and 72 °C for 10 s, with a melting curve over a temperature range of 60–95 °C. Each experiment used three biological replicates and three technical replicates. The relative expression levels of genes were calculated using the 2^−ΔΔCt method (Schmittgen, 2008). Gene-specific DNA primers for qPCR are listed in File S7.

Identification and characterization of the MADS-box gene family in Avena sativa L

Based on the oat genome (OT3098), a total of 16 MADS-box genes with the PF00319 domain were identified through the Hidden Markov Model, Blastp search and Pfam annotation, which were named AsMADS1 to AsMADS16. Their basic properties were systematically evaluated, including gene ID and location, protein size, number of exons and introns, molecular weight, isoelectric points, instability index, GRAVY, and subcellular location (Table 1). The protein length of the AsMADS genes ranged from 142 to 278 amino acids, with 1 to 10 exons. The protein molecular weight ranged from 16.11 kDa (AsMADS11) to 31.61 kDa (AsMADS2), and the isoelectric points varied from 5.52 (AsMADS3) to 9.28 (AsMADS16). The hydrophobicity of the 16 proteins was negative. The instability index ranged from 44.82 to 65.59. The subcellular location prediction indicated that all members were located in the nucleus. The above results indicated that AsMADS genes were unstable hydrophobic nucleoproteins.

Multiple sequence alignment and phylogenetic relationship analysis of MADS-box genes in Avena sativa L

A multiple sequence alignment analysis of the 16 AsMADS proteins was conducted using mafft. The N-terminus of the sequence was conserved (File S1). The 16 AsMADS protein sequences were analyzed for the presence of MADS-box and other domains using an NCBI batch CDD search. Our results showed that all the AsMADS proteins contained the M domain, and 13 out of 16 AsMADS proteins contained the K domain, except for AsMADS4, AsMADS9, and AsMADS11 (Fig. 1). Moreover, 14 out of the 16 members belonged to type II MADS, and the remaining two (AsMADS9 and AsMADS11) belonged to type I MADS (Fig. 1 and Files S2 and S3).

To evaluate the evolutionary relationships among the AsMADS proteins, a phylogenetic tree was constructed using a bootstrapped neighbor-joining (NJ) method with the amino acid sequences of 207 MADS genes (16 from oat, 108 from Arabidopsis, and 83 from rice). The phylogenetic analysis indicated that the 16 AsMADS members were grouped into nine subfamilies: SEP/AGL2, AGL6, SQUA, AGL12, AG, BS, SVP, Mα, and MIKC* (Fig. 2). No members were grouped into the TM3/SOC1, AGL17, AGL15, FLC or Mγ subfamilies. The results showed that the SVP subfamily contained 3 AsMADS members, the SEP/AGL2, SQUA, AG, BS and Mα subfamilies contained two members, and the remaining three subfamilies contained only one gene (Fig. 2 and File S3). In short, MADS genes in oats consisted of type I genes (AsMADS9 and AsMADS11), MIKC* (AsMADS4), and MIKC^c-type genes (13 other AsMADS members).

Structural analysis and conserved motif composition of the AsMADS gene family

Gene structure prediction was performed in GSDS website (http://gsds.gao-lab.org/), which showed that 15 AsMADS members contained multiple introns and exons, with the number of exons ranging from five to 10 and the number of introns ranging from 4 to 9; AsMADS9 only contained one exon (Table 1 and Fig. 3C). AsMADS3, AsMADS8, and AsMADS13, which were grouped in the SVP subfamily, had a similar intron/exon pattern, but AsMADS9 and AsMADS11 had radically different structures, even though they were grouped in the same subfamily.

The secondary structure of AsMADS proteins comprised an alpha helix, extended strand, beta turn, and random coil. The AsMADS proteins had a large proportion of alpha helix amino acids (>47%), followed by random coils (File S4). Except for AsMADS11, 15 of the 16 AsMADS proteins had similar 3D structures (Fig. 4), indicating similar functions. Compared with other genes, it is clear that the 3D predicted structure of AsMADS11 contains less α-helices and more β-sheets. This makes the 3D predicted structure of AsMADS11 look simpler. However, the secondary structural pattern of AsMADS11 was similar to that of the other members (File S4).

The full-length protein sequences of 16 AsMADS were used to investigate the conserved motifs. In total, 15 different motifs were identified, and motifs 1, 2, 3 and 5 were widely distributed, indicating that they might be more conserved. Interestingly, motif 11 only existed in the BS subfamily (AsMADS7 and AsMADS15), which may be related to its evolutionary history (Fig. 3B and File S3).

Chromosome localization and collinearity analysis of AsMADS genes

Chromosomal distribution of the AsMADS gene family was visualized using TBtools and the oat genome annotation information. The 16 AsMADS genes were located on chromosomes 1D, 2D, 3D, 4D, 6D, 7D, 1A, 3A, 6A, 7A, 2C and 6C (Table 1). Briefly, chromosomes 1D, 2C, 3A, and 4D harbored 2 AsMADS genes, whereas the other chromosomes possessed only one AsMADS gene. Additionally, most members were located on the ends of the chromosome, except for AsMADS8 and AsMADS10 (Fig. 5A).

The intraspecific and interspecific collinearity of MADS-box genes was visualized using TBTools. There was no collinearity among the 16 AsMADS genes on the oat chromosomes. The interspecific collinearity analysis among Arabidopsis, rice and oat showed that there were 16 pairs of collinear genes between oat and rice, and one gene was syntenic among oat (AsMADS2), rice (Os01g10504.1) and Arabidopsis (AT4G18960), which was annotated as AG (Fig. 5B). The results showed that there was a high degree of conservation and consistency in the linear relationship of MADS-box gene evolution between oat and rice.

Prediction of cis-elements in the promoters of AsMADS genes

The 2 kb upstream region of the AsMADS genes was extracted by TBtools, and submitted to the online website plantregmap to search for cis-elements. The top 12 putative cis-elements in the AsMADS promoters are shown (Fig. 6). In addition to typical promoter elements, such as TATA boxes and CAAT boxes, G-box, circadian, MYB, MYC, STRE, ABRE, as-1, CGTCA-motif and TGCAG-motif elements were predicted. The promoters of five AsMADS members contained circadian elements, namely, AsMADS1, AsMADS7, AsMADS8, AsMADS11 and AsMADS15. Fourteen of the 16 members contained the light-responsive element G-box, except AsMADS9 and AsMADS15, but AsMADS9 contained other light-responsive elements, such as the I-box and sp1 (File S5). These results indicated that the promoters of all AsMADS members contained light-responsive cis-acting elements.

GO enrichment and transcriptome expression analysis of AsMADS genes

For the GO classification, the 16 AsMADS genes were categorized into three main categories: biological processes, cellular components and molecular functions (Fig. 7). AsMADS genes categorized as flower development (GO:0009908), floral meristem determinacy (GO:0010582), floral whorl structural organization (GO:0048459), floral organ formation (GO:0048449), and rhythmic process may be related to the photoperiod insensitivity of oat (File S6).

Based on previous transcriptome data from oat, a heatmap of the AsMADS gene expression levels at the initial stage of oat panicle differentiation in HQ2 and MSY4 under short days was drawn. Ten of 16 AsMADS members were differentially expressed in the transcriptome data. No matter in panicles or leaves, AsMADS2, AsMADS11, and AsMADS16 were upregulated in MSY4, and AsMADS3 and AsMADS12 were downregulated in MSY4 (Fig. 8).

Expression analysis of AsMADS genes in HQ2 and MSY4 under short-day conditions

In order to explore the expression specificity, the expression of AsMADS gene among materials with different photoperiod sensitivities was carried out. In HQ2, the panicle can only develop to the branch differentiation stage under short days, whereas in MSY4, it can complete the whole panicle differentiation process (An et al., 2018). Gene expression data for qPCR are listed in File S8. The expression levels of MSY4 and HQ2 were compared in the first three panicle differentiation stages. Compared with HQ2, AsMADS3, AsMADS8, AsMADS11, AsMADS13, and AsMADS16 were upregulated at these three differentiation stages in MSY4, while AsMADS12 was downregulated. AsMADS6 and AsMADS9 were downregulated in the initial stage and elongation stage and then upregulated in the branch differentiation stage in MSY4 than in HQ2, while AsMADS2 was downregulated in the initial stage and showed increased expression in the following two stages. The AsMADS15 expression level remained basically unchanged throughout the panicle differentiation process between MSY4 and HQ2. With the development of the panicle in HQ2, the expression levels of the AsMADS3, AsMADS8, AsMADS12 and AsMADS13 genes gradually decreased (Fig. 9).

In order to study the expression profile of AsMADS genes in MSY4 under short-day conditions during the entire panicle differentiation stages, real-time fluorescence quantitative analysis was carried out. With the development of the panicle in MSY4, the expression levels of the AsMADS9 and AsMADS11 genes gradually decreased. The expression levels of AsMADS15 and AsMADS16 were basically unchanged. The expression levels of the other six members fluctuated at different stages of panicle differentiation (Fig. 10).