Вгрузим же
! wget https://kodomo.fbb.msu.ru/~golovin/bilayer/dppc.itp
! wget https://kodomo.fbb.msu.ru/~golovin/bilayer/lipid.itp
! wget https://kodomo.fbb.msu.ru/~golovin/bilayer/dppc.gro
! wget https://kodomo.fbb.msu.ru/~golovin/bilayer/b.top
! wget https://kodomo.fbb.msu.ru/~golovin/bilayer/em.mdp
! wget https://kodomo.fbb.msu.ru/~golovin/bilayer/pr.mdp
! wget https://kodomo.fbb.msu.ru/~golovin/bilayer/md.mdp
Перешли в рабочую директорию
! cd ../zhenia/
! source /home/preps/golovin/progs/bin/GMXRC.bash
На основе одного липида создали ячейку с 64 липидами
! gmx genconf -f dppc.gro -o b_64.gro -nbox 4 4 4
Создали pdb файлы
! editconf -f dppc.gro -o dppc.pdb
! editconf -f b_64.gro -o b_64.pdb
Теперь посмотрим на то, что у нас получилось
import __main__
from IPython.display import Image
__main__.pymol_argv = [ 'pymol', '-c' ] # открыть pymol без использования внешнего GUI
import pymol
pymol.finish_launching()
from pymol import cmd,stored
import numpy as np
import time
import pandas as pd
import matplotlib.pyplot as plt
cmd.load('dppc.pdb')
cmd.bg_color('white')
cmd.do('''set stick_radius, 0.5
show stick, all
orient dppc
zoom dppc, complete = 1
''')
cmd.do('png lipid.png, width=1080, height=1000, ray=1')
Image(filename='lipid.png')
Видим, что у меня истекла лицензия, а файл dppc содержит структуру фосфатидилхолина у которого одна из жирных кислот содержит циклопентановый элемент.
cmd.delete("all")
cmd.bg_color('white')
cmd.load('b_64.pdb')
cmd.do('''as cartoon, all
set stick_radius, 0.5
show stick, all
zoom b_64, complete = 1
turn x, 20
turn y, 20
turn z, 20'''
cmd.do('png lipids.png, width=1080, height=1000, ray=1')
Image(filename='lipids.png')
Получилась ячейка!
В текстовом редакторе в файле b.top установили правильное количество липидов в системе - 64. В файле с параметрами для молекулярной динамики установили шаг по времени 0.004 ps.
Сделали небольшой отступ в ячейке от липидов, что бы добавить примерно 2500 молекул воды.
! gmx editconf -f b_64.gro -o b_ec -d 0.5
Проведём оптимизацию геометрии системы, что бы удалить "плохие" контакты молекул
! gmx grompp -f em -c b_ec -p b -o b_em -maxwarn 2
! gmx mdrun -deffnm b_em -v
Потенциальная энергия в ходе 68 шагов уменьшилась с 4.63573e+05 Кдж/моль до -3.11687e+04 Кдж/моль.
Добавим в ячейку молекулы воды
! gmx solvate -cp b_em -p b -cs spc216 -o b_s
Проведём "утряску" воды:
%%bash
gmx grompp -f pr -c b_s -p b -o b_pr -maxwarn 1
gmx mdrun -deffnm b_pr -v
#если водяной взрыв, то
gmx grompp -f em -c b_s -p b -o b_empr -maxwarn 1
gmx mdrun -deffnm b_empr -v
gmx grompp -f pr -c b_empr -p b -o b_pr -maxwarn 1
gmx mdrun -deffnm b_pr -v
Теперь переведем файл с добавленной водой b_s.gro и файл с оптимизированной водой b_pr.gro в pdb формат и посмотрим, как это выглядит
! gmx editconf -f b_s.gro -o b_s.pdb
! gmx editconf -f b_pr.gro -o b_pr.pdb
Посмотрим на бислой до "утряски" воды
cmd.do("reinitialize")
cmd.delete("all")
cmd.bg_color('white')
cmd.load('b_s.pdb')
cmd.do('''as cartoon, all
set stick_radius, 0.5
show stick, all
show nonbonded
zoom b_s, complete = 1
turn x, 20
turn y, 20
turn z, 20''')
cmd.do('png before.png, width=1080, height=1000, ray=1')
Image(filename='before.png')
... и после "утряски" воды
cmd.do("reinitialize")
cmd.delete("all")
cmd.bg_color('white')
cmd.load('b_pr.pdb')
cmd.do('''as cartoon, all
set stick_radius, 0.5
show stick, all
show nonbonded
zoom b_pr, complete = 1
turn x, 20
turn y, 20
turn z, 20''')
cmd.do('png after.png, width=1080, height=1000, ray=1')
Image(filename='after.png')
В структуре после оптимизации молекулы воды расположены более равномерно и, соответственно, более оптимально для образования водородных связей. В структуре до оптимизации молекулы воды расположены чуть более "клочками".