Package rdkit :: Package Chem :: Package Fraggle :: Module FraggleSim
[hide private]
[frames] | no frames]

Module FraggleSim

source code

Functions [hide private]
 
delete_bonds(mol, bonds, ftype, hac) source code
 
is_ring_cut_valid(smi) source code
 
select_fragments(f_smi, ftype, hac) source code
 
generate_fraggle_fragmentation(mol)
>>> q = Chem.MolFromSmiles('COc1cc(CN2CCC(NC(=O)c3cncc(C)c3)CC2)c(OC)c2ccccc12')...
source code
 
atomContrib(subs, mol, tverskyThresh=0.8) source code
 
compute_fraggle_similarity_for_subs(inMol, qMol, qSmi, qSubs, tverskyThresh=0.8) source code
 
GetFraggleSimilarity(queryMol, refMol, tverskyThresh=0.8)
return the Fraggle similarity between two molecules
source code
 
_test() source code
Variables [hide private]
  rdkitFpParams = {'fpSize': 1024, 'maxPath': 5, 'nBitsPerHash': 2}
  acyc_smarts = Chem.MolFromSmarts("[*]!@!=!#[*]")
  cyc_smarts = Chem.MolFromSmarts("[R1,R2]@[r;!R1]")
  cSma1 = Chem.MolFromSmarts("[#0][r].[r][#0]")
  cSma2 = Chem.MolFromSmarts("[#0][r][#0]")
  modified_query_fps = {}
  __package__ = 'rdkit.Chem.Fraggle'

Imports: sys, Chem, DataStructs


Function Details [hide private]

generate_fraggle_fragmentation(mol)

source code 

>>> q = Chem.MolFromSmiles('COc1cc(CN2CCC(NC(=O)c3cncc(C)c3)CC2)c(OC)c2ccccc12')
>>> list(generate_fraggle_fragmentation(q))
 ['[*]C(=O)NC1CCN(Cc2cc(OC)c3ccccc3c2OC)CC1', '[*]C(=O)c1cncc(C)c1.[*]C1CCN(Cc2cc(OC)c3ccccc3c2OC)CC1', '[*]C(=O)c1cncc(C)c1.[*]Cc1cc(OC)c2ccccc2c1OC', '[*]C(=O)c1cncc(C)c1.[*]c1cc(OC)c2ccccc2c1OC', '[*]C1CCN(Cc2cc(OC)c3ccccc3c2OC)CC1', '[*]C1CCN(Cc2cc(OC)c3ccccc3c2OC)CC1.[*]c1cncc(C)c1', '[*]Cc1cc(OC)c2ccccc2c1OC.[*]NC(=O)c1cncc(C)c1', '[*]Cc1cc(OC)c2ccccc2c1OC.[*]c1cncc(C)c1', '[*]N1CCC(NC(=O)c2cncc(C)c2)CC1.[*]c1cc(OC)c2ccccc2c1OC', '[*]NC(=O)c1cncc(C)c1.[*]c1cc(OC)c2ccccc2c1OC', '[*]NC1CCN(Cc2cc(OC)c3ccccc3c2OC)CC1', '[*]NC1CCN(Cc2cc(OC)c3ccccc3c2OC)CC1.[*]c1cncc(C)c1', '[*]c1c(CN2CCC(NC(=O)c3cncc(C)c3)CC2)cc(OC)c2ccccc12', '[*]c1c(OC)cc(CN2CCC(NC(=O)c3cncc(C)c3)CC2)c(OC)c1[*]', '[*]c1cc(CN2CCC(NC(=O)c3cncc(C)c3)CC2)c(OC)c2ccccc12', '[*]c1cc(OC)c2ccccc2c1OC.[*]c1cncc(C)c1']

GetFraggleSimilarity(queryMol, refMol, tverskyThresh=0.8)

source code 
return the Fraggle similarity between two molecules

>>> q = Chem.MolFromSmiles('COc1cc(CN2CCC(NC(=O)c3cncc(C)c3)CC2)c(OC)c2ccccc12')
>>> m = Chem.MolFromSmiles('COc1cc(CN2CCC(NC(=O)c3ccccc3)CC2)c(OC)c2ccccc12')
>>> sim,match = GetFraggleSimilarity(q,m)
>>> sim
0.980...
>>> match
'[*]C1CCN(Cc2cc(OC)c3ccccc3c2OC)CC1'

>>> m = Chem.MolFromSmiles('COc1cc(CN2CCC(Nc3nc4ccccc4s3)CC2)c(OC)c2ccccc12')
>>> sim,match = GetFraggleSimilarity(q,m)
>>> sim
0.794...
>>> match
'[*]C1CCN(Cc2cc(OC)c3ccccc3c2OC)CC1'

>>> q = Chem.MolFromSmiles('COc1ccccc1')
>>> sim,match = GetFraggleSimilarity(q,m)
>>> sim
0.347...
>>> match
'[*]c1ccccc1'