hw 9

Моделирование самосборки липидного бислоя из случайной стартовой конформации

In [53]:
#import warnings
#warnings.filterwarnings("default")
import io

import ipywidgets
from IPython.display import display, display_svg, SVG, Image, HTML

import numpy as np

%matplotlib inline
import matplotlib.pyplot as plt
import matplotlib.patches as mpatches

from base64 import b64encode
In [2]:
import __main__
__main__.pymol_argv = [ 'pymol', '-qc' ]
 
import pymol
 
# Call the function below before using any PyMOL modules.
pymol.finish_launching()

Скачиваем файлы в рабочую директорию:

In [3]:
%%bash
wget http://kodomo.cmm.msu.ru/~golovin/bilayer/dppc.itp
wget http://kodomo.cmm.msu.ru/~golovin/bilayer/lipid.itp
wget http://kodomo.cmm.msu.ru/~golovin/bilayer/dppc.gro
wget http://kodomo.cmm.msu.ru/~golovin/bilayer/b.top
wget http://kodomo.cmm.msu.ru/~golovin/bilayer/em.mdp
wget http://kodomo.cmm.msu.ru/~golovin/bilayer/pr.mdp
wget http://kodomo.cmm.msu.ru/~golovin/bilayer/md.mdp

На основе одного липида созадим ячейку с 64 липидами:

In [4]:
%%bash
genconf -f dppc.gro -o b_64.gro -nbox 4 4 4

С помощью editconf преобразуйте dppc.gro и b_64.gro в pdb файлы:

In [5]:
%%bash
editconf -f dppc.gro -o dppc.pdb
In [6]:
%%bash
editconf -f b_64.gro -o b_64.pdb

Просмотрите результат в PyMol:

In [8]:
pymol.cmd.load("./dppc.pdb")
pymol.cmd.png("./dppc.png")
Image("./dppc.png")
Out[8]:
In [10]:
pymol.cmd.load("./b_64.pdb")
pymol.cmd.png("./b_64.png")
Image("./b_64.png")
Out[10]:

В текстовом редакторе в файле b.top установите правильное количество липидов в системе:

In [11]:
# nano b.top 
# DDPS    34 -> DDPS    64

Сделаем небольшой отступ в ячейке от липидов, что бы добавить примерно 2500 молекул воды:

In [12]:
%%bash
editconf -f b_64.gro -o b_ec -d 0.5

Проведём оптимизацию геометрии системы, что бы удалить "плохие" контакты молекул:

In [13]:
%%bash
grompp -f em -c b_ec -p b -o b_em -maxwarn 2
In [14]:
%%bash
mdrun -deffnm b_em -v

Начальное значение максимальной силы - 4.37970e+05

Конечное значение максимальной силы - 6.16887e+02

In [19]:
Image("./maxforce.png")
Out[19]:

Добавим в ячейку молекулы воды типа spc:

In [20]:
%%bash
genbox -cp b_em -p b -cs spc216 -o b_s

Проведём "утряску" воды:

In [21]:
%%bash
grompp -f pr -c b_s -p b -o b_pr -maxwarn 1
In [22]:
%%bash
mdrun -deffnm b_pr -v
In [23]:
%%bash
grompp -f em -c b_s -p b -o b_empr -maxwarn 1
mdrun -deffnm b_empr -v
In [24]:
%%bash
grompp -f pr -c b_empr -p b -o b_pr -maxwarn 1
mdrun -deffnm b_pr -v
In [25]:
%%bash
editconf -f b_pr.gro -o b_pr.pdb
editconf -f b_s.gro -o b_s.pdb
In [27]:
pymol.cmd.load("./b_pr.pdb")
pymol.cmd.png("./b_pr.png")
Image("./b_pr.png")
Out[27]:
In [29]:
pymol.cmd.load("./b_s.pdb")
pymol.cmd.png("./b_s.png")
Image("./b_s.png")
Out[29]:

Структура на втором изображении кажется менее разреженной (более плотной).

Запуск на lomonosov:

cd ..

scp -r ./hw9 lom:_scratch/chem/volosnikova

ssh lom

cd _scratch/chem/volosnikova

cp /home/users/golovin/progs/share/gromacs/top/residuetypes.dat .

grompp -f md -c b_pr -p b -o b_md -maxwarn 1

sbatch -n 4 -e error.log -o output.log -t 5 -p test impi /opt/ccoe/gromacs-5.0.4/build/bin/gmx_mpi mdrun -testverlet -deffnm b_md -v

1654607

sbatch -N1 --ntasks-per-node=2 -e error-gpu.log -o output.log -t 350 -p gpu impi /opt/ccoe/gromacs-5.0.4/build/bin/gmx_mpi mdrun -testverlet -deffnm b_md -v

1654623

Анализ результатов моделирование самосборки липидного бислоя

Любой анализ начинают с визуального анализа движений молекул. При вопросе о выводк групп выберите DPPC.

In [41]:
%%bash
trjconv -f b_md.xtc -s b_md.tpr -o b_pbc_11.pdb -skip 20 -pbc mol
DPPC
In [95]:
video = io.open('./b_pbc_1.mp4', 'r+b').read()
video_encoded = b64encode(video)
html_code = '<video controls alt="PyMol Movie" src="data:video/mp4;base64,{}" type="video/mp4">'.format(video_encoded)
display(HTML(html_code))

#на всякий случай приложила видео, файл под файлом hw9.ipynb

Занесите ваши наблюдения в журнал. Отметье в какой момент происходит образование бислоя. Для опредления соотвествия между номером модели и временем моделирования найдите в pdb файле выражение "MODEL 50" (цифра 50 это и есть номер модели) и двумя строчками выше будет упосинание о времени в моедлировании.

In [49]:
Image("./b_pbc_110067.png")
Out[49]:

Кажется, что бислой появляется +- с 67 кадра, что соответствует следующим параметрам:

TITLE bilayer in water t= 33000.00000

REMARK THIS IS A SIMULATION BOX

CRYST1 103.540 73.487 21.079 90.00 90.00 90.00 P 1 1

MODEL 67

Определите площадь занимаемую одним липидом. Для этого вам придется получить размеры ячейки из траектории.

In [61]:
%%bash
g_traj -f b_md.xtc -s b_md.tpr -ob box_1.xvg
DPPC

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

In [75]:
# time x y z
box = np.loadtxt("./box_1.xvg", skiprows=24)
plt.figure(1)
plt.plot(box[:,0], box[:,2] * box[:,3] / 64., 'r')
plt.xlabel('time')
plt.ylabel('Syz')
plt.show()

Нормалью к поверхности бислоя является ось X. Площадь уменьшается, что говорит о том, что система стремится принять оптимальное состояние.

Определите изменение гидрофобной и гидрофильной поверхности в ходе самосборки.

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

In [63]:
%%bash
g_sas -f b_md.xtc -s b_md.tpr -o sas_b.xvg -dt 100
DPPC
In [72]:
sas = np.loadtxt('sas_b.xvg', skiprows=22)
plt.figure(1)
plt.plot(sas[:,0], sas[:,1], 'r', label = 'hydrophob', )
plt.plot(sas[:,0], sas[:,2], 'b', label = 'end')
plt.xlabel('time')
plt.legend()
plt.show()

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

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

Для конца траектории и для начала траектории:

In [87]:
%%bash 
wget http://kodomo.cmm.msu.ru/~golovin/bilayer/sn1.ndx
g_order -s b_md -f b_md.xtc -o ord_end.xvg -n sn1.ndx -b 45000 -d X
DPPC
g_order -s b_md -f b_md.xtc -o ord_start.xvg -n sn1.ndx -e 5000 -d X
DPPC

Постройте зависимости и добавтье к отчёту. Сделайте вывод об изменении параметра порядка при образовании бислоя.

In [94]:
start = np.loadtxt('#ord_start.xvg.2#', skiprows=12)
end = np.loadtxt('#ord_end.xvg.2#', skiprows=12)

plt.figure(1)
plt.plot(start[:,0], start[:,1], 'r', label = 'start')
plt.plot(end[:,0], end[:,1], 'b', label = 'end')
plt.xlabel('atom')
plt.ylabel('order')
plt.legend()
plt.show()

В конце сборки бислоя порядок выше, чем в начале.