Blog Michael-Projekt

[0] C.M.: Projektplan

Es sollen statistische Fluktuationen in bestimmten biochemischen Reaktionsnetzwerken erforscht werden, insbesondere im Grenzfall kleiner Molekülzahlen, wo starke Abweichungen von analytischen Näherungen zu erwarten sind.

Die Dynamik der Reaktionen wird durch exakte Monte-Carlo-Simulation der Reaktionsprozesse untersucht. Die Resultate werden mit analytischen Näherungen (insbes. mit Lösungen entsprechender stochatischer Differentialgleichungen) verglichen.

Besonderes Ziel: Suche nach einem einfachen Netzwerk, welches statistisch schwankende Molekülzahlen $n$produziert, die gemäß einem Potenzgesetz $P(n)\propto n^{-\gamma}$ mit gebrochenzahligem Exponenten $\gamma$ verteilt sind.

Zunächst soll aber die korrekte Funktion des Simulationsprogrammes anhand verschiedener relativ einfacher Reaktions-Netzwerke getestet werden:

  • Ideale Produktion (Unerschöpflicher Rohstoff)

- (A=const)k>(B)
- 1 Parameter: k
- [B](t): Vergleich mit analyt. Erwart.
- Statistik der Zeitintervalle zwischen Ereignissen

  • Nichtideale Produktion (Begrenzter Rohstoff)

- (A)k>(B)
- 2 Parameter: [A], k
- [B](t): Vergleich mit analyt. Erwart.

  • Zyklus aus 2 Poisson-Prozessen

- (A)k1>(B)k2>(A)
- 3 Parameter: A_sum=[A]+[B], k1, k2
- P(A),C_AA(tau): Vergleich mit analyt. Erwart.

  • Enzymatische Umwandlung

- (A)+(e) <k1,k2> (eA) k3> (e)+(B)
- 5 Parameter: A_sum=[A]+[B][A]+[eA]+[B], [e], k1..k3
- Vergleich mit Michaelis Menten

  • Goldbeter-Koshland-Zyklus

- (X)+(a) <k1,k2> (aX) k3> (a)+(Y)
- (Y)+(d) <k4,k5> (dY) k6> (d)+(X)
- 9 Parameter: A_sum=[X]+[aX]+[Y]+[dY], [a], [d], k1..k6
- P(A),C_AA(tau): Numerisch erforschen, Parameter-Studien


[1] M.S.: Problem im Modell "Ideale Produktion"

Ich poste die Mail der Vollständigkeit halber auch einfach mal hier:

Ich sitze jetzt schon einige Zeit hier am Programm und ich glaube ich hab mich etwas verrannt:

in „idealproduction“ sollte doch alle Schwankung aus der Poissonverteilung der Zeit kommen …

wobei das „Poisson“ – Lambda = 1/k ist und x ein zufälliger integer (

*tmDelay_=PoissonRnd(&RNDINI,1.0/totRate);

)somit sollte nach Definition die Poissonverteilung aussehen wie :

(1)
\begin{align} $$ P_{\lambda}(X = x) = \frac{\lambda^x}{x!}exp(-\lambda) $$ $$ dt = P_{\frac{1}{k_{tot}}}(X = x) = (\frac{1}{k_{tot}})^{x} \frac{exp(-\frac{1}{k_{tot}})}{x!} $$ \end{align}

ich habe jetzt verschiedene k-werte eingesetzt und bekomme immer ein ähnliches bild … das passt aber nicht (siehe angehängte dateien)


[2] C.M.: Zu [1]

  • Die ideale Produktion ist ein sog. Poisson-Prozeß. Dies bedeutet, daß ein Ereignis (die Produktion eines neuen B-Moleküls) mit einer bestimmten Rate (=Wahrscheinlichkeit pro Zeiteinheit) k passiert. Die Zeitpunkte, an denen die Ereignisse tatsächlich passieren, seien mit $t_n$ bezeichnet.
  • Die Zeitdauern $\tau_n=t_{n+1}-t_n$ zwischen zwei aufeinanderfolgenden Ereignissen sollten in einem Poisson-Prozeß exponential-verteilt sein: $P(\tau)\propto e^{-\tau/\tau_c}$ mit $\tau_c=1/k$. Wenn ich mir Deine timedif_hist-Kurven anschaue, scheint dies ja auch zu funktionieren: Lauter e-Funktionen mit verschiedenen Abfallkonstanten. Für Deinen Bericht würde ich übrigens mehrere solche Kurven für verschiedene k in einfach-logarithmischer Auftragung in den gleichen Graphen zeichnen (—> Bündel abfallendender Geraden mit versch. Steigungen).
  • Die von Dir gesuchte Poisson-Verteilung kommt etwas anders in Spiel: Wenn man ein festes endliches Zeitintervall T betrachtet, dann geschehen in T im Mittel $\overline{n}=kT$ Produktions-Ereignisse, manchmal aber auch etwas mehr oder weniger viele. Die tatsächliche Anzahl sei $n$. Die Verteilung $P(n)$ für gegebenes $\overline{n}$ ist eine Poissonverteilung: $P_{\overline{n}}(n)=e^{-\overline{n}}\left(\frac{\overline{n}^n}{n!}\right)$. Um zu testen, ob das korrekt herauskommt, müßtest Du erst ein paar neue Routinen programmieren.

[3] M.S.: Laborbuch

Ich stelle hier auch gleich mein "Laborbuch" (oder wie man sowas in der Theoretischen Physik nennt) online (siehe files: Enzymatische Stoffumwandlung….pdf).


[4] C.M.: Zweistufige Produktion

Bei einem zweistufigen Produktionsprozeß (mit unerschöpflichem Substrat)

(2)
\begin{align} S+E \rightarrow ES \rightarrow E+P \end{align}

würde man für die Zeitdifferenzen zwischen zwei aufeinanderfolgenden Produktionsereignissen von P keine Exponential-Verteilung, sondern die Faltung zweier Exp.-Vert. erwarten. Dies schien sich zunächst numerisch nicht zu bestätigen. Dabei liegt der Fehler darin, daß im simulierten Molekülzahl-Verlauf DRecXY state[P_sub] die aufeinanderfolgenden Wertpaare solchen Zeitpunkten entsprechen, an denen IRGENDEIN Prozeß im System passiert ist. Hier interssieren aber nur die Zeitpunkte, in denen ein neues P-Molekül entstanden ist. Das Gillespie-Programm sollte entsprechend modiifiziert werden.


[5] C.M.: Zu [4]

  • Das Gillespie-Programm (siehe files: gillespie.h) wurde nun so modifiziert, daß für jeden Stoff die Molekülzahl nur dann in seinen state-Verlauf eingetragen wird, wenn sie sich wirklich geändert hat. Zeitdifferenzen zwischen bestimmten Molekülzahl-Änderungen lassen sich so leicht auswerten.
  • Das Programm wurde für die in [4] beschriebene Zweistufen-Produktion getestet (siehe files: prodChain.h). In der Simulation war [S]=1 und [E]=1. Die Raten für beide Stufen waren dabei k1=k2=1.0/sec:

twoStepProd1.gif

  • Interessanterweise bleibt die qualititive Form von P(dt) auch dann erhalten, wenn man z.B. zwei E-Moleküle im System hat. In diesem Fall arbeiten die individuellen Enzym-Moleküle unabhängig voneinander und die Produktions-Ereignisse von P-Molekülen erfolgen aysnchron. Trotzdem kommen offensichtlich sehr kleine dt nicht vor.

[6] C.M.: Einfacher Poisson-Prozeß als stochastische DGL

Betrachte Produktion eines Stoffes X aus dem Nichts mit Rate k:

(3)
\begin{align} \Phi \rightarrow k \rightarrow X \end{align}

In einer Kontinuumsnäherung würde dies auf die DGL $\dot{x}=\overline{k}$ führen. Zur näherungsweisen Berücksichtigung der stochastischen Effekte kann man neben der mittleren Rate $\overline{k}$ auch den fluktuierenden Teil $\Delta k(t)$ berücksichtigen:

(4)
\begin{align} \dot{x}=\overline{k}+\Delta k(t). \end{align}

Um den Theorie-Apparat der stochastischen Differentialgleichungen benutzen zu können, sollen die Ratenfluktuationen durch Gaussches weißes Rauschen angenähert werden, d.h.

(5)
\begin{align} \left< \Delta k(t) \Delta k(t^{\prime})\right>=\Gamma \delta(t-t^{\prime}). \end{align}

Wie groß muß die Fluktuationsstärke $\Gamma$ gewählt werden, damit das Verhalten eines Poisson-Prozesses möglichst gut reproduziert wird ?

Wir betrachten zur Beantwortung dieser Frage die Anzahl n(T) der X-Moleküle, die während eines Zeitintervalls der Länge T produziert werden:

(6)
\begin{align} n(T)=\int_0^T \dot{x}(t) dt = \overline{k}T +\int_0^T \Delta k(t) dt = \overline{n}+\Delta n. \end{align}

Für einen Poissonprozeß muß im Ensemble-Mittel gelten:

(7)
\begin{align} \left< (\Delta n)^2 \right> = \overline{n}, \end{align}

oder

(8)
\begin{align} \left< \left( \int_0^T \Delta k(t) dt\right)^2 \right> = \overline{k}T. \end{align}

Schreibt man das quadrierte Integral als Produkt zweier Integrale mit verschiedenen Zeitvariablen und vertauscht Zeit-Integration mit Ensemble-Mittelung, so kann die linke Seite obiger Gleichung in $\Gamma T$ umgeformt werden, so daß $\Gamma = \overline{k}$ folgt. Die Ratenfluktuationen müssen für einen Poissonprozeß also gerade die Bedingung

(9)
\begin{align} \left< \Delta k(t) \Delta k(t^{\prime})\right>=\overline{k} \delta(t-t^{\prime}) \end{align}

erfüllen. Wir dividieren diese Bedingungsgleichung auf beiden Seiten durch $\overline{k}$:

(10)
\begin{align} \left< \frac{\Delta k(t)}{\sqrt{\overline{k}}} \frac{\Delta k(t^{\prime})}{\sqrt{\overline{k}}}\right> = \delta(t-t^{\prime}) \end{align}

und definieren eine neue Zufallsgröße

(11)
\begin{align} \zeta(t)=\frac{\Delta k(t)}{\sqrt{\overline{k}}}. \end{align}

Für diese Größe hat der Korrelator die Standardform, wie sie in der Theorie stochastischer DGLen üblich ist:

(12)
\begin{align} \left< \zeta(t) \zeta(t^{\prime}) \right> = \delta(t-t^{\prime}). \end{align}

Resultat: Um einen einfachen Poissonprozeß näherungsweise als stochastische DGL mit weißem Rauschen zu beschreiben, müssen wir ansetzen:

(13)
\begin{align} \dot{x}=\overline{k} + \sqrt{\overline{k}} \; \zeta(t). \end{align}

[7] C.M.: 2-Poisson-Zyklus als stochastische DGL

Wir betrachten nun folgenden Zyklus:

(14)
\begin{align} X_0 \rightarrow a \rightarrow X \rightarrow d \rightarrow X_0. \end{align}

Dabei soll die Gesamtkonzentration der beiden Molekülformen X_0 und X erhalten sein:

(15)
\begin{equation} x_0+x=x_g. \end{equation}

Die Aktivierungs- bzw. Deaktivierungs-Reaktionen werden in diesem einfachen Modell jeweils durch Poisson-Prozesse mit mittleren Ratenkonstanten $\overline{a}$ bzw. $\overline{d}$ beschrieben.

Ohne Fluktuationen würde die Dynamik der X-Molekülzahl folgender DGL gehorchen:

(16)
\begin{align} \dot{x} = a (x_g-x) - d x. \end{align}

Im stationären Fall $\dot{x}=0$ ergibt sich die Gleichgewichtskonzentration

(17)
\begin{align} \overline{x}=\frac{x_g}{1+(d/a)}. \end{align}

In diesem stationären Zustand laufen ständig Aktivierungs- und Deaktivierungs-Prozesse mit gleicher mittlerer Rate $R_{eq}$ ab. Diese mittlere Rate ergibt sich durch Einsetzen der Gleichgewichtskonzentration $x=\overline{x}$ in den Term $a (x_g-x)$ oder $d x$. Das Ergebnis lautet in jedem Fall:

(18)
\begin{align} R_{eq}=\frac{a d}{a+d}x_g. \end{align}

Im nächsten Schritt sollen wieder die zeitlichen Fluktuationen der Raten berücksichtigt werden. Im Gegensatz zu [6] hat man es nun mit zwei stochastischen Prozessen zu tun, die gleichzeitig ablaufen. Eigentlich sind diese beiden Prozesse keine einfachen Poissonprozesse, weil ja z.B. nicht einmal die mittlere Rate $a (x_g-x)$ konstant ist, sondern vielmehr von der Molekülzahl x(t) abhängt. Außerdem sind der Aktivierungs- und Deaktivierungsprozeß nicht statistisch unabhängig. Wir werden diese Komplikationen aber näherungsweise vernachlässigen.

Sieht man den Aktivierungsprozeß als einfachen Poissonprozeß und verwendet die Überlegungen aus Post [6], so sollte für die Fluktuation der Aktivierungsrate um ihren Mittelwert gelten:

(19)
\begin{align} \Delta k_{act}(t) = \sqrt{R_{eq}}\; \zeta_{act}(t). \end{align}

Analog folgt für die Deaktivierung:

(20)
\begin{align} \Delta k_{dea}(t) = \sqrt{R_{eq}}\; \zeta_{dea}(t). \end{align}

Die Rauschfunktionen $\zeta_i(t)$ sind wieder $\delta$-korreliert mit Vorfaktor 1.

Somit lautet die stochastische DGL des Systems:

(21)
\begin{align} \dot{x} = [a (x_g-x) - d x] + \Delta k_{act}(t) + \Delta k_{dea}(t), \end{align}

oder

(22)
\begin{align} \dot{x} = [a (x_g-x) - d x] + \sqrt{R_{eq}}\;\zeta_{act}(t) + \sqrt{R_{eq}}\; \zeta_{dea}(t). \end{align}

Ignoriert man Kreuzkorrelationen zwischen den Prozessen und unterstellt $\alpha$-Stabilität, so können die beiden Rauschterme zu einem einzigen zusammengefaßt werden, wobei sich die Varianzen addieren:

(23)
\begin{align} \dot{x} = [a (x_g-x) - d x] + \sqrt{R_{eq}+R_{eq}}\;\zeta(t), \end{align}

oder

(24)
\begin{align} \dot{x} = [a (x_g-x) - d x] + [\sqrt{2 R_{eq}}]\;\zeta(t). \end{align}

Diese Gleichung hat die Standardform einer stochastischen DGL:

(25)
\begin{align} \dot{x} = f(x) + g(x)\;\zeta(t). \end{align}

Es folgt

(26)
\begin{equation} A(x)=f(x)=a (x_g-x) - d x = ax_g -(a+d)x \end{equation}

und

(27)
\begin{align} B(x)=\frac{1}{2}g^2(x)=R_{eq}. \end{align}

Somit:

(28)
\begin{align} \frac{A(s)}{B(s)}=\frac{ax_g}{R_{eq}}-\frac{(a+d)}{R_{eq}} s. \end{align}
(29)
\begin{align} \int_0^x \frac{A(s)}{B(s)} ds = \frac{ax_g}{R_{eq}} x - \frac{(a+d)}{R_{eq}} \frac{x^2}{2}. \end{align}
(30)
\begin{align} \int_0^x \frac{A(s)}{B(s)} ds = \frac{a+d}{d}x - \frac{(a+d)^2}{2adx_g}x^2. \end{align}

Für große x dominiert der quadratische Term und man erhält als stationäre Lösung der Fokker-Planck-Gleichung eine Gaußverteilung:

(31)
\begin{align} P(x) \propto \exp \left[ -\frac{1}{2} \frac{(a+d)^2}{adx_g} x^2 \right] \end{align}

mit Varianz

(32)
\begin{align} \sigma_x^2=\frac{adx_g}{(a+d)^2}, \end{align}

was sich im Speziallfall $a=b=k$ zu $\sigma_x^2=\frac{k^2 x_g}{(2k)^2}=\frac{x_g}{4}$ vereinfacht. Die Standardabweichung wäre dann also $\sigma_x = \sqrt{\frac{x_g}{4}}$. Bei 10000 Substrat-Molekülen ergibt sich z.B. eine Gaußverteilung der Breite $\sigma_x = 50$ um den Mittelwert 5000.


[8] C.M.: Zusatz zu [7]

Betrachte nochmals das Integral

(33)
\begin{align} \int_0^x \frac{A(s)}{B(s)} ds = \left[ \frac{a+d}{d}\right]x + \left[-\frac{(a+d)^2}{2adx_g}\right]x^2 = \beta x + \alpha x^2. \end{align}

Mittels quadratischer Ergänzung folgt

(34)
\begin{align} \alpha x^2 + \beta x = \alpha (x+\frac{\beta}{2\alpha})^2-\frac{\beta^2}{4\alpha}. \end{align}

Teilweises Einsetzen von $\alpha$ und $\beta$ liefert:

(35)
\begin{align} \int_0^x \frac{A(s)}{B(s)} ds = \frac{ax_g}{2d} + \alpha\left[ x - (\frac{ax_g}{a+d}) \right]^2 = \mbox{const} + \alpha\left[ x - \overline{x} \right]^2. \end{align}

Für den Faktor $\alpha$ vor der eckigen Klammer kann man auch $\frac{-1}{2\sigma_x^2}$ schreiben. Somit:

(36)
\begin{align} \int_0^x \frac{A(s)}{B(s)} ds = \mbox{const} + \frac{-(x - \overline{x})^2}{2\sigma_x^2}. \end{align}

Die stationäre PDF hat damit die erwartete Form:

(37)
\begin{align} P(x) \propto \exp \left[ - \frac{(x - \overline{x})^2}{2\sigma_x^2} \right]. \end{align}

[9] C.M.: Autokorrelation von x(t) im 2-Poisson-Zyklus

Wir betrachten im folgenden die Fluktuation $\Delta x(t)$ der Molekülzahl X um ihren Mittelwert. Es gilt also:

(38)
\begin{align} x(t)=\overline{x}+\Delta x(t) = \frac{ax_g}{a+d}+\Delta x(t). \end{align}

Einsetzen obiger Gleichung in den Driftterm $f(x)$ der stochastischen DGL liefert:

(39)
\begin{align} f(\overline{x}+\Delta x)=-(a+d)\Delta x. \end{align}

Da $\dot{x}=\frac{d \Delta x}{dt}$, kann somit die stoch. DGL geschrieben werden als

(40)
\begin{align} \frac{d \Delta x}{dt} = -(a+d)\Delta x + \sqrt{2R_{eq}}\;\zeta(t), \end{align}

oder

(41)
\begin{align} \left[ \frac{d}{dt} + (a+d) \right] \Delta x(t) = \sqrt{2R_{eq}}\;\zeta(t). \end{align}

Für den Differentialoperator in der eckigen Klammer kann man sich eine Greensfunktion definieren:

(42)
\begin{align} \left[ \frac{d}{dt} + (a+d) \right] G(t) = \delta(t-0). \end{align}

Sie lautet explizit

(43)
\begin{align} G(t)=\theta(t) e^{-(a+d)t}. \end{align}

Somit ist die allgemeine Lösung der stoch. DGL gegeben durch die Faltung

(44)
\begin{align} \Delta x(t) = \int_{-\infty}^t dt^{\prime} G(t-t^{\prime}) \sqrt{2R_{eq}}\;\zeta(t^{\prime}). \end{align}

Mit Hilfe obiger Gleichung läßt sich nun die Autokorrelationsfunktion $C_{xx}(\tau)$ der X-Molekülzahl bestimmen:

(45)
\begin{align} C_{xx}(\tau) = \left< \Delta x(\tau) \Delta x(0) \right> = 2R_{eq} \left< \int_{-\infty}^\tau dt^{\prime} G(\tau-t^{\prime}) \zeta(t^{\prime}) \; \int_{-\infty}^0 dt^{\prime\prime} G(0-t^{\prime\prime}) \zeta(t^{\prime\prime})\right>. \end{align}
(46)
\begin{align} C_{xx}(\tau) = 2R_{eq} \int_{-\infty}^\tau dt^{\prime} \int_{-\infty}^0 dt^{\prime\prime} G(\tau-t^{\prime}) G(0-t^{\prime\prime}) \left< \zeta(t^{\prime}) \zeta(t^{\prime\prime})\right>. \end{align}
(47)
\begin{align} C_{xx}(\tau) = 2R_{eq} \int_{-\infty}^\tau dt^{\prime} \int_{-\infty}^0 dt^{\prime\prime} e^{-(a+d)(\tau-t^{\prime})} e^{(a+d)t^{\prime\prime}} \delta(t^{\prime}-t^{\prime\prime}). \end{align}
(48)
\begin{align} C_{xx}(\tau) = 2R_{eq}e^{-(a+d)\tau} \int_{-\infty}^\tau dt^{\prime} \int_{-\infty}^0 dt^{\prime\prime} e^{(a+d)(t^{\prime}+t^{\prime\prime})} \delta(t^{\prime}-t^{\prime\prime}). \end{align}

Als nächstes eliminieren wir mit Hilfe der $\delta$-Funktion das $t^{\prime}$-Integral. Hierdurch kollabiert das Doppelintegral zu

(49)
\begin{align} \int_{-\infty}^0 e^{2(a+d)t^{\prime\prime}} dt^{\prime\prime}, \end{align}

was sich weiter zu $\frac{1}{2(a+d)}$ berechnet.

Insgesamt erhält man für die Autokorrelationsfunktion

(50)
\begin{align} C_{xx}(\tau) = \frac{R_{eq}}{a+d} e^{-(a+d)\tau}, \end{align}

oder

(51)
\begin{align} C_{xx}(\tau) = \left[\frac{adx_g}{(a+d)^2}\right] e^{-\tau/\left[(a+d)^{-1}\right]}. \end{align}

Die Autokorrelation zerfällt also exponentiell mit einer Zeitkonstanten

(52)
\begin{align} \tau_c=\frac{1}{a+d}. \end{align}

[10] M.S.: Anmerkungen zu [9]

ab Gleichung (51) wurde a+d mit a+b verwechselt

Gl. (42)

(53)
\begin{align} \left[ \frac{d}{dt} + (a+d) \right] G(t) = \delta(t-0). \end{align}

einfach $\delta(t-0) = \delta(t)$

in Gl. (45)

(54)
\begin{align} C_{xx}(\tau) = \left< \Delta x(\tau) \Delta x(0) \right> = 2R_{eq} \left< \int_{-\infty}^\tau dt^{\prime} G(\tau-t^{\prime}) \zeta(t^{\prime}) \; \int_{-\infty}^0 dt^{\prime\prime} G(0-t^{\prime\prime}) \zeta(t^{\prime\prime})\right>. \end{align}

sollte doch über t gemittelt werden (Autokorrelation bzgl. der Zeit t ). Auch wenn i.a. Zeit- und Ensemblemittelwert für stationäre Systeme gleich ist (war das nicht die Definition von stationär?) sollte man dasselbe Ergebnis kriegen, wenn man "sauber" nicht t=0 setzt, sondern über t mittelt, also:

(55)
\begin{align} C_{xx}(\tau) = \left< \Delta x(\tau +t) \Delta x(t) \right> = 2R_{eq} \left< \int_{-\infty}^{\tau +t} dt^{\prime} G(\tau +t -t^{\prime}) \zeta(t + t^{\prime}) \; \int_{-\infty}^0 dt^{\prime\prime} G(t-t^{\prime\prime}) \zeta(t+ t^{\prime\prime})\right>. \end{align}

dann kann man aber nicht einfach die Greensfunktionen aus dem Mittelwert herausziehen …


[11] C.M.: Zu [10]

  • Die b/d-Verwechslung habe ich ausgebessert.
  • Ich habe bisher noch keine Herleitung der Autokorr.fkt. hingekriegt, die ohne Verwendung des Ensemble-Mittels (statt Zeitmittel) auskommt. Aber man kann ja die Stationarität getrost voraussetzen.

[12] M.S.: Goldbeter in 2-Poisson-Näherung

Ich habe soeben ein bisschen gerechnet um den Goldbeter im 2 Poisson-Zyklus nachzubauen und bin dabei auf die Rate $k_1 = k_2 = 0.0001$ gekommen, was bei 100 Teilchen jeweils eine totale Rate von .01 bedeutet.

Simuliert man im Normalen 2 Poisson Zyklus nur "kurz" d.h. unter T=10000 (was beim Goldbeter ja bereits nicht mal mehr auf dem 16 GB Rechner funktionierte) so bekommt man dasselbe schlechte Ergebnis, bei Zeiten von T = 1000000 und höher wird das Ergebnis jedoch wieder der erwartete Gauß …

D.h. unsere Statistik beim Goldbeter reicht einfach nicht aus …

Ich habe die Beispieldateien für den 2 Poisson Zyklus angehängt und die Rechnung ist ganz am ende des Pdf vom 2009-2-12

Ich werde jetzt mal die Enzymzahl a und b erhöhen und damit die Rate erhöhen, was effektiv mehr Werte in derselben Gesamtzeit bedeuten sollte. Das wird aber den Rechenaufwand für den PC nicht vermindern.


[13] M.S.: Neue Werte für Goldbeter

Habe ein bisschen herumgespielt und die totale Rate erhöht, ohne unsere Prä-Equilibrium-Näherung zu verlassen. Dann zeigt sich doch ein ziemlich exakter Gauß. Eine Analyse der Autokorrelationsfunktion scheint mir nicht sinnvoll, da die möglichen Zeiten so gering sind.

Ich beschreibe die Gewählten Werte in einer extra pdf.


[14] C.M.: Goldbeter ohne Prä-Equilibrium

Habe den stochastischen Goldbeter-Koshland-Zyklus in einem Parameter-Regime außerhalb der Prä-Equilibrium-Näherung untersucht:

Reaktionen:
X0 + A <-> AX0 -> A + X
X + B <-> BX -> B + X0

Ratenkonstanten:
k_AX0_Binding = 1.0;
k_AX0_Decay = 0.0;
k_X_Production = 1.0; (wird als kX variiert)
k_BX_Binding = 1.0;
k_BX_Decay = 0.0;
k_X0_Production = 1.0;

Teilchenzahlen:
[X0]+[X]=200;
[A]=[B]=10;

Simulationsdauer:
100 000 sec;

Resultate:

  • Solange kX > kX0, also eine asymmetrische Situation vorliegt, gibt es im Mittel mehr X- als X0-Teilchen. Man erhält also 2 getrennte Peaks für P(X0) und P(X).
  • Die einander zugewandten Ausläufer der PDFs sind exponentiell verteilt, wobei die Abfallkonstante umso kleiner wird, je weiter sich kX von oben kX0 annähert.
pdfEinfLog.gif
  • In einem sehr engen Bereich um kX = kX0 verschmelzen die beiden PDFs (bzw. werden identisch). Die verschmolzene PDF hat auf der Seite kleiner Molekülzahlen die Form eines Potenzgesetzes (siehe doppelt-logarithmische Auftragung im unteren Bild). Der Exponent des Potenzgesetzes läßt sich durch die Parameter in weiten Grenzen einstellen (in den Bildern unten nicht gezeigt).
pdfLogLog.gif
  • Das Verhalten ändert sich nicht qualitativ, wenn die Zerfalls-Ratenkonstanten für die Enzym-Substrat-Komplexe nicht Null sondern nur klein gewählt werden.
  • Bei Erhöhung der Teilchenzahlen werden die PDFs sehr "rauschig". Warum ist unklar.
  • Die Autokorrelation—Funktion von [X](t) zerfällt exponentiell. Dabei nimmt die Korrelationszeit zu, je weiter sich kX an kX0 annähert.
acc1.gif

[15] C.M.: Protokoll letzter Projekt-Tag

  • Das Ziel des letzten Tages bestand darin, zu belegen, daß eine enzymatische Umwandlung vom Typ $A+e \leftrightarrow eA \rightarrow B+e$, wobei [A] als konstant angenommen wird, tatsächlich näherungsweise als einfacher Poissonprozeß $A \rightarrow B$ genähert werden kann (im Sinne eines Coarse Grainings).
  • Zunächst wurde die effektive Produktionsrate $k_{eff}$ der B-Moleküle als Funktion der mikroskopischen Parameter berechnet, und zwar sowohl für das Regime des Pre-Equilibriums (Michaelis-Menten-Theorie) als auch für den Grenzfall einer rein sequentiellen Zwei-Stufen-Reaktion (kein Rückzerfall von eA). Die analytische Produktionsrate stimmt mit der numerischen Simulation der enzymatischen Umwandlung in beiden Fällen überein. Für den sequentiellen Prozeß ergibt sich eine zu Michaelis-Menten äquivalente - und recht ähnliche - Formel, obwohl hier kein Gleichgewicht zwischen A+e und eA vorliegt.
  • Weiterhin wurde getestet, ob auch die zeitliche Statistik der diskreten Produktionsereignisse denen eines Poisson-Prozesses entspricht. Hierzu wurde die PDF der Zeitintervalle zwischen aufeinanderfolgenden Produktions-Ereignissen numerisch bestimmt. Es ergibt sich eine Exponentialverteilung mit der charakteristischen Abfallzeit $\tau_c = 1/k_{eff}$, in Übereinstimmung mit dem Poissonprozeß.
  • Somit darf der Goldbeter-Zyklus durch einen 2-Poisson-Zyklus genähert werden.

[16] M.S.: Frage zu Poissonverteilung

Ich habe die Herleitung der Wahrscheinlichkeit der Zeitabstände $\tau$zwischen zwei Ereignissen eines Poissonverteilten Prozesses mal grundsätzlich hergeleitet, da die Unterscheidung zwischen der Ereignisrate $\lambda$ welche die mittlere Anzahl der Ereignisse in einem festen Zeitraum T darstellt und der Ratenkonstante $k_{total}$ in "Max-Thesis" nicht ganz klar wurde - ich bitte um Korrekturen (evtl. habe ich wieder mal 1/k mit k vertauscht etc.)

Die Poissonverteilung habe die Form:

(56)
\begin{align} P_{\lambda}(X = x) = \frac{\lambda^x}{x!}exp(-\lambda)~~~mit~ \lambda~\in~\mathbb{R}_{>0}~und~x~\in~\mathbb{N}_0 \end{align}

Ich habe einfach die Ratenkonstante als $k_{total}*T = \lambda$ definiert wobei $\lambda$ die Anzahl der Ereignisse in dem festen Zeitraum T darstellt (ansonsten macht ja obige Definition der Poissonverteilung keinen Sinn) …

k erweitert quasi diese Definition nur auf beliebige Zeiträume k …

Der von uns vorher immer benutzte Zusammenhang $k = \frac{1}{\lambda}$ kann eigentlich nicht so einfach sein, da:

- in lambda der betrachtet Zeitraum T schon drinsteckt, k jedoch nicht
- bei Erhöhung von der Rate k der Mittelwert $1/ k = \lambda$ niedriger werden würde, was bedeuten würde höhere Rate k führt zu GERINGERER Mittlerer Häufigkeit des Auftretens eines des Ereignisses in einem festen Zeitraum T

meine ausführliche Herleitung siehe PDF "Zeitabstände_bei_Poissonverteilung.pdf"

Zeitabstände zwischen zwei Ereignissen einer Poissonverteilung


Könntest du mir bitte noch den genauen Titel des Buches aus dem du mir den Abschnitt über braunsche Bewegung und die FPG - Herleitung kopiert hast schicken - für die Referenzenliste.

Kann man die Langevin-Gleichung als Verallgemeinerung einer normalen deterministischen DGL verstehen? Ich habe gelesen sie wäre eine Verallgemeinerung der Newtonschen BGL, was aber meiner Meinung nach nicht passt, da eine BGL ja keine Differntialgleichung mehr sein sollte, sondern nur mehr von der From x(t) = …. (ohne d / dt x etc.)


[17] C.M.: Antwort auf [16]

Ich glaube, wir sind hier nur einer Sprach-Verwirrung zum Opfer gefallen. Meine persönliche Ansicht ist wie folgt:

  • Man sollte das $\lambda$ in der Poisson-Verteilungs-Formel nicht als "Ratenkonstante" bezeichnen. Eine echte Rate hat immer die Einheit $\frac{1}{sec}$, was ja bei der dimensionslosen Größe $\lambda$ nicht der Fall ist. Für mich wäre eine gute Bezeichnung für $\lambda$ z.B. "Mittlere Ereigniszahl".
  • Leider weiss ich nicht genau, welches $k_{total}$ in der "Max-Thesis" Du meinst. Normalerweise würde das Symbol auf eine echte Rate hindeuten. In diesem Fall wäre also $\lambda = k_{total}T$. Allerdings wird in der Max-Thesis z.B. in Formel (3.4) das Symbol k auch für eine mittlere Ereigniszahl benutzt. Man muß eben in jedem Einzelfall die Einheiten analysieren, um die korrekte Bedeutung jedes Symbols zu erschließen.
  • Den falschen Zusammenhang $k = \frac{1}{\lambda}$ haben wir sicher nicht benutzt. Das k muß in jedem Fall proportional zu $\lambda$ sein und niemals umgekehrt proportional.
  • Deine Herleitung der Poisson-Verteilung ist korrekt - bis auf die nach meiner Meinung falsche Verwendung des Begriffes "Ratenkonstante".
  • Die Referenz für die Herleitung der Fokker-Planck-Gleichung ist das Buch von Schweitzer.
  • "Kann man die Langevin-Gleichung als Verallgemeinerung einer normalen deterministischen DGL verstehen?". Ich würde schon sagen. Eine normale, nicht-autonome, deterministische DGL enthält irgendwelche externen Variablen (Stimulus), welche auf die Zeitentwicklung der inneren Systemvariablen (Response) einen Einfluß haben. Normalerweise haben die Stimulus-Variablen einen wohldefinierten Zeitverlauf. Bei der Langevin-Gleichung ist aber mindestens eine der Stimulus-Variablen eine in ihrem Zeitverlauf unbestimmte Zufallsvariable, über die nur statistische Eigenschaften (Mittelwerte, Verteilungsdichten, Autokorrelationsfunktionen, ..) bekannt sind.
  • Die Newtonsche BGL ist nach allgemeinem Sprachgebrauch einfach die berühmte $F=m\cdot a$-Formel, oder ausführlicher $\frac{\partial^2}{\partial t^2} x(t) = F_{total}(t)/m$. Bewegungsgleichungen sind in der Mechanik in der Regel Differentialgleichungen.

[18] M.S.: zu [16]/[17]

Ja anscheinend muss ich nur mit der Sprache aufpassen (… Formeln sind der Prosa im Punkt Exaktheit einfach überlegen ;-)

Laut abstrakter definition der Poisson-Verteilung wird der Parameter $\lambda$ für die mittlere Häufigkeit des Auftretens des Ereignisses in einem festen Zeitraum T mit abstrakt mit "Ereignisrate" bezeichnet - das werde ich einfach übernehmen.

Für $k_{total}$ also die Anzahl der Ereignisse in einem festen Zeitraum T geteilt durch T habe ich keine "offizielle" Bezeichnung gefunden … um verwirrung zu vermeiden habe ich jetzt einfach "… dazu wird jetzt die neue Konstante $k_{total}$ definiert " geschrieben.


[18] M.S.: $\alpha$-stabilität, Art der ensemblemittelwerte

Zunächst habe ich hier kurz die "Eckpunkte" der $\alpha$-stabilität aufgeschrieben. Mehr ist für die Projektarbeit auch nicht nötig denke ich:

$\alpha$-Stabilität einer Zufallsgröße $\zeta (t)$:
$a*\zeta_A (t) + b*\zeta_ (t) = \sqrt{a^2 + b^2} \zeta (t)$

wobei $\zeta_A (t), \zeta_B (t)~und~\zeta (t)$ unabhängige Zufallsgrößen des gleichen Typs sind - dies ist z.B. für die Poisson- und die Gaußverteilung erfüllt.

Wie sind die Mittelwerte in Gleichungen (5) bis (10) zu interpretieren?

Da hier und auch später (z.B. Gl. (45) bis (48)) öftmals verwendet wird, dass Mittelwert und Zeitliche Integration vertauschen, sollten diese Mittelwerte als Ensemblemittelwerte interpretiert werden (oder habe ich hier etwas übersehen?)


[19] C.M.: Zu [18]

  • $\alpha$-Stabilität: Ich würde darüber auch nicht zu viele Worte verlieren (Dein Bericht ist ohnehin schon lang genug !). Nähere Infos zum Thema finden sich bei Bedarf hier.
  • Mittelwerte: Alle $<>$-Mittelwerte in obigen Herleitungen sind im Ensemble-Sinne zu interpretieren, was also Stationarität voraussetzt.
  • Projekt-Bericht: Nur ein rein ästhetischer Vorschlag: Ich würde die Bilder kleiner machen und mit viel weißem Rand versehen - das sieht wesentlich edler aus. Insbesondere die Bild-Überschriften und Achsen-Beschriftungen finde ich zu groß geraten. Ansonsten finde ich Deinen Bericht beim ersten Überfliegen sehr gut. Ich werde ihn während der nächsten Tage genauer lesen.

[20] M.S.: neue Bezeichnungen der Konstanten

Ich schreibe hier nur kurz alle Bezeichnungen der Raten … usw. auf, auf welche wir uns so eben geeinigt haben, dass ab jetzt keine Uneinheitlichkeiten mehr auftreten bzw. man wieder nachsehen kann:

$\lambda$ = "mittlere Ereigniszahl" Einheit: [1]

$k_{total}$ = "mittlere Ereignisrate" Einheit: $[\frac{1}{Zeit}]$

$k_1 , k_2 , k_{act},usw.$ = "chemischer Ratenkoeffizient" Eineit: $[ \frac{1}{Zeit * (Moleküleinheit)^n} ]$


Unless otherwise stated, the content of this page is licensed under Creative Commons Attribution-ShareAlike 3.0 License