基于全基因組重測序的山羊產羔數性狀關鍵調控基因的篩選

李恒，字向東，王會，熊燕，呂明杰，劉宇，蔣旭東

基于全基因組重測序的山羊產羔數性狀關鍵調控基因的篩選

李恒1，字向東1，王會2，熊燕1，呂明杰1，劉宇1，蔣旭東1

1西南民族大學動物科學國家民委重點實驗室，成都 610041；2西南民族大學青藏高原動物遺傳資源保護與利用教育部重點實驗室，成都 610041

【目的】對產羔數不同的山羊進行全基因組重測序分析，挖掘參與調控川中黑山羊產羔數性狀關鍵調控基因，為解析山羊產羔數性狀遺傳機制及分子遺傳改良提供理論依據。【方法】選擇6只產4—6羔的川中黑山羊為高繁組（high fecundity, HF）和6只產單羔的川中黑山羊為低繁組（low fecundity, LF），采集頸靜脈血液樣本提取基因組DNA，構建350 bp雙末端測序文庫，利用Illumina HiSeq PE150平臺對12個文庫進行全基因組重測序。測序產出的凈數據經BWA軟件比對至山羊參考基因組ARS1，所獲得的高質量SNPs通過兩種全基因組掃描分析方法（、）的綜合分析確定候選區域，候選區域的注釋基因分別利用g:Profiler和KOBAS在線數據庫進行GO分析與KEGG通路分析，以篩選調節川中黑山羊產羔數性狀候選基因。為了進一步鑒定調節山羊產羔數目的關鍵遺傳標記，根據基因組重測序變異分析，對繁殖候選基因的同義與非同義SNPs進行定位篩選，后續將12個山羊樣本的擴增產物進行Sanger測序以驗證重測序結果。【結果】12只山羊共獲得431.50 Gb 凈數據，經變異檢測與注釋發現，HF組山羊共發現7 771 417個單核苷酸多態性(single nucleotide polymorphism, SNPs)，LF組山羊檢測到8 935 907個SNPs，且LF組各類SNPs 均多于HF組。設置同時達到top 5%最大值和top 5%最小值的窗口為候選區域，在低雜合性、高遺傳分化的區域共注釋130個強選擇信號，其中HF組、LF組以及共享窗口的注釋基因分別為84、59和13個，經GO富集與KEGG通路分析發現，19個候選基因參與川中黑山羊的繁殖、繁殖過程和胚胎發育等調控，包括11個HF組特異性候選基因（和），5個LF組特異性候選基因（、、、和）和3個HF組與LF組共享窗口基因（、和）。同時，大多數GO分析，如G蛋白偶聯受體活性、激素反應和神經肽信號通路等，都包含這19個候選基因。此外，14個HF候選基因有9個顯著富集在代謝途徑、神經活性配體-受體相互作用、糖胺聚糖-硫酸乙酰肝素/肝素的生物合成、鈣離子信號通路、cAMP信號通路和葉酸生物合成等KEGG通路中（＜0.05）。19個繁殖候選基因中共有2個同義突變（，10 A4662G）與2個非同義突變（G529A，），且僅定位于HF候選基因中。Sanger測序發現，、和突變位點均可檢測到多態性，與基因組重測序結果一致，其中G529A多態性導致第177位丙氨酸突變為蘇氨酸，多態性導致翻譯提前終止。【結論】本研究共發現11個HF組特異性候選基因，推測是川中黑山羊多羔性狀的關鍵調控基因，外顯子G529A與外顯子A281T突變可能是調控山羊多羔性狀的關鍵遺傳標記，在改良山羊繁殖性能方面具有較大的應用價值。

川中黑山羊；基因組重測序；多羔性狀；候選基因

0 引言

【研究意義】產羔性狀是山羊的重要經濟性狀，提高窩產羔數不僅可以提高出欄率，提高經濟效益，而且還可以提高選擇強度，加快山羊育種進程。產羔數是一個低遺傳力（0.08—0.18）的限性性狀，難以用常規育種方法進行快速改良，但適合用標記輔助選擇（marker-assisted selection, MAS）等分子育種技術來改良[1]。實現MAS的先決條件是找到與數量性狀基因座相連鎖的分子遺傳標記，而要找到輔助選擇的分子標記，則需要解析山羊高繁殖力的分子遺傳基礎和形成原因。【前人研究進展】在過去20多年里，在一些綿羊品種的突變體中已成功分離鑒定出控制卵泡發育和排卵數的多胎基因，主要包括骨形態發生蛋白受體-1B（BMP receptor-1B，）基因的1個突變（FecB）[2]、骨形態發生蛋白-15（bone morphogenetic protein-15，）基因的6個突變（FecX、FecX、FecX、FecX、FecX和FecX）[3-4]和生長分化因子-9（growth and differentiation factor-9，）基因的2個突變（FecG和）[5-6]，這些突變與綿羊的多羔性狀密切相關。但是，迄今為止，在已開展過檢測的印度[7]、伊朗[8]、泰國[9]和我國[10-11]的30多個山羊品種（類群）中均未檢測到這些突變，排除了這些突變作為山羊品種多羔性狀主基因的可能性，山羊多羔性狀的遺傳機制有待研究。山羊全基因組測序技術的發展與完善，為發掘參與調控川中黑山羊產羔數性狀關鍵調控基因提供了全新的視角[12]。對嶗山奶山羊高、低繁殖力的兩個極端種群進行基因組重測序，發現細胞周期蛋白B2（cyclin B2,2）雄激素受體（androgen receptor,）腺苷酸環化酶1（adenylate cyclase 1,）SMAD家庭成員2（SMAD family member 2,）等基因在高繁殖力群體中被特異性選擇[13]。絲氨酸/蘇氨酸激酶3（serine/threonine kinase 3,）蛋白磷酸酶3催化亞基α （protein phosphatase 3 catalytic subunit alpha,）等96個候選基因與大足黑山羊產羔數性狀顯著相關[14]，其中也與阿爾巴斯絨山羊產羔數性狀相關[15]。進一步分析67個單核苷酸多態性（SNP）與大足黑山羊產羔數的關聯，發現半胱氨酰轉運RNA合成酶2（cysteinyl-tRNA synthetase 2，）Ⅰ型血小板結合蛋白基序的解聚蛋白樣金屬蛋白酶（ADAM metallopeptidase with thrombospondin type 1 motif 14，）和甲基轉移酶 25（methyltransferase like 25，）基因編碼區SNP與頭胎產羔數目顯著相關[16]。對頭胎產單羔、雙羔、三羔的三組濟寧青山羊進行全基因組掃描，結果在雙羔和三羔群體中發現酪氨酸激酶受體（KIT proto-oncogene, receptor tyrosine kinase，）、P21-活化激酶1P21-activated kinases 1，腺苷活化蛋白激酶α1亞基（protein kinase AMP-activated catalytic subunit alpha 1，）等多個產羔數調節基因[17]。WANG等[18]在相同飼養管理條件下，對在5個連續分娩記錄中表現出穩定的產羔數差異的12只山羊進行基因組重測序，結果發現細胞分裂周期蛋白25同源蛋白C（cell division cycle 25C，）、核酸內切酶G（endonuclease G，）和Nanos同源基因3（nanos C2HC-type zinc finger 3，）的變異與產羔數性狀顯著相關，螺旋卷曲過程對生殖能力有潛在的調控作用。這些候選基因的發現拓展了人們對繁殖力遺傳基礎的認知，為山羊產羔數性狀的遺傳機制的解析提供重要線索。【本研究切入點】川中黑山羊是我國優良地方山羊品種，具有生長速度快，產肉性能高，產羔率高等優良特性[19-20]。我們在川中黑山羊和的前期研究中未檢測到與綿羊多羔性狀相關突變的存在，但高繁川中黑山羊與單胎藏山羊之間檢測到5個新的堿基突變，其中4個導致氨基酸改變，有2個堿基突變，其中1個導致氨基酸突變，兩個品種的15核苷酸序列則完全一致[21]。說明山羊和中氨基酸改變的突變與綿羊中不同，這些突變是否與山羊產羔數相關尚有待進一步研究。從卵巢轉錄組和miRNA水平也不能完全解析川中黑山羊多羔性狀的分子遺傳機制[22-23]。【擬解決的關鍵問題】本研究選擇高、低繁殖力的川中黑山羊為研究對象，利用基因組重測序技術掃描其低雜合性、高遺傳分化區域，挖掘高、低繁殖力山羊關鍵基因遺傳變異，并對候選靶基因進行GO富集與KEGG通路分析，以期篩選調節川中黑山羊產羔數的關鍵候選基因，為山羊產羔數性狀遺傳機制的闡釋及后續山羊分子育種工作提供理論基礎。

1 材料與方法

1.1 樣品的采集

本研究根據第一、二胎產羔記錄，于2019年9月選擇飼養于四川省樂至縣天龍育種場，胎產羔數差異顯著的高繁殖力（high fecundity，HF）和低繁殖力（low fecundity，LF）三歲齡川中黑山羊各6只，其中HF組山羊連續兩胎胎產羔數4—6只，LF組連續兩胎胎產羔數1只（表1）。每只山羊采集其頸靜脈血液2 mL于EDTA-K2抗凝管中，試劑盒提取基因組DNA后用于基因組重測序研究。

1.2 山羊血液基因組DNA的提取及質量檢測

使用1%瓊脂糖凝膠電泳與紫外分光光度計分別鑒定山羊基因組DNA的完整性、濃度和純度，12管質檢合格的HF組和LF組DNA樣品各抽取1.5 μg分別建庫。利用Agilent 2100檢測文庫中插入片段的長度，片段符合預期后，送至Illumina HiSeq PE150平臺進行雙端測序，整個流程委托北京諾禾致源生物信息科技有限公司完成。

表1 12只川中黑山羊的胎產羔數信息

編號HF1—6為高繁殖力山羊，編號LF1—6為低繁殖力山羊

The symbols from HF1 to HF6 were high fecundity goats, and the symbols from LF1 to LF6 were low fecundity goats

1.3 測序數據的變異檢測與注釋

原始數據經初步質控處理，將接頭序列、未測出堿基超過10%和低質量堿基數超過50%的paired reads去除，獲得凈數據。利用BWA[24]軟件將clean reads比對至山羊參考基因組ARS1[25]，比對結果經SAMTOOLS軟件去除重復[26]。利用GATKv3.2-2進行SNP鑒定與分型[27]，過濾后獲得的高質量SNPs用ANNOVAR[28]軟件工具進行功能注釋。

1.4 群體選擇分析

利用滑動窗口法（sliding windows）對HF、LF組川中黑山羊進行全基因組掃描，按照窗口內SNPs個數不超過20，確定選擇消除窗口大小（100 kb）。計算每個窗口的固定系數（Fst）值并轉換為值，選擇top 5% 最大t值窗口為選擇區域。通過計算窗口內SNPs位點的雜合性（Hp），進而對Selective sweep進行評估。分別計算窗口的Hp值并轉化為值，選擇top 5%最小值窗口作為選擇區域。在計算得到每個窗口Hp和Fst值的基礎上，選擇最小值和最大值均達到Top 5%的區域為候選區域，并對該區域進行基因注釋。

1.5 候選基因的功能注釋

基于Bio-Mart在線數據庫獲得山羊候選區域內的基因注釋信息（基因類型和ENSEMBL號）后，利用g:Profiler網站將山羊基因符號ID轉換為小鼠同源基因。利用g:Profiler[29]網站和KOBAS[30]在線數據庫分別進行Gene ontology（GO）富集分析和Kyoto Encyclopedia of Genes and Genomes（KEGG）通路分析。

1.6 SNP位點的驗證

為進一步鑒定調節山羊產羔數目的關鍵遺傳標記，據重測序變異分析報告對19個繁殖候選基因的同義與非同義SNPs進行定位篩選。候選位點經Primer Premier 5.0 軟件設計引物（表2），12個山羊樣本的擴增產物送至成都擎科生物有限公司進行Sanger測序。

2 結果

2.1 高、低繁殖力組山羊基因組測序數據質控

分別完成12只山羊血液樣本基因組重測序工作，共獲得431.50 Gb 凈數據，錯誤分布率均為0.03%，測序質量較高（Q20≥96.9%，Q30≥91.69%），所有樣本比對率在99.61%—99.73%之間，平均測序深度10.70×，1×覆蓋度在94.91%—95.32%之間，4×覆蓋度在88.70%—93.13%之間。所有原始序列數據已上傳至國家基因組數據中心基因組序列檔案庫（https://bigd.big.ac.cn/gsa.），注冊號為CRA003846。

表2 SNP位點驗證引物

2.2 變異檢測與注釋

基于嚴格的閾值，兩組山羊共檢測到16 707 324個高質量SNPs位點，其中LF組8 935 907個SNPs位點，HF組7 771 417個SNPs位點，且LF組的各類SNPs 均多于HF組（表3）。在這些SNPs中，HF組和LF組山羊外顯子SNPs平均各有53 425（0.68%）和60 588（0.69%）個，基因上游1 kb區域分別有44 949個和50 798個，基因下游1 kb區域分別分布44 201個和50 349個，兩組基因間區SNPs占比最多，分別占總SNPs的70.01%和70.09%；其次是分布在內含子區域SNPs，HF和LF組分別含有2 145 836個和2 473 891個；位于剪切位點的SNPs數目最少。在編碼區SNPs 中，同義突變和錯義突變分別占外顯子總SNPs的58.66%和40.77%，且兩組相差不大。此外，HF組和LF組在終止密碼子處分別檢測到300和323個SNPs，在剪切位點處分別檢測到214個和237個SNPs，這些位點可能會影響轉錄剪接，從而影響蛋白質產物及功能。

表3 川中黑山羊SNPs檢測及注釋結果統計

2.3 群體選擇分析

經不同大小窗口調試發現，窗口大小為100 kb時，SNPs數目小于20個并逐漸趨于穩定，因此，選擇100 kb作為選擇信號檢測最佳窗口長度，使用50%的重疊區為滑動尺，計算川中黑山羊HF組與LF組之間的Fst值，并標準化作圖，選取LF和HF之間的閾值為1.9，共篩選到579個窗口，注釋623個基因（圖1）。

分別計算川中黑山羊HF組與LF組染色體窗口內SNP位點的值，Z轉換后并標準化作圖，選擇top 5%最小值為閾值，HF群體的閾值為-2.05，篩選到579個窗口（圖2），共注釋617個候選基因。 LF群體的閾值為-2.08，篩選到579個窗口（圖3），注釋到610個基因。

圖1 高繁組和低繁組川中黑山羊1-29號常染色體平均固定系數ZFst的分布

圖2 高繁組川中黑山羊1-29號常染色體平均雜合度ZHp的分布

為進一步篩選川中黑山羊產羔性狀關鍵調控基因，本研究初步掃描了具有低雜合性、高遺傳分化的區域，即候選窗口值均達到top 5%最大值和top 5%最小值標準。其中HF組中鑒定65個候選窗口，共注釋84個候選基因（圖4-A），LF組山羊共鑒定了42個候選窗口，注釋了59個候選基因（圖4-B）。因兩組候選基因中存在、等13個交集基因，所以HF組和LF組共注釋了130個差異候選基因。

2.4 高、低繁殖力山羊候選基因基因GO富集分析

進一步揭示與產羔數性狀的候選靶基因功能，本研究對兩群體差異表達基因（130個）進行生物功能分析，-value經Bonferroni校正后，設置corrected-value ≤0.05為閾值，選擇候選靶基因中顯著富集的GO terms。

圖3 低繁組川中黑山羊1-29號常染色體平均雜合度ZHp的分布

圖4 川中黑山羊高繁組（A）與低繁組（B）受選擇區域

GO結果顯示，川中黑山羊群體中篩選的130個強選擇性候選基因，共富集到193條顯著性GO功能條目。其中歸類為生物學過程的GO條目最多，共113條，主要行使生長發育和代謝功能，其次是細胞組成，共50條GO條目，分子功能本體條目最少，其富集的30條GO條目主要參與跨膜信號受體活性、神經肽受體活性和蛋白多糖生物合成的過程等。

值得注意的是，19個基因與繁殖緊密相關，如繁殖、繁殖過程、胚胎發育和生殖過程調控（表4）。其中HF組特異性強選擇性基因為11個，包含腺苷酸環化酶10（adenylate cyclase 10，）、多巴胺受體D1（dopamine receptor D1，）、硫酸肝素6-鄰磺酸轉移酶1（heparan sulfate 6-O-sulfotransferase 1，）、胰島素樣生長因子結合蛋白7（insulin like growth factor binding protein 7，）、促黑激素同源框2（msh homeobox 2，）、頭蛋白（noggin，）、腎單位腎癆4（nephronophthisis 4 (juvenile) homolog，）、妊娠相關血漿蛋白A（pregnancy- associated plasma protein-A，、睪丸發育相關蛋白（testis development-Related Protein，）、木糖轉移酶1（xylosyltransferase 1，）和催乳素釋放激素受體（prolactin releasing hormone receptor，），HF和LF組共享窗口基因3個，如醛酮還原酶家族成員B3（aldo-keto reductase family 1，member B3，）、組蛋白脫乙酰酶4（histone deacetylase 4，）和mu阿片受體（opioid receptor mu 1，），以及LF組特異性候選基因5個，如膜聯蛋白5（annexin A5，）、內皮素A型受體(endothelin receptor type A，)、FA核心復合體連接酶（FA complementation group L，）、胰島素樣生長因子1（insulin-like growth factor1，）和速激肽前體1（tachykinin precursor 1，）。同時，大多數GO術語，如G蛋白偶聯受體活性、激素反應和神經肽信號通路等，都包含這19個候選基因。

表4 川中黑山羊繁殖相關GO條目

2.5 高、低繁殖力山羊候選基因KEGG分析

KEGG通路分析結果顯示，川中黑山羊差異候選基因集共富集112條代通路，選擇前20條顯著的代謝通路作為展示（圖5）。KEGG代謝通路結果表明，14個HF山羊候選基因中，9個基因顯著富集在代謝途徑（、）、糖胺聚糖-硫酸乙酰肝素/肝素的生物合成（、）、神經活性配體-受體相互作用（）、鈣離子信號通路（）、cAMP信號通路（）和葉酸生物合成（）通路（corrected-value ≤0.05）中，這些通路可能在山羊產羔數性狀的調控中極具關鍵作用。

2.6 SNP位點的驗證

為篩選多羔性狀關鍵遺傳標記，本研究掃描LF和HF組外顯子SNPs分布。19個繁殖候選基因中，共有2個同義突變（，）與2個非同義突變（G529A，），且僅定位于HF候選基因中。其中G529A多態性導致第177位丙氨酸突變為蘇氨酸，多態性導致翻譯提前終止。重要的是，4個SNPs僅有3個與高通量測序結果呈現相同趨勢（，G529A，），的外顯子SNP測序結果與重測序變異報告不一致，未發現該突變位點（圖6）。

3 討論

山羊產羔數是復雜的表型特征，其不僅由排卵率、胚胎或胎兒成活率、子宮容受性等多個次級性狀調控，也受遺傳、季節、降雨、氣溫和營養等因素影響[31-33]。GWAS雖是檢測復雜性狀的遺傳學基礎的最有力的工具[34]，但深入挖掘GWAS的統計能力需要巨大的樣本量和適當的對照設計[35]。因此，大樣本、家系信息詳盡等要求限制了GWAS在山羊群體中的應用。基因組重測序仍是目前分析山羊復雜性狀的強有力手段。該技術不僅可發掘出特定性狀的多個候選基因，其公布的變異數據也為其他性狀研究提供思考。例如最新發現的山羊產羔數候選基因、等[36-37]，遺傳標記信息均是來自先前的山羊基因組重測序報道。

圖5 川中黑山羊選擇信號KEGG通路分析

本研究對高產山羊品種川中黑山羊不同產羔數群體進行了全基因組測序，獲得了大量遺傳變異，為山羊高繁殖力遺傳機制解析提供了重要的遺傳資料。為闡明川中黑山羊產羔性狀的選擇特征，設置top 5%為閾值[13-14]，在川中黑山羊群體的低雜合性、高遺傳分化的區域共注釋130個候選基因。

通過GO富集進一步發現了19個候選基因可能參與川中黑山羊產羔性狀的調節，參與繁殖、繁殖過程、胚胎發育等生殖過程調控。在這19個候選基因中，11個候選基因在HF組中特異性選擇，可能與川中黑山羊的多羔性狀有關。這些基因的功能主要與生殖細胞發育及胚胎發育等生殖過程相關。其中編碼的可溶性腺苷酸環化酶（SAC），對男性生殖至關重要[38]。SAC調控cAMP/PKA通路與AMPK的活性，其中cAMP/PKA通路在顆粒細胞FSH的功能調節和E2產生過程中發揮重要作用[39]，AMPK的活性影響卵母細胞成熟[40]。推測SAC通過調節生殖激素分泌水平進而調控動物繁殖過程，其表達的豐度制約了雌性個體排卵率。與與家禽產卵相關[41-42]，可能在山羊群體的排卵機制中發揮特殊作用。可能是卵巢早衰病因學的候選基因[43]，該基因突變也被認為與特發性低促性腺激素減退癥有關[44]。牛繁殖候選基因下調促進了細胞凋亡，阻礙了細胞增殖，增加了孕酮和雌二醇的產生[45]，其還與[46]功能類似，對胚胎著床、妊娠的建立與維持至關重要[47]。抑制IGFBP7表達可顯著降低植入胚胎數與妊娠率[47]，這意味著除排卵數外，早期胚胎的丟失也是山羊窩產羔數目減少的誘因。目前并未發現功能突變與山羊早期胚胎丟失有關，其調控機制還需深入研究。在卵巢中表達并與BMP蛋白相互作用，參與正常胚胎的發育調控[48]。PAPAA蛋白的缺乏會下調小鼠卵巢類固醇生成和損害女性生育能力[49]，其mRNA水平與妊娠期間血清水平分別是識別卵母細胞能正常發育至胚胎和評估子宮內膜容受性的潛在重要標志[50-51]。PAPAA蛋白具有成為預測高繁殖力個體診斷工具的潛力。和是參與精子發育的關鍵基因[52-53]。是一種木糖基轉移酶，參與對著床期囊胚的附著和生長有直接作用的蛋白多糖肝素/硫酸肝素合成過程[54]。這些證據均表明了這11個繁殖調控基因可能是調控山羊多羔性狀的關鍵基因，在改良山羊繁殖性能方面具有較大的應用價值，但其在產羔調節網絡中的具體功能還需進一步研究。

A：PRLHR，B：MSX2，C：DRD1，D：ADCY10

LF組差異表達基因（和）在川中黑山羊排卵、植入后胚胎發育等生殖過程可能呈現負調控作用。研究表明，可抑制應激條件引發的牛卵丘細胞的凋亡[55]，對牛卵母細胞發育可能具有增益作用。其5′側翼區雜合SNPs雖不影響其在卵巢中的表達，但對山羊的多產性具有顯著影響[56]，且編碼區同義突變還可通過調控基因表達影響與受體結合的親和力[57]，進一步影響動物繁殖。因此，后續擬定分析部分SNPs是否與山羊排卵率相關，對于解析川中黑山羊產羔性狀遺傳基礎可能具有理論指導意義。啟動子區域內的SNPs可下調其基因的表達[58]，進而導致胎盤介導的妊娠并發癥及女性的反復流產[59-60]。EDNRA作為一種內皮素受體，不僅參與山羊毛色調控[61]，還在肥胖引起的小鼠排卵缺失過程中發揮重要作用[62]。突變可破壞DNA損傷修復過程，可能是引發人類卵巢早衰的潛在病因[63-64]，具體的生殖功能調控機制還有待探索。基因編碼速激肽家族成員P（SP）和神經激肽A（Neurokinin A，NKA）[65]，參與調控哺乳動物生殖激素分泌[66]。速激肽信號對排卵前促黃體素分泌激增具有促進作用，也被認為是調節排卵的一個負面因素，其功能障礙導致雌鼠排卵障礙[67]。此外，川中黑山羊LF組次優選擇基因與美姑山羊的繁殖性狀相關[68]。這些發現均為產羔數調控網絡解析提供了新視角。

HF與LF組間存在多個共享窗口基因，正如之前嶗山奶山羊與大足黑山羊的相關報道[13-14]，川中黑山羊19個繁殖候選基因中也存在3個HF和LF組的交集基因（）。其中調控前列腺素水平[69]，參與睪酮的生成[70]，存在于卵母細胞和顆粒細胞中，體外成熟液中其激動劑的添加可提高卵母細胞的囊胚率[71]，然用其激動劑孵育精子，囊胚率的結果相反[72]。這些基因的功能與動物繁殖緊密相關，其功能的深入挖掘對山羊繁殖學研究是極為重要的，但該相關報道并未對共享窗口基因功能定位做深入剖析，這3個共享基因在川中黑山羊產羔性狀調控機制中發揮主控或微效調控作用還有待探查。

進一步的功能分析表明，神經活性配體-受體相互作用、糖胺聚糖-硫酸乙酰肝素/肝素的生物合成、鈣離子信號通路、cAMP信號通路和葉酸生物合成通路，在動物繁殖活動中也發揮重要功能。例如硫酸肝素蛋白聚糖是在GDF9信號通路中起著重要作用，并參與排卵周期卵泡中卵母細胞信號和卵丘細胞功能的形成[73]。鈣離子信號通路參與調控哺乳動物卵母細胞發育[74-75]。cAMP通路調節卵巢類固醇P4、E2的產生[76]，參與卵泡顆粒細胞分化和成熟[77]。神經活性配體-受體相互作用途徑不僅調節家禽卵泡發育和卵子生產[78-79]，還參與軟體動物性腺發育和產卵[80-81]。推測這些通路可能參與川中黑山羊卵母細胞發育及排卵過程，并且對排卵率至關重要。此外，葉酸(維生素B9)是機體健康和發育的最佳營養必需物質。葉酸缺乏會引發貧血、生殖健康和胎兒發育受損[82-83]。多個HF特異性強選擇候選基因在這些通路中顯著富集，暗示這些通路可能是川中黑山羊產羔數性狀關鍵調控途徑，對其多羔性狀調節起促進作用。

本研究還特別分析19個候選基因外顯子區域的SNPs，以期進一步篩選山羊多羔機制關鍵分子標記。發現導致氨基酸序列變化的外顯子突變只在HF組中，并且僅定位于（G529A）與（A281T）候選基因中。測序結果表明，外顯子G529A與A281T突變在人工選擇條件下被強烈選擇，表現出較高的遺傳分化，改變了、的翻譯，可能在山羊繁殖力中發揮關鍵作用，有望作為高產標記用于山羊分子育種。后續在川中黑山羊群體中擴群驗證G529A與A281T突變位點，將成為本課題組下一步工作的方向與重點。

總而言之，本研究在川中黑山羊群體中發現了許多新的變異體。這些可能是調控生殖的關鍵的基因。其中、、、、、和可能是調控山羊多羔性狀的關鍵基因，同時這些候選基因大部分顯著富集在神經活性配體-受體相互作用、糖胺聚糖-硫酸乙酰肝素/肝素的生物合成、鈣離子信號通路、cAMP信號通路和葉酸生物合成等潛在的生殖相關通路中。外顯子G529A與突變可能是川中黑山羊多羔性狀的關鍵遺傳標記，在改良山羊繁殖性能方面可能具有較大的應用價值。

4 結論

本研究篩選到19個與產羔性狀相關的新候選基因，包括11個HF組特異性候選因（和），5個LF組特異性候選因（和）和3個LF組和HF組共享窗口基因（和）。

[1] DE LIMA L G, DE SOUZA N O B, RIOS R R, DE MELO B A, DOS SANTOS L T A, DE MORAES SILVA K, MURPHY T W, FRAGA A B. Advances in molecular genetic techniques applied to selection for litter size in goats (): a review. Journal of Applied Animal Research, 2020, 48(1): 38-44. doi:10.1080/09712119.2020.1717497.

[2] MULSANT P, LECERF F, FABRE S, SCHIBLER L, MONGET P, LANNELUC I, PISSELET C, RIQUET J, MONNIAUX D, CALLEBAUT I, CRIBIU E, THIMONIER J, TEYSSIER J, BODIN L, COGNIé Y, CHITOUR N, ELSEN J M. Mutation in bone morphogenetic protein receptor-IB is associated with increased ovulation rate in Booroola Mérino ewes. Proceedings of the National Academy of Sciences of the United States of America, 2001, 98(9): 5104-5109. doi:10.1073/pnas.091577598.

[3] GALLOWAY S M, MCNATTY K P, CAMBRIDGE L M, LAITINEN M P E, JUENGEL J L, JOKIRANTA T S, MCLAREN R J, LUIRO K, DODDS K G, MONTGOMERY G W, BEATTIE A E, DAVIS G H, RITVOS O. Mutations in an oocyte-derived growth factor gene () cause increased ovulation rate and infertility in a dosage- sensitive manner. Nature Genetics, 2000, 25(3): 279-283. doi:10. 1038/77033.

[4] MARTINEZ-ROYO A, JURADO J J, SMULDERS J P, MARTí J I, ALABART J L, ROCHE A, FANTOVA E, BODIN L, MULSANT P, SERRANO M, FOLCH J, CALVO J H. A deletion in the bone morphogenetic protein 15 gene causes sterility and increased prolificacy in Rasa Aragonesa sheep. Animal Genetics, 2008, 39(3): 294-297. doi:10.1111/j.1365-2052.2008.01707.x.

[5] HANRAHAN J P, GREGAN S M, MULSANT P, MULLEN M, DAVIS G H, POWELL R, GALLOWAY S M. Mutations in the genes for oocyte-derived growth factors GDF9 and BMP15 are associated with both increased ovulation rate and sterility in Cambridge and belclare sheep (). Biology of Reproduction, 2004, 70(4): 900-909. doi:10.1095/biolreprod.103.023093.

[6] NICOL L, BISHOP S C, PONG-WONG R, BENDIXEN C, HOLM L E, RHIND S M, MCNEILLY A S. Homozygosity for a single base-pair mutation in the oocyte-specific GDF9 gene results in sterility in Thoka sheep. Reproduction (Cambridge, England), 2009, 138(6): 921-933. doi:10.1530/rep-09-0193.

[7] AHLAWAT S, SHARMA R, MAITRA A. Screening of indigenous goats for prolificacy associated DNA markers of sheep. Gene, 2013, 517(1): 128-131. doi:10.1016/j.gene.2012.12.015.

[8] TEJANGOOKEH H M, SHAHNEH A Z, ZAMIRI M J, DALIRI M, KOHRAM H, JAVAREMI A N. Study of BMP15 gene polymorphism in Iranian goats. African Journal of Biotechnology, 2009, 8(13), 2929-2932.

[9] SUPAKORN C, PRALOMKARN W. Sheep FecB gene polymorphism role in Thai meat goat proliferation rate//Proceedings of 9th World Congress Genetics Applied to Livestock Production, Leipizig Germany, 2010.

[10] HUA G H, CHEN S L, AI J T, YANG L G. None of polymorphism of ovine fecundity major genes FecB and FecX was tested in goat. Animal Reproduction Science, 2008, 108(3/4): 279-286. doi:10.1016/ j.anireprosci.2007.08.013.

[11] HE Y Q, MA X K, LIU X Y, ZHANG C X, LI J. Candidate genes polymorphism and its association to prolificacy in Chinese goats. Journal of Agricultural Science, 2010, 2(1): 88–92. doi:10.5539/jas. v2n1p88.

[12] 李恒, 字向東. 全基因組測序在山羊上的研究進展. 中國畜牧雜志, 2021, 57(10): 29-34. doi:10.19556/j.0258-7033.20200910-01.

LI H, ZI X D. Research progress on whole-genome sequencing on goat. Chinese Journal of Animal Science, 2021, 57(10): 29-34. doi:10.19556/j.0258-7033.20200910-01. (in Chinese)

[13] LAI F N, ZHAI H L, CHENG M, MA J Y, CHENG S F, GE W, ZHANG G L, WANG J J, ZHANG R Q, WANG X, MIN L J, SONG J Z, SHEN W. Whole-genome scanning for the litter size trait associated genes and SNPs under selection in dairy goat (). Scientific Reports, 2016, 6: 38096. doi:10.1038/srep38096.

[14] GUANG-XIN E, ZHAO Y J, HUANG Y F. Selection signatures of litter size in Dazu black goats based on a whole genome sequencing mixed pools strategy. Molecular Biology Reports, 2019, 46(5): 5517-5523. doi:10.1007/s11033-019-04904-6.

[15] ISLAM R, LIU X X, GEBRESELASSIE G, ABIED A, MA Q, MA Y H. Genome-wide association analysis reveals the genetic locus for high reproduction trait in Chinese Arbas Cashmere goat. Genes & Genomics, 2020, 42(8): 893-899. doi:10.1007/s13258-020-00937-5.

[16] E G X, ZHOU D K, YANG B G, DUAN X H, NA R S, HAN Y G, ZENG Y. Association analysis of sixty-seven single nucleotide polymorphisms with litter size in Dazu Black goats. Animal Genetics, 2020, 51(1): 151-152. doi:10.1111/age.12879.

[17] WANG J J, ZHANG T, CHEN Q M, ZHANG R Q, LI L, CHENG S F, SHEN W, LEI C Z. Genomic signatures of selection associated with litter size trait in Jining gray goat. Frontiers in Genetics, 2020, 11: 286. doi:10.3389/fgene.2020.00286.

[18] WANG K, LIU X F, QI T, HUI Y Q, YAN H L, QU L, LAN X Y, PAN C Y. Whole-genome sequencing to identify candidate genes for litter size and to uncover the variant function in goats (). Genomics, 2021, 113(1): 142-150. doi:10.1016/j.ygeno.2020.11.024.

[19] ZI X D, MU X K, LU J Y, MA L, WANG Y. Polymorphisms of growth hormone(GH) and insulin-like growth factor I(IGF-I) genes in prolific Lezhi Black Goat: Possible association with litter size. Journal of Southwest University for Nationalities (Natural Science Edition), 2014, 40(3): 344-349.

[20] Lü M J, LI H, ZI X D. Assessment of estrous synchronization protocols and pregnancy specific protein B concentration for the prediction of kidding rate in Lezhi black goats. Small Ruminant Research, 2021, 195: 106299. doi:10.1016/j.smallrumres.2020.106299.

[21] YANG C X, ZI X D, WANG Y, YANG D Q, MA L, LU J Y, NIU H R, XIAO X. Cloning and mRNA expression levels of GDF9, BMP15, and BMPR1B genes in prolific and non-prolific goat breeds. Molecular Reproduction and Development, 2012, 79(1): 2. doi:10. 1002/mrd.21386.

[22] ZI X D, LU J Y, MA L. Identification and comparative analysis of the ovarian microRNAs of prolific and non-prolific goats during the follicular phase using high-throughput sequencing. Scientific Reports, 2017, 7: 1921. doi:10.1038/s41598-017-02225-x.

[23] ZI X D, LU J Y, ZHOU H, MA L, XIA W, XIONG X R, LAN D L, WU X H. Comparative analysis of ovarian transcriptomes between prolific and non-prolific goat breeds via high-throughput sequencing. Reproduction in Domestic Animals, 2018, 53(2): 344-351. doi:10. 1111/rda.13111.

[24] LI H, DURBIN R. Fast and accurate short read alignment with Burrows-Wheeler transform. Bioinformatics, 2009, 25(14): 1754- 1760. doi:10.1093/bioinformatics/btp324.

[25] BICKHART D M, ROSEN B D, KOREN S, SAYRE B L, HASTIE A R, CHAN S, LEE J, LAM E T, LIACHKO I, SULLIVAN S T, BURTON J N, HUSON H J, NYSTROM J C, KELLEY C M, HUTCHISON J L, ZHOU Y, SUN J J, CRISà A, PONCE DE LEóN F A, SCHWARTZ J C, HAMMOND J A, WALDBIESER G C, SCHROEDER S G, LIU G E, DUNHAM M J, SHENDURE J, SONSTEGARD T S, PHILLIPPY A M, VAN TASSELL C P, SMITH T P L. Single-molecule sequencing and chromatin conformation capture enable de novo reference assembly of the domestic goat genome. Nature Genetics, 2017, 49(4): 643-650. doi:10.1038/ng.3802.

[26] LI H, HANDSAKER B, WYSOKER A, FENNELL T, RUAN J, HOMER N, MARTH G, ABECASIS G, DURBIN R. 1000 genome project data processing subgroup. The sequence alignment/ map format and SAMtools. Microbiology Spectrum, 2009, 25(16): 2078-2079. doi:10.1093/bioinformatics/btp352.

[27] MCKENNA A, HANNA M, BANKS E, SIVACHENKO A, CIBULSKIS K, KERNYTSKY A, GARIMELLA K, ALTSHULER D, GABRIEL S, DALY M, DEPRISTO M A. The Genome Analysis Toolkit: a MapReduce framework for analyzing next-generation DNA sequencing data. Cell Reports, 2010, 20(9): 1297-1303. doi:10.1101/ gr.107524.110.

[28] YANG H, WANG K. Genomic variant annotation and prioritization with ANNOVAR and wANNOVAR. Nature Protocols, 2015, 10(10): 1556-1566. doi:10.1038/nprot.2015.105.

[29] RAUDVERE U, KOLBERG L, KUZMIN I, ARAK T, ADLER P, PETERSON H, VILO J. G: Profiler: a web server for functional enrichment analysis and conversions of gene lists (2019 update). Nucleic Acids Research, 2019, 47(W1): W191-W198. doi:10.1093/ nar/gkz369.

[30] BU D C, LUO H T, HUO P P, WANG Z H, ZHANG S, HE Z H, WU Y, ZHAO L H, LIU J J, GUO J C, FANG S S, CAO W C, YI L, ZHAO Y, KONG L. KOBAS-i: intelligent prioritization and exploratory visualization of biological functions for gene enrichment analysis. Nucleic Acids Research, 2021, 49(W1): W317-W325. doi:10.1093/nar/gkab447.

[31] MOKHTARI M S, ASADI FOZI M, GUTIERREZ J P, NOTTER D R. Genetic and phenotypic aspects of early reproductive performance in Raeini Cashmere goats. Tropical Animal Health and Production, 2019, 51(8): 2175-2180. doi:10.1007/s11250-019-01915-0.

[32] ?URI?I? D, BENI? M, ?AJA I ?, VALPOTI? H, SAMARD?IJA M. Influence of season, rainfall and air temperature on the reproductive efficiency in Romanov sheep in Croatia. International Journal of Biometeorology, 2019, 63(6): 817-824. doi:10.1007/s00484-019- 01696-z.

[33] ASTUTI D A, KHOTIJAH L, MAIDIN M S, NUGROHO P. Reproductive profile of etawah crossbred does fed Flushing diet containing different kinds of plant oil and animal fat. Pakistan Journal of Biological Sciences, 2020, 23(5): 650-657. doi:10.3923/pjbs.2020. 650.657.

[34] MCCARTHY M I, ABECASIS G R, CARDON L R, GOLDSTEIN D B, LITTLE J, IOANNIDIS J P A, HIRSCHHORN J N. Genome-wide association studies for complex traits: consensus, uncertainty and challenges. Nature Reviews Genetics, 2008, 9(5): 356-369. doi:10. 1038/nrg2344.

[35] HONG E P, PARK J W. Sample size and statistical power calculation in genetic association studies. Genomics & Informatics, 2012, 10(2): 117-122. doi:10.5808/gi.2012.10.2.117.

[36] LIU N, CUI W B, CHEN M Y, ZHANG X L, SONG X Y, PAN C Y. A 21-bp indel within thegene is significantly associated with litter size in goat. Animal Biotechnology, 2021, 32(2): 213-218. doi:10.1080/10495398.2019.1677682.

[37] JIANG E H, KANG Z H, WANG X Y, LIU Y, LIU X F, WANG Z, LI X C, LAN X Y. Detection of insertions/deletions (InDels) within the goat Runx2 gene and their association with litter size and growth traits. Animal Biotechnology, 2021, 32(2): 169-177. doi:10.1080/10495398. 2019.1671858.

[38] BALBACH M, FUSHIMI M, HUGGINS D J, STEEGBORN C, MEINKE P T, LEVIN L R, BUCK J. Optimization of lead compounds into on-demand, nonhormonal contraceptives: Leveraging a public- private drug discovery institute collaboration. Biology of Reproduction, 2020, 103(2): 176-182. doi:10.1093/biolre/ioaa052.

[39] CHEN H, CHAN H C. Amplification of FSH signalling by CFTR and nuclear soluble adenylyl cyclase in the ovary. Clinical and Experimental Pharmacology & Physiology, 2017, 44(Suppl 1): 78-85. doi:10.1111/1440-1681.12756.

[40] JAYARAJAN V, APPUKUTTAN A, ASLAM M, REUSCH P, REGITZ-ZAGROSEK V, LADILOV Y. Regulation of AMPK activity by type 10 adenylyl cyclase: Contribution to the mitochondrial biology, cellular redox and energy homeostasis. Cellular and Molecular Life Sciences, 2019, 76(24): 4945-4959. doi:10.1007/s00018-019-03152-y.

[41] WANG C, LI S J, LI C, FENG Y P, PENG X L, GONG Y Z. Molecular cloning, expression profile, polymorphism and the genetic effects of the dopamine D1 receptor gene on duck reproductive traits. Molecular Biology Reports, 2012, 39(9): 9239-9246. doi:10.1007/ s11033-012-1797-3.

[42] LIU Z, YANG N, YAN Y, LI G, LIU A, WU G, SUN C. Genome-wide association analysis of egg production performance in chickens across the whole laying period. BMC Genetics, 2019, 20(1): 67. doi:10.1186/ s12863-019-0771-7.

[43] BARONCHELLI S, VILLA N, REDAELLI S, LISSONI S, SACCHERI F, PANZERI E, CONCONI D, BENTIVEGNA A, CROSTI F, SALA E, BERTOLA F, MAROZZI A, PEDICINI A, VENTRUTO M, POLICE M A, DALPRà L. Investigating the role of X chromosome breakpoints in premature ovarian failure. Molecular Cytogenetics, 2012, 5(1): 32. doi:10.1186/1755-8166-5-32.

[44] FESTA A, UMANO G R, MIRAGLIA DEL GIUDICE E, GRANDONE A. Genetic evaluation of patients with delayed puberty and congenital hypogonadotropic hypogonadism: Is it worthy of consideration? Frontiers in Endocrinology, 2020, 11: 253. doi:10.3389/fendo.2020. 00253.

[45] LI J, LIU J, CAMPANILE G, PLASTOW G, ZHANG C, WANG Z, CASSANDRO M, GASPARRINI B, SALZANO A, HUA G, LIANG A, YANG L. Novel insights into the genetic basis of buffalo reproductive performance. BMC Genomics, 2018, 19(1): 814. doi:10. 1186/s12864-018-5208-6.

[46] NALLASAMY S, KAYA OKUR H S, BHURKE A, DAVILA J, LI Q X, YOUNG S L, TAYLOR R N, BAGCHI M K, BAGCHI I C. Msx homeobox genes act downstream of BMP2 to regulate endometrial decidualization in mice and in humans. Endocrinology, 2019, 160(7): 1631-1644. doi:10.1210/en.2019-00131.

[47] LIU Z K, WANG R C, HAN B C, YANG Y, PENG J P. A novel role of IGFBP7in mouse uterus: regulating uterine receptivity through Th1/Th2 lymphocyte balance and decidualization. PLoS ONE, 2012, 7(9): e45224. doi:10.1371/journal.pone.0045224.

[48] GERHART J, SCHEINFELD V L, MILITO T, PFAUTZ J, NEELY C, FISHER-VANCE D, SUTTER K, CRAWFORD M, KNUDSEN K, GEORGE-WEINSTEIN M. Myo/Nog cell regulation of bone morphogenetic protein signaling in the blastocyst is essential for normal morphogenesis and striated muscle lineage specification. Developmental Biology, 2011, 359(1): 12-25. doi:10.1016/j.ydbio. 2011.08.007.

[49] NYEGAARD M, OVERGAARD M T, SU Y Q, HAMILTON A E, KWINTKIEWICZ J, HSIEH M, NAYAK N R, CONTI M, CONOVER C A, GIUDICE L C. Lack of functional pregnancy- associated plasma protein-A (PAPPA) compromises mouse ovarian steroidogenesis and female Fertility1. Biology of Reproduction, 2010, 82(6): 1129-1138. doi:10.1095/biolreprod.109.079517.

[50] KORDUS R J, HOSSAIN A, CORSO M C, CHAKRABORTY H, WHITMAN-ELIA G F, LAVOIE H A. Cumulus cell pappalysin-1, luteinizing hormone/choriogonadotropin receptor, amphiregulin and hydroxy-delta-5-steroid dehydrogenase, 3 beta- and steroid delta- isomerase 1 mRNA levels associate with oocyte developmental competence and embryo outcomes. Journal of Assisted Reproduction and Genetics, 2019, 36(7): 1457-1469. doi:10.1007/s10815-019- 01489-8.

[51] YU M, WANG J, LIU S, WANG X Q, YAN Q. Novel function of pregnancy-associated plasma protein A: promotes endometrium receptivity by up-regulating N-fucosylation. Scientific Reports, 2017, 7: 5315. doi:10.1038/s41598-017-04735-0.

[52] WON J, DE EVSIKOVA C M, SMITH R S, HICKS W L, EDWARDS M M, LONGO-GUESS C, LI T S, NAGGERT J K, NISHINA P M. NPHP4is necessary for normal photoreceptor ribbon synapse maintenance and outer segment formation, and for sperm development. Human Molecular Genetics, 2010, 20(3): 482-496. doi:10.1093/hmg/ ddq494.

[53] MAO S, WU F, CAO X, HE M, LIU N, WU H, YANG Z, DING Q, WANG X. Tdrp deficiency contributes to low sperm motility and is a potential risk factor for male infertility. American Journal of Translational Research, 2016, 8(1): 177-187.

[54] FARACH M C, TANG J P, DECKER G L, CARSON D D. Heparin/heparan sulfate is involved in attachment and spreading of mouse embryos. Developmental Biology, 1987, 123(2): 401-410. doi:10.1016/0012-1606(87)90398-8.

[55] 潘陽陽, 王萌, 芮弦, 王立斌, 何翃閎, 王靖雷, 馬睿, 徐庚全, 崔燕, 樊江峰, 余四九. IGF-1調控RBM3表達抑制低溫應激誘導牦牛卵丘細胞凋亡. 中國農業科學, 2020, 53(11): 2285-2296.

PAN Y Y, WANG M, RUI X, WANG L B, HE H H, WANG J L, MA R, XU G Q, CUI Y, FAN J F, YU S J. RNA-binding motif protein 3(RBM3) expression is regulated by insulin-like growth factor(IGF-1) for protecting yak() cumulus cells from apoptosis during hypothermia stress. Scientia Agricultura Sinica, 2020, 53(11): 2285-2296. (in Chinese)

[56] THOMAS N, VENKATACHALAPATHY T, ARAVINDAKSHAN T, RAGHAVAN K C. Molecular cloning, SNP detection and association analysis of 5' flanking region of the goat IGF1gene with prolificacy. Animal Reproduction Science, 2016, 167: 8-15. doi:10.1016/j. anireprosci.2016.01.016.

[57] CHENG Y Y, LIU S C, WANG G, WEI W Z, HUANG S, YANG R, GENG H W, LI H Y, SONG J, SUN L D, YU H, HAO L L. Porcine IGF1synonymous mutation alter gene expression and protein binding affinity with IGF1R. International Journal of Biological Macromolecules, 2018, 116: 23-30. doi:10.1016/j.ijbiomac.2018.05.022.

[58] INAGAKI H, OTA S, NISHIZAWA H, MIYAMURA H, NAKAHIRA K, SUZUKI M, NISHIYAMA S, KATO T, YANAGIHARA I, KURAHASHI H. Obstetric complication-associatedpromoter polymorphisms may affect gene expression via DNA secondary structures. Journal of Human Genetics, 2019, 64(5): 459-466. doi:10.1038/s10038-019-0578-4.

[59] ARANDA F, UDRY S, PERéS WINGEYER S, AMSHOFF L C, BOGDANOVA N, WIEACKER P, LATINO J O, MARKOFF A, LARRA?AGA G. Maternal carriers of the ANXA5 M2 haplotype are exposed to a greater risk for placenta-mediated pregnancy complications. Journal of Assisted Reproduction and Genetics, 2018, 35(5): 921-928. doi:10.1007/s10815-018-1142-4.

[60] DRYLLIS G, GIANNOPOULOS A, ZOI C, POULIAKIS A, LOGOTHETIS E, VOULGARELIS M, ZOI K, KOUSKOUNI E, DINOU A, STAVROPOULOS-GIOKAS C, KREATSAS G, KONSTANTOPOULOS K, POLITOU M. Correlation of single nucleotide polymorphisms in the promoter region of the(annexin A5) gene with recurrent miscarriages in women of Greek origin. The Journal of Maternal-Fetal & Neonatal Medicine, 2020, 33(9): 1538-1543. doi:10.1080/14767058.2018.1521799.

[61] DI GERLANDO R, MASTRANGELO S, MOSCARELLI A, TOLONE M, SUTERA A M, PORTOLANO B, SARDINA M T. Genomic structural diversity in local goats: Analysis of copy-number variations. Animals, 2020, 10(6): E1040. doi:10.3390/ani10061040.

[62] HOHOS N M, ELLIOTT E M, GIORNAZI A, SILVA E, RICE J D, SKAZNIK-WIKIEL M E. High-fat diet induces an ovulatory defect associated with dysregulated endothelin-2 in mice. Reproduction (Cambridge, England), 2021, 161(3): 307-317. doi:10.1530/rep-20- 0290.

[63] YANG Y, GUO T, LIU R, KE H, XU W, ZHAO S, QIN Y.gene mutations in premature ovarian insufficiency. Human Mutation, 2020, 41(5): 1033-1041. doi:10.1002/humu.23997.

[64] YANG Y, ZHAO S, QIN Y. Response to “Should FANCL heterozygous pathogenic variants be considered as potentially causative of primary ovarian insufficiency? ”. Human Mutation, 2020, 41(9): 1700-1701. doi:10.1002/humu.24073.

[65] FERGANI C, NAVARRO V M. Expanding the role of tachykinins in the neuroendocrine control of reproduction. Reproduction (Cambridge, England), 2016, 153(1): R1-R14. doi:10.1530/rep-16-0378.

[66] LEóN S, FERGANI C, TALBI R, SIMAVLI S, MAGUIRE C A, GERUTSHANG A, NAVARRO V M. Characterization of the role of NKA in the control of puberty onset and gonadotropin release in the female mouse. Endocrinology, 2019, 160(10): 2453-2463. doi:10. 1210/en.2019-00195.

[67] LEóN S, FERGANI C, TALBI R, MAGUIRE C A, GERUTSHANG A, SEMINARA S B, NAVARRO V M. Tachykinin signaling is required for induction of the preovulatory luteinizing hormone surge and normal luteinizing hormone pulses. Neuroendocrinology, 2021, 111(6): 542-554. doi:10.1159/000509222.

[68] GUO J Z, TAO H X, LI P F, LI L, ZHONG T, WANG L J, MA J Y, CHEN X Y, SONG T Z, ZHANG H P. Whole-genome sequencing reveals selection signatures associated with important traits in six goat breeds. Scientific Reports, 2018, 8: 10405. doi:10.1038/s41598-018- 28719-w.

[69] PASTEL E, POINTUD J C, LOUBEAU G, DANI C, SLIM K, MARTIN G, VOLAT F, SAHUT-BARNOLA I, VAL P, MARTINEZ A, LEFRAN?OIS-MARTINEZ A M. Aldose reductases influence prostaglandin F2α levels and adipocyte differentiation in male mouse and human species. Endocrinology, 2015, 156(5): 1671-1684. doi:10. 1210/en.2014-1750.

[70] ZHANG Q, PEI L G, LIU M, LV F, CHEN G H, WANG H. Reduced testicular steroidogenesis in rat offspring by prenatal nicotine exposure: Epigenetic programming and heritability via nAChR/ HDAC4. Food and Chemical Toxicology, 2020, 135: 111057. doi:10.1016/j.fct.2019.111057.

[71] OLABARRIETA E, TOTORIKAGUENA L, AGIRREGOITIA N, AGIRREGOITIA E. Implication of mu opioid receptor in thematuration of oocytes and its effects on subsequent fertilization and embryo development in mice. Molecular Reproduction and Development, 2019, 86(9): 1236-1244. doi:10.1002/mrd.23248.

[72] OLABARRIETA E, TOTORIKAGUENA L, ROMERO- AGUIRREGOMEZCORTA J, AGIRREGOITIA N, AGIRREGOITIA E. Mu opioid receptor expression and localisation in murine spermatozoa and its role in IVF. Reproduction Fertility and Development, 2020, 32(4): 349-354. doi:10.1071/rd19176.

[73] WATSON L N, MOTTERSHEAD D G, DUNNING K R, ROBKER R L, GILCHRIST R B, RUSSELL D L. Heparan sulfate proteoglycans regulate responses to oocyte paracrine signals in ovarian follicle morphogenesis. Endocrinology, 2012, 153(9): 4544-4555. doi:10. 1210/en.2012-1181.

[74] TIWARI M, PRASAD S, SHRIVASTAV T G, CHAUBE S K. Calcium signaling during meiotic cell cycle regulation and apoptosis in mammalian oocytes. Journal of Cellular Physiology, 2017, 232(5): 976-981. doi:10.1002/jcp.25670.

[75] STEWART T A, DAVIS F M. An element for development: Calcium signaling in mammalian reproduction and development. Biochimica et Biophysica Acta (BBA) - Molecular Cell Research, 2019, 1866(7): 1230-1238. doi:10.1016/j.bbamcr.2019.02.016.

[76] LU N S, LI M J, LEI H L, JIANG X Y, TU W L, LU Y, XIA D. Butyric acid regulates progesterone and estradiol secretion via cAMP signaling pathway in porcine granulosa cells. The Journal of Steroid Biochemistry and Molecular Biology, 2017, 172: 89-97. doi:10. 1016/j.jsbmb.2017.06.004.

[77] JOZKOWIAK M, HUTCHINGS G, JANKOWSKI M, KULCENTY K, MOZDZIAK P, KEMPISTY B, SPACZYNSKI R Z, PIOTROWSKA- KEMPISTY H. The stemness of human ovarian granulosa cells and the role of resveratrol in the differentiation of MSCs-A review based on cellular and molecular knowledge. Cells, 2020, 9(6): E1418. doi:10.3390/cells9061418.

[78] ZHANG T, CHEN L, HAN K P, ZHANG X Q, ZHANG G X, DAI G J, WANG J Y, XIE K Z. Transcriptome analysis of ovary in relatively greater and lesser egg producing Jinghai Yellow Chicken. Animal Reproduction Science, 2019, 208: 106114. doi:10.1016/j.anireprosci. 2019.106114.

[79] CHEN X, SUN X, CHIMBAKA I M, QIN N, XU X, LISWANISO S, XU R, GONZALEZ J M. Transcriptome analysis of ovarian follicles reveals potential pivotal genes associated with increased and decreased rates of chicken egg production. Frontiers in Genetics, 2021, 12: 622751. doi:10.3389/fgene.2021.622751.

[80] XU R Y, PAN L Q, YANG Y Y, ZHOU Y Y. Characterizing transcriptome in female scallopfarreri provides new insights into the molecular mechanisms of reproductive regulation during ovarian development and spawn. Gene, 2020, 758: 144967. doi:10. 1016/j.gene.2020.144967.

[81] HUANG D X, ZHANG B, HAN T, LIU G B, CHEN X, ZHAO Z H, FENG J Q, YANG J W, WANG T M. Genome-wide prediction and comparative transcriptomic analysis reveals the G protein-coupled receptors involved in gonadal development of. Genomics, 2021, 113(1): 967-978. doi:10.1016/j.ygeno.2020.10.030.

[82] NADERI N, HOUSE J D. Recent developments in folate nutrition// Advances in Food and Nutrition Research. Amsterdam: Elsevier, 2018: 195-213. doi:10.1016/bs.afnr.2017.12.006.

[83] BROWN L L, COHEN B E, EDWARDS E, GUSTIN C E, NOREEN Z. Physiological need for calcium, iron, and folic acid for women of various subpopulations during pregnancy and beyond. Journal of Womens Health (Larchmt), 2021, 30(2): 207-211. doi:10.1089/jwh. 2020.8873.

Screening of Key Regulatory Genes for Litter Size Trait Based on Whole Genome Re-Sequencing in Goats ()

LI Heng1, ZI XiangDong1, WANG Hui2, XIONG Yan1, Lü MingJie1,LIU Yu1, JIANG XuDong1

1Key Laboratory of Animal Science of National Ethnic Affairs Commission, Southwest Minzu University, Chengdu 610041;2Key Laboratory of Qinghai-Tibetan Plateau Animal Genetic Resource Reservation and Utilization of Ministry of Education, Southwest Minzu University, Chengdu 610041

【Objective】The purpose of this study was to analyze the genome of different fecundity populations of goats (s) and to explore the key regulatory genes involved in the regulation of litter size traits of Chuanzhong black goats (CBGs), and to provide the theoretical reference for analyzing the genetic mechanism of litter size traits and molecular genetic improvement of fecundity in goats. 【Method】The high fecundity (HF) CBG does (= 6) that produced 4-6 kids per doe kidding and low fecundity (LF) does (= 6) that produced only one kid per doe kidding were chosen in this study. The jugular blood samples were collected to extract genomic DNA. The 350 bp double-terminal sequencing library was constructed, and then 12 whole genome libraries were resequenced by IlluminaHiSeqPE150 platform. The clean data from sequencing were mapped to goat reference genome ARS1 by using BWA software, and two whole-genome scanning analysis methods (and) were used to comprehensively analyze the high-quality SNPs obtained to identify candidate regions. GO analysis and KEGG pathway analysis were performed on the G:Profiler and KOBAS online databases, respectively, to screen candidate genes for regulating the number of kids in CBGs. To further identify the key genetic markers that regulate the number of kids, the synonymous and non-synonymous single nucleotide polymorphisms (SNPs) of reproductive candidate genes were mapped and screened according to the variation analysis report of genome resequencing. The amplified products of 12 goat samples were sequenced by Sanger sequencing to verify the resequencing results.【Result】A total of 431.50 Gb clean data were obtained from the genome resequencing study of 12 CBGs. Through mutation detection and annotation, 7 771 417 SNPs were detected in HF group and 8 935 907 SNPs were detected in LF group, and all types of the LF group SNPs were more than those in HF group. The windows that reach the maximumvalue of top 5% and the minimumvalue of top 5% were set as candidate regions. A total of 130 strong selection signals were annotated in the regions with low heterozygosity and high genetic differentiation, of which 84, 59 and 13 genes were annotated in HF group, LF group and shared window, respectively. GO enrichment analysis and KEGG pathway showed that 19 candidate genes were involved in the regulation of reproduction, reproduction and embryonic development of CBG, including 11 HF group-specific candidate genes (,,,,,,,,,, and), and five strong selection signal genes (,,,, and) in LF group, and three window genes (,and) in HF group shared with LF group. The most GO terms, such as G-protein-coupled receptor activity, hormone response and neuropeptide signal pathway, contained these 19 candidate genes. In addition, nine of the 14 HF candidate genes were significantly enriched in metabolic pathway, neuroactive ligand-receptor interaction, glycosaminoglycan-heparan sulfate/heparin biosynthesis, calcium signal pathway, cAMP signal pathway and folate biosynthesis KEGG pathways (＜0.05). Among the 19 reproductive candidate genes, there were two synonymous mutations (G771T,G529,,andgene mutations could be detected, and this result was consistent with the results of genome resequencing, in whichG529A polymorphism led to alanine mutation to threonine, and【Conclusion】A total of 11 HF group-specific candidate genes were found in this study, which were speculated to be the key regulatory genes for fecundity trait. The mutations ofgene exon G529A andexon A281T might be the key genetic markers for regulating prolificacy traits in goats, which had great application value in improving reproductive performance of goats.

Chuanzhong black goat; genome resequencing; fecundity; candidate genes

10.3864/j.issn.0578-1752.2022.23.015

2021-08-17；

2022-10-12

國家自然科學基金（31902154）、西南民族大學中央高校基本科研業務專項（2021PTJS26）

李恒，E-mail：lih199501 @sina.com。通信作者字向東，E-mail：zixd@sina.com

（責任編輯 林鑒非）