Genome-Wide Association Study Using a Multiparent Advanced Generation Intercross (MAGIC) Population Identified QTLs and Candidate Genes to Predict Shoot and Grain Zinc Contents in Rice

Zinc (Zn) is an essential trace element for the growth and development of both humans and plants. Increasing the accumulation of Zn in rice grains is important for the world’s nutrition and health. In this study, we used a multiparent advanced generation intercross (MAGIC) population constructed using four parental lines and genotyped using a 55 K rice SNP array to identify QTLs related to Zn2+ concentrations in shoots at the seedling stage and grains at the mature stage. Five QTLs were detected as being associated with shoot Zn2+ concentration at the seedling stage, which explained 3.7–5.7% of the phenotypic variation. Six QTLs were detected as associated with grain Zn2+ concentration at the mature stage, which explained 5.5–8.9% of the phenotypic variation. Among the QTLs, qSZn2-1/qGZn2 and qSZn3/qGZn3 were identified as being associated with both the shoot and grain contents. Based on gene annotation and literature information, 16 candidate genes were chosen in the regions of qSZn1, qSZn2-1/qGZn2, qSZn3/qGZn3, qGZn7, and qGZn8. Analysis of candidate genes through qRT-PCR, complementation assay using the yeast Zn-uptake-deficient double-mutant ZHY3, and sequencing of the four parental lines suggested that LOC_Os02g06010 may play an important role in Zn2+ accumulation in indica rice.


Introduction
Zinc (Zn) is an essential micronutrient for plant growth and development [1]. Zn is also a component of numerous enzymes [2,3]. For instance, alcohol dehydrogenase (ADH), carbonic anhydrase (CA), copper/zinc superoxide dismutase (CSD), and zinc finger domain-containing proteins (ZFPs) are related to Zn in plants [1]. Furthermore, Zn is an essential component of the zinc finger motif found in many DNA-binding transcription regulators (TRs) [4]. Therefore, Zn is necessary for many biochemical processes to function properly, and adequate Zn supply is critical to crop productivity [5]. However, Zn deficiency is the most widely occurring soil micronutrient deficiency worldwide and adversely affects crop production on millions of hectares of arable land, especially in alkaline soil [6].
Zn deficiency in the edible parts of crops is one of the main reasons for Zn deficiency in humans, since the primary source of Zn intake is through food and vegetables. Zn deficiency is related to health problems such as growth retardation, impaired immune function,

Measurement of Zinc Concentration
The shoot and grain samples were dried at 65 • C for 3 days. The dried samples were crushed, wet-digested in concentrated HNO 3 at 120 • C for 30 min, and further digested with HClO 4 at 180 • C until the samples became transparent. The samples were then diluted with ultrapure water. The zinc concentrations were determined by ICP-MS (inductively coupled plasma mass spectrometry).

SNP Genotyping and Association Analysis
Meng et al. [37] genotyped the MAGIC population with a 55 K SNP array. A three-step filtering strategy was used to select high-quality SNPs for QTL mapping. First, monomorphic markers among the four parents were removed. Second, markers with missing values higher than 10% were deleted. Finally, markers with a minor allele frequency of less than 3% were filtered out. After filtering, 22,160 high-quality markers were kept to be used for analysis. A mixed linear model (MLM) implemented in TASSEL version 5.2.3 [41] was used to analyze the associations between SNP markers and traits. A value of p < 10 −2.5 was used as the threshold to declare the significance of marker-trait associations. R 2 was used to evaluate the percentage of phenotypic variance explained.

RNA Extraction and Real-Time PCR
To examine the expression response of the candidate genes to Zn 2+ deficiency, root samples were taken at the three-leaf stage from seedlings of rice variety cv. Nipponbare or the four parental lines grown in 1/4 strength IRRI (International Rice Research Institute) solution with 0.77 µM Zn 2+ or without Zn 2+ . A randomized complete block design with three replicates and a plot size of 24 seedlings were used as the experimental layout. Samples were taken from all 24 seedlings in each plot and mixed for RNA extraction. The composition of the full-strength IRRI solution was as follows: 1 To investigate the expression pattern of LOC_Os02g06010 in different organs at different growth stages, the 3 week old seedlings of cv. Nipponbare precultured hydroponically were transplanted to the paddy field in the Shenzhen Experimental Farm of the Agricultural Genomics Institute in Shenzhen, Chinese Academy of Agricultural Sciences . Tissue samples collected included roots, basal stem, leaf sheath, and leaf blade at the vegetative stage,  and roots, basal stem, lower leaf sheath, lower leaf blade, flag leaf sheath, flag leaf blade,  node I, node II, inter node II, peduncle, rachis, spikelet, husk, and seed at the reproductive stages. A single plant was regarded as a biological replicate and three biological replicates were used.
The total RNA was extracted by Trizol (Vazyme Biotech Co. Ltd., Nanjing, China). The total RNA was then reverse-transcribed using the HiScript Q RT SuperMix for qPCR kit (Vazyme Biotech Co. Ltd., Nanjing, China). The AceQ Universal SYBR qPCR Master Mix kit (Vazyme Biotech Co. Ltd., Nanjing, China) was used for quantitative analysis [42]. The primers used for qRT-PCR are given in Supplementary Table S2.

Expression of Candidate Genes in Yeast
The Zn 2+ translocation ability of each candidate gene protein was examined via a yeast complementation assay. The yeast strain used was the Zn-uptake-deficient double-mutant ZHY3 (MATα ade6 can1 his3 leu2 trp1 ura3 zrt1::LEU2 zrt2::HIS3). The open reading frames of all candidate genes were amplified from the full-length cDNA of rice cv. Nipponbare, and the primer sequences used are shown in Supplementary Table S3. Each candidate gene was ligated into a pYES2 vector with the correct direction.
Using a yeast transformation kit (Coolaber Technology Co. Ltd., Beijing, China), the empty vector pYES2, positive control OsZIP5, or candidate gene vector was introduced into the ZHY3 yeast cells. The transformed yeast was selected on a SD medium without uracil (SD-U). Positive clones were cultured in the SD-U liquid medium until the early logarithmic phase, concentrated, and washed three times with sterile ultrapure water. The yeast cell suspensions of ZHY3 transformed with empty vector pYES2, positive control OsZIP5, or candidate genes were serially diluted (1:10), and 6 µL of the cell suspension was spotted onto solid media containing 1, 3, or 100 µmol/L ZnSO 4 . The growth phenotypes were evaluated by culturing the plates at 30 • C for three days.
The OD 600 of the overnight yeast cells in the SD-U liquid medium was adjusted to 1.0 with sterile distilled water. Subsequently, 20 µL cell suspensions were added to 20 mL liquid SD-U media containing 3 or 100 µmol/L ZnSO 4 in each bottle. The OD 600 was determined at the indicated time.

Sequence Analysis of LOC_Os02g06010 in Four Parental Lines and 30 Other Rice Varieties
DNA samples of the four parental lines were extracted from the seedling samples by CTAB. The DNA samples were then used as templates to amplify the full-length genomic sequence and promoter of LOC_Os02g06010 via the KOD-FX polymerase (Toyobo, Japan) using specific primers (Supplementary Table S4). The PCR products were sequenced by Sangon Biotech Co., Ltd. (Shanghai, China). LOC_Os02g06010 sequences of the four parental lines were aligned and analyzed using DNAMAN.
The accessions were genotyped using the Illumina HiSeq 2000 (PE150) (50X) by Berry Genomics Corporation; the average sequencing depth of each accession genome was 50× [40]. The sequence of LOC_Os02g06010 in different rice varieties was analyzed according to the sequencing data.

Distribution of Zn 2+ Concentration in Shoots and Grains of MAGIC Population
We evaluated the Zn 2+ accumulation in the 215 lines of the MAGIC population and its four parental lines. The zinc concentrations in shoots of the vegetative growth period (3 weeks after being transferred) and brown rice at the mature stage were determined, respectively. The Zn 2+ concentrations in shoots and brown rice were highest in the parental line B, and the lowest in the parental line D ( Figure 1A,B). The MAGIC population exhibited significant phenotypic variation in Zn 2+ concentrations in shoots and brown rice. The distributions of the two traits were approximately normal ( Figure 1A,B). respectively. The Zn concentrations in shoots and brown rice were highest in the paren-tal line B, and the lowest in the parental line D ( Figure 1A,B). The MAGIC population exhibited significant phenotypic variation in Zn 2+ concentrations in shoots and brown rice. The distributions of the two traits were approximately normal ( Figure 1A,B).
The correlation between Zn concentration in the shoots of seedlings and concentration in the mature brown rice was positive and moderate (correlation coefficient was 0.39; Figure 2).
We evaluated the Zn 2+ accumulation in the 215 lines of the MAGIC population and its four parental lines. The zinc concentrations in shoots of the vegetative growth period (3 weeks after being transferred) and brown rice at the mature stage were determined, respectively. The Zn 2+ concentrations in shoots and brown rice were highest in the parental line B, and the lowest in the parental line D ( Figure 1A,B). The MAGIC population exhibited significant phenotypic variation in Zn 2+ concentrations in shoots and brown rice. The distributions of the two traits were approximately normal ( Figure 1A,B).
The correlation between Zn concentration in the shoots of seedlings and concentration in the mature brown rice was positive and moderate (correlation coefficient was 0.39; Figure 2).

Identification of Candidate Genes Related to Zn 2+ Accumulation in Shoots and Grains
The QTLs qSZn1, qSZn2-1/qGZn2, qSZn3/qGZn3, qGZn7, and qGZn8 were relatively more relevant for Zn 2+ concentration in shoots, concentration in brown rice, or both (Table 1, Figure 3). Based on annotation (http://rice.plantbiology.msu.edu/index.shtml) and literature information, 16 genes were chosen as the candidate genes responsible for Zn 2+ accumulation ( Table 2). Among them, 15 encode proteins with transmembrane structure (Supplementary Figure S1). Six genes, namely LOC_Os02g06010 (qSZn2-1/qGZn2), LOC_Os07g37510 (qGZn7), LOC_Os08g41370, LOC_Os08g42020, LOC_Os08g42150, and LOC_Os08g42170 (qGZn8) were found to be induced by Zn 2+ deficiency (Figure 4).   To determine the Zn 2+ transport activity of the expressed proteins of six candidate genes, expression levels of the candidate genes in the yeast mutant strain ZHY3 were compared with the empty vector pYES2 and positive control OsZIP5. In terms of growth, there were no significant differences between the ZHY3 transformed by the six candidate genes, OsZIP5 and pYES2 under the control medium and the solid medium containing 100 µM Zn 2+ (Supplementary Figure S2). There was no significant difference between the LOC_Os07g37510-, LOC_Os08g41370-, LOC_Os08g42020-, LOC_Os08g42150-, LOC_Os08g42170-, and pYES2-transformed ZHY3 on the solid medium containing 1 µM or 3 µM Zn 2+ either (Supplementary Figure S2), suggesting that these proteins may not have Zn 2+ transport activity. The growth of LOC_Os02g06010-transformed ZHY3 was significantly higher than that of pYES2-transformed ZHY3 in the solid medium containing 1 µM or 3 µM Zn 2+ after 3 days (Supplementary Figure S2; Figure 5A), suggesting that LOC_Os02g06010 may have Zn 2+ transport activity. Experiment using a liquid medium with different Zn 2+ concentrations confirmed that the growth of ZHY3 under 3 µM Zn 2+ supply was increased by LOC_Os02g06010 ( Figure 5A,B).

Expression Pattern and Sequence Analysis of LOC_Os02g06010 in Parental Lines
The difference of LOC_Os02g06010 in the four parental lines was amplified by PCR and sequenced. The coding region sequence of LOC_Os02g06010 showed no difference between the four parental lines ( Figure 6A). The sequence in the 2000 bp region of the upstream of the ATG and the 3'-UTR sequence after the TGA of LOC_Os02g06010 were identical between A, B, and C ( Figure 6A). However, D had four single-nucleotide mutations in the first intron, three single-nucleotide mutations in the second intron, and one single-nucleotide mutation in the 3'-UTR ( Figure 6A).

Expression Pattern and Sequence Analysis of LOC_Os02g06010 in Parental Lines
The difference of LOC_Os02g06010 in the four parental lines was amplified by PCR and sequenced. The coding region sequence of LOC_Os02g06010 showed no difference between the four parental lines ( Figure 6A). The sequence in the 2000 bp region of the upstream of the ATG and the 3'-UTR sequence after the TGA of LOC_Os02g06010 were identical between A, B, and C ( Figure 6A). However, D had four single-nucleotide mutations in the first intron, three single-nucleotide mutations in the second intron, and one single-nucleotide mutation in the 3'-UTR ( Figure 6A).
The sequence of LOC_Os02g06010 in different rice varieties was analyzed according to the sequencing data, and SNP: −951 bp and SNP: −830 bp in the intron region of the upstream of the ATG were found in the accessions. The two SNPs of most accessions were the same as those of A, B, and C (ABC-TYPE), and only seven accessions were the same as D (D-TYPE). We selected 30 accessions of ABC-TYPE and seven accessions of D-TYPE to measure the Zn 2+ concentration in brown rice (Supplementary Table S5). It was found that the Zn 2+ concentration in brown rice of ABC-TYPE was 13.1% higher than that of D-TYPE (Supplementary Figure S3). that the Zn 2+ concentration in brown rice of ABC-TYPE was 13.1% higher than that of D-TYPE (Supplementary Figure S3).
The expression pattern of LOC_Os02g06010 in response to Zn 2+ deficiency was further analyzed in four parental lines. Compared with the normal Zn 2+ supply, the expression of LOC_Os02g06010 was significantly increased by Zn 2+ deficiency in A, B, and C, while it was significantly reduced in D ( Figure 6B).

Expression of LOC_Os02g06010 in Different Organs at Different Growth Stages
To further analyze the biological functions of LOC_Os02g06010, its expression levels in different organs at different stages were investigated. LOC_Os02g06010 was expressed in all the measured tissues at different growth stages, and was strongly expressed at the flowering stage (Figure 7). The expression pattern of LOC_Os02g06010 in response to Zn 2+ deficiency was further analyzed in four parental lines. Compared with the normal Zn 2+ supply, the expression of LOC_Os02g06010 was significantly increased by Zn 2+ deficiency in A, B, and C, while it was significantly reduced in D ( Figure 6B).

Expression of LOC_Os02g06010 in Different Organs at Different Growth Stages
To further analyze the biological functions of LOC_Os02g06010, its expression levels in different organs at different stages were investigated. LOC_Os02g06010 was expressed in all the measured tissues at different growth stages, and was strongly expressed at the flowering stage (Figure 7).

Discussion
A set of rice MAGIC populations has been developed by IRRI to better integrate QTL discovery and breeding. Their use in QTL mapping for a range of agronomic traits has been previously reported [36][37][38][39]44]. We previously identified QTLs associated with the toxicity tolerance of rice to three essential metals (Fe, Zn, and Al) by using three of the MAGIC populations, including DC1, DC2, and eight-way populations, genotyped using

Discussion
A set of rice MAGIC populations has been developed by IRRI to better integrate QTL discovery and breeding. Their use in QTL mapping for a range of agronomic traits has been previously reported [36][37][38][39]44]. We previously identified QTLs associated with the toxicity tolerance of rice to three essential metals (Fe, Zn, and Al) by using three of the MAGIC populations, including DC1, DC2, and eight-way populations, genotyped using a 55 K SNP array [37]. In this study, we found that the four parents of the DC1 population displayed substantial differences in Zn 2+ concentrations in shoots at the seedling stage and in grains at the mature stage ( Figure 1A,B). GWAS was conducted using a mixed linear model for the two traits.
Five QTLs for the concentration of Zn 2+ in shoots at the seedling stage, and six QTLs for the concentration of Zn 2+ in brown rice at the mature stage were identified (Table 1). These QTLs did not co-locate with those QTLs for Zn 2+ concentration in leaves or shoots at the mature stage reported by Norton et al. [35] and Yang et al. [29]. This surprising inconsistency may be attributed to the obvious differences in genetic populations, growth conditions, and growth stages for sample collection. However, QTLs qSZn2-2 and qGZn5-1 identified in this study are physically very close to the known zinc-affecting genes ( Table 1), suggesting that our study was of good quality. Functional analysis of some YSL members has shown that they are not only involved in Fe uptake but also in the transport and homeostasis of other transition metals such as Zn 2+ [45]. The expression level of OsYSL14 in flag leaves exhibited significant correlations with Fe and/or Zn concentrations in the seeds [43]. OsZIP7 expression in Arabidopsis resulted in a 25% increase in shoot Zn concentrations compared to nontransformed plants [25]. OsZIP7 functions in xylem loading in roots and intervascular transfer in nodes to deliver Zn to the grain in rice [24]. Two QTLs for the two traits were co-located; qSZn2-1 was co-located with qGZn2 and qSZn3 was co-located with qGZn3. This could explain the moderate and positive correlation observed between the two traits ( Figure 2).
Screening of 16 candidate genes for five important QTLs, chosen based on annotation and literature information for response to Zn 2+ deficiency, identified six responsive genes ( Figure 4). Among them, only LOC_Os02g06010 was found to have Zn 2+ transport activity, using the transformation of the yeast mutant strain ZHY3. In rice, many genes such as OsZIP1, OsZIP3, OsZIP5, and OsZIP9 have been proven to have Zn 2+ transport activity in yeast [23,27,46]. LOC_Os02g06010 codes for an integral membrane protein, with a domain of unknown function domain (DUF) and transmembrane structure (Supplementary Figure S1).
Further functional analysis of LOC_Os02g06010 in four parental lines showed that the expression of LOC_Os02g06010 was significantly induced in A, B, and C but inhibited in D under zinc-deficient conditions ( Figure 6B). The different expression patterns of LOC_Os02g06010 in the parental lines may have caused the difference in zinc accumulation among the four parental lines. To find the reasons for the different expression levels of the four parental lines, we analyzed the promoter sequence and genome sequence of LOC_Os02g06010. The coding region of LOC_Os02g06010 showed no difference between the four parental lines ( Figure 6A). In the 3'-UTR and the introns, parental lines A, B, and C were identical, while the parental line D showed one and seven single-nucleotide mutations, respectively. ( Figure 6A). The UTR region does not encode amino acids, but plays an important role in the regulation of gene expression. In post-transcriptional regulation, the noncoding regions of mRNA (including 5'-UTR and 3'-UTR, especially 3'-UTR) determine the translation speed, stability, subcellular localization, degradation, and other features of mRNA [47][48][49]. The mutation of the 3'-UTR in the rice OsCKX2-2 gene can affect CKX enzyme activity to regulate rice development and grain size [50]. Intron retention, as a component of regulated gene expression programs, could increase transcription levels by affecting the rate of transcription, nuclear export, and transcript stability [51,52]. An intron was converted into one that increased mRNA accumulation 24-fold and reporter enzyme activity 40-fold relative to the intronless control by introducing 11 copies of the more active TTNGATYTG motif in Arabidopsis [53]. An intron-derived motif strongly increased gene expression from transcribed sequences via a splicing independent mechanism in Arabidopsis thaliana [54]. An intron-containing plasmid increased the level of GUS enzyme activity 10to 40-fold and 80-to 90-fold compared with the intronless plasmid, pBI221, in transgenic rice protoplasts and transgenic rice tissues, respectively [55]. The intron elevated GUS gene expression mainly at the mRNA accumulation level, but also stimulated enhancement at the translational level in rice [56]. The 5' UTR intron of the rice rubi3 gene enhanced GUS reporter gene activity in transgenic lines about 29-fold [57]. Wu et al. [58] found that the addition of the maize Ubi intron 1 significantly enhanced the OsMT2b promoter activity in rice embryos, suggesting that the natural variations observed in the introns or UTRs may cause differences in expression. SNP: −951 bp and SNP: −830 bp were found in the sequenced varieties, and the Zn 2+ concentration in brown rice of ABC-TYPE was higher than that of D-TYPE (Supplementary Table S5; Supplementary Figure S3). It is speculated that SNP: −951 bp and SNP: −830 bp in the first intron region of the upstream of the ATG may be the key sites for regulating LOC_Os02g06010 expression. LOC_Os02g06010 was expressed in all the measured tissues of different growth stages of rice, and was more strongly expressed at the flowering stage ( Figure 7). Taken together, the present results suggest that LOC_Os02g06010 may play an important role in Zn 2+ accumulation in indica rice. Further characterization of LOC_Os02g06010 is needed to elucidate its functional significance in Zn 2+ uptake, translocation, and accumulation in rice.

Conclusions
In this study, we used a MAGIC population to locate QTLs related to Zn 2+ concentration in shoots at the seedling stage and in grains at the mature stage. Five QTLs (qSZn1, qSZn2-1, qSZn2-2, qSZn3, and qSZn8) were detected to be relevant for Zn 2+ concentration in shoots at the seedling stage, and explained 3.7-5.7% of the phenotypic variation. Six QTLs (qGZn2, qGZn3, qGZn5-1, qGZn5-2, qGZn7, and qGZn8) were detected to be relevant for Zn 2+ concentration in grains at the mature stage, and explained 5.5-8.9% of the phenotypic variation. Sixteen candidate genes were predicted, and LOC_Os02g06010 was selected as our candidate gene through qRT-PCR and yeast heterologous functional complementation verification test. The expression of LOC_Os02g06010 was significantly induced in A, B, and C, but inhibited in D by Zn 2+ deficiency. The coding region of LOC_Os02g06010 showed no difference between the four parental lines. In the 3'-UTR and introns, parental lines A, B, and C were identical, while the parental line D showed one and seven singlenucleotide mutations, respectively. SNP: −951 bp and SNP: −830 bp in the first intron may be the key sites for regulating LOC_Os02g06010 expression. LOC_Os02g06010 was strongly expressed at the flowering stage of rice. These results suggest that LOC_Os02g06010 may play important roles in Zn 2+ accumulation in indica rice.
Supplementary Materials: The following are available online at https://www.mdpi.com/2077-0 472/11/1/70/s1, Table S1: Description of the four parental lines used for developing the MAGIC DC1 population, Table S2: Primers used for qRT-PCR, Table S3: Primers used to amplify the open reading frames of six candidate genes and OsZIP5, Table S4: Primers used to amplify the full length of LOC_Os02g06010 and its promoter, Table S5: Detailed information of 37 rice accessions, Figure S1: Predicted topological models of proteins encoded by 16 candidate genes, Figure S2: Complementation of the yeast mutant strain ZHY3 by OsZIP5 and six candidate genes, Figure S3: Haplotype analyses of LOC_Os02g06010.