轉(zhuǎn)錄組分析鑒定腰椎椎間盤退行性變中潛在免疫相關(guān)生物標(biāo)志物

李志超，王 磊，薛景才，王文波，陳 煜，齊軍強(qiáng)，許 鵬

1.山東中醫(yī)藥大學(xué)第二附屬醫(yī)院骨科，濟(jì)南 250002

2.山東中醫(yī)藥大學(xué)研究生院，濟(jì)南 250001

3.海軍軍醫(yī)大學(xué)研究生院，上海 200433

4.海軍軍醫(yī)大學(xué)長征醫(yī)院骨科，上海 200003

以椎間盤退行性變（IDD）為特征的椎間盤疾病，給患者帶來了生理和心理上的痛苦，給家庭和社會造成巨大的經(jīng)濟(jì)負(fù)擔(dān)［1］。IDD 的致病因素包括營養(yǎng)物質(zhì)減少、基質(zhì)降解、細(xì)胞凋亡、炎性因子增加、免疫應(yīng)答、遺傳因素和機(jī)械負(fù)荷等［2］，但其具體機(jī)制目前仍不清楚。IDD 與免疫之間有很強(qiáng)的相關(guān)性［3］，椎間盤的駐留細(xì)胞（纖維環(huán)和髓核）和浸潤的免疫細(xì)胞分泌的炎性細(xì)胞因子是引發(fā)IDD 的關(guān)鍵因素［4］。發(fā)生退行性變的椎間盤使巨噬細(xì)胞傾向于促炎性極化，這似乎加劇了IDD［5］。

競爭性內(nèi)源RNA（ceRNA）是一種新的基因表達(dá)調(diào)控模式。在調(diào)控過程中，一些長鏈非編碼RNA（lncRNA）和微RNA（miRNA）結(jié)合，lncRNA 起到miRNA 海綿的作用，降低miRNA 對其靶基因的抑制作用，使靶基因在細(xì)胞中的表達(dá)量增加［6］。有研究［7-8］表明，IDD 的病理變化與lncRNA、miRNA 及其靶點(diǎn)的失衡以及許多不同細(xì)胞生物過程有關(guān)。失衡的ceRNA 相互作用軸可能是治療IDD 的關(guān)鍵目標(biāo)［9］。本研究通過加權(quán)基因共表達(dá)分析（WGCNA）構(gòu) 建IDD 中l(wèi)ncRNA-miRNA-mRNA 的ceRNA 網(wǎng) 絡(luò)，篩選與IDD 免疫相關(guān)的關(guān)鍵基因，為尋找IDD 的新靶點(diǎn)和治療策略提供參考。

1 材料和方法

1.1 數(shù)據(jù)來源

從基因表達(dá)匯編（GEO）數(shù)據(jù)庫（https：//www.ncbi.nlm.nih.gov/geo）GSE67567 數(shù) 據(jù) 集 下 載IDD 人類退行性髓核組織和對照的髓核組織微陣列芯片數(shù)據(jù)，GSE67567 數(shù)據(jù)集包含GSE56081（lncRNA 和mRNA 數(shù)據(jù)）和GSE63492（miRNA）2 個(gè)子集，各包括來自5名健康志愿者和5例IDD患者的髓核樣本。

1.2 數(shù)據(jù)預(yù)處理和差異表達(dá)分析

使用GEO 的在線分析工具GEO2R 獲得標(biāo)準(zhǔn)化的差異表達(dá)基因。根據(jù)GPL15314、GPL19449文件和GEO中的注釋信息，將探針轉(zhuǎn)換為相應(yīng)的基因符號。對于差異表達(dá)的miRNA，將閾值設(shè)置為|log2FC| ＞ 2、調(diào)整P＜ 0.05；對于差異表達(dá)的lncRNA，將閾值設(shè)置為|log2FC| ＞ 1、調(diào)整P＜ 0.05。

1.3 共表達(dá)模塊建立和基因集富集分析（GSEA）

使用Image GP工具（https：//www.ehbio.com/Cloud Platform/front/#/）分 析GEO 中的mRNA 表達(dá)矩陣。首先構(gòu)建差異表達(dá)基因的共表達(dá)網(wǎng)絡(luò)，將其轉(zhuǎn)換為拓?fù)渲丿B矩陣，通過分層聚類和動態(tài)樹切割算法，使用模塊特征基因和模塊關(guān)系來確定與表型相關(guān)的共表達(dá)模塊，最小模塊大小為30［10］。使用Enrichr數(shù)據(jù)庫（https：//maayanlab.cloud/Enrichr/）分析每個(gè)模塊中的mRNA，進(jìn)行京都基因與基因組百科全書（KEGG）通路和基因本體（GO）生物過程的GSEA，篩選出與免疫調(diào)節(jié)關(guān)系最為密切的模塊。

1.4 ceRNA 網(wǎng)絡(luò)建立和GSEA

通過Targetscan 數(shù)據(jù)庫（https：//www.targetscan.org）和miRTarbase 數(shù)據(jù)庫（http：//mirtarbase.mbc.nctu.edu.tw/php/index.php）篩 選miRNA 靶 向 的mRNA，將這些mRNA 與上述所選模塊中的差異表達(dá)mRNA 取交集獲得免疫相關(guān)mRNA。同理，通過StarBase 數(shù)據(jù)庫（http：//starbase.sysu.edu.cn/index.php）預(yù)測差異表達(dá)miRNA 調(diào)控的lncRNA。最后，使用Cytoscape 3.7.2建立lncRNA-miRNA-mRNA ceRNA網(wǎng) 絡(luò)。使 用Enrichr 對ceRNA 網(wǎng) 絡(luò) 中 的mRNA 的KEGG 通路和GO 生物過程進(jìn)行GSEA。

1.5 免疫浸潤分析

使用ImmuCellAI（http：//bioinfo.life.hust.edu.cn/web/ImmuCellAI/）［11］比較GEO 數(shù)據(jù)庫中IDD 組和對照組的免疫浸潤特征。ImmuCellAI 是一種基于轉(zhuǎn)錄組測序或芯片工具估計(jì)免疫細(xì)胞浸潤豐度數(shù)據(jù)的工具，可以判斷18 種T 細(xì)胞和6 種其他類型免疫細(xì)胞的比例。

1.6 與CD4 幼稚T 細(xì)胞（CD4 naive）相關(guān)的ceRNA子網(wǎng)絡(luò)的構(gòu)建

使用Pearson 相關(guān)分析確定ceRNA 網(wǎng)絡(luò)中CD4 naive 浸潤與mRNA 和lncRNA 表達(dá)的相關(guān)性，以相關(guān)系數(shù)絕對值＞ 0.8、P＜ 0.01 為條件篩選mRNA和lncRNA［12］。同法分析mRNA 和lncRNA 之間的相關(guān)性，篩選出關(guān)鍵lncRNA。利用篩選出的mRNA 和lncRNA，在之前的ceRNA網(wǎng)絡(luò)基礎(chǔ)上構(gòu)建與CD4 naive相關(guān)的ceRNA 子網(wǎng)絡(luò)。

1.7 確定關(guān)鍵靶點(diǎn)的通路

使用關(guān)鍵lncRNA 表達(dá)量的中位數(shù)對所有表達(dá)量進(jìn)行分組，并使用進(jìn)行GO 生物過程和KEGG 通路的GSEA。通過獲得的標(biāo)準(zhǔn)化富集評分（NES）及其相應(yīng)的P值（NOMP-val）篩選關(guān)鍵的生物過程和途徑。

1.8 驗(yàn)證

使用GEO 數(shù)據(jù)庫中的另一個(gè)IDD 數(shù)據(jù)集GSE124272（包括8 例IDD 患者和8 名健康志愿者），通過IDD 和對照組基因表達(dá)的比較，分析差異表達(dá)基因，對獲得的關(guān)鍵基因進(jìn)行驗(yàn)證。

2 結(jié) 果

2.1 IDD 中mRNA、lncRNA 和miRNA的差異表達(dá)

通過GEO2R 分析GSE56081 和GSE63492 數(shù)據(jù)集中的數(shù)據(jù)，經(jīng)過校正，得到15 220 個(gè)mRNA、370個(gè)miRNA 和623 個(gè)與IDD 相關(guān)的lncRNA。進(jìn)一步根據(jù)|log2FC| ＞ 2、調(diào)整P＜ 0.05 篩選出35 個(gè)差異表達(dá)miRNA，均為低表達(dá)；根據(jù)|log2FC| ＞ 1、調(diào)整P＜ 0.05篩選出87個(gè)差異表達(dá)lncRNA，其中45個(gè)為低表達(dá)，42 個(gè)為高表達(dá)（圖1）。

圖1 差異表達(dá)基因Fig.1 Differentially expressed genes

2.2 共表達(dá)模塊構(gòu)建及GSEA

通過WGCNA 完成IDD 和對照組的共表達(dá)模塊構(gòu)建，通過聚類分析和動態(tài)切樹算法，得到5 個(gè)不同顏色的模塊，并識別出與臨床特征相關(guān)的共表達(dá)模塊（圖2）。

圖2 WGCNA 結(jié)果Fig.2 Results of WGCNA

對每個(gè)模塊的KEGG 通路和GO 生物過程進(jìn)行GSEA，結(jié)果表明，與免疫調(diào)節(jié)最相關(guān)的是藍(lán)綠色（turquoise）模塊，總共包含1 808 個(gè)mRNA。藍(lán)綠色模塊的KEGG 通路主要包括剪接體、溶酶體、RNA降解、亨廷頓病、人T 細(xì)胞白血病病毒1 型感染、氧化磷酸化等，GO 生物學(xué)過程主要包括參與免疫反應(yīng)的中性粒細(xì)胞活化、纖溶酶原激活的調(diào)控、巨自噬調(diào)節(jié)、TORC1 信號的調(diào)節(jié)、中性粒細(xì)胞介導(dǎo)的免疫、中性粒細(xì)胞脫顆粒等（圖3）。

圖3 藍(lán)綠色模塊的KEGG通路和GO生物過程Fig.3 KEGG pathways and GO biological processes of turquoise module

2.3 ceRNA 網(wǎng)絡(luò)構(gòu)建和GSEA

使用Targetscan 和miRTarbase 數(shù)據(jù)庫預(yù)測受差異表達(dá)miRNA 調(diào)節(jié)的mRNA，Targetscan 數(shù)據(jù)庫中獲得2 617 個(gè)、miRTarbase 數(shù)據(jù)庫中獲得3 940 個(gè)mRNA，取交集后獲得160 個(gè)mRNA，再與免疫相關(guān)模塊中的mRNA 取交集得到25 個(gè)mRNA。使用StarBase 數(shù)據(jù)庫預(yù)測出差異表達(dá)miRNA 調(diào)節(jié)的797個(gè)lncRNA，將其與差異表達(dá)lncRNA 取交集后獲得8 個(gè)lncRNA。使 用Cytoscape 3.7.2 構(gòu) 建lncRNAmiRNA-mRNA ceRNA 網(wǎng)絡(luò)（圖4a），獲得ceRNA 網(wǎng)絡(luò)中的mRNA，并對其KEGG 通路和GO 生物過程進(jìn)行GSEA。結(jié)果表明，GO 生物過程主要包括T 細(xì)胞耐受誘導(dǎo)調(diào)控、細(xì)胞對粒細(xì)胞-巨噬細(xì)胞集落刺激因子刺激的反應(yīng)、CD4+α-β T細(xì)胞活化的正調(diào)節(jié)、細(xì)胞對轉(zhuǎn)化的反應(yīng)生長因子刺激的反應(yīng)、α-β T細(xì)胞分化的正調(diào)節(jié)等（圖4b），KEGG通路主要包括Fox O信號通路、細(xì)胞衰老、SNARE 在囊泡轉(zhuǎn)運(yùn)中的相互作用、MAPK 信號通路（圖4c）。

圖4 免疫相關(guān)ceRNA 網(wǎng)絡(luò)及GSEAFig.4 Immune-related ceRNA network and GSEA

2.4 與CD4 幼稚T 細(xì)胞相關(guān)的ceRNA 子網(wǎng)絡(luò)構(gòu)建及GSEA

使用ImmuCellAI 進(jìn)行免疫細(xì)胞浸潤分析，結(jié)果顯示，CD4 naive、CD8 幼稚T 細(xì)胞（CD8 naive）、衰竭T 細(xì)胞（Tex）、自然調(diào)節(jié)性T 細(xì)胞（nTreg）、誘導(dǎo)調(diào)節(jié)性T 細(xì)胞（iTreg）、輔助性T 細(xì)胞1（Th1）、輔助性T 細(xì)胞17（Th17）、濾泡輔助性T 細(xì)胞（Tfh）、中央記憶細(xì)胞（Tcm）、自然殺傷細(xì)胞（NKT）、黏膜相關(guān)不變T 細(xì)胞（MAIT）、樹突狀細(xì)胞（DC）、單核細(xì)胞（monocyte）、γ-σ T 細(xì)胞（Tgd）、CD4 T 細(xì)胞在IDD組與對照組間差異有統(tǒng)計(jì)學(xué)意義（P＜ 0.05，圖5a），其中以CD4 naive T 細(xì)胞差異最為顯著，IDD 組豐度顯著低于對照組（P＜ 0.05）。

圖5 免疫浸潤分析及l(fā)ncRNA 與mRNA 的Pearson 相關(guān)分析Fig.5 ImmuCellAI immune infiltration analysis and Pearson correlation analysis of lncRNAs and mRNAs

利 用ceRNA 網(wǎng) 絡(luò) 中mRNA 和lncRNA 的 表 達(dá)量篩選與CD4 naive T 細(xì)胞免疫浸潤相關(guān)的關(guān)鍵基因，經(jīng)Pearson 相關(guān)分析，發(fā)現(xiàn)8 個(gè)lncRNA 和25 個(gè)mRNA 與CD4 naive T細(xì)胞浸潤具有強(qiáng)相關(guān)性，其中1 個(gè)lncRNA（LINC00641）與10 個(gè)mRNA 呈強(qiáng) 正相關(guān)性（圖5b）。

在ceRNA網(wǎng)絡(luò)基礎(chǔ)上構(gòu)建一個(gè)與CD4 naive T細(xì)胞相關(guān)的ceRNA子網(wǎng)絡(luò)（圖6a），關(guān)鍵mRNA可能是UHMK1、ZFP36L2、ZCCHC3、ZBTB20。在GSE56081數(shù)據(jù)集的GSEA 中，使用LINC00641 的中位表達(dá)值進(jìn)行分組，明確關(guān)鍵因子的調(diào)控途徑。結(jié)果顯示，許多機(jī)制與免疫密切相關(guān)，這些相關(guān)的KEGG 通路包括嗅覺轉(zhuǎn)導(dǎo)、轉(zhuǎn)化生長因子信號通路、溶酶體、胞質(zhì)DNA 感應(yīng)通路、剪接體等（圖6b）；GO 生物過程主要包括負(fù)調(diào)控細(xì)胞基質(zhì)黏附、白細(xì)胞遷移、細(xì)胞對細(xì)胞因子刺激的反應(yīng)、炎性反應(yīng)、白細(xì)胞遷移的調(diào)節(jié)等（圖6c）。

圖6 與CD4 naive相關(guān)的ceRNA 子網(wǎng)絡(luò)及GSEAFig.6 ceRNA sub-network associated with CD4 naive and GSEA

2.6 驗(yàn)證結(jié)果

通過GEO2R 分析GSE124272 數(shù)據(jù)集的差異表達(dá)基因，結(jié)果顯示，LINC00641、UHMK1、ZFP36L2、ZCCHC3和ZBTB20這5個(gè)關(guān)鍵基因在IDD組和對照組中均存在差異表達(dá)（表1）。

表1 GSE124272數(shù)據(jù)集中關(guān)鍵基因的差異表達(dá)Tab.1 Differential expression of key genes in GSE124272 dataset

3 討 論

IDD 是臨床常見病，發(fā)生機(jī)制復(fù)雜［13］。一項(xiàng)薈萃分析［14］表明，突出的椎間盤會自發(fā)重吸收，即IDD 患者未經(jīng)手術(shù)治療或侵入性治療，脫出的髓核自發(fā)消失或明顯縮小，但現(xiàn)有文獻(xiàn)報(bào)道多為病例報(bào)告或回顧性研究，其機(jī)制可能與炎性反應(yīng)和免疫調(diào)節(jié)有關(guān)［15］。因此，尋找新的免疫靶點(diǎn)和新的治療方法具有重要意義。近年來，基因表達(dá)譜的研究為越來越多的疾病篩選新型生物標(biāo)志物帶來了新思路與新方法［16］，本研究旨在探索IDD 與免疫細(xì)胞的相關(guān)性并篩選相關(guān)的特征基因。

本研究通過基因表達(dá)譜鑒定差異表達(dá)基因，使用WGCNA 和GSEA 篩選免疫相關(guān)模塊，并構(gòu)建了 由8 個(gè)lncRNA、21 個(gè)miRNA 和25 個(gè)mRNA 組成的ceRNA 網(wǎng)絡(luò)。mRNA 的GO 生物過程與KEGG通路分析表明其與免疫關(guān)系密切，如T 細(xì)胞耐受誘導(dǎo)的調(diào)節(jié)、CD4+的正向調(diào)節(jié)、T 細(xì)胞活化與分化的正向調(diào)節(jié)等生物過程，表明IDD 與T 細(xì)胞的調(diào)節(jié)有關(guān)，T 細(xì)胞參與IDD 的病理過程。而結(jié)果中的其他相關(guān)通路，如Fox O 信號通路、MAPK 信號通路與炎性反應(yīng)及氧化應(yīng)激關(guān)系密切［17-18］，鑒于炎性反應(yīng)是IDD 或椎間盤再生機(jī)制的中心，其平衡是維持組織穩(wěn)態(tài)的關(guān)鍵［19］，抑制炎性反應(yīng)可能有助于治療IDD［20］。

本研究使用免疫細(xì)胞浸潤分析的方法探索與IDD 具有強(qiáng)相關(guān)性的免疫細(xì)胞，結(jié)果發(fā)現(xiàn)，IDD患者與對照組的CD4 naive 表達(dá)差異性最為顯著。Naive T 細(xì)胞也稱非致敏T 細(xì)胞，在暴露于抗原刺激之前處于相對靜止的狀態(tài)。CD4 T 細(xì)胞是人體免疫系統(tǒng)中重要的免疫細(xì)胞［21］。最近的研究［22］表明，CD4 naive 的功能比以前認(rèn)為的要多，并且它們在表型、分化階段、持久性、功能和解剖位置方面是多種多樣的。動物實(shí)驗(yàn)［19］表明，椎間盤內(nèi)促炎治療可使血液和淋巴結(jié)中的T 細(xì)胞亞群增加。另一項(xiàng)研究［23］表明，在髓核突出的椎間盤中存在CD4 和CD8 T 細(xì)胞。在本研究中，CD4 naive 和CD4 T 細(xì)胞在免疫細(xì)胞浸潤分析后均存在顯著差異。Pearson 相關(guān)分析發(fā)現(xiàn)，UHMK1、ZFP36L2、ZCCHC3、ZBTB20和LINC00641 可能是IDD 免疫相關(guān)的關(guān)鍵靶點(diǎn)，即IDD 的特征基因。本研究還篩選了1 個(gè)與CD4 naive相關(guān)的ceRNA 子網(wǎng)絡(luò)，并將LINC00641 識別為關(guān)鍵因素，可將其作為IDD的潛在治療靶點(diǎn)。GSEA結(jié)果顯示，轉(zhuǎn)化生長因子信號通路可能是IDD 潛在的信號通路。上述關(guān)鍵基因的表達(dá)在GSE124272 數(shù)據(jù)集中得到了驗(yàn)證。

目前大多數(shù)與LINC00641相關(guān)的研究都與惡性腫瘤有關(guān)［24］，也與糖尿病足的一些生物學(xué)過程有關(guān)［25］。在IDD 的研究中，Wang 等［26］發(fā)現(xiàn)LINC00641 可以直接與miR-153-3p 結(jié)合，作為miR-153-3p 的天然海綿，通過miRNA介導(dǎo)的沉默途徑靶向自噬相關(guān)基因5，并調(diào)節(jié)髓核細(xì)胞的自噬和死亡。本研究結(jié)果表明，LINC00641 可能通過免疫相關(guān)競爭性調(diào)節(jié)基因機(jī)制影響IDD 的發(fā)生與發(fā)展。

本研究發(fā)現(xiàn)的與CD4 naive 相關(guān)的mRNA 目前研究尚少，這些mRNA與腫瘤和免疫密切相關(guān)［27-29］。UHMK1 可能參與RNA 加工，也可能參與白細(xì)胞的形成［30］。ZFP36L2 在CD4 naive 中高表達(dá)，可抑制調(diào)節(jié)性T 細(xì)胞的功能［31］。本研究的免疫浸潤分析結(jié)果顯示，CD4 和CD4 naive 與IDD 密切相關(guān)，間接驗(yàn)證了ZFP36L2 是IDD 的重要靶點(diǎn)之一。同時(shí)，ZBTB20 可抑制CD8+T 細(xì)胞免疫代謝［32］，并能促進(jìn)巨噬細(xì)胞炎性反應(yīng)，啟動先天免疫反應(yīng)［33］。這些結(jié)果為IDD 的診斷和靶向治療提供了新思路。

本研究通過GSEA 來識別LINC00641 的關(guān)鍵調(diào)控途徑，發(fā)現(xiàn)與免疫相關(guān)的關(guān)鍵途徑是轉(zhuǎn)化生長因子信號通路、胞質(zhì)DNA 感應(yīng)通路、溶酶體通路等。胞質(zhì)DNA 感應(yīng)通路與先天免疫反應(yīng)密切相關(guān)［34］，溶酶體是真核細(xì)胞的主要分解代謝區(qū)［35］，許多炎性疾病的根本原因可能是自噬-溶酶體途徑的功能障礙［36］。最近有研究［37］表明，溶酶體與椎間盤細(xì)胞的死亡有關(guān)。近年來，有研究［38］發(fā)現(xiàn)轉(zhuǎn)化生長因子在調(diào)節(jié)細(xì)胞生長分化和免疫功能方面具有重要作用。轉(zhuǎn)化生長因子信號通路對椎間盤的發(fā)育和生長至關(guān)重要［39］，但過度激活轉(zhuǎn)化生長因子信號對椎間盤是有害的，抑制異常的轉(zhuǎn)化生長因子信號可以延緩椎間盤的退行性變［40］。因此，通過lncRNA 抑制轉(zhuǎn)化生長因子信號通路來減少炎性反應(yīng)可能是治療IDD 的有效方法。

本研究尚有一定的局限性，受現(xiàn)有數(shù)據(jù)庫中數(shù)據(jù)的限制，本研究選取的樣本量有限，以后將增加樣本量進(jìn)一步研究探索；此外，本研究組未開展細(xì)胞實(shí)驗(yàn)或動物實(shí)驗(yàn)，未來將通過基礎(chǔ)實(shí)驗(yàn)研究進(jìn)一步驗(yàn)證關(guān)鍵lncRNA 和mRNA 在IDD 中的作用。

總之，本研究在WGCNA 的基礎(chǔ)上，將ceRNA、IDD 和免疫調(diào)節(jié)聯(lián)系起來，構(gòu)建了一個(gè)與免疫調(diào)節(jié)相關(guān)的IDD ceRNA 網(wǎng)絡(luò)。在這個(gè)網(wǎng)絡(luò)中，lncRNA LINC00641 是一個(gè)關(guān)鍵靶點(diǎn)，它干預(yù)了mRNA（UHMK1、ZFP36L2、ZCCHC3、ZBTB20）的表達(dá)。同時(shí)，CD4 naive 表達(dá)減少，這可能與IDD 的發(fā)生機(jī)制密切相關(guān)。lncRNA LINC00641 可能通過轉(zhuǎn)化生長因子信號通路調(diào)節(jié)IDD。這些結(jié)果將有助于從ceRNA 的角度更好地探索免疫調(diào)節(jié)在IDD 發(fā)展中的作用機(jī)制，并為治療靶點(diǎn)的開發(fā)提供思路。