本发明涉及固态发酵技术领域，尤其涉及固态发酵过程解析与调控，具体涉及一种基于组学整合技术确定固态发酵核心功能模块的方法及系统。
背景技术：
固态发酵(ssf)是指一类使用不溶性固体基质培养微生物的工艺过程，既包括将固体悬浮在液体中的深层发酵，也包括在没有游离水或游离水很少的培养基上培养微生物的过程。ssf具有数千年悠久的应用历史，广泛应用于面包、酒曲、制酱、生产豆豉和奶酪等多种传统风味发酵食品。近年来，因其具有节水、节能、高得率、清洁等优势，逐渐受到各国的重视。
经过近几十年的迅速发展，我国已成为发酵工业大国，发酵工业占全国工业总产值的比例逐年提高。但由于人们对固态发酵本质认识不够深入，以及相关的发酵设备及其配套工艺研发滞后，造成了固态发酵在大规模应用中仍存在诸多问题。因此对于固态发酵工艺基本原理的探究变得极为重要。
由于菌落在酿造过程中的发挥主体作用，对于固态发酵工艺基本原理的研究，需要探明固态发酵过程中微生物群落结构和功能演替规律，阐明微生物代谢活动和产品成分、核心代谢模块之间的关系，即“合成产品成分的核心微生物及其代谢路径”所组成的核心功能模块。近些年，尤其是宏基因组技术的发展对固态发酵的研究产生了深远的影响。
高通量的宏基因组数据使得人们能够系统地分析微生物的代谢、微生物群落的相互作用及其对环境的反应机制。高通量的宏基因组数据的解析，一般包括：第一步—可视化多组学数据；第二步—挖掘数据。
第一步中，目前全基因尺度可视化方法有两类，第一类以反应物为端点，以反应为边生成代谢网络，包括cytoscape、3omics、genetrailexpress、david、webgestalt、paint4net。第一类只适用于反应数量较少的网络，展示全基因尺度网络时会重叠交织，无法辨识，解读较为困难。第二类是以kegg网络数据代谢图为底板，通过颜色渲染进行数据分析。如pathview具有基因组尺度可视化的功能，但功能比较单一，目前只能对全基因组尺度代谢途径中的代谢物进行颜色修订，无法实现多组学数据的有效整合。
第二步中，通过组学数据与目标研究对象的关联性进行确定核心模块。现在一般是基于统计学方法将组学数据与目标研究对象的关联性进行确定核心模块。
cn201310756418.0公开了一种刺糖多孢菌基因组尺度代谢网络模型及构建方法及应用。构建方法包括：根据kegg及ncbi数据库中刺糖多孢菌基因组序列的注释信息，添加多杀菌素生物合成的特征反应和菌体合成反应，并对网络反应进行手动精炼，获得刺糖多孢菌基因组尺度代谢网络模型。该申请利用刺糖多孢菌基因组尺度代谢网络模型，预测潜在基因靶点对多杀菌素产量提高的影响，最终确定改造方向，实现菌株的途径分子改造方法。
cn201210066912.x公开了一种氧化葡萄糖酸杆菌的基因组尺度代谢网络构建和分析方法。构建方法包括：利用kegg数据库中氧化葡萄糖酸杆菌的基因注释信息和对应酶的生物化学信息构建反应列表草图；对草图进行修改；添加生物量合成、运输和交换反应形成代谢网络；将代谢网络转化为sbml(系统生物学标记语言)格式，利用matlab和cobratoolbox对代谢网络进行调试，分析存在的缺口和无效循环，并根据调试结果进行修正；根据代谢网络，利用pajek绘制代谢网络图，并进行网络拓扑结构分析；根据代谢网络，利用matlab和cobratoolbox进行鲁棒性、基因必要性和通量可变性分析中的一种或多种分析。
cn201410814025.5公开了一种基于基因组尺度代谢网络模型指导胞内还原力调节的异丁醇合成菌株构建方法。构建方法包括：基于基因组尺度代谢网络模型，采用流量平衡分析和代谢最小调节分析，模拟胞内还原力代谢不同改造方式对菌株生长和异丁醇合成的作用规律，根据表型系数得出甘油醛-3-磷酸脱氢酶是异丁醇合成菌株胞内还原力调节的关键靶点。
cn201611072942.6公开了一种出芽短梗霉基因组尺度代谢网络模型及其应用。构建方法包括网络数据库的建立、数据库的精炼和校准、数学模型的建立和模型的验证与分析四个步骤；构建获得的出芽短梗霉基因组尺度代谢网络模型能够对潜在提高聚苹果酸的基因靶点进行预测。
技术实现要素：
本发明提供一种基于组学整合技术调控固态发酵的方法，实现固态发酵核心功能模块的准确挖掘，为解析固态发酵过程、构建理性调控策略提供了新方法。
本发明提供一种基于组学整合技术调控固态发酵的方法，包括以下步骤：
对固态发酵过程取样，对样品进行组学检测；
将检测所得组学数据通过kegg数据库比对，获得整个发酵过程的生化信息；
根据所得生化信息在kegg数据库中所属的代谢模块进行组装，构建全基因尺度代谢网络；
在所得全基因尺度代谢网络中，采用不同的标识表示基因组、化合物、酶活的含量水平，以固态发酵过程进行可视化；
根据代谢网络属性划分所得全基因尺度代谢网络；
基于所得划分结果，对组学数据与目标研究对象依据统计学方法进行关联性分析；
由所得分析结果确定固态发酵的核心功能模块。
其中，所述组学数据包括基因组、蛋白组和代谢组的数据。
进一步地，获得的生化信息包括微生物种类、基因、rna、酶、化合物及其对应关系。
进一步地，所述代谢网络属性包括拓扑学结构属性和网络统计学属性，所述拓扑学结构属性包括节点度、网络密度、网络直径、网络的连通度；所述网络统计学属性包括网络特征向量中心度和模块化划分。
进一步地，网络直径＜50、节点度>5、网络密度>0.4、网络的连通度>3、特征向量中心度>0.3、模块化划分为qmodularity、edge-betweenness、leadingeigenvector、fastgreedy、multilevel中的一种。
进一步地，所述统计学方法包括主因子成分分析、聚类分析、冗余度分析、偏最小二乘法回归分析、偏最小二乘判别分析和正交偏最小二乘法回归分析；
采用3种或3种以上方法对组学数据与目标研究对象进行关联性分析；
由所得3组或3组以上分析结果的交集确定固态发酵的核心功能模块。
又一方面，本发明提供一种利用上述任一项所述方法的系统，包括下述单元：
组学检测单元，所述组学检测单元对固态发酵过程取样，对样品进行组学检测；
数据比对单元，所述数据比对单元将检测所得组学数据通过kegg数据库比对，获得整个发酵过程的生化信息；
代谢网络构建单元，所述代谢网络构建单元根据所得生化信息在kegg数据库中所属的代谢模块进行组装，构建全基因尺度代谢网络；
发酵过程可视化单元，所述发酵过程可视化单元在所得全基因尺度代谢网络中，采用不同标识表示基因组、化合物、酶活的含量水平，以固态发酵过程进行可视化；
第一关联单元，所述第一关联单元根据代谢网络属性划分所得全基因尺度代谢网络；
第二关联单元，所述第二关联单元基于所得划分结果，对组学数据与目标研究对象依据统计学方法进行关联性分析；由所得分析结果确定固态发酵的核心功能模块。
又一方面，本发明提供一种上述任一项所述方法的的应用，其特征在于，根据确定的核心功能模块调控固态发酵过程用于筛选微生物、改造微生物、或者设计混菌发酵的菌群。
与现有技术相比，本发明具有如下优势：
1)本发明所述方法为解析固态发酵过程、构建理性调控策略提供了新方向。本发明所述方法构建了“基因-酶-代谢物”多尺度数据的整合，并将图论技术引入代谢网络分析，考虑到了网络自身的特性以及具体的生化反应关联，通过代谢网络的属性并结合组学数据与目标研究对象演替过程的关联性分析，确定了固态发酵的核心功能模块。生物学功能可能并不是单独一个基因所实现的，而是几个基因的共同发挥，可能是几个基因所连成的一个“小模块”，因此相比于确定单个节点，如单独一个基因或者一个酶，可准确挖掘更为有效的固态发酵数据。
2)本发明所述方法可以用于筛选微生物或者设计混菌发酵的菌群。
3)本发明所述方法将“静态”kegg数据库转为“动态”。本发明所述方法对组学数据通过kegg数据库比对，获得整个发酵过程的生化信息，并根据生化信息在kegg数据库中所属的代谢模块进行组装，构建全基因尺度的代谢网络。而后，用不同的颜色和图形表示基因组、化合物、酶活的含量水平，对固态发酵过程进行可视化，而后通过代谢网络的属性结合组学数据与目标研究对象演替过程的关联性分析，确定固态发酵的核心功能模块。核心功能模块即为kegg数据库中ko编号所标识的通路。由此可依据kegg数据库调控固态发酵，实现kegg数据库的“动态”应用。
具体实施方式
下面将对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例仅仅是本发明的一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。
本发明实施例提供一种基于组学整合技术调控固态发酵的方法，包括：
对固态发酵过程取样，对样品进行组学检测；
将检测所得组学数据通过kegg数据库比对，获得整个发酵过程的生化信息；
根据所得生化信息在kegg数据库中所属的代谢模块进行组装，构建全基因尺度代谢网络；
在所得全基因尺度代谢网络中，采用不同的标识分别表示基因组、化合物、酶活的含量水平，以固态发酵过程进行可视化；
根据代谢网络属性划分所得全基因尺度代谢网络；
基于所得划分结果，对组学数据与目标研究对象依据统计学方法进行关联性分析；
由所得分析结果确定固态发酵的核心功能模块。
kegg数据库((京都基因与基因组百科全书)是一个整合了基因组、化学和系统功能信息的数据库。与其他数据库相比，具有强大的图形功能。kegg是一个综合数据库，它们大致分为系统信息、基因组信息和化学信息三大类。进一步可细分为16个主要的数据库。可以通过不同的颜色编码来区分。其中keggpathway数据库是一个手工画的代谢通路的集合。ko编号表示一个通路，这个通路是不分物种的，相当于所有物种的这一通路的并集。k编号表示一个基因，是ko通路中的基本单位，某一k编号代表的不是某一具体物种的基因，而是所有物种的某一同源基因的统称。一旦基因被分配了ko标识或k编号，通过基因组注释流程，物种特异性的通路图、brite功能层次和kegg模块就可以自动产生了。
本发明中，“核心功能模块”指的是发酵过程中关联程度较为紧密、对发酵过程影响较大的模块，即代谢网络中的显著性通路(ko编号所标识的那个通路/模块)。
高通量的宏基因组数据使得人们能够系统地分析微生物的代谢、微生物群落的相互作用及其对环境的反应机制。但本申请发明人发现，现有的高通量的宏基因组数据解析中，一方面其代谢网络构建中，以反应物为端点，以反应为边进行构建，在展示全基因尺度代谢网络时会重叠交织，无法辨识，不易解读。另一方面，在数据挖掘中，只通过统计学方法对组学数据与目标研究对象进行关联性分析，忽略了网络自身的特性以及具体的生化反应关联，以致数据挖掘不够准确。在数据挖掘中，一般通过统计学方法分析得出某个基因或者某个酶具有显著性差异，进而将其视为关键点，由此忽略了代谢网络的整体性，例如忽略了基因与基因之间的关联性。
本发明实施例提供一种组学整合技术调控固态发酵的方法，包括基于组学检测数据的基因组尺度代谢网络的构建以及据此确定固态发酵核心功能模块，而非某个独立的基因。在代谢网络构建中，对固态发酵过程取样，对样品进行组学检测，对组学数据通过kegg数据库比对，获得整个发酵过程的生化信息，并根据生化信息在kegg数据库中所属的代谢模块进行组装，构建全基因尺度的代谢网络。而后，采用不同的标识(例如可以用不同的颜色和图形)表示基因组、化合物、酶活的含量水平，对固态发酵过程进行可视化，通过代谢网络的属性结合组学数据与目标研究对象演替过程的关联性分析，确定固态发酵的核心功能模块。在本发明实施例方法中，首次实现了“基因-酶-代谢物”多尺度数据的整合，并将图论技术引入代谢网络分析，考虑到了网络自身的特性以及具体的生化反应关联，通过代谢网络的属性并结合组学数据与目标研究对象演替过程的关联性分析，确定固态发酵的核心功能模块，实现固态发酵核心功能模块的准确挖掘，为解析固态发酵过程、构建理性调控策略提供了新的方法。
本发明实施例所述方法可用于筛选微生物、改造微生物、或者设计混菌发酵的菌群。
例如，某种产毒素的微生物，通过强化或破坏相应核心功能模块，得到高产某种毒素或者不产某种毒素的微生物。
例如，根据本发明实施例方法确定微生物的核心功能模块后，将模块互补的几个微生物构成小的菌群，进行混菌发酵。
在进一步实施例中，组学整合指的是将基因组、蛋白组和代谢组进行整合。换言之，组学数据包括基因组、蛋白组和代谢组的数据。
将这些组学数据通过kegg数据库比对，获得整个发酵过程的生化信息。获得的生化信息具体可以包括微生物种类、基因、rna、酶、化合物及其对应关系。根据生化信息在kegg数据库中所属的代谢模块进行组装，构建得到全基因尺度代谢网络。所谓“全基因尺度”指的是组学实际所能检测到的发酵过程中所有的微生物、基因、rna、酶和化合物全部考量在内，实现“基因-酶-代谢物”多尺度的整合。
将固态发酵过程可视化后应用统计学方法进行关联性分析，确定核心功能模块。具体地，在本发明实施例中，将图论技术引入代谢网络分析，考虑到了网络自身的特性以及具体的生化反应关联，通过代谢网络的属性并结合组学数据与目标研究对象演替过程的关联性分析，得以确定固态发酵的核心功能模块，也就是代谢网络中的“显著性模块”。
在进一步实施例中，所述代谢网络属性包括拓扑学结构属性和网络统计学属性，所述拓扑学结构属性包括节点度、网络密度、网络直径、网络的连通度；所述网络统计学属性包括网络特征向量中心度和模块化划分。
其中，节点度是指和该节点相关联的边的条数,又称关联度。本发明实施例所述方法中，节点度可以＞5。
网络密度是指实际有的边数与最大可能边数之比。本发明实施例所述方法中，网络密度可以＞0.4。
网络直径是指一个路径中，所有最短路径中的最大值称为网络直径。本发明实施例所述方法中，网络密度可以＜50(链路树)。
网络的连通度是网络节点中所有节点直接联通范围内的节点个数的算术平均值。本发明实施例所述方法中，网络的连通度可以>3。
特征向量中心度，设网络中每个节点的初始值为1，每个节点的度等于相邻节点的度，在一轮计算结束后对所有节点的值做归一化处理，经过多次迭代到达稳态时各个节点的值：
其中，c为比例常数，记做x＝[x1,x2,x3,…,xn]t,达到稳态时，各个节点的值为特征向量中心度。
本发明实施例所述方法中，特征向量中心度可以>0.3。
模块化划分指的是将代谢网络划分成若干小的独立、互相作用的模块。可以是包括一个或多个基因所连成的通路(模块)。本发明实施例所述方法中，根据图的连接关系对节点进行分类，进行qmodularity(q模块度)、edge-betweenness(连边紧密度)、leadingeigenvector(矩阵-特征向量)、fastgreedy(快速贪婪)、multilevel(多级)中任意一种的模块化划分。
在进一步实施例中，统计学方法包括主因子(成分)分析(principalcomponentsanalysis,pca)、聚类分析(clusteranalysis,ca)、冗余度分析(redundancyanalysis,rda)、偏最小二乘法回归分析(partialleastsquareregressionanalysis,pls)、偏最小二乘判别分析(pls-da)和正交偏最小二乘法回归分析(opls/o2pls)。并且采用3种或3种以上方法对组学数据与目标研究对象进行关联性分析；由所得3组或3组以上分析结果的交集确定固态发酵的核心功能模块。
具体地，每组分析结果可以选定前100位，3组或多组分析结果中均在前100位者，即为固态发酵的核心功能模块。
本文在此对根据统计学方法进行关联分析的具体过程不再赘述，依据现有技术的关联分析过程进行即可。需要注意的是，本申请旨在构建“基因-酶-代谢物”多尺度的整合，在关联分析中将网络自身的特性以及具体的生化反应关联与结合与组学数据生物物质之间演替过程的关联性分析相结合，得到固态发酵过程的核心功能模块，以此调控固态发酵。
相应的，本发明又一实施例提供一种利用上述任一项所述方法的系统，包括下述单元：
组学检测单元，所述组学检测单元对固态发酵过程取样，对样品进行组学检测；
数据比对单元，所述数据比对单元将检测所得组学数据通过kegg数据库比对，获得整个发酵过程的生化信息；
代谢网络构建单元，所述代谢网络构建单元根据所得生化信息在kegg数据库中所属的代谢模块进行组装，构建全基因尺度代谢网络；
发酵过程可视化单元，所述发酵过程可视化单元在所得全基因尺度代谢网络中，采用不同标识(例如可以用不同的颜色和图形)表示基因组、化合物、酶活的含量水平，以固态发酵过程进行可视化；
第一关联单元，所述第一关联单元根据代谢网络属性划分所得全基因尺度代谢网络；
第二关联单元，所述第二关联单元基于所得划分结果，对组学数据与目标研究对象依据统计学方法进行关联性分析；由所得分析结果确定固态发酵的核心功能模块。
需要说明的是，图论技术以及基于统计学方法将组学数据与目标研究对象的关联性分析是行业内本领域技术人员的公知，在此不再赘述。本发明在确定核心模块过程中也应用到了统计学方法。然而本发明不同的是，在此之前，均是仅通过统计学方法确定单个基因，而本发明发现生物学功能可能并不是单独一个基因所实现的，而是几个基因的共同发挥，可能是几个基因所连成的一个“小模块”。本发明通过构建“基因-酶-代谢物”多尺度数据的整合，将图论技术引入代谢网络分析，考虑到了网络自身的特性以及具体的生化反应关联，通过代谢网络的属性并结合组学数据与目标研究对象演替过程的关联性分析，确定了固态发酵的核心功能模块。相比于确定单个节点，如单独一个基因或者一个酶，可准确挖掘更为有效的固态发酵数据。将“静态”kegg数据库转为“动态”。
在本说明书的描述中，参考术语“一实施例”、“一些实施例”、“进一步实施例”、“具体示例”、或“一些示例”等的描述意指结合该实施例或示例描述的具体特征、结构、材料或者特点包含于本发明的至少一个实施例或示例中。在本说明书中，对上述术语的示意性表述不一定指的是相同的实施例或示例。而且，描述的具体特征、结构、材料或者特点可以在任何的一个或多个实施例或示例中以合适的方式结合。
下文以具体示例对本发明实施例所述基于组学整合技术确定固态发酵核心功能模块，调控固态发酵的方法进一步阐述。
实施例1
对山西老陈醋乙偶姻的代谢模块进行了分析。
对山西老陈醋醋酸发酵阶段1天、3天、5天、7天和9天的样本进行基因组、转录组和蛋白组的检测，通过kegg数据库比对，共获得18500种化合物，11144个反应，7524种酶，共有341微生物，绘制基因组尺度组学数据演替过程图谱。以乙偶姻为中心，筛选网络直径为50以内的，节点度>5、网络密度>0.4、网络的连通度>3、特征向量中心度>0.3、经qmodularity模块化划分包含乙偶姻在内的模块网络。
依据统计学方法筛选同时满足在偏最小二乘法回归分析、偏最小二乘判别分析和正交偏最小二乘法回归分析中与乙偶姻相关性前100的化合物、酶和基因。最终筛选出了合成乙偶姻的核心模块。偏最小二乘法回归分析、偏最小二乘判别分析和正交偏最小二乘法回归分析的具体实施方法为：以各个酶、基因、化合物的变化为自变量，以乙偶姻的变化为因变量，计算其相关性。最终确定的核心功能模块包括ko00500、ko00010、ko00620、ko00360、ko00030、ko00520、ko00400、ko00650、ko00061、ko00250、ko00330、ko00290、ko00910、ko00473、ko01040代谢路径及其包含酶类。
实施例2
对实验室酱油固态发酵工艺进行组学检测。共测得10144个反应，7524种酶，体系中优势细菌分别属于weissella(魏斯氏菌),staplococcus(葡萄球菌),pediococcus(小球菌),lactobacillus(乳酸杆菌),bacillus(芽孢杆菌)等菌属。在属水平上，酱油微生物包含184个细菌属以及116个真菌属。以这些化合物绘制基因组尺度组学数据演替过程图谱。以4-乙基愈创木酚为中心，筛选网络直径为50以内的，节点度>5、网络密度>0.4、网络的连通度>3、特征向量中心度>0.3、经qmodularity模块化划分包含4-乙基愈创木酚在内的模块网络。
同时满足在偏最小二乘法回归分析、偏最小二乘判别分析和正交偏最小二乘法回归分析分析中与乙偶姻相关性前100的化合物、酶和基因，最终筛选出了合成4-乙基愈创木酚的核心模块。偏最小二乘法回归分析、偏最小二乘判别分析和正交偏最小二乘法回归分析的具体实施方法为：以各个酶、基因、化合物的变化为自变量，以4-乙基愈创木酚的变化为因变量，计算其相关性。核心功能模块包括ko00090、ko02400、ko00620、ko00360、ko00030、ko00520、ko00480、ko00650、ko00240、ko00650、ko00230、ko00290、ko01910、ko00473、ko02040代谢路径及其包含酶类。
实施例3
对实验室固态发酵白酒酒醅中进行组学检测，共测得9986个反应，6548种酶，在属水平上，酿酒微生物包含168个细菌属以及104个真菌属。以这些化合物绘制基因组尺度组学数据演替过程图谱。以四甲基吡嗪为中心，筛选网络直径为50以内的，节点度>5、网络密度>0.2、网络的连通度>6、特征向量中心度>0.2、经qmodularity模块化划分包含苯乙醇在内的模块网络。
同时满足在偏最小二乘法回归分析、偏最小二乘判别分析和正交偏最小二乘法回归分析分析中与乙偶姻相关性前100的化合物、酶和基因，最终筛选出了合成苯乙醇的核心模块。偏最小二乘法回归分析、偏最小二乘判别分析和正交偏最小二乘法回归分析的具体实施方法为：以测得的各个酶、基因、化合物的变化为自变量，以苯乙醇的变化为因变量，计算其相关性。最终筛选出了四甲基吡嗪的核心模块。该模块的核心路径包括ko00058、ko02300、ko00119、ko00240、ko00031、ko00354、ko00442、ko00623、ko00161、ko00231、ko00645、ko00165、ko01921、ko00433、ko02020代谢路径及其包含酶类。
以上所述，仅为本发明的具体实施方式，但本发明的保护范围并不局限于此，任何熟悉本技术领域的技术人员在本发明揭露的技术范围内，可轻易想到变化或替换，都应涵盖在本发明的保护范围之内。因此，本发明的保护范围应所述以权利要求的保护范围为准。