Практикумы 1 и 2. Знакомство и работа в Pymol

# Запускаем pymol
import __main__
__main__.pymol_argv = [ 'pymol', '-x' ]

import pymol
pymol.finish_launching()
from pymol import cmd, stored

---

Отчетность по первому практикуму - прокомментировать приведенный в упражнениях код.

1.1.

cmd.fetch('1lmp')                                      # Загружаем структуру 1lmp

import time                                            # Становимся властелинами времени
for i in range (1,200):                                # 200 раз делаем следующее непотребство, варьируя i:
     cmd.zoom(str(i-1)+"+"+str(i)+"+"+str(i+1)+'/CA')  # Приблизим CA атомы i-1, i и i+1 остатков. 
                                                       # Для выделения атомов мы используем selection macro
                                                       # https://pymolwiki.org/index.php/Selection_Macros
     time.sleep(.5)                                    # Передышка полсекунды

1.2.

cmd.do('''             # передаем на выполнение паймолу набор pymol команд
fetch 1cll, async=0    # Загружаем структуру 1cll, async=0 говорит паймолу, что надо дождаться выполнения этой команды
                       # прежде чем переходить к следующей. 
as lines, n. C+O+N+CA  # Отобразить в виде линий атомы остова (= атомы с именами C, O, N, CA)
zoom i. 4+5            # Приблизить камеру к остаткам 4 и 5
mset 1 x1000           # Повторить первое (и единственное) состояние модели в течение 1000 кадров мультика
mview store''')        # Сделать текущий вид ключевым. Между ключевыми кадрами паймол может делать плавные переходы.
                       # Тут у нас такой кадр один, так что никакой анимации не получится.

1.3*

stored.r = []  # Заводим список r в специальной переменной pymol stored. Ее видим как мы из ipython ноутбука,
               # так и паймол во время исполнения команд.
cmd.iterate('1cll and n. CA','stored.r.append(resi)')  # Проитерируемся по всем CA атомам объекта 1сll,
                                                       # По мере итерации будем записывать номера остатков в stored.r


import numpy as np  # Импортируем нампай

length = len(stored.r)               # получаем длину списка номеров остатков
colors = np.linspace(1,0.5, length)  # организовываем массив чисел от 1 до 0.5 той же длины, что и stored.r
for k,i in enumerate(stored.r):      # итерируемся по номерам остатков, k - индекс номера остатка i в списке stored.r
    cmd.set_color('col%d' % k, [colors[k],0.5,0.75])  # Создаем цвет с названием col_k в RGB представлении
                                                      # В цикле мы варьируем красный от 1 до 0.5, зеленый и синий фиксирован.
                                                      # Получаются цвета от розового (много красного) до темно-сиреневого
                                                      # (мало красного)
    print([1,1,colors[k]])  # строчка, которая ничего особо не делает. Элемент обфускации? :D
    cmd.set('cartoon_color','col%d' % k ,'resi %d' % int(i))  # Устанавливаем цвет отображения cartoon для i-ого остатка
                                                              # равный col_k
cmd.show_as('cartoon','all')  # Показываем все в виде cartoon

Получилась такая красота. Если дама будет делать практикум в районе 8 марта, ей будет приятно, т.к. напоминает цифру 8.

1.4*

for i in range(100):     # ранее мы сказали паймолу, что хотим мультик из 1000 кадров. Пора этим воспользоваться
    cmd.frame((10*i)+1)  # Переключаемся на 10*i + 1 кадр (от 1 до 991 кадра). Между 991 и 1 кадром остается место, 
                         # чтобы камера плавно приплыла в начало мультика.
    cmd.zoom('n. CA and i. %d+%d' % (i,i+7))  # Зумим на CA атомы i и i+7 остатков. В итоге, в цикле, где i от 0 до 99 пробегает,
                                              # Мы увидим СА атомы остатков до 106, то есть не совсем все. Номера остатков в stored.r
                                              # у нас от 1 до 152 (см код ниже).
    cmd.mview('store')                        # Сохраняем вид в качестве ключевого кадра. Между ключевыми кадрами паймол делает плавный переход
                                              # насколько ему позволяет количество кадров между ключевыми кадрами.

min(map(lambda x : int(x), stored.r))

1

max(map(lambda x : int(x), stored.r))

152

---

Отчетность по второму практикуму - код и картинки того, что код делает. Иногда только картинки.

2.1. Sculpting

cmd.fetch("1LMP")

'1LMP'

Было.

Стало.
Согласно https://pymolwiki.org/index.php/Molecular_Sculpting, sculpting позволяет изменять геометрию молекул, при этом сохраняя длины связей, углы, хиральность и планарность химических групп. Однако настоящей минимизации энергии не происходит, как могло бы показаться. Pymol просто стремится сделать перечисленные выше параметры как можно ближе к исходным параметрам из pdb.

2.2. Mutagenesis

cmd.reinitialize()
cmd.fetch("1LMP")
cmd.bg_color("white")
cmd.select("ligand", "r. NAG or r. NDG")
cmd.select("site", "byres all within 6 of ligand")
cmd.show("licorice", "site")
cmd.color("green", "site")
cmd.select("hph_core", "byres all within 5 of i. 56")
cmd.show("licorice", "hph_core")
cmd.color("cyan", "hph_core")
cmd.color("red", "i. 56")
cmd.color("orange", "ligand")
cmd.set("cartoon_transparency", 0.8)
cmd.set("cartoon_side_chain_helper", 1)

Наверное, одним из вариантов замены было бы создать внутри гидрофобного ядра белка аргинин. Например, замена F56R. Это бы скорее всего нарушило фолд, потому что аргинин заряжен и может образовывать до 5 водородных связей. Внутри гидрофобного ядра ему трудно будет образовать все 5 водородных связей, да и чтобы переместить заряд из среды с высокой диэлектрической константой (поверхность белка) в среду с низкой диэлектрической константой (гидрофобное ядро белка) нужно затратить работу.
Однако, в случае такой замены не понятно, просто глядя на структуру, что станет с фолдом. Во что он свернется, если вообще свернется.

Окружение F56 (красный) показано циановым, окружение лиганда (оранжевый) показано зеленым.

Было.

Стало.

2.3. Мультик

cmd.set_name("1LMP", "mut")
cmd.fetch("1LMP", "wt")
cmd.hide("licorice", "all")
cmd.show("licorice", "i. 56")
cmd.color("cyan", "m. mut")
cmd.color("white", "m. wt")
cmd.color("red", "i. 56 and m. mut")
cmd.color("white", "i. 56 and m. wt")
cmd.orient("all")
cmd.translate([55,0, 0], object="m. wt")
cmd.orient("all")

cmd.mview("purge")
dur = 10
n_frames = 24 * dur
cmd.mset(f"1x{n_frames}")
cmd.mview("store", 1, object="m. mut")
cmd.mview("store", 1, object="m. wt")

cmd.translate([-55, 0, 0], object="m. wt")
cmd.mview("store", 0.15 * n_frames, object = "m. wt")
cmd.mview("store", 0.85 * n_frames, object = "m. wt")

cmd.orient("m. mut")
cmd.mview("store", 0.5 * n_frames, object="m. mut")
cmd.mview("store", 0.5 * n_frames, object="m. wt")
cmd.mview("store", 0.5 * n_frames)

cmd.translate([55, 0, 0], object="m. wt")
cmd.orient("all")
cmd.mview("store", n_frames, object="m. wt")
cmd.mview("store", n_frames, object="m. mut")
cmd.mview("store", n_frames)
cmd.frame(1)

cmd.movie.produce(r"D:\tmp\movie.mpg", mode="draw", preserve=True)

2.4. Пришиваем TAMRA

Заходим на Pubchem, скачиваем .sdf для 3D и картинку метки.

У TAMRA 2 карбоксильные группы, пришьем группой в мета положении в бензольном кольце относительно пиронинового фрагмента.

cmd.delete("all")
cmd.fetch("1LMP")
cmd.load(r"C:/Users/Zerg/Downloads/Conformer3D_CID_2762604.sdf", "tag")
cmd.orient("i. 81")  # Пришиваем к Ser 81
cmd.show("licorice", "i. 81")
cmd.fuse("r. UNK and idx. 5", "i. 81 and n. CB")
cmd.remove("tag")
cmd.remove("i. 81 and n. OG")
cmd.set_dihedral("i. 81 and n. N",
                "i. 81 and n. CA",
                "i. 81 and n. CB",
                "r. UNK and idx. 5", -120)  # Отвернем метку подальше от белка

cmd.edit("r. UNK and idx. 4", "r. UNK and idx. 32")  # Косметические манипуляции со сложноэфирной группой. 
cmd.cycle_valence()                                  # Двойную связь с одного кислорода на другой перетаскиваем
cmd.cycle_valence()
cmd.edit("r. UNK and idx. 5", "r. UNK and idx. 32")
cmd.cycle_valence()
cmd.cycle_valence()

Готово.

2.5. Генерируем альфа-спираль

Значения $\phi, \psi$ взяты из https://www.ncbi.nlm.nih.gov/pmc/articles/PMC3061398/: "...canonical α-helix peak of ϕ, ψ = (−63, −43)...".

def draw_seq(aa, length, phi, psi):
    cmd.delete("all")
    cmd.fragment(aa, "helix")  # Загрузим первый фрагмент
    cmd.alter("helix", "resi=1")  # Почему-то фрагмент загружается как остаток №2, исправим это
    for i in range(1, length):
        cmd.edit(f"i. {i} and n. C")
        cmd.editor.attach_amino_acid("pk1", aa)
        cmd.set_dihedral(f"i. {i} and n. C",
                         f"i. {i+1} and n. N",
                         f"i. {i+1} and n. CA",
                         f"i. {i+1} and n. C", phi)  # phi
        cmd.set_dihedral(f"i. {i} and n. N",
                         f"i. {i} and n. CA",
                         f"i. {i} and n. C",
                         f"i. {i+1} and n. N", psi)  # psi
    return

draw_seq("ala", 100, -63, -43)

Можно проще :D

cmd.fab(input="A"*100, name="helix", ss=1)

Вообще, обладая такой невероятной мощью, можно делать разные нефизические, но прикольные штуки.

draw_seq("lys", 100, 0, 180)
cmd.color("blue","sidechain and not n. CB+2HB")
cmd.color("white","n. CB+2HB")
cmd.color("cyan","bb. and not n. O")
cmd.color("black", "n. O")

$\phi, \psi$ = 0, 180

draw_seq("met", 50, 90, 90)
cmd.remove("n. H and i. 1")
cmd.remove("n. O and i. 50")
cmd.set("solvent_radius", 3.25)
cmd.show("surface", "all")
cmd.color("yellow", "all")
cmd.set("surface_quality", 2)
# Размножим объект и подвигаем руками

$\phi, \psi$ = 90, 90

2.6. Генерируем B-форму ДНК

Как и в предыдущем задании, тут можно считерить:

cmd.fnab("ATGC" * 25, dbl_helix = True)

Но давайте сделаем по-честному, как нас просят.

Возьмем геометрию из структуры Dickerson dodecamer-a: 1BNA.

cmd.delete("all")
cmd.fetch("1BNA")
cmd.select("p1", "i. 6 + i. 19")
cmd.select("p2", "i. 7 + i. 18")
cmd.extract("pair_1", "p1")
cmd.extract("pair_2", "p2")
cmd.pair_fit("m. pair_2 and bb.", "m. pair_1 and bb.")
p2_to_p1_transformation = cmd.get_object_matrix("pair_2")
cmd.delete("m. pair_2")
cmd.set_name("pair_1", "p1")
cmd.delete("1BNA")

length = 100
cmd.alter(f"/p1///DA", f"resi = {length}")
cmd.alter(f"/p1///DT", f"resi = {length + 1}")

for i in range(1, length):  # arrange bp in helical fashion
    cmd.create(f"p{i+1}", f"p{i}")
    cmd.transform_selection(f"p{i+1}", p2_to_p1_transformation)
    cmd.alter(f"/p{i+1}///DA", f"resi = {length - i}")
    cmd.alter(f"/p{i+1}///DT", f"resi = {length + i + 1}")
    
cmd.create("DNA", "all")
cmd.delete("p1")

for i in range(1, length):  # bond bp together
    cmd.delete(f"m. p{i+1}")  # cleanup
  
    cmd.select("pk1", f"i. {i} and n. O3'")
    cmd.select("pk2", f"i. {i+1} and n. P")
    cmd.bond()
  
    cmd.select("pk1", f"i. {length + i} and n. O3'")
    cmd.select("pk2", f"i. {length + i+1} and n. P")
    cmd.bond()

cmd.show("licorice", "all")

Было непросто, но получилось.