drrdom / crem Goto Github PK
View Code? Open in Web Editor NEWCReM: chemically reasonable mutations framework
License: BSD 3-Clause "New" or "Revised" License
CReM: chemically reasonable mutations framework
License: BSD 3-Clause "New" or "Revised" License
Hi, thanks for creating this useful repo which allows us to generate new molecules from a framework!
I have been playing around with the grow_mol function, and I want to specify 3 positions in my starting molecule to grow the molecule. I can do so with the parameter replace_ids and specify the 3 positions in a list. However, I realised that all of the generated molecules only had 1 replacement per molecule. What I am looking for is a way to possibly generate new molecules with 2 or even all 3 positions replaced at once.
I have read the documentation for crem, and as far as I'm aware, there is no way to do this currently. Do you have any advice regarding this issue?
Thank you!
Alan
Hi! Thanks for maintaining this repo! I was wondering which parameters were relevant for speeding up mutate_mol
? I've found better speeds by increasing n_cores
but if I wanted to speed it up even more would decreasing radius
help too?
Thank you for creating this valuable tool. I'm wondering how I can visualize the database I've created following the documentation. Also, could you please summarize the variations I can expect in the outcome database if I've created it with different radii (for example: 1, 2, 3, 4, 5)? Additionally, what if I only want to create fragmented databases with different heteroatoms and not couple them?
Hi! I followed the instruction to install the crem package. However, when I ran the demo code, I found the below issues:
---------------------------------------------------------------------------
OperationalError Traceback (most recent call last)
/tmp/ipykernel_1568311/2420725268.py in <module>
1 m = Chem.MolFromSmiles('c1cc(OC)ccc1C') # methoxytoluene
----> 2 mols = list(mutate_mol(m, db_name='replacements.db', max_size=1,radius=3))
/home/shuangjia/anaconda3/lib/python3.7/site-packages/crem/crem.py in mutate_mol(mol, db_name, radius, min_size, max_size, min_rel_size, max_rel_size, min_inc, max_inc, max_replacements, replace_cycles, replace_ids, protected_ids, symmetry_fixes, min_freq, return_rxn, return_rxn_freq, return_mol, ncores, **kwargs)
493 protected_ids_1=protected_ids, protected_ids_2=None,
494 min_freq=min_freq, symmetry_fixes=symmetry_fixes,
--> 495 **kwargs):
496 for smi, m, rxn in __frag_replace(mol, None, frag_sma, core_sma, radius, ids, None):
497 if max_replacements is None or len(products) < (max_replacements + 1): # +1 because we added source mol to output smiles
/home/shuangjia/anaconda3/lib/python3.7/site-packages/crem/crem.py in __gen_replacements(mol1, mol2, db_name, radius, dist, min_size, max_size, min_rel_size, max_rel_size, min_inc, max_inc, max_replacements, replace_cycles, protected_ids_1, protected_ids_2, min_freq, symmetry_fixes, **kwargs)
342 max_atoms = num_heavy_atoms + max_inc
343
--> 344 row_ids = __get_replacements_rowids(cur, env, dist, min_atoms, max_atoms, radius, min_freq, **kwargs)
345
346 if max_replacements is None:
/home/shuangjia/anaconda3/lib/python3.7/site-packages/crem/crem.py in __get_replacements_rowids(db_cur, env, dist, min_atoms, max_atoms, radius, min_freq, **kwargs)
284 else:
285 sql += f" AND {k} = {v}"
--> 286 db_cur.execute(sql)
287 return set(i[0] for i in db_cur.fetchall())
288
OperationalError: no such table: radius3
I'm not sure if it's a bug after the package update, whether I use pypi to install or manually install this problem exists, look forward to your answer, thank you!
Implement random seed to make results reproducible
Was a bit confused by the fragment_mol
method in fragmentation.py
, specifically these lines:
frags = rdMMPA.FragmentMol(mol, pattern="[!#1]!@!=!#[!#1]", maxCuts=4, resultsAsMols=False, maxCutBonds=30)
frags += rdMMPA.FragmentMol(mol, pattern="[!#1]!@!=!#[!#1]", maxCuts=3, resultsAsMols=False, maxCutBonds=30)
Shouldn't maxCuts=3 be included as a subset of all fragmentations for maxCuts=4?
What exactly is the difference between maxCutBonds and maxCuts? I tried looking for documentation on this and was unfortunately unable to find it.
Would it make more sense for the maxCuts parameter to be user-defined rather than preset?
Thanks for your help!
Hello,
Thank you for maintain CReM! I ran into an error while using the program, and would like to see if you have any insights into this problem.
I'm trying to grow a fragment generated using RGroupDecomposition in rdkit, like so:
from rdkit import Chem
from rdkit.Chem import rdRGroupDecomposition as rgd
from crem.crem import grow_mol
mol = Chem.MolFromSmiles("Cc4cccc(Cc3cnc(NC(=O)c1cnn2ccc(C)nc12)s3)c4")
core = Chem.MolFromSmiles("Cc1cnc(NC(=O))s1")
params = rgd.RGroupDecompositionParameters()
params.RGroupLabels = rgd.RGroupLabels.AtomIndexLabels
fragment = rgd.RGroupDecompose([core], [mol], options=params)[0][0]['R1']
display(fragment)
The R1:1 is a dummy atom generated by RGroupDecomposition. If I try to run grow_mol at this point, I receive this error:
derivatives = list(grow_mol(r, db_name='/home/javi/data/rgroup/replacements02_sa2.db', max_replacements=9, return_mol=True))
>>>
ValueError Traceback (most recent call last)
crem.ipynb Cell 3' in <cell line: 62>()
---> 66 derivatives = list(grow_mol(r, db_name='replacements02_sa2.db', max_replacements=9, return_mol=True))
File ~/workspace/sci-dev/sci-dev-venv/lib/python3.8/site-packages/crem/crem.py:487, in mutate_mol(mol, db_name, radius, min_size, max_size, min_rel_size, max_rel_size, min_inc, max_inc, max_replacements, replace_cycles, replace_ids, protected_ids, symmetry_fixes, min_freq, return_rxn, return_rxn_freq, return_mol, ncores, **kwargs)
483 protected_ids = sorted(protected_ids)
485 if ncores == 1:
--> 487 for frag_sma, core_sma, freq, ids in __gen_replacements(mol1=mol, mol2=None, db_name=db_name, radius=radius,
488 min_size=min_size, max_size=max_size,
489 min_rel_size=min_rel_size, max_rel_size=max_rel_size,
490 min_inc=min_inc, max_inc=max_inc,
491 max_replacements=max_replacements,
492 replace_cycles=replace_cycles,
493 protected_ids_1=protected_ids, protected_ids_2=None,
494 min_freq=min_freq, symmetry_fixes=symmetry_fixes,
495 **kwargs):
496 for smi, m, rxn in __frag_replace(mol, None, frag_sma, core_sma, radius, ids, None):
497 if max_replacements is None or (max_replacements is not None and len(products) < max_replacements):
File ~/workspace/sci-dev/sci-dev-venv/lib/python3.8/site-packages/crem/crem.py:339, in __gen_replacements(mol1, mol2, db_name, radius, dist, min_size, max_size, min_rel_size, max_rel_size, min_inc, max_inc, max_replacements, replace_cycles, protected_ids_1, protected_ids_2, min_freq, symmetry_fixes, **kwargs)
334 hac_ratio = num_heavy_atoms / mol_hac
336 if (min_size <= num_heavy_atoms <= max_size and min_rel_size <= hac_ratio <= max_rel_size) \
337 or (replace_cycles and cycle_pattern.search(core)):
--> 339 frag_sma = combine_core_env_to_rxn_smarts(core, env)
341 min_atoms = num_heavy_atoms + min_inc
342 max_atoms = num_heavy_atoms + max_inc
File ~/workspace/sci-dev/sci-dev-venv/lib/python3.8/site-packages/crem/mol_context.py:349, in combine_core_env_to_rxn_smarts(core, env, keep_h)
346 links[i].append(a.GetNeighbors()[0].GetIdx())
347 att_to_remove.append(a.GetIdx())
--> 349 for i, j in links.values():
350 m.AddBond(i, j, Chem.BondType.SINGLE)
352 for i in sorted(att_to_remove, reverse=True):
ValueError: too many values to unpack (expected 2)
However, when I try to replicate the problem later on in a notebook, I get a different error:
list(grow_mol(fragment, db_name='replacements02_sa2.db', max_replacements=9, return_mol=True))
>>>
RuntimeError Traceback (most recent call last)
crem.ipynb Cell 4' in <cell line: 10>()
---> 10 list(grow_mol(fragment, db_name='replacements02_sa2.db', max_replacements=9, protected_ids=[9,10], return_mol=True))
File ~/workspace/sci-dev/sci-dev-venv/lib/python3.8/site-packages/crem/crem.py:487, in mutate_mol(mol, db_name, radius, min_size, max_size, min_rel_size, max_rel_size, min_inc, max_inc, max_replacements, replace_cycles, replace_ids, protected_ids, symmetry_fixes, min_freq, return_rxn, return_rxn_freq, return_mol, ncores, **kwargs)
483 protected_ids = sorted(protected_ids)
485 if ncores == 1:
--> 487 for frag_sma, core_sma, freq, ids in __gen_replacements(mol1=mol, mol2=None, db_name=db_name, radius=radius,
488 min_size=min_size, max_size=max_size,
489 min_rel_size=min_rel_size, max_rel_size=max_rel_size,
490 min_inc=min_inc, max_inc=max_inc,
491 max_replacements=max_replacements,
492 replace_cycles=replace_cycles,
493 protected_ids_1=protected_ids, protected_ids_2=None,
494 min_freq=min_freq, symmetry_fixes=symmetry_fixes,
495 **kwargs):
496 for smi, m, rxn in __frag_replace(mol, None, frag_sma, core_sma, radius, ids, None):
497 if max_replacements is None or (max_replacements is not None and len(products) < max_replacements):
File ~/workspace/sci-dev/sci-dev-venv/lib/python3.8/site-packages/crem/crem.py:314, in __gen_replacements(mol1, mol2, db_name, radius, dist, min_size, max_size, min_rel_size, max_rel_size, min_inc, max_inc, max_replacements, replace_cycles, protected_ids_1, protected_ids_2, min_freq, symmetry_fixes, **kwargs)
312 else:
313 mol = mol1
--> 314 f = __fragment_mol(mol, radius, protected_ids=protected_ids_1, symmetry_fixes=symmetry_fixes)
316 if f:
317 mol_hac = mol.GetNumHeavyAtoms()
File ~/workspace/sci-dev/sci-dev-venv/lib/python3.8/site-packages/crem/crem.py:95, in __fragment_mol(mol, radius, return_ids, keep_stereo, protected_ids, symmetry_fixes)
92 atom.SetIntProp("Index", atom.GetIdx())
94 # heavy atoms
---> 95 frags = rdMMPA.FragmentMol(mol, pattern="[!#1]!@!=!#[!#1]", maxCuts=4, resultsAsMols=True, maxCutBonds=30)
96 frags += rdMMPA.FragmentMol(mol, pattern="[!#1]!@!=!#[!#1]", maxCuts=3, resultsAsMols=True, maxCutBonds=30)
97 # hydrogen atoms
RuntimeError: Pre-condition Violation
RingInfo not initialized
Violation occurred on line 83 in file Code/GraphMol/RingInfo.cpp
Failed Expression: df_init
RDKIT: 2021.09.5
BOOST: 1_67
This expression works fine when I use a molecule without a dummy atom, for example the core molecule defined above. Any help with this would be greatly appreciated.
Thank you!
Perhaps a novice question, but opening the pre-defined fragment databases given here: http://www.qsar4u.com/pages/crem.php.
using gzip, etc. gives an unexpected end of file error. Any advice on properly unarchiving these?
per the python 3.11 documentation:
The population parameter of random.sample() must be a sequence, and automatic conversion of sets to lists is no longer supported. Also, if the sample size is larger than the population size, a ValueError is raised. (Contributed by Raymond Hettinger in bpo-40465.)
This causes crem to crash when running, as it attempts to use random.sample
on a dict.keys object, which apparently is now a non-sequence type (and no longer gets automatically converted).
This is fixed by first converting to a list (or other sequence datatype) prior to calling random.sample
.
Hello, I am getting a strange error with this molecule.
mol = Chem.MolFromSmiles('C[C@H](C(:[O]):[O-])[C@H](c1ccc2c(c1)O[C@@H](C1CCN([C@@H](C)c3cc(F)ccc3OC(F)(F)F)CC1)CC2)C1CC1')
if mol: print(True)
a =list( mutate_mol(mol, db_name='replacements02_sc2.5.db', radius=3, min_size=0, max_size=2, min_inc=-5, max_inc=6))
The molecule looks like a correct RDKit mol. But i got the following error:
[16:29:26] non-ring atom 2 marked aromatic
---------------------------------------------------------------------------
AtomKekulizeException Traceback (most recent call last)
Input In [14], in <cell line: 1>()
----> 1 a =list( mutate_mol(mol, db_name='/home/ale/GITLAB/bi_crem_database/replacements02_sc2.5.db', radius=3, min_size=0, max_size=2, min_inc=-5, max_inc=6))
File ~/anaconda3/envs/moldrug/lib/python3.9/site-packages/crem/crem.py:487, in mutate_mol(mol, db_name, radius, min_size, max_size, min_rel_size, max_rel_size, min_inc, max_inc, max_replacements, replace_cycles, replace_ids, protected_ids, symmetry_fixes, min_freq, return_rxn, return_rxn_freq, return_mol, ncores, **kwargs)
483 protected_ids = sorted(protected_ids)
485 if ncores == 1:
--> 487 for frag_sma, core_sma, freq, ids in __gen_replacements(mol1=mol, mol2=None, db_name=db_name, radius=radius,
488 min_size=min_size, max_size=max_size,
489 min_rel_size=min_rel_size, max_rel_size=max_rel_size,
490 min_inc=min_inc, max_inc=max_inc,
491 max_replacements=max_replacements,
492 replace_cycles=replace_cycles,
493 protected_ids_1=protected_ids, protected_ids_2=None,
494 min_freq=min_freq, symmetry_fixes=symmetry_fixes,
495 **kwargs):
496 for smi, m, rxn in __frag_replace(mol, None, frag_sma, core_sma, radius, ids, None):
497 if max_replacements is None or len(products) < (max_replacements + 1): # +1 because we added source mol to output smiles
File ~/anaconda3/envs/moldrug/lib/python3.9/site-packages/crem/crem.py:314, in __gen_replacements(mol1, mol2, db_name, radius, dist, min_size, max_size, min_rel_size, max_rel_size, min_inc, max_inc, max_replacements, replace_cycles, protected_ids_1, protected_ids_2, min_freq, symmetry_fixes, **kwargs)
312 else:
313 mol = mol1
--> 314 f = __fragment_mol(mol, radius, protected_ids=protected_ids_1, symmetry_fixes=symmetry_fixes)
316 if f:
317 mol_hac = mol.GetNumHeavyAtoms()
File ~/anaconda3/envs/moldrug/lib/python3.9/site-packages/crem/crem.py:113, in __fragment_mol(mol, radius, return_ids, keep_stereo, protected_ids, symmetry_fixes)
110 output.add((env, frag, ids_0))
111 else: # multiple cuts
112 # there are no checks for H needed because H can be present only in single cuts
--> 113 env, frag = get_canon_context_core(chains, core, radius, keep_stereo)
114 output.add((env, frag, get_atom_prop(core) if return_ids else tuple()))
116 if symmetry_fixes:
File ~/anaconda3/envs/moldrug/lib/python3.9/site-packages/crem/mol_context.py:302, in get_canon_context_core(context, core, radius, keep_stereo)
299 def get_canon_context_core(context, core, radius, keep_stereo=False):
300 # context and core are Mols or SMILES
301 # returns SMILES by default
--> 302 res = get_std_context_core_permutations(context, core, radius, keep_stereo)
303 if res:
304 env, cores = res
File ~/anaconda3/envs/moldrug/lib/python3.9/site-packages/crem/mol_context.py:247, in get_std_context_core_permutations(context, core, radius, keep_stereo)
245 # remove Hs from context and core
246 if context: # context cannot be H (no check needed), if so the user will obtain meaningless output
--> 247 context = Chem.RemoveHs(context)
248 if core and Chem.MolToSmiles(core) != '[H][*:1]':
249 core = Chem.RemoveHs(core)
AtomKekulizeException: non-ring atom 2 marked aromatic
Any idea what is going on? A lot of thanks in advance!
If I draw the following molecule in your web app:
And try to use the "Scaffold" (e.g. grow) functionality at the selected index, I get a number of results (maybe a few hundred/thousand) at radius = 1.
However, if I iterate over the atoms in the molecule in a local version, crem reports it cant generate any results at the same position:
smiles = 'c1cc(OC)ccc1C'
base_mol = Chem.MolFromSmiles(smiles)
atoms = len(list(base_mol.GetAtoms()))
for atom in range(0, atoms):
gen_results = list(grow_mol(base_mol, return_mol=True, max_replacements=20000,
db_name='replacements02_sc2.5.db', radius=1, replace_ids=[atom], max_atoms=15, ncores=4))
num_results = len(gen_results)
print(atom, num_results)
The output of this is
[10:36:04] WARNING: not removing hydrogen atom with dummy atom neighbors
0 19999
[10:36:06] WARNING: not removing hydrogen atom with dummy atom neighbors
1 20000
2 0
3 0
[10:36:10] WARNING: not removing hydrogen atom with dummy atom neighbors
4 19999
5 0
6 0
7 0
[10:36:15] WARNING: not removing hydrogen atom with dummy atom neighbors
8 19995
[10:36:18] WARNING: not removing hydrogen atom with dummy atom neighbors
9 19999
[10:36:21] WARNING: not removing hydrogen atom with dummy atom neighbors
10 19998
[10:36:25] WARNING: not removing hydrogen atom with dummy atom neighbors
11 19999
12 0
13 0
14 0
15 0
[10:36:29] WARNING: not removing hydrogen atom with dummy atom neighbors
16 19994
17 0
18 0
And just as a sanity check, this is the numbering of the molecule:
I'm not sure why index 1 is able to produce results, but index 5 doesn't (even though it should by symmetry), and the same position (index 5) works in the web app. Do you know why these would produce different results?
Dear Pavel,
I'm so grateful for your powerful work. when I tried to generate the a custom fragment database as you described.
but after runing the py file fragmentation -i input.smi -o frags.txt -c 32 -v
, the error happened and indicated that
<Can't generate mol for: SMILES
import sys; print('Python %s on %s' % (sys.version, sys.platform))>
could you please give some suggestions how to fix it and what errors?
neigher did I use the wrong format input.smi file? could you please show me some examples how to use the right smi file?
many many thanks,
best,
Sh-Y
Fragmentation of certain molecules fails, it is because ',' as a delimiter doesn't always work . A fix might be to use '\t' as a delimiter as it doesn't appear in SMARTS patterns.
Error in process_line() of frag_to_env_mp.py:
original line CC[C@H]1C[C@@h]1C(=O)NC(C)C,|&1:2,4|,C1C@H[C@H]1[:2],CC(C)NC(=O)[:1].CC[:2]
too many values to unpack (expected 4) ['CC[C@H]1C[C@@h]1C(=O)NC(C)C', '|&1:2', '4|', 'C1C@H[C@H]1[:2]', 'CC(C)NC(=O)[:1].CC[:2]']
There are 5 unpacked values and not 4 as expected because '|&1:2,4|' has a ',' in the pattern
Hi, excellent package! I just wondering how to stop the warning message that pops up during link_mol and mutate_mol. I tried some
import warnings
warnings.filterwarnings("ignore", message='not removing hydrogen atom with dummy atom neighbors')
without succes.
Hi all
I was trying the script with guacamol test and got the warning for GPU from tensorflow (I don't have GPU set up), but what I had hard time to understand is this:
(my-rdkit-env) J:\esp\Personal\AndreaZaliani\pubchem_yojana\crem>python crem/guacamol_crem_test.py
2022-07-04 11:22:30.488812: W tensorflow/stream_executor/platform/default/dso_loader.cc:64] Could not load dynamic library 'cudart64_110.dll'; dlerror: cudart64_110.dll not found
2022-07-04 11:22:30.492783: I tensorflow/stream_executor/cuda/cudart_stub.cc:29] Ignore above cudart dlerror if you do not have a GPU set up on your machine.
Traceback (most recent call last):
File "J:\esp\Personal\AndreaZaliani\pubchem_yojana\crem\crem\guacamol_crem_test.py", line 30, in
from crem import mutate_mol2
File "J:\esp\Personal\AndreaZaliani\pubchem_yojana\crem\crem\crem.py", line 11, in
from crem.mol_context import get_canon_context_core, combine_core_env_to_rxn_smarts
ModuleNotFoundError: No module named 'crem.mol_context'; 'crem' is not a package
Any hint for me?
Best
andrea
Hi, Crem,
I tried the new function enumerate_compounds
in utils
module.
but there is an error, when I use the method to repeat the example you provied.
73dd185
TerminatedWorkerError: A worker process managed by the executor was unexpectedly terminated. This could be caused by a segmentation fault while calling the function or by an excessive memory usage causing the Operating System to kill the worker.
my PC memory is 32G, so I didnot why this error occurs
could you please provide some suggesions to fix it? many thanks,
Shengyang
how can we implement the function in the new version 0.2.11, auxiliary function get_replacements?
could you please provide some examples how to run the auxiliary function get_replacements? many thanks,
best,
Hi! I was wondering whether you might have some more jupyter tutorials that aren't uploaded on the current GitHub repository? Thanks for building this software, it's super cool!
Hi,
I am trying to build my own db according to the instructions, but when I try mols = list(mutate_mol(m, db_name='test.db', max_size=3))
in the tutorial with my own db, the result list of mols
is always empty, but when I try with the preproduced db replacements02_sa2.db
there is no problem; I used the CHEMBL231.smi
file in the example folder, and follow the instructions
fragmentation -i CHEMBL231.smi -o frags.txt -c 32 -v
frag_to_env -i frags.txt -o r3.txt -r 3 -c 32 -v
sort r3.txt | uniq -c > r3_c.txt
env_to_db -i r3_c.txt -o tert.db -r 3 -c -v
got the result
root@872fabd400c3:/home/crem# python test.py
[]
is there something iI missed?
hi, DrrDom,
could you please provide some example for using grow_mol2
and mut_mol2
functions?
In particular, how to set the parameters, such as replace_ids
?
Many many thanks
Hi,
I am using the enumerate_compounds
function to generate new molecules, with 3 positions for potential replacement. Hence, I set the n_iterations = 3
. As far as I'm aware, this does a stepwise growing of my base molecule mol
.
However, I am only interested in the generated molecules where at least 1 specified position is grown. Is there any way to filter out molecules which don't meet this criteria, or how do you suggest I go about this problem?
Thank you so much!
A declarative, efficient, and flexible JavaScript library for building user interfaces.
๐ Vue.js is a progressive, incrementally-adoptable JavaScript framework for building UI on the web.
TypeScript is a superset of JavaScript that compiles to clean JavaScript output.
An Open Source Machine Learning Framework for Everyone
The Web framework for perfectionists with deadlines.
A PHP framework for web artisans
Bring data to life with SVG, Canvas and HTML. ๐๐๐
JavaScript (JS) is a lightweight interpreted programming language with first-class functions.
Some thing interesting about web. New door for the world.
A server is a program made to process requests and deliver data to clients.
Machine learning is a way of modeling and interpreting data that allows a piece of software to respond intelligently.
Some thing interesting about visualization, use data art
Some thing interesting about game, make everyone happy.
We are working to build community through open source technology. NB: members must have two-factor auth.
Open source projects and samples from Microsoft.
Google โค๏ธ Open Source for everyone.
Alibaba Open Source for everyone
Data-Driven Documents codes.
China tencent open source team.