李 穎 宋 峰 蘇德其其格 AKONYANI Zaccheaus Pazamilala 烏蘭吐雅 吳江鴻

(內(nèi)蒙古民族大學(xué) 動(dòng)物科技學(xué)院，內(nèi)蒙古 通遼 028000)

內(nèi)蒙古絨山羊是一種兼具絨肉兼用性的地方優(yōu)良品種[1]，毛被潔白有光澤，是我國(guó)重要的創(chuàng)匯動(dòng)物[2]，是經(jīng)過(guò)長(zhǎng)期人工選育而成可適應(yīng)飼養(yǎng)地區(qū)寒冷、干旱等氣候環(huán)境的優(yōu)良選育能適應(yīng)惡劣環(huán)境變化的絨山羊品種、提高絨山羊優(yōu)良產(chǎn)絨性狀以及環(huán)境與基因互作已成為研究熱點(diǎn)。既使基因組序列已經(jīng)公布，但絨山羊參考基因組收錄的數(shù)據(jù)其特征和mRNA 結(jié)構(gòu)并沒(méi)有得到深入分析且基因信息收錄不全，因此不能滿足現(xiàn)階段絨山羊育種及優(yōu)良性狀挖掘等相關(guān)研究。伴隨著高通量測(cè)序技術(shù)的迅猛發(fā)展，轉(zhuǎn)錄組測(cè)序(RNA-Seq)技術(shù)在各個(gè)研究領(lǐng)域得到廣泛的推廣與運(yùn)用，成為研究基因表達(dá)調(diào)控的主要方法[3]。RNA-Seq可用于發(fā)現(xiàn)低豐度轉(zhuǎn)錄本和新轉(zhuǎn)錄本，并識(shí)別不同樣本之間轉(zhuǎn)錄本的差異表達(dá)[4-5]。第二代高通量測(cè)序技術(shù)可用于篩選差異表達(dá)基因[6]，第三代測(cè)序技術(shù)可以不間斷地直接讀取全長(zhǎng)轉(zhuǎn)錄本，可以有效地獲得單個(gè)RNA分子的完整高質(zhì)量序列，并準(zhǔn)確識(shí)別下一代測(cè)序(NGS)無(wú)法識(shí)別的差異表達(dá)異構(gòu)體、同源基因、超家族基因和等位基因的轉(zhuǎn)錄本[7]。它可以改善基因表達(dá)的實(shí)時(shí)定量聚合酶鏈反應(yīng)(RT-PCR)結(jié)果，識(shí)別選擇性剪接和基因融合現(xiàn)象，發(fā)現(xiàn)新基因和轉(zhuǎn)錄異構(gòu)體[8]。為能夠更加全面得到絨山羊皮膚中的基因及轉(zhuǎn)錄本，本研究以內(nèi)蒙古白絨山羊?yàn)樵囼?yàn)動(dòng)物，利用 SMRT 技術(shù)對(duì)絨山羊皮膚進(jìn)行全長(zhǎng)轉(zhuǎn)錄組測(cè)序及生物信息學(xué)分析，通過(guò)與絨山羊參考基因組進(jìn)行數(shù)據(jù)比對(duì)，鑒定新基因及轉(zhuǎn)錄本，分析絨山羊皮膚基因的可變剪接(Alternative splicing,AS)類型并統(tǒng)計(jì)不同類型的數(shù)量，預(yù)測(cè)絨山羊皮膚中 lncRNA 及簡(jiǎn)單重復(fù)序列(Simple sequence repeats,SSR)，為有效提高絨山羊基因組注釋水平，豐富基因序列及結(jié)構(gòu)信息，深入挖掘絨山羊品種資源，優(yōu)良性狀的提升及家畜基因組與環(huán)境互作機(jī)理研究提供理論依據(jù)。

1 材料與方法

1.1 試驗(yàn)材料

為更加全面的得到絨山羊皮膚基因的所有轉(zhuǎn)錄本，本研究在呼和浩特市土默特左旗分4個(gè)時(shí)間點(diǎn)(1、4、7和10月)采集了3只內(nèi)蒙古白絨山羊的體側(cè)皮膚和耳部皮膚樣品，使用Trizol分別從皮膚樣本中提取RNA，并將樣品RNA混合構(gòu)建文庫(kù)后進(jìn)行測(cè)序。

1.2 試驗(yàn)方法

1.2.1絨山羊皮膚樣本RNA提取

每100 mg新鮮皮膚組織樣品中加入1 mL Trizol試劑，使用無(wú)菌手術(shù)刀在冰上將樣品切碎，并用無(wú)菌勻漿器勻漿，將組織或細(xì)胞裂解液轉(zhuǎn)移到1.5 mL 無(wú)RNA酶EP管中。置冰上，靜置5 min。每管加入200 μL氯仿，充分混合后冰上靜置10 min，使核蛋白復(fù)合物完全解離。4 ℃ 13 000 r/min離心15 min。期間，取一新EP管，加入500 μL異丙醇，冰上預(yù)冷。離心后，將上層水相 (約500 μL)轉(zhuǎn)移到該新的EP管中。冰上靜置，酒精沉淀10 min，13 000 r/min 離心10 min，去除上清液。用1 mL 75%乙醇洗滌RNA沉淀1次，12 000 r/min離心5 min。去掉上清，真空干燥RNA沉淀5～10 min。將RNA溶解在30～50 μL DEPC處理過(guò)的去離子水中，最后進(jìn)行分光光度分析以確定樣品濃度和純度。

1.2.2絨山羊皮膚全長(zhǎng)轉(zhuǎn)錄組建庫(kù)及SMRT測(cè)序

在文庫(kù)構(gòu)建之前，從每個(gè)皮膚樣本中收集等量的總RNA，根據(jù)Iso-Seq協(xié)議構(gòu)建SMART測(cè)序文庫(kù)。使用SmarterTMPCR cDNA synthesis kit合成mRNA的全長(zhǎng)cDNA。PacBio Iso-Seq文庫(kù)按照PacBio標(biāo)準(zhǔn)文庫(kù)制備方案制備，在PacBio RsII儀器上進(jìn)行測(cè)序，PacBio測(cè)序由中國(guó)北京百邁客測(cè)序公司進(jìn)行。

1.2.3絨山羊皮膚全長(zhǎng)轉(zhuǎn)錄組數(shù)據(jù)處理

使用Iso-seq管道將原始數(shù)據(jù)進(jìn)行處理，運(yùn)行參數(shù)為minFullPass=3和minPredictedAccuracy=0.9。然后通過(guò)poly A尾部信號(hào)和5'及3' cDNA引物來(lái)確定全長(zhǎng)非嵌合的轉(zhuǎn)錄本。ICE(Iterative clustering for error correction)被用來(lái)獲得共識(shí)異構(gòu)體，使用Quiver對(duì)來(lái)自ICE的FL共識(shí)序列進(jìn)行拋光，高質(zhì)量的FL轉(zhuǎn)錄本被分類，其標(biāo)準(zhǔn)是校正后的準(zhǔn)確率高于99%。使用GMAP將FL共識(shí)序列映射到參考基因組上，用Pbtranscript-tofu軟件包對(duì)映射的讀數(shù)進(jìn)行進(jìn)一步折疊，最小覆蓋率(Min-coverage)為85%，最小一致性(Min-identity)為90%，在折疊多余的轉(zhuǎn)錄本時(shí)不考慮5′差異。

1.2.4絨山羊皮膚全長(zhǎng)轉(zhuǎn)錄組數(shù)據(jù)分析

新轉(zhuǎn)錄本和新基因功能注釋使用BLAST[9]軟件與NR、Swissprot、GO、COG、KOG、eggNOG、Pfam、KEGG數(shù)據(jù)庫(kù)進(jìn)行比對(duì)，獲取注釋信息；使用TransDecoder[10]軟件對(duì)其編碼區(qū)序列及對(duì)應(yīng)氨基酸序列進(jìn)行預(yù)測(cè)；通過(guò)Astalavista[11]軟件獲取每個(gè)樣品間存在的可變剪接類型；使用Cytoscape[12]插件ClueGO對(duì)可變剪接基因進(jìn)行功能富集分析；利用MEME對(duì)所有轉(zhuǎn)錄本poly A位點(diǎn)上游50 bp的序列進(jìn)行分析；從新轉(zhuǎn)錄本中篩選500 bp以上的轉(zhuǎn)錄本，利用MISA軟件做SSR分析；對(duì)新發(fā)現(xiàn)的轉(zhuǎn)錄本進(jìn)行l(wèi)ncRNA的預(yù)測(cè)，主要包括：CPC分析、CNCI[13]分析、pfam蛋白結(jié)構(gòu)域分析和CPAT[14]分析4種方法。

2 結(jié)果與分析

2.1 絨山羊皮膚SMRT測(cè)序數(shù)據(jù)結(jié)果統(tǒng)計(jì)

采用SMRT測(cè)序技術(shù)，完成樣品的全長(zhǎng)轉(zhuǎn)錄組測(cè)序，使用2個(gè)Cell獲得 42.41 Gb clean data的數(shù)據(jù)。測(cè)序獲得CCS 591 585 條，其中全長(zhǎng)非嵌合序列466 058條。對(duì)全長(zhǎng)非嵌合序列進(jìn)行聚類得到一致序列149 604條，對(duì)一致序列進(jìn)行Polish得到高質(zhì)量一致序列共137 211條，用二代轉(zhuǎn)錄組數(shù)據(jù)對(duì)低質(zhì)量一致序列進(jìn)行校正，校正后與高質(zhì)量一致序列合并進(jìn)行去冗余得到82 382條轉(zhuǎn)錄本序列。對(duì)去冗余后的82 382個(gè)轉(zhuǎn)錄本進(jìn)行可變剪接分析，檢測(cè)到基因18 919個(gè)，其中新基因6 166個(gè)，新發(fā)現(xiàn)轉(zhuǎn)錄本55 875 個(gè)。對(duì)新發(fā)現(xiàn)的轉(zhuǎn)錄本進(jìn)行序列結(jié)構(gòu)分析，共預(yù)測(cè)出 66 951個(gè)SSR、39 346個(gè)完整ORF序列和10 927個(gè)lncRNA，并完成了49 573個(gè)新轉(zhuǎn)錄本的功能注釋。

2.2 絨山羊皮膚全長(zhǎng)轉(zhuǎn)錄組測(cè)序完整性評(píng)估

利用BUSCO[15]對(duì)去冗余后的轉(zhuǎn)錄組進(jìn)行完整性評(píng)估，BUSCO使用OrthoDB作為參考數(shù)據(jù)庫(kù)，構(gòu)建了多個(gè)進(jìn)化分支的單拷貝基因集。評(píng)估結(jié)果如圖1所示，完整覆蓋基因數(shù)為3 145個(gè)，其中單拷貝基因數(shù)1 260個(gè)，多拷貝基因數(shù)1 885個(gè)，占總基因數(shù)的79.62%，測(cè)序獲得的絨山羊皮膚全長(zhǎng)轉(zhuǎn)錄組數(shù)據(jù)的完整性及準(zhǔn)確性是非?？煽康?。

2.3 絨山羊皮膚新轉(zhuǎn)錄本及新基因功能注釋

將得到的新轉(zhuǎn)錄本和新基因與NR、Swiss-prot、GO、COG、KOG、Pfam和KEGG等數(shù)據(jù)庫(kù)比對(duì)進(jìn)行功能注釋，各數(shù)據(jù)庫(kù)注釋到的轉(zhuǎn)錄本數(shù)量統(tǒng)計(jì)見(jiàn)表1。由表1可知，共對(duì)49 573個(gè)轉(zhuǎn)錄本進(jìn)行注釋，其中通過(guò)NR數(shù)據(jù)庫(kù)成功注釋的轉(zhuǎn)錄本數(shù)量最多，共有49 374條，在總注釋量中占比99.60%。新轉(zhuǎn)錄本NR數(shù)據(jù)庫(kù)注釋分布前幾位的分別是：山羊(Caprahircus，11 873，24.06%)、牛(Bostaurus，9 918，20.10%)、野牦牛(Bosmutus，6 109，12.38%)、綿羊(Ovisaries，4 737，9.60%)和藏羚羊(Pantholopshodgsonii，3 331，6.75%)(圖2(a))。新基因NR數(shù)據(jù)庫(kù)注釋分布前幾位的分別是：牛(Bostaurus，1 787，42.65%)、野牦牛(Bosmutus，1 318，31.46%)、山羊(Caprahircus，182，4.34%)、綿羊(Ovisaries，138，3.29%)和野牛(Bisonbison，135，3.22%)(圖2(b))，被注釋到的物種主要為牛科物種。KEGG分析結(jié)果顯示：33 200個(gè)新轉(zhuǎn)錄本被富集到302條通路中，2 252個(gè)新基因被富集到245條通路中。GO數(shù)據(jù)庫(kù)將43 582個(gè)新轉(zhuǎn)錄本及3 170個(gè)新基因進(jìn)行了注釋，結(jié)果可分為三大類，細(xì)胞組分、分子功能及生物學(xué)過(guò)程(圖2(c))。一條基因通過(guò)內(nèi)含子的不同剪接可構(gòu)成不同的轉(zhuǎn)錄本，SMRT測(cè)序技術(shù)可以獲得基因詳盡的轉(zhuǎn)錄本，本研究所預(yù)測(cè)的轉(zhuǎn)錄本數(shù)量遠(yuǎn)高于參考基因組，新轉(zhuǎn)錄本及新基因功能注釋將完善絨山羊轉(zhuǎn)錄組數(shù)據(jù)信息。

圖中n及對(duì)應(yīng)數(shù)字為近緣物種單拷貝基因集和該基因集內(nèi)基因數(shù)目。n and corresponding numbers in the figure are single-copy gene sets of closely related species and the number of genes in the gene set.圖1 絨山羊皮膚轉(zhuǎn)錄組完整性評(píng)估分析Fig.1 Analysis of transcriptome integrity assessment in the skin of cashmere goats

表1 絨山羊皮膚新轉(zhuǎn)錄本及基因功能注釋統(tǒng)計(jì)表Table 1 Statistical table of new transcripts and gene function annotations in the skin of cashmere goats

2.4 絨山羊皮膚可變剪接及l(fā)ncRNA分析

通過(guò)Astalavista軟件鑒定存在的可變剪接類型，主要的基因可變剪接類型如圖3(a)所示。全長(zhǎng)轉(zhuǎn)錄組測(cè)序共鑒定了25 717個(gè)可變剪接，其中可變外顯子(Mutually exclusive exon)631個(gè)、內(nèi)含子保留(Intron retention)9 899個(gè)、外顯子跳躍(Exon skipping)9 184個(gè)、可變轉(zhuǎn)錄起始位點(diǎn)(Alternative 5′splice site)2 325個(gè)、可變轉(zhuǎn)錄終止位點(diǎn)(Alternative 3′splice site)3 678個(gè)(圖3(b))。將AS分析得到的新轉(zhuǎn)錄本對(duì)其編碼區(qū)序列及對(duì)應(yīng)氨基酸序列進(jìn)行預(yù)測(cè)，共獲得ORF 49 294條，其中完整ORF 39 346條，預(yù)測(cè)得到的完整ORF區(qū)編碼蛋白序列長(zhǎng)度分布為：0～100 aa，15 039條(30.51%)、100～200 aa，11 650條(23.63%)、200～300 aa，6 879條(13.96%)、300～400 aa，4 761條(9.66%)及400～500 aa，3 282條(6.66%)。

使用Cytoscape 插件ClueGO對(duì)可變剪接基因進(jìn)行功能富集分析，1 052個(gè)可變剪接基因共富集到624個(gè)GO term。對(duì)富集到的term進(jìn)行分類發(fā)現(xiàn)主要涉及mRNA分解代謝(mRNA catabolic process)、mRNA代謝(mRNA metabolic process)、內(nèi)皮細(xì)胞發(fā)育(Endothelial cell development)、腎小球?yàn)V過(guò)調(diào)節(jié)(Regulation of glomerular filtration)等生物學(xué)過(guò)程，其中核糖體蛋白家族和泛素結(jié)合酶家族基因是富集到term最多的基因家族。同時(shí)找到了一條與皮膚形態(tài)發(fā)育相關(guān)的term(Wnt信號(hào)通路,Wnt signaling pathway)，并基于可變剪接位點(diǎn)信息發(fā)現(xiàn)基因黑素皮質(zhì)素1受體(Melanocortin 1 receptor,MC1R，NC_030825.1:16104986-16105935)與微管蛋白3類(Tubulin-β-III,TUBB3，NC_030825.1:16104986-16116602)共用TUBB3第一外顯子。TUBB3第一內(nèi)含子中含有終止密碼子(TGA)，如在翻譯過(guò)程中保留第一內(nèi)含子，則翻譯過(guò)程就會(huì)在第一內(nèi)含子終止密碼子處停止，產(chǎn)生一個(gè)新蛋白即MC1R所編碼蛋白。

通過(guò)CPC分析、CNCI分析、pfam蛋白結(jié)構(gòu)域分析和CPAT分析4種方法，共鑒定出lncRNA 10 927個(gè)(圖3(c))。根據(jù)基因組上的位置關(guān)系，可將lncRNA分為基因間lncRNA(Intergenic lncRNA，lincRNA)2 585條占總lncRNA的24.4%、內(nèi)含子lncRNA(Intronic lncRNA)5 356條占總lncRNA的50.6%、正義lncRNA(Sense lncRNA)2 125條占總lncRNA的20.1%和反義lncRNA(Antisense lncRNA)528條占總lncRNA的5.0%(圖3(d))?；谖恢藐P(guān)系及互補(bǔ)序列共預(yù)測(cè)出lncRNA 8 094條，其中基于位置關(guān)系預(yù)測(cè)出lncRNA 5 889條占比72.8%。通過(guò)與參考基因組注釋的lncRNA對(duì)比分析，共鑒定出8 251個(gè)新lncRNA，這將為絨山羊基因組數(shù)據(jù)做重要補(bǔ)充。

(a)新轉(zhuǎn)錄本NR同源物種分布圖；(b)新基因NR同源物種分布圖；(c)轉(zhuǎn)錄本GO功能注釋圖，其中藍(lán)色代表細(xì)胞成分，紅色代表分子功能，綠色代表生物過(guò)程。(a)與(b)圖例中數(shù)值表值表示該物種與NR數(shù)據(jù)庫(kù)的比對(duì)注釋數(shù)量，百分比為該物種比對(duì)數(shù)在總比對(duì)數(shù)中的占比。(a) Distribution map of new transcript NR homologous species;(b) New gene NR homologous species distribution map;(c) Annotated map of transcript GO function,in which blue represents cellular components,red represents molecular function,and green represents biological process.The values in (a) and (b) indicate the number of matches of the species with the NR database,and the percentage is the proportion of the number of matches of the species in the total number of matches.圖2 絨山羊皮膚新轉(zhuǎn)錄本及基因功能注釋圖Fig.2 Annotation map of new transcripts and gene functions in the skin of cashmere goats

(a)基因可變剪接類型(a:Exon skipping,外顯子跳躍；b:Alternative 3′ splice site,可變轉(zhuǎn)錄終止位點(diǎn)；c:Mutually exclusive exon,可變外顯子；d:Alternative 5′ splice site,可變轉(zhuǎn)錄起始位點(diǎn)；e:Intron retention,內(nèi)含子保留)；(b)可變剪接事件數(shù)量統(tǒng)計(jì)圖；(c)預(yù)測(cè)得到的非編碼RNA韋恩圖；(d)lncRNA分類圖。(a) Gene alternative splicing type (a:Exon skipping,exon skipping;b:Alternative 3′ splice site,variable transcription termination site;c:Mutually exclusive exon,variable exon;d:Alternative 5′ splice site,alternative transcription initiation site;e:Intron retention,intron retention);(b) Statistics of the number of alternative splicing events;(c) Venn diagram of predicted non-coding RNAs;(d) lncRNA Sequence source pie chart.圖3 絨山羊皮膚可變剪接及l(fā)ncRNA統(tǒng)計(jì)圖Fig.3 Alternative splicing and lncRNA statistics in the skin of cashmere goats

2.5 絨山羊皮膚可變多聚腺苷酸化(APA)及SSR預(yù)測(cè)

通過(guò)分析發(fā)現(xiàn)，SMRT測(cè)序鑒定的基因中9 485個(gè)基因至少存在1個(gè)polyA位點(diǎn)，151個(gè)基因存在至少5個(gè)polyA位點(diǎn)，平均每個(gè)基因含有1.68個(gè) polyA位點(diǎn)(圖4(a))。利用MEME對(duì)所有轉(zhuǎn)錄本polyA位點(diǎn)上下游50 bp的序列進(jìn)行分析，鑒定得到的上下游堿基分布如圖4(b)所示，絨山羊皮膚polyA剪切位點(diǎn)上游富含腺嘌呤，且存在保守元件AATAAA。

(a)ployA位點(diǎn)個(gè)數(shù)統(tǒng)計(jì)圖;(b)APA上下游堿基分布圖;(c)SSR類型分布圖，其中包括單堿基重復(fù)(p1)、雙堿基重復(fù)(p2)、三堿基重復(fù)(p3)、四堿基重復(fù)(p4)、五堿基重復(fù)(p5)、六堿基重復(fù)(p6)和混合SSR(c)。(a) Statistics of the number of ployA sites;(b) Distribution of upstream and downstream bases of APA;(c) SSR type distribution map,including single-base repeats (p1),double-base repeats (p2),three-base repeats (p3),four-base repeats (p4),five-base repeats (p5),six-base repeats (p6),and mixed SSRs (c).圖4 絨山羊皮膚APA及SSR分析統(tǒng)計(jì)圖Fig.4 APA and SSR analysis statistical chart in the skin of cashmere goats

從新轉(zhuǎn)錄本中篩選500 bp以上的轉(zhuǎn)錄本，利用MISA軟件做SSR分析。共評(píng)估序列78 601條，評(píng)估序列的總堿基數(shù)為253 230 057 bp，識(shí)別的SSR總數(shù)為66 951條。鑒定出7種類型SSR:Mono-nucleotide(單堿基)42 864條、Di-nucleotide(雙堿基)13 696條、Tri-nucleotide(三堿基)9 300條、Tetra-nucleotide(四堿基)706條、Penta-nucleotide(五堿基)327條、Hexa-nucleotide(六堿基)58條、Compound SSR(混合微衛(wèi)星，2個(gè)SSR距離小于100 bp)6 898條，對(duì)不同SSR類型的密度分布進(jìn)行統(tǒng)計(jì)如圖4(c)。

2.6 絨山羊皮膚TF預(yù)測(cè)分析

轉(zhuǎn)錄因子(Transcription factor,TF)是指能夠結(jié)合在某基因上游特異核苷酸序列上的蛋白質(zhì)，本研究分析TF使用動(dòng)物轉(zhuǎn)錄因子數(shù)據(jù)庫(kù)animalTFDB 2.0[16]，共預(yù)測(cè)得到轉(zhuǎn)錄因子5 641個(gè)，對(duì)不同類型的轉(zhuǎn)錄因子個(gè)數(shù)進(jìn)行統(tǒng)計(jì)如圖5。預(yù)測(cè)的TF類型主要是Zf-c2h2，預(yù)測(cè)數(shù)量是1 920，占總預(yù)測(cè)量的34.1%，其次是Zbtb、Homeobox以及Bhlh等，c2h2家族是已知TF家族中比較大的一類，能結(jié)合多種蛋白質(zhì)發(fā)揮其功能。

3 討 論

圖5 絨山羊皮膚轉(zhuǎn)錄因子類型預(yù)測(cè)Fig.5 Transcription factor type prediction in the skin of cashmere goats

近年來(lái)，隨著國(guó)內(nèi)外絨山羊遺傳育種研究的推進(jìn)，絨山羊在產(chǎn)絨、繁殖、生長(zhǎng)等多個(gè)重要經(jīng)濟(jì)性狀上得到不斷地改進(jìn)和提升[17]，然而選育絨山羊新品種、提高絨山羊優(yōu)良產(chǎn)絨性狀以及環(huán)境與基因互作仍然是現(xiàn)階段研究熱點(diǎn)。隨著測(cè)序技術(shù)的進(jìn)步，SMRT為分析轉(zhuǎn)錄組特征及mRNA 結(jié)構(gòu)等工作提供了新的方法，雖然第二代測(cè)序技術(shù)比Sanger測(cè)序提供了巨大的改進(jìn)，但它們的局限性不適用于某些特定的生物學(xué)問(wèn)題，由PacBio開(kāi)發(fā)的SMRT測(cè)序提供了一種替代方法來(lái)克服這些限制[18]，與第二代測(cè)序技術(shù)不同，PacBio測(cè)序是一種實(shí)時(shí)測(cè)序方法，在讀取步驟之間不需要暫停[19]。全長(zhǎng)轉(zhuǎn)錄組序列對(duì)于品種選育、性狀提升、基因組注釋和功能研究具有非常重要的作用。雖然現(xiàn)階段山羊以具有良好的遺傳和基因組資源，包括高度連續(xù)的參考基因組(ARS1)[20]。然而，Muriuki等[21]認(rèn)為與其他反芻動(dòng)物相比，山羊參考基因組基因表達(dá)信息是有限的，因此，Muriuki等[21]構(gòu)建了山羊基因表達(dá)圖譜，提供了一組功能信息來(lái)注釋當(dāng)前參考基因組(ARS1)中15 %未注釋的基因。本研究通過(guò) SMRT 對(duì)絨山羊皮膚進(jìn)行全長(zhǎng)轉(zhuǎn)錄組測(cè)序分析，共獲得42.41 Gb clean data的數(shù)據(jù)，獲得CCS 591 585條，平均長(zhǎng)度2 358 bp，其中FLNC序列466 058條占CCS序列的78.78%。SMRT測(cè)序技術(shù)對(duì)發(fā)現(xiàn)新基因與新轉(zhuǎn)錄本具有絕對(duì)的優(yōu)勢(shì)，通過(guò)SMRT測(cè)序在毛竹中發(fā)現(xiàn)了8 091個(gè)新轉(zhuǎn)錄本[22]，在穿山甲中鑒定出8 014個(gè)新基因[23]。本研究測(cè)序共檢測(cè)到基因18 919個(gè)，其中新基因6 166個(gè)，新發(fā)現(xiàn)轉(zhuǎn)錄本55 875個(gè)。對(duì)新發(fā)現(xiàn)的轉(zhuǎn)錄本進(jìn)行序列結(jié)構(gòu)分析，共預(yù)測(cè)出66 951 個(gè)SSR及10 927個(gè)lncRNA，并完成了49 573 個(gè)新轉(zhuǎn)錄本和4 237個(gè)新基因的功能注釋。本研究所預(yù)測(cè)到的轉(zhuǎn)錄本數(shù)量遠(yuǎn)高于參考基因組，該數(shù)據(jù)集補(bǔ)充了山羊現(xiàn)有的遺傳和基因組資源。具有高質(zhì)量功能注釋的高度連續(xù)參考基因組對(duì)于完善基因信息、開(kāi)發(fā)畜禽物種新資源及對(duì)品種優(yōu)良性狀的挖掘是非常必要的。

SSR作為一種共顯性遺傳標(biāo)記，具有特異性強(qiáng)、重復(fù)性好和多態(tài)性高等優(yōu)點(diǎn)，被廣泛用于遺傳多樣性等相關(guān)研究[24]。目前，轉(zhuǎn)錄組數(shù)據(jù)信息在?？苿?dòng)物[25]、山羊[26]等物種中成功地開(kāi)發(fā)了SSR標(biāo)記。本研究中共評(píng)估序列78 601條，識(shí)別的SSR總數(shù)為66 951條，其中單堿基SSR 42 864條、雙堿基SSR 13 696條、三堿基SSR 9 300條、四堿基SSR 706條、五堿基SSR 327條、六堿基SSR 58條、混合微衛(wèi)星6 898條。上述對(duì)SSR分析，為進(jìn)一步做絨山羊特異SSR標(biāo)記及遺傳圖譜構(gòu)建分析等研究打下了良好基礎(chǔ)。

SMRT 測(cè)序技術(shù)還能更好的分析AS，使用該技術(shù)在絨山羊皮膚中共鑒定了25 717個(gè)AS，其中內(nèi)含子保留9 899個(gè)為主要的剪接類型。Xu等[27]對(duì)山羊中的AS分析的數(shù)據(jù)結(jié)果顯示，總共有14 521 個(gè)基因經(jīng)歷了AS，內(nèi)含子保留最為普遍占AS事件總數(shù)的37.04%，與本研究AS預(yù)測(cè)結(jié)果相一致。本研究發(fā)現(xiàn)ClueGO將可變剪接基因富集到內(nèi)皮細(xì)胞發(fā)育、mRNA分解代謝、mRNA代謝等生物學(xué)過(guò)程，同一個(gè)基因轉(zhuǎn)錄而來(lái)的前體mRNA通過(guò)可變剪切可形成不同的剪接異構(gòu)體，最終形成不同的蛋白質(zhì)而發(fā)揮不同的功能[28]，提示在動(dòng)物體生長(zhǎng)發(fā)育過(guò)程中參與上述過(guò)程的基因容易發(fā)生可變剪接。結(jié)合功能富集結(jié)果及可變剪接發(fā)生位點(diǎn)信息發(fā)現(xiàn)了MC1R與TUBB3基因間存在可變剪接位點(diǎn)，共用TUBB3第一外顯子。Herraiz等[29]研究結(jié)果表明，UV會(huì)導(dǎo)致人體內(nèi)MC1R和TUBB3基因向MC1R-TUBB3嵌合體發(fā)生轉(zhuǎn)化。當(dāng)暴露于 UV 后，MC1R-TUBB3嵌合異構(gòu)體會(huì)阻止cAMP途徑的過(guò)度刺激[29],cAMP途徑的過(guò)度刺激會(huì)導(dǎo)致生物體皮膚變黑，在人體內(nèi)過(guò)度刺激該途徑可能會(huì)增加黑色素瘤患病風(fēng)險(xiǎn)[30]。MC1R與TUBB3基因間發(fā)生可變剪接可能是山羊及綿羊在日常接觸UV輻射而觸發(fā)的一種自我保護(hù)機(jī)制，阻止 cAMP 途徑的過(guò)度刺激減少了真黑素的產(chǎn)生使背負(fù)毛發(fā)顏色沒(méi)有發(fā)生改變的同時(shí)又防止了UV對(duì)其產(chǎn)生的危害。同時(shí)，lncRNA 作為目前研究的熱點(diǎn)，在絨山羊毛囊生長(zhǎng)發(fā)育中起著重要的調(diào)控作用[31]。Wang等[32]基于高通量測(cè)序和生物信息學(xué)分析，在絨山羊皮膚的生長(zhǎng)期和休止期鑒定了1 108個(gè) lncRNA。Ge等[33]通過(guò)lncRNA觀點(diǎn)綜合分析了褪黑激素促進(jìn)絨山羊次級(jí)毛囊生長(zhǎng)。上述研究表明了lncRNA在絨山羊毛囊發(fā)育及皮膚穩(wěn)態(tài)過(guò)程中發(fā)揮著不可或缺的作用。本研究通過(guò)CPC、CNCI、pfam蛋白結(jié)構(gòu)、CPAT四種分析方法，共鑒定出lncRNA 10 927條，Intronic lncRNA 5 356條為主要類型，而參考基因組(ARS1)注釋lncRNA數(shù)量?jī)H為2 676條平均長(zhǎng)度1 434 bp[20]。本研究測(cè)序獲得的lncRNA預(yù)測(cè)數(shù)據(jù)遠(yuǎn)高于參考基因組及其前人的試驗(yàn)發(fā)現(xiàn)。SMRT測(cè)序鑒定的基因中9 485個(gè)基因至少存在1個(gè)poly(A)位點(diǎn)，預(yù)測(cè)到轉(zhuǎn)錄因子5 641個(gè)主要類型是Zf-c2h2，這與在巨菌草[34]及北京鴨[35]上的試驗(yàn)研究結(jié)果相一致?？勺兗艚臃治鲋械玫降男罗D(zhuǎn)錄本對(duì)其編碼區(qū)序列及對(duì)應(yīng)氨基酸序列進(jìn)行預(yù)測(cè)獲得ORF 49 294條，其中完整ORF 39 346條，本研究獲得的數(shù)據(jù)將豐富絨山羊轉(zhuǎn)錄組信息，同時(shí)也為絨山羊基因組做出重要的數(shù)據(jù)補(bǔ)充。

4 結(jié) 論

本研究對(duì)內(nèi)蒙古白絨山羊皮膚樣品進(jìn)行了全長(zhǎng)轉(zhuǎn)錄組測(cè)序并對(duì)SSR、ORF、可變剪接、新基因及新轉(zhuǎn)錄本進(jìn)行了系統(tǒng)分析。試驗(yàn)結(jié)果為研究絨山羊皮膚mRNA全長(zhǎng)序列、掌握完整的基因結(jié)構(gòu)信息、選育絨山羊新品種和優(yōu)良性狀的挖掘提供了理論依據(jù)，為絨山羊基因組資源進(jìn)行數(shù)據(jù)補(bǔ)充及提供有價(jià)值的參考。