БИОЛОГИЧЕСКИЕ МЕМБРАНЫ, 2008, том 25, № 1, с. 66-75
УДК 539.1
СРАВНИТЕЛЬНОЕ ИЗУЧЕНИЕ МОЛЕКУЛЯРНОЙ ДИНАМИКИ, ДИФФУЗИИ И ПРОНИЦАЕМОСТИ ПО ОТНОШЕНИЮ К ЛИГАНДАМ ДЛЯ БИОМЕМБРАН С РАЗЛИЧНЫМ ЛИПИДНЫМ СОСТАВОМ
© 2008 г. К. В. Шайтан1, М. Ю. Антонов1' 2, Е. В. Турлей1, О. В. Левцова1, К. Б. Терёшкина1, И. Н. Николаев2, М. П. Кирпичников1
Московский государственный университет им. М.В. Ломоносова, 119991 Москва, Ленинские горы,
электронная почта: shaitan@moldyn.org, mikhail@moldyn.org, yegor@moldyn.org, 2Якутский государственный университет им. М.К. Аммосова, 677000 г. Якутск, ул. Белинского, 58,
электронная почта: n_ivan_n@mail.ru Поступила в редакцию 21.08.2007 г.
Методом молекулярной динамики проведено сравнительное изучение липидных бислоев различного состава. Использован метод управляемой молекулярной динамики для оценки кинетических параметров и проникновения малых молекул разной химической природы через биомембраны. В рамках единого протокола молекулярного моделирования обсуждается соответствие результатов численного эксперимента с данными по структуре биомембран, коэффициентам латеральной диффузии и сравнительной проницаемости мембран для различных лигандов.
Структурные и кинетические свойства биологических мембран играют важную роль в организации различных клеточных процессов [1]. Молекулярная динамика (МД) биомембранных структур представляет значительный интерес как в связи с развитием молекулярных и мембранных технологий [2], так и благодаря тем уникальным возможностям детализации и визуализации молекулярных процессов в достаточно сложно организованных структурах, которые предоставляет современный численный эксперимент [3-5]. Особенно это касается микроскопической картины массопереноса в сильно анизотропных структурированных средах, диффузионных процессов на границе водной и мембранной фаз, формирования и релаксации неравновесных гетерофазных структур [6, 7].
Следует отметить, что исследование методом МД в мембранах [8, 9] зачастую сталкивается со значительными сложностями. Молекулярные модели мембранных систем содержат не менее 104 атомов, и в полноатомном силовом поле расчет траектории даже длиной порядка 100 нс является весьма трудоемким. С другой стороны, характерные времена пассивного транспорта составляют порядка микросекунд. В связи с этим необходимо развитие таких методов численного эксперимента, которые позволяли бы получать информативный результат за разумное время. Имеется несколько подходов к этой проблеме. Относительно простые модели виртуальной гидрофобной среды [10, 11] позволяют достаточно быстро оценить влияние фактора гидрофобности на возможные структурные изменения молекул на границе фаз. При этом,
как правило, теряется информация о вкладе куло-новских взаимодействий мембранной структуры с лигандами. Использование огрубленных моделей липидного слоя (coarse-grained [12]) позволяет также сократить время получения результата. Однако даже более детальные тяжелоатомные модели липидного бислоя не дают надежных данных для кинетических коэффициентов [13, 14]. Конечно, почти любую модель для конкретной структуры можно откалибровать таким образом, что некоторые расчетные параметры будут соответствовать экспериментальным значениям. При этом остается открытым вопрос об универсальности разработанного протокола молекулярного моделирования при переходе к другим структурам. В ряде предыдущих работ по динамике гидратированного бислоя 1-пальмитоил-2-олеоил-,от-глицеро-3-фос-фатидилхолина (ПОФХ) был разработан протокол молекулярного моделирования на основе полноатомного силового поля Amber 1999 [15], который дает удовлетворительные результаты как для основных структурных, так и кинетических параметров. При этом важным условием является приготовление системы в состоянии, которое характеризуется равновесными распределениями флукту-аций основных макропараметров (температуры, объема, давления) [16, 17]. В данной работе расширяется набор мембранных структур, изучаемых с использованием этого подхода. Отметим, что используемые значения латерального давления согласуются с данными [18].
Ниже на примере трех типов распространенных липидов и образуемых ими бислоев прослеживаются закономерности в структуре и динамике
мембран, по составу близких к биомембранам высших организмов. При этом рассматривается также задача проницаемости биомембран для низкомолекулярных лигандов, встречающихся в тканевых жидкостях. Экспериментальной информации о физических деталях этих процессов относительно немного, хотя известно достаточно много данных в том числе об относительно легко и трудно проникающих молекулах [19]. Для сравнительного изучения проницаемости мембран по отношению к этим лигандам в данной работе используется метод управляемой (направленной) МД [16, 17, 20, 21]. Суть метода состоит во включении внешнего воздействия на лиганд и направлении эволюции системы по определенному сценарию. Это позволяет наблюдать трансмембранный транспорт даже на сравнительно коротких по длине траекториях и дает возможность оценить количественные характеристики параметров, характеризующих перенос лигандов.
ЭКСПЕРИМЕНТАЛЬНАЯ ЧАСТЬ
В работе изучались модели липидных бислой-ных мембран, состоящих из следующих липидов в различных пропорциях: 1-пальмитоил-2-олеоил-,от-глицеро-3 -фосфатидилхолин, 1,2-дипальмито-ил-.от-глицеро-3-фосфатидилхолин (ДПФХ), 1,3-(1-стеароил-2-пальмитоил-,от-глицеро-3-фосфати-дил)глицерин (кардиолипин, КЛ). Для расчета траекторий использовался пакет молекулярной динамики PUMA [22, 23], для решения классических уравнений - алгоритм Верле с потенциальным полем Amber 1999 [15].
Проводились расчеты для трех гидратирова-ных бислоев, содержащих: 8 молекул КЛ, 16 - ПО-ФХ, 16 - ДПФХ (система I); 30 молекул ПОФХ, 30 - ДПФХ (система II); 64 молекулы ПОФХ (система III). Использовались периодические граничные условия. Начальная структура бислоя соответствовала перпендикулярному положению оси наибольшей протяженности молекул липидов относительно плоскости мембраны. Молекулы липидов предварительно поворачивались вокруг своей длинной оси на случайный угол. На молекулу липида приходилось 34-43 молекулы растворителя (для большинства липидов величиной полной гидратации считается не менее 27 молекул воды на липид [24]). В стартовой конфигурации молекулы воды помещались на расстоянии не менее 2.3 А от крайних атомов мембраны. Система изначально приготовлялась в молекулярном конструкторе так, чтобы удельная площадь, приходящаяся на молекулу липида, соответствовала экспериментальным данным [25-30] и составляла в течение расчета для ПОФХ 62-68 А2, для ДПФХ 59-62 А2, для КЛ 100-120 А2. В случае бислоев из разных липидов удельная площадь вычислялась как среднее
Таблица 1. Величина латерального давления в системах
Система Средняя удельная площадь, А2 Латеральное давление,бар
I 78.0 ± 5.5 -242
II 65.4 ± 1.5 -300
III 66.8 ± 3.7 -265
арифметическое площадей участвующих липидов согласно их концентрации.
Парциальные заряды и силовые константы для молекул липидов вычислялись в соответствии с процедурой, описанной в [16]. Отрицательный заряд молекулы кардиолипина (-2) компенсировался добавлением в воду ионов натрия. Использовалась модель воды Т1Р3Р, в которой валентные связи и углы не фиксировались. Ван-дер-ваальсовы взаимодействия вычислялись с применением специальной сглаживающей (переключающей) функции, а кулоновский потенциал умножался на специальную экранирующую функцию, как это описано, например, в [16]. Радиус обрезания кулоновских взаимодействий составлял от 16 до 20 А в зависимости от типа мембраны. Константа диэлектрической проницаемости полагалась равной единице. Шаг численного интегрирования - 1 фс.
Расчет проводился с периодическими граничными условиями при постоянной температуре и при условии постоянного давления (ЛРГ-ансамбль). Ба-ростатирование осуществлялось баростатом Бе-рендсена с одинаковым для всех направлений временем релаксации, равным 100 пс. Для учета эффектов поверхностного натяжения бислоя и поддержания удельной площади, приходящейся на липид в составе бислоя, соответствующей экспериментальным данным [25-30], латеральная компонента давления баростата бралась отрицательной [31]. Температура поддерживалась равной 300 К с помощью столкновительной среды (столк-новительного термостата [23]). Средняя частота столкновений виртуальных частиц составляла 10 пс-1, масса частиц 1 а.е.м.
РЕЗУЛЬТАТЫ И ИХ ОБСУЖДЕНИЕ
На первом этапе эксперимента проводилась релаксация бислоя в течение 200-500 пс при температуре 500 К. Затем проводился набор рабочего участка МД траектории в течение не менее 2 нс.
При изучении динамики к системам прилагалось отрицательное латеральное давление (табл. 1), необходимое для поддержания средней удельной площади на начальном уровне.
На рис. 1а-в приведены плотности вероятности флуктуаций объема расчетной ячейки для рассматриваемых систем. Видно, что плотность вероятности флуктуаций объема расчетной ячейки
68
ШАЙТАН и др.
р х 103 1. 1.
9.70 9.75
V х 10"
9.80
4, А3
9.55 9.60 9.65 А2
81 80 79 78 77 76
1600 2000 2400 2800 3200 3600
г, пс
г, а
36.0
1600 2000 2400 2800 3200 3600
г, пс
1.2801.2841.288 1.2921.296
V х 10-5, А3
1.60
1.61
1.62 1.63
V х 10-5, А3
А2 66.2 66.0 65.8 65.6 65.4 65.2 65.0 64.8 64.6
1600 1800 2000
2200 г, пс
3200 3600 4000 4400 4800 г, пс
г, а
36.6 36.4 36.2 36.0 35.8 35.6 35.4 35.2 35.0 34.8 1500
1700 1900 2100 3200 3600 4000 4400 4800 г, пс г, пс
Рис. 1. Параметры расчетной ячейки, а-в - плотность вероятности значений объема расчетной ячейки, и гауссовская аппроксимация кривой; г-е - флуктуации удельной площади, приходящейся на молекулу липида; ж-и - флуктуации толщины бислоя для систем 1-111 соответственно.
имеет вид Гауссова распределения. Согласно теории термодинамических флуктуаций Эйнштейна, плотность вероятности равновесных флуктуаций объема р(АУ) задается гауссовым распределением:
(А V)2
/АТА л 2 <А V2)
р (А V) = Ave \
символом < ) здесь обозначена операция взятия среднего, А^^ - нормировочный множитель, <АV2) = = къ7У%т - дисперсия, %т - коэ
Для дальнейшего прочтения статьи необходимо приобрести полный текст. Статьи высылаются в формате PDF на указанную при оплате почту. Время доставки составляет менее 10 минут. Стоимость одной статьи — 150 рублей.