林雨濃，馬萬欣，包玲玲，何曙光，石順利，張秋生，趙澈勒格日，高會江，李俊雅*，王澤昭*

（1.中國農(nóng)業(yè)科學院北京畜牧獸醫(yī)研究所，北京 100193；2.內(nèi)蒙古自治區(qū)通遼市畜牧業(yè)發(fā)展中心，內(nèi)蒙古通遼 028000）

肉牛產(chǎn)業(yè)是內(nèi)蒙古自治區(qū)經(jīng)濟發(fā)展的支柱產(chǎn)業(yè)，是農(nóng)牧民增收的重要途徑?！翱茽柷吲！笔且晕鏖T塔爾牛為父本，蒙古牛及三河牛和蒙古牛的雜種牛為母本而培育成的乳肉兼用品種。作為通遼地區(qū)地方畜牧業(yè)的重要資源，具有適應(yīng)性強，宜牧，耐粗飼，耐寒，抗病力強等特點。然而，隨著肉牛市場需求的增加和養(yǎng)殖業(yè)的發(fā)展，“科爾沁?！痹谏a(chǎn)性能和經(jīng)濟效益方面遇到了瓶頸。為了提升群體肉用性能，當?shù)剞r(nóng)業(yè)部門利用北美肉用型西門塔爾牛為父本，本地“科爾沁牛”為母本進行雜交改良，歷經(jīng)30 年，經(jīng)過“雜交改良、橫交固定和選育提高”三個階段，培育出肉用型“科爾沁肉牛”群體。

近年來，群體遺傳背景分析成為評估和研究動物群體遺傳結(jié)構(gòu)的重要工具。隨著SNP 芯片和測序技術(shù)的發(fā)展，有助于人們更全面地了解不同群體之間的遺傳多樣性和遺傳聯(lián)系。Eva M.Strucken 等[1]采用777K 芯片對印度本地的15 個牛種進行遺傳多樣性與群體結(jié)構(gòu)研究，發(fā)現(xiàn)印度瘤牛比印度本地普通牛有更高的遺傳多樣性，他們推斷本地普通牛較小的有效群體導致了群體基因頻率受到遺傳漂變的影響更大，而低有效群體數(shù)量增加了近交衰退發(fā)生的概率，使其在選擇的過程中遺傳多樣性降低速度快于印度瘤牛。在動物育種過程中，保持一定的遺傳多樣性有利于提高群體的環(huán)境適應(yīng)性[2]。Alejandro Amaya 等[3]對約2.7 萬頭哥倫比亞西門塔爾牛（1975—2017年）的數(shù)據(jù)進行了分析，發(fā)現(xiàn)該種群的近交系數(shù)由5.06%下降到了2.25%，認為是得益于歐洲和美系西門塔爾牛之間的基因交流，才使該種群的遺傳多樣性得以提升。

我國疆域遼闊、資源豐富，但目前由于經(jīng)濟效益等問題致使國外種源占有率居高不下。物種作為戰(zhàn)略性資源、是農(nóng)業(yè)的基石，但是目前關(guān)于國內(nèi)地方品種的研究存在著許多空白和滯后。馬均等[4]研究發(fā)現(xiàn)，我國自主培育的專門化肉牛新品種—華西牛與蒙古牛、三河牛和其他家系西門塔爾牛出現(xiàn)不同的比例的群體結(jié)構(gòu)組成，形成了自己獨特的遺傳結(jié)構(gòu)特征。劉正喜等[5]在2022 年使用全基因組重測序數(shù)據(jù)對渤海黑牛的群體結(jié)構(gòu)、遺傳多樣性、親緣關(guān)系分析后發(fā)現(xiàn)，渤海黑牛與郟縣紅牛和魯西牛的關(guān)系密切，該研究為逐年減少的渤海黑牛群體的種質(zhì)資源保護和品種發(fā)展提供了關(guān)系的理論基礎(chǔ)。

本研究目的是分析“科爾沁肉?！?、華西牛和美系西門塔爾牛這三個群體之間的遺傳關(guān)聯(lián)，評估“科爾沁肉?！备牧加N在遺傳層面的影響，旨在為“科爾沁肉牛”新品種培育提供遺傳背景支持。

1 材料與方法

1.1 試驗群體及SNP 分型

“科爾沁肉牛”（Kerqin beefcattle，“KBC”），來自內(nèi)蒙古科爾沁牛業(yè)股份有限公司的437 頭核心群牛只；美系西門塔爾牛（American Simmental cattle，ASC）來自于通遼京緣種牛繁育有限責任公司的25 頭種牛；華西牛（HuaXi cattle，HXC），來自內(nèi)蒙古奧科斯牧業(yè)有限公司的55 頭核心群牛只。

每個樣本使用含有EDTA 抗凝劑的真空采血管采血5mL 送樣檢測，并取2mL 血樣進行DNA提取。使用Illumina BovineHD BeadChip 芯片進行變異檢測，該芯片含有覆蓋?；蚪M的777962個標記位點。

1.2 數(shù)據(jù)質(zhì)控

使用Genome Studio 軟件進行基因分型。通過PLINK v1.9[6]軟件提取出常染色體后進行數(shù)據(jù)質(zhì)量控制（Quality control，QC）。QC 順序如下：①刪除等位基因檢出率（Call rate）小于90%的位點；②刪除最小等位基因頻率（Minor allele frequency，MAF）小于0.01 的位點；③刪除哈代溫伯格平衡檢驗P 值小于1×10 箒6 的SNP 位點；④刪除SNP 缺失率大于10%的樣本個體。通過QC 后數(shù)據(jù)用于后續(xù)分析。單獨使用PLINK 中的--indep-pairwise 參數(shù)對所有個體的連鎖位點進行過濾，剩余位點進行群體結(jié)構(gòu)分析。過濾后的個體數(shù)見表1。

表1 試驗樣本數(shù)量

1.3 遺傳多樣性分析

對QC 后數(shù)據(jù)統(tǒng)計群體期望雜合度（Expected heterozygosity，He）、群體觀察雜合度（Observed heterozygosity，Ho）、MAF 分布情況。 MAF 分布情況是按照不同MAF 大小將SNP 標記分為5 組，分別是（MAFSNP<0.1）、（0.1≤MAFSNP＜0.2）、（0.2≤MAFSNP＜0.3）、（0.3≤MAFSNP＜0.4）和（0.4≤MAFSNP），統(tǒng)計3 個群體SNP 標記的分布情況。

1.4 群體分析

使用GCAT 軟件進行了群體主成分（Principal component analysis，PCA）計算，輸出數(shù)據(jù)導入到R 中，使用ggplot2 包[7]進行可視化。隨后使用鄰接法（Neighbor-Joining，NJ）構(gòu)建進化樹[8]：從成對FST 值的距離矩陣生成鄰居網(wǎng)絡(luò)樹，并在Splitstree 4.14.5[9]中實現(xiàn)建樹。使用Admixture v.1.3 軟件的默認設(shè)置進行群體結(jié)構(gòu)分析。設(shè)置了10 倍交叉驗證選擇最合適的K 值[10]。使用PopLDdecay 軟件[11]計算了所有數(shù)據(jù)的LD 衰減距離，并繪制了LD 衰減圖。

2 結(jié)果

2.1 群體遺傳多樣性

由表2 可知，經(jīng)過質(zhì)控后的“KBC”群體具有最高的平均MAF 值（0.28），ASC 群體的平均MAF 值最低（0.23）。平均He 的范圍取值從0.34（HXC）到0.37（“KBC”）。平 均Ho 最高為“KBC”（0.37），最低的是ASC（0.36），HXC和ASC 的平均Ho 略大于He，“KBC”的Ho 和He 幾乎相等。

表2 群體遺傳多樣性參數(shù)

由圖1 可知，“KBC”群體被過濾掉了最多的遺傳標記?！癒BC”群體過濾掉了81778 個MAF 小于0.01 的標記。QC 后，ASC 群體擁有最多的低頻標記（MAF＜0.1）。隨著MAF 的增大，HXC 和“KBC”群體內(nèi)中高MAF 的SNPs 多于ASC 群體。在三個群體中，“KBC”具有最多的MAF≥0.3 的SNP 個數(shù)，具體數(shù)量為（413400），其次是HXC（363989），最少的是ASC（344174）。

圖1 不同群體的MAF 分布

2.2 PCA 和進化樹

由圖2（a）可知，三個群體主要聚集成為兩個簇，將“KBC”群體和HXC 群體分開。此外，“KBC”群體與ASC 群體中部分個體出現(xiàn)明顯的聚集和重疊。根據(jù)群體在圖中的空間位置，可以判斷出“KBC”和HXC 的遺傳距離較遠，和ASC 之間遺傳距離較近。

圖2 三個牛群體的主成分分析結(jié)果和鄰居連接系統(tǒng)發(fā)育樹

為了觀察三個群體在進化距離上的相對位置，本研究構(gòu)建了環(huán)形進化樹（圖2b）。“KBC”群體和HXC 群體分別位于進化樹的起點位置和末端位置，表明這兩個群體的進化距離相對較遠。ASC 群體和HXC 群體處于兩個比較靠近的分支，暗示這兩個群體經(jīng)歷了不同的進化路徑，但是其分化路徑相對較近。ASC 群體處于“KBC”群體的一個分支下，表明這兩個群體可能具有共同的祖先，使得在基因?qū)用嫔暇邆湟欢ǖ南嗨菩浴?/p>

2.3 種群結(jié)構(gòu)分析

通過假設(shè)群體結(jié)構(gòu)數(shù)的先驗值K，推斷群體內(nèi)亞群個數(shù)。通過對比不同祖先群體遺傳成分比例觀察“KBC”、ASC 和HXC 群體間遺傳相似性。由圖3 可知，在K=2 時，“KBC”的祖先比例和HXC 相似，少數(shù)個體祖先比例與ASC 相似；當K=3 時，“KBC”群體和ASC 群體部分個體祖先比例相似，此時HXC 群體出現(xiàn)與“KBC”群體和ASC 群體截然不同的祖先成分比例；隨著K 值不斷增大，“KBC”群體和HXC群體差異的趨勢逐漸增大；當K=5 時，“KBC”出現(xiàn)了更為復(fù)雜的遺傳比例；K=2～4 時，ASC群體的祖先組成比例變化不明顯，“KBC”的祖先比例變化復(fù)雜。

圖3 當K=2～5 時，基于LD 過濾的三種牛種群的SNPs 的結(jié)構(gòu)分析

2.4 連鎖不平衡

如圖4 所示，三個群體的LD-decay 曲線中“KBC”群體衰減速度最快，ASC 的LD 衰減速度較慢，HXC 群體的LD 衰減速度位于兩者之間。三個群體的遠端SNP 之間的連鎖水平出現(xiàn)差異，在≥200Kb 之后并未衰減到相近的r2范圍。隨著距離逐漸增大，三個群體的LD 衰減曲線逐漸平緩，r2結(jié)果始終保持是r2ASC≥r2HXC≥r2KBC。

圖4 三個群體的連鎖不平衡衰減圖

3 討論

3.1 遺傳多樣性參數(shù)

遺傳多樣性是生物多樣性的核心維度之一，反映了群體內(nèi)部基因的多樣性水平。它在維持種群的健康、適應(yīng)性和進化潛力方面起著關(guān)鍵作用[12,13]。目前，通過基因芯片和測序技術(shù)輔助分析種群遺傳進展成為熱點研究內(nèi)容[14]，通過評估種群遺傳多樣性對于制定未來的育種目標至關(guān)重要[2]。本研究對“KBC”群體、HXC 群體和ASC群體的基因型數(shù)據(jù)進行分析，通過對群體遺傳多樣性的統(tǒng)計，以期發(fā)現(xiàn)“KBC”在育種過程中的遺傳變化。研究發(fā)現(xiàn)，“KBC”的MAF、He 和Ho 分別為0.28、0.37 和0.37，在三個群體中均為最高，并且“KBC”群體He 和Ho 的一致性在三個群體中為最高。這一結(jié)果反映出“KBC”群體內(nèi)的遺傳多樣性比ASC 和HXC 群體較高，同時暗示著近親交配的水平較低[15]。種群遺傳多樣性越高，意味著進化潛力越大，對環(huán)境變化的反應(yīng)能力越強[16]。此外，人工選擇也會極大的影響群體的多樣性[17]。“KBC”是以“科爾沁?！保ㄈ优！廖鏖T塔爾牛）為母本，以引進的美系西門塔爾牛為父本雜交培育而成的群體，與ASC 相比，其培育期間的選擇強度相對較低，使其還保持著一定的遺傳改良潛力。Maiorano AM 等[18]研究發(fā)現(xiàn)，親本自交系雜交產(chǎn)生的雜種優(yōu)勢賦予了后代更高的遺傳多樣性，并加強了雜種個體優(yōu)勢性狀表現(xiàn)[19,20]?！癒BC”群體在遺傳改良的過程中，保留了“科爾沁?！比后w的適應(yīng)性強的同時，借助ASC 群體的基因交流提升了肉用性能。在哈迪-溫伯格平衡下，“KBC”的Ho 大于He，說明雜交系統(tǒng)仍然是擴大品種遺傳變異的重要方法。

3.2 PCA 和進化樹

PCA 是利用解釋群體最大遺傳方差的特征值向量對群體進行分離[21]。本結(jié)果表明，“KBC”和HXC 群體各自聚集成了兩個簇并保持了一定的距離。同時，“KBC”群體中部分個體與ASC群體出現(xiàn)明顯的聚集，這些“重疊”，標志著個體之間的基因型具有一定的相似性。這種相似性是由于群體之間存在相似的遺傳背景或近親關(guān)系導致了其共享某些遺傳特征或具有相似的基因型組合。

通過構(gòu)建進化樹，發(fā)現(xiàn)“KBC”群體和ASC群體進化距離較近，與HXC 群體距離較遠。需要注意的是，環(huán)形進化樹是無根樹，僅提供了群體之間遺傳和進化距離的相對關(guān)系。群體在進化樹上的位置不代表起源和進化方向，而是反映了在進化距離上的相對位置。最初“KBC”群體在培育中使用了ASC 改良其肉用性能，在基因組上留下了選擇足跡（Selective Footprint）[22]。經(jīng)過長期的人工選擇，增加了“KBC”群體和ASC 群體的遺傳相似性，因此其在進化樹上的位置處于同一分支。后續(xù)研究如果要全面地理解群體間的關(guān)系和演化過程，需要分析表型特征信息或分子鐘估計[23]等進行更深入的研究。

3.3 群體結(jié)構(gòu)分析

群體遺傳結(jié)構(gòu)分析基于遺傳標記信息將種群分為若干亞群，對于研究群體和個體間的遺傳關(guān)聯(lián)具有重要意義[24]。揭示了群體間的遺傳交流程度、遷移歷史以及潛在的遺傳壁壘[25]。研究發(fā)現(xiàn)，在K=2～3 時“KBC”群體內(nèi)部分個體與ASC 祖先組成比例十分相似，在K=2～4 時，ASC 群體的組成比例比較穩(wěn)定，而“KBC”群體祖先比例變化明顯。Ocampo 等[26]發(fā)現(xiàn)不同品種之間的分化模式和遺傳關(guān)系與每個品種的進化歷史高度一致。“KBC”群體和HXC 群體在雜交改良階段都引入過西門塔爾的血統(tǒng)，因此，當K值比較小的時候，更容易地捕捉到群體之間的主要遺傳差距。在“KBC”群體內(nèi)，不同個體之間的祖先組成比例也在不斷變化，一方面是其原始群體本身的多態(tài)性所引起；另一方面，經(jīng)過不斷的雜交改良，致使“KBC”群體存在多祖先的基因滲入。

本研究中，在假設(shè)K 值范圍內(nèi)并未出現(xiàn)交叉驗證（Cross-Validation，CV）錯誤率的拐點，可能是研究中ASC 群體和HXC 群體使用的樣本數(shù)量較小，捕捉三個群體間的遺傳分化效果較差。當K 大于4 時，三個群體祖先比例出現(xiàn)明顯不同，且個體祖先比例變得相對混亂。這是因為，隨著K 值增大，模型捕捉群體之間細微的遺傳差異的靈敏度也隨之提高。

3.4 LD 衰減

通過對比群體之間的連鎖不平衡水平，對于合理利用種質(zhì)資源至關(guān)重要。三個群體中，“KBC”群體衰減的最快，ASC 群體的LD-decay衰減的最慢。而人工選擇的效果可以反映在每個群體的基因組連鎖不平衡水平上[27]。在物種遺傳改良過程中，在高強度的選擇下，家畜群體中有利的染色體片段頻率迅速增加[28]，并促進群體中LD 水平的增加[29,30]。當群體中LD 水平增加，提高后續(xù)選擇的準確性，可以更快地改善目標性狀。HXC 群體已經(jīng)于2021 年12月通過國家畜禽遺傳資源委員會審定，被認定為肉牛新品種。而“KBC”群體目前還處于新品種培育的過程中。因此，“KBC”群體內(nèi)遠端標記的連鎖效果比其他兩個群體較差。連鎖不平衡分析可以反映相應(yīng)的獨特育種選擇歷史和種群結(jié)構(gòu)[27,31,32]，同時，也對后續(xù)基因組信息挖掘提供指示。值得注意的是，在研究中隨著標記之間的距離逐漸增大，三個群體的LD 衰減逐漸平緩，但是一直存在差距。后續(xù)如要進行GWAS 分析時，為了保證準確性，“KBC”群體可能會需要比HXC 和ASC 更多的基因組信息。

4 結(jié)論

本研究使用遺傳信息數(shù)據(jù)評估了“KBC”、HXC 和ASC 群體的遺傳多樣性和種群結(jié)構(gòu)?！癒BC”群體內(nèi)中高頻（MAF≥0.3）的標記數(shù)量最多，且具有最高的遺傳多樣性和最低的連鎖不平衡衰減速度。在進化方向上，“KBC”群體與HXC 群體的距離更遠，但部分個體的種群組成與ASC 群體分化較差。通過進行群體遺傳背景分析，發(fā)現(xiàn)“科爾沁肉?！苯?jīng)過雜交改良、橫交固定、選育提高三個階段的持續(xù)選育，在遺傳方面（分子水平）與華西牛和美系西門塔爾牛品種存在明顯差異，已具備形成新培育品種的先決條件。