付言峰，趙為民，李 輝，李碧俠，程金花，徐小波，任守文

(江蘇省農(nóng)業(yè)科學(xué)院 畜牧研究所；江蘇現(xiàn)代農(nóng)業(yè)(生豬)產(chǎn)業(yè)技術(shù)體系阜寧推廣示范基地， 江蘇省農(nóng)業(yè)種質(zhì)資源保護(hù)與利用平臺(tái)，農(nóng)業(yè)部種養(yǎng)結(jié)合重點(diǎn)實(shí)驗(yàn)室，江蘇 南京 210014)

產(chǎn)仔數(shù)是養(yǎng)豬生產(chǎn)中一個(gè)重要的育種性狀、繁殖性狀和經(jīng)濟(jì)性狀[1-2]，而總產(chǎn)仔數(shù)和產(chǎn)活仔數(shù)一直是種豬選育中最常記錄的2個(gè)產(chǎn)仔數(shù)指標(biāo)。由于產(chǎn)仔數(shù)遺傳力較低且又是伴性性狀[3]，利用常規(guī)的育種方法選育難度較大，且進(jìn)展緩慢[4]。隨著分子生物學(xué)的發(fā)展，全基因組SNPs分子標(biāo)記逐漸應(yīng)用于豬產(chǎn)仔數(shù)的選擇育種中。

利用全基因組SNPs分子標(biāo)記，可以進(jìn)行全基因組關(guān)聯(lián)分析(Genome-wide association study，GWAS)[5]和全基因組選擇育種(Genomic selection)[6]。GWAS是利用全基因組范圍內(nèi)篩選出的高密度分子標(biāo)記(現(xiàn)大多為SNP)進(jìn)行試驗(yàn)群體的全基因組掃描，再通過統(tǒng)計(jì)學(xué)方法將掃描所得的分子標(biāo)記基因型與表型性狀關(guān)聯(lián)分析，最后鑒定出影響表型性狀的相關(guān)分子標(biāo)記和候選基因[7]，近年來已在奶牛和生豬育種生產(chǎn)上得到廣泛的應(yīng)用。其中在豬上，研究人員利用全基因組關(guān)聯(lián)分析技術(shù)在嵊縣花豬、金華豬群體分別篩選了160，124 個(gè)達(dá)到顯著性水平的SNP[8]，還利用全基因組SNPs芯片在杜洛克豬中篩選了乳腺癌抗雌激素耐藥基因3(BCAR3)、N-乙基馬來酰亞胺敏感因子(NSF)等初情期日齡潛在的候選基因[9]。

本研究將利用Illumina全基因組豬50K SNPs芯片檢測(cè)了186頭加系大白豬的基因組，并與這些豬的總產(chǎn)仔數(shù)、產(chǎn)活仔數(shù)性狀進(jìn)行了關(guān)聯(lián)分析，進(jìn)一步對(duì)獲得的顯著關(guān)聯(lián)SNPs染色體物理位置100 kb范圍的基因進(jìn)行篩選挖掘，最后對(duì)這些候選基因進(jìn)行GO功能和KEGG通路注釋，擬解決試驗(yàn)豬高產(chǎn)仔數(shù)關(guān)鍵位點(diǎn)和基因的篩選，研究結(jié)果將對(duì)進(jìn)一步揭示豬高產(chǎn)仔性狀的遺傳基礎(chǔ)，提高豬育種效率和經(jīng)濟(jì)效益具有較大的借鑒意義。

1 材料和方法

1.1 試驗(yàn)動(dòng)物和樣品采集

試驗(yàn)動(dòng)物來自江蘇天兆實(shí)業(yè)有限公司。首先，在1 800頭加系大白豬基礎(chǔ)群母豬中選取了300頭組成核心群，其次，又在核心群中挑選了186頭能繁母豬為研究對(duì)象，按豬場(chǎng)常規(guī)進(jìn)行飼喂和管理，同時(shí)在自由采食及飲水條件下，注意進(jìn)行疾病的防治和行為體況的觀察。分2次進(jìn)行耳組織樣的采集，然后置于干冰中冷凍保存，回到實(shí)驗(yàn)室后置于-20 ℃冰箱保存，用于提取基因組DNA。

1.2 DNA提取和全基因組SNPs芯片檢測(cè)

DNA提取方法：將約2 mg 豬組織樣剪碎，加入600 μL 組織DNA 提取液(組織裂解液)，上下顛倒10次混勻；加入蛋白酶K(20 mg/mL)6 μL，使其至終濃度200 μg/mL，50～55 ℃水浴鍋內(nèi)消化過夜；加入等體積(606 μL)苯酚抽提，12 000 r/min 離心10 min，取上清(≤700 μL)；加入等體積苯酚∶氯仿(1∶1)抽提，12 000 r/min 離心10 min，取上清(≤700 μL)；加入等體積氯仿抽提，12 000 r/min 離心10 min，取上清(≤500 μL)；加入2倍體積的無水乙醇沉淀DNA，12 000 r/min 離心10 min，棄上清。用至少500 μL 70%乙醇快速洗滌沉淀，重復(fù)2次(為防止DNA溶解，此步操作應(yīng)快速進(jìn)行)；常溫晾干，用30，100～200 μL，1 000～2 000 μL滅菌雙蒸水(或TE溶液)溶解DNA；0.7%的瓊脂糖電泳初步檢測(cè)濃度和純度，紫外分光光度儀測(cè)定DNA 含量。

DNA經(jīng)過質(zhì)控檢測(cè)合格之后，進(jìn)行Illumina全基因組豬50K SNPs芯片檢測(cè)，質(zhì)控后共得到1.75 Gb的有效數(shù)據(jù)，用于后期的變異檢測(cè)。

1.3 全基因組關(guān)聯(lián)分析(Genome-wide association studies，GWAS)

本研究使用GEMMA軟件[10]中的混合線性模型，進(jìn)行豬產(chǎn)仔表型性狀與全基因組SNPs之間的GWAS分析，GWAS分析模型為：Y=μ+Xb+u+e，其中，Y表示加系大白豬產(chǎn)仔性狀表型值；μ表示產(chǎn)羔數(shù)平均值；X表示固定效應(yīng)矩陣，b為固定效應(yīng)向量；u表示剩余多基因效應(yīng)；e表示產(chǎn)仔表型值的隨機(jī)殘差，u和e均服從正態(tài)分布。

本研究對(duì)湖羊GWAS 結(jié)果中的差異顯著性水平(P值)采用Bonferroni 校正，并使用R軟件中的”qqman”[11]進(jìn)行曼哈頓圖和Q-Q圖(Quantile-Quantile Plot)繪制。

1.4 候選基因篩選

利用Ensembl網(wǎng)站(http：//www.ensembl.org/index.html)的大白豬基因組序列Sscrofa 11.1(GCA_001700135.1)進(jìn)行每個(gè)SNP前后40 kb范圍內(nèi)定位的候選基因篩選，再利用其基因注釋工具和NCBI網(wǎng)站(http：//www.ncbi.nlm.nih.gov/)Pubmed文獻(xiàn)數(shù)據(jù)庫(kù)進(jìn)行基因的功能注釋。

1.5 基因型質(zhì)控

凡是符合以下3個(gè)條件中的任何一個(gè)，基因型質(zhì)控過程將刪除該個(gè)體或位點(diǎn)：分型缺失數(shù)>10%的位點(diǎn)、最小等位基因頻率低于0.01的位點(diǎn)、哈代-溫伯格檢驗(yàn)P值<0.000 001 的位點(diǎn)。

2 結(jié)果與分析

2.1 表型數(shù)據(jù)統(tǒng)計(jì)

選用的186頭母豬均有詳細(xì)的產(chǎn)仔數(shù)記錄，這些記錄包括：耳號(hào)、樣品號(hào)、豬只類型、父親耳號(hào)、出生日期、出生體質(zhì)量、乳頭數(shù)、達(dá)100 kg日齡、備注、配種日期、產(chǎn)仔日期、總產(chǎn)仔數(shù)、產(chǎn)活仔數(shù)、健仔數(shù)、斷奶仔豬數(shù)、出生窩質(zhì)量、出生個(gè)體均質(zhì)量。這些母豬均含有第1，2胎產(chǎn)仔數(shù)記錄，最高胎次為第8胎(表1)。

表1 試驗(yàn)豬產(chǎn)仔數(shù)記錄Tab.1 Litter size records of experimental pigs

2.2 全基因組關(guān)聯(lián)分析

在經(jīng)SNP質(zhì)控之后，在18對(duì)(36條)常染色體上一共得到36 867個(gè)SNP標(biāo)記，用于試驗(yàn)豬產(chǎn)仔數(shù)(總產(chǎn)仔數(shù)、產(chǎn)活仔數(shù))的全基因組關(guān)聯(lián)分析研究。這些加系大白豬所有胎次總產(chǎn)仔數(shù)共顯著關(guān)聯(lián)到2個(gè)SNP，分別為seq-rs323899658和seq-rs329781338，其中seq-rs323899658 SNP沒有注釋到基因，而seq-rs329781338 SNP又注釋到3個(gè)基因，分別為SSBP1(線粒體 DNA 單鏈結(jié)合蛋白基因)、WEE2和KIAA1147(表2)。

表2 豬總產(chǎn)仔數(shù)顯著關(guān)聯(lián)SNP的候選基因注釋結(jié)果Tab.2 Genes associated with total number born(TNB)of pigs

這些試驗(yàn)豬產(chǎn)活仔數(shù)共顯著關(guān)聯(lián)到7個(gè)SNP，分別為：seq-rs81238474、seq-rs321377412、seq-rs340736313、seq-rs80782154、seq-rs81244816、seq-rs81467772和seq-rs329781338。這7個(gè)SNP中有5個(gè)SNP有基因注釋，2個(gè)沒有基因注釋；這5個(gè)有注釋的SNP共注釋到22個(gè)基因，分別為TAS2R60、EPHA1、FAM131B、CLCN1、CASP2、TMEM139、TRBV21OR9-2、PRSS2、TRBV19、TRBV24-1、U6、SSBP1、WEE2、KIAA1147、NKX2-8、ALKBH3、HSD17B12、U6、HOMER2、WHAMM、FSD2和SCARNA15(表3)。

表3 豬產(chǎn)活仔數(shù)顯著關(guān)聯(lián)SNP的候選基因注釋結(jié)果Tab.3 Genes associated with number born alive(NBA)of pigs

加系大白豬所有胎次總產(chǎn)仔數(shù)和產(chǎn)活仔數(shù)的GWAS結(jié)果如圖1，2所示，從曼哈頓圖結(jié)果(圖1-A、圖2-A)來看，總產(chǎn)仔數(shù)有2個(gè)SNP達(dá)到Bonferroni校正的閾值線，產(chǎn)活仔數(shù)有7個(gè)SNP達(dá)到Bonferroni校正的閾值線，表明加系大白豬總仔數(shù)和產(chǎn)活仔數(shù)分別有2，7個(gè)顯著的SNP位點(diǎn)。從總產(chǎn)仔數(shù)和產(chǎn)活仔數(shù)的QQ圖結(jié)果來看(圖1-B、圖2-B)，總產(chǎn)仔數(shù)有2個(gè)P值的觀測(cè)值明顯超過期望值，剩下P值的觀測(cè)值都低于期望值或與其相當(dāng)。產(chǎn)活仔數(shù)有較多的P值觀測(cè)值都高于期望值，剩下的2個(gè)值相當(dāng)，基本沒有低于期望值的觀測(cè)值。由此可見，QQ圖與曼哈頓圖結(jié)果一致。

圖1 豬總產(chǎn)仔數(shù)曼哈頓圖(A)和QQ-Plot圖(B)Fig.1 Manhattan chart(A)and QQ plot chart(B)of total number born(TNB)in pigs

圖2 豬產(chǎn)活仔數(shù)曼哈頓圖(A)和QQ-Plot圖(B)Fig.2 Manhattan chart(A)and QQ plot chart(B)of number born alive(TNA)in pigs

2.3 候選基因功能富集分析

在ENSEMBL網(wǎng)站上下載家豬的參考基因組信息，將顯著SNP位點(diǎn)注釋到距離其100 kb以內(nèi)的基因上，總產(chǎn)仔數(shù)的2個(gè)顯著SNP注釋到3個(gè)基因，產(chǎn)活仔數(shù)的7個(gè)SNP注釋到21個(gè)基因。然后使用clusterProfiler包將注釋到的候選基因根據(jù)GO(Gene Ontology)和KEGG數(shù)據(jù)庫(kù)進(jìn)行基因功能富集分析。

豬總產(chǎn)仔數(shù)候選基因GO富集到的生物學(xué)過程一共42個(gè)(P<0.05)按顯著性由高到低(P值由低到高)分別為：多生物過程的調(diào)控、減數(shù)分裂細(xì)胞周期的負(fù)調(diào)控、減數(shù)分裂核分裂的調(diào)控、參與卵母細(xì)胞成熟的減數(shù)分裂、細(xì)胞周期過程的調(diào)控、參與卵母細(xì)胞成熟的減數(shù)分裂細(xì)胞周期過程的調(diào)控、卵母細(xì)胞成熟、卵母細(xì)胞發(fā)育、卵母細(xì)胞分化、生殖過程的負(fù)調(diào)控和減數(shù)分裂細(xì)胞周期的調(diào)控。這些生物學(xué)過程均來自WEE2基因(圖3-A)。

豬總產(chǎn)仔數(shù)候選基因GO富集到的分子功能一共5個(gè)(P<0.01)，按顯著性由高到低(P值由低到高)分別為：非跨膜蛋白酪氨酸激酶活性、跨膜受體蛋白酪氨酸激酶活性、蛋白酪氨酸激酶活性、跨膜受體蛋白激酶活性和鎂離子結(jié)合，這些分子功能均來自WEE2基因(圖3-B)。豬總產(chǎn)仔數(shù)候選基因生物通路一共5條(P<0.05)，按顯著性由高到低(P值由低到高)分別為：錯(cuò)配修復(fù)、DNA復(fù)制、同源重組、細(xì)胞周期和人類免疫缺陷病毒1型感染。這些生物通路均來自WEE2和SSBP1基因(圖3-C)。

豬產(chǎn)活仔數(shù)候選基因GO富集到的分子功能一共3個(gè)(P<0.01)，按顯著性由高到低(P值由低到高)分別為：跨膜受體蛋白酪氨酸激酶活性、蛋白酪氨酸激酶活性和跨膜受體蛋白激酶活性。這3種富集到的分子功能均同時(shí)來自EphA1和WEE2基因(圖4)。

3 結(jié)論與討論

3.1 群體分層分析

GWAS研究是基于全基因組范圍內(nèi)不同遺傳變異位點(diǎn)(多為SNPs)之間會(huì)出現(xiàn)連鎖不平衡的現(xiàn)象來確定影響某些表型性狀或數(shù)量性狀的基因[12]，而影響連鎖不平衡的因素很多，如遺傳漂變、群體分層等，群體分層會(huì)導(dǎo)致GWAS結(jié)果假陽性[12-13]。根據(jù)本研究得到的QQ圖結(jié)果可以看出，P值觀察值和預(yù)測(cè)值基本相同，說明試驗(yàn)所用的加系大白豬樣本未出現(xiàn)群體分層現(xiàn)象，本研究所得到的GWAS結(jié)果是可靠的，所使用的混合線性模型是合理的。

A.生物學(xué)過程；B.分子功能；C.KEGG通路。 A.Biological process；B.Molecular function；C.KEGG pathway.

圖4 豬產(chǎn)活仔數(shù)候選基因GO功能富集分析-分子功能Fig.4 Go function enrichment analysis of candidate genes associated with total number born(TNB) in pigs，molecular function

3.2 全基因組關(guān)聯(lián)分析

利用SNP 生物芯片技術(shù)和全基因組關(guān)聯(lián)分析(GWAS)技術(shù)能夠使用較低的成本和較高的準(zhǔn)確率篩選出與繁殖性狀關(guān)聯(lián)性強(qiáng)的遺傳標(biāo)記，從而進(jìn)一步篩選出相應(yīng)的候選基因，為后期的功能驗(yàn)證以及選種選育、培養(yǎng)更高繁殖力的商品豬種打下基礎(chǔ)[8]。GWAS曼哈頓圖顯示，本研究所篩選到的所有SNP都均勻分布，QQ圖也顯示，本研究選擇的混合線性模型是適合的，前面的研究結(jié)果也佐證了這一結(jié)論[12]。

3.2.1 總產(chǎn)仔數(shù) 本研究所用的加系大白豬的總產(chǎn)仔數(shù)性狀鑒定到2個(gè)顯著性的SNP位點(diǎn)，其中1個(gè)SNP又注釋到3個(gè)基因SSBP1、WEE2和KIAA1147，而這3個(gè)基因中只有WEE2基因GO功能富集到生物學(xué)過程和分子功能，只有WEE2和SSBP1基因KEGG通路能富集到，且前3位分別為錯(cuò)配修復(fù)、DNA復(fù)制和同源重組。有研究表明，WEE2基因(WEE1同源體，也稱為WEE1B)屬于母源基因，僅在卵巢的卵母細(xì)胞中高表達(dá)。在小鼠GV 期卵子中降調(diào)WEE2基因表達(dá)后，小鼠卵母細(xì)胞能夠完成減數(shù)第一次分裂后排出極體，由MI 卵子發(fā)育至PB1 卵子，但通過IVF受精后，卵子不能形成原核，表現(xiàn)為受精障礙[14]。由此可見，WEE2基因在卵母細(xì)胞受精過程中發(fā)揮至關(guān)重要的作用。WEE2基因純合突變c.G585C(p.K195N)純合錯(cuò)義突變將導(dǎo)致受精失敗和女性不孕[15]。本研究結(jié)果表明，WEE2基因與加系大白豬總產(chǎn)仔數(shù)相關(guān)。

另一個(gè)基因SSBP1(又稱mtSSB)，中文名稱為線粒體 DNA 單鏈結(jié)合蛋白，具有促進(jìn)線粒體 DNA復(fù)制和損傷修復(fù)進(jìn)而調(diào)控線粒體生成的重要作用[16]，近年來，SSBP1在腫瘤中作用受到越來越多的關(guān)注，有研究表明，SSBP1在卵巢癌組織中的表達(dá)顯著低于正常卵巢組織[17]。此基因的生物學(xué)過程、分子功能和信號(hào)通路與總產(chǎn)仔數(shù)的緊密相關(guān)。

3.2.2 產(chǎn)活仔數(shù) 本研究所用的加系大白豬的產(chǎn)活仔數(shù)性狀鑒定到7個(gè)顯著關(guān)聯(lián)SNP，其中5個(gè)SNP有基因注釋，共注釋到22個(gè)基因，這些基因GO功能富集到3種分子功能，跨膜受體蛋白酪氨酸激酶活性、蛋白酪氨酸激酶活性和跨膜受體蛋白激酶活性，對(duì)應(yīng)的基因?yàn)镋phA1和WEE2。GO功能中的生物學(xué)過程、細(xì)胞學(xué)分區(qū)和KEGG通路均未富集到。

目前已有許多研究表明，EphA1(促紅細(xì)胞生成素產(chǎn)生肝細(xì)胞受體A1)和ephrin A1(促紅細(xì)胞生成素產(chǎn)生肝細(xì)胞配體A1)在人[18]、小鼠和豬[19]的胚胎附植過程中發(fā)揮著重要的作用。EphA1基因第15 外顯子上有一個(gè)G→C的突變，且該SNP突變的GG基因型可作為大白豬和長(zhǎng)白豬中的優(yōu)勢(shì)基因型[20]。豬胚胎附植前、中、后期ephrin A1和Eph A1在其子宮內(nèi)膜附植點(diǎn)mRNA表達(dá)量均呈先升高后降低的趨勢(shì)，且妊娠豬表達(dá)量極顯著高于空懷豬(P<0.01)[21]，EphA1的CpG島甲基化很可能會(huì)影響此基因的表達(dá)[22]。梅山豬胚胎附植前、中、后期ephrin A1和Eph A1在其子宮內(nèi)膜附植點(diǎn)mRNA和蛋白表達(dá)量均呈先升高后降低的趨勢(shì)[23]。另外，胚胎附植期在人的子宮內(nèi)膜上皮細(xì)胞和子宮內(nèi)膜癌細(xì)胞中均能檢測(cè)到Eph A1的mRNA和蛋白表達(dá)[18]。

綜上所述，加系大白豬所有胎次母豬總產(chǎn)仔數(shù)和產(chǎn)活仔數(shù)的全基因組范圍內(nèi)的顯著關(guān)聯(lián)SNP分別為seq-rs329781338和seq-rs81244816，這2個(gè)SNP注釋到的功能候選基因分別為WEE2和EphA1。本研究結(jié)果可為下一步提高豬產(chǎn)仔數(shù)和產(chǎn)仔數(shù)的基因編輯精準(zhǔn)育種提供一定參考。