基于调控-代谢网络整合模型预测微生物生长表型的方法

摘要

本发明公开了一种基于调控-代谢网络整合模型预测微生物生长表型的方法。该方法首先构建基因调控网络，基于大量基因表达谱数据，运行多次线性回归推断每个目标基因随转录因子表达变化的线性方程，计算错误发现率FDR；然后取FDR<＝0.05的调控关系作为全局的调控网络，找出调控网络中调控代谢基因的转录因子，根据转录因子的类型，计算转录因子敲除时生长速率对应的反应流值F，并由原始代谢网络进行同样的流平衡分析得到最大细胞生长速率Fmax，通过计算F/Fmax预测微生物在转录因子敲除时的生长表型的变化。本发明方法使得分析的精确度提高，得以更好地预测微生物的生长表型。

说明书

技术领域

本发明属于微生物技术领域，具体的说，涉及一种基因调控网络与代谢网络联合建模的 方法，其可用于预测微生物的生长表型。

背景技术

目前，代谢网络建模是预测基因敲除对代谢表型影响的主要方法。代谢网络包含所有酶、 代谢物和生化反应之间的相互作用，全基因组范围代谢网络重构与分析能够发现基因敲除、 插入、异常表达、环境变化等对生物系统表型的影响。

动力学模拟是分析代谢网络最直接有效的手段,但由于很多动力学参数未知使得大规模网 络的动力学分析受到限制。而基于约束的代谢流分析、基元模式和极端途经分析、最优化方 法等得到更加广泛的应用,其中流平衡分析最为常用。流平衡分析(FluxBalanceAnalysis)是一 种基于约束的代谢网络分析方法，通过加入化学计量信息、热力学信息、代谢能力约束等,将 系统行为约束在一个封闭的解空间内,然后设置目标函数通过线性规划得到最优解。通常目标 函数设为细胞生物质合成最大、生长速率最快或某个目标代谢物产量最高。

以大肠杆菌和酵母为代表的微生物代谢模型研究中，流平衡分析能够比较准确地模拟基 因敲除对细胞生长的影响。流平衡分析可以不依赖于反应动力学参数等信息，快速识别对生 长或目标代谢物合成具有重要影响的必需基因，如果某个基因敲除造成生长速率为零，则属 于致死基因，不能改动；如果某个基因敲除或过表达能够提高目标代谢物产量，则成为候选 的改造基因。通过这种计算模拟的手段能够寻找可行的改造途径，对于代谢工程高产菌株设 计具有重要的指导意义。

如果能够将基因调控网络整合进来，采用合理方式设定约束，将有助于提高代谢流的预 测，发现生物在不同环境下的响应变化，使得现有的代谢模型能更加精确地模拟基因敲除对 生长表型的影响。基因调控网络描述调控因子与目标基因之间的相互作用关系，近年来有很 多研究分别对两种网络进行构建和分析，但是将二者整合起来揭示转录调控对代谢表型影响 的研究较少。这是由于调控网络是基于统计推断的，而代谢网络是基于生化知识的，如何实 现基因组范围的调控网络与代谢网络整合是极富挑战的问题。

发明内容

针对现有技术的不足，本发明的目的在于提供一种基于调控-代谢网络预测微生物在不同 环境下生长表型的方法。该方法简单、高效，预测结果准确。

本发明提供一种基于调控-代谢网络整合模型预测微生物生长表型的方法，具体步骤如下：

(1)首先构建基因调控网络，基于大量基因表达谱数据，用线性回归推断每个目标基 因随转录因子表达变化的线性方程，然后随机抽取部分表达谱数据进行若干次bootstrap线性 回归，根据转录因子和目标基因的调控关系在这若干个线性方程中出现与否的概率，计算错 误发现率FDR(FalseDiscoveryRate)，FDR表明该转录因子没有调控该基因的比例。

(2)然后取FDR<＝0.05的调控关系作为全局的调控网络；

(3)接着找出调控网络中调控代谢基因的转录因子，根据转录因子的类型，设定敲除 该转录因子时目标代谢基因表达的概率P(Gene＝ON|Factor＝OFF)，如果是激活因子,则P＝FDR, 如果是抑制因子，则P＝1-FDR；

(4)由代谢网络的流变分析FVA(FluxVariabilityAnalysis)可得每个反应的最大流值 Vmax，对每个转录因子进行模拟敲除时，对每个反应施加相应的约束为Vmax*P，接着进行 目标函数为生长速率最大的流平衡分析，得出细胞生长速率对应的反应流值F。

(5)对于没有任何基因敲除的野生型代谢网络进行同样的流平衡分析得出最大细胞 生长速率Fmax，然后计算突变体相对于野生型的生长速率之比F/Fmax，即代表生长表型的 变化。

本发明中，随机抽取部分表达谱数据进行150-300次bootstrap线性回归。

本发明的有益效果在于，对比不施加其他约束的流平衡分析，结合了基因调控网络的信 息，使得分析的精确度提高，得以更好地预测微生物的生长表型。

附图说明

图1是本发明方法预测生长表型的总体流程示意图。

图2是以葡萄糖－铵为细胞生长底物时的ROC曲线的结果图示。

图3是以半乳糖－铵为细胞生长底物时的ROC曲线的结果图示。

图4是以葡萄糖－尿素为细胞生长底物时的ROC曲线的结果图示。

具体实施方式

以下结合附图和实施例对本发明的技术方案作进一步描述。

本发明的实施例中，采用的算法流程图如图1所示。本实施例以酵母为例，首先基于已搜 集的2929组酵母的基因表达谱数据，用线性回归推断每个目标基因随转录因子表达变化的线 性方程，若某转录因子在方程中的系数为正，则代表有激活作用，若系数为负，则代表有抑 制作用，若系数为零，则代表没有调控作用。然后随机抽取2929组表达谱数据中的子集进行 200次bootstrap线性回归，根据转录因子和目标基因的调控关系在这200个线性方程中出现与 否的概率，计算错误发现率FDR；然后取FDR<＝0.05的可靠调控关系，形成包含2588对转录 因子与目标基因的相互作用的酵母全局调控网络。

接下来输入酵母的代谢网络，以Yeast6模型为例,该模型包含1888个反应,1458个代谢 物,900个编码酶的基因，用流平衡分析得出最大细胞生长速率Fmax。找出其中有调控关系的 代谢基因，根据调控其的转录因子类型，设定转录因子敲除时，该代谢基因表达的概率 P(Gene＝ON|Factor＝OFF),如果是激活因子,则P＝FDR,如果是抑制因子,则P＝1-FDR；再对完整的 Yeast6代谢网络进行流变分析，可得每个反应的最大流值Vmax，则对每个转录因子进行模拟 敲除时，将对每个反应施加新的约束为Vmax*P；接着进行目标函数为生长速率最大的流平衡 分析，得出细胞生长速率的反应流值F；最后计算突变体相对于野生型的生长速率之比 F/Fmax。

下面验证预测结果的有效性，实验数据来源于文献i，共有三种环境下的实验数据，分别 以葡萄糖－铵，半乳糖－铵，葡萄糖－尿素为细胞生长的底物，对此我们也预测了这三种情 况下敲除转录因子的生长表型。将我们预测的突变体相对于野生型的生长速率之比与实验中 对应的比值进行比较，与实验对应上的有51个转录因子敲除的突变体，对每个环境计算预测 值与实验值这两组数据的pearson相关性，结果如表1所示。

表1预测值与实验值这两组数据的pearson相关性结果

葡萄糖－铵 半乳糖－铵 葡萄糖－尿素 相关系数 0.4594311 0.3128724 0.3535231 p-value 0.0006953 0.02539 0.01093

结果显示p值均小于0.05，可以认为我们的预测值与实验值显著线性相关。另外，对实验 值设定阈值，将数据二元化，表示是否影响生长；以此为标准，通过ROC曲线判断预测结果 是否能正确反映转录因子敲除对突变体生长的影响。

图2、图3、图4分别显示了以葡萄糖－铵，半乳糖－铵，葡萄糖－尿素为细胞生长底 物的三种环境下ROC曲线的结果，其中对实验值设定了三个阈值分别为0.2、0.5和0.8；总 体来说，所有的ROC曲线都在对角线之上，曲线下面积值如各图所示，因此说明预测的结果 可以有效判断转录因子敲除对生长的影响。综上，我们通过DIRAM方法可以有效地预测酵 母在三种不同环境下的生长表型。

参考文献：

1、FendtSM,OliveiraAP,ChristenS,etal.Unravelingcondition‐dependentnetworksof transcriptionfactorsthatcontrolmetabolicpathwayactivityinyeast[J].Molecularsystemsbiology, 2010,6(1):432.