Transcript Általános reakciókoordináta használata QM és QM

Általános reakciókoordináta használata QM
és QM/MM felszínen
Fuxreiter Mónika, Petr Kulhanek, Alessandro Laio,
Simon István, Csányi Gábor és Mones Letif
Enzimológiai Intézet és University of Cambridge
[email protected]
2
Reakciók vizsgálata oldatfázisban és enzimatikus
környezetben
 1. Potenciális energiafelszín definiálása (hibrid erőtér)
 2. Mintavételezés molekuladinamika segítségével
 Reakciókoordináta kiválasztása  r 
 Mintavételezési technika
 3. Szabadenergia számítása
„Potential of Mean Force”
F  r    1 ln P r 
Összevetés a kísérleti adatokkal (pl. kcat)
3
Mi legyen a reakciókoordináta?
 Geometriai koordináták: - kötéshossz
- kötésszög
- torzió
- kötéshossz-kötéshossz különbség
- koordinációs szám
-…
 Gond: - rossz átmeneti állapot indikáció és
- rossz a valódi átmeneti állapot mintavételezése
 Létezik olyan univerzális reakciókoordináta, amelyről a priori
feltételezhető, hogy jól indikálja az átmeneti állapotot?
4
Empirikus vegyértékkötés módszer (EVB)
I. resonance state
BO
H

 (x)   D 1  e
q1
i
q
r3
PT
P
O
O
BO
O
(1)
1s  1o  gsol
aqi ( rq r0i , q )
  U
2
A
O
b1
b2
1. rezonancia állapot
Nq
O
Nk
j 1
i
kötés, j
N sz
 U
j 1
O
b1
i
szög , j
Nt
 U
j 1
 Általános reakciókoordináta:
Egapx  1s x   2s x
i
torzió , j
H
O
P
O
O
H
O
b2
2. rezonancia állapot
Szabadenergia
r3
o
i
H
AO
O
II. resonance state
( 2)
 2s   2o  gsol
N nk
i
 Unemkötő
,kl
k ,l 1
1
2
A. Warshel and R. M. Weiss, (1981) Ann. N Y Acad. Sci. 367: 370
Egap  1s   2s
5
Egap mint reakciókoordináta a QM/MM felszínen
Egap számításához klasszikus potenciálfüggvények szükségesek
Indirekt módszer
•
•
Direkt módszer
•
•
Egap vezérelt MD az MM felszínen
QM/MM számítás konfigurációkra
MD a QM/MM felszínen
Egap számítása külön lépésben
A
A
MM
QM/MM
QM/MM
!
Egap
Egap
DFQM /MM ( Egap ) = DFmapping ( Egap ) + DFmapping®QM /MM ( Egap )
(
)
{
(
)
(
)}
DFmapping®QM /MM ( Egap ) = -b -1 ln d Egap - Egap exp éë- EQM /MM Egap - em Egap b ùû
fi = m
¶ é QM /MM
MM
H
( r, Y ) +Vbias
( Egap ( r), t )ùû
ë
¶ri
6
Kiméra programok
Kvantumdinamikai
program
koordinátá k
fbias
1,  2 , Egap 1   2
fQM , fQM / MM , f MM
fi = -
•
¶ é QM /MM
MM
H
Egap ( r), t )ùû
( r, Y ) +Vbias
(
ë
¶ri
DFT alapú: CPMD/Gromos
R. Car and M. Parrinello, (1985) Phys. Rev. Lett. 55: 2471
•
XdynBP
PMF könyvtár
Szemiempirikus módszerek: AMBER
D.A. Case et al., (2005) J. Comp. Chem. 26: 1668
7
A modell

Modell rendszer:
-
Cl + MeCl
ClMe + Cl
-
• Gázfázisban 300 K-en
• Oldatfázisban 300 K-en (659 TIP3P vízmolekula)

Program: AMBER + XdynBP

Felszín: QM(PM3)/MM

Vizsgált koordináták: • Egap
• DD
Dist2
DD=Dist1-Dist2
Dist1
Dist2
Dist1
8
Alkalmazott mintavételi technikák
A

FEP/US
G. M. Torrie, J. P. Valleau, (1977) J. Comput. Phys. 23: 187
!
²0
²
A

Blue Moon
J. Schlitter et al., (2003) J. Chem. Phys. 118: 2057
!
!
!
²
A

Metadinamika
A. Laio, M. Parrinello, (2002) PNAS 99: 12562
!
!
²0
A
F²

Adaptive Biasing Force
E. Darve, A. Pohorille, (2001) J. Chem. Phys. 115: 9169
!
!
²
²
!
!
²
9
Indirekt vs. direkt módszer
Blue Moon
FEP/US_1
FEP/US_2
FEP/US_3
12
F / kcal mol
-1
10
• Indirekt módszer hatékonyabb, ha a
klasszikus felszín alakjában „közel”
áll a QM/MM felszínhez (< 2-3 RT)
8
6
• Ez azonban a priori nem tudható
4
• Ha nincs „közel”: rossz a mintavétel a
QM/MM felszín kritikus tartományaiban
2
0
-200
-150
-100
-50
0
50
-1
Egap / kcal mol
100
150
200
10
A direkt módszer alkalmazhatósága
Ebond / kcal mol
-1
300
= 1.510
= 1.510
= 1.510
= 1.510
= 1.510
= 3.020
= 0.755
200
-2
-4
-6
-1
De= 71.0 r0 = 1.8
De= 152.0 r0 = 1.8
De= 35.5 r0 = 1.8
De= 71.0 r0 = 1.6
De= 71.0 r0 = 2.0
De= 71.0  r0 = 1.8
De= 71.0 r0 = 1.8
F / kcal mol
400
-8
-10
-12
100
-14
-16
0
0
2
4
6
d/A
8
10
-300
-200
-100
0
100
Egap / kcal mol
200
300
-1
• Az erőterek paramétereinek bizonyos határon belüli variálására invariáns az
Egap-függő aktiválásia gát és szabadenergia-különbség
• Direkt módszer nagyobb különbség esetén is biztonságosabban használható
11
MTD
BM
ABF
-2 0 0
-1 0 0
0
E ga p / kc a l m ol
100
200
-1
26
24
22
20
18
16
14
12
10
8
6
4
2
0
-3 0 0
DA / kc a l m ol
DA / kc a l m ol
-1
Egap vs. DD profilok a QM(PM3)/MM felszínen
26
24
22
20
18
16
14
12
10
8
6
4
2
0
300
-2
-1
-1
0
DD /Å
• Egap esetén szimmetrikusabb profil
• Gátmagasságban 10-15% különbség
Különböző TS mintavétel
1
2
12
Reakciókoordináták megbízhatóságának vizsgálata
1. Reakciókoordináta
TS(CV)
 r 
30
Mintavételezési
technika
2. Szabadenergia-profil és
TS indikáció
Kényszerezett MD
 r    ‡
 r   
25
20
F / kcal/mol
MD
15
10
5
0
‡
-300
-100
0
100
200
300
100
200
300
CV
3. Konfigurációk a jelzett TS-ben
K1
K2
K3
K4
K5
… Kn
30
25
20
F / kcal/mol
véletlen sebességekkel
Rövid MD trajektóriák
minden konfigurációból
4. „Lecsurgás” az egyik vagy másik minimumba
K1
K2
K3
K4
K5
… Kn
-200
15
10
5
0
-300
-200
-100
0
CV
13
Egap és DD Geissler-tesztjének eredménye
I. konklúzió:
• DD még egy ilyen egyszerű rendszer
esetén sem ad megbízható eredményt
• Egap esetén a TS indikáció sokkal jobb
Egap
DD
14
Egap vs. DD profilok konvergenciája
Egap
DD
Hatékonyságvizsgálat
Hiszterézisvizsgálat
II. konklúzió:
• DD esetén nagy a hiszterézis, lassabb a konvergencia
• Egap esetén a mintavételezés sokkal jobb
15
Egap alkalmazhatóságának korlátai

Jelenlegi korlátok:
 elsősorban kémiai reakciók vizsgálhatók
 megfelelő minőségű reakció (vegyértékváltozás)
 végállapotok definiálása szükséges (kémiai intuíció!)

Néhány probléma és lehetséges megoldásuk:
Probléma
Lehetséges kutatási irány
• hiányzó erőtérparaméterek
• átmeneti fémek reakciói
• redoxreakciók
• charge constrained DFT
• többlépéses reakciók
• alternatív útvonalak vizsgálata
• multidimenziós szabadenergia-felszín
több Egap terében (metadinamika)
• automatikus Egap kiválasztás több
lehetségesből
16
Többállapotú rendszerek vizsgálata

PT reakció
Modell rendszer:
H
Od
oldatfázisban
Oa

Felszín: QM(PM3)/MM

Egyszerű reakciókoordináták: • DIS = d(Od-H) vagy DIS = d(Oa-H)
• DD = d(Oa-H) – d(Od-H)

Mellékreakciók:
17
Többállapotú rendszerek vizsgálata

Megoldás:
• DIS vagy DD + egyéb O-H kötésekre restraint/constraint
• más geometriai rekciókoordinák alkalmazása (CN, MINDIS)
DF / kcal mol-1
60
50
40
30
20
10
0
0.8
1
1.2
1.4
MINDIS / A
1.6
MINDIS = min {dij (Oi - H j )}
1.8
1
2
3
4
5
6
7
8
9
10
2
18
EGAP többállapotú rendszerekre
Reaktáns állapot(ok)
Termék állapot
6 lehetséges ekvivalens vegyértékállapot
(2 oxigén x 3 hidrogén)
e1R, e2R, e3R, e4R, e5R, e6R

Ekvivalens EGAPek:
EGAPk = ekR - e P
1 vegyértékállapot
e
P
19
Két- és többállapotú EGAPek
EGAPk = ekR - e P
Mellékreakciók!
1 n R P
EGAP = åek - e
n k=1
Nem diszkriminál, túl
sok energiaállapot keveredik!
MINEGAP = min {ekR } - e P
20
MINEGAP
DF / kcal mol-1
50
100
150
MINEGAP / kcal mol
Jól működik!
diszkriminál és irányít
-1
200
250
1
2
3
4
5
6
7
8
9
10
MINEGAP folytonos, de
a deriváltjai NEM!
35
30
25
20
15
10
5
0
0
MINEGAP = min {ekR } - e P
300
21
Az EGAP-család új gyermeke: EWEGAP
E
E
…
l
k
EGAP
…
l
k
EWEGAP
j
i
j
i
1 n R P
EGAP = åek - e
n k=1
n
EWEGAP =
åe
k=1
n
j=1
MINEGAP = min {ekR } - e P
- bekR R
k
åe
MINEGAP
lim ( EWEGAP) = EGAP
e
- be Rj
b ®0
-eP
lim ( EWEGAP) = MINEGAP
b ®¥
22
EWEGAP:  megválasztása
n
EWEGAP =
åe
k=1
n
- bekR R
k
åe
e
- be Rj
j=1
 ~ 1.0 mol kcal-1
jó választás
-eP
23
Köszönetnyilvánítás
Fuxreiter Mónika
Petr Kulhánek
Simon István
Csányi Gábor
Enzimológiai Intézet elméleti
csoportjának munkatársai
Noam Bernstein
MTA, SZBK Enzimológiai Intézet
University of Cambridge
Alessandro Laio
SISSA, Trieszt, Olaszország