Faculteit Wetenschappen Vakgroep Toegepaste Wiskunde en Informatica
Modellering van biochemische signaalnetwerken
Nele Delerue Promotor: Prof. Dr. W. Govaerts Begeleider: C. Sonck
Masterproef ingediend tot het behalen van de academische graad van master in de wiskunde, afstudeerrichting toegepaste wiskunde
Academiejaar 2010-2011
Voorwoord Deze pagina wil ik graag wijden aan alle mensen die me bijgestaan hebben tijdens het ontstaan van deze masterproef. Iedereen opnoemen, lijkt me wat te veel van het goede, maar ik wil toch enkele personen expliciet vermelden. Ik wil vooral een woord van dank richten aan mijn promotor en mijn begeleidster, namelijk Prof. Dr. W. Govaerts en C. Sonck. Zij reikten mij niet alleen een onderwerp aan, maar gaven ook sturing en ik had nooit het gevoel dat ik aan mijn lot overgelaten werd. Het was duidelijk dat ze er beide moeite in wilden steken en hun best deden om mijn redenering te volgen. Na een druk eerste semester vroegen ze zelf om elke twee weken af te spreken en ik weet zeker dat deze opgelegde druk heel belangrijk voor mij is om goed te presteren. Hiernaast getuigde dit ook van hun interesse in mijn vorderingen en het was aangenaam om te merken dat ze ook steeds effectief gelezen hadden wat ik gerealiseerd had. Dit ging van het controleren van de wiskunde tot het zelf implementeren van de verschillende modellen, het zelf opstellen van andere bruikbare redeneringen tot het corrigeren van spelling, taal en zinsbouw. Als ik het al eventjes minder zag zitten, dan kwam ik na onze tweewekelijkse afspraak vaak vol goede moed buiten, klaar om meteen verder te werken met de nieuw verworven inzichten. Het was heel aangenaam dat ze hielpen redeneren als ik dit nodig had. Ook Evgeni Nikolaev (Thomas Jefferson University, Jefferson Medical College, Philadelphia, PA, USA) wil ik bedanken om de artikels aan te reiken aan Prof. Dr. W. Govaerts als onderwerp voor een masterproef. Tot slot wil ik mijn zus, Veerle Delerue, bedanken. Hoewel ze niet steeds haar volledige enthousiasme kon behouden bij het lezen van een masterproef in de wiskunde, heeft ze doorgezet en mij zinnige opmerkingen bezorgd voor de afwerking van mijn masterproef.
De auteur geeft de toelating deze masterproef voor consultatie beschikbaar te stellen en delen van de masterproef te kopi¨eren voor persoonlijk gebruik. Elk ander gebruik valt onder de beperkingen van het auteursrecht, in het bijzonder met betrekking tot de verplichting de bron uitdrukkelijk te vermelden bij het aanhalen van resultaten uit deze masterproef. Gent, juni 2011
Nele Delerue
Inhoudsopgave Inleiding
1
1 Biochemische achtergrond 1.1 Activeringsenergie . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1.2 Hoe gaan enzymen te werk? . . . . . . . . . . . . . . . . . . . . . . . . . . . 1.3 Het proces van (de)fosforylatie . . . . . . . . . . . . . . . . . . . . . . . . . .
3 3 4 4
2 Signaalschakelaars en bistabiliteit door dubbele fosforylatie 2.1 M AP K . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2.1.1 Werking van een M AP K-cascade . . . . . . . . . . . . 2.1.2 Toepassingsgebied van de M AP K-cascade . . . . . . . 2.2 Biochemische verklaring van bistabiliteit . . . . . . . . . . . . 2.3 Beschrijving van een duale fosforylatie-defosforylatiecyclus . . 2.4 Reactievergelijkingen . . . . . . . . . . . . . . . . . . . . . . . 2.5 Uitdrukkingen voor de snelheid . . . . . . . . . . . . . . . . . 2.6 Interpretatie van de uitdrukkingen voor de snelheid . . . . . . 2.7 Eigenschappen van een duale fosforylatie-defosforylatiecyclus. 2.8 Implementatie . . . . . . . . . . . . . . . . . . . . . . . . . . . 2.8.1 Tijdsintegratie . . . . . . . . . . . . . . . . . . . . . . . 2.8.2 Vari¨eren van de MM-constanten Km1 en Km2 . . . . . 2.8.3 Vari¨eren van de parameters k1cat en k2cat . . . . . . . . . 2.8.4 Besluit voor bistabiliteit . . . . . . . . . . . . . . . . .
. . . . . . . . . . . . . .
. . . . . . . . . . . . . .
. . . . . . . . . . . . . .
. . . . . . . . . . . . . .
. . . . . . . . . . . . . .
. . . . . . . . . . . . . .
. . . . . . . . . . . . . .
. . . . . . . . . . . . . .
6 6 6 7 8 8 9 11 18 18 20 21 22 27 28
3 Signaalschakelaars en bistabiliteit door dubbele multisite fosforylatie 3.1 Beschrijving van een duale fosforylatie-defosforylatiecyclus met multisite-bindingen . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3.2 Reactievergelijkingen . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3.3 Uitdrukkingen voor de snelheid . . . . . . . . . . . . . . . . . . . . . . . . . 3.4 Implementatie . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3.4.1 Tijdsintegratie en continuatie . . . . . . . . . . . . . . . . . . . . . .
30
4 Analyse van bistabiliteit door dubbele fosforylatie 4.1 Twee-staps cycli: beschrijving van het model . . . . . . . . . . . . . . . . . . 4.1.1 Implementatie . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4.2 Bistabiliteit in twee-staps cycli: wiskundig . . . . . . . . . . . . . . . . . . .
38 38 41 42
i
30 31 32 34 35
Inhoudsopgave . . . . . . .
43 46 55 56 57 58 59
5 Bistabiliteit door driedubbele fosforylatie 5.1 Drie-staps cycli: beschrijving van het model . . . . . . . . . . . . . . . . . . 5.2 Bistabiliteit in drie-staps cycli: numeriek . . . . . . . . . . . . . . . . . . . .
61 61 65
6 Multistabiliteit door een cascade van dubbele fosforylaties 6.1 Twee-staps cycli in een cascade: beschrijving van het model . . . . . . . . . 6.2 Multistabiliteit in twee-staps cycli in een cascade: numeriek . . . . . . . . .
68 68 72
Algemeen besluit
80
Glossarium
84
A . . . . . .
86 86 89 90 90 90 92
. . . . . . . . . . .
92
. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
92 93 94
4.3
4.4
A.1 A.2 A.3 A.4 A.5 A.6
4.2.1 Bepalen van de evenwichtstoestand . . . . . . . . . . . . . 4.2.2 Bifurcatie-analyse . . . . . . . . . . . . . . . . . . . . . . . Bistabiliteit in twee-staps cycli: numeriek . . . . . . . . . . . . . . 4.3.1 Parameterrestricties voor een stabiele evenwichtstoestand . 4.3.2 Parameterrestricties die bistabiliteit toestaan . . . . . . . . 4.3.3 Overzicht van equilibria en de mogelijkheid tot bistabiliteit βmax . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
Appendix Appendix Appendix Appendix Appendix Appendix A.6.1 A.6.2 A.6.3 A.6.4
A B C D E F
. . . . . .
. . . . . .
. . . . . .
. . . . . .
. . . . . .
. . . . . .
. . . . . .
. . . . . .
. . . . . .
. . . . . .
. . . . . .
. . . . . .
. . . . . .
. . . . . .
. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . ∂β =0. . . . . . Oplossingen van de vergelijking ∂γ Bestaansvoorwaarden nagaan . . . . . . . . . . . Andere veranderlijken in βmax . . . . . . . . . . Berekende waarde voor β: minimaal of maximaal?
Bibliografie
. . . . . .
. . . . . .
. . . . . .
. . . . . .
. . . . . .
. . . . . .
. . . . . . .
. . . . . .
. . . . . . .
. . . . . .
. . . . . . .
. . . . . .
. . . . . . .
. . . . . .
. . . . . . .
. . . . . .
95
ii
Inleiding Als afstudeerrichting heb ik voor toegepaste wiskunde gekozen aangezien ik het heel interessant vind om te zien hoe de vaak abstracte wiskunde ook aangewend kan worden in het werkelijke leven, waar de functionaliteit ervan ligt. Ik hou ervan om te zien hoe wiskundigen een heel belangrijke rol kunnen spelen in aspecten van het leven waar de grote meerderheid van de bevolking dit niet zou verwachten. Het was dus belangrijk voor mij dat ik dit ook werkelijk kon toepassen in mijn masterproef en ’wiskunde met een nut’ kon produceren. Bijkomend moest er in de bacheloropleiding een minor gekozen worden en opteerde ik voor de minor biologie aangezien dit net als vele andere wetenschappen een vakgebied is dat me heel erg aanspreekt. Bij het overlopen van de mogelijke onderwerpen voor een masterproef was ik dan ook meteen ge¨ınteresseerd in het onderwerp van de modellering van biochemische signaalnetwerken. De toepassing ervan in ziektes zoals bijvoorbeeld kanker leek me heel interessant en was een extra stimulans om dit onderwerp te kiezen. Aangezien het vakgebied chemie al tot een ver verleden behoorde, schrok de term ’biochemisch’ me in eerste instantie een klein beetje af, maar aan de andere kant leek het me ook een uitdaging om me te verdiepen in deze combinatie van twee verschillende (en voor mij persoonlijk toch wat minder gekende) wetenschappen. Via deze masterproef willen we beter begrijpen hoe biologisch gedrag volgt uit de organisatie van eiwitten in cascades en netwerken. We weten namelijk dat chemische signaalnetwerken van regulerende eiwitten ervoor zorgen dat een cel informatie kan ontvangen van signalen van buitenaf en hierop kan reageren door zich te vermenigvuldigen, celdood te ondergaan,... Deze signaalnetwerken die via de eiwitten een signaal doorgeven, kunnen beschreven worden aan de hand van modellen bestaande uit een stelsel van differentiaalvergelijkingen. In het wiskundig computerprogramma Matcont kunnen we de modellen implementeren en zoeken naar evenwichtstoestanden van het model en eventuele bifurcaties van die evenwichtspunten bestuderen. Onze bedoeling is om eerst via numerieke experimenten, maar later ook analytisch, verschillende vereenvoudigde modellen van signaalnetwerken te bestuderen. Met behulp van dit onderzoek willen we onder andere aantonen dat bistabiliteit (waarbij het systeem kan afwisselen tussen twee stabiele toestanden, maar niet in een andere toestand kan blijven) en hysteresis (de stimulus moet een drempelwaarde overschrijden om het systeem naar een andere stabiele toestand te brengen, waar het kan blijven als de stimulus terug daalt) kunnen ontstaan uit een distributief kinetisch mechanisme van de twee-staps MAPK fosforylatie en defosforylatie. We voeren dit onderzoek uit op meerdere modellen en kijken ook of er eventuele multistabiliteit kan ontstaan (een systeem met meer dan twee stabiele evenwichtspunten). Op aanraden van E. Nikolaev zoeken we ook of de verschillende modellen
1
periodiek gedrag kunnen vertonen voor bepaalde parameterwaarden. In deze masterproef komen in de verschillende hoofdstukken de verschillende bestudeerde modellen aan bod. We hebben ervoor gekozen om de algemene structuur van elk hoofdstuk steeds in grote lijnen te behouden, zodat het onderscheid duidelijk aan het licht komt. Aangezien dit een masterproef is voor het behalen van het diploma van master in de wiskunde, ligt de nadruk meer op wiskunde dan op biochemie. Toch is het onmogelijk om het verloop en de opbouw te volgen zonder de achterliggende biochemie te begrijpen. Hier wordt dan ook aandacht aan besteed in de tekst en de biochemische termen worden nog eens samengevat in het glossarium dat teruggevonden kan worden op pagina 84. Dimensies bij de verschillende concentraties en parameterwaarden worden weggelaten in de tekst, maar staan vermeld in de tabellen met de waarden. Deze masterproef werd voornamelijk gebaseerd op [6] en [8]. Hierbij wil ik opmerken dat het bij het onder handen nemen van deze artikels duidelijk was dat de teksten niet door wiskundigen geschreven waren. De wiskunde was vaak onvolledig, weggestopt in de appendix en/of gewoon weggelaten, wat het voor mij zeker niet gemakkelijk maakte om de redeneringen steeds volledig te doorgronden. Resultaten die bekomen werden, stonden daarenboven vaak in een onlogische volgorde waarbij gevolgen uit bekomen uitdrukkingen vaak al voor de benodigde uitdrukkingen te vinden waren. Dit wijst er volgens mij op dat een wiskundige wel zorgde voor de berekeningen en deze bezorgde aan een biochemicus, maar deze de tekst schreef zonder verder overleg met de wiskundige. Het lijkt me dus van groot belang om een goede samenwerking te hebben tussen de verschillende vakspecialisten en je ook voldoende in te werken in het domein van de andere wetenschapper. Tijdens het opstellen van deze masterproef kwam ik ook meerdere problemen tegen. Zo werd er in de artikels vaak onvoldoende informatie gegeven over de gebruikte modellen en de parameterwaarden. Dit zorgde ervoor dat ik zelf modellen moest opstellen bij de kinetische diagrammen en ook zelf moest bepalen welke vereenvoudigingen ik invoerde. Samen met het feit dat er wel conclusies getrokken werden uit bifurcatiediagrammen horend bij deze modellen, maar deze figuren vaak niet expliciet gegeven werden, kon ik meermaals niet beschikken over controlemateriaal. Ook bij het gebruik van Matcont moest soms naar oplossingen gezocht worden. Zo detecteert dit programma sommige speciale (dichtgelegen) punten niet. Een oplossing hiervoor was om de continuatie van de evenwichtspunten uit te voeren met een kleinere tolerantie of kleinere staplengte. Als we te maken hadden met stijve problemen, konden we zien dat er zich kleine onregelmatigheden bevonden rond de equilibria. Dit werd opgelost door een numerieke methode te selecteren die speciaal gericht is op stijve problemen.
2
Hoofdstuk 1 Biochemische achtergrond Biochemische netwerken stellen de cel in staat om informatie te ontvangen, te verwerken en te reageren op die informatie. Dit is van levensbelang omdat de cel anders niet in staat is om zich te vermenigvuldigen, te differenti¨eren of celdood te ondergaan. De enzymen, een specifieke klasse van biologische katalysatoren, spelen een essenti¨ele rol in het reguleren van deze processen in cellen en organismen. Hoewel zij onmogelijke reacties niet mogelijk kunnen maken en geen invloed hebben op de positie van het equilibrium van een reactie, zijn ze nodig om reacties te versnellen en deze zo bruikbaar te maken voor het biologische systeem. Na het proces worden de enzymen onveranderd gerecupereerd.
1.1
Activeringsenergie
Voor alle reacties geldt dat er een energiebarri`ere moet doorbroken worden vooraleer zij kunnen uitgevoerd worden. De energie die nodig is om de barri`ere te doorbreken en een reactie in werking te stellen, wordt de activeringsenergie van de reactie genoemd.
Figuur 1.1: De invloed van enzymen op de activeringsenergie (uit [11]). De taak van de enzymen bestaat er uit de reactiesnelheid te be¨ınvloeden door de activeringsenergie te verlagen (zie Figuur 1.1). Daardoor verhogen ze het aantal moleculen dat 3
1.2. Hoe gaan enzymen te werk?
de overgangsenergie bereikt en dus de energiebarri`erre doorbreekt. Ze zijn niet in staat om onmogelijke reacties te katalyseren (het verschil in vrije energie tussen de substraten en het product blijft namelijk gelijk), maar ze versnellen een reactie in een bepaalde richting.
1.2
Hoe gaan enzymen te werk?
Reacties hangen vaak af van toevallige botsingen van substraten in de driedimensionale ruimte. Het is duidelijk dat dit sneller zou gebeuren als die botsing in de hand zou gewerkt worden en dit is de rol van de enzymen (zie Figuur 1.2): de enzymen binden de substraten op hun actieve site die zich meestal op het oppervlak van het enzym bevindt. Deze actieve sites hebben een complementaire vorm aan het substraat en dit betekent dat enzymen specifiek zijn. De substraten worden in een bepaalde richting vastgehouden en dicht bij elkaar gebracht, waardoor de botsingen niet meer toevallig zijn.
Figuur 1.2: Werking van een enzym via induced fit (uit [19]). Vooraleer de binding plaatsvindt, ondergaan zowel de substraten als de actieve sites een lichte verandering (induced fit). Dit heeft een gunstig effect op de reactiviteit van bepaalde chemische resten op de substraten zodat de reactie sneller kan doorgaan. Na het proces laat het enzym het product los en keert het terug naar zijn originele staat zodat het klaar is voor de katalyse van een nieuwe reactie.
1.3
Het proces van (de)fosforylatie
ATP is de afkorting van de Engelse term ’adenosine triphosphate’ waaruit meteen duidelijk is dat er zich in deze stof drie fosfaatgroepen bevinden. ATP is meer bepaald een nucleotide dat opgebouwd is uit meerdere nucle¨ınezuren (de purine adenine, de suiker ribose en de eerder vermelde drie fosfaatgroepen). ATP zorgt voor de uitwisseling van energie tussen cellen en is een energiemolecule voor biochemische reacties. Veel van deze reacties eisen namelijk een grote toevoeging van energie (zie activeringsenergie in paragraaf 1.1) en het grootste deel hiervan wordt voorzien door het uiteenvallen van ATP in ADP (adenosine diphosphate) en een vrije fosfaatgroep. Hierbij verkrijgen we heel wat vrije energie, maar dan is het ook meteen duidelijk dat ADP heel wat minder energie bevat dan ATP.
4
1.3. Het proces van (de)fosforylatie
De omgekeerde reactie is echter ook mogelijk: de hoogenergetische (en dus levensnoodzakelijke) stof ATP kan terug gevormd worden uit ADP. Dit gebeurt via een proces dat oxidatieve fosforylatie genoemd wordt en hierbij wordt opnieuw een binding gemaakt tussen ADP en een vrije fosfaatgroep. Het ligt voor de hand dat hiervoor chemische energie vereist is en deze ontstaat uit het afbreken van bepaalde moleculen (opslagmoleculen genoemd, bv. vetten, zetmeel,...). De hoge vrije energie die men bekomt uit ATP wordt bijvoorbeeld gebruikt in het biochemisch proces van de fosforylatie (zie Figuur 1.3). Bij het fosforyleren van een molecule (bv. een eiwit) wordt een fosfaatgroep (P O4 ) gebonden op bepaalde residuen van het molecule mits toevoeging van energie. Zoals beschreven in paragraaf 1.1 kunnen enzymen deze activeringsenergie verlagen en in het geval van de fosforylatie worden deze enzymen kinasen genoemd. Dit biochemisch proces is gemakkelijk reversibel: het eiwit kan zijn oorspronkelijke vorm terug aannemen nadat de fosfaatgroep verwijderd is. Dit wordt logischerwijze de defosforylatie genoemd en de enzymen die hier een handje toesteken worden fosfatasen genoemd.
Figuur 1.3: Het proces van de (de)fosforylatie van een eiwit. Dit proces van reversibele fosforylatie wordt heel vaak toegepast en is het meest algemene controlemechanisme in de cel. Het zorgt voor structurele veranderingen in het eiwit en is dan ook een belangrijk regulerend mechanisme. Eiwitenzymen worden op deze manier namelijk vaak geactiveerd of gedeactiveerd en door deze steeds veranderende cellulaire gedragingen ontstaan moleculaire wissels die ervoor zorgen dat de cel zich anders gaat gedragen. De activiteit van een eiwit wordt dus rechtstreeks bepaald door de activiteit van de kinasen en fosfatasen die erop inwerken. Zo vinden we vaak kinasen terug die op een specifieke manier georganiseerd zijn, namelijk als een fosforylatie-cascade. Hierbij wordt een eerste kinase geactiveerd door bijvoorbeeld het ontvangen van een bepaald signaal, waarna deze actieve kinase een volgende kinase in de cascade fosforyleert en activeert en dit gaat zo verder tot het signaal doorgestuurd wordt door de cel.
5
Hoofdstuk 2 Signaalschakelaars en bistabiliteit door dubbele fosforylatie Dit hoofdstuk is gebaseerd op [6].
2.1
M AP K
M AP K staat voor ’mitogen-activated protein kinase’. Een mitogen is een stof die de mitose stimuleert en de mitose is op zijn beurt een onderdeel van de celcyclus. De mitose is namelijk de eigenlijke celdeling. ERK (Extracellular signal-regulated kinase) wordt vaak als synoniem gebruikt voor M AP K, maar is in feite een specifieke, dierlijke M AP K. M AP K wordt geactiveerd binnen een eiwit-kinase-cascade die toepasselijk de naam ’M AP K-cascade’ draagt.
2.1.1
Werking van een M AP K-cascade
De M AP K-cascade is een ketting van eiwitfosforylaties in de cel die een signaal van een receptor, gebonden door een groeigen op het celoppervlak, doorgeeft aan het DNA in de celkern. Elk tussenliggend enzym fungeert door zijn werking als ’aan’- en ’af’-knop voor het nabijgelegen eiwit. Dit verloopt via een ingewikkeld proces, waarvan de belangrijkste zaken hier op een rijtje gezet worden (zie ook Figuur 2.1). Elke M AP K-cascade bestaat uit drie kinases, namelijk M AP kinase (M AP K of ERK), M AP kinase kinase (M AP KK, zoals bijvoorbeeld M EK) en M AP kinase kinase kinase (M AP KKK zoals bijvoorbeeld Raf, M KP 3), die geactiveerd worden in een reeks. Een geactiveerde M AP KKK fosforyleert M AP KK op diens serine- en threonine-residuen, dewelke op zijn beurt zorgt voor de activering van M AP K door fosforylatie op zijn threonineen tyrosine-residuen. Om de M AP KKK geactiveerd te krijgen, worden enkele stappen uitgevoerd. Het proces begint met een groeigen dat bindt op een receptor die zich op het oppervlak van de cel bevindt. Deze receptor wordt zo in staat gesteld om zichzelf te fosforyleren op zijn tyrosineresidu. In de gefosforyleerde staat interageert de receptor met een eiwit (GRB2) dat een complex vormt met het eiwit SOS (Son Of Sevenless), dat op zijn beurt het G-eiwit RAS
6
2.1. M AP K
Figuur 2.1: Een vereenvoudigde versie van de M AP K/ERK baan (uit [20]). activeert door de GDP-GTP uitwisseling te stimuleren. De geactiveerde RAS zorgt dan dat de activering van M AP KKK plaatsvindt, zodat de eerder besproken opeenvolging van fosforylaties kan uitgevoerd worden. De geactiveerde M AP K wordt verplaatst binnenin de kern waar ze transcriptiefactoren activeert en fosforyleert die uiteindelijk de genexpressie controleren.
2.1.2
Toepassingsgebied van de M AP K-cascade
M AP K-cascades zorgen voor het overbrengen van extracellulaire signalen naar de celkern waar dan op hun beurt specifieke genen geactiveerd worden die zorgen voor de celdifferentiatie, -deling en -groei. Mutaties van ´e´en van de eiwitten in de baan kan zorgen voor een permanent verblijf in de aan-/af-staat van de activering van een bepaald eiwit. Dit ontregelt het volledige systeem en onderzoek wijst erop dat dit vaak aan de basis ligt van ongecontroleerde celgroei en het ontstaan van oncogenen, wat op zijn beurt zorgt voor het ontstaan van kanker. Verder onderzoek naar deze cascades zou goede manieren aan het licht kunnen brengen om deze permanente staten om te keren. Het onderzoek naar deze cascades is bijgevolg volop aan de gang om geneesmiddelen tegen kanker te ontwikkelen.
7
2.2. Biochemische verklaring van bistabiliteit
2.2
Biochemische verklaring van bistabiliteit
In dit en het volgende hoofdstuk zullen we aantonen dat bistabiliteit volgt uit competitieve inhibitie en verzadiging, meer precies uit: • substraatverzadiging van minstens ´e´en van de twee enzymen, • competitieve inhibitie van de tweede modificatiestap door het substraat van de eerste stap.
2.3
Beschrijving van een duale fosforylatie-defosforylatiecyclus
Zoals beschreven bij de werking van een M AP K-cascade (zie paragraaf 2.1.1) hebben we in de M AP K-cascade te maken met de verschillende vormen M AP K, M AP KK en M AP KKK. Houden we rekening met het feit dat M AP K op twee verschillende residuen kan gefosforyleerd worden (namelijk het tyrosine- of het threonine-residu), dan zien we dat er vier vormen aangenomen worden. Als we M AP K afkorten tot M, dan hebben we de volgende vier gedaantes: • M, • M pY : de gefosforyleerde vorm van M waar de fosfaatgroep op het tyrosine-residu gehecht is, • M pT : de gefosforyleerde vorm van M waar de fosfaatgroep op het threonine-residu gehecht is, • M pp: de dubbel gefosforyleerde vorm van M (fosfaatgroep op beide residuen). In deze paragraaf vereenvoudigen we dit model echter en houden we geen rekening met de bindingsplaats van de fosfaatgroep. We beschouwen dus slechts ´e´en intermediair en werken met M , M p en M pp. Uit onderzoek blijkt dat de (de)fosforylatie via een distributief systeem gebeurt. Om de dubbel gefosforyleerde vorm te bekomen, moet er met andere woorden twee keer een binding met het enzym plaatsvinden. Dit in tegenstelling tot processieve (of niet-distributieve) systemen waarbij een enzym gebonden blijft aan zijn substraat nadat het ´e´en keer zijn werking heeft uitgevoerd en dus meerdere malen zijn functie kan vervullen terwijl er maar ´e´en binding van enzym en substraat vereist is.
8
2.4. Reactievergelijkingen
Figuur 2.2: Cyclus van M AP K. In Figuur 2.2 is de cyclus van M AP K weergegeven. We merken hierbij op dat de kinase M AP KK zorgt voor de fosforylatie van zowel M als van M p en de fosfatase M KP 3 voor de defosforylatie van M pp en M p. De snelheden van deze fosforylaties zijn respectievelijk v1 en v2 , die van de defosforylaties respectievelijk v3 en v4 . In dit model nemen we eveneens aan dat de som van de concentraties van alle vormen van M AP K constant is. We defini¨eren deze totale concentratie als Mtot en hebben dan Mtot = [M ] + [M p] + [M pp] = constant. Samen met Figuur 2.2 kunnen we het volgende stelsel van vergelijkingen opstellen: d[M ] = v4 − v1 , dt d[M pp] (2.1) = v2 − v3 , dt [M p] = Mtot − [M ] − [M pp], waarbij de afgeleiden naar de tijd beschouwd worden.
2.4
Reactievergelijkingen
Uit Figuur 2.2 en het proces waarbij het enzym inwerkt op zijn substraten (zie Hoofdstuk 1), volgen de stappen in de fosforylatie zoals weergegeven in Figuur 2.3 (met P een fosfaatgroep). We houden hierbij rekening met het feit dat dit een distributief proces is.
Figuur 2.3: Symbolische weergave van de fosforylatie van M . Een analoge figuur kan opgesteld worden voor de fosforylatie van M p. We nemen aan dat er een constante hoeveelheid ATP/ADP in overvloed aanwezig is en dat de eerste en laatste 9
2.4. Reactievergelijkingen
stap omkeerbaar zijn. Het equilibrium in die stappen ligt dus niet specifiek langs ´e´en van de twee kanten van de reactie. Voor de middelste stap nemen we aan dat de condities zo zijn dat de omgekeerde reactie verwaarloosbaar is. Dit wil zeggen dat het equilibrium zich duidelijk rechts van die pijl bevindt in de reactie. We hebben dan dus M + M AP KK ↔ M − M AP KK ∗ → M p − M AP KK ↔ M p + M AP KK, M p + M AP KK ↔ M p − M AP KK ∗ → M pp − M AP KK ↔ M pp + M AP KK. We gebruiken hier een ∗ als notatie voor een complex waar substraat en enzym gebonden zijn, maar het enzym zijn taak nog niet vervuld heeft en dus nog actief is. Doordat het proces distributief is, weten we dat het enzym zijn activiteit verliest na ´e´en inwerking, waardoor het sterretje in de notatie verdwijnt bij het nieuw gevormde complex. De defosforylatie van M pp en M p m.b.v. het enzym M KP 3 is analoog en is schematisch weergegeven in Figuur 2.4.
Figuur 2.4: Symbolische weergave van de defosforylatie van M p. We beschouwen dezelfde (on)omkeerbaarheid en bekomen M pp + M KP 3 ↔ M pp − M KP 3∗ → M p − M KP 3 ↔ M p + M KP 3, M p + M KP 3 ↔ M p − M KP 3∗ → M − M KP 3 ↔ M + M KP 3. We hebben nu vier reactievergelijkingen. Voor de fosforylatie zullen we het systeem echter vereenvoudigen. We zullen namelijk aannemen dat de tussenvormen M p − M AP KK en M pp − M AP KK niet expliciet gevormd worden, maar dat het actieve complex na inwerking van het enzym meteen uiteenvalt in het oorspronkelijk enzym en de gefosforyleerde vorm M p, respectievelijk M pp. We werken dus met de volgende vier reactievergelijkingen k1 k2 ∗
→ M + M AP KK M − M AP KK M p + M AP KK, k−1 (2.2) k3 k4 M p + M AP KK M p − M AP KK ∗ → M pp + M AP KK, k−3 h1 h2 h3 ∗ M pp + M KP 3 M pp − M KP 3 → M p − M KP 3 M p + M KP 3, h−1 h−3 (2.3) h4 h5 h6 M p + M KP 3 M p − M KP 3∗ → M − M KP 3 M + M KP 3, h−4 h−6 10
2.5. Uitdrukkingen voor de snelheid
waarbij de ki ’s en de hi ’s de snelheidsconstanten van de reacties voorstellen. Deze constanten zijn een maat voor hoe snel de reactie verloopt.
2.5
Uitdrukkingen voor de snelheid
Aan de hand van het model dat we tot nu toe opgesteld hebben en de gestelde voorwaarden, willen we de Michaelis-Menten (MM) vergelijkingen (zie Glossarium) afleiden voor de verschillende snelheden die weergegeven zijn in Figuur 2.2. Als eerste gebruiken we hiervoor de reactievergelijkingen in (2.2), afgeleid voor de fosforylatie van M AP K door het inwerken van de kinase M AP KK. We willen nu een MM-vergelijking opstellen voor de snelheid v1 . Uit Figuur 2.2 kunnen we duidelijk afleiden dat dit de snelheid is waarmee M omgezet wordt naar zijn gefosforyleerde vorm door de inwerking van M AP KK. Samen met (2.2) krijgen we de volgende uitdrukking: v1 = k2 [M − M AP KK ∗ ]. Het is dus duidelijk dat we een uitdrukking nodig hebben voor [M − M AP KK ∗ ]. De kinase M AP KK zit gedurende de cyclus verspreid in gebonden en ongebonden vormen. Voor de totale concentratie van deze kinase volgt uit Figuur 2.2 de volgende uitdrukking: [M AP KK]tot = [M AP KK] + [M − M AP KK ∗ ] + [M p − M AP KK ∗ ].
(2.4)
Hierbij nemen we aan dat deze totale concentratie constant is en uitdrukking (2.4) zullen we verder nodig hebben in de afleiding1 . We nemen nu aan dat er evenwicht is tussen de opbouw en de afbraak van de intermediairen M − M AP KK ∗ en M p − M AP KK ∗ in (2.2). Dit is analoog aan de gebruikelijke benadering die men maakt bij het opstellen van de klassieke MM-vergelijkingen, cf. [2]. Uit de eerste reactievergelijking krijgen we dan de volgende uitdrukking: k1 [M ][M AP KK] = k−1 [M − M AP KK ∗ ] + k2 [M − M AP KK ∗ ] k1 [M ][M AP KK] k−1 + k2 1 = [M ][M AP KK], Km1
⇒ [M − M AP KK ∗ ] =
(2.5)
k−1 + k2 gesteld hebben. Door het evenwicht in de tweede reactievergek1 lijking van (2.2) hebben we analoog waarbij we Km1 =
k3 [M p][M AP KK] = k−3 [M p − M AP KK ∗ ] + k4 [M p − M AP KK ∗ ] 1
We willen nog opmerken dat we gedurende het opstellen van de vier MM vergelijkingen bepaalde parameters invoeren die afhankelijk zijn van de oorspronkelijke parameters. Op het einde van deze paragraaf worden deze nog eens expliciet weergegeven.
11
2.5. Uitdrukkingen voor de snelheid k3 [M p][M AP KK] k−3 + k4 1 = [M p][M AP KK], Km2
⇒ [M p − M AP KK ∗ ] =
(2.6)
k−3 + k4 gesteld hebben. Aangezien we een uitdrukking voor [M − k3 M AP KK ∗ ] willen bekomen en we gebruik kunnen maken van de uitdrukking (2.4) voor [M AP KK]tot (waarbij deze laatste een constante is), drukken we [M AP KK] en [M p − M AP KK ∗ ] uit in termen van [M − M AP KK ∗ ]. Op deze manier krijgen we uit onze vergelijking (2.4) een uitdrukking voor [M − M AP KK ∗ ] en bijgevolg een uitdrukking voor v1 . Uit (2.5) krijgen we dus waarbij we Km2 =
[M AP KK] = Km1
[M − M AP KK ∗ ] , [M ]
wat voor (2.6) betekent dat [M p − M AP KK ∗ ] =
Km1 [M p][M − M AP KK ∗ ] · . Km2 [M ]
We brengen nu alles samen in (2.4) en krijgen dan Km1 [M − M AP KK ∗ ] + [M − M AP KK ∗ ] [M ] Km1 [M p][M − M AP KK ∗ ] + · Km2 [M ] Km1 [M p] Km1 +1+ = [M − M AP KK ∗ ] [M ] Km2 [M ] Km1 Km2 + Km2 [M ] + Km1 [M p] [M − M AP KK ∗ ] = Km2 [M ]
[M AP KK]tot =
⇒ [M − M AP KK ∗ ] =
Km2 [M ][M AP KK]tot . Km1 Km2 + Km2 [M ] + Km1 [M p]
We bekomen dus de volgende MM-vergelijking aan de hand van de invoering k1cat = k2 (dit wordt om notatieredenen ingevoerd naar analogie met latere afleidingen): v1 = k2 [M − M AP KK ∗ ] = k1cat [M − M AP KK ∗ ] k1cat Km2 [M ][M AP KK]tot = Km1 Km2 + Km2 [M ] + Km1 [M p] k cat [M ][M AP KK]tot /Km1 = 1 . 1 + [M ]/Km1 + [M p]/Km2
(2.7)
Een analoge redenering is gevolgd bij het opstellen van de MM-vergelijking voor v2 . Uit 12
2.5. Uitdrukkingen voor de snelheid
Figuur 2.2 blijkt duidelijk dat dit de snelheid is waarmee de gefosforyleerde vorm van M AP K (M p) omgezet wordt naar zijn dubbel-gefosforyleerde vorm (M pp) door de inwerking van M AP KK. Samen met (2.2) krijgen we dus de volgende uitdrukking: v2 = k4 [M p − M AP KK ∗ ]. Het is duidelijk dat we een uitdrukking nodig hebben voor [M p − M AP KK ∗ ]. We nemen opnieuw aan dat er evenwicht is tussen de opbouw en de afbraak van de intermediairen M − M AP KK ∗ en M p − M AP KK ∗ in (2.2) en we kunnen dus opnieuw met de uitdrukkingen (2.5) en (2.6) werken. Aangezien we nu een uitdrukking voor [M p−M AP KK ∗ ] willen bekomen, vormen we deze vergelijkingen om tot de volgende: [M p − M AP KK ∗ ] , [M p] Km2 [M ][M p − M AP KK ∗ ] . [M − M AP KK ∗ ] = Km1 [M p] [M AP KK] = Km2
Brengen we opnieuw alles samen in (2.4), dan krijgen we [M p − M AP KK ∗ ] Km2 [M ][M p − M AP KK ∗ ] + · [M p] Km1 [M p] ∗ +[M p − M AP KK ] Km2 Km2 [M ] = + + 1 [M p − M AP KK ∗ ] [M p] Km1 [M p] Km1 Km2 + Km2 [M ] + Km1 [M p] = [M p − M AP KK ∗ ] Km1 [M p]
[M AP KK]tot =Km2
⇒ [M p − M AP KK ∗ ] =
Km1 [M p][M AP KK]tot . Km1 Km2 + Km2 [M ] + Km1 [M p]
We bekomen dus de volgende MM-vergelijking (met de invoering van k2cat = k4 ): v2 = k4 [M p − M AP KK ∗ ] = k2cat [M p − M AP KK ∗ ] k2cat Km1 [M p][M AP KK]tot = Km1 Km2 + Km2 [M ] + Km1 [M p] k cat [M p][M AP KK]tot /Km2 = 2 . 1 + [M ]/Km1 + [M p]/Km2 Na het opstellen van het schema (2.2) voor de fosforylatie van M AP K door het inwerken van de kinase M AP KK en de bijhorende MM-vergelijkingen af te leiden, stellen we nu analoog de MM-vergelijkingen op horende bij (2.3) (de defosforylatie van M pp en M p door het inwerken van de fosfatase M KP 3). We willen een MM-vergelijking opstellen voor de snelheid v3 . Uit Figuur 2.2 kunnen we duidelijk zien dat dit de snelheid is waarmee M pp omgezet wordt naar zijn gedefosforyleerde vorm M p door de inwerking van M KP 3. Samen met (2.3) krijgen we de volgende uitdrukking: v3 = h2 [M pp − M KP 3∗ ], 13
2.5. Uitdrukkingen voor de snelheid
waarbij we aangenomen hebben dat de concentratie van [M p−M KP 3] steeds constant blijft. Het is duidelijk dat we ditmaal een uitdrukking nodig hebben voor [M pp − M KP 3∗ ]. De fosfatase M KP 3 zit gedurende de cyclus verspreid in zowel gebonden als ongebonden vormen. We hebben dus (zie (2.3)) [M KP 3]tot =[M KP 3] + [M − M KP 3] + [M p − M KP 3] + [M p − M KP 3∗ ] + [M pp − M KP 3∗ ],
(2.8)
waarbij we opnieuw veronderstellen dat de totale concentratie van deze fosfatase constant is. We nemen nu aan dat er evenwicht is tussen de opbouw en de afbraak van de intermediairen M pp − M KP 3∗ , M p − M KP 3, M p − M KP 3∗ en M − M KP 3 in (2.3). Uit de eerste reactievergelijking van (2.3) krijgen we de volgende uitdrukkingen: h1 [M pp][M KP 3] = h−1 [M pp − M KP 3∗ ] + h2 [M pp − M KP 3∗ ] h1 [M pp][M KP 3] ⇒ [M pp − M KP 3∗ ] = h−1 + h2
(2.9)
en h2 [M pp − M KP 3∗ ] + h−3 [M p][M KP 3] = h3 [M p − M KP 3].
(2.10)
Door het evenwicht in de tweede reactievergelijking van (2.3) hebben we analoog h4 [M p][M KP 3] = h−4 [M p − M KP 3∗ ] + h5 [M p − M KP 3∗ ] h4 ⇒ [M p − M KP 3∗ ] = [M p][M KP 3] h−4 + h5
(2.11)
en h5 [M p − M KP 3∗ ] + h−6 [M ][M KP 3] = h6 [M − M KP 3].
(2.12)
Aangezien we een uitdrukking voor [M pp − M KP 3∗ ] willen bekomen en we gebruik kunnen maken van de uitdrukking (2.8) voor [M KP 3]tot (waarbij deze laatste een constante is), drukken we [M KP 3], [M − M KP 3], [M p − M KP 3] en [M p − M KP 3∗ ] uit in termen van [M pp − M KP 3∗ ]. Op deze manier krijgen we uit onze vergelijking (2.8) een uitdrukking voor [M pp − M KP 3∗ ] en bijgevolg een uitdrukking voor v3 . Uit (2.9) krijgen we dus [M KP 3] =
h−1 + h2 [M pp − M KP 3∗ ] , h1 [M pp]
wat voor (2.11) betekent dat [M p − M KP 3∗ ] =
(h−1 + h2 )h4 [M p] [M pp − M KP 3∗ ]. (h−4 + h5 )h1 [M pp]
Uit (2.10) volgt dat h2 h−3 [M pp − M KP 3∗ ] + [M p][M KP 3] h h3 3 h2 (h−1 + h2 )h−3 [M p] = + [M pp − M KP 3∗ ] h3 h1 h3 [M pp]
[M p − M KP 3] =
14
2.5. Uitdrukkingen voor de snelheid
en uit (2.12) dat h5 h−6 [M p − M KP 3∗ ] + [M ][M KP 3] h6 h6 (h−1 + h2 )h4 h5 [M p] (h−1 + h2 )h−6 [M ] = + [M pp − M KP 3∗ ]. (h−4 + h5 )h1 h6 [M pp] h1 h6 [M pp]
[M − M KP 3] =
We brengen nu alles samen in (2.8) en krijgen dan (h−1 + h2 )h4 h5 [M p] (h−1 + h2 )h−6 [M ] h2 h−1 + h2 + + + [M KP 3]tot = h1 [M pp] (h−4 + h5 )h1 h6 [M pp] h1 h6 [M pp] h3 (h−1 + h2 )h−3 [M p] (h−1 + h2 )h4 [M p] = + + + 1 [M pp − M KP 3∗ ] h1 h3 [M pp] (h−4 + h5 )h1 [M pp] N3 = (h−4 + h5 )h1 h3 h6 [M pp] met N3 =(h−1 + h2 )(h−4 + h5 )h3 h6 + (h−1 + h2 )(h−4 + h5 )h3 h−6 [M ] + (h−1 + h2 ) h3 h4 h5 + (h−4 + h5 )h−3 h6 + h3 h4 h6 [M p] + (h−4 + h5 )h1 h6 (h2 + h3 )[M pp] =(h−1 + h2 )(h−4 + h5 )h3 h6 · h3 h4 h5 + (h−4 + h5 )h−3 h6 + h3 h4 h6 h1 (h2 + h3 ) h−6 [M ] + [M p] + [M pp] . 1+ h6 (h−4 + h5 )h3 h6 (h−1 + h2 )h3 We bekomen dus de MM-vergelijking v3 = h2 [M pp − M KP 3∗ ] (h−4 + h5 )h1 h2 h3 h6 [M pp][M KP 3]tot = N3 h1 h2 [M pp][M KP 3]tot (h−1 +h2 ) = +h5 )h−3 h6 +h3 h4 h6 1 + hh−6 [M ] + h3 h4 h5 +(h(h−4 [M p] + 6 −4 +h5 )h3 h6 =
h1 (h2 +h3 ) [M pp] (h−1 +h2 )h3
k3cat [M pp][M KP 3]tot /Km3 , 1 + [M ]/Km5 + [M p]/Km4 + [M pp]/Km3
waarbij we Km3 =
h−1 + h2 , Km4 = h1 + h1 h2 /h3
k−4 + k5 , Km5 = h4 1 + h5 /h6 + h−3 (h−4 + h5 )/(h3 h4 )
h6 h2 en k3cat = gesteld hebben. h−6 1 + h2 /h3 Analoog willen we nu een MM-vergelijking opstellen voor de snelheid v4 . Uit Figuur 2.2 kunnen we duidelijk zien dat dit de snelheid is waarmee M p omgezet wordt naar zijn gedefosforyleerde vorm M door de inwerking van M KP 3. Samen met (2.3) krijgen we de volgende uitdrukking: v4 = h5 [M p − M KP 3∗ ], 15
2.5. Uitdrukkingen voor de snelheid
waarbij we aangenomen hebben dat de concentratie van [M − M KP 3] steeds constant blijft. Het is duidelijk dat we een uitdrukking nodig hebben voor [M p − M KP 3∗ ]. We nemen opnieuw aan dat er evenwicht is tussen de opbouw en de afbraak van de intermediairen M pp − M KP 3∗ , M p − M KP 3, M p − M KP 3∗ en M − M KP 3 in (2.3) en we kunnen dus opnieuw met de vier uitdrukkingen (2.9) - (2.12) werken. Aangezien we nu een uitdrukking voor [M p − M KP 3∗ ] willen bekomen, vormen we deze vergelijkingen om tot de volgende: h−4 + h5 [M − M KP 3∗ ] [M KP 3] = , h4 [M p] (h−4 + h5 )h1 [M pp] [M p − M KP 3∗ ], [M pp − M KP 3∗ ] = (h−1 + h2 )h4 [M p] (h−4 + h5 )h1 h2 [M pp] (h−4 + h5 )h−3 [M p − M KP 3] = + [M p − M KP 3∗ ], (h−1 + h2 )h3 h4 [M p] h3 h4 h5 (h−4 + h5 )h−6 [M ] [M − M KP 3] = + [M p − M KP 3∗ ]. h6 h4 h6 [M p] We brengen alles samen in (2.8) en krijgen dan h−4 + h5 h5 (h−4 + h5 )h−6 [M ] (h−4 + h5 )h1 h2 [M pp] [M KP 3]tot = + + + h4 [M p] h6 h4 h6 [M p] (h−1 + h2 )h3 h4 [M p] (h−4 + h5 )h−3 (h−4 + h5 )h1 [M pp] = + +1+ [M p − M KP 3∗ ] h3 h4 (h−1 + h2 )h4 [M p] N4 = (h−1 + h2 )h3 h4 h6 [M p] met N4 = (h−1 + h2 )(h−4 + h5 )h3 h6 + (h−1 + h2 )(h−4 + h5 )h3 h−6 [M ] = + (h−1 + h2 ) h3 h4 h5 + (h−4 + h5 )h−3 h6 + h3 h4 h6 [M p] = + (h−4 + h5 )h1 h6 (h2 + h3 )[M pp] =(h−1 + h2 )(h−4 + h5 )h3 h6 · h−6 h3 h4 h5 + (h−4 + h5 )h−3 h6 + h3 h4 h6 h1 (h2 + h3 ) 1+ [M ] + [M p] + [M pp] . h6 (h−4 + h5 )h3 h6 (h−1 + h2 )h3 We bekomen dus de volgende MM-vergelijking v4 = h5 [M p − M KP 3∗ ] (h−1 + h2 )h3 h4 h5 h6 [M p][M KP 3]tot = N4 h4 h5 [M p][M KP 3]tot (h−4 +h5 ) = +h5 )h−3 h6 +h3 h4 h6 1 + hh−6 [M ] + h3 h4 h5 +(h(h−4 [M p] + 6 −4 +h5 )h3 h6 =
k4cat [M p][M KP 3]tot /Km4 , 1 + [M ]/Km5 + [M p]/Km4 + [M pp]/Km3 16
h1 (h2 +h3 ) [M pp] (h−1 +h2 )h3
2.5. Uitdrukkingen voor de snelheid
waarbij we k4cat =
h5 gesteld hebben. 1 + h5 /h6 + h−3 (h−4 + h5 )/(h3 h4 )
Samengevat hebben we de volgende vier MM-vergelijkingen: k1cat [M ][M AP KK]tot /Km1 , 1 + [M ]/Km1 + [M p]/Km2 k cat [M p][M AP KK]tot /Km2 , v2 = 2 1 + [M ]/Km1 + [M p]/Km2 k3cat [M pp][M KP 3]tot /Km3 v3 = , 1 + [M ]/Km5 + [M p]/Km4 + [M pp]/Km3 k4cat [M p][M KP 3]tot /Km4 v4 = , 1 + [M ]/Km5 + [M p]/Km4 + [M pp]/Km3 v1 =
met k1cat = k2 , k2cat = k4 , k−1 + k2 , Km1 = k1 k−3 + k4 Km2 = , k3 h2 k3cat = , 1 + h2 /h3 h5 , 1 + h5 /h6 + h−3 (h−4 + h5 )/(h3 h4 ) h−1 + h2 = , h1 + h1 h2 /h3 k−4 + k5 , = h4 1 + h5 /h6 + h−3 (h−4 + h5 )/(h3 h4 )
k4cat = Km3 Km4
Km5 =
h6 . h−6
Ons model dat hier besproken wordt, ziet er dan uit als volgt:
17
(2.13) (2.14) (2.15) (2.16)
2.6. Interpretatie van de uitdrukkingen voor de snelheid
Model 1. d[M ] = v4 − v1 , dt d[M pp] = v2 − v3 , dt [M p] = Mtot − [M ] − [M pp], k cat [M ][M AP KK]tot /Km1 , v1 = 1 1 + [M ]/Km1 + [M p]/Km2 k cat [M p][M AP KK]tot /Km2 v2 = 2 , 1 + [M ]/Km1 + [M p]/Km2 k3cat [M pp][M KP 3]tot /Km3 , v3 = 1 + [M ]/Km5 + [M p]/Km4 + [M pp]/Km3 k4cat [M p][M KP 3]tot /Km4 v4 = . 1 + [M ]/Km5 + [M p]/Km4 + [M pp]/Km3
2.6
Interpretatie van de uitdrukkingen voor de snelheid
We hebben in ons model bij het defosforylatieschema de werkelijke desfosforylatiestappen beschouwd, evenals de stap waar het nieuw gevormde product terug opgesplitst werd. Bij het fosforylatieschema hebben we dit vereenvoudigd door deze twee stappen samen te nemen als ´e´en stap en te veronderstellen dat het product meteen uiteenvalt nadat het enzym zijn werking uitgevoerd heeft. We merken uit onze MM-vergelijkingen dat dit niet echt zorgt voor een fundamentele verandering in de uitdrukkingen. We zien echter wel dat er een extra term tevoorschijn komt in de noemer van de MM-vergelijkingen als we de extra stap beschouwen in de defosforylatie en er aangenomen wordt dat de stap van het splitsen van product en enzym wederkerig is. De katalytische constanten kicat en de MM-constanten Kmi vertonen eveneens een wat ingewikkeldere vorm.
2.7
Eigenschappen van een duale fosforylatie-defosforylatiecyclus.
d[M ] = v4 − v1 = 0. We We bekijken eerst ´e´en evenwichtstoestand, namelijk het geval waar dt nemen aan dat zowel M AP KK als M KP 3 verzadigd zijn door hun substraten M , respectievelijk M pp. We willen onderzoeken wat er gebeurt met de concentraties van M en M p als we door verandering van parameters een stijging van [M pp] ondervinden. Door het feit dat Mtot een constante is en er geldt dat Mtot = [M ] + [M p] + [M pp], zien we al meteen in dat [M ] en [M p] niet beiden kunnen toenemen.
18
2.7. Eigenschappen van een duale fosforylatie-defosforylatiecyclus.
Aangezien we geen vrije M KP 3 meer hebben in ons systeem en zowel M pp als M p wil binden met dit enzym, is het duidelijk dat deze twee stoffen concurreren om de M KP 3. v4 ondervindt dus competitieve inhibitie van M pp (tegen M p). Dit laatste kan afgeleid worden uit de voorgaande intu¨ıtieve interpretatie samen met Figuur 2.2, maar kan ook wiskundig geverifieerd worden in (2.16): een stijging van [M pp] zorgt duidelijk voor een daling van v4 . Analoog kun je hieruit opmaken dat v4 ook tegengewerkt wordt door M . Een besluit dat we hieruit kunnen trekken is dat v1 ook competitieve inhibitie ondervindt van M p: ditmaal is de concurrentie tegen M gericht. Opnieuw is dit zowel uit Figuur 2.2 af te lezen als uit (2.13). Uit al deze resultaten kunnen we nu besluiten wat er gebeurt bij deze evenwichtstoestand als [M pp] stijgt. We willen dan namelijk v1 = v4 behouden, maar door een stijging van [M pp] hebben we een daling van v4 (door de competitieve inhibitie). We kunnen dit compenseren door ofwel opnieuw een stijging van v4 te bekomen, ofwel een daling van v1 (of een combinatie van beide). Willen we v4 doen stijgen, dan kunnen we de parameterwaarden zo aanpassen dat [M ] daalt (werkt v4 tegen). Om v1 te laten dalen, kunnen we het voorkomen van competitieve inhibitie uitbuiten en er dus voor zorgen dat [M p] stijgt. Samengevat hoort bij deze eerste evenwichtstoestand een daling van [M ] en/of een stijging van [M p] als [M pp] gaat stijgen. We willen nu ook rekening houden met de tweede evenwichtsvoorwaarde waarbij geldt dat d[M pp] =v2 −v3 = 0. Aangezien we al voldoen aan onze eerste evenwichtsvoorwaarde v4 = v1 , dt kunnen we reeds meer zeggen over wat er gebeurt met v2 bij een stijging van [M pp]. Uit de eerste evenwichtsvoorwaarde hadden we namelijk besloten dat [M ] & [M pp] %⇒ . [M p] % Aangezien M en M p onderlinge competitie voeren om het enzym M AP KK, houdt dit laatste in dat er een verhoging van v2 plaatsvindt (opnieuw zowel te zien in Figuur 2.2 als in (2.14)) en dat [M pp] toeneemt. We hebben dus te maken met een positieve feedback loop: een stijging van [M pp] zorgt ervoor dat de productie ervan nog verhoogd wordt. Om de evenwichtstoestand te behouden, moet v3 dus stijgen om gelijk te blijven aan de stijgende waarde van v2 . Door de symmetrie in de cyclus die weergegeven wordt in Figuur 2.2, kunnen we meteen besluiten dat we eveneens van een positieve feedback loop kunnen spreken als we [M ] bekijken. Ook hier impliceert een stijging van deze concentratie onder de evenwichtstoestanden dat er extra M aangemaakt wordt. Op deze manier voelen we reeds aan dat er mogelijks bistabiliteit voorkomt in dit model waarbij we ofwel te maken hebben met hoge [M ] en lage [M pp], ofwel omgekeerd.
19
2.8. Implementatie
2.8
Implementatie
We onderzoeken nu dit eerste model aan de hand van een implementatie in Matcont. Aangezien de drie veranderlijken [M ], [M p] en [M pp] samen gelijk zijn aan de constante Mtot , kunnen we uit het verloop van twee veranderlijken meteen het verloop van de derde afleiden. We kiezen er hier voor om het verloop van [M ] en [M pp] expliciet te bekijken. Ons model wordt dan in Matcont ingegeven als in Figuur 2.5:
Figuur 2.5: Implementatie van Model 1 (zie pagina 18) in Matcont. Voor de parameters kiezen we de waarden die na grondig onderzoek gekozen werden in [6]. Hierbij merken we op dat in het artikel met verschillende waarden voor [M AP KK]tot gewerkt wordt. In dit onderzoek willen we echter de continuatie van verschillende parameters bekijken en daarom is het van belang dat de andere parameterwaarden vast gekozen worden zodat een grondige studie opgezet kan worden. Dit geldt ook voor de parameter [M AP KK]tot . We hebben ervoor gekozen om deze parameter de vaste waarde 50 te geven gedurende de volledige implementatie. Ook merken we op dat er in het artikel geen parameterwaarde gegeven wordt voor Km5 . Aan de hand van tabel S1 in het supplementaire materiaal op [12] en de reeds h6 , kunnen we de afgeronde waarde Km5 = 78 afleiden. gedefinieerde uitdrukking Km5 = h−6 Dit werd achteraf bevestigd toen bleek dat de auteurs van [6] hun model elektronisch ter beschikking gesteld hebben in de SBML bibliotheek (Systems Biology Markup Language), [17] als de BIOMD0000000027.xml.origin file. Bij het inladen van deze file bleek daar Km5 inderdaad de waarde 78 te hebben. De gekozen parameterwaarden staan samengevat in Tabel 2.1.
20
2.8. Implementatie
Parameter Km1 Km2 Km3 Km4 Km5 k1cat k2cat k3cat k4cat Mtot [M AP KK]tot [M KP 3]tot
Waarde 50 500 22 18 78 0.01 15 0.084 0.06 500 50 100
Dimensie nM nM nM nM nM s−1 s−1 s−1 s−1 nM nM nM
Tabel 2.1: Parameterwaarden bij Model 1.
2.8.1
Tijdsintegratie
We kunnen nu de veranderingen van de concentraties in de tijd bekijken bij verschillende startwaarden voor [M ] en [M pp]. Hierbij kiezen we deze startwaarden willekeurig, maar houden we er wel rekening mee dat de som van deze concentraties nooit groter is dan 500. Met het oog op verder onderzoek hebben we ervoor gekozen om de tijdsintegratie te bekijken voor twee verschillende waarden van Km1 . We zien meteen een fundamenteel verschil in deze twee grafieken van [M pp] versus [M ] in Figuur 2.6.
(a) Km1 = 200
(b) Km1 = 50
Figuur 2.6: Tijdsintegratie van [M pp] versus [M ] met parameterwaarden zoals in Tabel 2.1, met uitzondering van Km1 (Model 1).
21
2.8. Implementatie
Als we tijdsintegratie uitvoeren over een voldoende grote tijdsspanne, dan merken we op dat we in het geval van Km1 = 200 (Figuur 2.6(a)) te maken hebben met ´e´en stabiel vast punt. Het systeem convergeert namelijk voor elke startwaarde naar ([M ], [M pp]) = (48.08893352; 445.5156715). In het geval waar Km1 = 50 (Figuur 2.6(b)) vinden we ´e´en onstabiel vast punt en twee stabiele vaste punten (bistabiliteit). In dit geval convergeert het systeem voor bepaalde startwaarden naar ([M ], [M pp]) = (11.53634286; 481.9477132) en voor andere naar ([M ], [M pp]) = (437.7299837; 49.4175348).
2.8.2
Vari¨ eren van de MM-constanten Km1 en Km2
We stellen dus vast dat het vari¨eren van een parameter het faseportret grondig kan veranderen. Daarom willen we nu het effect hiervan onderzoeken op het aantal equilibria en de stabiliteit ervan. Aangezien we op basis van Km1 verschillende gevallen hebben kunnen onderscheiden in Figuur 2.6, beginnen we met het effect van de parameter Km1 na te gaan. Hiervoor selecteren we de equilibria in Figuur 2.6(b) en nemen we de parameter Km1 als veranderlijke.
Figuur 2.7: E´en-parameter bifurcatiediagram bij het vari¨eren van Km1 met de andere parameterwaarden zoals in Tabel 2.1 (Model 1). Om de volledige grafiek te construeren, is het noodzakelijk dat we een bifurcatiediagram opstellen vanuit de twee verschillende equilibria. Als we namelijk de equilibriumbifurcaties tekenen in het (Km1 , [M pp])-vlak merken we op dat we te maken hebben met twee nietsamenvallende takken van equilibria. Het schijnbare snijpunt van deze verschillende takken is een transkritische bifurcatie. Hierbij kruisen twee takken van equilibria elkaar waar de ene tak een overgang inhoudt van stabiele naar onstabiele equilibria en de andere van onstabiele naar stabiele. We spreken hier van takken die elkaar schijnbaar snijden omdat we dit op deze figuur niet met zekerheid kunnen zeggen. We hebben namelijk een projectie van een vierdimensionaal systeem op een tweedimensionale ruimte. Om te confirmeren of de takken elkaar al dan niet 22
2.8. Implementatie
snijden, zetten we ook de ´e´en-parameter bifurcatiediagrammen met vari¨erende Km1 uit ten opzichte van [M ], respectievelijk [M p] in Figuur 2.8.
(a)
(b)
Figuur 2.8: E´en-parameter bifurcatiediagrammen bij het vari¨eren van Km1 met de andere parameterwaarden zoals in Tabel 2.1 (Model 1). Hier wordt aangetoond dat het snijpunt wel degelijk schijnbaar was, want in het (Km1 , [M p])vlak snijden de twee equilibriumtakken elkaar niet. De evenwichtscurve in Figuur 2.7 heeft een limietpunt in ([M ], [M pp], Km1 ) = (346.762093; 137.936921; 78.006341). In dit limietpunt hebben we twee zuiver re¨ele eigenwaarden, namelijk −0.331667 en 0. We konden op voorhand verwachten dat ´e´en van de twee eigenwaarden nul zou zijn aangezien dit een eigenschap is van limietpunten. Het negatief zijn van de andere eigenwaarde wil zeggen dat het limietpunt de grens vormt tussen stabiele en onstabiele equilibria. Dit wordt ook duidelijk als we de eigenwaarden onderzoeken van de andere punten op dit ´e´en-parameter bifurcatiediagram. We merken vooreerst op dat deze allemaal zuiver re¨eel zijn voor een positieve Km1 . Meer bepaald hebben we in de punten gelegen op de stippellijnen te maken met een negatieve en een positieve eigenwaarde en de punten op de volle lijnen hebben twee negatieve eigenwaarden. We weten dus dat de evenwichtspunten corresponderend met de stippellijnen onstabiel zijn en deze corresponderend met de volle lijnen stabiel. Als we dit model puur wiskundig verder ontleden, merken we op dat er zich een neutraal zadelpunt (aangeduid met H1) en twee Hopf-punten H2 en H3 in dit model bevinden (zie Figuur 2.7). Een neutraal zadelpunt wordt gekenmerkt door twee re¨ele eigenwaarden die op hun teken na gelijk zijn. H1 bevindt zich in ([M ], [M pp], Km1 ) = (−0.060936; 493.510391; −0.267592) en voldoet aan deze eigenschap aangezien ze de eigenwaarden ±1.07605 heeft. H2 (in ([M ], [M pp], Km1 ) = (−0.259535; 470.631445; −0.053128)) en H3 (in ([M ], [M pp], Km1 ) = (509.001208; −20.546097; −29.958857)) zijn Hopf-punten aangezien we in beide met twee zuiver imaginaire eigenwaarden te maken hebben met een tegengesteld teken (0.501583i en −0.501583i, respectievelijk 0.0531013i en −0.0531013i). In H2 is de Lyapunov co¨effici¨ent 1.081608 · 10−3 en dus positief, wat erop wijst dat H2 een subkritisch Hopf-punt is. Door het negatief zijn van de eerste Lyapunov co¨effici¨ent van H3 (= −1.701988 · 10−7 ) weten we dat 23
2.8. Implementatie
dit Hopf-punt superkritisch is. Aangezien we hier echter met een biochemisch model werken en concentraties steeds positief moeten zijn (wat hier niet het geval is bij [M ] en [M pp]) en door de definitie van Km1 ook geldt dat deze niet negatief kan zijn, hebben dit zadelpunt en deze Hopf-punten geen fysische betekenis. Toch zouden deze Hopf-punten voor ons van belang kunnen zijn, aangezien er periodieke banen kunnen ontstaan uit een Hopf-punt. Als deze periodieke banen zich zouden verderzetten tot bij positieve waarden voor Km1 , dan is dit van belang voor ons model. We kunnen echter reeds opmerken dat de periodieke banen zullen ontstaan voor dalende waarden van Km1 , aangezien we bij H3 weten dat dit Hopf-punt superkritisch is. Dit houdt in dat eventuele periodieke banen ontstaan aan de kant van de onstabiele equilibria en de banen zelf stabiel zijn. We onderzochten de stabiliteit van de equilibria en merkten reeds op dat we vanuit het Hopf-punt H3 bij stijgende waarden van de parameter Km1 te maken hadden met stabiele equilibria en bij dalende waarden met onstabiele, waardoor de eventuele periodieke banen zouden ontstaan aan de linkerkant van het Hopf-punt. Analoog ontstaan de periodieke banen aan de linkerkant van ons subkritisch Hopf-punt H2 aangezien dit gedeelte stabiel is. Deze vermoedens worden bevestigd met behulp van Matcont in Figuur 2.9.
(a)
(b)
Figuur 2.9: E´en-parameter bifurcatiediagrammen bij het vari¨eren van Km1 met de andere parameterwaarden zoals in Tabel 2.1 (Model 1), waarbij de periodieke banen weergegeven worden. Hierbij geeft 2.9(b) een gedetailleerder beeld weer van de periodieke banen die vanuit het Hopf-punt H3 ontstaan. Uit Figuur 2.7 kunnen we duidelijk afleiden dat we inderdaad een verschil in aantal equilibria moesten bekomen voor de tijdsintegratie als we verschillende waarden voor Km1 beschouwden. We kunnen meer bepaald zeggen dat we bij het constant houden van de andere parameterwaarden drie equilibria hebben als Km1 < 78 (waarvan twee stabiele) en ´e´en stabiel equilibrium als Km1 > 78. We kunnen analoog een bifurcatiediagram opstellen voor veranderende waarden van de parameter Km2 aan de hand van de equilibria in Figuur 2.6(b).
24
2.8. Implementatie
Figuur 2.10: E´en-parameter bifurcatiediagram bij het vari¨eren van Km2 met de andere parameterwaarden zoals in Tabel 2.1 (Model 1). In dit bifurcatiediagram (Figuur 2.10) vinden we nu twee limietpunten terug bij het vari¨eren van Km2 , namelijk in ([M ], [M pp], Km2 ) = (49.332649; 434.231791; 825.650688) en in ([M ], [M pp], Km2 ) = (354.249402; 129.630555; 370.509427) met twee zuiver re¨ele eigenwaarden −0.485888 en 0, respectievelijk −0.307781 en 0. Op dezelfde manier als in Figuur 2.7 kunnen we afleiden dat de equilibria gelegen op het gedeelte tussen de twee limietpunten onstabiel zijn en de overige equilibria stabiel. We hebben dus te maken met een bistabiel systeem en dit systeem vertoont hysteresis. Dit houdt in dat een parameter een bepaalde drempelwaarde moet overschrijden om het systeem in een andere stabiele toestand te krijgen en het systeem niet noodzakelijk terugkeert naar zijn oorspronkelijke stabiele toestand als de parameter de drempelwaarde terug in de andere richting overschrijdt. Dit kunnen we gemakkelijk aantonen op Figuur 2.10. Nemen we namelijk een lage waarde voor de parameter Km2 , dan zien we dat het systeem convergeert naar een stabiele toestand met hoge [M pp]. Laten we in deze stabiele toestand de parameter Km2 toenemen, dan neemt [M pp] geleidelijk af tot Km2 de drempelwaarde 826 aanneemt in het limietpunt. Bij het passeren van Km2 door deze drempelwaarde verspringt het systeem naar zijn andere stabiele toestand met lage [M pp]. Als we Km2 dan laten afnemen, blijven we nu echter in de stabiele toestand met lage [M pp] tot de andere drempelwaarde Km2 = 371 bereikt wordt in het andere limietpunt, waarna er terug versprongen wordt naar de andere stabiele toestand (met hoge [M pp]). Het onstabiele gebied tussen de twee limietpunten wordt dus nooit bereikt. We kunnen ook opmerken dat we bij hoge, respectievelijk lage waarden van Km2 meteen weten naar welke stabiele toestand het systeem zal convergeren. Bij tussenliggende waarden (tussen de twee drempelwaarden Km2 = 371 en Km2 = 826) kunnen we dit niet zomaar 25
2.8. Implementatie
bepalen. We selecteren nu ´e´en van de limietpunten uit Figuur 2.7 en 2.10, continueren deze met Km1 en Km2 variabel en stellen zo een bifurcatiediagram op in het (Km1 , Km2 )-vlak (zie Figuur 2.11).
Figuur 2.11: Continuatie van LP bij het vari¨eren van Km1 en Km2 met de andere parameterwaarden zoals in Tabel 2.1 (Model 1). Hier hebben we te maken met een Cusp-punt in ([M ], [M pp], Km1 , Km2 ) = (248.699454; 239.432353; 241.939473; 749.430954) met zuiver re¨ele eigenwaarden −0.535044 en 0. Opnieuw vinden we in dit wiskundig model een biologisch irrelevant speciaal punt, namelijk een Bogdanov-Takens punt dat enkel hypothetisch van tel is aangezien [M] en Km1 hier negatief zijn. Dit punt is terug te vinden op ([M ], [M pp], Km1 , Km2 ) = (−0.277373; 482.045704; −0.262908; 855.892385). In Figuur 2.11 zijn er voor parameterwaarden binnen het afgebakende gebied drie equilibria en erbuiten ´e´en equilibrium. Voor Km2 = 500 zien we dus inderdaad dat we 3 equilibria hebben bij Km1 < 78 en 1 equilibrium elders, zoals we uit Figuur 2.7 konden aflezen. Een tweede twee-parameter bifurcatiediagram dat we willen opstellen is deze met veranderende waarden voor Km1 en Mtot (zie Figuur 2.12). Opnieuw bakent de limietpuntcurve een gebied af waarvoor geldt dat er voor parameterwaarden binnen het gebied bistabiliteit is en voor waarden erbuiten slechts ´e´en stabiel equilibrium. Op deze grafiek vinden we een Cusp-punt terug in ([M ], [M pp], Km1 , Mtot ) = (148.521044; 151.171130; 136.584366; 307.614564) met eigenwaarden −0.78341 en 0. Verder hebben we ook te maken met een biologisch irrelevant Bogdanov-Takens punt in ([M ], [M pp], Km1 , Mtot ) = (−0.161634; 261.806726; −0.152999; 271.921490).
26
2.8. Implementatie
Figuur 2.12: Continuatie van LP bij het vari¨eren van Km1 en Mtot met de andere parameterwaarden zoals in Tabel 2.1 (Model 1).
2.8.3
Vari¨ eren van de parameters k1cat en k2cat
We kunnen verder nog experimenteren met het vari¨eren van andere parameters. Dit gebeurt analoog als in paragraaf 2.8.2, maar nu kiezen we ervoor om de invloed van de parameters k1cat en k2cat te onderzoeken door het opstellen van de bifurcatiediagrammen in Figuur 2.13 en 2.14.
Figuur 2.13: E´en-parameter bifurcatiediagram bij het vari¨eren van k1cat met de andere parameterwaarden zoals in Tabel 2.1 (Model 1).
27
2.8. Implementatie In Figuur 2.13 hebben we het ´e´en-parameter bifurcatiediagram opgesteld waarbij k1cat als vrije parameter beschouwd werd. We merken twee limietpunten op, namelijk in ([M ], [M pp], k1cat ) = (48.302882; 441.662799; 0.006181) en ([M ], [M pp], k1cat ) = (349.931191; 128.983760; 0.012899) met zuiver re¨ele eigenwaarden −0.784647 en 0, respectievelijk −0.244293 en 0. De equilibria gelegen tussen de twee limietpunten zijn onstabiele equilibria (ze hebben zowel een positieve als een negatieve eigenwaarde) en de overige equilibria zijn stabiel (ze hebben twee negatieve eigenwaarden).
Figuur 2.14: Continuatie van LP bij het vari¨eren van k1cat en k2cat met de andere parameterwaarden zoals in Tabel 2.1 (Model 1). In Figuur 2.14 hebben we de continuatie van zowel k1cat als k2cat uitgevoerd bij Model 1. Er bevindt zich een Cusp-punt op ([M ], [M pp], k1cat , k2cat ) = (137.675378; 199.319355; 0.077105; 0.912146) met eigenwaarden −0.0586781 en 0. Bij parameterwaarden die zich binnen het gebied begrensd door de LP-curves bevinden, is er bistabiliteit in het systeem. Bij parameterwaarden die zich buiten dit gebied bevinden, hebben we te maken met ´e´en stabiel equilibrium.
2.8.4
Besluit voor bistabiliteit
Uit Figuur 2.2 konden we intu¨ıtief reeds verwachten dat we met bistabiliteit zouden te maken hebben als we aannemen dat er substraatverzadiging van de enzymen optreedt. Door deze substraatverzadiging en het gegeven dat de twee fosforylatiestappen (respectievelijk defosforylatiestappen) door dezelfde kinase (respectievelijk fosfatase) uitgevoerd worden, treedt er competitieve inhibitie op tussen M en M p (om M AP KK) en tussen M pp en M p (om M KP 3) en dit zorgt voor de stabiele toestanden van hoge [M ] en lage [M pp] enerzijds en lage [M ] en hoge [M pp] anderzijds (zie ook paragraaf 2.7). We bekijken nu de opgestelde figuren in paragraaf 2.8.2 en 2.8.3 en proberen die te kaderen binnen deze intu¨ıtieve redenering. 28
2.8. Implementatie
Als eerste zien we in Figuur 2.11 dat we voor deze specifieke parameterwaarden hebben dat Km2 > Km1 nodig is voor bistabiliteit. Uit de definitie van Km1 en Km2 kunnen we opmaken dat dit de relatieve afbraak is van het complex M − M AP KK (respectievelijk M p − M AP KK). M p − M AP KK ondervindt dus een grotere relatieve afbraak dan M − M AP KK of er wordt met andere woorden relatief meer M − M AP KK gevormd dan M p − M AP KK. Dit betekent dat zowel M als M p willen binden met het enzym M AP KK en dat de ongefosforyleerde vorm M meer enzymen aantrekt om mee te binden dan M p. Deze figuur weerspiegelt bijgevolg de competitieve inhibitie van de tweede reactie door het substraat van de eerste reactie bij de fosforylatie. In Figuur 2.12 zien we de relatie Mtot > Km1 voor het verkrijgen van bistabiliteit waarbij deze beide grootheden uitgedrukt worden in nM . Er moet dus meer vrije M , M p en M pp zijn dan dat er M − M AP KK afgebroken wordt. Dit impliceert dat er niet genoeg enzymen aanwezig zijn om alle vrije substraten te binden of dat we met andere woorden te maken hebben met substraatverzadiging van de enzymen. Tot slot zien we in Figuur 2.14 voor het bekomen van bistabiliteit dat k2cat > k1cat bij deze specifieke parameterwaarden.
29
Hoofdstuk 3 Signaalschakelaars en bistabiliteit door dubbele multisite fosforylatie Dit hoofdstuk is gebaseerd op [6].
3.1
Beschrijving van een duale fosforylatie-defosforylatiecyclus met multisite-bindingen
Zoals reeds aangegeven, hebben we tot nu toe onderzoek gedaan naar een model dat vereenvoudigd was en waar geen rekening gehouden werd met de bindingsplaats van de fosfaatgroep bij de fosforylatie van M AP K. Om nu meer aan te sluiten bij de werkelijkheid breiden we ons onderzoek uit naar een model waar wel een onderscheid gemaakt wordt tussen de verschillende bindingsplaatsen (het theronine- (T) en tyrosine- (Y) residu). Opnieuw wordt dit systeem verondersteld een distributief mechanisme te zijn en moet er dus per (de)fosforylatie een nieuwe binding met een enzym plaatsvinden.
Figuur 3.1: Cyclus van MAPK met multisite-bindingen.
30
3.2. Reactievergelijkingen
In Figuur 3.1 wordt de cyclus van M AP K weergegeven als we rekening houden met deze verschillende bindingsplaatsen. Er wordt in deze cyclus aangenomen dat de kinase M EK (M AP K/ERK kinase) zowel in staat is om M tot M pY te fosforyleren als M tot M pT , evenals M pY tot M pp en M pT tot M pp. De defosforylatie in de omgekeerde richting wordt uitgevoerd door M KP 3 en zorgt voor de omzetting van M pp tot M pT , van M pT tot M en van M pY tot M . In deze cyclus is de fosfatase M KP 3 niet in staat om M pp te defosforyleren tot M pY . De fosforylaties worden uitgevoerd met snelheden v1 , v2 , v3 en v4 en de defosforylaties met snelheden v5 ,v6 en v7 zoals af te lezen is in Figuur 3.1. We nemen opnieuw aan dat de som van de concentraties van alle vormen van M AP K constant is, wat in dit model leidt tot de uitdrukking Mtot = [M ]+[M pY ]+[M pT ]+[M pp] = constant. Samen met Figuur 3.1 geeft dit dan het volgend stelsel van vergelijkingen: d[M ] = −v1 − v3 + v6 + v7 , dt d[M pY ] = v1 − v2 − v7 , (3.1) dt d[M pp] = v2 + v4 − v5 , dt [M pT ] = Mtot − [M ] − [M pY ] − [M pp], waarbij de afgeleiden naar de tijd beschouwd worden.
3.2
Reactievergelijkingen
We kunnen nu opnieuw reactievergelijkingen opstellen bij dit model waarbij de aanname van een constante hoeveelheid ATP/ADP geldt en de (on)omkeerbaarheid van bepaalde stappen. Eveneens maken we de vereenvoudiging bij de fosforylatie dat de tussenvormen niet expliciet gevormd worden en het uiteenvallen van het complex bestaande uit het enzym en de gefosforyleerde vorm onmiddellijk gebeurt. We bekomen dan de volgende reactievergelijkingen: k1 k2 Y∗
→ M + M EK M − M EK M pY + M EK, k−1 k3 k4 ∗
→ M pY + M EK M pY − M EK M pp + M EK, k−3 (3.2) k k6 5 M + M EK M − M EK T ∗ → M pT + M EK, k−5 k7 k8 ∗ M pT + M EK M pT − M EK → M pp + M EK. k−7
31
3.3. Uitdrukkingen voor de snelheid h1 h2 h3 ∗ M pp + M KP 3 M pp − M KP 3 → M pT − M KP 3 M pT + M KP 3, h−1 h−3 h4 h5 h6 ∗ T M pT + M KP 3 M pT − M KP 3 → M − M KP 3 M + M KP 3, h−4 h−6 h7 h8 h9 ∗ Y M pY + M KP 3 M pY − M KP 3 → M − M KP 3 M + M KP 3. h−7 h−9
(3.3)
Hierbij voorzien we de actieve complexen (waarbij de (de)fosforylatie nog niet tot stand gekomen is en het enzym zijn werking nog moet uitvoeren) opnieuw van een * in de notatie. In (3.2) en (3.3) is er echter nog een nieuwe notatie weergegeven. Bij complexen gevormd door het substraat M AP K en een enzym (zowel de kinase als de fosfatase) wordt er namelijk een letter weergegeven als bovenindex. Deze letter wijst op de plaats van de binding van het enzym en komt zowel voor bij actieve als niet-actieve complexen. Zo zien we bijvoorbeeld in de eerste regel van (3.2) het intermediair complex M − M EK Y ∗ . Dit is dan het intermediair waarbij de kinase gebonden is aan het substaat op het Y-residu en de kinase dus zal instaan voor de fosforylatie van M naar M pY . De ki ’s en hi ’s zijn constanten en meer bepaald de snelheidsconstanten van de reacties.
3.3
Uitdrukkingen voor de snelheid
Op basis van (3.2) en (3.3) kunnen we nu analoog als in paragraaf 2.5 MM-vergelijkingen opstellen voor de verschillende snelheden die weergegeven zijn in Figuur 3.1. Hierbij veronderstellen we dat [M EK]tot en [M KP 3]tot constant zijn, waarbij [M EK]tot =[M EK] + [M − M EK Y ∗ ] + [M − M EK T ∗ ] + [M pY − M EK ∗ ] + [M pT − M EK ∗ ], [M KP 3]tot =[M KP 3] + [M pp − M KP 3∗ ] + [M pT − M KP 3] + [M pT − M KP 3∗ ] + [M − M KP 3T ] + [M pY − M KP 3∗ ] + [M − M KP 3Y ]. We maken eveneens de aanname dat er evenwicht is tussen de opbouw en de afbraak van de intermediairen en bekomen zo de volgende uitdrukkingen voor de snelheid: k1cat [M EK]tot [M ]/Km1 , 1 + [M ](Km1 + Km3 )/(Km1 · Km3 ) + [M pY ]/Km2 + [M pT ]/Km4 k2cat [M EK]tot [M pY ]/Km2 v2 = , 1 + [M ](Km1 + Km3 )/(Km1 · Km3 ) + [M pY ]/Km2 + [M pT ]/Km4 k3cat [M EK]tot [M ]/Km3 v3 = , 1 + [M ](Km1 + Km3 )/(Km1 · Km3 ) + [M pY ]/Km2 + [M pT ]/Km4 k4cat [M EK]tot [M pT ]/Km4 v4 = , 1 + [M ](Km1 + Km3 )/(Km1 · Km3 ) + [M pY ]/Km2 + [M pT ]/Km4 v1 =
32
(3.4) (3.5) (3.6) (3.7)
3.3. Uitdrukkingen voor de snelheid k5cat [M KP 3]tot [M pp]/Km5 , 1 + [M pp]/Km5 + [M pT ]/Km6 + [M pY ]/Km7 + [M ]/Km8 k6cat [M KP 3]tot [M pT ]/Km6 v6 = , 1 + [M pp]/Km5 + [M pT ]/Km6 + [M pY ]/Km7 + [M ]/Km8 k7cat [M KP 3]tot [M pY ]/Km7 , v7 = 1 + [M pp]/Km5 + [M pT ]/Km6 + [M pY ]/Km7 + [M ]/Km8 v5 =
waarbij we k1cat k2cat k3cat k4cat
= k2 , = k4 , = k6 , = k8 ,
k5cat =
h2 , 1 + h2 /h3
k6cat =
h5 , 1 + h5 /h6 + h−3 (h−4 + h5 )/(h3 h4 ) h8 , 1 + h8 /h9 k−1 + k2 , k1 k−3 + k4 , k3 k−5 + k6 , k5 k−7 + k8 , k7 h−1 + h2 , h1 + h1 h2 /h3 k−4 + k5 , h4 1 + h5 /h6 + h−3 (h−4 + h5 )/(h3 h4 )
k7cat = Km1 = Km2 = Km3 = Km4 = Km5 = Km6 =
h−7 + h8 , h7 + h7 h8 /h9 h6 = , h−6 + h6 h−9 /h9
Km7 = Km8
gesteld hebben. Deze resultaten leiden tot Model 2:
33
(3.8) (3.9) (3.10)
3.4. Implementatie
Model 2. d[M ] = −v1 − v3 + v6 + v7 , dt d[M pY ] = v1 − v2 − v7 , dt d[M pp] = v2 + v4 − v5 , dt [M pT ] = Mtot − [M ] − [M pY ] − [M pp], k1cat [M EK]tot [M ]/Km1 v1 = , 1 + [M ](Km1 + Km3 )/(Km1 · Km3 ) + [M pY ]/Km2 + [M pT ]/Km4 k2cat [M EK]tot [M pY ]/Km2 v2 = , 1 + [M ](Km1 + Km3 )/(Km1 · Km3 ) + [M pY ]/Km2 + [M pT ]/Km4 k3cat [M EK]tot [M ]/Km3 , v3 = 1 + [M ](Km1 + Km3 )/(Km1 · Km3 ) + [M pY ]/Km2 + [M pT ]/Km4 k4cat [M EK]tot [M pT ]/Km4 v4 = , 1 + [M ](Km1 + Km3 )/(Km1 · Km3 ) + [M pY ]/Km2 + [M pT ]/Km4 k5cat [M KP 3]tot [M pp]/Km5 v5 = , 1 + [M pp]/Km5 + [M pT ]/Km6 + [M pY ]/Km7 + [M ]/Km8 k6cat [M KP 3]tot [M pT ]/Km6 , v6 = 1 + [M pp]/Km5 + [M pT ]/Km6 + [M pY ]/Km7 + [M ]/Km8 k7cat [M KP 3]tot [M pY ]/Km7 v7 = , 1 + [M pp]/Km5 + [M pT ]/Km6 + [M pY ]/Km7 + [M ]/Km8 waarbij we het verloop van [M pT ] niet expliciet beschouwen, maar dit afleiden uit het verloop van de andere concentraties en de relatie Mtot = [M ] + [M pT ] + [M pY ] + [M pp].
3.4
Implementatie
We onderzoeken dit model met behulp van een implementatie in Matcont en dit wordt ge¨ımplementeerd als in Figuur 3.2. Bij de implementatie van dit model moeten uiteraard opnieuw parameterwaarden gekozen worden. Deze worden weergegeven in Tabel 3.1. De parameters in de laatste kolom worden analoog gekozen als aangegeven in Figuur 5 van het artikel [6] en de eerste twee kolommen zijn gebaseerd op Tabel 3 in het extra materiaal bij dit artikel dat we vinden in [12].
34
3.4. Implementatie
Figuur 3.2: Implementatie van Model 2 (pagina 34) in Matcont. Parameter Km1 Km2 Km3 Km4 Km5 Km6 Km7 Km8
Waarde (nM ) 410 40 20 300 22 18 34 40
Parameter k1cat k2cat k3cat k4cat k5cat k6cat k7cat
Waarde (s−1 ) 1.08 0.007 0.008 0.45 0.084 0.06 0.108
Parameter Mtot [M EK]tot [M KP 3]tot
Waarde (nM ) 1000 400 180
Tabel 3.1: Parameterwaarden bij Model 2.
3.4.1
Tijdsintegratie en continuatie
Na het uitvoeren van tijdsintegratie met deze parameterwaarden en willekeurige startwaarden voor de veranderlijken (waarbij de som maximum 1000 is) kunnen we zien dat onze veranderlijken naar een equilibrium convergeren in ([M ], [M pY ], [M pp]) = (0.46681315; 5.7508542; 964.47253) (figuur niet weergegeven). We selecteren dit equilibrium en voeren continuatie uit met ofwel de parameter [M EK]tot vrij, ofwel de parameter [M KP 3]tot vrij, waarbij de overige parameterwaarden constant gehouden worden. Deze grafieken worden weergegeven in Figuur 3.3.
35
3.4. Implementatie
Figuur 3.3: E´en-parameter bifurcatiediagram bij het vari¨eren van [M EK]tot , respectievelijk [M KP 3]tot met de andere parameterwaarden zoals in Tabel 3.1 (Model 2). Beide grafieken kunnen op een analoge manier besproken worden. We merken namelijk op dat we twee limietpunten hebben, gelegen op ([M ], [M pY ], [M pp]) = (4.678733; 42.377537; 876.531333) en ([M ], [M pY ], [M pp]) = (106.657539; 348.221670; 388.863128). In die punten heeft [M EK]tot de waarde 272.003226, respectievelijk 355.541146 en [M KP 3]tot 264.702743, respectievelijk 202.508207. Als we de eigenwaarden in de verschillende punten van de equilibriumkromme bekijken, dan stellen we vast dat deze allemaal zuiver re¨eel zijn. In de punten die gelegen zijn tussen de twee limietpunten hebben we twee negatieve eigenwaarden en ´e´en positieve en elders hebben we drie negatieve eigenwaarden. Bij lage en hoge parameterwaarden voor zowel [M EK]tot als [M KP 3]tot hebben we dus te maken met ´e´en stabiel equilibirum. Bij waarden die daartussen liggen (tussen de twee limietpunten) kunnen we spreken van bistabiliteit waarbij er een onstabiel equilibrium ligt tussen twee stabiele equilibria. Bij de uitgevoerde tijdsintegratie lagen beide gekozen parameterwaarden buiten het gebied van bistabiliteit, dus hadden we bij andere startwaarden voor de verandelijken geen ander equilibrium kunnen bekomen aangezien er in dit geval maar ´e´en equilibrium in het systeem aanwezig is. Vertrekken we vanuit de limietpunten die vastgesteld werden in Figuur 3.3, dan kunnen we continuatie uitvoeren met de andere parameter vrij en dit wordt gedaan in Figuur 3.4. We merken hierbij op dat het niet vanzelfsprekend is om deze figuur te construeren. Matcont geeft de foutmelding dat de stapgrootte te klein is en de grafiek wordt niet getekend. De poging om de stapgrootte te vergroten, levert geen betere resultaten op. Dit probleem kan echter ook veroorzaakt worden doordat de limietpunten niet nauwkeurig genoeg berekend werden en we keren dus een stap terug om de limietpunten opnieuw te berekenen met behulp van een kleinere tolerantie. Dit kunnen we aanpassen in het kader met de naam ’Continuer’ in Matcont en dit zorgt ervoor dat de numerieke methode die gebruikt wordt om de equilibriumkromme te construeren kleinere integratiestappen gebruikt. De vergrote precisie in deze stap zorgt ervoor dat het twee-parameter bifurcatiediagram nu wel opgesteld kan worden.
36
3.4. Implementatie
Figuur 3.4: Continuatie van LP bij het vari¨eren van [M EK]tot en[M KP 3]tot met de andere parameterwaarden zoals in Tabel 3.1 (Model 2). Alhoewel negatieve waarden van onze twee veranderende parameters fysisch niet toegelaten zijn, worden ze toch weergegeven in Figuur 3.4. De bedoeling hiervan is dat het zeker duidelijk wordt dat we te maken hebben met twee verschillende limietpuntkrommen. Op beide bevindt zich een Zero Hopf-punt dat in feite een neutraal zadelpunt voorstelt dat tevens een limietpunt is. De twee neutrale zadelpunten hebben parameterwaarden ([M EK]tot , [M KP 3]tot ) = (0; 0), waardoor ze lijken samen te vallen op deze figuur. Als we echter de waarde van de veranderlijken nagaan in de neutrale zadelpunten, dan vinden we deze enerzijds in ([M ], [M pY ], [M pp]) = (106.657538; 348.221669; 388.863130) en anderzijds in ([M ], [M pY ], [M pp]) = (4.678732; 42.377535; 876.531338). Dit geeft meteen aan dat de twee limietpuntkrommen in Figuur 3.4 elkaar kruisen en niet snijden. Beide neutrale zadelpunten hebben drie eigenwaarden die 0 zijn. In het gebied dat afgebakend is door de twee limietpuntkrommen kan er zich bistabiliteit voordoen, terwijl de gebieden erbuiten enkel ´e´en equilibrium kunnen bevatten.
37
Hoofdstuk 4 Analyse van bistabiliteit door dubbele fosforylatie Dit hoofdstuk is gebaseerd op [8].
4.1
Twee-staps cycli: beschrijving van het model
In dit hoofdstuk werken we met het niet-specifieke eiwit W in een distributief mechanisme met covalente modificatie.
Figuur 4.1: Kinetisch diagram waarbij het eiwit W drie verschillende vormen Wα , Wβ en Wγ heeft, die gevormd worden door de inwerking van de kinase e1 en de fosfatase e2 . Het eiwit W kan drie vormen aannemen: Wα , Wβ en Wγ . Hierbij stelt Wα de niet-gefosforyleerde vorm van het eiwit voor, Wβ de gefosforyleerde vorm van Wα en Wγ de gefosforyleerde vorm van Wβ (of dus de dubbel gefosforyleerde vorm van Wα ). We nemen aan dat de twee fosforylatiestappen door hetzelfde enzym e1 uitgevoerd worden (de kinase) en de twee defosforylatiestappen door het enzym e2 (de fosfatase). We nemen eveneens aan dat elke (de)fosforylatiestap volgend distributief MM-mechanisme
38
4.1. Twee-staps cycli: beschrijving van het model
volgt: kai ki WS + e eWS → WP + e, kdi
(4.1)
waarbij kai , kdi en ki respectievelijk de samenstellings- (association), splitsings- (dissociation) en katalytische constante zijn van stap i (i ∈ {1, 2, 3, 4}). eWS is het MM-complex dat gevormd werd door het enzym e en zijn substraat WS om het product WP te produceren. Hierbij nemen we aan dat de nodige energie en fosfaatgroepen in constante hoeveelheid aanwezig zijn en hiermee reeds rekening gehouden wordt in de kinetische constanten. We merken op dat we hier (in tegenstelling tot in paragraaf 2.4) enkel rekening houden met het vormen van het substraat-enzym-complex en niet met het vormen van het product-enzymcomplex, wat de analyse zal vereenvoudigen. Nemen we nu vi als de snelheid waarmee stap i uitgevoerd wordt, dan zien we aan de hand van Figuur 4.1 dat de differentiaalvergelijkingen die de tijdsevolutie van dit systeem weergeven als volgt zijn: dWα = v2 − v1 , dt dWβ = v1 − v2 − v3 + v4 , dt dWγ = v3 − v4 . dt
(4.2)
We weten dat de enzymen zich zowel in gebonden als ongebonden vorm manifesteren en we hebben meer bepaald dat e1T = [e1 ] + [e1 Wα ] + [e1 Wβ ], e2T = [e2 ] + [e2 Wγ ] + [e2 Wβ ], waarbij eiT de totale concentratie van het enzym ei voorstelt. We nemen aan dat de totale concentratie van de enzymen e1 en e2 gedurende de cyclus constant is. Als we nu ook aannemen dat er evenwicht is tussen de opbouw en de afbraak van de intermediairen in (4.1), dan bekomen we de volgende MM-vergelijkingen (analoge uitwerking als in paragraaf 2.5): v1 =
v2 =
α] Vm1 [W Km1
[W ]
,
v3 =
[Wβ ] [Wα ] + Km3 Km1 [Wβ ] Vm2 Km2 , [Wβ ] [Wγ ] +K + Km2 m4
1+
1
v4 =
Hierbij hebben we de volgende notaties ingevoerd: kdi + ki Kmi = , kai 39
β Vm3 Km3
1
,
[Wβ ] [Wα ] + Km3 Km1 [Wγ ] Vm4 K m4 . [Wβ ] [Wγ ] +K + Km2 m4
1+
i = 1, 2, 3, 4
4.1. Twee-staps cycli: beschrijving van het model
en Vm1 = k1 e1T , Vm2 = k2 e2T ,
Vm3 = k3 e1T , Vm4 = k4 e2T .
Om de wiskundige analyse te vereenvoudigen, kunnen we de concentraties van de veranderlijken Wα , Wβ en Wγ dimensieloos maken. Hierbij beschouwen we de totale concentratie van het eiwit. We nemen aan dat deze constant is en door het combineren van de gegevens in Figuur 4.1 en (4.1) is deze van de volgende vorm: WT = [Wα ] + [Wβ ] + [Wγ ] + [e1 Wα ] + [e1 Wβ ] + [e2 Wγ ] + [e2 Wβ ]. De dimensieloze concentraties van de veranderlijken worden dan weergegeven als α=
[Wα ] , WT
β=
[Wβ ] , WT
γ=
[Wγ ] , WT
en deze liggen bijgevolg allemaal tussen 0 en 1. Samen met (4.2) bekomen we de volgende differentiaalvergelijkingen voor de dimensieloze veranderlijken voor de tijdsevolutie van dit systeem: dα v2 − v1 = , dt WT dβ v1 − v2 − v3 + v4 = , dt WT v3 − v4 dγ = . dt WT Voeren we nu analoog de volgende herschaling in: KS1 =
Km1 , WT
KS2 =
Km2 , WT
KS3 =
Km3 , WT
KS4 =
Km4 , WT
dan kunnen we onze MM-vergelijkingen herleiden tot v1 = v2 =
Vm1 KαS1 1+
α KS1
+
β KS3
Vm2 KβS2 1+
γ KS4
+
β KS2
,
v3 =
,
v4 =
Vm3 KβS3 1+
α KS1
+
β KS3
Vm4 KγS4 1+
γ KS4
+
β KS2
, .
(4.3)
Voor het verdere verloop van deze analyse defini¨eren we nog enkele dimensieloze grootheden: Vm3 k3 = , Vm1 k1 Vm2 k2 = = , Vm4 k4 Vm1 k1 e1T = = . Vm4 k4 e2T
r31 = r24 χ14
40
4.1. Twee-staps cycli: beschrijving van het model
Hierbij kunnen we duidelijk zien dat r31 en r24 onafhankelijk zijn van de enzym-concentraties. Tot slot wordt een nieuwe parameter Θ ingevoerd die we de asymmetrische factor noemen en die gedefinieerd wordt als Θ := r31 r24 =
Vm3 Vm2 k3 k2 = . Vm1 Vm4 k1 k4
De asymmetrische factor is met andere woorden een quoti¨ent met in de teller het product van de katalytische constanten van de stappen die Wβ consumeren (k2 en k3 ) en in de noemer het product van de katalytische constanten van de stappen die Wβ produceren (k1 en k4 ).
4.1.1
Implementatie
Nemen we de restrictie KS1 = KS2 = KS3 = KS4 = KS dan ziet het gebruikte model er als volgt uit: Model 3. v2 − v1 dα = , dt WT v1 − v2 − v3 + v4 dβ = , dt WT γ = 1 − α − β, Vm1 KαS , v1 = 1 + KαS + KβS v3 = v2 = v4 =
Vm3 KβS 1+
α KS
+
β KS
Vm2 KβS 1+
γ KS
+
β KS
Vm4 KγS 1+
γ KS
+
β KS
, , ,
Herschrijven we het model hier zodanig dat we de waarden voor de parameters die voor ons van belang lijken, kunnen aanpassen, dan wordt dit model in Matcont vertaald als in Figuur 4.2. We kiezen ervoor om het verloop van de veranderlijken α en β expliciet weer te geven en het verloop van de veranderlijke γ enkel te defini¨eren door de relatie α + β + γ = 1. Om grafieken te kunnen tekenen waarbij γ geplot wordt, wordt ook het model met α en γ als expliciete veranderlijken ingevoerd in Matcont.
41
4.2. Bistabiliteit in twee-staps cycli: wiskundig
Figuur 4.2: Implementatie van Model 3 (pagina 41) in Matcont. Voor het bepalen van de parameterwaarden baseren we ons op het Ortega-model in [13]. We zullen echter verschillende waarden voor Θ en r31 beschouwen, aangezien uit de wiskundige analyse het belang van Θ bleek bij het al dan niet voorkomen van bistabiliteit. De parameterwaarde voor de totale concentratie van het eiwit W werd gekozen naar analogie met [6] of Hoofdstuk 2. Meer bepaald zien onze vaste parameterwaarden eruit als in Tabel 4.1. Parameter KS χ14 Vm1 WT
Waarde 0.01 1.10 1 500
Dimensie 1 1 nM s−1 nM
Tabel 4.1: Parameterwaarden bij Model 3.
4.2
Bistabiliteit in twee-staps cycli: wiskundig
We willen opnieuw onderzoeken hoe dit systeem evolueert na een bepaalde tijdsspanne. Aangezien het model goed lijkt op Model 1 dat ge¨ımplementeerd werd in paragraaf 2.8, kunnen we reeds veronderstellen dat we bepaalde intervallen van parameterwaarden hebben waar bistabiliteit mogelijk is. In andere gebieden zullen we echter slechts ´e´en stabiel equilibrium vinden. Hoewel het numerieke gedeelte van het onderzoek naar dit model pas later gebeurt (in paragraaf 4.3), willen we dit al aantonen met de hulp van Matcont. Daartoe voeren we tijdsintegratie uit met het model ge¨ımplementeerd in Matcont zoals in Figuur 4.2 en de 42
4.2. Bistabiliteit in twee-staps cycli: wiskundig
parameterwaarden als in Tabel 4.1. We voeren vanuit een bekomen equilibrium continuatie uit met χ14 als vrije parameter en verschillende waarden van Θ (zie Figuur 4.3).
Figuur 4.3: E´en-parameter bifurcatiediagram bij het vari¨eren van χ14 voor verschillende waarden van Θ (Model 3) met vaste parameterwaarden r31 = 2, KS = 0.01, Vm1 = 1 en WT = 500. Het is meteen duidelijk dat de nieuw ingevoerde parameter Θ een invloed heeft op het aantal equilibria in het systeem. Aan de hand van deze figuur kunnen we aantonen dat we ergens een waarde voor Θ zouden moeten vinden die de overgang van ´e´en equilibrium naar drie equilibria (bistabiliteit) inhoudt en dus het bifurcatiepunt bepaalt. In Figuur 4.3 zien we namelijk dat we bij een kleine Θ (bv. Θ = 1) te maken hebben met ´e´en stabiel equilibrium voor alle waarden van χ14 . Voor een grotere Θ (bv. Θ = 4) is er een gebied van waarden van χ14 waarbij twee stabiele equilibria mogelijk zijn met ´e´en onstabiel equilibrium tussen de twee stabiele. Tussen Θ = 1 en Θ = 4 ligt er dus een waarde voor Θ waar we ons bifurcatiepunt terugvinden en dit willen we in deze paragraaf wiskundig onderzoeken.
4.2.1
Bepalen van de evenwichtstoestand
In een eerste stap zoeken we naar de evenwichtspunten. In een evenwichtstoestand geldt dat de rechterleden van (4.2) gelijkgesteld mogen worden aan 0 aangezien we hier te maken hebben met behoud van de verschillende vormen van het substraat. We hebben dan v1 = v2 en v3 = v4 . Kijken we nu naar de uitdrukkingen voor deze snelheden in (4.3), dan merken we op dat de noemers van v1 en v3 gelijk zijn, evenals de noemers van v2 en v4 . Als we dus v2 v1 en beschouwen, worden de termen uit de noemers weggedeeld. Uit de de quoti¨enten v3 v4 gelijkheden v1 = v2 en v3 = v4 , volgdend uit het behoud van de stoffen, kunnen we de relatie
43
4.2. Bistabiliteit in twee-staps cycli: wiskundig v1 v2 = halen. Als we de vergelijkingen uit (4.3) invullen in deze gelijkheid bekomen we v3 v4 Vm1 KαS1 Vm2 KβS2 v1 v2 = ⇔ = v3 v4 Vm4 KγS4 Vm3 KβS3 α KS3 γ KS2 Vm2 Vm3 ⇔ = KS1 β KS4 β V V | m1{z m4} Θ
KS2 KS3 αγ ⇔ = Θ. KS1 KS4 β 2
(4.4)
Om het model verder te vereenvoudigen, nemen we aan dat de concentratie van de verschillende vormen van het eiwit die voorkomen veel groter is dan de totale concentratie van de twee katalysatoren, m.a.w. WT e1T en WT e2T . Dit betekent dus dat alle enzymen verzadigd zijn door het eiwit en dat [WS ] [eWS ], waarbij de concentraties van deze complexen constant verondersteld werden in de MM-vergelijkingen. De laatste vier termen van WT = [Wα ] + [Wβ ] + [Wγ ] + [e1 Wα ] + [e1 Wβ ] + [e2 Wγ ] + [e2 Wβ ] zijn dan verwaarloosbaar t.o.v. de eerste drie termen. Laten we deze termen dus wegvallen en maken we de vergelijking dimensieloos door beide leden te delen door WT , dan krijgen we α + β + γ = 1.
(4.5)
We maken opnieuw een vereenvoudiging van dit model door aan te nemen dat de MMconstanten van de kinase, respectievelijk van de fosfatase gelijk zijn, d.i. KS1 = KS3 en KS2 = KS4 . Rekening houden met deze veronderstelling zorgt ervoor dat (4.4) equivalent is met αγ = Θ. β2
(4.6)
Uit (4.5) en (4.6) kunnen we nu een uitdrukking afleiden voor α en β in functie van γ en Θ. Zonderen we namelijk α af uit deze twee uitdrukkingen en stellen we deze gelijk aan elkaar, dan bekomen we β 2Θ γ 2 2 ⇔β Θ + βγ + γ − γ = 0 p −γ ± γ 2 − 4Θγ 2 + 4Θγ ⇒β = . 2Θ 1−β−γ =
We weten dat de concentraties α, β en γ door hun definities in het interval [0, 1] liggen en uit de definitie van Θ volgt ook dat deze positief is. Daaruit volgt dat γ 2 ≤ γ, of dus −4Θγ 2 + 4Θγ ≥ 0. De wortel van de berekende discriminant is bijgevolg minstens γ. Voor een positieve p β moet de berekende teller door het positief zijn van Θ dus ook positief zijn, wat door γ 2 − 4Θγ 2 + 4Θγ ≥ γ enkel kan als het plusteken gekozen wordt in de teller 44
4.2. Bistabiliteit in twee-staps cycli: wiskundig
van β. We hebben p γ 2 − 4Θγ 2 + 4Θγ β= 2Θ ⇒α=1−β−γ p 2Θ + γ − γ 2 − 4Θγ 2 + 4Θγ − 2Θγ 2γ = · 2Θ 2γ 2 p −γ + γ 2 − 4Θγ 2 + 4Θγ = . 4Θγ −γ +
Ook voor χ14 kunnen we een uitdrukking berekenen, hetzij in functie van γ, Θ, KS1 en KS2 . Vm1 en de voorwaarde voor de Hiervoor gebruiken we de definitie van χ14 , namelijk χ14 = Vm4 evenwichtstoestand v3 = v4 . Combineren van de aanname dat KS1 = KS3 en KS2 = KS4 met (4.3), geeft ons Vm3 β Vm4 γ = KS1 + α + β KS2 + γ + β Vm3 β(KS2 + γ + β) ⇒Vm4 = , γ(KS1 + α + β) en dus Vm1 Vm4 Vm1 γ(KS1 + α + β) = Vm3 β(KS2 + γ + β) αγ + βγ + γKS1 = r31 β(KS2 + β + γ) β 2 Θ + βγ + γKS1 = , r31 β(KS2 + β + γ)
χ14 =
waarbij de voorlaatste stap volgt uit de definitie van r31 en de laatste uit αγ = β 2 Θ (uit (4.6)). Samengevat bekomen we dus de volgende drie uitdrukkingen bij een evenwichtstoestand: 2 p −γ + γ 2 − 4Θγ 2 + 4Θγ α= , 4Θγ p −γ + γ 2 − 4Θγ 2 + 4Θγ β= , (4.7) 2Θ βγ + β 2 Θ + γKS1 χ14 = . r31 β(β + γ + KS2 )
45
4.2. Bistabiliteit in twee-staps cycli: wiskundig
4.2.2
Bifurcatie-analyse
βγ + β 2 Θ + γKS . r31 β(β + γ + KS ) Aangezien β uitgedrukt wordt in termen van γ en Θ (zie (4.7)) kunnen we χ14 weergeven als functie van γ door vaste waarden te kiezen voor Θ, KS en r31 (zie Figuur 4.4).
We nemen bijkomend aan dat KS1 = KS2 = KS en hebben dan χ14 =
Figuur 4.4: χ14 als functie van γ met KS = 0.01 en r31 = 2 voor drie verschillende parameterwaarden Θ: (A) Θ = 1/2, (B) Θ = 4, (C) Θ = 1.062 (m.b.v. Maple). Als we kijken naar Figuur 4.4, dan zien we dat we bij grafiek A geen bistabiliteit hebben. Bij grafiek B daarentegen merken we op dat er een interval van waarden van χ14 bestaat waarvoor het systeem wel bistabiliteit vertoont. Dit leiden we af door het aantal mogelijke waarden voor γ te bekijken bij een vaste waarde voor χ14 . In grafiek A bestaat er bij een vaste χ14 slechts ´e´en mogelijke waarde voor γ, terwijl er in grafiek B een gebied is waar er bij een vaste χ14 drie mogelijke waarden voor γ bestaan (twee stabiele equilibria met een onstabiel equilibrium ertussen). Grafiek C vormt hier de overgang tussen de twee mogelijkheden en dit is dus de grafiek waar we kunnen spreken van een bifurcatiepunt. Uit de figuur kunnen we nu duidelijk opmaken dat het fundamentele verschil tussen de ∂χ14 . Deze uitdrukking heeft namelijk geen drie grafieken te vinden is in de uitdrukking ∂γ ∂χ14 nulpunten in grafiek A (er geldt meer bepaald dat > 0), terwijl er in grafiek B twee ∂γ ∂χ14 nulpunten te vinden zijn van en in grafiek C ´e´en nulpunt. Voor het bifurcatiepunt ∂γ zoeken we dus een oplossing voor ∂χ14 = 0. ∂γ
46
(4.8)
4.2. Bistabiliteit in twee-staps cycli: wiskundig
Nulpunten van de vergelijking in (4.8) Met behulp van het wiskundig computerprogramma Maple kunnen we nu de nulpunten van deze afgeleide berekenen (zie A.1): s ! 1 (1 + KS )KS γ1,2 = 1+ ± A1 , (4.9) 2 Θ − 1 + KS (4Θ − 1) s ! (1 + KS )KS 1 1− ± A2 , (4.10) γ3,4 = 2 Θ − 1 + KS (4Θ − 1) met p √ √ KS 2 (−4Θ + 1) + Θ − 1 + 2 KS KS + 1 Θ − 1 + KS (4Θ − 1) A1 = , Θ − 1 + KS (4Θ − 1) s p √ √ KS 2 (−4Θ + 1) + Θ − 1 − 2 KS KS + 1 Θ − 1 + KS (4Θ − 1) . A2 = Θ − 1 + KS (4Θ − 1) s
Voorwaarde voor een re¨ ele oplossing Vooreerst merken we op dat er in het bifurcatiepunt wel degelijk een re¨ele wortel voor (4.8) moet bestaan, of dat minstens ´e´en van de voorgaande uitdrukkingen (4.9) en (4.10) voor γi , i ∈ {1, 2, 3, 4}, met andere woorden re¨eel moeten zijn. Uit de tweede term van de γ’s en de uitdrukking (1 + KS )KS > 0 volgt dan dat er moet gelden: Θ − 1 + KS (4Θ − 1) > 0 ⇔ (1 + 4KS )Θ > 1 + KS 1 + KS Θ > 1 + 4KS ⇔ 1 + KS Θ < 1 + 4KS
1 4 1 als KS < − 4 als KS > −
Aangezien we uit de definitie van KS weten dat deze positief is, is de voorwaarde van de tweede ongelijkheid onmogelijk. We hebben dus voor eender welke fysisch mogelijke waarde van KS dat 1 1 + KS Θ> > . (4.11) 1 + 4KS 4 Dit is meteen een voorwaarde die nodig is gedurende de volledige bespreking om te kunnen spreken van een bifurcatiepunt. Reduceren van oplossingen door ze te laten samenvallen ∂χ14 = 0. In het ∂γ bifurcatiepunt moeten we echter twee samenvallende waarden voor γ overhouden binnen het interval ]0, 1[ die voldoen aan deze uitdrukking.
Nu hadden we in (4.9) en (4.10) vier mogelijke waarden voor γ waar
47
4.2. Bistabiliteit in twee-staps cycli: wiskundig
Een eerste logische stap is om te kijken wanneer γ1 en γ2 , respectievelijk γ3 en γ4 samenvallen of met andere woorden wanneer A1 = 0, respectievelijk A2 = 0. Rekening houdend met de voorwaarde (4.11) kunnen we opmerken dat de noemers van A1 en A2 positief zijn. De voorwaarde impliceert ook dat de laatste wortel in de laatste term van de teller van A1 en A2 bestaat. Daaruit blijkt dat de laatste term in de teller van A1 (respectievelijk A2 ) positief (respectievelijk negatief) is. We hebben bijgevolg A1 = 0
⇒ KS 2 (−4Θ + 1) + Θ − 1 moet negatief zijn,
A2 = 0
⇒ KS 2 (−4Θ + 1) + Θ − 1 moet positief zijn.
We vinden dus het onderscheid in KS 2 (−4Θ + 1) + Θ − 1 = 0 ⇔ (1 − 4KS 2 )Θ = 1 − KS 2 1 − KS 2 . ⇔Θ= 1 − 4KS 2 Willen we een oplossing voor A1 = 0, dan hebben we de volgende voorwaarden: KS 2 (−4Θ + 1) + Θ − 1 < 0 1 − KS 2 Θ < 0 < KS < 12 2 1 − 4KS ⇔ 1 − KS 2 Θ > KS > 21 2 1 − 4KS Willen we een oplossing voor A2 = 0, dan hebben we de volgende voorwaarden: KS 2 (−4Θ + 1) + Θ − 1 > 0 2 Θ > 1 − K S 0 < KS < 1 − 4KS 2 ⇔ 2 1 − KS Θ < KS > 21 1 − 4KS 2
1 2
(4.12)
(4.13)
1 1 1 − KS 2 Voor 0 < KS < geldt: hebben dat 2 ∈ ]1, +∞[, terwijl we voor KS > 2 2 1 − 4KS 1 − KS 2 1 . Dit laatste leidt samen met het tweede geval in (4.13) tot de on2 ∈ −∞, 4 1 − 4KS 1 1 gelijkheid Θ < , wat in tegenstrijd is met de eerder gestelde voorwaarde Θ > uit (4.11). 4 4 1 We zien dus in dat A2 = 0 geen oplossing heeft als KS > . 2 We gaan nu de vergelijking A1 = 0 van naderbij bekijken: p p p A1 = 0 ⇔KS 2 (−4Θ + 1) + Θ − 1 = −2 KS KS + 1 Θ − 1 + KS (4Θ − 1) i
⇔KS 4 (−4Θ + 1)2 + 2KS 2 (−4Θ + 1)(Θ − 1) + (Θ − 1)2 = 4KS (KS + 1) (Θ − 1 + KS (4Θ − 1)) (1 + KS )2 (1 + KS )2 ii ⇔Θ1 = ∨ Θ = . 2 (2KS − 1)2 (2KS + 1)2 48
(4.14)
4.2. Bistabiliteit in twee-staps cycli: wiskundig
Voor ii werd Maple gebruikt (zie A.2). De equivalentie i werd bekomen door beide leden te kwadrateren en is enkel geldig mits er aan bepaalde bestaansvoorwaarden voldaan is: • BV 1 (volgt uit het noodzakelijk positief zijn van de term onder het wortelteken): Θ − 1 + KS (4Θ − 1) > 0. Deze voorwaarde is automatisch voldaan door de gestelde voorwaarde in (4.11). • BV 2 (volgt uit het negatief zijn van het rechterlid): KS 2 (−4Θ + 1) + Θ − 1 ≤ 0. We gaan deze voorwaarde na in de verschillende mogelijke gedaantes voor Θ en vullen met andere woorden Θ1 en Θ2 uit (4.14) in in deze uitdrukking en bekijken hun teken. We doen dit met behulp van Maple en geven de code weer in A.3: ◦ 6KS (KS + 1) KS 2 (−4Θ1 + 1) + Θ1 − 1 = − 2KS − 1
( < 0 ⇔ KS > > 0 ⇔ KS <
1 2 1 2
(4.15)
◦ KS 2 (−4Θ2 + 1) + Θ2 − 1 = −
2KS (KS + 1) <0 2KS + 1
(4.16)
Een analoge redenering kan opgebouwd worden voor het oplossen van de vergelijking A2 = 0. Hier hebben we echter als tweede bestaansvoorwaarde dat moet gelden: KS 2 (−4Θ + 1) + Θ − 1 > 0. Θ1 is dus enkel een oplossing van deze vergelijking als KS < 12 , terwijl Θ2 wegens deze tweede bestaansvoorwaarde nooit een oplossing is. We kunnen volgende besluiten samenvatten uit de berekeningen in deze deelparagraaf: 1 • A1 = 0 heeft 1 wortel Θ2 als KS < , 2 1 • A1 = 0 heeft 2 wortels Θ1 en Θ2 als KS > , 2 1 • A2 = 0 heeft 1 wortel Θ1 als KS < , 2 1 • A2 = 0 heeft geen wortels als KS > . 2
Consequentie voor de veranderlijken en fysische relevantie In het bifurcatiepunt hebben we dus de volgende waarden voor γ, waarbij gebruik gemaakt werd van Maple (zie A.4):
49
4.2. Bistabiliteit in twee-staps cycli: wiskundig
• in Θ1 : ◦ γ1,2
1 1 = + 2 2
p (2KS − 1)2 3
(uit Θ1 is enkel oplossing van A1 = 0 als KS >
1 volgt hier: 2KS − 1 > 0) 2
1 2KS − 1 + 2 6 1 + KS = . 3 =
◦ γ3,4
1 1 = − 2 2
p
(2KS − 1)2 3
(uit Θ1 is enkel oplossing van A2 = 0 als KS <
1 volgt hier: 2KS − 1 < 0) 2
1 1 − 2KS − 2 6 1 + KS . = 3 =
• in Θ2 : ◦ 1 1 + (2KS + 1) 2 2 = 1 + KS .
γ1,2 =
◦ γ3,4 : bestaat niet in Θ2 aangezien Θ2 geen wortel is van A2 . We merken op dat de twee overgebleven γ’s (γ1 en γ2 ) beschouwd in Θ2 geen rol spelen in ons model. We weten namelijk dat KS > 0 en door de definitie van γ geldt dat deze kleiner is dan 1, wat onmiddellijk tot een tegenspraak leidt, namelijk γ1,2 = 1 + KS > 1 in dit bifurcatiepunt. De enige fysiek mogelijke oplossing tot nu toe is dus Θ1 . We moeten geen rekening meer houden met de eerder gevonden oplossing Θ2 . Opnieuw door het gebruik van Maple vinden we in A.5 dat in het bifurcatiepunt geldt: 1 : 2 1 + KS (KS − 2)(2KS − 1) γ= en daaruit volgt dat β = − , 3 3(KS + 1) (KS − 2)2 2(KS2 + 2KS + 1)(KS + 1)2 α= en χ14 = − , 3(KS + 1) (2KS2 + 10KS − 1)(2KS − 1)(KS − 2)r31
• als KS >
50
4.2. Bistabiliteit in twee-staps cycli: wiskundig 1 • als KS < : 2 1 + KS 1 − 2KS 1 + KS γ= en daaruit volgt dat β = ,α= = γ en 3 3 3 KS + 1 χ14 = − . r31 (2KS − 1) 1 rekening met het feit dat α, β, γ ∈]0, 1[, dan kunnen we besluiten 2 dat voor een geldig bifurcatiepunt de restrictie KS < 2 moet genomen worden (uit α, γ < 1, β > 0).
Houden we voor KS >
Eliminatie van een laatste mogelijke oplossing 1 We bekijken dit eerste geval voor KS > opnieuw en gaan nu op een andere manier te 2 1 werk. We onderzoeken de voorwaarde A1 = 0 in het interval KS ∈ , 2 . Vooreerst bekijken 2 1 we wat dit wil zeggen voor γ en meer bepaald γ1,2 in KS = . 2 1 In (4.14) zien we meteen dat Θ1 niet bestaat in KS = aangezien de noemer dan 0 wordt. 2 1 1 We hebben dus Θ2 als enige wortel van A1 = 0 in KS = . In KS = geldt dan: 2 2 (1 + KS )2 (2KS + 1)2 3 2 9 2 = 2 = 2 16 s
Θ2 =
⇒ γ1,2
1 = 2
1+
3 2
7 − 16 √ 1 = 1+ 4 2 3 = 2
· 12 + 12 ·
! 5 4
1 In KS = geldt bijgevolg dat γ > 1, wat niet toegelaten is als we met een werkelijk biologisch 2 model werken. 1 Willen we nu een bifurcatiepunt voor KS binnen ons interval , 2 , dan is het noodzakelijk 2 dat we voor dit punt een γ vinden die kleiner is dan 1. Doordat γ1,2 continu is in zijn domein 1 en γ1,2 > 1 in KS = , hebben we een waarde voor KS en Θ nodig waarbij A1 = 0 en γ1,2 = 1, 2 zodat de overgang naar γ1,2 < 1 binnen ons interval mogelijk is. Door de uitdrukking van γ1,2 bij A1 = 0 (zie (4.9)) zien we dat de voorwaarde voor γ1,2 = 1 (1 + KS )KS te vinden is in de tweede term, namelijk F (KS ) := = 1. Voor deze net Θ − 1 + KS (4Θ − 1) 51
4.2. Bistabiliteit in twee-staps cycli: wiskundig 1 gedefinieerde uitdrukking F (KS ) wisten we al dat F = 4, wat er door het groter zijn 2 dan 1 voor zorgde dat ook γ1,2 > 1. ( A1 = 0, Onze opdracht bestaat er nu in om een oplossing KS te vinden met : F (KS ) = 1 p p p A1 = 0 ⇔ −KS 2 (−4Θ + 1) − Θ + 1 = 2 KS KS + 1 Θ − 1 + KS (4Θ − 1) (4.17) F (KS ) = 1 ⇔ (1 + KS )KS = Θ − 1 + KS (4Θ − 1) (4.18) • In het rechterlid van (4.17) vinden we onder de wortels van de middenste factoren het linkerlid van (4.18) terug. Vervangen we dit in (4.17), dan krijgen we de volgende uitdrukking: − KS 2 (−4Θ + 1) − Θ + 1 = 2 (Θ − 1 + KS (4Θ − 1)) ⇔ (4KS2 − 1 − 2 − 8KS )Θ = KS2 − 1 − 2 − 2KS K 2 − 2KS − 3 , ⇔ Θ = S2 4KS − 8KS − 3
(4.19)
mits de bestaansvoorwaarde 4KS2 − 8KS − 3 6= 0 √ 8±4 7 ⇔ KS 6= 8√ ⇔ KS 6= 1 ± 7 ⇔ KS 6= −1.64575 ∧ KS 6= 3.64575. 1 Aangezien we aan het werken zijn binnen KS ∈ , 2 moeten we verder geen rekening 2 houden met deze bestaansvoorwaarde. • Doordat we ook het rechterlid van (4.18) terugvinden in het rechterlid van (4.17) (onder de wortel van de laatste term), kunnen we dit op een analoge manier vervangen in (4.17): − KS 2 (−4Θ + 1) − Θ + 1 = 2KS (KS + 1) ⇔ (4KS2 − 1)Θ = 2KS (KS + 1) + (KS − 1)(KS + 1) (3KS − 1)(KS + 1) , ⇔Θ= 4KS2 − 1
(4.20)
mits de bestaansvoorwaarde 4KS2 − 1 6= 0 ⇔ (2KS − 1)(2KS + 1) 6= 0 1 1 ⇔ KS 6= − ∧ KS 6= 2 2 vervuld is. Opnieuw is dit een bestaansvoorwaarde die niet van belang is binnen ons interval. 52
4.2. Bistabiliteit in twee-staps cycli: wiskundig We hadden reeds eerder aangezien dat A1 = 0 twee oplossingen had voor Θ als KS > 12 . Aangezien we een stelsel aan het oplossen zijn waarbij zowel aan deze vergelijking (4.17) moet voldaan zijn als aan (4.18) (waarbij deze laatste lineair is in Θ) merken we op dat Θ maar ´e´en unieke oplossing kan hebben. Het is dus noodzakelijk dat onze berekende waarden voor Θ in (4.19) en (4.20) samenvallen en aangezien dit uitdrukkingen zijn in KS , kunnen we daaruit meteen de waarde voor KS halen waarbij aan dit stelsel voldaan is. 3KS2 + 2KS − 1 KS2 − 2KS − 3 = 4KS2 − 8KS − 3 4KS2 − 1 ⇔ 4KS4 − 8KS3 − 13KS2 + 2KS + 3 = 12KS4 − 16KS3 − 29KS2 + 2KS + 3 ⇔ KS2 (8KS2 − 8KS − 16) = 0 ⇔ KS2 − KS − 2 = 0 ⇔ KS = −1 ∨ KS = 2. Met dit eerste geval moet geen rekening gehouden worden door de voorwaarde KS > 0. 3 1 We hebben dus dat γ = 1 voor KS = 2, wat samen met γ = voor KS = ervoor zorgt 2 2 1 dat we geen biologisch toegestaan bifurcatiepunt bekomen voor KS ∈ , 2 . 2 Besluit voor het bifurcatiepunt Samen met het eerder uitgesloten gebied hebben we dus geen fysisch relevant bifurcatiepunt 1 voor KS > en in ons gezochte bifurcatiepunt geldt dat: 2 1 KS < , 2 (1 + KS )2 , Θ= (2KS − 1)2 1 + KS γ= = α, 3 1 − 2KS β= , 3 KS + 1 χ14 = − . r31 (2KS − 1)
(4.21)
Voorwaarden voor bistabiliteit We hebben eerder berekend dat α + β + γ = 1 en door het positief zijn van deze drie dimensieloze concentraties kunnen we de volgende voorwaarden berekenen: 1 + KS <1 3 ⇔ 0 < 1 + KS < 3 ⇔ −1 < KS < 2,
α, γ ∈]0, 1[ ⇔ 0 <
53
4.2. Bistabiliteit in twee-staps cycli: wiskundig 1 − 2KS <1 3 ⇔ 0 < 1 − 2KS < 3 ⇔ −1 < −2KS < 2 1 ⇔ > KS > −1. 2
β ∈]0, 1[ ⇔ 0 <
Uit de definitie van KS volgt dat deze noodzakelijk positief moet zijn. Een restrictie die hier 1 opgelegd wordt zodat we te maken hebben met een bifurcatiepunt, is KS < (en dit hadden 2 we eerder in ons onderzoek ook al bekomen). Samengevat hebben we dus te maken met bistabiliteit bij parameterwaarden waarbij aan de 1 (1 + KS )2 voorwaarde KS < voldaan is en ook geldt dat Θ > . Dit wordt ge¨ıllustreerd in 2 (1 − 2KS )2 Figuur 4.5.
Figuur 4.5: Verloop van de asymmetrische factor Θ als functie van de Michaelis constante KS in het bifurcatiepunt. Bereik veranderlijken en parameterwaarden in het bifurcatiepunt 1 en de vergelijkingen in (4.21) kunnen we Aan de hand van dit restrictie-interval KS ∈ 0, 2 het bereik van de andere parameterwaarden en veranderlijken bepalen in het bifurcatiepunt:
54
4.3. Bistabiliteit in twee-staps cycli: numeriek
• KS → 0: Θ→1 α=γ→ β→
1 3
1 3
χ14 →
1 . r31
1 • KS → : 2 Θ → +∞ α=γ→
1 2
β→0 χ14 → +∞. We hebben dus de volgende voorwaarden voor de parameterwaarden van het bifurcatiepunt: Θ ∈]1, +∞[ 1 1 α, γ ∈ , 3 2 1 β ∈ 0, 3 1 , +∞ . χ14 ∈ r31 Uit het besluit over Θ kunnen we als nodige voorwaarde voor bistabiliteit stellen dat Θ > 1 k2 k3 in dit geval (namelijk KS1 = KS2 = KS3 = KS4 = KS ). Doordat Θ = , is het noodzakelijk k1 k4 voor bistabiliteit dat k2 k3 > k1 k4 . Een voorwaarde voor bistabiliteit is dus dat het product van de katalytische constanten van de fosforylatie en de defosforylatie van Wβ groter moet zijn dan het product van de katalytische constanten van de fosforylatie van Wα en de defosforylatie van Wγ .
4.3
Bistabiliteit in twee-staps cycli: numeriek
In de voorgaande paragraaf bekeken we enkele voorwaarden die nodig waren voor het bereiken van bistabiliteit in een twee-staps cyclus. We willen deze voorwaarden nu numeriek nagaan en experimenteren daarom met verschillende parameterwaarden voor dit model in Matcont (zie ook paragraaf 4.1.1).
55
4.3. Bistabiliteit in twee-staps cycli: numeriek
4.3.1
Parameterrestricties voor een stabiele evenwichtstoestand
Aangezien we hier de vaste waarde KS = 0.01 gekozen hebben en deze kleiner is dan 1/2, weten we uit de wiskundige analyse in paragraaf 4.2 dat het systeem ´e´en evenwichtstoestand (1 + KS )2 heeft voor parameterwaarden die voldoen aan Θ < . De kritieke waarde van Θ in (1 − 2KS )2 het bifurcatiepunt is in dit geval (KS = 0.01) 1.06 en lagere waarden zouden dus moeten leiden tot ´e´en evenwichtstoestand. We tonen dit numeriek aan in Figuur 4.6 door verschillende bifurcatiediagrammen te tekenen bij het vari¨eren van χ14 , waarbij we verschillende waarden kiezen voor Θ en r31 . We houden er wel rekening mee dat Θ steeds onder de kritieke grens van 1.06 blijft.
Figuur 4.6: E´en-parameter bifurcatiediagram van respectievelijk α, β en γ bij het vari¨eren van χ14 met de andere parameterwaarden zoals in Tabel 4.1 en aangeduid op de grafiek (Model 3). In Figuur 4.6 zien we dat de ultrasensitiviteit niet enkel afhankelijk is van de waarde van KS , maar dit ook bepaald wordt door de afstand van de asymmetrische factor Θ tot zijn kritische waarde. Onder ultrasensitiviteit verstaan we het fenomeen waarbij een kleine storing op een bepaalde parameter een grote impact heeft op een veranderlijke, wat we hier kunnen aflezen door de stijlheid van de grafieken. Bij Θ = 1 zien we een meer abrupte overgang van hoge α naar lage α, die gepaard gaat met de abrupte overgang van lage γ naar hoge γ. 56
4.3. Bistabiliteit in twee-staps cycli: numeriek
4.3.2
Parameterrestricties die bistabiliteit toestaan
(1 + KS )2 een interval para(1 − 2KS )2 meterwaarden voor χ14 bestaat dat bistabiliteit impliceert bij deze vaste waarde KS = 0.01. We hebben dan meer bepaald twee stabiele evenwichtstoestanden en ´e´en onstabiele. We tonen dit opnieuw numeriek aan door de verschillende bifurcatiediagrammen te tekenen bij het vari¨eren van χ14 waarbij we verschillende waarden kiezen voor Θ en r31 . Θ is in dit geval groter dan de kritieke grens van 1.06 gekozen (Figuur 4.7). We weten dat er bij parameterwaarden die voldoen aan Θ >
Figuur 4.7: E´en-parameter bifurcatiediagram van respectievelijk α, β en γ bij het vari¨eren van χ14 met de andere parameterwaarden zoals in Tabel 4.1 en aangeduid op de grafiek (Model 3). In Figuur 4.7 hebben we steeds met twee limietpunten te maken en door het beschouwen van de eigenwaarden merken we inderdaad dat het deel van de equilibriumkromme tussen de twee limietpunten onstabiel is en de andere delen stabiel.
57
4.3. Bistabiliteit in twee-staps cycli: numeriek
4.3.3
Overzicht van equilibria en de mogelijkheid tot bistabiliteit
Zoals in de voorgaande paragrafen duidelijk aangetoond werd met behulp van verschillende implementaties en berekeningen, zijn er gebieden in de ruimte van de parameters waarin bistabiliteit mogelijk gemaakt wordt naast het voorkomen van ´e´en equilibrium. Zo zagen we in paragraaf 4.3.1 dat er bepaalde waarden waren voor Θ (in het geval van KS = 0.01 meer bepaald waarden voor Θ die kleiner zijn dan 1.06) waar we ongeacht van de waarde voor χ14 slechts ´e´en equilibrium terug vinden in ons systeem. In paragraaf 4.3.2 vonden we dan voor Θ > 1.06 dat het aantal equilibria afhankelijk was van de waarde van χ14 en dat er gebieden zijn waar bistabiliteit mogelijk is. We kunnen hieruit afleiden wat we verwachten te zien als we een twee-parameter bifurcatiediagram opstellen met Θ en χ14 als vrije parameters. Bij Θ > 1.06 (in het geval KS = 0.01) kunnen we verwachten dat er twee limietpuntkrommen zijn die het gebied afbakenen waar er bistabliliteit mogelijk is, terwijl we voor Θ < 1.06 geen limietpuntkrommen zullen opmerken. Deze vorm van limietpuntkrommen impliceert een Cusp-punt in het hoekpunt van de limietpuntkrommen en deze vermoedens worden bevestigd als we deze continuatie uitvoeren in Figuur 4.8, waarbij we r31 = 2 gesteld hebben.
Figuur 4.8: Continuatie van LP bij het vari¨eren van Θ en χ14 en de andere parameterwaarden zoals in Tabel 4.1 en r31 = 2 (Model 3). Het Cusp-punt bevindt zich in (α, β, Θ, χ14 ) = (0.336957; 0.326667; 1.062162; 0.515306), waar we inderdaad de berekende kritieke waarde voor Θ in terug vinden.
58
4.4. βmax
4.4
βmax
Uit de voorgaande numerieke experimenten blijkt dat β bij evenwicht een bovengrens heeft die lager is dan de door het model opgelegde waarde 1 (zie Figuur 4.6 en 4.7). Deze maximale waarde lijkt af te hangen van de gekozen waarde voor Θ. We willen nu expliciet berekenen wat deze maximale waarde is. Dit kunnen we doen door de afgeleide van β naar de veranderlijke χ14 te berekenen en te kijken wanneer deze nul is. Als we kijken naar (4.7), dan merken we meteen op dat deze afgeleide niet zo simpel te berekenen is aangezien de veranderlijke χ14 niet voorkomt in de uitdrukking voor β, maar dit wel omgekeerd het geval is. We herschrijven onze te berekenen uitdrukking dus als volgt: ∂β =0⇔ ∂χ14
∂β ∂γ ∂χ14 ∂γ
= 0.
Deze uitdrukking heeft een oplossing als de teller gelijk is aan nul en de noemer niet. In paragraaf 4.2 hebben we al berekend wanneer de noemer gelijk is aan nul. We houden dit nu in ons achterhoofd bij het berekenen van de uitdrukking waar de teller nul is. Opnieuw zullen we ons hier onnodig rekenwerk besparen door gebruik te maken van het wiskundig computerprogramma Maple. De resultaten worden weergegeven in A.6. Deze appendix wordt opgesplitst in verschillende stukjes, zodat de redenering nog steeds duidelijk te volgen is achter het programmeerwerk. We berekenen dus eerst oplossingen voor de uitdrukking ∂β =0 ∂γ
(4.22)
in A.6.1 en houden er > 0. We hebben dan twee oplossingen voor (4.22), √ rekening mee dat Θ√ 2Θ + Θ −2Θ + Θ namelijk γ = en γ = − . −1 + 4Θ −1 + 4Θ Nu willen we nagaan of er voldaan is aan de verschillende bestaansvoorwaarden en of de oplossingen dus wel mogelijk zijn. • BV 1: β moet re¨eel zijn als ze beschouwd wordt in de berekende γ’s. Door de aanwezigheid van een term onder de wortel in de uitdrukking voor β is deze voorwaarde equivalent met de voorwaarde dat de term onder de wortel niet-negatief is, of met andere woorden γ 2 − 4Θγ 2 + 4Θγ ≥ 0. Het is dezelfde term die voorkomt onder het ∂β . wortelteken in de uitdrukking voor ∂γ • BV 2: In A.6.1 worden de oplossingen van (4.22) gezocht door de vergelijking B = 0 te berkenenen. Als we de berekening van Maple phier volgen, dan komt dit neer op het zoeken van de oplossing van de vergelijking 2 γ 2 − 4Θγ 2 + 4Θγ = 2γ − 8Θγ + 4Θ. Maple zal hier beide leden kwadrateren om tot een oplossing te komen, maar hierbij wordt niet nagegaan of het rechterlid wel degelijk positief is, wat noodzakelijk is door het positief zijn van het linkerlid. Als tweede bestaansvoorwaarde hebben we dus nodig dat 2γ − 8Θγ + 4Θ > 0.
59
4.4. βmax Deze bestaansvoorwaarden worden nagegaan in A.6.2. We houden dus ´e´en oplossing γ over die voldoet aan (4.22): √ −2Θ + Θ γ=− −1 √ + 4Θ √ − Θ(2 Θ − 1) √ =− √ (2 Θ − 1)(2 Θ + 1) √ Θ = √ . 2 Θ+1 Met behulp van deze waarde voor γ kunnen we nu eveneens de andere veranderlijken berekenen in dit punt (zie A.6.3). 1 β= √ , 2 Θ+1 √ Θ = γ, α= √ 2 Θ+1 √ √ Θ + Θ + 2ΘKS + ΘKS χ14 = √ √ 1 + Θ + 2 ΘKS + KS r31 √ √ √ Θ 1 + Θ + 2 ΘKS + KS = √ √ 1 + Θ + 2 ΘKS + KS r31 √ Θ = r √31 r24 =√ , r31 waarbij de laatste gelijkheid volgt uit de definitie van de asymmetrische factor Θ = r31 r24 . We merken op dat dit een oplossing is voor (4.22), maar we nog niet berekend hebben of dit een oplossing is waarbij β minimaal of maximaal is. We kunnen wel al intu¨ıtief vermoeden dat dit een maximale waarde is als we kijken naar Figuur 4.6 en 4.7. Voor de effectieve berekening echter bekijken we de tweede afgeleide van β naar γ in A.6.4 en zien we dat we door het negatief zijn van deze tweede afgeleide te maken hebben met een maximale waarde 1 . Deze is enkel afhankelijk van de asymmetrische factor Θ voor β, namelijk βmax = 2√Θ+1 √ r24 en wordt bereikt in χ14 = √ . De overige concentraties in dit specifieke punt zijn dan r31 √ Θ α=γ= √ . 2 Θ+1
60
Hoofdstuk 5 Bistabiliteit door driedubbele fosforylatie Dit hoofdstuk is gebaseerd op [8].
5.1
Drie-staps cycli: beschrijving van het model
In dit hoofdstuk werken we met het niet-specifieke eiwit W in een distributief mechanisme met covalente modificatie.
Figuur 5.1: Kinetisch diagram waarbij het eiwit W vier verschillende vormen Wα , Wβ , Wγ en Wδ heeft die gevormd worden door de inwerking van de kinase e1 en de fosfatase e2 . We nemen nu aan dat het eiwit W vier vormen kan aannemen, namelijk Wα , Wβ , Wγ en Wδ . Hierbij stelt Wα de niet-gefosforyleerde vorm van het eiwit voor, Wβ de gefosforyleerde vorm van Wα , Wγ de gefosforyleerde vorm van Wβ (of dus de dubbel gefosforyleerde vorm van Wα ) en Wδ de gefosforyleerde vorm van Wγ (of dus de driedubbel gefosforyleerde vorm van Wα ). We nemen aan dat de drie fosforylatiestappen door hetzelfde enzym e1 gekatalyseerd worden (de kinase) en de drie defosforylatiestappen door het enzym e2 (de fosfatase). 61
5.1. Drie-staps cycli: beschrijving van het model
We veronderstellen dat elke (de)fosforylatiestap volgend distributief MM-mechanisme volgt: kai ki WS + e eWS → WP + e, kdi
(5.1)
waarbij kai , kdi en ki respectievelijk de samenstellings-, splitsings- en katalytische constante zijn van stap i, i ∈ {1, . . . , 6}. eWS is het MM-complex dat gevormd werd door het enzym e en zijn substraat WS om het product WP te produceren. Hierbij wordt verondersteld dat de nodige energie en fosfaatgroepen in constante hoeveelheid aanwezig zijn en ingerekend werden in de kinetische constanten. We merken op dat we hier (in tegenstelling tot in [6] of Hoofdstuk 2) enkel rekening houden met het vormen van het substraat-enzym-complex en niet met het vormen van het product-enzym-complex, wat de analyse zal vereenvoudigen. Nemen we nu vi als de snelheid waarmee stap i uitgevoerd wordt, dan zien we aan de hand van Figuur 5.1 dat de differentiaalvergelijkingen die de tijdsevolutie van dit systeem weergeven als volgt zijn: dWα dt dWβ dt dWγ dt dWδ dt
= v2 − v1 , = v1 − v2 + v4 − v3 ,
(5.2)
= v3 − v4 + v6 − v5 , = v5 − v6 .
We weten dat de enzymen zich zowel in gebonden als ongebonden vorm manifesteren en we hebben meer bepaald dat e1T = [e1 ] + [e1 Wα ] + [e1 Wβ ] + [e1 Wγ ], e2T = [e2 ] + [e2 Wδ ] + [e2 Wγ ] + [e2 Wβ ], waarbij eiT de totale concentratie van het enzym ei voorstelt. We nemen aan dat de totale concentratie van het enzym e1 en e2 gedurende de cyclus constant is. Als we nu ook veronderstellen dat er evenwicht is tussen de opbouw en de afbraak van de intermediairen in (5.1), dan bekomen we de volgende MM-vergelijkingen (analoge uitwerking als in paragraaf 2.5): v1 =
v3 =
v5 =
[W ]
α] Vm1 [W Km1
1+
1+
1+
[Wβ ] Km3 [Wβ ] Vm3 Km3 [Wβ ] [Wα ] + Km3 Km1 [Wγ ] Vm5 K m5 [Wβ ] [Wα ] + Km3 Km1 [Wα ] Km1
+
+
[Wγ ] Km5
+
[Wγ ] Km5
+
[Wγ ] Km5
,
v2 =
,
v4 =
,
v6 =
62
β Vm2 Km2
1+
1+
1+
[Wγ ] Km4 [Wγ ] Vm4 K m4 [Wγ ] [Wδ ] + Km4 Km6 [Wδ ] Vm6 K m6 [Wγ ] [Wδ ] + Km4 Km6 [Wδ ] Km6
+
+
[Wβ ] Km2
+
[Wβ ] Km2
+
[Wβ ] Km2
,
,
.
5.1. Drie-staps cycli: beschrijving van het model
Hierbij hebben we de volgende notaties ingevoerd: kdi + ki Kmi = , i = 1, 2, 3, 4, 5, 6 kai en Vm1 = k1 e1T , Vm3 = k3 e1T , Vm5 = k5 e1T ,
Vm2 = k2 e2T , Vm4 = k4 e2T , Vm6 = k6 e2T .
Om de wiskundige analyse te vereenvoudigen, kunnen we de concentratie van de veranderlijken Wα , Wβ , Wγ en Wδ dimensieloos maken. Hierbij beschouwen we de totale concentratie van het eiwit. We nemen aan dat deze constant is en door het combineren van de gegevens in Figuur 5.1 en (5.1) is deze van de volgende vorm: WT = [Wα ] + [Wβ ] + [Wγ ] + [Wδ ] + [e1 Wα ] + [e1 Wβ ] + [e1 Wγ ] + [e2 Wδ ] + [e2 Wγ ] + [e2 Wβ ]. De dimensieloze concentraties van de veranderlijken worden dan weergegeven door α=
[Wα ] , WT
[Wβ ] , WT
β=
γ=
[Wγ ] , WT
δ=
[Wδ ] , WT
en deze liggen bijgevolg allemaal tussen 0 en 1. Samen met (5.2) bekomen we de volgende differentiaalvergelijkingen voor de dimensieloze veranderlijken voor de tijdsevolutie van dit systeem: dα dt dβ dt dγ dt dδ dt
v2 − v1 , WT v1 − v2 + v4 − v3 = , WT v3 − v4 + v6 − v5 = , WT v5 − v6 = . WT =
Voeren we nu analoog de volgende herschaling in: KS1 =
Km1 , WT
KS2 =
Km2 , WT
KS3 =
Km3 , WT
KS4 =
Km4 , WT
KS5 =
Km5 , WT
KS6 =
Km6 , WT
dan kunnen we onze MM-vergelijkingen herleiden tot v1 = v3 = v5 =
Vm1 KαS1 1+ 1+ 1+
β KS3 Vm3 KβS3 α + KβS3 KS1 Vm5 KγS5 α + KβS3 KS1 α KS1
+
+
γ KS5
+
γ KS5
+
γ KS5
,
v2 =
,
v4 =
,
v6 =
63
Vm2 KβS2 1+ 1+ 1+
γ KS4 Vm4 KγS4 δ + KγS4 KS6 Vm6 KδS6 δ + KγS4 KS6 δ KS6
+
+
β KS2
+
β KS2
+
β KS2
, , .
(5.3)
5.1. Drie-staps cycli: beschrijving van het model
Voor het verdere verloop van deze analyse defini¨eren we nog enkele dimensieloze parameters: Vm3 k3 Vm2 k2 = , r24 = = , Vm1 k1 Vm4 k4 k5 Vm4 k4 Vm5 = , r46 = = , r53 = Vm3 k3 Vm6 k6 Vm1 k1 e1T χ16 = = . Vm6 k6 e2T r31 , r24 , r53 en r46 zijn hierbij duidelijk onafhankelijk van de enzym-concentraties. Ook hier maken we de vereenvoudiging dat WT e1T en WT e2T en dus dat [WS ] [eWS ]. De laatste zes termen van WT = [Wα ] + [Wβ ] + [Wγ ] + [Wδ ] + [e1 Wα ] + [e1 Wβ ] + [e1 Wγ ] + [e2 Wδ ] + [e2 Wγ ] + [e2 Wβ ] zijn bijgevolg verwaarloosbaar t.o.v. de eerste vier termen, dus kunnen we deze termen laten wegvallen. Door het dimensieloos maken van de vergelijking (beide leden delen door WT ) krijgen we r31 =
α + β + γ + δ = 1.
(5.4)
Model 4. v2 − v1 , WT v3 − v4 + v6 − v5 = , WT v5 − v6 = , WT = 1 − α − γ − δ, Vm1 KαS1 , v1 = 1 + KαS1 + KβS3 + KγS5
dα dt dγ dt dδ dt β
=
v2 = v3 = v4 = v5 = v6 =
Vm2 KβS2 1+ 1+ 1+ 1+ 1+
γ KS4 Vm3 KβS3 α + KβS3 KS1 Vm4 KγS4 δ + KγS4 KS6 Vm5 KγS5 α + KβS3 KS1 Vm6 KδS6 δ + KγS4 KS6 δ KS6
+
64
+
β KS2
+
γ KS5
+
β KS2
+
γ KS5
+
β KS2
, , , , .
5.2. Bistabiliteit in drie-staps cycli: numeriek
5.2
Bistabiliteit in drie-staps cycli: numeriek
De voorgaande resultaten leiden tot Model 4. De implementatie van dit model in Matcont is weergegeven in Figuur 5.2 en de bijhorende parameterwaarden in Tabel 5.1.
Figuur 5.2: Implementatie van Model 4 (pagina 64) in Matcont.
Parameter Vm1 KS χ16 WT
Waarde 1 0.01 1.5 500
Dimensie nM s−1 1 1 nM
Tabel 5.1: Parameterwaarden bij Model 4. In Figuur 5.3 bekijken we het bifurcatiediagram met χ16 als vrije parameter voor verschillende waarden van Θ.
65
5.2. Bistabiliteit in drie-staps cycli: numeriek
Figuur 5.3: E´en-parameter bifurcatiediagram bij het vari¨eren van χ14 met de andere parameterwaarden zoals in Tabel 5.1 en aangeduid op de grafieken (Model 4). Er is geen fundamenteel verschil met Model 3 dat numeriek onderzocht werd in paragraaf 4.3. We hebben opnieuw een grens voor Θ waaronder we ´e´en equilibrium hebben voor alle waarden van χ16 en boven die grens hebben we een interval van waarden voor χ16 waar bistabiliteit mogelijk is, terwijl er buiten dat interval terug slechts plaats is voor ´e´en equilibrium. Zitten we in het geval waar slechts ´e´en equilibrium mogelijk is (Θ beneden zijn grens) dan zijn alle equilibria stabiel. Boven die grens hebben we echter twee limietpunten waartussen onstabiele equilibria voorkomen, erbuiten zijn ze stabiel. Het verschil in aantal equilibria lezen we ook af als we een twee-parameter bifurcatiediagram tekenen waarbij de parameters Θ en χ16 vrij genomen worden (zie Figuur 5.4). Het Cusp-punt bevindt zich in (α, β, δ, Θ, χ16 ) = (0.252575; 0.247475; 0.252425; 1.020202; 0.510101). Hieruit kunnen we zonder de analytische berekening afleiden dat de kritieke grens voor Θ tussen ´e´en equilibrium of de mogelijkheid tot meerdere equilibria ligt bij Θ = 1.020202.
66
5.2. Bistabiliteit in drie-staps cycli: numeriek
Figuur 5.4: Continuatie van LP bij het vari¨eren van Θ en χ16 met de andere parameterwaarden zoals in Tabel 5.1 en r31 = 2 (Model 4). In Hoofdstuk 4 hadden we te maken met een twee-staps cyclus die leidde tot bistabiliteit. In Hoofdstuk 5 voegden we nu een extra stap toe aan het systeem met het idee in het achterhoofd om misschien multistabiliteit (en dus meer dan twee stabiele equilibria) te detecteren. Ondanks experimenten met verschillende parameterwaarden werd echter steeds ten hoogste bistabiliteit teruggevonden in dit model.
67
Hoofdstuk 6 Multistabiliteit door een cascade van dubbele fosforylaties Dit hoofdstuk is gebaseerd op [8].
6.1
Twee-staps cycli in een cascade: beschrijving van het model
In dit hoofdstuk wordt het derde model besproken uit [8]. Dit is meer bepaald het model dat in dit artikel weergegeven wordt als Figuur 8(B) (hier Figuur 6.1).
Figuur 6.1: Kinetisch diagram waarbij de eiwitten W en Z drie verschillende vormen Wα , Wβ en Wγ , respectievelijk Zα , Zβ en Zγ kunnen aannemen die gevormd worden door de inwerking van de kinase e1 en de fosfatase e2 , respectievelijk Wγ en e3 .
68
6.1. Twee-staps cycli in een cascade: beschrijving van het model
Uit deze figuur kunnen we volgend stelsel van differentiaalvergelijkingen opstellen: dWα dt dWβ dt dWγ dt dZα dt dZβ dt dZγ dt
= v2 − v1 , = v1 − v2 + v4 − v3 , = v3 − v4 ,
(6.1)
= v6 − v5 , = v5 − v6 + v8 − v7 , = v7 − v8 ,
waarbij vi de snelheid is van de reactie in stap i (i = 1, . . . , 8). Hoewel Wγ ook gebruikt wordt in stap 5 en 7 (als kinase) wordt het verbruik van Wγ in deze stappen niet opgenomen in de differentiaalvergelijking voor de verandering van Wγ . Dit is te verantwoorden door het eerder gestelde feit dat een enzym onveranderd gerecupereerd wordt eens het zijn werking vervuld heeft, zodat we niet kunnen spreken van verlies of winst van Wγ in stap 5 of 7. In de vergelijkingen voor Zα , Zβ en Zγ speelt Wγ de rol die e1 speelt in de vergelijkingen voor Wα , Wβ en Wγ . Opnieuw aan de hand van Figuur 6.1 en rekening houdend met het distributief gedrag verkrijgen we de volgende reactievergelijkingen voor de (de)fosforylatie van de verschillende vormen van de eiwitten W en Z: ka1 k1 Wα + e1 e1 Wα → Wβ + e1 , kd1 ka3 k3 Wβ + e1 e1 Wβ → Wγ + e1 , kd3 ka4 Wγ + e2 e2 Wγ kd4 ka2 Wβ + e2 e2 Wβ kd2
69
k4 → Wβ + e2 , k2 → Wα + e2 ,
6.1. Twee-staps cycli in een cascade: beschrijving van het model ka5 Zα + Wγ Wγ Zα kd5 ka7 Zβ + Wγ Wγ Zβ kd7 ka8 Zγ + e3 e3 Zγ kd8 ka6 Zβ + e3 e3 Zβ kd6
k5 → Zβ + Wγ , k7 → Zγ + Wγ ,
k8 → Zβ + e3 , k6 → Zα + e3 ,
waarbij kai , kdi en ki respectievelijk de samenstellings-, splitsings- en katalytische constante zijn van stap i, i ∈ {1, . . . , 8}. Bij deze reactievergelijkingen veronderstellen we dat de nodige energie en fosfaatgroepen in constante hoeveelheid aanwezig zijn en ingerekend werden in de kinetische constanten. We stellen de MM-vergelijkingen van dit wat uitgebreider model op en nemen hierbij aan dat er evenwicht is tussen de opbouw en de afbraak van de intermediairen in de reactievergelijkingen. Op die manier krijgen we op het eerste zicht heel vergelijkbare uitdrukkingen met deze in Model 3. v1 = v3 = v5 = v7 =
Vm1 KαW S1 1+
αW KS1
+
βW KS3
Vm3 KβW S3 1+
αW KS1
+
βW KS3
Z Vm5 KαS5
1+
αZ KS5
+
βZ KS7
Z Vm7 KβS7
1+
αZ KS5
+
βZ KS7
,
v2 =
,
v4 =
,
v6 =
,
v8 =
Vm2 KβW S2 1+
γW KS4
+
βW KS2
Vm4 KγWS4 1+
γW KS4
+
βW KS2
Z Vm6 KβS6
1+
γZ KS8
+
βZ KS6
Z Vm8 KγS8
1+
γZ KS8
+
βZ KS6
Hierbij hebben we de volgende notaties ingevoerd: kdi + ki Kmi = , i = 1, 2, 3, 4, 5, 6, 7, 8 kai en Vm1 Vm2 Vm3 Vm4
= k1 e1T , = k2 e2T , = k3 e1T , = k4 e2T ,
Vm5 Vm6 Vm7 Vm8 70
= k5 [Wγ ], = k6 e3T , = k7 [Wγ ], = k8 e3T .
, , , .
(6.2)
6.1. Twee-staps cycli in een cascade: beschrijving van het model
We merken op dat we hier wel een fundamenteel verschil hebben met de eerder ingevoerde modellen. Bij het constant houden van de verschillende parameterwaarden bleven de nieuw ingevoerde parameters Vmi constant onder de tijdsintegratie. Hier geldt dit opnieuw voor bepaalde Vmi ’s, maar andere (meer bepaald Vm5 en Vm7 ) zijn veranderlijk bij tijdsintegratie. Ze zijn namelijk beiden afhankelijk van [Wγ ] die enkel constant beschouwd wordt in de evenwichtsoestand v3 = v4 (zie differentiaalvergelijking in (6.1)). Om te werken met dimensieloze concentraties van de veranderlijken voeren we een herschaling door, meer bepaald werden volgende grootheden en parameters ingevoerd: [Wβ ] , WT [Zβ ] βZ = , ZT
[Wα ] , WT [Zα ] αZ = , ZT
[Wγ ] , WT [Zγ ] γZ = , ZT
βW =
αW =
Kmi , WT Kmi = , ZT
γW =
KSi =
i ∈ {1, . . . , 4},
KSi
i ∈ {5, . . . , 8},
waarbij WT de totale concentratie is van W en ZT de totale concentratie van Z. Beide totale concentraties worden constant verondersteld en door het combineren van de gegevens in Figuur 6.1 en de opgestelde reactievergelijkingen zien deze totale concentraties er als volgt uit: WT = [Wα ] + [Wβ ] + [Wγ ] + [e1 Wα ] + [e1 Wβ ] + [e2 Wγ ] + [e2 Wβ ] + [Wγ Zα ] + [Wγ Zβ ], (6.3) ZT = [Zα ] + [Zβ ] + [Zγ ] + [Wγ Zα ] + [Wγ Zβ ] + [e3 Zγ ] + [e3 Zβ ]. (6.4) Ook bij dit model worden dimensieloze grootheden ingevoerd voor het verdere verloop van de analyse, namelijk: Vm3 = Vm1 Vm7 = = Vm5 = r31 r24 , Vm1 = = Vm8
r31 = r75 ΘW χ18
k3 , k1 k7 , k5
Vm2 = Vm4 Vm6 r68 = = Vm8 ΘZ = r75 r68 , r24 =
k2 , k4 k6 , k8
(6.5)
k1 e1T . k8 e3T
Als laatste maken we de vereenvoudiging dat er veel meer eiwitten aanwezig zijn in ons model dan enzymen, waardoor er veel meer vrije eiwitten te vinden zijn dan eiwitten die gebonden zijn in een complex. Op deze manier zijn de laatste zes termen in de uitdrukking voor WT (zie (6.3)) verwaarloosbaar t.o.v. de eerste drie termen, analoog voor de laatste vier termen in de uitdrukking voor ZT (zie (6.4)). We kunnen deze termen dus laten wegvallen in (6.3) en (6.4) en na het dimensieloos maken van beide vergelijkingen krijgen we αW + βW + γW = 1, αZ + βZ + γZ = 1. 71
(6.6) (6.7)
6.2. Multistabiliteit in twee-staps cycli in een cascade: numeriek
6.2
Multistabiliteit in twee-staps cycli in een cascade: numeriek
Samengevat bekijken we het volgende model (waarbij we de herschaling doorgevoerd hebben in de differentiaalvergelijkingen uit (6.1)): Model 5. v1 − v2 + v4 − v3 , WT v3 − v4 = , WT = 1 − βW − γW , v5 − v6 + v8 − v7 = , ZT v7 − v8 = , ZT = 1 − βZ − γZ , Vm1 KαW S1 v1 = , αW 1 + KS1 + KβW S3
dβW dt dγW dt αW dβZ dt dγZ dt αZ
=
v2 = v3 = v4 = v5 = v6 = v7 = v8 =
Vm2 KβW S2 1+
γW KS4
+
βW KS2
Vm3 KβW S3 1+
αW KS1
+
βW KS3
Vm4 KγWS4 1 + KγWS4 + KβW S2 Z Vm5 KαS5 1+
αZ KS5
+
βZ KS7
Z Vm6 KβS6
1+
γZ KS8
+
βZ KS6
Z Vm7 KβS7
1+
αZ KS5
+
βZ KS7
Z Vm8 KγS8
1+
γZ KS8
+
βZ KS6
, , , , , , .
Samen met de invoering van bepaalde parameters die afhangen van de oorspronkelijk ingevoerde parameters (zie (6.5)) vertaalt dit zich in Matcont tot het systeem in Figuur 6.2.
72
6.2. Multistabiliteit in twee-staps cycli in een cascade: numeriek
Figuur 6.2: Implementatie van Model 5 in Matcont. Hierbij hebben we de parameter Vm5 ingevoerd als Vm5 = k5 · γW · WT . Dit volgt meteen uit [Wγ ] . de uitdrukkingen Vm5 = k5 [Wγ ] en γW = WT Voor de keuze van de parameters hebben we ons eerst gebaseerd op Figuur S1 uit het extra materiaal bij [8] die we terugvinden in [14]. Daaruit hebben we de gegeven waarden voor ΘW , ΘZ en KS overgenomen. Aan de hand van deze figuur kozen we ervoor om χ18 = 14 te nemen, aangezien we bij deze waarde meteen meerdere evenwichtspunten terugvinden. Doordat we echter niet weten of we met hetzelfde model werken doordat de vergelijkingen en eventuele vereenvoudigingen van het model niet weergegeven zijn in het artikel ([8]), is het niet zeker dat dit in ons model ook zo is. Voor de andere parameters kiezen we gelijkaardige waarden als in de voorgaande modellen, namelijk Vm1 = 1, Vm4 = 1, k5 = 0.006, √ WT = 500√en ZT = 500. Analoog met de voorgaande modellen kiezen we r31 = ΘW en r75 = ΘZ . Om een poging te doen om de parameters te optimaliseren, voeren we opeenvolgende continuaties uit, steeds naar verschillende parameters. Bij elk ´e´en-parameter 73
6.2. Multistabiliteit in twee-staps cycli in een cascade: numeriek
bifurcatiediagram kiezen we dan een nieuwe waarde voor de parameter door deze gelijk te stellen aan een waarde waarbij zich een maximum aantal evenwichtspunten bevinden. We voeren dit uit voor alle parameters, waarna we opnieuw beginnen met de continuatie en eventuele aanpassing van de eerste parameters totdat we parameterwaarden gekozen hebben die het maximum aantal evenwichtspunten bereiken in hun bifurcatiediagram. Uiteindelijk bekomen we op deze manier de parameterwaarden die weergegeven zijn in Tabel 6.1. Parameter ΘW r31 ΘZ r75 χ18 KS
Waarde 30 5.4772 50 7.0711 1 0.01
Dimensie 1 1 1 1 1 1
Parameter Vm1 Vm4 k5 WT ZT
Waarde 1.4 1 0.004 500 500
Dimensie nM s−1 nM s−1 s−1 nM nM
Tabel 6.1: Parameterwaarden bij Model 5. Zoals steeds beginnen we de numerieke analyse met het uitvoeren van een tijdsintegratie. Deze figuur wordt niet weergegeven in dit verslag, maar wordt wel gebruikt om equilibria te ontdekken en deze te selecteren voor de verdere continuatie. In Figuur 6.3 zullen we vooreerst de continuatie van de evenwichtspunten uitvoeren naar de parameter Vm1 . In Figuur 6.3(a) geven we een algemeen overzicht van de continuatie van equilibiria met vrije parameter Vm1 in Model 5. We merken op dat we te maken hebben met twee neutrale zadelpunten (op de figuur aangegeven als ’H’) en zes limietpunten. We overlopen de grafiek door rechts onderaan te beginnen en de equilibriumkromme te volgen. In Figuur 6.3 nummeren we alle speciale punten in de volgorde waarop ze bereikt worden en we stellen aan de hand van die nummering Tabel 6.2 op met de waarden van de variabelen en de parameter Vm1 in de speciale punten. In Figuur 6.3(d) wordt het deel van het bifurcatiediagram weergegeven met γZ ∈ [−0.001; 0.002] zodat het duidelijk is dat LP6 niet als eerste speciaal punt bereikt wordt als we rechts onderaan starten op de equilibriumkromme.
LP1 LP2 LP3 H1 LP4 LP5 H2 LP6
βW 0.003853 0.000313 0.018291 0.073895 0.073797 0.057449 0.057611 0.018291
γW 0.995700 0.999684 0.971377 0.238098 0.237084 0.822115 0.821129 0.010332
βZ 0.014074 0.014074 0.000278 0.012215 0.014074 0.014074 0.012671 0.000011
γZ 0.010150 0.975776 0.999719 0.980175 0.975776 0.010150 0.008199 0.000000
Vm1 0.668345 5.957370 0.374618 1.410262 1.412845 0.551830 0.553018 2.669412
Tabel 6.2: Waarden van de variabelen en de parameter Vm1 in de speciale punten uit Figuur 6.3. In Figuur 6.3(b) en Figuur 6.3(c) wordt een ingezoomd stuk weergegeven van Figuur 6.3(a). 74
6.2. Multistabiliteit in twee-staps cycli in een cascade: numeriek
(a)
(b)
(c)
(d)
Figuur 6.3: E´en-parameter bifurcatiediagram bij het vari¨eren van Vm1 met de andere parameterwaarden als in Tabel 6.1 (Model 5).
75
6.2. Multistabiliteit in twee-staps cycli in een cascade: numeriek
Het algemeen overzicht gaat bij deze figuren verloren, maar we kunnen op deze manier beter uitmaken met hoeveel equilibria we voor een bepaald interval parameterwaarden Vm1 te maken hebben en we zien dat we tot zes equilibria kunnen hebben. Dit vinden we bijvoorbeeld terug bij Vm1 = 0.68. We zien namelijk in Figuur 6.3(c) al vier snijpunten met de equilibriumkromme als we een rechte trekken waarvoor geldt dat Vm1 = 0.68, wat reeds vier equilibria impliceert. We mogen echter het algemeen overzicht in Figuur 6.3(a) of 6.3(b) niet vergeten, waarbij we nog een vijfde en zesde equilibrium vinden. We willen nu uiteraard ook de stabiliteit van de verschillende equilibria nagaan. Dit doen we aan de hand van de eigenwaarden in de verschillende punten en deze van de speciale punten vinden we terug in Tabel 6.3.
LP1 LP2 LP3 H1 LP4 LP5 H2 LP6
EW 1 −0.399961 −6.1504 −2.60738 −0.184398 −0.1656 −0.192722 −0.208601 −3.7673
EW 2 −0.233414 −1.15165 −0.376687 −0.0575021 −0.0576225 −0.0447309 −0.0446729 −0.533255
EW 3 −0.0803419 −0.698267 −0.0880907 −0.00332566 0 0 −0.00274853 −0.23515
EW 4 0 0 0 0.00332566 0.00333409 0.00275327 0.00274853 0
Tabel 6.3: Eigenwaarden van de speciale punten uit Figuur 6.3. Hieruit kunnen we meteen vaststellen dat de limietpunten een eigenwaarde hebben die gelijk is aan nul, wat een vereiste is om te kunnen spreken van een limietpunt. Bij de punten die aangeduid worden met ’H’ zien we telkens twee re¨ele eigenwaarden die op een teken na gelijk zijn, wat erop wijst dat dit neutrale zadelpunten zijn. Bekijken we nu ook de eigenwaarden van de andere punten in Figuur 6.3(a) en houden we er rekening mee dat de equilibria enkel stabiel zijn als ze vier negatieve eigenwaarden hebben, dan kunnen we de volledige stabiliteit van het bifurcatiediagram bespreken. Beginnen we dus rechts onderaan in Figuur 6.3(a) dan bevinden we ons op een stabiel stuk tot LP1, waarna er overgegaan wordt tot een onstabiel stuk en deze onstabiliteit terug verloren gaat in LP2. Tussen LP2 en LP3 zijn de equilibria opnieuw stabiel tot ze in LP3 terug onstabiel worden. Deze onstabiliteit blijft behouden in H1, LP4, LP5 en H2 en in LP6 is er terug overgang naar stabiliteit. De eigenwaarden van de verschillende limietpunten gaven aan wanneer de stabiliteit zou omslaan. Hebben we namelijk allemaal negatieve eigenwaarden als we kijken naar de nietnul-eigenwaarden, dan betekent dit dat er een overgang is van stabiliteit naar onstabiliteit (of omgekeerd) in dat limietpunt, terwijl we weten dat de (on)stabiliteit behouden blijft als het limietpunt ook strikt positieve eigenwaarden heeft. Hernemen we nu opnieuw het voorbeeld waar Vm1 = 0.68, dan vinden we van boven naar onder een stabiel equilibrium, twee onstabiele, opnieuw een stabiel en een onstabiel en een stabiel equilibrium. In tegenstelling tot de voorgaande modellen (zie Model 3 en 4) hebben 76
6.2. Multistabiliteit in twee-staps cycli in een cascade: numeriek
we hier dus te maken met multistabiliteit, we hebben namelijk drie stabiele equilibria voor bijvoorbeeld Vm1 = 0.68 We kunnen op een analoge manier een onderzoek uitvoeren bij de continuatie naar Vm4 en dit wordt weergegeven in Figuur 6.4.
Figuur 6.4: E´en-parameter bifurcatiediagram bij het vari¨eren van Vm4 met de andere parameterwaarden als in Tabel 6.1 (Model 5). Bij het vari¨eren van de parameter Vm4 hebben we twee neutrale zadelpunten (opnieuw aangegeven met ’H’) en drie limietpunten. Onderzoek naar de waarden van de variabelen en de parameter Vm4 in deze speciale punten geeft ons Tabel 6.4, waarbij we de nummering gekozen hebben als we links bovenaan starten in ons bifurcatiediagram (zie ook Figuur 6.4).
LP1 H1 LP2 LP3 H2
βW 0.018291 0.073727 0.073587 0.018291 0.016643
γW 0.971377 0.236368 0.234929 0.971377 0.974832
βZ 0.001128 0.012213 0.014074 0.059980 0.059855
γZ 0.998808 0.980179 0.975776 0.267454 0.265528
Vm4 3.737142 0.989626 0.987047 3.737142 3.728131
Tabel 6.4: Waarden van de variabelen en de parameter Vm1 in de speciale punten uit Figuur 6.4. Uit Figuur 6.4 zien we ook dat we maximaal vier equilibria kunnen onderscheiden bij een vrije parameter Vm4 . Voor het onderzoek naar de stabiliteit bekijken we de eigenwaarden (zie Tabel 6.5).
77
6.2. Multistabiliteit in twee-staps cycli in een cascade: numeriek
LP1 H1 LP2 LP3 H2
EW 1 −2.23861 −0.183078 −0.164095 −0.329207 −0.351877
EW 2 −0.338426 −0.0571093 −0.0571315 −0.100758 −0.101136
EW 3 −0.329207 −00330543 0 0 −0.00499283
EW 4 0 0.00330544 0.00330881 0.00497047 0.00499282
Tabel 6.5: Eigenwaarden van de speciale punten uit Figuur 6.4. Doordat we zowel in H1 als in H2 opnieuw twee eigenwaarden hebben die op een teken na gelijk zijn, zijn dit eveneens neutrale zadelpunten. Bekijken we de tekens van de eigenwaarden van de limietpunten van naderbij, dan merken we meteen op dat er zich maar ´e´en punt bevindt op deze grafiek waar er sprake is van een wisseling van stabiliteit, namelijk in LP1. Dit wordt bevestigd als we de overige eigenwaarden bekijken. Als we links bovenaan starten dan vinden we namelijk stabiele equilibria tot we in LP1 komen en we overgaan naar onstabiele equilibria. Voor de rest van de equilibriumkromme blijven de equilibria onstabiel. Een laatste onderzoek in dit model doen we aan de hand van de vrije parameter χ18 . In dit geval ziet het ´e´en-parameter bifurcatiediagram eruit als in Figuur 6.5.
Figuur 6.5: E´en-parameter bifurcatiediagram bij het vari¨eren van χ18 en de andere parameterwaarden als in Tabel 6.1 (Model 5). In Figuur 6.5 hebben we twee limietpunten LP1 en LP2, waarvan de corresponderende waarden van de variabelen en χ18 weergegeven zijn in Tabel 6.6 en de eigenwaarden in Tabel 6.7.
78
6.2. Multistabiliteit in twee-staps cycli in een cascade: numeriek
LP1 LP2
βW 0.001490 0.001490
γW 0.998443 0.998443
βZ 0.014074 0.014074
γZ 0.975776 0.010150
χ18 0.235295 2.088971
Tabel 6.6: Waarden van de variabelen en de parameter Vm1 in de speciale punten uit Figuur 6.5. We hebben dus maximaal drie equilibria en voor de eigenwaarden geldt er:
LP1 LP2
EW 1 −1.17557 −1.17557
EW 2 −0.6974 −0.235473
EW 3 −0.235473 −0.234057
EW 4 0 0
Tabel 6.7: Eigenwaarden van de speciale punten uit Figuur 6.5. Uit de eigenwaarden volgt dat de equilibria van stabiliteit veranderen in LP1 en LP2 en dat we in het gebied tussen deze twee limietpunten drie equilibria hebben, meer bepaald ´e´en onstabiel equilibrium tussen twee stabiele equilibria. Hier valt ons op dat niet voor alle parameters een interval van multistabiliteit aanwezig is: in het geval van het vari¨eren van χ18 hebben we opnieuw bistabiliteit, in tegenstelling tot het bifurcatiediagram dat opgesteld werd met een vrije parameter Vm1 . Als laatste wil ik hierbij opmerken dat ik duidelijk een verschillende grafiek bekom van deze die gegeven wordt in Figuur S1 in [14]. Dit kan te wijten zijn aan het feit dat ik misschien een ander model gebruik dan de makers van deze figuur (door het aannemen van andere vereenvoudigingen bij het opstellen ervan), maar het kan eveneens te wijten zijn aan een andere keuze van parameterwaarden.
79
Algemeen besluit Bij het schrijven van deze masterproef hebben we heel wat meer geleerd over de werking van biochemische signaalnetwerken. Uit de verschillende bestudeerde modellen kunnen we opmerken dat de organisatie en de structuur van de verschillende cycli heel belangrijk zijn voor de werking van het systeem. We vinden namelijk voor twee-staps cycli, twee-staps cycli met multisite bindingen en drie-staps cycli de mogelijkheid tot bistabiliteit, terwijl bij tweestaps cycli die in een cascade geordend zijn en waarbij een bepaalde stof zowel als substraat als enzym optreedt, zelfs multistabiliteit mogelijk is. In de eerste gevallen werken we dus op ´e´en niveau en merken we ongeacht het aantal cycli bistabiel gedrag op. Door het toevoegen van een tweede niveau, waarbij met andere woorden een cascade gevormd wordt, krijgen we multistabiliteit. Dit laatste is onafhankelijk van de aanwezigheid van feedbackloops tussen de verschillende niveaus in de cascade. Hierbij willen we wel opmerken dat er heel wat parameters meespelen in de verschillende modellen en de waarden en combinaties hiervan ook cruciaal zijn voor het aantal evenwichtspunten van elk model. Door ons wiskundig onderzoek in paragraaf 4.2 hebben we aangetoond dat een berekende grenswaarde bij bepaalde parameters zorgt voor het onderscheid tussen de mogelijkheid van bistabiliteit en het voorkomen van ´e´en enkel stabiel equilibrium. Door een invloed uit te oefenen op deze verschillende parameterwaarden in werkelijke biochemische systemen kunnen we de uitkomst van de cyclus dus in grote mate bepalen. Hierbij moeten we wel rekening houden met de vooropgestelde voorwaarden van de verzadiging van een enzym door zijn substraat en moeten de fosforylatiestappen of de defosforylatiestappen allemaal door dezelfde kinase of dezelfde fosfatase uitgevoerd worden. Door deze voorwaarden wordt de vereiste competitieve inhibitie en positieve feedback mogelijk, wat kan leiden tot bistabiliteit. Het voorkomen van deze multistabiliteit verklaart dat een bepaald biochemisch model vaak meerdere functies kan uitvoeren. De netwerken die zich in de cel bevinden zorgen er op deze manier namelijk voor dat het overbrengen van signalen via ´e´en en dezelfde baan kan resulteren in verschillende biologische processen. Hoewel E. Nikolaev in het verder onderzoek naar deze artikels een mogelijkheid zag tot het vinden van oscillaties in deze specifieke signaalnetwerken, vonden wij dit niet terug tijdens ons onderzoek als we enkel rekening houden met biologisch relevante situaties. We merken op dat B. N. Kholodenko dit wel detecteerde in [4] voor een MAPK cascade, maar hij voegde een negatieve feedback loop toe aan het model en vermeldt dat dit in combinatie met ultra-
80
sensitiviteit van de MAPK cascade kan leiden tot oscillaties. Aangezien dit in onze modellen niet aanwezig is, lijkt het niet onlogisch dat we geen biologisch relevante oscillaties ontdekken in deze masterproef. We hebben wel duidelijk het fenomeen van ultrasensitiviteit opgemerkt in onze modellen. Een kleine mutatie in ´e´en van de parameters kan een groot effect hebben op de uitkomst van een cascade. De voortplanting van het signaal doorheen de volledige baan leidt ertoe dat kleine veranderingen in de input een grote invloed hebben op het eindresultaat.
81
Verder onderzoek We merken op dat we heel wat vereenvoudigingen hebben doorgevoerd in onze modellen (gelijkstellen van de Michaelis-constanten, intermediairen in evenwicht, verwaarloosbare complexen, snelheid die bijvoorbeeld door de verzadiging niet afhankelijk is van de reactanten,...). Het spreekt voor zich dat we ook onderzoek zouden moeten voeren naar modellen waarbij deze vereenvoudigingen niet gemaakt zijn, om meer bij de werkelijkheid aan te sluiten. Ook complexere vormen van cascades, waarbij rekening gehouden wordt met de werkelijke voorafgaande signalen, lijken me noodzakelijk om het werkelijke biochemisch proces volledig te doorgronden en hier ook op in te kunnen spelen met behulp van geneesmiddelen. Zo hebben we bijvoorbeeld in paragraaf 3.4.1 duidelijk gemerkt dat de hoeveelheid van [M EK]tot en [M KP 3]tot een niet te onderschatten invloed heeft op de uitkomst van het systeem. Het aanpassen van deze hoeveelheden zou dus in de hand kunnen gewerkt worden en zorgen voor andere (en betere) effecten en uitkomsten. Het lijkt ons heel belangrijk om te blijven werken met numerieke simulaties en deze voor zoveel mogelijk verschillende parameterwaarden uit te voeren. Het is duidelijk dat er heel veel verschillende combinaties bestaan en het vaak onmogelijk lijkt om deze allemaal na te gaan. Toch kan er heel wat informatie gewonnen worden uit deze simulaties. We kunnen echter enkel een minimum aantal stabiele equilibria halen uit deze experimenten als niet alle combinaties van parameterwaarden uitgeprobeerd worden. Daarom lijkt het me van groot belang om samen met deze simulaties ook verder wiskundig onderzoek uit te voeren op meer ingewikkelde modellen, zodat we eveneens het maximale aantal stabiele equilibria kunnen berekenen. De numerieke simulaties kunnen dan als een belangrijk controlemiddel gebruikt worden. In [9] wordt er gewezen op het belang van verder onderzoek naar complexe interacties en eigenschappen van eiwitten die als een deel van grote biochemische netwerken binnen een kankercel functioneren. Het onderzoek naar deze modellen wordt dus gelinkt aan de zoektocht naar nieuwe effici¨ente geneesmiddelen tegen kanker. In dit artikel wordt nadruk gelegd op de grote hoeveelheid experimentele informatie die beschikbaar is over individuele eiwitten, maar het gebrek aan informatie over hun interacties in grote netwerken. Het doel in de kankerstudie is om de moleculaire basis van het ontstaan en het verderzetten van bepaalde kankersoorten beter te begrijpen en deze informatie te gebruiken voor behandelingen. Dit artikel handelt over kankersoorten die ontstaan door het voorkomen van pathologische RAS (cf. paragraaf 2.1.1), wat een gemuteerde vorm is van het wildtype-RAS. Een studie naar de reductie van RAS-signalering in netwerken met oncogenische RAS-mutanten waarbij de
82
reductie van wildtype-RAS beperkt blijft, wordt door de onderzoekers in deze studie beoogd. Hierbij leggen de auteurs de nadruk op het belang van een geneesmiddel dat werkt op het volledige netwerk en niet enkel op het specifieke eiwit. Ook in [10] wordt er gekeken naar het ontstaan van kanker, maar hierbij wordt eveneens de bistabiliteit betrokken. Er worden verscheidene experimenten uitgevoerd op menselijk beenmergcellen en de cellulaire differentiatie wordt onderzocht. Dit houdt een reeks beslissingen in over het lot van de cel en elk ervan wordt normaal als irreversibel beschouwd in hogere organismen. Nu blijkt uit deze onderzoeken dat er vaak toch een deel plasticiteit teruggevonden wordt en er dus invloed kan uitgeoefend worden op deze differentiatie. Dit besluit werd bekomen doordat de aandacht verschoven werd naar een kwantitatieve analyse van de interacties tussen de verschillende individuele componenten in een biologisch systeem en hoe deze interacties evolueren in de tijd. De studie van hoe menselijke stamcellen osteogenische (staat in voor de aanmaak van beendercellen) en myogenische (staat in voor de aanmaak van spiercellen) differentiatie ondergaan, toont aan dat een bistabiel wisselmechanisme samen met feedback regulatie dominerend kan zijn voor deze celdifferentiatie. Door het experimenteren met behandelde cellen werden verschillende dosissen van een bepaald eiwit (in dit onderzoek BMP2) toegevoegd en werd er opgemerkt dat er een groot verschil te zien was als de dosis BMP2 een bepaalde drempelwaarde overschreed (hysteresis). Dit hield meer bepaald de overgang in van een ongedifferentieerde toestand van de cel naar een gedifferentieerde toestand. Dit betekent dat er een keuze is tussen myogenische versus osteogenische differentiatie onder verschillende condities. Er werd opgemerkt dat MAPK-inhibitie zorgt voor de osteogenische differentiatie en MAPK-activatie noodzakelijk is voor de myogenische activatie. Dit wil meer bepaald zeggen dat het myogenisch vermogen onderdrukt wordt bij de inhibitie van MAPK, maar er wel osteogenische stimulatie is. Bij een blijvende MAPK-inhibitie zorgt dit voor een opstapeling van beendercellen en een tekort aan spiercellen. Ook brengt dit onderzoek het goede nieuws dat de myogenische activiteit herwonnen kan worden na osteogenese, wat van levensbelang is voor het vinden van geneesmiddelen tegen kanker. Het onderzoek naar deze interessante materie is nog volop aan de gang is (de artikels [9] en [10] dateren beide uit 2009) en de onderzoekers hebben een heel specifieke toepassing in gedachten bij dit onderwerp, namelijk het vinden van effici¨ente geneesmiddelen om een zware en veel voorkomende ziekte als kanker te genezen. Een samenwerking van verschillende onderzoeksgebieden zoals wiskunde, biochemie, farmaceutica,... is hierbij vereist en het resultaat kan van groot belang zijn binnen de gehele maatschappij.
83
Glossarium Actieve site: plaats aan het oppervlak van het enzym waarop het substraat bindt en een chemische reactie ondergaat. Activeringsenergie: de energie die een systeem nodig heeft om een reactie te starten (zie ook paragraaf 1.1). ADP: adenosinedifosfaat, ontstaat uit ATP door het wegnemen van een fosfaatgroep. ATP: adenosinetrifosfaat, de belangrijkste energiedragende verbinding in cellen van levende organismen. Bistabiliteit: een systeem met twee stabiele toestanden waarbij een grondige storing vereist is om van de ene evenwichtstoestand naar de andere over te gaan. Competitieve inhibitie: het fenomeen waarbij een inhibitor en een substraat streven naar een binding op de actieve site van een enzym en elkaar zo tegenwerken. Covalente binding = atoombinding: een binding tussen atomen waarbij deze atomen ´e´en of meerdere gemeenschappelijke elektronenparen hebben. Defosforylatie: het biochemisch proces van het verwijderen van een fosfaatgroep van een enzym of een andere organische molecule (zie ook paragraaf 1.3). Distributief mechanisme: het modificatieproces waarbij een enzym zijn substraat na ´e´en modificatie terug loslaat en er dus een nieuwe binding moet optreden vooraleer er een nieuwe modificatie kan uitgevoerd worden (niet-processief). Enzym: een eiwit dat zich als biologische katalysator gedraagt. Feedback mechanisme: een ketting van gebeurtenissen waarin een bepaalde stap een voorgaande stap controleert of be¨ınvloedt (positief of negatief). Fosfaatgroep: het negatief geladen ion P O43− bestaande uit een fosforatoom en vier zuurstofatomen. Fosfatase: een enzym dat de binding van fosfaatgroepen met organische moleculen verbreekt tijdens de defosforylatie. Fosforylatie: het biochemisch proces van het binden van een enzym of een andere organische molecule met een fosfaatgroep (zie ook paragraaf 1.3).
84
Hydrolyse: de splitsing van een chemische verbinding die veroorzaakt wordt door de opname van water. Induced fit: een verandering in de samenstelling van het enzym door substraatbinding. Inhibitor: een molecule die een enzym verhindert te werken op de gewone manier. Katalysator: een stof die de snelheid van een chemische reactie verhoogt, waarbij deze katalysator onveranderd uit de chemische reactie gerecupereerd wordt. Kinase: een enzym dat fosfaatgroepen van ATP kan overdragen op organische moleculen tijdens de fosforylatie. MAPK: mitogen-activated protein kinase (zie ook paragraaf 2.1). Michaelis-Menten vergelijking: de vergelijking die de wiskundige analyse van de snelheden van reacties beschrijft van een substraat S met een enzym E tot een enzym-substraatcomplex ES dat reageert tot een product P (E + S ES → E + P ). Standaard wordt dan aangenomen dat [S] >> [E] en dat de opbouw en de afbraak van de concentratie van het intermediair ES in evenwicht is. De Michaelis-Menten vergelijking van deze simpele reac[S] , waarbij v de snelheid van de tievergelijking ziet er dan uit als volgt: v = vmax [S] + Km reactie is, vmax de maximale snelheid en Km de MM-constante. MM: Michaelis-Menten. Multistabiliteit: een systeem met meer dan twee stabiele toestanden. Negatieve feedback regulatie: een reactiepad waarbij een product zijn eigen aanmaak tegenwerkt. Positieve feedback regulatie: een reactiepad waarbij een product zijn eigen aanmaak bevordert. Processief mechanisme: het enzym blijft gebonden aan zijn substraat zodat er herhalingen van de katalytische gebeurtenis kunnen optreden totdat het eindproduct gevormd is. Substraat: een molecule dat in een enzymreactie gebonden wordt met het enzym en omgezet wordt tot een product.
85
Bijlage A
A.1 > > > > > > >
Appendix A
restart: b:=(-g+sqrt(g^2-4*t*g^2+4*t*g))/(2*t): chi:=(b*g+b^2*t+g*K)/(r*b*(b+g+K)): A:=diff(chi,g): a:=solve(A=0,g): m:=allvalues(a): g1:=simplify(m[1]); / (1 + K) K \1/2 |-----------------| \4 K t - 1 + t - K/ / / g1 := 1/2 + ---------------------- + |- | 2 \ \ / (1 + K) K \1/2 2 2 4 |-----------------| t K - 2 K \4 K t - 1 + t - K/ / (1 + K) K \1/2 2 / (1 + K) K \1/2 - |-----------------| K - 2 K - |-----------------| t \4 K t - 1 + t - K/ \4 K t - 1 + t - K/ / (1 + K) K \1/2\ / // (1 + K) K \1/2 + |-----------------| | / ||-----------------| \4 K t - 1 + t - K/ / / \\4 K t - 1 + t - K/ \\1/2 (4 K t - 1 + t - K)|| /2 //
> g2:=simplify(m[2]);
86
A.1. Appendix A
/ (1 + K) K \1/2 |-----------------| \4 K t - 1 + t - K/ / / g2 := 1/2 + ---------------------- - |- | 2 \ \ / (1 + K) K \1/2 2 2 4 |-----------------| t K - 2 K \4 K t - 1 + t - K/ / (1 + K) K \1/2 2 / (1 + K) K \1/2 - |-----------------| K - 2 K - |-----------------| t \4 K t - 1 + t - K/ \4 K t - 1 + t - K/ / (1 + K) K \1/2\ / // (1 + K) K \1/2 + |-----------------| | / ||-----------------| \4 K t - 1 + t - K/ / / \\4 K t - 1 + t - K/ \\1/2 (4 K t - 1 + t - K)|| /2 // > g3:=simplify(m[3]); / (1 + K) K \1/2 |-----------------| \4 K t - 1 + t - K/ / / g3 := 1/2 - ---------------------- + |- | 2 \ \ / (1 + K) K \1/2 2 2 4 |-----------------| t K + 2 K \4 K t - 1 + t - K/ / (1 + K) K \1/2 2 / (1 + K) K \1/2 - |-----------------| K + 2 K - |-----------------| t \4 K t - 1 + t - K/ \4 K t - 1 + t - K/ / (1 + K) K \1/2\ / // (1 + K) K \1/2 + |-----------------| | / ||-----------------| \4 K t - 1 + t - K/ / / \\4 K t - 1 + t - K/ \\1/2 (4 K t - 1 + t - K)|| /2 // 87
A.1. Appendix A
> g4:=simplify(m[4]); / (1 + K) K \1/2 |-----------------| \4 K t - 1 + t - K/ / / g4 := 1/2 - ---------------------- - |- | 2 \ \ / (1 + K) K \1/2 2 2 4 |-----------------| t K + 2 K \4 K t - 1 + t - K/ / (1 + K) K \1/2 2 / (1 + K) K \1/2 - |-----------------| K + 2 K - |-----------------| t \4 K t - 1 + t - K/ \4 K t - 1 + t - K/ / (1 + K) K \1/2\ / // (1 + K) K \1/2 + |-----------------| | / ||-----------------| \4 K t - 1 + t - K/ / / \\4 K t - 1 + t - K/ \\1/2 (4 K t - 1 + t - K)|| /2 //
We hebben dus q 2 2 2 S −BKS −2KS −BΘ+B − 4BΘKS −2K B(4KS Θ−1+Θ−KS )
1 B + + , 2 2 2 q 2 2 2 S −BKS −2KS −BΘ+B − 4BΘKS −2K 1 B B(4KS Θ−1+Θ−KS ) γ2 = + − , 2 2 2 q γ1 =
2
2
2
S −BKS +2KS −BΘ+B − 4BΘKS +2K 1 B B(4KS Θ−1+Θ−KS ) γ3 = − + , 2 2 2 q 2 2 2 S −BKS +2KS −BΘ+B − 4BΘKS +2K 1 B B(4KS Θ−1+Θ−KS ) γ4 = − − , 2 2 2
met s B=
(1 + KS )KS . (4KS Θ − 1 + Θ − KS )
88
A.2. Appendix B
Na deling van teller en noemer van het quoti¨ent onder de wortel van de laatste term door B, bekomen we dan uiteindelijk s ! 1 (1 + KS )KS γ1 = 1+ + A1 , 2 Θ − 1 + KS (4Θ − 1) s ! (1 + KS )KS 1 γ2 = 1+ − A1 , 2 Θ − 1 + KS (4Θ − 1) s ! 1 (1 + KS )KS γ3 = 1− + A2 , 2 Θ − 1 + KS (4Θ − 1) s ! 1 (1 + KS )KS 1− − A2 , γ4 = 2 Θ − 1 + KS (4Θ − 1) waarbij we p √ √ KS 2 (−4Θ + 1) + Θ − 1 + 2 KS KS + 1 Θ − 1 + KS (4Θ − 1) A1 = , Θ − 1 + KS (4Θ − 1) s p √ √ KS 2 (−4Θ + 1) + Θ − 1 − 2 KS KS + 1 Θ − 1 + KS (4Θ − 1) A2 = Θ − 1 + KS (4Θ − 1) s
gesteld hebben.
A.2
Appendix B
> w:=K^2*(-4*t+1)+t-1+2*sqrt(K)*sqrt(K+1)*sqrt(t-1+K*(4*t-1)): > s:=solve(w=0,t): > o1:=factor(s[1]); 2 (1 + K) o1 := ---------2 (2 K - 1) > o2:=factor(s[2]); 2 (1 + K) o2 := ---------2 (2 K + 1)
89
A.3. Appendix C
A.3
Appendix C
> t:=(K+1)^2/(2*K-1)^2: > simplify(K^2*(-4*t+1)+t-1); 6 K (K + 1) - ----------2 K - 1 > t:=(K+1)^2/(2*K+1)^2: > simplify(K^2*(-4*t+1)+t-1); 2 K (K + 1) - ----------2 K + 1
A.4
Appendix D
> t:=(K+1)^2/(2*K-1)^2: > simplify((1+K)*K/(t-1+K*(4*t-1))); 2 (2 K - 1) ---------9 > t:=(K+1)^2/(2*K+1)^2: > simplify(((1+K)*K)/(t-1+K*(4*t-1))); 2 (2 K + 1)
A.5
Appendix E
> t:=(K+1)^2/(2*K-1)^2: > assume(K>1/2); > g:=simplify(1/2+1/2*sqrt((1+K)*K/(t-1+K*(4*t-1)))); K~ g := 1/3 + ---3 > b:=simplify((-g+sqrt(g^2-4*t*g^2+4*t*g))/(2*t));
90
A.5. Appendix E
(K~ - 2) (2 K~ - 1) b := - ------------------3 (K~ + 1) > a:=simplify(b^2*t/g); 2 (K~ - 2) a := ---------3 (K~ + 1) > chi:=simplify((b*g+b^2*t+g*K)/(r*b*(b+g+K))); 2 2 2 (K~ + 2 K~ + 1) (K~ + 1) chi := - ----------------------------------------2 (2 K~ + 10 K~ - 1) (2 K~ - 1) (K~ - 2) r > assume(K<1/2); > g:=simplify(1/2-1/2*sqrt((1+K)*K/(t-1+K*(4*t-1)))); K~ g := 1/3 + ---3 > b:=simplify((-g+sqrt(g^2-4*t*g^2+4*t*g))/(2*t)); (2 K~ + 3 signum(K~ + 1) - 1) (2 K~ - 1) b := -1/6 ---------------------------------------K~ + 1 > b:=(1-2*K)/3: > a:=simplify(b^2*t/g); K~ a := 1/3 + ---3 > chi:=simplify((b*g+b^2*t+g*K)/(r*b*(b+g+K))); K~ + 1 chi := - -----------r (2 K~ - 1)
91
A.6. Appendix F
A.6 A.6.1 > > > >
Appendix F Oplossingen van de vergelijking
∂β =0 ∂γ
restart: assume(t>0); b:=(-g+sqrt(g^2-4*t*g^2+4*t*g))/(2*t): B:=diff(b,g); 2 g - 8 t~ g + 4 t~ -1 + ---------------------------2 2 1/2 2 (g - 4 t~ g + 4 t~ g) B := --------------------------------2 t~
> g12:=solve(B=0,g); 1/2 1/2 2 t~ + t~ -2 t~ + t~ g12 := ------------, - -------------1 + 4 t~ -1 + 4 t~
A.6.2
Bestaansvoorwaarden nagaan
We stellen eerst g gelijk aan de eerste oplossing: > > g:=(2*t+t^(1/2))/(-1+4*t): We gaan BV 1 na voor g1: > > r1:=simplify(g^2-4*t*g^2+4*t*g): > is(r1>0); true Nu we weten dat voldaan is aan BV 1, gaan we ook BV 2 na voor g1: > > l1:=simplify(2*g-8*t*g+4*t): > is(l1>0); false We weten nu dat deze term niet overal positief is, maar we moeten nog
92
A.6. Appendix F
controleren of er eventueel gebieden zijn waarin deze term wel positief is en deze oplossing g1 dus geldig is in een deelgebied: > is(l1<0); true Er is geen deelgebied waar deze oplossing bestaat, want de term is overal negatief en dus onmogelijk gelijk aan de wortel in het linkerlid. > We doen hetzelfde voor de tweede oplossing g2: > g:=-(-2*t+t^(1/2))/(-1+4*t); 1/2 -2 t~ + t~ g := - -------------1 + 4 t~ > r2:=simplify(g^2-4*t*g^2+4*t*g): > is(r2>0); true b is dus overal gedefinieerd als ze beschouwd wordt in deze tweede oplossing g2. Nu rest ons nog BV 2: > l2:=simplify(2*g-8*t*g+4*t): > is(l2>0); true
A.6.3
Andere veranderlijken in βmax
> g:=sqrt(t)/(2*sqrt(t)+1): > b:=simplify(b); 1 b := ----------1/2 2 t~ + 1 > a=simplify(1-b-g);
93
A.6. Appendix F
1/2 t~ a = ----------1/2 2 t~ + 1 > chi:=simplify((b*g+b^2*t+g*K)/(r*b*(b+g+K))); 1/2 1/2 t~ + t~ + 2 t~ K + t~ K chi := ----------------------------1/2 1/2 (1 + t~ + 2 t~ K + K) r
A.6.4 > > > > > > >
Berekende waarde voor β: minimaal of maximaal?
restart: assume(t>0); b:=(-g+sqrt(g^2-4*t*g^2+4*t*g))/(2*t): B:=diff(b,g): C:=diff(B,g): g:=sqrt(t)/(2*sqrt(t)+1): simplify(C); 2 - ----1/2 t~
94
Bibliografie [1] G. Dejaegher, cursus Beginselen van de celbiologie en genetica, hoofdstuk 3: De energetische en chemische huishouding van de cel, Universiteit Gent (2007) [2] S.P. Ellner, J. Guckenheimer, Dynamic Models in Biology (Princeton University Press, 2006). [3] W. Govaerts, cursus Toegepaste Wiskundige Evolutiemodellen, universiteit Gent (2009) [4] B. N. Kholodenko, Negative feedback and ultrasensitivity can bring about oscillations in the mitogen-activated protein kinase cascades, FEBS Journal 267 (2000), 1583-1588. [5] Yu. A. Kuznetsov, Elements of Applied Bifurcation Theory, Applied Mathematical Sciences 112, Third Edition (Springer, 2004). [6] N.I. Markevich, J.B. Hoek, B.N. Kholodenko, Signaling switches and bistability arising from multisite phosphorylation in protein kinase cascades, The Journal of Cell Biology 164, 3 (2004), 353-359. [7] C.K. Mathews, K.E. van Holde, K.G. Ahern, Biochemistry (Third edition) (Addison Wesley Longman, 2000). [8] F. Ortega, J.L. Garc´es, F. Mas, B.N. Kholodenko, M. Cascante, Bistability from double phosphorylation in signal transduction: Kinetic and structural requirements, FEBS Journal 273 (2006), 3915-3926. [9] E. C. Stites, K. S. Ravichandran, A Systems Perspective of Ras Signaling in Cancer, Clin Cancer Res 15, 5 (2009), 1510-1513. [10] L. Wang, B. L. Walker, S. Iannaccone, D. Bhatt, P. J. Kennedy, W. T. Tse, Bistable switches control memory and plasticity in cellular differentiation, PNAS 106, 16 (2009), 6638-6643. [11] http://163.16.28.248/bio/activelearner/06/ch6c1.html [12] http://jcb.rupress.org/content/164/3/353 [13] http://jjj.biochem.sun.ac.za/database/ [14] http://onlinelibrary.wiley.com/doi/10.1111/j.1742-4658.2006.05394.x/ suppinfo 95
Bibliografie
[15] http://sourceforge.net/projects/matcont/ [16] http://www.blackwell-synergy.com [17] http://www.ebi.ac.uk [18] http://www.medicinenet.com [19] http://www.mun.ca/biology/scarr/Induced-Fit_Model.html [20] http://www.nature.com/ [21] http://www.scq.ubc.ca/ [22] http://www.wikipedia.org/
96