PA054: Formální modely v systémové biologii David Šafránek 23.4.2010 Obsah Spojitý vs. diskrétní model Uvažujme modely tvaru M = S, R, reanet, ∅, map, rates (tzv. reakční sítě), rozšířené o komponentu rates : R → R+ (ohodnocení reakcí konstantami). spojitá sémantika diskrétní sémantika doména s M s con M ∈ R+ 0 s dsc M ∈ N0 molární koncentrace [M] počet molekul det. kinetická konstanta parametr exp. rozložení (CTMC) význam rates(r) rychlost reakce frekvence provedení reakce [s−1], [M−1 · s−1] 1 rates(r) . . . prům. čas mezi reakcemi Spojitý vs. diskrétní model abstrahovan modelovan cas promenne diskretni spojite kvalitativni model stochasticky model spojity modelaproximace abstrakce abstrakce Spojitý vs. diskrétní model abstrahovan modelovan cas promenne diskretni spojite kvalitativni model stochasticky model spojity model kvalitativni Petriho site LTL/CTL model checking aproximace abstrakce abstrakce Spojitý vs. diskrétní model abstrahovan modelovan cas promenne diskretni spojite kvalitativni model stochasticky model spojity model kvalitativni Petriho site LTL/CTL model checking stochasticke Petriho site PCTL/CSL model checking aproximace abstrakce abstrakce Spojitý vs. diskrétní model abstrahovan modelovan cas promenne diskretni spojite kvalitativni model stochasticky model spojity model kvalitativni Petriho site LTL/CTL model checking stochasticke Petriho site PCTL/CSL model checking spojite Petriho site diferencialni rovnice aproximace abstrakce abstrakce Spojitá interpretace míst a diskrétní aproximace • tokeny lze interpretovat v R+ 0 jako molární koncentrace [M] • pro pravděpodobnostní model checking nutno aproximovat: • předpokládejme molární koncentraci všech látek omezenou (interval 0, max) ⊂ R) • rovnoměrné rozdělení do N intervalů: 0, (0, 1 · max N , (1 · max N , 2 · max N , . . . , (N − 1 · max N , N · max N Kalibrace aproximace • spojitý průběh koncentrace pro p ∈ P můžeme porovnat se stochastickým (diskrétním) modelem: s con(t) = σ · N i=1 (i · Pr{ s dsc(t) = i}) kde: • σ ∈ R+ určuje jemnost rozdělení (délku intervalu) [M] • N je počet úrovní • s dsc (t) (resp. s con(t)) určuje aktuální diskrétní (resp. spojitou) úroveň v čase t Kalibrace aproximace • stochastický model lze převést na deterministický (spojitý) • pro ordinární sítě dostáváme pro hazardní funkci reakce t ∈ T: ht = kt · N · p∈•t ( m(p) N ) kde kt je deterministická kinetická konstanta Konverze/rozklad látky v čase A k −→ B • předpokládejme [A] nádoba obsahující nA molekul • kolik molekul se rozpadne/zkonvertuje v čase t? • hodnota přímo úměrná hodnotě nA v daném okamžiku − dnA(t) dt = k · nA(t) • koeficient úměrnosti je konstanta k [s−1 ] tzv. reakční konstanta (koeficient) - determinuje rychlost reakce Exponenciální rozklad/konverze A k −→ B −dnA(t) dt = k · nA(t) ⇔ nA(t) = nA(0) · e−kt • lineární dif. rce 1. řádu • jednoznačné řešení • numericky aproximovatelné • jednotka k [s−1] Exponenciální rozklad/konverze Diferenciální rovnice reakcí 1. řádu – konverze A k −→ B d[A] dt = −k[A] d[B] dt = k[A] Diferenciální rovnice reakcí 1. řádu – inflow k −→ B d[B] dt = k Diferenciální rovnice reakcí 1. řádu – rozklad A k −→ d[A] dt = −k[A] Diferenciální rovnice reakcí 2. řádu – syntéza A + B k −→ AB d[A] dt = −k[A][B] d[B] dt = −k[A][B] d[AB] dt = k[A][B] • zákon o aktivním působení hmoty (Law of Mass Action) • A i B musí mít řádově srovnatelné koncentrace • jednotka k [M−1 · s−1] Diferenciální rovnice reakcí 2. řádu – homeosyntéza A + A k −→ AA d[A] dt = −k[A]2 d[B] dt = −k[A]2 d[AB] dt = k[A]2 • zákon o aktivním působení hmoty (Law of Mass Action) • jednotka k [M−1 · s−1] Diferenciální rovnice reakcí 2. řádu – rozklad AB k −→ A + B d[A] dt = k[AB] d[B] dt = k[AB] d[AB] dt = −k[AB] Stochastický vs. deterministický model • molární koncentrace [M]: c = n V kde n je množství látky [mol], V je objem roztoku [l] • vyjadřuje se pomocí Avogadrovy konstanty (počet částic v 1 molu): c = N NA · V kde NA Avogadrova konstanta [mol−1], V objem roztoku [l] a N je počet molekul. • převodní faktor γ = NA · V : N = c · γ Stochastický vs. deterministický model Mějme stochastický model M = S, R, reanet, ∅, map, rates . Vytvoříme deterministický model M = S, R, reanet, ∅, map, rates : Pro reakci r ∈ R definujeme převodní vztah mezi rates(r) a det. kinetickou konstantou rates (r): typ reakce r ∈ R M → M A → B rates (r) = rates(r) A + B → AB rates (r) = rates(r) · γ A + A → AA rates (r) = rates(r)·γ 2 Formální sémantika modelu Klasický spojitý deterministický model Uvažujme modely tvaru M = S, R, reanet, ∅, map, rates (tzv. reakční sítě), rozšířené o komponentu rates : R → R+ (ohodnocení reakcí konstantami). • SVal = R+ ... koncentrace substance v médiu [M] • RVal = R+ ... reakční tok [M · s−1] ri M(t) = rates(ri ) · sj ∈In(ri ) ( sj M(t))|map( sj ,ri )| • mass action kinetics: • rates(ri ) je tzv. kinetická (reakční) konstanta pro reakci ri Formální sémantika modelu Klasický spojitý deterministický model • výpočet efektu reakce na reagující substance: r1 P2P1 Pm S2 SnS1 i1 in j2 jm i2 j1 • označme fluxi (t) = rates(ri ) · sj ∈In(ri )( sj M(t))|map( sj ,ri )| • effekt pro proměnné: d si M(t) dt = i,j ∈reanet map( i, j ) · fluxj (t) Spojité Petriho sítě Mějme Petriho síť N = P, T, f , m0 . Síť nazveme spojitou Petriho sítí, pokud: • P je konečná neprázdná množina míst (places), • T je konečná neprázdná množina přechodů (transitions), • f : ((P × T) ∪ (T × P)) → R+ 0 je množina orientovaných hran vážených celými čísly, • v : T → H je přiřazení kinetické funkce ht každému přechodu t ∈ T H := t∈T {ht|ht : R|•t| → R} • m0 : P → R+ 0 je iniciální označkování (marking). Kinetiku proměnných definujeme: dm(p) dt = t∈•p f (t, p)v(t) − t∈p• f (p, t)v(t) Obsah Stoichiometrická matice Uvažme model (reakční síť) M = S, R, reanet, ∅, map, rates . Definujeme stoichiometrickou matici M modelu M po prvcích: ∀si ∈ S, rj ∈ R. Mij = map( si , rj ), je-li definováno, 0, jinak. Model jako síť reakčních komplexů Uvažme model (reakční síť) M = S, R, reanet, ∅, map, rates . Zavedeme alternativní definici komponenty reanet tak, aby zachycovala relaci mezi reakcemi a příslušnými reakčními komplexy. Relaci budeme značit recnet, recnet ⊆ |S| i=1 N0 × |S| i=1 N0, a definujeme: • pro každé r ∈ R definujeme vstupní a výstupní komplex (po složkách i ∈ {1, ..., |S|}): cin(r)i = ( −K, reanet( si , r ) ∧ K ≡ map( si , r ) < 0 0, jinak. cout (r)i = ( K, reanet( si , r ) ∧ K ≡ map( si , r ) > 0 0, jinak. • recnet = { cin(r), cout (r) |r ∈ R ∧ cin(r) + cout (r) = 0} Množinu všech komplexů modelu M značíme complexes(M): complexes(M) = {cin(r)|r ∈ R} ∪ {cout(r)|r ∈ R} Třídy souvislosti Množinu L ⊂ complexes(M) nazýváme třídou souvislosti pokud splňuje následující podmínky: 1. c ∈ L ⇒ (∀c ∈ complexes(M). c, c ∈ recnet ⇒ c ∈ L) 2. c ∈ L ⇒ (∀c ∈ complexes(M). c , c ∈ recnet ⇒ c ∈ L) Počet tříd souvislosti modelu M značíme λM. M. Feinberg: Lectures on Chemical Reaction Networks http://www.che.eng.ohio-state.edu/~FEINBERG/LecturesOnReactionNetworks/ Třídy souvislosti – příklad complexes(M) = L1 ∪ L2 : L1 = {A1 + A2, A3, A4 + A5, A6} L2 = {2A1, A2 + A7, A8} ⇒ λM = 2 Třídy souvislosti – příklad MAPK λM = 6 Komplexnost modelu Nechť M = S, R, recnet, ∅, map, rates model (síť reakčních komplexů). Uvažujme jeho stoichiometrickou matici M. Definujeme komplexnost modelu M jako nezáporné číslo κ ∈ N0 definované vztahem: κ = |complexes(M)| − λM − h(M) kde h(M) je hodnost stoichiometrické matice Komplexnost modelu Nechť M = S, R, recnet, ∅, map, rates model (síť reakčních komplexů). Uvažujme jeho stoichiometrickou matici M. Definujeme komplexnost modelu M jako nezáporné číslo κ ∈ N0 definované vztahem: κ = |complexes(M)| − λM − h(M) kde h(M) je hodnost stoichiometrické matice Interpretace: Čím větší je κ, tím méně stoichiometrická matice vystihuje dynamiku modelu (= tím složitější je model). Komplexnost modelu – příklad |complexes(M)| = 7 λM = 2 h(M) = 5 ⇒ κ = 7 − 2 − 5 = 0 Třídy souvislosti – příklad MAPK |complexes(M)| = 26 λM = 6 h(M) = 15 ⇒ κ = 26 − 6 − 15 = 5 Hodnost stoichiometrické matice • dimenze stoichiometrického prostoru je dána h(M) • uvažme následující značení: • NN ∈ Zh(M)×|R| matice lineárně nezávislých řádků M • ND ∈ Z(|S|−h(M))×|R| matice lineárně závislých řádků M M = NN ND Hodnost stoichiometrické matice • dimenze stoichiometrického prostoru je dána h(M) • uvažme následující značení: • NN ∈ Zh(M)×|R| matice lineárně nezávislých řádků M • ND ∈ Z(|S|−h(M))×|R| matice lineárně závislých řádků M M = NN ND Definujeme spojovací (link) matici L ∈ Zh(M)×(|S|−h(M)) jako matici splňující vztah: ND = L · NN Význam hodnosti stoichiometrické matice • konzervace mas/energií — moiety conservation • definováno jako v čase konstantní součet koncentrací substrátů • např. ATP + ADP • zachyceno lineární závislostí řádků v M • každý substrát v ND je konzervován lineární kombinací substrátů v NN • spojovací matice zachycuje právě tuto závislost Podprostory stoichiometrické matice Definujeme nulový podprostor M, značíme np(M), jako prostor generovaný vektory splňujícími M · v = 0 • dim(np(M)) = |S| − h(M) • nulový prostor zachycuje stabilní distribuci reakčního toku (flux) • báze tohoto prostoru určuje tzv. elementární módy, které vymezují podsítě modelu se specifickou dynamikou Podprostory stoichiometrické matice Definujeme levý nulový podprostor M, značíme lnp(M), jako prostor generovaný vektory splňujícími MT · v = 0 • dim(lnp(M)) = |R| − h(M) • levý nulový prostor zachycuje konzervační a časové invarianty • všechny reakce zahrnuté v tomto prostoru manipulují s konzervovanou masou/energií Základní vlastnosti ODEs • uvažujme nejprve skalární rovnici dx(t) dt = f (x(t), t) • pokud f nezávisí na čase hovoříme o autonomním systému: dx(t) dt = f (x(t)) zjednodušeně lze psát: dx dt = f (x) • nullcline: {x|f (x) = 0} Základní vlastnosti ODEs Funkce f : D → R, D ⊂ R, je nazývána Lipschitz-spojitá právě když existuje κ > 0 t.ž. pro všechna x, y ∈ D platí: |f (x) − f (y)| ≤ κ|x − y| Picard-Lindelövova věta Nechť f : D → R je Lipschitz-spojitá a nechť x0 ∈ D. Potom existuje > 0 t.ž. iniciální problém dx dt = f (x), x(0) = x0 má právě jedno řešení x(t) pro všechna 0 ≤ t ≤ . Základní vlastnosti ODEs • libovolná dvě řešení pro libovolné dva různé iniciální problémy mají prázdný průnik • výsledek lze přímo zobecnit na systémy rovnic • zejména každá funkce, která je spojitá a v každém bodě existuje její derivace, je Lipschitz-spojitá Systémy (soustavy) obecných diferenciálních rovnic • dimenze systému n, uvažujeme stavový prostor X = Rn + • stav je vektor x(t) = x1(t), x2(t), ..., xn(t) T ∈ X • f (x(t)) = f1(x(t)), f2(x(t)), ..., fn(x(t)) T • předpokládáme ∀i. fi (x(t)) spojitá a diferencovatelná na X • soustava ODEs: d dt       x1 . . . xn       =       f1(x1, ..., xn) . . . fn(x1, ..., xn)       Stabilní stavy Mějme systém dx dt = f (x(t)). Stav ¯x ∈ X je stabilním stavem (ekvilibriem) systému, právě když platí f (¯x) = 0. Stabilní stavy Mějme systém dx dt = f (x(t)). Stav ¯x ∈ X je stabilním stavem (ekvilibriem) systému, právě když platí f (¯x) = 0. Typy stabilních stavů: 1. ¯x je stabilní ekvilibrium pokud pro každé > 0 existuje δ > 0 t.ž. řešení pro x0 = x0 1 , x0 2 , | x0 1 , x0 2 − ¯x1, ¯x2 | < δ, splňují | x1(t), x2(t) − ¯x1, ¯x2 | < pro libovolné t > 0. Stabilní stavy Mějme systém dx dt = f (x(t)). Stav ¯x ∈ X je stabilním stavem (ekvilibriem) systému, právě když platí f (¯x) = 0. Typy stabilních stavů: 1. ¯x je stabilní ekvilibrium pokud pro každé > 0 existuje δ > 0 t.ž. řešení pro x0 = x0 1 , x0 2 , | x0 1 , x0 2 − ¯x1, ¯x2 | < δ, splňují | x1(t), x2(t) − ¯x1, ¯x2 | < pro libovolné t > 0. 2. ¯x je nestabilní ekvilibrium pokud není stabilní. Stabilní stavy Mějme systém dx dt = f (x(t)). Stav ¯x ∈ X je stabilním stavem (ekvilibriem) systému, právě když platí f (¯x) = 0. Typy stabilních stavů: 1. ¯x je stabilní ekvilibrium pokud pro každé > 0 existuje δ > 0 t.ž. řešení pro x0 = x0 1 , x0 2 , | x0 1 , x0 2 − ¯x1, ¯x2 | < δ, splňují | x1(t), x2(t) − ¯x1, ¯x2 | < pro libovolné t > 0. 2. ¯x je nestabilní ekvilibrium pokud není stabilní. 3. ¯x je asymptoticky stabilní ekvilibrium pokud ¯x je stabilní a existuje δ > 0 t.ž. všechna řešení pro x0 = x0 1 , x0 2 , | x0 1 , x0 2 − ¯x1, ¯x2 | < δ, splňují limt→∞| x1(t), x2(t) − ¯x1, ¯x2 | = 0. Charakterizace stability v rovině Uvažujme následující soustavu dimenze 2: x1 = k1x1 x2 = k2x2 Tuto soustavu lze zapsat maticí (zn. M): dx dt = k1 0 0 k2 x1 x2 Řešení této soustavy pro iniciální problém x0 = x(0) je následující: x1(t) = x1(0)ek1t x2(t) = x2(0)ek2t Charakterizace stability v rovině Označme M = k1 0 0 k2 . Pro stabilní stav ¯x platí f (¯x) = 0, zde máme ¯x = 0, 0 . Stabilní stav lze blíže charakterizovat vlastními čísly M. Charakterizace stability v rovině Označme M = k1 0 0 k2 . Pro stabilní stav ¯x platí f (¯x) = 0, zde máme ¯x = 0, 0 . Stabilní stav lze blíže charakterizovat vlastními čísly M. Výpočet vlastních čísel λ = λ1, λ2 : det(M − λE) = 0 Charakterizace stability v rovině Označme M = k1 0 0 k2 . Pro stabilní stav ¯x platí f (¯x) = 0, zde máme ¯x = 0, 0 . Stabilní stav lze blíže charakterizovat vlastními čísly M. Výpočet vlastních čísel λ = λ1, λ2 : det(M − λE) = 0 ⇒ λ1 = k1, λ2 = k2 Charakterizace stability v rovině Označme M = k1 0 0 k2 . Pro stabilní stav ¯x platí f (¯x) = 0, zde máme ¯x = 0, 0 . Stabilní stav lze blíže charakterizovat vlastními čísly M. Výpočet vlastních čísel λ = λ1, λ2 : det(M − λE) = 0 ⇒ λ1 = k1, λ2 = k2 Charakterizace stability v rovině Nyní uvažme soustavu: dx dt = α β −β α x1 x2 Matice této soustavy má komplexně sdružená vlastní čísla: λ1 = α + βi λ2 = α − βi Řěšení této soustavy mají následující tvar (pro a, φ parametry): x(t) = aeαt cos(βt + φ) −sin(βt + φ) Charakterizace stability v rovině Obecné lineární systémy Uvažujme systém: d dt x1 x2 = a b c d x1 x2 a označme M = a b c d . Existuje invertibilní matice přechodu P t.ž.: y1 y2 = P−1 x1 x2 Označíme-li B = P−1MP, můžeme psát: d dt y1 y2 = B y1 y2 Obecné lineární systémy • transformace daná maticí přechodu P zachovává vlastní čísla • pokud M má dvě různá vlastní čísla λ1 = λ2, pak lze zvolit B = λ1 0 0 λ2 • pokud M má dvě komplexně sdružená vlastní čísla λ1 = ¯λ2 = α ± βi, pak lze zvolit B = α β −β α Klasifikace ekvilibrií lineárního systému Uvažujme systém: d dt x1 x2 = a b c d x1 x2 a označme M = a b c d . Definujeme stopu matice M, značíme tr(M), vztahem tr(M) = a + d. Platí: • tr(M) = λ1 + λ2 • det(M) = λ1λ2 λ1,2 = tr(M) 2 ± 1 2 (tr(M)2 − 4det(M) Klasifikace ekvilibrií lineárního systému • tr(M) < 0 ... nutná podmínka pro asymptotickou stabilitu • 0, 0 asymptoticky stabilní ⇔ všechna vlastní čísla mají záporné reálné složky. • 0, 0 asymptoticky stabilní ⇔ det(M) > 0 a tr(M) < 0.