Stoffwechsel-Modellierung & -Simulation
(Elektronentransportkette von Escherichia coli)

(Stand: Juni 2016)

Im Stoffwechsel (Metabolismus) reagieren Stoffwechselpartner (Me­ta­bo­li­te) miteinander, ändern dadurch ihre Konzentration und bewirken bestimmte makroskopische Funktionen (beispielsweise Ener­gie­ge­win­nung oder Wachstum). Dabei erfüllen Proteine die Aufgabe als Enzyme, also als Katalysatoren von Reaktionen. Die Aktivität eines Enzyms, welche vor allem durch die Menge bzw. Konzentration gegeben ist, hat einen wesentlichen Einfluss auf die Geschwindigkeit einzelner Reak­tio­nen sowie dadurch auch auf andere Reaktionen und letztendlich auf makroskopische Phänomene. Die Enzymmenge kann in Ab­hän­gig­keit geänderter Be­din­gun­gen mithilfe der Genexpression reguliert werden. Dadurch erreichen Lebewesen insgesamt eine Anpassung des Stoffwechsels.

Zur Beantwortung biologischer Fra­ge­stel­lun­gen zu Teilaspekten der Stoffwechsel unterschiedlicher Lebewesen wird die systembiologische Modellierung eingesetzt. In diesem Fall liegt der Fokus weniger auf dem Finden der Reak­tionsstruk­tu­ren als vielmehr auf Mechanismen, ihren Pa­ra­me­tern und insbesondere auf den Teil­mo­del­len, die die Regulation beschreiben. Der Modellierung schließt sich in den meisten Fällen eine Identifikation der unbekannten Pa­ra­me­ter­wer­te an; gefolgt von Si­mu­la­tio­nen und ihrer Auswertung. Im Folgenden wird anhand eines realitätsnahen Beispiels die Anwendung dieser einzelnen Schritte gezeigt.

Fragestellung

---

In den Reaktionen der Elek­tro­nen­trans­port­ket­te (ETK, auch: At­mungs­ket­te) werden Elektronen des zentralen Stoff­wech­sels auf Elek­tro­nen­ak­zep­to­ren übertragen und gleichzeitig ein Pro­to­nen­po­ten­tial geschaffen, welches schließlich zur Bildung des biochemischen Energieträgers ATP verwendet wird. Zur Anpassung an unterschiedliche Be­din­gun­gen werden die Reak­tio­nen der ETK reguliert. Unterstützt durch ein kinetisches Modell wurde durch Henkel et al. die Anpassung des Bak­te­riums Escherichia coli an unterschiedliche Sauer­stoff­be­din­gun­gen untersucht:

[1] S. G. Henkel, A. Ter Beek, S. Steinsiek, S. Stagge, K. Bettenbrock, M. J. Teixeira de Mattos, T. Sauter, O. Sawodny, M. Ederer: "Basic regulatory principles of Escherichia coli's electron transport chain for varying oxygen conditions." In: PLoS ONE, 9(9):e107640, 2014. doi: 10.1371/journal.pone.0107640.

Dabei wurden ausschließlich stationäre Arbeitspunkte im Chemostat (kontinuierliche Kultivierung), siehe [2, 3], bei unterschiedlichen Sauerstoffbedingungen untersucht. Als wesentliche Grundlage für den biologischen und experimentellen Hin­ter­grund sei außerdem auf [4–6] verwiesen. Die vorliegenden Messdaten zeigen ein charakteristisches Verhalten in bestimmten Komponenten (sowohl Me­ta­bo­lit­kon­zen­tra­tio­nen als auch eine Tran­skrip­tions­fak­tor-Aktivität), welches den Aus­gangs­punkt für die kinetische Mo­del­lie­rung bildet. Im Folgenden soll in stark vereinfachter Form (und auf der Basis künstlicher Messdaten) nur ein Teilmodell der regulierten ETK hergeleitet werden (für die Formeln muss JavaScript aktiviert sein), während für eine ausführliche Be­schrei­bung auf die genannte Publikation verwiesen wird. Es soll gezeigt werden, wel­che regulatorische Strukturen die Mess­da­ten erklären können und ob dabei Ubichinon der wesentliche inhibierende Fak­tor für ArcA sein kann.

Modellierung

---

Den Ausgangspunkt der Mo­del­lie­rung bildet die folgende Abbildung der vergleichsweise verlässlichen strukturellen In­for­ma­tion über die Reak­tio­nen der ETK sowie die Existenz der zwei Tran­skrip­tions­fak­toren FNR und ArcA, deren regulatorische Wir­kung im Rahmen der Mo­del­lie­rung analysiert werden soll (da die vorhandenen In­for­ma­tio­nen entweder weniger verlässlich oder widersprüchlich sind bzw. Mess­daten nur unzureichend erklären):

(Im Folgenden werden bei der Beschreibung nicht alle Annahmen begründet) Sauerstoff wird aufgenommen und reagiert in Reaktionen vom Oxidase-Typ. Dadurch wird Ubichinol (QH$_2$) zu Ubichinon (Q) oxidiert. Dieses Redoxpaar ist - wie auch NADH|NAD - außerdem an De­hy­dro­ge­na­se-Reak­tio­nen beteiligt. Die verfügbaren Elektronen des zentralen Stoffwechsels werden durch NADH repräsentiert, wobei das Modell hier so entwickelt wurde, dass es aus den oxidativen Species aufgebaut ist, weshalb noch ein Abfluss von NAD benötigt wird ("Schnitt­stel­le" zum zentralen Koh­len­hy­drat-Me­ta­bo­lis­mus). Als mathematisches Modell dieser Struktur werden folgende Dif­fe­ren­tial­glei­chun­gen aufgestellt:

\begin{align} {\dot c}_\mathrm{O_2} &= v_\mathrm{O_{2},in}-v_\mathrm{Oxi} \\ {\dot c}_\mathrm{Q} &= v_\mathrm{Oxi}-v_\mathrm{Dh} \\ {\dot c}_\mathrm{NAD} &= v_\mathrm{Dh}-v_\mathrm{NADH} \end{align}

Es wird angenommen, dass erstens die Reaktionen quasi irreversibel sind, dass zweitens zusammengefasste Reaktionen der Isoenzyme verwendet werden können, weshalb im Folgenden nur noch von DER Oxidase und DER Dehydrogenase die Rede sein wird und dass drittens die Ge­samt­men­ge an oxidativen und reduzierten Species konstant ist, also Q+QH$_2$=const. und NAD+NADH=const. Für die Ra­ten­glei­chun­gen werden bekannte, das heißt physiologisch plausible, Kinetiken eingesetzt:

\begin{align} v_\mathrm{O_{2},in} &= \frac{v_\mathrm{O_{2},in,100\%}\times t}{t_\mathrm{end}}\\ v_\mathrm{Oxi} &= a_\mathrm{Oxi}\frac{c_\mathrm{O_{2}}}{(K_\mathrm{m,Oxi,O_{2}}+c_\mathrm{O_{2}})}\frac{c_\mathrm{QH_{2}}}{(K_\mathrm{m,Oxi,O_{2}}+c_\mathrm{QH_{2}})}\\ v_\mathrm{Dh} &= a_\mathrm{Dh}\frac{c_\mathrm{Q}}{(K_\mathrm{m,Dh,Q}+c_\mathrm{Q})}\frac{c_\mathrm{NADH}}{(K_\mathrm{m,Dh,NADH}+c_\mathrm{NADH})}\\ v_\mathrm{NADH} &= v_\mathrm{max,NADH}\frac{c_\mathrm{NAD}/c_\mathrm{NADH}}{K_\mathrm{m,NAD}+c_\mathrm{NAD}/c_\mathrm{NADH}} \end{align}

Bei der Oxidase und Dehydrogenase sind die Enzymaktivitäten mit den Maxi­mal­ge­schwin­dig­kei­ten gleichgesetzt. Über $v_\mathrm{O_{2},in}$ wird der Sauerstoff langsam verfahren, wodurch die resultierenden Simulationen mit den stationären(!) Mess­da­ten vergleichbar werden. Nun ist die Frage, ob und wie in diesem Setting die folgenden Messdaten von Sauerstoff, Ubi­chi­non und ArcA für verschiedene Sauerstoffbedingungen (anaerob...mikroaerob...aerob, repräsentiert durch prozentuale "Aerobiose" $a$, siehe [4–6]) erklärt werden können:

Im wesentlichen ist zu beobachten, dass sowohl Sauerstoff als auch Ubichinon im mikroaeroben Bereich ein nahezu konstantes Verhalten aufweisen (mit starken An­stie­gen um 100% Aerobiosis), während ArcA einen charakteristischen "Zickzack"-Verlauf besitzt. Zur sensitiven Regulation ("Stabilisierung", "Homöostase") der Sauer­stoff­kon­zen­tra­tion wird folgende Teil­struk­tur eingeführt (es handelt sich aufgrund der beiden Inhibierungen um eine positive Regulation, jedoch bezüglich Sauerstoff aufgrund der Oxidasereaktion um eine negative Rück­kopp­lung):

Für die beobachtete sensitive Regulation werden zu dieser Struktur noch folgende Gleichungen eingeführt:

\begin{align} R_\mathrm{FNR} &= \frac{c_\mathrm{O_{2}}^{n_\mathrm{FNR}}}{(K_\mathrm{A,FNR,O_{2}}^{n_\mathrm{FNR}}+c_\mathrm{O_{2}}^{n_\mathrm{FNR}})}\\ v_\mathrm{syn,Oxi} &= p_\mathrm{syn,Oxi,min} + p_\mathrm{syn,Oxi,FNR} \times (1-R_\mathrm{FNR})\\ {\dot a}_\mathrm{Oxi} &= v_\mathrm{syn,Oxi}-D \times a_\mathrm{Oxi} \end{align}

Mit der letzten Gleichung ist nun als weitere Variable eine veränderliche En­zym­ak­ti­vi­tät eingeführt worden, welche in der Oxi­da­se-Ra­ten­glei­chung die Funktion einer variablen Maximalaktivität erfüllt. Die Veränderung dieser Aktivität hängt von der Synthese und der "Ver­dün­nung durch Wachstum" ab. Die Syntheserate $v_\mathrm{syn,Oxi}$ wiederum ist eine algebraische Gleichung, in welcher der erste Term eine grundlegende Expression repräsentiert und es einen zweiten Term gibt, welcher negativ vom Tran­skrip­tions­fak­tor abhängt. Das entscheidende für eine sensitive Regulation ist nun die Gleichung für den Tran­skrip­tions­fak­tor, da es sich um eine algebraische Gleichung ähnlich einer Hill-Kinetik handelt, bei welcher ein betraglich hoher negativer Exponent eine inhibierende sensitive Funktionalität gewährleistet (Rand­be­mer­kung: ähnlich zur Ko­ope­ra­ti­vi­tät kann auch hier ein solches Verhalten durch sequentielle Bindungen, mehrere Bin­de­stel­len und andere Effekte erklärt werden). Diese sensitive Regulation ist in der folgenden linken Abbildung gezeigt:

Dabei besitzt FNR eine sinkende Aktivität bei steigender Sauer­stoff­kon­zen­tra­tion mit einem quasi-schaltenden Verhalten um $K_A$ herum. In der rechten Abbildung ist nun der Vergleich der Simulationen für den unregulierten (gestrichelt) und regulierten (durchgezogen) Fall dargestellt. Es ist gut erkennbar, wie ein nahezu konstanter Verlauf der Sauerstoffkonzentration (grün, mit den entsprechneden charakteristischen Parametern der Enzymkinetik und Regulation) bei unterschiedlichem Sauer­stoff­ein­trag in den Reaktor erreicht werden kann (nicht gezeigt: die resultierende, ansteigende Rate der Oxi­da­se­reak­tion).

Nun stellt sich die Frage, ob und wie auf ähnliche Weise auch die Ubichinon-Konzentration und die ArcA-Aktivität erklärt werden kann. Dazu werden die folgenden zwei Varianten (im Folgenden "M1" und "M2") von Regulationsstrukturen analysiert:

Die Gleichungen werden analog zur Re­gu­la­tion der Oxidase aufgestellt. Bei M2 kommt neben einer "lokalen" Regulation der Dehydrogenase durch ArcA noch eine vorsteuernde Regulation durch FNR hinzu (entsprechend ergeben sich dann drei additive Terme in der Syn­the­se­glei­chung).

Da die Bestimmung unbekannter Pa­ra­me­ter­wer­te mithilfe von Messdaten ein wesentlicher Schritt in der Stoff­wech­sel­mo­del­lie­rung darstellt, wird dies im Fol­gen­den anhand der Anpassung dieser beiden Varianten vorgestellt.

Parameteridentifikation

---

Bei der automatischen Bestimmung von Parametern mithilfe von Messdaten, also der Parameteridentifikation (oder Schät­zung), werden die Parameter $\underline{\theta}$ so gefunden, dass die parameterabhängigen, simulierten Verläufe $\hat{y}_i$ der Variablen ihren gemessenen Werten so gut wie möglich entsprechen. Eine Simulation von Mo­del­len, die durch gewöhnliche Dif­fe­ren­tial­glei­chun­gen beschrieben werden, erfolgt durch das automatische "Lösen eines An­fangs­wert­pro­blems" mit entsprechenden Verfahren, die wiederum in diversen Computer-Umgebungen oder -Sprachen implementiert sind. Die Abweichung zwischen Simulation und Messung wird durch eine Zielfunktion (Gütefunktion, u.a.) beschrieben. Unter der Annahme normalverteilter Unsicherheiten der Messwerte (mit Varianzen $\sigma_{i}^2$) kann folgende Zielfunktion aufgestellt werden, bei der die Pa­ra­me­ter­iden­ti­fi­ka­tion durch die Mi­ni­mie­rung der Summe der Feh­ler­qua­dra­te erreicht wird:

\begin{align} \chi^2(\underline{\theta}) &= \sum_{i=1}^{N}\frac{\left(y_{i}-\hat{y}_{i}(\underline{\theta})\right)^2}{\sigma_{i}^2} \end{align}

Die Minimierung erfolgt iterativ ebenfalls durch Implementierungen entsprechender (nichtlinearer) Optimierungsverfahren. Die identifizierten Parameterwerte der beiden Modellvarianten können der folgenden Ta­bel­le entnommen werden:

Es ist hervorzuheben, dass starke Un­ter­schie­de in den Parameterwerten der Re­gu­la­tion auftreten. Vor allem bei der Re­gu­la­tion der De­hy­dro­ge­na­se ist zu beobachten, dass einerseits bei Fehlen der Re­gu­la­tion durch FNR (M1) die ArcA-abhängige Re­gu­la­tion stärker ist als bei M2, bei welcher andererseits FNR den stärksten Einfluss auf die Regulation der De­hy­dro­ge­na­se ausübt.

Systembiologische Analyse und Interpretation

---

In der folgenden Abbildung werden die Simulationen der parameteridentifizierten Modellvariante M1 mit den Messdaten verglichen:

Es ist gut erkennbar, dass die sensitive Re­gu­la­tion tatsächlich auch eine vergleichsweise gute Anpassung für Ubichinon bedeutet. Es ist jedoch nicht möglich, den charakteristischen Verlauf der ArcA-Messdaten zu erklären. Anders sieht es da bei der angepassten Modellvariante M2 aus. Hier wird nicht nur der ArcA-Verlauf gut getroffen, sondern auch der zuvor kaum erkennbare schwache Abfall von Ubichinon im mikroaeroben Bereich gut wiedergegeben:

Dieses Verhalten wird letztendlich durch den zusätzlichen regulatorischen Einfluss von FNR auf die Dehydrogenase erklärbar: die In­for­ma­tion über eine gesteigerte Sauerstoffverfügbarkeit wirkt sich überproportional auf diese Enzymaktivität aus, wodurch der leichte Ubichinonabfall im oberen mikroaeroben Bereich einen An­stieg von ArcA mit sich bringt. Aus biologischer Sicht hat das zum einen den Vorteil, mit dem Transkriptionsfaktor ArcA eine differenziertere Regulation für andere Stoff­wech­sel­we­ge und Isoenzyme zu erreichen, und zum anderen besitzt dieses regulatorische Motiv bessere Eigenschaften bei dynamischen Änderungen, [7], hier also bei schnelleren Schwankungen der Sauer­stoff­ver­füg­bar­keit.

In der oben genannten Publikation wurden die hier beschriebenen Prinzipien der regulierten ETK innerhalb eines umfangreicheren Mo­dells (u.a. mit phämomenologischem Biomassemodell) angewendet und diskutiert.

Referenzen

---

1. [1] S. G. Henkel, A. Ter Beek, S. Steinsiek, S. Stagge, K. Bettenbrock, M. J. Teixeira de Mattos, T. Sauter, O. Sawodny, M. Ederer: Basic regulatory principles of Escherichia coli's electron transport chain for varying oxygen conditions. In: PLoS ONE, 9(9):e107640, 2014. doi: 10.1371/journal.pone.0107640
2. [2] A. Novick, L. Szilard: Experiments with the Chemostat on spontaneous mutations of bacteria. In: Proc. Natl. Acad. Sci. USA, 36(12):708–719, 1950. doi: 10.1073/pnas.36.12.708
3. [3] A. Novick, L. Szilard: Description of the Chemostat. In: Science, 112:715–716, 1950. doi: 10.1126/science.112.2920.715
4. [4] S. Alexeeva, B. de Kort, G. Sawers, K. J. Hellingwerf, M. J. Teixeira de Mattos: Effects of limited aeration and of the ArcAB system on intermediary pyruvate catabolism in Escherichia coli. In: J. Bacteriol., 182(17):4934–4940, 2000. doi: 10.1128/JB.182.17.4934-4940.2000
5. [5] S. Alexeeva, K. J. Hellingwerf, M. J. Teixeira de Mattos: Quantitative assessment of oxygen availability: perceived aerobiosis and its effect on flux distribution in the respiratory chain of Escherichia coli. In: J. Bacteriol., 184(5):1402–1406, 2002. doi: 10.1128/JB.184.5.1402-1406.2002
6. [6] S. Alexeeva, K. J. Hellingwerf, M. J. Teixeira de Mattos: Requirement of ArcA for redox regulation in Escherichia coli under microaerobic but not anaerobic or aerobic conditions. In: J. Bacteriol., 185(1):204–209, 2003. doi: 10.1128/JB.185.1.204-209.2003
7. [7] S. Mangan, U. Alon: Structure and function of the feed-forward loop network motif. In: Proc. Natl. Acad. Sci. USA, 100(21):11980–11985, 2003. doi: 10.1073/pnas.2133841100