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

Module QED

source code



QED stands for quantitative estimation of drug-likeness and the concept was for the first time
introduced by Richard Bickerton and coworkers [1]. The empirical rationale of the QED measure
reflects the underlying distribution of molecular properties including molecular weight, logP,
topological polar surface area, number of hydrogen bond donors and acceptors, the number of
aromatic rings and rotatable bonds, and the presence of unwanted chemical functionalities.

The QED results as generated by the RDKit-based implementation of Biscu-it(tm) are not completely
identical to those from the original publication [1]. These differences are a consequence of
differences within the underlying calculated property calculators used in both methods. For
example, discrepancies can be noted in the results from the logP calculations, nevertheless
despite the fact that both approaches (Pipeline Pilot in the original publication and RDKit
in our Biscu-it(tm) implementation) mention to use the Wildmann and Crippen methodology for the
calculation of their logP-values [2]. However, the differences in the resulting QED-values
are very small and are not compromising the usefulness of using Qed in your daily research.

[1] Bickerton, G.R.; Paolini, G.V.; Besnard, J.; Muresan, S.; Hopkins, A.L. (2012)
    'Quantifying the chemical beauty of drugs',
    Nature Chemistry, 4, 90-98 [http://dx.doi.org/10.1038/nchem.1243]

History:
  2012-04 Adapted to internal RDkit implementation
  2013-05 moved to rdkit.Chem.QED

Classes [hide private]
  QEDproperties
QEDproperties(MW, ALOGP, HBA, HBD, PSA, ROTB, AROM, ALERTS)
  ADSparameter
ADSparameter(A, B, C, D, E, F, DMAX)
Functions [hide private]
 
ads(x, adsParameter)
ADS function
source code
 
properties(mol)
Calculates the properties that are required to calculate the QED descriptor.
source code
 
qed(mol, w=QEDproperties(MW=0.66, ALOGP=0.46, HBA=0.05, HBD=0.61, PSA=0.0..., qedProperties=None)
Calculate the weighted sum of ADS mapped properties
source code
 
weights_max(mol)
Calculates the QED descriptor using maximal descriptor weights.
source code
 
weights_mean(mol)
Calculates the QED descriptor using average descriptor weights.
source code
 
weights_none(mol)
Calculates the QED descriptor using unit weights.
source code
 
default(mol)
Calculates the QED descriptor using average descriptor weights.
source code
 
_runDoctests(verbose=None) source code
Variables [hide private]
  WEIGHT_MAX = QEDproperties(MW=0.5, ALOGP=0.25, HBA=0.0, HBD=0....
  WEIGHT_MEAN = QEDproperties(MW=0.66, ALOGP=0.46, HBA=0.05, HBD...
  WEIGHT_NONE = QEDproperties(MW=1.0, ALOGP=1.0, HBA=1.0, HBD=1....
  AliphaticRings = Chem.MolFromSmarts('[$([A;R][!a])]')
  AcceptorSmarts = ['[oH0;X2]', '[OH1;X2;v2]', '[OH0;X2;v2]', '[...
  Acceptors = [Chem.MolFromSmarts(hba) for hba in AcceptorSmarts]
  StructuralAlertSmarts = ['*1[O,S,N]*1', '[S,C](=[O,S])[F,Br,Cl...
  StructuralAlerts = [Chem.MolFromSmarts(smarts) for smarts in S...
  adsParameters = {'ALERTS': ADSparameter(A=0.01, B=1199.094025,...
  __package__ = 'rdkit.Chem'
  hba = '[$([N;+0;X3;v3]);!$(N[C,S]=O)]'
  smarts = '[34S]'

Imports: namedtuple, math, Chem, MolSurf, Crippen, rdmd, setDescriptorVersion


Function Details [hide private]

qed(mol, w=QEDproperties(MW=0.66, ALOGP=0.46, HBA=0.05, HBD=0.61, PSA=0.0..., qedProperties=None)

source code 
Calculate the weighted sum of ADS mapped properties

some examples from the QED paper, reference values from Peter G's original implementation
>>> m = Chem.MolFromSmiles('N=C(CCSCc1csc(N=C(N)N)n1)NS(N)(=O)=O')
>>> qed(m)
0.253...
>>> m = Chem.MolFromSmiles('CNC(=NCCSCc1nc[nH]c1C)NC#N')
>>> qed(m)
0.234...
>>> m = Chem.MolFromSmiles('CCCCCNC(=N)NN=Cc1c[nH]c2ccc(CO)cc12')
>>> qed(m)
0.234...

Decorators:
  • @setDescriptorVersion(version= '1.1.0')

Variables Details [hide private]

WEIGHT_MAX

Value:
QEDproperties(MW=0.5, ALOGP=0.25, HBA=0.0, HBD=0.5, PSA=0.0, ROTB=0.5,\
 AROM=0.25, ALERTS=1.0)

WEIGHT_MEAN

Value:
QEDproperties(MW=0.66, ALOGP=0.46, HBA=0.05, HBD=0.61, PSA=0.06, ROTB=\
0.65, AROM=0.48, ALERTS=0.95)

WEIGHT_NONE

Value:
QEDproperties(MW=1.0, ALOGP=1.0, HBA=1.0, HBD=1.0, PSA=1.0, ROTB=1.0, \
AROM=1.0, ALERTS=1.0)

AcceptorSmarts

Value:
['[oH0;X2]',
 '[OH1;X2;v2]',
 '[OH0;X2;v2]',
 '[OH0;X1;v2]',
 '[O-;X1]',
 '[SH0;X2;v2]',
 '[SH0;X1;v2]',
 '[S-;X1]',
...

StructuralAlertSmarts

Value:
['*1[O,S,N]*1',
 '[S,C](=[O,S])[F,Br,Cl,I]',
 '[CX4][Cl,Br,I]',
 '[#6]S(=O)(=O)O[#6]',
 '[$([CH]),$(CC)]#CC(=O)[#6]',
 '[$([CH]),$(CC)]#CC(=O)O[#6]',
 'n[OH]',
 '[$([CH]),$(CC)]#CS(=O)(=O)[#6]',
...

StructuralAlerts

Value:
[Chem.MolFromSmarts(smarts) for smarts in StructuralAlertSmarts]

adsParameters

Value:
{'ALERTS': ADSparameter(A=0.01, B=1199.094025, C=-0.09002883, D=1e-09,\
 E=0.185904477, F=0.875193782, DMAX=417.725314),
 'ALOGP': ADSparameter(A=3.172690585, B=137.8624751, C=2.534937431, D=\
4.581497897, E=0.822739154, F=0.576295591, DMAX=131.3186604),
 'AROM': ADSparameter(A=3.21778897, B=957.7374108, C=2.274627939, D=1e\
-09, E=1.317690384, F=0.375760881, DMAX=312.337261),
 'HBA': ADSparameter(A=2.948620388, B=160.4605972, C=3.615294657, D=4.\
435986202, E=0.290141953, F=1.300669958, DMAX=148.7763046),
...