張成龍，鄭浪漫，劉春潔，高慶華，劉書東*

(1.塔里木大學(xué)動(dòng)物科學(xué)與技術(shù)學(xué)院，阿拉爾 843300；2.新疆生產(chǎn)建設(shè)兵團(tuán)塔里木畜牧科技重點(diǎn)實(shí)驗(yàn)室，阿拉爾 843300)

新疆南疆和田地區(qū)位于塔克拉瑪干沙漠南緣，昆侖山北麓。年降水量在21～146 mm，年均蒸發(fā)量達(dá)2 600 mm，每年沙塵天氣220 d以上，其中濃浮塵(沙塵暴)天氣在60 d以上，生態(tài)環(huán)境脆弱。該地區(qū)漢朝已有養(yǎng)羊業(yè)記載，地方品種綿羊與和田羊、策勒黑羊兩個(gè)品種的相同點(diǎn)：鼻梁微隆，耳大、下垂，體格較小，舍飼狀態(tài)常年發(fā)情。不同點(diǎn)：和田羊被毛白色，以產(chǎn)地毯毛為主；策勒黑羊被毛黑色，以產(chǎn)羔皮和多胎著稱。經(jīng)過長期選擇能抵抗多變的沙漠環(huán)境，是我國獨(dú)特的優(yōu)勢種質(zhì)資源綿羊群體。

家畜在歷史馴化和品種培育及改良過程中受到強(qiáng)烈的自然選擇和人工選擇，在家畜基因組上留下相應(yīng)的遺傳印記，具體表現(xiàn)為基因組上連鎖不平衡程度增加或降低和受選擇位點(diǎn)及其附近遺傳多樣性降低。2004年美國、英國和新西蘭的科學(xué)家公布了第三代高密度的綿羊遺傳連鎖圖[1]。綿羊基因組學(xué)的研究是開展綿羊分子遺傳標(biāo)記研究與應(yīng)用的基礎(chǔ)[2-3]，2014年，華大基因和中國科學(xué)院的研究者公布了綿羊基因組序列，以雌性陶賽特為綿羊材料，通過二代測序技術(shù)與基因組圖譜技術(shù)相結(jié)合，測序深度為75X，其中Contig N50為39 959 bp，Scaffold N50為2 231 873 bp，共包括20 908個(gè)基因，141個(gè)大結(jié)構(gòu)改變，近一萬個(gè)拷貝數(shù)變異，大約有1 000萬個(gè)SNPs，基因組全長2.61 Gb[4]。綿羊基因組測序的完善為功能基因定位、相關(guān)分子標(biāo)記的研究提供了遺傳信息。Li等[5]選擇99只Finnsheep不同顏色的綿羊群體，通過全基因組關(guān)聯(lián)分析和選擇性清掃相結(jié)合的方法，獲得了毛色相關(guān)的基因(TYRP、ASIP、MITF)。Kijas等[6]基于Ovine SNP50 BeadChip通過Fst方法對全世界74個(gè)綿羊品種進(jìn)行了品種間的選擇性清掃檢測，共獲得181個(gè)基因，其中KIT、MITF和色素沉積相關(guān)；ASIP和黑、白色或黑白相間的毛色相關(guān)；NPR2、HMGA2、BMP2和動(dòng)物的身體大小、肌肉肥大等生產(chǎn)性狀相關(guān)；FGF5和毛質(zhì)性狀相關(guān)等。Zhao等[7]用3個(gè)典型的中國本土綿羊品種，基于 Ovine Infinium HD SNP BeadChip 進(jìn)行選擇性清掃，其中25 個(gè)基因組區(qū)域被選擇，篩選出73 個(gè)候選基因，這些基因與尾部脂肪的形成、繁殖、生長性狀等有關(guān)。沙漠環(huán)境抗逆性是在沙漠極端環(huán)境下生物生長發(fā)育和生產(chǎn)性能表現(xiàn)正常的能力，外來引入品種對沙漠環(huán)境適應(yīng)性差[8]，生產(chǎn)性能下降[9]，Yang 等[10]使用重測序方法從不同的生態(tài)環(huán)境中獲得77個(gè)綿羊基因組，并確定綿羊與沙漠(ANXA6、GPX3、GPX7和PTGS2)和高原(IFNGR2、MAPK4、NOX4、SLC2A4和PDK1)等極端環(huán)境相適應(yīng)有關(guān)的相關(guān)基因。

和田地區(qū)存在和田羊和策勒黑羊兩個(gè)地方綿羊品種，不斷的自然選擇和人工選擇在這兩種綿羊基因組上留下了大量的遺傳印記。和田羊與策勒黑羊均具有適應(yīng)沙漠環(huán)境的能力，但是關(guān)于兩個(gè)品種的基因組差異研究較少。本研究旨在利用基因組選擇性清掃檢測和田羊和策勒黑羊的基因組差異，尋找適應(yīng)南疆沙漠環(huán)境的優(yōu)異分子標(biāo)記，為南疆土著綿羊保護(hù)和育種提供理論基礎(chǔ)。

1 材料與方法

1.1 樣品采集

本試驗(yàn)于和田地區(qū)策勒縣綠源種畜場隨機(jī)選取策勒黑羊41只，于和田羊保種場隨機(jī)選取和田羊84只，均為成年母羊，采樣時(shí)間為2018—2021年。藏綿羊37只，數(shù)據(jù)來自International Stroke Genetics Consortium(ISGC)(http://www.Sheephapmap.org)。

1.2 DNA提取和基因組芯片制備

利用耳缺鉗在每個(gè)個(gè)體的耳緣處剪下直徑3 mm的組織塊，放入己滅菌含有75%酒精的離心管。采用苯酚-氯仿法提取DNA，TAE電泳和核酸儀檢測分析后，進(jìn)行Ilumina Ovine SNP50 BeadChip制備。BeadScan軟件對圖像進(jìn)行轉(zhuǎn)換，利用Genome Studio軟件對數(shù)據(jù)結(jié)果進(jìn)行處理。

1.3 質(zhì)控分析

采用Plink1.07[11]軟件進(jìn)行數(shù)據(jù)質(zhì)控，質(zhì)控標(biāo)準(zhǔn)為：去掉樣本中檢出率小于95%的個(gè)體樣品，SNPs檢出率小于90%、最小等位基因頻率小于5%、哈代-溫伯格平衡檢驗(yàn)P值小于1×10-6。

1.4 主成分分析

使用GCTA軟件[12]對過濾后的SNPs數(shù)據(jù)進(jìn)行PCA分析，該主成分分析基于個(gè)體間的遺傳相關(guān)性識別，代表種群結(jié)構(gòu)的主要成分[13]。

1.5 基因組選擇信號的計(jì)算

本試驗(yàn)以和田羊和策勒黑羊?yàn)檠芯繉ο螅捎肏apFLK方法，選擇藏綿羊作為外群。HapFLK軟件[14]基于群體間單體型頻率來檢測選擇差異。使用Arlequin軟件[15]來構(gòu)建雷諾(Reynolds)親屬矩陣，在fastPHASE模型中，親緣關(guān)系矩陣假設(shè)有10個(gè)群集(cluster)，使用 20 次迭代的期望最大化算法(EM)對所有染色體進(jìn)行 FLK 分析。得出結(jié)果后，利用以下公式校正HapFLK統(tǒng)計(jì)量：

其中，mean(hapFLK)為hapFLK統(tǒng)計(jì)量的平均數(shù)，sd(hapFLK)為hapFLK統(tǒng)計(jì)量的標(biāo)準(zhǔn)偏差。校正以后的HapFLK統(tǒng)計(jì)量近似服從正態(tài)分布，利用Excel 2016的函數(shù)1-NORM.DIST(HapFLKadj,0,1,TRUE)得到每個(gè)位點(diǎn)的P值,若-lgP≥3認(rèn)為該位點(diǎn)受到了選擇[16]。

1.6 GO和KEGG富集分析

所有獲得SNPs基于Oar_v4.0(https://www.sheephapmap.org/)進(jìn)行注釋?；虻墓δ芊治鰠⒖糔CBI數(shù)據(jù)庫(http:/www.ncbi.nlm.nih.gov/gene)、OMIM數(shù)據(jù)庫(http://www.ncbi.nlm.nih.gov/omim)以及相應(yīng)的參考文獻(xiàn)。使用DAVID工具[17](http：//david.Abcc.ncifcrf.gov/)進(jìn)行GO 和 KEGG(Kyoto Encyclopedia of Genes and Genomes)分析。每個(gè)富集條目的重要程度由P值衡量，P值遵循超幾何分布，公式為：

其中，N為總基因數(shù)；n為N中差異表達(dá)基因的總數(shù)；M為N中屬于某個(gè)富集條目的基因個(gè)數(shù)；i為n中屬于某個(gè)富集條目的基因個(gè)數(shù)。為了避免結(jié)果假陽性，對富集分析結(jié)果P采用多重檢驗(yàn)(False Discovery Rate, FDR)進(jìn)行校正。FDR的公式為：

其中，P是原始P值，n是測試次數(shù)，rankP是特定原始P值的等級。FDR≤0.05時(shí)，符合此條件的GO術(shù)語和途徑被定義為候選基因顯著富集。

2 結(jié) 果

2.1 主成分分析

如圖1所示，藏羊與策勒黑羊、和田羊群體可以由第一主成分分為兩組，前兩個(gè)主成分(PC1和PC2)可以將所有樣本分為4個(gè)子組，PC2與PC3兩個(gè)主成分也大致將所有樣本分為4個(gè)組。其中藏羊群體與和田羊、策勒黑羊群體明顯區(qū)分為兩個(gè)群體，和田羊群體與策勒黑羊群體基因組差異很小。

A、B、C分別為PC1-PC2、PC1-PC3和PC2-PC3所對應(yīng)的區(qū)分結(jié)果；PCA結(jié)果將所有個(gè)體大致分為3組，其中和田羊與策勒黑羊遺傳關(guān)系較近

2.2 選擇性清掃結(jié)果

采用HapFLK結(jié)果的前1%作為有效結(jié)果，共計(jì)470個(gè)SNPs，結(jié)果如圖2所示，OAR1、OAR3、OAR6號染色體上的選擇信號較強(qiáng)烈。通過基因注釋(距離該位點(diǎn)最近的基因)找到這些位點(diǎn)相對應(yīng)的候選基因，并去除重復(fù)與無效的基因，共得到了167個(gè)基因，其中綿羊毛生長相關(guān)基因有KRT82、KRT80、KRT1、KRT79、KRT8、KRT7、KRT77、KRT5、KRT84、KRT18；毛色相關(guān)基因?yàn)镵IT，生長發(fā)育相關(guān)基因有AMHR2、SP1、OARDACVR1B、ADIPOQ。候選位點(diǎn)及注釋結(jié)果見表1。

黑色虛線表示顯著水平，黑色實(shí)線表示極顯著水平

表1 HapFLK候選位點(diǎn)及其基因

2.3 目的基因富集分析

將選擇信號突出的1、3、6號染色體的候選基因進(jìn)行富集分析，內(nèi)容包括細(xì)胞組分(cellular component)、分子功能(molecular function)、生物學(xué)過程(biological process)和KEGG信號通路分析。167個(gè)候選基因參與富集，對P值進(jìn)行FDR校正后，共得到7個(gè)顯著GO組分和6個(gè)顯著KEGG通路。候選基因GO分析顯著結(jié)果如表2所示，每個(gè)條目名稱及其數(shù)量見圖3；候選基因KEGG通路及其關(guān)系見圖4。

表2 候選基因顯著富集的GO條目

圖3 候選基因富集分析得到的GO條目及其基因數(shù)

3 討 論

3.1 基因組選擇性清掃檢測揭示兩個(gè)綿羊品種生產(chǎn)性狀的不同

通過基因組選擇性清掃檢測物種的起源進(jìn)化和基因交流[18]，探究家畜在環(huán)境變化過程中的遺傳機(jī)理。不同綿羊品種之間分子遺傳標(biāo)記不同，利用基因組選擇性清掃的方法獲得具有一定意義的分子標(biāo)記。本研究采用HapFLK方法，以藏綿羊作為外群，在剔除祖先成分之后，選擇性清掃結(jié)果集中體現(xiàn)了兩個(gè)品種間的基因組差異。兩個(gè)綿羊品種的主要差異體現(xiàn)在毛質(zhì)、毛色、尾型等性狀，這些結(jié)果在表型上也有體現(xiàn)。在表型上，和田羊與策勒黑羊在全年發(fā)情及多胎性上不一致，但在基因組清掃中并未體現(xiàn)其差異，認(rèn)為和田羊與策勒黑羊類似，均存在多胎潛力，可能因?yàn)轱曫B(yǎng)方式不同，而沒有表現(xiàn)出多胎性狀，這個(gè)結(jié)論與王瓊等[19]的研究結(jié)果相似。

3.2 角蛋白基因家族

基因組選擇信號分析獲得了KRT82、KRT80、KRT1、KRT79、KRT8、KRT7、KRT77、KRT5、KRT84、KRT18等角蛋白基因家族基因，富集分析結(jié)果顯示，KRT基因數(shù)量較多，且富集程度高。KRT基因主要為角蛋白基因和角蛋白相關(guān)基因[20]，角蛋白是羊毛纖維的主要成分，是形成皮膚毛囊細(xì)胞的主要結(jié)構(gòu)蛋白質(zhì)，是維持毛囊結(jié)構(gòu)并在毛囊中表達(dá)最豐富的蛋白質(zhì)，在很大程度上決定了羊毛的結(jié)構(gòu)特性，其編碼基因是毛囊基因表達(dá)和毛發(fā)生物學(xué)研究中重要的候選基因。該基因家族廣泛存在于綿羊的毛發(fā)和皮膚中，調(diào)節(jié)與控制綿羊毛囊及毛發(fā)的發(fā)育。同時(shí)，KRT基因參與動(dòng)物皮膚組織分化和綿羊毛長、彎曲度、卷曲方式[21]和纖維直徑[22-24]調(diào)控，是羊毛質(zhì)量[25-26]的主要調(diào)控基因。

經(jīng)過長時(shí)間的選育，兩個(gè)綿羊品種毛質(zhì)性狀存在很大差異，在生產(chǎn)性狀上策勒黑羊被毛為異質(zhì)毛，而和田羊?yàn)榈靥盒兔?。相關(guān)研究顯示，在綿羊機(jī)體中，KRT基因的突變會(huì)導(dǎo)致羊毛分叉，變細(xì)，變禿。在其他生物中，KRT基因也廣泛參與了其相關(guān)表型的表達(dá)與調(diào)控[20,27-28]，在人類的大部分皮膚病中都能找到KRT基因的影子，如人類KRT1突變可能會(huì)導(dǎo)致表皮溶解性魚鱗病(EI)的發(fā)生[29]，另外，人類的頭發(fā)相關(guān)表型也與KRT基因有關(guān)[30-32]。KRT18參與到小鼠肌肉的功能表達(dá)和小鼠表皮細(xì)胞的分化[33]。KRT7從頭表達(dá)導(dǎo)致上皮細(xì)胞可塑性和干細(xì)胞特性的改變是增加終末期腎腫瘤發(fā)生風(fēng)險(xiǎn)的關(guān)鍵因素[34]。

3.3 其他的基因及其功能

血小板衍生生長因子受體ɑ(platelet derived growth factor receptor alpha,PDGFRA)位于OAR6染色體，物理位置為69 735 803～69 771 932 bp，該基因在機(jī)體器官發(fā)育、傷口愈合中起著重要的作用，其突變位點(diǎn)與特發(fā)性嗜酸性粒細(xì)胞增多癥和家族性胃腸道間質(zhì)瘤以及其他多種癌癥有關(guān)，并參與皮膚及毛發(fā)色素沉著的表達(dá)[35-36]，是綿羊抗逆性和毛發(fā)顏色研究的候選基因。

KIT原癌基因，受體酪氨酸激酶(KIT proto-oncogene, receptor tyrosine kinase,KIT)位于OAR6染色體，物理位置為70 189 729～70 234 612 bp。該蛋白磷酸化多種細(xì)胞內(nèi)蛋白，這些蛋白與多種不同類型的細(xì)胞增殖、分化、遷移和凋亡相關(guān)，在造血、干細(xì)胞維持中起重要作用。KIT是調(diào)節(jié)黑色素生成、控制黑色素細(xì)胞增殖和凋亡的關(guān)鍵基因，是影響動(dòng)物毛皮色素沉著機(jī)制的重要基因[37]。和田羊被毛為白色，策勒黑羊被毛為黑色，該基因可能與兩個(gè)綿羊品種毛色相關(guān)。此外，KIT基因是生物環(huán)境適應(yīng)性相關(guān)基因，參與綿羊的環(huán)境適應(yīng)性及馴化[38-39]。

激酶插入域受體(kinase insert domain receptor,KDR)位于OAR6染色體，物理位置為70 607 356～70 562 533 bp。該基因?yàn)檠馨l(fā)育中至關(guān)重要的蛋白質(zhì)的編碼基因[40]，與類風(fēng)濕性關(guān)節(jié)炎(RA)的易感性[41]和原發(fā)性乳腺血管肉瘤的發(fā)生[42]緊密相關(guān)。

脂聯(lián)素，C1Q和膠原結(jié)構(gòu)域(adiponectin, C1Q and collagen domain containing,ADIPOQ)位于OAR1染色體，物理位置為198 631 126～19 861 993 bp，脂聯(lián)素蛋白基因完全在脂肪組織中表達(dá)[43-44]。編碼蛋白在血漿中循環(huán)，并參與代謝和荷爾蒙生成的過程。該基因極有可能與綿羊尾部脂肪的大小相關(guān)，綿羊尾部脂肪是綿羊能量的儲(chǔ)存庫，具有供能與保暖的作用，是綿羊適應(yīng)惡劣環(huán)境的一個(gè)重要特征。該基因可以作為培育適應(yīng)南疆環(huán)境綿羊的候選基因。

4 結(jié) 論

本研究通過基因組選擇性清掃的方法，解析和田羊和策勒黑羊的遺傳規(guī)律和基因組差異，發(fā)現(xiàn)兩個(gè)綿羊品種主要在毛質(zhì)性狀上存在差異，繁殖相關(guān)性狀差異較小，這表明和田羊經(jīng)過純繁能夠顯著提高繁殖性能。該結(jié)果可為南疆綿羊品種的遺傳改良提供理論基礎(chǔ)和分子選育方向。