2. 华中农业大学, 农业部猪遗传育种重点实验室, 武汉 430070
, MA Yunlong1,2, XIANG Tao1,2, ZHU Mengjin1,2, YU Mei1,2, LI Xinyun1,2, Liu Xaolei1,2*
, Zhao Shuhong1,2*
2. Key Laboratory of Swine Genetics and Breeding of Ministry of Agriculture, Huazhong Agricultural University, Wuhan 430070, China
动物育种从常规表型选育到利用最佳线性无偏预测(best linear unbiased prediction，BLUP)估计育种值(estimated breeding values, EBVs)进行选育，再到标记辅助选育(marker assisted selection, MAS)，使得育种的遗传进展不断加快，选育效果愈见明显[1]。但是，一方面，由于目前经过功能验证的主基因数目有限，解释的遗传变异较少，另一方面，部分经济性状遗传力较低、遗传机制复杂，少数几个育种标记的遗传改良效果较为有限[2]。随着高通量测序和芯片技术的发展，新一代的选育技术被逐渐推广并应用到动植物遗传改良中，即全基因组选择技术(genomic selection, GS)。该技术由Meuwissen等[3]于2001年首次提出，结合个体间基因组遗传信息的差异，从表型中剖分出真实遗传差异，能够有效改良低遗传力的性状，真正实现由分子技术指导育种实践的突破。近年来，全基因组选择的计算方法、影响因素、应用策略和育种方案在众多学者的不断开发与研究中大量涌现。
全基因组选择在奶牛育种中的应用效果尤为突出，自2009年开始奶牛的选育已完全由全基因组选择主导，相比于后裔测定，全基因组选择可以显著地缩短世代间隔，极大地降低奶牛选育成本，并且已取得较大遗传进展[4-6]。猪的全基因组选择育种于2009年在欧美初步尝试，2012-2014年数个国外大型育种集团，如丹育、PIC、海波尔等相继开展猪的全基因组选择育种[7]。目前，我国的猪育种企业也逐步启用全基因组选择育种策略，与各高校、研究所之间紧密合作，建立基因组选育方案和体系，全基因组选择成为了当今猪育种的一大热点[8]。相比奶牛及猪育种，绵羊和山羊的全基因组选择研究较少，有研究表明[9]，在父母本和子一代组成的群体家谱进行羊奶产量育种值评估中，一步法(single-step blup, SSBLUP)评估的群体育种值准确性优于GBLUP，远超过BLUP方法。禽类及水产类动物也在尝试全基因组选择育种策略[10-11]。因此，针对全基因组选择模型在动物遗传育种领域的研究进行综述，归纳全基因组选择应用的重点、难点、突破点，为全基因组选择技术研究提供新的视角尤为重要。本文综述了最新全基因组选择的统计模型及该技术的应用现状，同时对全基因组选择中存在的问题进行讨论，并对全基因组选择的未来进行展望。
1 全基因组选择概况及其原理动物重要经济性状多为复杂的数量性状，常规育种手段主要利用性状记录值、基于系谱计算的个体间亲缘关系，通过最佳线性无偏估计(best linear unbiased predication，BLUP)来估计各性状个体育种值(estimated breeding values, EBVs)[12]，通过加权获得个体综合选择指数，根据综合选择指数高低进行选留。随着微卫星(microsatellite)、单核苷酸多态性(single nucleotide polymorphisms，SNP)等遗传标记的开发和应用，将部分功能验证的候选标记联合BLUP计算育种值，即标记辅助选择(marker assisted selection, MAS)[13]，这样不仅可以提高育种值估计的准确性，而且可以在能够获得DNA时进行早期选择，缩短世代间隔，加快遗传进展。由于畜禽的重要经济性状大部分为微效多基因控制，难以找到大效应功能突变位点，并且解释的遗传变异有限，因此标记辅助选择存在一定的局限性。另外，基于系谱计算的个体亲缘关系中全同胞所有个体具有相同的育种值，而实际上，个体间性状表现存在差异，因此，系谱推断的亲缘关系存在一定缺陷。全基因组选择(genomic selection, GS)是指通过覆盖全基因组范围内的高密度标记进行育种值估计，继而进行排序、选择，简单可以理解为全基因组范围内的标记辅助选择[14]，主要方法是通过全基因组中大量的遗传标记估计出不同染色体片段或单个标记效应值，然后将个体全基因组范围内片段或标记效应值累加，获得基因组估计育种值(genomic estimated breeding value, GEBV)，其理论假设是在分布于全基因组的高密度SNP标记中，至少有一个SNP能够与影响该目标性状的数量遗传位点(quantitative trait loci, QTL)处于连锁不平衡(linkage disequilibrium, LD)状态，这样使得每个QTL的效应都可以通过SNP得到反映。相比BLUP方法，全基因组选择可以有效降低计算个体亲缘关系时孟德尔抽样误差的影响；相比MAS方法，全基因组选择模型中包括了覆盖于全基因组的标记，能更好地解释表型变异[15-16]。随着芯片和测序技术的快速发展，高通量SNP标记检测成本不断降低，使得全基因组选择技术应用于育种实践成为可能。
2 全基因组选择模型统计模型是全基因组选择的核心，极大地影响了基因组预测的准确度和效率。根据统计模型的不同，全基因组选择的模型大体可分为两大类：第一类是直接法，此方法把个体作为随机效应，参考群体和预测群体遗传信息构建的亲缘关系矩阵作为方差协方差矩阵，通过迭代法估计方差组分，进而求解混合模型获取待预测个体的估计育种值；第二类是间接法，此方法则首先在参考群中估计标记效应，然后结合预测群的基因型信息将标记效应累加，获得预测群的个体估计育种值[17-18]。
2.1 直接法模型直接法模型：
| $ y = Xb + Z\mu + e $ |
其中，y为性状向量；b为固定效应；μ为随机效应，且服从均值为0，方差为Gσa2的正态分布，可记作μ~N(0, Gσa2)，σa2为遗传方差；X和Z分别为b和μ的关联矩阵；e为残差效应，服从正态分布N(0, Iσe2)，G为个体间的亲缘关系矩阵。G的构造方法较多，目前Vanraden[19]提出的方法应用最广，其计算公式：
| $ G = M'M/2\sum\nolimits_{i = 1}^m {{p_i}\left( {1 - {p_i}} \right)} $ |
其中，M为m×n标准化的基因型矩阵，m为标记个数，n为分型个体数；pi为第i个位点最小等位基因频率，通过约束最大似然法(restricted maximum likelihood, REML)等估计得到方差组分，计算个体育种值：
| $ \mu = {\left( {Z'Z + {G^{ - 1}} \cdot \sigma _e^2/\sigma _\alpha ^2} \right)^{ - 1}}Z'\left( {y - Xb} \right) $ |
直接法与传统BLUP方法原理一致，只是用基于标记计算的G矩阵代替了基于系谱计算的A矩阵，为了与传统BLUP相区分，该方法被称之为GBLUP。GBLUP的实现只需要构建G矩阵，因此，计算速度快，相比于A阵，G阵能真实反映个体间遗传信息的差异，降低了孟德尔抽样造成的偏差，因此GBLUP相比于BLUP具有更高的准确性[20-21]。由于GBLUP计算高效，使之成为实际育种中应用最广的方法，但所有标记对G矩阵具有等同的贡献，且不同性状利用相同的G矩阵，而实际上不同性状遗传机制不同，复杂程度不同，因此GBLUP方法具有很大的改进空间[22-23]。
2.2 基于直接法的模型改进不同性状具有不同的遗传构建(genetic architecture)，直接法利用相同的关系矩阵，并且在构建矩阵的过程中，给予所有标记等同权重，针对直接法的这些不足，目前对其改进方法主要分为两大类：一类仍然在GBLUP模型中设置一个随机效应(不包含残差效应)，但是在构建G矩阵过程中，对不同标记给予权重，称之为性状特异关系矩阵；另一类将标记分类，按照不同染色体区域、与性状关联程度大小等条件，将标记分为不同的组别，在模型中设置两个或多个随机效应。
2.2.1 单随机效应Zhang等[24]提出对基因组关系矩阵进行优化，构建性状特异关系矩阵，称之为TABLUP(trait-specific relationship matrix BLUP)。TABLUP利用IBS (identity by state)构建关系矩阵，并给每一个标记分配不同权重，权重分两种：一种是直接用BayesB方法估计出来的标记方差，用TAB表示；另一种是利用RRBLUP估计的标记方差得到2pi(1－pi)gi2，用TAP表示。结果显示，性状特异关系矩阵能显著增加GBLUP/RRBLUP准确性，略差于BayesB，并且TAB优于TAP。然而BayesB方法计算复杂度较高，计算速度慢，在BayesB基础上优化GBLUP意义不大。2014年Zhang等[25]将标记权重与全基因组关联分析结果联合起来，利用公共数据库中全世界研究报道该性状候选位点的出现次数作为标记的权重，也在一定程度上提高了模型的预测准确性。
受TABLUP的启发，Zhang等[26]再次对关系矩阵构建方法进行优化，提出了GBLUP|GA方法。BLUP|GA关系矩阵构造方法：
| $ \begin{array}{*{20}{l}} {S = M_1'D{M_1}/2\sum\nolimits_{i = 1}^{{M_1}} {{p_i}\left( {1 - {p_i}} \right)} }\\ {T = \omega S + \left( {1 - \omega } \right)G} \end{array} $ |
其中，M1为通过RRBLUP计算效应值之后，按照效应值绝对值从大到小排序前n1%的标记，同时考虑LD的原因，在这些标记上下游附加n2个标记，D为对角矩阵，对角元素为标记对应的效应值，并标准化到均值为1，将构建的S矩阵与所有标记构建的G进行权重相加得到T，再进行预测。参数n1、n2、ω在参考群内采用交叉验证得到最佳组合。结果表明，BLUP|GA相比于GBLUP准确性显著提高，而且在部分性状的预测准确度达到或超越了BayesB方法，使得性状特异矩阵成为了研究热点。
任何与目标性状无关的标记被用于估计亲缘关系矩阵都会稀释QTL的作用，因此，有研究者提出，将大效应标记放入模型中作为固定效应[27]，解释主要的遗传方差，剩余遗传方差由随机效应部分获取，该策略可显著提高GBLUP准确性，但如何选择大效应标记则具有极大的难度。
随着多组学技术的发展，大量组学数据并未得到充分利用，整合不同物种组学数据进行个体育种值估计也成为全基因组选择的新热点[28-29]。Hao等[30]利用组织特异SNP功能注释信息构建关系矩阵，不仅可以提高全基因组关联分析统计效率，并且可用于全基因组选择育种值估计。Gao等[31]提出了一种新的思路，将标记与基因通路信息整合构建个体亲缘关系矩阵，通过与GBLUP对比发现，在拟南芥数据上准确性可提高0.56%~26.67%，在水稻数据上准确性可提高1.62%~16.53%。
GBLUP及其优化方法只利用基因组信息构建矩阵，这意味着被预测的所有个体具有相同的高密度标记，未分型的个体无法被预测。Aguilar等[32]以及Christensen和Lund[33]几乎同时提出一种新的构造亲缘关系矩阵的方法，同时使用系谱关系矩阵A和基因组关系矩阵G，获得新的矩阵H，该方法称之为SSBLUP (single-step BLUP)，其构造方法：
| $ \begin{array}{l} H = \\ \left[ \begin{array}{l} {A_{11}} + {A_{12}}A_{22}^{ - 1}\left( {{G_w} - {A_{22}}} \right)A_{22}^{ - 1}{A_{21}}\;{A_{12}}A_{22}^{ - 1}{G_w}\\ \;\;\;\;\;\;\;\;\;\;{G_w}A_{22}^{ - 1}{A_{21}}\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;{G_w}\;\;\; \end{array} \right] \end{array} $ |
其中，A11为未用高密度分型的个体间血缘关系矩阵；A22为用高密度分型的个体间血缘关系矩阵；Gw=(1－α)G*+αA22，G*为校正之后的基因组关系矩阵，是一种常用校正方法[34]：
Avg(diag(G))·a+b=Avg(diag(A22))
Avg(offdiag(G))·a+b=Avg(offdiag(A22))
其中，Avg为平均值；diag为矩阵对角元素；offdiag为矩阵非对角元素；a、b为校正系数。通过迭代求解MME的方法，可以获得所有分型以及未分型个体育种值。但在求解MME过程中涉及到对H矩阵的求逆计算，融合了G矩阵的H矩阵变得稠密，增加了其直接求逆的计算强度。Aguilar等[32]及Christensen和Lund[33]推导出H逆矩阵的构造方法：
| $ \begin{array}{l} {H^{ - 1}} = {A^{ - 1}} + \\ \left[ \begin{array}{l} 0\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;0\\ 0\;\;\;\;\;\;\;\;\;\;\left[ {\left( {1 - \alpha } \right)G* + \alpha {A_{22}}} \right]{\;^{ - 1}} - \omega A_{22}^{ - 1}\;\; \end{array} \right] \end{array} $ |
其中，ω一般为1，α为血缘关系权重，一般为0.05。构造H逆矩阵中A的逆矩阵可根据Henderson方法直接由系谱导出，不需要对A直接求逆。因此，只需要对G*进行求逆，而G*的矩阵维度远远小于H矩阵，因此直接求逆可节省大量时间。
SSBLUP不仅可以估计被基因分型个体的育种值，而且可以估计未进行基因分型个体的育种值。对于分型个体，由于囊括了更多有系谱记录个体的表型数据信息，理论上相比于GBLUP具有更高的准确性；对于未分型个体，由于亲缘关系被分型个体的基因组亲缘关系矩阵校正，理论上比未校正的BLUP有更高的准确性。研究表明，相比于GBLUP，SSBLUP对于猪不同性状预测准确性具有较大提升。Christensen等[34-35]分析了猪产仔数、仔猪成活率、日增重及饲料转化效率，对比不同方法发现SSBLUP表现最好。Miar[36]对肉质性状和胴体品质使用GBLUP和SSBLUP方法对比发现，SSBLUP准确性比GBLUP分别高出17%和16%。众多研究表明，SSBLUP在猪全基因组选择上存在较大优势，是当前猪全基因组选择中最为广泛使用的方法之一。另外，在实现SSBLUP的过程中，G矩阵的计算仍然采用标准Vanraden方法，并没有针对不同遗传构建的性状特异性优化。因此，若采用Zhang等[26]策略与SSBLUP进行整合，采用性状特异矩阵T代替G矩阵，可能进一步提高SSBLUP准确性。
2.2.2 多随机效应Sarup等[37]提出将已验证功能的突变位点区域单独作为一个随机效应，剩余标记作为第二个随机效应，称之为GFBLUP (genomic feature BLUP)，这样在一个模型中同时放入2个随机效应使得候选QTL区域捕获更多遗传方差，通过对杜洛克群体的日增重、饲料转化效率、眼肌面积3个性状分析发现，GFBLUP能显著增加预测准确性。GFBLUP问题在于大部分性状已验证的QTL区域较少，而且模型只限制2个随机效应。Speed和Balding[38]首次将GS与GWAS结合起来，根据GWAS结果信息，结合似然比检验(likelihood ratio test，LRT)，挑选多个(M)区域作为随机效应，不限制模型中随机效应个数，称之为MultiBLUP：
| $ y = Xb + \sum\limits_{i = 1}^M {Z{\mu _i}} + e $ |
其中，μi为第i个区域随机效应值，服从分布μi~N(0, Giσai2)，Gi为利用第i个区域SNP计算的基因组关系矩阵，σai2为通过多随机效应方差组分估计的第i个区域加性遗传方差。MultiBLUP相比于GFBLUP更加灵活，每一个QTL区域作为一个随机效应，避免放在一个随机效应中导致QTL区域间相互影响，尽可能使每一个QTL分配到准确的遗传方差。通过利用人数据对比发现，MultiBLUP准确性高于Bayes系列方法。Weissbrod等[39]在MultiBLUP基础上继续优化，增加不同区域间互作效应，称之为MKLMM (multi-kernel linear mixed model)，使得预测准确性进一步提高。多随机效应灵活多变，但是当群体不断增加，多随机效应的方差组分估计成为一大难题，也成为多随机效应模型受制约的关键因素。
2.3 间接法模型间接法模型：
| $ y = Xb + \sum\limits_{i = 1}^M {{Z_i}{g_i}} + e $ |
其中，y为表型向量；X为固定效应系数矩阵；b为固定效应；间接法中Z与直接法不同，Zi为第i个位点数字化基因型向量(如：0, 1, 2)；gi为第i个位点效应值；e为模型拟合残差，服从分布N~(0, Iσe2)。对于多元回归，有g=(Z′Z)－1Z′(y－Xb)，可得到标记效应值，但是由于在模型求解中，需要估计的模型效应变量数远远超过观察值的样本数，会导致多重共线性和过度参数化[40-41]，Z′Z无法求逆，故通过参数对角权重法解决，即：
| $ g = {\left( {Z\prime Z + I\cdot\sigma _e^2/\sigma _{{g_i}}^2} \right)^{ - 1}}Z\prime \left( {y - Xb} \right) $ |
其中，σgi2为第i个标记方差，直接与性状遗传构建相关。间接法重点和难点在于如何对超参的先验分布，即对gi及其方差服从的分布进行合理假设。Whittaker等[42]假设所有标记都具有效应，且来源于同一个分布，即σgi2相等，该方法被称之为RRBLUP (ridge regression BLUP)，Goddard[43]和Hayes等[44]从理论上证明了RRBLUP与GBLUP方法是等价的。然而，不同性状遗传机理不同，控制的基因数目及遗传效应大小不一样，因此，认为所有标记都具有效应，显然不合理，故假设所有标记方差不等更符合实际情况。Meuwissen等[3]提出了两种假设：一种认为所有标记都具有效应，gi服从正态分布N(0, σg2)，标记方差不等并服从逆卡方分布x－2(v, S)，自由度v和尺度参数S直接关系性状遗传构建，效应值gi的条件后验分布为t分布，即大部分标记效应较小，只有少部分大效应标记，由于参数求解过程结合Bayes理论，这种假设被称为BayesA；另一种假设不同于BayesA，认为大部分(π，百分比)标记无效应，只有少部分(1-π)标记具有效应，这小部分标记的条件后验分布采用与BayesA相同的t分布，基于该假设开发了一系列的Bayes方法，如Bayes B、C、Cpi等。各种模型的预测准确度较大程度的取决于其模型假设是否适合所预测表型的遗传构建，以下列举几种经典Bayes方法的先验假设区别[45-46]：
| $ \begin{array}{l} {g_i}|\left( {\pi ,\sigma _g^2} \right) \sim \left\{ {\begin{array}{*{20}{l}} {\sigma _g^2 = 0;\pi \sim dis{t_0}}\\ {{g_i}|\sigma _g^2 \sim dis{t_1};\sigma _g^2 \sim dis{t_2};\left( {1 - \pi } \right)} \end{array}} \right.\\ \sim \left[ \begin{array}{l} BayesA:dis{t_0} = 0\& {\rm{dis}}{{\rm{t}}_{\rm{1}}}{\rm{ = N}}({\rm{0}},\sigma _{\rm{g}}^{\rm{2}})\& {\rm{dis}}{{\rm{t}}_{\rm{2}}}{\rm{ = }}{{\rm{x}}^{{\rm{ - 2}}}}\left( {{\rm{v}},{\rm{S}}} \right)\\ {\rm{BayesB}}:{\rm{dis}}{{\rm{t}}_{\rm{0}}}{\rm{ = 0}}.{\rm{947}}\& {\rm{dis}}{{\rm{t}}_{\rm{1}}}{\rm{ = N}}({\rm{0}},\sigma _{\rm{g}}^{\rm{2}})\& {\rm{dis}}{{\rm{t}}_{\rm{2}}}{\rm{ = }}{{\rm{x}}^{{\rm{ - 2}}}}\left( {{\rm{v}},{\rm{S}}} \right)\\ {\rm{BayesC}}:{\rm{dis}}{{\rm{t}}_{\rm{0}}}{\rm{ = U}}\left( {{\rm{0}},{\rm{1}}} \right)\& {\rm{dis}}{{\rm{t}}_{\rm{1}}}{\rm{ = N}}({\rm{0}},\sigma _{\rm{g}}^{\rm{2}})\& {\rm{dis}}{{\rm{t}}_{\rm{2}}}{\rm{ = }}{{\rm{x}}^{{\rm{ - 2}}}}\left( {{\rm{v}},{\rm{S}}} \right)\\ {\rm{BayesC}}\pi :{\rm{dis}}{{\rm{t}}_{\rm{0}}}{\rm{ = U}}\left( {{\rm{0}},{\rm{1}}} \right)\& {\rm{dis}}{{\rm{t}}_{\rm{1}}}{\rm{ = N}}({\rm{0}},\sigma _{\rm{g}}^{\rm{2}})\& {\rm{dis}}{{\rm{t}}_{\rm{2}}}{\rm{ = C}}\\ {\rm{BayesD}}\pi :{\rm{dis}}{{\rm{t}}_{\rm{0}}}{\rm{ = U}}\left( {{\rm{0}},{\rm{1}}} \right)\& {\rm{dis}}{{\rm{t}}_{\rm{1}}}{\rm{ = N}}({\rm{0}},\sigma _{\rm{g}}^{\rm{2}})\& {\rm{dis}}{{\rm{t}}_{\rm{2}}}{\rm{ = }}{{\rm{x}}^{{\rm{ - 2}}}}\left( {{\rm{v}},{\rm{Gamma}}\left( {{\rm{1}},{\rm{1}}} \right)} \right)\\ {\rm{BayesLASSO}}:{\rm{dis}}{{\rm{t}}_{\rm{0}}}{\rm{ = 0}}\& {\rm{dis}}{{\rm{t}}_{\rm{1}}}{\rm{ = N}}({\rm{0}},\sigma _{\rm{g}}^{\rm{2}})\& {\rm{dis}}{{\rm{t}}_{\rm{2}}}{\rm{ = Expon}}({\lambda ^{\rm{2}}}{\rm{/2}}) \end{array} \right. \end{array} $ |
其中，dist0为无效应标记分布；dist1为标记效应值分布；dist2为标记效应方差分布。符合更复杂假设的模型，如Bayes模型，往往具有更多的待估参数，在提高预测准确度的同时带来了更大的计算量。参数求解过程一般先假设模型中待求变量的分布类型，即假设参数的先验分布，确定样本和参数的联合分布，并根据贝叶斯定理推导出参数的后验分布，构造一个蒙特卡洛马尔科夫链(markov chain monte carlo, MCMC)，采用吉普斯(gibbs)抽样或者MH(metropolis-hasting)抽样方法，设置充足的迭代次数直至收敛，达到平稳状态。MCMC往往需要数万次的迭代，每一次迭代需要重估所有标记效应值，该过程连续且不可并行，需消耗大量的计算时间，限制了其在时效性需求较强的动物育种实践中的应用。
2.4 基于间接法的模型改进Zhou等[47]结合Bayes和GBLUP思想，提出一种新的假设，认为所有标记都具有效应，大部分(π)标记效应可以通过GBLUP获得，服从正态分布N~(0, σa2)，另外少部分(1－π)标记除了σa2之外，具有额外的方差σb2无法被GBLUP捕获，这少部分标记服从正态分布N~(0, σa2+σb2)，区别于经典Bayes理论中每个标记具有独立的方差，该假设仅将标记方差分为两类，降低了假设复杂度，所有标记效应混合两部分正态分布：
| $ \pi N\left( {0,\sigma _a^2} \right) + \left( {1 - \pi } \right)N\left( {0,\sigma _a^2 + \sigma _b^2} \right) $ |
该假设称之为BSLMM(Bayesian sparse linear mixed model)。其中π也是待估参数，根据π的值不同，BSLMM具有灵活的可变性：当π等于1时，BSLMM为标准GBLUP模型；当σa等于0时，BSLMM为BVSR模型。Zhou等[47]对比了BSLMM与GBLUP、BVSR[48]、BayesCπ[49]、BayesLASSO[46]等，结果表明BSLMM表现更加优越。
Moser等[50]将标记效应假设更加复杂，对标记按照效应大小分为4个组，分别对应大、中、小及无效应4个水平，水平中具有相同的方差，水平间呈现相同的梯度比值，其标记效应分布：
| $ \begin{array}{l} {\pi _1}N\left( {0,0\sigma _g^2} \right) + {\pi _2}N\left( {0,{{10}^{ - 4}}\sigma _g^2} \right) + \\ {\pi _3}N\left( {0,{{10}^{ - 3}}\sigma _g^2} \right) + {\pi _4}N\left( {0,{{10}^{ - 2}}\sigma _g^2} \right),\\ {\pi _1} + {\pi _2} + {\pi _3} + {\pi _4} = 1, \end{array} $ |
该假设称之为BayesR。随着每个组内标记个数的分布不同，BayesR能够适用于由简单到复杂的所有性状，因此具有更高的灵活性和广泛性。大量的数据应用表明，相比其他间接法模型，BayesR往往具有较高的准确性[51-52]。
Zeng和Zhou[53]认为，从经典的BayesB模型到BSLMM模型，再到BayesR模型，虽然在不断的优化升级，但是每一个方法假设都是强制将标记分为有限的几个类别。对于一个未知性状，无法获知其遗传机制的情况下，固然无法知晓标记的分类数，因此提出了一种新的理论，该理论结合了狄利克雷过程(Latent Dirichlet Process)，利用狄利克雷回归分析将标记自动分为若干类，类别数目无限制，根据性状复杂程度自动优化分类的数目，该方法称之为DPR(dirichlet process regression)，假设：
| $ \sigma _g^2 \sim G,G \sim DP\left( {H,\lambda } \right) $ |
其中，H是一个基分布，可以看做G的期望；λ是系数，可以看做G的方差矩阵的“倒数”。Zeng和Zhou[53]对比DPR与其他Bayes方法发现，DPR相比于BayesR并无太大优势，在部分性状上甚至会略逊于BayesR。
Zeng等[54]认为，全基因组范围内的标记由于LD的存在使得标记间并不独立，而目前所有Bayes模型前提假设标记独立，因此其在BayesB模型基础上，将标记单独拟合升级为窗口式拟合，增加参数控制每一个标记在窗口中的贡献大小，同时考虑不同标记MAF(minor allele frequency)影响，称之为BayesN。通过基于奶牛群体的模拟数据发现，对于50K密度，BayesN相比BayesB可提高2%的准确性；对于600K密度，BayesN相比BayesB可提高7%的准确性。
间接法的假设更加复杂，但是更加符合性状的遗传构建，对于性状的遗传解析具有很好的理论研究价值，缺陷是计算速度较慢，难以用于育种实践。
3 全基因组选择不同模型预测准确性的比较不同模型假设存在差异，计算时间及准确度也不同，为了定量比较不同模型的优劣，展示模型的前提假设优化所带来的预测准确性提高，本研究利用第16届QTL-MAS Workshop公布的3个模拟性状，对文章中讨论的部分方法进行比较。该数据共包含4 100个个体，其中4 000(3 000个有表型信息，1 000个无表型信息)个个体具有基因型，因此需要预测的个体为1 000个具有基因型的个体及100个无基因型个体。结果如表 1所示，其中BLUP、GBLUP、Fixed GBLUP、Weighted GBLUP、SSBLUP为自编程序实现，Fixed GBLUP中作为协变量的QTNs为全基因组关联分析后超过显著水平的标记，Weighted GBLUP中利用性状特异基因组关系矩阵，由全基因组关联分析P值权重获得，权重方式与BLUP|GA类似，SSBLUP中亲缘关系矩阵权重为0.05。MultiBLUP模型由LDAK[38]软件实现，窗口大小设置为550 kb，关联显著水平为1×10-12，似然比检验显著水平为0.01。BayesB、BayesC、BayesLASSO由R软件中BGLR包实现，总迭代次数为2万。BSLMM由GEMMA[55]软件实现，BayesR及DPR[53]由对应软件实现，迭代次数为软件默认设置。
从表 1可以看出，利用系谱信息的BLUP模型准确性明显低于利用基因组信息的模型；间接法模型准确性优于GBLUP模型，计算时间更长，但与基于GBLUP改进的模型准确性基本等同；虽然间接法模型不断改进，但准确性并未明显提升；将大效应标记作为固定效应的Fixed GBLUP可提高GBLUP准确性，但是，由于加入的大效应标记较少，只能解释部分遗传变异，提升效果没有利用权重基因组关系矩阵的Weighted GBLUP高；对于分型个体而言SSBLUP相对于GBLUP并没有明显优势，而对于未分型的个体能够大大提高预测准确性。
|
表 1 不同模型在第16届QTL-MAS Workshop模拟数据的预测结果比较 Table 1 The comparison of predicted results of 16th QTL-MAS Workshop simulation data using different models |
性状遗传构建复杂多样，目前还没有一种模型能广泛适用于所有性状[56-57]。随着全基因组选择统计模型的不断改进优化，模型的稳定性及准确性不断提高，但是依然面临两个重要的挑战，即计算准确性和计算效率；直接法(GBLUP为代表)计算效率较高，但是计算准确性略差于间接法(BayesB为代表)，虽然学者对直接法进行了改进，但是由于改进的策略中人为设定参数较多，因此模型的预测准确性受主观因素影响较大；间接法计算准确性较高，但是由于参数求解过程中计算量庞大，且无法实现并行运算，而育种讲求时效性，所以难以高效指导育种实践；因此，如何优化模型，尽可能减少人为设定参数，与机器学习方法有效结合，并融入高效可并行运算，既能保证较高准确性的同时，大大提升计算效率，是未来全基因组选择模型优化的方向。
全基因组选择领域已经有大量的研究和应用成果，然而理论和应用方面面临着诸多挑战：1)全基因组选择主要考虑加性效应，对于显性效应及互作效应等未纳入到育种值估计模型中[58]；2)全基因组选择目前主要在品种内进行，品种间由于遗传背景不同，跨品种预测准确性难以保证；3)同品种间亲缘关系太远的个体育种值预测效果也不理想，如不同育种公司间由于育种策略不同，选择方向差异，导致同品种间遗传背景也不同，难以实现跨公司预测；4)全基因组选择只用到基因组信息，大量的多组学研究结果利用不够充分，如何将多组学信息进行整合，通过整合组学提高选择准确度也是目前待解决的问题；5)随着全基因组选择的逐渐应用，分型个体数目越来越大，相比传统BLUP的稀疏矩阵，利用基因组信息计算的稠密矩阵给混合模型参数估计及模型求解带来了巨大的挑战，通过数学或着计算机手段简化计算复杂度，才能更高效利用庞大的基因组数据甚至其他各组学数据；6)个体分型主要是芯片技术，如猪illumina 60K SNP芯片等，芯片分型具有良好的稳定性，但由于密度不足，使得全基因组选择对LD的依赖性强，通过测序手段可以得到较高密度SNP标记从而减少对LD的依赖，同时测序方法可以捕获不同品种间所有遗传变异，可能实现跨品种预测，并且测序能够得到更丰富的遗传信息，如CNV等，对于亲缘关系较近的群体，可以通过填充技术将芯片个体标记密度填充到测序水平。因此，测序技术的应用将成为全基因组选择新时代的转折点。尽管测序技术对于全基因组选择具有众多好处，但也存在一些问题，测序技术已经经历了3代技术革新，检测质量及完整性越来越高，高质量的测序结果需要更高的测序深度，意味着测序成本更昂贵，并且测序数据庞大，主流的分析软件处理速度较慢，使用复杂繁琐，对于计算资源的配置需求较高，因此如何快速、有效地储存、处理及分析数据是测序技术应用于全基因组育种的重要挑战，另外，测序只能检测参考基因组中已知的序列和基因信息，对于未知的基因序列和基因还不能进一步深入研究。当然，随着测序方法和芯片技术的不断成熟，未来个体分型费用将不断降低，分型准确性不断提高，全基因组选择将逐步替代传统育种方法，为动物育种改良带来一次新的技术革命。
| [1] | MEUWISSEN T, HAYES B, GODDARD M. Genomic selection:A paradigm shift in animal breeding[J]. Anim Front, 2016, 6(1): 6–14. |
| [2] | VISSCHER P M, WRAY N R, ZHANG Q, et al. 10 years of GWAS discovery:biology, function, and translation[J]. Am J Hum Genet, 2017, 101(1): 5–22. DOI: 10.1016/j.ajhg.2017.06.005 |
| [3] | MEUWISSEN T H E, HAYES B J, GODDARD M E. Prediction of total genetic value using genome-wide dense marker maps[J]. Genetics, 2001, 157(4): 1819–1829. |
| [4] | HAYES B J, BOWMAN P J, CHAMBERLAIN A J, et al. Invited review:Genomic selection in dairy cattle:progress and challenges[J]. J Dairy Sci, 2009, 92: 433–443. DOI: 10.3168/jds.2008-1646 |
| [5] | WELLER J I, EZRA E, RON M. Invited review:A perspective on the future of genomic selection in dairy cattle[J]. J Dairy Sci, 2017, 100(11): 8633–8644. DOI: 10.3168/jds.2017-12879 |
| [6] | GARCÍA-RUIZ A, COLE J B, VANRADEN P M, et al. Changes in genetic selection differentials and generation intervals in US Holstein dairy cattle as a result of genomic selection[J]. Proc Natl Acad Sci U S A, 2016, 113(28): E3995–E4004. DOI: 10.1073/pnas.1519061113 |
| [7] |
王晨, 秦珂, 薛明, 等. 全基因组选择在猪育种中的应用[J]. 畜牧兽医学报, 2016, 47(1): 1–9.
WANG C, QIN K, XUE M, et al. Application of genomic selection in swine breeding[J]. Acta Veterinaria et Zootechnica Sinica, 2016, 47(1): 1–9. (in Chinese) |
| [8] | SAMORÈ AB, FONTANESI L. Genomic selection in pigs:state of the art and perspectives[J]. Ital J Anim Sci, 2016, 15(2): 211–232. DOI: 10.1080/1828051X.2016.1172034 |
| [9] | MUCHA S, MRODE R, MACLAREN-LEE I, et al. Estimation of genomic breeding values for milk yield in UK dairy goats[J]. J Dairy Sci, 2015, 98(11): 8201–8208. DOI: 10.3168/jds.2015-9682 |
| [10] | WOLC A, ZHAO H H, ARANGO J, et al. Response and inbreeding from a genomic selection experiment in layer chickens[J]. Genet Sel Evol, 2015, 47: 59. DOI: 10.1186/s12711-015-0133-5 |
| [11] | LÍPEZ M E, NEIRA R, YÁÑEZ J M. Applications in the search for genomic selection signatures in fish[J]. Front Genet, 2015, 5: 458. |
| [12] | HENDERSON C R. Best linear unbiased estimation and prediction under a selection model[J]. Biometrics, 1975, 31(2): 423–447. |
| [13] | LANDE R, THOMPSON R. Efficiency of marker-assisted selection in the improvement of quantitative traits[J]. Genetics, 1990, 124(3): 743–756. |
| [14] | MEUWISSEN T. Genomic selection:marker assisted selection on a genome wide scale[J]. J Anim Breed Genet, 2007, 124(6): 321–322. DOI: 10.1111/j.1439-0388.2007.00708.x |
| [15] | YANG J A, BENYAMIN B, MCEVOY B P, et al. Common SNPs explain a large proportion of the heritability for human height[J]. Nat Genet, 2010, 42(7): 565–569. DOI: 10.1038/ng.608 |
| [16] | GODDARD M E, HAYES B J, MEUWISSEN T H E. Using the genomic relationship matrix to predict the accuracy of genomic selection[J]. J Anim Breed Genet, 2011, 128(6): 409–421. DOI: 10.1111/jbg.2011.128.issue-6 |
| [17] | ZHANG Z, ZHANG Q, DING X D. Advances in genomic selection in domestic animals[J]. Chin Sci Bull, 2011, 56(25): 2655–2663. DOI: 10.1007/s11434-011-4632-7 |
| [18] | MISZTAL I, LEGARRA A. Invited review:efficient computation strategies in genomic selection[J]. Animal, 2017, 11(5): 731–736. |
| [19] | VANRADEN P M. Efficient methods to compute genomic predictions[J]. J Dairy Sci, 2008, 91(11): 4414–4423. DOI: 10.3168/jds.2007-0980 |
| [20] | HBIERA D, FERNANDO R L, DEKKERS J C M. The impact of genetic relationship information on genome-assisted breeding values[J]. Genetics, 2007, 177(4): 2389–2397. |
| [21] | MUIR W M. Comparison of genomic and traditional BLUP-estimated breeding value accuracy and selection response under alternative trait and genomic parameters[J]. J Anim Breed Genet, 2007, 124(6): 342–355. DOI: 10.1111/j.1439-0388.2007.00700.x |
| [22] | GODDARD M E, KEMPER K E, MACLEOD I M, et al. Genetics of complex traits:prediction of phenotype, identification of causal polymorphisms and genetic architecture[J]. Proc Biol Sci, 2016, 283: 20160569. DOI: 10.1098/rspb.2016.0569 |
| [23] | YANG J, ZENG J, GODDARD M E, et al. Concepts, estimation and interpretation of SNP-based heritability[J]. Nat Genet, 2017, 49(9): 1304–1310. DOI: 10.1038/ng.3941 |
| [24] | ZHANG Z, LIU J F, DING X D, et al. Best linear unbiased prediction of genomic breeding values using a trait-specific marker-derived relationship matrix[J]. PLoS One, 2010, 5(9): e12648. DOI: 10.1371/journal.pone.0012648 |
| [25] | ZHANG Z, OBER U, ERBE M, et al. Improving the accuracy of whole genome prediction for complex traits using the results of genome wide association studies[J]. PLoS One, 2014, 9(3): e93017. DOI: 10.1371/journal.pone.0093017 |
| [26] | ZHANG Z, ERBE M, HE J L, et al. Accuracy of whole-genome prediction using a genetic architecture-enhanced variance-covariance matrix[J]. G3(Bethesda), 2015, 5(4): 615–627. |
| [27] | MOORE J K, MANMATHAN H K, ANDERSON V A, et al. Improving genomic prediction for pre-harvest sprouting tolerance in wheat by weighting large-effect quantitative trait loci[J]. Crop Sci, 2017, 57(3): 1315–1324. DOI: 10.2135/cropsci2016.06.0453 |
| [28] | SCHRAG T A, WESTHUES M, SCHIPPRACK W, et al. Beyond genomic prediction:combining different types of omics data can improve prediction of hybrid performance in maize[J]. Genetics, 2018, 208(4): 1373–1385. DOI: 10.1534/genetics.117.300374 |
| [29] | AN N R, LEE S S, PARK J E, et al. Current status of genomic prediction using Multi-omics data in livestock[J]. J Biomed Translat Res, 2017, 18(4): 151–156. DOI: 10.12729/jbtr. |
| [30] | HAO X J, ZENG P, ZHANG S J, et al. Identifying and exploiting trait-relevant tissues with multiple functional annotations in genome-wide association studies[J]. PLoS Genet, 2018, 14(1): e1007186. DOI: 10.1371/journal.pgen.1007186 |
| [31] | GAO N, TENG J Y, YE S P, et al. Genomic prediction of complex phenotypes using genic similarity based relatedness matrix[J]. Front Genet, 2018, 9: 364. DOI: 10.3389/fgene.2018.00364 |
| [32] | AGUILAR I, MISZTAL I, JOHNSON D L, et al. Hot topic:A unified approach to utilize phenotypic, full pedigree, and genomic information for genetic evaluation of Holstein final score[J]. J Dairy Sci, 2010, 93(2): 743–752. DOI: 10.3168/jds.2009-2730 |
| [33] | CHRISTENSEN O F, LUND M S. Genomic prediction when some animals are not genotyped[J]. Genet Sel Evol, 2010, 42: 2. DOI: 10.1186/1297-9686-42-2 |
| [34] | CHRISTENSEN O F, MADSEN P, NIELSEN B, et al. Single-step methods for genomic evaluation in pigs[J]. Animal, 2012, 6(10): 1565–1571. DOI: 10.1017/S1751731112000742 |
| [35] | GUO X, CHRISTENSEN O F, OSTERSEN T, et al. Improving genetic evaluation of litter size and piglet mortality for both genotyped and nongenotyped individuals using a single-step method[J]. J Anim Sci, 2015, 93(2): 503–512. DOI: 10.2527/jas.2014-8331 |
| [36] | MIAR Y.Genomic selection for pork quality and carcass traits in both cross-and pure-bred populations[D]. Edmonton, AB, Canada: University of Alberta, 2015. |
| [37] | SARUP P, JENSEN J, OSTERSEN T, et al. Increased prediction accuracy using a genomic feature model including prior information on quantitative trait locus regions in purebred Danish Duroc pigs[J]. BMC Genet, 2016, 17: 11. |
| [38] | SPEED D, BALDING D J. MultiBLUP:improved SNP-based prediction for complex traits[J]. Genome Res, 2014, 24(9): 1550–1557. DOI: 10.1101/gr.169375.113 |
| [39] | WEISSBROD O, GEIGER D, ROSSET S. Multikernel linear mixed models for complex phenotype prediction[J]. Genome Res, 2016, 26(7): 969–979. DOI: 10.1101/gr.201996.115 |
| [40] |
王重龙, 丁向东, 刘剑锋, 等. 基因组育种值估计的贝叶斯方法[J]. 遗传, 2014, 36(2): 111–118.
WANG C L, DING X D, LIU J F, et al. Bayesian methods for genomic breeding value estimation[J]. Hereditas, 2014, 36(2): 111–118. (in Chinese) |
| [41] | WANG C L, DING X D, WANG J Y, et al. Bayesian methods for estimating GEBVs of threshold traits[J]. Heredity, 2013, 110(3): 213–219. DOI: 10.1038/hdy.2012.65 |
| [42] | WHITTAKER J C, THOMPSON R, DENHAM M C. Marker-assisted selection using ridge regression[J]. Genet Res, 2000, 75(2): 249–252. |
| [43] | GODDARD M. Genomic selection:prediction of accuracy and maximisation of long term response[J]. Genetica, 2009, 136(2): 245–257. DOI: 10.1007/s10709-008-9308-0 |
| [44] | HAYES B J, VISSCHER P M, GODDARD M E. Increased accuracy of artificial selection by using the realized relationship matrix[J]. Genet Res, 2009, 91(1): 47–60. |
| [45] | HABIER D, FERNANDO R L, KIZILKAYA K, et al. Extension of the bayesian alphabet for genomic selection[J]. BMC Bioinformatics, 2011, 12: 186. DOI: 10.1186/1471-2105-12-186 |
| [46] | YI N J, XU S Z. Bayesian LASSO for quantitative trait loci mapping[J]. Genetics, 2008, 179(2): 1045–1055. DOI: 10.1534/genetics.107.085589 |
| [47] | ZHOU X, CARBONETTO P, STEPHENS M. Polygenic modeling with bayesian sparse linear mixed models[J]. PLoS Genet, 2013, 9(2): e1003264. DOI: 10.1371/journal.pgen.1003264 |
| [48] | LOGSDON B A, HOFFMAN G E, MEZEY J G. A variational Bayes algorithm for fast and accurate multiple locus genome-wide association analysis[J]. BMC Bioinformatics, 2010, 11: 58. DOI: 10.1186/1471-2105-11-58 |
| [49] | HABIER D, FERNANDO R L, KIZILKAYA K, et al. Extension of the bayesian alphabet for genomic selection[J]. BMC Bioinformatics, 2011, 12: 186. DOI: 10.1186/1471-2105-12-186 |
| [50] | MOSER G, LEE S H, HAYES B J, et al. Simultaneous discovery, estimation and prediction analysis of complex traits using a bayesian mixture model[J]. PLoS Genet, 2015, 11(4): e1004969. DOI: 10.1371/journal.pgen.1004969 |
| [51] | FIKERE M, BARBULESCU D M, MALMBERG M M, et al. Genomic prediction using prior quantitative trait loci information reveals a large reservoir of underutilised blackleg resistance in diverse canola (Brassica napus L.) lines[J]. Plant Genome, 2018, 11: 170100. |
| [52] | HAYES B J, CORBET N J, ALLEN J M, et al. Towards multi-breed genomic evaluations for female fertility of tropical beef cattle[J]. J Anim Sci, 2018. DOI: 10.1093/jas/sky417 |
| [53] | ZENG P, ZHOU X. Non-parametric genetic prediction of complex traits with latent dirichlet process regression models[J]. Nat Commun, 2017, 8: 456. DOI: 10.1038/s41467-017-00470-2 |
| [54] | ZENG J, GARRICK D, DEKKERS J, et al. A nested mixture model for genomic prediction using whole-genome SNP genotypes[J]. PLoS One, 2018, 13(3): e0194683. DOI: 10.1371/journal.pone.0194683 |
| [55] | ZHOU X, STEPHENS M. Efficient multivariate linear mixed model algorithms for genome-wide association studies[J]. Nat Methods, 2014, 11(4): 407–409. DOI: 10.1038/nmeth.2848 |
| [56] | WRAY N R, WIJMENGA C, SULLIVAN P F, et al. Common disease is more complex than implied by the core gene omnigenic model[J]. Cell, 2018, 173(7): 1573–1580. DOI: 10.1016/j.cell.2018.05.051 |
| [57] |
张巧霞, 张玲妮, 刘飞, 等. 基于GBLUP与惩罚类回归方法的猪血液性状基因组选择研究[J]. 畜牧兽医学报, 2017, 48(12): 2258–2267.
ZHANG Q X, ZHANG L N, LIU F, et al. A study of genomic selection on porcine hematological traits using GBLUP and penalized regression methods[J]. Acta Veterinaria et Zootechnica Sinica, 2017, 48(12): 2258–2267. (in Chinese) |
| [58] |
王延晖, 朱波, 李俊雅. 基于显性模型的基因组选择中贝叶斯方法研究[J]. 畜牧兽医学报, 2017, 48(1): 60–67.
WANG Y H, ZHU B, LI J Y. Bayesian models including dominant effects for genomic selection[J]. Acta Veterinaria et Zootechnica Sinica, 2017, 48(1): 60–67. (in Chinese) |