Докинг низкомолекулярных лигандов в структуру белка.
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
# organic
from rdkit.Chem import Draw
from rdkit.Chem.Draw import IPythonConsole
import pmx # pdb manipulation
# protein preparation
pdb=pmx.Model('alignment_with_ligand.B99990001.pdb')
for r in pdb.residues[135:]:
print r #show residues
# create protein and ligand objects
newpdb = pdb.copy()
for r in newpdb.residues[-3:]:
newpdb.remove_residue(r)
newpdb.writePDB(fname="protein_without_ligand.pdb")
ligand = pdb.copy()
for r in ligand.residues[:-3]:
ligand.remove_residue(r)
# find geometric ligand center
center = np.mean([a.x for a in ligand.atoms],axis=0)
# prepare protein for docking
prot = oddt.toolkit.readfile('pdb','protein_without_ligand.pdb').next()
prot.OBMol.AddPolarHydrogens()
prot.OBMol.AutomaticPartialCharge()
NAG содержит в себе CH3C(=O)NH группу
Создадим лиганды где метильный радикал этой группы будет заменён на:
OH
NH3+
H
Ph
COO-
# NAG SMILES notation: CC(=O)NC1C(C(C(OC1O)CO)O)O
# list ligands for docking
smiles = ['OC(=O)NC1C(C(C(OC1O)CO)O)O',
'[NH3+]C(=O)NC1C(C(C(OC1O)CO)O)O',
'C(=O)NC1C(C(C(OC1O)CO)O)O',
'C1=CC=C(C=C1)C(=O)NC1C(C(C(OC1O)CO)O)O',
'[O-]C(=O)C(=O)NC1C(C(C(OC1O)CO)O)O']
ligands = []
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()
ligands.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)
#create docking object
dock_obj= oddt.docking.AutodockVina.autodock_vina(protein=prot,size=(20,20,20),center=center,
executable='/usr/bin/vina',autocleanup=True, num_modes=9)
print dock_obj.tmp_dir
# vina parameters:
# --center - coordinates of the center
# --size - size by dimensions
# --cpu - the number of CPUs to use
# --exhaustiveness - exhaustiveness of the global search (roughly proportional to time)
# --num_modes - maximum number of binding modes to generate
# energy_range - maximum energy difference between the best binding mode and the worst one displayed (kcal/mol)
print " ".join(dock_obj.params)
# do docking
res = dock_obj.dock(ligands,prot)
import pandas as pd
# show results: affinity for every ligand
formula = []
affinity = []
for i,r in enumerate(res):
formula.append(r.formula)
affinity.append(r.data['vina_affinity'])
res_table = pd.DataFrame()
res_table['ligand'] = formula
res_table['affinity'] = affinity
# sort by affinity
res_table.sort_values('affinity', inplace=True)
res_table
# choose best ligand by its affinity to protein: C7H13NO6
print res_table.iloc[0].tolist()
# save pdbs for protein with best docked ligand
for i,r in enumerate(res):
if r.formula == 'C7H13NO6':
r.write(filename='docking-%s-%s.pdb' % (r.formula, i), format='pdb', overwrite=True)
# Launch PyMOL
import __main__
__main__.pymol_argv = [ 'pymol', 'x' ]
import pymol
from pymol import cmd
pymol.finish_launching()
# Set image
import time
from IPython.display import Image
defaultImage = 'pymolimg.png'
def prepareImage(width=500, height=500, sleep=2, filename=defaultImage):
cmd.ray(width, height)
cmd.png(filename)
time.sleep(sleep)
cmd.load('docking-C7H13NO6-26.pdb', 'ligand')
cmd.load('protein_without_ligand.pdb', 'protein')
cmd.show_as('sticks', 'ligand')
cmd.show_as('cartoon', 'protein')
cmd.color('white','protein')
cmd.center('ligand')
cmd.origin('ligand')
cmd.zoom('ligand', 4)
prepareImage()
Image(defaultImage)
%%bash
jupyter nbconvert --to html hw7_filippova.ipynb