胰腺癌中lncRNA和mRNA表达谱分析及调控网络研究

张美晨，魏文竹，卢星宇，汤文瀛，王迪，刘跃娟

（哈尔滨医科大学大庆校区医学信息学院，黑龙江 大庆 163319）

0 引言

胰腺癌（pancreaticcancer，PC）因其起病隐匿、早期诊断难、手术切除困难、放化疗不敏感以及预后差、死亡率高等原因，成为危害最大的恶性肿瘤之一。我国胰腺癌发病率近20年来增长了将近6倍，居恶性肿瘤死亡率的第5位，因此积极寻找对胰腺癌有更好诊断意义的临床指标尤为重要。目前，越来越多研究表明，癌症的发生与基因突变、缺失和表达障碍有关。本研究利用胰腺癌芯片数据，计算得到差异表达的mRNA和lncRNA，结合蛋白质互作数据，构建失调的lncRNA-mRNA网络，通过挖掘网络功能模块，初步探索其作用机制，为揭示胰腺癌发生发展过程中lncRNA与其靶基因的关系和作用机制提供新线索。

1 材料与方法

1.1 胰腺癌的lncRNA和mRNA表达谱

我们从基因表达数据库(Gene Expression Omnibus，GEO)中下载了编号为GSE61166的基因芯片数据。本研究从中选取4个胰腺癌患者和4个胰腺炎患者胰腺组织的lncRNA和mRNA表达谱数据作为研究对象。

1.2 BLAST序列比对

利用BLAST进行序列比对，在基因芯片GPL61166平台上下载序列信息，同时在GENECODE数据库下载人类基因序列，利用计算机软件比对，通过BioMart数据库，将Ensembl Gene ID转换为Gene Symbol。然后我们做了标准化和以2为底的对数转换，若多个探针集对应一个基因，则取其平均值作为该基因的表达量，共获得6530个mRNA和31个lncRNA，全部数据处理流程如图1所示。

1.3 差异基因筛选

差异基因通过基因芯片显著性分析方法（SAM）筛选，SAM方法根据芯片数据噪音大小与表达峰度进行修改，其优点是在得到较多特征基因的同时，错误发现率（FDR）还保持在较低水平。利用SAM方法，阀值设定为P＜0.05，得到差异表达的820个mRNA和16个lncRNA。

1.4 相关分析

图1 流程图

相关分析法是对具体有依存关系的现象探讨其相关方向以及程度的一种统计方法。本研究对差异表达的820个mRNA和16个lncRNA计算皮尔森相关系数，取系数大于0.8且P＜0.05的强相关对构建lncRNA-mRNA共表达网络。

1.5 构建lncRNA-mRNA互作网络

根据上一步得到的共表达网络，从HPRD（Human Protein Reference Database）数据库下载人类蛋白质和蛋白质相互作用（PPI）数据，得到人的PPI网络，并将lncRNA-mRNA共表达网络和PPI网络取交集，然后对交集的mRNAs在PPI中选取直接邻居结点，和原共表达网络一起形成lncRNA-mRNA互作网络。

1.6 网络模块挖掘

利用Cytoscape软件中的JActiveModules插件进行模块挖掘，获得5个网络功能模块。用ClusterOne和Netmine对模块进行聚类和进一步挖掘，寻找p-value＜0.05的显著性lncRNA-mRNA互作对。

1.7 Gene ontology(GO)分析

通过来自R/bioConductor的clusterProfiler包完成了对模块中编码基因的GO terms富集度统计学分析，计算出这些基因的GO term的p-value和p-value的FDR值，定位这些基因最可能相关的GO term。

2 结果

2.1 差异表达 lncRNAs、mRNAs

对胰腺癌表达谱数据进行处理，得到差异表达的mRNAs、lncRNAs数量分别为820、16个。通过比较差异mRNAs的前5%和全部的lncRNAs在不同样本中的表达量（见图2、图3），我们发现这些mRNAs与lncRNAs在胰腺癌患者中低表达。

图2 lncRNA和mRNA表达

2.2 网络的全局特征

构建的lncRNA—mRNA共表达网络共包含了826个结点和7176条边，涉及到16个lncRNA和810个mRNA，利用Cytoscape可视化网络（如图4）对共表达网络和PPI网络整合并得到一个包含2373个结点和8626条边，涉及到16个lncRNA和2357个mRNA的lncRNA—mRNA互作网络。通过拓扑属性分析(如图5)，我们发现MIR7-3HG具有较大的度585，平均邻居结点12，较小的聚类系数0.447，大多数的路径长度在2-5之间。

图3 lncRNA和mRNA直方图

图4 lncRNA-mRNA共表达网络

图5 网络度分布

2.3 网络功能模块

对lncRNA—mRNA互作网络进行模块挖掘，得到5个模块。我们发现LINC00670与CCL2关系对在四个模块中都出现，LINC00670是一个存在转录因子CTCF绑定位点的长链非编码RNA，CTCF蛋白在印记调控区域和分化甲基化区域1和MAR3结合可抑制胰岛素样生长因子2基因[14]。另外，我们对与LINC00670有互作关系的靶基因做GO功能富集分析发现，其靶基因以最显著的形式富集在了胰液分泌通路上。

在模块一中，LINC00670的介数达到0.021,远高于平均介数值0.008。CCL2趋化因子及其受体可加强肿瘤细胞的增殖和侵袭能力，并聚集免疫细胞等机制以促进肿瘤的发生发展，CCL2-CCR2信号通道在肿瘤细胞的转移及增殖过程中发挥着重要的作用，尤其是在肿瘤细胞转移的早晚期[15]。LINC00670同时还与IL6和NFASC相连。其中，IL6（白细胞介素-6），是炎症介质网络的关键成分。有研究表明IL6与癌细胞增殖以及许多胰腺癌病人的恶病质过程有关。在模块二中受体活性修饰蛋白1（RAMP1）可以与降钙素受体（CTR）直接作用产生胰淀粉样酶(AMY)受体，在胰腺肿瘤引起的胰腺导管阻塞、胰腺脓肿、胰腺损伤中AMY高表达。有研究表明ADM能够增强胰腺癌细胞的侵袭性。与LINC00671相联的IAPP在II型糖尿病中会引起β细胞凋亡。用一种具有中和抗体的细胞特异性促进剂明显减少了病变的活性巨噬细胞积累，从而减少纤维化和病灶的生长。在模块四中CKM,DCX,IL6,IAPP都在模块一、二、三中出现，CKM,DCX,IL6同时连接LINC00670。有研究表明在胰腺癌中从UCM中更快地释放出DOX，在凋亡的外围细胞中可能会停止内吞作用和胞吐作用，有研究表明组织和血液的MMP9增加反映了胰腺癌与糖尿病相关。

2.4 模块的功能分析

对模块中编码基因进行GO功能富集分析结果见图6。通过来自R/bioConductor的clusterProfiler数据包，我们完成了模块中编码基因的GO term富集度统计学分析，阀值选取为Bonferroni校正P值 ＜0.01。

图6第一排从左至右分别是模块一、三、五，第二排从左至右分别是模块二和四。得到的五个模块中都有与炎症显著性功能相关的注释；模块1中基因功能注释到炎症、增殖、凋亡；几个模块里还有关于侵袭转移的功能注释。其中多个功能有重复出现。如:taxis、regulation of developmental process、multi-organism cellular process。 通 过GO注释的模块功能与免疫反应有关。

3 结论

通过本项目的初步研究，找到与胰腺癌发生密切相关的的五个模块，而这些模块又分别与炎症、免疫、细胞增殖和凋亡等有关。表明mRNA和lncRNA通过多种途径参与胰腺癌的发生及发展过程，在多个模块中发现lncRNA00670与胰腺癌密切相关。随着mRNA和lncRNA研究技术的不断发展和成熟，利用生物信息学分析将会发掘更多的功能性mRNA和lncRNA，为诊断和治疗胰腺癌等恶性病疾病提供新方向。

图6 GO功能富集注释

[1] Liu R C, Traverso L W. Diagnostic laparoscopy improves staging of pancreatic cancer deemed locally unresectable by computed tomography[J]. Surgical Endoscopy And Other Interventional Techniques, 2005, 19(5): 638-642.

[2] Lin Q J, Yang F, Jin C, et al. Current status and progress of pancreatic cancer in China[J]. World Journal of Gastroenterology:WJG, 2015, 21(26): 7988.

[3] Li N, Yang M, Shi K, et al. Long non-coding RNA HOXA11-AS in human cancer: A meta -analysis[J]. Clinica Chimica Acta,2017,474:165-170.

[4] Bister V, Skoog T, Virolainen S, et al. Increased expression of matrix metalloproteinases-21 and-26 and TIMP-4 in pancreatic adenocarcinoma[J]. Modern Pathology, 2007, 20(11): 1128 -1140.

[5] Kayed H, Jiang X, Keleg S, et al. Regulation and functional role of the Runt-related transcription factor-2 in pancreatic cancer[J].British journal of cancer, 2007, 97(8): 1106-1115.

[6] Namkoong H, Shin S M, Kim H K, et al. The bone morphogenetic protein antagonist gremlin 1 is overexpressed in human cancers and interacts with YWHAH protein[J]. BMC cancer, 2006, 6(1): 74.

[7] Zhu Y, Xu H, Chen H, et al. Proteomic analysis of solid pseudopapillary tumor of the pancreas reveals dysfunction of the endoplasmic reticulum protein processing pathway[J]. Molecular &Cellular Proteomics, 2014, 13(10): 2593-2603.

[8] Wang Y, Yu S, Huang D, et al. Cellular Prion Protein Mediates Pancreatic Cancer Cell Survival and Invasion through Association with and Enhanced Signaling of Notch1[J]. The American journal of pathology, 2016, 186(11): 2945-2956.

[9] Qin C F, Zhao F L. Long non-coding RNA TUG1 can promote proliferation and migration of pancreatic cancer via EMT pathway[J].European Review for Medical and Pharmacological Sciences, 2017,21: 2377-2384.

[10] Tahira A C, Dazzani B, Reis E M, et al. Long noncoding intronic RNAs are differentially expressed in primary and metastatic pancreatic cancer[J]. Molecular cancer, 2011, 10(1): 141.

[11] Ma C, Nong K, Zhu H, et al. H19 promotes pancreatic cancer metastasis by derepressing let-7’s suppression on its target HMGA2-mediated EMT[J]. Tumor Biology, 2014, 35(9): 9163-9169.

[12] Saung M T, Zheng L. Current standards of chemotherapy for pancreatic cancer[J]. Clinical therapeutics, 2017.

[13] Yi X, Zai H, Long X, et al. Krüppel-like factor 8 induces epithelialto-mesenchymal transition and promotes invasion of pancreatic cancer cells through transcriptional activation of four and a half LIM-only protein 2[J]. Oncology letters, 2017, 14(4): 4883-4889.

[14] Huang X, Ta N, Zhang Y, et al. Microarray Analysis of the Expression Profile of Long Non-Coding RNAs Indicates lncRNA RP11-263F15. 1 as a Biomarker for Diagnosis and Prognostic Prediction of Pancreatic Ductal Adenocarcinoma[J]. Journal of Cancer, 2017, 8(14): 2740.

[15] Li H, Wang X, Fang Y, et al. Integrated expression profiles analysis reveals novel predictive biomarker in pancreatic ductal adenocarcinoma[J]. Oncotarget, 2017, 8(32): 52571.