Бесплатный автореферат и диссертация по биологии на тему

Молекулярная динамика пептидов и сравнительное изучение динамических свойств природных L-аминокислотных остатков

ВАК РФ 03.00.02, Биофизика

Текст научной работыДиссертация по биологии, кандидата физико-математических наук, Сарайкин, Сергей Сергеевич, Москва

р

(л ' ^ 10^ ~

Московский государственный университет им. М.В. Ломоносова Биологический факультет

Сарайкин Сергей Сергеевич

УДК 577.3

МОЛЕКУЛЯРНАЯ ДИНАМИКА ПЕПТИДОВ И СРАВНИТЕЛЬНОЕ ИЗУЧЕНИЕ ДИНАМИЧЕСКИХ СВОЙСТВ ПРИРОДНЫХ Ь-АМИНОКИСЛОТНЫХ

ОСТАТКОВ.

03.00.02 - биофизика

Диссертация на соискание ученой степени кандидата физико-математических наук

Научный руководитель доктор физ.-мат. наук профессор К.В. Шайтан

Москва, 1999

Оглавление.

Введение 4

Глава 1. Методы компьютерного моделирования динамики молекул.

1.1. Молекулярная динамика.

1.1.1. Расчёт ньютоновских траекторий движения

1.1.2. Методы ускорения расчётов молекулярной динамики.

1.1.3. Методы молекулярной динамики при постоянной температуре.

1.1.3.1. Методы, фиксирующие кинетическую энергию.

1.1.3.2. Методы, моделирующие динамику системы в тепловом резервуаре.

1.2. Метод Монте-Карло.

1.2.1. Метод Монте-Карло с критерием Метрополиса.

1.2.2. Модель роста цепи.

1.2.3. Модификации метода Монте-Карло с критерием Метрополиса.

1.2.4. Молекулярная динамика с Window Moves.

1.2.5. Монте-Карло в пространстве последовательностей аминокислот.

1.3. Метод минимизации потенциальной энергии.

1.4. Методы изучения молекулярной подвижности белка вблизи одной конформации.

1.4.1. Метод нормальных мод.

1.4.2. Поиск путей перехода из одного минимума энергии в другой.

1.4.3. Изучение гиперповерхностей потенциальной энергии и построение карт потенциальной энергии.

1.5. Упрощенные методы моделирования полипептидов.

1.5.1. Двумерные решетки.

1.5.2. Трехмерные решетки. Глава 2. Методы и выбор оптимальных условий расчетов для изучения молекулярной динамики олигопептидов. 33 2.1. Метод молекулярной динамики. 3 3

6

"б" --у

9 "12 12 ""14 20 21 "21" 21 "22" "23" "23

24 "24 "25 "25"

30 30

"зТ

2.2. Методы обработки траекторий молекулярной динамики.

35

2.2.1. Корреляционные функции. 35

2.2.2. Построение двумерных и одномерных распределений вероятностей значений торсионных углов. 42 2.3. Выбор оптимальных условий для расчета траекторий молекулярной динамики. 43 Глава 3. Изучение молекулярной динамики серий

модифицированных дипептидов. 48 Глава 4. Изучение свойств столкновительного термостата и его применение для расчётов коэффициентов диффузии

низкомолекулярных соединений в различных жидких средах. 58

4.1. Столкновительная динамика. 5 9

4.2. Диффузия в столкновительном термостате. 61 Глава 5. О влиянии амплитуды флуктуаций на коэффициент трения броуновского осциллятора в водной среде. 70

5.1. Введение. 70

5.2. Молекулярно-динамические модели для осциллятора в воде. 71

5.3. Модель броуновского осциллятора. 73

5.4. Обсуждение результатов молекулярно-динамического моделирования. 75

5.5. Заключение. 79 Заключение и выводы 81 Список литературы 84

Введение.

В последнее десятилетие взгляд на биологические макромолекулы в нативном состоянии как на объекты, имеющие статическую пространственную структуру, претерпел существенные изменения. Была обнаружена и понята большая роль пространственной подвижности макромолекул в их биологической функции.

В настоящее время имеется широкий набор экспериментальных методов для изучения динамики биомолекулярных систем: тушение флуоресценции, ЯМР, ЭПР, ИК и Рамановская спектроскопии, неупругое рассеяние нейтронов [1,2]. Компьютерное моделирование позволяет получать информацию о процессах, происходящих в системе на уровне движений отдельных атомов, недоступную экспериментальным методам, но необходимую для понимания механизмов, а также проверки существующих и развития новых теоретических моделей процессов играющих важную биологическую роль.

Метод МД успешно используется в теоретических исследованиях структуры и динамики биологических макромолекул. Однако, в данное время, представляется возможным моделирование внутренней подвижности белков лишь в субнаносекундных диапазонах времен. Это на несколько порядков меньше характерных времен даже быстрых ферментативных реакций. Для ускорения расчетов иногда применяют сильно упрощенные модели, однако, как будет показано ниже, они часто приводят к противоречивым выводам, и результаты, полученные при таком моделировании, не могут непосредственно использоваться для изучения динамики конкретных белков.

В настоящее время имеется реальная потребность в переходе от рекордных расчетов больших белков к планомерному сравнительному изучению динамических свойств аминокислот и олигопептидов. В литературе до сих пор нет надежных данных о фундаментальных принципах, лежащих в основе динамического поведения элементарных "кирпичиков", из которых построены белки.

Настоящая работа посвящена изучению динамических свойств природных биологических Ь-аминокислот в составе серий модифицированных дипептидов. Особое внимание уделяется моделированию влияния среды на конформационную подвижность дипептидов. В связи с этим разбираются свойства различных моделей виртуальных сред в методе МД.

В данной работе также рассматривается применение метода МД в двух задачах: диффузия молекул в столкновительном термостате и режимы движения гармонического осциллятора в плотной среде.

Глава 1. Методы компьютерного моделирования динамики молекул.

В данной работе мы используем определение потенциальной энергии биомакромолекулы как суммы атом-атомных потенциалов. В рамках данного метода молекула рассматривается как механическая система, и с помощью современной вычислительной техники и молекулярной графики можно изучать взаимное расположение молекулярных групп и отдельные конформационные превращения [3-8]. Основные методы компьютерного анализа конформационных изменений могут быть разбиты на два класса: статические и динамические [1]. В свою очередь динамические методы делятся на детерминистские и стохастические. Основная задача статических методов - поиск оптимальных с точки зрения потенциальной энергии конформаций. При этом переход от одной конформации к другой осуществляется только на основе правил численных методов минимизации. Динамические методы позволяют в той или иной мере исследовать траектории движения. Идея детерминистических динамических методов состоит в использовании собственной динамики для исследования движений системы, то есть в задании уравнений движения и их решении. В стохастических методах переход от одной конформации к другой осуществляется в марковском процессе, который является вероятностным аналогом собственной динамики, или путем введения случайной силы в уравнениях движения и переходу к системе уравнений Ланжевена. Рассмотрим имеющийся набор тесно связанных между собой методов и решаемых с их помощью проблем.

1.1. Молекулярная динамика.

Метод молекулярной динамики позволяет моделировать детальную микроскопическую картину внутренней подвижности макромолекулы. В его основе лежит расчет классических (ньютоновских) траекторий движения макромолекулы в фазовом пространстве координат и импульсов ее атомов [4]. Метод молекулярной динамики успешно используется в теоретических исследованиях структуры и динамики биологических макромолекул.

1.1.1. Расчёт ньютоновских траекторий движения.

В методе молекулярной динамики [9, 10] рассчитываются классические (ньютоновские) траектории движения атомов макромолекулы в силовом поле эмпирического атом-атомного потенциала, т. е. моделируется детальная микроскопическая картина внутренней тепловой подвижности макромолекулы в субнаносекундных интервалах времен. Основу метода составляет численное решение классических уравнений Ньютона для системы взаимодействующих частиц:

= ¿ = 1,2 „.„и, (1.1)

где г1 - радиус-вектор ьго атома, mi - его масса, В] суммарная сила, действующая на ьый

атом со стороны остальных частиц:

ВД-^Р (1.2)

Здесь: г = {гх,г2,...,гп}ш, и (г) -потенциальная энергия, зависящая от взаимного расположения всех атомов; п - число атомов.

Задав координаты и скорости всех частиц в начальный момент времени, числено решают уравнения движения, вычисляя на каждом шаге все силы и новые координаты и скорости частиц. Температура определяется как средняя кинетическая энергия, приходящаяся на одну степень свободы системы:

1 ^ , <3г.

7X0 =-Ущ*?, • (1-3)

Здесь N - полное число степеней свободы молекулы, кв - постоянная Больцмана. В случае изолированной системы N=311-6, поскольку сохраняется ее полный импульс и момент импульса. Кроме того, в этом случае сохраняется полная энергия системы, а температура получается усреднением ее мгновенных значений Т(0 по некоторому интервалу времени. Потенциальная энергия молекулы задается в виде:

Щг)= иь + иу+иф+ив + ии + иа + ии (1.4)

где слагаемые отвечают следующим типам взаимодействий: иь - химическим связям;

17у - валентным углам; IIф - торсионным углам; иа - плоским группам; ии - ван-дер-ваальсовым контактам; 11е1 - электростатике; IIкЬ - водородным связям. Указанные слагаемые имеют различный функциональный вид. Валентные длины поддерживаются за счет потенциала

иь=\^К„{г-Ь0)\ (1.5)

^ ь

где суммирование проводится по всем химическим связям, Ь0 - обозначение для равновесных валентных длин, г - текущие длины связей, Кь - соответствующие силовые константы.

Валентные углы задаются потенциалами

^^^Хо9-^)2, (1.6)

где <90 - равновесные значения углов, 3 - их текущие значения, К9 - силовые константы. Энергия торсионных взаимодействий и потенциалов, отвечающих плоским группам, записываются в одинаковом виде:

иф = £ Кф [соз(пФ - 5) +1], (1.7)

ф

где п - кратность торсионного барьера, 5 - сдвиг фазы, константы Кф определяют высоты потенциальных барьеров двугранных углов Ф.

Ван-дер-ваальсовые взаимодействия атомов, разделенных тремя и более валентными связями описываются с помощью потенциалов Леннард-Джонса:

в

г12 г6. ч ч

(1.8)

I ] г " г

Параметры потенциала А и В зависят от типов атомов 1 и участвующих во взаимодействий; - расстояние между этими атомами.

Электростатические взаимодействия задаются кулоновским потенциалом

(1.9)

где ql, qj - парциальные заряды на атомах, £ - диэлектрическая проницаемость среды.

Водородные связи возникают и исчезают в процессе движения атомов между теми из них, которые имеют донорно-акцепторный статус. Функциональный вид потенциала водородной связи аналогичен потенциалу ван-дер-ваальсовым взаимодействий:

Существуют различные наборы параметров для потенциалов взаимодействий. Их значения определяются из учета различных типов экспериментальных данных (спектральные, калориметрические, кристаллографические) и квантовомеханических расчетов.

зависимым координатам, что требует использование специальных алгоритмов [11] (в циклах и при фиксированных концах молекулы [12]).

1.1.2. Методы ускорения расчётов молекулярной динамики.

Время, необходимое для расчета траектории молекулы, можно значительно сократить, уменьшая число степеней свободы. Существует два способа ограничения движений длин валентных связей и углов. В одном случае длины валентных связей и значения валентных углов жестко фиксированы, в другом случае на них накладываются упругие ограничения с очень большой константой упругости. Статистические свойства жестко и упруго ограниченных систем, вообще говоря, различны. В работе [13] показано, что при разных формах упругого потенциала получаются статистически разные результаты, один из таких потенциалов соответствует жестко фиксированным валентным связям и углам. При этом, в общем случае, упругие потенциалы статистически предпочтительнее как для валентных связей, так и для валентных углов. Рассмотрим, например, следующий численный эксперимент: с помощью метода молекулярной динамики моделировали движение трехатомной (рис.1, а) и четырехатомной молекул (рис.1, б) в растворе со сферическими молекулами. Оказалось, что в случае жестких

А' В'

(1.10)

Некоторую трудность представляет вычисление производных энергии по

ограничений, в отличие от упругих, вектора, соединяющие первый и третий атомы в случае (а) и первый и четвертый атомы в случае (б), неравномерно распределены по сфере.

Рис. 1. Модели трехатомной и четырехатомной молекул. Обозначенные на рисунке вектора оказываются неравномерно распределенными по сфере при моделировании движения, если длины валентных связей и значения валентных углов жестко фиксированы.

В работах [14] и [15] сравниваются траектории молекулярной динамики BPTI (bovine pancreatic trypsin inhibitor) с жесткими и упругими ограничениями длин валентных связей и значений валентных углов. Расчеты с фиксированными связями дают почти те же результаты, что и с нефиксированными, особенно если используются специальные корректирующие алгоритмы. Однако, колебания валентных углов, по-видимому, связаны с коллективными движениями в молекуле и, из-за плотной упаковки атомов внутри белка небольшие флуктуации валентных углов (±4°) играют существенную роль в двйжениях, включающих другие степени свободы. При фиксации валентных углов амплитуда флуктуаций торсионных углов уменьшается в 2 раза, а конформационные переходы по торсионным углам из одного минимума энергии в другой исчезают совсем.

В некоторых случаях степени свободы, соответствующие изменениям значений валентных углов, учитываются неявно [16]. Этот учет валентных углов незначительно увеличивает время счета, но значительно увеличивает конформационную подвижность.

Иногда используют алгоритмы, в которых переменные, соответствующие медленным степеням свободы, постоянны на протяжении некоторого числа шагов.

а

Б

Однако, при таких расчетах, происходит довольно быстрое накопление ошибки. Этого недостатка лишены методы MTS (multiple-time-step methods). В этих методах для вычисления сил, соответствующих быстрым и медленным степеням свободы, используются разные временные интервалы [17]. Время счета при этом сокращается в 4 -5 раз.

Как правило, в методе молекулярной динамике для ускорения расчетов ван-дер-ваальсовые, водородные и электростатические взаимодействия рассчитываются только

о

между атомами, находящимися на расстоянии меньшем, чем радиус обрезания (10-15 А ). В работе [18] на примере изучения молекулярной динамики человеческого лизоцима в воде показано, что в результате обрезания электростатических взаимодействий увеличивается подвижность воды и изменяется динамика заряженных групп. В этой работе предложен также метод РРРС (particle-particle and particle-cell) расчета кулоновских взаимодействий. В методе РРРС каждый атом взаимодействует с ближними атомами путем обычных кулоновских взаимодействий, а с далеко отстоящими ячейками через общий заряд и дипольный момент этой ячейки. Размеры ячеек возрастают как функции расстояния от атома (рис.2).

Рис.2. Метод РРРС расчета кулоновских взаимодействий. Стрелками обозначены диполъные моменты ячеек, пунктирными линиями - кулоновские взаимодействия.

1.1.3. Методы молекулярной динамики при постоянной температуре.

Поскольку кинетическая энергия системы в молекулярной динамике определяется как 3/2*кьТ, то постоянную температуру можно поддерживать, зафиксировав общую кинетическую энергию.

1.1.3.1. Методы, фиксирующие кинетическую энергию.

Под изучаемой системой будет пониматься система из N атомов — материальных точек массы та, а=1, ...,К Это может быть атомная жидкость или молекула. Состояние системы описывается декартовыми координатами Х(х15х2,...х^)и скоростями У(\>1,л>2....ук,). Здесь ха и \асуть координаты атома номер а. Кинетическая энергия во время движения должна оставаться постоянной:

\^Ув-К = 0 (1.11)

здесь ^-фиксируемое значение кинетической энергии. Наложение связи добавляет в уравнения движения некоторую обобщённую силу ()(Х,¥). В работе [19] показано, что существует бесконечно много функций ()(Х, V), которые обеспечивают такой обмен энергией между системой и окружением, при котором выполняется соотношение (1.11).

Конкретное значение <2 естественно выбирать, исходя из физических принципов. Исходным принципом механики является принцип стационарности действия [20-21]. Он гласит: истинное движение Х($ есть экстремаль функционала действия. Вариация функционала равна нулю на множестве путей с закреплёнными концами

'1

'о

ЗДо) = г* > ха(ь) = ьа (1.12)

, 1 «

здесь саи Ьа- константы, а=1,...,И; 1(Х,Х) = — ^та\га ~и(х1,х2,...хм)- лагранжиан

2 а=1

системы. При наличии связи, множество путей, на котором определён функционал действия, сужается: допустимый путь кроме (1.12) должен удовлетворять наложенной связи. Естественно потребовать, чтобы и вариация пути удовлетворяла связи. Однако, если связь неголономна, то это требование приводит к трудностям [20,22]. Связь (1.12),

фиксирующая кинетическую энергию - неголономна. Существует два способа учёта неголономных связей: традиционная неголономная механика и вакономная механика [22]. Они отличаются друг от друга определением вариации допустимого пути.

В традиционной неголономной механике не требуется, чтобы вариация допустимого пути также удовлетворяла наложенной связи. Здесь лишь требуется, ч�