玉米主要株型性状与产量的全基因组关联分析

白明兴 陈奋奇 陆晏天 丁永福 姬祥卓 彭云玲,2,*

(1 甘肃农业大学农学院/甘肃省干旱生境作物学重点实验室, 甘肃 兰州 730070；2 甘肃省作物遗传改良与种质创新重点实验室, 甘肃 兰州 730070)

玉米(ZeamaysL.)是全球最重要的粮食、饲料和经济能源作物，其生产安全对保障国家粮食增产、促进畜牧业发展以及改善人民生活质量发挥着极其重要的作用。株高和穗位高等性状是构成玉米株型的主要因子，适当调整玉米株型可以提高玉米单产的形成[1-2]。因此，通过研究玉米株型相关性状不仅能有效改善玉米植株的空间分布形态，促进玉米增量，也可以为选育理想株型和分子标记辅助选择(marker-assisted selection, MAS)奠定基础。近年来，国内外诸多学者已经利用数量性状定位(quantitative trait loci, QTL)分析和全基因组关联分析(genome-wide association study, GWAS)等手段对玉米株型和产量相关性状进行了研究，吴建伟等[3]以A188×91黄15构建了1套F2∶3群体，在水分胁迫下检测到了10个QTL与株高和穗位高相关，并在8号染色体上检测到了1个同时调控株高和穗位高的主效QTL。郑德波等[4]以K22为共同父本，CI7和Dan340分别为母本开发了2套F2群体,分别构建了含237和218个家系的两套F2∶3群体，并采用1 536个单核苷酸多态性标记(single nucleotide polymorphism, SNP)的基因芯片在不同环境下对2个群体的株高、穗位高等性状进行了QTL定位分析，共定位到48个QTL，其中多数QTL的作用方式以加性与部分显性效应为主。Pan等[5]采用10个重组自交系(recombinant inbred lines, RIL)群构建了一套含有1 887个品系的ROAM(random-open-parent association mapping)群，并采用GWAS、SLM(single linkage mapping)和JLM(joint linkage mapping)等3种不同的方法对玉米的10个株型性状进行了定位，在3种不同分析方法下分别定位到38、83、86个与株高相关的QTL，找到1个控制株高的主效QTL(qPH3)，并将其定位至600 kb(165.73～166.37 Mb)。另外，赵小强等[6]也通过构建两套F2群体对玉米的单穗重、穗轴重、百粒重、穗长等穗部性状进行QTL定位，共检测到49个QTLs。上述结果表明不同研究所采用的群体及QTL定位方法不同，导致所检测到相关农艺性状的QTL信息也存在差异，同时也反映出了这些玉米性状遗传基础的复杂性，导致研究者还不能清晰剖析其遗传机理。

关联分析也称作连锁不平衡作图，是利用分子标记对供试群体基因进行检测或扫描，进而获得丰富的基因位点及其等位基因信息的技术，可分为全基因组关联分析和候选基因关联分析两类[7]。近年来，随着高通量测序技术的快速发展，许多植物的基因组测序已逐渐完成，同时也开发了大量的分子标记，为关联分析的应用提供了良好的条件。并且随着SNP标记的不断开发，全基因组关联分析已成为剖析一些作物复杂农艺性状分子机理的主要手段[8-9]，如在玉米[10-11]、水稻(OryzasativaL.)[12-13]、小麦(TriticumaestivumL.)[14-15]和油菜(BrassicacampestricL.)[16-17]等主要农作物上均已被广泛应用。本研究以204份玉米自交系为关联群体，并借助分布于玉米20条染色体上的558 529个SNP，对玉米株型与产量相关的5个性状进行GWAS分析，以期为玉米株型的改良和玉米理想株型MAS育种提供参考。

1 材料与方法

1.1 试验材料

从中国农业科学院提供的368份材料中选择在武威适应性较好的204份玉米自交系为试材。

1.2 田间数据采集

于2018年在甘肃省武威试验站(37.97°N, 102.63°E)进行试验。材料采用完全随机区组设计，3次重复，且供试材料为双行区种植，行长4 m，株距25 cm，行距40 cm。在玉米大喇叭口期从供试自交系中选择长势整体一致的植株3株，测定其株高(plant height，PH)、穗位高(ear height, EH)、总叶片数(leaf number，LN)、穗上叶数(leaf number above ear，LNAE)，并在2019年3月进行室内考种，采集供试材料的单穗重(ear weight, EW)，计算均值。

1.3 表型数据处理

采用Excel 2016对各性状进行数据统计分析，并采用SPSS 22.0软件对5个相关性状进行正态性检验及主成分分析(principal component analysis，PCA)，对各性状进行相关分析，热图层次聚类分析采用 MeV4.9.0处理[18]，距离度量标准采用Pearson相关系数距离。

1.4 全基因组关联分析

利用分布于玉米全基因组的最小等位基因频率≥0.05的558 529个SNP位点进行全基因组关联分析，基因型数据可在网站(http://www.Maizego.org/Resources.html)获得(n=368)。在TASSEL 5.0软件中，采用3种不同的模型(即只控制群体结构的Q模型，只控制亲缘关系的K模型以及同时控制群体结构和亲缘关系的Q+K模型)对玉米的株高、穗位高、总叶片数、穗上叶数及单穗重等5个性状进行关联分析，并根据其Quantile-Quantile(QQ)散点图选择最佳模型对各性状进行GWAS分析。QQ图利用R语言的“CMplot”包绘制。

1.5 相关位点候选基因分析

前人已利用约56万个SNP对该群体的连锁不平衡(linkage disequilibrium，LD)衰减程度进行了评价，发现该关联群体的LD衰减距离为50 kb(R2=0.1)[19]。因此，本研究将以此衰减距离作为该群体的LD衰减距离，基于显著关联的SNP标记在玉米B73基因组序列RefGen_v3上的物理位置上、下游共100 kb范围内，通过NCBI(https://www.ncbi.nlm.nih.gov/)和MaizeGDB(https://www.maizegdb.org/)搜索与各性状相关的所有候选基因，并根据其基因功能注释选择最佳基因作为候选基因。

2 结果与分析

2.1 玉米自交系相关性状间的表型鉴定和对产量的影响

由供试自交系各性状的表型鉴定结果可知(表1)，玉米各性状存在丰富的遗传变异，其变异系数在9.00%～50.00%之间。通过对各性状进行描述统计发现，各性状的峰度和偏度都基本介于-1.0～1.0 之间，呈典型的数量遗传特性，因此可以对各个性状进行关联分析。此外，各性状表型值的变幅波动范围较大，这可能是试验点环境因素复杂和材料自身差异共同作用的结果。

2.2 供试自交系各性状间的相关性分析和聚类分析

由表2可知，供试自交系的株高(PH)与穗位高(EH)、总叶片数(LV)、穗上叶数(LNAE)及单穗重(EW)性状间均呈显著正相关；穗位高与总叶片数、穗上叶数及单穗重呈显著正相关；总叶片数与穗上叶数呈显著正相关，与单穗重无显著相关性；穗上叶数与单穗重无显著相关性。且由各性状的聚类分析(图1)可知，5个性状共聚成两类，其中总叶片数和穗上叶数2个性状聚成第一类，株高、穗位高和单穗重3个性状聚成第二类，热图层次聚类分析也进一步表明玉米自交系的这5个性状彼此间相互关联，相互协同，其相互作用共同构成了不同玉米自交系间的差异。

表1 不同玉米自交系相关性状的表型鉴定Table 1 Phenotypic identification of maize inbred lines

表2 不同玉米自交系各性状的相关性分析Table 2 Correlation analysis of traits in different inbred lines of maize

2.3 主成分分析

通过对各性状进行主成分分析可知(表3、表4)，在关联群体的各性状间筛选出特征值大于1的主成分2个，且总体方差累计贡献率达到71.667%，表明这2个主成分可以代表供试群体5个相关性状的绝大多数信息。其中PC1的特征值是2.361，贡献率为47.223%，除单穗重外，其他性状载荷值都较高；PC2的特征值是1.222, 贡献率达到24.444%，其单穗重的载荷值较高。

表3 各性状特征值及贡献率Table 3 Eigenvalue and contribution rate of each traits

2.4 玉米相关性状的全基因组关联分析

2.4.1 阈值确定 利用R语言进行Bonferroni阈值的确定，本研究采用P≤1×10-5作为此次全基因组关联分析的显著阈值。

表4 各表型因子载荷矩阵Table 4 Matrix of phenotypic factor loading

2.4.2 模型选择 本研究采用3种模型(Q、K和Q+K)对各性状进行GWAS分析，由各性状的QQ图(图2)可知，PH、EH、LN、LNAE及EW等5个性状中，K模型和Q模型总体对假阳性控制较差，而Q+K模型虽然对部分性状假阳性的控制过于严格，但总体对假阳性控制效果较好，综合考虑，本研究采用Q+K模型的结果进行后续分析。

2.4.3 关联分析 由5个相关性状的GWAS分析结果可知，以P≤1.0×10-5为显著阙值，在关联群体的5个性状间共检测到13个显著SNP(附表1)，分布于玉米的第3、第5、第6、第7号染色体上。其中，在与株高有关的9个SNP位点中，8个位于5号染色体上，1个在6号染色体上；与穗位高相关的3个SNP位点分别位于第3、第6和第7号染色体上以及检测到1个与单穗重有关的SNP位点，位于第6号染色上；但在总叶片数和穗上叶数2个性状间未关联到显著的SNP位点。

2.5 玉米候选基因

本研究在50 kb(R2=0.1)LD衰减距离范围内对13个与玉米相关性状显著SNP位点进行相关基因的搜集，在所有显著位点内共搜集到39个相关基因，有12个基因可以搜索到相关的注释信息(附表1)。其中，在与株高相关的9个显著SNP位点间共搜索到12个基因，有注释的基因5个；在与穗位高相关的3个显著SNP位点间共搜索到25个基因，有注释的基因4个；在与单穗重相关的1个显著SNP位点间共搜索到2个基因，未搜到相关注释。

3 讨论

高产、稳产、优质一直是玉米育种的主要目标，而株高、穗位高、总叶数、穗上叶片数等性状与玉米单产的提高有着密不可分的关系。本研究以204份玉米自交系为关联群体，通过对其株高、穗位高、总叶片数、穗上叶片数及单穗重等5个性状进行相关分析和聚类分析发现，各性状间存在显著的正相关性，说明各性状紧密联系，相互协同，共同影响玉米产量。随着分子标记技术的不断发展，玉米株型及产量相关性状的研究报道也逐渐增多，何坤辉等[20]利用玉米自交系许178×K12衍生的150个F7重组自交系在6种不同环境下定位到10个株高QTL和8个穗位高QTL，并找到了3个主效QTL，可以在3种及以上环境中被定位到。赵小强等[21]通过对两套群体的7个农艺性状进行QTLs检测，共定位到了56个QTLs位点,并找到了12个在各性状间稳定表达的QTL位点。李凯等[22]以360份玉米自交系为试材，利用全基因组关联分析在4种环境下对玉米株高、穗位高进行了定位分析，共检测到6个与株高显著相关的SNP位点，18个与穗位高显著相关的SNP位点。马娟等[23]以一套BC3F4群体为材料，利用56K SNP芯片对玉米产量相关性状进行QTL定位分析，并通过Meta-QTL分析对多年来玉米产量相关性状的QTL位点进行了整合。以上研究表明，不同群体在不同环境下对同一性状检测到的QTL位点虽有差异，但均可为今后对株型性状和产量性状相关基因的挖掘和基因克隆提供良好的基础。本研究通过对3种模型进行比较分析发现，Q+K模型能较好地控制各性状的假阳性，因此采用Q+K模型对玉米自交系的5个相关性状进行全基因组关联分析，在株高、穗位高和单穗重等3个性状间共检测到13个显著的SNP位点。而在总叶片数和穗上叶片数间未检测到显著的SNP位点，这可能是由于单一位点对表型的贡献率不高，造成位点效应太弱，不能达到显著水平。

图1 PH、EH、LN、LNAE、EW 5个性状的热图层次聚类分析Fig.1 Heat map with hierarchical clustering analysis of PH, EH, LN, LNAE and EW

图2 5个相关性状的3种GWAS模型的QQ图Fig.2 Quantile-quantile (QQ) plots resulting from GWAS results using three methods of five maize related traits

本研究也得到一些与前人QTL定位结果相似的SNP位点，Yan等[24]以杂交种综3×87-1的F2衍生出的266份F2∶3家系为材料，在2种环境下对不同发育时期玉米的株高进行了QTL定位，并在Bin5.03处bnlg1879～bcd207标记间检测到株高相关QTL。汤继华等[25]利用一套“永久F2”群体在Bin5.03处的bnlg1879～umc1935以及Bin5.03～ Bin5.04多处标记间检测到与株高相关的QTL。许诚等[26]也利用一套单片段代换系群体对玉米5个株型相关性状进行了QTL定位，指出umc1587-mmc0081-umc2072(Bin5.02～ Bin5.06)标记间存在同时控制株高、穗位高、穗上叶片数和总叶片数的QTL。另外，李贤唐等[27]利用一套重组自交系对玉米株型性状进行了QTL分析，在Bin6.04处检测到1个与株高相关的QTL(mmc0523～umc2006)，同时也在Bin6.03～Bin6.04处检测到1个与穗位高相关的QTL(umc1887～umc1796)。李浩川等[28]也采用一套DH系对玉米株高和穗位高进行遗传分析，并在Bin3.01处定位到与穗位高相关的QTL。而本研究也在Bin5.03和Bin6.04处检测到了9个与株高显著相关的SNP，在Bin3.01、Bin6.03及Bin7.01处分别定位到1个与穗位高相关的SNP，结合前人研究进一步表明，Bin 5.03、Bin6.04及Bin 3.01附近是调控玉米株高、穗位高等重要性状的富集点，且对玉米生长发育发挥着重要的积极效应。

此外，本研究对玉米各性状显著关联的SNP位点进行了可能的候选基因预测，在株高、穗位高、总叶片数、穗上叶数及单穗重等5个性状中，仅株高、穗位高及单穗重间检测到了显著的SNP位点，并对每一性状的显著位点搜索其候选基因及其功能注释(附表1)。在与株高相关的12个基因中，发现基因GRMZM2G048363(Bin5.03)编码无机磷酸盐转运蛋白(inorganic phosphate cotransporter)，而磷元素作为植物生长发育所必需的三大元素之一，其对植物生长发育起着至关重要的作用，同时也参与植物体内三磷酸腺苷和磷脂等重要分子的构成，在能量转移、代谢调节和蛋白活性调控等多方面起作用[29]，因此推测该基因可能是通过对磷素的吸收来达到植物生长所需的最佳比例，从而促进植物茎秆的生长。另外，本研究还发现在与穗位高相关的25个基因中，基因GRMZM2G042752(Bin7.01)编码F-box蛋白At5g07610，Dieterle[30]和吴丹等[31]指出F-box蛋白在植物生长发育中发挥着重要作用，尤其在调节植物光形态建成方面有着极其重要的正向作用。秘彩莉等[32]也指出F-box蛋白在植物生长发育中至关重要，同时也参与植物激素(乙烯、生长素、赤霉素、茉莉酸)的信号传导及花器官发育等多种生物学过程。而生长素和赤霉素在促进植物生长方面有着积极的正向作用，且赤霉素最突出的生理效应是促进茎的伸长。因此，基因GRMZM2G042752可能通过调节植物体内生长素和赤霉素等多种激素的共同作用，调节植物茎的伸长生长。

4 结论

本研究采用3种模型(Q、K和Q+K)对玉米株高、穗位高、总叶片数、穗上叶片、单穗重等5个株型相关性状进行全基因组关联分析，发现Q+K模型对假阳性的控制效果最好。基于Q+K模型的 GWAS 结果表明，以P≤1×10-5为显著阙值在群体间共检测到13个与玉米株高、穗位高及单穗重显著相关的SNP，并搜素到39个相关基因，其中12个基因在NCBI和MaizeGDB中可以搜到相关的功能注释。此外，通过对株高、穗位高及单穗重等3个性状所有检测到的SNP位点进行候选基因挖掘及其功能注释，推测出控制株高和穗位高最佳的候选基因。

附表1 玉米相关性状候选基因及功能注释Schedule 1 Candidate genes and functional annotation of corn-ear trifoliate traits