Содержание к диссертации
Введение
Глава 1. Молекулярное моделирование биомембран 7
Молекулярно-динамические расчеты 7
Силовое поле в молекулярной динамике 8
Определение парциальных атомных зарядов методами квантовой химии . 16
Анализ заселенностей орбиталей 20
Уравнения движения 22
Выбор ансамбля 23
Контроль давления и температуры 26
Экспериментальные данные по бислоям 29
Монослои - модель мембраны в молекулярно-динамических расчетах 31
Модели неявно заданной мембраны 31
Изучение гидрофобных взаимодействий на мембранном интерфейсе 33
Изучение диффузии малых молекул через мембраны 36
Глава 2. Гидрофобная мембрана на основе тетрадекана 39
Простые модели гидрофобных сред 39
Протокол молекулярной динамики 39
Результаты и обсуждение 41
Глава 3. Динамическая гетерогенность липидного бислоя 55
Постановка задачи 55
Протокол молекулярной динамики 56
Релаксация бислоя 58
Флуктуации параметров бислоя 58
Поверхностное натяжение бислоя 62
Распределения атомных групп и характер укладки атомов бислоя 63
Параметры порядка 65
Профиль электростатического потенциала 67
Латеральная диффузия липидов и микровязкость мембраны 68
Обсуждение оптимальной длины траектории 76
Глава 4. Взаимодействие углеродных нан отру бок и их комплексов с биомембранами 81
Перспективы молекулярного моделирования в бионанотехнологии 81
Постановка численного эксперимента 83
Нанотрубка 83
Внутренне содержимое нанотрубки 84
Ван-дер-ваальсовы сферы 84
Липидный бислой 85
Расчеты с абсорбцией 85
Моделирование «нановыстрела» 85
Абсорбция холестерина 86
Взаимодействие полипептида с нанотрубкой 88
Взаимодействие углеродной нанотрубки с фосфолипидным бислоем 90
«Наношприц» в действии 91
Заключение 95
Выводы 100
Благодарности 102
Литература 103
- Определение парциальных атомных зарядов методами квантовой химии
- Протокол молекулярной динамики
- Распределения атомных групп и характер укладки атомов бислоя
- Внутренне содержимое нанотрубки
Введение к работе
МД мембранных структур в настоящее время представляет значительный интерес в связи с развитием молекулярных технологий и биотехнологий. Вместе с тем эти объекты достаточно трудны для экспериментальных исследований, и фундаментальные закономерности динамического поведения таких структур остаются все еще не вполне ясными, несмотря на огромный прогресс в этой области. Это касается и детальной микроскопической картины массо- и энергопереноса в сильно анизотропных структурированных гетерогенных средах, формирования и релаксации неравновесных надмолекулярных структур, особенностей распределения молекулярных групп с различной полярностью на границах раздела фаз. В данной работе метод МД с использованием полноатомных силовых полей, специальных процедур и достаточно длинных траекторий применяется для уточнения микроскопической картины диффузионных процессов на границе водной и мембранной фаз. Использование полноатомного приближения делает МД-расчеты столь больших систем достаточно трудоемкими и не позволяет надеяться на достижение термодинамического равновесия. Поэтому весьма актуальным является развитие таких методов МД и протоколов расчета, которые позволяли бы за разумные времена вычислять величины, сопоставимые с экспериментальными данными. Сравнительно новым подходом здесь является развиваемый ниже вариант метода управляемой динамики, который позволяет стимулировать молекулярные процессы по определенным степеням свободы. Использование этого метода при одновременном контроле над равновесным характером распределений значимых с физической точки зрения величин (скоростей, объема, площади и давления) на масштабе характерных времен интересующих процессов позволяет дать количественную оценку параметрам, характеризующим физические механизмы элементарных актов переноса в микрогетерогенных структурах. Важным с биологической точки зрения явля-
ется использование подобного подхода в дизайне и моделировании функциональных наноструктур, задействованных, например, при направленной доставке веществ через мембрану.
Определение парциальных атомных зарядов методами квантовой химии
Квантово-химические методы исследования отталкиваются от главного постулата квантовой механики - о существовании волновой функции ц/ (q\, #2, , 7n, t), которая полностью описывает поведение системы частиц во времени. Эта функция должна удовлетворять уравнению Шрёдингера где р - импульс частицы, q - ее координата, t - время, Н - гамильтониан системы.
Точное решение этого уравнения возможно только для водородоподоб-ных атомов, поэтому используют приближение Борна-Оппенгеймера, согласно которому энергетические термы и соответствующие волновые функции для разных типов движения, сильно отличающихся своими временами, в произведении могут давать общую характеристику системы. Так что волновую функцию системы можно представить как произведение электронной и ядерной волновой функций. Далее работают с электронной волновой функцией. При этом гамильтониан тоже приходится упрощать и представлять через разбиение на сумму одноэлектронных вкладов. Но невозможно разбить многочастичный гамильтониан на набор операторов без нарушения принципа неразличимости идентичных частиц. Все электроны, согласно уравнению Шрёдингера, вносят одинаковый вклад в электронную плотность в любой точке реального пространства, и при разбиении пространства все электроны трактуются эквивалентно. Таким образом, определение энергии атомов осу ществляется не через пространственное разбиение гамильтониана (или элементов гильбертова пространства, к которым оно прилагается), а скорее через разбиение гамильтониана на сумму эффективных одноэлектронных вкладов [35].
В современной квантовой химии межэлектронное отталкивание учитывается в форме взаимодействия между электроном на данной орбитали и усредненным полем остальных электронов молекулы. Этот подход известен как метод ССП, или метод Хартри.
Здесь первый член - оператор кинетической энергии электронов. Второй -оператор потенциальной энергии взаимодействия п электронов с ядром. Третий - оператор энергии межэлектронного отталкивания. В методе ССП мож 1 но заменить потенциал типа —, зависящий от координат двух электронов, выражением, описывающим межэлектронное взаимодействие как функцию координат каждого отдельного атома. Полная волновая функция записывается как произведение волновых функций отдельных электронов:
Полная энергия достигает минимума при 8Ф = 0. Отсюда находят функции щ. Метод ССП включает итерационную процедуру. На каждом ее шаге молекулярные орбитали изменяются таким образом, что полная электронная энергия стремится к своему минимальному значению. Этот процесс называется самосогласованием.
Часто в расчетах по методу МО ССП необходимо рассматривать молекулы и с замкнутыми, и с открытыми оболочками. То есть возникает необходимость проводить расчеты с неспаренными электронами и без них.
Молекулы со спаренными электронами обычно рассчитываются в приближении ограниченного метода Хартри-Фока (RHF - Restricted Hartree-Fock). Каждая молекулярная орбиталь в этом случае занята двумя электронами или свободна. Орбиталь радикала при расчете по методу RHF и определитель для волновой функции имеют вид [36 г 37]
Проблема спина электрона здесь не возникает, так как все электроны спарены. Метод расчета систем с открытыми оболочками - неограниченный метод Хартри-Фока (UHF - Unrestricted Hartree-Fock). Здесь возникает два набора молекулярных орбиталей, отдельно для каждой ориентации спина.
Протокол молекулярной динамики
Процесс проникновения малых молекул через липидный бислой представляет интересную тему для изучения МД [89], так как экспериментальные подходы оказываются неспособными детально выявить перемещения малых пенетрантов внутри мембранного окружения. Для большинства пенетрантов равновесная концентрация в мембране слишком мала, чтобы быть выявленной экспериментально. Присоединение специальных меток не кажется возможным, поскольку они меняют проницаемость весьма драматично, особенно у малых пенетрантов. Лучшие экспериментальные данные по скоростям проникновения могут быть получены с помощью осмотических, радиографических, ЭПР- и ЯМР-измерений [90-95].
Было проведено несколько исследований большого числа пенетрантов, от слабополярных до полярных. Главные выводы, сделанные в этих исследованиях, следующие: растворимость пенетрантов в мембране лучше всего коррелирует с их растворимостью в алканах, лимитирующий скорость этап определяется прохождением участка мембраны сразу за головками липидов, для малых пенетрантов существует сильная зависимость от размера, отражающая эффект, найденный в мягких полимерах. Главная цель МД-расчетов - понять эти явления на более детальном уровне.
Ранние расчеты фокусировались, главным образом, на том, как реалистично симулировать сами мембраны, тогда как изучение переноса молекул началось не так уж давно. Так как процесс проникновения большинства молекул достаточно медленный, по крайней мере, при временных масштабах, доступных МД (не более десятков наносекунд), то приходится идти на различные ухищрения, чтобы получить статистически значимую информацию. Например, несколько секунд МД-траектории необходимо для мембраны поперечным размером 4x4 нм, чтобы наблюдать спонтанный транспорт воды через мембрану. Стандартной процедурой в этой связи является искусственное проталкивание пенетранта в мембрану либо помещение его стартовой позиции на уровне мембраны. Как показано в результате расчетов липидных мембран, мембранное окружение весьма далеко от гомогенности. Поэтому, чтобы получить полное описание процесса проникновения, необходимо сделать выборку взаимодействий пенетрант-мембрана на протяжении всей мембраны.
Чтобы это сделать, можно применить какие-либо ограничения, чтобы произвести локально выборку взаимодействий пенетрант-мембрана, либо следить за траекторией пенетранта во времени. Чтобы полностью описать процесс проникновения, необходимо знание локальных скоростей диффузии и локальной растворимости. В настоящее время большинство МД-исследований фокусируется пока что только на диффузии.
Первое МД-исследование в этой области относится к 1992 г. [96]. McKinnon и др. изучили диффузию кислорода через легочные мембраны и зависимость диффузионной скорости от концентрации холестерола. Однако, модель мембраны была достаточно грубой, состоящей из монослоя молекул гексадекана, ограниченных с одной стороны гармоническим потенциалом. Другой гармонический потенциал был приложен к молекуле кислорода. Его минимум перемещался по нормали к монослою с заданной скоростью, чтобы ускорить проникновение кислорода через монослой. Коэффициент диффузии мог быть вычислен по усредненной силе трения, действующей на кислород. Хотя модель была достаточно простой, вычисленный коэффициент диффузии кислорода через гексадекановый монослой оказалась в хорошем соответствии с экспериментальными коэффициентами кислорода в ДПФХ, полученными в эксперименте с тушением флуоресценции пирена в бислое (2.6x10" 9 S 9 см /с в МД при 2.3x10" см /с в эксперименте). Как предполагается, такой коэффициент диффузии кислорода в мембране тесно связан со скоростью перемещения кинков [97].
В других исследованиях выяснилось, что мембрана ведет себя как ускоритель проникновения неполярных молекул (кислорода в т.ч.). Рис. 10 пока зывает поразительную разницу между сопротивлением мембраны проникновению сильно полярных пенетрантов (воды), умеренно полярных (аммиака) и неполярных (кислорода). Происхождение этого эффекта в большей растворимости неполярных пенетрантов в липидной мембране.
Исследования проникновения заряженнных пенетрантов (ионов) показало, что вычисляемая свободная энергия сильно зависит от состояния гидра-тированности иона по мере того, как он перемещается из водного слоя в неполярную среду. Экспериментальные исследования проникновения ионов через границу раздела фаз вода-масло показали также хорошее согласие полученной энергии иона в потенциальной яме водного слоя с борновской энергией сольватации [98].
Интенсивное изучение динамики биологических мембран [1-3, 3-11, 11, 99] и белок-мембранных комплексов [100-106]ставит также вопрос об оптимизации расчетных процедур и моделировании мембранных структур некоторой гидрофобной средой [107]. Имеются относительно простые модели виртуальной гидрофобной среды [80]. Параметры изменения энергии групп при переходе из воды в мембрану в этих моделях определяются из данных по константам равновесия в системе вода-октанол [80, 107]. Метод МД позволяет рассчитать функции распределения молекулярных групп в системе вода-мембрана и оценить возможные неаддитивные эффекты.
В работах, посвященных изучению диффузии малых молекул в мембранах, состоящих из липидов или даже из более простых органических соединений [10] обычно не удается полностью проследить пассивный трансмембранный транспорт даже очень маленьких молекул, например воды, кислорода, этанола, аммиака, ионов калия, хлора. Сравнительно новым подходом здесь является развиваемый ниже вариант метода управляемой (стимулированной) динамики SMD [108, 109], который позволяет за разумные времена расчета получить несколько полных проходов исследуемой молекулы-пенетранта через мембрану и оценить картину диффузионных процессов.
Ниже изучена диффузия нескольких типов молекул (как реалистических и биологически важных, так и модельных ван-дер-ваальсовых сфер) в гидра-тированной углеводородной мембране.
Распределения атомных групп и характер укладки атомов бислоя
Следует отметить, что согласование удельной поверхностной площади липидов и толщины мембраны с экспериментальными значениями часто является основным критерием реалистичности модели. Более тонкое согласование моделей осуществляется и по другим важным параметрам. В первую очередь это касается функций распределения по нормали к мембране усредненных электронной плотности и плотности атомных групп. Аналогичные зависимости вычисляли для ДПФХ- [152], ДОФХ- [153] и ПОФХ- [137] мембран. На Рис. 22 представлены распределения трех типов: усредненных атомной, массовой и электронной плотностей.
Характер нерегулярной укладки атомов фосфора и азота липидных голов в плоскости мембраны изображен на Рис. 23 а.
В ходе исследования также изучались радиальные функции распределения атомов g{r) в плоскости мембраны, которые определяют вероятность нахождения атома определенного типа на заданном расстоянии от другого атома в проекции на плоскость мембраны. Число атомов dN, находящихся в круговом слое площадью dS и толщиной dr на расстоянии г от центрального атома (N - общее число атомов данного типа) связано с g(r) формулой dN = Ng{r)— = —g{r)2nrdr. (45)
На Рис. 236 представлены радиальные функции распределения для атомов фосфора и атомов азота в плоскости мембраны (усреднение вели по обоим слоям мембраны). Два главных максимума находятся на5.8Аи9.іА для атомов фосфора и 6.5 А и 8.3 А для атомов азота. Радиальная функция распределения атомов азота демонстрирует возможность довольно близкого (на расстоянии менее 1.5 А) нахождения атомов азота в слое липидных голов. При этом толщина слоя гидрофильных голов липидов такова, что позволяет разместить на разной глубине в мембране относительно слабо заряженные и отталкивающиеся друг от друга атомы азота. Отметим, что сильно заряженные атомы фосфора на таком близком расстоянии располагаться не могут. Интегрирование радиальных функций распределения дает следующие координационные числа для первой координационной сферы: 1.92 для атомов фосфора и 1.47 для атомов азота. Во второй координационной сфере -2.38 для атомов фосфора и 3.17 для атомов азота. Характер кривых на Рис. 236 указывает на отсутствие дальнего порядка в упаковке голов липидов.
Характеристики расположения голов липидов в плоскости бислоя. а - распределение атомов азота и фосфора в одном из слоев мембраны. Ортогональная проекция. Черным отображены атомы фосфора, светлым - азота. Показаны также вектора, соединяющие атом фосфора и атом азота в пределах одной молекулы липида; б - радиальные функции распределения атомов азота и фосфора в плоскости мембраны. 1 - Р-Р, 2 - N-N. Параметры порядка
По результатам МД-расчетов определяли параметр порядка SQB (Рис. 24), который может быть также получен экспериментально из данных ЯМР-спектроскопии для дейтерированных липидов. Усредненные значения параметра порядка могут быть получены и по данным ИК-спектроскопии [154]. Параметр порядка С-Н-связей в алкильной части липида вычисляется следующим образом: 5CH(0 = (3cos2 -l). (46) где в[ - угол между С-Н-связью при /-том атоме углерода в алкильной цепи и нормалью к мембране, а угловые скобки обозначают усреднение по времени. Максимальное значение этой величины 1 (когда все связи параллельны нормали), минимальное -0.5 (все связи лежат в плоскости мембраны). Полученные средние значения -SCH ДЛЯ обеих углеводородных цепей липида превышают значения, полученные в экспериментах [155-162] и других МД-расчетах [4, 163, 164] на 0.02-0.03. Это обусловлено, по-видимому, конечным радиусом обрезания электростатических взаимодействий. Чем меньше радиус обрезания, тем более упорядоченными являются липидные цепи [27]. Рассчитанные значения -SCR для начального участка пальмитоиловой цепи уменьшаются в направлении липидных голов. Это соответствует, за исключением второго атома углерода, экспериментальным данным [156], а также МД-расчетам [4, 164]. Однако экспериментальные значения -SCR, согласно [157-159], возрастают по мере приближения к липидным головкам. Наилучшее согласование по профилю параметра порядка для пальмитоила с другими МД расчетами наблюдается с данными [4]. С другой стороны, наилучшее согласование с экспериментальными данными показано для результатов измерений [162]. В целом, картина изменения расчетных значений -SCR вдоль пальмитоиловой цепи в составе липида остается похожей на профиль -SCR для свободной молекулы н-гексадекана [165]. Полученный профиль -SCR для олеоиловой цепи в целом также согласуется с экспериментальными данными [156] и компьютерными расчетами других авторов [4, 5, 134, 151, 166]. Особенно хорошее согласование достигается для поведения -SCR на участке в центре цепи (атомы углерода 8-13) и для спада значений -SCR по мере удаления от центра цепи к свободному концу.
Внутренне содержимое нанотрубки
В ходе рассматриваемого процесса молекула полиаланина испытывает конформационные изменения. Начальная спиральная конформация наиболее сильно деформируется при выбросе полипептида в вакуум и менее всего - в мембрану. По-видимому, среда играет в этом процессе демпфирующую и структурирующую роль. На Рис. 40 видно также, что при уменьшении скорости «нановзрыва» конформационные изменения полипептида обычно уменьшаются. Небольшое уширение распределения по конформациям при выбросе полиаланина в воду может быть связано с дополнительным напряжением, возникающим при входе гидрофобной молекулы в водную среду.
Данный полноатомный МД-расчет подтвердил предположение о возмо-жости использовать углеродные нанотрубки как наноконтейнеры для биомолекул [198]. В самом деле, вне зависимости от того, нанотрубка предварительно гидратирована или нет, она легко может поглотить небольшие негидрофильные молекулы. Другой особенностью подобных наноконтейнеров является их логическое развитие в устройства на основе нанотрубок, такие как «наношприц». Компьютерные расчеты, пусть даже в слегка упрощенной форме, выявляют базовую возможность их должного функционирования [199-202]. Имея в распоряжении молекулы, которые могут быть активированы светом [182, 203], можно переходить к экспериментальной стадии разработки.
В будущем, модифицируя нанотрубку (путем добавления функциональных групп) можно достичь как селективной абсорбции, так и селективной посадки нанотрубки на клеточных мембранах. Результаты данного МД-исследования являют собой удобный инструмент в проектировании будущего реального устройства.
Разработанные методы и протоколы МД, в том числе SMD, позволяют получить новую информацию о динамике и функционировании молекулярных конструкций на основе мембран. Использованная простейшая полноатомная модель гидратированной мембранной структуры на основе тетраде-кана позволяет исследовать закономерности формирования распределения различных молекул между двумя фазами - гидрофобной и водной. Как было показано, лишь очень маленькие и сильно гидрофобные молекулы, например кислород, способны за времена меньше 10 не самостоятельно проникнуть в мембранный слой при нормальных условиях. При фиксированном объеме проникновение в структурированную мембранную среду сильно затруднено. Обнаружено накопление кислорода в середине бислоя, что обусловлено повышенной рыхлостью структуры на интерфейсе двух монослоев.
Изучение динамики трансмембранной диффузии за времена порядка 10 не требует проведения расчетов при повышенной температуре. Энергия переноса, вычисленная при температуре 1000 К для основных типов атомов Amber 1999, варьируется от 0.9 до 5.5 ккал/моль. Эти изменения энергии находятся в тесной корреляции с изменением соответствующей борновской энергией сольватации. Важно отметить отсутствие аддитивности вкладов атомов в свободную энергию переноса функциональных групп. Это затрудняет простую континуальную имитацию гидрофобной среды путем включения в потенциальную энергию дополнительного терма, описывающего взаимодействие атомов с этой средой.
Вычисленные распределения аминокислотных остатков между тетраде-кановым монослоем находятся в соответствии с их известными гидрофобными свойствами.
Рассмотрение процессов равновесной динамики даже относительно небольших молекул в рассматриваемых системах не позволяет за разумные времена численного эксперимента сделать окончательные выводы о балансе гидрофобных сил при комнатной температуре. Развиваемый метод неравновесной динамики (SMD) позволяет организовать направленный и более быстрый сценарий развития событий по определенным степеням свободы. При этом вместо накопления равновесных траекторий осуществляется контроль над характером распределений значимых с физической точки зрения величин (скоростей, объема, площади и давления), который позволяет судить о достижении локального равновесия. Вычисляемые таким способом кинетические коэффициенты термодинамическими и статистическими соотношениями связаны с равновесными параметрами системы. Модельные расчеты диффузии на молекулярных масштабах показывают уже при скоростях 1 А/пс существенные отклонения от движения сферических частиц в сплошной среде, описываемого гидродинамической формулой Стокса, что вполне естественно - приближение сплошной среды на молекулярном уровне работает плохо. Формула Стокса в этом случае может использоваться лишь для качественных оценок. Аналогично обстоит дело и с динамикой образования и релаксации надмолекулярных мембранных структур (в контексте мембраны, пор). Образование поры кардинально меняет динамику трансмембранного транспорта. Время релаксации поры в тетрадекановой мембране относительно велико (более 10 не), и эффекты памяти в динамике такой простейшей мембраны могут быть весьма существенными.