楊紅麗 劉楠 楊聯(lián)貴
(內(nèi)蒙古大學(xué)數(shù)學(xué)科學(xué)學(xué)院, 呼和浩特 010021)
轉(zhuǎn)錄因子p53是細胞應(yīng)激網(wǎng)絡(luò)的核心, 以動態(tài)響應(yīng)的方式控制基因毒性壓力下的細胞命運抉擇.Mdm2是種E3泛素連接酶, 既能破壞p53的穩(wěn)定性又能提高p53的生成效率.Mdm2對p53的抑制性功能在p53-Mdm2振子中扮演著構(gòu)建性角色, 而Mdm2對p53的促進性功能如何調(diào)控這個基因網(wǎng)絡(luò)的動力學(xué)仍缺少研究.因此, 本文利用數(shù)學(xué)模型, 全面探究了Mdm2上調(diào)p53的這條通路對p53動力學(xué)的影響.結(jié)果表明:Mdm2在Ser395位點的磷酸化作用對p53的振蕩必不可少; 之前報道的磷酸酶Wip1被p53振蕩所需要, 可能僅發(fā)生在Mdm2所介導(dǎo)的正反饋通路強度較高的情景下; Mdm2促進p53的失活以及泛素化降解也是p53反復(fù)振動動力學(xué)發(fā)生的關(guān)鍵因素, 與以往的結(jié)論一致.本文的結(jié)果可對今后p53動力學(xué)的相關(guān)實驗起到一定的指導(dǎo)作用.
癌癥是種細胞信號通路的疾病, 嚴重威脅著人類生命健康[1].抑癌因子p53作為細胞信號轉(zhuǎn)導(dǎo)網(wǎng)絡(luò)的樞紐, Puszyński等[2]斷言在人類癌癥中約有一半會出現(xiàn)p53自身的突變, 而另一半會出現(xiàn)p53調(diào)控子的突變.在靜息細胞中, E3泛素連接酶Mdm2 (murine double minute 2) 使得p53多泛素化, 進而被蛋白酶體識別并快速降解; 此時p53的濃度被抑制在一個基底水平[3].而在缺氧、紫外線照射、依托泊苷處理、電離輻射等細胞壓力下, p53被穩(wěn)定并激活[4?7].活化的p53與靶基因的綁定親和力增強, 誘導(dǎo)多種下游蛋白的表達.這些下游蛋白參與到細胞周期停滯、DNA修復(fù)、細胞程序性死亡等生理過程中[8].據(jù)實驗報道, 細胞命運的信息被p53動力學(xué)編碼[9].因此, 還有一些問題亟需解決.例如, 解碼p53動力學(xué)所關(guān)聯(lián)的細胞命運抉擇的細致機制是個極具應(yīng)用價值的問題;研究p53蛋白信號網(wǎng)絡(luò)中本身的動力學(xué)控制也是個很有意義的問題, 等等.本文主要關(guān)注的是有關(guān)p53蛋白信號網(wǎng)絡(luò)中的動力學(xué)控制問題, 即Mdm2在Ser395位點磷酸化這條信號通路對p53動力學(xué)調(diào)控的問題.
2000年, Lev Bar-Or等[10]做了一項關(guān)于p53動力學(xué)的實驗, 發(fā)現(xiàn)電離輻射下細胞群體中p53的衰減式脈沖, 也就是阻尼振蕩.后來的實驗在單細胞中觀測, 發(fā)現(xiàn)p53的振蕩是無阻尼的, 并且脈沖的個數(shù)與輻射劑量正相關(guān)[11,12].后續(xù)的很多理論模型圍繞p53的反復(fù)式脈沖動力學(xué)開展, 例如, 張麗娟等[13]用p53-Mdm2負反饋環(huán)附加時滯構(gòu)建的振子模型揭示了群體細胞p53的阻尼振蕩可能是單細胞p53振蕩持續(xù)時間不一、脈沖數(shù)量差異的結(jié)果;Bottani和Grammaticos[14]研究了帶時滯的p53-Mdm2振子在p53, Mdm2半衰期參數(shù)控制下的分岔現(xiàn)象; Ma等[15]建立了復(fù)雜的混合模型 (隨機過程對接微分方程) 來模擬p53的數(shù)字式脈沖;Xia和Jia[16]重建了包含p53振子模塊和DNA修復(fù)模塊的簡化模型, 用更簡單的方式解釋了p53脈沖數(shù)正比于輻射劑量的可行性機制;等[17]探索了單細胞中Wip1 (wild type p53-induced phosphatase 1) 參與下p53發(fā)生振蕩的常微分模型和基于反應(yīng)擴散的偏微分模型; 畢遠宏等[18]研究了PDCD5 (programmed cell death 5) 和Mdm2最大生成速率調(diào)控下的p53全局動力學(xué)與穩(wěn)定性等問題, 展示了p53動力系統(tǒng)中存在的豐富動力學(xué)現(xiàn)象; Hat等[19]研究了p53基因網(wǎng)絡(luò)中常見的正負反饋通路及反饋背后的雙穩(wěn)態(tài)、振蕩等分岔機制; 王道光等[20]提取了p53-Mdm2負反饋環(huán)和
p53-Wip1-ATM (ataxia telangiectasia-mutated)
負反饋環(huán)共存的模型, 分析了p53信號轉(zhuǎn)導(dǎo)網(wǎng)絡(luò)中的雙節(jié)率現(xiàn)象; Ochab和Puszynski[21]開發(fā)了一種分段線性方法, 并以p53振子為例對比了新方法與傳統(tǒng)的非線性方法; 劉楠等[22]重建了p53-Mdm2通路和p53-Wip1-ATM通路耦合的模型, 在Ser395位點磷酸化的Mdm2被考慮在內(nèi), 并探索了這種情景下Wip1管理的p53潛在動力學(xué); 等等.雖然p53-Mdm2之間的正反饋被提及, 但據(jù)我們所知,尚未有細致的研究.
基于上述考慮, 本文重新審視了Mdm2促進p53的作用對p53動態(tài)的影響.首先給出這個基因網(wǎng)絡(luò)模型對應(yīng)的動力學(xué)方程; 然后利用分岔分析的辦法, 找出了關(guān)鍵的分岔點, 系統(tǒng)的動力學(xué)特征由分岔點劃分的參數(shù)區(qū)間所決定; 最后根據(jù)分岔曲線, 在二維參數(shù)平面上找出了振蕩的區(qū)域.數(shù)值研究結(jié)果表明: p53-Mdm2之間正反饋強度的差異可能是p53面對兩種都能激活A(yù)TM的細胞壓力時出現(xiàn)不同動力學(xué)響應(yīng)方式的關(guān)鍵原因; Wip1對這個基因網(wǎng)絡(luò)的振蕩會時而有利時而不利, 取決于Mdm2對p53的正調(diào)控強度; 在這個模型中, p53-Mdm2的正負雙重反饋, 是振蕩出現(xiàn)的根基, 也就是說p53-Mdm2之間的正負反饋對p53振蕩而言是同等重要的.總之, 振蕩出現(xiàn)需要正負兩種反饋強度接近、相互偶聯(lián).
事實上, p53的動力學(xué)行為是壓力類型依賴的[23].振蕩的p53動力學(xué)出現(xiàn)在那些能引起DNA雙鏈斷鏈 (DSB) 的細胞壓力中, 例如電離輻射和依托泊苷刺激[12,24].ATM作為DSB的信號傳感器, 接收到DNA損傷信號并迅速傳導(dǎo)到p53-Mdm2模塊[25].一方面激活p53的轉(zhuǎn)錄活性, 另一方面通過催化Mdm2在Ser395位點的磷酸化促使Mdm2自我泛素化并被快速分解來增強p53的穩(wěn)定性[26].在Ser395位點磷酸化的Mdm2能結(jié)合p53的信使RNA, 使其翻譯效率被大幅提升[27].此外, Mdm2還有一種在Ser166或Ser168位點的磷酸化方式, 這類Mdm2的作用是促進p53失活和降解[28].Mdm2是p53的靶蛋白之一, 實驗表明另一種p53的靶蛋白Wip1也對維持形狀均勻的p53脈沖至關(guān)重要[29].作為一種蛋白激酶, Wip1能促進ATM的去磷酸化, 從而抑制DNA損傷誘導(dǎo)下的ATM激活.因此, p53-Mdm2之間既有正反饋環(huán)又有負反饋環(huán), 而p53-Wip1-ATM僅形成了負反饋環(huán), 如圖1所示, 其中紅色線形成的回路為負反饋環(huán), 綠色線形成的回路為正反饋環(huán).
圖1 模型示意圖.箭頭線表示促進或狀態(tài)轉(zhuǎn)移; 杠頭線表示抑制Fig.1.Schematic diagrams of the model.Arrow-headed lines indicate promotion or state transition; bar-head lines indicate inhibition.
要構(gòu)建基因網(wǎng)絡(luò)的動力學(xué)模型, 首先要明確網(wǎng)絡(luò)中的節(jié)點.我們的模型中共有6個節(jié)點: 活性的ATM, 活性的p53, Wip1, 在Ser395位點磷酸化的Mdm2 (Mdm2p1), 未修飾的Mdm2 (Mdm2i), 在Ser166或Ser168位點磷酸化的Mdm2 (Mdm2p2).根據(jù)模型所考慮的生化反應(yīng) (也就是ATM被DSB激活, 被Wip1滅活; p53的激活需要ATM,合成受到Mdm2p1影響; Wip1基因的轉(zhuǎn)錄依賴于p53; Mdm2基因的轉(zhuǎn)錄也依賴p53; 不考慮能終止p53振蕩的Akt通路[30], Mdm2可自由地在Ser166或Ser168位點被磷酸化以及去磷酸化, 導(dǎo)致Mdm2i與Mdm2p2相互轉(zhuǎn)化; Mdm2被ATM催化轉(zhuǎn)變?yōu)镸dm2p1; Mdm2p1又能自發(fā)地去磷酸化, 變回Mdm2i), 可給出如圖1 所示的模型框架.網(wǎng)絡(luò)中的每個節(jié)點對應(yīng)著一個方程.用方括號“[]” 表示無量綱的蛋白質(zhì)濃度, 在方括號里寫了蛋白質(zhì)的種類, 用一類 “S” 形函數(shù) (廣義的希爾函數(shù)) 來刻畫蛋白質(zhì)之間的相互調(diào)控, 則由這個基因網(wǎng)絡(luò)“翻譯”出的動力學(xué)方程組如下:
其中濃度單位為“C”, 時間單位為分鐘(min), 希爾系數(shù)s0=s1=s2=s3=s4=s6=4 ,s5=2 , 其他的參數(shù)值及描述如表1[22,30,31]所列.在標準參數(shù)下p53的振蕩周期在4—7個小時之間 (圖2), 與實驗數(shù)據(jù)一致.
表1 模型中的參數(shù)Table 1.Parameters of the model.
圖2 p53的時間歷程圖Fig.2.Time occurs of p53.
Mdm2對p53正調(diào)控的必要途徑是Mdm2在Ser395位點的磷酸化, 與參數(shù)kp對應(yīng).換言之,kp的值反應(yīng)了p53-Mdm2正反饋的強度.首先在圖3給出了p53平衡態(tài)與kp的函數(shù)關(guān)系, 即余維1分岔圖.在不動點分支上出現(xiàn)了三個分岔點: 鞍結(jié)分岔點 (SN)、鞍結(jié)同宿分岔點 (SNIC)、霍普夫分岔點 (HB).在SN和SNIC之間的藍線是不穩(wěn)定的鞍點, 在SN和HB之間的紅色線是不穩(wěn)定的結(jié)點和焦點, 在SNIC和HB之外的黑色線代表p53穩(wěn)定的穩(wěn)態(tài).在SNIC和HB之間, 這個動力系統(tǒng)中僅存在穩(wěn)定的極限環(huán), 也就是p53動態(tài)行為僅有穩(wěn)定振蕩.這里極限環(huán)的分支之所以沒有被延拓是因為極限環(huán)上的折疊分岔點靠近HB, 因此系統(tǒng)的動力學(xué)由不動點分支上的分岔點基本確定.當kp在(0, SNIC) 參數(shù)區(qū)間中時, p53最終會穩(wěn)定在低穩(wěn)態(tài); 當kp在(SNIC, HB) 參數(shù)區(qū)間中時, p53動力學(xué)呈現(xiàn)反復(fù)的、非衰減的脈沖; 當kp在(HB, +∞)參數(shù)區(qū)間中, p53最終會達到一個較高的穩(wěn)態(tài).
圖3 關(guān)于參數(shù) kp 的一維分岔圖Fig.3.One-dimensional bifurcation graph on parameter kp.
為了驗證一維分岔分析的結(jié)果, 取了3個kp值在([p53], [Mdm2])相平面畫出軌線(圖4), 這里的[Mdm2]為三種形式Mdm2的濃度之和.當kp=0.3時(圖3中的黑色虛線), 分岔圖顯示p53有穩(wěn)定的低穩(wěn)態(tài), 在圖4中對應(yīng)著黑色軌道, 從原點開始最終停在了不動點.當kp=0.65 時(圖3中的紅色虛線), 不動點是不穩(wěn)定的焦點, 在圖4中對應(yīng)著紅色軌道, 始于原點的軌跡最終收斂到極限環(huán).當kp=1 時(圖3中的藍色虛線), p53最終會穩(wěn)定在較高的濃度, 在圖4中對應(yīng)著藍色軌道,從原點出發(fā)沿著漸進的環(huán)狀軌跡最終到達不動點,并不再移動.相平面分析(圖4)的結(jié)果與分岔分析(圖3)的結(jié)果保持一致, 證實了分岔分析的正確性.在正反饋強度很弱的時候, 極限環(huán)消失, 說明了Mdm2介導(dǎo)的正反饋對p53振蕩發(fā)生的必要性.此外, 正反饋強度過強時候振蕩也會消失.實際上,正反饋過強時可通過調(diào)節(jié)一些負反饋的強度來抵消正反饋, 進而使得振蕩復(fù)現(xiàn), 詳情見3.2和3.3小節(jié).
圖4 相平面分析Fig.4.Phase plane analysis.
上述的分析表明, 要獲取二維參數(shù)平面上的振蕩區(qū)域, 只需延拓SNIC和HB兩個分岔點.由于ATM的激活參數(shù)vatm是這個基因網(wǎng)絡(luò)接收刺激信號的開端, 受到DSB的調(diào)控, 所以我們把vatm作為第二參數(shù)張成(kp,vatm)參數(shù)平面(圖5).通過延拓SNIC分岔點(圖5中的藍色線)和HB分岔點(圖5中的紅色線), 得到二維參數(shù)平面上的兩條分岔曲線, 由這兩條分岔線圍成的灰色區(qū)域為振蕩區(qū)域.圖5表明, 當正反饋強度較弱, 即kp較小時, 例如kp=0.6 , p53的振蕩行為對DNA損傷更魯棒,也就是說只要DNA損傷輸入的信號超過某個閾值就能觸發(fā)p53振蕩.這與Lahav等[12]的實驗結(jié)果一致, 即電離輻射強度(DNA損傷程度)不能改變p53脈沖式動力學(xué)的應(yīng)答模式.而正反饋強度較強, 例如kp=1 時, 過多或過少的DNA損傷都不能使p53振蕩, 與Chen等[24]的實驗一致, 即依托泊苷的劑量(DNA損傷程度)需要適中才能出現(xiàn)p53振蕩.因此我們推斷, 在引起DSB的兩類刺激下, p53出現(xiàn)對DSB相異的魯棒性有可能是p53-Mdm2之間正反饋強度的差異引起的, 這與Sun和Cui[32]的觀點(ATM誘導(dǎo)的Mdm2磷酸化通路有利于雙模式的p53動力學(xué))有相似之處.
圖5 二維參數(shù)平面( kp , v atm )上的振蕩區(qū)域Fig.5.Oscillation area on the ( kp , v atm ) two-dimensional parameter plane.
Batchelor等[29]的實驗發(fā)現(xiàn), 在細胞中加入Wip1信使RNA的抑制劑后, p53振蕩消失.此時的p53出現(xiàn)單個大幅度脈沖, 類似于在紫外線照射實驗中觀測到的p53脈沖.因此, Wip1通路作用強度的差異很可能是兩種刺激下出現(xiàn)不同的p53動力學(xué)應(yīng)答的根本原因[23].為了探究Wip1和振蕩發(fā)生的關(guān)系, 在(kp,vwip1)參數(shù)平面上延拓出SNIC分岔線和HB分岔線; 同樣, 灰色區(qū)域為振蕩區(qū)域(圖6).當vwip1增大到一定程度時, 兩條分岔曲線與x軸接近垂直, 此時Wip1的作用飽和,也就是Wip1含量雖然充足, 但是反應(yīng)底物活性ATM的含量較低, 使得ATM失活反應(yīng)的最大速率被限制.這點在模型中體現(xiàn)在關(guān)于Wip1的希爾函數(shù)中, 在Wip1的濃度很高的情形下, 繼續(xù)增加Wip1濃度不再對這個基因網(wǎng)絡(luò)的其他節(jié)點有影響.顯然, 在標準參數(shù)下(kp=0.65 ), Wip1是振蕩出現(xiàn)的重要條件, 與文獻[22]結(jié)論一致.
圖6 二維參數(shù)平面( kp , v wip1 )上的振蕩區(qū)域Fig.6.Oscillation area on the ( kp , v wip1 ) two-dimensional parameter plane.
有趣的是, 當kp較小時,vwip1=0 也會發(fā)生振蕩.這意味著Wip1被抑制的情況下, 通過削弱p53-Mdm2之間的正反饋也能使得p53出現(xiàn)振蕩,在生物實驗中應(yīng)注意這一點.為了進一步說明這些情況, 在兩種kp值下做了[p53] 關(guān)于參數(shù)vwip1的余維1分岔圖(圖7).當kp較小時, 如圖6黑色小三角所指(kp=0.5 ), 此時Wip1可能不利于p53振蕩; 在圖7(a)中, 隨著vwip1增大, 到超過SNIC分岔點時, 不動點為穩(wěn)定的結(jié)點, 極限環(huán)消失.而當kp較大時, 如圖6紅色小三角所 (kp=0.7 ),此時Wip1有利于p53振蕩; 在圖7(b)中, 隨著vwip1增大, 到超過HB分岔點時, 不動點為不穩(wěn)定的焦點, 極限環(huán)出現(xiàn).因此, 在這個模型中Wip1對系統(tǒng)振蕩的影響是由p53-Mdm2之間的正反饋強度決定的.非時滯模型的振蕩需要恰當?shù)恼答亸姸? 默認參數(shù)下(kp=0.65 )正反饋較強, 附加的p53-Wip1-ATM負反饋回路與p53-Mdm2p1正反饋回路拮抗, 所以能出現(xiàn)振蕩.
圖7 關(guān)于參數(shù) k wip1 的一維分岔圖Fig.7.One-dimensional bifurcation graph on parameter k wip1.
最后探究p53-Mdm2之間的正負反饋協(xié)作對系統(tǒng)振蕩出現(xiàn)的影響.模型中假設(shè)Mdm2要實現(xiàn)對p53的抑制作用, 即破壞p53的穩(wěn)定性和轉(zhuǎn)錄活性, 首先要經(jīng)過Akt的磷酸化.事實上, 在文獻[30]的模型中, Akt的活性也受到p53的間接調(diào)控, 并且這種作用很可能引起p53的雙穩(wěn)態(tài).在本文考慮的模型中,kp?被視為定常量, 是控制Mdm2抑制p53這條通路的關(guān)鍵參數(shù).如圖8所示, 在(kp,kp?)參數(shù)平面中的振蕩區(qū)域(灰色部分)與x軸和y軸都不相交, 這代表p53-Mdm2p1通路和p53-Mdm2p2通路都是振蕩的要素.負反饋回路是振蕩發(fā)生的網(wǎng)絡(luò)結(jié)構(gòu)基礎(chǔ), 在很多p53振子的理論模型中, p53-Mdm2的負反饋都起到了重要作用[13?16], 本文的模型中起主導(dǎo)作用的負反饋回路同樣是p53-Mdm2.參數(shù)kp的振蕩區(qū)間(SNIC,HB)隨著kp?的增大而擴寬, 進一步說明了p53-Mdm2p2負反饋環(huán)有利于p53的振蕩.隨著負反饋強度的增加, 振蕩所需要的正反饋強度也隨之增加, 再次說明了振蕩的出現(xiàn)需要兩者同時參與、相互制約.事實上, Zhang等[33]提出, 用正負反饋回路耦合構(gòu)建的p53振子模型能產(chǎn)生更具魯棒性的振蕩.
圖8 二維參數(shù)平面( kp , k p? )上的振蕩區(qū)域Fig.8.Oscillation area on the ( kp , k p? ) two-dimensional parameter plane.
總而言之, 本文全面探討了Mdm2對p53的正調(diào)控(即Mdm2p1促進p53的翻譯)和p53振蕩發(fā)生之間的關(guān)聯(lián).首先利用生物事實給出一個p53基因網(wǎng)絡(luò)的數(shù)學(xué)模型; 然后利用數(shù)值分岔分析的方法, 研究了刻畫p53-Mdm2正反饋強度的參數(shù)kp: 1)kp與DSB的協(xié)作關(guān)系, 提出了實驗發(fā)現(xiàn)的p53對兩類能引起DSB的刺激(電離輻射、依托泊苷處理)出現(xiàn)不同動力學(xué)響應(yīng)模式的原因可能是kp的不同; 2)kp與Wip1的協(xié)作關(guān)系, 指出了Wip1一方面在高強度正反饋下能促進振蕩, 另一方面在低強度正反饋下能抑制振蕩; 3)kp與Mdm2p2的協(xié)作關(guān)系, 再次證實了p53-Mdm2的負反饋回路是振蕩發(fā)生的結(jié)構(gòu)基礎(chǔ).本文強調(diào)的是正負反饋對p53基因網(wǎng)絡(luò)發(fā)生振蕩的同等重要性.值得注意的是, 希爾函數(shù)和微分方程刻畫的是宏觀層次的模型, 今后的研究還需建立更為細致的多態(tài)模型[34].p53的振蕩與細胞周期停滯或細胞凋亡等過程密不可分[35], 希望本研究能為今后的實驗設(shè)計和醫(yī)學(xué)應(yīng)用提供思路.
此外, 本文的結(jié)果是在非線性程度(協(xié)同性)很高的調(diào)節(jié)關(guān)系上得出的, 這里的非線性程度很高主要體現(xiàn)在數(shù)值較大的希爾系數(shù)上.為了使結(jié)論更加可靠, 檢測希爾系數(shù)較小時振蕩發(fā)生的條件很有必要.把s0調(diào)到1;s1,s2和s3下調(diào)到2;s4,s6的
值不變 (因為p53以4聚體的形式綁定靶基因[36],其希爾系數(shù)理論上應(yīng)取為4).為了使p53出現(xiàn)振蕩,k1的值改為0.09,k3的值改為0.05, 其他參數(shù)保持默認值, 在圖9中作出關(guān)于參數(shù)kp的分岔圖組.從圖9可以看出, 不同于強非線性的分岔現(xiàn)象是分岔點SN以及SNIC消失, 振蕩區(qū)間變?yōu)閮蓚€HB之間的參數(shù)區(qū)間(HB1, HB2); 其他的振蕩控制現(xiàn)象未發(fā)生變化.也就是(kp,vatm)參數(shù)平面上的振蕩區(qū)域仍由兩條形分岔曲線所圍成, 當kp較小時p53能在ATM過度激活的情況下保持振蕩; (kp,vwip1) 參數(shù)平面上的振蕩區(qū)域由兩條形分岔曲線所圍成, 存在一部分kp的值使得Wip1抑制振蕩; (kp,kp?)參數(shù)平面上的振蕩區(qū)域由兩條“ / ”形直線所圍成, p53-Mdm2的正負反饋偶聯(lián)依然是這個基因網(wǎng)絡(luò)振蕩的“引擎”.因此, 本文的定性結(jié)論可能僅取決于這個基因網(wǎng)絡(luò)的結(jié)構(gòu),對希爾函數(shù)的非線性程度并不敏感.結(jié)論的可靠性得到了進一步的驗證.
圖9 分岔圖組Fig.9.Bifurcation graph group.