Практикум 3: хемоинформатика

In [1]:
from rdkit import Chem
from rdkit.Chem import AllChem
from rdkit import RDConfig
from rdkit.Chem.Draw import IPythonConsole 
from rdkit.Chem import Draw
import numpy as np
from IPython.display import display,Image
import rdkit.Chem.Lipinski as Lipinksy
In [2]:
# Does molecule suits Lipinsky rule?
def True_Lipinsky(mol):
    don = Lipinksy.NumHDonors(mol) <= 5
    acc = Lipinksy.NumHAcceptors(mol) <= 10
    mw = Lipinksy.rdMolDescriptors.CalcExactMolWt(mol) <= 500
    coef = Lipinksy.rdMolDescriptors.CalcCrippenDescriptors(mol)[0] < 5
    if don and acc and mw and coef:
        return True

Построим ибупрофен по SMILES

In [3]:
ibu_form = 'CC(C)CC1=CC=C(C=C1)C(C)C(=O)O'
ibu=Chem.MolFromSmiles(ibu_form)
AllChem.Compute2DCoords(ibu)
display(ibu)

Проверим, что он соотвествует правилам Липински:

  1. Молекула имеет не более 5 доноров водородных связей
  2. Молекула имеет не более 10 акцепторов водородных связей
  3. Молекулярная масса менее 500 Да
  4. log P (мера липофильности, коэффициент распределения в системе октанол-1 - вода) меньше 5.
In [4]:
print(Lipinksy.NumHDonors(ibu))
print Lipinksy.NumHAcceptors(ibu)
print Lipinksy.rdMolDescriptors.CalcExactMolWt(ibu)
print Lipinksy.rdMolDescriptors.CalcCrippenDescriptors(ibu)[0]

1
1
206.130679816
3.0732

Теперь нарисуем производное ибупрофена, которое могло бы участвовать в азид-алкиновом циклоприсоединении. Для этого в молекуле должна быть концевая тройная связь.

In [5]:
ibu_m_f = 'C#CCC1=CC=C(C=C1)C(C)C(=O)O'
ibu_m=Chem.MolFromSmiles(ibu_m_f)
AllChem.Compute2DCoords(ibu_m)
display(ibu_m)

Теперь сделаем фрагмент, кторый будет получаться после реакции циклоприсоединения.

In [6]:
ibu_m_f = 'N1C=C(CC2=CC=C(C=C2)C(C)C(=O)O)N=N1'
ibu_m=Chem.MolFromSmiles(ibu_m_f)
AllChem.Compute2DCoords(ibu_m)
display(ibu_m)
In [7]:
import pubchempy as pcp

Теперь найдем радикалы с азидом, подходящие для клик-химии. Я попробовала скачать необходимые вещества программно, но не разобралась, как задавать дополнительные опции поиска. Поэтому здесь привожу не совсем удачный вариант кода, сам файл со SMILES'aми получила скачала вручную.

azide = '[N-]=[N+]=[N-]'
per_page = 10**2
list = []
for i in range(20):
    try:
        list1 = pcp.get_compounds(azide,
                             'smiles', 
                             searchtype='substructure',  
                             listkey_count = per_page,
                             listkey_start = i*per_page
                            )
    except:
        pass
    list.extend(list1)    
print(len(list))

Выделим из полученного файла SMILES'ы без разрывов, не слишком динные (и избавимся от дефектных)

In [8]:
strings=np.genfromtxt('pr3/smiles.txt',dtype=np.str)
smiles = []
for line in strings:
    if len(line[1]) < 30 and not '.' in line[1] and len(line[1])>5:
        smiles.append(line[1])
        
print(smiles[:10])
print(len(smiles))

['C=C(C(=O)NCC', 'COC1=NC(=NC=C1N=[N+]=[N-])OC', 'CC(C(=O)OCCC', 'C[C@@H](C(=O)NC1=CC=CC(=N1)C', 'CCOP(=O)(C(C', 'CC1=NN(C(=C1CC(C', 'C1CC=C(C(C1)N=[N+]=[N-])I', 'CCOCCOCCOCN=[N+]=[N-]', 'CC(C(=O)OCCN=[N+]=[N-])Cl', 'C1=CC=C(C=C1)C(C']
12544

Построим новые молекулы, используя полученные радикалы и измененный ибупрофен. Сохраним в отдельный список те, которые удовлетворяют правилам Липински.

In [9]:
azide = 'N=[N+]=[N-]'
right_mol = []
right_mol_Lipinsky = []
for smi in smiles[:1500]:
    if azide in smi:
        newsmi=smi.replace('N=[N+]=[N-]',ibu_m_f)
        right_mol.append(newsmi)
    else:
        continue
    try:
        newmol=Chem.MolFromSmiles(newsmi)        
        if True_Lipinsky(newmol):
            right_mol_Lipinsky.append(newmol)
    except:
        pass

RDKit ERROR: [03:40:06] SMILES Parse Error: extra open parentheses for input: 'CC(C)(CC(CN1C=C(CC2=CC=C(C=C2)C(C)C(=O)O)N=N1)C'
RDKit ERROR: [03:40:06] SMILES Parse Error: extra open parentheses for input: 'C(/C=C(/CN1C=C(CC2=CC=C(C=C2)C(C)C(=O)O)N=N1)\C'
RDKit ERROR: [03:40:06] SMILES Parse Error: extra open parentheses for input: 'C1=CC2=C(C=C1N(CN1C=C(CC2=CC=C(C=C2)C(C)C(=O)O)N=N1)C'
RDKit ERROR: [03:40:06] SMILES Parse Error: extra open parentheses for input: 'CC(CCC(=O)NCCCN1C=C(CC2=CC=C(C=C2)C(C)C(=O)O)N=N1)(C'
RDKit ERROR: [03:40:06] SMILES Parse Error: extra open parentheses for input: 'CC(CCC(=O)NCCCN1C=C(CC2=CC=C(C=C2)C(C)C(=O)O)N=N1)(C'
RDKit ERROR: [03:40:07] SMILES Parse Error: extra open parentheses for input: 'C/C=C/CC(CCCN1C=C(CC2=CC=C(C=C2)C(C)C(=O)O)N=N1)(C'
RDKit ERROR: [03:40:07] SMILES Parse Error: extra open parentheses for input: 'C=CCC(CCN1C=C(CC2=CC=C(C=C2)C(C)C(=O)O)N=N1)(C'
RDKit ERROR: [03:40:09] SMILES Parse Error: extra open parentheses for input: 'C1CCC(C(C1)N1C=C(CC2=CC=C(C=C2)C(C)C(=O)O)N=N1)(C'
RDKit ERROR: [03:40:09] SMILES Parse Error: extra open parentheses for input: 'C1=CC(=C(C=C1N1C=C(CC2=CC=C(C=C2)C(C)C(=O)O)N=N1)C'

Нарисуем 50 первых молекул, проверим, что все в порядке.

In [10]:
print(len(right_mol_Lipinsky))
Draw.MolsToGridImage(right_mol_Lipinsky[:50], molsPerRow=5)

1318
Out[10]:

Теперь построим карту схожести ибупрофена с веществами из массива.

In [30]:
from rdkit.Chem.Draw import SimilarityMaps
fig, maxweight = SimilarityMaps.GetSimilarityMapForFingerprint(ibu, right_mol_Lipinsky[5], SimilarityMaps.GetMorganFingerprint)
print(maxweight)

Построим 3D-структуру данной молекулы

In [12]:
m3d=Chem.AddHs(right_mol_Lipinsky[5])
Chem.AllChem.EmbedMolecule(m3d)
AllChem.MMFFOptimizeMolecule(m3d,maxIters=500,nonBondedThresh=200 )
Out[12]:
0

Попробовала посмотреть в nglview, но почему-то не получилось.

import nglview as nv
nv.show_rdkit(m3d)

Поэтому воспользуемся паймолом. Файл я получила с помощью программы obgen из пакета openbabel:

obgen mol_3d.smi > mol_3d_ready.mol
In [18]:
f = open('mol_3d.smi', 'w')
f.write(Chem.MolToSmiles(m3d))
f.close()
In [23]:
from xmlrpclib import ServerProxy
import __main__
__main__.pymol_argv = [ 'pymol', '-x' ]
import pymol
pymol.finish_launching()
from pymol import cmd,stored
In [28]:
cmd.load('mol_3d_ready.mol')
cmd.orient()
cmd.bg_color('white')
cmd.png(filename='1', ray=1)