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
RDKit WARNING: [12:02:42] Enabling RDKit 2019.09.3 jupyter extensions
Нарисуем ибупрофен
ibu=Chem.MolFromSmiles('CC(C)CC1=CC=C(C=C1)C(C)C(=O)O')
AllChem.Compute2DCoords(ibu)
display(ibu)
Посчитаем параметры для правила Лепински: соединение, чтобы быть похожим на лекарство, должно:
• Иметь менее 5 атомов-доноров водородной связи;
• Обладать молекулярным весом менее 500;
• Иметь липофильность (log P — коэффициент распределения вещества на границе раздела вода-октанол) менее 5;
• Иметь суммарно не более 10 атомов азота и кислорода (грубая оценка количества акцепторов водородной связи)
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.073200000000001
Модифицируем ибупрофен так, чтобы к нему можно было присоединить азид, то есть сделаем ему концевую алкильную группу.
ibuin=Chem.MolFromSmiles('CC(C#C)CC1=CC=C(C=C1)C(C)C(=O)O')
AllChem.Compute2DCoords(ibuin)
display(ibuin)
В результате клик-реакции получим такое соединение:
ibuzid=Chem.MolFromSmiles('N1C=C(N=N1)C(C)CC1=CC=C(C=C1)C(C)C(=O)O')
AllChem.Compute2DCoords(ibuzid)
display(ibuzid)
Скачаем как можно больше различных азидов.
import pubchempy as pcp
# Скачивать долго, контроль, чтобы случайно не запустить заново.
computed = True
if computed == False:
azides = []
num_per_page = 100000
#Ищем структуры, содержащие N=N=N, NN#N или N#NN в качестве подструктур
for smiles in ['N=[N+]=[N-]', '[N-]=[N+]=[N-]', '[N-][N+]#[N]', '[N]#[N+][N-]', '[N-][N]#[N+]', '[N+]#[N][N-]']:
for i in range(10000):
try:
res = pcp.get_properties(
properties="CanonicalSMILES",
identifier=smiles,
namespace="smiles",
searchtype="substructure",
RingsNotEmbedded=True,
listkey_count=num_per_page, listkey_start=i*num_per_page
)
except:
break
print(f"Fetched page {i+1} of {smiles}".format(i+1, smiles))
azides.extend(res)
Fetched page 1 of N=[N+]=[N-] Fetched page 2 of N=[N+]=[N-] Fetched page 3 of N=[N+]=[N-] Fetched page 4 of N=[N+]=[N-] Fetched page 5 of N=[N+]=[N-] Fetched page 6 of N=[N+]=[N-] Fetched page 7 of N=[N+]=[N-] Fetched page 8 of N=[N+]=[N-] Fetched page 9 of N=[N+]=[N-] Fetched page 10 of N=[N+]=[N-] Fetched page 11 of N=[N+]=[N-] Fetched page 12 of N=[N+]=[N-] Fetched page 13 of N=[N+]=[N-] Fetched page 14 of N=[N+]=[N-] Fetched page 15 of N=[N+]=[N-] Fetched page 16 of N=[N+]=[N-] Fetched page 17 of N=[N+]=[N-] Fetched page 1 of [N-]=[N+]=[N-] Fetched page 2 of [N-]=[N+]=[N-] Fetched page 3 of [N-]=[N+]=[N-] Fetched page 4 of [N-]=[N+]=[N-] Fetched page 5 of [N-]=[N+]=[N-] Fetched page 6 of [N-]=[N+]=[N-] Fetched page 7 of [N-]=[N+]=[N-] Fetched page 8 of [N-]=[N+]=[N-] Fetched page 9 of [N-]=[N+]=[N-] Fetched page 10 of [N-]=[N+]=[N-] Fetched page 11 of [N-]=[N+]=[N-] Fetched page 12 of [N-]=[N+]=[N-] Fetched page 13 of [N-]=[N+]=[N-] Fetched page 14 of [N-]=[N+]=[N-] Fetched page 15 of [N-]=[N+]=[N-] Fetched page 16 of [N-]=[N+]=[N-] Fetched page 1 of [N-][N+]#[N] Fetched page 2 of [N-][N+]#[N] Fetched page 3 of [N-][N+]#[N] Fetched page 4 of [N-][N+]#[N] Fetched page 5 of [N-][N+]#[N] Fetched page 6 of [N-][N+]#[N] Fetched page 7 of [N-][N+]#[N] Fetched page 8 of [N-][N+]#[N] Fetched page 9 of [N-][N+]#[N] Fetched page 10 of [N-][N+]#[N] Fetched page 11 of [N-][N+]#[N] Fetched page 12 of [N-][N+]#[N] Fetched page 13 of [N-][N+]#[N] Fetched page 14 of [N-][N+]#[N] Fetched page 15 of [N-][N+]#[N] Fetched page 16 of [N-][N+]#[N] Fetched page 1 of [N]#[N+][N-] Fetched page 2 of [N]#[N+][N-] Fetched page 3 of [N]#[N+][N-] Fetched page 4 of [N]#[N+][N-] Fetched page 5 of [N]#[N+][N-] Fetched page 6 of [N]#[N+][N-] Fetched page 7 of [N]#[N+][N-] Fetched page 8 of [N]#[N+][N-] Fetched page 9 of [N]#[N+][N-] Fetched page 10 of [N]#[N+][N-] Fetched page 11 of [N]#[N+][N-] Fetched page 12 of [N]#[N+][N-] Fetched page 13 of [N]#[N+][N-] Fetched page 14 of [N]#[N+][N-] Fetched page 15 of [N]#[N+][N-] Fetched page 16 of [N]#[N+][N-] Fetched page 1 of [N-][N]#[N+] Fetched page 1 of [N+]#[N][N-]
len(azides)
6502750
azides[0]
{'CID': 155859234, 'CanonicalSMILES': 'CC(=O)OCC1C(C(C(C(O1)OC2C(C(OC3C2OC(OC3)C4=CC=CC=C4)OC5=CC=C(C=C5)OC)N=[N+]=[N-])OC(=O)C)OC(=O)C)OC(=O)C'}
smiles_list = [x["CanonicalSMILES"] for x in azides]
print(len(smiles_list))
smiles_list = list(filter(lambda x: len(x) < 30 and not '.' in x,
smiles_list))
print(len(smiles_list))
6502750 869853
from random import shuffle, seed
groups = ['N=[N+]=[N-]', '[N-]=[N+]=[N-]', '[N-][N+]#[N]', '[N]#[N+][N-]', '[N-][N]#[N+]', '[N+]#[N][N-]']
template = 'N1C=C(N=N1)C(C)CC1=CC=C(C=C1)C(C)C(=O)O'
result = []
num_to_check = 1500
# Новую молекулу лучше создавать в try из-за битых Smiles
#При первом запуске ячейки убрать решетки, чтобы перемешать
#shuffle(smiles_list)
for smi in smiles_list[:num_to_check]:
for group in groups:
if group in smi:
newsmi = smi.replace(group, template)
break
else:
continue
# Если новая молекула удолтворяет правилу пяти, сохраним в массив и покажем
try:
new_mol=Chem.MolFromSmiles(newsmi)
if (Lipinksy.NumHDonors(new_mol) <= 5 and
Lipinksy.NumHAcceptors(new_mol) <= 10 and
Lipinksy.rdMolDescriptors.CalcExactMolWt(new_mol) <=500 and
Lipinksy.rdMolDescriptors.CalcCrippenDescriptors(new_mol)[0] <= 5):
result.append((new_mol, Lipinksy.rdMolDescriptors.CalcCrippenDescriptors(new_mol)[0]))
except Exception:
continue
print(len(result))
result.sort(key=lambda x: x[1])
mols = [x[0] for x in result]
1309
Посмотрим на топ-12 гидрофильных молекул. Вот они:
Draw.MolsToGridImage(mols[:12],useSVG=True, molsPerRow=3, subImgSize=(250, 200))
Посмотрим на 5-го представителя из наиболее гидрофильных молекул.
from rdkit.Chem.Draw import SimilarityMaps
somemol = mols[4]
SimilarityMaps.GetMorganFingerprint(somemol, fpType='bv')
SimilarityMaps.GetSimilarityMapForFingerprint(ibu,
somemol, SimilarityMaps.GetMorganFingerprint)
(<Figure size 180x180 with 1 Axes>, 0.17876913057635951)
В самом деле, остаток ибупрофена очень похож на ибупрофен, а присоединенный азид не похож.
Теперь получим 3D структуру этой молекулы и покрутим ее.
m3d=Chem.AddHs(somemol)
Chem.AllChem.EmbedMolecule(m3d)
AllChem.MMFFOptimizeMolecule(m3d,maxIters=500,nonBondedThresh=200 )
0
import nglview as nv
nv.show_rdkit(m3d)
Молекулу покрутили, но в HTML она не экспортируется, было очень красиво, всем советую.
Так что просто покажем структуру в виде снимка.
display(m3d)