牛曉艷，曹 亮，詹海杰，馮國亮，李燕平，明世清，唐 貴，梁茂文，任克良*

（1.山西農(nóng)業(yè)大學(xué)動物科學(xué)學(xué)院，山西太谷 030801；2.山西省畜禽繁育工作站，山西太原 030001；3.廣靈縣優(yōu)種驢場，山西廣靈 037500）

作為一種經(jīng)濟型草食動物，驢具有耐粗飼、抗逆性強、適應(yīng)性強、易于飼養(yǎng)且飼料來源廣的特點。與牛羊肉相比，驢肉脂肪含量明顯較低，且膽固醇含量低，不飽和脂肪酸含量豐富，可預(yù)防心腦血管疾病。驢皮是熬制阿膠的主要原料，具有補血、滋陰、養(yǎng)肝、益氣、潤肺等功效。養(yǎng)驢業(yè)是傳統(tǒng)畜牧業(yè)的一個重要組成部分，隨著社會的發(fā)展，驢的役用價值逐漸被其肉用、皮用、奶用等經(jīng)濟價值所替代。驢存欄量以每年3.5%的速度下降[1]，地方品種更是面臨種群數(shù)量縮減、種質(zhì)退化甚至瀕臨滅絕的危險境地，亟需進行保護性研究。

廣靈驢主要分布在山西省廣靈、靈邱兩縣及其周圍各縣的邊緣地區(qū)[2]，以廣靈縣南村鎮(zhèn)、壺泉鎮(zhèn)、加斗鄉(xiāng)分布最多，質(zhì)量最好，其品種的形成具有二百多年的歷史[3]。該地區(qū)屬于黃土高原東緣、太行山北部，全年無霜期120 d，年降水量400 mm 左右。廣靈驢俗稱廣靈畫眉驢，具有體格高大、肌肉豐滿、結(jié)構(gòu)勻稱、“五白一黑”等特點，是改良小型驢種的優(yōu)良品種之一，也是全國24 個優(yōu)秀地方驢品種之一。1983 年，廣靈驢被編入全國《馬驢品種志》，2006 年農(nóng)業(yè)部（現(xiàn)為農(nóng)業(yè)農(nóng)村部）第662 號公告確定廣靈驢為國家級遺傳資源保護品種[4]。根據(jù)2019 年底的統(tǒng)計數(shù)字，廣靈驢共存欄1 960頭，其中基礎(chǔ)母驢1 020 頭，公驢52 頭，青年驢和驢駒888 頭[3]。

目前，對于廣靈驢品種特征、生產(chǎn)性能等方面的描述較少，所報道的數(shù)據(jù)集中于20 世紀(jì)八九十年代[2,5]，不能適應(yīng)現(xiàn)代社會的發(fā)展要求，此外，對于其體重體尺之間主成分分析的研究尚未見報道?；诖?，本研究對廣靈縣優(yōu)種驢場現(xiàn)存欄的不同年齡階段的種驢進行了體重體尺的相關(guān)分析和主成分分析，并建立了體重與體尺的最優(yōu)回歸模型，為種驢的測量評估和進一步的選育研究提供了方法。

1 材料與方法

1.1 實驗動物 體重、體尺數(shù)據(jù)通過對山西省廣靈縣優(yōu)種驢場2008—2018 年現(xiàn)存驢群進行測定和數(shù)據(jù)采集而獲得。分析過程中將驢群按照年齡劃分為6 個月～2 歲的幼驢階段、3～4 歲的青年驢階段、5 歲及以上的成年驢3 個階段，其中幼驢29 頭、青年驢20 頭、成年驢51 頭，共計100 頭。

驢群的飼喂采用以谷草為基礎(chǔ)粗飼料，秋季加喂青玉米秸稈（帶穗），精料補充以玉米、黑豆、麩皮為主的傳統(tǒng)飼養(yǎng)方法。

1.2 體重、體尺測定 分別測定不同年齡階段驢群的體重、體高、體長、胸圍、管圍5 個指標(biāo)。其中體重由地秤測定活重，讀取儀表值。體高為鬐甲最高點到地面的垂直距離。體長為肩端到臀端的斜直線距離。胸圍是指從肩胛骨后端引一垂線，繞體軀一周的周長。管圍是指左前管上1/3 處至管骨最細之處的周長。同時記錄采樣個體的毛色、系譜、性別等指標(biāo)。

1.3 統(tǒng)計分析 不同年齡階段驢群體重、體尺的平均值、標(biāo)準(zhǔn)差、標(biāo)準(zhǔn)誤等描述性指標(biāo)采用SPSS 24.0 進行計算。

主成分分析的基本思想是把原有的多個指標(biāo)轉(zhuǎn)化成少數(shù)幾個代表性較好的綜合指標(biāo)，這些少數(shù)指標(biāo)能夠反映原來大部分信息（85%以上），并且各個指標(biāo)間保持獨立，可起到降維和簡化數(shù)據(jù)結(jié)構(gòu)的作用[6]。采用R 語言中cor_pearson 函數(shù)分別對不同年齡組種驢5 個體重、體尺性狀進行相關(guān)性分析，采用corr 函數(shù)對相關(guān)性分析結(jié)果進行顯著性檢驗。

采用R 語言中princomp 函數(shù)進行主成分分析，采用summary 函數(shù)提取主成分的具體信息，采用loadings函數(shù)分析主成分的載荷內(nèi)容，采用screeplot 函數(shù)繪制主成分分析的碎石圖。

采用R 語言中l(wèi)m 進行體高、體長、胸圍、管圍相對于體重的多元回歸方程構(gòu)建；采用step 函數(shù)進行方程的逐步回歸，并參考AIC（赤池信息量）準(zhǔn)則和update函數(shù)進行模型的優(yōu)化，最終得到最優(yōu)的擬合模型。

2 結(jié)果與分析

2.1 年齡組和性別對廣靈驢體重體尺的影響 廣靈驢不同年齡組和性別對其體重體尺的影響見表1～2。由表1 可知，不同年齡階段對廣靈驢體重（F=36.67，P＜0.001）、體高（F=37.85，P＜0.001）、體長（F=37.56，P＜0.001）、胸圍（F=50.41，P＜0.001）、管圍（F=32.31，P＜0.001）均具有極顯著影響。幾個體重體尺指標(biāo)都表現(xiàn)出幼年組極顯著低于青年組和成年組，但青年組和成年組間差異不顯著，表明年齡對驢的生長發(fā)育性狀有極顯著的效應(yīng)，但4 歲之后，驢的生長趨于穩(wěn)定。同時也說明在廣靈驢快速生長階段（4 歲以前）應(yīng)注意飼料中的營養(yǎng)補充，使其快速達到成年體重。

表1 不同年齡組對廣靈驢體重體尺的最小二乘分析結(jié)果（最小二乘均值±標(biāo)準(zhǔn)誤）

表2 不同性別對廣靈驢體重體尺的最小二乘分析結(jié)果（最小二乘均值±標(biāo)準(zhǔn)誤）

2.2 廣靈驢體重和體尺間的相關(guān)性分析 由表3 可知，5 周歲及以上的成年驢體重與體高、體長、胸圍、管圍均表現(xiàn)出極顯著的正相關(guān)關(guān)系，其中體重與體長的相關(guān)性最高（γ=0.648）。各指標(biāo)間均呈現(xiàn)出極顯著的正相關(guān)關(guān)系，其中體長和體高的正相關(guān)關(guān)系最高（γ=0.914），胸圍與管圍的正相關(guān)關(guān)系最低（γ=0.368）。

表3 廣靈驢成年驢體重和體尺的相關(guān)性分析

由表4 可知，3～4 歲的青年驢體重與體高（P＜0.01）、體長（P＜0.05）、胸圍（P＜0.05）之間存在顯著或極顯著正相關(guān)關(guān)系，但體重與管圍之間的相關(guān)關(guān)系不顯著。青年驢的體長和體高相關(guān)性最高（γ=0.908），體重和胸圍之間的相關(guān)性最低（γ=0.497），胸圍與體高、體長以及體重與管圍，胸圍與管圍之間均無顯著的相關(guān)關(guān)系。

表4 廣靈驢青年驢體重和體尺的相關(guān)性分析

由表5 可知，6 個月～2 歲幼驢體重與各體尺指標(biāo)間均存在極顯著的正相關(guān)關(guān)系，且相關(guān)系數(shù)很高。其中，體重與體高的相關(guān)性最高（γ=0.945），體重與管圍的相關(guān)性最低（γ=0.817）。各指標(biāo)間的相關(guān)關(guān)系中，體長和體高的相關(guān)性最高（γ=0.970），胸圍與體重（γ=0.908）、胸圍與體高（γ=0.947）、胸圍與體長（γ=0.931）的正相關(guān)關(guān)系均較高。從上述廣靈驢不同生長發(fā)育階段體重與各體尺指標(biāo)的相關(guān)分析結(jié)果來看，不同階段的相關(guān)性關(guān)系存在差異。

表5 廣靈驢幼驢體重和體尺的相關(guān)性分析

2.3 廣靈驢5 個性狀的主成分分析 由表6～7 可知，6 個月～2 歲階段廣靈驢的特征值為4.61，貢獻率為0.92。第一主成分標(biāo)準(zhǔn)變化量所表達的關(guān)系式為：Y1=-0.45BW–0.46BH–0.46BL–0.45CM–0.42CC。由于貢獻率代表每個特征值所能夠解釋方差的比例，特征向量代表變量預(yù)測某一主成分的回歸系數(shù)。從貢獻率和累積貢獻率來看，協(xié)方差矩陣第一個特征根的貢獻率達到92%，說明第一個主成分代表了原來5 個因素90%以上的信息，所以第一主成分又稱為體重因子。3～4 歲青年驢階段廣靈驢的前3 個特征值分別為3.12、0.75、0.64，前3 個主成分的貢獻率分別為0.62、0.15、0.13，前3 個主成分的累積貢獻率為0.90（＞0.80），所以分別取體重、體高和體長3 個主成分代表原來5 個因素的信息，第一主成分的關(guān)系式為：Y1=-0.42BW–0.52BH–0.49BL–0.38CM–0.41CC；第二主成分的關(guān)系式為：Y2=0.32BW–0.37BH–0.43BL+0.76CM–0.04CC；第 三主成分的關(guān)系式為：Y3=-0.59BW–0.15BH–0.11BL+0.16CM+0.77CC。5 歲以上成年驢階段廣靈驢的前2 個特征值分別為3.33 和0.76，前2 個主成分的貢獻率分別為0.67 和0.15，累積貢獻率為0.82（＞0.80），所以分別取體重和體高代表原來5 個因素的信息。第一主成分的關(guān)系式為：Y1=-0.46BW–0.48BH–0.50BL–0.36CW–0.43CC；第二主成分的關(guān)系式為：Y2=-0.28BW+0.36BH+0.35BL–0.80CW–0.43CC。

表6 廣靈驢不同生長發(fā)育階段相關(guān)矩陣的特征值、貢獻率和累計貢獻率

圖1 為廣靈驢不同年齡階段5 個體重體尺性狀的碎石圖，反映了不同年齡階段廣靈驢不同因子所涵蓋的信息量。從圖1 可以直觀地看出，對于6 個月～2 歲的幼驢，第一主成分涵蓋了絕大多數(shù)的變異，代表了90% 以上的信息；對于3～4 歲的青年驢，前3 個主成分涵蓋了90%以上的變異；而對于5 歲以上的成年驢，前2 個主成分代表了80%以上的變異。

圖1 不同月齡廣靈驢的碎石圖

表7 廣靈驢不同生長發(fā)育階段入選主成分的特征向量

2.4 廣靈驢體重和體尺間最優(yōu)回歸模型的建立 應(yīng)用R語言中l(wèi)m 函數(shù)和step 函數(shù)建立廣靈驢體尺相對于體重的多元線性回歸方程，結(jié)果列于表8。

根據(jù)AIC 準(zhǔn)則，從一組可供選擇的模型中選擇最佳模型時以AIC 最小時為最佳[6]。由表8 可知，在擬合的4 個多元線性方程中，Y2和Y4的AIC 含量較小，但Y2中所采用的2 個體尺性狀（體高、體長）的回歸系數(shù)檢驗不顯著，故該方程不成立，而Y4所采用的2個體尺性狀（體長、胸圍）回歸系數(shù)檢驗均極顯著，說明體長和胸圍對廣靈驢的體重具有決定性作用，擬合的方程可用于預(yù)測廣靈驢體重，最終獲得的擬合方程為：BW=4.35BL+2.18CM－579.39。進一步的，隨機選取30 頭不同年齡階段的廣靈驢體重體尺數(shù)據(jù)分別代入上述公式和經(jīng)驗公式：體重=+25 進行檢驗，結(jié)果如表9 所示。

表8 廣靈驢體重相對于不同體尺性狀的多元線性回歸方程

由表9 可知，通過對本研究得出的擬合方程和經(jīng)驗公式擬合度（R2）的比較，發(fā)現(xiàn)2 個公式的擬合度都隨著廣靈驢年齡的增加而提高，對于成年驢的估計擬合度均可達到0.99 以上，較為理想。但也不難看出，對于成年驢階段和幼驢階段，擬合方程的擬合度均比經(jīng)驗公式高，說明擬合方程的效果更好；對于青年驢階段，2種公式的擬合度接近。另外，由表9 還可以看出，擬合公式擬合效果隨著個體體重的增加而提高，對青年和成年驢階段個體的擬合效果普遍好于幼年階段。

表9 采用擬合方程與經(jīng)驗公式對廣靈驢體重的估計

3 討 論

3.1 廣靈驢體重、體尺分析 廣靈驢屬大型驢種，據(jù)《山西省家畜家禽品種志》記載，廣靈驢的體型外貌特征為：體型接近正方形、體格高大、骨骼粗壯、結(jié)構(gòu)勻稱、耐寒耐粗飼。全身被毛短而粗密，毛色以黑化眉為主，青化眉、灰色、純黑次之。文獻報道廣靈驢成年公驢平均體高為138.4 cm，體長138.5 cm，胸圍為147.2 cm，管圍為17.8 cm；成年母驢平均體高為134.1 cm，體長為131.6 cm，胸圍為146.9 cm，管圍為15.7 cm[7]。與本研究所獲得數(shù)據(jù)相比，體高和體長數(shù)據(jù)一致，但胸圍和管圍較小，這可能反映了不同飼養(yǎng)環(huán)境下廣靈驢生產(chǎn)性能的差異，也和不同種群的選擇強度有關(guān)。

與其他大型驢種相比，廣靈驢體型較晉南驢大，較德州驢小。據(jù)聶鴻瑤等[2]對1979—1980 年在山西省聞喜縣、夏縣測定的晉南驢的體重體尺數(shù)據(jù)，成年公驢體重為127.3 cm，母驢為124.3 cm；體長公驢為123.2 cm，母驢為125.2 cm；胸圍公驢為136.3 cm，母驢為135.6 cm；管圍公驢為15.2 cm，母驢為14.8 cm。而馮志華[8]對100 頭晉南驢（其中公驢28 頭，母驢72 頭）體重體尺進行測定，得到公驢平均體高為135.5 cm，母驢為134.5 cm；體長公驢為132.3 cm，母驢為128.4 cm；胸圍公驢為132.3 cm，母驢為128.4 cm；管圍公驢為17.2 cm，母驢為15.9 cm，均略小于廣靈驢，但差異不大。農(nóng)博[9]描述了德州驢的品種特征，成年公驢體高為142 cm，成年母驢體高為140 cm；體長公驢可達143.6 cm，母驢達137 cm；胸圍公驢為152.8 cm，母驢為160 cm；管圍公驢為18.7 cm，母驢為16.4 cm，均高于廣靈驢或與之相當(dāng)。

不同的年齡組因素對廣靈驢體重體尺性狀均具有極顯著影響，表現(xiàn)為青年組和成年組的個體極顯著高于幼年組個體，但青年組和成年組之間差異不顯著。這與楊莉等[10]研究結(jié)果一致，他們采用3 種曲線對德州驢的生長發(fā)育規(guī)律進行擬合分析，得到公、母驢達到成熟時的拐點月齡分別為2.43 和1.95，說明在2 歲以后，驢的生長速度趨于平穩(wěn)，接近成年水平。不同性別因素對廣靈驢體重、體高、體長、胸圍均無顯著影響，但公驢管圍極顯著大于母驢，這與其他一些大型驢品種相似。

3.2 廣靈驢體重、體尺的相關(guān)性分析 本研究結(jié)果表明，不同發(fā)育階段的廣靈驢體重體尺之間均存在顯著的正相關(guān)，且成年驢階段體重與體長相關(guān)性最高（γ=0.648，P＜0.001），而青年和幼年階段體重與體高的相關(guān)性最高（γ=0.606，P＜0.001；γ=0.945，P＜0.001）。這與楊莉等[10]報道的公驢和母驢體重與胸圍的相關(guān)性最高不一致（公驢：γ=0.963；母驢：γ=0.971）。Pierre 等[11]計算了西非布基納法索等4 個地區(qū)1 001 頭公驢和351 頭未懷孕母驢體重和體尺的相關(guān)關(guān)系，也得到體重與胸圍具有最強的相關(guān)性（公驢：γ=0.89；母驢：γ=0.90）。但陸會寧等[12]在對80 只蘭州大尾羊× 小尾寒羊F1代體重和體尺的相關(guān)性研究中發(fā)現(xiàn)，體重和體長（γ=0.875）的相關(guān)性與體重和胸圍的相關(guān)性接近（γ=0.907）。說明體重和各體尺指標(biāo)間的相關(guān)性與物種和樣本量有關(guān)。本研究與其他研究相比樣本量較少，可能會對不同指標(biāo)間的相關(guān)關(guān)系計算結(jié)果產(chǎn)生影響。

3.3 廣靈驢體重、體尺的主成分分析 主成分分析法可以在力求數(shù)據(jù)信息丟失最少的情況下，對高維的變量空間進行降維，研究指標(biāo)體系的少數(shù)幾個線性組合，而且這幾個線性組合所構(gòu)成的綜合指標(biāo)能夠盡可能的保留原來指標(biāo)的信息。目前，國內(nèi)尚未見對地方驢品種體重體尺性狀進行主成分分析的報道，通過本研究結(jié)果可以看出，不同年齡階段對廣靈驢生長發(fā)育起到主要作用的性狀不同。在幼年階段，體重的影響最大；青年階段，體重、體高和體長則代表了其主要體型特征；成年階段，體重和體高最具代表性。這說明在廣靈驢的選育過程中，應(yīng)根據(jù)不同階段性狀之間的關(guān)聯(lián)性，并兼顧其他體尺性狀，對育種目標(biāo)進行優(yōu)化，保留地方品種的體型外貌特征。

3.4 廣靈驢體重、體尺的回歸分析 采用多元回歸分析中的逐步回歸分析法結(jié)合AIC 準(zhǔn)則得出廣靈驢體重與體尺性狀的最優(yōu)回歸線性模型，發(fā)現(xiàn)廣靈驢的體重主要與其體長和胸圍密切相關(guān)。楊陽等[13]對3～9 歲178 頭德州驢進行了體重和體尺的回歸分析，建立了德州驢體重（Y）與其體長（X1）、胸圍（X2）和尻長（X3）的多元線性回歸方程：Y=2.6X1+1.8X2+3.0X3–500.2，與本研究結(jié)果相似，但本研究未對廣靈驢的尻長數(shù)據(jù)進行采集，所以最終的方程僅包括了體長和胸圍2 個性狀。肖海霞等[14]采用R 軟件包對62 頭吐魯番驢進行體重和體尺分析，并建立了最優(yōu)的回歸模型，得出回歸方程為：Y（體重）=0.71X1（體高）+1.54X2（體長）+2.84X3（胸圍）–447.36，與本研究結(jié)果相似。

在實際生產(chǎn)中，對驢的體重測定通常較為困難，一些參考資料給出了如下計算公式：體重=+25，此公式的運算較為繁瑣難記[7]，本研究給出的體重估計公式將體重作為體長和胸圍的線性函數(shù)進行擬合，校正擬合度達到0.85 以上，更加簡單實用。比較經(jīng)驗公式與線性擬合方程的擬合度可知，線性擬合方程的擬合度在廣靈驢生長發(fā)育的各個階段均高于經(jīng)驗公式或與之接近，這說明線性方程的擬合和預(yù)測更加準(zhǔn)確，可代替?zhèn)鹘y(tǒng)經(jīng)驗公式在生產(chǎn)中進行推廣使用。

本研究提出的擬合公式是通過對廣靈驢體重體尺性狀進行分析得出的，通過與其他地方品種比較可知，不同品種具有不同的體型特征和生長發(fā)育特點，因此此公式僅適用于廣靈驢體重的預(yù)測。另外，雖然此公式是通過對6 個月～5 歲（幼年～成年）不同生長階段個體數(shù)據(jù)的擬合而得到，但也不難發(fā)現(xiàn)，公式估計值的擬合度隨著個體年齡的降低而迅速下降，所以本公式用于估計青年和成年驢的體重更加準(zhǔn)確適用。最后，與其他研究者的估計方程相比，本研究未考慮尻長、十字部高等數(shù)據(jù)，在表現(xiàn)形式上略有差異，今后，可進一步對廣靈驢的體尺數(shù)據(jù)進行收集，并擴大采樣規(guī)模，得到更為完善的公式。

4 結(jié) 論

本研究基于R 語言對不同年齡階段廣靈驢的體重體尺進行了比較，計算了各體重體尺性狀間的相關(guān)性，并對各性狀進行了主成分分析，最后得出了體重與體尺間的最優(yōu)回歸線性方程。研究結(jié)果表明：①幼年和青年階段體重與體高的相關(guān)性最高（γ=0.945，γ=0.606），成年階段體重和體長的相關(guān)性最高（γ=0.648）；②主成分分析結(jié)果表明，幼驢階段第一主成分的特征值為4.61，貢獻率為0.92；青年驢階段第一、第二、第三主成分的特征值分別為3.12、0.75、0.64，累積貢獻率達到0.90；成年驢階段第一、第二主成分的特征值分別為3.33 和0.76，累積貢獻率為0.82；③通過逐步回歸分析得到廣靈驢體重和體尺的最優(yōu)回歸方程為：BW=4.35BL+2.18CM–579.39。