К.В. Шайтан, К.Б. Терёшкина
Молекулярная динамика белков и пептидов
Методическое пособие

2.3. Обработка результатов.
2.3.1. Построение графиков статистик, которые были рассчитаны параллельно с траекторией. Использование пакета Matlab.
Построение графиков. Использование drawstat.m
Рассчитанные параллельно с траекторией статистики, как правило, находятся в файлах с расширениями dat. Проверить, какое имя соответствует какой статистике, можно, открыв tsb-файл или файл *.dat, если при задании статистики указывалось её полное название.
Полученные результаты удобнее всего обрабатывать в программе Matlab. Для этого понадобится открыть Matlab и задать рабочей директорией (Current Directory) ту директорию, где сохранены статистики. Для простого построения графиков нужно скопировать в директорию со статистиками файл drawstat.m. Чтобы увидеть список всех файлов со статистиками, в командной строке Matlab нужно набрать:
dir_*.dat
Примечание: символом _ обозначен пробел. Если в команде не встречается этот символ, значит текст нужно вводить без пробела.
При выполнении этой команды на экран будут выведены все файлы с расширением dat. Для того чтобы вывести названия только определённых типов файлов, наряду со звёздочкой можно записать набор символов из названия файла, характерных для данного типа статистик. Например, для монопептидов в воде названия обычно задаются в следующем виде:
Название файла Статистика Примечания
ala2adv2000.dat Расширенная ala – код аминокислоты

2 – to

adv, max, db, 2d, 3d, af, cf – названия статистик

f, p, h – углы

2000 – температура
ala2max2000.dat Распределение Максвелла по скоростям
ala2db2ff2000.dat Одномерное распределение по углу φ
ala2db2pp2000.dat Одномерное распределение по углу ψ
ala2dbhh22000.dat Одномерное распределение по углу χ
ala22d2fp2000.dat Двумерное распределение по углам φ и ψ
ala22d2fh2000.dat Двумерное распределение по углам φ и χ
ala22d2ph2000.dat Двумерное распределение по углам ψ и χ
ala23d2000.dat Трёхмерное распределение по углам φ, ψ и χ
ala2af2ff2000.dat Автокорреляционная функция угла φ
ala2af2pp2000.dat Автокорреляционная функция угла ψ
ala2af2hh2000.dat Автокорреляционная функция угла χ
ala2cf2fp2000.dat Кросскорреляционная функция углов φ и χ
ala2cf2fh2000.dat Кросскорреляционная функция углов φ и ψ
ala2cf2ph2000.dat Кросскорреляционная функция углов χ и ψ
Для того чтобы вывести названия только двумерных распределений по углам φ и ψ, нужно записать:
dir_*2D2fp*.dat
Запускает построение графиков следующая команда:
drawstat('moya2D2fpstatistika.dat')
Здесь moya2D2fpstatistika.dat – имя файла со статистикой, которую нужно визуализировать. Вводить это имя не нужно, достаточно выделить его левой кнопкой мыши из рабочего окна Matlab и перенести его (также левой кнопкой) в строку с командой.
При построении графиков двумерных распределений плотностей вероятностей иногда может потребоваться изменить цветовую палитру на чёрно-белую. В таком случае необходимо переписать файл minusgray.m в рабочую директорию, а в файле drawstat.m дописать строку с палитрой для данного типа графика или исправить её на: colormap (minusgray(256)). В файле drawstat.m типы статистик разделяются на трёхмерные, двумерные и одномерные. При необходимости построить трёхмерное распределение плотности вероятности в подпространстве трёх торсионных углов, перед командой построения нужно ввести значение параметра fdiv, который определяет уровень обрезания энергии. Для этого ввести в командной строке команду:
Global fdiv
затем задать значение fdiv из командной строки. Точки для положений молекулы с энергией выше заданной (вероятностью ниже заданной) не выводятся. Чем меньше fdiv, тем меньше точек будет выведено на график.
Примечание: существует несколько версий файла drawstat.m. Необходимо следить, чтобы при простом построении графиков в первой строке содержалась надпись function DrawStat(name), то есть для выполнения программы нужно было ввести только имя файла как указано в примере выше. В более поздней версии необходимо указывать ещё имена файлов с графикой, рисунком, а также название графика. Эта версия drawstat.m предназначена для обработки нескольких графиков одновременно и работает совместно с loadstat.m.
Если требуется построить несколько графиков для ряда однотипных функций, подписать оси и автоматически сохранить их, то понадобятся файлы loadstat.m и nomes.dat. Так как в файле loadstat.m сразу содержится модуль для проведения дисперсионного анализа, описание этого файла будет приведено в следующем разделе. При этом нужны будут файлы readstat.m и соответствующий файл dendrogramXX.m, где XX – число объектов, для которых строится одномерный дисперсионный анализ.

Примеры графиков для различных функций
Рассмотрим характерный вид графиков для описанных выше статистик. Как уже указывалось во "Введении", необходимо следить за тем, чтобы траектория обладала эргодическими свойствами. Только в этом случае можно говорить о некоей достоверности результатов. Помимо теоретического расчёта, о равномерном посещении фигуративной точной конформационного пространства можно судить и по виду графика. Например, на Рис. 18 приведены графики одномерного распределения плотностей вероятностей по торсионному углу для монопептида аланина. Слева длина траектории составляла 10 пс, справа – 1 нс. Видно, что график слева не гладкий и имеет всего один максимум. Это означает, что бóльшая часть конформационного пространства не посещалась молекулой, а попадания в остальные области были случайными. Таким образом, результаты не могут говорить о какой либо закономерности нахождения молекулы в подпространстве данного торсионного угла. При посещении молекулой значительной части подпространства, график выглядит как на Рис. 18 справа.
Рис. 18. Графики одномерного распределения плотностей вероятностей по торсионному углу для монопептида аланина. Слева – результаты, полученные после расчёта 10пс траектории, справа – после 1нс. Т=300К.
По мере расчёта траектории, в программе Modyp можно следить за изменением вида графиков. Чтобы вывести графики на экран, нужно нажать Calculations –> Graphics. При этом надо обратить внимание на то, что в начале расчёта траектории все графики будут негладкими из-за недостатка статистических данных. По мере расчёта их вид будет сглаживаться. Здесь же следует помнить, что при расчёте корреляционных функций первый визуальный результат можно получить лишь по истечении времени траектории равному времени наблюдения корреляции. То есть если оно задано равным 70пс, то до достижения этого времени на графике будут отображаться только координатные оси.
Ниже приведены типичные графики функций, полученные после расчёта различных статистик. Подписи осей и название графика извлекаются из соответствующих файлов dat. Название графика соответствует задаваемому названию статистики в файле tsb (текст в кавычках). Подписи к осям различаются в зависимости от типа статистики.
1. tAdvanced
В файле со статистикой tAdvanced содержатся данные по параметрам последней точки расчёта траектории. Данные этой статистики не используются при построении графиков. Не следует пытаться строить их с помощью модуля drawstat.m!
2. tMaxwell
Типичный вид распределения Максвелла по скоростям приведён на Рис. 19. Отклонения от данного распределения говорят о серьёзных нарушениях, приводящих к значительным ошибкам.
Рис. 19. Распределение Максвелла по скоростям.
3. tDistDb
Рис. 20. Распределение по расстоянию между двумя атомами.
4. tProbDb
5. tProb2D
Рис. 21. Двумерное распределение плотности вероятности по торсионным углам.
График получен с использованием палитры minusgray. Наиболее часто посещаемые области имеют более тёмный цвет. Для аминокислотных остатков чаще всего рассматривают двумерные распределениям по торсионным углам φ, ψ и χ. Среди возможных вариантов двумерных распределений обычно уделяют особое внимание сечению по углам φ, ψ. Основные варианты вторичной структуры на фоне разрешённых и запрещённых областей представлены на карте Рамачандрана (Рис. 22).
Рис. 22. Карта Рамачандрана для аминокислотного остатка [21,22]. Конформации, которые могут быть достигнуты любым аминокислотным остатком, представлены тёмно-серым цветом. Большинство аминокислот может заселять области, обозначенные светло-серым цветом. Белым обозначены запрещённые конформации, которые, тем не менее, могут встречаться в некоторых белковых структурах.
Здесь:
1 – вторая спираль полипролина (коллагеновая спираль),
2 – антипараллельная β-конформация,
3 – параллельная β-конформация,
4 – левая π-спираль,
5 – правая 27-спираль,
6 – левая α-спираль,
7 – левая 310-спираль,
8 – правая 310-спираль,
9 – правая α-спираль,
10 – левая 27-спираль,
11 – правая π-спираль.
6. tProb3D
Рис. 23. Трёхмерное распределение плотности вероятности по торсионным углам φ, ψ, и χ.
На Рис. 23 приведены трёхменые сечения Пункаре в подпространстве торсионных углов φ, ψ, и χ для аланина. Границы поверхностей соответствуют уровню свободной энергии 2,74 ккал/моль. Наиболее заселёные области конформационного пространства выделены тёмным цветом. Параметр fdiv был взят равным 100. График построен с использованием палитры minusgray.
7. tAutoCf
Рис. 24. Автокорреляционная функция для торсионного угла.
При анализе автокорреляционных функций различают в общем случае три параметра – скорость выхода на асимптоту, характерное время затухания (τ) и величину остаточной корреляции. Скорость выхода на асимптоту на интервале времени [0, τ] говорит о скоррелированности движения по данному торсионному углу, чем быстрее спадает функция на данном интервале, тем менее скоррелировано вращение по торсионному углу. Ограниченность движения в потенциальной яме приводит к появлению остаточной корреляции. Характерное время затухания автокорреляционной функции позволяет судить о времени конформационного перехода по торсионному углу.
8. tAutoCfD
Рис. 25. Нормированная автокорреляционная функция для торсионного угла.
9. tCrossCf
Рис. 26. Кросскорреляционная функция для торсионного угла.
Значения действительной части кросскорреляционной функции могут быть как положительными, так и отрицательными. Это говорит о положительной или отрицательной корреляции двух торсионных углов. Положительная корреляция означает, что оба угла движутся в одном направлении, отрицательная – в противоположных. Корреляция может меняться скачкообразно. Если значение действительной части кросскорреляционной функции близко к нулю, движения двух торсионных углов не скоррелированы.
10. tDistDevCf
Рис. 27. Корреляционная функция отклонения от среднего.
11. tDist2AixCf
Рис. 28. Кросскорреляционная функция отклонения атомов от выбранной оси.

Построение группы графиков. Проведение одномерного дисперсионного анализа. Форматы файлов nomes.dat, loadstat.m и dendrogram.m
Если требуется построить несколько однотипных графиков, удобнее использовать приложение loadstat.m. В нём содержится также информация для проведения одномерного дисперсионного анализа и построения кластерного дерева. Имена файлов со статистиками записываются в nomes.dat. При этом надо следить, чтобы названия файлов содержали одинаковое количество символов и имели расширение dat. Кластерное дерево строится с помощью файла dendrogramXX.m, где XX – число анализируемых статистик.
Вид файла nomes.dat для анализа двумерных распределений по углам φ и ψ в ряду модифицированных тирозинов:

ty222d2fp2000.dat
ty322d2fp2000.dat
tyc22d2fp2000.dat
tyo22d2fp2000.dat
tyr22d2fp2000.dat
tys22d2fp2000.dat

В зависимости от количества и типа статистик в файл loadstat.m необходимо внести некоторые изменения.
Параметр, который необходимо изменить Название параметра и пример его задания Комментарий
Число файлов, подлежащих анализу numeroarq=6 Будут обрабатываться шесть первых файлов из nomes.dat
Названия графиков titlefig(i,:)=[figstat(i,1:3) ' in water TIP3P, 2D (\phi and \psi), 2000K, 10ns']; В названии графика будут первые три символа из названия статистики и далее – текст в апострофах
Название файла fig с кластерным деревом clustfig  
Название файла emf с кластерным деревом clustris  
Столбец с данными, которые будут анализироваться в процессе дисперсионного анализа b(:,j)=a(:,3);
pl(:,j) = sum (a(:,3));
В данном случае анализируется третий столбец. Для трёхмерных распределений это число должно быть равно 4, для двумерных – 3, для остальных – 2
Названия файла для построения кластерного дерева H = dendrogram6 (Z); Необходимо изменить номер (здесь 6) в зависимости от числа статистик. В директории также должен находиться соответствующий файл dendrogramXX.m
Название графика с кластерным деревом title ('Amino acids in TIP3P, 2D, angles \phi and \psi, 2000K, 10ns', 'FontSize', 16)  
В файле dendrogramXX.m нужно изменить часть, отвечающую за создание подписей к данным: 1) задать изменение переменной i в нужных пределах; 2) для каждой подписи создать/изменить две строки – case и label. Порядок названий должен строго соответствовать порядку в файле nomes.dat

for i=1:6
switch v(i)
case 1
    label(i,:) = 'ty2';
case 2
    label(i,:) = 'ty3';
case 3
    label(i,:) = 'tyc';
case 4
    label(i,:) = 'tyo';
case 5
    label(i,:) = 'tyr';
case 6
    label(i,:) = 'tys';
end
end

На Рис. 29 приведён пример кластерного дерева для двадцати природных монопептидов. Анализировались двумерные распределения плотностей вероятностей по углам φ и ψ.
Рис. 29. Кластерное дерево для двумерных распределений плотностей вероятности у двадцати природных монопептидов по углам φ и ψ.

Выделите орфографическую ошибку мышью и нажмите Ctrl+Enter

МГУ им. Ломоносова Rambler's Top100

 

Rambler's Top100