（青島科技大學(xué)化工學(xué)院，山東青島266042）

通量平衡分析被廣泛應(yīng)用于分析代謝網(wǎng)絡(luò)的功能［1］。基于穩(wěn)態(tài)的假設(shè)，該方法能夠準(zhǔn)確預(yù)測細(xì)胞的最大生長速率并給出一組滿足條件的通量分布。但是最大生長速率下的通量分布不止一組，要獲得全部的解空間可使用基元模式分析［2］，該方法所得的基元模式個數(shù)隨代謝網(wǎng)絡(luò)中反應(yīng)個數(shù)的增加呈指數(shù)增加，因此，基因組尺度的代謝網(wǎng)絡(luò)通過枚舉基元模式的方法獲得所有通量分布也是不可行的。通量可變性分析可以對基因組尺度的代謝網(wǎng)絡(luò)進(jìn)行計(jì)算得到目標(biāo)函數(shù)值下每個反應(yīng)通量值的范圍（即通量解空間），該通量解空間的體積可用來衡量微生物的代謝靈活性。

微生物為了生存會利用自身的代謝能力適應(yīng)不斷改變的外界環(huán)境，改變內(nèi)部代謝通量分布。該適應(yīng)性行為具有多樣性但必須滿足化學(xué)計(jì)量學(xué)約束，而當(dāng)菌體的適應(yīng)性行為只有化學(xué)計(jì)量學(xué)約束，而沒有動力學(xué)約束和代謝調(diào)控約束時，可以利用通量可變性分析對其解空間進(jìn)行計(jì)算，進(jìn)而衡量菌體的代謝靈活性。西班牙學(xué)者San等［3］將菌體代謝靈活性改變的原因分為3部分：內(nèi)部反應(yīng)的可變性、交換反應(yīng)的可變性和生長速率的可變性，并利用基因組尺度的大腸桿菌代謝網(wǎng)絡(luò)iJO1366，以葡萄糖為底物對這3 個因素的重要性進(jìn)行了分析。

克雷伯氏桿菌以生物柴油的副產(chǎn)物甘油為底物發(fā)酵生產(chǎn)重要化學(xué)原料1，3－丙二醇的生物方法［4］近年來已經(jīng)實(shí)現(xiàn)了工業(yè)化。許多學(xué)者對如何提高1，3－丙二醇產(chǎn)率進(jìn)行了多方面研究，并取得了相應(yīng)的成果。曾宏等［5］對甘油脫氫酶的生產(chǎn)工藝進(jìn)行模擬并進(jìn)行經(jīng)濟(jì)分析；朱曉麗等［6］使用大氣壓介質(zhì)阻擋放電等離子體與LiCl進(jìn)行復(fù)合誘變，結(jié)合FeCl3和產(chǎn)物耐受的篩選方法獲得1，3－丙二醇高產(chǎn)突變菌株；陳琳等［7］以木糖作為發(fā)酵過程中的輔助底物與甘油共發(fā)酵生產(chǎn)1，3－丙二醇，1，3－丙二醇產(chǎn)率和甘油轉(zhuǎn)化率均得到提高；徐衛(wèi)濤等［8］通過控制發(fā)酵過程中的pH 值提高1，3－丙二醇的產(chǎn)率。以上均是通過改變外部條件來提高1，3－丙二醇產(chǎn)率。目前，還沒有學(xué)者對克雷伯氏桿菌代謝甘油發(fā)酵生長應(yīng)對外界環(huán)境改變的代謝靈活性進(jìn)行理論分析。鑒于此，作者使用代謝流分析軟件COBRA Toolbox，對克雷伯氏桿菌基因組尺度代謝網(wǎng)絡(luò)進(jìn)行通量平衡分析與通量可變性分析，分別得到內(nèi)部反應(yīng)的可變性、交換反應(yīng)的可變性和生長速率的可變性等3個因素對菌體代謝靈活性的影響程度，通過進(jìn)一步研究生長速率對菌體代謝靈活性的影響，找到提高代謝靈活性的兩個必要條件，從而增強(qiáng)菌體應(yīng)對外界環(huán)境改變的能力，指導(dǎo)菌體實(shí)際生長。

1 模型與方法

1.1 克雷伯氏桿菌基因組尺度代謝網(wǎng)絡(luò)

本實(shí)驗(yàn)使用Liao等［9］重構(gòu)的克雷伯氏桿菌基因組尺度代謝網(wǎng)絡(luò)iYL1228，該網(wǎng)絡(luò)包括1 229個基因、1 658個代謝物、2 262個生物反應(yīng)，并按照功能將反應(yīng)劃分為62個子系統(tǒng)。在構(gòu)建該基因組尺度代謝網(wǎng)絡(luò)時，學(xué)者們對84%的底物進(jìn)行生長模擬，發(fā)現(xiàn)其模擬數(shù)據(jù)與實(shí)驗(yàn)數(shù)據(jù)的錯誤率低于31%，大部分的錯誤率集中在0.4%～26.3%之間，顯示出良好的一致性，因此可以用該網(wǎng)絡(luò)進(jìn)行模擬計(jì)算。克雷伯氏桿菌基因組尺度代謝網(wǎng)絡(luò)iYL1228的參考條件見表1。

表1 克雷伯氏桿菌基因組尺度代謝網(wǎng)絡(luò)iYL1228的參考條件Tab .1 Reference conditions of genome－scale metabolic network iYL1228for Klebsiella pneumonlae

1.2 算法

本實(shí)驗(yàn)使用的代謝流分析軟件COBRA Toolbox（Constraint－based Reconstruction and Analysis）［10－11］是運(yùn)行在MATLAB環(huán)境下的工具包，該工具包結(jié)合GUN 線性規(guī)劃工具（glpk）可對基因組尺度的代謝網(wǎng)絡(luò)進(jìn)行通量平衡分析與通量可變性分析。

1.2.1 通量平衡分析

通量平衡分析以生長速率最大為目標(biāo)函數(shù)，利用線性規(guī)劃找到穩(wěn)態(tài)條件下的通量分布。該線性規(guī)劃問題用公式描述為：

式中：A為化學(xué)計(jì)量矩陣，若一個系統(tǒng)有m種代謝物，r個反應(yīng)，則矩陣就有m行r列；反應(yīng)通量值v作為列向量；穩(wěn)態(tài)條件下，b為零向量；代表反應(yīng)i在參考條件下計(jì)算出的通量最大值與最小值。

1.2.2 通量可變性分析

通量可變性分析可將可行的解空間量化，在與通量平衡分析相同的約束條件下，通量可變性分析可以將每個反應(yīng)的最大值與最小值計(jì)算出來。該線性規(guī)劃問題用公式描述為：

其中，相對于參考條件的通量可變性用公式表達(dá)為：

有兩類反應(yīng)對通量可變性沒有影響，計(jì)算Δ時不需要計(jì)算，這兩類反應(yīng)是：內(nèi)部循環(huán)反應(yīng)、通量最大值與最小值相同的反應(yīng)。其中，后者又分為兩類：反應(yīng)通量值為零的反應(yīng)和不為零的反應(yīng)。

2 結(jié)果與討論

2.1 通量可變性的組成及計(jì)算

通量可變性來源于3部分：內(nèi)部反應(yīng)的通量可變性（Δint）、交換反應(yīng)的通量可變性（Δext）和生長速率的通量可變性（Δgro）。為計(jì)算這3部分?jǐn)?shù)值，分別用vint、vext與vgro表示內(nèi)部反應(yīng)的通量值、交換反應(yīng)的通量值與生長速率的通量值。給定條件下，分別表示最優(yōu)生長速率、最優(yōu)生長條件下的交換反應(yīng)速率。Δint就是在約束條件下的值（圖1中a 曲線）。僅固定，那么所計(jì)算的Δ（Δint＋ext）就是內(nèi)部反應(yīng)固有變化和外部反應(yīng)變化之和（圖1中b曲線）。Δint＋ext減去Δint就是外部反應(yīng)的變化（Δext）。將約束條件去除，所得到的Δtot就是3個原因?qū)е碌耐靠傋兓藭rΔtot＝Δgro＋Δint＋Δext（圖1中c曲線）。

圖1 通量可變性的組成Fig.1 Composition of flux variability

2.2 底物吸收速率對通量可變性的影響

外部環(huán)境的改變會影響通量可變性，進(jìn)而影響菌體的代謝靈活性。以下就底物（甘油、氧氣和氨）的吸收速率對通量可變性的影響進(jìn)行研究。

2.2.1 甘油吸收速率對通量可變性的影響

有氧條件與無氧條件下，甘油吸收速率對通量可變性的影響如圖2所示。

圖2 甘油吸收速率對通量可變性的影響Fig.2 Effect of glycerol uptake rate on flux variability

由圖2可看出：有氧條件下，甘油吸收速率為0.46 mmol·gDW－1·h－1時能夠滿足ATP所需的最小能量；無氧條件下，甘油吸收速率為5.44 mmol·gDW－1·h－1時才能滿足ATP所需的最小能量。

為了考察超出生理范圍甘油吸收速率（10.6 mmol·gDW－1·h－1）的代謝靈活性與生理范圍內(nèi)的規(guī)律是否一致，將甘油吸收速率增大到16 mmol·gDW－1·h－1，結(jié)果顯示，在甘油吸收速率生理范圍外，Δint與Δext對代謝靈活性的影響幾乎可以忽略不計(jì)，Δgro是影響通量可變性的唯一重要因素。

2.2.2 氧氣吸收速率對通量可變性的影響（圖3a）

圖3 氧氣吸收速率（a）和氨吸收速率（b）對通量可變性的影響Fig.3 Effect of oxygen uptake rate（a）and ammonia uptake rete（b）on flux variability

由圖3a可看出，細(xì)胞生長所帶來的可變性對代謝靈活性的影響最大，其它2個因素所帶來的影響可以忽略不計(jì)。

2.2.3 氨吸收速率對通量可變性的影響（圖3b）

由圖3b可看出：當(dāng)氨吸收速率大于5.6mmol·gDW－1·h－1時，Δgro依然是影響通量可變性的最主要因素；當(dāng)氨吸收速率小于5.6 mmol·gDW－1·h－1時，Δint與Δext對代謝靈活性也有重要的影響，其中Δint影響更大。

2.3 碳氮調(diào)節(jié)機(jī)制對通量可變性的影響

在大腸桿菌中，葡萄糖的吸收速率隨氮的浮動而被調(diào)控［12］，進(jìn)而影響通量可變性。為此，添加了碳氮調(diào)節(jié)機(jī)制，研究其對克雷伯氏桿菌通量可變性的影響。

圖4 碳氮調(diào)節(jié)機(jī)制下，氨吸收速率對通量可變性的影響Fig.4 Effect of ammonia uptake rate on flux variability under coordination of carbon and nitrogen utilization

由圖4可知，碳氮調(diào)節(jié)機(jī)制下，Δgro對通量可變性的影響很大，Δint與Δext幾乎為零，對通量可變性的影響可以忽略。因此，碳氮調(diào)節(jié)機(jī)制下，生長速率的變化仍然是影響通量可變性的唯一重要因素。

以上結(jié)果顯示，生長速率的變化是引起菌體代謝靈活性改變的主要原因。

2.4 生長速率對通量可變性的影響

由圖2可知，固定甘油吸收速率就能固定Δ和最大生長速率。因此固定菌體的生長速率，就可以得到Δ和的關(guān)系曲線（圖5中▼曲線），這里的是甘油吸收速率為10.6mmol·gDW－1·h－1時的最大生長速率，vgro為較小甘油吸收速率下的生長速率。甘油吸收速率能夠約束生物體的最大生長速率，因此能夠影響Δ。另一種約束生長速率的方法就是固定其通量值的上限和下限，而不對外部通量值加以限制。本實(shí)驗(yàn)在限制生長速率的范圍而沒有改變其可以得到的最大生長速率的情況下，考察生長速率對通量可變性的影響，結(jié)果如圖5所示。

圖5 生長速率對通量可變性的影響Fig.5 Effect of growth rate on flux variability

由圖5可看出：方形（■）中的vgro代表0），Δ值始終保持在1.0附近而沒有變化，說明通量可變性隨上限值確定的生長速率值變化不敏感；圓圈（○）中的vgro代表，即參考條件下獲得的最大生長速率）；實(shí)心圓（●）中生長速率的上下限值相同），這兩條曲線顯示出幾乎相同的變化趨勢，Δ值均隨著值的增大而減小。

由圖5還可看出：要使菌體獲得較大的代謝靈活性，首先要提供較大的甘油吸收速率，這樣菌體才會有大的生長速率；另外要將生長速率調(diào)控至比理想值小的值，即菌體在次優(yōu)模式下生長。這樣克雷伯氏桿菌就可以具備較大的代謝靈活性來應(yīng)對外界環(huán)境的改變。

3 結(jié)論

使用COBRA Toolbox軟件對基因組尺度的克雷伯氏桿菌代謝網(wǎng)絡(luò)進(jìn)行了通量平衡分析和通量可變性分析，將外界環(huán)境的改變對菌體代謝靈活性的影響量化，分別計(jì)算出通量可變性3個因素?cái)?shù)值。結(jié)果表明，生長速率的可變性是影響其代謝靈活性的唯一主要因素；并進(jìn)一步分析了生長速率如何影響代謝靈活性，得到兩個必要條件：（1）提供較大的甘油吸收速率供菌體生長；（2）調(diào)控菌體至次優(yōu)模式下生長。這兩點(diǎn)可使菌體獲得較大的代謝靈活性以應(yīng)對外界環(huán)境的改變。

［1］ORTH J D，THIELE I，PALSSON B O.What is flux balance analysis？［J］.Nature Biotechnology，2010，28（3）：245－248.

［2］SCHUSTER S，F(xiàn)ELL D A，DANDEKAR T.A general definition of metabolic pathways useful for systematic organization and analysis of complex metabolic networks［J］.Nature Biotechnology，2000，18（3）：326－332.

［3］SAN R M，CANCELA H，ACERENZA L.Source and regulation of flux variability inEscherichiacoli［J］.BMC Systems Biology，2014，8（1）：67.

［4］LEE S M，HONG W K，HEO S Y，et al.Enhancement of 1，3－propanediol production by expression of pyruvate decarboxylase and aldehyde dehydrogenase fromZymomonasmobilisin the acetolactate－synthase－deficient mutant ofKlebsiellapneumoniae［J］.Journal of Industrial Microbiology ＆ Biotechnology，2014，41（8）：1259－1266.

［5］曾宏，方柏山，汪排芬，等.甘油脫氫酶生產(chǎn)工藝仿真設(shè)計(jì)與技術(shù)經(jīng)濟(jì)分析［J］.化工學(xué)報(bào)，2013，64（6）：2169－2176.

［6］朱曉麗，南南，修志龍.以粗甘油為原料高產(chǎn)1，3－丙二醇菌株的等離子體復(fù)合誘變［J］.過程工程學(xué)報(bào)，2011，11（6）：1030－1037.

［7］陳琳，王領(lǐng)民，金平，等.甘油和木糖共發(fā)酵生產(chǎn)1，3－丙二醇的研究［J］.化學(xué)工程，2013，41（7）：5－8.

［8］徐衛(wèi)濤，鐘志輝，付水林，等.pH 值對克雷伯氏肺炎桿菌發(fā)酵甘油產(chǎn)1，3－丙二醇的影響與控制［J］.化學(xué)與生物工程，2011，28（8）：81－83.

［9］LIAO Y C，HUANG T W，CHEN F C，et al.An experimentally validated genome－scale metabolic reconstruction ofKlebsiellapneumoniaeMGH 78578，iYL1228［J］.Journal of Bacteriology，2011，193（7）：1710－1717.

［10］BECKER S A，F(xiàn)EIST A M，MO M L，et al.Quantitative prediction of cellular metabolism with constraint－based models：The COBRA Toolbox［J］.Nature Protocols，2007，2（3）：727－738.

［11］SCHELLENBERGER J，QUE R，F(xiàn)LEMING R M T，et al.Quantitative prediction of cellular metabolism with constraint－based models：The COBRA Toolbox V2.0［J］.Nature Protocols，2011，6（9）：1290－1307.

［12］DOUCETTE C D，SCHWAB D J，WINGREEN N S，et al.α－Ketoglutarate coordinates carbon and nitrogen utilizationviaenzyme I inhibition［J］.Nature Chemical Biology，2011，7（12）：894－901.