MD szimuláció

A generált modellszerkezeteken, azaz az apo hPGK-n és négy terner komplexén (D-ADP-, L-(D-ADP-, D-CDP- és L-CDP komplex) MD szimulációt hajtottunk végre. A szimulációk futtatásához az NAMD programcsomagot (129), és annak all-atom típusú CHARMM 22-es erőterét (130) használtuk. A szimulációkat egy SGI Altix 350 típusú szerver nyolc darab processzorán futtattuk. A szimulációk előkészítését és a trajektóriák analízisét a CHARMM programcsomaggal (26) végeztük. A MD futásokhoz az NAMD programot használtuk, mivel parallel környezetben az NAMD hatékonyabb, mint a CHARMM.

MD szimuláció előtt az egyes modellszerkezetek potenciális energiáját minimalizáltuk, hogy eltávolítsuk az esetlegesen jelen lévő, energetikailag kedvezőtlen kölcsönhatásokat, amik például a fiziológiailag nem természetes rácserőkből adódhatnak. A következő minimalizálási protokollt követtük mind az öt rendszer esetében. 1500 lépés steepest descent módszerrel (131) kezdtük az energia-minimalizációt. Ennek során a H-atomok kivételével minden atomra harmonikus kényszereket róttunk ki. Ezzel megakadályoztuk, hogy a nehéz atomok túl nagy amplitúdójú mozgásokat végezzenek. Így az energia-minimalizáció stabilabb, kisebb a valószínűsége a nagy fluktuációknak az energiában. Az alkalmazott harmonikus kényszer potenciális energiája:





c i i

nehéz atomok

E 



ahol a H-atomok kivételével az összes atomra történik az összegzés, K az erőállandó, r_i az i-edik atom helyvektora, r_i0 az i-edik atom referencia helyvektora. A referencia helyvektorokat mindig a kényszerek kirovása pillanatában érvényes konfiguráció adta.

A minimalizálás során fokozatosan csökkentettük a kényszererők nagyságát. A steepest descent első 500 lépésében 10 kcal/mol/Ų , a következő 500 lépésben 1 kcal/mol/Ų , míg az utolsó 500 lépésben 0,1 kcal/mol/Ų volt a K értéke. 0,1 kcal/mol/Ų értéken tartva az erőállandót, 200 lépés conjugate gradient módszerrel (131) folytattuk a minimalizálást. Majd a kényszereket megszüntettük, és 100 lépés steepest descent módszerrel tovább minimalizáltuk a rendszerek potenciális energiáját. Ezt 200 lépés conjugate gradient követte. Majd 1000 lépés adopted-basis Newton-Raphson módszerrel (131) fejeztük be a minimalizációs protokollt.

Ezt követően minden rendszert 99 Å X 74 Å X 63 Å térfogatú vízdobozba helyeztünk.

Ezzel definiáltuk a szimuláció elemi celláját. A vízdoboz felépítéséhez explicit vízmolekulákat, azaz a CHARMM TIP3 vízmodelljét (26) használtuk. Az így kapott elemi cellára periodikus határfeltételeket róttunk ki, azaz a vízdoboz szemközti oldalait azonosítottuk egymással. Tehát ha a szimuláció során egy vízmolekula kilépett a vízdoboz egyik oldalán, akkor ezzel egy időben ugyanez a vízmolekula megjelent a szemben lévő oldalon. A periodikus határfeltételek segítségével megakadályozhatjuk, hogy a fehérje környezetéből eldiffundáljon a víz. Egyszerű eltolási transzformációk bevezetésével az elemi cellából ortorombos kristályrácsot konstruáltunk. Ezzel nagy mennyiségű víz környezet válik modellezhetővé a fehérje körül. Úgy választottuk meg a vízdoboz méreteit, hogy az a fehérje felületétől számítva x, y és z irányban 12 Å vastagságú legyen. Ekkora vízdoboz vastagsággal elkerülhetők a szomszédos elemi cellák közti fehérje-fehérje nemkötő kölcsönhatások (pl. van der Waasls), ugyanis ezen kölcsönhatások hatótávolságát (cut-off távolság) 12 Å-nek választottuk. Véletlenszerű pozíciókba Na-ionokat vagy Cl-ionokat (rendszertől függően) helyeztünk a vízdobozba, hogy annak elektromos össztöltése zérus legyen. Az ilyen módon szolvatált rendszereket ismét energia-minimalizációnak vetettük alá, hogy az esetleges nagy (3)

energiájú fehérje-víz kölcsönhatásokat eltávolítsuk. Az energia-minimalizációt a fentebb tárgyalt vákuumbeli minimalizációval azonos módon végeztük.

A Particle Mesh Ewald (PME) (132,133) módszert használtuk az elektrosztatikus energiák számolására. A PME rácsállandót 1 Å-nek, az interpoláció rendjét 6-nak választottuk. A valós térbeli szummáció cut-off értékét 12 Å-nek, ennek megfelelően a PME Gauss függvény szélességét 0,34 Å^-1-nek állítottuk be. A dielektromos állandó értékét 1-nek vettük. A van der Waals kölcsönhatásokat ún. switching függvénnyel moduláltuk. Lényegében ez annyit jelent, hogy definiáltunk két cut-off értéket, egy belső és külső cut-offot. Ezen két érték között egy szigmoid függvény által, a távolság növekedésével fokozatosan zérusra redukáltuk a van der Waals energiát. Ezzel mintegy kikapcsoljuk a kölcsönhatást, innen a függvény elnevezése. Esetünkben a belső cut-off 10 Å, míg a külső cut-off értéke 12 Å volt.

A szolvatált rendszerek energia-minimalizálása után kezdetét vette maga a MD szimuláció. A MD szimuláció során a klasszikus mozgásegyenletek integrálása 1 femtoszekundumos időlépésekben történt. Az energia-minimalizált szerkezeteket (apo és a négy terner komplex) először 300 K hőmérsékletre, a szimulációs hőmérsékletre fűtöttök fel. Az eljárást a kezdeti sebességek kiosztásával kezdtük. A Maxwell-Boltzmann eloszlás alapján a rendszerek minden egyes atomjához véletlenszerűen sebességértékeket rendeltünk: újraosztásnál a hőmérséklet 10 K-nel nőjön, így 30000 lépés (30 pikoszekundum) után elértük a kívánt 300 K-t.

Miután felfűtöttük a rendszereket 300 K hőmérsékletre, ekvilibrációval folytattuk a szimulációt. Az ekvilibráció során az a cél, hogy a rendszert termodinamikai (4)

egyensúlyba hozzuk. Több termodinamikai mennyiség (pl.: energia, hőmérséklet, nyomás és térfogat) és az aktuális konformáció együttes vizsgálatával ellenőriztük, hogy mikor kerül egyensúlyba a rendszer. Elsőként kanonikus sokaságként kezelve a rendszereket, azaz a térfogatot és a hőmérsékletet állandó értéken tartva (N=áll., V=áll., T=áll.), hagytuk a rendszereinket a mozgás egyenletek szerint mozogni. A kanonikus ekvilibráció 50 pikoszekundumon (50000 időlépés) keresztül zajlott. Majd 500 pikoszekundum (500000 időlépés) NpT ekvilibrációval folytattuk (N=áll., p=áll., T=áll.) a szimulációt, azaz a Nosé-Hoover algoritmust (134,135) használva izobár-izoterm körülmények között vizsgáltuk a rendszert 300 K hőmérsékleten és 1 atm nyomáson.

Az ekvilibráció után kezdetét vette a konformációs adatgyűjtés-mintavételezés, az ún.

production run. 20 nanoszekundum izobár-izoterm (p=1 atm, T=300 K) MD szimulációt hajtottunk végre mind az öt rendszeren (apo és négy terner komplex). Az atomi koordinátákat 1 pikoszekundumonként mentettük, így mind az öt rendszerre megkaptuk az atomi mozgások trajektóriáit (hely-idő függvény).