Auto-HMM-LMF：基于特征选择的自动编码器和隐马尔可夫模型药物反应预测方法

摘要

背景

预测癌细胞株对特定药物的反应是个性化医学中的一个基本问题。由于药物反应与癌细胞中的基因组信息密切相关，因此数百个人类癌症细胞系的一些大型小组是用基因组和药物基因组数据组织起来的。虽然已经开发了几种方法来预测药物反应，但在实现准确预测方面存在许多挑战。本研究提出了一种新的基于特征选择的方法，称为Auto-HMM-LMF，用于准确预测细胞系药物相关性。由于用于预测药物反应的特征空间具有很大的维数，Auto-HMM-LMF将重点放在特征选择问题上，以利用具有重要贡献的输入子集。

结果

本研究提出了一种基于特征分配和隐马尔可夫模型的突变数据特征选择新方法。此外，我们使用自动编码器模型来选择基因表达和拷贝数变化数据的特征。在选择特征后，应用logistic矩阵分解模型预测药物反应值。此外，通过与最强大的特征选择方法之一集成特征选择方法（EFS）的比较，我们表明基于本文所介绍的所选特征的预测模型在药物反应预测方面的性能要好得多。使用两个数据集，癌症药物敏感性基因组学（GDSC）和癌症细胞系百科全书（CCLE）来表明该方法在未发现的患者细胞系中的效率。对该模型的评估表明，Auto-HMM-LMF可以提高最新算法结果的准确性，并且可以为logistic矩阵分解方法找到有用的特征。

结论

我们描述了Auto-HMM-LMF在探索头颈癌新候选药物中的应用，表明所提出的方法在药物重新定位和个性化药物方面是有用的。Auto-HMM-LMF方法的源代码可在https://github.com/emdadi/Auto-HMM-LMF网址.

背景

个性化医学的计算模型使得在基因组信息的基础上理解癌细胞系成为可能。这一知识使得通过测量药物反应向不同类型癌症患者推荐个性化治疗成为可能。针对每种癌症类型，如乳腺癌、肺癌、卵巢癌和脑癌，已经开发出许多有效的抗癌药物。例如，多西紫杉醇、紫杉醇、卡铂、顺铂、长春瑞滨和埃里布林只是用于治疗乳腺癌的药物的几个例子。由于癌症治疗的药物反应取决于多种因素，如患者的基因组特征，这一过程在癌症治疗中是一个复杂的问题。这些挑战引发了对人类癌症细胞系和各种抗癌药物的大规模实验。例如，两个数据集癌症药物敏感性基因组学[1]（GDSC）和癌症细胞系百科全书[2]（CCLE）是根据已建立的抗癌药物对不同癌细胞株的药物敏感性数据创建的。已经提供了癌症细胞系的各种遗传特征，如基因表达谱、拷贝数改变、单核苷酸突变和甲基化数据。通过使用这些数据库，机器学习算法通过以具有统计意义的方式整合不同来源的数据，越来越多地应用于药物反应预测。

提出了几种基于推荐系统的模型来预测药物反应。Wang等人采用相似正则矩阵分解（SRMF）方法，利用细胞系中的基因表达谱和药物的化学结构预测细胞系的抗癌药物反应。他们指出雷帕霉素（一种mTOR抑制剂）可能是一种治疗非小细胞肺癌的新药物[三]. Suphavilai等人开发了一个模型，称为使用推荐系统（CaDRReS）的癌症药物反应预测，以了解药物和细胞系对潜在空间的预测。此外，他们还演示了如何利用已获得的特征探索药物机制和药物途径关联[4]. Emdadi等人提出了基于逻辑矩阵分解方法的DSPLMF方法，用于预测抗癌药物反应。DSPLMF致力于发现细胞株和药物的显著特征和潜在载体，以计算细胞株对药物敏感的概率。他们使用获得的潜在载体来识别癌细胞系的亚型和药物途径相关性[5].

从许多遗传候选特征中确定最佳特征子集是预测药物反应分类模型的关键问题。因此，许多算法提出了使用不同的方法进行特征选择。Xu等人提出了基于特征选择的AutoBorutaRF方法来预测药物反应。该方法首先构建了一个自动编码器网络，并使用Boruta算法[6]选择应用RandomForest分类器预测药物反应的重要特征[7]. Dong等人提出了一种称为支持向量机递归特征消除（SVM-RFE）的模型，该模型使用递归特征选择和SVM分类器的包装方法来预测药物反应[8].

本研究提出了一种基于特征选择的药物反应预测方法，称为Auto-HMM-LMF，用于有效预测细胞系药物相关性。纳入基因表达谱、拷贝数改变、单核苷酸突变、细胞系组织类型信息和药物化学结构信息。采用基于自动编码器和隐马尔可夫模型-多项式混合模型的两种策略来选择输入信息的基本特征。自动编码网络应用于基因表达谱、拷贝数变化数据。同时，将隐马尔可夫模型和多项式混合模型应用于突变数据。使用十倍交叉验证对Auto-HMM-LMF方法进行适当评估，以将其与最先进的方法进行比较。结果表明，该算法对测试数据集的性能是优越的。此外，通过与集合特征选择方法（EFS）的比较，我们表明Auto-HMM-LMF方法中两种考虑的特征选择策略可以选择适当的特征，从而显著改善预测结果。

方法

本文提出了一种新的方法（Auto-HMM-LMF），通过结合并有效地使用特征选择方法来有效预测细胞系药物关联。Auto-HMM-LMF算法的主要方案如图所示1在第一步中，使用了两种策略来选择输入数据的重要特征。基于自动编码网络的特征选择方法被应用于细胞系的基因表达谱，并且相似矩阵（Sim_经验)使用选定的特征构造。类似地，相似矩阵（Sim_CNV公司)通过对拷贝号更改信息应用自动编码器特征选择方法，使用所选特征创建。在下一步中，相似矩阵（Sim_MUT公司)基于单核苷酸突变数据的隐马尔可夫模型和多项式混合模型，使用一种新的特征选择方法生成。两个相似矩阵（Sim_IC50型)和（Sim_组织)分别使用跨药物的细胞系的IC50值和每个细胞系的组织类型信息来实现。最后，为了构建每个细胞株和药物的潜在载体，并预测细胞株是否对药物敏感，应用了逻辑矩阵分解方法。为了将IC50值分配给灵敏度和电阻这两个标签，我们使用了先前研究中介绍的策略[5,7,9]，它使用单个药物IC50值的中位数作为分类模型的阈值。如果IC50值小于单个药物细胞系的中位数，则分配给敏感度或标签为1的类别的细胞系，否则分配给耐药性或标签为0的类别的电池系。在下一节中，我们首先描述了研究和数据预处理中使用的数据集，然后解释了上述每个步骤的细节。

图1

Auto-HMM-LMF示意图概述。两个基于自动编码器的特征选择模型用于基因表达和拷贝数变化信息。针对单核苷酸突变数据，提出了一种基于隐马尔可夫模型（HMM）和多项式混合模型（MMM）的特征选择方法。在下一步中，定义每对细胞系之间的五个相似性和每对药物之间的相似性。为了预测跨药物细胞系的药物反应，应用了logistic矩阵分解方法

数据集和预处理

在这项工作中，我们使用了由1001个癌细胞系和265种测试药物组成的GDSC数据集，以及分析了1457个癌细胞系及其对24种药物的基因组图谱的CCLE数据集。在这些数据集中，细胞系以基因表达谱、拷贝数改变和单核苷酸突变等基因组特征为特征。半最大抑制浓度（IC50）值用于跨药物细胞系的敏感性测量。我们重点关注了98种和24种药物，其中SDF格式（编码药物的化学结构）分别可从GDSC和CCLE的NCBI PubChem Repository获得。在这些数据集中，基因表达特征没有缺失值。然而，一些细胞系的响应值、单核苷酸突变特征和拷贝数改变特征缺失。在第一步中，删除了包含一半以上特征缺失值的细胞系。

k近邻的已知值对其余缺失值进行了插补。每对细胞系的欧氏距离c_我和c_j个基于他们的基因表达谱x_我和x_j个定义如下：

$${\text{Dis}}_{{text{E}}\左（{{\text}}_}{{text}i}}，{text{c}}__{text{j}}}\右）\，=\，\左|{left|{text{x}}_[{text{i}}}{-}{text{x}}{{{{x}{}}}|_{{2}}^{2}$$

(1)

然后使用细胞系c的k个最近细胞系的平均特征值来插补药物d的缺失药物反应值（IC50），如下所示：

$${\text{IC5}}0\左（{{text{c}}，{text{d}}\right）=\mathop\sum\limits_{{\text}i}}=1}}^{text{k}}\frac{Dis_{E}\left（{{\text{c}，}c{i}}\right）=1}}^{k}Dis_{E}\左（{{\text{c}}，{}c{i}}\右）$$

(2)

类似地，使用细胞系c的k最近细胞系的平均特征值来插补基因g的缺失拷贝数改变值（CNV），如下所示：

$${\text{CNV}}\left（{{text{c}}，{text{g}}\right）=\mathop\sum\limits_{{\text}i}}=1}}^{\text{k}}\frac{{Dis_{E}\left（{{\text{c}{，{}c{i}}\right）=1}^{k}Dis_{E}\左（{{\text{c}}，{}c{i}}\右）$$

(3)

单核苷酸突变特征的值是二值的，即1表示突变，0表示野生型。细胞系c的k最近细胞系的平均特征值被认为是对基因g缺失的MUT（单核苷酸突变）值的插补，如下所示：

$${\text{MUT}}\left（{{text{c}}，{text{g}}\right）=\left\{{begin{array}{*{20} 我}1\h将&{if\，\，\left（{\mathop\sum\limits_{{\text{i}}=1}}^{{\text{k}}{\text}MUT}}\left c{{text{i}}，{\text{g}}\right）}\rift）}\hfill\\0\hfil&{otherwill\\end{array}}\right$$

(4)

最后，GDSC数据集考虑555个细胞株和98种药物，CCLE数据集考虑363个细胞系和24种药物。

使用自动编码器选择特征

特征选择方法旨在通过识别数据集中信息性和非冗余特征的子集来降低数据维度。自动编码器是一种用于无监督学习的非递归神经网络，它将数据集从初始特征空间减少到更重要的特征空间。它有一个输入层、一个输出层和一个或多个隐藏层。输出层中的节点数（神经元）与输入层中的相同。自动编码器通过假设输出层向量为输入层向量来学习权重向量。为了构建用于基因表达谱和拷贝数改变信息特征选择的自动编码网络，Xu等人引入了该策略[7]使用。为了筛选出基因表达特征和拷贝数变化数据，考虑了两个具有单个隐层（含100个神经元）和双曲正切函数的自动编码网络。在选择特征子集后，使用Boruta算法将另一组重要特征确定为逻辑矩阵分解模型的两类输入[6]. 为了确定基本基因，由自动编码网络选择的一组特征，以及对应于导入Boruta算法的细胞系和药物的敏感性和耐药性标签。Boruta算法是一种基于随机森林分类构建的包装器，它通过统计测试迭代删除不太重要的特征。该算法添加了使用自动编码器获得的所有特征的副本，并对复制的特征的值进行混洗以构建阴影特征，并试图找到基本特征。在扩展信息系统上运行随机森林分类器，Z值计算所有属性的重要性。Boruta算法重复查找过程（查找属性中的最大Z分数），直到为所有属性分配了重要性[6].

第一个单隐藏层自动编码器和Boruta算法分别用于GDSC和CCLE两个数据集的11、712和19、389个基因的基因表达谱。GDSC和CCLE选择的必需基因数量分别为798和1189。此外，对两个GDSC和CCLE数据集的24、959和24960个基因，应用类似的自动编码器和Boruta算法进行拷贝数改变。分别为GDSC和CCLE数据集选择了67和127个特征。

基于隐马尔可夫模型和多项式混合模型的特征选择

了解突变过程的活性对于癌症治疗和个性化治疗至关重要。由于突变过程在癌症基因组中留下了其活性的特征，因此从患者的单碱基替换模式中表征其活性突变过程的特征非常重要。在本研究中，我们使用了Wojtowicz等人提出的策略，将已知特征分配给相应的个体突变，以选择癌症类型中的基本突变基因[10]. 在这项工作中，我们只考虑了癌症体细胞突变目录（COSMIC）中经验证的突变特征[11]，我们重点研究了以前在癌症类型中被确定为活性的特征[12]. 表1显示了与GDSC和CCLE数据集中的癌细胞系相对应的14种癌症类型的活性特征（仅考虑GDSC和CCLE数据集中具有至少15种细胞系的癌症类型）。

表1 GDSC和CCLE数据集中与细胞株相对应的14种癌症类型的活性特征

因为有六类碱基替换（C:G > A： T、C:G > T： A、C:G > G： C、A:T > C： G、A:T > T： A、A:T > G： C）和四种可能的5′，我们将癌症基因组中的突变分为96类，包括碱基替换和四种潜在的3′碱基[13,14]. 我们从国际癌症基因组联盟数据门户下载了癌症类型的单碱基替换[15]. 我们分析了几种癌症类型患者的单碱基替代，这些患者（患者组1）对应于每种癌症类型的数量如表所示2对于每种癌症类型，应用以下隐马尔可夫模型和多项式混合模型，并确定所考虑的癌症类型的重要基因。在这个模型中，每种癌症类型的状态数是根据表中所示的相应特征的数量来确定的1例如，BRCA癌症的状态数（t）为12。

表2学习HMM和MMM模型的患者数量（患者组1）以及其基因表达信息可在国际癌症基因组联盟数据门户网站（患者组2）中获得的14种癌症类型的患者数量

详细的逐步特征选择程序如下所述：

识别紧密突变和孤立突变

我们使用2000 bp的距离阈值将突变分为两类，即近突变和孤立突变（孤立突变与任何其他突变都很遥远）。我们将每个突变序列的第一个突变设置为关闭。对于其他突变，如果与前一个突变的对应距离大于2000 bp，则该突变标记为孤立突变，否则接近。因此，从患者的一个突变序列中，我们可以得到几个子序列，一些对应于闭合突变，而另一些对应于孤立突变。例如，与BRCA癌症患者的孤立突变相对应的两个子序列和与BRCA癌患者的闭合突变对应的三个子序列如下：

$$\开始{aligned}和\下大括号{T>G，T>C，G>A，C>G，T>C，}_{close}\下大括弧{G>T，T>C$$

(5)

模拟闭合突变

由于与闭合突变相对应的子序列彼此接近，因此可以假设它们之间存在依赖性。因此，使用隐马尔可夫模型（HMM）对这些子序列进行建模。

具有t（突变特征数）隐藏状态的HMM M表示为

• Σ = {抄送₁，…，c_秒}是所有序列的字母集。

• 问 = {问₁，…，q_t吨}是一组状态，每个状态都能够发出字母∑的符号。

• π_我,∀我 = 1，…，t是从第i个状态开始的概率。

• A类 = 【a】_{i、 j个}]_{i、 j=1，…。，t吨}哪一个_{i、 j个}是q的转移概率_我到q_j个.

• E类 = [英]_{i、 j个}]_{i=1…，t，j=1，…，s}其中e_{i、 j个}是状态q的概率_我发射c_j个.

该模型假设每个代表突变类别的观察值都是由t状态之一发出的。生成观察序列的状态序列是未知的，每个状态都取决于之前的状态。为了学习模型的参数π、A、E，在第一步中都获得了作为HMM参数估计训练集的闭子序列。在本研究中，AntMarkov算法（受蚁群优化启发的隐马尔可夫模型参数估计算法）[16]用于估计HMM参数。

孤立突变建模

由于孤立的子序列远离任何其他突变，因此它们之间的依赖性假设动机较小。因此，使用多项式混合模型（MMM）对孤立突变进行建模。MMM由t突变特征边缘概率向量g和t × s发射矩阵E，（s = 96). 第一步中获得的所有孤立子序列都被视为MMM参数估计的训练集。向量g和发射矩阵E是根据样本中观察到的孤立突变次数来估算的（实验分布）。我们考虑一个大小为96的向量T，其中T[i]是样本中观察到第i个突变类别的孤立突变的总次数。通过将向量T应用于初始化向量g和发射矩阵E，得到了它们的估计值。

计算突变序列发生率

在训练了上述两个模型后，序列发生的概率O₁、…、O_T型，它被分解为封闭和孤立的子序列，{C₁，我₁，C₂，我₂，…，C_k1型，我_{k2（平方公里）}}，公式如下：

$${\text{P}}=\左}2}}P{\text{（i}}_{text{j}}}{\text}|MMM}}）}\right）{}$$

(6)

维特比算法[17]应用于查找生成闭合子序列的最可能状态序列的路径。为了确定与孤立突变相对应的最可能路径，使用了g矢量和发射矩阵E的估计值。对于每个单独的突变类别（O_t吨)，概率值最大的状态（Q_t吨)由以下公式得出：

$${\text{max}}_{j=1\ldots.t}\左$$

(7)

最后，我们附加这两个最可能的状态路径，以构建与患者对应的最终路径。然后，将观察到的状态数（特征）计算为每个样本的特征频率或每个路径的特征活动。例如，BRCA癌症患者有一个大小为12的向量，对应于签名数，该向量的元素是根据在最终路径中观察到每个状态的次数来计算的。

识别重要基因

为了考虑每个患者的特征活动和基因表达谱之间的关系，我们从国际癌症基因组联盟数据门户下载了患者的基因表达文件[15]. 可获得基因表达信息的患者数量（患者组2）如表所示2此外，由于两个GDSC和CCLE数据集的54和1667个基因的单核苷酸突变信息是可获得的，我们分析了这些基因的表达，以计算Spearman相关系数。因此，计算了1721个基因的表达与样本间特征活性之间的Spearman相关系数。通过这种方法，我们在每种癌症类型的近距离和孤立区域中确定了具有高Spearman相关系数（大于0.2）的基本基因。附加文件中说明了具有高Spearman相关系数的14种癌症类型中某些基因与特征活性之间的相关性结果1：表S1。我们认为这些基因是GDSC和CCLE数据集中单核苷酸突变数据的基本特征。最后，分别基于GDSC和CCLE中的隐马尔可夫模型和多项式混合模型，采用上述策略选择22个和72个基因。这些基因的列表如附加文件所示1：表S1。

相似性定义

由于相似细胞系和相似药物可能具有相似的药物反应，因此细胞系和药物之间的相似性可以改进药物反应预测[5,18].

在logistic矩阵分解模型中，识别最近邻需要相似度矩阵。基因表达谱、拷贝数改变、单核苷酸突变和组织类型信息用于细胞系相似性，药物化学结构用于药物相似性。因此，每对细胞系之间的五个相似性和每对药物之间的相似性定义如下：

细胞系相似性

• （模拟_经验)是基于基因表达谱的选定特征的相似性，其中，对于GDSC和CCLE两个数据集，通过自动编码器识别的基因表达谱基本基因的数量分别为798和1189。模拟_经验定义为每对n个细胞系的基因表达载体之间的皮尔逊相关性，排列成n个 × n矩阵。

• （模拟_CNV公司)是基于拷贝数改变数据的选定特征的相似度，其中67个和127个有用特征分别由GDSC和CCLE中的自动编码器选择。模拟_CNV公司矩阵定义为n × n矩阵通过Pearson相关性分析每对细胞系的拷贝数变化矢量。

• （模拟_MUT公司)是基于隐马尔可夫模型和多项式混合模型选择的单核苷酸突变信息特征的相似度。通过该策略分别从GDSC和CCLE数据集中鉴定出22个和67个必需基因。然后，将Jaccard相似性应用于对应于n个细胞系的每对单核苷酸突变载体_MUT公司构造为n × n矩阵。

• （模拟_IC50型)是基于细胞系的IC50值的细胞系之间的相似性。刘提出的细胞系间相似性的定义是基于细胞系响应IC50值之间的相关性[19]. 模拟_IC50型定义为n个细胞系之间的皮尔逊相关性 × n矩阵。

• （模拟_组织)是基于组织类型的细胞系之间的相似性。全套样本由GDSC和CCLE数据集组成，癌细胞系起源于大约14个组织部位。模拟_组织是一个n × n二值矩阵，如果两个单元格线c，则对应于第i行和第j列的条目为1_我和c_j个具有相同的组织类型，否则为零。模拟人生_组织与GDSC和CCLE细胞系对应的矩阵在附加文件中表示2：表S2和附加文件三：表S3。

由于上述每对相似矩阵之间的相关系数很低，因此矩阵之间不存在共线性，可以线性组合。我们构建了一个集成的矩阵相似性Sim_氯 = [SC_ij公司]_n×n，使用Sim的组合_经验，模拟_CNV公司，模拟_MUT公司，模拟_IC50型和Sim_组织通过以下公式：

$$\frac{{\uplambda{\text{Sim}}_{{\text}EXP}}}}+\upgamma{\text[Sim}}{{{CNV}}}+\upphi{\text]Sim}}{\text{TISSUE}}}}{{\uplambda+\upgamma+\upphi+\uppsi+\uprho}}$$

(8)

其中γ、λ、、ψ和ρ是控制每个矩阵重要性并在模型中调整的参数。我们定义了集合N_k个（c）_我)表示与c最相似的k个细胞系_我（c除外_我)使用（Sim_氯)矩阵。我们构造了邻接矩阵A = 【a】_ij公司]_n×n表示细胞系邻域信息如下：

$${\text{a}}_{{text{ij}}}=\left\{{\begin{array}{*{20} 我}{SC_{{\text{ij}}}}\hfill&{c_{{text{N}}}\在{\text}{{N}}_{text{k}}\left（{c_}{i}}\right）}\hfirl\\0\hfil&{otherwill\\end{array}}\ right中$$

(9)

药物相似性

药物之间的相似性是基于化学亚结构（Sim_药物). 对于每种药物，考虑881大小的零-一向量，其中881是药物已知化学亚结构的数量。在这个向量中，1表示药物的亚结构，否则为0。模拟_药物 = [标准偏差_ij公司]_米×米构造为m × m矩阵与m药物对应的每个化学子结构向量之间的Jaccard相似性。对于药物d_我，集合N_k个（d）_我)表示与d最相似的k种药物_我（d除外_我)使用Sim_药物矩阵。邻接矩阵，B = 【b】_ij公司]_米×米，描述药物邻域信息如下：

$${\text{b}}_{{text{ij}}}=\left\{{\begin{array}{*{20} 我}{SD_{{\text{ij}}}}\hfill&{d_{{text{N}}}\在{\text}{{N}}_{text{k}}\左（{d_}{i}}\右）}\\hfill\\0\hfill&否则}\\end{array}}\right$$

(10)

Logistic矩阵分解

为了利用选定的特征从GDSC和CCLE数据集预测癌细胞株的药物反应，引入了基于logistic矩阵分解方法的DSPLMF方法[5]基于以下目标函数应用：

$$开始{对齐}和{\text{min}}_{{{text{U}}，{text{V}}、\upbeta^{c}，\upbeta ^{d}}\mathop\sum\limits_{{\text}}=1}}^{n}\matshop\sum\ limits_{{{\text{i}}=1{}}{}q{{\text{ij}}}\右）\log\左}}^{c} +\upbeta_{{text{j}}}^{d}}\右）}\右}}^{d}}\右）+\frac{1}{2}tr[{text{u}}^}（\uplambda_{{text{c}}i+\upalpha{\text{H}}^c}）{text{u}}]+{}\frac}{1}}{2{tr[{\text}v}}^{T} （\uplambda_{{text{d}}}I+\uppeta{\text{H}}^{d}）{text{V}}]\\end{aligned}$$

(11)

其中u_我和v_j个是与细胞系c相对应的L大小的潜在载体_我和药物d_j个所有细胞系和所有药物的潜在载体分别用U和V表示。阳性值\（{\upbeta}_{{\text{i}}}^{c}\）和\（{\upbeta}_{{\text{j}}}^{d}\）是符合细胞系c的偏差参数_我和药物d_j个和\（{\upbeta}_{}^{c}{\text{和}}{\uppeta}_}^{d}\）分别是细胞系和药物的偏置载体[20]. 两个参数λ_c（c） = \（\frac｛1｝｛｛upsigma｝_｛\text｛c｝｝｝^｛2｝｝｝\）, λ_d日 = \（\frac{1}{{{\upsigma}_{{\text{d}}}^{2}}\），其中\（{\upsigma}_{{\text{c}}}^{2}\）和\（{\upsigma}_{{\text{d}}}^{2}\）是控制细胞系和药物的先验分布方差的参数。参数α和β决定了DSPLMF方法中细胞系相似性和药物相似性的有效性。（r≥1）是控制观察到的相互作用的重要性水平的参数。由于敏感性和耐药性等级在药物反应预测问题中具有相同的重要性，我们将r设为1。此外，H^c（c） = （E）^c（c） + \（\widetilde{{\text{E}}^{c}}}\）) − （A） + A类^T型)，E^c（c）和\（\widetilde{{\text{E}}^{c}}}\）是两个对角矩阵\（{\text{E}}_{{text{ii}}}^{c}=\sum\nolimits_{{\text}}=1}}^}{{\text{n}}}{（{\text}a}{{text}}}}）}\）和\（\widetilde{{\text{E}}_{{text{jj}}}^{c}}=\sum\nolimits_{{{\text{i}}=1}}^{{text}}{（{text{a}}_}{text{ij}}）}\），小时^d日 = （E）^d日 + \（\widetilde{{{text{E}}^{d}}}\）) − （B） + B类^T型)作为对角线元素（n是细胞系的数量）。E类^d日和\（\widetilde{{{text{E}}^{d}}}\）是两个对角矩阵\（{\text{E}}_{{text{ii}}}^{d}=\sum\nolimits_{{\text}}=1}}^}{{\text{m}}}{（{\text}b}}{{text}ij}}}）}\）和\（\widetilde{{\text{E}}_{{text{jj}}}^{d}}=\sum\nolimits_{{{\text{i}}=1}}^{{text}m}}{（{text{b}}_}{text{ij}}）}\），作为对角元素（m是药物的数量）。在训练所提出的模型后，使用以下公式确定细胞株和药物的潜在载体11然后，为了预测给定新细胞系在所有药物中的IC50值，选择新细胞系的k个最近邻，并根据其相邻潜伏载体的平均值估计该新细胞株的潜伏载体。由于（Sim）的元素_IC50型)矩阵未知，（Sim_氯)矩阵不能用于查找新细胞系的k近邻。我们使用DSPLMF方法中引入的策略进行估算（Sim_IC50型)DSPLMF方法设计了一个用于估计的决策树分类器模型（Sim_IC50型)利用新细胞系的基因表达谱、拷贝数改变和单核苷酸突变信息构建矩阵[5]. 然后，通过类似的方法，我们估计了与新细胞系相对应的潜在载体，以预测新细胞系对等式所示药物敏感的概率12对于一组细胞系和药物，细胞系c的概率_我对药物d敏感_j个可以建模为逻辑函数，如下所示：

$${\text{p}}_{{{text{ij}}}=\frac{{\exp\左{1+{\text{exp}}\左$$

(12)

最后，对概率应用阈值，为每个新的细胞系-药物对分配敏感或耐药类别。

结果

Auto-HMM-LMF预测性能评估

在药物反应预测问题中，使用特征选择方法是降低特征维数的常用方法之一。在以前的一些预测方法中，例如AutoBorutaRF，自动编码器方法用于选择基因组信息的重要特征。最强大的特征选择方法之一是Neumann等人提出的EFS方法。EFS方法集成了八种不同的特征选择法，并将所有单个输出归一化为一个公共尺度，即0到1的区间[21,22]. 首先，为了评估Auto-HMM-LMF模型中特征选择策略的效率，我们使用EFS方法选择基因表达谱、拷贝数变化和单核苷酸突变数据中的重要特征。在该方法中，为每组数据选择的特征数量等于Auto-HMM-LMF方法选择的特征的数量。然后，我们将这些特征与Auto-HMM-LMF方法中由Autocoder和HMM-MMM选择的特征交替使用，并将所得结果与其他方法进行比较。该方法应用于两个CCLE和GDSC数据集，我们在表中用EFS-LMF的名称表示该方法的结果三和4在本研究中，十倍交叉验证重复30次，其平均值用作评估AutoHMM-LMF方法预测性能的标准。

表3基于GDSC数据集七个指标的不同算法结果的性能比较

表4基于CCLE数据集七个指标的不同算法结果的性能比较

我们将Auto-HMM-LMF方法与六种分类模型（DSPLMF、EFS-LMF、CaDRReS、SRMF、AutoBorutaRF和SVM-RFE）进行了比较，以确定不同的指标。DSPLMF和AutoBorutaRF被设计为分类模型，但CaDRReS和SRMF方法预测IC50值作为输出。因此，为了将这些模型与Auto HMM LMF和EFS-LMF方法进行比较，我们将每种药物的预测IC50值的中位数作为分类阈值。如果与细胞系药物对相对应的预测IC50值小于此阈值，则将敏感类别分配给它；否则，它被标记为抗性等级。使用了七个指标：准确度、召回率、精确度、特异性、F1Score、马修斯相关系数（MCC）和接收器工作特性曲线下面积（AUC）；这些标准制定如下：

$${\text{Accuracy}}=\frac{{text{TP}}+{text{TN}}}{{text}TP}}+{text}FP}}+}text{TN{}}++{text$$

$${\text{Recall}}=\frac{{text{TP}}}{{{text}TP}}+{text{FN}}}}$$

$${\text{Precision}}=\frac{{text{TP}}}{{{text}TP}}+{text{FP}}}}$$

$${\text{Specificity}}=\frac{{text{TN}}{{{text}TN}+{text{FP}}}}$$

$${\text{F1Score}}=\frac{{2{text{TP}}}{{2}\text{TP}{+{text{FP}}+{text}FN}}}$$

$${\text{MCC}}={}\frac{{text{TP*TN}}-{text{FP*FN}}}{{sqrt{left（{{text}}+{text{FP}}\right）\ left text{FN}}+{\text{TN}}}\right）}}}$$

(13)

哪里

• TP（真阳性）：标记为敏感性并预测为敏感性的细胞株数量。

• TN（真阴性）：标记为耐药性并预测为耐药性的细胞株数量。

• FP（假阳性）：标记为耐药性并预测为敏感性的细胞株数量。

• FN（假阴性）：标记为敏感性并预测为耐药性的细胞株数量。

桌子三和4显示了在GDSC和CCLE数据集上进行的比较实验的结果（粗体数字表示最佳结果）。如表所示三，与最佳算法DSPLMF相比，Accuracy、Recall、Precision、F1Score、MCC和AUC标准的值分别增加了0.02、0.03、0.01、0.03、0.02和0.02。在特异性标准中，AutoBorutaRF方法的性能明显优于其他方法。关于其他标准，Auto-HMM-LMF方法的结果与AutoBorutaRF方法的结果非常显著。在表中4与其他算法的结果相比，Auto-HMM-LMF的所有标准值都有所增加，并且Auto-HMM-LMF在该数据集中显著优于最先进的方法。如表所示三和4与其他算法的结果相比，Auto-HMM-LMF的所有准则的值都有所增加。这些观察结果表明，HMM和MMM策略选择的突变数据特征以及基因表达和拷贝数变异数据的自动编码技术是非常有效和必要的。此外，EFS方法选择的特征在预测药物反应方面不如Auto-HMM-LMF方法选择的功能强大。

细胞系类型的组织特异性

为了证明Auto-HMM-LMF方法在不同组织类型中的性能，我们检查了在考虑特定细胞系组织类型时，我们提出的方法是否能够实现良好的性能。通过这种方法，我们考虑了GDSC数据集中的73个造血细胞和淋巴细胞系，并用7个标准评估了Auto-HMM-LMF方法。我们用这些细胞株训练Auto-HMM-LMF方法，并应用十倍交叉验证方法预测考虑的细胞株的药物反应。如表所示5这些结果证明，Auto-HMM-LMF方法也可以在造血细胞和淋巴细胞系上获得一致的或在某些标准下更高的性能。

表5基于七个标准的Auto-HMM-LMF方法对GDSC数据集73个造血细胞系的预测性能

预测和观测响应值之间的相关性

我们绘制了CCLE数据集中24种药物的观察到的药物反应和预测值的皮尔逊相关系数条形图。如图所示2以上相关系数（70%），从24中取17，均大于0.5。其中四种药物（PD-0325901、托普替康、AZD6244和伊立替康）的相关系数大于0.65。这些图显示了Auto HMM LMF方法在预测药物反应值方面的优异性能。上述4种药物的Auto-HMM-LMF模型的观察和预测药物反应的散点图如图所示三，CCLE数据集中20种其他药物的散点图如附加文件所示4.

图2

使用Auto-HMM-LMF模型预测和观察CCLE中24种药物的疗效值之间的皮尔逊相关系数

图3

使用Auto-HMM-LMF方法对4种药物的CCLE细胞系的观察到的和预测到的活性区域之间的相关性。绿色点表示观测值，紫色点表示Auto-HMM-LMF方法的预测值。一PD-0325901的观察和预测药物反应散点图（皮尔逊相关系数 = 0.81).b条拓扑替康的观察和预测药物反应散点图（皮尔逊相关系数 = 0.75).c（c）AZD6244的观察和预测药物反应散点图（皮尔逊相关系数 = 0.72).d日伊立替康的观察和预测药物反应散点图（皮尔逊相关系数 = 0.65)

药物重新定位申请

药物重新定位是为替代药理学目的选择已知药物的过程。对于这个问题，我们考虑了37种美国食品和药物管理局（FDA）批准的药物，这些药物没有在Choi等人的研究中的GDSC数据集中进行测试[23]. 在GDSC数据集上训练Auto-HMM-LMF模型，并预测20个头颈癌细胞株（HNSC）对37种药物中的20种抗癌药物的敏感性概率，如图所示4.

图4

20株HNSC细胞对20种抗癌药物敏感性的概率。已对名称为绿色的药物进行了治疗HNSC的研究。大于0.5的概率意味着相应的行药物可能是用于治疗相应列细胞系的新的重新定位药物。使用发散颜色方案时，忽略所有0.5左右的值，并添加红色的高概率与蓝色的低概率的视觉重要性

如图所示4，以下11种药物被确定为治疗HNSC的有效药物：

• 氮杂胞苷：氮杂胞嘧啶是一种被称为低甲基化剂的药物，一项研究报告称氮杂胞甙和顺铂对头颈癌的治疗有效[24].

• 雷戈拉非尼：雷戈拉非尼是一种口服多激酶抑制剂，靶向受体酪氨酸激酶（RTK）。Klinghammer等人建立了一个包含65例头颈部鳞状细胞癌的小组，他们证明雷戈拉非尼和埃弗罗莫司联合治疗对这些患者是有用的[25].

• 依维莫司：在癌症治疗中，依维莫斯被用作免疫抑制剂，以防止器官移植排斥反应。最近，一项研究表明，TP53突变的患者在头颈癌中明显受益于Everolimus[26].

• Fulvestrant：Fulvestlant是一种用于治疗激素受体（HR）阳性转移的药物。Grünow等人[27]结果表明，Fulvestrant抑制辐射诱导的ESR2表达，他们的发现证明了Fulvestlant联合放疗对HNSC患者的疗效。

• 干扰素Alf a–2B：干扰素α-2B是一种抗病毒或抗肿瘤药物，可有效治疗头颈癌[28].

• 卟啉：卟啉是一种光敏剂，用于癌症的放射治疗。体内研究表明，该药物可用于治疗HNSC患者[29].

• 咪喹莫特：咪喹莫特（INN）是一种处方药，作为免疫反应调节剂，用于治疗基底细胞癌。研究表明，外用咪喹莫特可能为患者提供合理且耐受性良好的姑息治疗选择[30].

• 地西他滨或5-氮杂-2′-脱氧胞苷是一种用于癌症治疗的核酸合成抑制剂。头颈部鳞状细胞癌对顺铂的耐药性降低了生存率。Viet等人[31]结果表明，地西他滨治疗可恢复HNSC细胞系对顺铂的敏感性，并显著降低诱导凋亡所需的顺铂剂量。

• 羟基脲：羟基脲是一种抗癌药物，用于治疗耐药、复发和转移的恶性黑色素瘤。一项研究表明，羟基脲是治疗头颈癌的单一活性剂。它已被临床用作放射增强药物和放射治疗[32].

• V andetanib:V andetan ib作为多种细胞受体的激酶抑制剂，是治疗癌症细胞系的抗癌药物。Sano等人[33]批准将V andetanib添加到与顺铂的联合治疗中，并且辐射可以在体外和体内HNSC模型中克服耐药性。

• 卡氟替佐米：卡氟替佐米是一种抗癌药物，作为一种选择性蛋白酶体抑制剂。通过上调促凋亡Bik、Carf ilzomib和ONX0912的表达，有效诱导HNSC细胞凋亡[34].

这些结果表明，Auto-HMM-LMF模型可以用于药物重新定位。此外，五种药物（依西美坦、伊布替尼、T沙利度胺、罗美司坦和氟氧美司酮）可能是治疗HNSC的新型药物。

超参数设置

由于GDSC数据集中的细胞株和药物数量高于CCLE数据集，因此我们调整了GDSC数据集上的超参数，并使用了两个数据集中超参数的获得值。这样，将十倍交叉验证程序应用于GDSC，并通过最大化AUC准则来确定超参数。

使用GDSC数据集学习的超参数如表所示6应用于公式12为了确定每个新细胞系的类别标签，从{0.1，…，1}中选择，该参数设置为0.4。潜在空间维度L从{1，…，min（n，m）}中选择，对于GDSC数据集，参数L设置为95，对于CCLE数据集，设置为23（其中n和m分别是细胞系数量和药物数量）。

表6基于GDSC数据集的Auto HMM LMF方法的学习超参数

讨论

本文提出了基于特征选择方法和logistic矩阵分解策略的Auto-HMM-LMF方法来预测药物反应。所提出的预测模型显示出比现有计算模型更高的预测效率。此外，我们还证明了Auto-HMM-LMF模型在药物重新定位中可能有用。因此，我们确定了五种用于HNSC治疗的药物（依西美坦、伊布替尼、T沙利度胺、罗美司坦和氟西美斯特龙）。为了说明隐马尔可夫模型（HMM）和多项式混合模型（MMM）对Auto-HMMLMF方法中突变数据所选特征的生物学意义，我们进一步考虑了与乳腺癌相关的癌细胞株（BRCA）和两个重要过程，即MMR和HRD。本研究通过考虑基因表达文件与12个BRCA癌细胞株特征信号活性之间的Spearman相关系数，选择了30个重要基因。在这些基因中，PMS2、MLH1、MSH2和MSH6四个基因的基因表达与标记6、20和26的活性具有较高的Spearman相关系数。这些基因的表达与BRCA中三种突变特征活性之间的Spearman相关性结果如图所示5另一方面，最近的一项研究[35]研究表明，这三个6、20和26个特征与乳腺癌中MMR缺乏相关。缺陷DNA错配修复（MMR）在许多癌症类型中都会发生，PMS2、MLH1、MSH2和MSH6基因的突变是错配修复缺陷的最常见原因。上述基因被称为DNA错配修复（MMR）基因，这些基因参与修复DNA复制中的错误（在准备细胞分裂时复制DNA时发生的错误）[36].

图5

BRCA样本中PMS2、MLH1、MSH2和MSH6表达与特征6、20和26活性的Spearman相关性比较。对于这三个特征，分离和接近的突变计数与考虑基因的表达呈正相关

此外，通过该模型选择的其他10个基因，即BRCA1、BRCA2、ATM、CHEK2、MRE11A、NBN、FANCA、PALB2、RAD51C、RAD50的表达与突变特征活性3具有较高的Spearman相关系数。这些具有标志性3活性的基因表达之间的Spearman相关系数如图所示6同样，在最近的一项研究中[35]研究表明，同源重组缺陷（HRD）与乳腺癌患者的标志3有关。同源重组缺陷是指人类细胞不能修复双链断裂。一些导致HRD的基因改变包括基因的体细胞突变，例如选定的10个基因。通过Auto-HMM-LMF方法为BRCA癌细胞株选择的其他八个基因是APOBEC家族中的APOBEC3A、APOBEC3B、APOBEC 1、APOBECM 3C、APOBEC3D、APOBEC3F、APOBEG3G和APOBEC4。这些具有2和13个特征活性的基因表达之间的Spearman相关系数如图所示7.最近的一项研究[35]研究表明，APOBEC将胞嘧啶脱氨基为尿嘧啶被认为会引发特征码2和13的突变。因此，这些结果表明，为乳腺癌选择的基因在生物学上是必要的，并且Auto-HMM-LMF方法能够检测到单核苷酸突变数据的显著特征。与其他模型相比，除了提高药物反应预测的准确性之外，Auto-HMM-LMF算法的一个最重要的优点是，该算法的运行时间显著低于其他提及方法的运行时间。由于该方法基于特征的选择，其局限性之一是结果取决于所选特征。因此，可以通过使用更强大的特征选择方法来改进该方法的结果。该方法的以下局限性是，它是为了解决新细胞系的感冒问题而设计的，而一些拟议的方法也可以预测新药或新的细胞系药物对。

图6

BRCA样本中BRCA1、BRCA2、ATM、CHEK2、MRE11A、NBN、FANCA、PALB2、RAD51C、RAD50表达与签名3活性的Spearman相关性比较。对于这个特征，分离和接近的突变计数与考虑的基因的表达呈正相关

图7

BRCA样本中APOBEC3A、APOBEC3B、APOBEC1、APOBEC3C、APOBEC3D、APOBEC3F、APOBEC3G和APOBEC4表达与信号2和13活性的Spearman相关性比较。对于这个特征，分离和接近的突变计数与考虑的基因的表达呈正相关

结论

在这项研究中，我们开发了一种基于特征选择的方法Auto HMM LMF，用于预测GDSC和CCLE数据集中癌细胞系对药物的反应。为了对基因表达和拷贝数变化数据进行特征选择，设计了两个自动编码网络。对于单核苷酸突变信息的特征选择，应用了基于隐马尔可夫模型（HMM）和多项式混合模型（MMM）的新方法。与最先进的预测方法相比，Auto-HMM-LMF显示出更好的整体预测性能。此外，通过与最强大的特征选择方法之一EFS方法进行比较，我们表明基于本文介绍的所选特征的预测模型在药物反应预测方面的性能要好得多。此外，我们认为该模型可以用于许多治疗研究领域，如药物重新定位和个性化医疗。最后，我们发现大量证据表明，Auto-HMM-LMF选择的特征和预测的响应与许多先前的研究具有显著的一致性。