結(jié)合GWAS先驗(yàn)標(biāo)記信息的肉雞RFI性狀全基因組選擇研究

杜永旺，黃 超，王一東，李 森，文 杰，陳智武，趙桂蘋，鄭麥青*

(1.中國農(nóng)業(yè)科學(xué)院北京畜牧獸醫(yī)研究所，北京 100193； 2.廣西金陵農(nóng)牧集團(tuán)有限公司，南寧 530049)

飼料成本占據(jù)肉雞養(yǎng)殖成本的70%左右，而遺傳選育是改進(jìn)飼料利用效率的主要途徑。家禽育種由最初的表型選擇、家系選擇正逐漸過渡到分子育種。隨著單核苷酸多態(tài)性(single nucleotide polymorphism, SNP)等分子標(biāo)記的開發(fā)與應(yīng)用，將部分功能驗(yàn)證的候選標(biāo)記聯(lián)合BLUP計(jì)算育種值的分子標(biāo)記輔助選擇不僅可以提高選擇準(zhǔn)確性，還可以實(shí)現(xiàn)早期選擇，縮短世代間隔。但是大多數(shù)畜禽的經(jīng)濟(jì)性狀都是受到微效多基因控制的數(shù)量性狀，獲取具有大效應(yīng)的變異位點(diǎn)較為困難，這也是分子標(biāo)記輔助選擇存在的局限性。而基因組選擇(genomic selection, GS)覆蓋了全基因組范圍內(nèi)的分子標(biāo)記，能更好的解釋遺傳變異。2006年在奶牛上的研究發(fā)現(xiàn)了GS在育種工作中的巨大潛力。GS方法主要包括以GBLUP為代表的直接法和以貝葉斯類方法為代表的間接法。間接法雖具有計(jì)算準(zhǔn)確性較高的優(yōu)勢，但是計(jì)算時(shí)間長，占用計(jì)算資源大的問題導(dǎo)致其在生產(chǎn)中的應(yīng)用受限。而GBLUP憑借計(jì)算時(shí)間短，使用簡單的優(yōu)勢廣泛應(yīng)用在豬、牛、羊、雞等畜禽的實(shí)際選育當(dāng)中。在GBLUP模型中，默認(rèn)所有的SNPs位點(diǎn)對性狀均有效應(yīng)且效應(yīng)相同，即未考慮先驗(yàn)信息。但是實(shí)際上不同的數(shù)量性狀的遺傳機(jī)制不同，復(fù)雜程度也不同，因此GBLUP方法存在一定的改進(jìn)空間。全基因組關(guān)聯(lián)分析作為挖掘畜禽經(jīng)濟(jì)性狀候選位點(diǎn)的重要手段，是基因組水平上最常用的先驗(yàn)信息獲取方式。對人類復(fù)雜性狀、畜禽重要經(jīng)濟(jì)性狀以及模擬數(shù)據(jù)的研究結(jié)果均表明，整合GWAS先驗(yàn)信息可以提高GS的準(zhǔn)確性。

金陵花雞屬于快大型黃羽肉雞，對生長速度和飼料消耗的選育要求較高，其出欄日齡為56日齡。出欄前兩周為金陵花雞的快速生長期，選育快速生長期的剩余采食量對生產(chǎn)具有重要的經(jīng)濟(jì)價(jià)值。本研究以廣西金陵花雞為試驗(yàn)素材，以42～56日齡的剩余采食量為目標(biāo)性狀，將試驗(yàn)群體隨機(jī)分為兩組，其中一組作為先驗(yàn)標(biāo)記信息發(fā)現(xiàn)群體，用于GWAS分析并篩選出最顯著的top5%、top10%、top15%和top20%的SNPs作為先驗(yàn)標(biāo)記信息；另外一組分別結(jié)合不同合集的先驗(yàn)標(biāo)記信息進(jìn)行遺傳參數(shù)估計(jì)并比較基因組育種值的預(yù)測準(zhǔn)確性，使用重復(fù)10次的五倍交叉驗(yàn)證法獲取準(zhǔn)確性，隨后兩組群體再進(jìn)行交叉驗(yàn)證。比較兩種策略的預(yù)測準(zhǔn)確性，為在金陵花雞剩余采食量性狀選育中實(shí)施基因組選擇提供理論依據(jù)。

1 材料與方法

1.1 試驗(yàn)材料

本試驗(yàn)所使用的材料是由廣西金陵農(nóng)牧集團(tuán)提供的快速型黃羽肉雞金陵花雞E系。該群體在0～42日齡采用群體飼養(yǎng)，在42～56日齡進(jìn)行剩余采食量測定，測定期間使用獨(dú)立的飲水和采食裝置，所有雞只自由采食，按照正常程序免疫和飼養(yǎng)管理。記錄42～56日齡期間的采食量和體重增長量，并計(jì)算剩余采食量。剩余采食量的計(jì)算方法為：

=+++++

本研究使用3個(gè)世代共2 510個(gè)個(gè)體，其中公雞1 648只，母雞862只，均具有完整的系譜和表型記錄。將群體隨機(jī)分為兩組，分組1和分組2各包含公雞824個(gè)，母雞431個(gè)，兩組RFI的性狀描述性統(tǒng)計(jì)見表1。

表1 剩余采食量的描述性統(tǒng)計(jì)Table 1 Descriptive statistics of RFI

1.2 基因分型與質(zhì)控填充

試驗(yàn)群體在56日齡翅下靜脈采血用于基因組DNA提取，并基于Illumina測序平臺使用雞55K SNP芯片對2 510個(gè)個(gè)體進(jìn)行基因型檢測，共得到44 561個(gè)SNPs，使用PLINK v1.9軟件對數(shù)據(jù)進(jìn)行質(zhì)量控制處理，質(zhì)控標(biāo)準(zhǔn)為：1)個(gè)體檢出率大于0.95；2)最小等位基因頻率(MAF)大于0.01；3)位點(diǎn)檢出率大于0.95，缺失的SNPs使用Beagle v5.0軟件進(jìn)行填充，最終保留2 510個(gè)個(gè)體和40 492個(gè)SNPs用于后續(xù)研究。

1.3 全基因組關(guān)聯(lián)分析

兩組試驗(yàn)群體互相作為先驗(yàn)信息發(fā)現(xiàn)群體和先驗(yàn)信息驗(yàn)證群體。使用先驗(yàn)信息發(fā)現(xiàn)群體進(jìn)行GWAS分析。所用模型為一般線性模型，具體為：

=++

其中，為所測定的性狀，即RFI，為基因型數(shù)據(jù)，這里指每一個(gè)SNP，為SNP的系數(shù)，為世代、性別等因素，為性別的系數(shù)，為殘差。將每個(gè)SNP作為固定因子進(jìn)行回歸分析，進(jìn)行顯著性檢驗(yàn)，獲取每一個(gè)SNP的值，并分別提取值較小的top5%、top10%、top15%、top20%的SNPs作為另外一組的先驗(yàn)標(biāo)記信息。

1.4 基因組選擇

獲取的先驗(yàn)標(biāo)記信息在先驗(yàn)信息驗(yàn)證群體上進(jìn)行驗(yàn)證。GS使用的GBLUP模型為：

=b+μ+

GBLUP使用的親緣關(guān)系矩陣是由VanRaden提出的，計(jì)算公式為：

其中，為×標(biāo)準(zhǔn)化的基因型矩陣，為標(biāo)記個(gè)數(shù)，為分型個(gè)體數(shù)，為第個(gè)位點(diǎn)的最小等位基因頻率。

結(jié)合GWAS先驗(yàn)信息的GS所用的親緣關(guān)系矩陣構(gòu)建方式為：使用GWAS結(jié)果中的先驗(yàn)標(biāo)記信息構(gòu)建親緣關(guān)系矩陣，剩余的位點(diǎn)構(gòu)建親緣關(guān)系矩陣，分別計(jì)算和對RFI解釋的遺傳方差，并以解釋遺傳方差的比例為權(quán)重?cái)M合出親緣關(guān)系矩陣：

=λ+(1-)

1.5 預(yù)測準(zhǔn)確性評估標(biāo)準(zhǔn)

本研究使用重復(fù)10次的5倍交叉驗(yàn)證來評估基因組育種值的預(yù)測準(zhǔn)確性，將GS分析群體隨機(jī)分為5組，其中4組合并為參考群，剩下一組作為驗(yàn)證群，將得到的基因組育種值與表型值的皮爾遜相關(guān)系數(shù)除以遺傳力開方作為預(yù)測準(zhǔn)確性，最后將重復(fù)10次得到的均值作為評價(jià)預(yù)測準(zhǔn)確性的標(biāo)準(zhǔn)。

2 結(jié) 果

2.1 GWAS結(jié)果

基于測定的RFI性狀，將3個(gè)世代2 510個(gè)個(gè)體隨機(jī)分為兩組，分組1和分組2各包含824只公雞和431只母雞，分別對兩個(gè)分組做GWAS分析，曼哈頓圖見圖1。獲取最顯著的top5%、top10%、top15%、top20%的SNPs作為另外一組的先驗(yàn)標(biāo)記信息。結(jié)果中最顯著的top5%有2 000個(gè)SNPs；top10%有4 000個(gè)SNPs；top15%有6 000個(gè)SNPs；top20%有8 000個(gè)SNPs，均散亂的分布在每條染色體上。

圖1 曼哈頓圖Fig.1 Manhattan plot

2.2 GS結(jié)果

分組1的遺傳參數(shù)見表2，其中矩陣計(jì)算的遺傳方差和環(huán)境方差分別為19.859和109.692，遺傳力為0.153。五倍交叉驗(yàn)證的預(yù)測準(zhǔn)確性結(jié)果為0.387。

表2 分組1遺傳參數(shù)Table 2 Genetic parameters of group 1

基于分組2的GWAS結(jié)果，篩選最顯著的top5%的SNPs構(gòu)建親緣關(guān)系矩陣，計(jì)算得到的遺傳方差和環(huán)境方差分別為14.297和115.721，遺傳力為0.110；剩余95%的SNPs構(gòu)建矩陣，解釋的遺傳方差和環(huán)境方差分別為19.506和110.061，遺傳力為0.151。計(jì)算權(quán)重為0.42，權(quán)重為0.58。使用結(jié)合GWAS先驗(yàn)標(biāo)記信息的親緣關(guān)系矩陣計(jì)算遺傳方差和環(huán)境方差分別為20.046和109.667，計(jì)算遺傳力為0.155，預(yù)測準(zhǔn)確性為0.401，準(zhǔn)確性提升3.62%。

篩選最顯著的top10%的SNPs構(gòu)建親緣關(guān)系矩陣，計(jì)算得到的遺傳方差和環(huán)境方差分別為16.038和113.559，遺傳力為0.124；剩余90%的SNPs構(gòu)建矩陣，解釋的遺傳方差和環(huán)境方差分別為19.366和110.247，遺傳力為0.149。計(jì)算權(quán)重為0.45，權(quán)重為0.55。使用結(jié)合GWAS先驗(yàn)標(biāo)記信息的親緣關(guān)系矩陣計(jì)算遺傳方差和環(huán)境方差分別為19.761和109.726，計(jì)算遺傳力為0.153，預(yù)測準(zhǔn)確性為0.406，準(zhǔn)確性提升4.91%。

篩選最顯著的top15%的SNPs構(gòu)建親緣關(guān)系矩陣，計(jì)算得到的遺傳方差和環(huán)境方差分別為18.173和111.394，遺傳力為0.140；剩余85%的SNPs構(gòu)建矩陣，解釋的遺傳方差和環(huán)境方差分別為19.012和110.619，遺傳力為0.147。計(jì)算權(quán)重為0.49，權(quán)重為0.51。使用結(jié)合GWAS先驗(yàn)標(biāo)記信息的親緣關(guān)系矩陣計(jì)算遺傳方差和環(huán)境方差分別為20.286和109.198，計(jì)算遺傳力為0.157，預(yù)測準(zhǔn)確性為0.407，準(zhǔn)確性提升5.17%。

篩選最顯著的top20%的SNPs構(gòu)建親緣關(guān)系矩陣，計(jì)算得到的遺傳方差和環(huán)境方差分別為18.618和110.809，遺傳力為0.144；剩余80%的SNPs構(gòu)建矩陣，解釋的遺傳方差和環(huán)境方差分別為18.630和111.021，遺傳力為0.144。計(jì)算權(quán)重為0.50，權(quán)重為0.50。使用結(jié)合GWAS先驗(yàn)標(biāo)記信息的親緣關(guān)系矩陣計(jì)算遺傳方差和環(huán)境方差分別為20.341和109.114，計(jì)算遺傳力為0.157，預(yù)測準(zhǔn)確性為0.404，準(zhǔn)確性提升4.39%。

分組2的遺傳參數(shù)見表3，其中矩陣計(jì)算的遺傳方差和環(huán)境方差分別為18.754和103.975，計(jì)算遺傳力為0.153。五倍交叉驗(yàn)證的預(yù)測準(zhǔn)確性結(jié)果為0.429。

表3 分組2遺傳參數(shù)Table 3 Genetic parameters of group 2

基于分組1的GWAS結(jié)果，篩選最顯著的top5%的SNPs構(gòu)建親緣關(guān)系矩陣，計(jì)算得到的遺傳方差和環(huán)境方差分別為10.506和112.212，遺傳力為0.086；剩余95%的SNPs構(gòu)建矩陣，解釋的遺傳方差和環(huán)境方差分別為18.694和104.087，遺傳力為0.152。計(jì)算權(quán)重為0.36，權(quán)重為0.64。使用結(jié)合GWAS先驗(yàn)標(biāo)記信息的親緣關(guān)系矩陣計(jì)算遺傳方差和環(huán)境方差分別為17.019和105.565，計(jì)算遺傳力為0.139，預(yù)測準(zhǔn)確性為0.446，準(zhǔn)確性提升3.96%。

篩選最顯著的top10%的SNPs構(gòu)建親緣關(guān)系矩陣，計(jì)算得到的遺傳方差和環(huán)境方差分別為13.400和109.125，遺傳力為0.109；剩余90%的SNPs構(gòu)建矩陣，解釋的遺傳方差和環(huán)境方差分別為18.590和104.235，遺傳力為0.151。計(jì)算權(quán)重為0.42，權(quán)重為0.58。使用結(jié)合GWAS先驗(yàn)標(biāo)記信息的親緣關(guān)系矩陣計(jì)算遺傳方差和環(huán)境方差分別為17.621和104.912，計(jì)算遺傳力為0.144，預(yù)測準(zhǔn)確性為0.448，準(zhǔn)確性提升4.43%。

篩選最顯著的top15%的SNPs構(gòu)建親緣關(guān)系矩陣，計(jì)算得到的遺傳方差和環(huán)境方差分別為14.160和108.323，遺傳力為0.116；剩余85%的SNPs構(gòu)建矩陣，解釋的遺傳方差和環(huán)境方差分別為18.640和104.216，遺傳力為0.152。計(jì)算權(quán)重為0.43，權(quán)重為0.57。使用結(jié)合GWAS先驗(yàn)標(biāo)記信息的親緣關(guān)系矩陣計(jì)算遺傳方差和環(huán)境方差分別為17.852和104.710，計(jì)算遺傳力為0.146，預(yù)測準(zhǔn)確性為0.441，準(zhǔn)確性提升2.80%。

篩選最顯著的top20%的SNPs構(gòu)建親緣關(guān)系矩陣，計(jì)算得到的遺傳方差和環(huán)境方差分別為14.670和107.813，遺傳力為0.120；剩余80%的SNPs構(gòu)建矩陣，解釋的遺傳方差和環(huán)境方差分別為18.685和104.195，遺傳力為0.152。計(jì)算權(quán)重為0.44，權(quán)重為0.56。使用結(jié)合GWAS先驗(yàn)標(biāo)記信息的親緣關(guān)系矩陣計(jì)算遺傳方差和環(huán)境方差分別為18.010和104.585，計(jì)算遺傳力為0.147，預(yù)測準(zhǔn)確性為0.438，準(zhǔn)確性提升2.10%。

3 討 論

在本研究中，基于矩陣和矩陣計(jì)算的RFI遺傳力范圍處于0.139～0.157之間，為低遺傳力性狀，與李森等和Li等的研究結(jié)果相近，而在劉天飛等的研究結(jié)果中，RFI為中高遺傳力，可能是由于本群體和其他研究群體的品種、飼養(yǎng)方式以及遺傳背景不同導(dǎo)致。與傳統(tǒng)的育種方法相比，基因組選擇能夠提高低遺傳力性狀的選育效果，在金陵花雞RFI性狀中應(yīng)用基因組選擇的選育效果較好。

在GBLUP模型中，默認(rèn)所有的SNPs具有相同的效應(yīng)且服從相同的分布，這顯然是不合理的。有學(xué)者使用整合先驗(yàn)信息的GS模型進(jìn)行了不同的改進(jìn)，如基因組特征BLUP、BayesRC等。Zhang等提出對基因組關(guān)系矩陣進(jìn)行優(yōu)化，構(gòu)建性狀特異關(guān)系矩陣，即TABLUP，其基于間接法Bayes及嶺回歸最佳線性無偏預(yù)測方法給所有標(biāo)記分配不同的權(quán)重，相較于GBLUP有顯著的提升，但是BayesB計(jì)算時(shí)間較長、復(fù)雜程度高，預(yù)測準(zhǔn)確性較高，基于BayesB提升GBLUP的準(zhǔn)確性意義不大。2015年，Zhang等基于rrBLUP計(jì)算每個(gè)標(biāo)記的效應(yīng)值，并按照絕對值大小進(jìn)行排序，分別構(gòu)建兩個(gè)不同的矩陣并相加，顯著提升了預(yù)測準(zhǔn)確性。同時(shí)，使用全基因組關(guān)聯(lián)分析對標(biāo)記分配權(quán)重也能提高預(yù)測準(zhǔn)確性。而本研究基于GWAS分析，將SNP位點(diǎn)按照值大小進(jìn)行排序，篩選對表型影響顯著的SNPs，并將其作為先驗(yàn)信息整合至GBLUP模型中。

RFI作為受到微效多基因調(diào)控的低遺傳力性狀，在GWAS分析中沒有定位到顯著的SNPs。本研究從最顯著的top5%的SNPs到top20%的SNPs中依次篩選并作為先驗(yàn)信息。研究結(jié)果中，隨著top合集的擴(kuò)大，矩陣解釋的遺傳方差逐漸擴(kuò)大；在將top20%的SNPs作為先驗(yàn)信息的策略中，矩陣解釋的遺傳方差逐漸接近矩陣解釋的遺傳方差。這與Romé等認(rèn)為20%的SNPs足以解釋群體的連鎖不平衡結(jié)果相近。

在本研究中，結(jié)合GWAS先驗(yàn)標(biāo)記信息的GS相較于GBLUP提升了2.10%到5.17%。隨著先驗(yàn)標(biāo)記數(shù)量的增加，GWAS捕獲到更多與RFI相關(guān)的SNPs位點(diǎn)，先驗(yàn)標(biāo)記信息解釋的遺傳方差和遺傳力也在不斷增大，構(gòu)建的矩陣權(quán)重也在逐漸增加。GS的預(yù)測準(zhǔn)確性逐漸上升。當(dāng)捕獲到所有與RFI性狀有關(guān)的位點(diǎn)時(shí)，再增加SNPs數(shù)量不會(huì)擴(kuò)大遺傳方差和遺傳力，但是GS預(yù)測準(zhǔn)確性會(huì)隨著無關(guān)位點(diǎn)的干擾呈現(xiàn)下降趨勢。整體來看，隨著先驗(yàn)標(biāo)記數(shù)量的增加，解釋的遺傳方差和遺傳力逐漸增加，而預(yù)測準(zhǔn)確性則呈現(xiàn)先上升再逐漸下降的趨勢。在分組1中，先驗(yàn)SNP標(biāo)記為top15%時(shí)取得最高預(yù)測準(zhǔn)確性，提升5.17%；分組2中是在top10%取得最大提升，為4.43%。Gao等在黃羽肉雞群體中，通過GWAS分析構(gòu)建單倍型并整合基因注釋信息到GS模型中，將RFI的預(yù)測準(zhǔn)確性從0.464提升到了0.468，提升了0.86%。Zhang等利用GWAS檢測到的QTL，構(gòu)建性狀特異性親緣關(guān)系矩陣，即BLUP|GA，并在奶牛和水稻的真實(shí)數(shù)據(jù)中證明BLUP|GA提升了GS的預(yù)測準(zhǔn)確性。袁澤湖在綿羊肉用性狀中整合GWAS先驗(yàn)信息的GS準(zhǔn)確性比基于芯片數(shù)據(jù)的GS準(zhǔn)確性更高。Schrag等在玉米上的研究表明，結(jié)合基因組、轉(zhuǎn)錄組和代謝組等其他組學(xué)的先驗(yàn)信息也可以提高GS的預(yù)測準(zhǔn)確性。目前，基因組選擇已經(jīng)廣泛應(yīng)用在實(shí)際生產(chǎn)中，對模型的改進(jìn)也在不斷進(jìn)行，但是不同群體的選育目標(biāo)和遺傳背景均不相同。因此，在生產(chǎn)中應(yīng)用基因組選擇應(yīng)該選擇使用整合本群體生物學(xué)先驗(yàn)信息的GS模型并調(diào)整相關(guān)參數(shù)。

4 結(jié) 論

本研究以金陵花雞為試驗(yàn)素材，基于GWAS分析結(jié)果中的值較小的top5%～20%的SNPs作為先驗(yàn)信息，并結(jié)合至GS模型中，對RFI進(jìn)行了遺傳參數(shù)估計(jì)和預(yù)測準(zhǔn)確性分析。金陵花雞RFI為低遺傳力性狀，將GWAS結(jié)果中最顯著的top10%～top15%的SNPs作為先驗(yàn)信息，并結(jié)合至GS的策略可以將基因組育種值的預(yù)測準(zhǔn)確性提高2.10%～5.17%，更適合本群體RFI性狀的選育。