韓若冰 李和平

(東北林業(yè)大學(xué)，哈爾濱，150040)

鹿茸的快速生長及其每年的周期性再生現(xiàn)象在哺乳動物中極為罕見。在梅花鹿鹿茸的快速生長期，其生長速度可達到12.5 mm·d-1[1]，不僅如此，鹿茸的生長速度堪比癌細胞，卻極少發(fā)生癌變，這些獨特的生物學(xué)特性使得鹿茸受到研究者的廣泛關(guān)注。了解鹿茸生長發(fā)育的分子調(diào)控機制對于鹿茸生物學(xué)的研究甚至是癌癥治療的進一步研究均具有重要的科學(xué)意義。

迄今為止，已有眾多研究關(guān)注鹿茸生長發(fā)育的分子調(diào)控機制。隨著高通量測序技術(shù)的迅速發(fā)展，選用梅花鹿鹿茸作為試驗材料，從轉(zhuǎn)錄組水平上探究鹿茸生長發(fā)育的研究應(yīng)運而生。通過對不同生長時期的梅花鹿鹿茸進行轉(zhuǎn)錄組測序，得到了幾種可能對鹿茸生長發(fā)揮重要調(diào)控作用的轉(zhuǎn)錄因子：PER1、EGFR1、GAS1[2]。通過對骨化期梅花鹿鹿茸的轉(zhuǎn)錄組測序得到了Sp7、SOX4、MMP9、MMP13等30種與鹿茸軟骨內(nèi)骨化相關(guān)的基因[3]。然而關(guān)于非編碼RNA在鹿茸生長發(fā)育過程中的研究還較為匱乏，主要集中于miRNA的研究，Yao et al.[4]通過高通量測序技術(shù)鑒定得到一些參與鹿茸軟骨發(fā)生和軟骨內(nèi)骨化的miRNA(miR-3072-5p、miR-1600、miR-34-5p、miR-6889-5p、miR-6729-5p)。

環(huán)狀RNA是一類在轉(zhuǎn)錄后水平具有重要基因表達調(diào)控作用的內(nèi)源性非編碼RNA[5-6]，不具有5’末端及3’末端。近些年來已經(jīng)成為1個新的研究熱點。然而，對于鹿科動物中環(huán)狀RNA的研究還較少。本研究運用高通量測序技術(shù)鑒定鹿茸2個生長發(fā)育時期(鹿茸生長前期、快速生長期)的環(huán)狀RNA并進行生物信息學(xué)分析，試圖挖掘?qū)β谷咨L發(fā)育發(fā)揮調(diào)控作用的環(huán)狀RNA，旨在為梅花鹿的多組學(xué)研究及其生長發(fā)育機理的探究提供參考。

1 材料與方法

1.1 樣品采集

試驗材料采集自哈爾濱金地鹿業(yè)有限公司的3只健康成年雄性梅花鹿鹿茸的頂端間充質(zhì)組織。在其鹿茸生長至小鞍子茸型期，取左茸頂端組織的間充質(zhì)組織為鹿茸生長前期樣品；在其生長至二杠茸型時，采集右茸頂端組織的間充質(zhì)組織為快速生長期樣品。樣品采集后，將間充質(zhì)組織切成小塊，迅速放入液氮中保存?zhèn)溆谩?/p>

1.2 RNA提取及文庫構(gòu)建、測序

按照TRIzol試劑盒的操作說明書提取6個間充質(zhì)組織樣品的總RNA，使用分光光度計檢測提取總RNA的濃度、純度并使用Agilent Bioanalyzer檢測總RNA的完整性。提取得到的總RNA放置于冰箱-80 ℃保存以防止RNA降解。得到質(zhì)檢合格的RNA后，去除其中的核糖體RNA，合成cDNA的第一鏈和第二鏈并且進行PCR擴增，利用Illumina HiseqTM 4000對構(gòu)建的文庫進行測序。

1.3 環(huán)狀RNA的鑒定與特征分析

利用TopHat[7]將過濾后的數(shù)據(jù)與馬鹿參考基因組進行比對，在比對結(jié)果中提取未比對上的數(shù)據(jù)，截取每一條未比對上的數(shù)據(jù)的兩端以得到短序列讀段，再次比對到參考基因組，并利用find_circ[5]軟件鑒定得到環(huán)狀RNA。鑒定條件為：只保留有且只有1個清晰斷點的環(huán)狀RNA；每條讀長的2個短序列讀段比對到基因組上的位置重疊不能超過2 bp且只允許2 bp錯配；獨特的讀長需大于2條；環(huán)狀RNA的長度小于100 kbp等。鑒定得到環(huán)狀RNA后，對環(huán)狀RNA類型、長度及環(huán)狀RNA在參考基因組上的分布進行統(tǒng)計。

1.4 表達量計算與差異分析

采用每千個堿基的轉(zhuǎn)錄每百萬映射讀取的轉(zhuǎn)錄本(RPM)計算環(huán)狀RNA的表達量。根據(jù)得到的環(huán)狀RNA的表達量對鹿茸2個生長發(fā)育時期的環(huán)狀RNA進行差異分析，從而得到顯著差異表達的環(huán)狀RNA，顯著差異表達的環(huán)狀RNA的篩選條件設(shè)置為P<0.05，|log2FC|>1。

1.5 富集分析

為探究環(huán)狀RNA在鹿茸生長發(fā)育過程中發(fā)揮的作用，對環(huán)狀RNA的來源基因進行GO富集分析和KEGG功能注釋分類。在差異表達基因中顯著性富集的通路以及差異表達基因中顯著富集的GO條目篩選閾值均設(shè)置為q≤0.05。

1.6 circRNA-miRNA-mRNA調(diào)控網(wǎng)絡(luò)的構(gòu)建

運用RNAhybrid (v2.1.2) 、svm_light (v6.01)，Miranda (v3.3)、TargetScan (v7.0)對最顯著的10個環(huán)狀RNA進行靶向miRNA的預(yù)測，結(jié)合本研究前期測得的miRNA數(shù)據(jù)[8]，從中篩選差異表達的miRNA。除此之外，利用mirTarBase預(yù)測得到與這些差異表達的miRNA具有靶向關(guān)系的mRNA，并從中篩選出差異表達的mRNA，進而得到miRNA-mRNA靶向作用關(guān)系，最終利用cytoscape[9]構(gòu)建與鹿茸生長發(fā)育相關(guān)的circRNA-miRNA-mRNA調(diào)控網(wǎng)絡(luò)。

2 結(jié)果與分析

2.1 環(huán)狀RNA的鑒定

通過對測序得到的原始數(shù)據(jù)進行過濾，得到的Q20為98.93%～99.05%，Q30為96.25%～96.58%，GC含量為51.35%～53.34%。從比對結(jié)果中提取未比對上的數(shù)據(jù)并截取其兩端，得到的短序列讀段為45865210～60425688，比對到參考基因組上的比對率為63.83%～67.32%。

在鹿茸發(fā)育的前期和快速生長期共鑒定得到7 790個環(huán)狀RNA。環(huán)狀RNA主要來源于蛋白質(zhì)編碼基因的已注釋的外顯子、單一外顯子、基因反義鏈、內(nèi)含子、外顯子和內(nèi)含子和基因間區(qū)[5，10]。對鑒定得到的環(huán)狀RNA類型進行了統(tǒng)計(表1)，其中，大多數(shù)環(huán)狀RNA來自于蛋白質(zhì)編碼基因的基因間區(qū)及外顯子、內(nèi)含子區(qū)，來自于內(nèi)含子區(qū)的環(huán)狀RNA數(shù)目最少，而環(huán)狀RNA的長度主要分布于1～1 000之間(表2)。對環(huán)狀RNA在染色體上的分布進行統(tǒng)計(表3)，結(jié)果表明7 790個環(huán)狀RNA分布在參考基因組的35條染色體上，其中分布在CM008042.1上的環(huán)狀RNA數(shù)量最少，而分布在CM008012.1上的環(huán)狀RNA數(shù)量最多。

表1 環(huán)狀RNA的類型分布

表2 環(huán)狀RNA的長度分布

表3 環(huán)狀RNA的染色體分布

2.2 環(huán)狀RNA的差異表達分析

通過對鑒定得到的2個發(fā)育時期的環(huán)狀RNA表達量進行差異基因的篩選，共得到208個顯著差異表達的環(huán)狀RNA，包括顯著上調(diào)的39個差異表達環(huán)狀RNA以及顯著下調(diào)的169個差異表達環(huán)狀RNA(圖1)。差異最顯著的10個環(huán)狀RNA的詳細信息見表4。

紅色表示表達量上調(diào)的環(huán)狀RNA；綠色表示表達量下調(diào)的環(huán)狀RNA；黑色表示表達量沒有差異。

表4 差異最顯著的前10個環(huán)狀RNA

2.3 差異表達基因的富集分析

對差異環(huán)狀RNA的來源基因進行富集分析可用于環(huán)狀RNA的功能預(yù)測。GO富集分析結(jié)果表明，環(huán)狀RNA的來源基因富集到166條與分子功能相關(guān)的GO條目上，顯著富集到生物調(diào)控，細胞過程、代謝過程等GO條目上；富集到156條與細胞組分相關(guān)的GO條目上，顯著富集到細胞、細胞部分等GO條目上；富集到980條與生物過程相關(guān)的GO條目上，包括結(jié)合與催化活性等相關(guān)的GO條目(表5)。

表5 差異表達的circRNA的來源基因的GO富集分析

環(huán)狀RNA來源基因的KEGG結(jié)果顯示共富集到81條信號通路上，包括粘著斑、PI3K-Akt信號通路、MAPK信號通路、TNF信號通路、FoxO信號通路(圖2)。其中，有13個差異顯著的環(huán)狀RNA的6個來源基因顯著富集到與細胞增殖、細胞周期等相關(guān)的PI3K-Akt信號通路上，分別是novel_circ_004050、novel_circ_004052、novel_circ_004059、novel_circ_004060、novel_circ_004062、novel_circ_004064、novel_circ_004069、novel_circ_004074、novel_circ_004863、novel_circ_002034、novel_circ_004646、novel_circ_003420、novel_circ_006430。

圖2 差異表達環(huán)狀RNA的來源基因的KEGG富集分析

2.4 控網(wǎng)絡(luò)的構(gòu)建

通過預(yù)測得到的環(huán)狀RNA、miRNA、mRNA之間的靶向互作關(guān)系構(gòu)建最終的調(diào)控網(wǎng)絡(luò)，我們發(fā)現(xiàn)在10個差異最顯著的環(huán)狀RNA中，有2個環(huán)狀RNA預(yù)測得到的miRNA在鹿茸發(fā)育的前期和快速生長期的表達量差異不顯著。因此，最終的調(diào)控網(wǎng)絡(luò)由8個差異表達的環(huán)狀RNA和與其具有靶向互作關(guān)系的25個差異表達的miRNA所構(gòu)成，通過進一步預(yù)測，得到這些差異表達的miRNA有120個存在靶向關(guān)系的差異表達的mRNA(圖3)。

圖3 差異表達的circRNA -miRNA-mRNA調(diào)控網(wǎng)絡(luò)

3 討論

自從環(huán)狀RNA在植物類RNA病毒中首次被發(fā)現(xiàn)后[11]，在多種生物體內(nèi)均發(fā)現(xiàn)有環(huán)狀RNA的存在[12]，其表達具有時空特異性。根據(jù)競爭性內(nèi)源RNA機制(ceRNA)，環(huán)狀RNA通過與miRNA競爭性結(jié)合進而調(diào)控下游靶基因的表達[13]，在多種腫瘤發(fā)生及組織器官的生長發(fā)育過程中均發(fā)揮著重要的調(diào)控作用，但在動物體中關(guān)于環(huán)狀RNA的研究還較少，且主要集中于環(huán)狀RNA的鑒定。

本研究首次鑒定了鹿茸生長發(fā)育過程中的環(huán)狀RNA，得到鹿茸發(fā)育前期及快速生長期的間充質(zhì)組織中共有7 790個環(huán)狀RNA，其中1 901個環(huán)狀RNA僅在鹿茸發(fā)育前期具有特異性表達，1 478個環(huán)狀RNA僅在鹿茸的快速生長期有特異性表達，有208個環(huán)狀RNA為顯著差異表達。我們猜測這些在不同生長發(fā)育時期有特異性表達的環(huán)狀RNA在鹿茸的某個生長發(fā)育時期發(fā)揮著獨特的調(diào)控作用。除此以外，我們發(fā)現(xiàn)一些環(huán)狀RNAs在鹿茸發(fā)育的2個時期均為顯著高表達(表達量高于其他環(huán)狀RNA 100倍以上，例如novel_circ_002441、novel_circ_007423、novel_circ_001328、novel_circ_007778、novel_circ_001447、novel_circ_005051、novel_circ_000884、novel_circ_005048、novel_circ_005975)，這些高表達的環(huán)狀RNA通過與miRNA及mRNA之間的靶向作用關(guān)系對鹿茸的生長發(fā)育發(fā)揮著某些重要的調(diào)控作用，但具體的作用機制還需進一步的功能試驗來探究。

環(huán)狀RNA的功能分析主要通過其來源基因富集分析及與其具有靶向作用關(guān)系的miRNA和mRNA來進行，顯著差異表達的環(huán)狀RNA來源基因的GO富集分析主要集中于生物調(diào)控、細胞過程、代謝過程、細胞代謝過程，而KEGG富集到代謝通路、粘著斑、PI3K-Akt信號通路、MAPK信號通路、FoxO信號通路上，說明差異顯著的環(huán)狀RNA主要涉及鹿茸發(fā)育的生物調(diào)控以及代謝等過程。在構(gòu)建的調(diào)控網(wǎng)絡(luò)中，1個環(huán)狀RNA可能會同時靶向多個miRNA，并且1個miRNA也可能會與多個環(huán)狀RNA具有靶向互作關(guān)系。其中，novel_circ_004646的靶向能力最強，靶向了13個差異表達的miRNA，而novel_circ_005407、novel_circ_000593只靶向了1個差異表達的miRNA。在與環(huán)狀RNA具有靶向關(guān)系的miRNA中，novel_circ_004646的其中1個靶向miR-457-x的靶向能力最強，與29個差異表達的mRNA均具有互作關(guān)系。但環(huán)狀RNA、miRNA、mRNA三者之間是否具有真實的靶向關(guān)系，其對鹿茸發(fā)育發(fā)揮著怎樣的調(diào)控作用以及如何發(fā)揮這些調(diào)控作用都需要進一步的研究。