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

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

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

ibu=Chem.MolFromSmiles('CC(C)CC1=CC=C(C=C1)C(C)C(=O)O')
AllChem.Compute2DCoords(ibu)
display(ibu)

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

import rdkit.Chem.Lipinski as Lipinksi
print(Lipinksi.NumHDonors(ibu))
print(Lipinksi.NumHAcceptors(ibu))
print(Lipinksi.rdMolDescriptors.CalcExactMolWt(ibu))
print(Lipinksi.rdMolDescriptors.CalcCrippenDescriptors(ibu)[0])

1
1
206.130679816
3.073200000000001

Мы видим, что у молекулы ибупрофена меньше пяти атомов - доноров водородной связи, меньше пяти акцепторов, его молекулярная масса меньше 500, а коэффициент гидрофобности меньше 5, так что эта молекула удовлетворяет требованиям правил Липнински.

Теперь скачаем из базы данных PubChem смайлсы всех соединений, содержащих в своём составе азидную группу, так как она необходима нам для клик-реакции.

import pubchempy as pcp

compounds = []
per_page = 10**5
for smiles in ["N=N=N", "NN#N",'N=[N+]=[N-]']:
    for i in range(200):
        try:
            a = pcp.get_properties(
              properties="CanonicalSMILES", 
              identifier=smiles, namespace="smiles", 
              searchtype="substructure",
              RingsNotEmbedded=True,
              listkey_count=per_page, listkey_start=i*per_page
            )
        except:
            break
        print("Retrieved page {} of {} search".format(i+1, smiles))
        compounds.extend(a)

Retrieved page 1 of N=N=N search
Retrieved page 2 of N=N=N search
Retrieved page 3 of N=N=N search
Retrieved page 1 of NN#N search
Retrieved page 1 of N=[N+]=[N-] search
Retrieved page 2 of N=[N+]=[N-] search
Retrieved page 3 of N=[N+]=[N-] search

len(a)

66906

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

neibu=Chem.MolFromSmiles('C#CCC1=CC=C(C=C1)C(C)C(=O)O')
AllChem.Compute2DCoords(neibu)
display(neibu)

В результате клик-реакции атомы азода замкнутся в цикл вокруг бывшей тройной связи. Прикинем, как будет выглядеть продукт такой циклизации, чтобы затем, используя этот SMILES, автоматически строить продукты клик-реакции.

Chem.MolFromSmiles('N1N=NC=C1(CC1=CC=C(C=C1)C(C)C(=O)O)')

neibu='N1N=NC=C1(CC1=CC=C(C=C1)C(C)C(=O)O)'

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

compounds=[]

for smi in a[:1500]:
    
    if 'N=[N+]=[N-]' in smi['CanonicalSMILES']:
        newsmi=smi['CanonicalSMILES'].replace('N=[N+]=[N-]', neibu)
    elif 'N=N=N' in smi['CanonicalSMILES']:
        newsmi=smi['CanonicalSMILES'].replace('N=N=N', neibu)
    elif 'NN#N' in smi['CanonicalSMILES']:
        newsmi=smi['CanonicalSMILES'].replace('NN#N', neibu)
    else:
        continue
    try:
        if newsmi.count('.')<2:
            newmol=Chem.MolFromSmiles(newsmi)
        else:
            pass

        
        if (Lipinksi.NumHDonors(newmol)<=5) and (Lipinksi.NumHAcceptors(newmol)<=10) \
        and (Lipinksi.rdMolDescriptors.CalcExactMolWt(newmol)<500) \
        and (Lipinksi.rdMolDescriptors.CalcCrippenDescriptors(newmol)[0]<=5):
            compounds.append(newmol)
    except:
        pass

len(compounds)

490

Итого мы получили 490 соединения - потенциальных лекартсв на основе ибупрофена. Посмотрим на первые 15 из них.

from IPython.display import SVG
display(Draw.MolsToGridImage(compounds[:15],useSVG=False, molsPerRow=3, subImgSize=(300, 300)))

Теперь для пятого соединения из нашего списка посмотрим карту сходства его с ибупрофеном.

from rdkit.Chem.Draw import SimilarityMaps

fp = SimilarityMaps.GetMorganFingerprint(compounds[4], fpType='bv')

fig, maxweight = SimilarityMaps.GetSimilarityMapForFingerprint(ibu, compounds[4], SimilarityMaps.GetMorganFingerprint)

Когда-нибудь на кодомо будет установлен модуль nglview, и мы сможем посмотреть на это соединение в 3D...

m3d=Chem.AddHs(compounds[4])
Chem.AllChem.EmbedMolecule(m3d)
AllChem.MMFFOptimizeMolecule(m3d,maxIters=500,nonBondedThresh=200 )

0

import nglview as nv

nv.show_rdkit(m3d)