A maize association panel consisting of 155 inbred lines (Yang et al., 2010) was used to detect associations between the nucleotide polymorphisms of ZmVDE1 and carotenoid content as well as composition in maize kernels. Eighty-nine inbred lines from the 155 lines in the association panel and a set of 44 teosinte accessions from the CIMMYT gene bank (Supplementary Table S1) were used to resequence ZmVDE1 in this study. The maize panel was planted in one-row plots with two replications in a randomized complete block design in Beijing, China during the summer of 2005, 2006, and 2007 and the winter of 2007 in Hainan, China. Each plot was 4 m long. Rows were spaced 0.67 m apart, and plants were grown at a planting density of 45,000 plants/ha; seeds were produced via self-pollination of each plant in the plot (Fu Z.Y. et al., 2013). The maize kernel carotenoids were measured after harvesting. The best linear unbiased predictors (BLUPs) for individual traits in each line were calculated using SAS 8.02 (SAS Institute, 1999) in the association mapping populations. BLUPs for each line across environments were used for the overall analysis (Fu Z.Y. et al., 2013).

According to the zeaxanthin, antheraxanthin and violaxanthin inter conversion pathway (EMP, Enzymes and Metabolic Pathways database¹), the protein sequence of the VDE homologous gene in Arabidopsis thaliana, AVDE1, was retrieved from TAIR². The sequence homologous to AVDE1 in maize was obtained via BLAST using AVDE1 protein sequence against the maize genomic sequence in the Maize Genetics Genomics Database (maizeGDB³). The related information for ZmVDE1, including nucleotide sequence, transcript and gene structure, was obtained from the Gramene database⁴. In addition, DNA motif was predicted using the PLACE⁵ and PlantCARE⁶ databases.

To further investigate the genetic relationship of the VDE genes among different species, reference sequences were retrieved from the non-redundant protein sequence database using protein tool BLAST at NCBI⁷. And then the alignment file was constructed using ClustalW in MEGA 5.0. Phylogenetic tree was constructed using neighbor-joining criterion in MEGA 5.0 with 1,000 bootstrap tests for every node (Tamura et al., 2011). The uniform rates was applied for rates among sites. Gaps or missing data was treated as complete deletion.

Multiple self-pollinated ears from each line for each replication were combined for high-performance liquid chromatography (HPLC) analysis. Measured metabolites included lutein, zeaxanthin, β-cryptoxanthin, α-carotene and β-carotene. Carotenoids were extracted as described in Chander et al. (2008) and quantified by standard regression against external standards (Kurilich and Juvik, 1999). External standard curves were constructed with eight serial dilutions and with three repeats for each dilution (R²≥ 0.99). The five carotenoids were separated on a reversed-phase C30 column (YMC Carotenoid, CT99S05-2546WT C30, 4.6 nm × 25 cm, 5 μm; Waters) at 30°C at 1.8 ml/min for the mobile phase, 75:20:5 (v/v/v) acetonitrile/methanol/dichloromethane, by scanning at 450 nm with a reference wave of 360 nm and were identified by the retention time of each standard. The peak times for lutein, zeaxanthin, β-cryptoxanthin, α-carotene and β-carotene were 8.37 min, 10.71, 19.94, 27.61, and 37.61 min, respectively. All phenotypic data were generated with ChemStation software (Agilent Technologies) (Xu et al., 2012).

Two primers that spanned most of the gene were used in conjunction with two additional sets of primers to sequence ZmVDE1 in 89 maize lines and 44 teosinte entries (Supplementary Table S2). A 30-μl aliquot of the resulting PCR product from each line was sequenced directly using an ABI3730 sequencer. The sequences were assembled using ContigExpress in Vector NTI Advance 10 (Invitrogen), aligned using MUSCLE (Edgar, 2004) and manually corrected with BioEdit (Hall, 1999). Additionally, to determine the expression level of ZmVDE1, 557,955 polymorphic sites with a minor allele frequency (MAF) ≥0.05 and missing rate <25% were generated by RNA-Sequencing maize kernels collected 15 days after pollination (DAP) from 368 maize inbred lines (Fu J.J. et al., 2013). Using the expression levels of ZmVDE1, we looked for expression QTL and analyzed their correlations with carotenoid content. Subsequently, >120,000 single nucleotide polymorphisms (SNPs) distributed across chromosome 2 in the 368 maize lines were used for testing whether there were selective sweeps around ZmVDE1.

The polymorphic sites including SNPs and insertions or deletions (InDels) of ZmVDE1 with MAF ≥0.05 in 89 maize lines from the 155-line association panel were extracted using TASSEL 2.1 (Bradbury et al., 2007). Associations between the polymorphic sites and five carotenoid-related traits were carried out using a mixed linear model (Yu et al., 2006), which incorporated population structure and kinship (Yang et al., 2010) in TASSEL 2.1.

Nucleotide diversity (π) and Tajima’s D statistic were calculated for ZmVDE1 and chromosome 2 using DNaSP version 5.0 (Librado and Rozas, 2009). The scaled per-nucleotide recombination parameter Rn (Hudson, 1987) is the length-weighted mean of Rn across sweep region estimated from 368 maize lines.

The 89 maize inbred lines and 44 teosinte entries for the two subspecies Z. mays ssp. mexicana and Z. mays ssp. parviglumis (Supplementary Table S1) were used to construct a minimum spanning tree for ZmVDE1. Arlequin version 3.5 (Excoffier and Lischer, 2010) was used to calculate the minimum spanning tree nucleotide polymorphisms. Arlequin’s distance matrix output was used in Hapstar-0.6 (Teacher and Griffiths, 2011) to draw the minimum spanning tree.

A BLAST search with AVDE1 protein sequences from A. thaliana (AT1G08550) against the maize genomic sequences in the maizeGDB database resulted in three hits for homologous genes, GRMZM2G027219 on chromosome 2, GRMZM2G701673 on chromosome 10 and GRMZM2G408706 on chromosome 8. The maize homologous gene GRMZM2G027219 had the highest identity with AVDE1 (E-value = 4.734 × 10^-134, 81.95% amino acid identity) and was referred to as ZmVDE1. ZmVDE1 from maize inbred line B73 has five exons, and the full-length cDNA sequence is 1,761 base pairs (bp), encoding 446 amino acids. The length of the 5′untranslated region (UTR) and 3′UTR are 148 and 272 bp, respectively (Supplementary Figure S1). Motif Scan analysis (Sigrist et al., 2010) showed that the VDE domain (amino acid 163–360) is present in the mature protein of ZmVDE1. Using the PLACE and PlantCARE databases, we also found the following putative DNA motifs, which are probably involved in light responsiveness, in the regulatory region and introns of ZmVDE1 (Supplementary Figure S1): I-box, G-box, GT-1 motif, TCT-motif, and MNF1-motif.

Functional domain analysis of this genomic fragment in InterProScan and EBI QuickGO confirmed that the VDE domain in Z. mays is highly consistent with those of the graminaceous plants, such as Triticum aestivum, Sorghum bicolor, and Oryza sativa (Figure 1A). To further investigate the genetic relationship of the VDE proteins among different species, phylogenetic analysis of these homologs was conducted in graminaceous plants, herbaceous plants and oil crops (Figure 1B). ZmVDE1 is grouped with all the graminaceous plant VDEs and is more closely related to the VDE of Sorghum bicolor than to other graminaceous crops. VDEs from other plants included oil crops and herbaceous plants in the another cluster.

Based on 53 samples ranging in developmental age from fertilization to maturity and consisting of the embryo, endosperm and whole-seed tissues from a previous study (Chen et al., 2014), we detected the expression pattern of ZmVDE1 during maize kernel development. ZmVDE1 was mainly expressed in the late stages of the embryo, endosperm and whole seed. Moreover, ZmVDE1 at most time points was expressed at a higher level in the endosperm than in the embryo samples (Supplementary Figure S2). Additionally, Stelpflug et al. (2015) had enhanced their previously published microarray-based gene atlas of maize inbred line B73 (Sekhon et al., 2011, 2012) to include 79 distinct replicated samples that had been analyzed using RNA-sequencing. The expression pattern of ZmVDE1 in the endosperm and whole-seed tissues was similar to the findings of Chen et al. (2014). Therefore, we inferred that ZmVDE1 is active in the late stages of seed development, in accordance with its function as an important synthetase involved in nutrient accumulation in maize kernels.

To investigate the nucleotide polymorphisms of ZmVDE1 in maize, the full-length sequence including the promoter, 5′UTR, exons, 3′UTR and most of the introns was obtained from 89 inbred lines. Our re-sequencing results identified five nucleotide polymorphisms with MAF ≥0.05, two InDels and three SNPs, in a 3,724 bp genomic region (Figure 2A). On average, one in every 745 bp was polymorphic. Additionally, the distribution of the nucleotide diversity was not even across ZmVDE1. Among these polymorphic sites, only one SNP (SNP283) was located in the coding region, which affected residue 65 (which encodes a p.Val65Ile alteration) with no change in charge or hydrophobicity. The limited number of nucleotide polymorphisms within ZmVDE1 across 89 maize inbred lines illustrates that ZmVDE1 has an important role in the maize kernel with a conserved function.

Association analysis with a mixed linear model that controlled for both population structure and individual relatedness revealed that no polymorphic sites in ZmVDE1 showed a significant association with any of the analyzed carotenoid-related traits at the P = 0.05 level (Figure 2A). Previous studies have shown the importance of transcript abundance in the control of carotenoid profiles (Aluru et al., 2008; Harjes et al., 2008; Bai et al., 2009). Thus, we tested the correlation between the polymorphisms identified in the DNA sequence and the ZmVDE1 expression level to investigate whether ZmVDE1 can regulate any carotenoid-related traits at the level of expression. At P < 1.8 × 10^-6, there was no statistical correlation between DNA sequence polymorphisms of ZmVDE1 and its expression level (Figure 2B). In addition, expression of ZmVDE1 was not correlated with phenotypic variation for any measured carotenoid trait (Figure 2C, Supplementary Figure S3). Thus the limited natural variation of ZmVDE1 makes no significant contribution to carotenoid variation.

To examine the levels of genetic diversity in ZmVDE1, we sequenced the promoter, 5′UTR, exons and introns of ZmVDE1 in 44 teosinte entries consisting of two subspecies, Zea mays ssp. parviglumis and Zea mays ssp. mexicana (Supplementary Table S1). Together with the ZmVDE1 sequences from 89 maize lines, 15 polymorphic sites were identified (Supplementary Table S3). When we used a sliding window of 100 bp with a step size of 25 bp to calculate nucleotide diversity, we noted that nucleotide diversity was not equally distributed in either maize or teosinte, with exons containing less genetic variation (Figure 3A). In most parts of ZmVDE1, the average π over all maize lines was lower than that in teosinte, especially for the promoter and exon 5 regions, which had dramatically reduced diversity in maize (π_M/π_T= 0.087 and π_M/π_T= 0.115, respectively) (Figure 3B). This finding suggests that the promoter and exon 5 in ZmVDE1 might have been under selection during the domestication of maize from teosinte.

To further investigate the evidence of selection during maize domestication, Tajima’s D test was conducted for the promoter, 5′UTR and exon regions of ZmVDE1 in maize and teosinte. Tajima’s D test showed that promoter, 5′UTR and exon 4 regions of ZmVDE1 had significantly negative parameter deviations (P < 0.01) (Figure 3B), which is consistent with directional selection (Jiao et al., 2012). The results suggest that a selective sweep might have occurred in those regions, leading to reduced nucleotide diversity across the ZmVDE1 locus in maize.

To further investigate the genetic signature of positive selection in the ZmVDE1 region in more detail, we selected 120,204 high-quality SNPs on chromosome 2 with MAF ≥0.05 from a panel of 368 diverse maize inbred lines (Fu J.J. et al., 2013). Using a sliding window of 500 SNPs with a step size of 100 SNPs to calculate the nucleotide diversity and conduct Tajima’s D test, we found that the maize inbred lines had severely reduced nucleotide diversity (π = 0.123) across a 4.2 Mb region surrounding the ZmVDE1 locus relative to other segments on chromosome 2 (Figure 4A). In addition, the statistic for Tajima’s D test was notably negative for a selective sweep in the 368 inbred lines, suggesting an excess of low-frequency polymorphisms relative to the expectation that might have resulted from the selection on the selective sweep in this population (Figure 4A). The large interval (4.2 Mb) associated with this selective sweep might be related to the low SNP or gene density and the recombination rate. As expected, no significant reduction in SNP or gene density was observed in the 4.2 Mb region surrounding the ZmVDE1 locus. There are 478 SNPs in the 4.2 Mb selective sweep with one in every 8,786 bp being polymorphic, and there are 45 genes in the large interval, with one every 93.3 kb (Figure 4B). In contrast, there is one polymorphic nucleotide every 1,951 bp and one gene every 49.5 kb in other regions of chromosome 2 (Figure 4B). The nucleotide estimate for the population recombination rate indicated that the recombination rate across this sweep region (Rn = 0.0757) was 1.59-fold lower than it was at the known selection target tga1 (Rn = 0.1205). These results suggested that few recombination events can be detected in this selective sweep because of the relatively low levels of polymorphism and that the low rate of recombination has contributed to the size of the sweep.

To gain deeper insights into the evolution of ZmVDE1, its full-length sequence including the promoter, 5′UTR, exons, 3′UTR and most of the introns from 89 Chinese elite inbred lines and 44 additional teosinte entries (Supplementary Table S1) was used to construct a minimum-spanning tree. In total, we identified five haplotypes for the 89 diverse maize lines, of which Hap1 predominated, and 15 haplotypes for 44 teosinte entries (Supplementary Table S4). With these haplotypes, a minimum-spanning tree was built (Figure 3C). The 17 haplotypes formed two distinguishable clusters: a maize haplotype cluster and a teosinte haplotype cluster. The maize haplotype cluster contained almost all the maize lines and two teosinte entries, whereas the teosinte haplotype cluster was composed of teosinte entries with one additional maize line. These results revealed that VDE1 was domesticated with fewer haplotypes in maize comparing with which in teosinte. (Figures 3A,B), the nucleotide diversity was much higher in teosinte than in maize.

In maize, ZmVDE1 contains five putative circadian DNA motifs, I-box (GATA motif), G-box, GT-1 motif, TCT-motif and MNF1-motif (Supplementary Figure S1). These motifs are thought to be involved in light regulation and circadian rhythm expression (Borello et al., 1993; Piechulla et al., 1998; Xu et al., 2003; Kidokoro et al., 2009). In addition to the light-regulated motifs, the other most important motifs in introns were cis-acting elements involved in plant responses to environmental stresses, including MYB2CONSENSUSAT(YAACKG), CBFHV (RYCGAC), RAV1AAT (CAACA), CCAATBOX1 (CCAAT), WBOXATNPR1(TTGAC), GT1GMSCAM4(GAAAAA), MYB1 AT (WAACCA), MYBCORE (CNGTTR), MYCCONSENSUSAT (CANNTG) and WRKY71OS (TGAC), which indicated that the expression of ZmVDE1 might also be regulated by transcription factors related to stress tolerance. All of these motifs showed high identities with other gramineous VDEs such as those from T. aestivum and Oryza sativa, especially in the regions corresponding to the mature proteins (Gao et al., 2013).

Violaxanthin de-epoxidase has a conserved function converting violaxanthin to antheraxanthin and then to zeaxanthin in the xanthophyll cycle, which is a photo-protection mechanism conserved among photosynthetic plants and algae to protect the photosynthetic apparatus from excess light (North et al., 2005). Meanwhile, downstream of the xanthophyll cycle, violaxanthin and zeaxanthin are the intermediate products in the biosynthesis pathway of ABA, which regulates the development process, seed maturation, stomatal closure and stress responses (Finkelstein et al., 2002; Zhu, 2002; Wang and Song, 2008). Therefore, VDE and zeaxanthin epoxidase (ZEP) have been the targets of genetic manipulation to better understand the carotenoid synthetic pathway. Römera et al. (2002) found that the amount of violaxanthin was diminished dramatically, and in some cases the amount of the monoepoxy intermediate antheraxanthin was increased, as a consequence of the genetic manipulation of transformation with sense and antisense constructs encoding zeaxanthin epoxidase, even though the total carotenoid content and ABA content were both increased. Thus VDE might regulate the synthesis of ABA (Owen et al., 1992; Marin et al., 1996; Pogson et al., 1996; Audran et al., 1998). The enhanced activity of VDE under excess light could reduce the content of ABA, perhaps by acting on the xanthophyll cycle and regulating the synthesis of the ABA precursor.

Candidate gene association mapping in DNA from kernels grown in multiple environments indicated that none of the five polymorphisms in ZmVDE1 were significantly associated with carotenoid-related traits in maize kernels. To confirm this result, we also extracted 220 SNPs of ZmVDE1 in 368 maize inbred lines (Fu J.J. et al., 2013) to perform a regional association analysis for carotenoid-related traits. No significant SNPs were identified, consistent with the association analysis in the 89 maize inbred lines. In addition, ZmVDE1 expression was not correlated with any of the measured carotenoid-related traits. One possible explanation for the absence of an association may be the small size of our association panel. It is well known that association mapping studies may be adversely affected by two major interconnected factors, the dominance of a SNP with a lower MAF and a small population size (Manolio et al., 2009). Simulation results of the effect on QTL detection power with increasing population size revealed that a population with 500 lines provides an 80% probability of detecting a gene that explains 3% or more of the phenotypic variation (Yan et al., 2011). So some important genetic variations may not have been detected in our association panel with its relatively small population size. Another explanation for this absence of an association is that ZmVDE1 has an important and conserved role in plants. Additionally, ZmVDE1 was expressed throughout the leaf with highest expression level at maturing and mature zones (Supplementary Figure S4). Moreover, expression of ZmVDE1 in seedling leaves was statistically higher than that in maize kernels (Chen et al., 2014; Stelpflug et al., 2015), similar to that of AtVDE in A. thaliana (North et al., 2005), suggesting that ZmVDE1 also has an important protective role in chloroplasts (Audran et al., 2001).

Previous studies have shown that association mapping results of candidate genes cloned by homology-based cloning is quite different among maize. For example, ZmGS3 and ZmGW2 have co-orthologous relationships with the rice region containing GS3 and GW2, respectively. ZmGS3 was slightly associated with maize kernel length and one-hundred kernel weight (HKW), with the P-values ranging between 0.05 and 0.01, whereas ZmGW2 was significantly associated with kernel size (P < 10^-3) (Li et al., 2010a,b). These results for ZmVDE1, ZmGS3 and ZmGW2 imply that the power of an association study is dependent on various factors, such as the gene effect, the genetic architecture of the quantitative trait and the size of association population and so on. So, in addition to the association analysis, we still need other genetic methods to validate gene function.

Because of the role of the VDE in the xanthophyll cycle, an important photo-protection mechanism, it is easy to understand why ZmVDE1 was a target of the domestication or improvement process in maize. In our study, ZmVDE1 was less diverse at the nucleotide level in maize than in teosinte, especially in the promoter and 5′UTR regions, consistent with a significantly negative Tajima’s D parameter, which may indicate that this region has undergone a recent and strong positive selection. While, significant reductions in diversity were detected in the VDE domain region (π_M/π_T = 0.079) with negative Tajima’s D parameter (Tajima’s D = -1.576), which may be because that VDE domain variants within these regions might have experienced positive selection during maize domestication and improvement. After this selective force on ZmVDE1, the VDE domain became fixed with no polymorphic sites, as indicated by widespread Chinese maize inbred lines (Figure 2A).

Positive directional selection results in reduced variability in particular regions, which are referred to as selective sweep regions (Tian et al., 2009). According to our results, ZmVDE1 is located in a genomic region of 4.2 Mb that was affected by a selective sweep based on the 368-line maize association population. A larger selected region (5.2 Mb) located on chromosome 2 has been reported, and there was a strong signal for ZmVDE1 in this region (Hufford et al., 2012), which supports our findings. Similar evidence for a domestication-related selection sweep has also been reported in maize. A selective sweep extending 60–90 kb upstream of tb1 indicated that tb1 had a role in the domestication of maize from teosinte between 6,000 and 10,000 years ago (Clark et al., 2004). Maize PSY1 on chromosome 6, which controls kernel color, has undergone very strong selection for the yellow kernel color phenotype. Patterns of diversity in this region are consistent with the occurrence of a large selective sweep (>0.5 Mb) at PSY1 in maize breeding lines (Palaisa et al., 2004). Maize also shows severely reduced nucleotide diversity relative to teosinte across a large selective sweep (1.1 Mb) on chromosome 10 that has been the target of strong selection during maize domestication (Tian et al., 2009). In addition, a study on chromosome 3, which contains the largest selected region (2.2 Mb), using 278 temperate maize inbred lines from different stages of breeding history showed that modern breeding has introduced highly dynamic genetic changes into the maize genome (Jiao et al., 2012). At the same time, artificial selection leads to a reduction in nucleotide diversity and an increase in the proportion of rare alleles (Jiao et al., 2012). Collectively, these results imply that ZmVDE1 is a major domestication- and improvement-related gene involved in the carotenoid biosynthesis pathway. The challenge for the future is to figure out how the evolution and selection of this gene have occurred.