于思?jí)?，郭占寶，?健，周正奎，侯水生

(中國(guó)農(nóng)業(yè)科學(xué)院北京畜牧獸醫(yī)研究所，農(nóng)業(yè)農(nóng)村部動(dòng)物遺傳育種與繁殖(家禽)重點(diǎn)實(shí)驗(yàn)室，北京 100193)

肉鴨產(chǎn)業(yè)是中國(guó)第三大肉類產(chǎn)業(yè)，是中國(guó)畜牧業(yè)的重要組成部分。中國(guó)鴨品種資源豐富，其中，北京鴨以其體型大、生長(zhǎng)速度快、飼料轉(zhuǎn)化率高、產(chǎn)蛋性能好等特點(diǎn)成為培育肉鴨新品種的理想材料[1]。英國(guó)櫻桃谷農(nóng)場(chǎng)以中國(guó)北京鴨品種資源為基礎(chǔ)，經(jīng)科學(xué)選育后，育成了瘦肉型白羽肉鴨品種——櫻桃谷鴨，該品種生長(zhǎng)速度快、飼料轉(zhuǎn)化率高，但其皮脂率較高，風(fēng)味物質(zhì)沉積能力差[2-3]，近年來(lái)，為了滿足中國(guó)居民肉鴨消費(fèi)需求，豐富中國(guó)肉鴨品種結(jié)構(gòu)，中國(guó)農(nóng)業(yè)科學(xué)院北京畜牧獸醫(yī)研究所自主培育了瘦肉型肉鴨品種中畜草原白羽肉鴨。中畜草原白羽肉鴨與櫻桃谷鴨均是以北京鴨資源群體為基礎(chǔ)選育，但又形成了各自特色。研究2個(gè)品種的遺傳差異，有望揭示人工選擇的遺傳變異機(jī)制。

鴨具有生長(zhǎng)快、世代間隔短，繁殖能力強(qiáng)等生物學(xué)特性[4-5]；且基因組小、重測(cè)序成本相對(duì)較低，是研究家養(yǎng)動(dòng)物人工選擇信號(hào)的理想模型。基因組測(cè)序技術(shù)的發(fā)展和應(yīng)用為遺傳育種研究提供了強(qiáng)有力的技術(shù)支撐[5-6]。目前，選擇信號(hào)分析已成為解析物種馴化和人工選擇印跡的重要方法，且在多種家養(yǎng)動(dòng)物中得到了有效利用，鑒定出了一批重要性狀候選基因。如Li等[7]通過(guò)對(duì)藏豬和中國(guó)地方豬種進(jìn)行全基因組重測(cè)序，結(jié)合群體遺傳分化指數(shù)(Fst)與群體核酸多樣性(θπ)選擇信號(hào)分析方法鑒定了一系列與缺氧、嗅覺(jué)、能量代謝和藥物反應(yīng)相關(guān)的基因，揭示了地方豬種在幾千年的馴化與人工選擇過(guò)程中基因組重要經(jīng)濟(jì)性狀相關(guān)基因的進(jìn)化方向。Zhang等[8]對(duì)19個(gè)世代的肉雞腹脂雙向選擇系進(jìn)行群體重測(cè)序，結(jié)合Fst分析確定了參與脂肪形成重要信號(hào)通路的基因，包括非經(jīng)典前折疊素RPB5相互作用因子1(URI1)、甘露糖結(jié)合凝集素2(MBL2)、乙酰膽堿轉(zhuǎn)移酶(CHAT)、前蛋白轉(zhuǎn)化酶枯草溶菌素1(PCSK1)處于強(qiáng)選擇狀態(tài)。在鴨中，F(xiàn)st結(jié)合全基因組關(guān)聯(lián)分析(GWAS)分析同樣鑒定出多個(gè)重要經(jīng)濟(jì)性狀相關(guān)基因，如胰島素樣生長(zhǎng)因子2信使核糖核酸結(jié)合蛋白1(IGF2BP1)可使北京鴨體格增大15%，并提高6%的飼料轉(zhuǎn)化率[5]；蘇氨酸天冬氨酸酶1(TASP1)有助于增加肌纖維直徑[9]；三磷酸腺苷結(jié)合盒轉(zhuǎn)運(yùn)蛋白G2(ABCG2)啟動(dòng)子區(qū)域的1個(gè)SNP變異影響鴨蛋殼顏色[10]。選擇信號(hào)分析與重要經(jīng)濟(jì)性狀候選基因的鑒定大大推動(dòng)了品種改良進(jìn)程，為畜牧業(yè)帶來(lái)了巨大效益。因此，利用選擇信號(hào)法對(duì)不同白羽肉鴨品種進(jìn)行分析，可以鑒定品種間遺傳分化的基因組區(qū)域，追溯品種在不同人工選擇下的基因組變異機(jī)制，以此闡明優(yōu)異性狀形成的遺傳基礎(chǔ)。

本研究旨在對(duì)2個(gè)白羽肉鴨品種的選擇信號(hào)進(jìn)行篩選，基于全基因組SNP變異，使用Fst和θπ2種分析方法發(fā)掘品種間遺傳分化信號(hào)，以期為瘦肉型白羽肉鴨品種的持續(xù)選育提供有效參考。

1 材料與方法

1.1 材料

于山東新希望六和集團(tuán)有限公司養(yǎng)殖基地和內(nèi)蒙古塞飛亞集團(tuán)有限公司養(yǎng)殖基地分別隨機(jī)挑選飼養(yǎng)至6周齡的櫻桃谷鴨商品代和中畜草原白羽肉鴨商品代各16只。翅靜脈采血2～5 mL，并使用肝素鈉抗凝，-20 ℃保存用于后續(xù)全基因組DNA提取與重測(cè)序。

1.2 方法

1.2.1 DNA提取與重測(cè)序 使用酚-仿法提取32個(gè)樣本的全基因組DNA?；贗llumina測(cè)序平臺(tái)進(jìn)行建庫(kù)和雙端測(cè)序。其中，櫻桃谷鴨(n=16)測(cè)序深度為5×，數(shù)據(jù)量5G/樣。中畜草原白羽肉鴨(n=16)測(cè)序深度為10×，數(shù)據(jù)量10G/樣。

1.2.2 SNP檢出 對(duì)raw reads去除接頭和低質(zhì)量reads后，得到clean reads。使用BWA[11]軟件將其與北京鴨參考基因組GCF_003850225.1(https:∥www.ncbi.nlm.nih.gov/assembly/GCF_003850225.1)進(jìn)行比對(duì)，使用GenomeAnalysisTK(GATK，v3.5)[12]確定需要重新比對(duì)的區(qū)域并進(jìn)行再比對(duì)，以增強(qiáng)后續(xù)SNP檢出的準(zhǔn)確性。使用GenomeAnalysisTK軟件進(jìn)行SNP變異信息的檢測(cè)與分型，個(gè)體基因組變異信息被合并后保存至VCF格式的文件中以供后續(xù)分析。使用VCFtools(v0.1.16)[13]對(duì)群體基因組變異信息進(jìn)行質(zhì)控并用于后續(xù)分析。質(zhì)控標(biāo)準(zhǔn)如下：①保留3×<平均測(cè)序深度<30×的SNP位點(diǎn)；②保留0.05<最小等位基因頻率<0.99的SNP位點(diǎn)；③保留檢出率>90%的SNP位點(diǎn)；④保留二等位的SNP位點(diǎn)。

1.2.3 主成分分析 使用PLINK(v1.90)[14]軟件對(duì)質(zhì)控后的SNP進(jìn)行主成分分析[15]，并用R語(yǔ)言plot函數(shù)對(duì)分析結(jié)果可視化。

1.2.4 選擇信號(hào)分析 本試驗(yàn)采用Fst來(lái)檢測(cè)櫻桃谷鴨與中畜草原白羽肉鴨群體間分化程度。為降低選擇信號(hào)的假陽(yáng)性，選擇以滑動(dòng)窗口的方式計(jì)算Fst以提高信號(hào)的靈敏性。利用VCFTools，以10 kb窗口、5 kb步長(zhǎng)滑動(dòng)計(jì)算2個(gè)群體間的Fst值。同時(shí)基于基因組雜合度計(jì)算兩個(gè)群體的θπ。每個(gè)群體的θπ值同樣以10 kb為窗口、5 kb為步長(zhǎng)在基因組上滑動(dòng)計(jì)算。對(duì)2個(gè)群體分別計(jì)算出的θπ值進(jìn)行標(biāo)準(zhǔn)化處理，每個(gè)θπ值取指數(shù)對(duì)數(shù)。然后以櫻桃谷鴨作為參考群，中畜草原白羽肉鴨作為試驗(yàn)群計(jì)算2個(gè)群體間標(biāo)準(zhǔn)化后的群體核酸多樣性比值(Pi)用于圖形展示。計(jì)算公式為：

Pi=ln(θπ)櫻桃谷鴨-ln(θπ)中畜草原白羽肉鴨

區(qū)域內(nèi)的Pi值越大說(shuō)明2個(gè)群體在該基因組區(qū)域的分化越明顯。將Fst計(jì)算值與Pi值的前1%重疊區(qū)域作為候選區(qū)域。

1.2.5 候選區(qū)域基因注釋及富集分析 使用NCBI數(shù)據(jù)庫(kù)下載的基因組注釋文件對(duì)候選區(qū)域進(jìn)行基因注釋。使用DAVID[16]與KOBAS[17]生物信息網(wǎng)站進(jìn)一步對(duì)候選基因進(jìn)行GO功能和KEGG通路富集分析，篩選候選基因。

2 結(jié) 果

2.1 櫻桃谷鴨與中畜草原白羽肉鴨的主成分分析

基于2個(gè)群體的全基因組SNP變異對(duì)2個(gè)白羽肉鴨品種進(jìn)行主成分分析。選取前2個(gè)主成分進(jìn)行繪圖，結(jié)果如圖1。2個(gè)白羽肉鴨品種在PC2上明顯分離，說(shuō)明2個(gè)白羽肉鴨品種雖然均源自北京鴨保種群，但經(jīng)過(guò)持續(xù)、定向化選育已經(jīng)分化成2個(gè)在基因組層面有顯著差異的群體。

2.2 櫻桃谷鴨與中畜草原白羽肉鴨的選擇信號(hào)檢測(cè)

以櫻桃谷鴨作為參考群，中畜草原白羽肉鴨作為試驗(yàn)群，根據(jù)10 kb滑動(dòng)窗口、5 kb步長(zhǎng)計(jì)算2個(gè)白羽肉鴨群體的Fst，所得Fst平均值為0.038，F(xiàn)st前1%閾值為0.177(圖2A)。高于閾值線的窗口共有1 372個(gè)，占據(jù)6.86 Mb的基因組區(qū)域。其中3號(hào)染色體100 265 001—100 270 000 bp區(qū)間內(nèi)Fst值最高，平均值為0.543。 在第5號(hào)染色體上受到選擇的顯著性區(qū)域最多，有344個(gè)窗口區(qū)間超過(guò)閾值線。

以相同窗口、步長(zhǎng)計(jì)算櫻桃谷鴨與中畜草原白羽肉鴨和核酸多樣性。中畜草原白羽肉鴨θπ平均值(0.00222)高于櫻桃谷鴨(0.00155)，說(shuō)明前者保留了更多的遺傳多樣性。所得群體間Pi的平均值為0.0377，Pi的前1%閾值為0.885(圖2B)。有1 487個(gè)窗口區(qū)間超過(guò)閾值線，共占據(jù)7.435 Mb大小的基因組區(qū)域。Pi值最高點(diǎn)位于5號(hào)染色體上，其中5號(hào)染色體52 520 001-52 525 000 bp區(qū)間內(nèi)選擇信號(hào)最強(qiáng)，平均值為4.97。 此外，在第1、2、9、11、13和16號(hào)染色體上也有較強(qiáng)的選擇信號(hào)。篩選Fst計(jì)算值與Pi比值前1%的重疊區(qū)域(圖2)，共獲得410 kb大小的候選區(qū)域。

圖1 2個(gè)白羽肉鴨品種的主成分分析Fig.1 Principal component analysis of two White-feather duck breeds

A，F(xiàn)st選擇信號(hào)分析；B，Pi比值選擇信號(hào)分析A，Selection signature analysis by Fst；B，Selection signature analysis by Pi圖2 2個(gè)白羽肉鴨品種常染色體上的選擇信號(hào)分布分析Fig.2 Distribution of selection signals identified on autosomes of two White-feather duck breeds

2.3 候選區(qū)域基因注釋與富集分析

使用參考基因組GCF_003850225.1對(duì)應(yīng)的gff注釋文件對(duì)候選區(qū)域進(jìn)行基因注釋，區(qū)域內(nèi)共注釋到21個(gè)基因。 對(duì)21個(gè)候選基因進(jìn)行GO功能富集分析，結(jié)果顯示共12條GO條目顯著富集(P<0.05)(圖3)，包含12個(gè)候選基因：PDE3A、PRKAR2B、GJD2、SEMA5A、LOC101803508(GOLGB1)、TMEM114、PTPRO、STXBP6、HS6ST1、SLC25A22、CNTNAP5和SHANK2。對(duì)21個(gè)候選基因進(jìn)行KEGG通路富集分析，結(jié)果表明，候選基因共富集到10條KEGG通路，基因功能涉及黏多糖生物合成、丙酸代謝、乙醛酸和二羧酸代謝、纈氨酸、亮氨酸和異亮氨酸的降解、細(xì)胞間隙連接、碳代謝、胰島素信號(hào)通路、嘌呤代謝、細(xì)胞周期等方面，其中有4條通路顯著富集(P<0.05)，涉及到2個(gè)候選基因：HS6ST1與PCCB(表1)。

基于Fst>0.250對(duì)GO與KEGG富集到的13個(gè)候選基因進(jìn)一步篩選，發(fā)現(xiàn)6個(gè)強(qiáng)烈受選擇的基因：PDE3A、PRKAR2B、SEMA5A、SHANK2、STXBP6和LOC101803508(GOLGB1)(Fst>0.250，Pi>0.936)(表2)，主要富集在脂質(zhì)代謝、氨基酸代謝、繁殖能力等通路。

圖3 21個(gè)候選基因的GO功能富集分析Fig.3 GO function enrichment analysis of 21 candidate genes

表1 21個(gè)候選基因的KEGG通路富集分析

表2 強(qiáng)選擇信號(hào)下的基因及其功能

3 討 論

為了滿足國(guó)內(nèi)對(duì)肉鴨多元化的消費(fèi)需求，提升中國(guó)肉鴨行業(yè)的國(guó)際競(jìng)爭(zhēng)力,本團(tuán)隊(duì)利用中國(guó)北京鴨品種資源培育了中畜草原白羽肉鴨新品種。該品種不僅更符合國(guó)內(nèi)市場(chǎng)需求，且其關(guān)鍵生產(chǎn)性能指標(biāo)如料重比、胸肉率、皮脂率等較引進(jìn)品種而言更具優(yōu)勢(shì)[18]。本研究中，基于重測(cè)序數(shù)據(jù)對(duì)櫻桃谷鴨與中畜草原白羽肉鴨基因組進(jìn)行遺傳分化分析，PCA結(jié)果證明2個(gè)品種已經(jīng)發(fā)生了明顯分化，F(xiàn)st與Pi平均值同樣證明2個(gè)群體有中等程度的遺傳分化,且中畜草原白羽肉鴨相較于櫻桃谷鴨具有更高的遺傳多態(tài)性。品種遺傳多態(tài)性越高，群體的遺傳變異程度就越高，群體適應(yīng)潛在的環(huán)境變化的能力越強(qiáng)[19]，同時(shí)為滿足多元化市場(chǎng)需求提供了變異基礎(chǔ)。

結(jié)合Fst與Pi綜合分析2個(gè)品種間選擇信號(hào)并對(duì)候選區(qū)域進(jìn)行基因注釋，共注釋到了21個(gè)基因,這些基因參與調(diào)節(jié)代謝、繁殖、免疫、神經(jīng)發(fā)育等。其中，KEGG顯著富集到4條通路，涉及2個(gè)候選基因HS6ST1與PCCB。研究表明，HS6ST1基因參與硫酸乙酰肝素的代謝，同時(shí)也參與牛肌肉脂質(zhì)代謝的調(diào)控[20]。PCCB基因是生物素依賴性丙酰輔酶A羧化酶(PCC)的2個(gè)亞基之一[21]，參與脂肪酸生物合成和短鏈脂肪酸分解代謝過(guò)程，為影響肉雞肌內(nèi)脂肪沉積的差異表達(dá)基因[22]。肉的多汁性、風(fēng)味和嫩度等品質(zhì)特征主要取決于肌肉總脂肪含量，推測(cè)中畜草原白羽肉鴨與櫻桃谷鴨在肌內(nèi)脂肪沉積能力上發(fā)生了分化，因而HS6ST1與PCCB基因受到了選擇。遺傳分化指數(shù)的提出者Wright認(rèn)為，當(dāng)Fst>0.250時(shí)，群體間有很大的遺傳分化[23]。因此，基于Fst>0.250對(duì)GO與KEGG富集到的共計(jì)13個(gè)候選基因進(jìn)一步篩選發(fā)現(xiàn)，有6個(gè)基因PDE3A、PRKAR2B、SEMA5A、SHANK2、STXBP6與LOC101803508(GOLGB1)受到了較強(qiáng)的選擇(Fst>0.250，Pi>0.936)。其中，PDE3A是在心肌細(xì)胞和血小板環(huán)磷酸腺苷(cAMP)介導(dǎo)的細(xì)胞內(nèi)信號(hào)調(diào)節(jié)中起重要作用的酶[24]。結(jié)合前期北京鴨不同組織器官時(shí)空轉(zhuǎn)錄組基因表達(dá)數(shù)據(jù)[5]，發(fā)現(xiàn)PDE3A在8周齡北京鴨的心臟與腹脂中高表達(dá)。已有研究表明，PDE3A參與代謝，可以促進(jìn)小鼠胃腸道蠕動(dòng)[25],還作用于動(dòng)物的卵泡，影響卵泡成熟[26]，從而影響動(dòng)物繁殖性能。PRKAR2B被證實(shí)為精子發(fā)生的熱敏基因[27]，且被鑒定為蛋重的候選基因[28]。此外，PRKAR2B基因影響脂質(zhì)儲(chǔ)存和代謝[29]，是小鼠和人類脂肪細(xì)胞分化所必需的[30],該基因在8周齡北京鴨腹脂中高表達(dá)[5]，推測(cè)PRKAR2B基因與中畜草原白羽肉鴨的繁殖與優(yōu)良肉質(zhì)性狀有關(guān)而受到選擇。研究證實(shí)，SEMA5A基因參與肉雞的免疫調(diào)控[31]，并在北京鴨腹脂中高表達(dá)[5]。SHANK2基因調(diào)節(jié)晝夜節(jié)律通路。大鼠節(jié)律基因的差異表達(dá)與脂質(zhì)代謝基因的差異表達(dá)呈正相關(guān)，晝夜節(jié)律影響甘油酯類和甘油磷脂類脂質(zhì)物質(zhì)的代謝[32]。甘油酯與甘油磷脂是肉制品風(fēng)味物質(zhì)或風(fēng)味前體物質(zhì)的重要來(lái)源[33]。北京鴨全譜轉(zhuǎn)錄組數(shù)據(jù)顯示，SHANK2基因在北京鴨腦中有較高的表達(dá)水平，且該基因在中畜草原白羽肉鴨群體中受到強(qiáng)選擇，說(shuō)明經(jīng)過(guò)長(zhǎng)期人工選育，中畜草原鴨可能在晝夜節(jié)律上同櫻桃谷鴨發(fā)生了較大分化，影響了甘油酯類和甘油磷脂類物質(zhì)的代謝，使其具有優(yōu)質(zhì)風(fēng)味。SHANK2基因影響肉雞神經(jīng)發(fā)育、細(xì)胞增殖和分化以及采食量[34]。 STXBP6是從腦中分離出來(lái)的突觸結(jié)合蛋白，在北京鴨腦、脾臟中高表達(dá)[5]。小鼠STXBP6表達(dá)量上調(diào)會(huì)間接抑制胰島素分泌而增加糖尿病的風(fēng)險(xiǎn)[35]，敲除STXBP6可能會(huì)引起小鼠代謝異常從而降低小鼠體重[36]。LOC101803508基因在人類中又名GOLGB1基因，參與蛋白質(zhì)代謝與囊泡運(yùn)輸，在北京鴨腹脂中高表達(dá)[5]。

脂肪酸與氨基酸的種類及含量是影響肉類風(fēng)味品質(zhì)的重要因素[37]，相對(duì)于櫻桃谷鴨，PDE3A、PRKAR2B、SEMA5A、SHANK2、STXBP6和LOC101803508(GOLGB1)基因在中畜草原白羽肉鴨基因組上受到較強(qiáng)的選擇，這可能影響中畜草原白羽肉鴨繁殖力與風(fēng)味肉品質(zhì)。但中畜草原白羽肉鴨與櫻桃谷鴨遺傳分化候選基因的功能及其精細(xì)定位仍需進(jìn)一步研究。隨著人民生活水平的提高，對(duì)肉類的要求逐漸從產(chǎn)量向品質(zhì)轉(zhuǎn)變，尤其是對(duì)禽肉的要求更為突出。定位影響肉品質(zhì)的主效基因，探索生長(zhǎng)速度、產(chǎn)量與肉品質(zhì)之間的平衡關(guān)系并將其用于家禽生產(chǎn)將是未來(lái)家禽育種的新方向。本研究篩選到的脂質(zhì)代謝、氨基酸代謝、免疫調(diào)控相關(guān)的基因，如PDE3A、PRKAR2B、SEMA5A、SHANK2、STXBP6與LOC101803508(GOLGB1),為后續(xù)篩選白羽肉鴨品種特異性分子標(biāo)記提供了參考。

4 結(jié) 論

本研究利用二代重測(cè)序數(shù)據(jù)，聯(lián)合Fst與Pi分析方法對(duì)2個(gè)源于北京鴨資源群體的瘦肉型白羽肉鴨品種進(jìn)行遺傳差異分析，發(fā)現(xiàn)中畜草原白羽肉鴨保留了更多的遺傳多樣性。鑒定到白羽肉鴨在選育過(guò)程中的候選基因21個(gè)，其中與繁殖力和風(fēng)味肉品質(zhì)相關(guān)的6個(gè)基因PDE3A、PRKAR2B、SEMA5A、SHANK2、STXBP6與LOC101803508(GOLGB1)受到了強(qiáng)選擇。