Genome-wide analysis of maize OSCA family members and their involvement in drought stress

Identification of OSCA protein-coding genes in the maize genome

Conserved OSCA domain DUF221 (Pfam accession: 02714) from the Pfam database (Finn et al., 2014) was used to build the Hidden Markov Model-based searches (http://hmmer.janelia.org/) and scanned against the maize genome (genome assembly: AGPv3) and the sorghum genome (genome assembly: V3). All the retrieved sequences were curated using the NCBI Conserved Domain Database (www.ncbi.nlm.nih.gov/Structure/cdd/wrpsb.cgi) (Marchlerbauer et al., 2011) to determine whether the proteins harbor the DUF221 domain or not.

Genic structure and phylogenetic relationships analysis

The protein sequences of the identified maize and sorghum OSCAs were downloaded from Phytozome v10.0. The protein sequence of previously identified OSCA protein-coding genes in the Arabidopsis and rice genome (Yuan et al., 2014; Li et al., 2015) were also downloaded form Phytozome v10.0. To show the exon/intron structure, the coding sequence of each OSCA gene were aligned to its corresponding genomic sequence and then a schematic representation was generated using GSDS 2.0 (http://gsds.cbi.pku.edu.cn) (Hu et al., 2015). To construct a phylogenetic tree of the identified OSCA proteins in rice, Arabidopsis, and sorghum, and maize genome, ClustalW (Thompson, Higgins & Gibson, 1994) was used to align multiple protein sequences of OSCA. The phylogenetic tree was constructed using this alignment output based on a neighbor-joining method in MEGA7 (Kumar, Stecher & Tamura, 2016) using the following parameters: Poisson correction, pairwise deletion, uniform rates, and bootstrap (1,000 replicates).

Prediction of transmembrane region

TMHMM Server V2 was used to predict the transmembrane region (TMs) of ZmOSCAs, and the prediction of TMs was manually edited according to the TMs of OSCA1.1, OSCA1.2, and OsOSCA1.2 (Jojoa-Cruz et al., 2018; Liu, Wang & Sun, 2018; Maity et al., 2018; Murthy et al., 2018; Zhang et al., 2018).

Gene duplication analysis

Chromosomal locations of ZmOSCAs were collected from the maize genome (assemble version: AGPv3) in Phytozome v10.0. The duplicated segmental blocks and gene pairs were analyzed on the Plant Genome Duplication Database (available online: http://chibba.agtec.uga.edu) with a display distance of 100 kilobases.

Expression profile analysis

Seedlings growth conditions and drought treatments of B73 (maize inbred line) were conducted according to Wang et al. (2016). Hydraulically cultured three-leaf stage seedlings were put on a plate and subjected to dehydration (40–60% relative humidity and 28 °C). As for the effect of drought implied on adult leaf, drought stress was applied by withholding water after the eight-leaf stage (V8), with well-watered plants (soil water content 40%) as control. The middle section of flag leaf which came from three replicates, were collected at the 12-leaf stage (V12), the 14-leaf stage (V14), the 16-leaf stage (V16), and the silking stage (R1) for both drought stressed and well-watered plants as control. Leaf samples from at least three replicates were frozen by liquid nitrogen and then stored at −80 °C before RNA isolation.

Raw RNA was extracted from leaf samples using TRI Reagent (Invitrogen, Carlsbad, CA, USA) according to the product manuals. The relative expression of ZmOSCAs was quantified using quantitative real time-PCR (qRT-PCR). qRT-PCR was tested in 96-well plates using an ABI7500 Real-Time PCR Systems (Applied Biosystems, Foster City, CA, USA). The PCR reaction system consists of one μl cDNA, 200 nM primers, and five μl SYBR Premix Ex Taq II (Takara, Beijing, China), and the reaction volume was 10 μl. The PCR reaction was performed with the following conditions: 10 min at 94 °C, 40 cycles of 15 s at 94 °C, and 30 s at 60 °C. The internal control was the expression of ZmUbi-2 (UniProtKB/TrEMBL; ACC: Q42415). The quantification method used was 2^−ΔCT, and the variation in the expression was derived from three biological replicates.

Association analysis

Association analysis for ZmOSCAs was conducted by using a mapping population containing 367 maize inbred lines and corresponding drought resistance phenotypes reported previously (Wang et al., 2016). The mapping population contains 556k single nucleotide polymorphism (SNP) markers, and the minor allele frequency of each marker was greater than or equal to 0.05. All identified ZmOSCAs harbored 168 SNPs in the coding region and both the 5′-, and 3′-untranslated region. Three statistical models including the general linear model (GLM) model, adjusting the first two principal components (PC₂), and the mixed linear model (MLM) model (incorporating PC₂ and a Kinship matrix) were selected to identify the SNPs significantly associated with drought resistance by using the TASSEL4.0 program (Yu et al., 2006; Bradbury et al., 2007).

Identification of the OSCA family members in maize

We carried out a systematic genome-wide screen of putative OSCA genes in maize. Initially, a Hidden Markov Model search was performed against the maize genome (genome version: AGPv3.0) utilizing the DUF221 domain (Pfam accession: 02714). Ultimately, 12 genes were found to be ZmOSCAs in the maize genome and named according to the Arabidopsis orthologs (Table 1). The physical location of each ZmOSCA in the genome was identified according to physical coordination provided by MaizeGDB (https://www.maizegdb.org). A total of 12 ZmOSCA genes were distributed unevenly on all the 10 chromosomes except chromosomes 2, 4, 7, and 10, without any clustering. Chromosomes 1 and 3 possessed as many as three ZmOSCA genes (the largest number of OSCA genes on a chromosome), and chromosomes 5 and 8 equally contained two genes. In contrast, chromosomes 6 and 9 each only harbored one ZmOSCA. The exon number of ZmOSCAs varied from 1 to 11. Approximately 75% (8/12) of the ZmOSCA genes contained more than eight exons, and only 17% (2/12) of the genes had less than two exons.

Phylogenetic relationship analysis of OSCA genes

In order to elucidate the phylogenetic relationships among OSCAs, a neighbor-joining tree of ZmOSCA proteins and the corresponding orthologs from rice, sorghum, and Arabidopsis was built, and the tree was based on the alignment of full-length OSCA proteins (Fig. 1; Table S1). As shown on the phylogenetic tree, 49 OSCA genes can be classified into four main classes: clades 1, 2, 3, and 4 (Fig. 1). All Arabidopsis and rice OSCAs fell in the same class or clade as previously reported, which is in agreement with previous work (Yuan et al., 2014). Interestingly, proteins derived from dicotyledonous Arabidopsis clustered separately from those of three monocot plants. Also, we found that some proteins from rice, sorghum, and maize displayed pairwise correspondence (blue box in Fig. 1), not only indicating that these genes are phylogenetically conserved among these species, but demonstrating that maize and sorghum possessed a closer phylogenetic relationship compared to maize and rice, conforming to the perspective that sorghum is a closer relative to maize than rice.

Throughout the phylogenetic tree, clades 1, 2, 3, and 4 contain 21, 20, 4, and 4 OSCA proteins, respectively (Fig. 1). Amongst the 49 OSCAs, there were 5, 5, 1, and 1 ZmOSCAs in clades 1, 2, 3, and 4, respectively. However, for the rice and sorghum genomes, there were 4, 5, 1, and 1 OSCAs in clades 1, 2, 3, and 4, respectively (Table 2). Obviously, compared with rice and sorghum, maize possessed one more OSCA in clade 1 (Fig. 1). Meanwhile, we detected the nodes that lead to maize-, sorghum-, and rice-specific clades (red circles in Fig. 1). These nodes denote the divergence point between maize, sorghum, and rice, and therefore reveal the most recent common ancestral (MRCA) genes before the split. After MRCA analysis, there were 11 clades (indicated by red circles in Fig. 1), and nearly all the clades were constituted by one gene from the rice, maize, and sorghum genomes. However, there was just one clade that contained two genes from the maize genome (indicated by the red arrow in Fig. 1). The results showed that between the maize, sorghum, and rice genomes, the OSCA gene family shared a similar evolutionary history, whilst the maize genome possessed a duplication event leading to a gain in one gene (indicated by red arrow in Fig. 1) after the maize, sorghum, and rice genomes split.

Synteny of OSCAs in maize, sorghum, and rice genomes

Gene collinearity comparative genomics indicates homologous gene function and phylogenetic relationships amongst several species and even the genome organization of extinct ancestral species. Thus, we investigated the colinearity of OSCA genes in rice, sorghum, and maize genomes. Initially, gene colinearity data were retrieved from the Plant Genome Duplication Database using ZmOSCAs as anchors, then the chromosomal segment containing multiple homologous genes across species was defined as a genomic syntenic block. After this analysis, 10 ZmOSCAs were found to have collinearity, and the syntenic members or collinear genes in rice and sorghum were shown in Fig. 2. As a result, eight chromosomal segments containing OSCAs were identified as being evolutionally conserved between rice, maize, and sorghum; these include ZmOSCA1.2, ZmOSCA1.3, ZmOSCA1.4, ZmOSCA2.2, ZmOSCA2.3, ZmOSCA2.4, ZmOSCA2.5, and ZmOSCA4.1 (Fig. 2; Table S2). In addition, genes within the eight syntenic block amongst the maize, sorghum, and rice genomes were identical to the phylogenetic analysis (Fig. 1), confirming the accuracy of our analysis. Syntenic blocks of ZmOSCA1.1b and ZmOSCA2.1b were observed between maize and sorghum genomes respectively, but no syntenic blocks of these genes were found between maize and rice genomes. Meanwhile, no collinear segments of ZmOSCA1.1a and ZmOSCA3.1 were found in the genome of maize, rice, and sorghum. The results suggest that most of the ZmOSCAs existed before the divergence of species, but some ZmOSCAs may have originated from duplication of the maize genome.

The conserved domain of ZmOSCAs

Conserved Domain Database was used to identify protein domains contained in ZmOSCAs. It was found that most ZmOSCAs harbored three domains except that ZmOSCA1.1a, ZmOSCA3.1, and ZmOSCA4.1 contained two domains (Fig. 3). It is noteworthy that DUF221 is located at the C-terminal of all ZmOSCAs, indicating that this domain is not only necessary but also relatively conservative compared with other domains in the gene family. It was reported that there were 11 TMs in the OSCA protein sequences (Jojoa-Cruz et al., 2018; Liu, Wang & Sun, 2018; Maity et al., 2018; Murthy et al., 2018; Zhang et al., 2018). We found that half of the 12 ZmOSCAs had exactly 11 transmembrane regions (Fig. 3). A multiple sequence alignment was used to show the presence of TMs in DUF221 of ZmOSCAs (Fig. 4). It was found that protein domain DUF221 contained a different number of TMs, while protein domain pfam13967 contained a fixed number of TMs. However, no TMs were detected in the protein domain pfam14703 in ZmOSCAs. Different ZmOSCA members contained two to seven TMs in the DUF221 region. Most of the ZmOSCAs contained at least eight TMs with two exceptions. ZmOSCA1.1a and ZmOSCA3.1 had fewer TMs than others, suggesting that a deletion event occurred for ZmOSCA1.1a and ZmOSCA3.1.

Expression profiles of ZmOSCA genes under drought stress

We studied the expression profiles of the 12 ZmOSCA genes in 15 different tissues using the reported transcriptomic data (Sekhon et al., 2011). As shown in Fig. 5A, expression patterns of ZmOSCA genes varied greatly. Expression values of ZmOSCA1.1a, ZmOSCA1.1b, ZmOSCA1.2, and ZmOSCA4.1 were relatively higher than the rest among the tissues tested. Previous studies have shown that OSCAs had key roles in different aspects of plant development, especially in stress responses (Li et al., 2015). In order to obtain further insight into the roles of ZmOSCAs in drought tolerance, their expression profiles were monitored by quantitative real-time PCR analyses using 3-week-old leaves of maize seedlings treated with drought stress for 5 or 24 h. As shown in Fig. 5B, in the genotype of B73, more than half of expression level of the ZmOSCAs were regulated by drought stress, including six ZmOSCAs that were up-regulated significantly and one down-regulated ZmOSCA; the others were hardly responsive to drought stress. In addition, genes showing up- or down-regulation were found unfavorably in each of the four major clades. The results showed that proteins from different clades might have overlapping and/or antagonistic functions in the regulation of drought responses in plants. Amongst the six up-regulated ZmOSCA genes, in comparison to normal growth conditions, the expression of ZmOSCA2.4 was induced more than sixfold, that being the largest folds change in relative expression levels (Fig. 5B). Five ZmOSCAs (ZmOSCA1.1b, ZmOSCA1.2, ZmOSCA1.4, ZmOSCA2.1, and ZMOSCA4.1) were up-regulated over twofold in response to drought (Fig. 5B). Additionally, a slightly up-regulated expression was observed in some ZmOSCAs following drought treatment for five hours or 24 h; these include ZmOSCA1.1a, ZmOSCA2.2, ZmOSCA2.4, ZmOSCA2.5, and ZmOSCA3.1. In contrast to the other ZmOSCAs, the expression level of ZmOSCA2.3 decreased about threefold after drought treatment (Fig. 5B). Relative expression of ZmOSCAs were also detected in adult leaves at four growth stages, V12, V14, V16, and R1. Drought stress for the adult leaves was applied by withholding water at the V8 stage, while the corresponding control was well-watered (soil water content 40%) plants. As shown in Fig. 6, five ZmOSCAs were up-regulated at least one stage, and one ZmOSCA was down-regulated at two stages under drought. Across the four stages, the relative expression levels of ZmOSCA2.1, ZmOSCA2.2, and ZmOSCA4.1 were up-regulated consistently. Taken together, the results clearly show the functional divergence of ZmOSCAs responding to drought stress in maize seedlings and adult leaves. Collectively, the data demonstrated that different ZmOSCAs showed variable expression patterns under drought stress, while the expression value of ZmOSCA4.1 was highly expressed in 15 tissues and up-regulated in both seedling and adult leaf by drought stress.

Association analysis of genetic variations in ZmOSCAs with drought tolerance

In order to further inquire whether the genetic variation in ZmOSCAs was associated with drought tolerance, a family-based association analysis was performed for these genes. At the seedling stage, the drought tolerance of the mapping population was assessed by evaluating seedling survival rate under severe drought stress. Based on previously reported methods and data (Liu et al., 2013; Wang et al., 2016), genetic polymorphism was characterized as the presence and number of SNP markers in each of the 12 ZmOSCAs. All the 12 identified ZmOSCAs were identified as polymorphic, with 14 SNPs on average for each identified gene (Table 3). As a consequence, ZmOSCA4.1 was identified as the second most polymorphic gene, with 34 SNPs in this mapping population. Subsequently, three statistical models were used to find out significant associations between genotype and phenotype. In brief, a GLM, with the first PC₂, and MLM incorporating both PC₂ and a kinship matrix correcting the effect of cryptic relatedness were used in the associations. GLM method was used to conduct the single-marker analysis. PC₂, via the first PC₂ on SNP data, was utilized to adjust spurious associations resulted from population structure. MLM was regarded as an effective method for controlling false positives in association analysis (Yu & Buckler, 2006; Yu et al., 2006). Subsequently, candidate gene association analysis identified significant associations between genetic variation of ZmOSCA4.1 with drought tolerance under different models with a P-value ≤ 0.01 (Table 3). The results appeared to show that ZmOSCA4.1 might be an important candidate gene for maize drought tolerance.