Учебный сайт Ксении Березиной

Назад к семестру

Хемоинформатика

Необходимо, используя пакет модулей RDkit, предложить аналог ибупрофена. Для этого нужно:

  • найти формулу ибупрофена и предложить способ изменения его SMILES для эмуляции продукта Click Chemistry;
  • на сайте PubChem найти все радикалы c азидом для Click Chemistry и скачать их SMILES нотации;
  • заменить в найденых радикалах азидную группу на модифицированный ибупрофен;
  • изобразить структуры модифицированных азидов;
  • отобрать те молекулы, которые удовлетворяют правилу пяти Липински.

Загрузим RDkit — пакет для работы с хемоинформатикой.

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

Нарисуем ибупрофен:

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

Правило пяти Липински говорит о том, какими свойствами должна обладать молекула, чтобы быть (скорее всего) активной как лекарство при оральном приеме:

  1. иметь не более 5 доноров водородных связей;
  2. иметь не более 10 акцепторов водородных связей;
  3. молекулярная масса не более 500 Да;
  4. log P (мера липофильности, коэффициент распределения в системе октанол-вода) меньше 5.

Посчитаем параметры правила пяти Липински для молекулы ибупрофена.

In [12]:
import rdkit.Chem.Lipinski as Lipinksy
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 [7]:
ibu_mod=Chem.MolFromSmiles('N1C=C(N=N1)C1=CC=C(C=C1)C(C)C(=O)O')
AllChem.Compute2DCoords(ibu_mod)
display(ibu_mod)

Скачали структуры всех радикалов с азидом в формате SMILES из базы PubChem ('azides_smiles.txt'). Извлечем все SMILES азидов:

In [2]:
strings=np.genfromtxt('azides_smiles.txt',dtype=np.str)
smiles = []
for line in strings:
    #formula length <30
    if len(line[1]) < 30 and not '.' in line[1]:
        smiles.append(line[1])

print smiles[0:100]
['C1CCC2=NN=NN2CC1', 'C1=CC2=NNN=C2C=C1', 'C1CC2=C(C1)NN=C2C3=NNN=N3', 'C1=CC=C2C(=C1)N=NN2N', 'C1=CC=C2C(=C1)N=NN2O', 'C1(=NNN=N1)C(C(=O)O)N', 'C1=CC(=CC=C1C', 'CCN1N=C(N=N1)C2=CCCN(C2)C', 'CC1=CC2=NNN=C2C=C1', 'C1(=NNN=N1)N', 'C1=NNN=N1', 'C1=CC=C(C=C1)N=[N+]=[N-]', 'C1=NC2=NNN=C2C(=N1)N', 'C1=CC2=NNN=C2C=C1C(=O)O', 'C1=CC(=CC=C1N)N=[N+]=[N-]', 'C1=CC=C2C(=C1)N=NN2Cl', 'CSC1=NNN=N1', 'C1=CC2=CC3=NNN=C3C=C2C=C1', 'C1=NN(N=N1)CC(=C(I)I)I', 'C1=CC(=CN=C1)C2=NNN=N2', 'C(CC1=NNN=N1)CN', 'C1=CC(=CC=C1N=[N+]=[N-])I', 'CC1=CC=CC(=N1)C(=O)NC2=NNN=N2', 'C1=NC2=C(N1)C(=NN=N2)N', 'CN1C2=C(C(=NC=N2)Cl)N=N1', 'C1[C@@H]2[C@@H](C2(C', 'CCOC1=NN=NN1C2=CC=CC=C2', 'CC1=CC=C(C=C1)N2C=C(N=N2)C=O', 'C(C(C(=O)O)N)N=[N+]=[N-]', '[N-]=[N+]=N[Pb]N=[N+]=[N-]', 'C1=NC=NC2=NNN=C21', 'C[Si](C)(C)N=[N+]=[N-]', 'C1=CC(=CC=C1C=O)N=[N+]=[N-]', 'CN1C(=NN=N1)SSC2=NN=NN2C', 'C(', 'C(C(CO)O)N=[N+]=[N-]', 'CN1C(=NN=N1)SSC', 'C1=CC=C2C(=C1)C=NN3C2=NN=N3', 'C1(=NNN=N1)[N+](=O)[O-]', 'C1=CC2=C(N=C1)N(N=N2)O', 'C1(=NNN=N1)[N+]', 'C1=CC(=CC=C1N=C=S)N=[N+]=[N-]', 'C1=CN2C(=NN=N2)N=C1', 'C1=CSC=C1C2=NN=NN2CC(=O)O', 'C1=CC2=NNN=C2C(=C1)C(=O)O', 'C1=CC=NC(=C1)N=[N+]=[N-]', 'C1=CN=CC=C1N=[N+]=[N-]', 'C1=CC(=CN=C1)N=[N+]=[N-]', 'CCN1C2=CC(=N[N+]', 'C1=CSC2=C3C(=NN=N3)N=CN21', 'CCN(CCC', 'CN(C)CCON1C2=CC=CC=C2N=N1', 'CN(C)CCCON1C2=CC=CC=C2N=N1', 'CN1C(=NN=N1)S(=O)(=O)CC(=O)OC', '[C-]', 'C1=CC(=C(C=C1N=[N+]=[N-])I)O', 'C1=CC(=CC=C1/C(=C/[N+]', '[3H]C([3H])([3H])N(CCCC(C', 'C1[C@@H]2[C@@H](C2(C', 'C1=CC2=NNN=C2C=C1[N+](=O)[O-]', 'C1=CC=C2C(=C1)N=NN2CCC(=O)O', 'CC1=NN(N=C1C(=O)N)C2=CC=CC=C2', 'CC1=CC=CC(=C1)C2=NNN=C2', 'CN1N=CC(=N1)C2=CCCNC2', 'C1=CC=C(C=C1)CC', 'CN1CCC=C(C1)C2=NN(N=C2)C', 'CCCCCC', 'C1=C(C=C(C2=NNN=C21)Cl)Cl', 'CN1C2=CC(=NC(=C2N=N1)C', 'C1=CC=C(C=C1)C2=NNN=C2', 'C1(=NNN=C1C(=O)O)C(=O)O', 'C1=CC(=CC=C1N2C=NN=N2)Cl', 'C1=CC=C(C=C1)CN2C(=C(N=N2)C', 'C1=CC=C(C(=C1)C2=NNN=N2)Cl', 'CC(CC(=O)O)N1C2=CC=CC=C2N=N1', 'C1=CC(=C(C=C1C2=NNN=N2)N)F', 'CC1=CC=C(C=C1)C2=NNN=C2', 'C1=CC(=CC=C1C2=NNN=C2)Cl', 'C1=CC=C(C=C1)C', 'CCC1=CC=C(C=C1)C2=NNN=N2', 'CC1=CC=C(O1)C2=NNN=N2', 'C1=CC(=CC=C1C2=NNN=N2)Br', 'CN1N=C(N=N1)NC(=O)C2=CC=CC=C2', 'CCCN1N=C(N=N1)NC(=O)C(C)C', 'CCN1C2=C(C=C(C=C2)C(=O)O)N=N1', 'CC(=O)NC1=CC=CC(=C1)C2=NNN=N2', 'C1=CC=C(C=C1)C2=NNN=C2C(=O)O', 'C(C1=NNN=N1)[C@@H](C(=O)O)N', 'CCCCOC1=CC=C(C=C1)C2=NNN=N2', 'C1=C(C=NC=C1Br)C2=NNN=N2', 'C(CC1=NNN=N1)CC2=NNN=N2', 'C1CCC(CC1)NC(=O)NC2=NNN=N2', 'CCC1=C(NN(C1=O)C2=NNN=N2)C', 'COC1=C(C=C2C(=C1)C(=C(C=N2)C', 'COC1=C(C=C2C(=C1)C(=C(C=N2)C', 'C1=CC(=CC(=C1)Br)C2=NNN=C2', 'C1=CC(=CC(=C1)Cl)C2=NNN=C2', 'C1=CC(=CC=C1C2=NNN=C2)Br', 'C1=CC=NC(=C1)C2=NNN=C2', 'CC1=CC=CC=C1C2=NNN=C2']

Теперь в скачанных SMILES азидов заменяем азидную группу на модифицированный ибупрофен (будем работать только с 1500 из азидов). Полученные структуры фильтруем по правилу Липински. А также отбираем молекулы более растворимые в воде, чем ибупрофен. Запоминаем самую водорастворимую и показываем ее в конце.

In [22]:
good_new_smiles=[] 
n=0
min = Lipinksy.rdMolDescriptors.CalcCrippenDescriptors(ibu)[0] #logP
min_index = 0
#Новую молекулу лучше создавать в try из-за возможных битых Smiles
for smi in smiles[:1500]:    
    if 'N=[N+]=[N-]' in smi:
        newsmi=smi.replace('N=[N+]=[N-]','N1C=C(N=N1)C1=CC=C(C=C1)C(C)C(=O)O') #replace with modified ibuprofen
    else:
        continue
    
    try:
        newmol=Chem.MolFromSmiles(newsmi)
        
        #see if Lipinski rule is followed AND molecule is more water-soluble than ibuprofen
        if Lipinksy.NumHDonors(newmol) <= 5 and Lipinksy.NumHAcceptors(newmol) <= 10 and Lipinksy.rdMolDescriptors.CalcExactMolWt(newmol) <=500 and Lipinksy.rdMolDescriptors.CalcCrippenDescriptors(newmol)[0] <= Lipinksy.rdMolDescriptors.CalcCrippenDescriptors(ibu)[0]:
            good_new_smiles.append(newmol)
            AllChem.Compute2DCoords(newmol)
            display(newmol)
            
            print n
            n=n+1
            if Lipinksy.rdMolDescriptors.CalcCrippenDescriptors(newmol)[0] < min: 
                min = Lipinksy.rdMolDescriptors.CalcCrippenDescriptors(newmol)[0]
                min_index = n-1

    except:
        pass
print 'good molecules: ', float(n)/float(1500)*100,'%'   
print 'molecule with min log P: molecule #%i with log P=%i' %(min_index,min)
display(good_new_smiles[min_index])
0
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
good molecules:  1.66666666667 %
molecule with min log P: molecule #21 with log P=0

Лучшая молекула - под номером 21. Так она выглядит в 3D:

In [24]:
m3d=Chem.AddHs(good_new_smiles[21])
Chem.AllChem.EmbedMolecule(m3d)
AllChem.MMFFOptimizeMolecule(m3d,maxIters=500,nonBondedThresh=200 )
Out[24]:
1
In [25]:
import nglview as nv
nv.show_rdkit(m3d)

Назад к семестру