朱墨, 鄭麥青, 崔煥先, 趙桂蘋?, 劉楊

1南京農(nóng)業(yè)大學(xué)動物科技學(xué)院動物繁育系，南京 210095；2中國農(nóng)業(yè)科學(xué)院北京畜牧獸醫(yī)研究所/動物營養(yǎng)學(xué)國家重點(diǎn)實(shí)驗(yàn)室，北京 100193

0 引言

【研究意義】畜禽的重要經(jīng)濟(jì)性狀大多是遺傳結(jié)構(gòu)復(fù)雜的數(shù)量性狀。常規(guī)的畜禽育種工作中，利用性狀的表型記錄值和系譜信息來計(jì)算個體之間的親緣關(guān)系，通過最佳線性無偏預(yù)測（best linear unbiased prediction, BLUP）進(jìn)行個體育種值（estimated breeding value, EBV）的估計(jì)，然后再進(jìn)行排序、選擇[1]。隨著分子遺傳學(xué)的發(fā)展，微衛(wèi)星、限制性片段長度多態(tài)性（restriction fragment length polymorphism, RFLP）、單核苷酸多態(tài)性（single nucleotide polymorphism,SNP）等分子標(biāo)記被開發(fā)出來，并用于動物育種，稱為標(biāo)記輔助選擇（maker assisted selection, MAS）[2]。然而，畜禽的大部分重要經(jīng)濟(jì)性狀受微效多基因控制，少數(shù)的幾個標(biāo)記解釋的遺傳變異十分有限，限制了標(biāo)記輔助選擇在畜禽育種中的應(yīng)用[3]。2001年，MEUWISSEN等提出了基因組選擇的概念[4]，其核心思想是利用覆蓋全基因組的標(biāo)記進(jìn)行基因組育種值估計(jì)，結(jié)合表型記錄和系譜信息進(jìn)行選擇。相比標(biāo)記輔助選擇，基因組選擇利用了覆蓋全基因組的標(biāo)記，可以更好地解釋表型變異[5-6]?；蚪M育種值估計(jì)是基因組選擇的核心內(nèi)容。根據(jù)使用的統(tǒng)計(jì)模型的不同，可以將基因組育種值的計(jì)算模型分為直接法和間接法，直接法模型利用參考群體和預(yù)測群體的標(biāo)記信息構(gòu)建個體間的親緣關(guān)系矩陣，并將關(guān)系矩陣代入混合模型方程組，通過求解混合模型方程組的方式獲得個體的基因組育種值，如GBLUP[7]、SSGBLUP[8]等方法；間接法模型利用參考群體估計(jì)標(biāo)記效應(yīng)，在預(yù)測群體中將標(biāo)記效應(yīng)累加，進(jìn)而直接獲得個體的基因組育種值，間接法模型的代表主要是貝葉斯類模型，MEUWISSEN在提出基因組選擇理論的同時也提出了兩種貝葉斯方法：BayesA[4]和 BayesB[4]，隨后，一系列改進(jìn)的貝葉斯法，如BayesC[9]、BayesCπ[9]、BayesDπ[9]和 BayesLASSO[10]等陸續(xù)被提了出來?！厩叭搜芯窟M(jìn)展】作為重要的農(nóng)業(yè)動物，雞的基因組草圖于2004年率先發(fā)布[11]。隨后，雞的60 K和600 K芯片相繼開發(fā)問世[12]。中國農(nóng)業(yè)大學(xué)和中國農(nóng)業(yè)科學(xué)院也推出了適合我國地方品種蛋雞50 K和肉雞55 K芯片[13]，極大地促進(jìn)了基因組選擇技術(shù)在家禽上的應(yīng)用。ZHANG等利用BLUP、GBLUP和BLUP|GA對中國地方雞生長和屠宰性狀的基因組育種值預(yù)測，結(jié)果發(fā)現(xiàn)，BLUP|GA的表現(xiàn)最佳[14]。LIU等在三黃雞的生長和屠宰性狀的研究中發(fā)現(xiàn)，在GBLUP、BayesLASSO和BayesMix4模型之間并沒有顯著的差異[15]。WOLC等在對蛋殼質(zhì)量的基因組選擇研究中發(fā)現(xiàn)，多性狀的一步法模型的表現(xiàn)最優(yōu)[16]?！颈狙芯壳腥朦c(diǎn)】屠宰性狀是肉雞的重要經(jīng)濟(jì)性狀。但是，屠宰性狀不能直接活體測量，只能在肉雞屠宰后才能測量，或者通過間接測量其他性狀的手段測量。使用間接測量的表型數(shù)據(jù)進(jìn)行育種值的估計(jì)，會導(dǎo)致估計(jì)的結(jié)果不夠準(zhǔn)確，進(jìn)而影響到個體的準(zhǔn)確選擇。近些年來，已經(jīng)出現(xiàn)了許多關(guān)于肉雞的基因組選擇研究，但是，關(guān)于屠宰性狀基因組選擇研究鮮有報(bào)道[17-19]。【擬解決的關(guān)鍵問題】在家禽育種中，基因組選擇主要是為了提高育種值估計(jì)的準(zhǔn)確性。本研究為探討肉雞屠宰性狀基因組選擇的準(zhǔn)確性，基于不同的模型和假設(shè)估計(jì)白羽肉雞屠宰性狀的基因組估計(jì)育種值（genomic estimated breeding value, GEBV），采用世代驗(yàn)證的方法評估預(yù)測結(jié)果的準(zhǔn)確性，第5、6世代作為參考群，第 7世代作為驗(yàn)證群，為白羽肉雞的育種中育種值估計(jì)的策略提供依據(jù)。

1 材料與方法

1.1 試驗(yàn)動物

試驗(yàn)動物群體來自廣東新廣農(nóng)牧股份有限公司的白羽肉雞祖代父系群體（B系）第5到7世代，共計(jì)3 362只，來源于227只公雞和1 305只母雞的后代。在42日齡時屠宰，記錄胸肌重（breast muscle weight,BrW）、屠體重（carcass weight, CW）、腿肌重（thigh muscle weight, ThW），計(jì)算胸肌率（breast muscle rate,BrR）、腿肌率（thigh muscle rate, ThR）。對表型數(shù)據(jù)進(jìn)行預(yù)處理，刪除表型的缺失值和異常值（平均值±3倍標(biāo)準(zhǔn)差）。利用R語言中的GLM函數(shù)校正表型數(shù)據(jù)的世代、性別和批次效應(yīng)。以上5個性狀的描述性統(tǒng)計(jì)量見表1。

表1 各屠宰性狀的描述性統(tǒng)計(jì)Table 1 Descriptive statistics for each carcass trait

1.2 基因型數(shù)據(jù)的獲取與質(zhì)控

1.2.1 基因組 DNA的提取及分型 采用常規(guī)酚-氯仿抽提法提取血樣基因組DNA，使用NanoDrop 2000核酸分析儀檢測 DNA的濃度和質(zhì)量。質(zhì)檢合格后的DNA樣品送至北京康普森生物技術(shù)有限公司，使用與中國農(nóng)業(yè)科學(xué)院北京畜牧獸醫(yī)研究所聯(lián)合研發(fā)的“京芯一號”雞 55 K SNP芯片進(jìn)行基因分型[13]。

1.2.2 基因型數(shù)據(jù)的質(zhì)控 采用 PLINK (V1.90) 軟件對芯片的基因型數(shù)據(jù)進(jìn)行質(zhì)量控制[20]。質(zhì)量控制的條件如下：1）保留樣本檢出率大于90%的個體；2）保留SNP檢出率大于90%的位點(diǎn)；3）保留次要等位基因頻率大于5%的SNP位點(diǎn)。芯片經(jīng)過質(zhì)控后，保留3 314個樣本和42 104個SNP用于后續(xù)分析。

1.3 基因組選擇的統(tǒng)計(jì)模型

1.3.1 GBLUP模型

式中，y是性狀的表型值向量；b是固定效應(yīng)的向量；u是加性遺傳效應(yīng)向量，服從正態(tài)分布：u～N(0,Gσ2 u)；e是隨機(jī)殘差效應(yīng)向量，服從正態(tài)分布：e～N(0,Iσ2 e)；X和Z分別為對應(yīng)的設(shè)計(jì)矩陣。

1.3.2 BayesB模型

式中，y是性狀的表型值向量；b是固定效應(yīng)的向量；Zi是第i個位點(diǎn)的基因型（0/1/2）；gi是第i個位點(diǎn)的效應(yīng)值；e是隨機(jī)殘差效應(yīng)向量；X是對應(yīng)的設(shè)計(jì)矩陣。BayesB 方法引入了一個SNP效應(yīng)指示變量π，假設(shè)大部分無效應(yīng)（比例為 π），只有一小部分標(biāo)記有效應(yīng)（比例為 1-π），且這一部分有效應(yīng)的方差服從逆卡方分布[4]。本研究中，使用基于 R 語言的BGLR 包進(jìn)行 BayesB 預(yù)測，設(shè)置蒙特卡洛馬爾科夫鏈（Markov chain monte carlo, MCMC）長為20 000，burn in 為 10 000，π 值為 0.95。

1.4 基因組選擇準(zhǔn)確性的分析方法

本研究采用世代驗(yàn)證法來評估基因組選擇的準(zhǔn)確性。采用第5和6世代的群體作為參考群體，第7世代的群體作為驗(yàn)證群體。本研究使用基因組估計(jì)育種值（GEBV）與校正了世代、性別、批次后的表型（y*）的皮爾遜相關(guān)系數(shù)作為基因組選擇準(zhǔn)確性的評價指標(biāo)。本研究中，對每個性狀的驗(yàn)證進(jìn)行20個重復(fù)。

2 結(jié)果

2.1 白羽肉雞B系屠體性狀的遺傳力估計(jì)

胸肌率、胸肌重、屠體重、腿肌率和腿肌重的描述性統(tǒng)計(jì)量匯總于表 1。分別使用基于系譜構(gòu)建的親緣關(guān)系A(chǔ)矩陣和基于全基因組SNP信息構(gòu)建的親緣關(guān)系G矩陣，采用平均信息約束最大似然算法（average information restricted maximum likelihood, AIREML）估計(jì)加性遺傳方差和殘差方差，本研究中的方差組分使用ASReml 4.1.0 軟件進(jìn)行估計(jì)[21]?；谙底V構(gòu)建的A陣估計(jì)BrR、BrW、CW、ThR和ThW的遺傳力均高于基于基因組SNP信息構(gòu)建的G陣估計(jì)BrR、BrW、CW、ThR和ThW的遺傳力。各性狀的遺傳力估計(jì)結(jié)果見表2。

表2 各屠宰性狀的遺傳力估計(jì)Table 2 Results of heritability for each carcass trait

2.2 白羽肉雞 B系屠宰性狀基因組選擇不同模型間比較

2.2.1 白羽肉雞B系胸肌率、胸肌重的世代驗(yàn)證 采用世代驗(yàn)證法，使用 GBLUP 和 BayesB 方法對BrR 和 BrW 進(jìn)行基因組預(yù)測分析，采用GEBV與y*的皮爾遜相關(guān)系數(shù)作為基因組選擇準(zhǔn)確性的評價指標(biāo)。兩種方法對 BrR 的預(yù)測準(zhǔn)確性分別為 0.3262、0.3765。兩種方法對 BrW 的預(yù)測準(zhǔn)確性分別為0.2871、0.2257。準(zhǔn)確性的結(jié)果見表3。由結(jié)果可以看出，在對 BrR 的基因組育種值估計(jì)準(zhǔn)確性上，BayesB方法估計(jì)的準(zhǔn)確性要高于 GBLUP 方法估計(jì)的準(zhǔn)確性；而在對 BrW 的基因組育種值估計(jì)準(zhǔn)確性上，GBLUP 方法估計(jì)的準(zhǔn)確性要高于 BayesB 方法估計(jì)的準(zhǔn)確性。

2.2.2 白羽肉雞 B系屠體重的世代驗(yàn)證 采用世代驗(yàn)證法，使用 GBLUP 和 BayesB 方法對 CW 進(jìn)行基因組預(yù)測分析，采用GEBV與y*的皮爾遜相關(guān)系數(shù)作為基因組選擇準(zhǔn)確性的評價指標(biāo)。兩種方法對 CW的預(yù)測準(zhǔn)確性分別為 0.2780、0.1376。準(zhǔn)確性的結(jié)果見表3。由結(jié)果可以看出，在對 CW 的基因組育種值估計(jì)準(zhǔn)確性上，GBLUP 方法估計(jì)的準(zhǔn)確性要高于BayesB 方法估計(jì)的準(zhǔn)確性。

2.2.3 白羽肉雞B系腿肌率、腿肌重的世代驗(yàn)證 采用世代驗(yàn)證法，使用 GBLUP 和 BayesB 方法對ThR和ThW進(jìn)行基因組預(yù)測分析，采用GEBV與y*的皮爾遜相關(guān)系數(shù)作為基因組選擇準(zhǔn)確性的評價指標(biāo)。兩種方法對ThR的預(yù)測準(zhǔn)確性分別為0.2153、0.2525。兩種方法對 ThW 的預(yù)測準(zhǔn)確性分別為0.2126、0.2844。準(zhǔn)確性的結(jié)果見表 3。由結(jié)果可以看出，在對ThR和ThW的基因組育種值估計(jì)準(zhǔn)確性上，BayesB方法估計(jì)的準(zhǔn)確性要高于 GBLUP方法估計(jì)的準(zhǔn)確性。

表3 各屠宰性狀基于GBLUP和BayesB方法的世代驗(yàn)證結(jié)果Table 3 Results of generation validation based on GBLUP and BayesB method for each carcass trait

3 討論

本研究系統(tǒng)探討了基于不同模型對白羽肉雞屠宰性狀的遺傳評估，屠宰性狀的遺傳力估計(jì)結(jié)果與先前的研究基本一致[14-15,22]?；谙底V構(gòu)建的親緣關(guān)系A(chǔ)矩陣估計(jì)的遺傳力高于基于全基因組SNP信息構(gòu)建的親緣關(guān)系 G矩陣估計(jì)的遺傳力?？赡艿脑蚴腔蚪M信息可以反映系譜記錄中無法體現(xiàn)的基因的同源狀態(tài)，更真實(shí)地反映個體之間的遺傳關(guān)系[7]。有研究發(fā)現(xiàn)，由于個體之間相似的環(huán)境組分，使得加性遺傳方差的估計(jì)值升高，從而導(dǎo)致遺傳力的過高估計(jì)[23]。因此，使用基于全基因組標(biāo)記進(jìn)行遺傳評估相較于系譜記錄，可以排除環(huán)境效應(yīng)的影響，估計(jì)結(jié)果更加準(zhǔn)確[24-25]。

本研究系統(tǒng)探討了GBLUP和BayesB方法對白羽肉雞屠宰性狀進(jìn)行基因組預(yù)測的效果。研究發(fā)現(xiàn)，基因組預(yù)測的準(zhǔn)確性與性狀的遺傳力大致呈正相關(guān)。使用GBLUP方法和BayesB方法時，預(yù)測準(zhǔn)確性最高的性狀均是胸肌率。對于胸肌率、腿肌率和腿肌重，BayesB方法的基因組預(yù)測準(zhǔn)確性高于GBLUP方法；對于屠體重和胸肌重，GBLUP方法的預(yù)測準(zhǔn)確性高于BayesB。BayesB 方法引入了一個SNP效應(yīng)指示變量π，假設(shè)大部分無效應(yīng)（比例為π），只有一小部分標(biāo)記有效應(yīng)（比例為1-π），且這一部分有效應(yīng)的方差服從逆卡方分布[4]。本研究中，除了在對屠體重的預(yù)測外，BayesB方法的預(yù)測準(zhǔn)確性均高于GBLUP方法的結(jié)果，可能是由于BayesB方法的模型假設(shè)更加符合數(shù)量性狀的微效多基因理論，但BayesB方法預(yù)測的偏差與GBLUP方法相比較大，這需要在后續(xù)研究中繼續(xù)探討。ZENG等研究發(fā)現(xiàn)，貝葉斯方法在性狀的遺傳結(jié)構(gòu)未知時優(yōu)于 GBLUP[26]。TENG等在對一個雜交雞群體生長性狀的基因組預(yù)測發(fā)現(xiàn)，BayesB方法優(yōu)于GBLUP方法[27]。因此，基因組選擇的計(jì)算模型沒有絕對的標(biāo)準(zhǔn)，要根據(jù)目標(biāo)性狀的遺傳結(jié)構(gòu)來具體選擇。

畜禽的育種工作注重時效性，計(jì)算效率是基因組選擇在實(shí)際生產(chǎn)中應(yīng)用時需要考慮的一個重要因素。貝葉斯模型的參數(shù)求解過程通過蒙特卡洛馬爾科夫鏈（Markov chain monte carlo, MCMC）過程，采用高強(qiáng)度的吉布斯（Gibbs）抽樣，往往需要經(jīng)過上萬次迭代，且無法并行計(jì)算，這常常限制了貝葉斯類的方法在育種中的應(yīng)用，使得基于混合線性模型的方法（GBLUP為代表）在基因組選擇中的應(yīng)用最為廣泛。但是，基于混合線性模型的方法的預(yù)測準(zhǔn)確性要低于貝葉斯類模型[9]。在一些研究中，將最顯著的SNP或驗(yàn)證過的QTL加入混合模型作為協(xié)變量來提高預(yù)測準(zhǔn)確性，然而這些SNP或QTL解釋的表型方差有限，可能會導(dǎo)致假陽性的結(jié)果[28-29]。還有研究者通過對SNP進(jìn)行加權(quán)，構(gòu)建性狀特異性的親緣關(guān)系矩陣，也能在一定程度上提高混合線性模型基因組預(yù)測的準(zhǔn)確性[30-32]。然而，這些SNP的效應(yīng)很容易被干擾。本研究中，BayesB方法的計(jì)算時長遠(yuǎn)遠(yuǎn)大于GBLUP方法（圖 1），但是BayesB方法的預(yù)測準(zhǔn)確性顯著高于GBLUP方法。而且，對于3 000的樣本量和55 K 的標(biāo)記密度，BayesB的計(jì)算時間大約為7 h，尚在可接受的范圍。在實(shí)際的育種中，使用低密度芯片數(shù)據(jù)，在樣本和標(biāo)記數(shù)量都不是十分龐大的情況下，為了提高基因組育種值預(yù)測的準(zhǔn)確性，使用貝葉斯方法進(jìn)行基因組育種值的估計(jì)是可行的方案。

4 結(jié)論

本研究使用GBLUP和BayesB兩種方法對白羽肉雞胸肌率、胸肌重、屠體重、腿肌率和腿肌重 5種屠宰性狀進(jìn)行了基因組選擇分析。研究發(fā)現(xiàn)，基因組選擇的準(zhǔn)確性與性狀的遺傳力大致呈正相關(guān)。使用GBLUP和BayesB方法時，預(yù)測準(zhǔn)確性最高的是胸肌率。對于胸肌率、腿肌率和腿肌重，BayesB方法的預(yù)測準(zhǔn)確性均高于GBLUP方法；對于屠體重和胸肌重，GBLUP方法的預(yù)測準(zhǔn)確性高于BayesB。但是，BayesB方法的計(jì)算時間要長于GBLUP方法。在實(shí)際的育種工作中，需要綜合考慮育種值估計(jì)的準(zhǔn)確性和育種的時效性來決定用何種方式估計(jì)基因組育種值。