Анализ молекулярной динамики биологических молекул в GROMACS (Анализ результатов моделирование перехода А-формы ДНК в В-форму в воде)

Скачиваем с "Ломоносов" нужные файлы

In [3]: 
# ssh skif
# Enter passphrase for key '/home/students/y12/alibek/.ssh/fbb-prac':
#  Last login: Thu Apr  7 15:47:45 2016 from 93.180.63.127
#┌─[fbbstudent @ access2 ~ ]
#└─> cd _scratch/fbb/Alibek/
#┌─[fbbstudent @ access2 ~/_scratch/fbb/Alibek ]
#└─>
In [4]: 
# На выходе получили папку с файлами и скопировали её в свою рабочую директорию
In [7]: 
trjconv -f dna_md.xtc -s dna_md.tpr -o dna_pbc_1.pdb -skip 20 -pbc mol
# Получился PDB-файл в котором есть все состояние перехода из А-формы в Б-форму. Нужно найти ссотояние, при котором происходит переход.
Найденное состояние - MODEL 13 (t= 2400.00000)
In [8]: 
%matplotlib inline
import numpy as np
import matplotlib.pyplot as plt
from scipy import optimize
In [37]: 
# Находим среднее квадратичное отклонение
# Были получены 2 файла командами:
## g_rms -f dna_md.xtc -s dna_md.tpr -o rms_1
## g_rms -f dna_md.xtc -s dna_md.tpr -o rms_2 -prev 400
%matplotlib inline
a = np.loadtxt("rms_1.xvg")
x_o=a[:,0]
y_o=a[:,1]
#Plot it
plt.scatter(x_o, y_o)
plt.xlim(1.3,10000)
plt.show()
print x_o
print y_o

[     0.     10.     20. ...,   9980.   9990.  10000.]
[  1.43900000e-04   1.55344900e-01   1.61314300e-01 ...,   3.21253800e-01
   3.54390300e-01   3.24048100e-01]
Видно, что имеет место разрыв в промежутке от 4500 до 5800, можно предположить, что именно в тот момент и происходит переход одной формы ДНК в другую
In [49]: 
%matplotlib inline
a = np.loadtxt("rms_2.xvg")
x_1=a[:,0]
y_1=a[:,1]
#Plot it
plt.plot(x_o, y_o,"ro")
plt.plot(x_o, y_1,"bo")
plt.xlim(1.3,10000)
plt.ylim(0,0.5)
plt.show()
print x_o
print y_o

[     0.     10.     20. ...,   9980.   9990.  10000.]
[  1.43900000e-04   1.55344900e-01   1.61314300e-01 ...,   3.21253800e-01
   3.54390300e-01   3.24048100e-01]
In [52]: 
# Определим изменение гидрофобной и гидрофильной поверхности в ходе конформационного перехода командой:
##  g_sas -f dna_md.xtc -s dna_md.tpr -o sas_dna.xvg 
%matplotlib inline
a = np.loadtxt("sas_dna.xvg")
x_o=a[:,0]
y_o=a[:,1]
#Plot it
plt.scatter(x_o, y_o)
plt.xlim(1.3,10000)
plt.ylim(0,15)
plt.show()
print x_o
print y_o

[     0.     10.     20. ...,   9980.   9990.  10000.]
[ 8.98307  8.48768  8.81794 ...,  9.11518  8.68584  8.98307]

После моделирования, площадь гидрофобной поверхности уменьшилась

In [60]: 
# Определим водородные связи между ДНК и водой:
##  g_hbond  -f dna_md.xtc -s dna_md.tpr -num hbond_dna_sol  
%matplotlib inline
a = np.loadtxt("hbond_dna_sol.xvg")
x_o=a[:,0]
y_o=a[:,1]
#Plot it
plt.scatter(x_o, y_o)
plt.xlim(1.3,10000)
plt.ylim(90,130)
plt.show()
print x_o
print y_o

[     0.     10.     20. ...,   9980.   9990.  10000.]
[  56.  115.  108. ...,  105.  114.  112.]

Большое количество водородных связей в системе ДНК-вода

In [58]: 
# Определим водородные связи между ДНК и ДНК:
##  g_hbond  -f dna_md.xtc -s dna_md.tpr -num hbond_dna 
%matplotlib inline
a = np.loadtxt("hbond_dna.xvg")
x_o=a[:,0]
y_o=a[:,1]
#Plot it
plt.scatter(x_o, y_o)
plt.xlim(1.3,10000)
plt.ylim(0,20)
plt.show()
print x_o
print y_o

[     0.     10.     20. ...,   9980.   9990.  10000.]
[ 14.  15.  15. ...,  14.  13.  13.]

Количество водородных связей уменьшилоась, чего следовало ожидать