GEO数据库分析环状RNA在乳腺癌组织中的表达

张鹤达, 李伟, 钟山亮, 姜林宏, 李建, 唐金海

根据2012年全球癌症统计显示,2012年估计有170万乳腺癌新发病例和52万的乳腺癌死亡人数,分别占女性癌症总发病率的25%和女性癌症死亡总数的15%,成为女性最常见的癌症,也是癌症导致死亡的主要原因[1]。随着医学的不断发展进步,乳腺癌的治疗手段也逐渐增加,但是仍不能彻底阻断乳腺癌的转移,因此我们需要迫切深入了解乳腺癌侵袭转移机制,为乳腺癌更有效地治疗提供新的方向。

环状RNA是一类新颖的单链、内源性RNA分子,通过与miRNA或其他分子结合,在转录或转录后水平介导基因的表达。环状RNA在多种肿瘤中异常调节,并参与重要的生物学功能的调控,包括细胞凋亡、增殖、迁移侵袭和肿瘤形成等[2-4]。因此环状RNA将成为新的肿瘤标志物。

乳腺癌常见的分子分型主要分为四类,分别是Luminal A型、Luminal B型、HER2阳性型和三阴性,其中Luminal A型乳腺癌占的比重较大,其次是三阴性乳腺癌,而且相对于Luminal A型乳腺癌,三阴性乳腺癌患者的转移复发率高,预后差,生存率低[5-6]。因此我们将通过分析GEO数据中Luminal A型和三阴性乳腺癌组织样本环状RNA的表达水平,进一步探索环状RNA在乳腺癌侵袭迁移中的作用,为乳腺癌治疗和转移预防提供新的前景。

1 材料和方法

1.1 微阵列数据 原始基因表达谱数据GSE101124是从GEO数据库(http://www.ncbi.nlm.nih.gov/geo)下载的。从GSE101124数据集中共获得4对Luminal A型乳腺癌组织和4对三阴性乳腺癌组织,并通过GPL 19978 Agraystar HumanCircRNA芯片V1平台进行分析,共分析了4 229个环状RNA。

1.2 GEO数据分析 差异表达基因是通过R软件中的limma软件包(版本3.4.4;http://www.r-project.org)确定,读取原始数据,对背景进行校正和归一化。与阴性对照相比,至少50%的基因芯片的探针的信号强度比阴性对照高10%。然后对微阵列权值进行估测,并将其应用到线性模型分析中,采用t检验方法鉴定Luminal A型和三阴性乳腺癌之间差异表达基因。采用Benjamini-Hochberg方法进行多次比较,以倍数差异≥2和调整后的P<0.05作为标准进行后续筛选。

1.3 环状RNA/miRNAs相互作用 利用miRanda和RNAhybrid 软件预测环状RNA与miRNA可能的结合位点。利用Cytoscape软件构建差异表达的环状RNA/miRNAs的相互作用网络。

1.4 基因在乳腺癌中的表达 通过来源于TCGA数据库的UALCAN在线软件分析环状RNA的亲本基因在乳腺癌组织中的表达,推测环状RNA与亲本基因的关系。

2 结果

2.1 筛选差异表达的环状RNA 通过对Luminal A型和三阴性乳腺癌组织进行差异分析,发现共有27个环状RNA差异表达,与Luminal A型相比,在三阴性乳腺癌组织中有6个环状RNA上调,21个环状RNA下调,见图1A,表1、表2。通过绘制火山图,更直接的显示环状RNA的差异倍数以及统计意义,见图1B。然后对差异的环状RNA进行聚类分析,见图1C,每行表示差异的环状RNA,每个垂直列对应于乳腺癌样本。

2.2 预测环状RNA与miRNA的相互作用 根据miRNA种子区碱基互补配对,通过miRanda 软件和RNAhybrid软件对环状RNA可能结合的miRNAs进行预测,构建环状RNA/miRNAs之间的联络网。通过对在三阴性乳腺癌组织中上调的6个环状RNA进行预测,发现hsa_circ_0001666吸附的miRNAs最多,达135个,其次是hsa_circ_0003611,有69个可结合的miRNAs,而且hsa_circ_0001666与其他环状RNA共同作用的miRNAs的数量也最多,其中与hsa_circ_0003511共同作用的miRNAs有12个(图2A)。对在三阴性乳腺癌组织中下调的前10个环状RNA进行预测,发现hsa_circ_0000965所结合的miRNAs最多,高达255个,其次是hsa_circ_0001907,可与211个miRNAs结合,而且hsa_circ_0000965和hsa_circ_0001907共同的作用miRNAs也最多,包含75个miRNAs。同时发现miR-5787可以被6个环状RNA共同作用(hsa_circ_0000965,hsa_circ_0000288,hsa_circ_0001907,hsa_circ_0005664,hsa_circ_0038645,hsa_circ_0020594);miR-1254,miR-4656,miR-4722-5p,miR-612,miR-6756-5p,miR-6803-5p,miR-6882-3p,miR-762可以被5个环状RNA共同作用。因此构建环状RNA/miRNAs网络将会更直接的观察其相互联系,为环状RNA通路机制的研究提供了重要的信息。

1A:Luminal A型和三阴性乳腺癌组织之间环状RNA微阵列信号值的散点图,与Luminal A型乳腺癌组织相比,红点代表上调的环状RNA,绿色点代表下调的环状RNA,黑点代表无明显意义的圆环RNA;1B:火山图,垂直线对应于差异倍数为2倍或1/2的分割线,水平线代表P值为0.05,图中的红点代表差异表达的环状RNA;1C:差异表达的环状RNA的热图,根据Pearson相关系数对环状RNA进行分类,颜色等级从绿色(低强度)到黑色(中等强度),再到红色(高强度)图1 筛选差异表达的环状RNA

表1 在三阴性乳腺癌组织中下调的环状RNA

环状RNA亲本基因P值倍数差异 hsa_circ_0000086ST6GALNAC30.0005914.158422 hsa_circ_0000965C19orf180.0045413.722668 hsa_circ_0001907LOC1005074120.0036373.51351 hsa_circ_0005664CAPZA10.0049072.875769 hsa_circ_0005516RGS100.0006442.651565 hsa_circ_0020594-0.0010832.500509 hsa_circ_0000288CAPRIN10.0282242.442410 hsa_circ_0000523METTL30.0041912.408066 hsa_circ_0033010LGMN0.0029222.400354 hsa_circ_0038645PRKCB0.0059162.364240 hsa_circ_0006220TADA2A0.0426382.347079 hsa_circ_0044234CDC270.0130022.325888 hsa_circ_0005221NPLOC40.0145892.293956 hsa_circ_0000849DYM0.0016412.291298 hsa_circ_0050205MEF2BNB-MEF2B0.0356912.262669 hsa_circ_0000977NOL100.0152082.158794 hsa_circ_0054254EML40.0363722.138124 hsa_circ_0058451ACSL30.0031392.102548 hsa_circ_0001421SEC31A0.0005522.097744 hsa_circ_0007380NSUN20.0162282.068703 hsa_circ_0008086CDKAL10.0434282.025271

表2 在三阴性乳腺癌中组织上调的环状RNA

环状RNA亲本基因P值倍数差异 hsa_circ_0078153STXBP50.0125742.01165 hsa_circ_0007762STXBP50.0159782.205221 hsa_circ_0001666FAM120B0.0013342.246033 hsa_circ_0003511LMBR10.0381662.502945 hsa_circ_0003611LPAR10.0441922.503140 hsa_circ_0091994FLNA0.0327722.931248

2.3 环状RNA与亲本基因的关系 通过UALCAN分析发现,与luminal型乳腺癌组织相比,C19orf18、METTL3、CDC27和ACSL3在三阴性乳腺癌组织中表达下调,见图3A。FLNA在三阴性乳腺癌组织中表达上调,见图3B,均与其相对应的环状RNA的表达呈正相关。

3 讨论

随着对环状RNA的不断研究,发现环状RNA不仅在多种肿瘤中异常调节,而且广泛参与肿瘤的发生发展[7],尤其在肿瘤的侵袭转移中起到了重要的作用。如hsa_circRNA_103809通过海绵miR-4302,上调了miR-4302的靶基因ZNF121的表达,从而促进肺癌细胞的侵袭[8];hsa_circ_0016788通过调控miR-486/CDK4途径促进肝癌细胞的侵袭[9];circHIPK3通过抑制miR-7促进结肠癌细胞的侵袭迁移[10]。肿瘤细胞的侵袭转移是癌症患者最主要的死亡原因,也是目前治疗的难点。因此本研究主要是通过对GEO数据库中4对Luminal A型乳腺癌组织和4对三阴性乳腺癌组织中的环状RNA进行相关分析,为环状RNA在乳腺癌侵袭转移的研究中提供新观念。通过分析发现与Luminal A型乳腺癌组织相比,在三阴性乳腺癌组织中共有6个环状RNA上调,21个环状RNA下调。然后我们对异常表达的环状RNA可结合的下游miRNAs进行预测,发现在上调的环状RNA中,hsa_circ_0001666吸附的miRNAs最多,而且hsa_circ_0001666均可以与hsa_circ_0078153、hsa_circ_0007762、hsa_circ_0003511、hsa_circ_0003611和hsa_circ_0091994共同调控部分miRNAs。同时在下调的环状RNA中,hsa_circ_0000965海绵的miRNAs最多,且同hsa_circ_0001907共同作用的miRNAs也最多。因此我们推测hsa_circ_0001666和hsa_circ_0000965具有较强的吸附作用。除此之外,我们还发现miR-5787可以被6个环状RNA共同作用,miR-1254,miR-4656,miR-4722-5p,miR-612,miR-6756-5p,miR-6803-5p,miR-6882-3p,miR-762可以被5个环状RNA共同作用,预示着这些miRNAs可能是侵袭力较强的因子,通过调控下游靶基因促进乳腺癌的迁移侵袭。而且通过文献查找发现miR-762可以通过抑制IRF7的表达促进乳腺癌细胞的增殖和侵袭[11];miR-1254通过靶向RASSF9促进乳腺癌细胞的增殖,抑制凋亡[12],然而miR-1254以及其他miRNAs在乳腺癌细胞中的迁移侵袭作用并未研究。因此在后续的研究中不仅需要双荧光素酶实验验证环状RNA与miRNA的相互作用,还需要证实miRNAs在乳腺癌中的侵袭迁移作用。

2A:预测三阴性乳腺癌组织中上调的6个环状RNA与miRNAs的相互作用;2B:预测三阴性乳腺癌组织中下调的前10个环状RNA与miRNAs的相互作用图2 预测环状RNA与miRNA的相互作用

3A:C19orf18、METTL3、CDC27和ACSL3在三阴性乳腺癌组织中表达;3B:FLNA在三阴性乳腺癌组织中表达图3 环状RNA与亲本基因的关系

最近研究发现环状RNA可以调控其亲本基因的表达,如在食管鳞状细胞癌中,随着circRNA ITCH的表达增加,ITCH的表达也增加[13];circGFRA1在乳腺癌中沉默后,GFRA1的表达也随之下调[14]。此外,circITGA7,作为一个抑癌基因,不仅可以通过海绵miR-370-3p上调NF1的表达,也可以通过抑制与RAS信号相关的剪接转录因子RREB1来上调其亲本基因ITGA7的表达,从而抑制结直肠癌的增殖和转移[18],并且通过组织分析发现circITGA7 和 ITGA7的表达呈正相关,因此环状RNA也可以通过上调亲本基因调控细胞的生物学功能[15]。在本研究中我们发现C19orf18、METTL3、CDC27和ACSL3在三阴性乳腺癌组织中表达下调,FLNA在三阴性乳腺癌组织中表达上调,与hsa_circ_0000965、hsa_circ_0000523、hsa_circ_0044234、hsa_circ_0058451和hsa_circ_0091994在乳腺癌组织中表达一致,预示着这些环状RNA可能正向调控其亲本基因,但是其具体发生机制仍待后续研究。

综上所述,用生物信息学的方法分析Luminal A型和三阴性乳腺癌组织中差异表达的环状RNA,并建立了环状RNA/miRNAs相互作用网络,同时发现部分环状RNA的表达与其亲本基因呈正相关,这为乳腺癌侵袭转移的治疗提供新的思路。但是其分子机制仍待进一步研究证实。