Докинг низкомолекулярных лигандов в структуру белка¶
Цель данного занятия - ознакомится с возможностями докинга низкомолекулярного лиганда в структуру белка.
В этом занятии использовали пакет Autodock Vina и oddt. Это программное обеспечение распространяется бесплатно для академических пользователей.
Работали с белком лизоцимом восточного моллюска Meretrix lusoria (идентификатор Uniprot - LYS_MERLU или P86383) структуру которого построили на основе гомологичного моделирования на прошлом практикуме.
Нужно было определить место докинга, удалить лиганд, запустить докинг и провести анализ.
# Подгружаем модули
import numpy as np
import copy
# Отображение структур
import IPython.display
import ipywidgets
from IPython.display import display,display_svg,SVG,Image
# Open Drug Discovery Toolkit
import oddt
import oddt.docking
import oddt.interactions
# Органика
from rdkit.Chem import Draw
from rdkit.Chem.Draw import IPythonConsole
import pmx # Модуль для манипулирования pdb
Подготавливаем белок для докинга¶
#подгружаем структура белка, полученную при помощи гомологического моделирования и сохранённую в данный pdb-файл
pdb=pmx.Model('new_seq.B99990001.pdb')
#посмотрим остатки
for r in pdb.residues[:125]:
print r
# создаём отдельные объекты для белка и для лиганда
newpdb = pdb.copy()
for r in newpdb.residues[122:]:
newpdb.remove_residue(r)
lig = pdb.copy()
for r in lig.residues[:122]:
lig.remove_residue(r)
# Ищем геометрический центр лиганда
count=0
x=0
y=0
z=0
for a in lig.atoms:
x=x+a.x[0]
y=y+a.x[1]
z=z+a.x[2]
print a.x
count=count+1
x_mean=x/count
y_mean=y/count
z_mean=z/count
print x_mean, y_mean, z_mean
Т. е. координаты геометрического центра лиганда - 38.9964186047 45.1291395349 27.5343023256. Эта информация нам пригодится далее для запуска докинга.¶
# Сохраняем белок без лиганда в отдельный pdb-файл
newpdb.writePDB("prot_without_lig.pdb")
# Готовим белок для докинга
prot = oddt.toolkit.readfile('pdb','prot_without_lig.pdb').next()
prot.OBMol.AddPolarHydrogens()
prot.OBMol.AutomaticPartialCharge()
print 'is it the first mol in 1lmp is protein?',prot.protein,':) and MW of this mol is:', prot.molwt
Молекулярная масса нашего белка - 13375.39512
Далее подготовим лиганды для докинга (NAG и его модификации, в которых вместо метильного радикала СH3C(=O)NH группы - OH, NH3+, H, Ph, COO-)¶
smiles = ['CC(=O)N[C@H]1[C@H](O)O[C@H](CO)[C@@H](O)[C@@H]1O','OC(=O)N[C@H]1[C@H](O)O[C@H](CO)[C@@H](O)[C@@H]1O','C(=O)N[C@H]1[C@H](O)O[C@H](CO)[C@@H](O)[C@@H]1O','[NH3+]C(=O)N[C@H]1[C@H](O)O[C@H](CO)[C@@H](O)[C@@H]1O','[O-]C(=O)C(=O)N[C@H]1[C@H](O)O[C@H](CO)[C@@H](O)[C@@H]1O','C2=CC=C(C=C2)C(=O)N[C@H]1[C@H](O)O[C@H](CO)[C@@H](O)[C@@H]1O']
mols= []
images =[]
for s in smiles:
m = oddt.toolkit.readstring('smi', s)
if not m.OBMol.Has3D():
m.make3D(forcefield='mmff94', steps=150)
m.removeh()
m.OBMol.AddPolarHydrogens()
mols.append(m)
###with print m.OBMol.Has3D() was found that:
### deep copy needed to keep 3D , write svg make mols flat
images.append((SVG(copy.deepcopy(m).write('svg'))))
display_svg(*images)
Проводим сам докинг¶
# Создаём объект для докинга
dock_obj= oddt.docking.AutodockVina.autodock_vina(
protein=prot,size=(20,20,20),center=[38.9964186047, 45.1291395349, 27.5343023256],
executable='/usr/bin/vina',autocleanup=True, num_modes=20)
print dock_obj.tmp_dir
print " ".join(dock_obj.params)
Выше можно видеть параметры, которые мы будем использовать для запуска докинга:
- --center_x (y,z): геометрический центр лиганда
- ---size_x (y,z): размер области, в которой будет производиться поиск
- ---cpu: сколько процессоров будет использоваться для выполнения нашей задачи
- --exhaustiveness: параметр, который говорит программе, насколько тщательно производить поиск
- --num_modes: количество конформаций лиганда в активном центре, которые будут найдены программой
- --energy_range: максимальная разница по энергии между лучшей и худшей конформацией из найденных
# запускаем сам докинг
res = dock_obj.dock(mols,prot)
# проанализируем результаты докинга
from operator import itemgetter
dock_res=[]
for i,r in enumerate(res):
hbs = oddt.interactions.hbonds(prot,r)#водородные связи
stack= oddt.interactions.pi_stacking(prot,r)#стэкинг-взаимодействия
phob = oddt.interactions.hydrophobic_contacts(prot,r)#гидрофобные контакты
dock_res.append([i,r.formula, float(r.data['vina_affinity']), float(r.data['vina_rmsd_ub']), r.residues[0].name, len(hbs[1]), len(stack[1]), len(phob[1])])
sort_res = sorted(dock_res, key=itemgetter(2))
for resi in sort_res:
print "%4d%10s%8s%8s%8s%4d%4d%4d" %(resi[0],resi[1],resi[2],resi[3],resi[4],resi[5],resi[6],resi[7])
Наилучшей оказалась структура под номером 45 (C13H17NO6), содержащая в виде заместителя фенил-радикал. Ещё можно посмотреть на следующую по энергии структуру - 36 (C8H13NO8), так как в ней лиганд имеет большее количество водородных связей с белком.
#запишем результаты в файлы pdb
for i,r in enumerate(res):
r.write(filename='r%s.pdb' %i, format='pdb')
from IPython.display import Image
Image(filename='docking.png')
Выше представлено визуальное отображение результатов докинга. Исходные лиганды структуры выделены зелёным цветом. Сиреневым цветов выделен лиганд из наихудшей структуры (44) по результатам докинга, синим (36) и красным (45) - двух наилучших. Две лучшие структуры попадают глубоко в карман белка в отличие от исходных лигандов и лиганда из худшей структуры, которые связываются ближе к поверхности. Вероятно, для некоторых модифицированных лигандов есть более выгодное положение в структуре, чем для самого NAG, однако даже самый худший по энергии результат докинга также кажется имеющим право на жизнь.