基于深度学习方法预测抗癌药物敏感性

Application of Deep Learning Methods in Prediction of Anticancer Drug Susceptibility

摘要: 近年来，个性化医疗引起研究者们的广泛关注，抗癌药物敏感性预测便是个性化医疗的一个主要挑战。本文将CCLE作为抗癌药物敏感性研究的数据集，选取了不同细胞系上的基因表达数据以及药物敏感性数据。同时我们设计了一种名为PCA Transformer (PCAT)的混合深度学习与机器学习的方法来对抗癌药物敏感性进行预测。首先构造一个PCA模型来提取在不同细胞系上的基因表达数据中的重要变量，使得约5万的基因维度降至500；随后基于降维后的基因表达值建立了一个神经网络Transformer模型来预测药物敏感性，通过均方根误差(RMSE)来评估我们模型的性能，以结果最优的潜变量数量建立的模型作为最终模型。为了验证PCA Transformer的性能，本文将Transformer模型与预测模型随机森林(RF)和支持向量回归(SVR)来进行对比，为了排除降维方法的影响，统一使用PCA进行降维。具体组合包括：PCA Transformer、PCA + SVR、PCA + RF。最后与前人研究方法(ISIRS)的结果进行比较并优化。最终的预测结果看出，对于CCLE中的24种药物，本方法预测得到的平均RMSE为0.7564，有6种药物的RMSE小于0.5 (L-685458、PF2341066等)，有18种药物的RMSE小于1。与其比较的预测方法的平均RMSE分别为：0.8284 (PCA + SVR)、0.8757 (PCA + RF)、ISIRS (0.9258)，体现出本方法有着更强的泛化能力。

Abstract: In recent years, personalized medicine has attracted extensive attention from researchers, and the prediction of anticancer drug susceptibility is a major challenge for personalized medicine. In this paper, CCLE was used as a dataset for anticancer drug susceptibility studies, and gene expression data and drug sensitivity data on different cell lines were selected. At the same time, we designed a hybrid deep learning and machine learning method called PCA Transformer (PCAT) to predict the susceptibility of anticancer drugs. Firstly, a PCA model was constructed to extract important variables in gene expression data on different cell lines, so that the gene dimension of about 50,000 was reduced to 500. Then, a neural network Transformer model was established based on the dimensionality reduction gene expression value to predict drug sensitivity, the performance of our model was evaluated by root mean square error (RMSE), and the model established with the optimal number of latent variables was used as the final model. In order to verify the performance of PCA Transformer, this paper compares the Transformer model with the prediction model random forest (RF) and support vector regression (SVR). Specific combinations include: PCA Transformer, PCA + SVR, PCA + RF. Finally, the results were compared and optimized with the results of previous research methods (ISIRS). The final prediction results showed that for the 24 drugs in CCLE, the average RMSE predicted by this method was 0.7564, 6 drugs had RMSE less than 0.5 (L-685458, PF2341066, etc.), and 18 drugs had RMSE less than 1. The average RMSE of the prediction method is 0.8284 (PCA + SVR), 0.8757 (PCA + RF) and ISIRS (0.9258), respectively, indicating that the proposed method has stronger generalization ability.

1. 引言

癌症是基因引起的疾病，当调控细胞生长的基因发生突变或损坏时，使得细胞失去控制，持续的生长及分裂而产生肿瘤。基因学的发展驱使许多研究者去发掘不同人所得癌症的不同之处，大量研究现象表明，即使是患有相同癌症类型的不同病人，对于同一种药物或者治疗方法，可能也会出现不同的药物反应症状。传统的治疗方式常常对患有相同类型，处于相同阶段癌症的病人采取一些类似的治疗，而由于某一类型癌症的异质性以及遗传多样性都会使得这些治疗并不那么有效，甚至对于不同癌症的子类也存在效果不明显的情况[1]，基于患者特异性分析的抗癌治疗引起研究者们的广泛关注。

国内外有许多学者投入到预测抗癌药物的反应能力的研究中，尽管原始组织或肿瘤样品的基因组有差异[2] [3]，具有异质基因组背景和基因表达的培养癌细胞系仍是研究药物活性的分子基础[4]并发现癌症生物学中新抗癌药物的基础材料。几项大规模高通量筛选工作已对一组体外细胞系的基因组信息及数百种化合物的药物敏感性概况进行了分类。在过去几十年，大量研究对基于这些工作得到的基因信息进行建模，并产生了各种预测模型用来预测药物敏感性。美国国家癌症研究所(the National Cancer Institute)在20世纪80年代发布的NCI-60数据即为几项高通量筛选工作之一，NCI-60数据在很长一段时间内都被用来测试候选药物抑制癌细胞增殖的能力，该基因组数据包含了人类跨越9个不同组织的60株癌症细胞系[5]。GDSC (The Genomics of Drug Sensitivity in Cancer)是由Yang等人在2012年开发的一组更大的数据集，描述了639种人类癌细胞系的基因组概况及对138种药物的药物反应数据进行分类，目的是为了鉴定癌细胞中药物敏感性的基因组生物标志物[6]。GDSC数据集将癌细胞系的抗癌药物敏感性数据与来自于COSMIC (the Catalogue of Somatic Mutations in Cancer)数据集结合。CCLE和GDSC数据集中的细胞系源自多种人类癌症组织，例如肺，乳腺和肾脏。CCLE和GDSC数据集都具有丰富的基因组信息数据，包括基因表达，DNA拷贝数，癌基因图谱突变等，其中CCLE包含了1000多个样本的基因特征，以及24种药物在504个样本上的敏感性值。基于这些数据集已有非常多的模型被提出，Staunton等人基于NCI-60的基因表达谱建立了预测模型[7]，他们使用加权投票方法将每个细胞系分类为在某些药物治疗下敏感或耐药；John等人利用结肠癌细胞的基因谱来预测和区分对多种化疗药物的反应[8]；Aben等人发现很难解释建立的模型，于是提出了一种名为TANDEM的算法，该算法使用由“上游”类型即突变、拷贝数、甲基化和癌症类型组成的数据和“下游”数据即基因表达数据来预测药物反应，最大限度地提高药物反应模型的可解释性[9]随着数据的更新，在样本量一定的情况下，基因数据维度在不断的变大，这也为预测造成了较大的阻碍。Riddick等人在2011年提出建立随机森林计算变量重要性值来筛选变量[10]，而在变量维数不断提升的情况下，随机森林计算速度与其他前沿模型有比较大的差距。Jordi等人在选择变量时采用了范剑青提出的Sure Independence Screening (SIS) [11]方法，该方法通过计算变量与药物敏感性之间的相关性系数去筛选变量。为解决皮尔逊相关系数选择变量受离群点的影响，Sure Independent Ranking and Screening (SIRS) [12]也被用于选择变量的使用。相关性筛选与随机森林筛选始终依赖于药物反应，且不能最大限定的筛选出重要的变量，也很难对最终预测有很大的提升。为了更好的提升预测及降维能力，Chiu等人采用了深度自编码以及深度神经网络来对基因进行降维预测[13]。而在一些小样本数据集上，深度神经网络较传统的机器学习模型更容易产生过拟合的情况，泛化能力并不突出。

本文以前人的研究为基础，为最大程度的提升预测能力和计算速度，将PCA模型与深度学习相结合，应用至CCLE数据集预测药物敏感性。

2. 研究方法与数据

2.1. 抗癌药物敏感性数据

在本文中，我们提出了一个深度学习与机器学习混合的方法，称为PCA Transformer，用于预测在细胞系上的药物敏感性。我们从CCLE数据集获取细胞系的基因表达谱和药物敏感性数据，其中的基因表达信息与药物敏感性数据作为模型的输入。

CCLE数据集包含1037个样本的细胞系基因表达数据、基因突变和拷贝数信息数据等。此外数据集中包含了24种抗癌药物在504个细胞系上测得的敏感性数据。研究发现，Koras等人于2020年对药物敏感性预测问题做了全面的评估[14]，评估认为基因表达是最有力的预测特征，因此我们选择了基因表达数据作为药物敏感性的预测变量。CCLE中抗癌药物敏感性主要包含三种数据：1) 半抑制浓度，是表示抗癌药物敏感性的一种数据类型，其表示抑制50%的细胞生长所需要的药物浓度；2) 半最大效应浓度(Concentration for 50% of Maximal Effect)是指能引起50%最大效应的浓度；3) Activity Area，药物剂量曲线上方的面积。根据指标意义我们可以知道，Activity Area越大，药物的敏感性越高，反之则越低。本文选取了Activity Area作为药物敏感性的表示指标。CCLE数据集数据可在www.broadinstitute.org/ccle/获取。

2.2. 主成分分析(PCA)

主成分分析(Principal Component Analysis, PCA)是一种常用的数据降维技术，本方法的思想可通过图1反映，它可以通过线性变换将原始数据转换为一组各维度线性无关的表示，转换后的这组变量称为主成分(Principal Components)。主成分分析的目标是减少数据集的维数，同时保留数据集中对方差贡献最大的特征，从而帮助理解数据的内在结构，同时减少计算量和噪声的影响。主成分分析的主要步骤包括：

1. 数据标准化：由于PCA对数据的尺度非常敏感，因此首先需要对原始数据进行标准化处理，使得每个特征均值为0，方差为1。

2. 计算协方差矩阵：标准化后的数据，计算其协方差矩阵。协方差矩阵的每一个元素是任意两个特征之间的协方差，代表了它们之间的线性相关程度。

3. 计算协方差矩阵的特征值和特征向量：特征值的大小表示了对应特征向量方向上数据变化的程度(即方差)，特征值越大，说明在该方向上的数据变化越大，即信息量越多。

4. 选择主成分：根据特征值的大小，选择前k个最大的特征值对应的特征向量。这些特征向量就是新的特征空间的方向，也就是主成分。通常选择主成分的数量k会远小于原始数据的维度，以达到降维的目的。

5. 将原始数据转换到新的特征空间：使用选定的主成分(特征向量)作为新的基，将原始数据投影到这个新的特征空间，得到降维后的数据。

Figure 1. Principal component analysis visualization

图1. 主成分分析可视化

2.3. Transformer模型

Transformer模型是一种深度学习模型，最初由Vaswani等人在2017年的论文《Attention Is All You Need》中提出，用于处理序列到序列的任务，如机器翻译。它的核心思想是使用自注意力机制(Self-Attention)来处理数据，该方法原理如图2所示，这使得模型能够并行处理序列中的所有元素，而不是像循环神经网络(RNN)那样按顺序处理。在自注意力机制中，模型会计算序列中每个位置与其他位置之间的关联程度(即注意力权重)，然后根据这些权重加权组合序列中的信息，从而得到每个位置的新表示。Transformer模型在自然语言处理(NLP)领域取得了巨大的成功，特别是在机器翻译、文本摘要、问答系统等任务中。此外，它也被扩展到其他领域，如图像处理和语音识别。随着模型的不断发展，出现了许多变体，如BERT、GPT、T5等，它们在各自的任务中都取得了显著的成就。

Figure 2. Transformer model

图2. Transformer模型

2.4. PCA Transformer设计

Figure 3. PCA Transformer flow chart

图3. PCA Transformer流程图

为了满足个性化医疗的需要，我们要利用基因信息数据来对抗癌药物敏感性进行预测。对于不同的抗癌药物，在不同的病人上显现出的效果是不相同的，如何选择合适的模型对其建模才能有好的预测效果是主要的问题。另外表达数据中存在上万种基因，如何从中找到较为重要的基因以减轻我们训练的压力，也是本研究需要解决的另一个问题。在本文中，我们基于本节所介绍的模型设计了名为PCA Transformer的方法来解决降维与预测这两个在抗癌药物敏感性预测建模中所需解决的问题。

在2.1节中我们指出了选择的数据，即在基因信息中我们选取基因表达数据，同时建模之前我们对基因表达数据与药物敏感性数据进行了预处理。CCLE中共有24种抗癌药物，对于每种药物，每个样本所测得的Activity Area即为一个细胞系。我们将抗癌药物数据和基因表达数据进行组合，选取二者共有的细胞系数据，每个药物平均包含五百多个样本，其中每个细胞系包含4万到5万的基因表达数据。PCA Transformer的第一个工作就是对维度大约5万的基因表达数据进行降维提取信息。我们利用PCA方法对这些基因表达数据进行一个无监督式的学习降维，并得到相应的潜变量，这些潜变量最大限度的涵括了原先高维数据中绝大部分信息。在借助PCA降维之后，我们用Transformer模型对提取出来的潜变量关于药物敏感性值建立深度学习模型，以RMSE (见3.1节)作为评价指标，利用模型中的注意力机制来对潜变量进行学习，最终通过学习好的模型对药物敏感性进行最终预测。

3. 实验结果

3.1. 预测效果的评价指标

有许多评价指标来评估回归模型的模型误差，例如均方误差(MSE)，均方根误差(RMSE)，绝对均方误差(MAE)等等。在此之中，RMSE对异常大的误差更加敏感，通常我们更加倾向于避免重大错误的发生，所以RMSE更适合用于度量药物敏感性模型。因此我们选择均方根误差(RMSE)作为我们模型的误差评价指标，具体见公式(1)中所写。

$\text{RMSE}=\sqrt{{\displaystyle \sum {\left(y_i-{\widehat{y}}_i\right)}^2/N}}$ (1)

其中N是测试集的大小，和是药物敏感性数据样本和相应第i个输入数据的预测值。为了对模型进行评估，我们需要用交叉验证来对模型进行检验，检验结果表示我们方法的好坏。留一法和K折交叉验证是两种主流的交叉验证方法，留一法相对于K折交叉验证需要更长的运算时间，因此本文所提出的方法以及比较方法均采用K折交叉验证来进行评估。K折交叉验证随机地将样本平均分为k组，其中的k − 1组样本作为训练集，通过对k − 1组样本建立模型来预测余下的1组样本得到该组的均方根误差。重复k次后我们就可以得到k个测试误差。最终K折交叉验证的测试误差为：

$C{V}_{(k)}=\frac{1}{k}{\displaystyle \sum _{i=1}^k{\text{RMSE}}_i}$ (2)

使用K折交叉验证能有效地提高模型的稳定性及泛化能力，避免数据的过拟合，且更为准确的评估模型的预测能力。本文实验将k设为10，即十折交叉验证来对模型进行检验。

3.2. PCA Transformer性能的评估

在统计学习中，有许多性能卓越的预测模型在过去几十年中被提出，其中最为流行且有效的模型当属支持向量回归(SVR)与随机森林(RF) [15] [16]。在药物敏感性预测领域中，二者也是研究及应用者最为常用的模型，过去大多药物预测方法基于这两种模型，并被不断被验证是非常有效的预测药物敏感性的方法。

本文将Transformer与上述两种模型进行预测能力的比较，以此来验证Transformer高效的预测能力，此外我们将降维方法固定为PCA降维，仅针对预测能力进行比较。对于三种模型，我们都对最后的结果计算RMSE，此外我们都用十折交叉验证来进行模型的评估。三种模型关于24种药物建模得到的预测结果见表1，图4用直方图直观的展示了三者的差别。由结果可知SVR与RF最终预测得到的24种药物平均RMSE为0.8284和0.8757。与Transformer结果比较可以看到，Transformer模型下有23种药物的预测结果都要好于另外两种模型，仅有1种药物略微低于另外两种模型，不论从具体药物还是总体来说，Transformer都有一个很好的性能提升。

本文使用的数据集相较于前人使用的数据集在变量个数上提升了一倍，这也为我们的结果预测加大了难度。本节将PCA Transformer得到的结果与安彪使用的ISIRS方法[17]的结果进行了比较。图4展示了三者的结果，我们可以看到，PCA Transformer相比ISIRS在RMSE结果上降低非常大，平均结果PCA Transformer比ISIRS低0.17，这对于RMSE来说是个非常大的提升。同时我们的结果是在新数据集中得到的，在变量数量翻倍的情况下仍然有着较大的提高，可以说PCA Transformer是个非常好的预测方法。此外同样在SIS情况下的预测，利用PCA Transformer进行预测很好的弥补SIS相对于ISIRS的劣势，由此可见PCA Transformer对于药物敏感性预测具有非常好的效果。

Table 1. Comparison of PCA Transformer, ISIRS, PCA + RF and PCA + SVR for 24 drug modeling predictions

表1. PCA Transformer、ISIRS、PCA + RF和PCA + SVR对24种药物建模预测的结果比较

药物	PCAT	ISIRS	PCA + RF	PCA + SVR	药物	PCAT	ISIRS	PCA + RF	PCA + SVR
AEW541	0.5770	0.9849	0.6941	0.6593	Irinotecan	1.0887	0.7036	1.1609	1.1445
Nilotinib	0.4584	0.8602	0.5679	0.4623	Topotecan	1.1198	0.8468	1.2008	1.2098
17-AAG	1.0263	1.0276	1.126	1.0558	LBW242	0.5292	0.9808	0.7692	0.7230
PHA-665752	0.4110	1.0531	0.5669	0.5159	PD-0325901	1.4513	0.7556	1.7337	1.7119
Lapatinib	0.8734	0.9466	0.9361	0.9099	Paclitaxel	1.4694	0.8832	1.8994	1.5238
Nutlin-3	0.4229	0.9586	0.5351	0.5236	AZD6244	1.1553	0.8390	1.2388	1.2249
AZD0530	0.7775	0.9607	0.8274	0.7796	PLX4720	0.8344	0.8649	0.8759	0.9216
PF2341066	0.4485	0.9230	0.506	0.4740	RAF265	0.8796	1.0183	0.9913	0.8906
L-685458	0.4039	0.8373	0.5846	0.5730	TAE684	0.7689	0.9818	0.9756	0.8572
ZD-6474	0.6543	1.0149	0.6904	0.6688	TKI258	0.5207	1.0291	0.6519	0.5626
Panobinostat	0.7027	0.9597	0.7505	0.7622	Erlotinib	0.5794	0.9476	0.5358	0.5481
Sorafenib	0.4674	0.9397	0.5279	0.5005	PD-0332991	0.5329	0.9022	0.6715	0.6777

Figure 4. Comparison of drug sensitivity prediction outcomes in PCA Transformer, ISIRS, PCA + RF, and PCA + SVR

图4. PCA Transformer、ISIRS、PCA + RF和PCA + SVR药物敏感性预测结果比较

4. 结论与展望

预测特定癌症类型对治疗的反应能力是个性化医学的主要目标之一。基于NCI-60、CCLE和GDSC等基因组数据，研究者们提出了各种模型方法来提高药物反应的预测能力或是提高模型解释能力。本文提出了一个名为PCA Transformer的方法，利用CCLE数据，将基因表达作为预测变量来预测药物反应值(Activity Area)，极大的提升了预测能力。从结果来看，CCLE中24种药物，PCA Transformer建模评估得到的结果较其它模型有很大的提升。另外在模型上，PCA虽然能很好地提取数据中的重要信息，提升预测能力，但因为其无监督的学习方式，使得模型缺乏了一定的解释性，另一方面在于本文训练过程中没有融合更多的生物信息，包括药物结构信息、基因突变信息等，单纯利用基因表达信息作为特征略显单薄，从这一方面来看PCA Transformer预测结果也有一定的提升空间。进一步的研究可以从两个方面着手。从数据上，我们可以利用更多能提供药物预测能力的信息，如基因突变、拷贝数等基因数据，显然当有了更多的信息之后，自然就能在预测能力上更进一步，另外我们也可以将药物结构信息包括进来，这样不再单一的对基因数据进行建模，而是与药物信息结合，也有可能对预测能力有较大的提升。从方法上看，对于降维方法，若是包含其余数据可能会有稀疏性的影响，可以考虑稀疏自编码等针对稀疏数据的降维方法，此外也可以考虑VAE等其它自编码器来提升基因表达数据的降维效果。对于如何利用其他信息，如结构信息，可以考虑在进一步的研究中使用基于注意力机制的模型，目前流行于NLP领域并被广泛用在药物结构翻译领域的BERT模型即为注意力机制模型的一种。此外，对于整体的考量，如考量基因之间的关系以及药物与基因的关系来进行药物敏感性的预测，类似图神经网络等图模型都可以考虑应用在将来的研究中。