1. 引言
前列腺癌(Prostate Cancer)是一种严重危害男性健康的恶性肿瘤[1]。从流行病学数据来看,前列腺癌的发病率和病死率在全球范围内均处于较高水平[2]。前列腺癌的发生受多种因素的影响,包括遗传因素、年龄、种族、雄激素水平及生活方式等[3]。尽管近年来,前列腺特异性抗原筛查的广泛应用提高了早期诊断率,并促进了手术、放疗及内分泌治疗等治疗手段的发展,但由于肿瘤的异质性和耐药性的存在,部分患者仍可能进展为去势抵抗性前列腺癌,导致预后较差[4]。因此,探索前列腺癌的分子机制并寻找有效的生物标志物和治疗靶点对于提高患者生存率具有重要意义[5]。
近年来,越来越多的研究发现,前列腺癌的发生和进展与特定的差异表达基因(DEGs)密切相关,这些基因可能在肿瘤的发生、侵袭和转移过程中发挥重要作用[6]。高通量基因测序技术的快速发展,使研究者能够更全面地解析前列腺癌的基因表达谱,并为个性化治疗提供新的可能性。随着机器学习和生物信息学的融合应用,研究者们利用特征选择算法和分类模型对关键基因进行筛选,从而构建更精准的前列腺癌预测和分类模型[7]。这些研究不仅有助于前列腺癌的早期诊断,还为潜在的治疗策略提供了新的方向。
2. 数据来源与处理
本文所用数据均来源于公开的UCSC Xena数据库(https://xena.ucsc.edu/)。首先从UCSC Xena数据库中的TCGA板块选取了203例样本,其中分为152例前列腺癌样本和51例正常样本。这里为了弥补正常样本的不足,又从GTEX板块中选取了122例样本,并将两数据集进行合并。将合并后的数据集进行数据预处理,这个过程中去除了基因表达量较低的基因以及基因表达异常的样本。再根据数据库特征,对样本进行分类,最终得到151个前列腺癌样本和152个正常样本。
3. 方法
3.1. 差异基因筛选
考虑到TCGA和GTEX数据库之间具有批次效应,本文首先对合并后的数据进行了PCA和数据归一化来去除批次效应,如图1,展示了部分样本去除批次效应后的变化。接着利用R语言中的Deseq2软件包对去除批次效应后的数据进行了差异基因筛选。这里,本文设置域值(P < 0.05,
)来筛序差异表达基因,然后使用火山图对结果进行可视化。并通过R语言程序对DEGs进行了功能富集分析(GO与KEGG),得到了与差异基因相关的生物学过程(BP)、细胞组分(CC)、分子功能(MF)以及相关的信号通路。
Figure 1. Distribution of samples after removal of batch effects
图1. 去除批次效应后样本分布图
3.2. 机器学习筛选关键基因
在创建基因筛选模型时,首先通过差异表达分析(如DESeq2)进行初步筛选。然而,单独依赖DESeq2包进行差异表达分析可能无法充分捕捉到所有关键基因,因此研究进一步应用了三种机器学习中的特征选择算法:随机森林(RF)、LASSO算法和GBM (Gradient Boosting Machine)算法。这些算法通过不断训练和调整模型参数,有效地优化了筛选基因的范围,为后面建立预测模型提供了有效帮助。
3.2.1. 随机森林
随机森林(如图2)是一种基于集成学习的算法,通过构建多个决策树并结合投票(分类)或平均(回归)结果,提高模型的准确性和鲁棒性。它具有良好的抗过拟合能力和容错性,能有效处理高维数据与变量间的非线性关系。在基因筛选和疾病分类中,随机森林不仅能识别出具有重要诊断价值的关键特征,还能提供特征的重要性排序,有助于发现潜在的生物标志物和构建稳定的预测模型。
3.2.2. LASSO回归
LASSO (Least Absolute Shrinkage and Selection Operator)是一种改进的线性回归方法,通过引入L1正则化项,实现变量筛选和模型压缩。其目标函数为:
Figure 2. Random forest flowchart
图2. 随机森林流程图
(1)
其中,
是L1正则化项,λ是正则化强度控制参数。当λ较大时,LASSO会
将部分系数
收缩至0,从而实现特征的自动选择。这种机制使LASSO不仅可以降低模型复杂度、防止过拟合,还能提升模型的可解释性,尤其适用于高维数据。在基因筛选或疾病分类等生物信息学应用中,LASSO能够从海量特征中识别出与疾病高度相关的核心基因,为生物标志物的发现和精准医学研究提供了有力支持。
3.2.3. GBM算法
GBM (梯度提升机)算法是一种集成学习算法,通过逐步构建多个弱学习器(通常是决策树),每一步利用前一轮的残差进行学习,逐步提升模型的预测精度。其核心思想是最小化损失函数,通过以下公式更新模型:
(2)
上式中,
是第m轮的预测模型,
是学习率,
是当前训练的基学习器。GBM通过优化损失函数(如均方误差)来训练模型,在每轮迭代中减小预测误差。此外,GBM能够计算特征的重要性,帮助识别关键特征或生物标志物,特别适用于基因筛选和疾病分类等任务。
3.3. 预测模型构建与评估
使用R语言程序构建了基于广义线性模型(GLM)、随机森林(RF)、支持向量机(SVM)、朴素贝叶斯(NBM)、K近邻(KNN)、LASSO回归(LASSO)、极致梯度提升(XGBoost)、多层感知机(MLP) 8种方法的前列腺癌诊断预测模型,同时采用了十折交叉验证,以减少模型训练过程中因数据划分不同而导致的偏差,提升模型的泛化能力和稳健性。
在训练不同的分类器模型后,需要对其性能进行评估,研究主要采用了ROC曲线和AUC值来评估各分类模型的预测能力。通过绘制ROC曲线,可以观察模型在不同判别阈值下的灵敏度(Sensitivity)与特异性(Specificity)变化,并通过AUC值定量评估模型的整体分类能力。AUC值越大,说明模型在区分阳性和阴性样本方面的性能越优。
4. 结果
4.1. 基因筛选
本研究对303例样本进行了分析,包含151个tumor样本和152个normal样本。在进行数据预处理之后,我们使用“DESeq2”包对数据集进行分析,同时设置域值(P < 0.05,
),最后共获得16,833个 差异表达基因,其中1499个差异表达基因下调,1087个差异表达基因上调,DEGs的结果显示在如下火山图和热图中(如图3、图4)。
Figure 3. Differential gene volcano map
图3. 差异基因火山图
注:方格中红色代表正相关,蓝色代表负相关。颜色越深,相关性越强。
Figure 4. Heat map of differential genes (top 50)
图4. 差异基因(前50个)热图
4.2. 富集分析
我们通过R语言程序对DEGs进行了功能富集分析,筛选出了校正后概率P < 0.05的富集途径。从图中可以看出,基因本体论(GO)包含869条目:生物过程(BP)为643条,细胞组分(CC)为89条,分子功能(MF)为137条。我们按照特定顺序排列条目的P值,选择每个过程中的前条记录,并在图5中展示这8条记录。如图5(A)~(C),生物过程主要和肌肉收缩、细胞间黏附、嗜同性细胞粘附等有关,细胞组分主要和细胞外基质、离子通道复合物、肌浆等密切相关,分子功能主要和被动跨膜转运体活性、通道活性、阳离子通道活性等密切相关。
(a) (b)
(c) (d)
Figure 5. Enrichment analysis. (a) BP in GO enrichment analysis; (b) CC in GO enrichment analysis; (c) MF in GO enrichment analysis; (d) KEGG enrichment analysis
图5. 富集分析。(a) GO富集分析中的BP方面;(b) GO富集分析中的CC方面;(c) GO富集分析中的MF方面;(d) KEGG富集分析
KEGG富集结果则包含335个条目,关键差异基因主要在神经活性配体与受体的相互作用、钙信号通路、肌肉细胞的细胞骨架和cAMP信号通路等通路富集明显。如图5(D),我们按照校正后概率P的升序顺序展示了前15个通路。
4.3. 关键基因筛选
在进行差异分析后,我们又分别使用了RF、LASSO回归和GBM算法重新对上面得到的2586个差异表达基因进行了筛选。随机森林通过生成多个决策树并计算每个基因的特征重要性,最终选出最具预测性的特征基因。在模型训练过程中,RF的参数设置为
,通过调整树的数量以优化模型的性能。如图6(a)和图6(b),最终筛选出21个较为显著的基因;LASSO算法则通过调整正则化参数来获得最佳的基因子集,从而提高模型的预测能力并减少过拟合。如图7(a)和图7(b),通过选择
,筛选出32个关键基因;如图8,GBM算法通过十折交叉验证选择出评分最高的13个基因。如图9,RF与LASSO交集为5个标志基因,RF与GBM交集为8个标志基因,LASSO与GBM的交集为3个标志基因。这里为防止筛选过度,研究将上述三个交集合并,最终得到12个关键基因,分别为DLX1、SLIT1、SIM2、TDRD1、HOXC6、RPS17、DHPS、ARL6IP4、RPS10、SEC31B、GSTP1、GOLM1。
(a) (b)
Figure 6. Random forest screening for biomarkers
图6. 随机森林筛选生物标志物
(a) (b)
Figure 7. LASSO screening biomarkers
图7. LASSO筛选生物标志物
Figure 8. GBM screening biomarkers
图8. GBM筛选生物标志物
Figure 9. Wayne’s map
图9. 韦恩图
4.4. 前列腺癌预测模型构建
利用上文筛选出的12个关键标志基因,研究构建了前列腺癌早期预测模型。构建步骤如下:将来源于UCSC Xena数据库的数据划分为测试集和训练集,其中,选取了
的数据作为测试集,用于验证模型的鲁棒性和泛化能力,剩下的数据则作为模型训练集,流程如图10。
Figure 10. Flowchart of ten-fold cross-validation
图10. 十折交叉验证流程图
采用十折交叉验证法对训练集构建了基于广义线性模型(GLM)、随机森林(RF)、支持向量机(SVM)、朴素贝叶斯(NBM)、K近邻(KNN)、LASSO回归(LASSO)、极限梯度提升(XGBoost)、多层感知机(MLP) 8种方法的前列腺癌诊断预测模型。接下来,我们将使用独立的测试集验证八种前列腺癌诊断预测模型的性能。同时,我们将通过计算混淆矩阵(如图11)和AUC值来评估这些模型的准确性和可靠性。根据图12可知,测试集在8种模型下表现效果良好,AUC均在0.7以上。由表1、图12知,测试集在XGBoost模型下表现最优,AUC值高达97%,准确度高达93%。结果表明,基于XGBoost构建的前列腺癌诊断预测模型性能最好,鲁棒性最强。
Table 1. Evaluation results of 8 models
表1. 8种模型评价结果
模型 |
准确度 |
精确度 |
召回率 |
F1分数 |
GLM |
0.762 |
0.758 |
0.768 |
0.763 |
RF |
0.925 |
0.932 |
0.914 |
0.923 |
SVM |
0.911 |
0.925 |
0.894 |
0.909 |
KNN |
0.881 |
0.857 |
0.914 |
0.885 |
Bayesian |
0.918 |
0.909 |
0.927 |
0.918 |
LASSO |
0.733 |
0.756 |
0.724 |
0.739 |
XGBoost |
0.931 |
0.915 |
0.927 |
0.921 |
MLP |
0.776 |
0.833 |
0.690 |
0.755 |
Figure 11. Confusion matrix
图11. 混淆矩阵
Figure 12. ROC curve of test set under 8 models
图12. 8种模型下测试集ROC曲线
5. 讨论
前列腺癌作为男性中最常见的恶性肿瘤之一,近年来其发病率和死亡率的不断上升引发了全球范围内的广泛关注。本研究旨在通过生物信息学和机器学习技术融合,探究前列腺癌与基因表达的关系并建立良好的预测疾病模型,从而为前列腺癌的诊断与治疗提供理论依据。
本研究首先使用生物信息学方法初步筛选出2586个差异基因,包含1499上调基因和1087个下调差异基因,然后对这些差异基因进行GO富集分析和KEGG富集分析。通过GO富集分析发现前列腺癌生物过程主要和肌肉收缩、细胞间黏附、嗜同性细胞粘附等有关等密切相关;通过KEGG富集分析发现其与神经活性配体与受体的相互作用、钙信号通路、肌肉细胞的细胞骨架和cAMP信号通路等密切相关。
研究基于三种机器学习方法(RF,LASSO和BGM)对上述所得的差异基因进行特征选择,最终得到了12个差异基因。同时,对上述12种关键基因构建了8种前列腺癌诊断预测模型。通过测试集验证并绘制ROC曲线、计算混淆矩阵,结果表明8种机器学习模型在诊断上表现出色,准确性和可靠性均较高,其中,XGBoost模型准确率和AUC值分别为93%、97%。
研究中所筛选出的这些基因,在已有文献中已有一定的相关性。Goel等[8]表明DLX1在前列腺癌组织中显著高表达,并且在晚期和转移性前列腺癌患者中呈上调趋势;其高水平与肿瘤侵袭性增强和预后不良密切相关,可作为前列腺癌诊断的潜在非侵入性生物标志物。Wu等[9]的研究表明,SLIT1基因在前列腺癌中的表达升高。SLIT1的高表达可能与前列腺癌的侵袭性及较差预后相关,在激素抵抗性的前列腺癌病例中尤为明显。这一发现提示SLIT1信号通路在前列腺癌进展中具有重要作用。通过基因敲降实验发现,Lu等[10]发现转录因子SIM2在前列腺癌细胞中发挥促肿瘤作用。SIM2在前列腺癌组织中呈上调表达;抑制SIM2可引发广泛的基因表达改变,影响细胞增殖和代谢通路,暗示SIM2的高表达有助于前列腺肿瘤细胞的增殖和分化异常,从而促进前列腺癌的发生发展。Chen等[11]研究指出,TDRD1在前列腺癌患者中异常高表达。TDRD1属于胚系特异基因,在多达68%的前列腺肿瘤中错误表达;其表达水平与雄激素受体驱动的TMPRSS2-ERG基因融合高度相关,是ERG转录因子的直接下游靶基因[11]。由于TDRD1过表达与前列腺癌早期复发风险相关,因此该基因有望用作前列腺癌的诊断及预后指标。通过基因组分析,Luo等[12]发现HOXC6在前列腺癌中异常高表达,其高水平与肿瘤侵袭性增加及治疗后复发密切相关。HOXC6被认为是侵袭性前列腺癌的临床相关生物标志物之一,可能通过与雄激素受体通路相关的转录调控机制促进前列腺癌的进展。Nasr等[13]的综述指出,多种核糖体蛋白在前列腺癌组织中呈现高表达趋势,是驱动肿瘤发生的重要因素。RPS17作为核糖体小亚基蛋白,其过表达可能通过增强肿瘤细胞的蛋白质合成能力,在前列腺癌的发生发展过程中起关键作用,提示其有潜力成为新的肿瘤标志物。Connell等[14]发现RPS10基因在前列腺癌组织和泌尿生殖液中均高水平表达。Proteomics分析显示RPS10在前列腺癌中蛋白水平显著上调。尽管其作为尿液标志物的研究尚不多见,RPS10的过表达可能促进肿瘤细胞的蛋白合成和增殖,在前列腺癌中发挥促进作用。Zhao等[15]发现DHPS (脱氧高胺合成酶)通过催化真核翻译因子eIF5A的假亮氨酸化,在肿瘤细胞生长中扮演角色。虽缺乏针对前列腺癌的直接研究,其他肿瘤研究表明抑制DHPS可减少eIF5A活化,从而抑制癌细胞增殖和迁移[15]。这暗示DHPS可能影响前列腺癌细胞的生长调控,对肿瘤发生发展产生作用。
ARL6IP4被Javier-DesLoges等[16]初步发现参与调控内质网应激(ER stress),而内质网应激是肿瘤细胞在化疗或缺氧条件下存活的重要机制。在前列腺癌中,内质网应激的激活可能通过未折叠蛋白反应(UPR)促进癌细胞存活和耐药性,这对前列腺癌的研究有着一定价值。Stankewich等[17]发现SEC31B作为COPII囊泡蛋白的组成部分,参与内质网和高尔基体之间的蛋白质转运。吴等[18]研究表明,内质网–高尔基体转运过程对肿瘤细胞的生长、迁移和侵袭至关重要。SEC31B可能通过调节细胞内蛋白的分泌和运输,在肿瘤细胞的适应性和生长中起到关键作用,对前列腺癌有一定研究价值。GSTP1是前列腺癌研究中最常见的分子标志物之一。几乎90%以上的前列腺癌标本中检测到GSTP1基因启动子区域的高水平异常甲基化[19],导致基因沉默和蛋白表达缺失。这种启动子高甲基化在良性前列腺组织中罕见,却从前列腺上皮内瘤变(PIN)阶段即出现并贯穿癌变全过程[19]。因此,GSTP1的甲基化状态被广泛用于前列腺癌的早期诊断,是公认的前列腺癌发生早期事件和诊断生物标志物。Varambally等[20]基于表达谱分析发现,GOLM1在前列腺癌组织中显著高表达,且主要由前列腺上皮肿瘤细胞产生。随后研究通过尿液检测证实,前列腺癌患者尿液中GOLM1 mRNA水平明显升高,并可作为前列腺癌的诊断指标,其诊断效能优于血清PSA。此外,近期研究显示GOLM1过表达还能通过激活TGF-β1/Smad2通路促进上皮–间质转化(EMT),在前列腺癌进展和转移中发挥促癌作用,体现了其作为诊断和预后生物标志物以及治疗靶点的潜力[21]。
6. 结论
本研究基于机器学习的方法,成功筛选出了关于前列腺癌的12个关键基因,这些基因在前列腺癌样本与正常样本中呈差异表达,可能作为未来诊断前列腺癌的潜在标志物。同时,构建了8种基于机器学习的前列腺癌预测诊断模型,其中XGBoost模型在测试集中表现最好,AUC高达97%,准确率高达93%,这也进一步表明这些基因与前列腺癌的产生和发展有着密切的联系,可作为前列腺癌早期诊断及研究的潜在靶点。
综上所述,本文结合生物信息学与机器学习技术,完成了对前列腺癌关键生物标志物的筛选,同时构建了对前列腺癌的早期诊断预测的模型,AUC值和准确率分别为97%和93%。本研究虽构建了基于12个关键基因的前列腺癌预测模型,但仍存在数据集规模有限、生物学机制缺乏实验验证、模型可解释性不足及临床因素未纳入等不足。未来研究将扩大数据集,提高模型泛化能力;通过细胞与分子实验验证关键基因功能;引入深度学习及特征选择算法优化模型性能;采用SHAP、LIME等方法增强模型可解释性;整合多组学数据提升诊断精准度;并推进模型临床转化,与PSA等生物标志物结合,开发精准检测工具,以促进前列腺癌的早期诊断与个性化治疗。
NOTES
*通讯作者。