Биология Метод Монте-Карло

редактировать

Биология Методы Монте-Карло (BioMOCA) были разработаны в университете из Иллинойса в Урбана-Шампейн для моделирования переноса ионов в среде электролита через ионные каналы или нанопоры, встроенные в мембраны. Это трехмерный симулятор на основе частиц Монте-Карло для анализа и изучения проблемы переноса ионов в системах ионных каналов или подобных нанопорах во влажных / биологических средах. Смоделированная система состоит из белка, образующего ионный канал (или искусственные нанопоры, такие как углеродная нанотрубка, CNT), с мембраной (то есть липидным бислоем), разделяющей две ионные ванны с обеих сторон. BioMOCA основана на двух методологиях, а именно: перенос Больцмана Монте-Карло (BTMC) и частицы-частицы-частицы-сетка (PM). Первый использует метод Монте-Карло для решения уравнения Больцмана, в то время как второй разделяет электростатические силы на ближнюю и дальнюю составляющие.

Содержание

1 Предпосылки
- 1.1 Решение для электростатического поля
- 1.2 Диэлектрический коэффициент
- 1.3 Анизотропная диэлектрическая проницаемость
2 Расчеты
- 2.1 Дискретизация интегрирования ящика
- 2.2 Размер иона
- 2.3 Ионно-белковое взаимодействие
- 2.4 Ионно-водное взаимодействие
- 2.5 Позиционно-зависимая диффузия
- 2.6 Гидратные оболочки
3 Условия и методы
- 3.1 Граничные условия
- 3.2 Метод многосеточной фокусировки
- 3.3 EMF и DBF
- 3.4 Визуализация с использованием VMD
- 3.5 Запись траекторий в двоичном формате
4 См. Также
5 Ссылки

Фон

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

Ансамбль ионов в области моделирования распространяется синхронно во времени и в трехмерном пространстве путем интегрирования уравнений движения с использованием схемы "чехарда" второго порядка точности. Положения ионов r и силы F определены на временных шагах t и t + dt. Скорости ионов определены как t - dt / 2, t + dt / 2. Основные разностные уравнения движения:

v → (t + dt 2) = v → (t - dt 2) + F → (t) dt {\ displaystyle {\ vec {v}} (t + {\ frac {dt} {2}}) = {\ vec {v}} (t - {\ frac {dt} {2}}) + {\ vec {F}} (t) \, dt}

{\ vec {v}} (t + {\ frac {dt} {2}}) = {\ vec {v}} (t - {\ frac {dt} {2}}) + {\ vec {F}} (t) \, dt

r → (T + dt) = r → (t - dt) + v → (t + dt 2) dt {\ displaystyle {\ vec {r}} (t + dt) = {\ vec {r}} (t-dt) + {\ vec {v}} (t + {\ frac {dt} {2}}) \, dt}

{\vec {r}}(t+dt)={\vec {r}}(t-dt)+{\vec {v}}(t+{\frac {dt}{2}})\,dt

где F - сумма электростатических и парных сил ион-ионного взаимодействия.

Решение для электростатического поля

электростатический потенциал вычисляется через регулярные интервалы времени путем решения уравнения Пуассона

∇ (ε (r) ∇ ϕ ( р, t)) знак равно - (ρ ионов (г, T) + ρ перм (г)) {\ Displaystyle \ набла (\ varepsilon (г) \ набла \ фи (г, т)) = - (\ rho _ { \ text {ion}} (r, t) + \ rho _ {\ text {perm}} (r))}

\nabla (\varepsilon (r)\nabla \phi (r,t))=-(\rho _{{\text{ions}}}(r,t)+\rho _{{\text{perm}}}(r))

где $ρ ion (r, t) {\ displaystyle \ rho _ {\ text {ion}} (r, t)}$ $\rho _{{\text{ions}}}(r,t)$ и $ρ perm (r) {\ displaystyle \ rho _ {\ text {perm}} (r)}$ $\rho _{{\text{perm}}}(r)$ являются плотность заряда ионов и постоянные заряды на белке соответственно. $ϵ (r) {\ displaystyle \ epsilon (r)}$ $\epsilon (r)$ - локальная диэлектрическая постоянная или диэлектрическая проницаемость, а $ϕ (r, t) {\ displaystyle \ phi (r, t)}$ $\phi (r,t)$ - локальный электростатический потенциал. Решение этого уравнения обеспечивает самосогласованный способ включения приложенного смещения и эффектов зарядов изображения, индуцированных на диэлектрических границах.

Ионы и частичные заряды на остатках белка назначаются конечной прямоугольной сетке с использованием схемы «облако в ячейке» (CIC). Решение уравнения Пуассона на сетке учитывает компонент сетки частиц схемы PM. Однако такая дискретизация приводит к неизбежному усечению короткодействующей составляющей электростатической силы, которую можно скорректировать, вычислив ближний заряд-заряд Кулоновские взаимодействия.

Диэлектрический коэффициент

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

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

Анизотропная диэлектрическая проницаемость

Стало очевидно, что макроскопические свойства системы не обязательно распространяются на масштабы молекулярных длин. В недавнем исследовании, проведенном Резой Тогри, Р. Джеем Машлом и Эриком Якобссоном из Университета Иллинойса, Урбана-Шампейн, они использовали моделирование молекулярной динамики для изучения свойств воды в безликих гидрофобных цилиндрах диаметром от 1 до 12. нм. Это исследование показало, что вода претерпевает определенные изменения в структуре, диэлектрических свойствах и подвижности при изменении диаметра трубки. В частности, они обнаружили, что диэлектрические свойства в диапазоне от 1 до 10 нм сильно отличаются от свойств объемной воды и фактически являются анизотропными по своей природе. Хотя такие невыразительные гидрофобные каналы не представляют собой настоящие ионные каналы, и необходимо провести дополнительные исследования в этой области, прежде чем можно будет использовать такие данные для ионных каналов, очевидно, что такие свойства воды, как диэлектрическая проницаемость внутри ионного канала или нанопоры может быть намного сложнее, чем считалось ранее. В то время как высокая осевая диэлектрическая проницаемость экранирует электростатические заряды иона в осевом направлении (вдоль канала), низкая радиальная диэлектрическая проницаемость увеличивает взаимодействие между подвижным ионом и частичными зарядами или изображениями диэлектрического заряда на канале, обеспечивая более высокую селективность в ионной среде. каналы.

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

Вычисления

Дискретизация блочного интегрирования

Чтобы использовать блочное интегрирование для дискретизации D-мерного уравнения Пуассона

∇ (ε ∇ φ) = ρ {\ displaystyle \ nabla (\ varepsilon \ nabla \ varphi) = \ rho}

\nabla (\varepsilon \nabla \varphi)=\rho

, где $ε {\ displaystyle \ varepsilon}$ $\varepsilon$ является диагональным тензором D × D, это дифференциальное уравнение переформулируется как интеграл уравнение. Интегрируя приведенное выше уравнение по D-мерной области $Ω {\ displaystyle \ Omega}$ $\Omega$ , и используя теорему Гаусса, получаем интегральную формулировку

∮ ∂ Ω n ^ (ε ∇ φ) = - ∫ Ω ρ {\ Displaystyle \ oint _ {\ partial \ Omega} {\ hat {n}} (\ varepsilon \ nabla \ varphi) = - \ int _ {\ Omega} \ rho}

\oint _{{\partial \Omega }}{\hat {n}}(\varepsilon \nabla \varphi)=-\int _{\Omega }\rho

В этом Приложение предполагается, что это двумерный случай. Обновление до трехмерной системы было бы простым и законным, поскольку теорема Гаусса также верна для одного и трех измерений. $ϵ {\ displaystyle \ epsilon}$ $\epsilon$ предполагается заданным в прямоугольных областях между узлами, а $φ {\ displaystyle \ varphi}$ $\varphi$ задано в сетке узлы (как показано на рисунке справа).

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

. Области интегрирования $Ω {\ displaystyle \ Omega}$ $\Omega$ затем выбираются в виде прямоугольников с центром вокруг узла и простирающихся до 4 ближайших соседних узлов.. Затем градиент $∇ φ {\ displaystyle \ nabla \ varphi}$ $\nabla \varphi$ аппроксимируется с использованием центрированной разности, нормальной к границе области интегрирования $Ω {\ displaystyle \ Omega}$ $\Omega$ и среднее значение $ϵ {\ displaystyle \ epsilon}$ $\epsilon$ по поверхности интегрирования $∂ Ω {\ displaystyle \ partial \ Omega}$ $\partial \Omega$ . Этот подход позволяет аппроксимировать левую часть приведенного выше уравнения Пуассона в первом порядке следующим образом:

∮ ∂ Ω n ^ (ε ∇ φ) = φ i + 1, j - φ i, jhix (hjy 2 ϵ i, jx + hj - 1 Y 2 ε я, j - 1 x) {\ displaystyle \ oint _ {\ partial \ Omega} {\ hat {n}} (\ varepsilon \ nabla \ varphi) = {\ frac {\ varphi _ {i + 1, j} - \ varphi _ {i, j}} {h_ {i} ^ {x}}} \ left ({\ frac {h_ {j} ^ {y}} {2}} \ epsilon _ {i, j} ^ {x} + {\ frac {h_ {j-1} ^ {y}} {2}} \ varepsilon _ {i, j-1} ^ {x} \ right)}

\oint _{{\partial \Omega }}{\hat {n}}(\varepsilon \nabla \varphi)={\frac {\ varphi _{{i+1,j}}-\varphi _{{i,j}}}{h_{i}^{x}}}\left({\frac {h_{j}^{y}} {2}}\epsilon _{{i,j}}^{x}+{\frac {h_{{j-1}}^{y}}{2}}\varepsilon _{{i,j-1 }}^{x}\right)

- φ я, j - φ я - 1, jhi - 1 Икс (HJY 2 ϵ я - 1, JX + HJ - 1 Y 2 ε я - 1, J - 1 х) {\ Displaystyle {} - {\ гидроразрыва {\ varphi _ {i, j} - \ varphi _ {i-1, j}} {h_ {i-1} ^ {x}}} \ left ({\ frac {h_ {j} ^ {y}} {2}} \ epsilon _ {i-1, j} ^ {x} + {\ frac {h_ {j-1} ^ {y}} {2}} \ varepsilon _ {i-1, j-1} ^ {x} \ right)}

{}-{\frac {\varphi _{{i,j}}-\varphi _{{i-1,j}}}{h_{{i-1}}^{x}}}\left({\frac {h_{j}^{y}}{2}}\epsilon _{{i-1,j}}^{x}+{\frac {h_{{j-1}}^{y}}{2}}\varepsilon _{{i-1,j-1}}^{x}\right)

+ φ i, j + 1 - φ i, jhjy (hix 2 ε i, jy + hi - 1 x 2 ε i - 1, jy) {\ displaystyle {} + { \ frac {\ varphi _ {i, j + 1} - \ varphi _ {i, j}} {h_ {j} ^ {y}}} \ left ({\ frac {h_ {i} ^ {x}} {2}} \ varepsilon _ {i, j} ^ {y} + {\ frac {h_ {i-1} ^ {x}} {2}} \ varepsilon _ {i-1, j} ^ {y} \ right)}

{}+{\frac {\varphi _{{i,j+1}}-\varphi _{{i,j}}}{h_{{j}}^{y}}}\left({\frac {h_{i}^{x}}{2}}\varepsilon _{{i,j}}^{y}+{\frac {h_{{i-1}}^{x}}{2}}\varepsilon _{{i-1,j}}^{y}\right)

- φ i, j - φ я, j - 1 hj - 1 Y (hix 2 ε i, j - 1 y + hi - 1 x 2 ε я - 1, j - 1 y) {\ displaystyle {} - {\ frac {\ varphi _ {i, j} - \ varphi _ {i, j-1}} {h_ {j-1} ^ {y}}} \ left ({\ frac {h_ {i} ^ {x}} {2 }} \ varepsilon _ {i, j-1} ^ {y} + {\ frac {h_ {i-1} ^ {x}} {2}} \ varepsilon _ {i-1, j-1} ^ { y} \ right)}

{}-{\frac {\varphi _{{i,j}}-\varphi _ {{i,j-1}}}{h_{{j-1}}^{y}}}\left({\frac {h_{i}^{x}}{2}}\varepsilon _{{ i,j-1}}^{y}+{\frac {h_{{i-1}}^{x}}{2}}\varepsilon _{{i-1,j-1}}^{y }\right)

где $ε x {\ displaystyle \ varepsilon ^ {x}}$ $\varepsilon ^{x}$ и $ε y {\ displaystyle \ varepsilon ^ {y}}$ $\varepsilon ^{y}$ - две компоненты диагонали тензора $ϵ {\ displaystyle \ epsilon}$ $\epsilon$ . Дискретизировать правую часть уравнения Пуассона довольно просто. $ρ {\ displaystyle \ rho}$ $\rho$ дискретизируется на тех же узлах сетки, как это было сделано для $φ {\ displaystyle \ varphi}$ $\varphi$ .

∫ Ω i ρ = ​​Volume (Ω я) ρ я {\ displaystyle \ int _ {\ Omega _ {i}} \ rho = {\ text {Volume}} (\ Omega _ {i}) \ rho _ {i}}

\int _{{\Omega _{i}}}\rho ={\text{Volume}}(\Omega _{i})\rho _{i}

Размер иона

Конечный размер ионов учитывается в BioMOCA с использованием парных сил отталкивания, полученных из 6–12 потенциала Леннарда-Джонса. В симуляторе используется усеченно-сдвинутая форма потенциала Леннарда-Джонса для имитации отталкивания ионного ядра. Модифицированная форма парного потенциала Леннарда-Джонса, сохраняющая только отталкивающую составляющую, имеет вид

ULJ (rij) = {4 ϵ LJ ((σ ijrij) 12 - (σ ijrij) 6) + ϵ LJ rij < 2 1 / 6 σ i j 0 r i j>2 1/6 σ ij {\ displaystyle U_ {LJ} (r_ {ij}) = {\ begin {cases} 4 \ epsilon _ {LJ} \ left (\ left ({\ frac {\ sigma _ {ij}) } {r_ {ij}}} \ right) ^ {12} - \ left ({\ frac {\ sigma _ {ij}} {r_ {ij}}} \ right) ^ {6} \ right) + \ epsilon _ {LJ} r_ {ij} <2^{1/6}\sigma _{ij}\\0r_{ij}>2 ^ {1/6} \ sigma _ {ij} \ end {cases}}}

U_{{LJ}}(r_{{ij}})={\begin{cases}4\epsilon _{{LJ}}\left(\left({\frac {\sigma _{{ij}}}{r_{{ij}}}}\right)^{{12}}-\left({\frac {\sigma _{{ij}}}{r_{{ij}}}}\right)^{6}\right)+\epsilon _{{LJ}}r_{{ij}}<2^{{1/6}}\sigma _{{ij}}\\0r_{{ij}}>2 ^ {{1/6}} \ sigma _ { {ij}} \ end {cases}}

Здесь $ϵ LJ {\ displaystyle \ epsilon _ {LJ}}$ $\epsilon _{{LJ}}$ - параметр энергии Леннарда-Джонса, а $σ ij = (σ i + σ j) / 2 {\ displaystyle \ sigma _ {ij} = (\ sigma _ {i} + \ sigma _ {j}) / 2}$ $\sigma _{{ij}}=(\sigma _{i}+\sigma _{j})/2$ - среднее значение отдельных параметров расстояния Леннарда-Джонса для частиц i и j.. Использование усеченной формы потенциала вычислительно эффективно, поскольку Это препятствует наложению или слиянию ионов, что было бы явно нефизическим.

Ионно-белковое взаимодействие

Доступность рентгеновских кристаллографических измерений с высоким разрешением полных молекулярных структур позволяет получить информацию о типе и местоположении всех атомов, образующих белок. В BioMOCA атомы белка моделируются как статические точечные заряды, встроенные в конечный объем, недоступный для ионов, и связанный с определяемым пользователем диэлектрическим коэффициентом. Кроме того, доступен ряд параметров силового поля, которые предоставляют информацию о заряде и радиусах атомов в различных аминокислотных группах. Сочетание молекулярной структуры и силовых полей обеспечивает координаты, радиусы и заряд каждого атома в белковом канале. BioMOCA использует такую информацию в стандартном формате PQR (Position-Charge-Radius) для отображения белковой системы на прямоугольной сетке.

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

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

Взаимодействия ионов с водой

В качестве подхода с уменьшенными частицами BioMOCA заменяет явные молекулы воды континуальным фоном и обрабатывает взаимодействия ионов с водой с использованием метода BTMC, в котором должны быть соответствующие скорости рассеяния. выбрал. Другими словами, траектории ионов случайным образом прерываются событиями рассеяния, которые объясняют диффузное движение ионов в воде. Между этими событиями рассеяния ионы следуют за силами Ньютона. Время свободного полета, T f, генерируется статистически из общей скорости рассеяния согласно

- ln ⁡ (r) = ∫ 0 T f λ (p → (t)) dt {\ displaystyle - \ ln (r) = \ int _ {0} ^ {T_ {f}} \ lambda ({\ vec {p}} (t)) \, dt}

- \ ln (r) = \ int _ {0} ^ {{ T_{f}}}\lambda ({\vec {p}}(t))\,dt

где r - случайное число, равномерно распределенное на единичный интервал. $λ {\ displaystyle \ lambda}$ $\lambda$ , функция от импульса, представляет собой общий коэффициент рассеяния для всех механизмов столкновений. В конце каждого свободного полета скорость иона выбирается случайным образом из максвелловского распределения. Поскольку правильный механизм рассеяния для взаимодействий ионов с водой в растворах негабаритных электролитов еще не разработан, в нашей модели используется позиционно-зависимая скорость рассеяния, связанная с локальной диффузией. Эта зависимость от положения проистекает из того факта, что молекулы воды могут иметь разный порядок организации в разных областях, что влияет на скорость рассеяния.

Позиционно-зависимый коэффициент диффузии

Широко признано, что ионы и молекулы воды не обладают такой же подвижностью или коэффициентом диффузии в ограниченных областях, как в объеме. Фактически, более вероятно, что будет иметь место уменьшение эффективной подвижности ионов в ионных каналах. В методах редуцированных частиц, где вода в канале рассматривается как неявный фон континуума, требуется средняя подвижность ионов, чтобы показать, как ионы могут диффундировать из-за локальных электростатических сил и случайных событий. В симуляциях Транспортного Монте-Карло предполагается, что полная скорость рассеяния ( $λ {\ displaystyle \ lambda}$ $\lambda$ ) является результатом только взаимодействия ионов с водой; это связано с диффузией ионов выражением

λ = k T m D {\ displaystyle \ lambda = {\ frac {kT} {mD}}}

\lambda ={\frac { kT}{mD}}

где m - масса иона, а D - его постоянная диффузии. Как показывает уравнение, пониженный коэффициент диффузии ионов внутри просвета канала приводит к увеличению числа случаев рассеяния.

Гидратные оболочки

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

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

U hy знак равно c 0 ехр ⁡ (c 1 - rc 2) cos (c 3 (c 1 - r) π) + c 4 (c 1 r) 6 {\ displaystyle U_ {hy} = c_ {0} \ exp \ left ({\ frac {c_ {1} -r} {c_ {2}}} \ right) cos (c_ {3} (c_ {1} -r) \ pi) + c_ {4} \ left ({\ frac {c_ {1}} {r}} \ right) ^ {6}}

U_{{hy}}=c_{0}\exp \left({\frac {c_{1}-r}{c_{2}}}\right)cos(c_{3}(c_{1}-r)\pi)+c_{4}\left({\frac {c_{1}}{r}}\right)^{6}

Коэффициенты c i были определены эмпирически для раствора 1 M KCl, используя моделирование методом МД для сравнения функций радиального распределения ионов с расчетами Равновесия Монте-Карло. Эффект гидратных оболочек оказался важным при моделировании при более высоких концентрациях соли, когда наблюдается насыщение проводимости многих ионных каналов, в том числе пориновых, по мере дальнейшего увеличения концентрации соли в электролитных ваннах. Более ранние симуляции, которые не включали модель гидратных оболочек, не воспроизводили поведение насыщения проводимости. Это предполагает дополнительный отталкивающий потенциал, действующий для предотвращения скопления ионов и, следовательно, ограничения концентрации ионов и плотности тока в ограниченном пространстве поры даже при высокой концентрации соли в ванне. Когда был включен потенциал отталкивания, наблюдалась умеренная проводимость канала .

Условия и методы

Граничные условия

Электрические и физиологические свойства ионных каналов экспериментально измеряются путем введения канала в липидную мембрану, разделяющую две ванны, содержащие растворы определенных концентраций.. Постоянное электростатическое смещение прикладывается к каналу за счет погружения электродов в две ванны. Формулировка граничных условий, которые точно представляют эти контактные области, может потребовать чрезвычайно больших областей ванны и является сложной задачей. За пределами дебаевского расстояния от мембраны электростатический потенциал и плотности ионов существенно не меняются. Это предположение было подтверждено результатами континуальных результатов, представленных ранее. Для типичных концентраций соли, используемых при моделировании ионных каналов, длина Дебая имеет порядок 10 Å. Используя предположение, граничные условия Дирихле накладываются на потенциал на двух плоскостях доменной границы, которые поперечны каналу, с учетом того, чтобы эти плоскости были достаточно далеко от мембраны.

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

Многосеточный метод фокусировки

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

В настоящее время в BioMOCA реализованы два решателя Пуассона на основе метода конечных разностей. Один из них использует предварительно подготовленную схему сопряженного градиента (pCG) и используется по умолчанию. Последний заимствован из решателя APBS, который использует схему V-multi-grid. Помимо численного подхода к решению уравнения Пуассона, основное различие между двумя решающими программами заключается в том, как они обращаются к диэлектрической проницаемости в системе. В первом решателе значение диэлектрической проницаемости присваивается каждой ячейке в сетке, в то время как в решающей программе APBS диэлектрические коэффициенты определяются в узлах сетки. Как обсуждалось ранее, в решателе pCG используется метод блочного интегрирования, который наиболее точно обрабатывает уравнение Пуассона. Несмотря на то, что полноценный многосеточный решатель, основанный на методе блочной интеграции, находится в стадии разработки, есть удобный способ повторно использовать уже существующий код и обрабатывать системы ионных каналов.

Моделирование ионных каналов требует наличия больших областей ванны для точной обработки экрана. Наличие таких областей ванны делает область сетки уравнения Пуассона большой и приводит либо к большому количеству узлов сетки с точным разрешением сетки, либо к небольшому количеству узлов сетки с очень грубой дискретизацией. Из объемного моделирования достаточно крупной сетки для описания ванн с использованием схемы PM. Однако в канальной области требуется высокое разрешение из-за сильно заряженной природы этих областей и наличия пространственно изменяющихся диэлектрических областей. Кроме того, основной интерес представляет изучение поведения канала с точки зрения ионной проницаемости, селективности, стробирования, плотности и т. Д. Другими словами, лучше разместить больше вычислительных ресурсов в области канала., и минимум ванн, чтобы снизить общие вычислительные затраты и ускорить моделирование с недель до, возможно, дней. Была разработана схема, основанная на методе фокусировки сетки, которая позволяет удовлетворить требования большой области ванны и точного разрешения сетки в канале одновременно эффективным с вычислительной точки зрения способом. Эта методология способна иметь несколько доменов с мелкими ячейками, которые могут потребоваться для описания нескольких каналов пор, таких как порин OmpF, или массива ионных каналов, разделяющих одни и те же области ванны, или даже с еще более мелкими ячейками внутри мелкой ячейки для относительно больших каналов с узкие проходы для ионов, такие как канал рецептора никотина.

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

ЭДС и DBF

Электродвижущая сила (ЭДС) - это измерение энергии, необходимой заряженной частице, такой как ион, для пересечения ионного канала, встроенного в мембрана. Частично этот потенциальный энергетический барьер обусловлен взаимодействием между пересекающим ионом и постоянными / частичными зарядами на остатках белка. Другая часть происходит от индуцированных диполей в диэлектрической среде белок / мембрана и называется диэлектрической граничной силой (DBF). Чтобы вычислить только DBF, можно отключить все статические заряды на остатках белка и перетащить ион через пору и вычислить энергетический барьер, используя

P D B F = ∫ - d z ^. E → {\ displaystyle P_ {DBF} = \ int -d {\ hat {z}}. {\ Vec {E}}}

P_{{DBF}}=\int -d{\hat {z}}.{\vec {E}}

Важно отметить, что измерения EMF или DBF - это просто качественные измерения, как ион не обязательно пересекает канал через центр своего просвета по прямой линии, и он часто сопровождается другими ионами, движущимися в том же или противоположном направлении, что резко меняет динамику системы. Более того, в отличие от расчетов управляемой МД, где остатки белка динамически перемещаются в качестве иона или ионы, подпрыгивающие через канал, в наших расчетах ЭМП или ДБФ белок моделируется как статический континуум, что в дальнейшем влияет на расчеты энергии более количественным образом. Другой проблемой, которая дополнительно влияет на измерения, является отсутствие молекул гидратации воды, которые движутся вместе с ионом и экранируют часть его заряда. Сказав все вышесказанное, все же вычисление ЭДС или DBF важно для решения проблемы избирательности канала или стробирования. Вычисление любого из этих двух энергетических барьеров доступно в качестве опции в BioMOCA.

Визуализация с использованием VMD

Визуализация VMD молекулы грамицидина 1MAG вместе со структурой, генерируемой BioMOCA, где зеленый цвет представляет белок, красный относится к мембране (т.е. липиду), а фиолетовый - канал, а также левая и правая ванны

VMD был оборудован возможностью загрузки структур BioMOCA. Это очень полезная функция, так как для сравнения можно загрузить как структуру белка (например, файл PDB или PQR), так и структуры, созданные BioMOCA. На рисунке справа показано, как BioMOCA создал структуру для канала грамицидина с обернутой вокруг него мембраной. Кроме того, BioMOCA также выгружает траектории ионов в стандартных форматах, чтобы их можно было позже загрузить в инструменты молекулярной визуализации, такие как VMD, и просмотреть кадр за кадром в формате фильма.

Запись траекторий в двоичном формате

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

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

См. Также

Ссылки