BIOMOLEKULÁRNÍ SIMULACE
Proč 3D struktury? - “anatomie” of biomolekulárních systémů (interakce, konformace, katalysa atd.) - vývoj léčiv založený na struktuře - proteinové a enzymové inženýrství
Proč biomolekulární simulace - “4D struktury” - “moleculární mikroscopie” - vysvětlení strukturních změn, interpretace experimentů - predikce
BIOMOLEKULÁRNÍ SIMULACE -
stabilita proteinů v různých prostředích vliv kovalentních modifkací na strukturu proteinů sbalování proteinů validace předpovězených struktur proteinů interakce protein-ligand, protein-protein membránové struktury struktury dalších biopolymerů mechanika biopolymerů a jejich komplexů přenos signálu v biomolekulárních systémech virtuální experimenty servis experimentálních metod
Newtonovy pohybové zákony
2
mi
∂ ri ∂t
2
=Fi
∂V Fi =− ∂ ri
Simulace molekulové dynamiky 1957
“hard spheres”
-
Alder & Waiinwright
1964
argon
?
Rahman
1971
voda
2.2 ps Rahman & Stillinger
1977
BPTI
8 ps
1988
fosfolipidová dvojvrstva 200 ps Egberts & Berendsen
1993
biotin-streptavidin
108 ps Myiamoto & Kollman
1995
bacteriorhodopsin
300 ps Edholm et al.
1998
porin
1 ns
Tieleman & Berendsen
1998
sbalování
1 μs
Duan & Kollman
2011
sbalování
ms
Lindorff-Larsen et al.
2013
simulace virové kapsidy
ns
Zhao et al.
McCammon et al.
Software pro biomolekulární simulace zdarma GROMACS – http://www.gromacs.org levné AMBER - http://ambermd.org/ GROMOS - http://www.gromos.net/ CHARMM - http://www.charmm.org/ drahé
S-peptide demo
S-peptide GROMACS 1. počáteční souřadnice 2. topologie 3. instrukce pro program
S-peptide GROMACS převedení speptide.pdb do topologie a souřadnic, přidání vodíků $ pdb2gmx -f speptide -o speptide -p speptide + vyber správné silové pole vytvoření boxu s proteinem uprostřed $ editconf -f speptide -o box -c -d 1 naplnění boxu vodou $ genbox -cp box -cs -p speptide -o solvated přidání proti-iontů, pokud třeba
S-peptide GROMACS 1. počáteční souřadnice (speptide.gro) Go Rough, Oppose Many Angry Chinese Serial killers 286 1LYS N 1 2.497 -0.065 2.231 1LYS H1 2 2.581 -0.048 2.180 1LYS H2 3 2.519 -0.086 2.326 1LYS H3 4 2.448 -0.142 2.190 ... 19ALA C 284 2.846 3.022 2.056 19ALA OC1 285 2.919 3.015 1.954 19ALA OC2 286 2.713 3.020 2.055 1.79949 3.37953 1.37997
S-peptide (19 aminokyselin, 286 atomů, C86H140N27O32S, 1 Cl-, 859 H2O)
S-peptide GROMACS 2. topologie (speptide.top) 22 typů atomů 286 atomů 287 vazeb, 733 1–4 interakcí 513 valenčních úhlů 798 torsí + topologie vody a iontů
S-peptide GROMACS 3. instrukce pro program (md.mdp) integrator = constraints = constraint_algorithm = dt = nsteps = nstcomm = nstxout = nstvout = nstfout = nstlog = nstenergy = nstlist = ns_type = coulombtype = rlist = rcoulomb = rvdw =
použij molecular dynamics
md fixované délky vazeb all-bonds lincs 0.002 ; ps ! 500000 ; total 1 ns. 1 simulovaný čas 250 (500 000 krát 2 fs = 1 ns) 1000 0 100 100 frekvence ukládání 10 grid PME 1.0 nastavení modelování 1.0 nekovalentních interakcí 1.0
S-peptide GROMACS 3. instrukce pro program (speptide.top) ; Berendsen temperature coupling is on in two groups Tcoupl = berendsen tc-grps = Protein SOL udržování teploty tau_t = 0.1 0.1 ref_t = 300 300 ; Energy monitoring energygrps = Protein SOL ; Isotropic pressure coupling is now on Pcoupl = berendsen Pcoupltype = isotropic tau_p = 0.5 udržování tlaku compressibility = 4.5e-5 ref_p = 1.0 ; Generate velocites is off at 300 K. gen_vel = no gen_temp = 300.0 teplota při startu gen_seed = 173529
S-peptide minimizace energie $ grompp -f em -c solvated -p speptide -o em1 $ mdrun -s em1 -o em1 -e em1 -g em1 -c after_em1 simulace molekulární dynamiky $ grompp -f md -c after_em1 -p speptide -o md1 $ mdrun -s md1 -o md1 -e md1 -g md1 -c after_md1 Step Time Lambda 2800 5.60000 0.00000 Rel. Constraint Deviation: Max between atoms RMS Before LINCS 0.058424 247 248 0.007393 After LINCS 0.000082 180 182 0.000029 Energies (kJ/mol) Angle Proper Dih. RyckaertBell. LJ14 Coulomb14 6.22621e+02 5.32697e+01 7.29416e+02 2.94892e+02 3.86087e+03 LJ (SR) Coulomb (SR) Potential Kinetic En. Total Energy 4.62848e+03 4.71919e+04 3.70024e+04 7.47176e+03 2.95306e+04 Temperature Pressure (bar) 3.14265e+02 2.02309e+02
S-peptide
NODE (s) Real (s) (%) Time: 573.400 580.000 98.9 9:33 (Mnbf/s) (GFlops) (ns/day) (hour/ns) Performance: 11.327 1.597 15.068 1.593 Finished mdrun on node 0 Sun Sep 20 11:21:17 2011
Příklady
V. Spiwok, P. Lipovová, T. Skálová, J. Dušková, J. Dohnálek, J. Hašek, N.J. Russell, B. Králová: J. Mol. Model. (2007) 13:485-497.
Příklady
V. Spiwok, P. Lipovová, T. Skálová, J. Dušková, J. Dohnálek, J. Hašek, N.J. Russell, B. Králová: J. Mol. Model. (2007) 13:485-497.
Příklady
V. Spiwok, P. Lipovová, T. Skálová, J. Dušková, J. Dohnálek, J. Hašek, N.J. Russell, B. Králová: J. Mol. Model. (2007) 13:485-497.
Uvnitř simulace molekulární dynamiky
Newtonovy pohybové zákony
síly
2
hmotnosti
mi
∂ ri ∂t
2
=Fi
∂V Fi =− ∂ ri
potenciální energie
Simulace molekulární dynamiky vs Minimalizace energie
r
r
Potenciální energie
“dobrá” struktura nízká energie
“špatná” struktura vysoká energie
Silová pole Proteiny, nukleové kyseliny, lipidy: AMBER GROMOS OPLS CHARMM obecné molekuly: GAFF MM2 MM3 MMFF speciální: Glycam (sacharidy) Martini (coarse grained)
Force fields
all atom
united atom
coarse grained
Silové pole
1 1 2 2 V = ∑ k r r −r 0 ∑ k −0 bonds 2 angles 2 ∑ k 1cos n −s torsions
∑
pairs
[
12 ij 12 ij
qi q j C 1 4 0 r r ij r
−
C
6 ij 6 ij
r
]
Silové pole
1 1 2 2 V = ∑ k r r −r 0 ∑ k −0 bonds 2 angles 2 ∑ k 1cos n −s torsions
∑
pairs
[
12 ij 12 ij
qi q j C 1 4 0 r r ij r bonds
−
C
6 ij 6 ij
r
]
Harmonický potenciál potenciál:
1 2 V = k r r−r 0 2 síla:
r
∥F∥=−k r r −r 0 r0
r
Silové pole
1 1 2 2 V = ∑ k r r −r 0 ∑ k −0 bonds 2 angles 2 ∑ k 1cos n −s torsions
∑
pairs
[
12 ij 12 ij
qi q j C 1 4 0 r r ij r valence angles
−
C
6 ij 6 ij
r
]
Silové pole
1 1 2 2 V = ∑ k r r −r 0 ∑ k −0 bonds 2 angles 2 ∑ k 1cos n −s torsions
∑
pairs
[
12 ij 12 ij
qi q j C 1 4 0 r r ij r torsions
−
C
6 ij 6 ij
r
]
Torse
vlastní torse
nevlastní torse
Silové pole
1 1 2 2 V = ∑ k r r −r 0 ∑ k −0 bonds 2 angles 2 ∑ k 1cos n −s torsions
∑
pairs
[
12 ij 12 ij
qi q j C 1 4 0 r r ij r
non-covalent electrostatic interactions
−
C
6 ij 6 ij
r
]
Parciální náboje
Silové pole
1 1 2 2 V = ∑ k r r −r 0 ∑ k −0 bonds 2 angles 2 ∑ k 1cos n −s torsions
∑
pairs
[
12 ij 12 ij
qi q j C 1 4 0 r r ij r
non-covalent van der Waals interactions
−
C
6 ij 6 ij
r
]
Lennardův-Jonesův potenciál
12
6
C C V = 12 − 6 r r
r
rvdW r0
1–2, 1–3, 1–4 interakce
1
4 2
3
1–2 – pouze kovalentní 1–3 – pouze kovalentní 1–4 – nekovalentní interakce sníženy 1–5, 1–6 atd. – nekovalentní interakce jako obvykle
[ atomtypes ] ;name bond_type mass charge ptype sigma epsilon opls_111 OW 8 15.99940 0.834 A 3.15061e01 6.36386e01 opls_112 HW 1 1.00800 0.417 A 0.00000e+00 0.00000e+00 [ moleculetype ] ; molname nrexcl SOL 2 [ atoms ] ; id at type res nr residu name at name cg nr charge 1 opls_111 1 SOL OW 1 0.834 2 opls_112 1 SOL HW1 1 0.417 3 opls_112 1 SOL HW2 1 0.417 [ bonds ] ; i j funct length force.c. 1 2 1 0.09572 502416.0 1 3 1 0.09572 502416.0 [ angles ] ; i j k funct angle force.c. 2 1 3 1 104.52 628.02
Voda (TIP3P model) [ atomtypes ] ;name bond_type mass charge ptype sigma epsilon opls_111 OW 8 15.99940 0.834 A 3.15061e01 6.36386e01 opls_112 HW 1 1.00800 0.417 A 0.00000e+00 0.00000e+00 [ moleculetype ] ; molname nrexcl SOL 2 [ atoms ] ; id at type res nr residu name at name cg nr charge 1 opls_111 1 SOL OW 1 0.834 2 opls_112 1 SOL HW1 1 0.417 3 opls_112 1 SOL HW2 1 0.417 [ bonds ] ; i j funct length force.c. 1 2 1 0.09572 502416.0 1 3 1 0.09572 502416.0 [ angles ] ; i j k funct angle force.c. 2 1 3 1 104.52 628.02
Jak ziskat (chybějící) parametry silového pole? 1. Jiná silová pole (s opatrností) 2. Experiment IČ, krystalografie, ... 3. Molekulární modelování kvantová chemie
Chemické reakce - kvantová chemie - kombinace moleculární a kvantové mechaniky (QM/MM) - molekulární mechanika “trénovaná” kvantovou chemií (empirical valence bond) - speciální “reaktivní” silová pole
QM/MM
M. Krupička, I. Tvaroška: J Phys Chem B (2009) 113(32):11314-11319.
Constraints
Periodická okrajová podmínka
Udržování teploty Termostat Berendsenův, Noseův-Hooverův, V-rescale Počáteční rychlost Barostat Berendsenův, Parrinellův-Rahmanův Udržování povrchového napětí
Analýza výstupů Časový vývoj strukturních parameterů: - energie, teplota - vzdálenosti, úhly, torse - počet nativních kontaktů - sekundární struktura - radius of gyration - root mean square deviations (RMSD) RMSD
RMSD
time
time
Analýza výstupů - root mean square fluctuations (RMSF)
RMSF
residue number
Analýza výstupů - esenciální dynamka trajektorie
molekula CV2
kolektivní pohyby CV1
Speciální otázky 1. voda 2. proteiny 3. nukleové kyseliny 4. sacharidy 5. jiné molekuly
Voda Proč je voda důležitá?
+
+
Voda TIP (Transferable intermolecular potential) TIP3P Vyladěno tak, aby přesně modelovaly experimentální parametry TIP4P
TIP5P
- radiální distribuční funkci - self-diffusion coefficient - tepelná kapacita - bod varu a tání - dielektrická konstanta ...
Force field customization
5
6
55 56 1 - 50
10-391
54 57
51 52
60 59 58
53
Remove atoms 59 and 60 Remove bonds 57–59 and 57–60 Remove angles 54-57-59, 54-57-60 58-57-59, 58-57-60, 59-57-60 Remove torsions and 1–4 interactions involving 59, 60 55 56
4
1
7
2
8 9 3 Remove atom 1 Remove bond 1–2 Remove angles 1-2-3, 1-2-4 Remove torsions 1-2-4-5, 1-2-4-6, 1-2-4-7, Remove 1–4 interactions 1-5, 1-6 and 1-7, renumber 62 63 67-448
1 - 50
61
54 51 52
53
57
59
58
60
64
65
66
Force field customization 55 56
62 63 67-448
1 - 50
54 57
51 52
61 59
53
58
Add bond 57–59 Modify bond 57–58 Add angle 54-57-59, 58-57-59, 57-59-60, 57-59-61 Modify angle 54-57-58 Add torsion 51-54-57-59, 55-54-57-59, 56-54-57-59, 54-57-59-60, 54-57-59-61, 58-57-59-60, 58-57-59-61, 57-59-61-62, 57-59-61-63, 57-59-61-64 Modify torsion 51-54-57-58, 55-54-57-58, 56-54-57-58, 60-59-61-62, 60-59-61-63, 60-59-61-64 Add 1–4 interactions 51-59, 55-59, 56-59, 54-60, 54-61, 58-60, 58-61, 57-62, 57-63, 57-64
60
64
65
66
Force field customization
62 63
60
55 56
67-448 1 - 50
59
54 57
51 52
53
58
61 64
65
66
Vzorkování
Vzorkování A
B
Vzorkování A
Vpot,A porovnání Vpot je (většinou) k ničemu: - spousta stupňů volnosti - voda - teplota, entropie
B
Vpot,B
Vzorkování A
B
A
B
time
Vzorkování A
B
A
B
time můžeme simulovat
Počítače klastry
Superpočítače http://www.top500.org
National Super Computer Center in Guangzhou 3 120 000 jader
Výpočetní centra v ČR Metacentrum http://www.metacentrum.cz/ CERIT-SC http://www.cerit-sc.cz/ CESNET http://www.cesnet.cz/ IT4Innovation http://www.it4i.cz/
Speciální hardware
Příklady simulací – 2-adrenergní receptor
Dror et al. (2009) Proc Natl Acad Sci USA, 106, 4689–4694
Distribuované výpočty http://folding.stanford.edu/
Posílené vzorkování Metadynamika, Umbrella sampling, ... A
B
čas Můžeme simulovat
Metadynamika
Spiwok & Tvaroška (2009) J Phys Chem B, 113, 9589–9594
Metadynamika
Spiwok et al. (2015) J Chem Phys, 113, 9589–9594
Paralelní temperování
12
teplota
10
8
6
4
2
0 0.00
2.00
4.00
6.00
8.00
10.00
12.00
čas Nguyen et al. (2003) Proteins, 61, 795–808