CN113588847A - 一种生物代谢组学数据处理方法、分析方法及装置和应用 - Google Patents

Abstract

本发明公开了一种生物代谢组学数据处理方法、分析方法及装置和应用。生物代谢组学数据处理方法包括筛选出多个待测样品中的差异代谢物的步骤，生物代谢组学数据包括待测样品的液相色谱‑质谱数据，筛选的方法包括：（1）对生物代谢组学数据进行预处理和标准化处理，得到二维矩阵数据；（2）缺失峰补齐处理：对所述二维矩阵数据进行稀疏重构处理，再进行降维处理；（3）对步骤（2）得到的数据进行化学计量分析，筛选出待测样品中的差异代谢物。本发明的方法可用于分析不断有新样品增加，样品收集时间长的代谢组学数据，识别和确定差异代谢物的结构，准确度和精度均较高。

Description

一种生物代谢组学数据处理方法、分析方法及装置和应用

技术领域

本发明属于代谢组学技术领域，尤其涉及一种生物代谢组学数据处理方法、分析方法及装置和应用。

背景技术

代谢组学的应用极为广泛，几乎已经渗透到生活的各个方面。尤其是在现代生命科学和医药科学领域，应用更为普遍，其应用范围包括临床医学、药物研发、生物农业、食品安全、环境科学以及生物能源等方面。

在疾病诊断方面，由于机体的病理变化，代谢产物也产生了某种相应的变化，应用代谢组学分析对这些由疾病引起的代谢产物的响应进行分析，可以促进疾病生物标志物的发现以及辅助临床诊断。如代谢组学技术被应用于胃癌、肝癌、前列腺癌、肠癌和心血管疾病等重大疾病方面的研究。

代谢组学的数据处理主要包括原始数据的预处理和多变量数据分析两个步骤，其目的在于揭示生物样品的关键性生物标志物、代谢节点和代谢途径，最终帮助揭示代谢组学的数据分析结果。数据处理是挖掘过程中的重要一步，如果存在许多不相关、冗余的信息、噪声或不可靠的数据，则使得数据训练期间的特征识别更加困难。

LC/MS成为复杂生物样品中分析代谢产物的首选方法，但是基于LC/MS的代谢组学数据，都具有信息量大、维度高以及数据复杂性强等特点，如何处理以及更好地整合代谢组学数据仍然面临着巨大的挑战，比如临床数据收集时间长，分析样本量巨大；为了确定各种生物系统中的生物扰动、分析代谢变化的时间依懒性、以及评估个性化医疗中的治疗效果和高危人群分层，来自大规模代谢学研究中的数据通常收集期可达数月甚至数年之久，并且必须分成多个批次，需要综合考虑各实验批次研究中所得到的全部数据，另外，相对于小样本而言巨大数量的变量所造成的维度灾难，以及当数据分析过程中需要评估大量参数以及测试大量统计假设时，如果样本量过小在估计参数或者测试统计假设上带来的困难，也需要慎重考虑数据整合的重要性和必要性。

目前，通常使用开源程序处理LC/MS代谢组学数据，但是由于这些软件不能对仪器产生的原始数据进行识别，所以需要使用者准备符合软件格式要求的数据，也因此需要使用与仪器配套的软件对原始数据进行处理，直至成为可供分析使用的数据集。当前的一些程序对于整合不同样本和不同批次的代谢组学数据，存在一些问题和局限性，数据处理时间和难度会随样本数增加而增加，不适用于不断有新样本需要测试的代谢组学数据分析，结果的准确度和精度相对较低。

因此，本领域亟需开发一种准确度和精度均较高，且可用于处理不断有新样品增加，样品收集周期较长的代谢组学数据的分析方法。

发明内容

本发明所要解决的技术问题是克服现有技术中在整合不同样本和不同批次的代谢组学数据时，数据处理时间和难度会随样品增加而增加，不适用于不断有新样品需要进行处理的代谢组学数据分析，且现有方法的准确度和精度都较低等缺陷，而提供一种生物代谢组学数据处理方法、分析方法及装置和应用。本发明的生物代谢组学数据处理方法可用于分析不断有新样品增加，样品收集时间长的代谢组学数据，例如结直肠癌高危人群血清成分的代谢组学数据，通过采用液质谱联用（LC-MS）方式检测正常人群和高危人群（腺瘤、癌前病变、癌症）血清样本的代谢组学数据，比较和分析正常人群和患者血清中不同代谢产物的质谱峰，识别和确定差异代谢物的结构，为对差异代谢物的功能和代谢途径的确定提供研究基础，便于进一步判断出不同个体（患者）的病理、生理状态，为后期临床诊断和治疗提供依据。本发明的生物代谢组学数据处理方法的准确度和精度更高。

本发明采用以下技术方案解决上述技术问题：

本发明提供一种生物代谢组学数据处理方法，所述生物代谢组学数据包括待测样品的液相色谱-质谱数据，所述生物代谢组学数据处理方法包括筛选出多个所述待测样品中的差异代谢物的步骤，所述筛选出多个所述待测样品中的差异代谢物的步骤包括：

（1）对所述生物代谢组学数据进行预处理和标准化处理，得到二维矩阵数据；

（2）缺失峰补齐处理：对所述二维矩阵数据进行稀疏重构处理，再进行降维处理；

（3）对步骤（2）得到的数据进行化学计量分析，筛选出所述待测样品中的所述差异代谢物。

一些实施例中，所述待测样品包括正常人群待测样品、高危人群待测样品和炎症人群待测样品中任意一种或多种。

一些实施例中，所述待测样品的液相色谱-质谱数据的检测过程中还可包括校验的操作，所述校验包括如下步骤：所述检测的过程中插入对QC样品的测试，所述QC样品为所述待测样品的混合物，基于所述QC样品的液相色谱-质谱数据对所述待测样品的液相色谱-质谱数据进行系统误差校正。

较佳地，所述校验过程中，删除所述QC样品的液相色谱-质谱数据中变异系数大于30%的特征峰。通过删除重复性差的特征峰有助于提供高质量数据集。所述变异系数可为本领域技术人员常规认为的标准偏差与测量结果平均值的比值，即相对标准偏差。

采用所述QC样品模拟数据采集过程中信号的变化，将所述QC样品作为训练集，建立预测模型，预测所述待测样品的信号变化，从而对所述待测样品的信号进行校正，其还可用于评价仪器的稳定性和重复性以确保实验结果的可靠性。

本发明一较佳实施例中，所述QC样品的液相色谱-质谱数据中变异系数小于30%的特征峰的比例为大于70%。此时，证明检测仪器的稳定性和重复性佳，确保实验结果具有可靠性。

一些实施例中，所述预处理可按照本领域常规包括基线校正、峰识别、峰过滤和峰对齐中的任意一种或多种。

一些实施例中，所述预处理可在AMIX tool-kit软件、MetaXCMS软件、MZmine2软件或Mnova软件中进行。

一些实施例中，所述预处理的操作前还可进一步包括将所述液相色谱-质谱数据原始数据转化为mzXML格式数据的操作。

一些实施例中，所述标准化处理可在MetaXCMS软件中进行。

一些实施例中，所述标准化处理的方法可为本领域常规，一般可包括权重转换，较佳地，包括如下步骤：对所述预处理后的生物代谢组学数据进行中心化，再除以列变量标准差。通过所述标准化处理让数据无量纲化，使不同性质的变量具有可比性，将不同数量级的变量数据经过不同的转换至合适范围，避免大值变量掩盖小值变量的波动。

一些实施例中，所述稀疏重构处理的方法包括二维离散余弦变换法。

其中，采用所述二维离散余弦变换法进行所述稀疏重构处理时，具体包括如下步骤：设步骤（1）得到的所述二维矩阵数据为

，其中，n为代谢组学特征数；形成二维观测矩阵y=Ax+N，其中，

，A为采用所述二维离散余弦变换进行稀疏构成的字典矩阵，x为稀疏矩阵，N为加性高斯白噪声，m为所述二维观测矩阵的维度，即输入信号x经过所述二维离散余弦变换得到的维度数。

一些实施例中，所述降维处理的方法包括压缩感知算法。

其中，采用所述压缩感知算法进行所述降维处理时，具体包括如下步骤：采用正交匹配追踪算法计算所述稀疏重构处理得到的数据的频率f_t和残差R_t；

，

，

其中，

为经过t次迭代的索引积核选出的所述字典矩阵A的列集合，f_t为t次迭代后的频率，R_t为t次迭代后的残差，上标H为所述字典矩阵的转置操作；

根据f_t，采用二维离散反余弦变换恢复压缩后的数据，得数据X’，其中，

，且m<<n，取m=20；

判断迭代完成后R_t是否满足允许的最大训练残差值，如不满足则迭代结束，否则重复步骤a、步骤b和步骤c，直到满足要求。

步骤（2）中，对所述二维矩阵数据进行稀疏重构处理，再进行降维处理，促使得到降维数据，能从压缩后的数据无损的恢复原始数据，方便后续的数据分析，有效提高了差异代谢物筛选的准确度和精度。

一些实施例中，所述化学计量分析包括核主成分分析KPCA。

一些实施例中，所述KPCA可按照本领域常规包括如下步骤：

a1：利用高斯径向基核函数计算核矩阵K：

b1：使所述核矩阵K进行中心化处理，得修正矩阵K _c；

c1：利用jacobi迭代计算K_c的特征值

和特征向量

；

d1：按所述特征值从大到小排序，得

，并对所述特征向量进行相应排序，取前L个所述特征值和前L个所述特征向量；

e1：通过施密特正交化方法处理步骤（d1）得到的前L个所述特征值和前L个所述特征向量，得单位正交化特征向量

；

f1：计算降维后的数据：

，其中，Y为L维的降维处理结果。

一较佳实施例中，

，其中，

为经步骤（2）处理后得到的数据，

，

为方差。

一较佳实施例中，所述K _c：

，其中，

为m行m列，主对角线上值为1/m，其他位置为0的矩阵。

一些实施例中，所述KPCA通过映射函数把步骤（2）得到的数据非线性映射到特征空间，在特征空间上进行PCA分析，对各组数据进行归类，抽取更多数据特征信息，去除重复性差的待测样品和异常待测样品，增强了对代谢组学数据处理能力。在生物代谢组学数据处理方法中使用所述KPCA后的分类效果明显好于PCA。采用所述KPCA可挖掘到数据集中蕴含的非线性关系，PCA无法达到该目的。

一些实施例中，所述KPCA的操作后还可进一步包括对所述KPCA得到的KPCA模型进行K折交叉验证的操作。

一些实施例中，所述差异代谢物的筛选条件包括VIP≥1；较佳地，所述差异代谢物的筛选条件包括：显著性差异p-value≤0.05，VIP≥1，变化倍数Fold change≥1.5或变化倍数≤0.667，单变量方差分析p-value≤0.05，两因素方差分析p-value≤0.05。

本发明还提供一种生物代谢组学数据的分析方法，所述方法包括：获取差异代谢物，所述差异代谢物根据如上所述的生物代谢组学数据处理方法处理后得到；采用正交偏最小二乘判断分析法OPLS-DA对所述差异代谢物进行分析。

一些实施例中，所述OPLS-DA的条件和方法可为本领域常规，可建立代谢物样品分组之间的关系模型，即OPLS-DA模型。所述OPLS-DA使用正交信号的校正技术，将X矩阵分解成与Y矩阵相关和不相关的两类信息，然后过滤掉所述不相关的信息，所述相关的信息用于筛选贡献较大的所述差异代谢物。

一些实施例中，所述OPLS-DA的操作后还可进一步包括对所述OPLS-DA得到的OPLS-DA模型进行K折交叉验证和/或根据MS/MS碎片，采用数据库确认所述差异代谢物的结构的操作。

其中，所述MS/MS碎片分子量的误差可按照本领域常规小于30ppm。

其中，所述数据库可为本领域常规使用的数据库，较佳地为Metlin数据库、MoNA数据库和KEGG数据库中的任意一种或多种。

本发明还提供一种如上所述的生物代谢组学数据处理方法或如上所述的生物代谢组学数据的分析方法在确定不同待测样品中差异代谢物的应用。

一些实施例中，所述应用为用于确定高危人群血清中异常代谢物的应用，较佳地为用于确定直肠癌高危人群血清中异常代谢物的应用。

本发明还提供一种生物代谢组学数据处理装置，所述生物代谢组学数据包括待测样品的液相色谱-质谱数据，所述生物代谢组学数据处理装置包括：

二维矩阵数据生成模块，用于对所述生物代谢组学数据进行预处理和标准化处理，得到二维矩阵数据；

缺失峰补齐模块，用于对所述二维矩阵数据进行稀疏重构处理，再进行降维处理；

差异代谢物筛选模块：用于对所述缺失峰补齐模块处理后的数据进行化学计量分析，筛选出不同待测样品之间差异代谢物。

一些实施例中，所述缺失峰补齐模块包括二维离散余弦变换单元和压缩感知计算单元。

一些实施例中，所述差异代谢物筛选模块包括KPCA处理单元。

较佳地，所述生物代谢组学数据处理装置还可进一步包括对所述KPCA得到的KPCA模型进行K折交叉验证模块。

本发明还提供一种生物代谢组学数据的分析装置，其包括数据获取模块和正交偏最小二乘判断分析模块；所述数据获取模块用于获取差异代谢物，所述差异代谢物根据如上所述的生物代谢组学数据处理装置处理后得到。

一些实施例中，所述生物代谢组学数据的分析装置中还可进一步包括对所述正交偏最小二乘判断分析模块得到的OPLS-DA模型进行K折交叉验证模块和/或确定所述差异代谢物结构模块。

在符合本领域常识的基础上，上述各优选条件，可任意组合，即得本发明各较佳实例。

本发明所用试剂和原料均市售可得。

本发明的积极进步效果在于：本发明的生物代谢组学数据处理方法可用于分析不断有新样品增加，样品收集时间长的代谢组学数据，例如结直肠癌高危人群血清成分的代谢组学数据，通过采用液质谱联用（LC-MS）方式检测正常人群和高危人群（腺瘤、癌前病变、癌症）血清样本的代谢组学数据，比较和分析正常人群和高危人群血清中不同代谢产物的质谱峰，识别和确定差异代谢物（即异常代谢物）的结构，结果的准确度和精确度更高，为对差异代谢物的功能和代谢途径的确定提供研究基础，便于进一步判断出不同个体（患者）的病理、生理状态，为后期临床诊断和治疗提供依据。

附图说明

本公开可以通过参考下文中结合附图所给出的描述而得到更好的理解。所述附图连同下面的详细说明一起包含在本说明书中并且形成本说明书的一部分，而且用来进一步举例说明本公开的优选实施例和解释本公开的原理和优点。其中：

图1为本发明一实施例提供的生物代谢组学数据处理方法的流程图；

图2为本发明另一实施例提供的生物代谢组学数据处理方法的流程图

图3为本发明一实施例提供的生物代谢组学数据的分析方法的流程图；

图4为本发明另一实施例提供的生物代谢组学数据的分析方法的流程图；

图5为本发明一实施例中QC样品制备模拟图；

图6为本发明一实施例中测试液相色谱-质谱数据时的样品测试顺序图；

图7为本发明一实施例中QC样品的不同变异系数的特征峰所占比例的结果图；

图8为本发明一实施例中采用MetaXCMS软件进行峰识别时，示例样品提取离子基峰的液相色谱图；

图9为本发明一实施例中上调代谢物和下调代谢物个数分布图。

具体实施方式

下面通过实施例的方式进一步说明本发明，但并不因此将本发明限制在所述的实施例范围之中。下列实施例中未注明具体条件的实验方法，按照常规方法和条件，或按照商品说明书选择。

下述实施例中所使用的实验方法如无特殊说明，均为常规方法。

本发明实施方式的第一方面提供一种确定结直肠癌高危人群血清中异常代谢物的生物代谢组学数据处理方法，如图1；

其中，生物代谢组学数据处理方法，所述生物代谢组学数据包括待测样品的液相色谱-质谱数据，所述生物代谢组学数据处理方法包括筛选出多个所述待测样品中的差异代谢物的步骤，所述筛选出多个所述待测样品中的差异代谢物的步骤包括：

步骤S110：对所述生物代谢组学数据进行预处理和标准化处理，得到可用于统计分析的二维矩阵数据；

步骤S120：缺失峰补齐处理：对所述二维矩阵数据进行稀疏重构处理，再进行降维处理；

步骤S130：对步骤S120得到的数据进行化学计量分析，筛选出所述待测样品中的所述差异代谢物。

在本实施方式的一个实施例中，所述待测样品包括119份样品，其中，健康人群49人，腺瘤人群（高危人群）45人，其余为炎症人群。本实施例中的样本是从2017年下半年至2019年采集，一般采集后按规定条件保存在医院的样品柜中，待累计到一定数量再进行测试。

在本实施方式的一个实施例中，生物代谢组学数据的测试方法包括如下步骤：

A：待测样品的制备，将所有待测样本在4℃下融化，样品量不足按照等比例缩减；自每个样本中取100µL于2mL离心管中；每个离心管加入400µL温度为-20℃的甲醇，振荡60s，充分混匀；在12000rpm，4℃条件下离心10 min，取全部上清液，转移至新的2mL离心管中，真空浓缩干燥；加入150µL的含2-氯苯丙氨酸的甲醇水溶液复溶，2-氯苯丙氨酸的浓度为4ppm，甲醇水溶液中甲醇的体积百分比为80%，再使用0.22µm滤膜过滤，得到待测样品；

B：QC样品的制备，取上述方法制得的待测样品，每个待测样品取20µL混合成QC样品，制备模拟图见图5，QC：quality control；

C：采用液相质谱联用设备分别测试上述制得的各待测样品和QC样品的液相色谱-质谱数据，每10个待测样品中插入一个QC样品，具体测试顺序见图6；

在本实施方式的一个实施例中，生物代谢组学数据的测试条件如下：

液相色谱在超高效液相色谱仪中进行测定，测试条件为：采用ACQUITY UPLC®HSS T3 1.8 µm（2.1×150 mm）色谱柱，自动进样器温度为8℃，流动相的流速为0.25mL/min，柱温为40℃，进样2μL后进行梯度洗脱；正离子模式下，流动相0.1%甲酸水（B₂）-0.1%甲酸乙腈（A₂）；正离子模式下，梯度洗脱程序为0~1min，2% A₂+98%B₂；1~9min，2% A₂+98%B₂→50% A₂+50%B₂；9~12min，50% A₂+50%B₂→98%A₂+2%B₂；12~13.5min，98% A₂+2% B₂；13.5~14min，98% A₂+2%B₂→2% A₂+98%B₂；14~20min，2% A2+98% B₂；

负离子模式下，流动相为5mM甲酸铵水（B₁）-乙腈（A₁）；梯度洗脱程序为0~1min，2%A₁+98%B₁；1~9min，2% A₁+98%B₁→50% A1+50%B₁；9~12min，50% A₁+50%B₁→98% A₁+2% B₁；12~13.5min，98%A₁+2% B₁；13.5~14min，98% A₁+2% B₁→2%A₁+98% B₁； 14~17min，2% A₁+98% B₁，百分比为体积百分比；

质谱在质谱仪中进行测试，测试条件为：采用电喷雾离子源（ESI），正负离子电离模式，正离子喷雾电压为3.50kV，负离子喷雾电压为2.50kV，鞘气流速为30arb，辅助气流速为10arb，毛细管的温度为325℃，以分辨率70000进行全扫描，扫描范围为81~1000，并采用HCD进行二级裂解，碰撞电压为30eV，同时采用动态排除去无必要的MS/MS信息；

本实施方式的另一个实施例中，见图2，所述待测样品的液相色谱-质谱数据的检测过程中还包括校验的操作S140，所述校验包括如下步骤：所述检测的过程中插入对QC样品的测试，所述QC样品为所述待测样品的混合物，基于所述QC样品的液相色谱-质谱数据对所述待测样品的液相色谱-质谱数据进行系统误差校正。

在所述校验过程中，删除所述QC样品的液相色谱-质谱数据中变异系数大于30%的特征峰。通过删除重复性差的特征峰有助于提供高质量数据集。所述变异系数为本领域技术人员常规认为的标准偏差与测量结果平均值的比值，即相对标准偏差。

采用所述QC样品模拟数据采集过程中信号的变化，将所述QC样品作为训练集，建立预测模型，预测所述待测样品的信号变化，从而对所述待测样品的信号进行校正，其还可用于评价仪器的稳定性和重复性以确保实验结果的可靠性。

本实施方式的一个实施例中，所述QC样品的液相色谱-质谱数据中变异系数小于30%的特征峰的比例为大于70%，如图7所示。可见，本实施例中，检测仪器的稳定性和重复性佳，确保实验结果具有可靠性。

本实施方式的一个实施例中，所述预处理包括基线校正、峰识别、峰过滤和峰对齐。

本实施方式的一个实施例中，所述预处理在MetaXCMS软件中进行，其他实施例中，所述预处理还可在AMIX tool-kit软件、MZmine2软件或Mnova软件中进行。

本实施方式的一个实施例中，所述预处理时，所述MetaXCMS软件的主要参数包括：bw=5，ppm=15，peakwidth=c(5,30)，mzwid=0.015，mzdiff=0.01。在其他实施例中，所述MetaXCMS软件的主要参数可包括：bw=4~6，ppm=13~17，peakwidth=c(5,30)，mzwid=0.013~0.017，mzdiff=0.005~0.015，可根据处理LC-MS数据时的具体情况，通过尝试选择最优参数，对数据峰进峰识别、峰过滤、峰对齐等操作，得到最优数据。

峰识别：将生物代谢组学数据分割成质量单位切片，根据每个切片的时间点对应的最大峰强度确定峰信号，提取离子基峰色谱图，见图8；

检测峰值之前，使用二阶导数高斯模型峰对每个切片匹配过滤，过滤后，使用信噪比截止值选择峰，信噪比截止值取未经过二阶导数高斯模型峰过滤的数据平均值，其中信噪比截止值为10最优；

峰宽大于切片宽度，峰由多个切片组成，切片易丢失，按峰强度罗列，消除临近次高强度峰；峰宽小于切片宽度，相邻切片最大强度峰组合成重叠的EIBPC（即100.0/100.1、100.1/100.2等），产生平滑和锯齿状色谱图，用于过滤和峰检测；在邻近消除处理步骤，将锯齿状轮廓（从不完整信号中积分）中检测到的峰由从平滑轮廓（从完整信号中积分）检测到的峰替代；

峰过滤：确定单个样品峰后，在样品之间进行峰匹配，计算保留时间偏差；使用峰匹配算法，采用0.25m/z宽的固定间隔区间来匹配质量域中的峰；为了避免由于任意切片边界而将组分开，本实施例使用重叠切片，其中相邻切片重叠一半（比如：100.0-100.25、100.125-100.375）；在分割切片时，每个峰在两个重叠的切片中被计数两次；与峰值拾取类似，后处理步骤用于去除源自重叠切片的峰值组；

峰切片后，每个切片中具有不同保留时间的峰组，也可以采用固定间隔对切片峰进行匹配，再计算色谱时间内峰的整体分布，并动态识别许多峰具有相似保留时间的区域的边界；计算峰分布方法是核密度估计器，从该分布中识别出所谓的主峰，它们代表了许多具有相似保留时间的峰；从分布中的最高峰开始，下降到主峰的任一侧，直到分布再次增加；该过程设置一个固定间隔，在该间隔中所有峰都被放入一个组，对分布中的所有主峰重复该过程，同时调节峰匹配兼容性；

去除不重要的峰组，考虑峰重现性，从少于一半的样本中去除峰重复出现差的峰组，若已知样本根据不同的条件分组，可以从少于一半样本组中去除重复出现差的峰组；

峰对齐：校正所有样品的保留时间，分组后，自动识别出数百个“表现良好”的峰组，将峰值匹配到“表现良好”的峰组作为临时标准；同时对于每组，计算该组中每个样品的中值保留时间和与中值的偏差；性能良好的峰组均匀分布在色谱图的重要部分，为每个样品构建详细的非线性保留时间偏差轮廓；

保留时间偏差，使用局部回归拟合方法loess拟合数据，采用分段拟合方式，纠正保留时间的局部扰动，loess拟合方法从数据中自动去除残余异常值峰；使用产生的偏差曲线来校正原始峰列表的保留时间，将它们再次匹配到组中；匹配/比对过程以迭代方式重复，连续识别越来越多的表现良好的峰组，以实现越来越精确的比对；

本实施方式的一个实施例中，所述预处理的操作前还进一步包括将液相色谱-质谱数据原始数据转化为mzXML格式数据的操作。

本实施方式的一个实施例中，所述标准化处理在MetaXCMS软件中进行。

本实施方式的一个实施例中，所述标准化处理的方法包括权重转换，具体包括如下步骤：对所述预处理后的生物代谢组学数据进行中心化，再除以列变量标准差。通过所述标准化处理让数据无量纲化，使不同性质的变量具有可比性，将不同数量级的变量数据经过不同的转换至合适范围，避免大值变量掩盖小值变量的波动。

本实施方式的一个实施例中，所述稀疏重构处理的方法包括二维离散余弦变换法。

采用所述二维离散余弦变换法进行所述稀疏重构处理时，具体包括如下步骤：设步骤（1）得到的所述二维矩阵数据为

，其中，n为代谢组学特征数；形成二维观测矩阵y=Ax+N，其中，

，A为采用所述二维离散余弦变换进行稀疏构成的字典矩阵，x为稀疏矩阵，N为加性高斯白噪声，m为所述二维观测矩阵的维度，即输入信号x经过所述二维离散余弦变换得到的维度数。

本实施方式的一个实施例中，所述降维处理的方法包括压缩感知算法。

其中，采用所述压缩感知算法进行所述降维处理时，具体包括如下步骤：采用正交匹配追踪算法计算所述稀疏重构处理得到的数据的频率f_t和残差R_t；

，

，

其中，

为经过t次迭代的索引积核选出的所述字典矩阵A的列集合，f_t为t次迭代后的频率，R_t为t次迭代后的残差，上标H为所述字典矩阵的转置操作；

根据f_t，采用二维离散反余弦变换恢复压缩后的数据，得数据X’，其中，

，且m<<n，取m=20；

判断迭代完成后R_t是否满足允许的最大训练残差值，如不满足则迭代结束，否则重复步骤a、步骤b和步骤c，直到满足要求。

本实施方式的一个实施例中，对所述二维矩阵数据进行稀疏重构处理，再进行降维处理，促使得到降维数据，能从压缩后的数据无损的恢复原始数据，方便后续的数据分析，有效提高了差异代谢物筛选的准确度和精度。

本实施方式的一个实施例中，所述差异代谢物的筛选条件包括：p-value≤0.05，VIP≥1，变化倍数（Fold change）≥1.5或变化倍数≤0.667，单变量方差分析p-value≤0.05，两因素方差分析p-value≤0.05。

本实施方式的另一个实施例中，所述化学计量分析包括核主成分分析（KPCA），见图2中S130，采用所述KPCA挖掘到数据集中蕴含的非线性关系，PCA无法达到该目的；所述KPCA包括如下步骤：

a1：利用高斯径向基核函数计算核矩阵K：

，其中，

为经步骤（2）处理后得到的数据，

，

为方差；

b1：使所述核矩阵K进行中心化处理，得修正矩阵K _c：

，其中，

为m行m列，主对角线上值为1/m，其他位置为0的矩阵；

c1：利用jacobi迭代计算K_c的特征值

和特征向量

；

d1：按所述特征值从大到小排序，得

，并对所述特征向量进行相应排序，取前L个所述特征值和前L个所述特征向量；

e1：通过施密特正交化方法处理步骤（d1）得到的前L个所述特征值和前L个所述特征向量，得单位正交化特征向量

；

f1：计算降维后的数据：

，其中，Y为L维的降维处理结果。

本实施方式的另一个实施例中，所述KPCA通过映射函数把步骤（2）得到的数据非线性映射到特征空间，在特征空间上进行PCA分析，对各组数据进行归类，抽取更多数据特征信息，去除重复性差的待测样品和异常待测样品，增强了对代谢组学数据处理能力。在生物代谢组学数据处理方法中使用所述KPCA后的分类效果明显好于PCA。

本实施方式的一个实施例中，对KPCA模型进行K折交叉验证时，参考的参数包括R2X，R2X为KPCA模型的可释放度，结果见表1，可见R2高于0.4，证明拟合度理想，KPCA模型的有效性佳；

表1

	pre	R2X(cum)
			AD&HC	10	0.434

表1中，pre为主成分数；R2X为模型（对X变量数据集）可解释度；Q2为模型可预测度；AD为腺癌人群；HC为高危人群。

经过上述一系列处理过程，得到包括质核比（m/z）和保留时间及峰面积等信息的数据矩阵，正离子模式获得10977个前体分子，负离子模式获得10551个前体分子，导出数据Excel进行后续分析。

本发明实施方式的第二方面提供了一种生物代谢组学数据分析方法，其包括获取差异代谢物，所述差异代谢物根据如上所述的生物代谢组学数据处理方法处理后得到；采用正交偏最小二乘判断分析法OPLS-DA对所述差异代谢物进行分析，见图3。

本实施方式的一个实施例中，采用所述OPLS-DA建立代谢物样品分组之间的关系模型，即OPLS-DA模型。所述OPLS-DA使用正交信号的校正技术，将X矩阵分解成与Y矩阵相关和不相关的两类信息，然后过滤掉所述不相关的信息，所述相关的信息用于筛选贡献较大的所述差异代谢物。

本实施方式的另一个实施例中，所述OPLS-DA的操作后还进一步包括对所述OPLS-DA得到的OPLS-DA模型进行K折交叉验证和根据MS/MS碎片，采用Metlin数据库确认所述差异代谢物的结构的操作，见图4中S220。MS/MS碎片分子量的误差小于30ppm。其他实施例中可使用MoNA数据库或KEGG数据库。

对OPLS-DA模型进行K折交叉验证结果见表2；R2高于0.4，证明拟合度理想，OPLS-DA模型的有效性佳。

表2

	pre	R2X(cum)	R2Y(cum)	Q2(cum)
					AD&HC	1+3+0	0.583	0.983	0.917

表2中，pre为主成分数；R2X为模型（对X变量数据集）可解释度；R2Y为模型（对Y变量数据集）可解释度；Q2为模型可预测度；AD为腺癌人群；HC为高危人群。

结果发现，异常代谢物有2300多种，其中，上调代谢物有1456种，下调代谢物有888种，结果见图9，其中确定结果的异常代谢物包括：（2S）-甘草素、（2S,3S）-2-羟基十三烷-1,2,3-三羧酸盐、（2Z，4S，5R）-2-氨基-4,5,6-三羟基己基-2-烯酸、（3R）-β-亮氨酸，各个异常代谢物性质见表3。

表3

	m/z(质核比)	rt（保留时间）	分子量	分子式
					（2S）-甘草素	255.07055409	719.4865	256.0736	C15H12O4
（2S,3S）-2-羟基十三烷-1,2,3-三羧酸盐	331.1761844	599.827	332.1835	C16H28O7
					（2Z，4S，5R）-2-氨基-4,5,6-三羟基己基-2-烯酸	178.0863185	345.952	177.0637	C6H11NO5
（3R）-β-亮氨酸	130.0851783	515.5005	131.0946	C6H13NO2

最后，还需要说明的是，在本发明中术语“包括”、“包含”或者其任何其他变体意在涵盖非排他性的包含，从而使得包括一系列要素的过程、方法、物品或者设备不仅包括那些要素，而且还包括没有明确列出的其他要素，或者是还包括为这种过程、方法、物品或者设备所固有的要素。

尽管上面已经通过本公开的具体实施例的描述对本公开进行了披露，但是，应该理解，本领域技术人员可在所附方案的精神和范围内设计对本公开的各种修改、改进或者等同物。这些修改、改进或者等同物也应当被认为包括在本公开所要求保护的范围内。

Claims (10)

1.一种生物代谢组学数据处理方法，其特征在于，所述生物代谢组学数据包括待测样品的液相色谱-质谱数据，所述生物代谢组学数据处理方法包括筛选出多个所述待测样品中的差异代谢物的步骤，所述筛选出多个所述待测样品中的差异代谢物的步骤包括：

（1）对所述生物代谢组学数据进行预处理和标准化处理，得到二维矩阵数据；

（2）缺失峰补齐处理：对所述二维矩阵数据进行稀疏重构处理，再进行降维处理；

（3）对步骤（2）得到的数据进行化学计量分析，筛选出所述待测样品中的所述差异代谢物。

2.如权利要求1所述的生物代谢组学数据处理方法，其特征在于，所述待测样品包括正常人群待测样品、高危人群待测样品和炎症人群待测样品中的任意一种或多种；

和/或，所述待测样品的液相色谱-质谱数据的检测过程中还包括校验的操作，所述校验包括如下步骤：所述检测的过程中插入对QC样品的测试，所述QC样品为所述待测样品的混合物，基于所述QC样品的液相色谱-质谱数据对所述待测样品的液相色谱-质谱数据进行系统误差校正；

较佳地，所述校验过程中，删除所述QC样品的液相色谱-质谱数据中变异系数大于30%的特征；

较佳地，所述QC样品的液相色谱-质谱数据中变异系数小于30%的特征峰的比例为大于70%。

3.如权利要求1所述的生物代谢组学数据处理方法，其特征在于，所述预处理包括基线校正、峰识别、峰过滤和峰对齐中的任意一种或多种；

和/或，所述预处理在AMIX tool-kit软件、MetaXCMS软件、MZmine2软件或Mnova软件中进行；

和/或，所述预处理的操作前还进一步包括将所述液相色谱-质谱数据原始数据转化为mzXML格式数据的操作；

和/或，所述标准化处理在MetaXCMS软件中进行；

和/或，所述标准化处理的方法包括权重转换，较佳地，所述标准化处理包括如下步骤：对所述预处理后的生物代谢组学数据进行中心化，再除以列变量标准差；

和/或，所述稀疏重构处理的方法包括二维离散余弦变换法；

和/或，所述降维处理的方法包括压缩感知算法。

4.如权利要求3所述的生物代谢组学数据处理方法，其特征在于，采用所述二维离散余弦变换法进行所述稀疏重构处理时，具体包括如下步骤：设步骤（1）得到的所述二维矩阵数据为

，其中，n为代谢组学特征数；形成二维观测矩阵y=Ax+N，其中，

，A为采用所述二维离散余弦变换进行稀疏构成的字典矩阵，x为稀疏矩阵，N为加性高斯白噪声，m为所述二维观测矩阵的维度，即输入信号x经过所述二维离散余弦变换得到的维度数；

和/或，采用所述压缩感知算法进行所述降维处理时，具体包括如下步骤：采用正交匹配追踪算法计算所述稀疏重构处理得到的数据的频率f_t和残差R_t；

，

，

其中，

为经过t次迭代的索引积核选出的所述字典矩阵A的列集合，f_t为t次迭代后的频率，R_t为t次迭代后的残差，上标H为所述字典矩阵的转置操作；

根据f_t，采用二维离散反余弦变换恢复压缩后的数据，得数据X’，其中，

，且m<<n，取m=20；

判断迭代完成后R_t是否满足允许的最大训练残差值，如不满足则迭代结束，否则重复步骤a、步骤b和步骤c，直到满足要求。

5.如权利要求1所述的生物代谢组学数据处理方法，其特征在于，所述化学计量分析包括核主成分分析KPCA；

较佳地，所述KPCA的操作后包括对所述KPCA得到的KPCA模型进行K折交叉验证的操作；

较佳地，所述KPCA包括如下步骤：

a1：利用高斯径向基核函数计算核矩阵K；

b1：使所述核矩阵K进行中心化处理，得修正矩阵K _c；

c1：利用jacobi迭代计算K_c的特征值

和特征向量

；

d1：按所述特征值从大到小排序，得

，并对所述特征向量进行相应排序，取前L个所述特征值和前L个所述特征向量；

e1：通过施密特正交化方法处理步骤（d1）得到的前L个所述特征值和前L个所述特征向量，得单位正交化特征向量

；

f1：计算降维后的数据：

，其中，Y为L维的降维处理结果；

更佳地，所述

，其中，

为经步骤（2）处理后得到的数据，

，

为方差；

更佳地，所述

，其中，

为m行m列，主对角线上值为1/m，其他位置为0的矩阵；

和/或，所述差异代谢物的筛选条件包括VIP≥1；较佳地，所述差异代谢物的筛选条件包括：显著性差异p-value≤0.05，VIP≥1，变化倍数Fold change≥1.5或变化倍数≤0.667，单变量方差分析p-value≤0.05，两因素方差分析p-value≤0.05。

6.一种生物代谢组学数据的分析方法，其特征在于，所述方法包括：获取差异代谢物，所述差异代谢物根据权利要求1-5中任一项所述的方法处理后得到；采用正交偏最小二乘判断分析法OPLS-DA对所述差异代谢物进行分析。

7.如权利要求6所述的生物代谢组学数据的分析方法，其特征在于，所述OPLS-DA建立代谢物样品分组之间的关系模型，得OPLS-DA模型；

和/或，所述OPLS-DA的操作后还进一步包括对所述OPLS-DA得到的OPLS-DA模型进行K折交叉验证和/或根据MS/MS碎片，采用数据库确认所述差异代谢物的结构的操作；较佳地，所述MS/MS碎片分子量的误差小于30ppm；较佳地，所述数据库为Metlin数据库、MoNA数据库和KEGG数据库中的任意一种或多种。

8.一种如权利要求1~5中任意一项所述的生物代谢组学数据处理方法，或如权利要求6或7所述的生物代谢组学数据的分析方法在确定不同待测样品中差异代谢物的应用；

较佳地，所述应用为用于确定高危人群血清中异常代谢物的应用，更佳地为用于确定直肠癌高危人群血清中异常代谢物的应用。

9.一种生物代谢组学数据处理装置，其特征在于，所述生物代谢组学数据包括待测样品的液相色谱-质谱数据，所述生物代谢组学数据处理装置包括：

二维矩阵数据生成模块，用于对所述生物代谢组学数据进行预处理和标准化处理，得到二维矩阵数据；

缺失峰补齐模块，用于对所述二维矩阵数据进行稀疏重构处理，再进行降维处理；

差异代谢物筛选模块：用于对所述缺失峰补齐模块处理后的数据进行化学计量分析，筛选出不同待测样品之间差异代谢物；

较佳地，所述缺失峰补齐模块包括二维离散余弦变换单元和压缩感知计算单元；

较佳地，所述差异代谢物筛选模块包括KPCA处理单元；较佳地，所述生物代谢组学数据处理装置还进一步包括对所述KPCA得到的KPCA模型进行K折交叉验证模块。

10.一种生物代谢组学数据的分析装置，其特征在于，其包括数据获取模块和正交偏最小二乘判断分析模块；所述数据获取模块用于获取差异代谢物，所述差异代谢物根据权利要求9所述的装置处理后得到；

较佳地，所述生物代谢组学数据的分析装置中还进一步包括对所述正交偏最小二乘判断分析模块得到的OPLS-DA模型进行K折交叉验证模块和/或确定所述差异代谢物结构模块。

Images