檢測(cè)控制動(dòng)態(tài)性狀的QTL隨機(jī)回歸模型

周曉晶　何倩毓　張巧生

摘 要：表型值隨著時(shí)間（生命時(shí)期、年齡、胎次等）或其他可以量化的因素（生理狀態(tài)、生產(chǎn)水平、代謝率和環(huán)境條件等）變化的性狀稱為動(dòng)態(tài)性狀，如產(chǎn)蛋、泌乳和生長(zhǎng)等性狀皆屬此類。該研究基于隨機(jī)回歸模型思想，提出了檢測(cè)控制動(dòng)態(tài)性狀的多個(gè)數(shù)量性狀位點(diǎn)（QTL）的貝葉斯定位分析方法，在這里視B樣條（B-spline）函數(shù)作為子模型。計(jì)算機(jī)模擬試驗(yàn)驗(yàn)證所建議分析方法的統(tǒng)計(jì)強(qiáng)度和QTL參數(shù)估計(jì)精度。相比于Legendre多項(xiàng)式，由低階多項(xiàng)式構(gòu)成的樣條函數(shù)更廣泛的應(yīng)用于非參數(shù)的數(shù)據(jù)分析，并且能很好的回避龍格現(xiàn)象。實(shí)際數(shù)據(jù)集分析了包含有420頭中國(guó)荷斯坦奶牛的胸圍數(shù)據(jù)，并且得出了具有生物學(xué)意義的結(jié)論。

關(guān)鍵詞：動(dòng)態(tài)性狀；隨機(jī)回歸模型；B-樣條；模擬

中圖分類號(hào) Q348 文獻(xiàn)標(biāo)識(shí)碼 A 文章編號(hào) 1007-7731（2017）14-0030-05

Abstract：Dynamic traits are those phenotypic values change with time and other quantifiable factors such as age，parities，physiological status，performance level and environment etc.On the basis of random regression models（RRM），Bayesian mapping measure has been developed to detect the QTLs controlling dynamic traits，in which the B-spline is imbedded into the RRM.The computer simulation experiments were conducted to demonstrate the statistical powers of QTL detection and the precision of QTL parameter estimation.Compared to the Legendre polynomial，the splines that construct curves from pieces of lower degree polynomials smoothed at selected pointed（knots）are more commonly used in nonparametric data analysis，which seldom have oscillations at two extremes of the curve.A real dataset with approximately 420 China Holstein cows，which contains the records of the chest from the local dairy farm was analyzed and the biological conclusions were derived.

Key words：Dynamic Traits；Random regression model；B-spline；Simulation

動(dòng)物體的很多經(jīng)濟(jì)性狀是隨著時(shí)間（生命時(shí)期、年齡、胎次等）或其他可以量化的因素（生理狀態(tài)、生產(chǎn)水平、代謝率和環(huán)境條件等）變化的，稱為動(dòng)態(tài)性狀。這些性狀廣泛存在于生物和醫(yī)學(xué)領(lǐng)域，如生長(zhǎng)和發(fā)育性狀、產(chǎn)奶量、產(chǎn)蛋量、泌乳和生長(zhǎng)等性狀。任何動(dòng)物和植物的發(fā)育同時(shí)經(jīng)歷系統(tǒng)的和個(gè)體特異性發(fā)育過(guò)程，QTL是整個(gè)基因組中控制該發(fā)育過(guò)程的系統(tǒng)組成部分的基因。由于動(dòng)態(tài)性狀在育種和生產(chǎn)中特殊的經(jīng)濟(jì)地位，揭示這類性狀遺傳規(guī)律及其改良提高的研究工作尤顯重要。

關(guān)于控制動(dòng)態(tài)性狀的QTL檢測(cè)方法較多，其中隨機(jī)回歸模型能夠擬合每個(gè)QTL的遺傳效應(yīng)的動(dòng)態(tài)模式，該方法不僅能夠檢測(cè)控制動(dòng)態(tài)軌跡的QTL，而且能夠推斷任意時(shí)間點(diǎn)的QTL。隨機(jī)回歸模型已是分析縱向數(shù)據(jù)或者重復(fù)記錄數(shù)據(jù)常用方法，有各種形式的隨機(jī)回歸模型，如基于重復(fù)記錄的隨機(jī)回歸模型，基于母體效應(yīng)的隨機(jī)回歸模型，二維隨機(jī)回歸模型等等。隨機(jī)回歸模型因其能靈活地描述個(gè)體在不同時(shí)間點(diǎn)的基因表達(dá)而備受關(guān)注。隨機(jī)回歸模型在動(dòng)物育種中的應(yīng)用范圍已覆蓋到奶牛性狀的遺傳評(píng)估、體質(zhì)評(píng)分、采食量、生長(zhǎng)性狀、生存性狀、環(huán)境互作的表型值等。隨機(jī)回歸子模型是建立隨機(jī)回歸模型的關(guān)鍵因素，不同的子模型將導(dǎo)致動(dòng)態(tài)性狀遺傳參數(shù)估計(jì)值相當(dāng)大的差異。標(biāo)準(zhǔn)化的多項(xiàng)式是正交多項(xiàng)式的一種，廣泛應(yīng)用于大型農(nóng)場(chǎng)動(dòng)物的各種動(dòng)態(tài)性狀育種值的遺傳評(píng)估。

非參數(shù)的Legendre多項(xiàng)式廣泛應(yīng)用于定位動(dòng)態(tài)性狀基因位點(diǎn)。高階的Legendre多項(xiàng)式能夠很好的擬合方差和均值的變化，但是，這樣的多項(xiàng)式經(jīng)常在極值點(diǎn)高估觀測(cè)值，導(dǎo)致龍格現(xiàn)象，也就是說(shuō)，由于曲線上極值點(diǎn)處的振蕩，曲線的擬合度隨著階數(shù)的升高而顯著降低。并且關(guān)于動(dòng)態(tài)性狀的基因定位分析方法幾乎都是建立在試驗(yàn)研究基礎(chǔ)上的。相反，B樣條函數(shù)（B-spline）能夠更廣泛的應(yīng)用于非參數(shù)的數(shù)據(jù)分析去推斷變量的經(jīng)驗(yàn)分布。1974年，Gordon和Riesenfeld用B樣條基函數(shù)代替了Bernstein基函數(shù)，構(gòu)造了B樣條曲線。B樣條曲線分段組成。每一段的參數(shù)t的區(qū)間為[0，1]。這樣就克服了Bezier曲線的缺點(diǎn)：改變Berier曲線任意一個(gè)控制點(diǎn)，曲線上的所有點(diǎn)都變換。BSpline曲線的優(yōu)點(diǎn)：修改某一控制點(diǎn)只引起與該控制點(diǎn)相鄰的曲線形狀發(fā)生變化，遠(yuǎn)處的曲線形狀不受影響。spline曲線段由幾段低階多項(xiàng)式構(gòu)成，他們?cè)谶x擇點(diǎn)相連接，這些選擇點(diǎn)就是所謂的節(jié)點(diǎn)。同使用多項(xiàng)式分析比較，B-spline可以提供更加準(zhǔn)確的遺傳特征函數(shù)的估計(jì)值，而且不易出現(xiàn)端點(diǎn)異常的現(xiàn)象，這是由于在每個(gè)區(qū)段內(nèi)使用了低階多項(xiàng)式的原因。一個(gè)好的隨機(jī)回歸模型應(yīng)該包含一組固定曲線，這些曲線可以描述不同群體變量之間的表型相關(guān)和加性遺傳效應(yīng)與永久環(huán)境效應(yīng)的隨機(jī)回歸部分；還可能存在其他的隨機(jī)回歸因素，例如群體效應(yīng)和母體遺傳效應(yīng)。

Meyer采用牛奶產(chǎn)量測(cè)定日記錄模型分析澳大利亞安格斯牛的生長(zhǎng)性狀，在該文中，作者比較了隨機(jī)回歸模型的11個(gè)不同組合，他們分別是二次和三次的B樣條曲線，結(jié)點(diǎn)個(gè)數(shù)達(dá)到了6個(gè)。通過(guò)計(jì)算，擬合生長(zhǎng)性狀的最優(yōu)隨機(jī)回歸模型為帶有3個(gè)結(jié)點(diǎn)的三次B樣條曲線。

本研究基于隨機(jī)回歸模型提出了貝葉斯檢測(cè)方法，在這里，B樣條函數(shù)作為子模型來(lái)表征群體均值的動(dòng)態(tài)性、QTL效應(yīng)和指定個(gè)體的隨時(shí)間獨(dú)立的環(huán)境效應(yīng)。計(jì)算機(jī)模擬說(shuō)明方法的有效性，實(shí)際數(shù)據(jù)集來(lái)源于當(dāng)?shù)啬翀?chǎng)的奶牛數(shù)據(jù)，包括420個(gè)個(gè)體的胸圍記錄，試圖得出具有生物學(xué)意義的結(jié)論。

1 研究方法

2 模擬

進(jìn)行模擬實(shí)驗(yàn)驗(yàn)證基于貝葉斯B樣條（B-spline）的檢測(cè)動(dòng)態(tài)性狀的定位分析方法的有效性和適應(yīng)性。在貝葉斯分析的框架內(nèi)，視B樣條函數(shù)為回歸模型的子模型，能夠同時(shí)描述群體均值、多個(gè)QTL的遺傳效應(yīng)，以及其他環(huán)境因子隨時(shí)間的變化。我們知道，在功能定位中采用的Logisitic曲線由于其不可加性，不能實(shí)現(xiàn)同時(shí)檢測(cè)多個(gè)QTL。

以含有150個(gè)獨(dú)立個(gè)體的回交群體為例。每一基因位點(diǎn)包含兩種基因型。在長(zhǎng)度為600cM的染色體片段上等距離設(shè)置61個(gè)共顯性分子標(biāo)記。放置10個(gè)控制動(dòng)態(tài)性狀軌跡的QTL。假定性狀表型值隨時(shí)間的變化及加性效應(yīng)隨時(shí)間的變化服從帶有4個(gè)節(jié)點(diǎn)及2階多項(xiàng)式的B-spline。每個(gè)個(gè)體的初始時(shí)間點(diǎn)指定為10，20，30，40，50，60，70，80，90，100。4個(gè)節(jié)點(diǎn)選取為10，40，70，100。節(jié)點(diǎn)的選擇是至關(guān)重要的，節(jié)點(diǎn)過(guò)多導(dǎo)致過(guò)高擬合，節(jié)點(diǎn)太少導(dǎo)致擬合效果差。

對(duì)于模擬數(shù)據(jù)集，用階數(shù)分別為4，5和6階的Legendre多項(xiàng)式替換B-spline，來(lái)檢測(cè)模擬的10個(gè)QTL。當(dāng)然，4階Legendre多項(xiàng)式和B-spline具有相同個(gè)數(shù)的回歸系數(shù)。同時(shí)，我們希望5階和6階Legendre多項(xiàng)式相比于4階Legendre多項(xiàng)式具有較高的擬合度。接下來(lái)，在相同的實(shí)驗(yàn)設(shè)計(jì)下，利用Legendre多項(xiàng)式生成模擬數(shù)據(jù)，然后利用貝葉斯B-spline檢測(cè)QTL。

3 實(shí)例分析

動(dòng)物資源群體來(lái)源于大慶本地牧場(chǎng)——紅驥牧場(chǎng)。奶牛品種為中國(guó)荷斯坦奶牛，重要經(jīng)濟(jì)性狀為胸圍動(dòng)態(tài)性狀，測(cè)定時(shí)間為2010年1月至2013年13月。具體為2010年出生的60頭奶牛、2011年出生的70頭奶牛及2012年出生的80頭奶牛。由于奶牛早期生長(zhǎng)發(fā)育狀況與日后的產(chǎn)奶量、繁殖問(wèn)題和健康問(wèn)題有著密切關(guān)系，因此對(duì)奶牛生長(zhǎng)性狀的研究主要針對(duì)的是從出生到第一個(gè)泌乳期結(jié)束這段時(shí)間，分為初生、斷奶、周歲和頭胎分娩這4個(gè)時(shí)間點(diǎn)。在實(shí)際牛場(chǎng)管理過(guò)程中，這幾個(gè)時(shí)間點(diǎn)正是牛只轉(zhuǎn)群的時(shí)間，飼料成分、管理方式、飼養(yǎng)環(huán)境都發(fā)生了明顯變化，這樣的劃分方式也與生產(chǎn)實(shí)際相契合。所以特別關(guān)注初生、斷奶、周歲和頭胎分娩這四個(gè)時(shí)間點(diǎn)的體重記錄，將它們視為一組。共4組數(shù)據(jù)。處理數(shù)據(jù)，去掉異常值后剩余143個(gè)個(gè)體的4組數(shù)據(jù)。

帶有3個(gè)結(jié)點(diǎn)和2個(gè)多項(xiàng)式片段的模型擬合動(dòng)態(tài)性狀的群體均值和每個(gè)QTL的遺傳效應(yīng)。當(dāng)結(jié)點(diǎn)數(shù)選為3，6，和12個(gè)月時(shí)，該方法檢測(cè)到了控制胸圍的生長(zhǎng)軌跡。參數(shù)估計(jì)見表5。另外，當(dāng)內(nèi)部結(jié)點(diǎn)選為4月和9月時(shí)，擬合結(jié)果幾乎一致。這說(shuō)明該方法的精度較強(qiáng)的依賴于內(nèi)部結(jié)點(diǎn)的選取。

4 結(jié)論

本文提出了視B-樣條函數(shù)為子模型的隨機(jī)回歸模型，目的是檢測(cè)控制動(dòng)態(tài)性狀的QTL，實(shí)現(xiàn)了群體均值、QTL遺傳效應(yīng)和其他隨時(shí)間變化的環(huán)境因素的估計(jì)。B-樣條函數(shù)相比于勒讓德多項(xiàng)式的優(yōu)勢(shì)是能夠避免端點(diǎn)異常值現(xiàn)象。不同結(jié)點(diǎn)的選擇能夠?qū)е螺^大變異，因此在基于B-樣條函數(shù)的QTL定位中結(jié)點(diǎn)的選取是關(guān)鍵因素。理論上結(jié)點(diǎn)越多越容易產(chǎn)生數(shù)據(jù)的過(guò)高擬合，結(jié)點(diǎn)越少越容易產(chǎn)生數(shù)據(jù)的過(guò)低擬合。

目前，奶牛育種的總體趨勢(shì)是在保持奶牛產(chǎn)奶量以及乳成分等優(yōu)良遺傳性狀的同時(shí)，兼顧奶牛的軀體結(jié)構(gòu)、趾蹄健康、使用壽命、繁殖性能等綜合遺傳性能的選育，以獲得奶牛養(yǎng)殖的最大經(jīng)濟(jì)效益。有研究表明-加強(qiáng)體型性狀的選擇對(duì)奶牛產(chǎn)奶量以及乳成分的提高有利，也有利于降低體細(xì)胞數(shù)，增強(qiáng)個(gè)體乳房炎的抗性，鑒于此，有必要研究奶牛除了產(chǎn)奶量之外的動(dòng)態(tài)性狀的表征。

參考文獻(xiàn)

[1]C.R Henderson.Analysis of Covariance in the Mixed Model：Higher-Level，Non-homogeneous，and Random Regressions[J].Biometrics，1992，38（3）：623-640.

[2]J.Jamrozik，L.R.Schaeffer and J.C.M.Dekkers.“Genetic evaluation of dairy cattle using test day yields and random regression model [J].Journal of Dairy Science，1997，80（6）：1217-1226.

[3]L.R.Schaeffer.Application of random regression models in animal breeding [J].Livestock Production Science，2004，86（1-3）：35-45.

[4]S.Macgregor，S.A.Knott，I.White and P.M.Visscher.Quantitative trait locus analysis of longitudinal quantitative trait data in complex pedigrees [J].Genetics，2005，171：1365-1376.

[5]R.Yang，Q.Tian，S.Xu.Mapping quantitative trait loci for longitudinal traits in line crosses [J].Genetics 2006，173：2339-2356.

[6]R.Yang，H.Gao，X.Wang，J.Zhang，Z.B.Zeng and R.Wu.A semiparametric approach for composite functional mapping of dynamic quantitative traits [J].Genetics 2007（177）：1859-1870.

[7]R.Yang，S.Xu.Bayesian shrinkage analysis of quantitative trait Loci for dynamic traits [J].Genetics，2007，176：1169-1185.

[8]C.de Boor.A Practical Guide to Splines [M].New York： Springer Verlag，2001.

[9]K.Meyer，W.G.Hill.Estimation of genetic and phenotypic covariance functions for longitudinal or “repeated” records by restricted maximum likelihood [J] Livest Prod Sci，1997，47：185200.

[10]B.J.De Groot，J.F.Keown，L.D.van Vleck，et al.Genetic parameters and responses of linear type，yield traits，and somatic cell scores to divergent selection for predicted transmitting ability for type in Holsteins [J].Journal of Dairy Science，2002，85（6）：1578-1585.

[11]Y.de Haas，L.L.G.Janss and H.N.Kadarmideen.Genetic and phenotypic parameters for conformation and yield traits in three Swiss dairy cattle breeds [J].Journal of Animal Breeding and Genetics，2007，124（01）：12-19.

[12]M.J.Sillanp??，E.Arjas.Bayesian mapping of multiple quantitative trait loci from incomplete inbred line cross data [J]Genetics，1998，148：1373-1388.

[13]M.J.Sillanp??，E.Arjas.Bayesian mapping of multiple quantitative trait loci from incomplete outbred offspring data [J].Genetics，1999，151：1605-1619.

[14]J.M.Satagopan，B.S.Yandell，M.A.Newton，T.C.Osborn.A bayesian approach to detect quantitative trait loci using Markov chain Monte Carlo [J].Genetics，1996，144（02）：805-816.

[15]N.Yi.A unified Markov chain Monte Carlo framework for mapping multiple quantitative trait loci [J].Genetics，2004，167（02）：967-975.

（責(zé)編：張宏民）