Бесплатный автореферат и диссертация по биологии на тему
Интегральные уравнения теории жидкостей в теоретическом изучении биологических макромолекул и их взаимодействий в растворах
ВАК РФ 03.01.02, Биофизика

Автореферат диссертации по теме "Интегральные уравнения теории жидкостей в теоретическом изучении биологических макромолекул и их взаимодействий в растворах"

На правах рукописи

005061495

Соболев Егор Васильевич

Интегральные уравнения теории жидкостей в теоретическом изучении биологических макромолекул и их взаимодействий в растворах

03.01.02 - Биофизика

АВТОРЕФЕРАТ диссертации на соискание ученой степени кандидата физико-математических наук

3 Ш) 2013

Пущино - 2013

Л

005061495

Работа выполнена в Федеральном государственном бюджетном учреждении науки Институте математических проблем биологии Российской академии наук

Научный руководитель:

кандидат физико-математических наук Тихонов Дмитрий Анатольевич

Официальные оппоненты: доктор физико-математических наук

Цыганкова Ирина Глебовна (в.н.с. ИТЭБ РАН, г. Пущино)

Защита состоится «26» июня 2013 г. в 16 часов на заседании совета Д 002.093.01 по защите диссертаций на соискание ученой степени кандидата наук, на соискание ученой степени доктора наук при Федеральном государственном бюджетном учреждении науки Институте теоретической и экспериментальной биофизики Российской академии наук по адресу: 142290, Московская область, г. Пущино, ул. Институтская, 3.

С диссертацией можно, ознакомиться в Центральной библиотеке НЦБИ РАН по адресу: 142290, Московская область, г. Пущино, ул. Институтская, д. 3.

Автореферат разослан «24» мая 2013 г.

Ученый секретарь 1' диссертационного совета, ,;

кандидат физкко:йатек!атйческих наук -//¿¿¿¿.і- ЛанинаН. Ф.

доктор физико-математических наук, профессор

Аграфонов Юрий Васильевич (декан ФФ ИГУ, г. Иркутск)

Ведущая организация:

Федеральное государственное бюджетное учреждение науки Институт биофизики клетки Российской академии наук, г. Пущино

Общая характеристика работы

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

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

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

Важный класс интегральных уравнений, описывающих равновесную структуру молекулярных жидкостей в терминах атом-атом парных корреляционных функций, предложили Чандлер и Андерсон в 1972 году. Эти уравнения основаны на модели связанных силовых центров (RISM, от англ. Reference Interaction Site Model).

Цели и задачи диссертационной работы. Цель работы заключается в исследовании особенностей модели связанных силовых центров при описании структуры и термодинамики сольватации биологических макромолекул и в создании методики учета растворителя на основе интегральных уравнений теории жидкостей.

Были поставлены следующие задачи:

1. Выполнить параметрический анализ уравнений ШБМ для плотного метана и провести сравнение поведения системы при использовании различных уравнений замыкания.

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

3. Сравнить метод ШБМ и континуальные подходы в задаче анализа термодинамики трех модельных состояний пептида окситоцина.

4. Сравнить решения усредненных по выборке из канонического ансамбля уравнения ЫЭМ и усредненных решений жестких уравнений ШБМ для каждого состояния из этой выборки.

5. Выполнить анализ термодинамики связывания 4',С-диамидино-2-фенилин-дола в малом желобе ДНК с целью выявить энергетически предпочтительный сайт связывания.

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

7. Разработать вычислительно-информационный Интернет-сервис для теоретического изучения гидратации биологических макромолекул.

Научная новизна. Предложены новые алгоритмы, позволившие решать уравнения ШБМ для множества мгновенных состояний биологической макромолекулы, состоящей из тысяч атомов. Алгоритмы реализованы в виде компьютерных программ и Интернет-сервиса.

Модель связанных силовых центров впервые применена для учета влияния растворителя в задаче о нековалентном связывании в растворе биологической молекулы (ДНК) и активного вещества (4',6-диамидино-2-фенилиндола) и показано преимущество интегральных уравнений перед континуальными подходами.

Впервые показано, что эффекты сольватации и термодинамику макромолекул можно описать с помощью уравнений ШБМ с усредненными по состояниям макромолекулы внутримолекулярными корреляционными функциями.

Впервые выполнен параметрический анализ уравнений ШЭМ.

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

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

Разработанными вычислительными программами можно пользоваться через вычислительно-информационный Интернет-сервис, который позволяет найти решения уравнений ЫБМ для любой молекулы и ознакомиться с результатами, представленными в наглядном виде. Сервис доступен через глобальную сеть Интернет по адресу http://www.rismproteins.org/online.html.

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

Положения, выносимые на защиту:

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

2. Оценка химического потенциала на основе численного термодинамического интегрирования решений уравнений с отталкивателыюй поправкой мости-кового функционала наилучшим образом описывает термодинамику растворения пептидов.

3. Уравнения RISM с использованием усредненных по выборке состояний макромолекулы внутримолекулярных корреляционных функций позволяют получить термодинамические характеристики сольватации макромолекул, подверженных температурному движению.

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

Степень достоверности и апробация результатов. Основные результаты диссертации докладывались на следующих конференциях: 4th International Symposium on Computational Methods in Toxicology and Pharmacology Integrating Internet Resources (Moscow, Russia, 2007); 40th International School of crystallography (Erice, Italy, 2008); 2-я и 3-я Международная конференция «Математическая биология и бионформатика» (Пущино, Россия, 2008, 2009); 9-я, 10-я, 11-я Пущинская международная школа-конференция молодых ученых (Пущино, Россия, 2005, 200G, 2007)

Исследование выполнено при финансовой поддержке РФФИ. Проект № 12-07-31085-мол_а «Система моделирования гидратации биологических макромолекул с учетом их конфигурационной подвижности» в 2012-1013 гг.

Проект № 10-07-00112-а «Система суперкомпыотерной поддержки научных исследований Пущинского научного центра РАН по физико-химической биологии и нанобиоэлектронике» в 2010-2012 гг.

Публикации. Материалы диссертации опубликованы в 15 печатных работах, из них 6 статей в рецензируемых журналах [1-6] и 9 тезисов докладов.

Личный вклад автора. Представленные в диссертации результаты получены лично соискателем.

Структура и объем диссертации. Диссертация состоит из введения, обзора литературы, 5 глав, заключения, списка сокращений и условных обозначений и библиографии. Общий объем диссертации 123 страницы, из них 102 страницы текста, включая 30 рисунков и 7 таблиц. Библиография включает 113 наименований на 15 страницах.

Содержание работы

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

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

Уравнение ШЗМ впервые предложили Чандлер и Андерсен в 1972 году. Оно является подобным уравнению Орнштейна-Цернике уравнением, связывающим

атом-атомные корреляционные функции, и записывается как интегральное уравнение в матричной форме:

Н = WCW + pWCH. (1)

Символы Н, W и С обозначают матрицы корреляционных функций. Элементами матрицы Н являются полные корреляционные функции (hai), которые суть сдвинутые на 1 радиальные функции распределения h,n (г) = дпу (г) — 1 . Уравнение (1) по сути есть определение прямых корреляционных функций са7, которые являются элементами матрицы С. Элементы матрицы W являются внутримолекулярные корреляционные функции, которые задают связи между атомами одной молекулы. Для жестких молекул с длинами связи Ьа1 они задаются через дельта-функции:

w„7 (г) = 6ау + (1 - <5а7) <5 (г - Lay).

Если корреляционные функции заданы в координатном пространстве, то умножение матриц подразумевает операцию свертки. Однако обычно уравнение (1) записывают в импульсном пространстве. В этом случае элементами матриц Н и С являются Фурье-образы соответствующих корреляционных функций.

Чтобы найти полные и прямые корреляционные функции, требуется другая связь между ними кроме (1), которую называют уравнением замыкания. Поскольку диаграммное разложение для RISM не известно, используют уравнения из теории простых жидкостей:

hay = ехр (-/Зиау + hay -Сау + В [hay, са7]) - 1. (2)

Конкретный вид уравнения замыкания зависит от выбора функционала В. Широко применяются замыкания Перкуса Иевика (PY), гиперцепное замыкание (HNC) и частично-линеаризованное гиперцепное замыкание (PLHNC).

Чтобы составить представление о поведении модели связанных силовых центров, необходимо получить зависимость решений от параметров:

Г = Г(р,Т) |r=const,

ю' г 6) РШИС

О 0.002 0.004 0.006 0.008 0.01 0.012 0.014 0.016 Р.А '

0 0.002 0.004 0.006 0.008 0.01 0.012 0.014 0.016 р,А 1

Рис. 1. Параметрические зависимости решений уравнений ЯШМ для метана: а) в замыкании РУ, б) в РЬНМС. Вверху изображены изотермические зависимости обратной сжимаемости от плотности. Внизу спинодаль (красная) и линии разрыва первых производных обратной сжимаемости по плотности или температуре (черные).

которая является решением уравнения движения по параметру р:

(1Г _ гдЪ

1р ~ ~ 2 V

(3)

где Ъ нелинейный оператор обозначающий систему уравнений (1) и (2).

Специфика (3) заключается в том, что в точках бифуркации решений системы Z (Г, р) = 0 матрица Якоби 3%, становится особенной. Для того чтобы преодолеть трудности, связанные с бифуркацией решений при обращении в нуль производных ¿¿Г/с?р, перейдем к зависимости уравнений от (х^1,?1), и будем численно интегрировать ту систему, для которой на очередном шаге приращение параметра меньше.

Изложенный метод был применен для расчетов изотерм обратной нормированной сжимаемости метана. Расчетами показано, что наличие или отсутствие бифуркаций решений существенно зависит от выбора уравнения замыкания.

С уравнением РЫШС система не имеет множественных решений на всей

плоскости параметров (рис. 1, б). Такое поведение аналогично поведению модели простых жидкостей и качественно описывает обратную сжимаемость метана. По-видимому, эта особенность позволяет применять уравнение PLHNC для описания биологических молекул любых размеров, в то время как при использовании других замыканий, например HNC, с ростом размера системы возникают трудности со сходимостью численных методов.

В случае других замыканий выявлены дополнительные особенности в поведении уравнений RISM. В частности с замыканием PY в области малых температур и плотностей появляются нефизические решения (рис. 1, а). Следовательно, это замыкание непригодно для изучения молекулярных жидкостей методом RISM.

Наиболее сложная картина возникает в случае HNC. Бифуркация решений возникает при температуре выше критической. При этом в точке бифуркации происходит излом изотермы. Такое поведение не наблюдаются для простых жидкостей, и является особенностью модели RISM в гиперцепном замыкании. Дополнительно, функции дин (г) могут быть отличны от нуля в области малых радиусов, что нехарактерно для парных функций распределения.

Результаты первой главы опубликованы в работах [3, 7 9].

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

Huv = WUCUV (Wv + pHvv-) _ (4)

где матрица внутримолекулярных корреляционных функций Wu = {uJa¡¡ (fc)} задает жесткие связи длиной L„/3 между атомами а и /? растворенной молекулы в виде дельта-функций в пространстве Фурье:

п\ A" -Í1 А N sinfc¿Q/3 (к) = д„р + (•!■— <wJ —ту-■

KLag

Поэтому, для описания равновесного состояния растворения подвижной молекулы необходимо усреднить решения уравнений ШБМ по ансамблю.

На практике можно усреднять корреляционные функции, решив уравнения ШЯМ в пределе бесконечного растворения (4) для каждой конфигурации некоторой репрезентативной выборки из конфигурационного пространства макромолекулы при заданной температуре. Выборку состояний в свою очередь можно получить с помощью имитационного моделирования, например, рассчитав методом молекулярной динамики отрезок равновесной траектории.

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

Прежде всего, следует воспользоваться симметрией молекулы растворителя. Так, для случая чистой воды, этот подход позволяет уменьшить размер вектора неизвестных на 1/3.

Значительной экономии вычислительных ресурсов можно добиться, если построить численную схему на основе неточного метода Ныотона-Рафсона и решать линеаризованную систему одним из итерационных методов подпространств Крылова. Такая схема позволяет реализовать безматричный метод, значительно сокращающий требования к оперативной памяти, потому что не требует хранения матрицы Якоби в явном виде. При этом итерационные методы решения линейных систем позволяют находить линейное приращение значительно быстрее, так как вдали от решения достаточно даже грубой оценки линейной прибавки.

Также было замечено, что в итерационном процессе численного решения системы уравнений ШБМ Фурье-образы корреляционных функций сильно меняются только вблизи нуля. Это позволяет применить одновременно неточный метод Ньютона-Рафсона на «грубой» сетке и метод простых итераций для точек дополняющих «грубую» сетку до «точной».

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

Ч,- (к) = ¿у + (п (кт - к)+ Н (RT - lij) ~

- Н (кт - к) Н (RT - /,,)) (1 -

где Я (ж) - функция Хевисайда.

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

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

В качестве объекта исследования выбран короткий пептид окситоцин. Всего было рассчитано три разных равновесных траектории пептида, каждая длиной 250 пс. Первая траектория (А) есть движение пептида в вакууме, вторая (В) -движение пептида в явном водном окружении. Третья траектория (С) - движение пептида в неявной водной среде, когда электростатические взаимодействия моделируются обобщенным методом Ворна и неполярные - методом доступной растворителю поверхности (GBSA).

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

и

кулы, взятых из соответствующей траектории с равномерным шагом по времени.

Каждое состояние характеризуется энергией Гиббса, которая является суммой молекулярно-механической энергии (и) и избыточного химического потенциала {А/л), усредненной по мгновенным конфигурациям:

В теории интегральных уравнений предложено несколько выражений для А/1, которые сделаны исходя из разных вариантов и приближений модели связанных силовых центров: в теории Гауссовых флуктуаций (вР), формула Сингера-Чанд-лера (БС), в теории парциальных волн (Р\У) и приближенные формулы полученные с помощью термодинамической теории возмущений для отталкивательных поправок бридж-функционала в виде отталкивателыюго члена потенциала Ле-нарда-Джонса (ТРТ/г12) и потенциала Викса Чандлера Андерсена (ТРТ/ШСА). Сравним результаты для этих выражений свободной энергии и разных замыканий, используя в качестве контроля метод СВЭА, параметризованный по расчетам малых пептидов.

Все выражения избыточного химического потенциала разделяются на две группы, в которых наблюдается схожее качественное поведение (рис. 2). Первую группу образуют выражения СР и Р\У. Особенностью этих выражений является отсутствие (вР) или модификация (Р\¥) слагаемого, пропорционального квадрату прямой корреляционной функции. Эти выражений дают качественно неверные результаты, показывая, что минимальной свободной энергией обладает пептид в вакууме. При этом результаты, полученные по ним, существенно зависят от выбора замыкания.

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

Рис. 2. Средние значения свободной энергии Гиббса трех состояний окситоцина. О ~ НИС замыкание, □ - РЬНМС замыкание.

Результаты второй главы опубликованы в работах [2, 10, 11].

В третьей главе рассмотрены уравнения, которые связывают полные и прямые корреляционные функции, описывающие среднее распределение молекул растворителя, вокруг молекулы растворенного вещества, подверженной температурным флуктуациям. Идея заключается в том, чтобы решать обычные уравнения ШЭМ в пределе бесконечного растворения (4), но при этом использовать матрицу внутримолекулярных корреляционных функций, усредненную по всем мгновенным конфигурациям макромолекулы:

= + (1 - 5ар)

(5 {г

., эт/сг

Уравнения ШБМ являются формальным определением полных парных корреляционных функций через прямые функции. Уравнение с усредненной матрицей внутримолекулярных корреляционных функций является другим подобным

-180

ТРТ/\УСА

-100

Таблица 1. Разность истинной средней энергии Гиббса и энергии Гиббса, полученной по решениям усредненного уравнения (ккал/моль). В скобках приведена дисперсия энергии Гиббса, полученной по мгновенным геометрическим конфигурациям

БС вР Р\¥ ТРТ/г"12 ТРТ/ШСА

А -20.4 (17.0) 18.0 (8.4) 12.6 (7.6) -20.0 (17.5) -19.4 (17.0)

В -12.2 (10.9) 12.5 (9.5) 9.0 (8.8) -12.0 (12.5) -11.7 (11.7)

С -12.2 (12.8) 14.7 (14.3) 11.2 (13.2) -11.3 (14.0) -11.2 (13.4)

определением.

В исходной модели ШБМ внутримолекулярные парные корреляционные функции являются дельта-функциями, полная матрица которых задает жесткую конфигурацию молекулы. Усреднение таких матриц по каноническому ансамблю по определению дает парные внутримолекулярные корреляционные функции подвижной молекулы. Поэтому можно ожидать, что новые уравнения будут связывать полные и прямые корреляционные функции, описывающие среднее распределение молекул растворителя вокруг подвижной молекулы. В пользу нового определения также говорит то, что формальный вывод уравнений МЭМ сделан без учета конкретной формы внутримолекулярных корреляционных функций.

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

Сравнение подходов выполним на пептиде окситоципе, воспользовавшись траекториями и расчетами по методу ШЭМ для жестких молекул, которые описаны в главе 2. : ■ ■

Разности энергии Гиббса, полученной усреднением мгновенных значений, и энергии, полученной но решения усредненного уравнения, сопоставимы с диспер-

Рис. 3. Энергии Гиббса трех состояний окситоцина, полученные по решениям усредненного уравнения ЫЭМ. Результаты термодинамического интегрирования отмечены пунктиром. О

ШС, □ РЦШС.

сией мгновенных значений (табл. 1). Различные выражения избыточного химического потенциала разделились на две группы точно так же, как и в главе 2. Так разности, вычисленные по формулам БС, ТРТ/У"12 и ТРТ/\¥СА, отрицательны, а по формулам СР и Р\¥ - положительны. Разница непосредственно корреляционных функций также сопоставима с дисперсией в районе первых пиков даже для наиболее отличающихся функций.

Качественная картина энергии Гиббса трех состояний пептида окситоцнна полностью аналогична картине, полученной в главе 2. Дополнительно, вычислены значения энергии Гиббса с отталкивательными поправками, полученными методом термодинамического интегрирования (ШТ/г 12 и ШТ/\¥СА). Если для варианта с отталкивательной частью потенциала Ленарда-Джонса приближенное и точное значения сопоставимы, то для варианта потенциала Викса- Чанд-лера-Андерсена величина изменилась примерно на 70 100 ккал/моль.

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

Результаты третьей главы опубликованы в работах [5, 12].

В четвертой главе выполнены оценки относительных энергий связывания препарата 4',6-диамидино-2-фенилиндол (БАР1) с последовательностью с!(ССА АТТСС^СС в специфичных сайтах связывания.

БАР1 обычно связывается с АТ-богатой последовательностью в малом желобе, где стесненное окружение стабилизирует связывание. Поэтому можно было бы ожидать, что в додекамере <1(ССССААТТСС)2 он будет локализован на ААТТ последовательности. Однако экспериментально в кристаллической структуре обнаружено, что БАР1 взаимодействует с последовательностью АТТС в додекамере

(1(ССССААТТСС)2.

Методом МЭМ были исследованы равновесные молекулярно-динамические траектории двух комплексов ДНК/ОАР1. Комплексы образованы последовательностью (ССААТТСС)2СС и ОАР1, локализованном на сайтах ААТТ и АТТв.

Полученные энергии Гиббса перехода от сайта ААТТ к сайту АТГС различны для различных способов оценки свободной энергии (табл. 2). Видно, что оценки разделись на две группы. В первой группе находятся формула ЭС и выражения для учета отталкивательных поправок по термодинамической теории возмущений. Оценки с помощью этих функционалов дают отрицательную разницу энергий и согласуются с экспериментальными данными. Вторую группу составляют функционалы СР и Р\У, которые дают положительные оценки. Положительная

Таблица 2. Свободная энергия связывания ОАР1 с сайтами ААТС и ААТТ в малом желобе ДНК, ккал/моль

Функционал ЙС СР РШ ТРТ/г"12 ТРТ/ШСА

без ионов

-2289.68 -5742.04 -4892.90 -3934.82 -3337.68

Сдатт -2284.92 -5751.85 -4901.72 -3931.28 -3333.66

Д<3 -4.76 9.81 8.82 -3.54 -4.02

с ионами

блАТС -3410.05 -6969.10 -6359.92 -5095.10 -4486.57

С ААТТ -3403.05 -6986.13 -6378.89 -5092.43 -3333.66

АС -7.00 17.03 18.96 -2.67 -2.32

разность означает, что предпочтительным сайтом связывания является ААТТ, и противоречит эксперименту.

Картина не меняется принципиально с учетом противоионов, кроме того, что разница между отдельными методами оценки свободной энергии Гиббса становится еще более выраженной.

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

Результаты четвертой главы опубликованы в работах [4, 13].

В пятой главе описана структура и возможности разработанного Интернет-сервиса для теоретического исследования гидратации биополимеров методом ин-

тегральных уравнений теории жидкостей (http: //www. rismproteins. org/online. html).

Сервис осуществляет полный цикл изучения структуры биополимера с заданной геометрией и силовым полем: готовит необходимые входные данные, решает уравнения RISM, по различным выражениям избыточного химического потенциала вычисляет атомные парциальные вклады в свободную энергию гидратации и заносит результаты в базу данных (БД). Данные в БД доступны через веб-интерфейс. Представление данных разработано с учетом опыта работы с результатами расчетов гидратации молекул и сделано максимально наглядным.

Результаты пятой главы опубликованы в работах [1, 14, 15].

Выводы:

1. Впервые выполнен параметрический анализ уравнений RISM и найдены точки бифуркации для различных замыканий на примере молекулы метана. Показано, что наиболее стабильным является частично-линеаризованное гиперцепное замыкание, что позволяет выбрать его для изучения макромолекул в растворе.

2. Показано, что химический потенциал пептида окситоцина в методе RISM полученный путем термодинамического интегрирования с учетом отталки-вательной поправки наилучшим образом согласуется с феноменологической моделью, параметризованной по расчетам малых пептидов.

3. Предложен новый метод, который позволяет расчитывать структурные и термодинамические свойства сольватации макромолекул, усредненные по ее состояниям. Метод основан на включении данных о флуктуациях межатомных расстояний внутри макромолекулы в уравнения RISM. Предложенный метод имеет большие перспективы для исследования макромолекул со значительной конформационной подвижностью,- таких как белки, так как позволяет работать с единственным решением усредненных уравнений.

4. Интегральные уравнения теории жидкостей впервые применены для изучения термодинамики связывания биологически активного вещества (4',6-ди-амидино-2-фенилиндола) с макромолекулой (ДНК). Уравнения RISM позволили учесть влияние растворителя точнее континуальных подходов и получить наименьшее значение свободной энергии для наблюдаемого в эксперименте места связывания.

5. Предложены параллельные алгоритмы решения уравнений RISM большой размерности, использующие физико-химические особенности макромолекул и сочетающие численные методы, эффективно работающие в разных интервалах пространственной размерности. Алгоритмы реализованы в виде вычислительных программ и Интернет-сервиса.

Список публикаций

1. Sobolev Е. V., Sobolev О. V., Tikhonov D. A. Online resource for theoretical study of hydration of biopolymers // SAR and QSAR in Environmental Research. 2008. Vol. 19, no. 3-4. P. 303-315.

2. Тихонов Д. А., Соболев E. В. Оценки энергии Гиббса гидратации по молеку-лярно-динамическим траекториям методом интегральных уравнений теории жидкостей в приближении RISM // Журнал физической химии. 2011. Т. 85, № 4. С. 732 738.

3. Соболев Е. В., Тихонов Д. А. Численное исследование сингулярности интегральных уравнений теории жидкостей в приближении RISM // Компьютерные исследования и моделирование. 2010. Т. 2, № 1. С. 51-G2.

4. Соболев Е. В., Тихонов Д. А., Фридман X., Труонг Т. Применение метода RISM для оценки свободной энергии связывания 4',6-диами-дино-2-фенилиндола в малом желобе ДНК по молекулярно-динамической тра-

ектории // Математическая биология и биоинформатика. 2010. Т. 5, № 2. С. 98-113. • . ,

5. Тихонов Д. А., Соболев Е. В. Усредненный по молекулярным траекториям метод интегральных уравнений в приближении RISM // Математическая биология и биоинформатика. 2010. Т. 5, № 2. С. 188-201.

6. Тихонов Д. А., Соболев Е. В. Метод псевдосредних функций в теории RISM. Температурная зависимость гидратации пептида окситоцина // Математическая биология и биоинформатика. 2010. Т. 5, № 2. С. 202-214.

7. Sobolev Е. V., Tihonov D. A. Solutions of integral equation in the theory of molecular liquids near phase transition region // Математическая биология и бионформатика: II Международная конференция, г. Пущино, 7-13 сентября 2008 г.: Доклады / Под ред. В. Д. Лахно. М.: МАКС Пресс, 2008. С. 50-51.

8. Соболев Е. В., Тихонов Д. А. Множественность решений уравнений RISM/HNC в пределе малой плотности // III Международная конференция «Математическая биология и бионформатика». 10-15 октября 2010 г., Пущино: материалы конференции / Под ред. В. Д. Лахно. М.; МАКС Пресс, 2010. С. 77-78.

9. Соболев Е. В., Тихонов Д. А. Расчет водного растворителя в методе RISM // БИОЛОГИЯ — НАУКА XXI ВЕКА: 10-я Пущинская международная школа-конференция молодых ученых, посвященная 50-летию Пущинского научного центра РАН (Пущино, 17-21 апреля 2006 года). Сборник тезисов. Пущино: 2006.-17-21 апреля. С. 348-349.

10. Соболев Е. В., Тихонов Д. А. Параллельный алгоритм численного решения уравнений RISM в пределе бесконечного разбавления // . III Международная конференция «Математическая биология и бионформатика». 10-15 октября

2010 г., Пущино: материалы конференции / Под ред. В. Д. Лахно. М.: МАКС Пресс, 2010. С. 79-80.

11. Lucka S., Sobolev Е. V., Tihonov D. A. Software for study of macromolecule hydration by the methods of integral equations of theory of liquids // Математическая биология и бионформатика: II Международная конференция, г. Пущино, 7 13 сентября 2008 г.: Доклады / Под ред. В. Д. Лахно. М.: МАКС Пресс, 2008. С. 162-163.

12. Тихонов Д. А., Соболев Е. В. Псевдосредние в методе интегральных уравнений теории жидкостей в приближении RISM // III Международная конференция «Математическая биология и бионформатика». 10-15 октября 2010 г., Пущино: материалы конференции / Под ред. В. Д. Лахно. М.: МАКС Пресс, 2010. С. 53-54.

13. Соболев Е. В., Тихонов Д. А. Теория жидкостей в изучении свойств растворов биологических макромолекул // БИОЛОГИЯ - НАУКА XXI ВЕКА: 10-я Пущинская международная школа-конференция молодых ученых, посвященная 50-летию Пущинского научного центра РАН (Пущино, 17-21 апреля 2006 года). Сборник тезисов. Пущино: 2006.— 17-21 апреля. С. 349.

14. Sobolev Е. V., Sobolev О. V., Tikhonov D. A. Online resource for theoretical study of hydration of biopolyiners // Fourth International Symposium on Computational Methods in Toxicology and Pharmacology Integrating Internet Resources. Book of Abstracts. Moscow: 2007. — September 1-5. P. 180.

15. Соболев О. В., Соболев Е. В., Тихонов Д. А. Гидратационный микроскоп // БИОЛОГИЯ — НАУКА XXI ВЕКА: 11-я Пущинская международная школа-конференция молодых ученых (Пущино, 29 октября-2 ноября 2007 года). Сборник тезисов. Пущино: 2007.— 29 октября-2 ноября. С. 61.

Подписано в печать: 20.05.2013

Заказ №8538 Тираж-100 экз. Печать трафаретная. Объем: 1 усл.п.л. Типография «11-й ФОРМАТ» ИНН 7726330900 115230, Москва, Варшавское 36 (499)788-78-56 . www.autoreferat.ru

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

Федеральное государственное бюджетное учреждение науки Институт математических проблем биологии Российской академии наук

На правах рукописи

04201358131

Соболев Егор Васильевич

Интегральные уравнения теории жидкостей в теоретическом изучении биологических макромолекул и их взаимодействий в растворах

03.01.02 - Биофизика

ДИССЕРТАЦИЯ на соискание ученой степени кандидата физико-математических наук

Научный руководитель

кандидат физико-математических наук

Тихонов Дмитрий Анатольевич

Пущине - 2013

Содержание

Введение ......................................................................5

Обзор литературы..........................................................10

1. Задачи и методы вычислительной молекулярной биологии ... 10

2. Интегральные уравнения теории жидкостей в биофизических задачах................................................................15

3. Парная корреляционная функция и уравнение Орнштейна-Цер-нике....................................................................18

4. Формализм связанных силовых центров ..........................28

5. Предел бесконечного растворения..................................34

6. Определение термодинамических свойств по парным корреляционным функциям ....................................................37

Глава 1. Численное исследование сингулярности интегральных

уравнений теории жидкостей в приближении ШБМ ..........42

1.1. Введение..............................................................42

1.2. Редукция по симметрии..............................................43

1.3. Частный случай двухцентровых жидкостей ......................44

1.4. Метод продолжения по параметру..................................46

1.5. Переход к зависимости от обратной сжимаемости и температуры 48

1.6. Схема численного интегрирования..................................49

1.7. Зависимость обратной сжимаемости метана от температуры и плотности..............................................................51

1.8. Выводы к первой главе..............................................54

Глава 2. Оценки энергии Гиббса методом ШЭМ по молекулярно-

динамическим траекториям..........................................55

2.1. Введение..............................................................55

2.2. Уравнения гидратации макромолекулы............................56

2.3. Переход к системе нелинейных алгебраических уравнений ... 59

2.4. Алгоритм численного решения дискретной системы..............61

2.5. Ослабление дальних внутримолекулярных корреляций..........63

2.6. Параллельный алгоритм............................................65

2.7. Оценки энергии Гиббса пептда окситоцина в рамках континуального подхода......................................................68

2.8. Оценки энергии Гиббса пептида окситоцина методом RISM . . 69

2.9. Выводы ко второй главе ............................................73

Глава 3. Уравнения сольватации молекул, подверженных тепловому движению..........................................................74

3.1. Введение..............................................................74

3.2. Уравнения RISM, усредненные по конфигурациям растворенной молекулы..............................................................75

3.3. Обрезание гистограмм внутримолекулярных корреляционных функций ....................................................................77

3.4. Гидратация пептида окситоцина по решениям усредненного уравнения RISM ..........................................................79

3.5. Выводы к третьей главе..............................................83

Глава 4. Анализ связывания 4',6-диамидино-2-фенилиндола в

малом желобе ДНК ....................................................85

4.1. Введение..............................................................85

4.2. Процедура и параметры моделирования ..........................87

4.3. Свободная энергия связывания DAPI и DNA в воде..............89

4.4. Анализ корреляций различных вкладов в энергию Гиббса ... 91

4.5. Выводы к четвертой главе..........................................93

Глава 5. Интернет-сервис для теоретического изучения гидрата-

ции биополимеров........................................................95

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

5.2. Описание работы сервиса............................................96

5.3. Использование ресурса для расчета гидратации пептида ТИР-Cage....................................................................100

5.4. Выводы к пятой главе................................................104

Заключение....................................................................105

Список сокращений и условных обозначений ......................107

Литература....................................................................109

Введение

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

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

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

Важный класс интегральных уравнений, описывающих равновесную структуру молекулярных жидкостей в терминах атом-атом парных корреляционных функций, предложили Чандлер и Андерсон в 1972 году. Эти уравнения основаны на модели связанных силовых центров (RISM, от англ. Reference Interaction Site Model).

Цели и задачи диссертационной работы. Цель работы заключается в исследовании особенностей модели связанных силовых центров при описании структуры и термодинамики сольватации биологических макромолекул и в создании методики учета растворителя на основе интегральных уравнений теории жидкостей.

Были поставлены следующие задачи:

1. Выполнить параметрический анализ уравнений ШЭМ для плотного метана и провести сравнение поведения системы при использовании различных уравнений замыкания.

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

3. Сравнить метод ШБМ и континуальные подходы в задаче анализа термодинамики трех модельных состояний пептида окситоципа.

4. Сравнить решения усредненных по выборке из канонического ансамбля уравнения ШБМ и усредненных решений жестких уравнений ШБМ для каждого состояния из этой выборки.

5. Выполнить анализ термодинамики связывания 4',6-диампдино-2-фенил-индола в малом желобе ДНК с целью выявить энергетически предпочтительный сайт связывания.

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

7. Разработать вычислительно-информационный Интернет-сервис для теоретического изучения гидратации биологических макромолекул.

Научная новизна. Предложены новые алгоритмы, позволившие решать уравнения ШБМ для множества мгновенных состояний биологической макромолекулы, состоящей из тысяч атомов. Алгоритмы реализованы в виде компьютерных программ и Интернет-сервиса.

Модель связанных силовых центров впервые применена для учета влияния растворителя в задаче о нековалентном связывании в растворе биологической молекулы (ДНК) и активного вещества (4',6-диамидино-2-фенилиндола) и показано преимущество интегральных уравнений перед континуальными подходами.

Впервые показано, что эффекты сольватации и термодинамику макромолекул можно описать с помощью уравнений ШБМ с усредненными по состояниям макромолекулы внутримолекулярными корреляционными функциями.

Впервые выполнен параметрический анализ уравнений ШБМ.

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

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

Разработанными вычислительными программами можно пользоваться через вычислительно-информационный Интернет-сервис, который позволяет найти решения уравнении ШБМ для любой молекулы и ознакомиться с результатами, представленными в наглядном виде. Сервис доступен через глобальную сеть Интернет по адресу http://www.rismproteins.org/online.html.

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

Положения, выносимые на защиту:

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

2. Оценка химического потенциала на основе численного термодинамического интегрирования решении уравнений с отталкивательной поправкой мостикового функционала наилучшим образом описывает термодинамику растворения пептидов.

3. Уравнения RISM с использованием усредненных по выборке состояний макромолекулы внутримолекулярных корреляционных функций позволяют получить термодинамические характеристики сольватации макромолекул, подверженных температурному движению.

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

Степень достоверности и апробация результатов. Основные результаты диссертации докладывались на следующих конференциях: 4th International Symposium on Computational Methods in Toxicology and Pharmacology Integrating Internet Resources (Moscow, Russia, 2007); 40th International School of crystallography (Erice, Italy, 2008); 2-я и 3-я Международная конференция «Математическая биология и бионформатика» (Пущино, Россия, 2008, 2009); 9-я, 10-я, 11-я Пущинская международная школа-конференция молодых ученых (Пущино, Россия, 2005, 2006, 2007)

Исследование выполнено при финансовой поддержке РФФИ.

Проект № 12-07-31085-мол_а «Система моделирования гидратации биологических макромолекул с учетом их конфигурационной подвижности» в 2012-2013 гг.

Проект № 10-07-00112—а «Система сунеркомпьютерной поддержки научных исследований Пугцинского научного центра РАН по физико-химической биологии и нанобиоэлектронике» в 2010-2012 гг.

Публикации. Материалы диссертации опубликованы в 15 печатных работах, из них 6 статей в рецензируемых журналах [1-6] и 9 тезисов докладов.

Личный вклад автора. Представленные в диссертации результаты получены лично соискателем.

Структура и объем диссертации. Диссертация состоит из введения, обзора литературы, 5 глав, заключения, списка сокращений и условных обозначений и библиографии. Общий объем диссертации 123 страницы, из них 102 страницы текста, включая 30 рисунков и 7 таблиц. Библиография включает 113 наименований на 15 страницах.

Обзор литературы

1. Задачи и методы вычислительной молекулярной биологии

Фундаментальная проблема молекулярного подхода в биологии заключатся в том, чтобы описать живые системы в терминах физики и химии. Компьютерные исследования мезоскопических систем, представляющих интерес с точки зрения биологии, начали развиваться относительно недавно. Модели классической механики адекватно описывают большинство свойств таких систем, а моделирование молекулярной динамики (МД) является наиболее важным теоретическим методом, применяемым в таких исследованиях. Первая статья о моделировании динамики белка ингибитора - трипсина поджелудочной железы быка (ВРТ1) - была опубликована около 30 лет назад [7].

Метод МД применяется для решения трех типов задач. В первом случае имитационное моделирование используют как средство получения точек в конфигурационном пространстве. Прежде всего, это задачи определения и уточнения структур биологических молекул, полученных в эксперименте, например, методами рептгепоструктурной кристаллографии. Во втором случае вычисляют средние в термодинамическом равновесии, в том числе структурные и динамические свойства (например, среднеквадратичные отклонения атомов) и термодинамику системы. Для этих задач, как и для первого случая, необходимо получать адекватные точки из конфигурационного пространства, но теперь с дополнительным условием, что бы они представляли репрезентативную выборку для заданной температуры. Третье направление использует имитационное моделирование для изучения непосредственно динамики. Здесь уже не достаточно репрезентативной выборки из конфигурационного пространства, необходимо получать верное представление о развитии системы во времени. Для первых двух приложений равнозначные результаты дают метод МД [8] и

и

метод Монте-Карло [9]. Для третьего случая подходит только метод МД.

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

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

Альтернативный подход заключается в том, чтобы неявно рассматривать растворитель, учитывая только его влияние. Неявные методы позволяют избавиться от имитационного моделирования динамики большого числа молекул растворителя, опираясь на теоретические знания об их поведении в целом. Такие схемы оценки свободной энергии сольватации могут быть полезными, оставаясь при этом нетребовательными в вычислительном плане. Неявные модели растворителя свободны от статистических ошибок связанных с нерепре-зентатпвностью выборок. Наконец, неявные представления растворителя могут быть очень полезными концептуальные инструментами, чтобы лучше по-

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

Однако можно ли строго выразить влияние растворителя, не учитывая явно его степени свободы? Проинтегрируем функцию распределения некоторой системы в термодинамическом равновесии по всем состояниям растворителя. Редуцированная функция распределения не будет зависеть координат растворителя, хотя будет включать среднее влияние растворителя на растворенное вещество. Это среднее влияние называется потенциалом средней силы. Фундаментальную концепцию потенциала средней силы ввел Кирквуд [10], чтобы описать равновесную структуру жидкостей. Потенциал средней силы равен свободной энергии растворения вещества в заданной конфигурации и не равен средней потенциальной энергии. Он не зависит от степеней свободы растворителя и никакой информации о влиянии растворителя на равновесные свойства системы не теряется.

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