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

Моделирование проводимости ионных каналов на основе методов молекулярной и броуновской динамики

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

Автореферат диссертации по теме "Моделирование проводимости ионных каналов на основе методов молекулярной и броуновской динамики"

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

Турченков Дмитрий Александрович

Моделирование проводимости ионных каналов на основе методов молекулярной и броуновской

динамики

Специальность 03.01.02 - Биофизика

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

11 ФЬЙ 2015

Москва - 2014

005559044

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

Научный руководитель: Быстрое Владимир Сергеевич, доктор физико-

математических наук

Официальные оппоненты: Алиев Рубин Ренатович, доктор физико-математи-

ческих наук, заведующий лабораторией электрофизиологии Федерального государственного образовательного учреждения высшего профессионального образования «Московский физико-технический институт»

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

Ведущая организация: Федеральное государственное бюджетное учреждение

науки Институт биофизики клетки Российской академии наук

Защита состоится 5 марта 2015 г. в 17:00 на заседании диссертационного совета Д 501.002.11 при Московском государственном университете имени М.В.Ломоносова, расположенном по адресу: 119991, ГСП-1, Москва, Ленинские горы, д.1, стр. 2, физический факультет МГУ, ЦФА.

С диссертацией можно ознакомиться в библиотеке Московского государственного университета имени М.В.Ломоносова.

Автореферат разослан « 2 » 2015 г.

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

кандидат технических наук

Сидорова Алла Эдуардовна

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

Актуальность темы исследования. Ионные каналы играют одну из ключевых ролей в жизнедеятельности клетки, принимая участие в процессах генерации и распространения нервного импульса, опосредования мышечных сокращений, регуляции ионного обмена, сигнальной функции и др. Любые незначительные структурные нарушения ионного канала могут привести к изменению его основных физиологических характеристик, определяющих его биологическую роль в клетке — уровне проводимости и селективности. К настоящему времени уже выявлено более 27 каналопатий (Roulesau, 2011), к которым относятся гнпо- и гиперкалиемическип периодические параличи, различные формы эпилепсии и миастенические синдромы (Kass, 2005). Более того, недавние ис-следоваши семейства пуринергических рецепторов показали не только нх исключительную важность в опосредовании сократительных ответов клетки, но и выявили их участие в механизмах пролиферации, дифференцировки и апо-птоза. Уже известно, что нарушения в работе отдельных каналов данного типа явлются причиной различных нейродегенративных расстройств, таких как болезнь Паркинсона, Альцгеймера и др (Burnstock, 2008). Поэтому совершенно очевидно, что изучение механизмов работы ионных каналов и опосредуемых ими биологических эффектов представляет не только научный интерес, но и является важным аспектом на пути к созданию новых узкоспециализированных лекарственных форм.

С этой точки зрения исследование структуры и функций ионных каналов, а так же происходящих в них процессов, имеет смысл только вместе с рассмотрением осуществляемых ими макроскопических функций, что требует привлечения современных, и во многом междисциплинарных научных дисциплин, методов и подходов. К настоящему времени разработано достаточно большое количество различных методологий изучения ионных каналов. Основным экспериментальным методом на данный момент является методика локальной фиксации потенциала «patch clamp», обладающая наибольшей чувствительностью, позволяющей получать характерные величины проводимости одиночных каналов. Однако, проведение подобных экспериментов задействует достаточно дорогостоящее оборудование и применяется исключительно в научных исследованиях в силу специфики подготовки объекта исследования.

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

можно выделить методы квантовой химии, молекулярной и броуновской динамики, а также элементы электродиффузионной теории (МаЯео, 2012). Несмотря на значительные достижения, полученные с использованием данных методов, каждый из них имеет ряд существенных ограничений, что делает практически невозможным их применение для макроскопического описания такой системы, как ионный канал на мембране. Это вызвано несколькими факторами. Первая причина в том, что для проведения моделирования методами квантовой химии или молекулярной динамики нам необходимо знать пространственную структуру ионного канала. Учитывая, что мембранные белки тяжело поддаются кристаллизации (а значит проведение рентгеноструктурного анализа проблематично), применение данных методов возможно только на некоторых модельных, хорошо изученных системах, таких как ацетилхолиновый рецептор, калиевый КсвА канал. Вторым сдерживающим фактором является размер системы: так, только КсвА состоит из более чем 15 ООО атомов, а интегрированный в фосфоли-пидный бислой с явным учетом растворителя размер системы будет составлять более 100 ООО атомов. Для таких систем вычисления методами квантовой химии даже на современных суперкомпьютерах невозможны. Применение же методов молекулярной динамики способно только смоделировать единичные акты прохода иона, и говорить о биологически значимых временных эволюциях системы в данном случае не приходится. Если говорить о применении электродиффузи-оннон теории и уравнений Пуассона-Нернста-Планка, то, несмотря на потенциально существенно большую временную эволюцию системы, которую можно получить данной методологией, ее применение ограничивается заложенным диффузионным механизмом переноса, который не наблюдается для большего числа ионных каналов.

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

Цели и задачи диссертационной работы: Разработка комбинированного алгоритма моделирования проводимости ионного канала и его применение для оценки величины трансмембранных токов на примере ионных каналов Р2Х2, Р2Х4 и Р2Х7 типа.

Для достижения указанной цели были поставлены следующие задачи:

• Создать математическую модель описания движения ионов и нейромеди-аторов в вязкой среде.

• Учесть зависимость вязкости и диэлектрической проницаемости раствора от концентрации растворенных электролитов.

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

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

• На основе разработанной методологии получить характерные величины ионных токов, селективности и проводимости ионных каналов на примере пуринергических рецепторов Р2Х2, Р2Х4 и Р2Х7 типа, сравнить полученные результаты с «patch clamp» экспериментами.

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

Данные алгоритмы реализованы в интерактивной студии моделирования, с возможностью ЗО-визуалнзации и вычислениями на графических видеокартах с использованием технологии NVIDIA CUDA. Интерактивность и высокая

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

На основе разработанных алгоритмов и ПО была впервые предсказана конфигурация селективного фильтра пуринергического рецептора Р2Хг типа с неизвестной на данный момент пространственной структурой, объясняющая его высокую проводимость и слабую селективность к катионам (Ding, 1999), обнаружен диффузионный механизм проводимости данного канала. Теоретическая и практическая значимость.

Разработанное программное обеспечение позволяет проводить компьютерные «patch clamp» эксперименты различной конфигурации и получать интересующие исследователя величины селективности и проводимости ионного канала без использования дорогостоящей аппаратуры и трудоемких экспериментов. Благодаря высокой производительности возможна длительная временная эволюция системы, что позволяет отслеживать изменения концентраций ионов и нейромедиаторов в отдельных компартментах рассматриваемой системы. Эти данные могут представлять высокую диагностическую значимость в клинической практике и существенно упрощать разработку новых лекарственных форм, целью которых являются различные ионные каналы. Кроме того, комбинированные алгоритмы моделирования, представленные в работе, могут лечь в основу разработки высокопроизводительных методов для суперкомпьютеров, которые позволят смоделировать целый участок нейрона и механизм ионной проводимости на молекулярном уровне, что необходимо для создания компьютерных моделей как отдельных функционирующих элементов, так и всего мозга в целом.

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

• Разработана новая математическая модель, описывающая движение ионов и нейромедиаторов в вязкой среде на основе комбинирования методов молекулярной и броуновской динамики.

• Установлено, что значение среднего квадрата скорости и перемещения броуновской частицы для разностных схем Эйлера и Хейна зависит от размера шага, что приводит к отклонениям от температуры термостата и изменению коэффициента диффузии. Область применимости существующих на данный момент основных разностных схем численного интегрирования уравнения Ланжевена лежит в диапазоне 7At < 1.

• Показано, что для концентраций электролита и растворе С > 0.5 М/л в классической модели Ланжевена необходимо учитывать изменения диэлектрической проницаемости и вязкости раствора.

• Установлено, что при размере частиц радиусом R < 1.5 А необходимо вводить поправки к закону Стокса на днэлектричекое трение.

• Предложена новая разностная схема численного интегрирования уравнения Ланжсисна в пространстве координат и скоростей с учетом скоррели-ровапности стохастических приращений на каждом итерационном шаге, не имеющая ограничений на шаг интегрирования, с асимптотическими значениями среднего квадрата скорости и перемещения, соответствующими точному решению.

• На основе разработанного метода удалось объяснить високос значение проводимости Р2Хг капала наличием Акрзда в области селективного фильтра. Механизм проводимости данного канала диффузионный. Показано, что наличие остатков Ser;«n и Ser3i2 в области селективного фильтра Р2Х7 канала (согласно данным моделирования по гомологии) приводит к величинам проводимости, отличным от экспериментальных данных.

Степень достоверности и апробация результатов. Основные результаты диссертации докладывались па научной сессии НИЯУ МИФИ (Москва. 2012); 15-ом научном симпозиуме международной исследовательской группы по системной биологии (Амслапд, 2012); двух семинарах IV-oro съезда Биофизиков России (Нижний Новгород, 2012); семинаре кафедры биофизики биологического факультета МГУ (Москва. 2012) XII и XlII-ofl ежегодной международной молодежной конференции ИБХФ РАН-Вузы (Москва, 2012.2013); международном симпозиуме вычислительного и теоретического моделирования межмолекулярных взаимодействий (Дубна, 2013); международной конференции актуальных вопросов современных физико-математических наук (2014); V международной конференции по математической биологии и биоипформатикс (Пущино, 2014); международной конференции вопросов современной биологии, физики и химии (2014).

Публикации. Материалы диссертации опубликованы в 19 печатных работах, из них G статей в рецензируемых журналах [1-6], из них 4 статьи в российских и зарубежных научных журналах, рекомендованных ВАК РФ для публикации материалов кандидатских и докторских диссертаций; 5 статей в сбор-

никах трудов [7 11] н 0 тезисов докладов на всероссийских и международных конференциях [12 17]; 2 авторских свидетельства РФ на разработку программного обеспечения для ЭВМ [18, 19].

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

Структура и объем диссертации. Диссертация состоит из введения, 6 глав, включающих обзор литературы и 5 глав авторских исследований, заключения, 3 приложений и списка литературы. Работа изложена на 150 страницах машинописного текста, включает 43 рисунка и 16 таблиц. Библиография включает 306 наименований на 19 страницах.

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

В последнее время именно разработка гибридных методов компьютерного моделирования находит все большую популярность среди научного сообщества. Начало разработке комбинированной методологии моделирования проводимости ионных каналов было заложено в конце 90-х годов XX века года в работах (Allen и Chung, 1999), которые впервые применили методы молекулярной и броуновской динамики для моделирования ионной проводимости калиевого канала. В данной работе методами МД был подсчитан эффективный коэффициент диффузии ионов Na+ и К+ в канале, который затем использовался в БД моделировании ионной проводимости на .основе упрощенной модели ионного капала. Несколько ранее в конце 70-х годов была предложена (Warshel, Karplus и Levitt) гибридная методология QM/MM моделирования, которая широко применяется для изучения отдельных элементов данной системы (например, возможность денротонирования аминокислотных остатков аснартата Asp«o н глутаминовой кислоты Glu71 в селективном фильтре К es А ионного канала), однако применение данной методологии для полноценного моделирования ионного канала ограничено современными вычислительными способностями. Из трудов российских ученых необходимо отметить работу Вороиовского С.Е., который впервые совместил методы БД и стохастическую активацию иониых каналов па основе известных констант скоростей ферментативных реакций. Более подробно существующие на данный момент методы изучения рассмотрены в Главе 1 диссертации.

В основе используемой в работе модели ионного капала па мембране лежит представление данной системы в виде трех комнартментов (Рис. 1): два из них представляют собой вне- (I компартмент) и внутриклеточную (II компарт-мент) среду, в которой растворены ионы и нейромедиаторы с соответствующей концентрацией, а второй компартмент представляет собой гидрофобный участок мембраны с ионным каналом. Необходимо отметить, что в данной работе элементы II компартмента считаются неподвижными и их характеристики не меняются в ходе моделирования (латеральная диффузия липидов, флин-флон эффект, а так же температурные флуктуации ионного канала не учитываются). Это приводит к тому, что заряженные частицы, формирующие данную биологическую систему разбиваются па типа:

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

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

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

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

где 7 = ^ удельный коэффициент трения, внешняя сила, действую-

щая па частицу, а стохастическая сила представлена через амплитуду силы а и винеровекпй процесс IV(Ь).

Для определения внешнего поля сил РГ(Ь) в уравнении (1), действующего на ьую частицу со стороны всей системы используется суммарный потенциал, состоящий из электростатического и парного потенциала взаимодействия типа Леннард-Джонса:

F (t)

dV = -7 Vdt + + adW (t)

(1)

11 -V 1 Vi-Я ,

и tot - У ----Mejj

^ 4я-£со Гц

Jf' L J

(2)

где e,j и <7;j рассчитываются по правилу Лоренца-Бертло:

gj + Qj

и

I II III

Рис. 1. Схематическое изображение используемой в работе комиартментной модели ионного канала, подчиняющейся комбинированному моделированию методами молекулярной (ком-партмент II) и броуновской динамики (компартменты I и III).

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

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

ныи коэффициент трении £ в таком случае состоит из двух составляющих: гидродинамического и диэлектрического (¡р трения. На основании закона Стокса и теории диэлектрического трения Адельмана коэффициент трения £ записывается в виде:

Су '---- '

Си

где </ — заряд частицы, Л радиус, £ — диэлектрическая проницаемость среды, е^ — высокочастотная диэлектрическая проницаемость воды, 7д — время релаксации (Дебая), 5 параметр, определяемый законом Стокса и зависящий от размера частицы (6 для частиц, размер которых превышает размер молекулы растворителя, и 4, если их размеры соизмеримы).

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

11 = щТТвх (4)

где X мольная доля электролита в растворе, А и В некоторые константы, подбираемые эмпирически.

__^ и; , г *

Здесь х мольная доля /-ого иона в растворе, 1Х ионная сила, выраженная через мольную долю, а

С,- = сп + саТ

Д = <1л + (¡¡2Т

где с,], с,-2 определяются из экспериментальных данных, а йа и ¿¡2 для широкого спектра электролитов принимают постоянные значения 1.44 • 10® и —1389 К-1.

Следующим важным аспектом является тот факт, что для большинства разностных схем интегрирования уравнения Ланжевена, применяемых на данный момент (это методы Эйлера (Euler), Хейна (Heun) и неявной средней точки (implicit midpoint) была обнаружена зависимость среднего квадрата скорости (и2) и перемещения (х2) частицы [3] от размера временного шага интегрирования. Подробнее эти вопросы обсуждаются в Главе 2 диссертации. Это является причиной отклонения от равновесной температуры термостата и изменения коэффициентов диффузии описываемых частиц. Для устранения данного недостатка мы использовали полученную нами разностную схему с учетом скорре-лированности стохастических приращений:

Ц = + Fe— --+ N ( 0, -

ту

1 - e~lAt Xi = Xi-1 + Vi-1-+ Fe

2(1 - е_2">л')

27

7

my

+N 0,

my*

27 д t

(6)

где нормировочная константа а находится из предположения того, что система находится в термодинамическом (тепловом) равновесии:

2укТ

а коэффициент корреляции р^):

(7)

туохог,

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

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

В качестве модели растворителя была выбрана универсальная 3-х частичная модель воды SPC, движение описывалось классической схемой Эйлера, а взаимодействие между частицами определялось потенциалом, состоящим из электростатического и парного потенциала взаимодействия типа Леннард-Джонса (2).

Помимо описания движения заряженных частиц в области ионной поры методы молекулярной динамики применялись для разрешения следующей проблемы. Как отмечалось выше, использование методов броуновской динамики для задач биологического моделирования параметризуется выбором усредненного действия растворителя на частицу. В классической модели Ланжевена таким параметром раствора является коэффициент трения вид которого зависит от используемого приближения. Однако, если в случае простых частиц (таких, как ионы Na+, К+, С1~) подобные приближения могут дать приемлемый результат, то в случае многоатомных молекул (NH£, ситуация заметно усложняется. Это связано с неоднозначностью выбора аппроксимирующей геометрии, распределения заряда для данных теорий. Поэтому единственной возможностью определить коэффициент диффузии D и трения £ таких сложных частиц—это использование прямого МД моделирования.

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

о)р} <8)

U0=oj=o )

где N — число исследуемых частиц в системе, и — число шагов моделирования,

причем t = п ■ At.

Другой метод предполагает построение автокорреляционной функции скорости (VACF) и использование соотношения Грина-Кубо :

<4

(V(0)-V(t))dt= 1

Шп

п N

+ (9)

fo=0 J=0

В свою очередь, для определения коэффициента трения £ броуновской используется автокорреляционная функция стохастической силы (ГАСР) :

с = 1

3 кьТ

(Fs(t0)F3(t))dt (10)

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

Разработанная методология комбинированного моделирования применялась для исследования ионной проводимости пуринергических рецепторов Р2Х2, Р2Х4 и Р2Х7 рецепторов. Выбор данных ионных каналов в качестве объекта исследования не случаен: среди семейства Р2Х каналов на сегодняшний день только для Р2Х4 получена пространственная структура (Hattorf, 2012) . Для Р2Х2 ионного канала только недавно была получена 3D реконструкция данных электронной микроскопии (Mio, 2009) , а пространственную организацию Р2Х7 рецептора подобрали по гомологии (Jiang, 2013). Поэтому проведенные в ходе данной работы численные эксперименты призваны не только подтвердить работоспособность методологии, но и внести ясность в пространственную организацию селективных фильтров Р2Х2 и Р2Х7 каналов. Ко всему прочему, с данными ионными каналами было проведено достаточно большое количество «patch clamp» измерений (Evans, 1996; Wong, 2000), что представляет собой идеальную экспериментальную базу для имитационного моделирования. Обсуждению данных вопросов и полученных результатов посвящена Глава 5 диссертации.

Для начала были построены упрощенные геометрии данных каналов. В качестве входных данных мы использовали PDB структуру Р2Х4 в открытом

1G

состоянии (PDB: 4DW1), идентичную геометрию для Р2Х7 (с отличием только в аминокислотных остатках в селективном фильтре), а геометрия Р2Х2 рецептора была получена на основе пространственной реконструкции электронной микроскопии.

Для проверки чувствительности модели к конфигурации селективного фильтра, для каждого ионного канала были определены заряженные или полярные аминокислотные остатки, которые могут быть повернуты в область ионной поры. Поскольку для Р2Х4 рецептора известна пространственная организация, то моделирование «нативной» и измененной геометрии канала позволит ответить на данный вопрос. Согласно PDB структуре данного канала, в области селективного фильтра отсутствуют какие-либо заряженные или полярные аминокислотные остатки. Однако поворот второго трансмембранного домена (ТМ2) по часовой стрелке на 35° развернет остаток аспарагиновой кислоты ASP357 в область поры, что приведет к изменению потенциала вдоль оси канала и к возможному увеличению уровня проводимости. Поскольку для Р2Х7 конфигурация селективного фильтра получена моделированием по гомологии, то вопрос о его пространственной организации остается открытым, в связи с чем были рассмотрены все 7 различных конфигураций фильтра, образованного остатками Эегззэ, Ser342 и Asp352. Анализ аминокислотной последовательности ТМ2 Р2Х2 рецептора выявил возможное участие остатков ТЬгззэ, Ser345 и ASP349 в формировании 7 различных конфигураций.

Согласно предлагаемой методологии, аминокислотные остатки и фосфоли-пиды формируют систему фиксированных зарядов II компартмента. В отдельных работах, посвященных изучению влияния распределения заряда вдоль оси канала на величины проводимости (Boda at al.) заряд аминокислотных остатков аппроксимируют величиной в ±0.5 или ±1 заряда электрона, однако предпосылки, заложенные в таком предположении, остаются далеко не очевидными. Поэтому в данной работе для вычисления распределения заряда по данным молекулам мы использовали один из наиболее популярных методов квантовой химии — метод функционала плотности DFT. Данные аминокислоты рассматривались в депротонированном состоянии. Сначала была проведена предварительная геометрическая оптимизация молекулы в вакууме. Затем электронная плотность основного состояния для данной геометрии была спроецирована на положения атомов методом Хиршфельда. Использовался локальный функционал в форме PWC, неограниченные по спину волновые функции, двойной численный базнс с учетом валентных орбиталеи для атомов водорода. Ради-

ус обрезания 3.7 А. Проводился полноэлектроннын расчет без использования псендопотенциалов. Порог сходимости для геометрической оптимизации составил 10_,> Хартрн, 0.005 Á по расстоянию. Относительная точность сходимости па этапе расчета самосогласованного поля 10~°. Мультипольное разложение до октуполя. Для улучшения сходимости применялся DIIS. Вычисления проводились с использованием пакета Accelrys Materials Studio, модуль DMol. Полученный по нашим расчетам заряд рассматриваемых в работе аминокислотных остатков составил: Asp (q = -0.73е, R = 1.7 Ä), Ser (q = -0.43е, R = 1.2 А) и Thr (q = —0.50c, R = 2.1 Á). Для каждой молекулы определялся размер bounding-сферы, аппроксимирующий полярный/ заряженный участок молекулы, используемый в нашей модели (Рис. 1).

Аналогичные вычисления были проведены для молекул основных мем-бранообразующнх фосфолинидов в модели растворителя COSMO. Результаты представлены в Таблице 1. Учитывая, что «patch clamp» измерения проводились на реальных биологических тканях, состав фосфолипидного бислоя, используемый и моделировании соответствовал фосфолипидному составу тела нейрона крысы: PC 28 %, РЕ 20%, SPH 4%, PI 0% и PS 4 %, что соответствует поверхностной плотности заряда на внутренней стороне мембраны а и 40 мКл/'м2.

Состав комиартментов I и III соответствовал «patch clamp» про токолам измерении одиночных каналов данного типа. Для P2X,i: [NaCl\¡ = 147, [NaCl]nl = 5, [NaF}m = 140 (мМ/л), внешний потенциал -150 мВ. Р2Х2: [NaCl\¡ = 145, [NaF]n] = 5, внешний потенциал -100 мВ. Р2Х7: [CsC6Hu07]r = 150, [CsCeHn07]

ш = внешний потенциал —110 мВ. Параметры моделирования системы были: размер ячейки 300 А, температура 298° К, шаг броуновской динамики Дtßp = 0.1 не, шаг молекулярной динамики At^o = 1 не, время моделирования 1 мкс.

Зная число попов, прошедших через ионный канал за время моделирования системы, легко определить величины проводимости G и амплитуды тока I. Для чистоты модельного эксперимента концентрации растворенных веществ в компартмептах I и III были фиксированы: после прохождения через канал нон заново помещался в свой стартовый комиартмеит со случайными координатами. В Таблице 2 представлены результаты численного моделирования проводимости ионных каналов Р2Хо, Р2Х i и Р2Х7 типа при различных конфигурациях селективных фильтров. Как видно, конфигурация P2Xi канала, при которой ни один заряженный или полярный аминокислотный остаток не повернут в

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

Липид Заряд, qe v, A3 <R>, A

Phosphatidylcholine (PC) +0.76 -0.78 286 ± 40

Phosphatidylethanolamine (РЕ) +0.37 -0.37 205 ± 31

Phosphatidylserine (PS) -0.84 247 ± 34

Phosphatidylinositol (PI) -0.92 261 ± 37

Phosphatidylglycerol (PG) -0.93 205 ± 30 4

Cardiolipin (CL, DPG) -1.6 391 ± 50

Sphingomyelin (SPH) +0.85 -0.83 230 ± 37

Lysophosphatidylcholine (LPC) +0.78 -0.78 286 ± 40

область ионной поры («без фильтра», БФ) наилучшим образом (9.0 пСм) согласуется с экспериментальными данными («patch clamp» измерениями и PDB структурой), тогда как наличие остатка аспарагиповой кислоты ASP357 увеличивает величину проводимости канала почти в два раза (16.3 пСм). Согласно структуре Р2Х7 канала, предсказанной моделированием по гомологии, селективный фильтр формируют только два остатка серина Эегззд и Ser^, обращенных в область ионной поры. Однако, величина проводимости для данной конфигурации (14.3 пСм), полученным на основе нашей модели, будет значительно отличаться от экспериментальных данных (8.3 пСм), тогда как конфигурация БФ находится в согласии с данными. Отсюда мы делаем вывод о возможном отличии реальной пространственной организации Р2Х7 канала от предсказанной поворотом ТМ2 домена па угол в 15-25° по часовой стрелки, во избежание наличия остатков данных аминокислот в области ионной поры. Моделирование Р2Хг рецептора позволило объяснить его высокий уровень проводимости (32 пСм) наличием остатка аспарагиповой кислоты Asp3ig в области ионной поры.

Поскольку атомистическая структура Р2Х2 и Р2Х7 рецепторов на настоя-

Таблица 2. Величины проводимости ионных каналов Р2Х2, Р2Х4 и Р2Х7 тина при различных конфигурациях селективных фильтров.

Канал Проводимость, G (пСм) Эксп.

1'2Х4 9.0 БФ ± 1.2 1-Я 10.5 ± 1.5

Р2Х7 Р2Х2 8.0 ВФ ± 0.9 Sn*» 11.1 ± 1.7 12.1 SL4\I42 ± 1.5 15.9 ' ± Акрам 2.1 Sern и 11.3 ± 1.9 Scrras 18..'i ± Aspt.vj 2.1 19.7 Smm ± Asp4'^ 2.7 Scr.4.4!) 2 I.Ii Sana ± Asp;K2 2-Я 8.3 ± (l.li

1 1.1 Пф ± 1.3 ТЬгзлэ 27.0 ± 2.1 25.1 Ксгл., ± 2.4 л о.'2 ± ASP349 2.1 ТЬгззо :ш.2 2.7 TIlTOQ II ± Asp-jio 2.4 :t9.2 ■^--r.i 15 ± Asp.li() 2.4 Tlll'339 17 Аярм» 2.8 32 2.0

щий момент отсутствует, для проверки пашей гипотезы были использованы их вольт-амиерные характеристики. Моделирование проводилось с шагом внешнего потенциала п 20 мВ в диапазоне —20 —120 мВ для Р2Х7 и —20 -.—100 мВ для Р2Х2 канала при различных конфигурациях селективных фильтров. Параметры моделирования не менялись. Как видно (Рис. 2) наилучшим образом экспериментальные зависимости описывают состояния «без фильтра» для Р2Хт и Aspaig для Р2Х2, что согласуется с нашими предположениями.

Исследование избирательности (селективности) капала к отдельным ионам позволяет охарактеризовать его механизм проводимости. Для этих целен мы провели моделирование Р2Х2 с Asp319 конфигурацией селективного фильтра при внешнем потенциале в —120 мВ для случая различных электролитов одинаковой концентрации (120 мМ/'л) в I и III компартментах: KCl, NaCl, LiCl, RbCl и CsCl. Полученные величины ионных токов и проводимости представлены в Таблице 3.

Таблица 3. Амплитуда тока I и проводимость G ионного канала P2Xtj тина для различных катионов

Ion Li+ Na+ K+ Rb+ C.s+

I, nA -3.52 -3.92 -4.52 -4.3G -4.38

G, HCM 29 33 38 36 37

Как видно, величины проводимости данных ионов возрастают согласно коэффициентам их свободной диффузии: К+ « Св+ и > Ма+ > Ы+. Это свидетельствует о том, что ярко выраженной селективности к катионам у данного капала нет, а отношение уровней проводимости для отдельных катионов позволяет говорить о диффузионном механизме ионного транспорта данного

Мембранный потенциал, мВ

Рис. 2. Вольт-амперные характеристики Р2Хг (сверху) и Р2Х7 ионных каналов (снизу) для различных конфигураций селективных фильтров.

канала. Это косвенно подтверждается и достаточно большим (11 А в диаметре) размером ионной поры.

В заключении отметим, что все компьютерные эксперименты в данной работе были проведены с использованием авторского программного пакета интерактивной студии моделирования Patch Clamp Simulation, реализующей описываемую в данной работе методологию мезоскопического моделирования проводимости ионных каналов (Рис. 3).

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

• полная интерактивноеть добавление компонентов системы может осу-

ществляется прямо во время моделирования;

• высокая производительность достигается благодаря использованию комбинированной методологии моделирования и реализации вычислительного ядра на графических видеокартах;

• длительная эволюция системы моделируемое время может достигать нескольких миллисекунд;

• масштабируемость системы позволяет с легкостью проводить моделирование системы из десятков различных ионных каналов;

• возможность проводить «patch clamp» эксперименты без использования дорогостоящего лабораторного оборудования;

Рис. 3. Patch Clamp Simulation. Вид моделируемой системы: растворенные ионы во внутри-п внеклеточных компартментах, фосфолипиды на внутренней стороне мембраны, желтым цветом выделена область ионный норы (молекулы воды не отображены).

Подробные характеристики, возможности разработанного ПО, и используемые методологии рассмотрены в Главе 6.

Выводы

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

• Показано, что асимптотические значения среднего квадрата скорости и перемещения броуновской частицы для разностных схем Эйлера и Хейна зависят от размера шага. Выявлена область применимости (jAt < 1) существующих на данный момент основных разностных схем численного интегрирования уравнения Ланжевена.

• Показана необходимость учета зависимостей диэлектрической проницаемости и вязкости раствора от концентраций растворенных электролитов при С > 0.5 М/л и введения поправок к закону Стокса на диэлектрическое трение для частиц, радиусом R < 1.5 Á.

• Предложена новая разностная схема численного интегрирования уравнения Ланжевена в пространстве координат и скоростей с учетом скоррели-ровапностн стохастических приращений па каждом итерационном шаге, не имеющая ограничений на шаг интегрирования, с асимптотическими значениями среднего квадрата скорости и перемещения, соответствующими точному решению.

• Показано, что существующая пространственная структура Р2Х7 рецептора, полученная моделированием по гомологии, приводит к величинам проводимости, отличным от экспериментальных данных, что свидетельствует об отличии реальной пространственной организации Р2Х7 канала от предсказанной поворотом ТМ2 домена на угол в 15-25° по часовой стрелки, во избежание наличия остатков Ser339 и Ser3i2 в области ионной поры.

• Объяснено высокое (30.2 ± 2.0 pS) значение проводимости Р2Х2 канала наличием Aspa» в области селективного фильтра, выявлена избирательность к катионам: К+ « Cs+ « Rb+ > Na+ > Lí+, свидетельствующая о диффузионном характере проводимости данного канала.

Список работ, опубликованных по теме диссертации

1. Turchenkov D. A., Bystrov V. S. Conductance Simulation of the Purinergic P2X2, P2X4, and P2X7 Ionic Channels Using a Combined Brownian Dynamics and Molecular Dynamics Approach // The Journal of Physical Chemistry B. 2014. Vol. 118, no. 31. P. 9119-9127.

2. Турченков Д. А., Бороновсюш С. E., Нарциссов Яр. Р. Моделирование диффузии ионов в синаитической щели с использованием стохастической модели Ланжевена в приближении диэлектрического трения // Биофизика. 2013. Т. 58, № б. С. 1013-1021.

3. Турченков Д. А., Турченков М. А. Анализ упрощения разностных схем для уравнения Ланжевена, влияние учета корреляции приращений // Компьютерные исследования и моделирование. 2012. Т. 4, № 2. С. 325-338.

4. Турченков Д. А., Быстров В. С. Экспериментальные и теоретические методы изучения ионных каналов // Математическая биология и биоинформатика. 2014. Т. 9, № 1. С. 112-148.

5. Шайтан К. В., Шайтан А. К., Турченков Д. А. и др. Алгоритмы и методы исследования трехмерных атомистических моделей молекул белков на основе анализа картины рассеяния мощного рентгеновского лазерного излучения // Наноструктуры. Математическая физика и моделирование. 2013. Т. 9, № 2. С. 33-74.

6. Turchenkov D. A. Determination of System Diffusion Coefficient // Journal of Nature Science and Sustainable Technology. 2013. Vol. 7, no. 4. P. 327-333.

7. Turchenkov D. A., Turchenkov M. A. Simulation of conductivity of ionic channel purinergic P2X2 receptor using combined Brownian and molecular dynamic approach // Computational and Theoretical Modeling of Biomolecular Interactions / Ed. by у Rubin A. B. et al.; Moscow-Izhevsk: Institute of Computer Science. 2013. P. 77-79.

8. Turchenkov D. A. Determination diffusion coefficient of Brownian particles using velocity or force autocorrelation function in molecular dynamic simulations // Quantitative Chemistry, Biochemistry and Biology. Steps Ahead / Ed. by у Gennady E. Zaikov et al.; New-York: Nova Publishers. 2013. P. 271-277.

9. Турченков Д. А., Быстров В. С., Турченков М. А. Анализ конфигурации селективного фильтра ионных каналов методами квантовой химии // Актуальные вопросы естественных наук: физика, химия, биология / Под ред. проф. Обухова А.Г.; Москва. 2014. С. 35-38.

10. Турченков Д. А., Быстров В. С., Турченков М. А. Применение методов квантовой химии для анализа молекул фосфолипидов в биологических мембранах // Актуальные вопросы современных физико-математических наук / Под ред. проф. Сынзыныса Б. И.; Москва. 2014. С. 63-66.

11. Turchenkov D. A., Bystrov V. S. Implementation of multiscale modeling technique of ionic channel conductance using combined QM/MD/BD dynamics approach // Доклады V Международной конференции «Математическая биология и биоинформатика» / Под ред. проф. Лахно В .Д.; Москва: МАКС Пресс. 2014. С. 118-119.

12. Турченков Д. А., Бороиовскин С. Е., Нарциссов Яр. Р. Моделирование диффузии ионов с использованием стохастического интегрирования уравнения Ланжевена / / Сборник трудов научной сессии НИЯУ МИФИ, Москва. 2012. Т. 2. С. 137.

13. Boronovskiy S. Е., Turchenkov D. A., Nartsissov Y. R. Modeling of ion diffusion in water solutions based on methods of Langevin dy-namics // Abstract Book 15th Workshop of the International Study Group for Systems Biology, Ameland. 2012. P. 42.

14. Турченков Д. А., Бороновский С. E., Нарциссов Яр. Р. Моделирование диффузии ионов в стохастической модели Ланжевена с использованием приближения диэлектрического трения // Материалы докладов IV съезда биофизиков России, Нижний Новгород. 2012. Т. 1. С. 293.

15. Турченков Д. А. Метод вычисления коэффициента трения броуновских частиц на основе автокорреляционной функции скорости или силы при моделировании биомембран методами молекулярной динамики // Материалы докладов IV съезда биофизиков России, Нижний Новгород. 2012. Т. 1. С. 294.

16. Турченков Д. А., Турченков М. А. Определение коэффициента диффузии броуновских частиц на основе автокорреляционной функции скорости и силы методами молекулярной динамики // Материалы XII ежегодной международной молодежной конференции ИБХФ РАН-Вузы, Москва. 2012.

17. Турченков Д. А., Быстров В. С. Моделирование проводимости ионного канала на примере пуринергического рецептора Р2Х2 типа комбинированными методами броуновской и молекулярной динамики // Материалы XIII ежегодной международной молодежной конференции ИБХФ РАН-Вузы, Москва. 2013.

18. Турченков Д. А., Турченков М. А. Интерактивная студия моделирования Patch Clamp Simulation (PCS) // Свидетельство государственной регистрации программы для ЭВМ: № 2013618398. 2013.

19. Шуров Д. Л., Шайтан А. К., Турченков Д. А. и др. Программный комплекс реконструкции пространственной структуры белков и комплексов на основе карт электронной плотности низкого разрешения // Свидетельство государственной регистрации программы для ЭВМ: № 2013618398. 2013.

Подписано в печать 23.12.2014 г. Бумага офсетная. Печать цифровая. Формат А4/2. Усл. печ. л.1. Заказ № 214. Тираж 100 экз. Типография «КОПИЦЕНТР» 119234, г. Москва, Ломоносовский пр-т, д.20 Тел. 8(495)213-88-17 \vw\v.autoreferat 1 .ги