基于TCGA数据库前列腺癌预后相关mRNA风险模型的构建

张汉荣，林联拯，赵 力，于忠英，李金雨

（中国人民解放军联勤保障部队第909医院（厦门大学附属东南医院）泌尿外科，福建 漳州363000）

0 引言

前列腺癌是一种常见而复杂的恶性肿瘤，根据最新的全球癌症统计数据，全球每年约有130万例前列腺癌新发病例，死亡病例达36万例，占所有肿瘤新发病例的15%，是男性第二常见的肿瘤类型和第五大癌症死亡原因［1-2］.目前，前列腺癌的治疗方法主要有手术治疗、放疗、化疗和激素治疗［3-4］.选择何种治疗手段与肿瘤的分期、分级及危险度等因素密切相关［5］.前列腺癌的预后模型的构建、或利用构建的预测模型指导前列腺癌患者分层治疗是前列腺癌的研究热点之一.目前的预测模型大多数是以患者的年龄、肿瘤大小、TNM分期等临床资料为影响因素，预测模型的效率和精度上仍存在不足［6-8］.随着高通量基因测序技术的发展和成熟，利用癌症基因组或转录组数据构建的癌症预后模型显示出良好的预测性能［9-10］.本研究的目的是，利用公共数据库TCGA中提供的前列腺癌RNA-Seq数据，筛选潜在的前列腺癌相关基因，并利用这些候选基因构建预测模型，为临床前列腺癌的风险评估、分层治疗提供参考依据.

1. 资料与方法

1.1 研究资料

从the cancer genome atlas（TCGA）数据库（http：//www.tcga.org/）下载前列腺癌及对应癌旁组织的mRNA表达数据与临床资料.TCGA数据库共收录前列腺癌癌组织492例，癌旁组织52例.

1.2 前列腺癌组与对应癌旁组织差异基因的筛选与富集分析

根据TCGA下载的数据，利用R软件的“edge”包筛选出癌组织与癌旁组织间差异表达的mRNA.筛选标准为：|Log2FC|≥1，Padj≤0.05.利用Volcano Plot包对差异基因数据绘制火山图.利用clusterProfiler包对前列腺癌组织的上调基因分别进行GO富集和KEGG富集分析，以了解前列腺癌与癌旁组织差异基因的基因功能和信号通路的富集情况.

1.3 差异基因对前列腺癌患者预后的影响

以上调最明显的10个基因为候选基因，分别以上述基因mRNA表达量中值为分界线，将492例前列腺癌患者分为高表达组和低表达组，并做生存分析.比较差异基因对前列腺癌预后的影响.生存分析采用Log-rank检验方法，P＜0.05为差异有统计学意义.

1.4 COX风险回归模型的构建

将上述生存分析有统计学意义的基因纳入COX回归分析，以各基因的表达量为自变量，患者的生存时间为因变量.COX回归模型计算各患者的风险得分（Risk score），最后根据患者的风险得分将492例患者分为高危组和低危组，进一步进行生存分析，比较高危组和低危组的生存差异.ROC曲线用于验证该风险模型的预测精度.

2. 结果

2.1 前列腺癌对应癌旁组织的差异基因及富集分析

基因表达差异分析结果如图1A所示，共筛选表达上调基因1978个，下调基因1644个.其中，上调最明显的基因为：PCA3、AMACR、MTND4P12、RNY3P8、DLX1、OR51E2、PCAT14、GOLM1、HPN、GLYATL1.下调最明显的基因为：ADAM33、MEG3、GSTM2、GABRE、CYP3A5、CORO6、GOLGA8A、DIO3OS、SLC26A10、HIF3A.

对上调基因做GO功能富集分析，结果如图1B所示，上调的基因主要与carboxylic acid biosynthetic process，organic acid biosynthetic process及sulfur compound metabolic process等代谢过程有关.KEGG富集分析（图1C）显示，上调基因主要集中在Phagosome、Cell adhesion molecules和Hematopopoietic aureus infection等信号通路.

图1 差异基因火山图及GO功能、KEGG通路富集分析结果

2.2 差异基因对前列腺癌患者预后的影响

上调最明显的10个基因对预后的影响分别如图2所示，生存分析有显著差异的基因有6个，分别为PCA3、MTND4P12、RNY3P8、OR51E2、PCAT14、GOLM1.且上述基因的生存分析中，均表现为低表达组的总生存率较高表达组高.因此，考虑上述6个基因可能是前列腺癌发生、发展的关键基因.

图2 10个差异基因对492例前列腺癌患者预后的影响

2.3 COX风险回归模型的构建

以各基因的表达量为自变量，生存时间为因变量，得到风险得分公式为：Risk score=0.5503×（PCA3）+0.6777×（MTND4P12）+0.5605×（RNY3P8）-0.365×（OR51E2）+0.5358×（PCAT14）-0.1972×（GOLM1）.计算492例前列腺癌患者的风险得分，并以风险得分的中值将患者分为高风险组与低分险组（图3A），该模型的ROC曲线（图3B）曲线下面积为0.801，故认为该模型的预测精度较好.高、低风险组的生存分析结果如图3C所示，结果提示，高风险组的总生存时间明显较低风险组低，P＜0.0001.从高、低风险组的生存曲线可看出，多个基因表达水平构建的预测模型，较单一基因表达水平的生存分析差异更大，提示多基因表达构建的预测模型比单一基因构建的模型预测精度更高.

图3 COX风险回归模型的构建

3. 讨论

与其他肿瘤相比，前列腺癌具有更明显的肿瘤异质性，不同个体间的差异大［11］.良好的预后模型可以促进临床咨询和指导医生制定治疗和随访计划.目前临床上用于评估前列腺癌患者预后的模型主要有前列腺特异抗原（PAS）水平，Gleason评分肿瘤临床分期等［12-14］.然而，上述用于评估前列腺癌患者危险度的模型的精确度仍不能满足准确的患者分类和分层的需求.因此，临床上亟需一种可精确预测患者危险度的预测模型，以指导前列腺癌患者的指导治疗，实现个体化治疗.目前已有大量基因单个基因表达量构建的前列腺癌预后模型，如Mee YK等［15］的研究提示，高表达FOXM1基因与前列腺癌预后差相关.除FOXM1基因外，近年来，ARHGAP10［16］，SNW1［17］，NCAPH［18］等基因也有类似报道.但是基于单基因构建的预测模型预测效能及模型的稳定性相对较差.L Schmidt等［19］通过对接受根治性前列腺切除术的123名男性的前列腺癌组织样本进行全基因组miRNA表达谱分析，筛选出6个与前列腺癌预后最相关的miRNA，并利用其中的4个miRNA的表达量比值构建预测模型.结果显示该预测模型对前列腺癌具有良好的预测作用.

本研究通过对TCGA数据库收录的492例前列腺癌组织和52例癌旁组织的全基因组测序结果分析，筛选出在癌组织中高表达的基因.以上调最明显的10个基因为候选基因，分别以各基因的表达水平中值为界，比较高表达组与低表达组间患者预后的差异.结果提示，10个候选基因，有6个（PCA3、MTND4P12、RNY3P8、OR51E2、PCAT14、GOLM1）基因的表达量与预后相关，均表现为高表达提示预后不良.进一步对联合上述6个基因做COX风险回归模型.以上述6个基因的表达量为自变量，患者生存资料为因变量.结果提示，联合6个基因构建的风险预测模型具有良好的预测精度（AUC=0.801）.根据模型可利用各患者上述6个基因的表达量计算该患者的风险得分.本研究进一步将492例患者风险得分的中值为分界线，比较高、低分险组生存差异.结果提示，高风险组的生存时间明显较低风险组短（P＜0.00001）.联合6个基因构建的风险模型做的生存分析的P值明显小于单一基因生存分析的P值，提示，该风险模型具有更好的预测精度.

综上所述，本研究利用生物信息学方法，成功构建了一个以6基因（PCA3、MTND4P12、RNY3P8、OR51E2、PCAT14、GOLM1）mRNA表达水平为自变量的风险模型.该模型显示出良好的预测精度.