0 引言
【研究意义】氮是植物生长发育中至关重要的营养元素,合理的氮肥施用将显著的增加作物产量[1],氮营养供应不充足则会导致作物减产,而过量施用氮肥会导致环境污染等问题。因而氮素的高效利用是当前水稻育种关注的重要问题,通过选育氮高效吸收利用的水稻品种,从而减少田间氮肥施用以提高水稻产量已成为水稻育种的重要目标之一[2]。【前人研究进展】水稻的氮素利用效率（nitrogen use efficiency,NUE）受品种、氮肥用量、氮肥施用方法和土壤含氮量等因素的影响。近年来,通过对水稻氮响应相关突变体的研究和QTL定位等方式克隆与鉴定了多个与NUE相关的基因及位点。ZHANG等[3]鉴定了一个位于第12染色体上有关氮耐受性的QTL位点,并克隆了TOND1。HU等[4]通过在高世代回交群体分离鉴定了一个硝酸盐转运基因NRT1.1B,该基因通过增强硝酸盐的摄取、根向茎转运和氮响应基因的表达从而提高氮素利用率。SUN等[5]通过分析RIL群体在不同氮浓度下株高和分蘖数的变化鉴定了一个水稻主效的NUE数量性状位点qNGR9,实为已报道的直立密穗基因DEP1,携带dep1-1的植株表现出对氮吸收和同化的增强,从而增加水稻产量。在长期的进化过程中,植物形成了适应不同生存环境的生理机制和分子机理来响应环境的变化。在应对缺氮环境时,植物提高NUE的主要方式是通过提高氮的吸收、转运、同化和再活化等方式,且与环境之间存在复杂的互作及遗传网络关系[6-9]。TONG等[10]通过在2种氮处理条件下的株高、穗数、叶绿素含量、干重、湿重和产量变化等指标鉴定了31个氮响应相关的QTL位点。CHO等[11]在两种氮浓度下,通过水稻发育后期籽粒中氮含量、秸秆中氮含量、收获指数等多个指标鉴定了多个氮响应相关的QTL位点。同样,在不同的氮水平下水稻叶片呈现明显的叶形变化,而水稻叶片的长、宽和厚的改变也直接影响水稻的光合作用[12]。殷春渊等[13]研究表明,氮肥处理对水稻叶片的影响主要表现在倒二叶叶长和剑叶叶宽2个特性上。因此,叶片形态结构的变化也是水稻对不同水平氮响应的指标。近年来,全基因组关联分析（genome-wide association study,GWAS）成功的应用在一些模式生物研究和农作物的育种当中[14~18]。【本研究切入点】在水稻中,利用GWAS分析了包括产量、品质、抽穗期、株形等多种重要的农艺性状,同样也对病虫害抗性等生物胁迫和盐碱耐受性等非生物胁迫性状进行了全基因组关联分析,并鉴定了许多重要的SNP位点[14,16-17,19-22]。但是,与氮素利用相关的GWAS研究尚未见报道。【拟解决的关键问题】本研究通过选取134份来自不同地区的水稻品种作为关联分析的群体,利用在不同氮肥条件下叶宽的变化,开展水稻叶宽及氮响应相关SNP位点分析,旨在为进一步揭示水稻叶片形态调控和氮素高效利用的生理机制提供基础。1 材料与方法
1.1 试验材料
选择保存于中国水稻所（CNRRI）,国际水稻所（IRRI）种质资源库及来自多个国家和地区的134个地方品种,其中包括41份来自国内各省农科院的地方品种;36份来自20个国家的地方品种;57份国内大面积推广种植品种。1.2 试验设置及性状考察
试验于2015年在中国水稻研究所实验基地（浙江,富阳）进行,试验田土壤为水稻土（土壤含有机质36.9 g·kg-1、全氮2.73 g·kg-1、速效钾104 mg·kg-1、速效磷52 mg·kg-1、碱解氮156 mg·kg-1,pH 6.24）。试验设3个氮肥水平：低氮（N0,不另施氮肥）、正常氮（N1,纯氮96 kg·hm-2）和高氮（N2,纯氮192 kg·hm-2）,以施氮量为主区,品种为裂区,3次重复。主区面积5.0 m×23.0 m,播种前将其做成3条畦宽1.5 m,沟宽25 cm,长23 m的畦。水稻种于畦上,每品种种植4行,每行6株,株行距为20 cm×20 cm。不同氮处理间作埂隔离,并用塑料薄膜覆盖埂体,保证单独排灌。试验用关联群体于2015年5月10日采用点直播栽培,每穴播3—5粒种子,2叶1心时定苗2株/穴。以尿素为氮肥,按基肥﹕分蘖肥﹕穗肥 = 5﹕3﹕2分3次施入,其中基肥（播种前1 d）与分蘖肥（播种后25 d）统一施用,穗肥由于品种间生育期的差异,于孕穗始期分品种单独施用。6月中旬至9月中旬每隔3d检查品种的生育期,标记孕穗始期相近的品种,将穗肥小心撒于每个标记品种叶片上。为防止水分流动导致的串肥,穗肥施用期间（6月中旬至9月中旬）保持畦面湿润但无水层,仅沟中有水。磷、钾施用量为P2O5 75 kg·hm-2、K2O 150 kg·hm-2,其中磷肥为底肥一次施入,钙、镁、钾肥按基、蘖肥各50%分2次施入。其他栽培管理措施同一般大田。
在抽穗后10 d对剑叶叶宽进行调查,每品种选取中间的6株无边际效应的植株进行考察。参考王兰等[23]方法量取每片叶的最宽处为叶宽,每株取主茎分别记录剑叶的叶宽。
1.3 数据分析
1.3.1 性状相关性分析 使用Excel 2010进行数据的录入和整理,并计算其平均值和标准差。使用SAS 9.2进行进行不同氮处理下叶宽变化的相关分析和方差分析。1.3.2 基因型测定 采用CTAB抽提法提取134份品种幼叶基因组DNA[24],检测合格后用于高通量测序。每品种测序reads数据量在5Gb以上,相当于10×以上的基因组测序。从reads到SNP的检测参考ZHOU等[17]方法进行。SNP使用ANNOVAR （Version: 2013-08-23）软件包进行注释[25]。
1.3.3 群体结构与亲缘关系分析 根据SNP信息,使用PHYLIP 3.68软件包基于邻接法构建品种间系统进化树[26],使用在线工具iTOL（http://itol.embl.de/）进行进化树的优化。使用EIGENSOFT软件中的smartpca程序包进行主成分分析（PCA）[27]。
1.3.4 叶宽性状的全基因组关联分析 采用基于混合线性模型（mixed linear model,MLM）的EMMAX（Efficient Mixed-Model Association eXpedited）软件包进行目标性状与SNP标记间的关联分析[20, 28],将EMMAX软件计算的亲缘关系K矩阵及基于EIGENSOFT软件计算的群体结构Q纳入分析模型中,关联分析参数设定为：emmax -v,-d 10,-t,-p,-k,-c,-o。在考虑关联分析群体大小和SNP数量的情况下,对于关联分析的显著性,采取P<10-6作为宽松阈值,而进一步的严格阈值采取Bonferroni校验阈值[29]。
2 结果
2.1 3种氮处理条件下叶宽性状多样性及其变化
在134份材料中,剑叶叶宽表现丰富的变异,正常氮处理下变幅为9.93—29.00 mm。如木邦谷、罗大穗等地方品种表现为株型较大,分蘖少,成熟期叶片少,剑叶叶宽大于20 mm;而地方品种英德大叶拍、合江23等表现为株型小,成熟期叶片多,叶片呈细长状,叶宽小于10 mm。由此表明,此关联群体具有充分的表型变异,能够进行基于叶宽的全基因组关联分析。对裂区设计整体进行方差分析,变异来源于氮肥使用量、品种及品种与氮肥交互作用的P值均小于0.0001,而区组与氮处理间的交互作用P为0.2444,这说明叶宽变异主要受氮处理及品种差异引起,不同品种在不同氮浓度下展现不同的叶宽变化。对134份品种成熟期剑叶宽调查结果表明,在低、中、高3种氮处理下,剑叶叶宽均呈正态分布（图1-a—图1-c）,从频率分布来看,中氮和高氮处理下剑叶宽为15—18 mm的品种明显增多。在反映氮响应变化中,正常氮处理相对于低氮处理表现为叶宽变大,而相对于高氮处理变化要小（图1-d）。相关性分析亦表明,氮处理与剑叶叶宽呈显著正相关（表1）,这表明,叶宽受外界氮水平的影响,且叶宽的变化能够反映不同氮处理。
Table 1
表1
表1不同氮浓度下叶宽统计及其相关性分析
Table 1Leaf widths at different nitrogen concentrations and their correlation analysis
| 氮水平 N level | 剑叶叶宽 Width of flag leaf(mm) | 相关性分析 Correlationship analysis | ||
|---|---|---|---|---|
| 均值 Mean | 变幅 Rang | N0 | N1 | |
| N0 | 15.55&#x000b1;2.91 | 10.17&#x02014;28.53 | ||
| N1 | 16.45&#x000b1;3.29 | 9.93&#x02014;29.00 | 0.8362** | |
| N2 | 16.55&#x000b1;3.15 | 11.00&#x02014;27.33 | 0.8263** | 0.8896** |
新窗口打开
显示原图|下载原图ZIP|生成PPT图1不同氮处理下水稻剑叶叶宽分布
-->Fig. 1The distribution of flag leaf width under different nitrogen treatments
-->
根据不同氮处理下叶宽相对变化,将叶宽变化分为四种响应模式（图2-e）。有35.2%的品种表现为在高氮条件下叶宽增加而低氮条件下叶宽减小,有18.1%的品种在低氮和高氮条件下剑叶宽均增大,另有36.2%的品种在低氮和高氮条件下叶宽均减小,还有5.7%的品种在低氮条件下剑叶变宽而在高氮条件下反而变窄。由此可知,该关联群体在不同氮处理下叶宽相对变化呈现丰富的变异,可以利用叶宽的相对变化来反映氮素利用效率。
2.2 关联分析SNP标记
基于134份水稻的测序数据,结合使用BWA和GATK软件进行SNP的提取,使用最小等位基因频率大于0.03（MAF>0.03）进行SNP的过滤,最终提取到3 356 591个SNP位点,平均每kb有8.6个SNP位点。为进一步分析SNP位点信息,使用ANNOVAR软件对SNP数据注释（图2）,其中,有2 814个位点注释为剪接位点相关,此类位点通过影响基因转录剪接的变化导致基因功能改变;409 387个SNP位点注释到外显子区,引起核酸编码信息的改变,其中59%的SNP为非同义突变,造成编码信息的改变,从而引起蛋白序列的变化;另外有12 221个位点注释与终止密码子相关,终止密码子的缺失或增加会导致编码蛋白移码突变或提前终止。
显示原图|下载原图ZIP|生成PPT图2关联分析SNP标记分布及其功能注释
-->Fig. 2The distribution and functional annotation of association analysis of SNP marker
-->
2.3 群体的遗传结构分析
基于134品种重测序的SNP数据,开展亲缘关系和系统进化分析。结果表明,此关联群体主要分为籼稻和粳稻两个亚群。由GAPIT软件基于SNP计算得到的亲缘关系系数表明,此134份水稻品种主要分为71份籼稻品种和63份粳稻品种,亚种内部分品种的亲缘关系较近,形成小的亚群（图3-a和图3-b）。基于SNP数据计算品种间遗传距离和NJ（Neighbor- joining）法构建系统发育树表明,此关联群体主要为籼粳之间的分歧（图3-c）。主成分分析（principal component analysis,PCA）同样也说明了此关联群体籼粳分化的群体结构。此外,基于PCA三点图可知,籼稻亚种间分化较为集中,粳稻亚种间存在不同程度的分化,较为分散,部分品种间存在小亚群结构（图3-d）。
显示原图|下载原图ZIP|生成PPT图3关联群体的群体结构
-->Fig. 3The population structures of association population
-->
2.4 叶宽性状的全基因组关联分析
通过对不同氮处理条件下134份水稻成熟期剑叶叶宽的全基因组关联分析表明,在不同氮处理条件下,与叶宽相关联的SNP位点各不相同（表2）。以-log10(p)大于6为阈值,共检测到14个与叶宽性状显著相关的SNP位点（图4）。在低氮和中氮水平下,均检测到6个显著的SNP位点与叶宽性状相关。而高氮水平只检测到位于第12染色体上的2个显著位点,其中一个位点与中氮水平一致,推测该调控叶宽的位点受氮水平影响较小。Table 2
表2
表23种氮处理下剑叶宽的关联位点
Table 2Association sites of flag leaf width in three nitrogen treatments
| 性状 Trait | 染色体 Chr. | 物理位置 Position (bp) | 主要等位基因 Major allele | 次要等位基因 Minor allele | 最小等位基因频率 MAF | P值(-lg p) P value |
|---|---|---|---|---|---|---|
| N0 | 9 | 19489161 | C | T | 0.50 | 7.46 |
| 4 | 19717235 | T | A | 0.50 | 7.06 | |
| 4 | 20263943 | C | T | 0.46 | 6.69 | |
| 4 | 20706884 | G | A | 0.49 | 6.43 | |
| 6 | 27244159 | T | C | 0.47 | 6.08 | |
| 1 | 35634562 | C | T | 0.48 | 6.04 | |
| N1 | 12 | 15066507 | C | T | 0.05 | 8.45 |
| 6 | 5564126 | G | A | 0.08 | 7.13 | |
| 12 | 14729445 | C | T | 0.04 | 6.32 | |
| 12 | 14869227 | G | T | 0.06 | 6.15 | |
| 9 | 11074359 | G | A | 0.08 | 6.15 | |
| 12 | 14955889 | G | A | 0.05 | 6.01 | |
| N2 | 12 | 15066507 | C | T | 0.05 | 6.90 |
| 12 | 20822979 | T | C | 0.26 | 6.08 |
新窗口打开
显示原图|下载原图ZIP|生成PPT图4不同氮水平下剑叶宽的全基因组关联分析
-->Fig. 4Genome-wide association studies of flag leaf width at different nitrogen levels
-->
在低氮水平下检测到的显著性SNP位点最小等位基因频率均大于0.46,此类SNP在关联群体中分布较为广泛。而在中氮和高氮水平下检测到的SNP位点最小等位基因频率均较小,属于一类较为稀有的突变位点。在中氮水平下,第12染色体上显著性SNP位点的P值（-lgP）为8.45,高于Bonferroni校正后的阈值水平。对中氮水平下叶宽性状显著关联的位点进行候选基因分析时发现,位于第12染色体上的Chr12_14869227显著位点±200 kb的候选区间内,包含一注释为Cytochrome P450功能相关蛋白的候选基因LOC_Os12g25660,其与业已报道的叶宽性状相关基因OsBR6ox同属于细胞色素P450家族（图5）。OsBR6ox编码水稻BR生物合成途径中的一个关键的BR-6氧化酶,该酶活下降的植株表现为BR缺乏症状,而使叶宽受到显著的影响[30-31]。在高氮水平下,检测到另一个位于第12染色体上显著性关联的峰,该区段包含已知的OsCSLD4,OsCSLD4编码一个细胞壁生物合成和植株生长必须的糖基转移酶家族的纤维素合酶类似蛋白（cellulose synthase-like,CSL）,OsCSLD4的突变导致水稻叶片变窄、半卷和矮化等表型[32-33]。
显示原图|下载原图ZIP|生成PPT图5第12染色体上在N1水平下的关联显著位点候选区间分析
-->Fig. 5Candidate interval analysis of significant loci on chromosome 12 at N1 level
-->
2.5 不同氮水平剑叶叶宽的全基因组关联分析
进一步使用不同氮水平下叶宽的相对变化值来分析叶宽对氮水平的响应,共检测到28个显著的SNP位点,其中20个与低氮响应相关,8个与高氮响应相关（表3、图6）。在低氮环境下,叶宽的响应位点比高氮环境下多,这与叶宽变化在低氮环境下表型变异丰富,而在高氮水平下叶宽变化不明显有关。在高氮水平下,位于第7染色体上的显著性位点Chr7_20639350,其-lg(P)为9.01,超过了Bonferroni检验阈值;另外,还检测到位于第1染色体上的显著性峰,该区段存在多个显著性峰值点（图7）,在此区段内存在一个已知的OsATG7,该基因与氮素利用相关,osatg7表现为叶片明显早衰,氮的再活化作用被抑制[34-35],推测该基因为此检测区间的主效基因。Table 3
表3
表3低氮和高氮水平下剑叶叶宽变化的关联位点
Table 3Association sites of flag leaf width variations under low N and high N stress
| 性状 Trait | 染色体 Chr. | 物理位置 Position (bp) | 主要等位基因 Major allele | 次要等位基因 Minor allele | 最小等位基因频率 MAF | P值(-lg P) P value |
|---|---|---|---|---|---|---|
| N0-N1 | 5 | 19246184 | C | T | 0.04 | 7.99 |
| 8 | 21630810 | C | A | 0.44 | 7.48 | |
| 1 | 5834252 | C | T | 0.04 | 7.39 | |
| 1 | 7340161 | C | T | 0.03 | 7.22 | |
| 9 | 6720357 | C | A | 0.03 | 7.15 | |
| 7 | 7512991 | A | G | 0.04 | 7.07 | |
| 4 | 2943821 | G | A | 0.04 | 6.99 | |
| 6 | 11045452 | G | A | 0.04 | 6.66 | |
| 2 | 17789413 | C | T | 0.45 | 6.60 | |
| 4 | 8563109 | G | A | 0.03 | 6.52 | |
| 3 | 21613067 | G | A | 0.41 | 6.43 | |
| 10 | 2300766 | C | T | 0.03 | 6.29 | |
| 7 | 17179329 | C | T | 0.05 | 6.28 | |
| 5 | 26409615 | T | C | 0.06 | 6.27 | |
| 8 | 21424507 | C | T | 0.06 | 6.22 | |
| 10 | 15030259 | T | C | 0.46 | 6.14 | |
| 10 | 14971679 | A | G | 0.03 | 6.12 | |
| 2 | 17573081 | C | T | 0.44 | 6.08 | |
| 3 | 22432529 | G | A | 0.03 | 6.05 | |
| 3 | 22169582 | G | A | 0.42 | 6.01 | |
| N2-N1 | 7 | 20639350 | G | A | 0.06 | 9.01 |
| 1 | 26146954 | T | C | 0.03 | 8.04 | |
| 1 | 26512077 | T | C | 0.06 | 7.20 | |
| 1 | 25812567 | G | A | 0.04 | 7.02 | |
| 7 | 20656204 | C | T | 0.03 | 6.93 | |
| 1 | 23639246 | G | T | 0.04 | 6.41 | |
| 5 | 4349907 | C | T | 0.04 | 6.19 | |
| 9 | 18747007 | T | C | 0.49 | 6.12 |
新窗口打开
图6剑叶叶宽对低氮和高氮差异响应的全基因组关联分析
-->Fig. 6Genome-wide association analysis of flag leaf width for response to low nitrogen and high nitrogen
-->
显示原图|下载原图ZIP|生成PPT图7第1染色体上在N2水平的关联显著位点候选区间分析
-->Fig. 7Candidate interval analysis of significant loci on chromosome 1 at high nitrogen level
-->
3 讨论
水稻成熟期剑叶形态是决定植株光合过程的受光面积和株型结构的关键因素,并进而影响灌浆期光合生物同化量的大小及水稻产量[36-37],因而改善水稻成熟期剑叶形态成为水稻保持优势株型和高产的重要手段。但水稻叶片形态是复杂的受微效多基因控制的数量性状,易受光、温、水、肥等环境影响[38]。近年来,随着GWAS在植物领域研究的应用,为发掘重要种质资源和快速鉴定新基因提供了新的分析工具[14,16]。从GWAS的结果来看,水稻剑叶叶宽对氮肥胁迫的响应更为显著。在低氮或高氮水平下检测到的关联位点远多于正常氮肥氮肥水平下剑叶叶宽的关联位点,尤其以低氮响应为甚,这可能与剑叶叶宽易受氮肥水平的影响以及不同种质资源对低氮或高氮的响应方式不同有关。因而,对于易受环境影响的性状在开展GWAS分析时,最好能结合多环境下的分析或严格控制其他环境因素的影响。本研究采用3次重复的裂区设计,但没有作多年/异地的重复,可能会对部分关联位点的精准选择产生影响。在进行籼粳型的群体进行关联分析时,将BN （Balding-Nichols）矩阵纳入基于MLM模型的关联分析中,进行群体结构的校正,可以有效规避由于群体结构造成的影响,从而减少假阳性位点的产生。
氮素的同化和利用对植物至关重要,前人的研究结果来看,分别在第1、2、3、5、7、10、11、12染色体上存在已报道的与氮素相关QTL位点。其中,位于第2染色体的RM5812 标记与本研究的低氮胁迫下的检测到的显著性位点区间相近[39]。同样,在低氮胁迫下检测到的位于第3、第5染色体上的显著性位点所在区间也有QTL位点报道过（http://archive.gramene.org/）。因此,基于显著性的关联位点,需要进一步对候选区间内的基因进行鉴定、分析和克隆。后期的进一步表达水平变化或基因敲除研究对关联位点进行功能验证十分重要,深入展开关联候选基因的功能研究将有助于准确揭示水稻氮素利用遗传机制奠定基础。
4 结论
在不同氮肥水平下对水稻剑叶叶宽进行全基因组关联分析,共检测到的42个显著的关联位点中,其中14个与剑叶叶宽相关,28个与氮肥的响应相关。对所检测到的位点进行候选区间分析,发现2个位于候选区间内分别调控水稻剑叶叶宽及氮响应的已知基因;检测到多条染色体上的显著位点与已报道与氮素利用相关QTL区间相近。此外,基于候选区间基因的功能注释,检测到1个与已知调控叶宽的同源基因。致谢：中国科学院遗传与发育生物学研究所田志喜研究员和李云海研究员为关联分析提供指导及关联群体材料,在此表示感谢。
The authors have declared that no competing interests exist.
