刘坤, 张雪海, 孙高阳, 闫鹏帅, 郭海平, 陈思远, 薛亚东, 郭战勇, 谢惠玲, 汤继华, 李卫华. 玉米株型相关性状的全基因组关联分析[J]. 中国农业科学, 2018, 51(5): 821-834 https://doi.org/10.3864/j.issn.0578-1752.2018.05.002
LIU Kun, ZHANG XueHai, SUN GaoYang, YAN PengShuai, GUO HaiPing, CHEN SiYuan, XUE YaDong, GUO ZhanYong, XIE HuiLing, TANG JiHua, LI WeiHua. Genome-Wide Association Studies of Plant Type Traits in Maize[J]. Scientia Acricultura Sinica, 2018, 51(5): 821-834 https://doi.org/10.3864/j.issn.0578-1752.2018.05.002

0 引言

【研究意义】玉米的理想株型^[1]具有较高的光能利用效率,具备制造更多干物质的生理潜力。实现生育期内源、流和库的协调发展,合理分配光合产物,可使植株间的竞争达到最小状态,从而提高产量。株型相关性状如株高、穗位高、叶片数等与植株产量、光合效率、抗倒性等密切相关,是理想株型设计育种的基础。因此,对其遗传基础进行研究,不但可以了解玉米理想株型建成的遗传机制,而且还可帮助在育种实践中有效降低植株高度与穗位高度、合理调控叶片数目及分布,从而提高作物产量。【前人研究进展】玉米株型性状是典型的数量性状,由多基因控制并受多种因素尤其是环境的影响十分显著。因此,科研工作者常利用数量性状位点（quantitative trait locus,QTL）定位的方法分析玉米株型性状的遗传机制,解析基因的遗传效应、基因之间的互作及基因与环境之间的互作关系。近年来,随着测序技术和统计方法的发展,研究者借助高密度标记分别对株高、穗位高、叶夹角、叶长、叶宽、节间数等玉米株型性状的遗传规律进行解析。玉米株高与穗位高是株型的重要构成因子,并且是玉米数量遗传学研究的模式性状,前人分别利用不同的分离群体（如重组自交系群体）和关联群体等^{[2,3,4,5,6,7,8]},对玉米株高和穗位高的遗传机制进行了大量研究,检测到多个与株高和穗位高相关的QTL,分布在玉米的10条染色体上^[9,10]。其中位于玉米第1染色体的Brachytic2所在的玉米株高效应位点qph1,可以解释20%的表型变异。TENG等^[11]克隆的株高基因ZmGA3ox2,被认为是调控玉米株高的主效QTL即qPH3.1的候选基因。YANG等^[12]对玉米株高进行了关联分析,在混合线性模型下仅找到一个位于第3染色体的SNP（Chr.3:162709488）,该SNP可解释9.61%的表型变异;PAN等^[13]利用3种不同的定位方法（SLM、JLM和GWAS）对10个重组自交系群体组成的ROAM群体的株型性状进行定位,在SLM、JLM和GWAS下分别鉴定到83、86和38个株高QTL,并精细定位了一个位于第3染色体的株高主效QTL到600 kb（165.73—166.37 Mb）。国内外许多****对玉米叶片数与产量之间的关系也进行了大量研究^[14,15,16],发现虽然多叶品种可以增加穗上部分的叶片面积,增加了受光面积指数,使得植株生物量增高,但是多叶品种并没有增加籽粒产量。而玉米品种最终产量的提高是在种植密度不变的条件下增加了玉米的叶面积指数,进而增加受光面积,加速灌浆速率^{[17,18,19,20]},其中穗上叶片数对单产有较大的正向直接效应^[21]。【本研究切入点】尽管前人利用不同的分离群体定位并克隆了一些与株高、穗位高、叶片数等株型性状相关的QTL或基因,但由于这些性状属于复杂的数量性状,人们对其遗传机理还知之甚少。【拟解决的关键问题】本研究以284份遗传多样性丰富的温带、亚热带和热带材料组成的关联群体为研究对象,借助覆盖玉米全基因组的约56万个SNP标记,对玉米多个株型相关性状进行全基因组关联分析,旨在为玉米理想株型的构建,为抗倒伏、耐密植的玉米新品种选育提供理论基础。

1 材料与方法

1.1 试验材料

所用材料为华中农业大学严建兵教授提供,是从关联群体（n=513）中选择的在郑州及鹤壁适应性较好的284份玉米自交系构成,主要来源为热带（104份）、亚热带（42份）及温带材料（138份）。

1.2 田间试验设计

关联群体于2012年分别种植在郑州河南农业大学试验站（郑州试验站）与鹤壁浚县的鹤壁农科院试验站（浚县试验站）。田间设计采用完全随机区组设计,3次重复,每个材料单行区种植,行长4 m,行距0.67 m,株距0.25 m,密度63 000株/hm²,田间管理采用常规管理。为了叶片数目调查的准确性,在小喇叭口和大喇叭口时期分别对玉米第5和第10片叶进行标记,并在玉米灌浆末期,每一行挑选10株长势均匀的植株对4个株型相关性状株高（plant height,PH）、穗位高（ear height,EH）、总叶片数（leaf number,LN）和穗上叶片数（leaf number above ear,LNAN）进行调查。

1.3 数据处理与分析

1.3.1 一般统计分析 采用Excel、SPSS19.0统计软件对2个环境下关联群体的株高（PH）、穗位高（EH）、总叶片数（LN）、穗上叶片数（LNAN）以及由此衍生的穗位高/株高（EH/PH）和穗上叶片数/总叶片数（LNAN/LN）等6个性状进行统计分析,并分别计算各性状的平均值、标准差、标准误、频率分布等,并进行相关性分析。
同时根据R语言的lme4包^[22]的混合线性模型进行最佳线性无偏预测（best linear unbiased prediction,Blup）,评估2个环境下每一个材料每一个性状的育种值,最终Blup值也用于后续的全基因组关联分析。
1.3.2 全基因组关联分析 所使用的基因型为华中农业大学严建兵教授提供的覆盖玉米全基因组的最小等位基因频率大于0.05的558 629个SNPs组成^[12],基因型数据可在网站http://www.Maizego.org/Resources.html下载。考虑到很多标记间可能存在高度连锁不平衡,因此,采用GEC软件^[23]计算了有效标记数（effective number,Ne）,并根据软件给出的建议P值（1/Ne）作为确定各性状-SNP间是否显著关联的依据。
在全基因组关联分析（GWAS）中,由于不同性状对GWAS不同模型控制假阳性和假阴性的敏感程度不同,因此统计功效是首先要考虑的因素。为了选择最优的GWAS模型,使用了3种常用的模型,即只控制群体结构的Q模型（Q）、只控制亲缘关系的K模型（K）及同时控制群体结构和亲缘关系的Q+K模型（Q+K）,对每一个性状进行了GWAS分析。通过对3种模型的QQ Plot的比较来选择最优模型并对最优模型的GWAS结果进行进一步解析,以上GWAS分析均在TASSEL3.0（http://www.maizegenetics.net/）软件中完成。先前研究表明,该群体的连锁不平衡（linkage disequilibrium,LD）衰减距离为50 kb^[24,25],利用558 629个SNP标记对所用材料的LD进行了重新评价,发现平均LD为50—100 kb（R²=0.1）,简便起见,选择50 kb作为LD衰减距离。因此定义峰值SNP上下游各50 kb的区间范围（即共100 kb）为一个位点;同时定义不同位点间,若存在重叠部分,则这些位点构成一个非重复的QTL。另外,基于B73参考基因组（B73 RefGen_v2）,从玉米数据库MaizeGDB（http://www.maizegdb.org）下载玉米基因组的基因列表搜索每一个位点内的所有可能的候选基因。根据InterProScan数据库（http://www.ebi.ac.uk/interpro/scan.html）并结合美国国立生物技术信息中心数据库（https://www.ncbi.nlm.nih.gov/）对每一个候选基因进行功能注释。最后根据每一个位点内所有基因的功能注释及其在B73的各组织的表达情况,选择一个最可能的候选基因。

2 结果

2.1 株型相关农艺性状表型数据分析

对此关联群体在不同环境下的6个株型相关性状LN、LNAN、LNAN/LN、PH、EH、EH/PH的表型进行统计分析（表1）,得知玉米株型相关性状存在广泛变异。同时,也发现浚县试验站的株高、穗位高、总叶片数及穗上叶片数均高于郑州试验站（图1）,主要是由于浚县位于更高纬度,日照长度要高于郑州,导致营养生产期延长所致,同时发现2个环境中株型相关性状的峰度和偏度的绝对值均小于1（表1）,暗示这些性状均遵从正态分布,表明这些性状均为典型的数量性状,由微效多基因控制。
Table 1
表1
表1玉米不同环境下株型相关性状描述统计分析
Table 1Descriptive statistics for plant-type traits of the association population under different environments

环境 Environment	性状 Traits	变异范围 Range	均值 Mean	标准差 SD	偏度 Ske.	峰度 Kur.
郑州Zhengzhou	LN	14.90—23.20	18.37	1.30	0.62	0.71
	LNAN	3.50—6.85	5.18	0.62	0.17	-0.01
	LNAN/LN	0.20—0.38	0.28	0.03	0.01	0.05
	PH(cm)	107.00—240.30	170.89	23.92	0.03	-0.01
	EH(cm)	27.00—112.38	63.61	15.21	0.47	0.55
	EH/PH	0.19—0.60	0.37	0.07	0.50	0.71
浚县 Xunxian	LN	11.50—24.70	20.10	1.78	-0.12	0.71
	LNAN	3.50—7.80	5.62	0.72	0.13	0.18
	LNAN/LN	0.19—0.35	0.28	0.03	-0.10	-0.51
	PH(cm)	118.50—270.40	197.47	30.02	-0.04	-0.43
	EH(cm)	35.58—140.50	81.66	17.95	0.28	0.21
	EH/PH	0.25—0.60	0.41	0.06	0.19	0.27
综合 Blup	LN	13.56—24.00	19.16	1.51	0.37	0.66
	LNAN	3.85—7.20	5.37	0.62	0.23	0.12
	LNAN/LN	0.19—0.38	0.28	0.03	-0.10	0.03
	PH(cm)	108.9—257.00	182.64	26.15	0.10	-0.05
	EH(cm)	27.00—124.05	71.79	16.81	0.36	0.40
	EH/PH	0.19—0.60	0.39	0.06	0.31	0.75

新窗口打开
显示原图|下载原图ZIP|生成PPT
图1株型相关性状在关联群体的分布情况 a：叶片数相关性状;b：株高相关性状;XX表示浚县试验站;ZZ表示郑州试验站;LN：总叶片数;LNAN：穗上叶片数;LNAN/LN表示穗上叶片数与总叶片数比值;EH：穗位高;PH：株高;EH/PH：穗位高与株高的比值。下同
-->Fig. 1Distribution of plant type traits in association mapping populations (AMP) a: Traits related to leaf number; b: Traits related to plant height related. Data from two locations, Zhengzhou (ZZ) and Xunxian (XX), for AMP are shown. LN: Total number of leaves; LNAN: Leaf number above ear; LNAN/LN: The ratio of LNAN to LN; EH: Ear height; PH: Plant height; EH/PH: The ratio of EH to PH. The same as below
-->

2.2 株型相关性状的相关性分析和方差分析

通过对2个环境的6个玉米株型性状进行相关性分析,发现除了浚县试验站穗位高与穗上叶片数及2个环境下株高与穗位高/株高的比值不存在显著相关之外。其余各性状之间均表现出显著正相关或显著负相关（表2）。2个环境中的穗上叶片数与总叶片数的比值和穗上叶片数之间、以及穗位高与株高比值和穗位高之间、株高和穗位高之间的相关系数较大,并处于极显著水平。说明各性状的发育存在着相互协同促进的作用。
Table 2
表2
表2关联群体株型性状的相关性分析
Table 2Correlation coefficients among plant type traits in the association population

性状 Trait	总叶片数 LN	穗上叶片数 LNAN	穗上叶片数/总叶片数LNAN/LN	株高 PH	穗位高 EH	穗位高/株高 EH/PH
LN		0.472**	-0.147*	0.234**	0.473**	0.482**
LNAN	0.513**		0.795**	0.356**	0.087	-0.188**
LNAN/LN	-0.205**	0.728**		0.240**	-0.222**	-0.533**
PH	0.389**	0.500**	0.260**		0.577**	0.105
EH	0.558**	0.225**	-0.187**	0.770**		0.693**
EH/PH	0.450**	-0.188**	-0.571**	0.112	0.711**

**: Correlation is significant at the 0.01 level; *: Correlation is significant at the 0.05 level. The same as below. Correlation coefficients in Xunxian are above the diagonal, in Zhengzhou are below the diagonal**：在0.01水平上显著相关;*：在0.05水平上显著相关。下同。对角线右上角表示浚县点,对角线左下角表示郑州点
新窗口打开
方差分析结果表明,关联群体所有6个株型性状的基因型效应都达到极显著水平,说明关联群体各材料间的株型性状存在显著的遗传变异（表3）。同时,环境间、基因型与环境互作的差异也达到极显著水平,说明玉米株型不同性状受环境影响较大,同时也存在着基因型与环境之间的互作关系。
Table 3
表3
表3玉米不同环境下株型相关性状的方差分析
Table 3Analysis of variance (ANOVA) for maize plant type traits in different locations

性状 Trait	变异来源 Source	离均差平方和 Type sum of squares	自由度 Df	均方 Mean square	F值 F
总叶片数 LN	环境E	231.359	1	231.359	157.332**
	基因型G	1684.621	247	6.82	4.638**
	基因型×环境G×E	461.096	247	1.867	1.269*
	误差Error	729.375	496	1.471
	总变异Total	381381.855	992
穗上叶片数 LNAN	环境E	51.138	1	51.138	1613.663**
	基因型G	370.545	247	1.5	47.339**
	基因型×环境G×E	94.819	247	0.384	12.114**
	误差Error	15.718	496	0.032
	总变异Total	29516.316	992
穗上叶片数/总叶片数LNAN/LN	环境E	0.024	1	0.024	60.086**
	基因型G	0.753	247	0.003	7.510**
	基因型×环境G×E	0.284	247	0.001	2.836**
	误差Error	0.201	496	0
	总变异Total	77.731	992
株高 PH	环境E	154293.307	1	154293.307	21896.670**
	基因型G	624125.793	247	2526.825	358.597**
	基因型×环境G×E	111389.229	247	450.969	64.000**
	误差Error	3495.028	496	7.046
	总变异Total	34834662.07	992
穗位高 EH	环境E	66081.315	1	66081.315	34206.174**
	基因型G	239554.291	247	969.855	502.034**
	基因型×环境G×E	52668.276	247	213.232	110.377**
	误差Error	958.199	496	1.932
	总变异Total	5702103.988	992
穗位高/株高 EH/PH	环境E	0.289	1	0.289	5239.943**
	基因型G	3.972	247	0.016	291.642**
	基因型×环境G×E	1.321	247	0.005	96.957**
	误差Error	0.027	496	0
	总变异Total	160.568	992

E: Environment; G: Genotype; G×E: Genotype-by-environment interactionE：环境;G：基因型;G×E：基因与环境互作
新窗口打开

2.3 株型相关性状全基因组关联分析

2.3.1 阈值确定 通过GEC软件^[23]计算了558 629个SNPs的有效标记数,结果表明,有效标记数Ne为250 345,软件给出的建议P值为3.99E-6（1/Ne）,因此,该值将作为后续GWAS分析的标记-性状间显著关联的阈值。
2.3.2 模型选择 以2个环境下6个性状的Blup值（作为一个新的环境）为例,用3种模型（Q、K和Q+K）分别进行了GWAS分析,通过不同性状3种模型的比较,从QQ图中不难看出（图2）,对这6个性状来说,Q模型（图中红色点）对假阳性的控制效果最差,而K（图中黑色点）和Q+K模型（图中蓝色点）对假阳性的控制较好,但是Q+K模型控制的过于严格,综合考虑,在后续分析中,对这6个性状,均选择K模型的结果进行进一步阐述。
显示原图|下载原图ZIP|生成PPT
图2不同环境下6个玉米株型相关性状Blup值的3种GWAS模型的QQ图 a：株高_2012;b：穗位高_2012;c：穗位高/株高_2012;d：穗上叶片数_2012;e：总叶片数_2012;f：穗上叶片数/总叶片数_2012。下同
-->Fig. 2Quantile-quantile (QQ) plots resulting from GWAS results using three methods for Blup data of six maize plant type traits in different environment a: Plant height (PH_2012), b: Ear height (EH_2012), c: The ratio of EH to PH (EH/PH_2012), d: Leaf number above ear (LNAN_2012), e: Total number of leaves (LN_2012), f: The ratio of LNAN to LN (LNAN/LN_2012). The same as below
-->

2.3.3 关联分析结果解析 利用覆盖全基因组的558 629个SNP标记对关联群体在2个环境下的6个株型性状进行关联分析,基于K模型,在P≤3.99E-6下,6个株型性状一共检测到56个SNPs-性状关联。2个环境联合分析条件下（Blup）,6个株型相关性状的曼哈顿图见图3。
显示原图|下载原图ZIP|生成PPT
图3浚县与郑州2个环境结合分析（Blup）下株型相关性状的曼哈顿图
-->Fig. 3The Manhattan plots of six maize plant type traits combine Xunxian and Zhengzhou location
-->

用55万SNPs标记对该群体LD的评估表明,该群体的衰减距离为50 kb^[24],基于此衰减距离,56个SNPs-性状关联涉及到的SNP一共分布于17个位点内,该结果与数量性状位点的本质相一致,即由微效多基因控制;其中在郑州试验站、浚县实验站及Blup下分别检测到8个、5个和9个位点。具体来说（表4和图4）,郑州环境下一共检测到10个SNP-性状关联,涉及到8个位点,分别是位于第5和第8染色体与株高有关的2个位点,分别解释8.73%和9.54%的株高表型变异;位于玉米第3、5、6、8染色体与穗位高有关的4个位点,每个位点解释的表型变异在8.19%—10.10%不等,平均为8.83%;及位于第4和第8染色体与穗位高与株高比值相关的2个位点,分别解释8.83%和10.42%的表型变异。并没有检测到与叶片数相关性状的显著位点。而在浚县环境下,一共检测到5个显著位点,分别为与株高有关的位于第1和第7染色体的2个位点,分别可解释9.69%和9.66%的表型变异;及与穗位高,穗位高与株高比值及总叶片数有关的分别位于第8、7、10染色体的3个位点,这3个位点分别能解释10.31%、9.34%和10.46%表型变异。2个环境联合分析（Blup）下,在第9染色体共检测到与穗位高相关的8个显著位点（涉及到16个显著的SNP）,每个位点平均可解释8.76%的表型变异;同时发现位于第5、8、9染色体的3个与株高有关的位点,分别解释8.12%、9.20%和8.58%的表型变异,其中位于第9染色体的位点里有11个显著的SNPs,并且该位点与穗位高的位点为同一个位点。还发现1个位于第8染色体的与穗位高与株高比值的位点,可解释10.12%的表型变异;2个与总叶片数有关的位点,位于第7和第10染色体,分别解释8.12%和8.22%的表型变异;另外还发现2个分别位于第3和第5与穗上叶片数有关的位点,分别可解释8.06%和8.88%的表型变异。
显示原图|下载原图ZIP|生成PPT
图4各个环境下株型相关性状显著位点的分布
-->Fig. 4Comparison of significant association loci for six plant type traits under different locations
-->

Table 4
表4
表4全基因组关联分析鉴定到的显著位点及有关信息以及每个位点里的最可能的候选基因及其功能注释
Table 4List of significant loci and their detailed information and a plausible candidate gene and its annotation within the significant loci identified by GWAS

位点a Loci	染色体 Chromosome	性状 Trait	试验点 Location	峰值SNP Lead SNP	物理位置b Position (bp)	P值c P value	观察数 N	R2 d (%)	QTL 区间 QTL region		候选基因e Candidate gene	基因功能注释 Gene annotation
									-50 kb	+50 kb
1	2	PH	XX	Chr.2.S_5839810	5839810	2.53E-06	241	9.69	5789810	5889810	GRMZM2G023921	DNA结合WRKY转录因子 DNA-binding WRKY factor
1	2	PH	XX	Chr2.S_5840172	5840172	2.53E-06	241	9.69	5790172	5890172	GRMZM2G023921	DNA结合WRKY转录因子 DNA-binding WRKY factor
1	2	PH	XX	Chr2.S_5840587	5840587	2.53E-06	241	9.69	5790587	5890587	GRMZM2G023921	DNA结合WRKY转录因子 DNA-binding WRKY factor
1	2	PH	XX	Chr2.S_5840714	5840714	2.53E-06	241	9.69	5790714	5890714	GRMZM2G023921	DNA结合WRKY转录因子 DNA-binding WRKY factor
1	2	PH	XX	Chr2.S_5840726	5840726	2.53E-06	241	9.69	5790726	5890726	GRMZM2G023921	DNA结合WRKY转录因子 DNA-binding WRKY factor
2	3	EH	ZZ	Chr3.S_173741502	173741502	3.41E-06	276	8.19	173691502	173791502	GRMZM2G167247	未知Unknown
3	3	LNAN	Blup	Chr3.S_199660028	199660028	2.72E-06	282	8.06	199610028	199710028	GRMZM2G005082	水解酶活性Hydrolase activity
3	3	LNAN	Blup	Chr3.S_199660031	199660031	2.72E-06	282	8.06	199610031	199710031	GRMZM2G005082	水解酶活性Hydrolase activity
4	4	EH/PH	ZZ	Chr4.S_198372985	198372985	2.16E-06	272	8.83	198322985	198422985	GRMZM2G071196	催化活性Catalytic activity
5	5	PH	ZZ	PZE-105103501	156842382	2.68E-06	269	8.73	156792382	156892382	GRMZM2G021998	F-box结构域蛋白 F-box domain containing protein
6	5	LNAN	Blup	Chr5.S_174323150	174323150	1.12E-06	284	8.88	174273150	174373150	GRMZM2G113409	离子通道活性Ion channel activity
7	5	PH	Blup	Chr5.S_203402792	203402792	2.61E-06	279	8.12	203352792	203452792	GRMZM2G542272	未知Unknown
7	5	PH	XX	Chr5.S_203402792	203402792	2.92E-06	235	9.66	203352792	203452792	GRMZM2G542272	未知Unknown
8	5	EH	ZZ	Chr5.S_208899879	208899879	5.35E-07	264	10.10	208849879	208949879	GRMZM2G045927	铜离子结合Copper ion binding
9	6	EH	ZZ	PZE-106067109	119876370	2.65E-06	264	8.35	119826370	119926370	GRMZM2G162755	花青素3-O-葡糖基转移酶 Anthocyanidin 3-O-glucosyltransferase
10	7	LN	Blup	Chr7.S_115182194	115182194	3.34E-06	276	8.12	115132194	115232194	GRMZM2G071059	CCR4-NOT转录复合物 CCR4-NOT transcription complex
10	7	LN	Blup	Chr7.S_115182505	115182505	1.81E-06	277	8.46	115132505	115232505	GRMZM2G071059	CCR4-NOT转录复合物 CCR4-NOT transcription complex
11	7	EH/PH	XX	Chr7.S_143823672	143823672	3.31E-06	235	9.34	143773672	143873672	GRMZM2G052288	多胺的生物合成过程 Polyamine biosynthetic process
12	8	EH	ZZ	Chr8.S_151906818	151906818	1.24E-06	272	9.09	151856818	151956818	GRMZM2G144021	硫酯酶/脂肪酶/酯酶 Esterase/lipase/thioesterase
12	8	EH	ZZ	Chr8.S_151907047	151907047	2.99E-06	271	8.41	151857047	151957047	GRMZM2G144021	硫酯酶/脂肪酶/酯酶 Esterase/lipase/thioesterase
13	8	EH/PH	Blup	Chr8.S_169624375	169624375	5.50E-07	269	10.12	169574375	169674375	GRMZM5G897342	细胞的氧化还原平衡 Cell redox homeostasis
位点a Loci	染色体 Chromosome	性状 Trait	试验点 Location	峰值SNP Lead SNP	物理位置b Position (bp)	P值c P value	观察数 N	R2 d (%)	QTL 区间 QTL region		候选基因e Candidate gene	基因功能注释 Gene annotation
									-50 kb	+50 kb
13	8	EH/PH	ZZ	Chr8.S_169624375	169624375	7.10E-07	257	10.28	169574375	169674375	GRMZM5G897342	细胞的氧化还原平衡 Cell redox homeostasis
13	8	EH/PH	Blup	Chr8.S_169624445	169624445	5.46E-07	269	10.13	169574445	169674445	GRMZM5G897342	细胞的氧化还原平衡 Cell redox homeostasis
13	8	EH/PH	ZZ	Chr8.S_169624445	169624445	5.14E-07	257	10.56	169574445	169674445	GRMZM5G897342	细胞的氧化还原平衡 Cell redox homeostasis
14	8	PH	Blup	Chr8.S_170276048	170276048	1.70E-06	262	9.20	170226048	170326048	GRMZM2G440125	锌离子结合Zinc ion binding
14	8	PH	ZZ	Chr8.S_170276048	170276048	1.80E-06	251	9.54	170226048	170326048	GRMZM2G440125	细胞的氧化还原平衡 Cell redox homeostasis
15	8	EH	XX	Chr8.S_174252562	174252562	3.90E-06	222	10.31	174202562	174302562	GRMZM2G148985	热、酸稳定磷酸化蛋白 Heat- and acid-stable Phosphoprotein
16	9	EH	Blup	Chr9.S_139516295	139516295	2.71E-07	273	10.20	139466295	139566295	GRMZM2G161293	乙酰葡糖氨基转移酶活性 Acetylglucosaminyltransferase activity
16	9	PH	Blup	Chr9.S_139516295	139516295	7.35E-07	273	9.29	139466295	139566295	GRMZM2G161293	乙酰葡糖氨基转移酶活性 Acetylglucosaminyltransferase activity
14	8	PH	ZZ	Chr8.S_170276048	170276048	1.80E-06	251	9.54	170226048	170326048	GRMZM2G440125	乙酰葡糖氨基转移酶活性 Acetylglucosaminyltransferase activity
15	8	EH	XX	Chr8.S_174252562	174252562	3.90E-06	222	10.31	174202562	174302562	GRMZM2G148985	热、酸稳定磷酸化蛋白 Heat- and acid-stable Phosphoprotein
16	9	EH	Blup	Chr9.S_139516425	139516425	2.37E-06	273	8.52	139466425	139566425	GRMZM2G161293	热、酸稳定磷酸化蛋白 Heat- and acid-stable Phosphoprotein
16	9	PH	Blup	Chr9.S_139516425	139516425	2.27E-06	273	8.49	139466425	139566425	GRMZM2G161293	热、酸稳定磷酸化蛋白 Heat- and acid-stable Phosphoprotein
16	9	EH	Blup	Chr9.S_139516431	139516431	2.37E-06	273	8.52	139466431	139566431	GRMZM2G161293	热、酸稳定磷酸化蛋白 Heat- and acid-stable Phosphoprotein
16	9	PH	Blup	Chr9.S_139516431	139516431	2.27E-06	273	8.49	139466431	139566431	GRMZM2G161293	热、酸稳定磷酸化蛋白 Heat- and acid-stable Phosphoprotein
16	9	EH	Blup	Chr9.S_139516436	139516436	2.37E-06	273	8.52	139466436	139566436	GRMZM2G161293	热、酸稳定磷酸化蛋白 Heat- and acid-stable Phosphoprotein
16	9	PH	Blup	Chr9.S_139516436	139516436	2.27E-06	273	8.49	139466436	139566436	GRMZM2G161293	热、酸稳定磷酸化蛋白 Heat- and acid-stable Phosphoprotein
16	9	EH	Blup	Chr9.S_139516462	139516462	2.37E-06	273	8.52	139466462	139566462	GRMZM2G161293	热、酸稳定磷酸化蛋白 Heat- and acid-stable Phosphoprotein
16	9	PH	Blup	Chr9.S_139516462	139516462	2.27E-06	273	8.49	139466462	139566462	GRMZM2G161293	热、酸稳定磷酸化蛋白 Heat- and acid-stable Phosphoprotein
16	9	EH	Blup	Chr9.S_139516723	139516723	1.33E-06	274	8.92	139466723	139566723	GRMZM2G161293	热、酸稳定磷酸化蛋白 Heat- and acid-stable Phosphoprotein
位点a Loci	染色体 Chromosome	性状 Trait	试验点 Location	峰值SNP Lead SNP	物理位置b Position (bp)	P值c P value	观察数 N	R2 d (%)	QTL 区间 QTL region		候选基因e Candidate gene	基因功能注释 Gene annotation
16	9	EH	Blup	SYN27192	139516741	3.84E-06	246	8.86	139466741	139566741	GRMZM2G161293	热、酸稳定磷酸化蛋白 Heat- and acid-stable Phosphoprotein
16	9	PH	Blup	Chr9.S_139516752	139516752	2.30E-06	273	8.48	139466752	139566752	GRMZM2G161293	热、酸稳定磷酸化蛋白 Heat- and acid-stable Phosphoprotein
16	9	EH	Blup	Chr9.S_139516765	139516765	1.42E-06	274	8.90	139466765	139566765	GRMZM2G161293	热、酸稳定磷酸化蛋白 Heat- and acid-stable Phosphoprotein
16	9	EH	Blup	Chr9.S_139516778	139516778	1.42E-06	274	8.90	139466778	139566778	GRMZM2G161293	热、酸稳定磷酸化蛋白 Heat- and acid-stable Phosphoprotein
16	9	EH	Blup	Chr9.S_139516803	139516803	1.90E-06	273	8.69	139466803	139566803	GRMZM2G161293	热、酸稳定磷酸化蛋白 Heat- and acid-stable Phosphoprotein
16	9	EH	Blup	Chr9.S_139516842	139516842	2.67E-06	273	8.42	139466842	139566842	GRMZM2G161293	热、酸稳定磷酸化蛋白 Heat- and acid-stable Phosphoprotein
16	9	EH	Blup	Chr9.S_139516886	139516886	1.24E-06	273	9.02	139466886	139566886	GRMZM2G161293	热、酸稳定磷酸化蛋白 Heat- and acid-stable Phosphoprotein
16	9	PH	Blup	Chr9.S_139516886	139516886	1.35E-06	273	8.89	139466886	139566886	GRMZM2G161293	热、酸稳定磷酸化蛋白 Heat- and acid-stable Phosphoprotein
16	9	EH	Blup	Chr9.S_139516887	139516887	2.37E-06	273	8.52	139466887	139566887	GRMZM2G161293	热、酸稳定磷酸化蛋白 Heat- and acid-stable Phosphoprotein
16	9	PH	Blup	Chr9.S_139516887	139516887	2.27E-06	273	8.49	139466887	139566887	GRMZM2G161293	热、酸稳定磷酸化蛋白 Heat- and acid-stable Phosphoprotein
16	9	EH	Blup	Chr9.S_139516888	139516888	2.37E-06	273	8.52	139466888	139566888	GRMZM2G161293	热、酸稳定磷酸化蛋白 Heat- and acid-stable Phosphoprotein
16	9	PH	Blup	Chr9.S_139516888	139516888	2.27E-06	273	8.49	139466888	139566888	GRMZM2G161293	热、酸稳定磷酸化蛋白 Heat- and acid-stable Phosphoprotein
16	9	EH	Blup	Chr9.S_139517017	139517017	1.61E-06	273	8.81	139467017	139567017	GRMZM2G161293	热、酸稳定磷酸化蛋白 Heat- and acid-stable Phosphoprotein
16	9	PH	Blup	Chr9.S_139517017	139517017	1.47E-06	273	8.75	139467017	139567017	GRMZM2G161293	热、酸稳定磷酸化蛋白 Heat- and acid-stable Phosphoprotein
16	9	EH	Blup	Chr9.S_139517032	139517032	2.92E-06	276	8.27	139467032	139567032	GRMZM2G161293	热、酸稳定磷酸化蛋白 Heat- and acid-stable Phosphoprotein
16	9	PH	Blup	Chr9.S_139517032	139517032	3.47E-06	276	8.06	139467032	139567032	GRMZM2G161293	热、酸稳定磷酸化蛋白 Heat- and acid-stable Phosphoprotein
17	10	LN	Blup	Chr10.S_128868215	128868215	3.95E-06	283	7.98	128818215	128918215	GRMZM2G074689	钙离子结合Calcium ion binding
17	10	LN	XX	Chr10.S_128868215	128868215	1.24E-06	239	10.46	128818215	128918215	GRMZM2G074689	钙离子结合Calcium ion binding

^a All QTLs with overlapping QTL regions were categorized as a non-redundant QTL or loci, ID of loci was numbered according to chromosome and position by a ascending order method; ^b Physical position of each SNP based on B73 RefGen_v2. ^c P value of the corresponding trait calculated by K model. ^d The phenotypic variance explained by the corresponding locus. ^e A plausible biological candidate gene in the locus or the nearest annotated gene to the lead SNP. ZZ: Zhengzhou; XX: Xunxian^a所有有重叠区间的QTL被归为一个非重复的QTL或位点,位点的ID根据染色体及物理位置升序排序命名;^b每个SNP的物理位置是基于B73参考基因组v2版本得到;^c相应性状的P值是通过K模型计算得到;^d相应位点解释的表型变异;^e位点内一个可能的生物学候选基因或离峰值SNP最近的候选基因。ZZ表示郑州,XX表示浚县,Blup表示2个环境的综合
新窗口打开
同时发现有4个位点,分别为表4中的第7、13、14和17位点可同时在2个以上环境中被检测到,表明这些位点受环境影响较小,可在不同环境下稳定遗传,这些信息为将来继续深入分析这些位点提供了有价值的信息。值得注意的是,在2个环境及Blup下均没检测到与穗上叶片数与总叶片数比值相关的位点（图3、图4和表4）,表明这一性状受环境影响较大或在该群体中不存在控制这一性状的主效位点。

2.4 株型相关候选基因确定

基于该群体的LD衰减距离,在峰值SNP上下游各50 kb范围内搜索候选基因。6个株型性状的17个显著关联位点上下游各50 kb范围内一共搜到80个候选基因,其中有功能注释的基因有42个。根据每一个位点内所有候选基因的功能注释及其在B73的各组织的表达情况,选择一个最可能的候选基因及其注释列于表4。

3 讨论

由于不同性状对GWAS不同模型在控制Ⅰ型错误（假阳性）的敏感性不同,即存在性状依赖性,因此,在进行GWAS分析前首要考虑的是GWAS的统计功效,即最大可能的降低假阳性。先前研究表明,对于玉米的开花期、穗位高和穗粗性状来说,K模型和Q+K模型均能很好的控制假阳性,但是Q+K模型比K模型更好^[26]。本研究发现Q+K模型对假阳性的控制较严格（图2）,若采用Q+K模型则会产生更多的假阴性结果,而K模型对假阳性的控制效果最好,因此,最终选择表现较好的K模型对GWAS结果进行进一步解析。
玉米的农艺及产量相关性状的研究近年来报道较多,陈荣江等^[27]通过研究玉米的穗长、穗行数、行粒数、百粒重、株高、穗位高、穗重及穗粗在内的8个性状与产量之间的关系,发现这8个性状与产量存在正相关;汤继华等^[2]以豫玉22构建的260多个F_2:3家系为材料,研究了不同环境条件下的穗位高、雄穗分枝数、茎粗与产量之间的关系,发现与产量之间也存在显著正相关。本研究以284份遗传性丰富的温带、热带、亚热带材料为研究对象,在郑州和浚县2个环境下对玉米株高、穗位高、总叶片数和穗上叶片数以及株高/穗位高、总叶片数/穗上叶片数等6个株型性状进行分析,发现在2个环境下的不同性状之间均存在显著正相关或者负相关,说明各性状的发育存在着相互协同促进作用。通过对郑州和浚县这6个株型性状的全基因组关联分析,在K模型下,一共检测到56个SNP-性状关联（株高、穗位高、穗位高与株高比值、总叶片数、穗上叶片数）,这些SNP涉及到17个位点。而穗上叶片数与总叶片数比值这一性状没有检测到显著的SNPs,表明该性状受环境影响较大或者在该关联群体中不存在控制这一性状的主效位点。在这17个位点里,同时发现有4个位点能够在2个环境中被检测到,表明这4个位点受环境影响较小,在不同环境下可以稳定存在。同时也发现,有的性状和多个位点显著关联,表明这些性状受微效多基因控制,符合数量性状的本质^[28]。
基于B73参考基因组及每一个显著位点的LD衰减距离,在这17个位点里一共发现80个候选基因,其中42个基因有功能注释,其余38个基因的功能未知;42个已知功能的基因涉及到转导蛋白、还原蛋白、金属离子转运蛋白、脂类蛋白、纤维素合酶等,功能涉及能量代谢、物质运输、生物调控、信号传导等生化代谢途径。
17个位点里,发现多个位点里存在多个显著的SNP,比如位点16里含有17个显著SNP,表明这些SNP间可能存在很强的连锁不平衡关系,这些成簇SNP所在区域,是控制株型性状的热点区域,有必要对其加以深入挖掘,为培育理想株型品种的选育提供重要研究价值。郑克志等^{[25, 29]}以玉米自交系T319与9406为亲本构建的242个重组自交系（F₈）对玉米株高和穗位高进行QTL分析,结果发现6个株高相关的QTL,其中有2个QTL同时控制穗位高。本研究也发现了4个位点（位点7、13、14和17）可同时在2个以上环境中被检测到（表4）,表明这些位点受环境影响较小,可在不同环境下稳定遗传,同时也发现有些位点同时控制多个性状,比如位点16同时控制株高和穗位高,这些信息为将来继续深入利用和克隆这些位点提供了有价值的信息。
本研究通过玉米株型性状显著关联的SNP所在位点对可能的候选基因进行预测。6个性状,有5个性状均能检测到显著位点,随后对每一个性状的显著位点的候选基因及其功能进行了预测（表4）,值得注意的是,在郑州和浚县2个环境的综合条件下,位于第9染色体的位点16（Chr.9：139516295—Chr.9：139517032）同时与株高和穗位高显著关联,并且该位点在先前的两个研究中也被报道与株高相关（假定2个研究定位的显著SNP物理位置小于5 Mb被认为同一个位点）。比如,YANG等^[12]通过A—D test共检测到株高相关的107个显著SNP,其中第9染色体的4个SNP（Chr9..S_140356597、Chr.9.S_140435314、Chr.9.S_142397931和Chr.9.S_142495942）与该位点共定位;SAMAYOA等^[30]通过关联分析检测到18个株高显著SNP,其中第9染色体的一个SNP, S9_142752809,与该位点共定位;该位点（位点16）内一共存在4个候选基因（GRMZM5G872568、GRMZM2G063875、GRMZM2G161293和GRMZM2G366532）,除GRMZM2G161293外,其余基因的功能未知。而GRMZM2G161293编码一个BGGP Beta-1-3-半乳糖-O-糖基-糖蛋白,具有乙酰葡糖转移酶活性,该酶催化UDP-N-乙酰氨基葡萄糖生成糖过程中的N-乙酰葡糖氨基残基的转移,猜测可能通过影响玉米籽粒可溶性糖含量进而影响产量,推测其为最可能的候选基因;另外发现在2个环境下能同时检测到的4个位点（位点7、13、14和17）,其中基因GRMZM2G005082、GRMZM5G897342、GRMZM2G440125和GRMZM2G074689分别位于位点7、13、14和17中,这4个基因分别编码水解酶,常温蛋白,锌离子结合蛋白及钙离子结合蛋白;这些同时控制多个性状的SNP或位点及其在2个环境中同时被检测到的SNP或位点可能与株型相关性状的候选基因紧密连锁,可以作为下一步研究的重点区域。

4 结论

本研究借助覆盖玉米全基因组的约56万个SNP标记,利用3种模型（Q、K和Q+K）对玉米总叶片数（LN）、穗上叶片数（LNAN）、穗上叶片数与总叶片数比值（LNAN/LN）、株高（PH）、穗位高（EH）、穗位高与株高的比值（EH/PH）等6个玉米株型相关性状进行全基因组关联分析,发现K模型对假阳性的控制效果最好。基于K模型的GWAS结果,一共检测到56个SNP-性状关联,涉及17个位点,对显著关联的SNP上下游各50 kb范围内候选基因进行搜索,共发现80个候选基因,其中42个具有功能注释。
致谢：华中农业大学作物遗传改良国家重点实验室严建兵教授为本试验提供关联群体材料及基因型,在此表示感谢。
（责任编辑 李莉）
The authors have declared that no competing interests exist.

---

参考文献文献选项 原文顺序
文献年度倒序
文中引用次数倒序
被引期刊影响因子