利用细胞代谢网络监测生物制造过程的方法

1.本发明涉及工业生物过程控制技术领域，尤其是指一种利用细胞代谢网络监测生物制造过程的方法。


背景技术：

2.随着全球经济发展和社会进步，为了满足人们日益增长的生活需求，以工业微生物技术为代表的小批量、高利润、多品种生产方式受到重视。生物制造过程一般在生物反应器中进行，其可以通过批次生产的方式在一定周期时间内获得经济效益高的产物，比如常用抗生素青霉素的生产，其外源分泌蛋白等一般由生物发酵反应而得。然而在发酵生产过程中经常会出现各类生产状况的失控或运行故障，如果不能及早监测、预警并及时调控，最终会影响生产质量和目标产物的产量。因此，生物制造过程运行状况评估与故障监测，具有重要的实际工程价值。
3.生物发酵过程中同时存在物化和生化反应，细胞内部机理复杂，而实际生产过程中一般只能检测到温度、ph、通风率、搅拌功率等基本生产工艺变量及一些可观测的细胞外物质浓度，对于细胞内的微观的代谢活动则无法观测，导致现有的发酵生产运行状况评估技术主要依据生产过程中各工艺变量和物质浓度进行综合分析，其监测的实效性和及时性有待提高。


技术实现要素：

4.为此，本发明所要解决的技术问题在于克服现有技术中于细胞内微观的代谢活动无法观测，导致现有的发酵生产运行状况评估不足的问题，提供一种利用细胞代谢网络监测生物制造过程的方法，该方法既利用微生物制造过程可检测或可观测信息，又充分挖掘细胞生长繁殖过程中不可检测的胞内代谢信息，能够更加有效的实施对生产制造过程运行状况评估及故障监测。
5.为解决上述技术问题，本发明提供了一种利用细胞代谢网络监测生物制造过程的方法，包括如下步骤：基于代谢网络进行通量均衡分析，得到反映细胞生长与繁殖过程的代谢通量；将细胞内代谢通量与生产过程变量一起构成训练数据集，建立生物发酵过程监测模型及统计量控制限；基于监测模型，利用实时采集的生产过程变量和代谢通量计算值进行运行状况在线监测。
6.在本发明的一个实施例中，所述通量均衡分析过程包括：建立通量均衡方程；建立胞外代谢物浓度与代谢通量的动态关系；根据多批次胞外代谢物浓度采样值估计代谢通量。
7.在本发明的一个实施例中，所述代谢通量的随时间变化，采集相邻两个时间点之间的代谢通量变化规律对代谢通量进行瞬态特征描述，包括常数、线性函数、二次函数。
8.在本发明的一个实施例中，多批次胞外代谢物浓度采样值根据生产实际来选取，在计算时间间隔内采样一次或多次采样，再进行曲线拟合。
9.在本发明的一个实施例中，在估计代谢通量时，通过求导数、线性规划、二次规划的优化算法求得最优代谢通量。
10.在本发明的一个实施例中，建立生物发酵过程监测模型包括：将生产过程变量和细胞内代谢通量合并构成扩展样品数据集；对样品数据集进行标准化处理生成标准样品数据集；利用主元分析建立监测模型；确定监测统计指标的控制限。
11.在本发明的一个实施例中，在扩展样品数据集时，所述细胞内代谢通量选取全部代谢通量或胞内外交换代谢通量或内部代谢通量。
12.在本发明的一个实施例中，所述主元分析包括一般主元分析、多阶段主元分析、核主元分析以及支持向量机方法。
13.在本发明的一个实施例中，所述运行状况在线监测包括：采集当前批次当前时刻的胞外代谢物浓度和生产过程变量；由胞外代谢物浓度估计获得的代谢通量；将生产过程变量与代谢通量合并构成测试数据集；对测试数据集进行标准化处理生成测试标准数据集；预估填充未来时刻的测试数据；计算统计标准，并与统计量控制限比较。
14.在本发明的一个实施例中，预估填充未来时刻的测试数据的方法包括直接用前若干批次对应时刻的均值数据、或考虑本批次当前及过去时刻实际数据与对应均值数据之差，预估填充未来数据。
15.本发明的上述技术方案相比现有技术具有以下优点：
16.本发明所述的利用细胞代谢网络监测生物制造过程的方法，在多批次生产过程数据基础上，利用细胞代谢网络和胞外可观测代谢物浓度，通过估计细胞代谢通量，将代谢通量与生产过程变量结合，构建扩展样本数据集，给出一种在线监测过程运行状况的方法，该方法既利用微生物制造过程可检测或可观测信息，又充分挖掘细胞生长繁殖过程中不可检测的胞内代谢信息，能够更加有效的实施对生产制造过程运行状况评估及故障监测，能够及时发现过程运行中的异常，以便做出相应的控制措施，对于减少安全隐患、提高企业综合效益具有明显的作用。
附图说明
17.为了使本发明的内容更容易被清楚的理解，下面根据本发明的具体实施例并结合附图，对本发明作进一步详细的说明，其中
18.图1是本发明利用细胞代谢网络监测生物制造过程的方法步骤图；
19.图2是本发明利用细胞代谢网络监测生物制造过程的方法流程示意图；
20.图3青霉素发酵过程部分可观测代谢物浓度计算结果。
21.图4青霉素发酵过程部分代谢通量动态估计。
具体实施方式
22.下面结合附图和具体实施例对本发明作进一步说明，以使本领域的技术人员可以更好地理解本发明并能予以实施，但所举实施例不作为对本发明的限定。
23.参照图1所示，本发明的利用细胞代谢网络监测生物制造过程的方法，包括如下步骤：基于代谢网络进行通量均衡分析，得到反映细胞生长与繁殖过程的代谢通量；将细胞内代谢通量与生产过程变量一起构成训练数据集，建立生物发酵过程监测模型及统计量控制
限；基于监测模型，利用实时采集的生产过程变量和代谢通量计算值进行运行状况在线监测；在多批次生产过程数据基础上，利用细胞代谢网络和胞外可观测代谢物浓度，通过估计细胞代谢通量，将代谢通量与生产过程变量结合，构建扩展样本数据集，给出一种在线监测过程运行状况的方法，该方法既利用微生物制造过程可检测或可观测信息，又充分挖掘细胞生长繁殖过程中不可检测的胞内代谢信息，能够更加有效的实施对生产制造过程运行状况评估及故障监测，能够及时发现过程运行中的异常，以便做出相应的控制措施，对于减少安全隐患、提高企业综合效益具有明显的作用。
24.参照图2所示，所述通量均衡分析过程包括：建立通量均衡方程；建立胞外代谢物浓度与代谢通量的动态关系；根据多批次胞外代谢物浓度采样值估计代谢通量；建立生物发酵过程监测模型包括：将生产过程变量和细胞内代谢通量合并构成扩展样品数据集；对样品数据集进行标准化处理生成标准样品数据集；利用主元分析建立监测模型；确定监测统计指标的控制限；所述运行状况在线监测包括：采集当前批次当前时刻的胞外代谢物浓度和生产过程变量；由胞外代谢物浓度估计获得的代谢通量；将生产过程变量与代谢通量合并构成测试数据集；对测试数据集进行标准化处理生成测试标准数据集；预估填充未来时刻的测试数据；计算统计标准，并与统计量控制限比较。
25.建立通量均衡方程：
26.微生物发酵过程中胞内代谢物的净积累或消耗速率远小于胞内外代谢物之间的交换反应速率，对于一个总化学反应个数为f并且包含了m个胞内代谢物和r个胞外可测代谢物的代谢网络，其通量均衡方程如下：
[0027][0028]
其中，c
r，t
∈r
r
表示r个可测的胞外代谢物在t时刻浓度向量，s
r
∈r
r
×
f
表示胞外可测代谢物的化学计量矩阵，v
t
∈r
f
表示代谢网络内代谢通量向量，包括胞内代谢通量及胞内外交换通量，s
m
∈r
m
×
f
表示m个胞内代谢物的化学计量矩阵。
[0029]
微生物代谢网络中化学反应个数f多于细胞内代谢物数量m，其自由度 d＝f
‑
rank(s
m
)，采用一组独立通量u
t
∈r
d
的线性组合表示代谢通量，即 v
t
＝k
·
u
t
，其中k∈r
f
×
d
表示s
m
零空间的一组标准正交基，u
t
称为自由通量。
[0030]
代谢通量的瞬态特征描述：
[0031]
考虑代谢通量v
t
随时间变化，引入计算时间序列点t
d
，d＝1，...，z，用以描述分析v
t
的动态变化，即将整个生产过程划分为z
‑
1个时间间隔，在每个时间间隔内，即在相邻两个时间点t
d
和t
d+1
之间，假设代谢通量按一定函数规律变化，可以是常数、线性函数、二次型函数。选用线性函数后的代谢通量v
t
表示如下：
[0032][0033]
其中，u
d
表示第d个计算时刻的自由通量，k(t，t
d
)是各时间间隔内代谢通量线性变化的系数函数，按表1选取。
[0034][0035]
表1
[0036]
基于代谢通量的胞外代谢物浓度计算：
[0037]
将式(2)代入式(1)，得：
[0038][0039]
对式(3)两边积分得：
[0040][0041]
其中，函数而c
r0
＝[c
1，0
，c
2，0
，
…
，c
r，0
]
t
表示r个胞外可观测代谢物的初始浓度向量。
[0042]
生产过程中实际采样观测时间点为t
′
n
，n＝1，...，n，根据式(4)，第r个胞外代谢物，r＝1，...，r，在实际采样时刻t
n
的浓度计算值为：
[0043][0044]
其中，c
r，0
表示第r个胞外代谢物初始浓度值。
[0045]
采样点可以根据生产实际来选取，简单操作是在代谢通量v
t
的每个计算时间间隔内进行一次采样观测，即z＝n
‑
1，也可以在计算时间间隔内进行多次采样；一般情形下代谢通量v
t
的第一个和最后一个计算时间点应分别对应于实际采样观测的初始和最终时刻。
[0046]
各采样时刻的代谢通量估计：
[0047]
在采样观测时间点t
′
n
，各胞外代谢物的实际采样值为c
r，n
，则浓度计算值与采样值的偏差为：
[0048][0049]
其中，为胞外代谢物r在采样时刻n的浓度误差项方差。
[0050]
以表示各个计算时刻的d维自由通量，通过优化算法寻找满足式(5)的u，使得偏差φ最小，即求解获得最优的自由通量u，也即各t
d
时刻的自由通量u
d
，再利用式(2)估计代谢通量v
t
。
[0051]
简单的优化算法是直接由式(6)对u求导数，即得到优化变量 u＝h
‑1·
j，其中矩阵h和矩阵j由下式给定：
[0052][0053][0054]
式中w
r
是由构成的对角矩阵。
[0055]
工程实际中可增加初始浓度c
r0
、代谢通量计算时刻t
d
为优化变量，利用优化算法求解获得最优的自由通量u、初始浓度c0以及时刻t
d
。
[0056]
建立运行状况监测模型及指标：
[0057]
扩展样本数据集
[0058]
生产中已有的i个正常批次数据，其中第i批次生产变量(包括温度、 ph通风率、搅拌功率、底物流加速率等)数据矩阵为x
i
∈r
n
×
o
，i＝1，...，l，n 为采样观测点数，0为生产变量数。根据估计所得的代谢通量v
t
，选取其中全部或部分代谢通量，若选取l个胞内外交换代谢通量在第i批次的数据构成集合y
i
∈r
n
×
l
，则该批次的生产变量数据与胞内外交换代谢通量合并构成扩展样本数据集g
i
＝[x
i
，y
i
]
t
∈r
n
×
h
，其中h＝o+l。全部i个批次的发酵过程扩展样本数据集则为按批次方向展开得 g∈r
i
×
n
·
h
。
[0059]
数据集的标准化处理
[0060]
将扩展样本数据集标准化变换为均值为0、方差为1的标准训练数据集，近似服从多维正态分布：
[0061][0062]
其中，mean(g)为g中同一时刻各个批次同一变量的均值，即多批次正常操作下的平均运行轨迹，std(g)是g中同一时刻各个批次同一变量的标准差。
[0063]
监测模型及控制限的求取
[0064]
构造标准数据集g
train
的协方差矩阵利用主元分析建立监测模型，先求取其特征值及其特征向量，按累计方差贡献率从大到小，选取a个主元特征值所对应的特征向量q1，
…
，q
a
，构成载荷矩阵q＝「q1，
…
，q
a
]，进而计算g
train
的各批次得分向量和残差向量，全部i批次得分向量和残差向量分别构成得分矩阵和残差矩阵，据此求取常用的t2控制限和spe控制限。
[0065]
监测模型也可以采用多阶段多向主元分析、多向核主元分析和多类支持向量机等方法。
[0066]
在线实时运行状况监测：
[0067]
当前时刻数据处理
[0068]
采集生产过程当前时刻n的生产变量数据x
test，n
∈r
n
×
o
，通过估计获得的代谢通量y
test，n
∈r
n
×
l
，将生产变量和代谢通量合并得到测试数据集 g
test，n
∈r
n
×
h
；以训练数据集中相应时刻的平均值和标准差对g
test，n
∈r
n
×
h
进行标准化处理。
[0069]
未来时刻数据估计
[0070]
当前时刻之后的n+1至批次结束n时刻的数据，采用预估的方法进行填充，用标准化的g
test，n
中第n时刻数据值进行填充，得到g
′
test，n
∈r
n
×
h
，再按时间方向展开为g
″
test，n
∈r1×
n
·
h
。也可以直接用前若干批次对应时刻的均值数据，或考虑当前及过去时刻实际数据与对应均值数据之差，预估填充未来数据。
[0071]
计算当前时刻t2统计量和spe统计量
[0072]
根据标准训练数据集g
train
的载荷矩阵q计算当前n时刻得分向量t
test，n
和残差向量e
test，n
：
[0073]
t
test，n
＝g
test，n
·
q
[0074]
g
new，n
＝t
test，n
·
q
y
ꢀꢀꢀꢀꢀ
(8)
[0075]
e
tesst，n
＝g
test，n
‑
g
new，n
[0076]
按下式计算n时刻t2统计量和spe统计量：
[0077][0078]
其中，ψ∈r
a
×
a
是由g
train
的第i批次得分向量中各个元素构成的对角矩阵。
[0079]
过程运行状况监测
[0080]
若t2和spe
n
低于统计量控制限，则过程运行正常；否则预警过程有异常或故障。
[0081]
为了验证本实施例的利用细胞代谢网络监测生物制造过程的方法的效果，将方法应用于青霉素生产过程运行监测：
[0082]
青霉素发酵代谢网络共包含66个代谢通量(其中胞内外交换通量5个) 和49种胞内代谢物，葡萄糖、青霉素、菌体生物质、氧气及二氧化碳等5 种胞外代谢物可检测，考虑反应温度、ph值、通风率、搅拌功率、底物流加速率和冷却水流加速率等生产工艺变量和胞内外交换代谢通量共11个，选取30批次、每批次800采样点，形成扩展样本数据集。
[0083]
按本发明方法完成代谢通量的胞外代谢物浓度计算和各采样时刻的代谢通量估计，图3给出胞内外交换代谢通量的估计值，由此扩展样本数据集不仅包含了生物发酵过程可检测的生产工艺变量信息，也包含了细胞内部不可测的代谢信息，更全面地反映生产过程特征。
[0084]
针对青霉素发酵过程由通风、搅拌及底物流加引起的典型故障，图4给出了本发明在线监测结果：
[0085]
其中(a)表示的过程监测只采用了生产工艺变量，而(b)则同时结合代谢网络信息。从t2监测图中可以发现，当故障发生时，本发明所提出的方法能够及时准确地监测故障，而基于生产变量监测有一定时间滞后。
[0086]
本领域内的技术人员应明白，本技术的实施例可提供为方法、系统、或计算机程序产品。因此，本技术可采用完全硬件实施例、完全软件实施例、或结合软件和硬件方面的实施例的形式。而且，本技术可采用在一个或多个其中包含有计算机可用程序代码的计算机可用存储介质(包括但不限于磁盘存储器、cd
‑
rom、光学存储器等)上实施的计算机程序产品的形式。
[0087]
本技术是参照根据本技术实施例的方法、设备(系统)、和计算机程序产品的流程图和/或方框图来描述的。应理解可由计算机程序指令实现流程图和/或方框图中的每一流
程和/或方框、以及流程图和/或方框图中的流程和/或方框的结合。可提供这些计算机程序指令到通用计算机、专用计算机、嵌入式处理机或其他可编程数据处理设备的处理器以产生一个机器，使得通过计算机或其他可编程数据处理设备的处理器执行的指令产生用于实现在流程图一个流程或多个流程和/或方框图一个方框或多个方框中指定的功能的装置。
[0088]
这些计算机程序指令也可存储在能引导计算机或其他可编程数据处理设备以特定方式工作的计算机可读存储器中，使得存储在该计算机可读存储器中的指令产生包括指令装置的制造品，该指令装置实现在流程图一个流程或多个流程和/或方框图一个方框或多个方框中指定的功能。
[0089]
这些计算机程序指令也可装载到计算机或其他可编程数据处理设备上，使得在计算机或其他可编程设备上执行一系列操作步骤以产生计算机实现的处理，从而在计算机或其他可编程设备上执行的指令提供用于实现在流程图一个流程或多个流程和/或方框图一个方框或多个方框中指定的功能的步骤。
[0090]
显然，上述实施例仅仅是为清楚地说明所作的举例，并非对实施方式的限定。对于所属领域的普通技术人员来说，在上述说明的基础上还可以做出其它不同形式变化或变动。这里无需也无法对所有的实施方式予以穷举。而由此所引伸出的显而易见的变化或变动仍处于本发明创造的保护范围之中。