葉緯明 呂彬彬 趙琛 狄增如

1)(北京師范大學管理學院系統(tǒng)科學系，北京 100875)

2)(北京有色金屬研究總院稀土材料國家工程研究中心，有研稀土新材料股份有限公司，北京 100088)

(2012年7月18日收到;2012年7月31日收到修改稿)

1 引言

隨著生物技術(shù)的不斷發(fā)展，基因工程被不斷向前推進，海量的數(shù)據(jù)為科學家探索基因的奧秘提供了有力的支持，大量生物組織的基因被成功測序.對基因調(diào)控網(wǎng)絡(genetic regulatory network，GRN)的研究也越來越受到各個領(lǐng)域的關(guān)注.其中，基因調(diào)控網(wǎng)絡的振蕩問題一直是大家關(guān)注的熱點之一[1-13].早在30年前，Thomas就曾提出：基因調(diào)控網(wǎng)絡要想形成振蕩，就必需至少要有一個負反饋環(huán)(環(huán)上負反饋的個數(shù)為奇數(shù)).該論斷于1995年被嚴格證明是正確的[14,15].近年來，該領(lǐng)域涌現(xiàn)出大量的科研成果.Goh等[16]認為基因調(diào)控網(wǎng)絡是否能夠形成振蕩主要取決于網(wǎng)絡中正負反饋環(huán)之間的競爭，如果負反饋環(huán)占優(yōu)，網(wǎng)絡呈現(xiàn)振蕩態(tài);如果正反饋環(huán)占優(yōu)，網(wǎng)絡則會呈現(xiàn)穩(wěn)定態(tài).2000年，Elowitz和Leibler等[17]成功地在實驗室中制造出了首個振蕩的3節(jié)點生物網(wǎng)絡.它包括3個轉(zhuǎn)錄抑制蛋白，每個抑制蛋白都抑制下一個基因的轉(zhuǎn)錄，當某種特殊蛋白質(zhì)含量發(fā)生變化時，細胞能在發(fā)光狀態(tài)和非發(fā)光狀態(tài)之間轉(zhuǎn)換，形成一個環(huán)狀回路，起到有機振蕩器的作用.2002年，U.Alon發(fā)現(xiàn)在真實的生物網(wǎng)絡中存在一系列非常小的結(jié)構(gòu)單元(被命名為“模體”)，它們在網(wǎng)絡中出現(xiàn)的頻率很大，并且對網(wǎng)絡的功能起著非常大的作用.至此，大家開始把目光轉(zhuǎn)移到這些小的結(jié)構(gòu)單元上來，希望通過對這些模體的研究來逐漸認識整個基因網(wǎng)絡.2006年，Nochomovitz和Li[18]運用布爾代數(shù)的方法研究了所有3節(jié)點和4節(jié)點的網(wǎng)絡，發(fā)現(xiàn)這種少節(jié)點的網(wǎng)絡中有許多是可以形成振蕩的.

目前我們發(fā)現(xiàn)的大多數(shù)基因調(diào)控網(wǎng)絡在大部分時間里都是處于靜息態(tài)的，只有當網(wǎng)絡受到外界刺激時，它們才會被激發(fā)到振蕩態(tài)，完成各種生物功能，然后再回到靜息態(tài)等待下一次刺激.對于一個原本應該處于靜息態(tài)的基因調(diào)控網(wǎng)絡，如果因為某種原因跳到了振蕩態(tài)上，那么它就會不斷地重復某種生物過程，導致機體功能紊亂，同時它對外界的刺激也可能不再敏感，使生物體處于嚴重失控狀態(tài).因此對基因調(diào)控網(wǎng)絡的控制就顯得非常重要了.

目前對生物網(wǎng)絡的控制研究主要集中在神經(jīng)網(wǎng)絡[19,20]、蛋白質(zhì)相互作用網(wǎng)絡，對基因調(diào)控網(wǎng)絡方面的控制還不是很多，特別是如何將一個振蕩的基因調(diào)控網(wǎng)絡控制到靜息態(tài)的研究還非常少[21].本文運用多相位超前驅(qū)動方法對基因調(diào)控網(wǎng)絡的控制進行了研究.通過大量的數(shù)值模擬，發(fā)現(xiàn)對于參數(shù)已經(jīng)確定的基因調(diào)控網(wǎng)絡，此方法的控制效率可以達到95%以上(10節(jié)點網(wǎng)絡)，對于參數(shù)不確定的網(wǎng)絡，此方法的控制效率也非常高.

2 基因調(diào)控網(wǎng)絡動力學模型

這里我們運用目前廣泛采用的描述基因聯(lián)合調(diào)控的動力學模型[21-25].它在抑制調(diào)控和促進調(diào)控中采用“AND”門的形式.具體的動力學模型如下：

這里xi表示第i個節(jié)點的濃度，acti(repi)代表促進(抑制)轉(zhuǎn)錄因子的總和.動力學(1)中一共有5個系統(tǒng)參量.μ表示在沒有促進作用和抑制作用時的漏轉(zhuǎn)錄率;K是網(wǎng)絡中促進子(抑制子)調(diào)控時的半高濃度;h是希爾系數(shù);γ代表衰減率;表示第j個節(jié)點對第i個節(jié)點促進(抑制)作用的權(quán)重.

動力學(1)是一種強抑制模型[26]，對所有的抑制調(diào)控采用“AND”門的形式[25]，而對所有的促進調(diào)控采用“OR”門的形式.為了使系統(tǒng)具有一般性，我們讓網(wǎng)絡中所有節(jié)點的參數(shù)相同.(μi=μ=0，γi=0.1，hi=2，Ki=0.3，αi=βi=1，i=1,2,···,N[27]).

3 數(shù)值模擬

3.1 主相位超前驅(qū)動方法

主相位超前驅(qū)動(dominant phase advanced driving，DPAD)思想最先是由Hu等[19,20]提出，并用于對神經(jīng)網(wǎng)絡的控制中.其主體思想是：對于一個隨機的網(wǎng)絡，網(wǎng)絡當中每一個節(jié)點都是不能自發(fā)產(chǎn)生 振蕩的，它們?nèi)绻胄纬烧袷幘捅仨毷艿絼e的節(jié)點的驅(qū)動，當然這里的驅(qū)動可以來自一個節(jié)點也可以來自多個節(jié)點，但必需至少要有一個來自外界的驅(qū)動.可以通過一系列的定義找到網(wǎng)絡中每個節(jié)點的驅(qū)動點(作為這個節(jié)點的上游節(jié)點).再通過一定的定義即可以從每個節(jié)點的所有上游節(jié)點中找出對其影響最大的節(jié)點，該節(jié)點的驅(qū)動即為主超前相位驅(qū)動(其他上游節(jié)點的驅(qū)動此時可以去掉).這樣網(wǎng)絡中的每個節(jié)點都將有且只有一個主超前相位驅(qū)動，即每個節(jié)點都只剩下一個上游節(jié)點，這樣整個網(wǎng)絡必定會變成一個一維環(huán)狀結(jié)構(gòu).運用這種方法找到的環(huán)狀結(jié)構(gòu)對于神經(jīng)網(wǎng)絡的調(diào)控具有非常重要的意義.而將這種方法應用到基因調(diào)控網(wǎng)絡中，找到的環(huán)狀結(jié)構(gòu)也意義非凡[21].如圖1所示，圖1(a)為一個10節(jié)點的隨機振蕩基因調(diào)控網(wǎng)絡，圖1(b)為運用主相位超前驅(qū)動方法找到的環(huán)狀結(jié)構(gòu).從這種環(huán)狀結(jié)構(gòu)中不僅可以清楚地找到振蕩基因調(diào)控網(wǎng)絡的振蕩源，還可以看到網(wǎng)絡信號在網(wǎng)絡中的傳播路徑，同時還能粗略地看出網(wǎng)絡中重要的和不重要的節(jié)點(分別用紅色和白色表示)，有利于對網(wǎng)絡進行調(diào)控.

圖1 少節(jié)點基因調(diào)控網(wǎng)絡 (a)10節(jié)點隨機振蕩基因調(diào)控網(wǎng)絡，其中綠色實線表示促進作用，紅色虛線表示抑制作用;(b)運用主相位超前驅(qū)動方法找到的一維環(huán)狀結(jié)構(gòu)，其中重要節(jié)點用紅色表示，不重要的節(jié)點用白色表示

圖2 不同節(jié)點數(shù)的基因調(diào)控網(wǎng)絡，主相位超前驅(qū)動環(huán)上的節(jié)點重要性

這里，網(wǎng)絡中重要的節(jié)點指的是去掉這個節(jié)點網(wǎng)絡就會從振蕩態(tài)到達穩(wěn)定態(tài);不重要的節(jié)點指的是去掉這個節(jié)點，網(wǎng)絡的振蕩態(tài)不會發(fā)生改變.運用主相位超前驅(qū)動方法對原始網(wǎng)絡進行簡化，可以得到：對于少節(jié)點的網(wǎng)絡，重要的節(jié)點幾乎都在環(huán)上，不重要的節(jié)點幾乎都在環(huán)外.但是隨著網(wǎng)絡節(jié)點數(shù)的增多，主相位超前驅(qū)動方法在區(qū)分節(jié)點的重要性上就會弱化.環(huán)上的節(jié)點不再重要的概率會逐漸增高.如圖2所示，經(jīng)過大量的數(shù)值模擬可以看出，當網(wǎng)絡節(jié)點數(shù)少于10個時，環(huán)上的節(jié)點是重要節(jié)點的概率超過了90%;而當網(wǎng)絡上的節(jié)點多于10個時，環(huán)上節(jié)點是重要節(jié)點的概率近似以線性比例下降;當網(wǎng)絡節(jié)點數(shù)達到40個時，環(huán)上節(jié)點為重要節(jié)點的概率下降到只有57.9%.

3.2 多相位超前驅(qū)動方法

我們對主相位超前驅(qū)動方法做了適當?shù)男薷模瑢⒅優(yōu)槎嘞辔怀膀?qū)動方法.運用這種方法可以將找到重要節(jié)點的概率大大提高.具體的做法如下：1)給網(wǎng)絡中每個節(jié)點定義一個相位[21];2)找出每個節(jié)點的上游節(jié)點(將上游節(jié)點指向下游節(jié)點的連邊定義為重要連邊);3)簡化原始網(wǎng)絡，將網(wǎng)絡中的重要連邊留下，組成相位超前驅(qū)動網(wǎng)絡;4)找出相位超前驅(qū)動網(wǎng)絡中所有的負反饋環(huán);5)統(tǒng)計相位超前驅(qū)動網(wǎng)絡中每個節(jié)點參與的負反饋環(huán)數(shù)并排序.如圖3所示，圖3(a)為圖1(a)的相位超前驅(qū)動網(wǎng)絡，圖3(b)為其相位超前驅(qū)動網(wǎng)絡中所有的負反饋環(huán)(共11個).各個節(jié)點的數(shù)據(jù)統(tǒng)計見表1.從表1可以看出：在相位超前驅(qū)動網(wǎng)絡中，10號節(jié)點參與的負反饋環(huán)數(shù)最多，它是一個重要節(jié)點，控制住它，就可以控制整個網(wǎng)絡.同時還可以看出節(jié)點參與負反饋環(huán)數(shù)排名中的前6位都是重要節(jié)點.而后4位都是不重要的節(jié)點.可見這一方法對于判斷節(jié)點的重要性具有很好的效果.

圖3 相位超前驅(qū)動網(wǎng)絡 (a)圖1(a)的相位超前驅(qū)動網(wǎng)絡;(b)(a)圖中所有的負反饋環(huán)

數(shù)值模擬結(jié)果表明：隨機選取10000個振蕩網(wǎng)絡，運用多相位超前驅(qū)動方法對節(jié)點參與負反饋環(huán)的排序，排在第一位的節(jié)點(以下稱為“首要節(jié)點”)是重要節(jié)點的概率比用主相位超前驅(qū)動方法尋找重要節(jié)點的概率要高很多.如圖4所示，可以看出運用多相位超前驅(qū)動方法找到的首要節(jié)點是重要點的概率并沒有隨著網(wǎng)絡節(jié)點數(shù)的增加而大幅度下降.當網(wǎng)絡節(jié)點數(shù)達到40個節(jié)點時，首要節(jié)點為重要節(jié)點的概率依然可以達到90%.

表1 多相位超前驅(qū)動法判斷節(jié)點的重要性(參數(shù)確定)

圖4 不同節(jié)點數(shù)的基因調(diào)控網(wǎng)絡，運用多相位超前驅(qū)動方法找到的首要節(jié)點為重要節(jié)點的概率

以上是在系統(tǒng)參數(shù)確定的情況下，運用多相位超前驅(qū)動方法尋找網(wǎng)絡重要節(jié)點的結(jié)果.那么當系統(tǒng)的參數(shù)不確定時，如何來尋找網(wǎng)絡的重要節(jié)點呢?數(shù)值模擬表明運用多相位超前驅(qū)動方法依然可以有效地找到網(wǎng)絡中的重要節(jié)點.具體方法是：將系統(tǒng)參數(shù)在其各自有效的范圍內(nèi)[27]劃分為若干個點，然后在每個參數(shù)點上統(tǒng)計出各個節(jié)點在相位超前驅(qū)動網(wǎng)絡上參與的負反饋環(huán)的個數(shù)，最后統(tǒng)計在所有參數(shù)點上各個節(jié)點在相位超前驅(qū)動網(wǎng)絡上參與的負反饋環(huán)的總個數(shù)，并依此排序.

還是以圖1(a)所示的網(wǎng)絡為例.通過數(shù)值模擬，我們發(fā)現(xiàn)：運用多相位超前驅(qū)動方法找到的在所有參數(shù)點上的相位超前驅(qū)動網(wǎng)絡中參與負反饋環(huán)最多的節(jié)點(10號節(jié)點，以下稱為“頭號節(jié)點”)，其平均重要性達到了100%(這里的平均重要性指的是：對所有的參數(shù)點[27]，控制目標節(jié)點后，網(wǎng)絡能夠被有效控制到穩(wěn)定態(tài)的平均概率)，其他節(jié)點的平均重要性見表2.從表2可以看出：網(wǎng)絡節(jié)點的平均重要性與其參與的負反饋環(huán)的總排名關(guān)系很大.

表2 多相位超前驅(qū)動法判斷節(jié)點的重要性(參數(shù)不確定)

我們用數(shù)值模擬的方法對10000個10節(jié)點的網(wǎng)絡進行了統(tǒng)計，結(jié)果表明：在系統(tǒng)參數(shù)不確定的情況下，運用多相位超前驅(qū)動方法找到的頭號節(jié)點，其平均重要性為100%的概率為92%.對于其他節(jié)點數(shù)的網(wǎng)絡，我們也做了類似的統(tǒng)計，結(jié)果如圖5所示.

圖5 不同節(jié)點數(shù)的基因調(diào)控網(wǎng)絡，頭號節(jié)點平均重要性為100%的概率

我們認為多相位超前驅(qū)動方法之所以比主相位超前驅(qū)動方法能更好地尋找到網(wǎng)絡當中的重要點，是因為在主相位超前驅(qū)動中只保留了對被驅(qū)動節(jié)點貢獻最大的上游節(jié)點，而忽視了其他上游節(jié)點的作用.這樣就導致了網(wǎng)絡簡化過程中只保留下來惟一一個負反饋環(huán)，并被認為是網(wǎng)絡振蕩的源.從數(shù)值統(tǒng)計上看，這惟一一個負反饋環(huán)上的節(jié)點是重要節(jié)點的概率比環(huán)外的節(jié)點要高很多.但這其中還是有節(jié)點是不重要的，它們出現(xiàn)在環(huán)上也許只是因為其所控制的下游節(jié)點只有它這一個上游節(jié)點，或者有多個上游節(jié)點但是每個上游節(jié)點的貢獻都很小.而運用多相位超前驅(qū)動方法就可以很好地避免上面的問題.在多相位超前驅(qū)動中，我們將原始網(wǎng)絡中所有的重要連邊都保留下來，形成相位超前驅(qū)動網(wǎng)絡，該網(wǎng)絡是對原始網(wǎng)絡的簡化，同時又保留著原始網(wǎng)絡中的重要信息.根據(jù)Goh等[16]的理論，以及我們之前的一些研究結(jié)果表明：負反饋環(huán)對于網(wǎng)絡的振蕩起著促進作用(當然這其中各種負反饋環(huán)所起的作用大小不同，但是定性上來說還是促進的[28])，那么簡化后的網(wǎng)絡中，每個負反饋環(huán)上的節(jié)點應該都對網(wǎng)絡的振蕩起著一定的作用.如果把這些作用累加起來，就可以看作是該節(jié)點對網(wǎng)絡振蕩所起的總貢獻，這樣總貢獻的排序即可以在一定程度上反映網(wǎng)絡中節(jié)點的重要性.

4 結(jié)論

基因調(diào)控網(wǎng)絡是目前大家研究的一個熱點，而它的自組織振蕩又是其中的核心問題.如何來控制自組織振蕩是問題的關(guān)鍵.主相位超前驅(qū)動方法對于少節(jié)點的網(wǎng)絡可以較好地找到網(wǎng)絡中的重要節(jié)點，但是對于節(jié)點數(shù)較多的網(wǎng)絡，這種方法找到重要節(jié)點的概率會呈線性下降.運用多相位超前驅(qū)動方法則可以很好地解決上述問題.對于系統(tǒng)參數(shù)已經(jīng)給定的網(wǎng)絡，運用該方法找到重要節(jié)點的成功率明顯提高.對于系統(tǒng)參數(shù)不確定的網(wǎng)絡，運用該方法也可以高效率地找到網(wǎng)絡中的重要節(jié)點.

[1］Tavazoie S，Hughes J D，Campbell M J，Cho G M 1999Nat.Genet.22 281

[2］Xu P X，Zhang X，Heaney S，Yoon A，Michelson A M，Maas R L 1999Development126 383

[3］Peterson K J，Davidson E H 2000P.Natl.Acad.Sci.USA97 4430

[4］ Lee T I，Rinaldi N J，Robert F，Odom D T 2002Science298 799

[5］Ihmels J，F(xiàn)riedlander G，Bergmann S，Sarig O 2002Nat.Genet.31 370

[6］Ma H W，Kumar B，Ditges U，Gunzer F 2004Nucleic Acids Res.32 6643

[7］Posfai G，Plunkett G 2006Science312 1044

[8］Gilchrist M，Thorsson V，Li B，Rust A G 2006Nature441 173

[9］Kurland C G，Collins L J，Penny D 2006Science312 1011

[10］Mcadams H H，Srinivasan B，Arkin A P 2004Nature Rev.Genet.5 1

[11］Shen-Orr S S，Milo R，Mangan S，Alon U 2002Nat.Genet.31 64

[12］Lipshtat A，Purushothaman S P，Lyengar R，Maayan A 2008Biophysical Journal94 2566

[13］Kim J R，Yoon Y，Cho K H 2008Biophysical Journal94 359

[14］Snoussi E H 1998J.Biol.Sys.6 3

[15］Plahte E H，Mestl T，Omholt W S 1995J.Biol.Sys.3 409

[16］Goh K G，Kahng B，Cho K H 2008Biophysical Journal94 4270

卓別林有一句很幽默的名言：“人們?yōu)槲覛g呼，是因為他們對我了如指掌；人們?yōu)閻垡蛩固箽g呼，是因為沒人弄得懂他?！钡拇_，電影不同于高深的物理理論，卓別林用最直接的方式向不同層次的觀眾們傳達著自己想要表達的觀點，這是人們喜歡他的原因。

[17］Elowitz M B，Leibler S 2000Nature403 335

[18］Nochomovitz Y D，Li H 2006P.Natl.Acad.Sci.USA103 4180

[19］Qian Y，Liao X，Huang X，Hu G 2010Phys.Rev.E 81 036101

[20］Liao X，Hu G 2011Front.Phys.6 124

[21］Ye W，Huang X，Li P，Hu G 2010Phys.Lett.A 374 2521

[22］Li C，Chen L，Aihara L S，Kazuyuki L 2006IEEE53 2451

[23］Mangan S，Alon U 2001P.Natl.Acad.Sci.USA100 11980

[24］Setty Y，Mayo A E，Surette M G，Alon U 2003P.Natl.Acad.Sci.USA100 7702

[25］Ishilara S，F(xiàn)ujimoto K，Shibata T 2005Genes to Cells10 1025

[26］Wu Y，Zhang X，Yu J，Ouyang Q 2009PLoS Computational Biology5 e1000442

[27］Ye W，Li P，Huang X，Hu G 2010Phys.Lett.A 374 4749

[28］Ye W，Zhang Z，L¨u B，Di Z，Hu G 2012Chin.Phys.B 21 060203