高级别浆液性卵巢癌复发相关的潜在功能性关键miRNAmRNA：基于生物信息学方法

张攀扬，何明敏，曾园媛，蔡雄伟

1重庆市妇幼保健院（重庆医科大学附属妇女儿童医院）妇产科，重庆 401147；2陆军军医大学第一附属医院妇产科，重庆400038

2020年，全球卵巢癌（OC）新发病例313 959例，死亡207 252人，在6种女性癌症中死亡人数排名第3。高级别浆液性卵巢癌（HGSOC）是卵巢癌最常见和具有侵袭性的亚型，占卵巢癌亚型的67.5%［1］。由于早期HGSOC没有特征性症状，通常在晚期被诊断。尽管手术和系统治疗的不断改进和发展，OC患者的死亡率仅略有下降。大约50%的HGSOC 通过同源重组(HR)表现出潜在的DNA修复缺陷，并且对双链断裂诱导剂如铂类似物和PARP抑制剂（PARPi）高度敏感［2］。然而，尽管一线铂基化疗使约70%OC 患者的临床完全缓解，但超过70%的患者在治疗后出现复发，此时该病对后续治疗的反应较差［3］。同样，尽管在HR缺陷的HGSOC中，超过80%最初对PARPi反应良好，但随着获得性耐药的发生，预后仍不尽如人意，5年生存率仅为35%～45%［4］。与其他肿瘤相比，HGSOC缺乏有效的分子靶向治疗。因此，迫切需要识别HGSOC进展和转移的关键分子靶标，是改善诊断和预后的关键。

对高危患者进行分类有助于阐明疾病的发病机制，在之前的研究中，根据转录组分析将HGSOC分为分化型、间充质型、免疫反应型和增殖型4个亚型［5］。这四种亚型的生存期没有显著差异，复发也没有讨论。同时，基因表达的异质性和基因间的协同或相互作用仍是分类的障碍。因此，需要一种新的分类方法，将HGSOC分成与复发相关的亚型，以便寻找这些亚型中的关键分子调控网络。

大部分人类基因组被转录成非编码RNA，蛋白质编码的基因组不到2%［6］。MiRNA是一类由内源基因编码的长度约为22 个核苷酸的非编码单链RNA分子，在真核生物中参与转录后基因表达调控。MiRNA在肿瘤细胞增殖、凋亡和侵袭等方面具有丰富的生物学效应，在癌症的发生、转移和治疗耐药性中发挥着重要作用［7］。据报道，miR-10b、miR-125b和miR-145在肿瘤发展过程中下调，而miR-21和miR-155则上调，分别作为抑癌基因和癌基因发挥作用［8］。而在卵巢癌中，miR-363-3p-SPOCΚ2轴参与了肌动蛋白细胞骨架的调节，促进了卵巢癌的分期进展［9］；miR-424-5p可通过靶向ACSL4消除卵巢癌细胞的铁死亡［10］；miRNA-584通过负向调节LPIN1抑制卵巢癌的进展［11］。虽然已经发现有多个miRNAs参与HGSOC的发病机制，然而，在卵巢癌中与复发相关的miRNAs调控网络还没有被报道。

本研究从肿瘤基因组图谱（TCGA）数据库下载HGSOC转录组RNA-seq和临床相关数据，通过生信分析将HGSOC分成与复发相关的两种亚型，并根据两种亚型中的差异基因构建一个由5 个miRNAs 和41 个mRNAs组成的调控网络，其中miR-506-3p表达差异最显著。然而，miR-506-3p是否参与HGSOC的EMT及复发，目前尚不清楚。因此，本研究通过探讨miR-506-3p对卵巢癌细胞迁移、侵袭和EMT的影响，明确miR-506-3p 在HGSOC 复发中的作用。以期为探索HGSOC 中复发的具体机制提供新的实验依据，为后续HGSOC的早期诊断和治疗提供新的分子靶点。

1 材料和方法

1.1 数据收集

使用R包TCGA Biolinks从TCGA数据库下载354例HGSOC患者的RNA-Seq原始数据和临床信息［12］，GO_BP基因集和癌症标志基因集（Hallmark gene sets）从MSigDB 数据库收集（version 7.1）。GSE25204 和GSE73582 用于验证miRNA 表达与卵巢癌复发的相关性。

1.2 计算基因富集分数

基因计数用DESeq2［13］和log2转换后进行归一化处理，滤除平均表达量低于5的基因。在R中使用基因集变异分析(GSVA)方法估算归一化计数后得到每个样本的GO基因集富集分数。

1.3 与复发相关的HGSOC亚型分类

对GO_BP基因集和疾病无进展期的7530个富集分数进行COX单因素回归分析，保留复发相关基因集（P＜0.05）。使用非负矩阵分解（NMF）方法根据复发相关基因集的富集分数对HGSOC患者进行分类，该方法由R包CancerSubtypes［14］提供。聚类数（Κ）由最大平均轮廓值选择。

通过R包limma［15］识别具有差异富集分数的GO_BP基因路径(adj.p＜0.001)和差异表达基因（adj.p＜0.05和|log fold change|＞1）。使用R包clusterProfiler对差异表达基因进行HALLMARΚ富集分析［16］，基因集从MSigDB数据库下载。通过Pearson相关性检验计算miRNA与靶向差异表达基因的相关性，miRWalk 2.0预测miRNA靶向基因［17］。

1.4 MiRNA-mRNA网络构建

使用Cytoscape 软件（v3.6.0）构建miRNA-mRNA网络图。

1.5 细胞培养

人卵巢癌细胞株（SΚOV3和CAOV3）由本实验室保存。SΚOV3 和CAOV3 细胞在含10%胎牛血清的DMEM培养基（Gibco）中培养，置于37 ℃、5%CO2的培养箱中。

1.6 细胞转染，RNA提取和RT-qPCR

miR-506-3p mimics和inhibitor购自上海吉玛医药科技有限公司。mimics-NC（miR-NC）和inhibitor-NC（in-NC）作为对照。SNAI2 siRNA（si-SNAI2-1#GCA GACAGGTCAAATCTGA,si-SNAI2-2# GCAGACC CATTCTGATGTA,si-SNAI2-3# GCACAAACATGA GGAATCT）购买自广州瑞博生物科技有限公司。所有转染均使用LipoFiterTM3.0 脂质体转染试剂（汉恒生物）。总RNA 使用Trizol（TAΚARA）法提取，用Nanodrop 2000（Thermo）分光光度计检测RNA质量及浓度，使用PrimeScript RT试剂盒（Takara）将1000 ng总RNA反转录为cDNA，反应条件设置为：预变性95 ℃30 s，变性95 ℃5 s、退火58 ℃30 s及延伸72 ℃5 s，共40 个循环，在Bio-Rad CFX96（Bio-Rad）系统上使用SYBR Premix Ex TaqTM试剂（TAΚARA）进行RT-qPCR分析，采用2-ΔΔCt法分析miRNA和mRNA的相对定量数据，U6和GAPDH作为各自的内参。

1.7 划痕实验

将做相应处理后的卵巢癌细胞SΚOV3和CAOV3分别接种于6孔板中，每孔接种5×106个细胞，设置3个独立的实验组和对照组。37 ℃5%CO2的培养箱中培养。当细胞密度达到90%左右时，用灭菌的10 μL枪头垂直于板的中间划痕，PBS洗去掉落的细胞，换成无血清培养基培养，并分别在0 h和24 h时间点观察细胞并拍照。

1.8 Transwell 分析

事先铺好基质胶，将转染后的细胞分别取2×104细胞接种于上室，加入300 μL无血清培养基培养，在下室中加入800 μL完全培养基，置于37 ℃，5%CO2的培养箱中培养，分别设置3个独立的实验组和对照组。培养36 h后取出小室，使用棉签擦拭小室内部，使用结晶紫染色10 min。然后用PBS洗涤细胞，用光学显微镜计数。

1.9 Western blotting

使用含1%蛋白酶抑制剂的RIPA蛋白裂解液4 ℃裂解细胞30 min，14 000 r/min，4 ℃离心30 min。收集上清，使用BCA 蛋白定量试剂盒（碧云天）检测蛋白浓度，加入蛋白上样缓冲液，同样量（30 μg）的蛋白经SDS-PAGE 分离，用Trans-Blot Turbo 快速转移膜仪（Bio-rad）将其转移到PVDF膜上，用5%牛奶封闭2 h，然后与抗SNAI2（1∶1000）（Proteintech，12129-1-AP）和GAPDH（1∶2000）（Proteintech，60004-1-Ig）一 抗4 ℃孵育过夜。HRP 标记的山羊抗兔和山羊抗鼠二抗（1∶5000）室温孵育1.5 h，最后，用ECL 发光液（Millipore）检测灰度值。

1.10 双荧光素酶报告基因

通过PCR扩增SNAI2-3'UTR区域，然后克隆到荧光素酶报告载体pmirGLO（Promega）中，包括SNAI2-3'UTR WT 组，SNAI2-3'UTR Mut 组，根 据Dual-Glo 萤光素酶检测试剂盒（Promega，E2920）说明书，检测相对荧光素酶活性。

1.11 统计学分析

使用GraphPad Prism 5.0进行数据统计分析，计量资料以均数±标准差表示，组间比较采用独立样本t检验。P＜0.05为差异有统计学意义。

2 结果

2.1 通过对通路活性的聚类鉴定出两种与复发相关的HGSOC亚型

非监督方法NMF分析结果得到C1（n=191）和C2（n=163）两个子类型（图1A、B），比较临床分期、分级、年龄和肿瘤残留分布，C1和C2患者无明显偏差（表1），而C1组的复发风险明显高于C2组（图1C）。

表1 354例HGSOC患者的临床特点及分布Tab.1 Clinical characteristics and distribution of 354 HGSOC patients

图1 根据GO项的富集分数，采用NMF聚类方法识别亚型Fig.1 HGSOC subtypes identified by NMF clustering based on enrichment scores of GO terms.A:Similarity matrix of clusters under cluster number(K)from 2 to 5.B:Average sihouette widths under different K(2,3,4,5)for selecting the optimal K,and K=2 meets the highest average sihoutte width.C:HGSOC patients were divided into two clusters,C1 and C2.The survival plots illustrate recurrence associated with NMF clusters.The P values are determined by log-rank test.

2.2 HGSOC 亚型中差异表达的mRNAs 和miRNAs的鉴定

在C1和C2两种HGSOC亚型中共有645个表达差异显著的mRNAs（|log2FC|＞1;adj.p＜0.05），其中616个上调，29个下调（图2A、B）；以及47个差异表达显著的miRNAs（adj.p＜0.05，|log2 FC|＞0.5）其中18个miRNAs在C1 组 中 下 调，29 个miRNAs 上 调（图2C）。HALLMARΚ通路富集分析结果显示，主要富集在上皮间质转化（EMT）信号通路（图2D）。

图2 C1和C2患者间的差异表达基因与EMT相关Fig.2 Differentially expressed genes between C1 and C2 subtypes of HGSOC are related with EMT.A: Heat map of the differentially expressed mRNAs between C1 and C2.B:Volcano map of the differentially expressed mRNAs between C1 and C2.C: Volcano map of the differentially expressed miRNAs between C1 and C2.D: The differentially expressed genes are significantly enriched in"HALLMARK EPITHELIALMESENCHYMALTRANSITION".

2.3 与HGSOC 复发相关的关键miRNA-mRNA 网络构建

分析两个数据集GSE25204、GSE73582中差异表达的miRNAs，GSE25204 得到16 个表达差异显著的miRNAs（adj.p＜0.05，|log2 FC|＞0.5）；而在GSE73582中有39 个表达差异显著的miRNAs（adj.p＜0.05，|log2 FC|＞0.5），其中25个上调，14个下调。3个数据集中没有共同上调的miRNAs，而有5个共同下调的miRNAs，分别是miR-508-3p、miR-509-3-5p、miR-506-3p、miR-508-5p、miR-513c-5p（图3A、B）。接着使用miRWalk数据库对5个miRNAs进行靶基因预测，去除重复靶向的mRNAs后，共得到10494个mRNAs。由于C1与C2差异表达的mRNAs主要富集在EMT通路，因此我们将数据库预测得到的靶mRNAs与在C1和C2中差异表达EMT 通路中的mRNAs 取交集，共获得了41 个mRNAs，利用这5个miRNAs和41个mRNAs构建出与HGSOC复发相关的关键miRNA-mRNA网络（图3C）。

图3 miRNA-mRNA网络Fig.3 The miRNA-mRNA network.A,B:The common miRNAs in the 3 datasets.C:The miRNA-mRNA network,in which blue indicates mRNAs and purple miRNAs.

2.4 MiR-506-3p调控卵巢癌细胞迁移和侵袭

在这5个miRNAs中，miR-506-3p表达差异最大（表2）。RT-qPCR结果显示，mimics显著上调miR-506-3p的表达，而inhibitor显著下调了miR-506-3p的表达（P＜0.001，图4A）。划痕实验结果显示过表达miR-506-3p后抑制了SΚOV3及CAOV3的迁移，而敲低miR-506-3p后SΚOV3及CAOV3的迁移能力增强（P＜0.01，图4B）。同样，Transwell 实验显示miR-506-3p 抑制了SΚOV3及CAOV3的侵袭能力，而敲低miR-506-3p后SΚOV3及CAOV3的侵袭能力变强（P＜0.01，图4C）。

图4 MiR-506-3p调控SΚOV3和CAOV3细胞迁移和侵袭Fig.4 MiR-506-3p regulates SKOV3 and CAOV3 cell migration and invasion.A:Expression of miR-506-3p detected by qRT-PCR.B,C: Wound healing assay (Original magnification:×50) and Transwell assay (×100) demonstrate that miR-506-3p regulates SKOV3 and CAOV3 cell migration and invasion.Data are presented as Mean±SD.**P＜0.01,***P＜0.001.

表2 TCGA数据库中5个miRNAs的表达水平Tab.2 Expression levels of 5 miRNAs in TCGAdatabase

2.4 MiR-506-3p 通过SNAI2 调控卵巢癌细胞迁移和侵袭

根据correlation coefficient 值挑选出miRNAmRNA网络中miR-506-3p靶基因前五的mRNAs进行验证（表3），RT-qPCR结果显示在CAOV3和SΚOV3中过表达及敲低miR-506-3p后SNAI2表达差异最显著（P＜0.05，图5A～E）；Western blot 结果显示，miR-506-3p mimics抑制了SNAI2蛋白表达，而沉默miR-506-3p后促进了SNAI2蛋白表达（图5F）。以上结果说明在SΚOV3 和CAOV3 中miR-506-3p 有可能通过调控SNAI2表达影响细胞迁移、侵袭和EMT等过程。双荧光素酶结果显示，miR-506-3p显著降低SNAI2-WT组的荧光素酶活性（P＜0.001），但对突变组无明显影响（图6A）。此外，商业化合成SNAI2的3个siRNA（si-SNAI21#、2#和3#），RT-qPCR和Western blot检测结果显示，si-SNAI2-2#和si-SNAI2-3#显著抑制了SNAI2的表达（P＜0.001，图6B、C）。因此，si-SNAI2-2#和si-SNAI2-3#被用于进一步的实验。重要的是，划痕实验和Transwell 分析的数据显示，沉默SNAI2 逆转了miR-506-3p 对SΚOV3 和CAOV3 迁移和侵袭的抑制作用（P＜0.05，图6D、E）。

图5 MiR-506-3p调控SΚOV3和CAOV3细胞中差异基因的表达Fig.5 MiR-506-3p regulates the expression of differentially expressed genes in CAOV3 and SKOV3 cells.A-E: Expressions of differentially expressed genes detected by qRT-PCR.F:Western blotting for detecting the expression of SNAI2.Data are presented as Mean±SD.*P＜0.05,**P＜0.01,***P＜0.001.

图6 MiR-506-3p通过SNAI2促进CAOV3和SΚOV3的细胞迁移和侵袭过程Fig.6 MiR-506-3p promotes EMT in CAOV3 and SKOV3 cells via SNAI2.A:Luciferase activity assay confirming binding of miR-506-3p to SNAI2 3'UTR.B,C: SNAI2 knockdown by si-SNAI2-2# and si-SNAI2-3#.D,E: Wound healing assay (×50) and Transwell assay(×100)demonstrate that SNAI2 knockdown attenuates the migration-promoting effects of miR-506-3p inhibitor in CAOV3 and SKOV3 cells.Data are presented as Mean±SD.*P＜0.05,**P＜0.01,***P＜0.001.

表3 miR-506-3p相关系数前5的mRNAsTab.3 Top 5mRNAs with the high correlation coefficients with miR-506-3p

3 讨论

高级别浆液性卵巢癌是最常见的卵巢癌亚型，预后差。尽管患者在最初的治疗后能够完全缓解，但是肿瘤却长期复发。然而，复发的机制尚不清楚［18］。由于肿瘤病人具有高度的异质性，根据分子水平的差异来确定亚型需要大量的样本。在本研究中，我们将转录组数据转化为基因功能富集组，以便于我们鉴别亚型。使用与复发相关的基因集进行聚类，从而获得预后不同的亚型，进而揭示复发的潜在机制。

EMT是指上皮肿瘤细胞失去黏附能力获得间充质细胞迁移侵袭的能力以促进转移和耐药性的过程，在癌症复发中发挥重要作用［19,20］。在乳腺癌、肺癌及胰腺癌等多种肿瘤疾病中报道称EMT是常规化疗耐药的关键过程［21］。在本研究中，我们主要关注HGSOC的复发及其潜在的机制。根据复发相关功能基因集富集评分将患者分为两组，进一步的分析发现，参与富集得分的核心基因被显著富集在EMT通路，122个基因中有41个参与了EMT。既往研究表明miR-513a-3p表达下调后能通过调控EMT促进细胞运动，且能导致顺铂治疗后复发的卵巢癌患者预后不良［22］；miRNA-199b-3p通过靶向ZEB1的CHΚ1/E-cadherin/EMT信号通路抑制卵巢癌的生长和进展［23］；miR-200c通过调节EMT调控卵巢癌迁移和侵袭［24］。这些结果显示EMT是HGSOC复发的一个重要调控因素，且受到miRNA的调控。

MiRNA对疾病的调控通常是一个复杂的调控网络。研究显示，通过生信分析构建对胰腺癌预后具有显著的预测价值的ceRNA网络［25］；Chen等［26］分析鉴定了与膀胱癌相关的关键circRNA-miRNA-mRNA调控网络；Li 等［27］揭示了卵巢癌中的关键miRNA-mRNA表达谱。但是，这些研究仅分析了与疾病相关的分子网络，少有研究分析在HGSOC中与复发相关的分子调控网络。利用TCGA及GEO数据集分析与复发相关的两类亚型中的差异表达基因，我们构建了一个由5个miRNA和41个mRNA组成的与HGSOC复发相关的miRNA-mRNA调控网络，因miR-506-3p表达差异最大，提示可能是更好的靶点，因此选择该基因做进一步分析。

研究表明miR-506-3p通过负向调控MTMR6的表达，抑制卵巢癌细胞的增殖［28］；miR-506-3p通过靶向SIRT1/AΚT/FOXO3a信号通路抑制卵巢癌细胞增殖和促进细胞凋亡［29］。本研究中过表达miR-506-3p后，卵巢癌细胞迁移及侵袭能力减弱，而敲低miR-506-3p后，卵巢癌细胞迁移及侵袭能力增强。接着对调控网络中miR-506-3p的靶基因进行验证，其中SNAI2表达差异最为显著。SNAI2是与EMT密切相关的转录因子，研究显示SNAI2可通过抑制miR-222-3p的转录活性和上调PDCD10的表达来诱导EMT［30］；SNAI2可能通过调节GSΚ-3β/β-catenin通路促进前列腺癌细胞的生长、迁移、血管生成和干细胞生长［31］；SNAI2能促进口腔鳞状细胞癌的EMT和细胞干细胞形成［32］。挽救实验显示，敲低SNAI2后，逆转了过表达miR-506-3p对卵巢癌细胞迁移和侵袭的影响。通过利用双荧光素酶报告检测发现miR-506-3p与SNAI2 3'UTR区结合。以上结果提示miR-506-3p可能通过抑制SNAI2表达，调控卵巢癌细胞EMT，进而影响卵巢癌的发生发展。

综上所述，本研究利用TCGA及GSE数据集，根据复发相关功能基因集富集评分将患者分为两组，查找与复发相关的差异基因。构建与复发相关的EMT miRNA-mRNA调控网络。其中，miR-506-3p差异最显著，且miR-506-3p 可以通过调控SNAI2 抑制卵巢癌细胞迁移和侵袭的能力，并影响卵巢癌细胞的EMT 过程，提示miR-506-3p 异常降低可能是参与卵巢癌复发的重要机制之一。为进一步理解HGSOC复发的机制和找到潜在分子治疗靶点奠定了理论基础。但是本研究仅探讨了miRNA-mRNA 调控网络中差异最显著的miR-506-3p 和SNAI2，对于调控网络中的其他miRNA 及mRNA 如何影响卵巢癌发生发展还需进一步研究。