雞繁殖性能近交衰退相關(guān)CpG島差異甲基化基因的篩選

薛 倩，李國輝，殷建玫，張會永，周成浩，朱云芬，邢偉杰，蘇一軍，鄒劍敏，韓 威*

(1.江蘇省家禽科學研究所，揚州 225125; 2.江蘇省家禽科學研究所 科技創(chuàng)新中心，揚州 225125)

由于近交而引起的后代群體適合度下降的現(xiàn)象被稱為近交衰退[1]，近交衰退是物種資源保護尤其是瀕危物種保護過程中備受關(guān)注的問題之一。群體遺傳學提出了兩個解釋近交衰退遺傳基礎(chǔ)的經(jīng)典假說: 顯性假說和超顯性假說。顯性假說認為，近交過程中，基因組純合性的增加導致了有害隱性等位基因的暴露；超顯性假說提出，在整個基因組中，雜合子一般優(yōu)于純合子，近交增加了純合子出現(xiàn)的概率，因此會表現(xiàn)出衰退[2]。近些年研究發(fā)現(xiàn)，這兩種假說對于某些近交衰退現(xiàn)象卻無法完全解釋[2-6]，比如近交衰退程度會隨著環(huán)境壓力的變化而變化的不確定現(xiàn)象等[7-8]。DNA表觀修飾也會受到環(huán)境的影響，進而影響到基因的表達和表型的可塑性[9-11]。并且已有研究表明，近交衰退與DNA表觀遺傳修飾密切相關(guān)[12-14]。表觀遺傳變化是發(fā)生在基因組DNA中的重要變異，它在不直接改變DNA編碼序列的情況下影響著基因的表達[15-16]。

DNA甲基化是目前研究最廣泛的表觀遺傳修飾方式之一[17]，其參與多種生物學過程[18-19]，在調(diào)控基因表達和維持基因組穩(wěn)定性方面起著關(guān)鍵的作用[20]。脊椎動物基因組DNA甲基化主要發(fā)生在 CpG胞嘧啶上，并且它的調(diào)控多為通過 CpG 島(CpG island)多個 CpG 位點的甲基化共同作用來完成，而不是通過對某單個 CpG 位點的甲基化或去甲基化來激活或失活基因[21]。與CpG 島鄰近(～2 kb) 的低CpG 密度區(qū)域被定義為CpG 島岸(CpG island shore)，這些區(qū)域的甲基化與轉(zhuǎn)錄的激活也密切相關(guān)。匡穎等[22]研究表明，近交系小鼠生長發(fā)育緩慢可能與H19基因上游區(qū)域被甲基化進而引起該基因的異常表達有關(guān)；蔣曹德[23]以梅山豬、大白豬群體為親本構(gòu)建了兩個自繁群體與兩個雜交群體，發(fā)現(xiàn)群體內(nèi)自繁的子代甲基化水平略高于雜交子代群體，且兩個群體間某些性狀差異可能與某特定位點甲基化差異相關(guān)[23]；Venney等[24]研究表明，三文魚的近交衰退與其CK-1、GTIIBS和hsp70三個特定基因的甲基化變化有關(guān)；楊晶淼等[25]研究報告了馬氏珠母貝的近交與雜交家系的生長性狀和甲基化水平均存在顯著性差異。植物上，關(guān)于山蘿卜[14]和茄科植物[26]的兩項研究均表明，DNA甲基化與近交衰退之間存在著密切的關(guān)系，并且近交個體比遠交個體存在更為廣泛的全基因組DNA甲基化現(xiàn)象。盡管有這些研究，但是人們對DNA甲基化，尤其是脊椎動物CpG 島及其附近區(qū)域甲基化在近交衰退中的調(diào)控作用還沒有完全深刻地理解。

家禽繁殖性能在家禽生產(chǎn)、育種和品種保護中具有極其重要的地位，尤其在地方雞遺傳資源活體保護過程中，群體繁殖性能近交衰退現(xiàn)象時有發(fā)生，嚴重影響著保種群的規(guī)模和世代的延續(xù)。本課題組在前期的研究中，已成功組建了狼山雞高近交組(分子近交系數(shù)FIS> 0.15)和低近交組(分子近交系數(shù)FIS< 0.04)，且觀察到狼山雞高近交組繁殖性能發(fā)生了顯著衰退現(xiàn)象[27]，然而，關(guān)于繁殖性能近交衰退個體DNA甲基化變化，尤其是CpG 島及其附近重要調(diào)控區(qū)域甲基化變化情況尚未進行深入探索。基于此，本研究利用新一代測序技術(shù)—全基因組重亞硫酸鹽測序(whole-genome bisulfite sequencing，WGBS)技術(shù)檢測和分析狼山雞高、低近交組性腺軸組織(包括卵巢和下丘腦)全基因組甲基化水平變化情況，著重篩選CpG 島及其附近區(qū)域內(nèi)存在的差異甲基化區(qū)域(differential methylation region，DMR)，并探討其在狼山雞繁殖性能近交衰退中可能發(fā)揮的調(diào)控作用，研究結(jié)果將為今后家禽育種和物種資源保護工作提供理論參考依據(jù)。

1 材料與方法

1.1 試驗動物和樣品采集

選取江蘇省家禽科學研究所國家級地方雞種基因庫保存的地方雞品種—狼山雞為試驗素材，基于前期狼山雞高、低近交組母雞繁殖性狀(包括開產(chǎn)日齡、300天產(chǎn)蛋數(shù)、開產(chǎn)體重和開產(chǎn)蛋重)的表型值記錄，選取高近交組中繁殖性能顯著衰退的個體3只 和低近交組中繁殖性能處于群體均值附近的個體3只，快速屠宰并采集其下丘腦和卵巢組織，迅速投入液氮中，然后轉(zhuǎn)移至-80 ℃冰箱中保存?zhèn)溆谩?/p>

1.2 DNA提取和WGBS建庫測序

使用德國凱杰公司提供的組織DNA快速提取試劑盒提取各樣本基因組DNA,并使用安捷倫2100分光光度計進行檢測。DNA樣品經(jīng)超聲片段化后進行亞硫酸氫鹽轉(zhuǎn)化，使用DNA甲基化測序文庫構(gòu)建試劑盒(斯威夫特生物科技，美國)，將單鏈DNA片段進行接頭連接，對DNA片段進行8個循環(huán)PCR擴增：95 ℃預變性2 min；8個循環(huán)(95 ℃變性30 s，65 ℃退火20 s，72 ℃延伸45 s)；最后72 ℃ 延伸性7 min。擴增反應體系包括預混反應體系25 μL (低濃度EDTA TE 10 μL, Buffer R1 10 μL, 反應物R2 4 μL，R3酶1 μL)，DNA純化產(chǎn)物和探針引物25 μL，共計50 μL。對擴增產(chǎn)物進行純化和完整性檢測，最后使用Illumina HiSeq 4000 測序平臺對測序文庫進行雙末端測序，每個文庫測序深度達30×。

1.3 數(shù)據(jù)過濾質(zhì)控和甲基化胞嘧啶(mC)位點識別

使用Cutadapt[28]軟件去除測序Reads中含有測序接頭、污染、低質(zhì)量堿基和未知堿基的Reads，序列質(zhì)量進一步使用FastQC (http://www.bioinformatics.babraham.ac.uk/projects/fastqc/)進行驗證。通過質(zhì)控的Reads被比對到雞的參考基因組中，使用Samtool[29]軟件對比對上的Reads進一步去重復。對于每個胞嘧啶位點，其甲基化水平取決于覆蓋到mC位點的Reads數(shù)與總Reads數(shù)的比值。使用二項分布檢驗鑒別錯誤檢出率(false-discovery rate，F(xiàn)DR)≤0.01的mC位點，所有mC位點測序覆蓋度需≥3。根據(jù)mC位點上、下游堿基種類，將mC位點分為3種類型：mCG、mCHG和mCHH，其中H代表腺苷(A)、胞嘧啶(C)或胸腺嘧啶(T)。

1.4 DMR檢測與CpG島區(qū)DMR相關(guān)基因功能富集分析

使用R語言包MethylKit[30]進行DMR檢測，采用滑窗方法(1 000 bp窗口，500 bp步長)掃描基因組，檢測高、低近交組中存在差異甲基化的區(qū)域，對每個窗口進行Fisher精確檢驗，并使用FDR對多次檢驗的P值進行校正。FDR<0.05且甲基化水平差異倍數(shù)至少為2倍的窗口被鑒定作為DMR。采用emboss[31-32]進行CpG 島鑒定，將CpG島定義為長度至少100 bp、GC含量大于50%，且CpG二核苷酸的出現(xiàn)率(觀測值/期望值)大于60%的基因組區(qū)域。對CpG 島區(qū)(CpG 島及其上、下游2 kb區(qū)域)甲基化信息進行統(tǒng)計分析，篩選位于CpG 島區(qū)的DMR。為探索這些DMR的功能，對包含這些DMR的基因(稱為DMR相關(guān)基因或差異甲基化基因)進行GO和KEGG數(shù)據(jù)庫功能注釋和富集分析，P<0.05的GO條目和Pathway信號通路為顯著富集的條目和通路。

1.5 亞硫酸氫鹽測序PCR(bisulfite sequencing PCR,BSP)方法驗證WGBS分析結(jié)果

為驗證WGBS測序分析結(jié)果的可靠性，從卵巢中隨機選取2個DMR相關(guān)基因(CDC27、SRD5A1)，對其CpG島區(qū)域甲基化水平進行驗證。利用Methylation Primer Express V1.0軟件設(shè)計擴增引物，引物信息見表1。選用WGBS相同的DNA樣品，每個樣品取1 mg，用EZ DNA Methylation-GoldTM試劑盒進行亞硫酸氫鹽轉(zhuǎn)化。利用設(shè)計的引物擴增轉(zhuǎn)化后的DNA目標區(qū)域，PCR產(chǎn)物純化后進行TA載體克隆，每個基因挑選10～15個克隆子進行Sanger測序，對測序結(jié)果采用BiQ Analyzer軟件分析。

表1 BSP驗證試驗引物信息

1.6 數(shù)據(jù)處理分析

采用IBM SPSS Statistics 20.0軟件統(tǒng)計分析高、低近交組各繁殖性狀，通過獨立樣本t(two-sample t-test)檢驗方法估計兩組間各性狀差異顯著性，當P<0.05時差異顯著，當P<0.01時差異極顯著；根據(jù)測序Reads中覆蓋到mC位點的Reads數(shù)與總Reads數(shù)的比值，計算全基因組甲基化水平；根據(jù)clean reads中唯一比對的reads，運用R語言psych包對各樣本進行主成分(principle components analysis, PCA)分析。

2 結(jié) 果

2.1 狼山雞高、低近交組繁殖性能差異

比較分析狼山雞高、低近交組間各繁殖性狀的差異，發(fā)現(xiàn)高、低近交組在開產(chǎn)日齡、300天產(chǎn)蛋數(shù)和開產(chǎn)體重3個性狀上均存在顯著或極顯著差異(P<0.05或P<0.01)(表2)，而在開產(chǎn)蛋重上差異不顯著(P>0.05)。

2.2 WGBS數(shù)據(jù)概述以及狼山雞全基因組DNA甲基化水平分析

測序樣本總共12個，包含4個組，每組3個生物學重復，每個樣本測序深度≥30×，產(chǎn)生約154.9 Mb clean reads，平均69.14% 的clean reads可以唯一比對到參考基因組中，這些clean reads用于后續(xù)的分析。為證明樣本采集的準確性和可靠性，對其進行主成分分析(PCA)，結(jié)果表明，來自同一組的3個生物學重復相距較近(圖1)，表明試驗樣本采集以及測序數(shù)據(jù)較為準確、可靠。

根據(jù)測序獲得reads中覆蓋到mC位點的reads數(shù)以及reads總數(shù)，計算并分析了狼山雞高近交組和低近交組中卵巢和下丘腦全基因組DNA甲基化水平，結(jié)果發(fā)現(xiàn)，狼山雞高近交組卵巢和下丘腦全基因組甲基化水平分別為3.49%和4.31%，低近交組卵巢和下丘腦全基因組甲基化水平分別為3.48%和4.35%。相同組織，狼山雞高近交組與低近交組全基因組甲基化水平差異不顯著(P>0.05)，然而兩種組織間比較發(fā)現(xiàn)，下丘腦基因組整體甲基化水平顯著高于卵巢組織(P<0.05，圖2)。

表2 高、低近交組狼山雞各繁殖性狀差異

2.3 狼山雞高、低近交組間差異甲基化區(qū)域檢測

差異甲基化區(qū)域(DMR)是表觀遺傳變化的一個重要標志，涉及到許多生物學過程調(diào)控。為探索與狼山雞繁殖性能近交衰退相關(guān)的差異甲基化區(qū)域, 本研究篩選高近交組和低近交組間甲基化水平差異大于等于2倍，且 FDR<0.05的區(qū)域作為差異甲基化區(qū)域。結(jié)果在兩組狼山雞的下丘腦和卵巢中分別檢測到5 948和4 593個差異甲基化區(qū)域，其中分別有1 798和995個差異甲基化區(qū)域位于CpG 島區(qū)。下丘腦CpG 島區(qū)檢測到的差異甲基化區(qū)域中，有809個上調(diào)，989個下調(diào)；卵巢CpG 島區(qū)的995個差異甲基化區(qū)域中，有374個上調(diào)，621個下調(diào)(圖3)。

脊椎動物CpG 島區(qū)DNA甲基化在基因轉(zhuǎn)錄、沉默等過程中起重要調(diào)控作用，因此，對下丘腦和卵巢基因組CpG 島區(qū)差異甲基化區(qū)域進行注釋分析發(fā)現(xiàn)，下丘腦和卵巢基因組CpG 島區(qū)差異甲基化區(qū)域分別注釋到了1 020 和 552個基因，這些基因稱為DMR相關(guān)基因或差異甲基化基因。

2.4 下丘腦基因組CpG島區(qū)差異甲基化基因功能注釋分析

對下丘腦基因組CpG 島區(qū)差異甲基化基因進行GO和KEGG注釋和功能富集分析發(fā)現(xiàn)，分別有920和686個基因被注釋到了GO和KEGG數(shù)據(jù)庫。GO富集分析結(jié)果顯示，有141個GO條目被這些基因顯著富集(P<0.05)，這些GO條目主要涉及到信號轉(zhuǎn)導、神經(jīng)系統(tǒng)發(fā)育、轉(zhuǎn)錄調(diào)控、蛋白激酶活性等生物學過程，另外，卵母細胞成熟調(diào)控以及生殖系統(tǒng)發(fā)育相關(guān)的GO條目也被顯著富集到。KEGG富集分析表明，這些差異甲基化基因顯著富集于19條KEGG信號通路中(P<0.05)，包括轉(zhuǎn)化生長因子β信號通路、乙型肝炎、脂肪酸代謝、胰島素信號通路、調(diào)控多能性干細胞信號通路、糖酵解/糖異生、氨基酸代謝等信號通路(圖4)，其中Smad7a、SMAD6、BMP2、TGFB2、MAPK6等繁殖相關(guān)基因包含在這些通路中。

2.5 卵巢基因組CpG島區(qū)差異甲基化基因功能注釋分析

對卵巢基因組CpG 島區(qū)差異甲基化基因進行GO和KEGG注釋和功能富集分析發(fā)現(xiàn)，分別有500和369個基因被注釋到了GO和KEGG數(shù)據(jù)庫。GO富集分析結(jié)果表明，有105個GO條目被這些基因顯著富集(P<0.05)，其中涉及到神經(jīng)系統(tǒng)發(fā)育、細胞增殖、離子轉(zhuǎn)運、轉(zhuǎn)錄調(diào)控等(P<0.05)。KEGG富集分析結(jié)果顯示，12條信號通路被顯著富集(P<0.05)，包括慢性骨髓白血病、流感A、精氨酸和脯氨酸代謝、黏著連接、調(diào)控多能性干細胞信號通路等(圖5)。另外，還發(fā)現(xiàn)一些與卵子發(fā)育和性激素分泌相關(guān)的通路被富集到，如黃體酮介導的卵母細胞成熟、卵母細胞減數(shù)分裂、GnRH信號通路、雌激素信號通路等，其中包含CDC27、ADCY8、AKT3等10個差異甲基化基因(表3)。

2.6 差異甲基化區(qū)域BSP驗證

為驗證全基因組重亞硫酸鹽測序分析結(jié)果的可靠性，隨機挑選2個CpG島區(qū)差異甲基化區(qū)域(Chr27:2838112-2838584和Chr2:79746219-79746715)，這兩個差異甲基化區(qū)域分別位于CDC27和SRD5A1基因的CpG島中。采用BSP—克隆測序法對這兩個區(qū)域甲基化情況進行檢測，結(jié)果顯示，在SRD5A1和CDC27基因的驗證區(qū)域分別檢測到了6和14個甲基化CpG位點(圖6B)，且這些位點的平均甲基化水平在高近交組與低近交組間的變化趨勢與全基因組重亞硫酸鹽測序結(jié)果(圖6A)相一致，表明了全基因組重亞硫酸鹽測序分析結(jié)果的準確性和可靠性。

表3 卵巢中富集的繁殖相關(guān)通路和CpG 島區(qū)差異甲基化基因

3 討 論

DNA甲基化尤其是脊椎動物CpG島區(qū)的甲基化在調(diào)控基因轉(zhuǎn)錄以及轉(zhuǎn)錄后生物學過程中發(fā)揮重要作用[21]，大量研究表明，DNA甲基化與動植物的近交衰退現(xiàn)象密切相關(guān)[10-12]。本研究基于課題組前期狼山雞高、低近交組的成功組建，以及狼山雞高近交組個體出現(xiàn)的繁殖性能近交衰退現(xiàn)象，通過全基因組重亞硫酸鹽測序技術(shù)開展狼山雞高近交組與低近交組間性腺軸組織DNA甲基化差異比較分析，試圖揭示DNA甲基化在狼山雞繁殖性能近交衰退中的調(diào)控作用。

通過比較分析狼山雞高、低近交組卵巢和下丘腦基因組整體甲基化水平差異，發(fā)現(xiàn)高近交組和低近交組間同一組織基因組整體甲基化水平差異不顯著(P>0.05)，然而，將基因組劃分成許多個小的區(qū)域，分析各區(qū)域高、低近交組間的甲基化水平差異時，卻發(fā)現(xiàn)了大量的上調(diào)或下調(diào)的差異甲基化區(qū)域，可能是由于這些區(qū)域上、下調(diào)趨勢正好相反，導致整個基因組甲基化水平在高、低近交組間差異不顯著。差異甲基化區(qū)域是表觀遺傳變化的一個重要標志，被看作為參與基因轉(zhuǎn)錄水平調(diào)控的功能性區(qū)域，尤其是脊椎動物CpG島區(qū)DNA甲基化的變化一直以來受到研究者們的極大關(guān)注[33-36]。本研究在下丘腦和卵巢基因組CpG島區(qū)分別檢測到了1 798和995個差異甲基化區(qū)域，這些差異甲基化區(qū)域分別注釋到了1 020 和 552個基因，為了進一步探索這些差異甲基化區(qū)域的功能，對差異甲基化區(qū)域相關(guān)基因進行了功能注釋和富集分析，發(fā)現(xiàn)了大量差異甲基化基因富集在神經(jīng)系統(tǒng)發(fā)育、信號轉(zhuǎn)導和離子轉(zhuǎn)運等神經(jīng)內(nèi)分泌調(diào)節(jié)相關(guān)生物學過程GO條目，這些生物學過程與雞的繁殖性能密切相關(guān)。KEGG功能富集分析顯示，卵巢中一些差異甲基化基因富集在了卵子發(fā)育和性激素分泌相關(guān)的通路，如黃體酮介導的卵母細胞成熟、卵母細胞減數(shù)分裂、GnRH信號通路、雌激素信號通路等，包含CDC27、ADCY8、AKT3等10個差異甲基化基因，推測這些基因CpG島區(qū)的差異甲基化在狼山雞繁殖性能近交衰退中可能發(fā)揮重要調(diào)控作用[37]。其中，CDC27基因CpG島區(qū)甲基化差異情況進一步得到了BSP-克隆測序方法的驗證，CDC27是編碼細胞周期進程所必需的一種蛋白質(zhì)[38-39]，被認為與細胞增殖和細胞分裂相關(guān)[40]，本研究中，該基因CpG島區(qū)甲基化水平在狼山雞高近交組中顯著上調(diào)，推測該基因CpG島區(qū)甲基化上調(diào)影響了其轉(zhuǎn)錄表達，進而影響了卵母細胞分裂和成熟，導致狼山雞高近交組繁殖性能衰退。

另外，在下丘腦和卵巢差異甲基化基因顯著富集的通路中，發(fā)現(xiàn)多條通路與氨基酸、脂肪和糖代謝相關(guān)，推測狼山雞的高度近交可能影響了其基礎(chǔ)代謝，進而可能間接影響到其生殖系統(tǒng)發(fā)育，導致其繁殖性能下降。一些與疾病相關(guān)的信號通路也被下丘腦和卵巢差異甲基化基因顯著富集到，如乙型肝炎、慢性骨髓白血病、流感A等，表明狼山雞高度近交后，其抗病力或免疫力方面也受到了影響，可能需要后續(xù)進一步關(guān)注和研究。

4 結(jié) 論

狼山雞高、低近交組間基因組整體DNA甲基化水平差異不顯著(P>0.05)，但全基因組范圍內(nèi)仍檢測到大量差異甲基化區(qū)域，且發(fā)現(xiàn)大量差異甲基化基因富集到繁殖相關(guān)的生物學過程和通路中，如CDC27、ADCY8、AKT3等差異甲基化基因，推測這些基因CpG島區(qū)甲基化在狼山雞繁殖性能近交衰退過程中發(fā)揮重要作用。本研究結(jié)果為今后深入研究雞繁殖性能近交衰退調(diào)控機制提供了有力線索，為物種資源保護工作提供了科學依據(jù)。