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

Загрузим модули RDkit:

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

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

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

Модифицируем ибупрофен для клик-химии:

ibu_mod=Chem.MolFromSmiles('CC(C)CC1=CC=C(C=C1)C(C1=CN=NN1)(C)C(=O)O')
AllChem.Compute2DCoords(ibu_mod)
display(ibu_mod)

И перевернем его:

ibu_template = Chem.MolFromSmiles('N1N=NC(=C1)C(C)(C(=O)O)C1=CC=C(C=C1)CC(C)C')
AllChem.Compute2DCoords(ibu_template)
display(ibu_template)

Найдем соединения, содержащие азидную группу, для клик-химии:

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(compounds)

517961

Соединений нашлось очень много... Проверим какое-то количество из них на удовлетворение правилам Липнински:

final_compounds = []

for smi in compounds[:2000]:
    
    if 'N=[N+]=[N-]' in smi['CanonicalSMILES']:
        newsmi=smi['CanonicalSMILES'].replace('N=[N+]=[N-]', 'N1N=NC(=C1)C(C)(C(=O)O)C1=CC=C(C=C1)CC(C)C')
    elif 'N=N=N' in smi['CanonicalSMILES']:
        newsmi=smi['CanonicalSMILES'].replace('N=N=N', 'N1N=NC(=C1)C(C)(C(=O)O)C1=CC=C(C=C1)CC(C)C')
    elif 'NN#N' in smi['CanonicalSMILES']:
        newsmi=smi['CanonicalSMILES'].replace('NN#N', 'N1N=NC(=C1)C(C)(C(=O)O)C1=CC=C(C=C1)CC(C)C')
    else:
        continue
        
    try:
        if newsmi.count('.')==0:
            newmol = Chem.MolFromSmiles(newsmi)
            if (Lipinksi.NumHDonors(newmol)<5) and (Lipinksi.NumHAcceptors(newmol)<10) and (Lipinksi.rdMolDescriptors.CalcExactMolWt(newmol)<500) and (Lipinksi.rdMolDescriptors.CalcCrippenDescriptors(newmol)[0]<5):
                final_compounds.append(newmol)
    except:
        pass

len(final_compounds)

126

Т.е. выбрали 126 соединений. Посмотрим на некоторые из них.

from IPython.display import SVG
display(Draw.MolsToGridImage(final_compounds[:6],useSVG=True, molsPerRow=3, subImgSize=(200, 200)))

Построим Similiraty Map с пятым веществом из этого массива.

from rdkit.Chem.Draw import SimilarityMaps

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

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

Интересно, что то, что раньше было ибупрофеном, теперь покрашено розовым - т.е. не похоже на него.

Теперь хотим посмотреть на 3D структуру.

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

1

import nglview as nv

nv.show_rdkit(m3d)

NGLWidget()

Молкула, в принципе, открылась, но в html оно в итоге не отображается. Покажем хотя бы картинкой.

display(Image(filename=f'pr3_3d.PNG', width=200))