王江翔, 肖清泉
(贵州大学大数据与信息工程学院, 贵阳 550025)
国内人口的乳腺癌发病率不容小觑,平均每年有大约30 万的女性乳腺癌患者和少数男性患者[1]。乳腺癌是乳腺上皮细胞在多种因素作用下,发生增殖失控的现象。 患者早期的症状为乳房肿块、乳头溢液、腋窝淋巴结肿大等症状,晚期可能发生癌细胞扩散转移,感染其他器官,直接威胁患者的生命。 随着医疗新兴技术的出现与发展,全球乳腺癌的发病率得到了初步的控制,但是乳腺癌发病机制和抗乳腺癌药物的关系依然是近年来研究的重点和难题。目前,在药物研发中,为了节约时间和成本,通常采用建立化合物活性预测模型的方法来筛选潜在活性化合物。
本文针对传统的回归预测模型对化合物活性预测的效果较差的问题,采用粒子群算法分别对4 种传统预测模型进行优化。 通过对比分析模型优化前后的均方误差、平均绝对误差等指标,验证了粒子群算法带来的模型预测性能的提升,可为回归预测模型的优化提供思路。
本文采用粒子群算法分别改进了SVR 算法、Random Forest 算法、XGBoost 算法和LightGBM 模型,用于构建生物活性的定量预测模型。 同时,以验证集MSE、MAE、拟合度作为评价指标,对比分析各个模型的预测效果,以验证所改进模型的性能。 具体流程如图1 所示。
图1 整体设计流程图Fig. 1 Overall design flow chart
粒子群算法[2]来自对鸟群捕食行为的研究,通过群体中个体之间的协作和信息共享来寻找最优解,该方法简单、且易于实现。 粒子群算法中每个粒子都关注2 个点—位置和速度,在寻找到自己的最优位置和速度时,将其与种群中的其它粒子进行信息共享,其它粒子共享到信息后会以此调整自身的最优位置与速度。 这样随着越来越多的位置与速度信息的出现,就能通过计算得到一个全局最优解。
其中,相较于位置属性,粒子对速度属性更加依赖,粒子根据自身前一时刻的速度来思考自身能够到达的位置,c1就是调节此次位移的步长,c2调节粒子靠近群体的位置的步长。
支持向量回归(SVR)模型[3]是以支持向量机(SVM)为基础的,两者间的区别在于SVR 是使得最远距离的样本点尽量靠近拟合线,SVM 是使得最近距离的样本点尽量疏远拟合线。 SVR 和SVM 对比示意如图2 所示。
图2 SVR 与SVM 对比示意图Fig. 2 Schematic diagram of the comparison between SVR and SVM
SVR 模型将变量的输入映射到特征空间,使得非线性的数据可以通过回归进行线性拟合,以此得到的构造决策函数如下:
SVR 在SVM 的中心线性函数设置了一个宽度为ε的容忍偏差,落到这个宽度内的样本不做损失计算。ε的损失函数计算公式如下:
此外,引入松弛变量ε、∊以适应数据在容忍偏差之外的松弛程度不同的需求,公式如下:
其中,C表示惩罚因子。 引入拉格朗日函数,将一般问题转化为对偶问题,公式如下:
对于非线性回归问题需要对映射空间进行一个转换,转化为一个可以进行线性回归的更高维度的特征空间,这个过程主要依赖于核函数,其高斯核函数公式如下:
映射到高维度的空间以后就能使用回归函数进行回归,回归函数如下:
粒子群算法的引入主要是为了对SVR 模型中惩罚参数C和高斯核函数的宽度系数g进行优化,寻找C、g的最优解。 SVR 对样本数据采用惩罚机制,惩罚力度的大小就由参数C控制,在容忍偏差之外的惩罚力度要适当,因为惩罚力度过大,会使得模型泛化能力较弱,而惩罚力度过小,会使得模型误差大,精度不高。高斯核函数的宽度系数g的大小表示样本之间的相关程度,g的值越大、样本的相关性就越大,g的值越小、样本的相关性就越小。g的值如果太大、模型精度会被影响,g的值如果太小、模型泛化性能力就不足。
因此,优化后的SVR 模型通过使用粒子群算法对惩罚参数C和高斯核函数的宽度系数g进行算法寻优,找到最优解,以此大幅度增加模型的精度以及泛化能力。
Random Forest 算法[4]是为了整合多种决策树算法而提出来的一种集成学习算法,其中的每棵决策树都是使用bootstrap 技术随机抽取样本构建的,每棵决策树的分裂点是从随机特征子空间中随机选择而来。 这样,每棵决策树都是随机生成,涵盖多种可能性,Random Forest 算法在其中通过投票的方式选择最好的预测结果。
决策树是一种分类器,通过模仿树木枝干分枝,从树的根部进行分裂,选取适当的分裂准则,在满足条件时进行分裂,直至生成一颗完整的树状,多个分枝就代表了多个路径规则,在预测数据时,依靠这些规则对数据进行分类和预测。 决策树的生成算法常见的有ID3、C4.5、CART,其中ID3、C4.5 基于香农熵进行分裂,CART 是基于基尼系数分裂。 香农熵的公式表达式如下:
香农熵可以用来计算信息增益的大小,ID3 算法就是通过选择信息增益最大的属性进行决策树生成的。 信息增益的计算公式如下:
但是,ID3 算法对于连续变量的选择性很差,C4.5 算法针对这一算法缺陷引入分裂信息比率因子,数学公式如下:
引入分裂信息比率因子后的分裂准则公式如下:
这一准则计算方法将连续变量进行离散化处理,使得决策树对连续变量的选择性变好。 而CART 采用一种全新的方法,把基尼系数衡量的信息纯度作为分裂准则,基尼系数公式如下:
决策树在建树时根据基尼系数最小进行分裂。
XGBoost 算法[5]是通过串联不同的知识点,合理高效地解决回归、排序等复杂问题,是近年来机器学习领域应用最广泛的算法之一。 boosting 是数学建模中运用最广泛的方法,就是通过集成不同的模型来得到一个高精度的模型,但是boosting 每次计算需要迭代几千上万次,因此效率低下。 XGBoost是在boosting 的基础上进行改进,XGBoost 是一种向前特征的算法,其模型公式为:
其中,k为树的棵树;f表示回归树,并且f(x)=wq(x),q:Rm→T,w∈RT;q为每一颗树的结构;T表示树中的叶子树;一个q对应一个叶子权重w;τ表示全部回归树组成的函数空间。 在实际计算当中,将规则安排到叶子里,然后由w进行求和, 得到最终的预测值。XGBoost的整个建模过程如下:
设定要预测模型:
设置目标函数为:
每个样本落在叶子节点上:
其中,T为t棵树叶子节点的的总数,mj是j个叶子的权重。
进一步,得到新的目标函数:
令:
对函数进行优化:
增加的切分点的增益:
LightGBM[6]是以GBTD 为基础的一种算法,而GBTD(Gradient Boosting Decision Tree)是一种通过迭代构建模型的算法,GBDT 具有梯度提升和决策树的特点,具有精度高、速度快、不容易过拟合的优点。 LightGBM 弥补了GBDT 处理海量数据时速度较慢的问题。 LightGBM 是一种以直方图决策树为基础的算法,其基本过程是:将特征划分成s个整数,然后生成一个宽度为s的直方图,在输入数据时,离散值被累积起来作为直方图的索引,在遍历一次数据后,所有累计的量全被直方图统计起来,然后根据生成的索引值进行遍历找到最佳的切分点。 这个算法的优点在于减小了存储量,使得运算速度加快,并且降低了运算成本,XGBoost 需要计算一次效益后,再遍历每一个特征,然而采用直方图算法的LightGBM只需要计算s次。 直方图算法如图3 所示。
图3 直方图算法Fig. 3 Histogram algorithm
迄今为止,LightGBM 已经普遍用于机器学习领域。 大量的数据表明,LightGBM 已比目前应用的Boosting 工具拥有更高的准确率和计算效率,LightGBM 比XGBoost 的计算速度更快,效率和精度也更高,LightGBM 采用直方图算法可以有效地解决XGBoost 算法把所有特征加入内存进行排序所造成的内存消耗。
本文数据主要源于第十八届华为杯数学建模D题,数据包括1 974 个化合物对ERα 的生物活性数据,以及729 个分子描述符信息(变量)。 由于数据中存在部分对生物活性影响不大的变量,本文通过相关性分析筛选出20 个对生物活性最具有显著影响的变量。 20 个变量之间相关度分布图如图4(a)所示。 由图4(a)可见,20 个主要变量的相关系数矩阵中,有部分变量的相关性过高,对模型的预测效果具有一定的影响。 因此,为了使得模型预测性能更佳,本文进一步研究了变量之间是否具有高度的线性关系或者非线性关系,去除了部分相关系数较大的变量,最后得到14 个变量,如图4(b)所示。 这14 个变量之间的整体独立性更高,对模型预测更加有利,故将这14 个变量的值作为建立生物活性定量预测模型。
图4 变量间的相关度分布图Fig. 4 Correlation distribution between variables
3.2.1 均方误差
均方误差(MSE) 是观测值与真值偏差的平方和观测次数n比值的平方根,当对某一变量进行多次的测量时,取这一测量结果误差的均方误差,称均方误差,公式如下:
3.2.2 平均绝对误差
平均绝对误差(MAE),又叫平均绝对离差,是所有单个观测值与算术平均值的偏差的绝对值的平均。 平均绝对误差可以避免误差相互抵消的问题,因而可以准确反映实际预测误差的大小。 绝对误差的平均值公式如下:
3.2.3 拟合优度
拟合优度(R2) 是指回归线对预测值的拟合程度。 度量拟合优度的统计量是可决系数(亦称确定系数)R2。R2最大值为1。R2的值越接近1,说明回归直线对观测值的拟合程度越好。 反之越差。
本文首先对1 974 个化合物进行样本划分,为了在保证精度的同时让模型具有更好的泛化能力,采用随机抽取样本的方法对数据进行划分,其中以90%的化合物为训练集,10%的化合物为验证集。构建了4 种算法模型,分别是SVR 算法、Random Forest 算法、XGBoost 算法、LightGBM 算法[6]。 对这4 种算法模型采用粒子群算法进行优化,对比分析各模型的最终实验结果,选择预测性能最好的模型作为本文的最终预测模型进行测试集的预测。
3.3.1 模型预测值对比分析
将划分好的训练集样本数据送入优化模型进行训练回归,得到训练好参数的优化模型,使用验证集对模型进行验证。 为了验证粒子群算法给各个原始模型带来的性能提升效果,本文对原始模型进行相同的训练与验证,将得到的预测值与样本真实值做对比,如图5 所示。
图5 4 个模型优化前后的预测结果与真实值对比Fig. 5 Comparison of the prediction results with the actual value before and after optimization of the four models
从图5(a)中可以看到,原始SVR 模型的预测结果拟合准确度并不高,预测值与真实值的误差较大。 这在很大程度上是由于惩罚参数C和高斯核函数的宽度系数g的优化程度不高所导致。 而用粒子群算法改进后的SVR 模型,C和g数值经过优化后,将明显提升模型预测准确率。 优化后SVR 模型的预测值与样本真实值对比如图5(b)所示。 由图5(b)可看出,粒子群算法优化后的SVR 模型的预测结果与真实值拟合程度更高,这也表明优化后的SVR 模型性能更好,误差更小。 同时,图5(c)至图5(h)也表明了经过粒子群算法优化后的模型预测效果更好,整体性能更佳。
3.3.2 拟合度对比分析
本文对模型拟合度进行了可视化分析,如图6所示。 从图6 中可以看出,优化后SVR 模型的预测拟合度在迭代回归5 次以后基本收敛,预测值与真实值的拟合度维持在79.52%。 而优化后的Random Forest 模型的误差较小,预测值拟合度稳定在81.29%。同时,优化后XGBoost 模型在训练过程中出现“局部最优”现象,拟合度维持在了77.20%,通过粒子群算法对参数进行寻优,使其训练能够自主跳出局部最优,拟合度最终达到84.52%。
图6 各个模型优化后的拟合度Fig. 6 Fitting degree of each model after optimization
相比之下,优化后的LightGBM 模型拟合误差更小,不断通过参数优化参数,跳出局部最优,拟合度最终稳定在84.68%。
3.3.3 预测误差对比分析
为了进一步验证粒子群算法给预测模型带来的影响,本文分别绘制了4 个优化后的模型的预测值与真实值的误差曲线图,如图7 所示。 由图7 可见,优化后的模型预测误差大多数都低于改进前的模型,预测效果得到了较好的改善。 其中,SVR 模型和XGBoost 模型的性能提升较大,整体误差基本都有所降低。 而相比之下,LightGBM 模型的整体预测误差更低,预测效果更好。
图7 各个模型的预测值与真实值的误差曲线Fig. 7 The error curve between the predicted value and the true value of each model
3.3.4 评价指标对比分析
本文给出了4 个模型评价指标MSE(均方误差)、MAE(平均绝对误差)、R2(拟合度),见表1。 由表1 可看到,优化后LightGBM 模型的MSE和MAE值分别为0.283 3和0.334 9,均低于其它3 个模型。 同时,优化后LightGBM 模型的R2值达到了最高的84.68%。
表1 优化后的模型的评价指标Tab. 1 Evaluation indicators of the optimized model
本文针对传统回归预测模型对化合物活性预测的效果较差的问题,分别采用粒子群算法优化了SVR 模型、Random Forest 模型、XGBoost 模型和LightGBM 模型。 对比分析各个模型的拟合度、MSE、MAE等评价指标的结果表明,优化后的LightGBM 模型对化合物活性的预测效果更好。 同时,通过可视化优化前后模型的预测结果和预测误差,直观地展示了粒子群算法优化带来的模型效果的提升,可为回归预测模型的优化方法提供思路。