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 import Chem
from rdkit.Chem.Draw import IPythonConsole
import mdtraj as md
u = md.load('TRINI_L.B99990001.pdb')
pdb = u.topology
for i,r in enumerate(pdb.atoms):
print(i,r, end=("\t" if i % 8 else "\n"))
#посмотрим на атомы
Последний атом белка #1134.
Последний атом лиганда #1177.
# сохраним pdb без лиганда
u2 = u.atom_slice(range(1134 + 1))
u2.save_pdb('TRINI_L.noligand.pdb')
Дальше найдем геометрический центр лиганда.
ligand = u.atom_slice(range(1134 + 1, 1177 + 1))
center = np.mean(ligand.xyz, axis=1) * 10 # To angstroms
center = center.reshape(3,)
center
prot = next(oddt.toolkits.rdk.readfile('pdb', 'TRINI_L.noligand.pdb'))
prot.protein, prot.molwt
smiles = ["CC(=O)NC1C(O)OC(CO)C(O)C1O", # NAG
"OC(=O)NC1C(O)OC(CO)C(O)C1O", # -OH
"[NH3+]C(=O)NC1C(O)OC(CO)C(O)C1O", # -NH3+
"C(=O)NC1C(O)OC(CO)C(O)C1O", # -H
"C=1(C=CC=CC1)C(=O)NC1C(O)OC(CO)C(O)C1O", # -Ph
"C(C([O-])=O)(=O)NC1C(O)OC(CO)C(O)C1O" # -COO-
]
mols= []
images =[]
for s in smiles:
m = oddt.toolkits.rdk.readstring('smi', s)
m.make3D(forcefield='mmff94', steps=150)
mols.append(m)
print(Chem.rdMolDescriptors.CalcMolFormula(m.Mol))
images.append((SVG(copy.deepcopy(m).write('svg'))))
display(*images)
Подготовили для докинга лиганды и белок. Отлично! Теперь... дочим!
#create docking object
dock_obj= oddt.docking.AutodockVina.autodock_vina(
protein=prot,size=(20,20,20),center=center,
executable='./bin/vina',autocleanup=True, num_modes=9,
prefix_dir='.')
print(dock_obj.tmp_dir)
print(" ".join(dock_obj.params))
# do it! *do it* do it. do-it-do-it-do-it-do-it!
res = dock_obj.dock(mols, prot)
И смотрим на результат:
for i,r in enumerate(res):
### Узнал дескриптор "формула" =)
formula = Chem.rdMolDescriptors.CalcMolFormula(r.Mol)
print("%4d%10s%8s%8s%8s" % (i, formula, r.data['vina_affinity'], r.data['vina_rmsd_ub'], r.residues[0].name))
import pandas as pd
results_df = pd.DataFrame([[Chem.rdMolDescriptors.CalcMolFormula(r.Mol), r.data['vina_affinity'], r.data['vina_rmsd_ub'], r.residues[0].name] for r in res],
columns=['brutto', 'affinity', 'ub', 'name'])
results_df.sort_values(by=['affinity'], ascending=False)
Ну и, наконец, смоделируем несколько лигандов!
indices = [36, 37, 38, 9, 0]
for i,r in [(i, res[i]) for i in indices]:
r.write(filename='r%s.pdb' % i, format='pdb')
import __main__
__main__.pymol_argv = [ 'pymol', '-x' ]
import pymol
pymol.finish_launching()
from pymol import cmd
cmd.reinitialize()
cmd.bg_color('white')
cmd.set('ray_trace_mode', 1)
cmd.load('TRINI_L.noligand.pdb', 'protein')
cmd.load('ligand.pdb', 'lig_x')
cmd.h_add('protein')
cmd.h_add('lig_x')
cmd.show('surface', 'protein')
for i in indices:
cmd.load(f'r{i}.pdb', f'r_{i}')
cmd.zoom('lig_x')
cmd.do('ray; png lig_x.png')
(Исходный лиганд)
Image('lig_x.png')
cmd.center('r_36')
cmd.do('ray; png r36.png')
Image('r36.png')
cmd.center('r_37')
cmd.do('ray; png r37.png')
Image('r37.png')
cmd.center('r_38')
cmd.do('ray; png r38.png')
Image('r38.png')
find_hbonds('r_0', 'protein')
cmd.center('r_0')
cmd.do('ray; png r0.png')
Image('r0.png')
cmd.center('r_9')
cmd.do('ray; png r9.png')
Image('r9.png')
Судя по всему, рассмотренные нами лиганды не очень хорошо подходят к этому белку: водородных связей немного (но я не рассмотрел гидрофобные взаимодействия), поэтому, скорее всего, в реальности взаимодействия происходить не будет.