Документ взят из кэша поисковой машины. Адрес оригинального документа : http://kodomo.cmm.msu.ru/~golikov_v/67.html
Дата изменения: Sun May 29 14:42:48 2011
Дата индексирования: Mon Oct 1 23:17:59 2012
Кодировка: Windows-1251
Молекулярная динамика биологических молекул в GROMACS

 

Цель данного занятия - используя пакет молекулярной динамики Gromacs смоделировать самосборку липидного бислоя из случайной стартовой конформации.

Даны файлы:
дополнительной топологии для липида DPPC, dppc.itp;
параметры для липидов lipid.itp;
координаты одного липида dppc.gro;
файл-заготовка тополгии системы b.top;
файл параметров для минимизации энергии em.mdp;
файл параметров для "утряски" воды pr.mdp;
файл параметров для молекулярной динамики md.mdp.

На основе одного липида создадим ячейку с 64 липидами: 
genconf -f dppc.gro -o b_64.gro -nbox 4 4 4
С помощью editconf преобразуем dppc.gro и b_64.gro в pdb файлы: 
editconf -f dppc.gro -o dppc.pdb
editconf -f b_64.gro -o b_64.pdb
В одном из хвостов жирная кислота образует цикл - позднее этот "плохой" контакт будет исправлен оптимизацией геометрии. Установим в файле b.top правильное количество липидов в системе - 64.
Сделаем небольшой отступ в ячейке от липидов, чтобы добавить примерно 2500 молекул воды: 
editconf -f b_64.gro -o b_ec -d 0.5
Проведем оптимизацию геометрии системы, чтобы удалить "плохие" контакты молекул. 
grompp -f em -c b_ec -p b -o b_em -maxwarn 2
mdrun -deffnm b_em -v
Начальное значение максимальной силы = 4.37970e+05, конечное = 6.4541919e+02.
Добавим в ячейку молекулы воды типа spc: 
genbox -cp b_em -p b -cs spc216 -o b_s
Проведем "утряску" воды: 
grompp -f pr -c b_s -p b -o b_pr -maxwarn 1
mdrun -deffnm b_pr -v
Переформатируем b_pr.gro и b_s.gro в pdb-формат: b_pr.pdb и b_s.pdb. Изменение в системах - исчез лишний цикл, образованный жирной кислотой; "утрясенная" вода распределена по объему более хаотично.

Скопируем эти файлы на суперкомпьютер. 
ssh skif
mkdir Golikov
exit
cd md
scp -r md/* skif:Golikov/
Запустим тестовое моделирование на суперкомпьютере: 
ssh skif
cd Golikov
grompp -f md -c b_pr -p b -o b_md -maxwarn 1
mpirun  -np 16 -q test -maxtime 5  /home/golovin/progs/bin/mdrun_mpi -deffnm b_md -v
Файл не содержит ошибок.

Запускаем основное моделирование на суперкомпьютере: 
mpirun  -np 16 -maxtime 1200  /home/golovin/progs/bin/mdrun_mpi -deffnm b_md -v

Анализ результатов

  • Силовое поле используемое при построении топологии топологии - ffgmx
  • Размер и форму ячейки - параллелепипед (6.26000х4.44300х5.77800 нм)
  • Минимизация энергии:
    • Алогритм минимизации энергии integrator = l-bfgs
    • Алгоритм расчета электростатики coulombtype = Cut-off
    • и Ван-дер-Ваальсовых взаимодействий vdw-type = Cut-off
  • Модель, которой описывался растворитель - flexspc
  • Утряска растворителя:
    • Число шагов - 10000
    • Длина шага - 0.001 пс
    • Алгоритм расчета электростатики coulombtype = pme
    • и Ван-дер-Ваальсовых взаимодействий. vdw-type = Cut-off
    • Алгоритмы термостата Tcoupl = Berendsen
    • и баростата. Pcoupl = no
  • Основной расчет МД:
    • Время моделирования 7 hours 15 minutes 03 seconds, количество процессоров - 16.
    • Длина траектории - 50нс
    • Число шагов - 10000000
    • Длина шага - 0.05 пс
    • Алгоритм интегратора md
    • Алгоритм расчета электростатики coulombtype = pme
    • и Ван-дер-Ваальсовых взаимодействий. vdw-type = Cut-off
    • Алгоритмы термостата Tcoupl = v-rescale
    • и баростата Pcoupl = berendsen
1. Визуальный анализ движения молекул. 
trjconv -f b_md.xtc -s b_md.tpr -o b_pbc_1.pdb -skip 20
b_pbc_1.pdb. Проблема - длинные связи, пересекающие ячейку от одного края до другого. Исправляем: 
trjconv -f b_md.xtc -s b_md.tpr -o b_pbc_2.pdb -skip 20 -pbc mol
b_pbc_2.pdb. Теперь все нормально, на анимации видно, как из бислоя периодически вылезают два хвоста.



Бислой начинает образовываться через 10500 фс после начала моделирования.

2. Определим площадь, занимаемую одним липидом. Получим размеры ячейки из траектории. 
   g_traj -f b_md.xtc -s b_md.tpr -ob box_1.xvg
В файле box_1.xvg содержатся размеры ячейки. В первой колонке время, во 2-4 колонках - размер ячейки в виде длин по трем осям. Ось Х является нормалью к поверхности бислоя. Зависимость площади по осям Y и Z от времени, нормированнвя на один липид в слое (т.е. делим на 32), представлена на графике:


Площадь, занимаемая одним липидом до формирования бислоя - 0.8 нм, после формирования - чуть меньше 0.7 нм.

3. Определим изменение гидрофобной и гидрофильной поверхностей в ходе самосборки: 
g_sas -f b_md.xtc -s b_md.tpr -o sas_b.xvg
Зависимость изменения гидрофобной (синим цветом) гидрофильной (красной) поверхностей доступных растворителю от времени:


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

4. Традиционной мерой оценки фазового состояния бифильных молекул является мера порядка. Для анализа нам понадобится специальный индекс файл. Запустим анализ для конца траектории 
   g_order -s b_md -f b_md.xtc -o ord_end.xvg -n sn1.ndx -b 45000 -d P


и для начала траектории 
   g_order -s b_md -f b_md.xtc -o ord_start.xvg -n sn1.ndx -e 5000 -d P


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

Назад