基于结肠癌癌变基因识别诊断分子标志

何 明，蔡 浩，罗 芸，郭 有

（赣南医学院第一附属医院医药大数据与生物信息研究中心，江西 赣州 341000）

结肠癌是我国发病率第四、死亡率第五的恶性肿瘤［1］。随着生活方式和饮食结构变化，结肠癌发病率逐年上升并呈年轻化趋势，患者早期无明显症状，大多数确诊已处于晚期，该肿瘤病死率高且预后差。1975 年研究者提出结肠癌发展遵循“正常黏膜—腺瘤—腺癌”的顺序［2］，结肠腺瘤是结肠癌最主要的癌前疾病［3］。超过80%结肠癌是由结肠腺瘤癌变而来，76%～90%的结肠癌可通过肠镜下腺瘤摘除进而降低发病率［4］。因此精准识别癌前疾病可有效抑制结肠癌发展，探索结肠癌的发病机制以及结肠癌变前诊断的分子标志迫在眉睫，以便患者尽早诊断、开展治疗，降低结肠癌发病率和死亡率。

普遍认为，结肠癌癌变是正常结肠黏膜在基因突变、表观遗传和基因表达变化逐渐积累发展的过程［5-9］。随着分子生物学技术发展，微阵列技术因可同时对样本的数万个基因表达水平进行检测，极大丰富了相关肿瘤的数据，已广泛应用于肿瘤致癌过程、诊断治疗、预后预测等研究。FEARON E R 等［5］研究认为APC、KRAS、TP53和DCC基因在结肠癌发展中承担着重要角色。KITA H 等［6］首次通过微阵列技术确定了在正常结肠黏膜发展为结肠腺瘤过程中显著差异表达的180 个基因，这些基因证实与结肠病变有密切关联。FONSECA A S 等［7］利用微阵列技术研究10 例配对结肠腺瘤和结肠腺癌基因表达谱，发现ETV4 在腺癌中显著高表达并参与细胞增殖和迁移过程。研究表明COX-2 从结肠腺瘤到结肠癌的发展中表达显著增加［8-9］。KANAOKA S等［10］研究发现，COX-2 mRNA 表达水平可有效区分结肠癌和无肿瘤患者，揭示可作为结肠癌检测的潜在标志物。但临床实践中活检组织才是结肠癌诊断治疗的重要组织类型，又因活检组织样本取材的限制，基于结肠活检组织的高通量基因表达研究极少，至今尚未发现应用于临床结肠癌分子诊断的生物标志物或特征。因此，基于活检组织的基因表达谱研究结肠癌癌变过程的关键基因及分子机制将有助于结肠癌诊断分子标志的挖掘。

基于此，本研究采用GEO 数据库的基因表达谱数据，分析活检组织中正常结肠黏膜—结肠腺瘤—结肠癌转变中表达差异基因，筛选从正常组织到肿瘤组织逐渐高表达或低表达的基因，通过生物信息学方法挖掘结肠癌癌变特征基因（Colorectal Cancer Carcinogenesis Genes，CRCCG）和分子机制，利用另一组织类型——手术切除组织和在线数据库验证癌变特征基因表达水平，GSVA 计算单个样本CRCCG-GSVA评分，ROC曲线评价CRCCG-GSVA评分的诊断价值，最后对CRCCG 进行生存分析，为开发临床诊断治疗结肠癌生物标志物提供新思路。

1 材料与方法

1.1 数据集来源 本研究在GEO 数据库中下载两组基于GPL570 平台检测的基因表达谱数据集：（1）GSE37364 为活检组织的数据包含38 例正常结肠黏膜、29 例结肠腺瘤和27 例结肠癌［11］；（2）GSE20916为手术切除组织的数据包含24例正常结肠黏膜、45例结肠腺瘤和36例结肠癌［12］。

1.2 数据预处理 对于affymetrix 平台下载的CEL格式的原始表达谱数据，在R 语言（版本3.6.1）中采用affy 软件包的RMA 算法对其背景校正、分位数校正和以2 为底log 转换，下载GPL570 平台的CDF文件，将探针ID对应至基因ID。去除数据中没有基因ID 对应或同一探针ID 对应多个基因ID 的探针ID，对于多个探针ID 对应同一个基因ID，则取它们的信号均值作为此基因的表达值。

1.3 筛选差异基因 使用R 语言中limma 软件包设置筛选条件为FDR＜0.05 且差异变化倍数（Fold Change， FC）＞1.5，获得正常结肠组织与结肠腺瘤、结肠腺瘤与结肠癌的差异表达基因，差异倍数≥1.5为上调基因，差异倍数＜0.67 为下调基因。利用VennDiagram 软件包分别对上调基因和下调基因绘制韦恩图，上调基因的交集和下调基因的交集定义为结肠癌癌变的差异基因。

1.4 分析差异基因参与的生物功能和通路 利用clusterProfiler、enrichplot 软件包对结肠癌癌变的差异基因进行GO 功能注释［13-14］和KEGG 通路富集分析 并 作 图［15-16］。GO 功 能 注 释 包 括 生 物 过 程（Biological process，BP）、细 胞 成 分（Cellular component，CC）和分子功能（Molecular function，MF）3 部分。显著富集筛选标准：P＜0.05 且FDR＜0.05。

1.5 构建蛋白质互作网络及筛选特征基因 将结肠癌癌变的差异基因导入STRING 在线数据库（https：//cn.string-db.org/），选择multiple proteins，物种选择homo sapiens，构建差异基因的蛋白质互作网络（Protein-protein interaction networks，PPI），置信分数设置为0.7，下载TSV 文件将其输入Cytoscape 软件（版本3.9.0）绘制可视化的网络图，通过cytoHubba插件［17］，采用MCC、DMNC、MNC、Degree、EPC、BottleNeck、EcCentricity、Closeness、Radiality、Betweenness、Stress、ClusteringCoefficient 12 种分析方法筛选前30个基因，挑选均出现的基因构建蛋白质互作网络并定义为CRCCG。

1.6 计算GSVA 评分与绘制ROC 曲线 基因集变异分析（Gene set variation analysis， GSVA）可以计算单个样本中特定基因集的富集评分。为评估CRCCG 的诊断性能，使用GSVA 软件包对每个样本计算CRCCG 的GSVA 评分。基于ROC 曲线，采用pROC 软件包根据GSVA 评分计算ROC 曲线下面积（Area Under Curve，AUC）分析CRCCG-GSVA 评分的诊断价值，AUC＞0.9 表明CRCCG-GSVA 评分具有很高的诊断价值，AUC＞0.8 表明CRCCG-GSVA 评分具有较高的诊断价值。

1.7 基因表达水平的验证 GSE20916 验证CRCCG 在3种结肠组织类型中的表达情况。GEPIA数据库（http：//gepia. cancer-pku. cn/）包括TCGA 数据库和GTEx 数据库中349 例癌旁正常结肠样本和275 例结肠癌组织样本［18］，验证CRCCG 在癌旁正常结肠样本和结肠癌样本的表达情况，并对其进行生存分析。生存分析根据结肠癌患者基因表达值的中位数将样本分为高表达组和低表达组。Log-rank检验基因高表达组与低表达组的总体生存率是否存在差异。COX 比例风险计算风险比（Hazard Ratio，HR）。P＜0.05 代表基因高表达组与低表达组的结肠癌患者生存率差异有统计学意义。

2 结果

2.1 差异表达基因 本研究利用limma 算法分析了数据集GSE37364 的基因表达谱，结果显示正常结肠黏膜与结肠腺瘤有1 845 个差异表达基因，838个上调基因，1 007 个下调基因；结肠腺瘤与结肠癌有1 061个差异表达基因，714个上调基因，347个下调基因。韦恩图发现，正常结肠黏膜和结肠腺瘤、结肠腺瘤和结肠癌的两组上调基因交叠了76 个基因（图1A），正常结肠黏膜和结肠腺瘤、结肠腺瘤和结肠癌的两组下调基因交叠了107 个基因（图1B），得到183 个结肠癌癌变的差异基因。这183 个基因表达从正常结肠黏膜—结肠腺瘤—结肠癌发展过程中，呈逐渐升高或降低趋势。热图聚类显示，183个基因可以将正常结肠黏膜、结肠腺瘤、结肠癌组织进行分类（图1C）。

图1 GSE37364中差异表达基因的韦恩图及热图

2.2 差异基因的功能分析 为进一步研究183 个结肠癌癌变差异基因的分子机制，基于GO 富集分析，差异基因主要参与白细胞迁移、类固醇激素反应、细胞趋化、白细胞趋化、细胞外基质分解、锌离子反应、白细胞聚集等生物过程；主要富集于受体配体活性、受体调控活性、G 蛋白偶联受体结合、细胞因子活性、趋化因子受体结合、趋化因子活性等分子功能（图2A）。通过KEGG 信号通路富集分析发现，差异基因主要参与10条信号通路：细胞因子受体相互作用、白细胞介素-17（Interleukin-17，IL-17）信号通路、病毒蛋白与细胞因子和细胞因子受体的相互作用、趋化因子信号通路、TNF 信号通路、NF-κB信号通路等信号通路（图2B）。

图2 差异基因GO注释和KEGG信号通路分析

2.3 PPI 网络的建立及关键基因的挖掘 基于STRING 在线数据库利用Cytoscape 软件对183 个结肠癌癌变相关的差异基因进行PPI 网络分析，得到差异基因之间有着显著的交互作用（图3A）。采用12 种分析方法筛选前30 个基因，取共有基因作为CRCCG，最 后 得 到11 个 基 因：ANLN、CXCL11、CXCL3、CXCL8、GZMB、IL1B、MMP3、SHCBP1、TIMP1、TRIP13、UBE2C。这些基因重新构建PPI 网络，结果11个CRCCG之间联系密切（图3B）。

2.4 CRCCG 的表达验证 在活检组织样本数据集GSE37364数据集中，CRCCG 在正常结肠黏膜、结肠腺瘤、结肠癌中的表达值逐渐增大（图4A）。在手术组织样本数据集GSE20916中进一步验证CRCCG的表达情况，结果显示CRCCG 在正常结肠黏膜、结肠腺瘤、结肠癌中的表达值也逐渐增大（图4B）。另外，应用GEPIA 在线数据库研究CRCCG 的表达情况，结果相比于正常结肠，CRCCG 在结肠癌中均显著高表达（图5）。

图5 GEPIA数据库验证CRCCG在癌旁正常结肠（灰色）和结肠癌（红色）表达水平

2.5 CRCCG-GSVA评分的分类效果 利用GSVA计算单个样本中CRCCG 的GSVA 评分。结果显示GSE37364 数据集中CRCCG-GSVA 评分随着肿瘤进展逐渐增大（图6A），ROC 曲线得到CRCCG-GSVA评分对正常结肠、结肠腺瘤和结肠癌进行分类的AUC 值为0.922，其中区分结肠癌与正常组织的AUC 为0.984，结肠腺瘤与正常组织的AUC 为0.929，结肠腺瘤与结肠癌的AUC 为0.823（图6B）。基于ROC 曲线选择最佳CRCCG-GSVA 评分分类阈值：＜-0.56判定为正常结肠，-0.56～0.3之间判定为结肠腺瘤，＞0.3 判定为结肠癌。CRCCG-GSVA评分的分类效果在GSE20916 数据集中也得到了验证，CRCCG-GSVA 评分也随着结肠癌发展逐渐增加（图6C），基于最佳分类阈值识别结肠腺瘤与正常结肠的准确率为84.06%，敏感性为79.17%，特异性为86.67%；识别结肠癌与正常结肠的准确率为93.33%，敏感性为91.67%，特异性为95.83%；识别结肠腺瘤与结肠癌的准确率为87.65%，敏感性为84.44%，特异性为91.67%。CRCCG-GSVA 评分可以很好区分正常结肠黏膜、结肠腺瘤、结肠癌。

图6 两组数据集的CRCCG-GSVA评分分布与ROC曲线

2.6 基因生存分析 对CRCCG进行生存分析，发现3个基因的表达与结肠癌患者的生存率相关。CXCL3和IL1B高表达组的总体生存期（Overall Survival，OS）显著高于低表达组（CXCL3：HR=0.62，P=0.047；IL1B：HR=0.61，P=0.043），TIMP1高表达组的OS 显著低于低表达组（HR=1.7，P=0.034）（图7）。

图7 GEPIA数据库中CXCL3、IL-1β、TIMP1基因的生存分析

3 讨论

结肠癌癌变的过程是个涉及多种基因变化的复杂过程。当下，临床上对结肠癌仍缺乏有效的治疗手段，因而深入探索结肠癌病变过程的分子机制尤为重要，对早期诊断、靶向治疗和预后评估具有指导作用。典型的结肠癌癌变过程就是从正常结肠黏膜—结肠腺瘤—结肠癌，根据这一发展特点，本研究关注结肠癌随着肿瘤进展而逐渐上调或者下调的差异基因。这些差异基因在正常结肠黏膜转变为结肠癌过程中表达水平是逐渐积累变化的，在结肠癌癌变过程的作用尤为显著。分析活检组织数据集GSE37364 获得从正常结肠组织到结肠癌进展中表达水平逐渐升高或者降低的183个差异基因，差异基因与多条肿瘤相关的生物过程和信号通路有关。基于PPI 网络利用cytoHubba 插件挖掘CRCCG，为避免单一种算法的计算偏差，采用了12种分析算法筛选在所有算法中均出现的11 个基因为CRCCG（ANLN、CXCL11、CXCL3、CXCL8、GZMB、IL1B、MMP3、SHCBP1、TIMP1、TRIP13、UBE2C）。手术切除组织样本和GEPIA 数据库验证了这11 个基因的表达趋势与活检组织样本的表达趋势一致，随着肿瘤进展表达水平显著升高。基于ROC 曲线确定了CRCCG-GSVA 评分最佳分类阈值：＜-0.56判定为正常结肠，-0.56～0.3 之间判定为结肠腺瘤，＞0.3 判定为结肠癌，独立数据集也验证CRCCG-GSVA 评分具有诊断价值。进一步挖掘CRCCG 与结肠癌预后关系发现IL1B、CXCL3和TIMP1这3个基因与结肠癌患者预后显著相关。

已有研究表明，部分结肠癌癌变基因与结肠癌发展密切相关。白细胞介素-1β（Interleukin-1 beta，IL1B）是癌症相关炎症的重要介质，在包括结肠癌等多种癌症中表达水平升高［19-20］。IL1B 通过Zeb1 激活结肠癌干细胞自我更新和上皮间质转化促进结肠癌生长和侵袭［21］。结肠癌细胞可诱导巨噬细胞释放IL1B，调控GSK3β磷酸化，稳定β-连环蛋白，增强T 细胞因子依赖性的基因活化，诱导肿瘤细胞中Wnt靶基因表达，促进肿瘤细胞生长［22］。趋化因子3（Chemokines 3，CXCL3）和趋化因子8（Chemokines 8，CXCL8）属于CXC 趋化因子家族的单链蛋白质，对炎症和抗肿瘤免疫至关重要。CXCL3 表达受KRAS-IRF2 调控，介导结直肠癌免疫治疗和免疫治疗耐药性［23］。在结肠癌组织中CXCL3 高表达，与肿瘤大小、肿瘤血栓、DNA 修复、细胞周期、细胞凋亡和P53 调控通路密切相关，致使结肠癌患者总体生存期较差［24-25］。CXCL8 在早期到晚期结肠癌以及转移和非转移病例中表达水平显著高于正常组织样本。另外发现血清样本中CXCL8 浓度在结肠癌患者中也显著高于健康对照组。与没有转移的结肠癌患者相比，有远处转移的结肠癌患者CXCL8 水平也显著升高［26］。组织金属蛋白酶抑制剂（Metalloproteinase inhibitor 1，TIMP1）具有抑制MMPs 能力，有效控制细胞外基质成分分解，影响肿瘤侵袭和转移。与正常结肠组织相比，TIMP1 在结肠肿瘤组织中表达水平升高，可特异性调节FAK-PI3K/AKT 和MAPK 通路，表达抑制导致结肠肿瘤增殖和转移能力减弱，细胞凋亡能力增强［27］。研究发现TIMP1 既在结肠癌组织样本中高表达，也在血液样本中高表达。在健康人的血小板中TIMP1 mRNA几乎检测不到，但在结肠癌患者的血小板中显著增加，检测血小板中TIMP1 mRNA 比检测CEA 或CA199 在结肠癌诊断具有更高的敏感性和特异性［28］。还有研究揭示血浆中TIMP1 表达水平是结肠癌患者生存的重要预后因素，术前血浆中高表达的TIMP1 预示结肠癌患者预后不良［29］。综合已有研究和本研究结果，CRCCG 不仅可以作为结肠癌组织检测的分子标志，或许还可以通过血液检测其表达值诊断结肠癌，早期对患者进行治疗，提高生存率。

本研究结合了公共数据库和多种组织类型分析结肠癌癌变过程中随肿瘤进展显著差异表达的基因，确定了CRCCG。CRCCG 是潜在的结肠癌诊断分子标志，其中IL1B、CXCL3和TIMP13个基因还与结肠癌的预后相关。这些结果为结肠癌的诊断和预后提供参考。肠镜活检组织检查是结肠癌诊断的金标准，但具有创伤性、需充分肠道准备、操作复杂、操作过程患者痛苦、患者依从性较差等缺点，大多数患者首诊已处于晚期。因此，急需开发敏感、特异、安全、无创的诊断方法便于诊断，从而降低结肠癌的死亡率。血液中含有与肿瘤相关的许多生物分子，同时血液又具有采样简便、风险较小等优点，因此基于血液的生物标志物可能是结肠癌诊断的更佳基质。本研究筛选得到的CRCCG 与分泌蛋白、免疫反应密切相关，其中CXCL8 和TIMP1已有研究分析其在结肠癌血液中的表达，CRCCG 或许可作为结肠癌诊断的血液分子标志。后续将基于血液样本进一步研究这个诊断标志，以期得到非侵入性的诊断标志便于临床应用。