Mechanika Molekularna

MECHANIKA MOLEKULARNA

Z chemii kwantowej wiemy, że energia całkowita układu (suma energii kinetycznej i potencjalnej) jest określona przez przybliżone rozwiązanie zależnego od czasu równania Schrödingera. Rozwiązanie takiego równania umożliwiło dopiero, ogłoszone rok później (1927), przybliżenie Borna-Oppenheimera. Przybliżenie zakłada, że ruch jąder można scharakteryzować energią potencjalną elektronów zależną od ustalonych położeń jąder. Czyli jądra atomowe znajdują się w potencjale pochodzącym od gęstości elektronowej i od innych jąder, przy czym potencjał elektronowy zależy oczywiście od położeń jąder. Założenie to uwzględnia fakt, że jądra poruszają się o wiele wolniej niż elektrony. Rozwiązaniem w tym przybliżeniu jest elektronowa funkcja falowa dla ustalonych położeń jąder. Przez wiele lat rozwiązanie takiego równania nawet dla bardzo małych cząsteczek było niemożliwe. Dopiero wprowadzenie komputerów umożliwiło symulacje zachowania prawdziwej cząsteczki na podstawie znanych nam praw chemii kwantowej. Niestety, tylko właściwości niewielkich cząsteczek można przewidywać przy użyciu wzorów mechaniki kwantowej, gdyż konieczne do tego typu przewidywań obliczenia wymagają ogromnych mocy obliczeniowych, nawet dla niewielkich cząsteczek. Dlatego dla bardzo dużych układów stosuje się techniki uproszczone, takie jak metody semiempiryczne czy metody mechaniki molekularnej.

Mechanika molekularna to metoda oparta na prostym mechanicznym modelu struktury cząsteczki. Taki model mechaniczny możemy uzyskać gdy w opisywanym wyżej przybliżeniu Borna-Oppenheimera zastąpimy potencjał pochodzący od gęstości elektronowej odpowiednio sparametryzowanym potencjałem klasycznym. Wtedy to cząsteczka jest traktowana jako zbiór atomów, który podlega zbiorowi klasycznych i mechanicznych funkcji potencjału.

Energia potencjalna jest zazwyczaj reprezentowana w postaci sumy łatwo interpretowalnych i weryfikowalnych eksperymentalnie członów: deformacji wiązań, deformacji kątów między wiązaniami, obrotów wokół wiązań oraz oddziaływań elektrostatycznych i van der Waalsa. W wielu wypadkach funkcja energii potencjalnej posiada dodatkowe człony równania, które dokładniej opisuję pewne właściwości molekuły (np.: oddziaływanie pochodzące od wiązań wodorowych czy człon opisujący jednoczesną deformację wiązań i kątów ). Parametryzacja funkcji energii potencjalnej jest wykonywana dla stosunkowo małego zbioru niewielkich molekuł zbudowanych z typowych fragmentów molekularnych i używana dla wszystkich układów zbudowanych z tych fragmentów. Takie podejście oczywiście zakłada, że atomy wbudowane w pewne otoczenie w jednej molekule powinny mieć takie same właściwości jak atomy w takim samym otoczeniu w innej molekule. Taką parametryzację nazywamy polem siłowym (ang. force field).

Funkcja energii potencjalnej może być parametryzowana na wiele sposobów. Różne parametryzacje są stosowane do różnych celów, występują również konkurencyjne parametryzacje dla tych samych zastosowań. W badaniach obiektów biologicznych stosowane są pola siłowe:

• AMBER (Kollman i współpracownicy, UCSF) [1-4]
• CHARMM (Karplus i współpracownicy, Harvard) [5-7]
• Gromos (Berendsen i współpracownicy, Nijmegen) [8, 9]
• CVFF i CFF9X (Accelrys Inc,) [10]
• Tripos (Tripos Inc.) [11]

Mechanika molekularna jest oparta na minimalizacji energii układu. W oparciu o fizykę statystyczną można przyjąć, że w niskich temperaturach, właściwości układu są określone przez właściwości stanu podstawowego - czyli stanu o najniższej energii. Niestety nie jest znany algorytm znajdowania globalnego minimum dla skomplikowanych funkcji wielu zmiennych, w związku z czym w praktyce wykorzystuje się algorytmy lokalnej minimalizacji, znajdujące najbliższe minimum lokalne. Najczęściej używane do tego celu algorytmy to algorytm najszybszego spadku, sprzężonych gradientów i Newtona-Raphsona.

Algorytm najszybszego spadku [12] polega na przesuwaniu atomów wzdłuż linii po hiperpowierzchni, w kierunku najszybciej malejącego gradientu (Rys. 1a). Kiedy funkcja przejdzie przez minimum, zmienia się gradient, a nowy jest prostopadły do poprzedniego, więc poszczególne kierunki minimalizacji będą prostopadłe względem siebie. Algorytm ten jest prosty i niezawodny, ale niestety w okolicy minimum mało efektywny. Algorytm sprzężonych gradientów [13-15] jest podobny, za wyjątkiem sposobu wybierania kierunku poszukiwania minimum. Kierunek gradientu został zastąpiony kierunkiem, który jest pewną kombinacja bieżącego i poprzednich gradientów (Rys. 1b). Algorytm ten jest szybszy od algorytmu najszybszego spadku.

Obydwa algorytmy są algorytmami używającymi jedynie informacji o pierwszych pochodnych. Dzięki temu są niezawodne ale i wolno zbieżne, zwłaszcza w pobliżu minimum. Algorytm Newtona wykorzystuje informacje o drugich pochodnych, dlatego jest szybko zbieżny w pobliżu minimum. Ale kiedy molekuła znajduje się daleko od minimum jest niestabilny. Wymaga on również obliczania i przechowywania macierzy drugich pochodnych rozmiaru NxN, co dla dużych molekuł jest niepraktyczne.

Dlatego też w praktyce używa się kombinacji tych trzech algorytmów. Najlepiej najpierw wykonać wstępną minimalizację za pomocą algorytmu najszybszego spadku, a następnie użyć algorytmu sprzężonych gradientów. Dla mniejszych układów, ewentualnie możemy zakończyć cykl minimalizacją z użyciem algorytmu Newtona dla szybkiego znalezienia minimum.

Rysunek 1. a) Ilustracja graficzna metody najszybszego spadku. b) Ilustracja graficzna metody sprzężonych gradientów.

Literatura

1. Ponder, J. W., and Case, D. A. (2003) Force fields for protein simulations. Adv Protein Chem 66, 27-85.
2. Cheatham, T. E., 3rd, and Young, M. A. (2000) Molecular dynamics simulation of nucleic acids: successes, limitations, and promise. Biopolymers 56, 232-56.
3. Cornell, W. D., Cieplak, P., Bayly, C. I., Gould, I. R., Merz, K. M. Jr., Ferguson, D. M. Spellmeyer, D. C., Fox, T., Caldwell, and J. W., Kollman, P. A. . (1995) A second generation force field for the simulation of proteins, nucleic acids and organic molecules. J. Am. Chem. Soc. 117, 5179-5197.
4. Pearlman, D. A., Case, D. A., Caldwell, J. C., Seibel, G. L., Singh, U. C., Weiner, P., & Kollman, P. A.,. (1991) AMBER 4.0. University of California, San Francisco.
5. MacKerell, A. D., Bashford, D., Bellott, M., Dunbrack, R. L., Evanseck, J. D., Field, M. J., Fischer, S., Gao, J., Guo, H., Ha, S., Joseph-McCarthy, D., Kuchnir, L., Kuczera, K., Lau, F. T. K., Mattos, C., Michnick, S., Ngo, T., Nguyen, D. T., Prodhom, B., Reiher, W. E., Roux, B., Schlenkrich, M., Smith, J. C., Stote, R., Straub, J., Watanabe, M., Wiorkiewicz-Kuczera, J., Yin, D., and Karplus, M. (1998) All-atom empirical potential for molecular modeling and dynamics studies of proteins. Journal of Physical Chemistry B 102, 3586-3616.
6. MacKerell, A. D., Wiorkiewicz-Kuczera, J., and Karplus, M. (1995) An all-atom empirical energy function for the simulation of nucleic acids. Journal of the American Chemical Society 117, 11946-11975.
7. Pavelites, J. J., Bash, P. A., Gao, J., and MacKerell, A. D. (1997) A Molecular Mechanics Force Field for NAD+, NADH, and the Pyrophosphate Groups of Nucleotides. Journal of Computational Chemistry 18, 221-239.
8. Hermans, J., Berendsen, H. J. C., van Gunsteren, W. F., and Postma, J. P. M.,. (1984) A Consistent Empirical Potential for Water-Protein Interactions. Biopolymers 23, 1.
9. Ott, K.-H., B. Meyer. (1996) Parametrization of GROMOS force field for oligosaccharides and assessment of efficiency of molecular dynamics simulations. J Comp. Chem. 17, 1068-1084.
10. (1994), Biosym Technologies, San Diego.
11. Clark, M., Cramer III, R. D., and van Opdenhosch, N. (1989) Validation of the General Purpose Tripose 5.2 Force Field. J. Comp. Chem. 10, 982-1012.
12. Levitt, M., and Lifson, S. (1969) Refinement of protein conformations using a macromolecular energy minimization procedure Journal of Molecular Biology 46, 269-279.
13. Gunsteren, W. F., and Karplus, M. J. (1980) A method for constrained energy minimization of macromolecules. Journal of Computational Chemistry 1, 266-274.
14. Fletcher, R., and Reeves, C. M. (1964) Function minimization by conjugate gradients. The Computer Journal 7, 149-154.
15. Fletcher, R. (1980) Practical Methods of Optimization, Vol. 1, John Wiley & Sons, New York.
16. Phillips, J. C., Braun, R., Wang, W., Gumbart, J., Tajkhorshid, E., Villa, E., Chipot, C., Skeel, R. D., Kale, L., and Schulten, K. (2005) Scalable molecular dynamics with NAMD. J Comput Chem 26, 1781-802.
17. Nilsson, L., and Karplus, M. (1986) Empirical Energy Functions for Energy Minimizations and Dynamics of Nucleic Acids. Journal of Computational Chemistry 7, 591-616.
18. MacKerell, A. D., Feig, M., and Brooks, C. L. (2004) Extending the treatment of backbone energetics in protein force fields: limitations of gas-phase quantum mechanics in reproducing protein conformational distributions in molecular dynamics simulations. Journal of Computational Chemistry 25, 1400-1415.
19. Feller, S. E., and MacKerell, A. D. (2000) An Improved Empirical Potential Energy Function for Molecular Simulations of Phospholipids. Journal of Physical Chemistry B 104, 7510-7515.