CN108509771A - 一种基于稀疏匹配的多组学数据关联关系发现方法 - Google Patents

Abstract

本发明公开了一种基于稀疏匹配的多组学数据关联关系发现方法，包括：对输入数据进行预处理，提高数据的质量；根据数据特点，选择合适的相似性度量，计算数据特征之间的相似性矩阵；基于特征之间的相似性网络，融入先验信息，挖掘数据特征之间的潜在关联关系。本发明方法可以充分利用现有已被证实的组学数据的特征的先验信息，减少噪音对结果的影响、降低数据误差所带来的不确定性，提高结果的精确性与鲁棒性。

Description

一种基于稀疏匹配的多组学数据关联关系发现方法

技术领域

本发明涉及生物信息技术领域，特别涉及一种基于稀疏匹配的多组学数据关联关系发现方法。

背景技术

生物组学由基因组学(Genomics)、转录组学(Transcriptomics)、蛋白质组学(Proteomics)和代谢组学(Metabolomics)共同构成，旨在从整体的角度出发，研究人类基因、核糖核酸、蛋白质及其代谢物等的相互的作用，通过整合分析人体各个层次组织的内在关联，为探索人类疾病的发病机制提供更加科学全面的方法。

随着科学技术的发展，高通量测序技术的出现极大地降低了测序的成本，提高了测序的性能，使得高效全面地测定同一样本不同层次的组学数据成为可能。TCGA(TheCancer Genome Atlas)数据库整合了33种癌症、26种组织类型、超过11000个肿瘤病例、大小约2.5PB的组学测序数据，其中绝大部分病例均同时具有单核苷酸变异数据、RNA-seq测序数据、基因拷贝数变异(CNV)、甲基化数据(DNA methylation)和临床数据等多个组学的数据。这些海量多源异构数据能在基因组、转录组、蛋白质组、代谢组水平多层次多方位提供描述肿瘤发生发展的重要信息，坚实地奠定了整合分析多组学数据关联关系发现的数据来源基础，为定性定量地揭示肿瘤发生发展的机理提供了的前所未有的机遇。

相关科学理论的不断深入以及相关的实践研究不断增多，使得研究人员发现基于单一组学数据的研究分析方法只能从单一的角度理解疾病发生发展在分子水平上的变化特征，已经无法满足进一步解释复杂的疾病机理过程的需求。整合多个组学数据能够补偿单一组学数据中的数据丢失及纠正部分错误信息，从整体的视角全面多方位地描述疾病发生发展过程中在各个组学层次上的变化特征，有效地降低了基于单一组学数据研究实验中出现假阳性的概率。鉴于此，研究学者致力于尝试整合分析多个组学层面的数据，全面研究病变过程中多组学数据之间内在的相互关联关系及其协同变化规律，提高对生物体系统机制的认识。例如，研究学者基于贝叶斯网络推理模型，对黑色素瘤病人的染色体拷贝数变异数据与基因表达数据进行整合分析，以检测诱发黑色素瘤的产生因素，最终证实了抗体TBC1D16及蛋白质RAB27A的异常调控会导致黑色素瘤的增殖。研究人员使用酵母种群的基因表达数据、转录因子结合位点(TFBS)和蛋白质-蛋白质相互作用(PPI)数据重构了基因调控网络，并通过预测常见的酵母群体基因表达活性调控证明其有效性。生物学家将单核苷酸多态性数据和mRNA基因表达做串联整合，随机选择变量作为层次贝叶斯模型的数据输入，检测单核苷酸多态性、基因表达变异和表型变化之间的关联，等等。虽然现有模型都通过整合多个组学数据以发掘它们之间的潜在关联关系，降低了单一组学数据随机性误差所造成的假阳性率，但大多数模型都存在一些不足之处，如整合研究方法过于简单、数据来源没有统一性、样本量不足等等。

同时，绝大部分的研究都仅仅聚焦于组学数据本身，很少考虑在模型中加入其他重要的先验信息，甚至完全没有利用现有的已证实的信息，而对先验信息的合理利用对模型准确性、健壮性和执行效率的提高有明显的积极影响。随着生物信息学的发展，越来越多的组学数据之间的潜在关联已被证实，如基因调控网络、蛋白质-蛋白质相互作用网络和代谢路径网络等等。充分地利用已被证实及广泛接受的先验信息能够降低因生物数据噪音大而导致结果误差大的影响，提高结果的准确性及可靠性。

发明内容

本发明的目的在于克服现有技术的缺点与不足，提供一种基于稀疏匹配的多组学数据关联关系发现方法，可以有效充分地利用被证实的先验组学信息，减少噪音对结果的影响、降低数据误差所带来的不确定性和提高结果的精确性与鲁棒性。

本发明的目的通过以下的技术方案实现：一种基于稀疏匹配的多组学数据关联关系发现方法，包括：

S1、对输入数据进行预处理；

S2、根据数据特点，选择合适的相似性度量，计算数据特征之间的相似性矩阵；

S3、基于特征之间的相似性网络，融入先验信息，挖掘数据特征之间的潜在关联关系。

优选的，预处理步骤包括：

(1)数据清洗：对于数据中缺失值，采取三次样条插值法对其进行补插；对于数据中的异常值，采用平均值对其进行替换修正；

(2)数据变换：对数据进行最小-最大规范化，使得所有数据集的取值范围均为[0,1]，消除指标间量纲影响，对数据进行变换的公式为：

其中，D_new为进行变换后所得的数据，D_old为原始数据，D_min为数据中的最小值，D_max为数据中的最大值；

(3)数据归约：利用主成分分析对数据进行降维处理：将原始数据按行排列组成矩阵X；对X进行数据标准化，使其均值变为零；求X的协方差矩阵C；将特征向量按特征值由大到小排列，取前k个按行组成矩阵P；通过计算Y＝PX，得到降维后数据Y。

优选的，选取最大信息系数、皮尔逊相关系数和互信息量中的一个或者多个作为特征相关性的度量方法。

具体的，通过下面的公式计算组学数据X和组学数据Y之间的互信息熵：

其中，p(x)为X的概率分布，p(y)为X的概率分布，p(x,y)为X和Y的联合概率分布；

通过下面的公式计算组学数据X和组学数据Y之间的皮尔逊相关系数：

其中，是基因X的均值，是药物Y的均值；

通过下面的公式计算组学数据X和组学数据Y之间的最大信息系数：

其中，max I(D|_G)是指在x列和y行的网格分割G上具有最大互信息量的分布D|_G。

进一步的，通过下面的公式计算组学数据X和组学数据Y之间最终的相似性度量：

优选的，建立模型发掘组学数据之间特征关联关系：令M₁和M₂表示利用高通量测序技术对同一批病人测序所得的两个组学数据矩阵，两个矩阵的维度分别m×n₁和m×n₂，其中m表示病人样本个数，n₁表示M₁中的特征个数，n₂表示M₂中的特征数量；在两个组学数据中找出相互关联的特征，此问题在数学上可以等价描述为如何确定一个n₁×n₂维的0-1矩阵X：

在模型中引入现有已被证实的先验信息，使用表示组学数据M₁的特征的相互关联网络，表示组学数据M₂的特征的相互关联网络，它们的元素的取值范围为[0,1]，表示组学数据自身特征之间的相关程度，如0则表示完全不相关；基于此，提出基于稀疏匹配的组学数据关联关系发现模型如下：

s.t.Xl₁≤b₁

X^Tl₂≤b₁

X_ij∈{0,1}

其中，矩阵X为要求解的目标变量；l₁和b₁均为n₂×1维的列向量，其中l₁的元素值全部为1，而b₁的元素值均为一常量参数c₁；类似地，l₂和b₂都是维数为n₁×1的列向量，其中l₂的元素值全部为1，而b₂的元素值均为c₂，c₂是一常量参数；矩阵H为组学数据M₁中的所有特征与组学数据中M₂的全部特征之间的相似性度量矩阵，维度为n₁×n₂；λ₁、λ₂和β为对应正则项的参数，用于调节各个正则项对最终结果的影响；为通过组学数据M₁的先验关联权重矩阵计算得出的拉普拉斯矩阵，其计算方式如下面公式所示；相似地，为利用组学数据M₂的先验关联权重矩阵计算得出的标准化拉普拉斯矩阵，其计算过程与类似；

其中，是一个n₁×n₁维的对角矩阵，其对角元素的值等于中该元素所在列的列和。

进一步的，求解组学数据关联关系发现模型的方法如下：

8)把模型最优化问题转化为一般形式：

s.t.g_k(X)≤0,(k＝1,2,...,N)

其中,N为把矩阵形式的不等式约束改写为普通形式后不等式约束的个数，即N＝p+2q，且p＝n₁+n₂和q＝n₁×n₂；g_k(X)的表达式为：

9)构造最优化问题的增广拉格朗日乘子罚函数如下：

其中，γ＝(γ₁，γ₂，...，γ_N)^T是不等式约束的拉格朗日乘子向量，M为罚因子；于是将1)中带约束的优化问题转化为无约束的优化问题：

10)以X^(n-1)为初始点，利用梯度下降法求解的极小值，其中，关于X_ij的偏导数为：

11)判断步骤3)中所求的的结果X⁽ⁿ⁾是否满足收敛条件：若||g(X⁽ⁿ⁾)||＜ε，表示X⁽ⁿ⁾满足算法的收敛条件，则停止计算，以X⁽ⁿ⁾为最终结果；否则转步骤5)；

12)判断迭代过程中每一步的收敛速度是否过慢：若则表示当前迭代收敛速度过慢，需要增大罚因子从而加快算法的收敛速度，令M＝ρM,转步骤6)；否则直接转步骤6)；

13)更新拉格朗日乘子向量

14)更新迭代次数n＝n+1，返回步骤1)。

优选的，得到结果矩阵X后，利用矩阵全部元素的均值作为阈值对其进行二值化，大于矩阵均值的设置为1,小于矩阵均值的设置为0；若X_ij的值为1,说明第i个基因和第j个药物之间存在关联关系。

本发明与现有技术相比，具有如下优点和有益效果：

本发明提供了一套完整的包括数据预处理，特征相似性度量，发掘特征关联关系的流程方法。本发明在模型中引入现有已被证实的先验信息，减少了噪音对结果的影响、降低数据误差所带来的不确定性和提高结果的精确性与鲁棒性，使得本方法的准确率得到了相比传统方法得到了显著提高。

附图说明

图1是实施例方法的基本流程图。

图2是实施例方法的具体示意图。

具体实施方式

下面结合实施例及附图对本发明作进一步详细的描述，但本发明的实施方式不限于此。

实施例1

一种基于稀疏匹配的多组学数据关联关系发现方法，基本流程如图1所示，包括：读入数据、数据清洗、数据变换、数据归约、计算特征相似矩阵、挖掘数据关联关系、得出结果。

基于稀疏匹配的多组学数据关联关系发现的方案如下：

1、对输入数据进行预处理，提高数据的质量，让数据更好地适应本方法。

海量的原始数据中存在着大量异常的数据，维度过大等多个问题，严重影响到模型的执行效率，甚至可能导致模型结果的偏差，所以进行数据预处理就显得尤为重要，本方法通过以下步骤对数据进行预处理：

1.1数据清洗：

对于数据中缺失值，采取三次样条插值法对其进行补插；对于数据中的异常值，采用平均值对其进行替换修正。

1.2数据变换：

对数据进行最小-最大规范化，使得所有数据集的取值范围均为[0,1]，消除指标间量纲影响。对数据进行变换的公式为：

其中，D_new为进行变换后所得的数据，D_old为原始数据，D_min为数据中的最小值，D_max为数据中的最大值。

1.3数据归约：

在保证数据信息丢失最少的原则下，用较少的综合特征代替原来全部的特征，利用主成分分析对数据进行降维处理：(1)将原始数据按行排列组成矩阵X；(2)对X进行数据标准化，使其均值变为零；(3)求X的协方差矩阵C；(4)将特征向量按特征值由大到小排列，取前k个按行组成矩阵P；(5)通过计算Y＝PX，得到降维后数据Y。

2、根据数据特点，选择合适的相似性度量，计算数据特征之间的相似性矩阵。

生物组学数据特征之间的相似性度量存在多种方法，互信息量是大多数研究者的选择。但是互信息量只能够有效地度量特征之间的线性和抛物线性，而无法有效地度量更广泛的统计特性(如指数性、正弦性和周期性等)，导致基于互信息量的生物组学数据特征之间的相似性度量存在一定局限性。Reshef于2011年提出了比互信息量能够有效度量更多种类的变量相关性的最大信息系数。为了综合利用多种相似性度量的特点，本方法选取综合最大信息系数、皮尔逊相关系数和互信息量作为特征相关性的度量方法。

3、基于特征之间的相似性网络，融入先验信息，挖掘数据特征之间的潜在关联关系。

组学数据特征之间的关联具有以下两个特点：第一，组学数据的特征关联具有稀疏性的特点；第二，组学数据特征之间的关联是“多对多”的关系。基于这些特性，建立模型发掘组学数据之间特征关联关系。令M₁和M₂表示利用高通量测序技术对同一批病人测序所得的两个组学数据矩阵，如基因表达数据、RNA-seq测序数据、甲基化数据、药物反应数据等。两个矩阵的维度分别m×n₁和m×n₂，其中m表示病人样本个数，n₁表示M₁中的特征个数，n₂表示M₂中的特征数量。问题是如何在两个组学数据中找出相互关联的特征，此问题在数学上可以等价描述为如何确定一个n₁×n₂维的0-1矩阵X：

此外，为了减少噪音对结果的影响、降低数据误差所带来的不确定性和提高结果的精确性与鲁棒性，在模型中引入现有已被证实的先验信息，如基因调控网络、蛋白质-蛋白质相互作用网络和代谢路径网络等等。使用表示组学数据M₁的特征的相互关联网络，表示组学数据M₂的特征的相互关联网络，它们的元素的取值范围为[0,1]，表示组学数据自身特征之间的相关程度，如0则表示完全不相关。基于此，提出基于稀疏匹配的组学数据关联关系发现模型(Network-based Sparse Binary Matching Model，NSBM)如下：

其中，矩阵X为要求解的目标变量；l₁和b₁均为n₂×1维的列向量，其中l₁的元素值全部为1，而b₁的元素值均为一常量参数c₁；类似地，l₂和b₂都是维数为n₁×1的列向量，其中l₂的元素值全部为1，而b₂的元素值均为c₂，c₂是一常量参数；矩阵H为组学数据M₁中的所有特征与组学数据中M₂的全部特征之间的相似性度量矩阵，维度为n₁×n₂。本方法采用两特征之间皮尔森相关系数的负绝对值作为度量标准，取值范围为[-1,0]，值越趋于-1，说明两特征之间的相关性越大，反之亦然；λ₁、λ₂和β为对应正则项的参数，用于调节各个正则项对最终结果的影响；为通过组学数据M₁的先验关联权重矩阵计算得出的拉普拉斯矩阵，其计算方式如公式(1-2)所示；相似地，为利用组学数据M₂的先验关联权重矩阵计算得出的标准化拉普拉斯矩阵，其计算过程与类似。

其中，是一个n1×n1维的对角矩阵，其对角元素的值等于中该元素所在列的列和。

在最优化模型(1-1)的目标函数中，第一项tr(HX^T)根据组学数据M₁中第i个特征与组学数据M₂中第j个特征的相似性来影响结果X_ij的值，组学数据M₁中第i个特征与组学数据M₂中第j个特征的相关性越大，则X_ij的值越趋于1，反之亦然；第二项的作用为根据组学数据M₁的特征的相互关联网络来影响结果X，先验相互关联网络中越相似的M₁的两个特征，在结果X中对于M₂中的同一个特征的值越相近。λ₁则控制调节此正则项对结果影响程度的大小；第三项的作用与第二项类似，融入组学数据M₂的特征的先验关联网络来提高结果X的精度，其影响程度的大小由参数λ₂控制调节；第四项β||X||₀为L₀-范数，定义为矩阵X中非0元素的个数，用来描述结果的稀疏性，非0元素越少，说明结果越稀疏，并由参数β来控制结果的稀疏程度；第一个约束条件Xl₁≤b₁用于约束与组学数据M₁中的某个特征关联的M₂的特征的个数不超过一定的数目。b₁的元素值可以设为不同的值，表示对M₁中的不同特征进行不同的约束。本方法中为简化模型，降低复杂度，对b₁的元素值均设为某一常量值c₂；第二个约束条件X^Tl₂≤b₁的作用与第一个约束条件的作用相类似。

式(1-1)是一个带线性约束的二次0-1规划问题，它的求解是一个NP难问题。为解决此问题，选择放松对X的取值范围的约束，不再限定它的取值只能是0或1，而假设它是一个取值范围为[0，1]的连续变量。其次，目标函数的第四项为L₀-范数，通过对其最小化来对结果的稀疏性来进行约束。而L₀-范数的最小化问题已被证实为NP难问题，现在有的计算技术没法对其进行有效求解。由于L₁-范数是L₁-范数的最优凸近似，所以实际中通常的解决办法用L₁-范数对其进行近似替换，利用L₁-范数对结果的稀疏性进行约束，本发明也采取这一策略。另外，由于前面假设X的元素的取值范围为[0，1]，所以β||X||₀可以等价地改写为βtr(UX^T)，其中U是一个n₁×n₂，元素值全为1的矩阵。所以，式(1-1)被改写为：

s.t.Xl₁≤b₁

X^Tl₂≤b₁

0≤X≤U (1-3)

式(1-3)是带线性不等式约束的凸二次规划问题，所以可以通过凸优化理论的方法来对结果进行求解，以保证其收敛性。在本发明中，通过利用增广拉格朗日乘子法(Augmented Lagrangian Multiplier Method，ALM)来对模型进行快速求解。增广拉格朗日乘子法具有完善的数学基础，是一种被广泛应用的优化算法，具有收敛速度快、普适性高等特点。它在罚函数法的基础上引入了拉格朗日乘子项，在迭代求解的过程中自适应地修正拉格朗日乘子与惩罚因子的值，并防止罚因子趋于无穷大，有效地避免了罚函数法中由于替代函数的海塞矩阵条件数无限增大而引起的病态问题。

最优化问题(1-3)的一般化形式为：

其中，N为把式(1-3)中矩阵形式的不等式约束改写为普通形式后不等式约束的个数，即N＝p+2q，且p＝n₁+n₂和q＝n₁×n₂。g_k(X)的表达式为：

通过引入松弛变量y_k，把带不等式约束的优化问题(1-4)转化为只有等式约束的最小化问题：

于是可以构造优化问题(1-5)的增广拉格朗日函数如下：

其中，γ＝(γ₁，γ₂，...，γ_N)^T是不等式约束的拉格朗日乘子向量；M为罚因子。于是，带约束的优化问题(1-4)就转化为无约束的优化问题：

首先，利用配方法将化为：

然后，令关于y_k的偏导数为0：

由式子(1-7)可知，当对于y_k取极小值时，y_k的取值如下：

当γ_k-Mg_k(X)≥0时，y_k的取值为：

当γ_k-Mg_k(X)＜0时，y_k的取值为：

y_k＝0

可将上述两种情况的结果整合到一个表达式，即有：

将式(1-8)代入(1-6)中，即可得到不等式约束问题(1-4)的乘子罚函数：

在每一步迭代求解的过程中，选择使用无约束最小优化方法中的梯度下降法求的极小值。关于X_ij的偏导数为：

另外，由经过n次迭代的γ⁽ⁿ⁾代得到第n+1次迭代的γ⁽ⁿ⁺¹⁾的修正公式为：

最后，由于问题是凸优化问题，所以在经过有限次的迭代步骤后，可以得到原问题的最优解。

实施例2

一种基于稀疏匹配的多组学数据关联关系发现方法，如图2所示，是方法的具体示意图，包括：

步骤1：获取基因表达数据M₁和药物反应数据M₂、基因-基因关联网络W₁和药物-药物关联网络W₂；

步骤2：分别对基因表达数据M₁和药物反应数据M₂进行数据预处理。

步骤3：利用基因表达数据M₁和药物反应数据M₂计算出基因-药物的相似性度量矩阵H。

步骤4：把基因-药物的相似性度量矩阵H、基因-基因关联网络W₁和药物-药物关联网络W₂，作为模型输入，挖掘基因和药物之间的关联关系。

步骤2具体包括：

步骤2.1：对数据进行清洗。对于数据中缺失值，采取三次样条插值法对其进行补插；对于数据中的异常值，采用平均值对其进行替换修正。

步骤2.2：对数据进行变换。对数据进行最小-最大规范化，使得所有数据集的取值范围均为[0,1]，消除指标间量纲影响。对数据进行变换的公式为：

其中，D_new为进行变换后所得的数据，D_old为原始数据，D_min为数据中的最小值，D_max为数据中的最大值。

步骤2.3：对数据进行规约。为了在保证数据信息丢失最少的原则下，用较少的综合特征代替原来全部的特征，利用主成分分析对基因表达数据M₁和药物反应数据M₂进行降维处理：(1)将原始数据按行排列组成矩阵X；(2)对X进行数据标准化，使其均值变为零；(3)求X的协方差矩阵C；(4)将特征向量按特征值由大到小排列，取前k个按行组成矩阵P；(5)通过计算Y＝PX，得到降维后数据Y。

步骤3具体包括：

步骤3.1：通过下面的公式计算基因X和药物Y之间的互信息熵：

其中，p(x)为X的概率分布，p(y)为X的概率分布，p(x,y)为X和Y的联合概率分布。

步骤3.2：通过下面的公式计算基因X和药物Y之间的皮尔逊相关系数：

其中，是基因X的均值，是药物Y的均值。

步骤3.3：通过下面的公式计算基因X和药物Y之间的最大信息系数：

其中，maxI(D|_G)是指在x列和y行的网格分割G上具有最大互信息量的分布D|_G。

步骤3.4：通过下面的公式计算基因X和药物Y之间最终的相似性度量：

步骤4具体包括：

步骤4.1：对模型的参数进行设置：X⁽⁰⁾＝0，γ⁽¹⁾＝0，M＝20，ρ＝2，δ＝0.1，ε＝0.001，n＝1。其中，X是求解结果；γ是拉格朗日乘子向量；M为罚因子；δ用于判断迭代求解过程中算法的收敛速度是否过慢；ρ是当算法的收敛速度过慢时增大罚因子的倍数；ε是用于判断最终结果是否收敛的精度。

步骤4.2：对模型进行求解：

1)把最优化问题转化为一般形式：

s.t.g_k(X)≤0,(k＝1,2,...,N)

其中,N为把矩阵形式的不等式约束改写为普通形式后不等式约束的个数，即N＝p+2q，且p＝n₁+n₂和q＝n₁×n₂。g_k(X)的表达式为：

2)构造最优化问题的增广拉格朗日乘子罚函数如下：

3)以X^(n-1)为初始点，利用梯度下降法求解由上面步骤所得的无约束问题，记结果为X⁽ⁿ⁾，其中，关于X_ij的偏导数为：

4)判断步骤3)中所求的的结果X⁽ⁿ⁾是否满足收敛条件：若||g(X⁽ⁿ⁾)||＜ε，表示X⁽ⁿ⁾满足算法的收敛条件，则停止计算，以X⁽ⁿ⁾为最终结果；否则转步骤5)；

5)判断迭代过程中每一步的收敛速度是否过慢：若则表示当前迭代收敛速度过慢，需要增大罚因子从而加快算法的收敛速度，令M＝ρM,转步骤6)；否则直接转步骤6)；

6)更新拉格朗日乘子向量

7)更新迭代次数n＝n+1，返回步骤1)；

步骤4.3：经过步骤4.2，得到结果矩阵X。利用矩阵全部元素的均值作为阈值对其进行二值化，大于矩阵均值的设置为1,小于矩阵均值的设置为0。若X_ij的值为1,说明第i个基因和第j个药物之间存在关联关系。

上述实施例为本发明较佳的实施方式，但本发明的实施方式并不受上述实施例的限制，其他的任何未背离本发明的精神实质与原理下所作的改变、修饰、替代、组合、简化，均应为等效的置换方式，都包含在本发明的保护范围之内。

Claims (8)

1.一种基于稀疏匹配的多组学数据关联关系发现方法，其特征在于，包括步骤：

S1、对输入数据进行预处理；

S2、根据数据特点，选择合适的相似性度量，计算数据特征之间的相似性矩阵；

S3、基于特征之间的相似性网络，融入先验信息，挖掘数据特征之间的潜在关联关系。

2.根据权利要求1所述的基于稀疏匹配的多组学数据关联关系发现方法，其特征在于，S1预处理步骤包括：

(1)数据清洗：对于数据中缺失值，采取三次样条插值法对其进行补插；对于数据中的异常值，采用平均值对其进行替换修正；

(2)数据变换：对数据进行最小-最大规范化，使得所有数据集的取值范围均为[0,1]，消除指标间量纲影响，对数据进行变换的公式为：

其中，D_new为进行变换后所得的数据，D_old为原始数据，D_min为数据中的最小值，D_max为数据中的最大值；

(3)数据归约：利用主成分分析对数据进行降维处理：将原始数据按行排列组成矩阵X；对X进行数据标准化，使其均值变为零；求X的协方差矩阵C；将特征向量按特征值由大到小排列，取前k个按行组成矩阵P；通过计算Y＝PX，得到降维后数据Y。

3.根据权利要求1所述的基于稀疏匹配的多组学数据关联关系发现方法，其特征在于，步骤S2中，选取最大信息系数、皮尔逊相关系数和互信息量中的一个或者多个作为特征相关性的度量方法。

4.根据权利要求3所述的基于稀疏匹配的多组学数据关联关系发现方法，其特征在于，通过下面的公式计算组学数据X和组学数据Y之间的互信息熵：

其中，p(x)为X的概率分布，p(y)为X的概率分布，p(x,y)为X和Y的联合概率分布；

通过下面的公式计算组学数据X和组学数据Y之间的皮尔逊相关系数：

其中，是基因X的均值，是药物Y的均值；

通过下面的公式计算组学数据X和组学数据Y之间的最大信息系数：

其中，maxI(D|_G)是指在x列和y行的网格分割G上具有最大互信息量的分布D|_G。

5.根据权利要求4所述的基于稀疏匹配的多组学数据关联关系发现方法，其特征在于，通过下面的公式计算组学数据X和组学数据Y之间最终的相似性度量：

6.根据权利要求1所述的基于稀疏匹配的多组学数据关联关系发现方法，其特征在于，步骤S3中，建立模型发掘组学数据之间特征关联关系：令M₁和M₂表示利用高通量测序技术对同一批病人测序所得的两个组学数据矩阵，两个矩阵的维度分别m×n₁和m×n₂，其中m表示病人样本个数，n₁表示M₁中的特征个数，n₂表示M₂中的特征数量；在两个组学数据中找出相互关联的特征，此问题在数学上可以等价描述为如何确定一个n₁×n₂维的0-1矩阵X：

在模型中引入现有已被证实的先验信息，使用表示组学数据M₁的特征的相互关联网络，表示组学数据M₂的特征的相互关联网络，它们的元素的取值范围为[0,1]，表示组学数据自身特征之间的相关程度，如0则表示完全不相关；基于此，提出基于稀疏匹配的组学数据关联关系发现模型如下：

s.t.Xl₁≤b₁

X^Tl₂≤b₁

X_ij∈{0,1}

其中，矩阵X为要求解的目标变量；l₁和b₁均为n₂×1维的列向量，其中l₁的元素值全部为1，而b₁的元素值均为一常量参数c₁；类似地，l₂和b₂都是维数为n₁×1的列向量，其中l₂的元素值全部为1，而b₂的元素值均为c₂，c₂是一常量参数；矩阵H为组学数据M₁中的所有特征与组学数据中M₂的全部特征之间的相似性度量矩阵，维度为n₁×n₂；λ₁、λ₂和β为对应正则项的参数，用于调节各个正则项对最终结果的影响；为通过组学数据M₁的先验关联权重矩阵计算得出的拉普拉斯矩阵，其计算方式如下面公式所示；相似地，为利用组学数据M₂的先验关联权重矩阵计算得出的标准化拉普拉斯矩阵，其计算过程与类似；

其中，是一个n₁×n₁维的对角矩阵，其对角元素的值等于中该元素所在列的列和。

7.根据权利要求6所述的基于稀疏匹配的多组学数据关联关系发现方法，其特征在于，求解组学数据关联关系发现模型的方法如下：

1)把模型最优化问题转化为一般形式：

s.t.g_k(X)≤0,k＝1,2,...,N

其中,N为把矩阵形式的不等式约束改写为普通形式后不等式约束的个数，即N＝p+2q，且p＝n₁+n₂和q＝n₁×n₂；g_k(X)的表达式为：

2)构造最优化问题的增广拉格朗日乘子罚函数如下：

其中，γ＝(γ₁，γ₂，...，γ_N)^T是不等式约束的拉格朗日乘子向量，M为罚因子；于是将1)中带约束的优化问题转化为无约束的优化问题：

3)以X^(n-1)为初始点，利用梯度下降法求解的极小值，其中，关于X_ij的偏导数为：

4)判断步骤3)中所求的的结果X⁽ⁿ⁾是否满足收敛条件：若||g(X⁽ⁿ⁾)||＜ε，表示X⁽ⁿ⁾满足算法的收敛条件，则停止计算，以X⁽ⁿ⁾为最终结果；否则转步骤5)；

5)判断迭代过程中每一步的收敛速度是否过慢：若则表示当前迭代收敛速度过慢，需要增大罚因子从而加快算法的收敛速度，令M＝ρM,转步骤6)；否则直接转步骤6)；

6)更新拉格朗日乘子向量

7)更新迭代次数n＝n+1，返回步骤1)。

8.根据权利要求6所述的基于稀疏匹配的多组学数据关联关系发现方法，其特征在于，得到结果矩阵X后，利用矩阵全部元素的均值作为阈值对其进行二值化，大于矩阵均值的设置为1,小于矩阵均值的设置为0；若X_ij的值为1,说明第i个基因和第j个药物之间存在关联关系。