Transcript Introducción a la Dinámica Molecular

M.Orozco
J.L.Gelpí
M.Rueda
J.R.Blas
Clase 2: Introducción
a la Dinámica
Molecular
The long range goal of molecular
approaches to biology is to describe
living systems in terms of chemistry
and physics
Prof. Dirac. Proc.Roy.Soc., 123,
714 (1929)
Dinámica molecular
• Técnica para obtener una visión en el
tiempo de un sistema molecular.
– Permite seguir sistemas que evolucionan en el
tiempo (transiciones)
– Permite guiar sistemas a zonas de energía
óptima
– Permite obtener visiones promediadas de
sistemas en equilibrio: Promedios de
Boltzman (< >)
Que es la dinámica molecular
• Técnica para muestrear
espacio conformacional
en zonas delimitadas
• Técnica para representar
interacciones moleculares
• Herramienta para seguir
transiciones
conformacionales.
• Técnica que permite
analizar el efecto de
perturbaciones en el
sistema
• Una caja negra
• Una técnica que permite
•
•
encontrar transiciones
dramáticas
Una técnica competitiva
con NMR o X-Ray en
resolución de estructura
Una herramienta de uso
universal
No olvidar!
• Experimento mal  No hay resultado
• Experimento correcto  Resultado correcto
• Cálculo MD mal  Hay resultado!
• Cálculo MD mal  Resultado incorrecto
• Cálculo MD correcto  ?
Dinámica molecular: Ideas
antiguas, métodos nuevos
Apuntes históricos
Isaac Newton desarrolló el
esqueleto formal de la
dinámica molecular
siglo XVIII
Apuntes históricos (2)
• Primer intento de trayectoria H+H2H2+H: Hirschfelder
1936. (J.A. Hirschfelder, H.Eyring and B.Topley. J.Chem.Phys, 4,
170 (1936)).
• Trayectoria líquidos ideales Alder & Wainright 1959.
(B.J.Alder and T.E.Wainwright, J.Chem.Phys, 31, 459 (1959)).
• Trayectoria argón liquida Rahman 1964. (A.Rahman,
Phys.Rev., A136, 405 (1964)).
• Trayectoria agua líquida Stillinger & Rahman 1974.
(F.H.Stillinger and A.Rahman, J.Chem.Phys., 60,1545 (1974).
• Trayectoria de la primera proteína (10 ps BPTI).
McCammon & Karplus 1977. (J.A.McCammon, B.R.Gelin and
M.Karplus, Nature, 267, 585 (1976))
Apuntes históricos (3)
• Van Gunsteren & Berendsen introducen constrains en
dinámica proteínas (BPTI) 1977. (W.F.van Gunsteren, and
H.J.C.Hermans, Molecular Physics, 34, 1311 (1977)).
• Introducción del efecto del entorno van Gunsteren (20 ps
BPTI) 1983. (W.F.van Gunsteren, H.J.C.Berendsen, J.Hermans and
W.G.J. Postma PNAS, 80, 4315 (1983)).
• Desarrollo de los primeros programas de MD (CHARMM,
•
GROMOS) 1976-1978.
Berendsen desarrolla métodos para mantener P y T
constante 1984. (H.J.C.Berendsen, J.P.M.Postma and W.F.van
Gunsteren, J.Chem.Phys., 81,3684 (1984).
• Primer cálculo MD/FEP. McCammon 1985. (T.P.Lybrand,
I.Ghosh and J.A.McCammon, J.Am.Chem.Soc, 107, 7793 (1985))
Dinámica molecular
• Su base física es la dinámica de Newton
en cualquiera de sus formulaciones.

E i
fi   
r i


fi  miai


v i   a i dt

ri 

 v i dt
El Hamiltoniano
(La función de energía potencial V)
• Puramente cuántico:
– H(x)= H(R,r)
– Los electrones están en equilibrio. Solo se
integra explícitamente el movimiento de
núcleos (BO).
– Se integran los movimientos de núcleos y
electrones.
Dinámica de Carr-Parrinello
• Se integra el movimiento de núcleos y de
•
•
•
•
“orbitales”.
Funciona con “ondas planas”.
Permite estudiar sistemas reactivos, efectos de
transferencia de carga, conductividad,...
Simulan procesos en escalas del picosegundo
Su uso en biología es restringido.
Acc. Chem. Res., 35, 455, 2002
El Hamiltoniano:
• Mixto QM/MM:
– H(x)= H(Rext,r) + H(Rext)+ HQM/MM(Rint,Rext)
– Se integran los movimientos de núcleos.
– Los electrones están en equilibrio. Solo se
integra el movimiento de núcleos (BO).
H
H
H
O
N
N
H
H
N
N
N
R
N
O
N
R
O
N
R
Colominas et al., JACS, 118, 6811, 1996
O
Métodos QM/MM
• Se acoplan a cálculos de MD o métodos de
•
•
Monte Carlo.
Permiten estudiar reactividad en
macromoléculas donde se pueda trazar la
frontera entre la parte clásica y la cuántica.
Han tenido un uso moderado:
– La partición QM/MM no es siempre fácil
– Lentos
– El cálculo QM de referencia es bajo y a menudo
requiere reparametrización
– El acoplamiento QM/MM no es trivial
El Hamiltoniano:
• Clásico:
– H(x)= H(R)
– Solo se integran los movimientos
nucleares.
– Los electrones están en equilibrio y su
efecto sobre el sistema se introduce
paramétricamente: Force-Field
La MD permite a priori simular
cualquier sistema de impacto
biológico, aún los no descritos
estructuralmente e inabordables
por cualquier otra técnica
(teórica o experimental)
La mayoría de los datos que se obtienen
en proyectos de genómica y proteómica
se realizan en medios muy diferentes
de los fisiológicos
Tienen estos resultados valor?
MS-ELECTROSPRAY
Determination of mass
for large macromolecules
with high accuracy
Does it provide structural
information?
“Although DNA duplex ions have been observed for
seven years, the nature of interaction between
strands in the gas phase remains a major unanswered
question”
S.A.Hofstadler and R.H. Griffey
Chem. Rev. 101, 377 (2001)
16-mer LC 298
16-mer LC 448
12-mer LC 298
16-mer water
16-mer DC 298
12-mer LC 448
12-mer water
16-mer DC 448
12-mer DC 298
12-mer DC 448
PROBLEMAS INTRÍNSECOS DE
LA DINÁMICA
MOLECULAR (Clásica)
• El force-field clásico está limitado a estudiar
•
•
•
•
cambios conformacionales e interacciones
débiles.
Típicamente la información estructural de
partida debe ser de alta calidad.
El sistema de simulación se suele limitar a
unos miles o decenas de miles de átomos.
Escala de simulación limitada (state of the art
10 ns, sistemas modelos hasta 1 ms). Ello
conduce a problemas de muestreos
insuficientes.
El cálculo inunda de información difícil de
desencriptar (data mining).
E(AMBER; kcal/mol)
E(stab) DNA pairs
35
30
25
20
15
10
5
0
-5
-5 0
5
10
15
20
25
30
35
E(MP2; kcal/mol)
Hobza et al., JCC,18,1136, 1997
MUESTREO INCOMPLETO
Mejorar general de la
capacidad de muestreo
• Aumento de temperatura
• Suavización del force-field
• Introducir “penalizaciones de visita”
– Van Gunsteren J.Comput., 8, 695, 1994
– Parrinello PNAS, 99, 12562, 2002
• Cálculos paralelos (M ns  n·m ns)
– Pande JMB., 323, 927, 2002
 ( t )  1  exp(  kt )
Mejorar de la capacidad de
muestreo en una dirección
• Modificación función potencial. Brooks,
•
•
•
•
Ann.Rev. Phys.Chem., 52, 499, 2001.
Multiple-copies. LES Elber & Karplus JACS, 112,
9161, 1990
Método intercambio de réplicas a distintas
temperaturas. Berne PNAS, 98, 14931, 2001.
Sesgo por replica seleccionada (M-evil)
Introducción de fuerzas (steered-MD),
simulando experimento AFM. Schulten PNAS,
99, 1351, 1999. Karplus PNAS, 97, 6521, 1997.
MEJORA DEL MUESTREO
ASOCIADO A CÁLCULOS DE
ENERGÍA LIBRE
CÁLCULO DE ENERGÍA LIBRE
• Una trayectoria de equilibrio permite
calcular
 G AB   k b T ln
 N B
 N A
• Alternativamente la Mec. Estadística nos
dice:
 G AB   
1
ln
 exp(  E   B
 exp(  E   A
PROBLEMAS CALCULO G DIRECTO
• Aunque en algún caso es aplicable en
general es difícil
– Problemas para tener muestreo simultáneo de
los dos estados A y B
– Problemas para construir la función partición
global para A y B.
• Es necesario sesgar la trayectoria para
mejorar el muestreo y buscar fórmulas
alternativas:
MD/TI y MD/FEP calculations
1  
 G   RT

 0
ln  exp(  E 
 
RT
 E )
1  
G    
 0
E
 d

OTRAS TÉCNICAS RELACIONADAS
• Umbrella sampling. Torrie&Valleau
J.Comp.Phys., 23, 187, 1977
• Adaptative umbrella sampling. Karplus & cow.
J.Chem.Phys., 111, 8048, 1999.
• Weigthed hystogram analysis (WHAM).
Rosenberg J.Comp.Chem., 13, 1011, 1992.
ALGUNAS IDEAS RECIENTES
• Steered-MD acoplada a cálculos energía
libre.
– Hummer-Szabo PNAS, 98, 3658, 2001.
(EFM)
– Schulten Science, 296, 525, 2002.
• Método de penalizaciones de visita
acoplado a cálculos de energía libre
– Laio-Parrinello. PNAS, 99, 12562, 2002
OTRAS APROXIMACIONES
PARA MEJORAR EL
MUESTREO
Mejorar la velocidad de cálculo
• Uno de los cambios mayores en la
aplicación útil de la MD va a venir de la
mejora en los ordenadores:
– Procesadores más rápidos
– Procesadores dirigidos únicamente a MD
– Mejora en los algoritmos de paralelización
– Aumento capacidad de almacenamiento
En nuestra experiencia,...
• En
•
•
•
Tirado et al., Biochemistry, 36,7313, 1997 presentamos
trayectorias de 2 ns de una proteína < 100 residuos. En
nuestra mejor máquina recogíamos 5 ps x día.
En Cubero et al., JACS, 121, 8653, 1999 presentamos 1
trayectoria de 5 ns de un 12-mer (PME). En nuestra mejor
máquina recogíamos 25 ps x día.
Nuestro standard actual de trabajo son 2-6 simulaciones de
10-15 ns sobre DNAs o proteínas (PME). En nuestra mejor
máquina recogemos > 100 ps x día.
Nuestras máquinas más rápidas eran WS (10M<precio<5M),
ahora son AMD(2 GHz) (0.5M<precio<0.3M).
90’s
100 ps
95’s
500 ps -1 ns
02’s
10-25 ns
00’s
5 ns
05’s
100 ns – 1 ms
Almacenamiento
• Las trayectorias generan una gran cantidad de
datos, pero la densidad de información es
pobre:
– Problemas de disco.
– El coste de los post-análisis puede ser mayor que el
costo de la simulación.
– El “data minning” es muy complejo.
EL PROBLEMA DEL ANÁLISIS
DE LAS TRAYECTORIAS
Buscar una aguja en billón de
pajares
Análisis de las trayectorias
• Geometría/conformación
• Determinación pautas de interacción con
solvente.
• Energética del sistema
• Cálculos de energía libre
• Cálculos de flexibilidad y entropía
Análisis de las trayectorias
• Datos de geometría del sistema.
– RMS, Rg, size-descriptors, Y, F,..., secondary
structure, helical analysis, interaction maps,
H-bond population,...
• Datos interacción sistema-solvente.
– Solvent maps, rdf, solvent-population,...
Surco reactivo
Shields et al., JACS, 119, 7463, 1997
SOLVATION MAPS
5'-CCAA(G·C)CTTGG-3'
5'-CCAA(G·T)CTTGG-3’
5'-CCAA(O·C)CTTGG-3'
5'-CCAA(O·T)CTTGG-3'
Hernández et al., NAR 28, 4873, 2000.
Analysis of the trajectories
• Flexibility
– Principal component analysis Essential
dynamics (Berendsen)
• Comparison of frequencies
• Comparison of eigenvalues
– Entropy calculations  Karplus-van Gunsteren
• Schlitter’s method (van Gunsteren)
• Extrapolation to infinite (Laughton-Orozco)
Principal component analysis(PCA)
Trayectoria
Matriz
covarianza
Diagonalización
Val. propios
frecuencias
Vect. propios
Nat. movimiento
Naturaleza eigenvector:
 AB   i  
Comparación de eigenvectors:
A

Para un conjunto:
Medida relativa:
AB


1
n
AB
n
n
 
j 1
 2
( i  
A
1
2
3
4
5
6
7
8
9
10
2
0,03
0,51
0,35
0,61
0,11
0,28
0,01
0,02
0,19
0,01
3
0,27
0,52
0,30
0,49
0,21
0,09
0,14
0,03
0,17
0,04
4
0,09
0,06
0,41
0,33
0,09
0,57
0,08
0,13
0,48
0,04
5
0,10
0,10
0,62
0,29
0,15
0,28
0,40
0,12
0,26
0,14
0.85
0.62
0.81
0.81
0.14
6
0,05
0,18
0,10
0,08
0,14
0,10
0,73
0,16
0,07
0,05
7
0,17
0,23
0,23
0,02
0,48
0,22
0,07
0,45
0,30
0,12
8
0,14
0,26
0,04
0,21
0,03
0,26
0,24
0,18
0,29
0,31
9
0,16
0,26
0,13
0,04
0,21
0,04
0,29
0,15
0,03
0,10
10
0,04
0,03
0,10
0,07
0,09
0,01
0,04
0,22
0,17
0,62
B
j
)
i 1

(
T
AA
AB

5’GATTAATTAATTAATC3’
1
0,87
0,28
0,22
0,00
0,24
0,06
0,14
0,05
0,00
0,03
B
j
0.94
0.77
0.72
0.30
0.51
 =2.72/5= 0.54
T
BB
)
2
Entropías configuracionales
• Karplus:
– S = 0.5 k ln(det s)
– s: covariance matrix of set of “important”
coordinates
• Schlitter:
– S = 0.5 k S ln[1 + (kTe2/h2) <qi2>]
– <qi2>: eigenvalues of diagonalized mass-weighted
cartesian coordinate covariance matrix
Entropías configuracionales
• Los valores de S calculados son sensibles al
En t r o p y K C al M o l
tiempo de simulación
775
750
Free
725
1:1 Complex
700
2:1 Complex
675
650
625
1000
2000
3000
Length of SamplingWindowps
4000
• Ajuste de una función 1/x para estimar S para
una simulación infinita S=So-a/t (Harris et
al., 2001)
Cálculos de energía libre total
MD-PB/SA (Kollman)
Gtot = Eintra + Gsolv – T Sintra
GA-Btot = (EAintra – EBintra) + (GAsolv – GBsolv)
Gtot.
Gsolv.
Eintra.
Total free energy
Solvation free energy
Intramolecular energy
Cubero et al., JACS. 123, 12018. 2001
GB/SA o PB/SA
Calc. solvatación
<X>
Cálculo directo
<Eint>
Trayectoria de
equilibrio (MD)
Cálculo S a partir
de Slitter
<X>
Parallel duplexes
Poly d(A·T)
• 5,6,7,8,9,10,11,12,15-mer
parallel duplexes
 Reverse WC *
 Hoogsteen
• 5,7,9,11,15-mer antiparallel
duplexes
* Pattabiraman, N. Biopolymers, 25, 1603 (1986)
Cubero et al., JACS. 123, 12018. 2001
Cubero et al., JACS. 124, 3133, 2002
Total free energy (kcal/mol)
WC
G= -343.2(0.5)x +84.6(4.9); R2=1.0000
rWC G= -342.0(0.4)x +81.3(2.9); R2=1.0000
H
G= -341.9(0.2)x +93.5(1.7); R2=1.0000
95
Estabilization
WC
85
E
-343.2(0.5)
rWC
-342.0(0.4)
H
-341.9(0.2)
Enucleation
75
65
0
1
2
3
4
5
(n-2) residues
6
7
(N-2) residues
8
9
10
11
12
13
WC
84.6(4.9)
rWC
81.3(2.9)
H
93.5(1.7)
Cubero et al., JACS. 123, 12018. 2001
EJEMPLOS DE APLICACIONES
RECIENTES DE LA DINÁMICA
MOLECULAR
Campos de aplicación más
activos
• DNA MD puede ser tan potente como las
•
•
•
•
técnicas experimentales.
Estudio de flexibilidad conformacional en
proteínas y ácidos nucleicos.
Folding/unfolding de péptidos/proteínas modelo.
Estudios de termodinámica de unión y reacción.
Dinámicas en condiciones límites: grandes
sistemas o grandes escalas temporales.
Dinámicas en grandes sistemas
• Se estudian sistemas > 100000 átomos con
•
•
•
•
condiciones de simulación realistas.
Típicamente sistemas de alto impacto
biológico (bombas, canales, receptores,...)
Trayectorias en el rango 2-10 ns
Requieren recursos computacional enormes
En algunos casos se utilizan “perturbaciones”
o dinámicas sesgadas para simular procesos
en escalas de tiempo ms.
Dinámicas en grandes sistemas (1)
• Simulaciones de aquaporinas (PNAS, 99,
6731, 2002; Science, 296, 525, 2002)
– Estudian el mecanismo de paso de glicerol y
agua por estos canales.
– Sistema proteico con membrana y solvente
(106000 átomos).
– Trayectorias libres y “forzadas” de 7-12 ns
Dinámicas en grandes sistemas (2)
• Simulaciones de F1-ATP sintasa (Nature
Struct. Biol., 9, 198, 2002)
– Mecanismo de transferencia de energía por
cambio conformacional.
– Sistema con 260000 átomos (proteína +
agua)
– Trayectorias “forzadas” de 7-ns
– 280000 horas CPU Power 3 (32 años CPU)
Dinámicas en grandes sistemas (3)
• Simulaciones de unfolding de
fibronectina FN-III. Gao et al., JMB, 323,
939, 2002.
– Simulaciones de steered molecular
dynamics
– Simulaciones 10 x 12 ns
– Sistema con 126000 átomos (proteína +
agua)
– 652000 horas CPU Athlon 1.3 GHz (74 años)
Dinámicas en grandes escalas
temporales
• Simulaciones que cubren escalas cercanas
al microsegundo con solvente explicito.
• Típicamente son dinámicas no sesgadas.
• Normalmente se enmarcan en estudios de
folding de proteínas.
• Son de utilidad dudosa como técnica
general “probe of concept”.
Dinámicas en grandes escalas
temporales (ejemplos)
• Simulaciones de folding reversible de péptidos
(20-200 ns) (ej Daura et al., Angew. Chem., 38,
236, 1999)
– Coleccionan “eventos” de folding reversible de
pequeños péptidos  Gfolding
• Simulaciones de 1 ms de “folding” de una
proteína modelo (Duan-Kollman: Science, 277,
1793, 1998)
– Consiguen representar usando “primeros principios” el
inicio del folding de una pequeña proteína
Dinámicas en pseudo-grandes
escalas temporales (ejemplos)
• Simulaciones distribuidas de folding de
Villin (36-residues) Zagrovic et al., JMB,
323, 927, 2002.
– GB/SA (no solvente explicito).
– “Equivalente” a 300 ms.
– RMSd 1.7 Angs estructura experimental.
– 10000 procesadores: 2.5 x 1011 integraciones
– C.a 1000 años CPU (PIII-500 MHz)
y ahora las ecuaciones,....