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

Module Pairs

source code

Contains an implementation of Atom-pair fingerprints, as
described in:

R.E. Carhart, D.H. Smith, R. Venkataraghavan;
"Atom Pairs as Molecular Features in Structure-Activity Studies:
Definition and Applications" JCICS 25, 64-73 (1985).

Functions [hide private]
 
pyScorePair(at1, at2, dist, atomCodes=None)
Returns a score for an individual atom pair.
source code
 
ExplainPairScore(score)
>>> m = Chem.MolFromSmiles('C=CC')...
source code
 
GetAtomPairFingerprintAsBitVect(mol)
Returns the Atom-pair fingerprint for a molecule as a SparseBitVect.
source code
 
_test() source code
Variables [hide private]
  numPathBits = 5
  _maxPathLen = 31
  numFpBits = 23
  fpLen = 8388608
  __package__ = 'rdkit.Chem.AtomPairs'

Imports: IntSparseIntVect, Chem, rdMolDescriptors, Utils, DataStructs, GetAtomPairFingerprint, GetHashedAtomPairFingerprint, GetAtomPairFingerprintAsIntVect


Function Details [hide private]

pyScorePair(at1, at2, dist, atomCodes=None)

source code 
Returns a score for an individual atom pair.

>>> m = Chem.MolFromSmiles('CCCCC')
>>> c1 = Utils.GetAtomCode(m.GetAtomWithIdx(0))
>>> c2 = Utils.GetAtomCode(m.GetAtomWithIdx(1))
>>> c3 = Utils.GetAtomCode(m.GetAtomWithIdx(2))
>>> t = 1 | min(c1,c2)<<numPathBits | max(c1,c2)<<(rdMolDescriptors.AtomPairsParameters.codeSize+numPathBits)
>>> pyScorePair(m.GetAtomWithIdx(0),m.GetAtomWithIdx(1),1)==t
1
>>> pyScorePair(m.GetAtomWithIdx(1),m.GetAtomWithIdx(0),1)==t
1
>>> t = 2 | min(c1,c3)<<numPathBits | max(c1,c3)<<(rdMolDescriptors.AtomPairsParameters.codeSize+numPathBits)
>>> pyScorePair(m.GetAtomWithIdx(0),m.GetAtomWithIdx(2),2)==t
1
>>> pyScorePair(m.GetAtomWithIdx(0),m.GetAtomWithIdx(2),2,
...  atomCodes=(Utils.GetAtomCode(m.GetAtomWithIdx(0)),Utils.GetAtomCode(m.GetAtomWithIdx(2))))==t
1

ExplainPairScore(score)

source code 

>>> m = Chem.MolFromSmiles('C=CC')
>>> score = pyScorePair(m.GetAtomWithIdx(0),m.GetAtomWithIdx(1),1)
>>> ExplainPairScore(score)
(('C', 1, 1), 1, ('C', 2, 1))
>>> score = pyScorePair(m.GetAtomWithIdx(0),m.GetAtomWithIdx(2),2)
>>> ExplainPairScore(score)
(('C', 1, 0), 2, ('C', 1, 1))
>>> score = pyScorePair(m.GetAtomWithIdx(1),m.GetAtomWithIdx(2),1)
>>> ExplainPairScore(score)
(('C', 1, 0), 1, ('C', 2, 1))
>>> score = pyScorePair(m.GetAtomWithIdx(2),m.GetAtomWithIdx(1),1)
>>> ExplainPairScore(score)
(('C', 1, 0), 1, ('C', 2, 1))

GetAtomPairFingerprintAsBitVect(mol)

source code 
Returns the Atom-pair fingerprint for a molecule as
a SparseBitVect. Note that this doesn't match the standard
definition of atom pairs, which uses counts of the
pairs, not just their presence.

**Arguments**:

  - mol: a molecule

**Returns**: a SparseBitVect

>>> m = Chem.MolFromSmiles('CCC')
>>> v = [ pyScorePair(m.GetAtomWithIdx(0),m.GetAtomWithIdx(1),1),
...       pyScorePair(m.GetAtomWithIdx(0),m.GetAtomWithIdx(2),2),
...     ]
>>> v.sort()
>>> fp = GetAtomPairFingerprintAsBitVect(m)
>>> list(fp.GetOnBits())==v
True