Моделирование перехода А-формы ДНК в В-форму в воде

После работы на сервере были получены следующие файлы:

dna_pbc_1.pdb
rms_1.xvg
rms_2.xvg
sas_dna.xvg
hbond_dna.xvg
hbond_dna_sol.xvg

Затем они были скачаны на домашний компьютер и дальнейший анализ проводился на нем.

Прежде всего файл с фрагментом ДНК был модифицирован для возможности дальнейшей работы. Для последующего моделирования (задача №1647686) была проведена оптимизация геометрии системы для удаления "плохих" контактов в молекуле. Начальное и конечное значение максимальной силы равно соответственно: 3045 для атома 50 и 593 для атома 32.

Также необходимо было добавить в систему воду и произвести утряску. После утряски вода располагается более равномерно и хаотично (справа), мы не наблюдаем никаких структур.

In [51]:
display(Image('/Users/anna/Desktop/with_water.png'))

Визуальный анализ движения молекулы ДНК показывает, что примерно на 5 кадре происходит переход в В-форму (справа): можно заметить, как пары азотистых оснований образуют параллельные "стопки", характерные для В-формы, меняется изгиб цепи.
Примерное время - 800 пс
(Улучшение визуализации не сильно помогло, все равно удобнее просматривать по кадрам)

In [56]:
display(Image('/Users/anna/Desktop/transformation.png'))
In [53]:
def video(fname, mimetype):
    from IPython.display import HTML
    video_encoded = open(fname, "rb").read().encode("base64")
    video_tag = '<video controls alt="test" src="data:video/{0};base64,{1}">'.format(mimetype, video_encoded)
    return HTML(data=video_tag)
video("/Users/anna/Downloads/pymol.mov", "mp4")
Out[53]:
In [2]:
%matplotlib inline
In [58]:
import numpy as np
import matplotlib.pyplot as plt
from scipy import optimize 
import IPython.display
import ipywidgets
from IPython.display import display,display_svg,SVG,Image, FileLink, FileLinks

Изменение среднеквадратичного отклонения от начального положения в зависимости от времени.

In [28]:
x, y = [], []
with open('/Users/anna/Downloads/struct/rms_1.txt') as f:
    for line in f:
        cols = line.split()
        x.append(float(cols[0]))
        y.append(float(cols[1]))
    plt.plot(x, y, 'b.')
    plt.xlabel('time')
    plt.ylabel('RMSD')
    plt.show()        

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

In [30]:
x, y = [], []
with open('/Users/anna/Downloads/struct/rms_2.txt') as f:
    for line in f:
        cols = line.split()
        x.append(float(cols[0]))
        y.append(float(cols[1]))
    plt.plot(x, y, 'b.')
    plt.xlabel('time')
    plt.ylabel('RMSD')
    plt.show()   

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

In [57]:
x, y1, y2 = [], [], []
with open('/Users/anna/Downloads/struct/sas_dna.txt') as f:
    for line in f:
        cols = line.split()
        x.append(float(cols[0]))
        y1.append(float(cols[1]))
        y2.append(float(cols[2]))
        

    plt.plot(x, y1,label='hydrophobic')
    plt.plot(x,y2,label='hydrophilic')
    plt.xlabel('time')
    plt.legend()
    fitfunc = lambda p, x: p[0]*pow(p[1]-x,1) + p[2] # Target function
    errfunc = lambda p, x, y: fitfunc(p, x) - y # Error function
    p0 = [1,1,0] # Initial guess for the parameters
    p1, success = optimize.leastsq(errfunc, p0[:], args=(x[100:], y1[100:]))
    p2, success = optimize.leastsq(errfunc, p0[:], args=(x[100:], y2[100:]))
    plt.plot( x,fitfunc(p1,x),"r-",c='blue',alpha=0.5)
    plt.plot( x,fitfunc(p2,x),"r-",c='red',alpha=0.5)

    plt.show()   

В целом гидрофобная и гидрофильная поверхности, доступные растворителю, изменяются в некоторых пределах на всем времени симуляции, но в среднем остаются примерно постоянными. Возможно, это связано с малым размером участка ДНК, т.е. на таком масштабе изменение не видно. Рассмотрим подробнее начало: здесь наблюдается локальное увеличение гидрофильной поверхности и уменьшение гидрофобной. Вспомним, что конформационный переход приходится примерно на 800. Это хорошо соотносится с изменением гидрофильной/гидрофобной поверхности. Т.о. возможно, переход из А-формы в B происходит из-за того, что уменьшение свободной энергии соответствует уменьшению гидрофобной поверхности, доступной растворителю.

Водородные связи внутри ДНК. Их количество меняется на протяжении моделирования, преобладают значения 13-14-15. Т.е. в среднем 14, что хорошо, поскольку совпадает со стандартным количеством водородных связей для утсон-криковских пар в нашей молекуле.

In [32]:
x, y1 =  [], []
with open('/Users/anna/Downloads/struct/hbond_dna.txt') as f:
    for line in f:
        cols = line.split()
        x.append(float(cols[0]))
        y1.append(float(cols[1]))
    plt.plot(x, y1,'ro')
    plt.xlabel('time')
    plt.ylabel('# hhbond DNA-DNA')
    plt.show()   

Зависимость количества водородных связей между водой и ДНК не была построена, так как файл выдавал нули либо скопление связей в одной временной точке, изменение параметра cut-off не помогло. Видимо, где-то закралась ошибка!

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

In [ ]: