易奇志 杜 焰 周天壽
1)(江西師范大學數學與信息科學學院,南昌 330022)
2)(江西師范大學化學與化工學院,南昌 330022)
3)(中山大學數學與計算科學學院,廣州 510275)
(2012年11月29日收到;2013年2月6日收到修改稿)
細胞生存于復雜的環(huán)境中,能夠感應來自其他細胞的生物信號,與周圍的細胞或環(huán)境形成各種相互作用.對于由多個相互作用的細胞組成的群體系統(tǒng),其群體行為并不僅僅依賴于單個細胞水平的動力學,而更多的是依賴于細胞之間相互作用的方式與效果[1-5].事實上,細胞之間的相互作用可以有效地調節(jié)和影響細胞的群體行為,如引起完全相同的細胞分化[6];又如P19細胞允許聚集時可以有效地分化成骨骼的肌肉(這是因為細胞聚集形成的相互作用可以調節(jié)轉錄因子myoD誘導的骨骼肌漿蛋白)[7].在細胞間有通訊的多細胞系統(tǒng)中,細胞之間的信息交換與整合會引起細胞群體協作響應,導致多細胞系統(tǒng)展示出各種有趣的現象,如生物節(jié)律[8],細胞斑圖[9-15],以及同步、聚類、多穩(wěn)態(tài)共存、多節(jié)律性等[14-21].
多細胞現象(表現為群體行為),除了與細胞之間的通訊方式、信號整合構件(又稱為順式調控構建)有關外,還可能與相互作用的細胞數目有關.例如,在生物體發(fā)育過程中,細胞的分裂會引起細胞數目的增加;在活的有機體內,細菌群體的增殖與擴張會引起細菌數目的增加.反過來,細胞數目的增加會引起細胞之間通信概率的增加,從而有可能影響和改變原有多細胞系統(tǒng)的群體動力學行為.例如,Kaneko等[3]模擬了一群有生化網絡相互作用的細胞,包含細胞分裂過程,發(fā)現當細胞數目到達某個閾值時,分裂導致幾乎完全相同的細胞同步振動.隨著數目的進一步增加,同步振動逐漸消失,導致細胞分成幾個具有不同振動相位的子類,從而導致不同振幅的聚類振動.他們由此提出了細胞分化的同構異素多元化理論(isologous diversification).又如,Nakajima等[1]設計了一個核心的基因調控模塊,并設計了細胞之間的3種不同方式的相互作用,分別對應于3個不同模型.他們展示出在這3種作用方式下細胞狀態(tài)是如何隨細胞總數的變化而分化.更細化地,他們顯示出:對第一種作用方式,細胞命運由細胞總數決定,即隨著總數的增加,細胞群體從1聚類切換到2聚類,再切換到另一個1聚類;第二種作用方式能夠使得細胞狀態(tài)出現多元化,即細胞數目增加到適當范圍時,1聚類與2聚類可以共存.這一理論模型所描述的現象與發(fā)育生物學上所觀察的現象[22,23]基本一致.此外,Koseska[20]考察了一個由群體感應機制耦合的壓制振動子組成的多細胞系統(tǒng),發(fā)現系統(tǒng)規(guī)模的增大不僅使異質平衡態(tài)的穩(wěn)定性區(qū)域擴大,而且使其吸引域也擴大.他由此提出了多細胞群體通過聚類進行協作分化的觀點.總之,系統(tǒng)規(guī)模的大小對于細胞采用何種策略取得群體行為起著非常關鍵的作用.研究細胞數目對于多細胞系統(tǒng)群體動力學行為的影響不僅具有重要的實際意義,而且有助于我們理解種群在擴張情形時的協作行為.
在文獻[21]中,我們研究了細胞通訊在由群體感應機制耦合的組合振子的多細胞系統(tǒng)中的作用,發(fā)現細胞通訊能夠誘導耦合系統(tǒng)的多穩(wěn)性和多節(jié)律性.盡管這一重要的定性結論并不受振子數目的影響,但鑒于細胞數目的增加可能會導致細胞群體動力學發(fā)生各種定性的或定量的改變,因此在前文工作的基礎上,這里我們將進一步調查細胞數目的改變是如何影響這個多細胞系統(tǒng)的群體行為的,特別關注細胞數目(或系統(tǒng)規(guī)模)的增加如何影響聚類(clustering)行為.這里的工作能夠看成是前一工作的延伸,但獲得了更多有趣且有意義的結果.
為理解起見,這里我們簡單地解釋聚類.聚類是一種穩(wěn)定的動力學狀態(tài),它可能包含若干個子群體,每個子群體具有相同或幾乎相同的行為.聚類在由大量相互作用單元的生物、化學、物理和社會系統(tǒng)中普遍存在,并吸引了眾多研究者的廣泛關注[1-3,14,15,19,20,24-27].在生物學上,聚類與細胞分化密切相關.多細胞生物系統(tǒng)的細胞分化機制極其復雜,一種普遍的觀點是細胞根據成形素濃度梯度設定其位置信息進行分化[28].然而,即使在相同的環(huán)境條件下,相同的細胞仍可以展示不同的表型(phenotype),從而由相同細胞組成的群體也可呈現出某種異質性(heterogeneity),主要表現在可能存在幾個子類,且每一個子類中的細胞表現出有組織的群體行為.正如文獻[2,20]所指出的,動力學聚類為細胞分化提供了另一種可能的機制或解釋.
本文主要關注系統(tǒng)規(guī)模對于組合振子多細胞系統(tǒng)兩種聚類行為——平衡態(tài)聚類和振動聚類的影響.利用分岔分析和數值模擬,我們發(fā)現:細胞數目的增加可以改變異質平衡態(tài)的穩(wěn)定性區(qū)間的大小,并誘導新的聚類行為;對于適當的細胞密度,細胞數目的增加有利于擴大異質平衡態(tài)的吸引域,特別是,當細胞數目充分大時,異質平衡態(tài)占主導地位,幾乎以概率1出現,表明細胞分化可能與細胞數目有密切關系;細胞數目的增加極大地豐富了平衡態(tài)聚類和振動聚類的表現形式和共存方式,為生物體的適應性和多功能性提供了良好的基礎.
這里,我們仍采用我們最近的工作[21]中的合成基因調控網絡模型。在這一模型中,首先,單個振子是由壓制振動子和松弛振子通過整合而成的組合振子 (見圖 1),其中 U(TetR),V(CI)和 W(LacI)分別是三個基因u(tetR),v(cI)和w(lacI)的產物,自體誘導物(一種小的信號分子)W的動力學可以通過在相應方程中引入一個小參數而減慢;其次,振子之間通過群體感應機制耦合而構成多細胞系統(tǒng),其中信號分子可以擴散到群體中每一個細胞的內部,并調控目標基因的表達(關于模型的詳細描述,見文獻[21]).類似于文獻[21]的簡化,我們合并轉錄、翻譯過程為單步過程,并采用擬平衡態(tài)假設,來建立數學模型,由此得到如下多細胞系統(tǒng)模型:
圖1 組合振子的示意圖
在下面的分析中,我們將主要從平衡態(tài)聚類(steady state clustering)、振動聚類(oscillatory clustering)兩個角度來刻畫系統(tǒng)規(guī)模N對多細胞系統(tǒng)聚類行為的影響.我們指出:平衡態(tài)聚類實際上是系統(tǒng)形成異質平衡態(tài)(inhomogeneous steady state,IHSS),更確切地說,在所有細胞中,對應成分(如W等)取兩種不同的常數值.通常,N個耦合振子在平衡態(tài)聚類中有N-1種不同的分布.特別地,當考慮耦合系統(tǒng)中的某個成分時,例如,Wi(i=1,2,···,N),假設m個W處于高的濃度狀態(tài)(記為mU),而其余的N-m個W 處于低的濃度狀態(tài)(記為(N-m)L),其中m可以取1到N-1之間的值,那么我們說系統(tǒng)具有異質平衡態(tài),并把這種異質平衡態(tài)模式記為mU|(N-m)L.與異質平衡態(tài)相對的是同質平衡態(tài)(homogeneous steady state,HSS),是指所有細胞中對應成分取相同的常數值.振動聚類是指系統(tǒng)中的所有單元分成若干組,每一組內的單元都有完全相同的狀態(tài),即以相同的振幅完全同步地振動.對于振動聚類,聚類模式更為復雜,可以按振幅和相位的不同分成不同的類(類的數目可以超過2,我們將在后面給出詳細的解釋).這些動力學聚類是理解細胞分化的重要基礎.在下面,我們將采用分岔分析和數值模擬來展示和解釋我們的結果,并用到縮寫符號:PB代表叉形分岔,LP代表極限點分岔或鞍結分岔,HB代表Hopf分岔.
本小節(jié)考察細胞數目的增加對于平衡態(tài)聚類的影響.為此,我們對耦合系統(tǒng)的規(guī)模做一個微小的改變,把系統(tǒng)從包含2個相同細胞擴充為包含3個或4個相同細胞.以細胞密度Q作為分岔參數,作W1關于Q的分岔圖.為了便于比較,我們只展示與穩(wěn)定的異質平衡態(tài)(IHSS,即平衡態(tài)聚類)有關的分岔圖,見圖2(a)—(c).對于N=2的情形,當Q介于LP1和LP2之間時,耦合系統(tǒng)處于穩(wěn)定的IHSS狀態(tài).此時,IHSS狀態(tài)只有一種可能的分布,即1U|1L或兩個振子中W的濃度一個處于高狀態(tài),一個處于低狀態(tài)(見圖2(a)中的粗虛線).對于N=3的情形,存在兩種不同的穩(wěn)定聚類模式(見圖2(b)中的粗虛線):2U|1L位于LP1與PB1之間(Q∈[0.1414,0.3679]),1U|2L位于PB2與LP2之間(Q∈[0.159,0.4671]).注意到這兩種不同模式的聚類位于不同的分岔分支上,因此IHSS的穩(wěn)定性區(qū)間為這兩種模式的存在區(qū)間的并集Q∈[0.1414,0.4671],相比N=2的情形,其IHSS的穩(wěn)定性區(qū)間Q∈[0.1493,0.4119](見圖2(a)中的粗虛線),有明顯的擴大(≈24.03%).當N=4時,有三種不同模式的穩(wěn)定聚類(見圖2(c)中的粗虛線),從左到右依次是:3U|1L位于LP1與LP2之間(Q∈[0.1378,0.35]),2U|2L位于PB1與PB2之間(Q∈[0.1489,0.4115]),及1U|3L位于LP3與HB1之間(Q∈[0.1624,0.446]).此時,IHSS的穩(wěn)定性區(qū)間為Q∈[0.1378,0.446],它的長度比N=3時的長度稍微小一點(≈5.37%),但比N=2時的長度要大.然而,對于N=4的情形,有一種新奇的現象發(fā)生在聚類狀態(tài)1U|3L所在的分支上.具體說來,發(fā)生在HB1處的Hopf分岔引起了一個簡單的周期1振動(見圖2(c)中的粗點線),它存在于HB1(Q=0.446)到LP4(Q=0.4807)之間的區(qū)間.此時,4個振子聚成兩類:1個振子處于高的濃度狀態(tài),以較小的振幅振動;3個振子處于低的濃度狀態(tài),以較大的振幅同步振動.
注意,對于N≥3的情形,有N-1種不同模式的平衡態(tài)聚類,這些聚類模式在一定的范圍內可以共存,系統(tǒng)最終演化到哪種IHSS狀態(tài),取決于初始條件.例如,在N=4時,三種IHSS狀態(tài):3U|1L,2U|2L和1U|3L可以在LP3和LP2之間共存(見圖2(c)).這三種分布的蛋白質濃度水平顯示了輕微的差別,這種差別可以看成是細胞群體的一種典型適應性,意味著當環(huán)境條件發(fā)生變化時,細胞群體容易調節(jié)其狀態(tài)分布,以便最優(yōu)地適應環(huán)境.
圖2 描繪穩(wěn)定的異質平衡態(tài)(IHSS)的分岔圖,細胞數目分別為:(a)N=2;(b)N=3;(c)N=4.這里,粗虛線代表各種穩(wěn)定模式的IHSS,細實線代表不穩(wěn)定的平衡態(tài),(c)中的粗點線代表簡單周期1的穩(wěn)定振動,圖中4個振子按成分W的濃度聚成兩類:1個振子在高狀態(tài)以小振幅振動,3個完全同步的振子處于低狀態(tài),以較大的振幅振動(這個分岔圖并沒有顯示出所有的分岔信息,只顯示了與當前討論相關的部分.符號的縮寫見正文)
由上可見,細胞數目小的改變也能引起耦合系統(tǒng)聚類行為的定量或定性變化,表現為細胞數目的變化可以調節(jié)平衡態(tài)聚類的穩(wěn)定性區(qū)間的大小,還可以導致新的聚類現象即振動聚類的發(fā)生.
現在,我們對系統(tǒng)規(guī)模做更大幅度的增加,將系統(tǒng)擴充為中等規(guī)模,即系統(tǒng)由20個相同的細胞組成.我們將探討和比較N=2和N=20這兩種情況下平衡聚類狀態(tài)的吸引域的異同.為此,我們讓細胞密度于[0,1]之間變化,其余系統(tǒng)參數值如前所述設置.對于每一種細胞數目N,隨機產生三個蛋白質U,V,W的初始條件,使之服從區(qū)間[0,2]上的均勻分布,通過直接的數值計算,可得到200組耦合振子所有成分的時間序列,然后根據每種情況下的200組時間序列數據統(tǒng)計各種穩(wěn)定吸引子的數目:同質平衡態(tài)(HSS),異質平衡態(tài)(IHSS)及振動(oscillation).通過這一方法,我們可以探測每一個有統(tǒng)計意義的穩(wěn)定吸引子,進一步得到穩(wěn)定性區(qū)域關于Q的依賴性(見圖3).此圖表明:系統(tǒng)規(guī)模的增加擴大了異質平衡態(tài)的穩(wěn)定性區(qū)間,并且明顯地擴大了異質平衡態(tài)的吸引域,但吸引域的擴大以同質平衡態(tài)吸引域的縮小為代價(比較圖3(a)和(b)中異質平衡態(tài)區(qū)域的邊界).
圖3 本圖表明系統(tǒng)規(guī)模對于動力學區(qū)域分布的影響(當細胞密度在[0,1]區(qū)間內變化時,對于N=2(a)和N=20(b)中的每一種情況都有三種動力學區(qū)域:同質平衡態(tài)(HSS),異質平衡態(tài)(IHSS)和振動(oscillation),但是在前兩種情況下吸引域的邊界有些不同,表現為:系統(tǒng)規(guī)模的增加擴大了異質平衡態(tài)的穩(wěn)定性區(qū)間,并且擴大了異質平衡態(tài)的吸引域)
下一步,我們對于N=2時HSS與IHSS共存的區(qū)間上分別任選一個Q的值(見圖3(a)),例如,取Q=0.25,0.4,對于這兩個Q值,我們分別研究系統(tǒng)規(guī)模的增加對于HSS和IHSS等動力學狀態(tài)的分布有何影響.為此,對于每一種細胞數目N,我們采用前面的方法,隨機產生200組蛋白質濃度的初始條件進行計算,可得到200組耦合振子所有成分的時間序列,根據這200組時間序列,可以計算出穩(wěn)定吸引子HSS和IHSS發(fā)生的次數.這樣,我們可以探測出每一種N對應的具有統(tǒng)計意義的穩(wěn)定吸引子.
圖4詳細地表明:對于兩個給定的細胞密度(Q=0.25和Q=0.40),系統(tǒng)規(guī)模對動力學狀態(tài)分布的效果.兩種情況擁有某些共同的特征:當細胞數目少量增加時,通過隨機初始條件到達IHSS的概率明顯單調地增加,且是以減少HSS的概率為代價的;當細胞數目繼續(xù)增加時,IHSS的概率有增大的趨勢,但不是單調地增加而是有所漲落;當細胞數目增加到充分大時,IHSS完全占主導地位,HSS幾乎不出現.兩者的區(qū)別在于:隨著細胞數目的增加,異質平衡態(tài)趨于概率1的速度稍有差別.具體說來,Q=0.25的情形要比Q=0.40的情形趨于1的速度快一點.這些結果說明細胞數目的增加有利于擴大平衡態(tài)聚類的吸引域,這也意味著細胞分化可能與細胞數目的大小有密切關系.
圖4 當細胞密度Q被固定為 (a)Q=0.25;(b)Q=0.40時,細胞數目的增加對動力學狀態(tài)(同質平衡態(tài)(HSS)和異質平衡態(tài)(IHSS))的分布產生的效果(這兩種情況具有某些共同的特征:當細胞數目少量增加時,通過隨機初始條件到達IHSS的概率明顯單調地增加,但是以減少HSS的概率為代價的;當細胞數目繼續(xù)增加時,IHSS的概率有增大的趨勢,但會出現漲落;當細胞數目增加到充分大時,IHSS完全占主導地位,HSS幾乎不出現)
通過數值計算,我們不難發(fā)現細胞數目的增加可以導致更豐富的平衡態(tài)聚類和振動聚類的模式及其共存方式.現以N=2和N=20的情況進行比較來說明.圖5展示了20個耦合振子隨細胞密度變化而出現的平衡態(tài)聚類和振動聚類.首先,隨著細胞數目(N)的增加,平衡態(tài)聚類的模式(N-1種)增加,且不同的IHSS共存的方式也更加多樣化.對于N=20的情形,當Q=0.36時,四種不同的IHSS狀態(tài)可以共存(見圖5(a)—(d)).其次,振動聚類共存的方式也更多樣化.對于N=20的情形,我們觀察到幾種周期振動的共存,如Q=0.9時反相周期2振動和不對稱周期2振動的共存(見圖5(f),(g)),以及Q=0.82時不對稱1111混合模式振動和不對稱1211混合模式振動的共存(見圖5(i),(j)),這些共存的方式在N=2的情況下對于所有的Q∈[0,1]都沒有觀察到;同時,也觀察到了與N=2類似的幾種振動聚類模式,如Q=0.43時的同相周期1振動,Q=0.76時的對稱12111213MMO,Q=0.51時的擬周期陣發(fā)振動以及Q=0.6時的混沌陣發(fā)振動(分別見圖5(e),(h),(k)和(l)).此外,系統(tǒng)規(guī)模的適當增加也使得異質平衡態(tài)可以和周期解共存(見圖3(b)),這種情況在N=2時不發(fā)生.
事實上,通過對振動聚類的仔細研究,我們發(fā)現,除了圖5(f)—(k)中展示的各種周期和擬周期振動聚類外,還有很多種其他的振動聚類模式,這是由于耦合系統(tǒng)中的各個振子可以按照振幅和相位使振動聚成不同的類.振幅不同可以導致聚類,振幅相同而相位有差異也可以導致聚類,從而對于N≥3的系統(tǒng),在振動中可以形成超過兩類的多個子類,這與平衡態(tài)聚類只有兩個子類不同.此外,振動方式的差異(如周期1振動、周期2振動、混合模式振動等方式)和細胞在不同類之間的分布情況使得振動聚類的表現形式極為豐富,從而也使得不同模式的聚類的共存方式更加多樣化.下面我們以N=20,Q=0.48的情況為例進行說明.
圖5 20個耦合的組合振子出現的平衡態(tài)聚類及振動聚類:(a)—(d)Q=0.36時有四種具有不同分布的異質平衡態(tài)共存:(a)10U|10L;(b)9U|11L;(c)8U|12L;(d)7U|13L;(e)—(j)多模式的振動聚類:(e)Q=0.43時的同相周期1振動,(f),(g)Q=0.9時的反相周期2振動和不對稱周期2振動,(h)—(j)三種類型的混合模式振動:(h)Q=0.76時的對稱12111213MMO,(i)Q=0.82時的不對稱1111MMO,(j)Q=0.82時的不對稱1211MMO;(k),(l)兩種模式的陣發(fā)振動,其中(k)對應于Q=0.51時的擬周期陣發(fā)振動,(l)對應于Q=0.6時的混沌陣發(fā)振動
圖6 當Q=0.48時,20個振子形成三種不同模式的不對稱周期1的聚類振動,在上下兩類中細胞數目的分布分別為(a)3:17;(b)4:16;(c)5:15
圖7 當Q=0.48時,20個振子形成三種不同模式的不對稱周期2的聚類振動,在第1類和第2類中細胞數目的分布分別為(a)7:13;(b)9:11;(c)8:12
圖6展示了20個振子以三種不同的不對稱周期1的方式振動,此時振子聚成兩類,上下兩類之間細胞數目的分布分別是 3:17,4:16 和 5:15. 圖 7則展示了不同于圖6的其他形式的振動聚類,盡管20個振子分別聚成兩類,兩類之間的細胞數目之比為 7:13,9:11 和 8:12,但是振動是以不對稱周期 2的方式進行(在圖7的每個子圖中兩類細胞的振幅有細微的區(qū)別).圖8展示了四種更為復雜的振動聚類,此時20個振子聚成三類,細胞在第1類、第2類和第 3 類之間的分布分別為 9:10:1,8:10:2,8:11:1和9:9:2.在這里,前三種情形中第1類和第2類的振幅有細微的差別,但是第4種情形中第1類和第2類的振幅相等但相位不同.注意:對于N=20,Q=0.48的系統(tǒng),按照3.2節(jié)所述的方法隨機產生初始濃度進行200次的計算,我們發(fā)現除了上述展示的振動聚類之外,還有其他的振動聚類模式存在,但是經統(tǒng)計,圖6中的三種振動聚類模式發(fā)生的頻率最高,共出現了90次,其中又以圖6(b)中的方式出現最為頻繁(200次中出現44次),圖8中的四種振動聚類的發(fā)生也較為頻繁.圖6至圖8中的這些狀態(tài)和其他未顯示的狀態(tài)在N=20,Q=0.48的系統(tǒng)中共存,系統(tǒng)最終演化到哪個狀態(tài)取決于初始條件,某種狀態(tài)出現頻率的高低與這種狀態(tài)吸引域的大小有關.
系統(tǒng)規(guī)模的增加大大豐富了平衡態(tài)聚類與振動聚類的表現形式,為細胞的分化和多功能性提供了良好的基礎,多種穩(wěn)定的動力學狀態(tài)的共存有效地提高了生物體適應環(huán)境的能力.
圖8 當Q=0.48時,20個振子形成四種不同模式的復雜的振動聚類,在第1類、第2類和第3類中細胞數目的分布分別為(a)9:10:1;(b)8:10:2;(c)8:11:1;(d)9:9:2
對于一個由群體感應機制耦合而成的組合振子多細胞系統(tǒng),我們已經考察了細胞數目的增加對于系統(tǒng)動力學行為的影響,聚焦于細胞數目對平衡態(tài)聚類和振動聚類的影響.通過細致的分析,我們發(fā)現:首先,細胞數目即使少量的增加也會引起異質平衡態(tài)的穩(wěn)定性區(qū)間的大小的明顯改變,并誘導出新的聚類現象;其次,細胞數目的增加有利于擴大異質平衡態(tài)的吸引域,且當細胞數目增加到充分大時,IHSS占主導地位,同質平衡態(tài)幾乎不出現,表明細胞分化依賴于系統(tǒng)規(guī)模的大小;最后,細胞數目的增加極大地豐富了平衡態(tài)聚類和振動聚類的表現形式和共存方式,為細胞分化和生物體的多功能性提供了基礎,也提高了系統(tǒng)對環(huán)境的適應性.
動力學聚類是細胞內部動力學與細胞之間相互作用的共同結果,為細胞分化提供了一種新的機制.如何把這種機制和成形素梯度的細胞分化機制結合起來考慮,對于理解多細胞生物體的發(fā)育可能是非常重要的.保持細胞之間的相互作用,僅由于細胞的增殖引起系統(tǒng)規(guī)模的增加,通過動力學聚類導致細胞分化的可能性提高,這種機制是容易達到的,因為它對成形素沒有要求,對于基因調控網絡的參數也不需要精細的調節(jié),從進化論的觀點來看,這是非常有用的.然而,在更真實的情形,何種機制導致細胞分化還不清楚.
最后指出:我們的模型在幾個方面進行了簡化.首先,沒有考慮詳細的生物過程,把轉錄、翻譯、啟動子結合等過程合并成了一步,并且沒有考慮組合調控的效果[14,15].其次,忽略了一些生物因素,如細胞的多樣性,基因表達過程中固有的隨機漲落,信息傳輸過程中的時間滯后,轉錄因子協作性的效果等等.這些方面可能對于規(guī)模增長的群體動力學行為產生重要的影響[20,34-39],值得進一步的探究.
[1]Nakajima A,KanekoK 2008 J.Theor.Biol.253 779
[2]KanenoK,YomoT 1994 Physica D 75 89
[3]KanekoK,YomoT 1997 B.Math.Biol.59 139
[4]Wang J W,Chen A M,Zhang J J,Yuan Z J,Zhou T S 2009 Chin.Phys.B 18 1294
[5]Wang B H,Lu Q S,Lv S J,Lang X F 2009 Chin.Phys.B 18 872
[6]Greenwald I,Rubin G M 1992 Cell68 271
[7]Armour C,Garson K,McBurney M W 1999 Exp.CellRes.251 79
[8]Garcia-OjalvoJ,Elowitz M B,Strogatz S H 2004 Proc.Natl.Acad.Sci.U.S.A.101 10955
[9]Blair S S 2007 Annu.Rev.CellDev.Biol.23 293
[10]Fields R D,Burnstock G 2006 Nat.Rev.Neurosci.7 423
[11]Scherrer R,ShullV 1986 Can.J.Microbiol.32 607
[12]Ben-Jacob E,Cohen I,Shochet O,Tenenbaum A,Czirok A,Vicsek T 1995 Phys.Rev.Lett.75 2899
[13]Basu S,Gerchman Y,Collins C H,Arnold F H,Weiss R 2005 Nature 434 1130
[14]Zhang J J,Yuan Z J,Zhou T S 2009 Phys.Rev.E 79 041903
[15]YiQ Z,Zhang J J,Yuan Z J,Zhou T S 2010 Eur.Phys.J.B 75 365
[16]McMillen D,KopellN,Hasty J,Collins J 2002 Proc.Natl.Acad.Sci.U.S.A.99 679
[17]Ullner E,Koseska A,Kurths J,Volkov E,Kantz H,Garca-OjalvoJ 2008 Phys.Rev.E 78 031904
[18]Ullner E,Zaikin A,Volkov E I,Garca-OjalvoJ 2007 Phys.Rev.Lett.99 148103
[19]Koseska A,Volkov E,Kurths J 2007 Phys.Rev.E 75 031916
[20]Koseska A,Ullner E,Volkov E,Kurths J,Garca-OjalvoJ 2010 J.Theor.Biol.263 189
[21]YiQ Z,Zhou T S 2011 Phys.Rev.E 83 051907
[22]Gurdon J,Lemaire P,KatoK 1993 Cell75 831
[23]Shimuta K,NakajoN,UtoK,HayanoY,OkazakiK,Sagata N 2002 EMBO J.21 3694
[24]Okuda K 1993 Physica D(Amsterdam)63 424
[25]Golomb D,HanselD,Shraiman B,Sompolinsky H 1992 Phys.Rev.A 45 3516
[26]Taylor A F,Kapetanopoulos P,Whitaker B J,Toth R,BullL,Tinsley M R 2008 Phys.Rev.Lett.100 214101
[27]KuramotoY 1984 Chemicaloscillations,waves and turbulence(Berlin:Springer-Verlag)
[28]Tabata T,TakeiY 2004 Development 131 703
[29]Zhou T S,Zhang J J,Yuan Z J,Chen L N 2008 Chaos 18 037126
[30]http://www.math.pitt.edu/bard/xpp/xpp.html.
[31]Ermentrout B 2002 Simulating,analyzing and animating dynamicalsystems:A guide toXppaut for researchers and students(software,environment and tools)1st ed.(Philadephia,PA:SIAM Press)
[32]Krupa M,Popovic N,KopellN,Rotstein H G 2008 Chaos 18 015106
[33]Izhikevich E M 2000 Int.J.Bifurcat.Chaos 10 1171
[34]Zhou T S,Chen L N,Aihara K 2005 Phys.Rev.Lett.95 178103
[35]Koseska A,Volkov E,Kurths J 2009 EuroPhys.Lett.85 28002
[36]Koseska A,Volkov E,Kurths J 2010 Chaos 20 023132
[37]Smolen P,Baxter D A,Byrne J H 2002 Biophys.J.83 2349
[38]Song H,Smolen P,Av-Ron E,Douglas D A,Byrne J H 2007 Biophys.J.92 3407
[39]Potapov I,Volkov E,Kuznetsov A 2011 Phys.Rev.E 83 031901