Illumina高通量測(cè)序技術(shù)研究肌少-骨質(zhì)疏松癥與骨質(zhì)疏松癥患者骨組織microRNAs差異表達(dá)譜及差異分析

陳錦成 朱國(guó)濤 秦曉飛 陳彥丞 羅駿 劉洪文 余博飛 徐杰*

1.福建中醫(yī)藥大學(xué)中醫(yī)骨傷及運(yùn)動(dòng)康復(fù)教育部重點(diǎn)實(shí)驗(yàn)室，福建 福州 350122 2.福建醫(yī)科大學(xué)省立臨床醫(yī)學(xué)院(福建省立醫(yī)院)骨科，福建 福州 350001

肌少癥是指與自然增齡相關(guān)的，其特征是肌肉強(qiáng)度的下降和(或)全身肌量的減少或人體肌肉的生理活動(dòng)功能衰退狀況呈進(jìn)行性發(fā)展[1-2]。骨質(zhì)疏松癥(osteoporosis，OP)屬于全身性骨病范疇，其特征是骨量逐步減少、骨的微觀結(jié)構(gòu)退化及脆性增加，其系統(tǒng)性退化伴隨骨折發(fā)生率的增加[2-3]。肌少-骨質(zhì)疏松癥(sarco-osteoporosis，SO)的概念是由Binkley 等[4]提出的，目前臨床上尚未有行業(yè)專(zhuān)家一致公認(rèn)的診治標(biāo)準(zhǔn)，本課題組前期的臨床實(shí)驗(yàn)研究規(guī)定的SO患者納入標(biāo)準(zhǔn)是合并骨質(zhì)疏松(腰椎BMD T評(píng)分<-2.5)及肌肉功能異常(肌力0～3 級(jí)或出現(xiàn)明顯肌肉萎縮)的全髖關(guān)節(jié)置換患者；OP患者納入標(biāo)準(zhǔn)是合并骨質(zhì)疏松但不伴肌肉功能異常的全髖關(guān)節(jié)置換患者。排除標(biāo)準(zhǔn)：已知體內(nèi)存在影響骨代謝的疾病，包括但不限于：糖尿病、甲狀腺疾病、腫瘤性疾病患者、免疫系統(tǒng)疾病等；近12個(gè)月內(nèi)服用影響骨代謝的藥物，包括但不限于：雙膦酸鹽制劑、利尿劑、腎上腺和甲狀旁腺藥物、糖皮質(zhì)激素、維生素D和鈣劑。符合以上任1條件可排除試驗(yàn)病例。根據(jù)前期課題組研究了SO與OP的肌肉組織與骨組織差異蛋白質(zhì)譜分析結(jié)果[5]，證實(shí)了SO組與OP組存在相互作用的顯著差異蛋白，那么SO與OP在miRNA水平是否也同樣存在差異基因調(diào)控其mRNA及下游靶基因。目前研究miRNA調(diào)控下游靶基因以揭示其潛在的調(diào)控關(guān)系，可以通過(guò)單純的生物信息學(xué)多網(wǎng)站和數(shù)據(jù)庫(kù)進(jìn)行預(yù)測(cè)分析靶基因[6]，但往往會(huì)得到龐大的靶向調(diào)控關(guān)系結(jié)果的籠統(tǒng)數(shù)據(jù)，其中不可避免存在大量的假陽(yáng)性結(jié)果[7]。因此，本課題組從臨床收集病例，術(shù)中分別取SO組與OP組的骨組織樣品，并運(yùn)用Illumina高通量測(cè)序技術(shù)鑒定、篩選肌少-骨質(zhì)疏松癥與骨質(zhì)疏松癥患者骨組織microRNAs的差異表達(dá)譜[8]，以期篩查出在肌少-骨質(zhì)疏松癥疾病演變過(guò)程中發(fā)揮重要的調(diào)控作用的候選miRNA。

1 材料和方法

1.1 一般資料與取材

選取2017年至2019年期間在福建省老年醫(yī)院和福建省立醫(yī)院骨二科行手術(shù)治療的股骨頸骨折患者中，隨機(jī)選擇并自愿參與本臨床研究6例患者的股骨頭組織(2～5 g)，SO組與OP組各選取3例；術(shù)中取一側(cè)開(kāi)路器切除的粗隆部位松質(zhì)骨組織并置于液氮保存。兩個(gè)組別的納入對(duì)象均為漢族女性。本課題已通過(guò)福建省立醫(yī)院的倫理委員會(huì)批準(zhǔn)，明確患者及其近親家屬知情同意且簽署了知情同意書(shū)，倫理審批號(hào)為K2019-03-034。

1.2 Illumina高通量測(cè)序建庫(kù)

首先提取SO組與OP組骨組織樣品的total RNA，利用small RNA的3'及5’端特殊結(jié)構(gòu)直接在small RNA兩端加上接頭，然后反轉(zhuǎn)錄為c-DNA并進(jìn)行PCR擴(kuò)增，再通過(guò)凝膠電泳分離目的片段，驗(yàn)證此次實(shí)驗(yàn)樣品RNA完整性后構(gòu)建small RNA文庫(kù)，隨后可進(jìn)行二代Illumina高通量測(cè)序。根據(jù)諾禾致源SE50測(cè)序策略進(jìn)行，使用TruSeq SR Cluster Kit v3-cBot-HS(Illumia)在cBot群集生成系統(tǒng)上對(duì)索引編碼的樣本進(jìn)行群集，產(chǎn)生簇后將文庫(kù)制備物在Illumina Hiseq 2500平臺(tái)和50 bp上測(cè)序。進(jìn)行已知的miRNA比對(duì)檢測(cè)，使用miRBase 20.0作為參考，使用改良軟件mirdeep2和srna-tools-cli獲取潛在的miRNA并繪制二級(jí)結(jié)構(gòu)；并運(yùn)用miRNA前體發(fā)夾結(jié)構(gòu)的特征來(lái)預(yù)測(cè)新型miRNA，將軟件miREvo和mirdeep2集成在一起，通過(guò)探索二級(jí)結(jié)構(gòu)、Dicer切割位點(diǎn)和未注釋的小RNA標(biāo)簽的最小自由能來(lái)預(yù)測(cè)新型miRNA。預(yù)測(cè)miRNA的靶基因是在miRanda軟件中進(jìn)行。對(duì)于SO組與OP組具有生物學(xué)重復(fù)的骨組織樣品進(jìn)行篩選miRNA的差異表達(dá)譜，使用DESeq R軟件(1.8.3)對(duì)這兩組進(jìn)行差異表達(dá)分析，并使用Benjamini&Hochberg方法將校正的P值0.05設(shè)置為顯著差異表達(dá)的閾值。對(duì)差異表達(dá)的miRNA的目標(biāo)基因候選物進(jìn)行基因本體論(GO)富集分析，并對(duì)目標(biāo)基因候選物進(jìn)行KEGG富集分析(http:// www. genome. jp/kegg/)，旨在從分子水平理解生物系統(tǒng)的高級(jí)功能和信號(hào)通路信息[9]。

1.3 統(tǒng)計(jì)學(xué)處理

2 結(jié)果

2.1 SO組與OP組比較組合的miRNA差異基因表達(dá)水平分析

利用DESeq軟件對(duì)差異表達(dá)基因進(jìn)行篩選，統(tǒng)計(jì)SO組與OP組各樣本中已知的和新篩選得到miRNA的表達(dá)量，并對(duì)兩組的骨組織樣本間表達(dá)量相關(guān)性進(jìn)行分析。SO組與OP組比較組合的差異基因統(tǒng)計(jì)結(jié)果顯示，共篩選到3 064(包括上調(diào)和下調(diào))個(gè)差異基因，其中上調(diào)的有1 432個(gè)，下調(diào)的有1 632，設(shè)置Pvalue <0.05且 ∣log2 Fold Change∣>0.0。①火山圖(volcano plot)分析結(jié)果：通過(guò)火山圖的分析可以更直觀的觀察到SO組與OP組之間的差異基因分布情況，橫坐標(biāo)X軸表示log2 Fold Change值，縱坐標(biāo)Y軸則表示為-log10(Pvalue)，綠色點(diǎn)表示有差異且下調(diào)miRNA，紅色的點(diǎn)則表示為有差異且上調(diào)miRNA。差異基因火山圖結(jié)果提示，肌少-骨質(zhì)疏松癥骨組織存在差異顯著表達(dá)的miRNA(圖1)。②差異miRNA層次聚類(lèi)分析分析結(jié)果：通過(guò)聚類(lèi)分析將表達(dá)模式相近的miRNA聚集成類(lèi)，進(jìn)而可以識(shí)別未知miRNA的功能或已知miRNA的未知功能；此外，通過(guò)聚類(lèi)分析還可以更好的認(rèn)識(shí)這些同類(lèi)的miRNA共同參與同一代謝過(guò)程或細(xì)胞信號(hào)通路。層次聚類(lèi)分析橫坐標(biāo)為樣品名，縱坐標(biāo)為差異miRNA的每百萬(wàn)個(gè)映射讀取的每千克堿基中的碎片歸一化后的數(shù)值，分析結(jié)果為SO組與OP組組間的層次聚類(lèi)分析(圖2)。③GO富集結(jié)果分析：對(duì)兩組差異表達(dá) miRNA的集合進(jìn)行g(shù)ene ontology富集分析。其中對(duì)差異顯著的miRNA進(jìn)行GO富集分析，分析結(jié)果提示差異miRNA不僅在參與生物學(xué)過(guò)程(biological process，BP)、還參與了細(xì)胞組分(cellular component，CC)和分子功能(molecular function，MF)調(diào)控機(jī)體生命功能的過(guò)程。研究miRNA差異基因在 gene ontology 中的分布狀況將闡明實(shí)驗(yàn)中SO組與OP組骨組織樣本差異在基因功能上的體現(xiàn)。在GO富集分析結(jié)果中，將選取最顯著的30個(gè)Term(keg-20)繪制柱狀圖(圖3)、散點(diǎn)圖(圖4)。根據(jù)分析結(jié)果得到差異表達(dá)miRNA，預(yù)測(cè)其靶向的mRNA，并對(duì)差異表達(dá)miRNA的靶基因進(jìn)行KEGG富集功能分析(圖5)，有助于構(gòu)建mRNA與miRNA調(diào)控網(wǎng)絡(luò)，為深入研究miRNA對(duì)肌少-骨質(zhì)疏松癥的調(diào)控機(jī)制以及后續(xù)研究何種治療方案對(duì)肌少-骨質(zhì)疏松癥患者機(jī)體細(xì)胞生物功能產(chǎn)生何種影響提供更多研究思路。

圖1 橫坐標(biāo)為log2 Fold Change值，縱坐標(biāo)為-log10(P value)Fig.1 The abscissa is log2FoldChange value, the ordinate is -log10 (P value)

圖3 橫坐標(biāo)是GO Term，縱坐標(biāo)則是GO Term富集的顯著性水平Fig.3 The abscissa is GO Term, and the ordinate is the significance level of GO Term enrichment

圖4 橫坐標(biāo)為注釋到GO Term上的差異基因數(shù)與差異基因總數(shù)的比值，縱坐標(biāo)為GO TermFig.4 The abscissa is the ratio of the number of differential genes annotated to GO Term to the total number of differential genes, and the ordinate is GO Term

圖5 KEGG富集散點(diǎn)圖Fig.5 KEGG enrichment scatter plot注：Gene Ratio 這一列表示是差異基因集；Padj 多重假設(shè)檢驗(yàn)校正后的P值；Count注釋到KEGG通路編號(hào)上的差異基因數(shù)。

2.2 qRT-PCR驗(yàn)證結(jié)果

內(nèi)參基因U6和目的基因miR-382-3p溶解曲線呈單峰狀曲線，表明樣品無(wú)污染并且設(shè)計(jì)的引物特異性足夠穩(wěn)定，實(shí)驗(yàn)數(shù)據(jù)與預(yù)期結(jié)果較為吻合，miR-382-3p的實(shí)驗(yàn)組與對(duì)照組相對(duì)表達(dá)量見(jiàn)圖6所示。

圖6 qRT-PCR驗(yàn)證miR-382-3p 的表達(dá)情況與測(cè)序結(jié)果一致Fig.6 qRT-PCR verifies that the expression of miR-382-3p is consistent with the sequencing results

3 討論

隨著全球老齡化進(jìn)程的加速，世界各國(guó)面臨老年人群的健康問(wèn)題所帶來(lái)的壓力越來(lái)越大，其中SO的發(fā)病率逐年增加，引起了各國(guó)醫(yī)療衛(wèi)生工作者的廣泛關(guān)注，SO隨著年齡的增加，患者往往出現(xiàn)全身肌肉的減少和骨密度的下降，并且伴隨而來(lái)的跌倒和骨折風(fēng)險(xiǎn)顯著增加。

近年來(lái)，隨著轉(zhuǎn)錄調(diào)控測(cè)序技術(shù)的發(fā)展與廣泛運(yùn)用，特別是small RNA測(cè)序技術(shù)在采用先進(jìn)的Illumina測(cè)序平臺(tái)后，可針對(duì)各種疾病的樣本中的miRNA、siRNA、piRNA展開(kāi)全面分析，既能鑒定已知sRNA、也能預(yù)測(cè)新的sRNA并預(yù)測(cè)sRNA的靶基因，為深入研究miRNA在SO疾病過(guò)程中的功能及調(diào)控機(jī)制提供有力手段[10]。

miRNAs屬于一類(lèi)內(nèi)源性的small RNA，其可以通過(guò)降解mRNA或者抑制翻譯過(guò)程來(lái)抑制靶基因的表達(dá)[11]。在肌少-骨質(zhì)疏松癥的研究中miRNA發(fā)揮的作用在很大程度上還是未知的。本課題組利用高通量測(cè)序手段研究SO中起關(guān)鍵調(diào)節(jié)作用的 miRNA及其靶基因，有望從表觀遺傳學(xué)層面揭示肌肉減少癥患者容易罹患OP的新機(jī)制。miRNA在骨病、骨性關(guān)節(jié)炎、肌少癥和OP等老年性疾病的組織樣品中存在差異表達(dá)，并且miRNA在老年性疾病中的調(diào)控作用已有廣泛研究[12]。根據(jù)Illumina高通量測(cè)序結(jié)果篩選出可能參與調(diào)節(jié)SO疾病發(fā)生演變過(guò)程的候選miRNA，并充分利用GO和KEGG數(shù)據(jù)庫(kù)的生物信息可以評(píng)估顯著差異的miRNA在SO發(fā)生、發(fā)展中所發(fā)揮的調(diào)控作用。

本課題組在研究中利用Illumina高通量測(cè)序技術(shù)對(duì)SO組與對(duì)照組OP進(jìn)行對(duì)比，差異表達(dá)顯著的miRNA有22個(gè)，其中15個(gè)發(fā)生了上調(diào)(表1)，7個(gè)(表2)發(fā)生下調(diào)；與OP對(duì)照組相比，在SO組骨組織中存在顯著高表達(dá)的miRNA，其中hsa-miR-382-3p、hsa-miR-27a-5p、hsa-miR-1226-5p和hsa-miR-3934-5p是發(fā)生上調(diào)較顯著的miRNA，顯著低表達(dá)的是hsa-miR-451a。研究結(jié)果提示在SO患者骨組織中存在部分差異表達(dá)的miRNA，其可能參與調(diào)控肌少癥罹患OP的疾病過(guò)程發(fā)生發(fā)展過(guò)程。

表1 顯著上調(diào)的miRNATable 1 Significantly up-regulated miRNAs

表2 顯著下調(diào)的miRNATable 2 Significantly down-regulated miRNAs

本研究中發(fā)現(xiàn)的特征性候選基因hsa-miR-382-3p在膝骨性關(guān)節(jié)炎存在差異表達(dá)，目前國(guó)內(nèi)外研究發(fā)現(xiàn)miR-382-3p可能通過(guò)直接靶向CX43通過(guò)TLR4/MyD88/NF-κB信號(hào)通路參與OA[13]；Heilmeier等[14]在糖尿病性骨病和絕經(jīng)后骨質(zhì)疏松癥研究中發(fā)現(xiàn)miR-382-3p的差異表達(dá)對(duì)機(jī)體成骨、脂肪形成和細(xì)胞增殖的存在調(diào)控作用。本研究中，肌少-骨質(zhì)疏松癥患者的miR-27a-5p表達(dá)水平較顯著上調(diào)，有研究表明miR-27a-5p通過(guò)抑制Atg7可減輕缺氧誘導(dǎo)的大鼠心肌細(xì)胞損傷[15]，其也可能在肌少-骨質(zhì)疏松癥的病理生理過(guò)程中發(fā)揮重要作用。hsa-miR-451a在肌少-骨質(zhì)疏松癥組中表達(dá)顯著下調(diào)，Lu等[16]研究發(fā)現(xiàn)miR-451a在骨質(zhì)疏松癥中存在差異表達(dá)，其可能通過(guò)BMP6信號(hào)傳導(dǎo)促進(jìn)成骨分化并抑制骨量丟失進(jìn)程；Karvande等[17]研究報(bào)道甲狀旁腺素促進(jìn)葡萄糖依賴(lài)性miR-451a表達(dá)，刺激成骨細(xì)胞分化，其機(jī)制是甲狀旁腺素可以通過(guò)PI3K-mTOR-AKT軸抑制AMPK磷酸化，從而防止作用于miR-451a啟動(dòng)子區(qū)域的八聚體結(jié)合轉(zhuǎn)錄因子1(OCT-1)的磷酸化和失活，AMPK活性的調(diào)節(jié)會(huì)影響成骨細(xì)胞中miR-451a的水平。

本課題組在后續(xù)基礎(chǔ)研究中將篩選出可能調(diào)控肌少-骨質(zhì)疏松癥的候選miRNA，并通過(guò)預(yù)測(cè)分析以確定候選miRNA所作用的靶基因，在細(xì)胞與動(dòng)物水平上深入研究候選的miRNA與靶基因的表達(dá)水平及其調(diào)控機(jī)制。