本發(fā)明屬于系統(tǒng)生物學技術領域，主要涉及生物信息學和生物數據挖掘，具體涉及一種基于模塊化因子圖的信號通路機制確認方法。

背景技術：

多發(fā)性骨髓瘤(multiple myeloma，MM)是一種以惡性漿細胞克隆性增殖為特點的惡性腫瘤，是目前血液系統(tǒng)第二大惡性腫瘤。目前MM的治療主要包括傳統(tǒng)化療、新藥靶向治療及免疫治療等。雖然新的靶向治療明顯提高了MM的療效,但患者中位生存時間仍在3～5年,其發(fā)病機制尚不明確。因此,進一步研究影響MM細胞生長的相關機制,尋建立模型研究方法，是亟需解決的問題。

無論從生物學還是從臨床表現看,MM細胞的特性不僅僅決定于其遺傳學特性(如染色體重排,缺失,擴增或某些特定基因的突變)。相反,該疾病的病理生理學表現明顯受MM細胞與其所處的骨髓微環(huán)境間雙向相互作用的影響。Virginia Hughes指出,骨髓微環(huán)境對MM細胞的存活、生長以及耐藥等重要環(huán)節(jié)有著息息相關的作用。骨髓基質細胞(bone marrow stromal cells,BMSCs)作為骨髓微環(huán)境的主要成員,與MM的發(fā)生、發(fā)展有著密切的關聯。

隨著對MM生物學研究的不斷深入,人們發(fā)現在腫瘤發(fā)生發(fā)展過程中，信號通路控制著眾多至關重要的細胞生物學過程。目前，對MM信號通路靶點的研究十分廣泛，主要靶向的通路包括這些信號通路包括PI3K/Akt/mTOR/P70S6K信號通路,IKK-αF/NF-κB信號通路,Ras/Raf/MAPK信號通路和JAK/STAT3通路,它們都可以通過以下途徑被激活:上游的細胞因子與相應受體的結合,或通過粘附啟動的激酶途徑直接由細胞粘附誘導增殖、抗凋亡信號通路的激活。

眾所周知，傳統(tǒng)的生物實驗非常昂貴并且要花費大量的時間，所以近年來越來越多的人在用生物模型去模擬生物生長狀況，從模擬的層面上去分析藥物影響或者提取關鍵蛋白質。Huiming Peng等人用系統(tǒng)生物學的方法研究p38MAPK異型的抗藥性，確定生物模型之后利用設置參數值的方式分別去探索p38的五種異型的抗藥性；Xiaoqiang Sun等人基于細胞內的信號通路利用微分方程建模的方法研究在組織骨再生的過程中細胞因子的組合預測，對人體組織骨再生成時不同細胞因子對成骨細胞和破骨細胞的刺激作用進行了探索，并且篩選出較好的細胞因子組合。但是，上面的研究方法存在一定的局限性，并沒有指出如何預測未知的致病因子的影響。

遺傳學改變和。生化條件引起的通路激活經常發(fā)生在腫瘤惡變早期和進展期，同時它們也是患者預后的重要指標。系統(tǒng)生物學方法期望通過建立細胞信號轉導過程的模型,找到參與此過程的各種分子之間相互作用的網絡,闡明其在基因調控、疾病發(fā)生中的作用。近幾年,對信號轉導網絡的定量分析逐日升溫,通常采用一系列的方程模型描述信號轉導通路的內部變化過程。

現有技術中，采用的主要方案包括以下幾種：

(1)臨床實驗。在對于多發(fā)性骨髓瘤的研究中，目前絕大多數還是處于利用實驗方法去觀察腫瘤細胞的生長發(fā)育，尤其治療期間，基本上是靠醫(yī)生的經驗去判斷。實驗成本比較昂貴。

(2)常微分方程(ordinary differential equations,ODE)。這是是描述動力學系統(tǒng)的常用方法，應用微分方程組可以構建一個復雜的數學模型，用以代表一系列生化反應的相互作用模式，并且模擬生物系統(tǒng)中各組分的時序性動態(tài)變化。常微分方程(ODE)是質量反應動力學過程的數學代表，可以用來描述連續(xù)時間范圍內生物系統(tǒng)各組分的動態(tài)變化。對于那些不考慮空間大小，并且反應速度和反應底物的濃度成一定的比例關系的生化反應系統(tǒng)比較適用。Chen利用ODEs來描述ErbB信號通路的輸入輸出對細胞分化和增值的影響,文中用299個ODE方程表示828個級聯反映，共有229個參數，計算規(guī)模很大。

(3)Petri Net的不確定性、并行性、異步性，以及對分布式系統(tǒng)的描述和分析能力使其在描述生物系統(tǒng)特性時有很大的優(yōu)勢。Chen和等用Petri nets建立了鳥氨酸循環(huán)的代謝模型。

上述方案主要存在的缺陷有以下幾點：

(1)人為經驗判斷，準確率不高。

(2)ODE等數學方法僅僅是描述生物量的變化，并不能直接以圖形的方式展現生物系統(tǒng)的結構特性。若可實現生物定量數據與圖形的結合與自動轉化，將能更好地刻畫生物系統(tǒng)結構與動態(tài)性質之間的關系。

(3)標準的Petri Net常用來定性分析生物網絡的結構性質，不能用于生物計算。

以下對本發(fā)明所涉及到的技術詞匯/技術術語注釋如下：

1、多發(fā)性骨髓瘤(multiple myeloma，MM)

2、骨髓基質細胞(bone marrow stromal cells,BMSCs)

3、常微分方程(ordinary differential equations,ODE)

4、反向蛋白質陣列(reverse phase protein arrays,RPPA)

技術實現要素：

有鑒于此，本發(fā)明在總結前人的研究基礎上，提出建立一個多層次的計算系統(tǒng)生物學模型來研究多發(fā)性骨髓瘤細胞的生長機制，利用現有RPPA數據，結合常微分方程組和Petri網來描述信號通路，并且對骨髓基質細胞和腫瘤細胞通路之間相互作用的組合影響進行量化的探索，并且使用模塊化的計算模型，降低計算成本。

具體而言，本發(fā)明所提出的技術方案如下：

本發(fā)明提供了一種基于模塊化因子圖的骨髓瘤信號通路機制確認方法，該方法包括：

步驟1、獲取RPPA數據；

步驟2、對所述RPPA數據進行預處理，粗粒度篩選關鍵蛋白質；

步驟3、基于所述粗粒度篩選的關鍵蛋白質，構建在細胞剛性環(huán)境下，細胞內蛋白質的相互作用通路，形成信號通路；

步驟4、采用常微分方程描述所述信號通路，并將所述信號通路分解成多個小模塊，針對每個所述小模塊，進行參數優(yōu)化，建立系統(tǒng)生物學模型；

步驟5、對所述系統(tǒng)生物學模型進行參數分析，所述參數分析包括穩(wěn)定性分析和敏感性分析。

優(yōu)選地，所述步驟1中，RPPA數據的獲取，通過以下方式：

步驟1.1、用壓強為100pa和400pa的細胞膠體模擬正常細胞和腫瘤細胞，記錄不同時間點細胞內蛋白質的濃度；

步驟1.2、利用蛋白質芯片，獲得正常細胞和腫瘤細胞的兩組RPPA數據。

優(yōu)選地，所述步驟2具體包括：

步驟2.1、對粗粒度篩選出的全部關鍵蛋白質數據，以t＝0min為標準，進行規(guī)范化，所述規(guī)范化方法為：

其中t₀表示t＝0min，表示第i個蛋白質在t_j時刻的濃度，表示第i個蛋白質在t₀時刻的濃度，為規(guī)范化后的蛋白質濃度；

步驟2.2、計算正常細胞、腫瘤細胞內蛋白質濃度變化率，具體方式為：

將其中濃度變化大于50％的蛋白質作為有意義的表達的蛋白質，作為粗粒度篩選出的關鍵蛋白質。

優(yōu)選地，所述步驟3具體包括：

步驟3.1、基于粗粒度篩選出的關鍵蛋白質，通過IPA數據庫，搜索相互作用的通路；

步驟3.2、在所述步驟3.1中搜索出的通路中，選擇p≤0.05的通路，作為信號通路，其中p表示某蛋白質在該通路中出現的誤差率。

優(yōu)選地，所述步驟4具體包括：

步驟4.1、使用常微分方程組描述信號通路，并使用RPPA數據中，高水平表達的蛋白質在不同時間點的采樣數據，確定信號通路中的關鍵參數；使用Petri網描述整個信號通路；

步驟4.2、基于所述Petri網描述的整個信號通路，將整個信號通路分解成多個小模塊；使用粒子群優(yōu)化方法優(yōu)化各個所述小模塊參數，獲得相對較小的參數范圍，所述參數優(yōu)化的目標函數是：

其中，表示蛋白質濃度時間序列數據，表示通過常微分方程獲得的模擬的蛋白質濃度時間序列，i表示蛋白質索引,t_j表示時間點，Θ表示常微分方程中的參數，M表示蛋白質數量，N表示時間點數量；

步驟4.3、把整個信號通路分解成兩個子網通路，并使用因子圖表示每個子網通路，為因子圖中的每個因子節(jié)點構造適應函數，所述適應函數為：

步驟4.4、使用置信度傳播方法，調和所述兩個子網通路中共享的蛋白質參數，并以所述步驟4.2中獲得的相對較小的參數范圍中的參數，作為置信度傳播方法的輸入參數，得到一個更小的參數范圍；

步驟4.5、以所述更小的參數范圍中的參數，作為粒子群優(yōu)化方法的輸入，對所述系統(tǒng)生物學模型進行參數優(yōu)化，獲得最終的系統(tǒng)生物學模型。

優(yōu)選地，所述步驟4.2中，整個信號通路分解成多個小模塊，可依據如下分解原則：a.每個子模塊中蛋白質數據盡可能少；b.每個子模塊中至少有一個蛋白質濃度是有實際數據依據。

優(yōu)選地，所述步驟4.3中，分解成子網通路可依據如下規(guī)則：a.從細胞表型出發(fā)，依次找出促進細胞增長或細胞死亡的蛋白質；b.如果按兩條表型找出的兩類蛋白質中共享蛋白質數量超過90％，則把兩條子網合并為一個大網；c.如果其中一個大的子網的蛋白質數量是另外一個子網蛋白質數量的2倍或以上，則重新分解這個大的子網。

優(yōu)選地，所述粒子群優(yōu)化方法的具體方式如下：

v_i(t+1)＝wv_i(t)+c₁·rand()·(p_i(t)-x_i(t))+c₂·Rand()·(p_g(t)-x_i(t))

x_i(t+1)＝x_i+v_i(t+1)。

優(yōu)選地，所述穩(wěn)定性分析通過計算變異系數對參數進行分析，所述變異系數的計算方法如下：

C·V＝(標準偏差SD/平均值Mean)×100％。

優(yōu)選地，所述敏感性分析中，參數的敏感性計算方法如下：

其中，ΔP_i是第i個參數值的變化量,表示一個很小的變化，例如可以是增加或減少1％，[ProteinName]表示系統(tǒng)輸出蛋白質，例如Casp3、p90RSK、CyclinD1、p21、p7056k，P_i表示優(yōu)化的參數，s_i表示敏感性參數。

與現有技術相比，本發(fā)明技術方案具有以下的有益效果：

(1)基于腫瘤細胞信號通路的互相作用，用系統(tǒng)生物學的方法建立了計算模型模擬腫瘤細胞的增值和凋亡。

(2)用常微分方程去描述信號通路中的反應，用petri網描述生物信號通路網絡結構，實現生物定量數據與圖形的結合與自動轉化，更好地刻畫生物系統(tǒng)結構與動態(tài)性質之間的關系。然后將整個信號通路基于規(guī)則進行模塊性的劃分，用尋優(yōu)算法獲取最優(yōu)值，降低了計算成本，提高了計算效率。

(3)通過模擬手段提高對腫瘤的生長過程中的分子機制的認識，實現對轉移精確地進行預測。

附圖說明

為了更清楚地說明本發(fā)明實施例或現有技術中的技術方案，下面將對實施例或現有技術描述中所需要使用的附圖作簡單地介紹，顯而易見地，下面描述中的附圖僅僅是本發(fā)明的一些實施例，對于本領域普通技術人員來講，在不付出創(chuàng)造性勞動的前提下，還可以根據這些附圖獲得其它的附圖。

圖1為本發(fā)明實施例的方法流程圖；

圖2為本發(fā)明實施例的RPPA粗粒度篩選蛋白質計算結果示例圖；

圖3為本發(fā)明實施例的構建與細胞剛性相關的信號通路示例圖；

圖4為本發(fā)明實施例的最大限度簡化通路示例圖；

圖5為本發(fā)明實施例的用混合Petri網去描述信號通路示例圖；

圖6為本發(fā)明實施例的將信號通路分解成n個小模塊示例圖；

圖7為本發(fā)明實施例的信號通路分解成的子網通路一示例圖；

圖8為本發(fā)明實施例的信號通路分解成的子網通路二示例圖；

圖9為本發(fā)明實施例的參數穩(wěn)定性分析結果示例圖；

圖10為本發(fā)明實施例的參數敏感性分析結果示例圖；

圖11為本發(fā)明實施例的不同條件下蛋白質對應的參數變化示例圖。

具體實施方式

下面結合附圖對本發(fā)明實施例進行詳細描述。應當明確，所描述的實施例僅是本發(fā)明一部分實施例，而不是全部的實施例?；诒景l(fā)明中的實施例，本領域普通技術人員在沒有作出創(chuàng)造性勞動前提下所獲得的所有其它實施例，都屬于本發(fā)明保護的范圍。

本領域技術人員應當知曉，下述具體實施例或具體實施方式，是本發(fā)明為進一步解釋具體的發(fā)明內容而列舉的一系列優(yōu)化的設置方式，而該些設置方式之間均是可以相互結合或者相互關聯使用的，除非在本發(fā)明明確提出了其中某些或某一具體實施例或實施方式無法與其他的實施例或實施方式進行關聯設置或共同使用。同時，下述的具體實施例或實施方式僅作為最優(yōu)化的設置方式，而不作為限定本發(fā)明的保護范圍的理解。

本發(fā)明是通過找出影響骨髓瘤細胞生存的信號通路中的關鍵蛋白質，輔助對于腫瘤細胞增生的模型建立方法研究。根據骨髓瘤細胞的剛性環(huán)境(骨髓瘤細胞壓強高)有助于腫瘤細胞的增生的結論，用壓強為100pa和400pa的細胞分別建立正常和腫瘤細胞的模型，測得得到兩組反相蛋白質陣列數據，通過實驗數據(100pa VS 400pa)粗粒度篩選出可能影響腫瘤細胞增長的蛋白質，并使用常微分方程組和Petri網描述信號相關的細胞信號通路，利用模塊化因子圖算法優(yōu)化模型中的關鍵參數，在最后通過合理和準確性驗證，最終通過該模型得出敏感度較高的蛋白質FAK、NFκB、mTOR1。圖1是本發(fā)明的總體流程圖，以下結合圖1對本發(fā)明的模型建立和計算方法進行詳細闡述。

(1)獲取RPPA(Reverse Phase Protein Arrays)數據，建立模型基礎數據庫

獲取RPPA數據的方式可以是多樣的，例如，可以采用現有的公用數據庫或資料庫中提供的已有的骨髓瘤細胞比對RPPA數據，即反向蛋白質陣列數據，其中對比數據模擬生長壓強需設置為100pa和400pa，以此作為模型建立的對比數據，從而建立數據庫?，F有技術中有多種途徑可以獲得上述RPPA公開數據，這里不再贅述。

也可以通過做模擬對比實驗，建立細胞生長的模型，設置例如壓強為100pa和400pa的細胞膠體模型模擬正常和腫瘤細胞，記錄在不同時間點細胞內蛋白質的濃度，在一個具體的實施方式中，該不同的時間點例如可以設置為0min,30min,60min,overnight，共4個時間點，當然，也可以根據具體的算法需要，設定不同數量的采樣時間點。利用蛋白質芯片得到兩組RPPA數據。該數據中包括在不同壓強下173種蛋白質在不同時間點下的濃度。

當然，此處也可以采用其他的已有技術或途徑，獲取可以作為后續(xù)模型算法所需要的骨髓瘤細胞的RPPA數據集。

(2)數據預處理，粗粒度篩選關鍵蛋白質

對所有蛋白質數據以t＝0min為標準對其規(guī)范化，其計算公式為：

其中t₀表示t＝0min，表示第i個蛋白質在t_j時刻的濃度，表示第i個蛋白質在t₀時刻的濃度，為規(guī)范化后的蛋白質濃度。

計算兩種細胞內蛋白質濃度變化率，計算公式為：

根據此公式找出其中濃度變化大于50％的蛋白質作為有意義的表達的蛋白質。圖2是根據RPPA粗粒度篩選蛋白質計算結果示例圖。

(3)構建多發(fā)性骨髓瘤信號通路

信號通路是指能將細胞外的分子信號經細胞膜傳入細胞內發(fā)揮效應的一系列酶促反應通路。這些細胞外的分子信號(稱為配體，ligand)包括激素、生長因子、細胞因子、神經遞質以及其它小分子化合物等。細胞內各種不同的生化反應途徑都是由一系列不同的蛋白組成的，執(zhí)行著不同的生理生化功能。各個信號通路中上游蛋白對下游蛋白活性的調節(jié)(包括激活或抑制作用)主要是通過添加或去除磷酸基團，從而改變下游蛋白的立體構象完成的。所以，構成信號通路的主要成員是蛋白激酶和磷酸酶，它們能夠快速改變和恢復下游蛋白的構象。

在本發(fā)明中，根據步驟(2)粗粒度篩選的蛋白質來構造在細胞剛性環(huán)境下MIC和MSC細胞內蛋白質的相互作用通路。首先將篩選出的蛋白質通過IPA(Ingenuity Pathway Analysis)數據庫尋找相互作用的重要的通路，IPA是基于云計算的一體化應用軟件，它可以分析來源于基因組、microRNA、SNP、芯片、代謝組、蛋白組、RNA-Seq實驗以及各類小規(guī)模實驗數據。利用IPA，可以搜索基因、蛋白、化學分子、藥物的各類信息，并且?guī)椭鷺嫿ㄏ嗷プ饔媚Ｐ汀Ｔ诟叨冉Y構化、集成豐富詳實生物化學知識的Ingenuity Knowledge Base支持下，IPA的分析和搜索可以幫助所獲得數據在生物體系中重要性。然后在搜索出的通路中選擇p≤0.05的通路，其中p值表示某蛋白質在該pathway(即通路)中出現的誤差率。由于在找出的通路中有一些其他的蛋白質，通過共享蛋白質來整合所有細胞通路并參考MM細胞的相關文獻或現有的公開數據，來構建與細胞剛性相關的信號通路，如圖3所示。由于估計由整個信號通路建立的常微分方程組中的所有參數是不可能的，因此我們需要重新定義信號通路，在保證通路結構完整的情況下最大限度簡化通路，如圖4所示。

(4)建立系統(tǒng)生物學模型

系統(tǒng)生物學模型所包含生化反應的動力學參數能夠明顯影響數學模擬結果，可以通過現有技術或相關文獻獲得，也可以通過體外的生理和生化實驗測定。但是，某些生化反應的動力學參數無法直接獲得，需要通過實驗數據分析加以估計幾乎所有描述真實生物系統(tǒng)的模型都過于巨大，其內部的動態(tài)變化過程無法被實驗數據全部描述。因此，系統(tǒng)生物學模型在不斷涌現的多實驗數據的整合分析中能夠起到關鍵的作用。

(4.1)在本發(fā)明中，我們使用常微分方程組描述信號通路，并采用RPPA試驗數據中高水平表達的蛋白質在各個時間點的采樣數據來確定信號通路中的關鍵參數，接上例，設定該采樣時間點位4個。在一個具體的實施方式中，結合圖5至圖8，常微分方程可以通過以下方式進行構建：以pCylinD為例，結合圖5至圖8中可見，GSK3β促進CylinD的磷酸化，而pp21抑制CylinD的磷酸化，根據這些相互作用可以得到：

其中[CyclinD]、[pp21]、[pGSK3β]分別表示蛋白質CyclinD、pp21、pGSK3β的濃度；kCyclinD_pp21表示CylinD被pp21抑制的磷酸率、kCyclinD_pGSK3β表示CylinD被pGSK3β激活的磷酸率。

遵循這樣的規(guī)律，根據圖6,我們得到用來描述其它蛋白質的40個常微分方程，該其它蛋白質即如附圖6信號通路中涉及到的蛋白質，由于其微分方程的表達式相同，此處不再一一贅述。由于這些常微分方程包含45未知參數并且蛋白質之間有相互作用，使用一般的智能算法來確定參數值是相當復雜的。本發(fā)明使用調整后的馬爾科夫鏈模塊分解方法把整個通路分解成多個小模塊，然后用智能算法來優(yōu)化每個小模塊的參數。該方法不僅可以減少計算機的壓力，也可以快速獲得更優(yōu)結果。

(4.2)采用信號通路的微分方程模型可以表征細胞在不同時間內分子濃度的變化,觀察在不同外界條件下細胞的動態(tài)響應過程,深入的揭示信號轉導的作用機制。同時，用混合Petri網去描述信號通路中的每一個生化反應及級聯反應，如圖5所示，把整個信號傳導通路基于Petri網和分解規(guī)則分解成n個小模塊，如圖6所示，應用PSO算法(即粒子群優(yōu)化算法)對每個小模塊進行估參(即參數估計)，從而在初始范圍搜索空間中獲得一個相對較小的參數范圍。

在一個具體的實施方式中，具體的分解規(guī)則可以設置為：a.每個子模塊中蛋白質數據盡可能少；b.每個子模塊中至少有一個蛋白質濃度是有實驗數據的，或者是有確定的數據支撐的，這樣才能保證計算的高效性。并使用粒子群算法來優(yōu)化各個小模塊的參數。

得到的蛋白質濃度時間序列數據用表示，通過ODE方程(即常微分方程)獲得的模擬的蛋白質濃度時間序列其中，i表示蛋白質索引,t^j表示時間點，Θ表示ODE方程中的參數，參數優(yōu)化的目標函數是：

其中，M表示蛋白質數量，N表示時間點數量，最終得到一個參數初始范圍。

(4.3)然后把整個信號通路分解成兩個子網通路，如圖7、圖8所示，為每一個子通路構造因子圖，并為每一個因子節(jié)點構造適應函數，即公式(1)。

在一個具體的實施方式中，可以為分解成子網通路設定規(guī)則：a.從細胞表型出發(fā)(細胞增長、細胞死亡)，依次找出促進細胞增長或細胞死亡的蛋白質；b.如果按兩條表型找出的兩類蛋白質中共享蛋白質數量超過90％，則把兩條子網合并為一個大網；c.如果其中一個大的子網的蛋白質數量是另外一個子網蛋白質數量的2倍或以上，則重新分解這個大的子網。另外，因子圖是指將一個具有多變量的全局函數因子分解，得到幾個局部函數的乘積，以此為基礎得到的一個雙向圖，叫做因子圖，因子圖是由變量節(jié)點和因子節(jié)點，以及連接兩個節(jié)點的邊構成。

(4.4)應用置信度傳播(Belief Propagation,BP)方法調和兩個子通路中共享的蛋白質參數，解決兩個子網通路相同蛋白質對應的沖突的參數，從而得到一更優(yōu)的參數范圍，以(4.2)步中輸出的相對較小的參數作為BP算法的輸入參數，然后得到一個更小的參數范圍，最后把這個范圍作為最后應用PSO算法進行系統(tǒng)估參的輸入。從而實現整個模型形成一個不斷縮小搜索范圍的過程。

其中表示適應度函數，表示模擬結果與實驗結果的誤差，分別表示參數集合和與因子節(jié)點對應的分子濃度水平的集合。表示通過包括參數集Θ的ODEs方程計算得到的蛋白質m在時間點t_j模擬濃度。表示蛋白質m在時間點t_j的實驗水平的濃度。

(4.5)估參方法中，本發(fā)明采用粒子群優(yōu)化算法PSO，以步驟(4.4)中得到的參數范圍作為初始搜索空間繼續(xù)對參數進行優(yōu)化，從而得到最優(yōu)值。

在一具體的實施方式中，可在該步中，重復運行例如5遍或更多遍，保證算法的穩(wěn)定性。

PSO中，每個優(yōu)化問題的解都是搜索空間中的一個粒子。所有的粒子都有一個由被優(yōu)化的函數決定的適應值(fitness value)，每個粒子還有一個速度決定他們飛翔的方向和距離。在本發(fā)明中，PSO初始化為一群隨機粒子(隨機解)。然后通過迭代找到最優(yōu)解。在每一次迭代中，粒子通過跟蹤兩個"極值"來更新自己。第一個就是粒子本身所找到的最優(yōu)解，這個解叫做個體極值。另一個極值是整個種群目前找到的最優(yōu)解，這個極值是全局極值。這個過程一直迭代，直到找到比較滿意的解為止。此解作為整個模型估參的最優(yōu)解。

(5)參數分析

a.穩(wěn)定性分析

某些生物系統(tǒng)(如代謝網絡和信號通路等)具有特殊的潛能，能夠聚合于特定的平衡態(tài)(或稱為穩(wěn)定狀態(tài))。當生物系統(tǒng)達到穩(wěn)定狀態(tài)時，系統(tǒng)內所有組分濃度改變速率為零。穩(wěn)態(tài)分析可以鑒定并合理描述這些穩(wěn)定狀態(tài)，這有助于闡釋生物系統(tǒng)內部的穩(wěn)態(tài)維持機制，以及細胞凋亡響應等生物學過程所涉及的狀態(tài)快速轉換機制。在本發(fā)明中，穩(wěn)定性分析指分析參數是否穩(wěn)定，對于某個參數，如果其變異系數大于1，則認為不穩(wěn)定，否則就是穩(wěn)定的。變異系數的計算公式為：C·V＝(標準偏差SD/平均值Mean)×100％

根據此公式，對本文中的所有參數進行穩(wěn)定性分析。如圖9所示。其中有10個參數的變異系數大于1，即所占比例為22％，意味著在這個模型里將近80％的參數是穩(wěn)定的。

b.敏感性分析

在系統(tǒng)生物學模型構建完成之后，簡單的結構化和參數化并不能很好地解釋生物系統(tǒng)的內在調控機制，因此對模型特性的整體分析是非常重要的。敏感性分析(Sensitivity Analysis)是指從定量分析的角度研究在一定范圍內擾動系統(tǒng)組分的初始濃度或反應動力學參數，對模型輸出結果影響的大小(特定組分的濃度變化或模型整體狀態(tài)的改變)。生物學實驗上可以理解為敲除或過表達特定基因產物，之后測定生物體或細胞特定表型的變化。敏感性分析能夠闡釋系統(tǒng)輸出對特定參數值的依賴性，這通常是研究者最為感興趣的。值得注意的是，從理論分析的角度進行敏感性分析時，初始濃度或反應參數的擾動范圍可以隨意選取。但是從實踐的視角出發(fā)，如果期望獲得的敏感性分析結果是可靠的，那么擾動的選取范圍必須符合真實的生理生化狀態(tài)。此外，敏感性分析方法對于尋找生物網絡的關鍵調控位點具有重要的指導意義。

在本發(fā)明中，敏感性分析是一種衡量某一特定參數對輸出的影響，分析參數的敏感性，即參數變化是否會引起整個系統(tǒng)輸出值的變化，簡單來說，給每一個參數增加或減少1％，觀察對系統(tǒng)輸出的影響，從而找出敏感的蛋白質。根據信號通路可知輸出為蛋白質Casp3、p90RSK、CyclinD1、p21、p7056k濃度。參數的敏感性計算公式：

其中，ΔP_i是第i個參數值的變化連,表示一個很小的變化(增加或減少1％)，[ProteinName]表示系統(tǒng)輸出蛋白質Casp3、p90RSK、CyclinD1、p21、p7056k，P_i表示優(yōu)化的參數，s_i表示敏感性參數。其中，每個參數通過增加或減少1％，觀察Casp3、p90RSK、CyclinD1、p21、p7056k變化百分比，示例如圖10，通過分析可得系統(tǒng)輸出變化較大的參數占所有參數的比例不超過2％，該結果表示所有參數中只有5-7個參數比較敏感，這說明系統(tǒng)輸出蛋白質的濃度變化受這幾個比較敏感的參數影響。

綜合上述結果，對比正常細胞和腫瘤細胞對應的參數變化值，如圖11所示，由上圖可以看出幾種蛋白質濃度變化范圍比較大，則認為FAK、NFκB、mTOR1這幾種蛋白質可能是通路中對腫瘤細胞生長或增生影響比較大的，從而實現對蛋白質的篩選。

以上所述，僅為本發(fā)明的具體實施方式，但本發(fā)明的保護范圍并不局限于此，任何熟悉本技術領域的技術人員在本發(fā)明揭露的技術范圍內，可輕易想到的變化或替換，都應涵蓋在本發(fā)明的保護范圍之內。因此，本發(fā)明的保護范圍應以權利要求的保護范圍為準。