基于circRNA芯片数据的肺癌生物信息学分析*

翟西华，李霞，孙秋文，白津，马凯

221000江苏 徐州，徐州医科大学 医学信息与工程学院(翟西华、马凯)；221000江苏 徐州，徐州医科大学附属医院 肿瘤科(李霞)；221000江苏 徐州，徐州医科大学 医学影像学院(孙秋文)；221000江苏 徐州，徐州医科大学 肿瘤防治研究所(白津)

肺癌是对人类生命威胁最大的恶性肿瘤之一。肺癌分为非小细胞肺癌和小细胞肺癌两种[1]，其中非小细胞肺癌包括腺癌、大细胞癌和鳞状细胞癌。近年来，尽管在肺癌的诊断、分类和治疗方面的研究逐步有了新的进展，但肺癌患者的总体生存率仍然很低，晚期或转移性肺癌患者的预后更差[2]。肺癌已成为全球公共卫生的挑战。因此，进一步了解肺癌的分子机制，寻找新的治疗靶标对于肺癌的治疗非常重要[3]。

CircRNA是一类特殊的非编码RNA，具有闭合环状结构，类似于长链非编码RNA[4]，被认为是转录物异常剪接产生的环状结构中的RNA分子。越来越多的数据表明，circRNA在许多生物学过程中都可以发挥重要作用。例如：细胞增殖[5]、侵袭[6]和分化[7]。相关报道表明，circRNA在肺癌、结直肠癌及肝癌等多种肿瘤中具有重要功能。Zong等[8]研究表明，circRNA_102231在肺癌组织表达上调，并可能成为肺癌患者的潜在生物标志物和治疗靶标。Yao等[9]报道在非小细胞肺癌组织中circRNA_100876表达显著上调，并与淋巴结转移和肿瘤分期密切相关。Wan等[10]提出circRNA-ITCH通过抑制Wnt/β-catenin途径抑制肺癌细胞的增殖。本文利用生物信息学方法分析肺癌的circRNA表达谱芯片，对比肺癌组织与正常组织，寻找差异表达circRNA并预测其潜在靶基因，为circRNA在肺癌中的作用机制与临床治疗提供新方向。

1 材料与方法

1.1 获取肺癌芯片数据

从美国生物技术信息中心的基因表达数据库(Gene Expression Omnibus,GEO)[11]中下载肺癌组织的circRNA 芯片数据，GSE101684包含4例肺癌患者的组织样本与4例肺癌癌旁正常组织的样本。GSE112214包含3例肺癌患者的组织样本与3例肺癌癌旁正常组织的样本。共有7例肺癌患者的组织样本(癌症组)与7例肺癌癌旁正常组织的样本(正常组)。

1.2 芯片数据处理

利用Perl脚本文件将芯片矩阵数据中安捷伦(Agilent)私有的circRNA名称转换为标准的circRNA 标志号(identity document,ID)(本文转化为：hsa_circ_0039908)。接着，利用R软件读取芯片矩阵数据，对其中circRNA的表达值进行log2对数转换[12]。表达矩阵中相同的circRNA通过取平均表达值的方法合并数据。CircRNA探针信息均通过平台文件GPL21825(芯片数据集的平台文件)注释。

1.3 差异表达circRNA

利用R语言中“limma”算法包对上述处理后的芯片数据进行差异表达circRNA的筛选。首先，对circRNA芯片数据进行归一化处理[13]，并在二维坐标系中绘制火山图。癌症组和正常组之间的circRNA表达值通过t检验计算P值，BH法校正P值[14]。差异circRNA的筛选条件为差异倍数|log2FoldChange|>1且校正后的P值<0.05。根据火山图找到差异表达的circRNAs，并展示分析结果。同时，使用R语言中的“pheatmap”绘图方法绘制热图，展现筛选出的circRNA表达量的差异性。

1.4 miRNA靶基因预测

环状RNA数据库(CircBase)[15]统一了成套的circRNA，能够在基因组环境中访问、下载和浏览支持它们表达的证据。数据库中收录了人类、小鼠等多个物种的circRNA序列信息[16]。根据circRNA的ID可以在该数据库中得到该circRNA所在的基因名称和染色体的位置，circRNA所在的样品名以及fasta碱基序列，便于寻找与circRNA相结合的miRNA。

CircRNA可通过吸附miRNA，导致miRNA不能与其他线性RNA结合，其他线性RNA表达变高，从而间接调控基因的表达[17]。根据circBase数据库得到的circRNA所在染色体的位置信息，在CSCD在线分析软件中，可得到circRNA相对应部分结合的miRNA ID及结合位置。利用Perl脚本文件在软件预测的哺乳动物miRNA靶基因数据库(Search for Predicted MicroRNA Targets in Mammals,TargetScan)、实验验证的miRNA靶基因数据库(Experimentally Validated MicroRNA-target Interactions Database,miRTarBase)与miRNA靶基因数据库(MicroRNA Target Prediction Database,miRDB)中进一步筛选，得到与circRNA结合的miRNA靶基因。

1.5 靶基因富集分析

首先利用Perl脚本文件将上述得到的miRNA的靶基因名转化为靶基因ID，利用R软件中“Cluster Profiler”包对其靶基因分别进行基因本体论(gene ontology,GO)生物学功能富集分析与京都基因与基因组大百科全书(Kyoto Encyclopedia for Genes and Genomes,KEGG)信号通路富集分析。GO富集包括分子功能、生物过程与细胞组成3个类别[18]。P<0.05为富集有统计学意义。KEGG由系统信息、基因组信息、化学信息与健康信息四类数据库组成，利用KEGG数据库[19]对差异转录本进行通路分析(结合KEGG注释结果)，得到肺癌发展过程中关键的通路。

2 结 果

2.1 差异circRNA筛选

如图1所示，对GSE101684与GSE112214进行数据处理，聚类分析表明两组样本中circRNA的表达存在差异。对比癌症组与正常组筛选出350个差异circRNA，其中上调的circRNA有169个，下调的circRNA有181个，如图2差异表达circRNA火山图所示。此外，如表1所示，在差异circRNA上调与下调最显著的前10个circRNA中，hsa_circ_0039908上调最明显。

图1 GSE101684与GSE112214数据集中差异表达基因热图Figure 1. Differentially Expressed Genes in Dataset GSE101684 and GSE112214

图2 差异表达circRNA火山图Figure 2. Differentially Expressed CircRNAs in the Subjects

2.2 靶基因预测结果

对上述差异最显著的circRNA，即hsa_circ_0039908，进一步分析。首先，根据hsa_circ_0039908的ID信息在CSCD在线软件中确定与该circRNA结合的miRNA。经软件分析，共有35个miRNA与hsa_circ_0039908结合，包括：hsa-miR-136-5p、hsa-miR-3119与hsa-miR-375等。基于miRDB、miRTarBase与TargetScan三个数据库，利用Perl脚本，得到这35个miRNA的靶基因。表2所示为与hsa_circ_0039908结合的35个miRNA中的10个代表性miRNA及其靶基因。通过构建竞争内源性RNA调控网络(ceRNA regulation network，ceRNA)与进一步分析，得到了关键的几对“circRNA-miRNA-mRNA”信息，具体包括：“hsa_circ_0039908”-“hsa-miR-136-5p”-“SOCS7”，“hsa_circ_0039908”-“hsa-miR-3119”-“BTG2”与“hsa_circ_0039908”-“hsa-miR-375”-“RLF”(图3)。

图3 ceRNA调控网络图Figure 3. ceRNA Regulatory Network

2.3 靶基因GO富集分析

对与hsa_circ_0039908结合的miRNA靶基因进行GO富集分析。如图4靶基因GO富集分析图所示，靶基因主要参与肌肉组织发育、对类固醇激素的反应和细胞酰胺代谢过程的负调控等生物学过程。

图4 靶基因GO富集分析Figure 4. GO Enrichment Analysis of the Target GenesGO: Gene Ontology.

DUS2基因与肺癌发生有关，DUS2L基因由DUS2编码，是hsa_circ_0039908的来源基因。DUS2编码的蛋白可以通过抑制干扰素诱导的蛋白激酶来提高翻译速率。由GO富集结果分析图可知，hsa_circ_0039908的靶基因主要富集于细胞质应激颗粒和蛋白激酶复合物。从而进一步说明，hsa_circ_0039908可能影响蛋白激酶参与的生物学过程，调控肺癌的发生发展。

2.4 靶基因KEGG富集分析

对与hsa_circ_0039908结合的miRNA靶基因进行KEGG富集分析。如图5靶基因KEGG富集分析图所示，共得到7条KEGG统计结果，靶基因主要与调节干细胞多能性的信号通路、FoxO信号通路、AMPK信号通路和催乳素信号通路等有关。

图5 靶基因KEGG富集分析Figure 5. KEGG Enrichment Analysis of the Target GenesKEGG: Kyoto Encyclopedia of Genes and Genomes.

FoxO信号通路包含circ-FoxO3基因。实验表明，circ-FoxO3在肺癌细胞中高表达。P53充当诱导细胞周期停滞并促进由于DNA损伤引起的细胞凋亡修复的介质。Circ-FoxO3可以促进癌基因MDM2与p53之间的结合，促进MDM2诱导的p53泛素化和降解，抑制p53的表达，进而抑制肿瘤细胞的增殖。由此说明，hsa_circ_0039908可以通过影响FoxO信号通路，影响肺癌的进展。

3 讨 论

传统分子生物学方法对circRNA数量和丰度的检测能效有限。因此，一直以来circRNA被认为是异常剪切的副产物。近年随着RNA测序技术和生物信息分析的发展，通过高通量测序和生物信息学方法已经鉴定了大量的circRNA，同时发现circRNA能够间接调控基因的表达。随着基因芯片技术的发展，越来越多地将生物信息学方法用于寻找癌症治疗靶点和诊断标志物，以实现对恶性肿瘤的诊断与治疗。研究发现大多数circRNA包含有miRNA结合位点，能吸附miRNA从而调控miRNA的靶基因[20]。

本研究选取GEO数据库的肺癌芯片数据，包括7例肺癌患者样本与7例癌旁正常组织样本。通过比较2组circRNA表达水平，共找出350个差异表达的circRNA，其中上调的circRNA有169个，下调的circRNA有181个。结果显示，与癌旁组织相比，hsa_circ_0039908在肺癌组患者中上调最显著，提示该基因可能在肺癌中具有一定特异性。

利用miRDB、miRTarBase与TargetScan三个数据库对hsa_circ_0039908进行miRNA靶基因预测。CircRNA富含miRNA结合位点，能作为miRNA海绵，间接调控miRNA靶基因的表达。因此，为了研究hsa_circ_0039908的功能，利用“Cluster Profiler”软件包对hsa_circ_0039908 的miRNA靶基因进行GO与KEGG富集分析。结果表明，靶基因参与肌肉组织发育、对类固醇激素的反应与细胞酰胺代谢过程的负调控等生物学过程。神经酰胺[21]与肿瘤的发生发展有密切的关系。研究表明[22]，神经酰胺通过不同信号途径，能促进或者抑制细胞增殖。“细胞酰胺代谢过程的负调控”是hsa_circ_0039908的miRNA靶基因参与的生物学过程。因此，hsa_circ_0039908的表达失常可能与癌细胞代谢异常高度相关。KEGG富集分析结果表明，miRNA靶基因富集的信号通路有调节干细胞多能性的信号通路、FoxO信号通路与AMPK信号通路。AMPK信号通路可以调控癌细胞自噬过程,在肺癌发生发展过程中自噬有两种作用方式，其一，“自噬”能清除细胞内有毒物质，预防肺癌；其二，“自噬”能为压力状态下的肺癌细胞提供能量，促进肺癌细胞存活。因此推测hsa_circ_0039908的上调可能会间接激活AMPK信号通路从而影响癌细胞“自噬过程”。

综上所述，利用生物信息学方法确定可能参与肺癌发生发展的circRNA，通过GO与KEGG分析确定了细胞酰胺代谢过程的负调控、AMPK信号通路与调节干细胞多能性的信号通路可能参与肺癌的发生发展。当前，circRNA在癌症发展过程中的作用机制与生物学功能尚不明确。需要更多深层次的研究探索circRNA在肺癌中的功能机制。由于样本数据的限制，无法根据肺癌的类型，针对性地研究circRNA对鳞状上皮细胞癌、小细胞未分化癌、大细胞未分化癌与腺癌等的作用机制。所以本研究存在一定局限性，仍需融合更多类型数据，进一步研究circRNA在不同类型肺癌发生发展中的作用机制。

作者声明：本文全部作者对于研究和撰写的论文出现的不端行为承担相应责任；并承诺论文中涉及的原始图片、数据资料等已按照有关规定保存，可接受核查。

学术不端：本文在初审、返修及出版前均通过中国知网(CNKI)科技期刊学术不端文献检测系统的学术不端检测。

同行评议：经同行专家双盲外审，达到刊发要求。

利益冲突：所有作者均声明不存在利益冲突。

文章版权：本文出版前已与全体作者签署了论文授权书等协议。