全基因組關(guān)聯(lián)分析鑒定法系大白豬繁殖性狀基因組區(qū)域和候選基因

李秀領(lǐng)，陳健梅，吳姿儀，韓雪蕾，喬瑞敏，王克君，楊 峰，王獻(xiàn)偉，李新建*

（1.河南農(nóng)業(yè)大學(xué)動(dòng)物科技學(xué)院，河南鄭州 450046；2.河南省畜牧總站，河南鄭州 450008）

母豬的繁殖能力在養(yǎng)豬業(yè)的經(jīng)濟(jì)效益構(gòu)成中起著重要作用。為了獲取更高的收益，提高母豬的繁殖性能尤為重要?？偖a(chǎn)仔數(shù)（Total Number Born，TNB）、產(chǎn)活仔數(shù)（Number Born Alive，NBA）和出生窩重（Litter Weight Born Alive，LWB）是評(píng)估母豬繁殖性能的關(guān)鍵指標(biāo)。然而，繁殖性狀的遺傳結(jié)構(gòu)較為復(fù)雜，由于其低遺傳力特性（約0.1），傳統(tǒng)的育種計(jì)劃獲得的遺傳改良非常緩慢。因此，標(biāo)記輔助選擇（Marker-Assisted Selection，MAS）和基因組選擇（Genomic Selection，GS）被廣泛應(yīng)用于繁殖性狀的遺傳改良。

近幾十年來，隨著分子育種技術(shù)的發(fā)展，全基因組關(guān)聯(lián)分析（Genome-Wide Association Study，GWAS）已被廣泛用于識(shí)別復(fù)雜性狀中的數(shù)量性狀基因座（Quan titative Trait Loci，QTL）和數(shù)量性狀核苷酸（Quantita tive Trait Nucleotides，QTN）。在豬相關(guān)的研究中，GWAS 已經(jīng)成為定位與復(fù)雜性狀相關(guān)的SNPs 和候選基因的有力方法之一。在GWAS 的統(tǒng)計(jì)方法中，單步全基因組關(guān)聯(lián)分析（Single-Step GWAS，ssGWAS）已被用于研究豬的經(jīng)濟(jì)性狀。與傳統(tǒng)GWAS 相比，ssGWAS 可以將個(gè)體的基因型、表型和系譜信息都整合進(jìn)混合線性模型中。因此，ssGWAS 可以提高低遺傳力性狀候選基因定位的精度。

此外，繁殖性狀受許多微效基因控制，但到目前為止，只有少數(shù)候選基因被確定，而且只是解釋了一小部分遺傳變異。因此，繁殖性狀的具體遺傳機(jī)制仍然需要進(jìn)一步研究。本研究的目的是通過在法系大白豬群體中進(jìn)行ssGWAS 檢測影響繁殖性狀（TNB、NBA 和LWB）的遺傳變異并確定候選基因。

1 材料與方法

1.1 動(dòng)物和表型數(shù)據(jù) 在本研究中，收集了河南某核心養(yǎng)豬場392 頭法系大白豬的系譜信息以及2020 年1—10 月的表型信息。其中，表型包括TNB、NBA 和LWB，所有表型記錄均為母豬第1 胎次產(chǎn)仔記錄。所有個(gè)體均被分組在標(biāo)準(zhǔn)的商業(yè)豬舍中并且每個(gè)個(gè)體都是在相同的飼料和飼養(yǎng)管理?xiàng)l件下飼養(yǎng)的。TNB、NBA和LWB 的表型記錄的平均值、標(biāo)準(zhǔn)差、最小值和最大值見表1。

表1 繁殖性狀的描述性統(tǒng)計(jì)

1.2 基因型數(shù)據(jù)和質(zhì)控 使用組織基因組DNA 提取試劑盒從豬耳組織樣中提取DNA，并使用紐勤公司GeneSeek Genomic Profiler（GGP）豬50K SNP 芯 片進(jìn)行基因分型。數(shù)據(jù)質(zhì)量控制由PLINK v1.90進(jìn)行。剔除基因型檢出率<0.90、次等位基因頻率<0.05 和Hardy-Weinberg 平衡檢驗(yàn)<10的SNP 標(biāo)記，刪 除SNP 標(biāo)記檢出率<0.90 的個(gè)體。所有位于未映射區(qū)域的SNPs 也被排除。最終保留了43 098 個(gè)SNPs 和392 個(gè)個(gè)體用于后續(xù)分析。

1.3 遺傳參數(shù)估計(jì) 使用GCTA v1.93.2進(jìn)行限制性最大似然算法來估計(jì)所有性狀的遺傳參數(shù)。所有的常染色體SNPs 被用來計(jì)算基因組關(guān)系矩陣（GRM），然后計(jì)算主成分。用單變量動(dòng)物模型來估計(jì)遺傳力，用雙變量GREML 分析來估計(jì)性狀間的遺傳相關(guān)。統(tǒng)計(jì)模型為：

其中，是表型觀測值的向量；代表固定效應(yīng)，包括前3 個(gè)主成分、產(chǎn)仔季節(jié)和妊娠天數(shù)，是的相關(guān)矩陣；是隨機(jī)多基因效應(yīng)的向量；是殘差效應(yīng)的向量。

1.4 ssGWAS GWAS 是使用BLUPF90 軟件家族進(jìn)行的。RENUMF90 模塊用于準(zhǔn)備表型和系譜數(shù)據(jù)，BLUPF90 模塊用于計(jì)算基因組估計(jì)育種值（Genomic Estimated Breeding Values，GEBV），POSTGSF90 模塊用于將GEBVs 轉(zhuǎn)換為SNPs 效應(yīng)，并計(jì)算相鄰窗口所解釋的方差。

在本研究中，單性狀動(dòng)物模型被用于GWAS 分析。統(tǒng)計(jì)學(xué)模型可以寫成：

其中，是基于所有動(dòng)物系譜的分子關(guān)系矩陣，是基因分型個(gè)體的分子關(guān)系矩陣的子集，是基因組關(guān)系矩陣。上標(biāo)-1 表示它是相關(guān)矩陣的逆數(shù)。VanRaden（2008）構(gòu)建G 矩陣的方法如下：

其中，Z 是根據(jù)等位基因頻率編碼的標(biāo)記矩陣，D 是SNP 效應(yīng)方差的權(quán)重對(duì)角矩陣（最初D=I）。ssGWAS的計(jì)算方法如下：

①構(gòu)建G 矩陣；

②計(jì)算所有動(dòng)物的GEBVs；

⑤歸一化SNP 權(quán)重以保持總方差不變；

⑥根據(jù)SNP 權(quán)重重新計(jì)算矩陣；='；

⑦循環(huán)步驟②。

從第②步到第⑦步重復(fù)一次迭代過程，按照Wang等的描述，計(jì)算第i 個(gè)連續(xù)SNP（SNP 窗口）解釋的遺傳變異百分比：

Beissinger 等報(bào)告說，與其他較大的窗口尺寸相比，5 個(gè)或10 個(gè)SNP 組成的滑動(dòng)窗口的檢測率與假陽性率都是最有利的。因此，本研究計(jì)算了由10 個(gè)相鄰SNP組成的滑動(dòng)窗口所解釋的加性遺傳變異的比例。

1.5 候選基因組區(qū)域注釋 解釋的遺傳變異大于1%的基因組窗口被認(rèn)為是與性狀相關(guān)的。在11.1基因組（http://www.ensembl.org/）上搜索了候選窗口內(nèi)的基因或部分重疊的基因。使用大鼠基因組數(shù)據(jù)庫（RGD，https://rgd.mcw.edu/rgdweb/homepage/）進(jìn)行人類表型本體論（Human Phenotype Ontology，HPO）分析，選擇Bonferroni 校正后<0.05 的條目進(jìn)行后續(xù)分析。隨后，使用KOBAS 在線網(wǎng)站對(duì)與生殖系統(tǒng)相關(guān)條目中涉及的候選基因進(jìn)行分析，以確定基因參與的KEGG 信號(hào)通路。<0.05 的KEGG 通路被選中。使用PigQTLdb 數(shù)據(jù)庫（http://www.animalgenome.org/cgi-bin/QTLdb/SS/index）對(duì)豬基因組中先前確定的QTL 進(jìn)行評(píng)估。

2 結(jié)果

2.1 表型統(tǒng)計(jì) 表1 列出了3 個(gè)繁殖性狀的描述性統(tǒng)計(jì)。其中，TNB 和NBA 的范圍是0～24，平均值分別為14.91 和13.11。LWB 平均值為15.14。此外，所有性狀的變異系數(shù)值都較高，這表明數(shù)據(jù)在某種程度上是離散的。另外，3 個(gè)性狀之間高度正相關(guān)，表型相關(guān)系數(shù)分別為0.86、0.67 和0.81（表2）。

2.2 遺傳力和遺傳相關(guān) 本研究對(duì)3 個(gè)繁殖性狀的遺傳力和遺傳相關(guān)都進(jìn)行了估計(jì)（表2）。3 個(gè)性狀的估計(jì)遺傳力在0.06～0.12 之間，其中，NBA 遺傳力最高，為0.12，TNB 遺傳力最低，為0.06。3 個(gè)性狀兩兩間的遺傳相關(guān)從0.62 到0.83，都呈現(xiàn)高度正相關(guān)。

表2 繁殖性狀的遺傳力（對(duì)角線，粗體）、遺傳相關(guān)（上三角）和表型相關(guān)系數(shù)（下三角）

2.3 GWAS 本研究將連續(xù)10 個(gè)SNPs 劃分為一個(gè)窗口后，共有43 098 個(gè)窗口用于后續(xù)的GWAS 分析，其中最后9 個(gè)窗口包含的SNPs 少于10 個(gè)。

本研究中3 個(gè)性狀分別鑒定到95、72、109 個(gè)相關(guān)的窗口，它們分別解釋了TNB、NBA 和LWB 總遺傳變異的1% 或以上（圖1）。共有40 個(gè)基因組區(qū)域與繁殖性狀相關(guān)，其中TNB、NBA 和LWB 分別包括15、14、11 個(gè)區(qū)域。對(duì)于TNB，共有15 個(gè)基因組區(qū)域，分布在染色體（SSC）1、3、5、7、13、15、16、17 和X 上，其中最高的區(qū)域解釋了4.794%的遺傳變異。有14 個(gè)基因組區(qū)域與NBA 顯著相關(guān)，它們分布在SSC 1、2、3、6、7、16、17 和X 上，最高的區(qū)域解釋了2.214%的遺傳變異。對(duì)于LWB，109 個(gè)窗口串聯(lián)成11 個(gè)基因組區(qū)域，這些基因組區(qū)域分布在SSC 2、3、6、10、13、16 和X 上，它們所解釋的最大遺傳變異為2.364%（表3）。

表3 繁殖性狀相關(guān)基因組區(qū)域

圖1 繁殖性狀ssGWAS 曼哈頓圖

2.4 候選基因鑒定 經(jīng)過注釋，共得到532 個(gè)基因，包括381 個(gè)編碼蛋白的基因。隨后，381 個(gè)基因被用來進(jìn)行HPO 分析，結(jié)果富集在以下3 個(gè)條目：泌尿生殖系統(tǒng)表型、心血管系統(tǒng)表型和神經(jīng)系統(tǒng)表型（圖2A）。在泌尿生殖系統(tǒng)表型中，共得到81 個(gè)候選基因（圖2B，2C）。有趣的是，大多數(shù)候選基因都位于2 號(hào)染色體上。隨后，對(duì)候選基因進(jìn)行了進(jìn)一步的KEGG 分析，得到了82 條顯著的通路。在82 條顯著的通路中共發(fā)現(xiàn)7 條與生殖相關(guān)的KEGG 通路，包括5 條與生殖激素相關(guān)的信號(hào)通路和2 條卵細(xì)胞減數(shù)分裂及成熟相關(guān)的通路（表4）。對(duì)于所有的候選基因，本研究查詢了以前的研究是否報(bào)道過與生殖性狀有關(guān)的文獻(xiàn)，結(jié)果列于表5。

表4 候選基因參與的顯著的繁殖相關(guān)KEGG 通路

表5 候選基因中先前報(bào)道過的與繁殖相關(guān)的基因

圖2 候選基因鑒定

3 討 論

產(chǎn)仔性狀包含仔豬的TNB、NBA 和LWB，與母豬的繁殖性能密切相關(guān)，影響豬的生產(chǎn)效率和經(jīng)濟(jì)利益。因此，探索豬產(chǎn)仔性狀的遺傳結(jié)構(gòu)是極為重要的。此外，探索產(chǎn)仔性狀的遺傳機(jī)制非常復(fù)雜，整個(gè)繁殖過程受多種激素和生理途徑的調(diào)節(jié)。GWAS 是檢測與目標(biāo)性狀相關(guān)的候選基因最有效方法之一。本研究利用ssGWAS方法探索了392 頭大白豬的TNB、NBA 和LWB 的遺傳變異，通過HPO 和KEGG 分析進(jìn)一步挖掘了與生殖相關(guān)的基因以及這些基因涉及的通路。

繁殖性狀的遺傳力通常在0.1 左右。在本研究中，3個(gè)繁殖性狀的估計(jì)遺傳力為0.06～0.12，與前人研究結(jié)果一致。此外，3 個(gè)性狀之間的遺傳和表型相關(guān)都呈現(xiàn)出高度正相關(guān)，這些結(jié)果與Wolf 等的研究一致。

本研究發(fā)現(xiàn)了40 個(gè)解釋1%以上遺傳變異的基因組區(qū)域。在這些區(qū)域中，還發(fā)現(xiàn)了一些以前報(bào)道的與繁殖有關(guān)的QTL。在3 個(gè)性狀中SSC13 上的48.65～52.85 Mb 區(qū)域包含了最多解釋超過1% 的遺傳變異的窗口，并且在這個(gè)基因組區(qū)域中有9 個(gè)先前報(bào)道的相關(guān)QTLs。這些QTLs 包括3 個(gè)黃體數(shù)量QTLs（493、8 826、24 285），2 個(gè)卵巢重量QTLs（18 354、18 356），1 個(gè)青春期年齡QTLs（8 824），1 個(gè)子宮角重量QTLs（18 409），1 個(gè)總產(chǎn)仔數(shù)QTLs（15 1333）和1 個(gè)死胎數(shù)QTLs（57 525）。此外，SSC14 的104.44～106.90 Mb 區(qū)域解釋了2.87%的TNB 性狀的遺傳變異，并且在此區(qū)域注釋到了基因。據(jù)報(bào)道，基因在胚胎發(fā)育中起著不可或缺的作用。

基于關(guān)聯(lián)分析鑒定到的顯著基因組區(qū)域內(nèi)的蛋白編碼基因，本研究進(jìn)行了HPO 分析，發(fā)現(xiàn)這些基因聚成了3 類：泌尿生殖系統(tǒng)、心血管系統(tǒng)和神經(jīng)系統(tǒng)。共有7 種表型涉及泌尿生殖系統(tǒng)，包括原發(fā)性閉經(jīng)、閉經(jīng)、月經(jīng)周期異常、女性生殖系統(tǒng)生理功能異常、泌尿生殖系統(tǒng)異常、生殖系統(tǒng)異常以及生殖系統(tǒng)生理學(xué)異常。這7 個(gè)表型共包含81 個(gè)候選基因，且大多聚集在SSC 2 上。經(jīng)過搜索文獻(xiàn)，在SSC2 上獲得了10 個(gè)先前報(bào)道的與生殖有關(guān)的基因。其中，畸形表皮自動(dòng)調(diào)節(jié)因子1（DEAF1）編碼一個(gè)多域轉(zhuǎn)錄因子。相關(guān)結(jié)合位點(diǎn)的變化可能與山羊的產(chǎn)仔數(shù)有關(guān)?；騾⑴c磷酸戊糖途徑（PPP）的第2 階段，并在糖酵解和PPP 之間提供一個(gè)鏈接。缺陷的患者在新生兒期出現(xiàn)低出生體重，甚至在嬰兒期早期死亡。胰島素樣生長因子2（）基因是一個(gè)父系印記基因，已知它能增加豬的肌肉生長并減少脂肪沉積。然而，也有報(bào)道稱對(duì)母豬的繁殖力有影響。此外，基因已被確認(rèn)為人類、小鼠和牛的一個(gè)印記基因。印記基因與早期胚胎和胎盤的發(fā)育有關(guān)系。這些基因的表達(dá)水平的變化可以導(dǎo)致異常的生長和胚胎死亡。據(jù)報(bào)道，在豬的心臟、卵巢等部位有母體表達(dá)，在肺和腎有雙側(cè)表達(dá)。因此，基因在豬身上的印記還有待驗(yàn)證?；蚴莝irtuin家族的一個(gè)重要成員，在調(diào)節(jié)細(xì)胞生理和衰老方面發(fā)揮了廣泛的作用。Gorczyca 等發(fā)現(xiàn)基因在豬卵巢卵泡細(xì)胞中表達(dá)。

為進(jìn)一步了解候選基因的功能和這些基因之間的聯(lián)系，本研究還進(jìn)行了KEGG 分析，在所有顯著的KEGG 信號(hào)通路中得到了7 條與繁殖性狀相關(guān)的通路。這7 條通路分別是催產(chǎn)素信號(hào)通路、催乳素信號(hào)通路、GnRH 信號(hào)通路、松弛素信號(hào)通路、雌激素信號(hào)通路、母細(xì)胞減數(shù)分裂和孕激素介導(dǎo)卵細(xì)胞成熟通路。其中催產(chǎn)素、催乳素、促性腺激素釋放激素（GnRH）、松弛素和雌激素都屬于生殖激素，與生殖活動(dòng)如排卵、懷孕、胚胎發(fā)育、分娩和泌乳直接相關(guān)。參與生殖激素相關(guān)途徑的共有7 個(gè)基因，即和CREB3L3。其中，MAP2K2編碼的蛋白質(zhì)屬于MAPK 激酶家族，基因?qū)μケP形成期間合體滋養(yǎng)層的正常發(fā)育有貢獻(xiàn)?；?qū)儆赗as 腫瘤基因家族，其編碼的產(chǎn)物在信號(hào)轉(zhuǎn)導(dǎo)途徑中發(fā)揮作用。有報(bào)道稱miR-210 通過和調(diào)節(jié)水牛排卵前的顆粒細(xì)胞功能?；騾⑴c調(diào)控血管生長和基礎(chǔ)血壓。Zou 等研究表明基因在胎盤血管的形成中有重要作用?；蚓幋aAMP 依賴性轉(zhuǎn)錄因子家族，它可以調(diào)節(jié)脂質(zhì)代謝并介導(dǎo)對(duì)感染的反應(yīng)。相關(guān)研究表明通過控制GH 信號(hào)傳導(dǎo)，對(duì)身體的生長有負(fù)向調(diào)節(jié)作用。此外，有2 個(gè)基因參與卵母細(xì)胞減數(shù)分裂和成熟相關(guān)通路，分別是和。Aurora-A（AURKA）是一種絲氨酸/蘇氨酸蛋白激酶。Yao 等研究結(jié)果顯示基因在豬卵母細(xì)胞減數(shù)分裂成熟、受精和早期胚胎有絲分裂過程中對(duì)微管的組裝起著關(guān)鍵的調(diào)節(jié)作用?；蚓幋a胰島素，這是一種多肽激素，可調(diào)節(jié)生理學(xué)的各個(gè)方面，是決定后生動(dòng)物雌性生殖的關(guān)鍵因素。據(jù)報(bào)道，胰島素參與了奶牛的卵泡成熟。

4 結(jié) 論

本研究使用50K SNP 芯片對(duì)法系大白豬群體中的TNB、NBA 和LWB 表型進(jìn)行了ssGWAS 分析，鑒定到了40 個(gè)基因組區(qū)域。區(qū)域內(nèi)注釋到的所有編碼蛋白的基因被用于進(jìn)行HPO 分析，發(fā)現(xiàn)了81 個(gè)與繁殖有關(guān)的候選基因，這些基因主要涉及GnRH 信號(hào)通路、雌性激素信號(hào)通路和卵細(xì)胞成熟通路等。其中，本研究還發(fā)現(xiàn)一些與產(chǎn)仔數(shù)有關(guān)的基因，包括、、、、等。盡管篩選出了一些與繁殖相關(guān)的候選基因，但部分基因在畜禽中還沒有被研究過，因此這些候選基因的功能還需要進(jìn)一步探討。下一步，將利用高密度測序數(shù)據(jù)對(duì)這些區(qū)域進(jìn)行精細(xì)定位，以進(jìn)一步確定與產(chǎn)仔數(shù)相關(guān)的基因并篩選出致因突變，為豬的遺傳機(jī)制分析和育種提供參考。