CN114913916A - 预测新冠病毒适应药物的药物重定位方法 - Google Patents

Abstract

本发明实施例公开了一种预测新冠病毒适应药物的药物重定位方法，包括：步骤1：采用WKNKN算法处理人类药物‑病毒关联数据集；步骤2：分别计算药物和病毒的杰卡德相似性和高斯核相似性；步骤3：利用中心核对称多核学习方式进行整合，从而分别获得药物和病毒的综合相似性矩阵；步骤4：结合处理完的人类药物‑病毒关联数据集信息，通过压缩感知算法补全原人类药物‑病毒关联数据集中缺失的关联，得到预测概率矩阵；步骤5：作为药物‑病毒的预测结果矩阵。本发明结合中心核对称的多核学习模型的压缩感知方法去预测药物与疾病的潜在关联概率，并进行实验验证该方法作为寻找能与新冠病毒作用的药物研发辅助工具的可行性。

Description

预测新冠病毒适应药物的药物重定位方法

技术领域

本发明涉及药物-疾病关联预测领域，尤其涉及一种预测新冠病毒适应药物的药物重定位方法。

背景技术

2019年被发现的COVID-19病毒，因外形类同于以往发现的冠状病毒而被命名为“新型冠状病毒”，是继人类发现并研究的六种冠状病毒外的(229E(HCoV-229E)、HCoV-OC43、HCoVNL63、HCoV-HKU1、严重急性呼吸综合征(SARS)相关冠状病毒(SARS-Cov)和中东呼吸综合征(MERS)相关冠状病毒(MERS-Cov))第7种冠状类病毒。纵观全球疫情情况，新冠病毒正在严重影响人类的正常生活，并且随时威胁着人类的生命安全。尽管现在已经有团队研发了针对性的新冠疫苗，但放眼全球，疫苗还是属于紧缺资源，目前对于新冠肺炎的治疗药物在紧急研发中。更严重的是，新冠病毒的遗传物质是单链RNA，这种遗传物质会使新冠病毒在转录复制阶段很容易发生变异，这种变异存在改变新冠病毒的体外特征的可能，随时威胁着现有新冠病毒疫苗的有效性。目前研究人员已经发现变异的新冠病毒毒株(如德尔塔毒株、奥密克戎毒株)，因此，加速研发治疗新冠的有效药物迫在眉睫。

药物开发通常包括三个阶段:发现、临床前和临床开发。每个阶段都需要花费大量的金钱和时间。药物重定位技术，是一种为现有药物寻找其他适应症状的技术(简称旧药新用)。得益于大数据时代的发展集成了很多生物临床试验结果，人们对于利用计算机算法模型作为辅助工具的药物重定位技术越来越感兴趣。与传统的药物临床试验相比，药物重定位技术能有效减少药物研发过程中的金钱和时间的花费。近年来不少学者已经提出了以计算机算法模型为辅助工具的药物重定位方法，并且均得到了有效的生物实验数据的验证。因此，采用药物重定位技术去筛选适应新冠病毒的治疗药物是可行的。

目前以计算机算法模型为辅助工具的药物重定位方法可分为三类：机器学习方法类、基于网络传播的方法和基于矩阵分解与补全的方法。这些方法是基于一种假设，即类似的药物与类似的疾病有关，反之亦然。

在基于机器学习方法中，预测潜在的药物-疾病关联可以被视为一个二元分类问题。药物-疾病关联被认为是样本，而药物和疾病先前的相似性被认为是特征。目前的机器学习方法类有基于有监督学习和半监督学习之分，它们的区别在于是否对训练集的样本标签存在依赖作用。有监督学习的方法严重依赖于数据集中的已知样本标签，但这些标签数据在实际情况下很难获得，从而限制了这些预测方法的应用范围。半监督学习方法不再需要负样本，但它们的性能并不稳定，这在一定程度上影响了半监督学习方法的应用。机器学习方法在药物重定位领域显示了巨大的潜力但当前方法依然存在着不足之处，限制了该类方法的实际应用。

药物和适应症状的关系结构利用网络结构可以很容易的表达出来，目前很多基于网络的药物重定位方法都是利用药物和适应症状的关联矩阵建立网络结构，基于各自的相似性信息进行权值的初始赋值，并使用游走或扩散的方式将关联资源从药物端传播到适应症状端，从而推断网络中的缺失边，在计算效率上有很大的优势。

对于一个部分仅知的信息矩阵，矩阵补全是用一个低秩的模拟矩阵近似该信息矩阵，以模拟矩阵对应的元素作为目标矩阵中缺失的信息的估计值。而矩阵分解通常是指用两个低秩的子矩阵的内积来近似目标矩阵，那么子矩阵的内积矩阵的元素就可以用于估计目标矩阵中对应不可见位置的元素值，而子矩阵可以看作是目标矩阵的分解。显然，通过以上描述可知，矩阵分解可以用于矩阵补全任务，但是，矩阵补全也并不总是采用分解的方法。

在上述方法中应用于药物重定位的关联数据集仅有已验证的“正”数据，压缩感知作为矩阵补全的一类方法，因其使用“子矩阵模拟”技术来预测药物-疾病的潜在相互作用而不将所有缺失的数据作为负数据处理，在应用于药物重定位的方面与机器学习的监督类方法相比存在先天优势。另一方面，由于现有的数据集的关联信息受生物实验的因素影响而变得稀疏，基于网络或图的这类方法大多依靠关联数据作为二分图的资源传播路径往往会出现冷启动问题而影响模型的预测性能，部分采用游走方式的方法甚至会因游走方式偏好权值大的路径而陷入无意义循环，增大模型的时间复杂性。相比之下，以现有的矩阵信息出发，压缩感知利用子矩阵对原矩阵信息进行采集，生成一个低秩的模拟矩阵去补全原矩阵的缺失部分，与基于网络或图的方法相比在一定程度上避免了模型的冗余。

药物重新定位旨在为现有药物确定新的适应症，这技术将大大降低了药物开发的成本和时间。但是，到目前为止的所有应用于药物重定位的方法都直接使用未经处理的原始的药物-疾病关联矩阵的信息，由于在实际情况下这些相似关联信息只有小部分被证明并记录下来，还有很多相似关联信息亟待证明，导致用于预测方法的原始关联矩阵数据稀疏性大。另外，对于药物之间和疾病之间的相似性特征的提取，目前还缺乏充分挖掘方法。

发明内容

本发明实施例所要解决的技术问题在于，提供一种预测新冠病毒适应药物的药物重定位方法，以使能够作为寻找能与新冠病毒作用的药物研发辅助工具。

为了解决上述技术问题，本发明实施例提出了一种预测新冠病毒适应药物的药物重定位方法，包括：

步骤1：采用WKNKN算法处理人类药物-病毒关联数据集HDVD，所述数据集包括药物-疾病关联邻接矩阵Y、药物化学结构相似度矩阵、疾病语义相似度矩阵；挖掘潜在相互作用的似然分数，降低药物-疾病关联矩阵Y的稀疏性，得到处理后的人类药物-病毒关联数据集HDVD；

步骤2：根据处理后的人类药物-病毒关联数据集HDVD分别计算药物和病毒的杰卡德相似性和高斯核相似性，得到药物和疾病的杰卡德相似矩阵DS_jac和dS_jac，得到疾病的高斯核相似矩阵dS_Gaus以及药物的高斯核相似矩阵DS_Gaus；

步骤3：分别结合药物的化学相似性和病毒的序列相似性利用中心核对称多核学习方式进行整合，从而分别获得药物和病毒的综合相似性矩阵；

步骤4：结合处理后的人类药物-病毒关联数据集HDVD信息，通过压缩感知算法补全原人类药物-病毒关联数据集中缺失的关联，得到预测概率矩阵；

步骤5：将得到预测概率矩阵作为药物-病毒的预测结果矩阵。

进一步地，步骤1包括以下子步骤：

子步骤1、对于每个药物D_i，采用与其最接近的K种已知药物的化学相似度矩阵DS_chem及其对应的关联作用来估计D_i的相互作用似然分数，推导公式如下：

其中D₁到

表示药物D_i的K个近邻药物，按降序排列；

是权重系数，其中η是衰减项并满足η≤1，

是正则化项；

子步骤2、相似地，对于每个疾病d_j，采用与其最接近的K种已知疾病的语义相似性及其对应的交互作用来估计d_j的交互似然分数，推导公式如下：

其中d₁ to

表示d_j的K个近邻药物，按降序排列；

是权重系数，其中η是衰减项并满足η≤1，

是正则化项；

子步骤3、采用下式计算得到处理后的药物-疾病关联邻接矩阵Y_F：

进一步地，步骤2中，采用下式计算药物的杰卡德相似矩阵DS_jac：

其中Γ(x)是与药物D_x相关联的疾病集合，Γ(y)是与药物D_y相关联的疾病集合；

同理，按照上述方式计算出疾病的杰卡德相似矩阵dS_jac。

进一步地，步骤2中，将疾病d_i和d_j之间的高斯核相似性表示为dS_Gaus，药物D_i和D_j之间的高斯核相似性表示为DS_Gaus，计算方法如下：

dS_Gaus(d_i，d_j)＝exp(-α_d||Q(d_i)-Q(d_j)||²)；

DS_Gaus(D_i，D_j)＝exp(-α_D||Q(D_i)-Q(D_j)||²)；

在药物-疾病关联邻接矩阵Y_F中，第i行表示药物D_i是否与每种疾病有关，第j列表示疾病d_j是否与每种药物有关；向量Q(D_i)和Q(d_j)分别将第i行向量和第j列向量表示为高斯核的特征向量；核带宽α_d和α_D定义为：

其中，初始核带宽系数α′_d和α′_D均设置为1。

进一步地，步骤3中，使用CKA-MKL算法分别在两个空间组合三个相似核矩阵，其中，特征核K^*采用下式计算：

其中k为核数，K_i分别在药物相似核集DS_chem、DS_jac、DS_Gaus和疾病相似核集dS_sem、dS_jac、dS_Gaus中取值，ω_i为核K_i的权值，N是节点数；

内核对齐的值定义如下：

其中P，Q∈R^N×N，N∈{n_D，n_d}，<P，Q>_F＝Trace(P^TQ)是Frobenius内积，

是Frobenius范数；

根据上式定义经过中心核对齐后的目标函数如下：

式中

表示一个中心核矩阵，I_N∈R^N×N是一个单位矩阵，I_N是单位向量，再将目标函数转变为：

其中，a∈R^k×1和M∈R^k×k分别由下式计算而得：

将最终目标函数表示为：

将上式最小化计算疾病和药物相似矩阵组合的最优权值ω_D和ω_d，并分别将它们与药物DS_chem、DS_jac、

和疾病dS_sem、dS_jac、

相似核按下式组合：

其中

和

和

分别是疾病核和药物核的最优权值；

分别计算药物和疾病的最佳组合相似矩阵。

进一步地，步骤4中，压缩感知算法的计算过程为：

构建初始子矩阵

和

根据下式计算初始概率矩阵：

其中p_i，j是药物D_i和疾病d_j的相互作用概率，f_i是药物子矩阵F的第i行，g_j是疾病子矩阵G的第j行，w_i，j是药物D_i和疾病d_j的初始权值；

构建损失函数式：

引入两个额外的可调参数λ_M和λ_N，转化得到损失函数为：

其中，

根据下式计算子矩阵F和矩阵G的偏导数：

其中

代表Hadamard乘积；

根据下式分别更新子矩阵F和矩阵G：

其中，k为迭代过程的学习率；上标n是当前迭代次数；子矩阵F和G更新的结束条件为max(ΔF，ΔG)＜10^-5。

进一步地，步骤4中，根据下式，用矩阵Fⁿ⁺¹和Gⁿ⁺¹计算概率矩阵Pⁿ⁺¹：

本发明的有益效果为：本发明结合中心核对称的多核学习模型的压缩感知方法去预测药物与疾病的潜在关联概率，并进行实验验证该方法作为寻找能与新冠病毒作用的药物研发辅助工具的可行性；本发明用药物-疾病关联数据集Fdataset和Cdataset来验证，最后迁移应用到人类药物-病毒关联数据集HDVD上。

附图说明

图1是本发明实施例的预测新冠病毒适应药物的药物重定位方法的流程示意图。

图2(a)和(b)分别是本发明实施例的方法与其他现有的方法在Fdataset数据集和Cdataset数据集的ROC图。

图3(a)和(b)是本发明实施例的方法与其他矩阵补全方法在新冠病毒数据集HDVD的ROC图和PR图。

图4是本发明实施例消融实验中的各模块组成的roc图，其中(a)是基于数据集Fdataset进行实验的roc图；(b)是基于数据集Cdataset进行实验的roc图。

具体实施方式

需要说明的是，在不冲突的情况下，本申请中的实施例及实施例中的特征可以相互结合，下面结合附图和具体实施例对本发明作进一步详细说明。

本发明先引入WKNKN算法处理药物-疾病的关联邻接矩阵，以降低其稀疏性。利用CKA-MKL处理多种相似性矩阵可以得到每个相似性矩阵(核)的权重，从而选择最佳的矩阵组合。本发明采用CKA-MKL算法将不同的相似信息矩阵进行融合，得到综合相似性矩阵。最后，利用压缩感知的矩阵补全算法计算药物-疾病关联预测评分。本发明的总体流程图如图1所示，图中右下方框内为压缩感知模型的迭代过程。

本发明的算法过程如下：

输入：邻接矩阵

药物化学相似矩阵

以及疾病语义相似矩阵

表示邻居个数的参数K和衰退参数η，学习率k，子矩阵F和G的Frobenius范数参数项系数λ_F和λ_G，子矩阵F和G的正则项参数λ_M和λ_N。

输出：关联预测矩阵Pⁿ⁺¹。

Step1使用WKNKN算法计算关联概率矩阵Y_F：

fori＝1→n_D do

根据式(1)计算Y_D；

end for

for j＝1→n_d do

根据式(2)计算Y_d；

end for

Step2根据式(4)分别计算药物和疾病的杰卡德相似矩阵DS_jac和dS_jac；

Step3根据式(5)计算疾病的高斯核相似矩阵dS_Gaus；

Step4同理，根据式(6)计算药物的高斯核相似矩阵DS_Gaus；

Step5使用CKA-MKL算法分别在两个空间组合三个相似核矩阵：

根据式(9)和式(10)定义目标函数式(11)和式(12)；

根据式(13)和式(14)将目标函数式(12)转化为式(15)；

通过使函数(15)最小化计算疾病和药物相似矩阵组合的最优权值ω_D和ω_d；

根据式(16)和式(17)分别计算药物和疾病的最佳组合相似矩阵；

Step6压缩传感模型的计算过程：

构建初始子矩阵

和

根据式(19)计算初始概率矩阵；

构建损失函数式(21)；

根据式(22)和式(23)将函数式(21)转化为式(24)；

Do

通过求解式(24)，从式(25)和式(26)计算矩阵F和矩阵G的偏导数；

根据式(27)和式(28)分别更新矩阵F和矩阵G；

Until max(ΔF，ΔG)＜10^-5；

Step7根据式(19)，用矩阵Fⁿ⁺¹和Gⁿ⁺¹计算概率矩阵Pⁿ⁺¹。

输出：Pⁿ⁺¹。

本发明将压缩感知作为寻找新冠病毒的潜在适应药物的方法。由于矩阵补全方法生成的模拟矩阵是经过低秩优化处理的，而低秩优化受原矩阵的秩影响，原则上模拟矩阵的秩越大，所能补全的信息越多。因此，引入权值K最近邻算法(Weight K Nearest KnownNeighbors algorithm，WKNKN)算法对原矩阵进行预处理，增大原矩阵的秩。此外，压缩感知允许本发明加入药物和疾病的相似性特征来辅助修正子矩阵。由于单一相似性信息对模型性能的稳定性影响较大，目前方法大多采用多种相似性的组合，但它们一般仅作求均值处理，缺乏系统的组合方式。为了更好的整合多种类型的相似性信息，利用中心核对称的多核学习模型(Centered Kernel Alignment-based Multiple Kernel Learning，CKA-MKL)对多种类型的药物和疾病相似性矩阵进行融合。经过CKA-MKL的处理，本发明可以得到每个相似矩阵(核)的权重，并选择最佳的矩阵组合，从而提高各种相似信息的互补性。

以下对本发明的各实施例原理及过程分别进行阐述：

1、验证数据集

在本发明中，采用Fdataset和Cdataset作为金标准数据集来测试本发明提出的模型的性能。

1.1、药物-疾病数据集

其中Fdataset是从不同的数据源中收集了593种药物、313种疾病和1,933种已验证的药物-疾病关联。而另一个数据集Cdataset包括409种疾病，663种药物和2532种疾病与药物的关联。两个数据集都各自包含三个信息矩阵：

(1)药物-疾病关联邻接矩阵

(n_D表示药物的数量，n_d表示疾病的数量)。若药物D(i)与疾病d(j)存在关联，则元素Y(D(i)，d(j))的值为1；否则，其值为0；

(2)药物化学结构相似度矩阵

该矩阵由从化学开发试剂盒导出的相应药物的化学结构信息构成，其中药物-药物对DS_chem(i，j)表示为药物D_i和D_j的二维化学指纹分数；

(3)疾病语义相似度矩阵

通过获取对应疾病对的相似度信息，经过文本挖掘计算出疾病之间的语义相似度得到的疾病语义相似度矩阵dS_sem。

1.2、药物-病毒数据集

HDVD是利用文本挖掘技术收集了大量经实验验证的药物-病毒相互作用词条构建的实验支持的人类药物-病毒关联数据库。本发明采用HDVD作为寻找新冠病毒的适应药物的数据集。HDVD可在(https：//github.com/luckymengmeng/HDVD)免费下载。HDVD收录了34种病毒、219种药物和455种证实的人类药物-病毒相互作用。与上述金标准数据集类同，HDVD包含三个信息矩阵：

(1)药物-病毒关联邻接矩阵

(mD表示药物的数量，mv表示疾病的数量)。若药物D(i)与病毒V(j)存在关联，则元素A(D_i，V_j)值为1；否则，将其设值为0；

(2)药物化学结构相似矩阵

SMILES(Simplified molecularInput Line Entry System)是一种以一维表示形式描述分子结构的信息格式。本发明采用SMILES格式从DrugBank数据库中下载相应药物的化学结构信息，通过Open Babel v2.3.1计算每种药物的分子访问系统(Molecular Access System，MACCS)指纹，并利用谷本系数测量两个分子之间的绝对相似性来构建药物化学结构相似矩阵DS_chem；

(3)病毒基因序列相似矩阵

MAFFT是一种基于相似性的多序列比对方法。本发明从国家生物技术信息中心(NCBI)下载了智人体内病毒的基因组核苷酸序列，并利用第7版MAFFT计算病毒之间的序列相似性来构建病毒基因序列相似矩阵VS_gen。

本发明各步骤的原理和流程细节如下：

2.1、WKNKN

压缩感知算法原理上是利用子矩阵对关联矩阵的信息进行采样，并用子矩阵还原的模拟矩阵来近似关联矩阵。目标矩阵的信息稀疏性会影响子矩阵中的元素对目标矩阵信息的描述性，从而影响模拟矩阵的恢复精度。Cdataset和Fdataset中的药物-疾病关联矩阵Y都是稀疏的(HDVD同理)，这是因为目前的研究进展中仍有很多药物-疾病关联尚未得到证实。因此，本发明通过引入WKNKN来挖掘潜在相互作用的似然分数，在一定程度上降低药物-疾病关联矩阵Y的稀疏性。即在邻接矩阵Y中，使用WKNKN算法来估算药物-疾病对之间的相互作用似然值。具体可分为三个步骤进行：

Step1、对于每个药物Di，采用与其最接近的K种已知药物的化学相似度矩阵DS_chem及其对应的关联作用来估计Di的相互作用似然分数。推导公式如下：

其中D₁到

表示药物D_i的K个近邻药物，按降序排列；

是权重系数，其中η是衰减项并满足η≤1，

是正则化项。

Step2、相似地，对于每个疾病d_j，采用与其最接近的K种已知疾病的语义相似性及其对应的交互作用来估计d_j的交互似然分数。推导公式如下：

其中d₁ to

表示d_j的K个近邻药物，按降序排列；

是权重系数，其中η是衰减项并满足η≤1，

是正则化项。

Step3、最后，如果Y_ij＝0，则将其替换为取Y_D和Y_d的平均值，如下所示：

2.2、药物和疾病的杰卡德相似性

Jaccard相似性是向量之间的公共邻域测度的标准化形式。以药物D_x和D_y为例，杰卡德相似性证明了从D_x和D_y节点对的疾病关联集合中随机选择是D_x和D_y节点对的共同作用疾病之一的概率。若D_x和D_y节点的共同疾病的数量越多，这个度量就越大，D_x和D_y就越相似。公式如下：

其中Γ(x)是与药物D_x相关联的疾病集合。

相似地，可以按照上述方式计算出疾病的杰卡德相似矩阵dS_jac。

2.3、药物和疾病的高斯核相似性

高斯核相似性是计算不同类型节点之间相似性的常用方法之一。通过径向基将数据进行高维投影，可计算不同节点向量之间的距离，从而获得节点之间的相似性权重，因此高斯核相似性也被称为径向基函数(RBF)核相似性。在邻接矩阵Y_F中，第i行表示药物D_i是否与每种疾病有关，第j列表示疾病d_j是否与每种药物有关。向量Q(D_i)和Q(d_j)分别将第i行向量和第j列向量表示为高斯核的特征向量。因此，本发明将疾病d_i和d_j之间的高斯核相似性表示为dS_Gaus，药物D_i和D_j之间的高斯核相似性表示为DS_Gaus。计算方法如下：

dS_Gaus(d_i，d_j)＝exp(-α_d||Q(d_i)-Q(d_j)||²) (5)

DS_Gaus(D_i，D_j)＝exp(-α_D||Q(D_i)-Q(D_j)||²) (6)

这里核带宽α_d和α_D被定义为：

其中初始核带宽系数α′_d和α′_D均设置为1。

2.4、基于中心核对齐的多核学习算法CKA-MKL

目前已有不少挖掘药物和疾病的相似性方法被提出，它们能从不同的角度去挖掘药物和疾病的相似特征。对于药物重定位方法而言，采用不同的相似性信息对模型预测性能的影响不同，因此目前的方法大多采用多种相似性信息的组合，但它们一般仅作求均值或互补等处理，缺乏系统的组合方式。为了合理地整合多种类型的相似性信息，本发明通过CKA-MKL处理多种相似性矩阵并得到每个相似性矩阵(核)的权重，从而选择最佳的相似矩阵组合，提高各种相似信息的互补性。具体地，由上述工作可以得到的药物和疾病的相似核集中各有三种相似核矩阵

然后，本发明使用基于中心核对齐的多核学习(CKA-MKL)方法分别在两个空间中组合三个核矩阵。最优内核计算如下：

其中k为核数，K_i分别在药物相似核集(DS_chem，DS_jac，DS_Gaus)和疾病相似核集(dS_sem，dS_jac，dS_Gaus)中取值，ω_i为核K_i的权值。N是节点数。

核对齐的得分可以通过计算两个核之间的余弦关联度来描述。核之间的关联度越大，核之间对齐度越高。因此内核对齐的值定义如下：

其中P，Q∈R^N×N，N∈{n_D，n_d}，<P，Q>_F＝Trace(P^TQ)是Frobenius内积，

是Frobenius范数。实际上，核对齐得分可以看作是两个核(特征核K^*和理想核矩阵Y_F

)之间的相关性。为了得到核的最优权值，应该让K^*，Y_F

之间的对齐得分最大。因此，经过中心核对齐后的目标函数如下：

式中

表示一个中心核矩阵。I_N∈R^N×N是一个单位矩阵。I_N是单位向量。式(11)也可以转变为：

其中，a∈R^k×1和M∈R^k×k分别由公式(13)和(14)计算而得：

最终目标函数也可以表示为：

为了获得每个相似核的重构权值，本发明采用标准二次规划来求解Eq.(15)。因此，分别得到药物和疾病的权重(ω_D，ω_d∈R^3×1)，并分别将它们与药物

和疾病

相似核按下式组合：

其中

和

和

分别是疾病核和药物核的最优权值。

2.5、压缩感知模型

压缩感知是一种矩阵补全类方法，其原理是通过子矩阵的内积得到的模拟矩阵近似目标矩阵，然后用模拟矩阵中的元素作为目标矩阵中不可观察部分的元素的估计。换句话说，给定一个相互作用信息存在缺失的矩阵

nD代表药物数量，n_d代表疾病的数量，本发明在WKNKN算法的基础上将Y_F看作是由两个低维的子矩阵

和

的内积

模拟所得，子矩阵分别对应于药物和疾病的潜在特征。通过将药物和疾病的关联信息映射到低维度的共同潜在的空间来近似于药物与疾病之间的相互作用概率。在补全过程中，通过减少模拟矩阵和目标矩阵的差距来驱使子矩阵的优化。为了避免模型出现过度优化的情况，在计算模拟矩阵和目标矩阵的差距时会计算子矩阵的复杂程度作为子矩阵优化的惩罚项。此外，在优化过程中，压缩感知模型允许引入药物和疾病的相似性信息来不断修正子矩阵模型。

根据上述描述，可以知道，给定药物D_i和疾病d_j，它们相互作用事件的概率是可以通过下式计算：

或者按矩阵方式表示：

其中p_i，j是药物D_i和疾病d_j的相互作用概率，f_i是药物子矩阵F的第i行，g_j是疾病子矩阵G的第j行，w_i，j是药物D_i和疾病d_j的初始权值。通过贝叶斯推理，可以推导出p(F，G|Y_F)的概率：

p(F，G|Y_F)∝p(Y_F|F，G)p(F)p(G)(20)

进而得到模型的损失函数：

为了提高方法预测的准确性，本发明在相似的药物可能与相似的疾病相互作用这一假设的基础上进一步扩展了损失函数。具体地，设

为药物相似度矩阵，其中每个条目DS_i，j表示药物D_i和D_j之间的相似度，设

为疾病相似度矩阵。相似的疾病与相似的药物结合是通过最小化药物之间的性质的距离来解释：

同样还有疾病之间的相似特性：

将上述正则化项(22)和(23)合并到(21)中，并引入两个额外的可调参数λ_M和λ_N，本发明的损失函数为：

有几种方法可以使损失函数(24)[46-48]最小化。与Liu等人[48]相似，压缩感知使用AdaGrad-一种迭代梯度下降方法，在迭代过程中，损失函数的偏导数可以写成：

其中

代表Hadamard乘积。子矩阵F和G按以下公式更新：

ΔF＝||Fⁿ⁺¹||_F-||Fⁿ||_F (29)

ΔG＝||Gⁿ⁺¹||_F-||Gⁿ||_F (30)

其中k为迭代过程的学习率；上标n是当前迭代次数；子矩阵F和G更新的结束条件为max(ΔF，ΔG)＜10^-5；||BF||_F是BF的Frobenius范式，||BG||_F同理，分别定义为：

本发明针对新冠病毒的药物研发构建一个性能优异、结果可靠的药物重定位辅助计算方法，因此本发明采用了在药物重定位研究邻域中被广泛应用的Fdataset和Cdataset作为金标准数据集，并将本发明提出的方法与近年来研究学者提出的方法进行比较，以此验证本发明的方法的性能表现。其中Fdataset数据集包括593种药物、313种疾病和1933种经过验证的药物-疾病关联，该数据集的药物信息从DrugBank数据库中收集，疾病信息收集自孟德尔人类遗传学(OMIM)在线数据库，该数据库专注于遗传疾病，包括文本信息和相关参考信息、序列记录、地图等相关数据库；而Cdataset数据集包括663种药物、409种疾病以及它们之间的2532种关联。然后，本发明在HDVD数据集上将本发明的方法与同类型的方法进行比较，进一步验证本发明方法在HDVD上的性能。药物-病毒数据集(human drug-virusassociation database，HDVD)收录了34种病毒(其中包括新冠病毒)、219种药物和455种证实的人类药物-病毒相互作用。

1、金标准数据集上的方法性能对比

为了测试本发明的方法的性能，本发明采用了10次10折CV在Fdataset和Cdataset两个金标准数据集上与其他经典的药物重定位方法进行比较。在测试了所有交互作用后，计算了真阳性率(TPR)和假阳性率(FPR)，如下式所示:

其中TP为正确识别的阳性样本数，FN为不正确识别的阴性样本数。TPR是阳性样本在阳性样本中被正确识别的比例。

其中FP为正确识别的阳性样本数，TN为正确识别的阴性样本数。FPR是错误识别的阴性样本占所有阴性样本的比例。

AUC是ROC曲线下的面积，它被广泛用于描述全局预测性能。AUC值越接近1表示性能越好，而AUC值越接近0.5表示性能越接近随机。本发明通过改变预设的排名阈值，分别计算真阳性率(TPR)和假阳性率(FPR)，并绘制ROC曲线。如图2所示，与其他方法相比，本发明的方法得到了更好的结果。本发明提出的方法在Cdataset中的AUC为0.955，比BNNR(0.950)高0.5％。DrugNet的AUC值最低，比本发明的方法低15.1％。本发明的方法在Fdataset中得分最高。本发明的方法的AUC值为0.937，比BNNR(0.933)高0.8％。DrugNet的AUC值最低，比本发明的方法低16.3％。综上所述，本发明的方法在金标准数据集上的性能表现优于图2中的其他方法。

2、HDVD数据集上的方法性能对比

为了验证本发明提出的方法在HDVD数据集上的性能表现，将本方法和其他近年来的同类型的预测方法在HDVD上进行10次10折CV验证，并绘画出每种对比方法的ROC曲线图并计算对应AUC值。如图3(a)所示，发现本发明的方法在10次10折CV实验结果中性能表现为AUC＝0.8771。其他方法分别为：BNNR(AUC＝0.8772)、HGIMC(AUC＝0.7851)、GMCLDA(AUC＝0.7350)和SIMCLDA(AUC＝0.7051)。其中本发明方法的AUC值与BNNR的AUC值接近。

为了进一步评价本发明的方法与BNNR及其他矩阵补全方法在HDVD数据集上的性能对比，本发明在相同的验证框架上计算了Precision指标，并绘制PR曲线和计算对应的AUPR值。由于Precision指标是正确识别的阳性样本在检索样本中所占的百分比。其精度值越大，该模型预测性能越好。因此本发明在HDVD数据集上计算了各方法的精确度(Precision)，并绘制Precision-recall曲线来综合评价本发明的方法在HDVD上的性能。Precision计算方式如下：

如图3(b)所示，本发明的方法AUPR值最高，为0.5060，比其他方法要高，分别是BNNR(AUPR＝0.4888)、HGIMC(AUPR＝0.4174)、GMCLDA(AUPR＝0.1771)和SIMCLDA(AUPR＝0.1518)。综上所述，本发明的方法在HDVD数据集的性能表现优于图3中的其他方法。

3、消融实验

为验证本发明提出的模型各部分的有效性，在Fdataset和Cdataset两个金标准数据集进行了消融实验。各种模块的组合对比如表1所示，其中模型1为本文所提出的方法(即WKNKN+CKA-MKL+压缩感知)，模型2为WKNKN算法、相似性均值融合和压缩感知的组合，模型3为CKA-MKL与压缩感知的组合，模型4为相似性均值融合与压缩感知的组合，模型5为WKNKN、单一相似性(没有引入药物和疾病的杰卡德相似性和高斯核相似性)和压缩感知的组合，模型6为单一相似性(没有引入药物和疾病的杰卡德相似性和高斯核相似性)和压缩感知方法。从图4可知，模型1与模型2相比略有提升，这是因为利用CKA-MKL处理多种相似性矩阵可以得到每个相似性矩阵(核)的权重，并选择最佳的矩阵组合，提高各种相似信息的互补性。模型1与模型3、模型2与模型4和模型5与模型6这3组对比模型均说明引入WKNKN对药物-疾病关联矩阵进行预处理能提升矩阵的奇异值采样数据对原矩阵数据的代表性，提高关联预测的精度。模型3与模型6的对比说明了通过富集多种相似性信息能在一定程度上提高药物重定位的效果。

表1

4、案例分析

为了确定预测的关联对是否真实或者是否具有参考价值，本实施例在HDVD数据集中选择新冠病毒作为应用对象用来验证本方法预测的潜在药物，并列出前15名候选药物的确认信息。在得到预测结果之后，通过权威公共数据库如DrugBank、CTD、KEGG和部分研究文献对候选药物对新冠病毒的适用性进行验证。预测结果和验证汇总见表2。在根据最终预测的关联分数排名的前15名候选药物中，有12种药物(80％的成功率)得到了各种证据的验证。这样的结果表明，本发明的方法在实际应用中预测新冠病毒的适应药物是有效可行的。

表2

尽管已经示出和描述了本发明的实施例，对于本领域的普通技术人员而言，可以理解在不脱离本发明的原理和精神的情况下可以对这些实施例进行多种变化、修改、替换和变型，本发明的范围由所附权利要求及其等同范围限定。

Claims (7)

1.一种预测新冠病毒适应药物的药物重定位方法，其特征在于，包括：

步骤1：采用WKNKN算法处理人类药物-病毒关联数据集HDVD，所述数据集包括药物-疾病关联邻接矩阵Y、药物化学结构相似度矩阵、疾病语义相似度矩阵；挖掘潜在相互作用的似然分数，降低药物-疾病关联矩阵Y的稀疏性，得到处理后的人类药物-病毒关联数据集HDVD；

步骤2：根据处理后的人类药物-病毒关联数据集HDVD分别计算药物和病毒的杰卡德相似性和高斯核相似性，得到药物和疾病的杰卡德相似矩阵DS_jac和dS_jac，得到疾病的高斯核相似矩阵dS_Gaus以及药物的高斯核相似矩阵DS_Gaus；

步骤3：分别结合药物的化学相似性和病毒的序列相似性利用中心核对称多核学习方式进行整合，从而分别获得药物和病毒的综合相似性矩阵；

步骤4：结合处理后的人类药物-病毒关联数据集HDVD信息，通过压缩感知算法补全原人类药物-病毒关联数据集中缺失的关联，得到预测概率矩阵；

步骤5：将得到预测概率矩阵作为药物-病毒的预测结果矩阵。

2.如权利要求1所述的预测新冠病毒适应药物的药物重定位方法，其特征在于，步骤1包括以下子步骤：

子步骤1、对于每个药物D_i，采用与其最接近的K种已知药物的化学相似度矩阵DS_chem及其对应的关联作用来估计D_i的相互作用似然分数，推导公式如下:

其中D₁到

表示药物D_i的K个近邻药物，按降序排列；

是权重系数，其中η是衰减项并满足η≤1，

是正则化项；

子步骤2、相似地，对于每个疾病d_j，采用与其最接近的K种已知疾病的语义相似性及其对应的交互作用来估计d_j的交互似然分数，推导公式如下:

其中d₁ to

表示d_j的K个近邻药物，按降序排列；

是权重系数，其中η是衰减项并满足η≤1，

是正则化项；

子步骤3、采用下式计算得到处理后的药物-疾病关联邻接矩阵Y_F：

3.如权利要求2所述的预测新冠病毒适应药物的药物重定位方法，其特征在于，步骤2中，采用下式计算药物的杰卡德相似矩阵DS_jac：

其中Γ(x)是与药物D_x相关联的疾病集合，Γ(y)是与药物D_y相关联的疾病集合；

同理，按照上述方式计算出疾病的杰卡德相似矩阵dS_jac。

4.如权利要求3所述的预测新冠病毒适应药物的药物重定位方法，其特征在于，步骤2中，将疾病d_i和d_j之间的高斯核相似性表示为dS_Gaus，药物D_i和D_j之间的高斯核相似性表示为DS_Gaus，计算方法如下：

dS_Gaus(d_i，d_j)＝exp(-α_d||Q(d_i)-Q(d_j)||²)；

DS_Gaus(D_i，D_j)＝exp(-α_D||Q(D_i)-Q(D_j)||²)；

在药物-疾病关联邻接矩阵Y_F中，第i行表示药物D_i是否与每种疾病有关，第j列表示疾病d_j是否与每种药物有关；向量Q(D_i)和Q(d_j)分别将第i行向量和第j列向量表示为高斯核的特征向量；核带宽α_d和α_D定义为：

其中，初始核带宽系数α′_d和α′_D均设置为1。

5.如权利要求4所述的预测新冠病毒适应药物的药物重定位方法，其特征在于，步骤3中，使用CKA-MKL算法分别在两个空间组合三个相似核矩阵，其中，特征核K^*采用下式计算：

K_i∈R^N×N，

其中k为核数，K_i分别在药物相似核集DS_chem、DS_jac、DS_Gaus和疾病相似核集dS_sem、dS_jac、dS_Gaus中取值，ω_i为核K_i的权值，N是节点数；

内核对齐的值定义如下：

其中P，Q∈R^N×N，N∈{n_D，n_d}，<P，Q>_F＝Trace(P^TQ)是Frobenius内积，

是Frobenius范数；

根据上式定义经过中心核对齐后的目标函数如下：

式中

表示一个中心核矩阵，I_N∈R^N×N是一个单位矩阵，I_N是单位向量，再将目标函数转变为：

其中，a∈R^k×1和M∈R^k×k分别由下式计算而得：

将最终目标函数表示为：

将上式最小化计算疾病和药物相似矩阵组合的最优权值ω_D和ω_d，并分别将它们与药物DS_chem、DS_jac、

和疾病dS_sem、dS_jac、

相似核按下式组合：

其中

和

和

分别是疾病核和药物核的最优权值；

分别计算药物和疾病的最佳组合相似矩阵。

6.如权利要求5所述的预测新冠病毒适应药物的药物重定位方法，其特征在于，步骤4中，压缩感知算法的计算过程为：

构建初始子矩阵

和

根据下式计算初始概率矩阵：

其中p_i，j是药物D_i和疾病d_j的相互作用概率，f_i是药物子矩阵F的第i行，g_j是疾病子矩阵G的第j行，w_i，j是药物D_i和疾病d_j的初始权值；

构建损失函数式：

引入两个额外的可调参数λ_M和λ_N，转化得到损失函数为：

其中，

根据下式计算子矩阵F和矩阵G的偏导数：

其中

代表Hadamard乘积；

根据下式分别更新子矩阵F和矩阵G：

其中，k为迭代过程的学习率；上标n是当前迭代次数；子矩阵F和G更新的结束条件为max(ΔF，ΔG)＜10^-5。

7.如权利要求6所述的预测新冠病毒适应药物的药物重定位方法，其特征在于，步骤4中，根据下式，用矩阵Fⁿ⁺¹和Gⁿ⁺¹计算概率矩阵Pⁿ⁺¹：

Images