基于BLUP的單性狀動物育種模型的研究及實現(xiàn)

胡世洋，楊慧贊，張永德，林勇，孫全亮

（1.賀州學(xué)院，計算機(jī)科學(xué)與信息工程學(xué)院，廣西賀州542899；2.廣西水產(chǎn)科學(xué)研究院，廣西南寧530021）

基于BLUP的單性狀動物育種模型的研究及實現(xiàn)

胡世洋1，楊慧贊2，張永德2，林勇2，孫全亮1

（1.賀州學(xué)院，計算機(jī)科學(xué)與信息工程學(xué)院，廣西賀州542899；2.廣西水產(chǎn)科學(xué)研究院，廣西南寧530021）

BLUP算法可以顯著提高動物遺傳進(jìn)展速度，但是它的計算過程比較復(fù)雜，數(shù)據(jù)量也相當(dāng)浩繁，若將這種復(fù)雜、抽象的數(shù)學(xué)算法轉(zhuǎn)換成可視化的計算機(jī)系統(tǒng)將大大有利于生產(chǎn)力的提高，針對這一問題，以羅非魚的遺傳育種為研究切入點，將BLUP算法應(yīng)用于羅非魚的育種值估計當(dāng)中，構(gòu)建出系統(tǒng)的模型，并采用面向?qū)ο蟮乃季S方法將BLUP的算法轉(zhuǎn)化為可執(zhí)行程序。實驗結(jié)果表明，開發(fā)出的軟件可以方便、高效地進(jìn)行羅非魚育種值的估計，當(dāng)影響羅非魚的外部因素發(fā)生改變時，最終的計算結(jié)果都能較準(zhǔn)確地反應(yīng)出相關(guān)因素對羅非魚的影響效應(yīng)。

BLUP；算法；羅非魚；軟件

1　引言

動物的選擇育種問題一直是養(yǎng)殖業(yè)的重要科題之一，隨著數(shù)理統(tǒng)計學(xué)與線性模型理論、計算機(jī)科學(xué)與互聯(lián)網(wǎng)絡(luò)技術(shù)的迅速發(fā)展，動物的育種值估計的方法發(fā)生了根本的變化［1］。最佳線性無偏預(yù)測法（Best Linear Unbiased Prediction，BLUP）由于能顯著提高遺傳進(jìn)展速度，特別是對于低遺傳力性狀和顯性性狀，因而被廣泛運用到動物的遺傳改良當(dāng)中。BLUP法可以通過對稱矩陣求解的方式把動物的種類效應(yīng)、動物自身效應(yīng)和機(jī)誤誤差區(qū)分開來，可以有效提高動物遺傳值的準(zhǔn)確性［2］。由于BLUP法的過程比較復(fù)雜，計算量也相當(dāng)浩繁，如何將BLUP算法的數(shù)學(xué)思想轉(zhuǎn)化成計算機(jī)思維，用程序設(shè)計語言將其運行并且實現(xiàn)代碼的可重用性是一個值得深入研究的問題［3］。文章以羅非魚的養(yǎng)殖育種為例研究BLUP算法的應(yīng)用及實現(xiàn)方法。

2　BLUP算法的工作原理

一般混合模型可表示為：

y是所有觀察值構(gòu)成的向量；b是所有固定效應(yīng)（包括）構(gòu)成的向量；X是固定效應(yīng)的關(guān)聯(lián)矩陣；u是所有隨機(jī)效應(yīng)構(gòu)成的向量；Z是隨機(jī)效應(yīng)的關(guān)聯(lián)矩陣；e是隨機(jī)殘差向量［4］。

隨機(jī)變量的數(shù)學(xué)期望：

方差-協(xié)方差矩陣結(jié)構(gòu)：

BLUP估計一般方程：

混合模型方程組的簡化形式：

3　BLUP算法在羅非魚研究過程中的具體研究及分析

某種魚場有如下性能測定資料，測定性狀為達(dá)500g日齡，已知該性狀的遺傳力為0.33，已知某群達(dá)500g的羅非魚相關(guān)記錄數(shù)據(jù)如表1所示，下面將對該群羅非魚的性狀資料進(jìn)行個體育種值估計。

表1　羅非魚達(dá)500g日齡記錄

根據(jù)資料性質(zhì)，可對羅非魚達(dá)500g日齡寫出如下動物模型：

用矩陣形式表示，則對于該資料有：

構(gòu)建混合模型方程組：

求解混合模型方程組有：

4　Blup算法在計算機(jī)中的算法實現(xiàn)方案

4.1系統(tǒng)框架的構(gòu)建

系統(tǒng)采用MVC的設(shè)計框架，由ODBC接口、相關(guān)數(shù)據(jù)庫構(gòu)件以及文件流操作模塊組成Model層［5］，由矩陣生成模塊、矩陣運算模塊、方程求解模塊以及相關(guān)賦值運算模塊組成Controller層，由于相關(guān)顯示窗體模塊組成View層［6］。系統(tǒng)可獲取的外部數(shù)據(jù)類型包括數(shù)據(jù)庫文件，XML文件、txt文件或xls文件，系統(tǒng)框架如果圖1所示。

圖1　系統(tǒng)框架圖

4.2系統(tǒng)相關(guān)操作類的實現(xiàn)

本系統(tǒng)是采用C#程序語言進(jìn)行開發(fā)，在系統(tǒng)實現(xiàn)過程中需要大量的矩陣運算，因此在系統(tǒng)中必須建立Matrix類來管理矩陣數(shù)據(jù)，同時需要有一個MatrixOperator類來操作Matrix類，從而實現(xiàn)矩陣的混合運算及相關(guān)信息的輸出［7］。系統(tǒng)要從數(shù)據(jù)庫文件中獲取外部數(shù)據(jù)文件并存取在內(nèi)存當(dāng)中，因此需要有DataModel類來負(fù)責(zé)外部數(shù)據(jù)的獲取來組織［8］。最后內(nèi)存中的數(shù)據(jù)需要在AnalysisOperator類中完成矩陣的構(gòu)建和運算，并將最終結(jié)果輸出到顯示窗體當(dāng)中。系統(tǒng)的相關(guān)類型如圖2所示。

圖2　系統(tǒng)類圖

根據(jù)上一節(jié)基于BLUP算法羅非魚育種值計算規(guī)則，可以構(gòu)建出其在計算機(jī)中的實現(xiàn)算法為：

步驟1：讀取數(shù)據(jù)，基礎(chǔ)數(shù)據(jù)可以是來自數(shù)據(jù)，或者是其它外部文件。

步驟2：構(gòu)建系數(shù)矩陣。針對羅非魚的數(shù)據(jù)特征，可以對體重達(dá)500g的日齡構(gòu)建動物模型，并由此得出固定因素和非固定因素的系數(shù)據(jù)矩陣。

步驟3：構(gòu)建親緣矩陣。根據(jù)羅非魚的個體系譜列表來構(gòu)建分子親緣矩陣。

步驟4：構(gòu)建混合模型矩陣。根據(jù)BLUP估計一般方程，混合模型方程組的簡化形式。

步驟5：混合矩陣方程的求解。得到了矩陣向量元素將是各個影響因素的育種估計值。

4.3算法實現(xiàn)的幾個具體問題

4.3.1矩陣類的定義

矩陣的存儲方式是采用二維數(shù)組存儲。數(shù)組類型為double型，Matrix類中不僅有行和列的相關(guān)屬性，而且有其對應(yīng)的構(gòu)造方法以及實例化矩陣的相關(guān)方法。

4.3.2矩陣操作類的定義算法

BLUP算法在運算過程中涉及矩陣相乘、矩陣數(shù)乘、矩陣轉(zhuǎn)置和矩陣求逆等相關(guān)操作。因此我們在MatrixOperator類必須要聲明這些方法，其中MatrixPrint（）可以實現(xiàn)矩陣內(nèi)部元素的輸出。

4.3.3系數(shù)矩陣的構(gòu)建算法

由上一節(jié)的動物模型可以得知，在單性狀遺傳值估計中，固定效應(yīng)的關(guān)聯(lián)矩陣X的維數(shù)是由測量記錄數(shù)和魚池這一固定效應(yīng)的數(shù)量決定的，而矩陣X的元素值是因固定效應(yīng)的量值決定。關(guān)聯(lián)矩陣X的構(gòu)建算法則是先確定魚場的個數(shù)，再遍歷測量記錄條目，若當(dāng)前記錄的魚場數(shù)與先前設(shè)置的魚池編號相同則將其在X矩陣中所對應(yīng)的無素值置為1，其它元素置為0。而隨機(jī)效應(yīng)的關(guān)聯(lián)矩陣Z則是與測量記錄數(shù)相關(guān)的單位矩陣。顯性向量矩陣y則是由各條記錄數(shù)中的日齡字段決定。這兩個矩陣的構(gòu)建較為簡單，在此不再贅述。

4.3.4親緣矩陣的構(gòu)建算法

構(gòu)建親緣矩陣的基本規(guī)則：（1）當(dāng)編號為t的個體的父母均未知時，att=1，ati=ait=0，i=1、2…t-1。（2）當(dāng)編號為t的個體的父母中的一個已知為p，而另一個未知時，att=1，ati=ait=0.5aip，i=1、2…t-1。（3）當(dāng)編號為t的個體的父母均已知（假定為p和q）時，att=1+0.5apq，ati=ait=0.5（aip+aiq），i=1、2…t-1。因此確定親緣矩陣的構(gòu)建算法為：（1）將羅非魚的觀測記錄讀入到數(shù)組arr［］［］當(dāng)中，創(chuàng)建空矩陣A，并使用數(shù)組a［］［］來存儲矩陣的元素。（2）循環(huán)數(shù)組arr［］［］中行記錄即逐個確定觀測記錄中的每一條記錄。在第t（1<=t<=arr.length）條記錄中，假定這條記錄的第a列是父本編號，第b列是其母本編號，若arr［t，a］和arr［t，b］的值均為空，那么判定t是不是第一條記錄，如果是那么a［t，t］=1。否則建立t的循環(huán)變量i（1<=i<=t），a［t，t］=1；a［t，i］=a［i，t］=0。（3）若arr［t，a］和arr［t，b］的值只有一個不為空記為p，那么判定t是不是第一條記錄，如果是那么a［t，t］= 1。否則建立t的循環(huán)變量i（1<=i<=t），a［t，t］=1，a［t，i］=a［i，t］=0.5*a［i，p］。（4）若arr［t，a］和arr［t，b］的值分別為p和q時，那么判定t是不是第一條記錄，如果是那么a［t，t］=1。否則建立t的循環(huán)變量i（1<=i<=t），a［t，t］=1+0.5*a［p，q］，a［t，i］=a［i，t］= 0.5*（a［i，p］+a［i，q］）。

4.3.5混合矩陣的構(gòu)建算法

5　系統(tǒng)的運行結(jié)果及應(yīng)用分析

系統(tǒng)測試環(huán)境：硬件環(huán)境：（1）硬件環(huán)境：PC機(jī)一臺。（2）軟件環(huán)境：Window7操作系統(tǒng)，Microsoft. NET Framework 4.0組件，SQL server2005數(shù)據(jù)庫。

實例測試1：提供2個魚場的5條記錄，并輸入性狀遺傳力值0.33，最終計算的育種值向量值為：（h1 h2 a1 a2 a3 a4 a5）=（114.44 112.72-0.46 2.57 -3.404.70-3.47）。測量記錄數(shù)據(jù)如圖3所示，系統(tǒng)輸出結(jié)果如圖4所示。

圖3　第一次測量記錄數(shù)據(jù)

圖4　第一次實驗輸出結(jié)果

實例測試2：提供3個魚場的6條記錄，并輸入性狀遺傳力值0.33，最終計算的育種值向量值為：（h1h2h3a1a2a3a4a5）=（142.84152.12169.62-2. 463.05-2.12-1.392.150.38）。測量記錄數(shù)據(jù)如圖5所示，實現(xiàn)輸出結(jié)果如圖6所示。

圖5　第二次測量記錄數(shù)據(jù)

圖6　第二次實驗輸出結(jié)果

6　結(jié)語

文章分析了BLUP算法的基本原理，并研究了其在羅非魚單性狀遺傳育種研究中的實現(xiàn)方法。在系統(tǒng)建模型的基礎(chǔ)上，開發(fā)了一套與之相配套的軟件系統(tǒng)，并對軟件實現(xiàn)過程的具體算法進(jìn)行了研究和描述。經(jīng)過軟件的測試和運行表明，對于單性態(tài)的育種案例，本軟件可以順利計算出羅非魚的育種值，因此本軟件有著一定的現(xiàn)實意義。

BLUP算法不僅可以用于單性狀育種值估計，而且同樣可以應(yīng)用于多性狀的育種值估計?，F(xiàn)實生活中影響魚類遺傳育種的性狀常常不僅僅一個。因此多性狀的育種值估計的研究和實現(xiàn)將是下一步工作的主要內(nèi)容。

［1］曾聰，曹小娟，高澤霞，等.團(tuán)頭魴生長性狀的遺傳力和育種值估計［J］.華中農(nóng)業(yè)大學(xué)學(xué)報，2014，33（2）：89-95．

［2］楊瑞飛，張凱，汪平，等.獺兔育種與生產(chǎn)信息管理系統(tǒng)的研發(fā)［J］.四川農(nóng)業(yè)大學(xué)學(xué)報，2014，32（1）：82-86．

［3］劉峰，陳松林，王磊，等.不同牙鲆群體遺傳力和育種值分析［J］.中國水產(chǎn)科學(xué)，2013，20（4）：691-697.

［4］劉瑩瑩，肖兵南，燕海峰，等.利用BLUP法估計荷斯坦種公牛的抗熱應(yīng)激育種值［J］.家畜生態(tài)學(xué)報，2012，33（1）：29-32．

［5］王光瓊，杜天行，未永慶，等.基于SOA構(gòu)架的分布式租車公司管理系統(tǒng)設(shè)計與實現(xiàn)［J］.軟件導(dǎo)刊，2014，13（6）：52-54．

［6］岳紹敏，王國春，董亞則，等.基于Spring MVC框架的電商平臺的研究與設(shè)計［J］.軟件工程師，2014，17（7）：3-7．

［7］鄒燕燕，安虹，崔慧敏，等.矩陣形狀無關(guān)的并行編譯優(yōu)化［J］.小型微型計算機(jī)系統(tǒng)，2014，35（7）：1518-1522．

［8］Seo-Young NOH，Shashi K.GADIA，Haengjin JANG．Comparisons of three data storage models in parametric temporal databases［J］.Journal of Central South University，2013，20（7）：1919-1927．

［責(zé)任編輯］張琴芳

Research and Implementation of Single Trait Animal Breeding Models Based BLUP

HuShi-yang1，YANG Hui-zan2，ZHANG Yong-de2，LINYong2，SUNQuan-liang1
（1.School ofComputerScienceandInformationEngineering，Hezhouuniversity，HezhouGuangxi 542899；2.Guangxi AcademyofFisherySciences，NanningGuangxi 530021）

BLUP algorithm can significantly improve the speed of animal genetic progress，but the calculations are more complex and the data size is quite voluminous.If this complex abstract mathematical algorithm changes into a visualization computer system，it will greatly help to improve productivity for this problem.This paper applies the BLUP algorithm in tilapia breeding value estimation to build the model of the system and use an object-oriented way of thinking to turn the BLUP into executable program.Experimental results show that the development of the software can easily and efficiently estimate breeding values of tilapia.When external factors change，the final result can more accurately reflect the effect of related factors affecting tilapia.

BLUP；algorithm；tilapia；software

TP311.11

A

1673—8861（2016）01—0144—05

2015-02-23

胡世洋（1987-），男，河南南陽人，賀州學(xué)院講師，碩士。主要研究方向：計算機(jī)網(wǎng)絡(luò)與信息管理系統(tǒng)。

廣西水產(chǎn)遺傳育種與健康養(yǎng)殖重點實驗室資助項目（GXKL-AQUA-2013-C-05）、賀州學(xué)院2014年度校級科研項目（2014YBZK09）、賀州學(xué)院2013年度校級科研項目（2013ZXZK14）。