Modellierung enzymatischer Reaktionen – Rechnungen zur Biokatalyse

Einleitung

Enzyme sind Biokatalysatoren, welche die meisten chemischen Reaktionen in lebenden Organismen steuern. Sie sind durch die Evolution über Millionen von Jahren optimiert worden und zeigen daher in der Regel eine sehr hohe Aktivität und Spezifität für die von ihnen katalysierten Umsetzungen. Die Aufklärung der Reaktionsmechanismen in Enzymen ist natürlich zunächst einmal von unmittelbarem biochemischen, medizinischen und pharmazeutischen Interesse. Darüber hinaus kann ein genaues Verständnis dieser Mechanismen und der zugrunde liegenden Prinzipien auch in der Chemie hilfreich sein, wenn es darum geht, komplexe Systeme im Nanobereich mit ganz spezifischen Eigenschaften zu entwerfen und zu synthetisieren.

Die Modellierung enzymatischer Reaktionen war lange ein Traumziel in der Theoretischen Chemie. Hierfür sind die etablierten klassischen Kraftfelder für biomolekulare Simulationen nicht geeignet, weil bei Reaktionen chemische Bindungen gebrochen und neu geknüpft werden – solche elektronischen Prozesse müssen explizit quantenmechanisch behandelt werden. Allerdings ist die komplette quantenmechanische Beschreibung eines gesamten Enzyms mit seinen Tausenden von Atomen sehr aufwändig und kaum praktikabel. Einen Ausweg bietet die Tatsache, dass das Reaktionsgeschehen in einem Enzym in der Regel in einem kleinen Raumbereich („aktives Zentrum“) lokalisiert ist. Dies legt den Ansatz nahe, das Gesamtsystem durch ein kombiniertes QM/MM-Verfahren zu beschreiben, bei dem der elektronisch relevante Teil quantenmechanisch (QM) und der Rest molekülmechanisch (MM) behandelt wird. Im Folgenden wird zunächst ein Überblick über aktuelle Entwicklungen bei QM/MM-Methoden gegeben, danach werden exemplarisch einige QM/MM-Untersuchungen zu den Reaktionsmechanismen in Molybdän-Enzymen vorgestellt.

QM/MM-Verfahren

Bei der Kombination von QM- und MM-Methoden sind eine Reihe konzeptioneller Probleme zu lösen, insbesondere hinsichtlich der Kopplung von QM- und MM-Region (mechanische oder elektronische Einbettung, mit oder ohne MM-Polarisation) und hinsichtlich der QM/MM-Grenzregion (Verwendung von Linkatomen oder Pseudopotenzialen, Minimierung elektrostatischer Artefakte). Hierfür sind in den letzten Jahren allgemein akzeptierte Standardprozeduren etabliert worden, sodass QM/MM-Rechnungen mittlerweile zur Methode der Wahl bei der Modellierung enzymatischer Reaktionen geworden sind [1]. Auf der Software-Seite werden QM/MM-Methoden am besten modular implementiert. Dabei kontrollieren die zentralen Module den Ablauf der Rechnung und erledigen die generischen QM/MM-Aufgaben, während die eigentlichen QM- und MM-Rechnungen in separaten Programmen erfolgen, die über geeignete Schnittstellen angebunden sind. Ein Beispiel für eine derartige Software ist das von uns mit entwickelte ChemShell-Paket, dessen Aufbau in Abbildung 1 skizziert ist.

Bei QM/MM-Anwendungen an Enzymen wird die QM-Region heutzutage meist mit der Dichtefunktionaltheorie (DFT) beschrieben, z.B. unter Verwendung des B3LYP-Funktionals, während die MM-Region durch eines der etablierten Protein-Kraftfelder mit festen Atomladungen (z.B. CHARMM, AMBER oder GROMOS) modelliert wird. Da die Proteinumgebung sehr polar ist, wird generell eine elektronische Einbettung der QM-Region bevorzugt, während offene Valenzen an der QM/MM-Grenze meist hinreichend gut durch Wasserstoff-Linkatome abgesättigt werden können. Zur Untersuchung von enzymatischen Reaktionen werden typischerweise auf der QM/MM-Potenzialfläche des Grundzustands Reaktionspfade berechnet und die relevanten stationären Punkte (Minima und Übergangszustände) lokalisiert, um die zugehörigen Aktivierungsbarrieren zu bestimmen. Wegen der konformationellen Vielfalt von Enzymen geschieht dies meist für eine Reihe von Startstrukturen, die aus Schnappschüssen einer klassischen Molekulardynamik(MD)-Simulation gewonnen werden. Bei dem anfänglichen Setup des Systems folgt man in der Regel den Protokollen, die für klassische MD-Simulationen etabliert worden sind.

Die hier skizzierte Standard-Vorgehensweise kann natürlich in vielerlei Hinsicht verbessert werden. Zur Erhöhung der Genauigkeit können zum einen korrelierte ab initio Methoden als QM-Komponente zum Einsatz kommen (zumindest in single-point Energieberechnungen), zum anderen kann auch die QM-Region sukzessive vergrößert werden (z.B. unter Verwendung von linear skalierenden QM-Algorithmen), um den Einfluss des MM-Teils auf die Ergebnisse zu minimieren. Alternativ kann man versuchen, durch den Einsatz von polarisierbaren Kraftfeldern mit variablen Atomladungen die MM-Region besser zu beschreiben, wobei man dann auch bei der Einbettungsprozedur die MM-Polarisation berücksichtigen muss. Um entropische und andere temperaturabhängige Effekte auf die Reaktivität zu erfassen, muss man die freie Energie entlang der Reaktionspfade berechnen. Hierzu stehen auf QM/MM-Niveau eine Reihe von Verfahren zur Verfügung (z.B. thermodynamische Integration, umbrella sampling, free energy perturbation theory), bei denen QM/MM MD-Simulationen entlang vorher optimierter Reaktionswege durchgeführt werden, um ein angemessenes Sampling aller konformationellen Freiheitsgrade zu gewährleisten. Weitere mögliche Varianten sind in der Literatur dokumentiert [1].

Auch wenn die QM/MM-Verfahren heute eine gewisse Reife und Standardisierung erreicht haben, gibt es noch viel Raum für methodische Entwicklungen. Ein Beispiel hierfür ist die Erweiterung des QM/MM-Ansatzes auf ein Dreischichten-Modell (QM/MM/Kontinuum), bei dem Randpotenziale benutzt werden, um den äußeren Teil der MM-Region und das umgebende Lösungsmittel zu beschreiben [2]. Ein solches Multiskalen-Verfahren stellt sicher, dass die langreichweitigen elektrostatischen Wechselwirkungen in einem solvatisierten Enzym konzeptionell sauber beschrieben werden, und es verringert gleichzeitig den Rechenaufwand relativ zu QM/MM signifikant, weil die Zahl der explizit behandelten Atome ohne merkliche Einschränkung der Genauigkeit um etwa eine Größenordnung verringert werden kann. Wir haben diesen Ansatz in ChemShell in modularer Form implementiert, sodass er für alle Typen von QM-Methoden effizient genutzt werden kann [2]. Ein zweites Beispiel betrifft den Einsatz von QM/MM-Energien als Randbedingungen bei der kristallographischen Bestimmung von Proteinstrukturen (quantum refinement). Hierfür werden derzeit standardmäßig klassische Kraftfelder benutzt, die jedoch gerade in den kritischen Bereichen (z.B. Inhibitoren, Chromophore oder Metallzentren) oft relativ ungenau sind, sodass die Verwendung von QM/MM-Methoden deutliche Vorteile bringen kann [3].

QM/MM-Rechnungen an Molybdän-Enzymen

Ein umfassender Überblick über aktuelle QM/MM-Studien zu enzymatischen Reaktionen findet sich in einem kürzlich erschienenen Übersichtsartikel [1]. Im Fokus unserer eigenen mechanistischen Arbeiten standen in den vergangenen Jahren die Reaktionen, die durch Cytochrom P450 [4] und durch Molybdän-Enzyme [5-8] katalysiert werden. Im Folgenden berichten wir exemplarisch über unsere QM/MM-Untersuchungen zu der Oxidation von Acetaldehyd zu Essigsäure durch Aldehyd-Oxidoreductase (AOR) [5, 6] und der Umwandlung von Xanthin in Harnsäure durch Xanthin-Oxidase (XO) [7, 8]. Beide Enzyme sind von medizinischem Interesse: AOR katalysiert den Abbau bestimmter Medikamente im Körper und kann hierbei Cytochrom P450-Enzyme ersetzen, während eine zu hohe XO-Aktivität in Niere und Leber zu einem erhöhten Harnsäurespiegel und damit zu Gicht führen kann. Vor diesem Hintergrund ist es sicher wünschenswert, die Wirkungsweise dieser beiden Enzyme in allen mechanistischen Einzelheiten zu verstehen.

Bei der theoretischen Modellierung von enzymatischen Reaktionen geht man von einer experimentell bestimmten Kristallstruktur aus, die nach den üblichen Protokollen bei klassischen Simulationen für die QM/MM-Rechnungen präpariert und in eine Wasserkugel eingebettet wird, um physiologischen Bedingungen nahe zu kommen. Im Falle von AOR starten wir mit einer publizierten Struktur von Desulfovibrio gigas (pdb 1VLB) und erhalten nach dem Setup das in Abbildung 2 gezeigte solvatisierte Modell mit insgesamt 23.357 Atomen (davon 9.711 Atome in den 3.237 Wassermolekülen). Im aktiven Zentrum befinden sich das Substrat Acetaldehyd und der Cofaktor, ein Molybdän(VI)-Komplex mit je einem Oxo-, Sulfido-, Hydroxy- und Dithiolen-Liganden. Die QM-Region (Abb. 3) umfasst neben dem Substrat und dem Cofaktor-Modell auch die Seitenkette der Aminosäure Glu869, welche in dieser Enzymfamilie konserviert ist und somit mechanistisch wichtig sein könnte (insgesamt 25 QM-Atome). Ungeachtet dieses experimentellen Befundes hatten frühere theoretische Modellstudien nur die Reaktion zwischen Substrat und Cofaktor in der Gasphase untersucht und hierfür sowohl einen konzertierten einstufigen Mechanismus als auch einen zweistufigen Mechanismus mit anfänglicher Substratkoordination am Metall in Betracht gezogen. Bei den QM/MM-Rechnungen finden wir diese beiden Reaktionspfade auch im Enzym, allerdings ist ein alternativer dritter Reaktionsweg deutlich günstiger (Abb. 3). Hierbei wird der Cofaktor in einem von der Lewis-Base katalysierten Mechanismus zunächst durch einen Protonentransfer von der Hydroxygruppe zur konservierten Aminosäure Glu869 aktiviert (Intermediat IM1), das entstehende Oxyanion kann dann sehr leicht durch einen nucleophilen Angriff auf das Carbonyl-Kohlenstoffatom des Acetaldehyds ein tetrahedrales Intermediat (IM2) bilden, welches schließlich durch einen formalen Hydridtransfer zum Produkt umlagert. Dieser letzte Schritt ist geschwindigkeitsbestimmend, er hat aber immer noch eine recht niedrige Barriere [5].

Die QM/MM-Rechnungen am AOR ergeben somit aufgrund der Einbeziehung von Glu869 einen qualitativ anderen Mechanismus als die früheren Gasphasen-Modellrechnungen. Ist dieses B3LYP/CHARMM-Ergebnis glaubhaft? Dies kann durch genauere Rechnungen geprüft werden. Zu diesem Zweck haben wir in einer weiteren QM/MM-Studie [6] hochgenaue korrelierte ab initio Methoden als QM-Komponente eingesetzt und zusätzlich mittels ausgedehnter QM/MM MD-Simulationen auch freie Aktivierungsenergien ermittelt (unter Einschluss entropischer Effekte). Die besten so erhaltenen freien Aktivierungsbarrieren unterscheiden sich für die einzelnen Reaktionsschritte (Abb. 3) zwar um typischerweise 3 kcal/mol von den ursprünglichen B3LYP/CHARMM-Werten, die qualitativen Schlussfolgerungen zum bevorzugten Mechanismus bleiben hiervon jedoch unberührt.

Die zweite untersuchte Reaktion, die Oxidation von Xanthin zu Harnsäure durch XO, ist komplizierter. Zum einen deuten Mutagenese-Experimente darauf hin, dass mindestens drei Aminosäuren in der unmittelbaren Umgebung des aktiven Zentrums für die Reaktion wesentlich sind (Glu802, Glu1261, Arg880), weil deren Substitution die enzymatische Aktivität beeinträchtigt oder sogar ganz unterbindet. Zum anderen ist das Substrat Xanthin (Abb. 4) sehr wandlungsfähig, da es in verschiedenen tautomeren Formen und außerdem auch als deprotoniertes Anion vorliegen kann – darüber hinaus kann das planare Xanthinmolekül in der Bindungstasche des Enzyms zwei verschiedene Orientierungen einnehmen (upside und upside-down). Wir mussten daher in unserer QM/MM-Studie sieben Varianten durchspielen, um alle plausiblen Reaktionswege auszuloten [7]. In dem besten gefundenen Mechanismus (Schema 2) reagiert Xanthin in der upside-down Orientierung. Zunächst erfolgen unter Beteiligung von Glu1261 und einem Wassermolekül drei sukzessive Protonentransfers (vom Reaktand zu IM1), bei denen sowohl der Cofaktor als auch das Substrat aktiviert werden. Bei dem Cofaktor geschieht dies – wie im Falle des AOR – durch eine Deprotonierung der Hydroxygruppe. Beim Xanthin-Substrat wird das ursprünglich vorliegende stabilste Xanthin-Tautomer (Abb. 4 links) in eine reaktive Form überführt, die bestens für die folgenden Reaktionsschritte präpariert ist (Abb. 4 rechts): Das Target-Kohlenstoffatom für den nucleophilen Angriff trägt eine große positive Ladung (rot), während im Sechsring die beiden unten stehenden Heteroatome stark negativ werden (blau) und somit im letzten Schritt der Reaktionssequenz (IM2 und folgender Übergangszustand) durch das benachbarte, positiv geladene Arg880 eine signifikante Stabilisierung erfahren. Auf diese Weise wird die Barriere für den geschwindigkeitsbestimmenden letzten Schritt im Enzym deutlich abgesenkt, während die Barrieren für die einfacheren ersten Schritte im Vergleich zur Gasphase sogar etwas ansteigen. Das Energieprofil für die bevorzugte Reaktionssequenz (Schema 1) ist im Enzym ziemlich „ausgewogen“, mit vergleichbar hohen Barrieren für mehrere Schritte, während es in der Gasphase viel ungleichmäßiger ist und am Ende eine wesentlich höhere Barriere aufweist. Man kann eigentlich nur staunen, wie geschickt im Enzym die polare Proteinumgebung angeordnet ist, sodass das „richtige“ Tautomer für die Reaktion selektiert und der strategisch wichtige Übergangszustand stabilisiert wird. Die katalytische Wirkung von XO ist somit auf die elektrostatische Präorganisation des aktiven Zentrums zurückzuführen.

Gibt es experimentelle Befunde, welche diesen plausiblen Mechanismus in XO stützen? Zur Beantwortung dieser Frage muss man die Literatur sichten. So ist beispielsweise bekannt, dass XO die Oxidation nicht nur bei Xanthin katalysiert, sondern auch bei 1-Methyl-2,6-dioxopurin, nicht jedoch bei 1-Methyl-6-oxopurin. Eine genauere Analyse zeigt, dass man dies bei unserem bevorzugten Mechanismus mit einer upside-down Orientierung des Substrats zwanglos erklären kann, während eine upside Orientierung zu analogem Verhalten bei den verschiedenen Substraten führen müsste [7]. Experimentell ist ferner bekannt, dass der Ersatz des Sulfido- durch einen Oxo-Liganden im Cofaktor von XO dessen Aktivität abschaltet, dass das Substrat 2-Oxo-6-methylpurin etwas leichter reagiert als Xanthin und dass die Mutation Glu802→Gln802 die Aktivität von XO etwas absenkt. Alle diese experimentellen Befunde stehen im Einklang mit QM/MM-Rechnungen [8], die zur Prüfung und Validierung des vorgeschlagenen Mechanismus durchgeführt wurden. Für einen Cofaktor mit Oxo- statt Sulfido-Ligand finden wir einen prohibitiv starken Anstieg der Barriere für den letzten Schritt (Hydridtransfer), bei dem Substrat 2-Oxo-6-methylpurin wird dieselbe Reaktionssequenz wie beim Xanthin realisiert (mit etwas niedrigeren Barrieren), und durch die Mutation Glu802→Gln802 wird der günstigste Reaktionsweg (Schema 1) gestört, sodass man mit der zweitbesten und daher etwas langsameren Variante vorlieb nehmen muss. Last but not least sei noch erwähnt, dass nach Publikation unserer QM/MM-Studien [5-8] zwei Kristallstrukturen in nahe verwandten Enzymen aus der XO-Familie gelöst wurden, bei denen das Substrat wie vorhergesagt in der upside-down Orientierung gefunden wurde [9].

Schlussbemerkungen

Die vorgestellten QM/MM-Rechnungen an den beiden Molybdän-Enzymen haben zu plausiblen Reaktionsmechanismen geführt, die mit den bekannten experimentellen Fakten im Einklang stehen. Sie bieten detaillierte mechanistische Informationen, die in dieser Form experimentell kaum zugänglich sind, und ergänzen somit die experimentellen Arbeiten auf diesem Gebiet durch komplementäre Einsichten. Man darf dabei natürlich nicht vergessen, dass die QM/MM-Rechnungen an Enzymen viele Annahmen und Fehlerquellen enthalten – im Gegensatz zu ab initio Rechnungen an kleinen isolierten Molekülen können wir die Ergebnisse nicht konvergieren und daher auch nicht mit Zuversicht die richtige Antwort vorhersagen. Es ist bei der Modellierung von enzymatischen Reaktionen somit essentiell, den engen Kontakt zum Experiment zu suchen und durch das Zusammenführen aller verfügbaren theoretischen und experimentellen Befunde ein in sich konsistentes Bild vom Mechanismus dieser Reaktionen zu entwickeln.