Задание 2. Работа с ProDy¶

Алимова Альфия

In [1]:
import prody
import seaborn as sns
import matplotlib.pyplot as plt
import numpy as np
import pandas as pd
In [2]:
P1 = prody.parsePDB('6l7k.pdb') # ЯМР-структура
P2 = prody.parsePDB('3akm.pdb') # РСА-структура

@> 2109 atoms and 20 coordinate set(s) were parsed in 0.17s.
@> 4463 atoms and 1 coordinate set(s) were parsed in 0.08s.
In [3]:
# средние RMSF для остатков в модели ЯМР
mean_RMSFs = [np.mean(prody.calcRMSF(res1)) for res1 in P1['A'].iterResidues()]
In [4]:
# средние B-факторы для каждого остатка в РСА-модели
mean_betas = []

for res2 in P2['A'].iterResidues():
    mean_beta = np.mean(res2.getBetas())
    mean_betas.append(mean_beta)

mean_betas_slice = mean_betas[0:131]
In [5]:
# проверим, что у нас получились массивы одинакового размера

print(len(mean_RMSFs), len(mean_betas_slice))

131 131
In [6]:
# точечная диаграмма по RMSF и B-фактору

plt.figure(1)
sns.set_style("darkgrid")
sns.relplot(x=mean_RMSFs, y=mean_betas_slice, color='#DE3163')
plt.xlabel("Mean RMSF")
plt.ylabel("Mean B-factor")
plt.show()

<Figure size 640x480 with 0 Axes>
No description has been provided for this image

Рис.5. Точечная диаграмма зависимости средних B-факторов от средних значений RMSF

Видно по графику, что при увеличении показателя RMSF увеличивается B-фактор, однако не во всех случаях. Большая часть данных лежит в диапазоне от 0 до ~ 1.5 (по оси Х), у которой средний В-фактор ниже 30 находится, и немного разбросана относительно друг друга. Скорее всего, это связано с неполнотой данных, так как мало моделей (всего лишь 20 (поданные в PDB) из 100 (всего сделано)) было рассмотрено. Поэтому положительная корреляция наблюдается, но немного слабая.

Посчитаем коррелляцию Пирсона:

In [7]:
import scipy.stats
In [8]:
scipy.stats.pearsonr(mean_RMSFs, mean_betas_slice)
Out[8]:
PearsonRResult(statistic=0.3172832886950378, pvalue=0.00022201881803534006)

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

Задание 3. Точная оценка RMSF¶

In [9]:
import math
In [10]:
selection_P2 = P2.select('chain A and resnum 1 to 130') # выбираем из PCA-модели только цепь А и выравненные остатки на ЯМР-модель
b_factors = selection_P2.getBetas()
In [11]:
RMSFs = prody.calcRMSF(P1.select('not hydrogen and not OXT')) # рассчитываем RMSF из ЯМР-модели
In [12]:
# расчёт RMSF^2 = 3*B-factor/(8 * pi**2)
# построение точечной диаграммы: RMSF^2(ожидаемое, через ProDy) = RMSF^2 (рассчитано по формуле выше)

RMSF_obs = [value ** 2 for value in RMSFs]
RMSF_calc = [3*B/(8 * math.pi**2) for B in b_factors]
In [13]:
plt.figure(1)
sns.set_style("darkgrid")
sns.relplot(x=RMSF_calc, y=RMSF_obs, color='#DE3163')
plt.xlabel("RMSF^2 (ProDy)")
plt.ylabel("RMSF^2 (calc)")
plt.show()

<Figure size 640x480 with 0 Axes>
No description has been provided for this image

Рис.6. Сравнение значений квадрата RMSF, полученных в результате расчёта при помощи пакета ProDy (RMSF^2 (ProDy)) и формулы (RMSF^2 (calc))

In [14]:
scipy.stats.pearsonr(RMSF_calc, RMSF_obs)
Out[14]:
PearsonRResult(statistic=0.1580265702981056, pvalue=2.480225257497689e-07)

Коэффициент корреляции в данном случае уменьшился, как и p-value (он равен 2.48 * 10-7, что меньше 0.05; то есть корреляция существует, но очень слабая), однако можно заметить, что на данном графике точки сконцентрированы больше всего ниже значения RMSF, рассчитанного по формуле и равного 5. Соответствие формуле плохое, потому что есть разброс в данных.

Итог для всех наших вычислений таков: хотя корреляция в обоих случаях существует, это не позволяет точно сказать, что полученный ансамбль из моделей, полученных в ходе ЯМР-эксперимента, представляет собой эволюцию структуры во времени. Тогда и нельзя сказать, что соответстуют B-факторы значениям RMSF.