三羧酸循環(huán)Petri網(wǎng)模型的改進與狀態(tài)輸出算法

王 龍，李孝忠，李 勇，劉曉琴，趙 朋

（天津科技大學計算機科學與信息工程學院，天津 300222）

三羧酸循環(huán)Petri網(wǎng)模型的改進與狀態(tài)輸出算法

王 龍，李孝忠，李 勇，劉曉琴，趙 朋

（天津科技大學計算機科學與信息工程學院，天津 300222）

針對已建立的三羧酸循環(huán)混合Petri網(wǎng)模型無法獲取循環(huán)反應過程中各個階段模型狀態(tài)的缺陷，使用混合函數(shù)Petri網(wǎng)對原混合Petri網(wǎng)模型進行改進，根據(jù)改進后的模型提出狀態(tài)算法，并采用簡化模型驗證了模型和算法的有效性.結果表明：使用該算法能夠依順序獲取反應循環(huán)過程的前K個階段各個庫所中物質（即反應物）的量，對監(jiān)控、研究三羧酸循環(huán)等生化循環(huán)的反應過程提供了極大便利.

混合函數(shù)Petri網(wǎng)；三羧酸循環(huán)；狀態(tài)輸出算法

生物體內所有的化學反應總稱為代謝.代謝反應進而形成生物的代謝網(wǎng)絡，它起到為生物體提供能量、調控、構建等諸多重要的作用.因此，生物代謝的研究一直是生命科學研究的基礎領域［1］.三羧酸循環(huán)是生物代謝過程中十分重要的一環(huán)，是三種重要生物大分子物質代謝的最終共同途徑，也是這幾種營養(yǎng)物質在體內相互聯(lián)系的樞紐.對其進行建模分析能夠為生物體內大分子物質的代謝循環(huán)的研究提供極大便利.

Hofest?dt［2］和Reddy等［3］最早提出使用Petri網(wǎng)進行生物網(wǎng)絡的建模.Hofest?dt的研究主要側重于用基本Petri網(wǎng)描述代謝過程，并對Petri網(wǎng)（Petri net，PN）的活性、可達性、P不變量、T不變量等結構特性進行了生物解釋，這為后續(xù)研究奠定了基礎. Reddy等則側重于分析，提出了對P不變量、T不變量、有界性等的定性分析.Matsuno等［4］介紹了一種擴展的基本Petri網(wǎng)——混合函數(shù)Petri網(wǎng)（hybrid functional Petri net，HFPN），以糖酵解過程為例進行HFPN建模.該模型的特性十分適于生物過程建模. Hardy等［5］提出了使用連續(xù)Petri網(wǎng)（continuous Petri net，CPN）對代謝過程進行建模的方法，使用不變量分析對由CPN模型模擬所獲取的數(shù)據(jù)進行分析，連續(xù)Petri網(wǎng)的連續(xù)概念同混合Petri網(wǎng)中離散元與連續(xù)元的應用有共通處.

文獻［6］用混合Petri網(wǎng)（hybrid Petri net，HPN）對三羧酸循環(huán)過程進行了建模，通過仿真能夠獲取反應過程的相關數(shù)據(jù).但是，該建模仿真方法存在無法獲取反應過程中各個時間點（或者說反應階段）的反應物狀態(tài)的缺點.而按需要提取各個階段反應物的消耗量、剩余量等信息，對代謝反應的研究有著重要的意義.為了彌補這一缺陷，本文在混合Petri網(wǎng)模型的離散變遷中加入兩個函數(shù)，將模型改進為混合函數(shù)Petri網(wǎng)模型.并根據(jù)循環(huán)過程中各個反應狀態(tài)的數(shù)據(jù)需求，提出一種優(yōu)化的狀態(tài)輸出算法.

1 建模背景及改進

1.1 三羧酸循環(huán)

三羧酸循環(huán)是需氧生物體內普遍存在的代謝途徑，又稱為檸檬酸循環(huán).三羧酸循環(huán)是生物代謝中一項十分重要的過程，它是三種重要生物大分子物質代謝的最終共同途徑，也是這幾種營養(yǎng)物質在體內相互聯(lián)系的樞紐［1］.它的發(fā)現(xiàn)是代謝化學最重要的成就之一，然而，代謝網(wǎng)絡相關研究大多集中在糖酵解和尿素循環(huán)等過程中，與它們擁有同等重要性的三羧酸循環(huán)的建模研究相對較少.因此，本文以三羧酸循環(huán)為對象進行建模.

1.2 三羧酸循環(huán)HPN模型改進

從HPN定義Σ=（P，T，Pre，Post，h，M0）［7］可知，混合Petri網(wǎng)由庫所、變遷、入射弧、出射弧、混合函數(shù)及初始標識六個元組成.僅靠這六元必然無法實現(xiàn)對三羧酸循環(huán)分階段輸出各個時間點反應物狀態(tài)的要求.

針對獲取代謝循環(huán)系統(tǒng)各狀態(tài)實時數(shù)據(jù)的需求，可以對該六元組進行相應改進.首先，為了輸出各階段狀態(tài)，就需要一個元組來保存各個階段的狀態(tài)值，添加庫所狀態(tài)集M，存儲待輸出的反應物狀態(tài).其次，由于是輸出的狀態(tài)以時間點進行劃分，故添加時間元tP（Petri time）［8］，表示一個基本時間單位，該時間單位的具體值可根據(jù)實際需要賦值；有了基本的時間單位之后，還需要添加抽樣間隔IS（interval of sampling），表示抽樣輸出的間隔所包含的時間單位的個數(shù)，該值同樣可按實驗需求自行設定；為了將反應過程分段，加入延遲元D，記錄每個變遷的發(fā)生延遲.除此之外，還要為每個變遷加入兩個函數(shù)：vI和vD［9］，分別表示該變遷如果發(fā)生將會使得相連的庫所中物質增加和減少的量，該函數(shù)還用于判斷變遷發(fā)生的可能性.在添加了增量、減量函數(shù)以及延遲元之后，需要一個集合保存每個變遷的發(fā)生條件，因此加入發(fā)生條件元C，由布爾函數(shù)來區(qū)分兩種發(fā)生條件.

使用混合函數(shù)Petri網(wǎng)［10］，依照上述對HPN的改進，就形成了三羧酸循環(huán)的混合函數(shù)Petri網(wǎng)模型.

2 混合函數(shù)Petri網(wǎng)

2.1 定義

為實現(xiàn)算法，經(jīng)過改進，將HPN擴展成為十三元組混合函數(shù)Petri網(wǎng)，定義為

其中：

P為庫所的非空有限集合，P=｛P1，P2，…，Pn｝；

T為變遷的非空有限集合，T=｛T1，T2，…，Tm｝；

Pre為輸入關聯(lián)映射，代表PN圖中由庫所指向變遷的有向弧的集合；

Post為輸出關聯(lián)映射，代表PN圖中由變遷指向庫所的有向弧的集合；

h為混合函數(shù)，表示系統(tǒng)中結點是離散的（D）還是連續(xù)的（C），系統(tǒng)中包含離散變遷T，D、連續(xù)變遷T，C、離散庫所PD、連續(xù)庫所PC；

M為模型的狀態(tài)集合，M=｛M1，M2，…，Mi｝，其中每個元素Mi表示在階段i時，系統(tǒng)中每個庫所的標識或者物質量，Mi=（m1（i），m2（i），…，mn（i））；

M0為系統(tǒng)初始狀態(tài)，表示系統(tǒng)初始狀態(tài)下各個庫所的值；

D為變遷延遲的有限集合，D=｛d1，d2，…，dm｝，表示相應的變遷在發(fā)生前必須等待的延遲時間，只有在延遲為零時，變遷才有可能發(fā)生；

C為變遷發(fā)生條件的有限集合，C=｛c1，c2，…，cm｝，由布爾函數(shù)給出值cm=｛0，1｝，該元控制相應變遷的發(fā)生條件，控制機制在后文介紹；

tP為系統(tǒng)的基本時間單元，由于生化反應極其迅速，常規(guī)的時間單位表述太過繁冗，故在系統(tǒng)中統(tǒng)一使用tP，表示一個基本時間單元；

IS為取樣間隔，表示兩步計算步驟的時間間隔，由tP的個數(shù)表示.當取樣間隔最小時，IS=tP；

vD為消耗速率的有限非負實數(shù)集合，vD=｛f1，f2，…，fx｝，其中x=（Pβ，tj），β=1，2，…，n，j=1，2，…，m，x是由庫所Pβ指向變遷tj的弧，vD表示當變遷tj發(fā)生時，庫所Pβ中消耗的物質的量；

vI為生成速率的有限非負實數(shù)集合，vI=｛g1，g2，…，gy｝，其中y=（tj，Pβ），j=1，2，…，m，β=1，2，…，n，y是由變遷tj指向庫所Pβ的弧，vI表示當變遷tj發(fā)生時，庫所Pβ中因生成而增加的物質的量.

2.2 變遷發(fā)生的控制機制

生化反應歸根結底是變化的發(fā)生，體現(xiàn)在PN系統(tǒng)中就是變遷，因此變遷是非常重要的一環(huán).在以上HFPN系統(tǒng)十三元中，T、Pre、Post、D、C、vD、vI七個元均與變遷有關，其中T、Pre、Post控制變遷在系統(tǒng)中的位置、結構，D、C控制變遷發(fā)生的條件，vD、vI控制變遷發(fā)生后對系統(tǒng)的改變.為了便于之后運用算法輸出HFPN系統(tǒng)各個階段的狀態(tài)，在這里先對變遷發(fā)生的控制機制進行討論.

從HFPN定義可以看出，延遲D及發(fā)生條件C對變遷的發(fā)生進行控制，根據(jù)兩者不同的值，系統(tǒng)中變遷的狀態(tài)可分為3種，在此引入函數(shù)lj代表變遷的狀態(tài)值.

（1）當dj＞0或cj=0時，lj=0.此時相應變遷等待延遲歸零或發(fā)生條件為真；

（2）當dj=0且cj=1時，lj=1.此時相應變遷立即發(fā)生；

（3）當dj＞0且cj=1時，lj=2.此時相應變遷已可以發(fā)生，只需等待延遲歸零.

2.3 系統(tǒng)模型改進

根據(jù)HFPN系統(tǒng)的定義，對原有的三羧酸循環(huán)HPN系統(tǒng)進行相應改進，系統(tǒng)模型如圖1所示.

圖1 三羧酸循環(huán)HFPN模型Fig. 1 The HFPN model of TCA cycle

3 HFPN系統(tǒng)狀態(tài)輸出算法

為了實現(xiàn)對循環(huán)系統(tǒng)狀態(tài)的監(jiān)控，并實時獲取各階段的反應數(shù)據(jù)，在已建立的三羧酸循環(huán)HFPN模型的基礎上，提出如下混合函數(shù)Petri網(wǎng)系統(tǒng)狀態(tài)輸出的優(yōu)化算法，使用本算法能夠獲取PN系統(tǒng)前K個狀態(tài)中各個庫所包含的物質量（即反應參與物的剩余量）.

3.1 算法

1. 聲明數(shù)組T′［z］，賦初值p=0，q=0

3.2 算法內容

該狀態(tài)輸出算法可分為三部分.首先對整個模型初始化，為后續(xù)計算做準備.由于模型中擁有離散和連續(xù)兩種變遷，而這兩種變遷在計算中算法并不相同，故首先對模型的變遷遍歷，將離散變遷與連續(xù)變遷區(qū)分開并按順序放入新的數(shù)組T′中.這一操作使變遷編號發(fā)生改變，所以還應對入射弧集合Pre和出射弧集合Post進行相應的調整.然后將變遷函數(shù)lj初始化.

為了獲取模型前K個階段的狀態(tài)M，接下來對1至K進行for循環(huán)計算.mβ（i+1）∶=mβ（i）表示由于還未進行新的點火（fire），每次循環(huán)開始各庫所的狀態(tài)值應與上一個階段相同.第一部分中已將變遷進行整合排序，這里根據(jù)變遷種類將整個循環(huán)分為兩個次級for循環(huán).在已建立的三羧酸循環(huán)系統(tǒng)中，規(guī)定連續(xù)變遷不消耗庫所中的托肯值.因此對連續(xù)庫所進行for循環(huán)，令其狀態(tài)值恒等于初始狀態(tài).

接著對離散變遷計算.循環(huán)開始討論發(fā)生條件cj（i）的值：如果與變遷j相連的前集庫所中物質量大于該變遷將要消耗的值（即fj），則cj（i）=1，反之cj（i）=0.這一操作避免了狀態(tài)值出現(xiàn)負數(shù)的情況，因為這在生化反應中是不可能發(fā)生的.若此時變遷狀態(tài)函數(shù)lj（i）為零且cj（i）=1，表明變遷的延遲dj仍大于零，變遷需持續(xù)等待.已知第i階段變遷延遲dj（M［i］）（其中M［i］表示系統(tǒng)在i階段的狀態(tài)集合）和采樣間隔IS，可知延遲過后系統(tǒng)處于第i+｛dj（M［i］）/IS｝-1階段，這樣便可根據(jù)發(fā)生條件cj（i）的值，對該階段的lj函數(shù)賦值（cj=1，lj=1；cj=0，lj= 0）.而在等待延遲過程中各階段（i到i+｛dj（M［i］）/IS｝-2）的lj值應為2.另一方面，如果第i階段lj的值為1，說明此時延遲已歸零，變遷tj已經(jīng)能夠發(fā)生，檢查發(fā)生條件是否為真，若cj（i）=1，則變遷發(fā)生.變遷tj發(fā)生自然對與其相連的庫所的狀態(tài)值產(chǎn)生增減變化，為了輸出正確的狀態(tài)值，算法接著對M［i+1］進行相應更改.依據(jù)定義中變遷tj的消耗速率集vD中的數(shù)據(jù)fx，將tj的前驅庫所的狀態(tài)值減去fx，即mβ（i+1）=mβ（i+1）-fx（M［i］）.同樣，依照變遷tj的生成速率集vI中的數(shù)值gy，將tj的后驅庫所的狀態(tài)值增加gy，即mβ（i+1）=mβ（i+1）+ gy（M［i］）.至此，算法核心部分，即計算第i階段各庫所的狀態(tài)值已經(jīng)完成.為了對第i+1階段循環(huán)繼續(xù)進行計算，還需要根據(jù)i階段系統(tǒng)的發(fā)生對模型狀態(tài)的改變，對變遷狀態(tài)函數(shù)lj（i+1）賦值：若此時延遲dj不為零，則lj（i+1）=0，反之lj（i+1）=1.

在主體for循環(huán)的最后將每階段的計算結果，即第i階段過后所有庫所的狀態(tài)值M［i+1］輸出，表示為M［i+1］=（m1（i+1），…，mn（i+1））.同時該數(shù)據(jù)也作為初始狀態(tài)進入下一輪循環(huán).

4 算法應用

為便于檢驗算法的正確性，根據(jù)三羧酸循環(huán)HFPN模型，取其中核心的八步循環(huán)反應構建如圖2所示簡化模型.其中

圖2 三羧酸循環(huán)簡化模型Fig. 2 Simplified model of Krebs cycle

由于三羧酸循環(huán)涉及的反應物量的準確數(shù)據(jù)缺乏，且這里僅為驗證算法，故將反應參與物的量設定為易于計算的整數(shù).

由算法可知，離散庫所中物質的量不會因離散變遷的發(fā)生而改變，因此僅取模型中連續(xù)庫所與連續(xù)變遷部分進行算法的檢驗.模型中混合函數(shù)h皆為C.

此外，為節(jié)約計算時間，將采樣區(qū)間設為最小，即IS=tP；將所有變遷的發(fā)生條件C的初值都設為1，這樣可使系統(tǒng)盡快進入計算階段.

將該模型數(shù)據(jù)帶入狀態(tài)輸出算法，經(jīng)過計算，輸出的模型前8個階段狀態(tài)如下：

5 結 語

前期工作中針對生物代謝網(wǎng)絡的特性提出了新生物代謝網(wǎng)絡建模步驟，并據(jù)此對三羧酸循環(huán)建立了混合Petri網(wǎng)模型.為了實現(xiàn)對代謝循環(huán)反應過程的監(jiān)控，獲取生化反應各個階段中各個反應參與物的狀態(tài)，本文針對這一需求，對原HPN模型進行了改進，添加數(shù)個函數(shù)及元，將原有的六元HPN擴展稱為十三元混合函數(shù)Petri網(wǎng)模型，并提出相應的HFPN系統(tǒng)狀態(tài)輸出算法.代入簡化模型可知算法是可行的，通過生化反應方程式的計算驗證了算法的正確性，實現(xiàn)了輸出反應各階段物質狀態(tài)的目標.

獲取代謝循環(huán)反應過程中各階段的反應物狀態(tài)數(shù)據(jù)，能夠為生化反應的實時監(jiān)控、反應調控、結果預測等提供極大便利.同時，使用計算機對循環(huán)反應過程進行建模以及仿真分析，相對于傳統(tǒng)生物實驗對于實驗環(huán)境、實驗器材、實驗時長的嚴格要求，節(jié)省了大量的物質成本，而且計算機對反應的快速仿真與計算也節(jié)約了可觀的時間成本.

文中的HFPN系統(tǒng)狀態(tài)輸出算法雖然是基于三羧酸循環(huán)模型提出的，但由于生化反應方程的相似性，該算法在其他生物反應網(wǎng)絡、生物代謝循環(huán)的建模分析中同樣可以發(fā)揮作用，但也存在連續(xù)性數(shù)據(jù)輸出表述方面的不足，以后可在此方面進一步研究.

［1］ 查錫良. 生物化學［M］. 2版. 上海：復旦大學出版社，2011：109-217.

［2］ Hofest?dt R A. Petri net application to model metabolic processes［J］. Systems Analysis Modeling Simulation，1994，16（2）：113-122.

［3］ Reddy V N，Liebman M N，Mavrovouniotis M L. Qualitative analysis of biochemical reaction systems［J］. Computers in Biology and Medicine，1996，26（1）：9-24.

［4］ Matsuno H，F(xiàn)ujita S，Doi A，et al. Towards biopathway modeling and simulation［J］. Lecture Notes in Computer Science，2003，2679：3-22.

［5］ Hardy S，Robillard P N. Petri net-based method for the analysis of the dynamics of signal propagation in signaling pathways［J］. Bioinformatics，2008，24（2）：209-217.

［6］ Li X Z，Wang L，Liu Y，et al. The modeling of the Krebs cycle based on hybrid Petri net［C］//Fang X W. Proceeding of the 11th International Conference on Natu-ral Computation（ICNC'15）. Changsha：University of Hunan Press，2015：1187-1191.

［7］ 吳哲輝. Petri網(wǎng)導論［M］. 北京：機械工業(yè)出版社，2006：181-209.

［8］ Castellini A，F(xiàn)ranco G，Manca V. Hybrid functional Petri nets as MP systems［J］. Natural Computing，2010，9（1）：61-81.

［9］ 丁德武，李文澤. 代謝網(wǎng)絡結構與功能的Petri網(wǎng)分析［J］. 計算機工程與應用，2012，48（28）：42-46.

［10］ Baldan P，Cocco N，Marin A，et al. Petri nets for modeling metabolic pathways：A survey［J］. Natural Computing，2010，9（4）：955-989.

責任編輯：常濤

State Algorithm and Improvement of the Petri Net Model of the Krebs Cycle

WANG Long，LI Xiaozhong，LI Yong，LIU Xiaoqin，ZHAO Peng
（College of Computer Science and Information Engineering，Tianjin University of Science & Technology，Tianjin 300222，China）

The established hybrid Petri net model of Krebs cycle in our preliminary work has a defect，so it is unable to get the model state of each period in the cyclic process.Aiming at remedying it，the original model has been improved with the hybrid functional Petri net.According to the improved model，a state algorithm is put forward and verified by a simplified model.Thus，the capacity of each place，i.e the quality of the reactants，in the first K phases of the cyclic reaction can be obtained and put out in sequence，which helps a lot in monitoring and studying the course of reaction of the Krebs cycle.

hybrid functional Petri net；the Krebs cycle；state algorithm

TP399

A

1672-6510（2016）05-0063-06

10.13364/j.issn.1672-6510.20150210

2015-11-16；

2016-02-04

國家自然科學基金資助項目（61070021，11301382）

王 龍（1990-），男，河南鄭州人，碩士研究生；通信作者：李孝忠，教授，lixz@tust.edu.cn.

數(shù)字出版日期：2016-07-11；數(shù)字出版網(wǎng)址：http://www.cnki.net/kcms/detail/12.1355.N.20160711.1612.008.html.