絨山羊皮膚全長轉錄組測序及生物信息學分析

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

(內蒙古民族大學 動物科技學院，內蒙古 通遼 028000)

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

1 材料與方法

1.1 試驗材料

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

1.2 試驗方法

1.2.1絨山羊皮膚樣本RNA提取

每100 mg新鮮皮膚組織樣品中加入1 mL Trizol試劑，使用無菌手術刀在冰上將樣品切碎，并用無菌勻漿器勻漿，將組織或細胞裂解液轉移到1.5 mL 無RNA酶EP管中。置冰上，靜置5 min。每管加入200 μL氯仿，充分混合后冰上靜置10 min，使核蛋白復合物完全解離。4 ℃ 13 000 r/min離心15 min。期間，取一新EP管，加入500 μL異丙醇，冰上預冷。離心后，將上層水相 (約500 μL)轉移到該新的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處理過的去離子水中，最后進行分光光度分析以確定樣品濃度和純度。

1.2.2絨山羊皮膚全長轉錄組建庫及SMRT測序

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

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

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

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

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

2 結果與分析

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

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

2.2 絨山羊皮膚全長轉錄組測序完整性評估

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

2.3 絨山羊皮膚新轉錄本及新基因功能注釋

將得到的新轉錄本和新基因與NR、Swiss-prot、GO、COG、KOG、Pfam和KEGG等數(shù)據(jù)庫比對進行功能注釋，各數(shù)據(jù)庫注釋到的轉錄本數(shù)量統(tǒng)計見表1。由表1可知，共對49 573個轉錄本進行注釋，其中通過NR數(shù)據(jù)庫成功注釋的轉錄本數(shù)量最多，共有49 374條，在總注釋量中占比99.60%。新轉錄本NR數(shù)據(jù)庫注釋分布前幾位的分別是：山羊(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ù)庫注釋分布前幾位的分別是：牛(Bostaurus，1 787，42.65%)、野牦牛(Bosmutus，1 318，31.46%)、山羊(Caprahircus，182，4.34%)、綿羊(Ovisaries，138，3.29%)和野牛(Bisonbison，135，3.22%)(圖2(b))，被注釋到的物種主要為?？莆锓N。KEGG分析結果顯示：33 200個新轉錄本被富集到302條通路中，2 252個新基因被富集到245條通路中。GO數(shù)據(jù)庫將43 582個新轉錄本及3 170個新基因進行了注釋，結果可分為三大類，細胞組分、分子功能及生物學過程(圖2(c))。一條基因通過內含子的不同剪接可構成不同的轉錄本，SMRT測序技術可以獲得基因詳盡的轉錄本，本研究所預測的轉錄本數(shù)量遠高于參考基因組，新轉錄本及新基因功能注釋將完善絨山羊轉錄組數(shù)據(jù)信息。

圖中n及對應數(shù)字為近緣物種單拷貝基因集和該基因集內基因數(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 絨山羊皮膚轉錄組完整性評估分析Fig.1 Analysis of transcriptome integrity assessment in the skin of cashmere goats

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

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

通過Astalavista軟件鑒定存在的可變剪接類型，主要的基因可變剪接類型如圖3(a)所示。全長轉錄組測序共鑒定了25 717個可變剪接，其中可變外顯子(Mutually exclusive exon)631個、內含子保留(Intron retention)9 899個、外顯子跳躍(Exon skipping)9 184個、可變轉錄起始位點(Alternative 5′splice site)2 325個、可變轉錄終止位點(Alternative 3′splice site)3 678個(圖3(b))。將AS分析得到的新轉錄本對其編碼區(qū)序列及對應氨基酸序列進行預測，共獲得ORF 49 294條，其中完整ORF 39 346條，預測得到的完整ORF區(qū)編碼蛋白序列長度分布為：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對可變剪接基因進行功能富集分析，1 052個可變剪接基因共富集到624個GO term。對富集到的term進行分類發(fā)現(xiàn)主要涉及mRNA分解代謝(mRNA catabolic process)、mRNA代謝(mRNA metabolic process)、內皮細胞發(fā)育(Endothelial cell development)、腎小球濾過調節(jié)(Regulation of glomerular filtration)等生物學過程，其中核糖體蛋白家族和泛素結合酶家族基因是富集到term最多的基因家族。同時找到了一條與皮膚形態(tài)發(fā)育相關的term(Wnt信號通路,Wnt signaling pathway)，并基于可變剪接位點信息發(fā)現(xiàn)基因黑素皮質素1受體(Melanocortin 1 receptor,MC1R，NC_030825.1:16104986-16105935)與微管蛋白3類(Tubulin-β-III,TUBB3，NC_030825.1:16104986-16116602)共用TUBB3第一外顯子。TUBB3第一內含子中含有終止密碼子(TGA)，如在翻譯過程中保留第一內含子，則翻譯過程就會在第一內含子終止密碼子處停止，產(chǎn)生一個新蛋白即MC1R所編碼蛋白。

通過CPC分析、CNCI分析、pfam蛋白結構域分析和CPAT分析4種方法，共鑒定出lncRNA 10 927個(圖3(c))。根據(jù)基因組上的位置關系，可將lncRNA分為基因間lncRNA(Intergenic lncRNA，lincRNA)2 585條占總lncRNA的24.4%、內含子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))。基于位置關系及互補序列共預測出lncRNA 8 094條，其中基于位置關系預測出lncRNA 5 889條占比72.8%。通過與參考基因組注釋的lncRNA對比分析，共鑒定出8 251個新lncRNA，這將為絨山羊基因組數(shù)據(jù)做重要補充。

(a)新轉錄本NR同源物種分布圖；(b)新基因NR同源物種分布圖；(c)轉錄本GO功能注釋圖，其中藍色代表細胞成分，紅色代表分子功能，綠色代表生物過程。(a)與(b)圖例中數(shù)值表值表示該物種與NR數(shù)據(jù)庫的比對注釋數(shù)量，百分比為該物種比對數(shù)在總比對數(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 絨山羊皮膚新轉錄本及基因功能注釋圖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,可變轉錄終止位點；c:Mutually exclusive exon,可變外顯子；d:Alternative 5′ splice site,可變轉錄起始位點；e:Intron retention,內含子保留)；(b)可變剪接事件數(shù)量統(tǒng)計圖；(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)計圖Fig.3 Alternative splicing and lncRNA statistics in the skin of cashmere goats

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

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

(a)ployA位點個數(shù)統(tǒng)計圖;(b)APA上下游堿基分布圖;(c)SSR類型分布圖，其中包括單堿基重復(p1)、雙堿基重復(p2)、三堿基重復(p3)、四堿基重復(p4)、五堿基重復(p5)、六堿基重復(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)計圖Fig.4 APA and SSR analysis statistical chart in the skin of cashmere goats

從新轉錄本中篩選500 bp以上的轉錄本，利用MISA軟件做SSR分析。共評估序列78 601條，評估序列的總堿基數(shù)為253 230 057 bp，識別的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個SSR距離小于100 bp)6 898條，對不同SSR類型的密度分布進行統(tǒng)計如圖4(c)。

2.6 絨山羊皮膚TF預測分析

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

3 討 論

圖5 絨山羊皮膚轉錄因子類型預測Fig.5 Transcription factor type prediction in the skin of cashmere goats

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

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

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

4 結 論

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