December 31, 2011

Chemdoodle web components tricks #3: Using mouseover functions to do multiple repaints on a canvas

Another step on the way to integrating ChemDoodle web components into SMARTCyp.

Here I show how you can use mouseover, and mouseout functions to redraw the canvas over and over again, and infinite number of times. The trick is to minimize the mouseover and mouseout functions, and use the drawChildExtras function to do all the drawings. Also, the first four lines and the last line in the drawChildExtras function are required to make this work (thanks to Kevin Theisen for helping me with those).

Here is the canvas, and below is the annotated code.



And here is the code:

 var tutorial3_testmol = new ChemDoodle.ViewerCanvas('tutorial3_testmol', 300, 300);
 tutorial3_testmol.specs.atoms_useJMOLColors = true;
 var caffeineMolFile = 'Molecule Name\n  CHEMDOOD08070920033D 0   0.00000     0.00000     0\n[Insert Comment Here]\n 14 15  0  0  0  0  0  0  0  0  1 V2000\n   -0.3318    2.0000    0.0000   O 0  0  0  1  0  0  0  0  0  0  0  0\n   -0.3318    1.0000    0.0000   C 0  0  0  1  0  0  0  0  0  0  0  0\n   -1.1980    0.5000    0.0000   N 0  0  0  1  0  0  0  0  0  0  0  0\n    0.5342    0.5000    0.0000   C 0  0  0  1  0  0  0  0  0  0  0  0\n   -1.1980   -0.5000    0.0000   C 0  0  0  1  0  0  0  0  0  0  0  0\n   -2.0640    1.0000    0.0000   C 0  0  0  4  0  0  0  0  0  0  0  0\n    1.4804    0.8047    0.0000   N 0  0  0  1  0  0  0  0  0  0  0  0\n    0.5342   -0.5000    0.0000   C 0  0  0  1  0  0  0  0  0  0  0  0\n   -2.0640   -1.0000    0.0000   O 0  0  0  1  0  0  0  0  0  0  0  0\n   -0.3318   -1.0000    0.0000   N 0  0  0  1  0  0  0  0  0  0  0  0\n    2.0640   -0.0000    0.0000   C 0  0  0  2  0  0  0  0  0  0  0  0\n    1.7910    1.7553    0.0000   C 0  0  0  4  0  0  0  0  0  0  0  0\n    1.4804   -0.8047    0.0000   N 0  0  0  1  0  0  0  0  0  0  0  0\n   -0.3318   -2.0000    0.0000   C 0  0  0  4  0  0  0  0  0  0  0  0\n  1  2  2  0  0  0  0\n  3  2  1  0  0  0  0\n  4  2  1  0  0  0  0\n  3  5  1  0  0  0  0\n  3  6  1  0  0  0  0\n  7  4  1  0  0  0  0\n  4  8  2  0  0  0  0\n  9  5  2  0  0  0  0\n 10  5  1  0  0  0  0\n 10  8  1  0  0  0  0\n  7 11  1  0  0  0  0\n  7 12  1  0  0  0  0\n 13  8  1  0  0  0  0\n 13 11  2  0  0  0  0\n 10 14  1  0  0  0  0\nM  END\n> \n07-08-2009\n';
 var caffeine = ChemDoodle.readMOL(caffeineMolFile);
 // get the dimension of the molecule
 var size = caffeine.getDimension();
 // find the scale by taking the minimum of the canvas/size ratios
 var scale = Math.min(tutorial3_testmol.width/size.x, tutorial3_testmol.height/size.y);
 // load the molecule first (this function automatically sets scale, so we need to change specs after)
 tutorial3_testmol.loadMolecule(caffeine);
 // change the specs.scale value to the scale calculated, shrinking it slightly so that text is not cut off
 tutorial3_testmol.specs.scale = scale*.8;
 //we are not showing atom numbers right away
 var showatomnumbers = false;
 
 tutorial3_testmol.mouseover = function(){
  //set atom numbers to be displayed
  showatomnumbers = true;
  this.repaint();
 }
 
 tutorial3_testmol.mouseout = function(){
  //set atom numbers to be hidden
  showatomnumbers = false;
  this.repaint();
 } 
 
 tutorial3_testmol.drawChildExtras = function(ctx){
  //four lines to make sure atom numbering, coordinates and scaling work correctly
  ctx.save();
  ctx.translate(this.width/2, this.height/2);
  ctx.rotate(this.specs.rotateAngle);
  ctx.scale(this.specs.scale, this.specs.scale);
  ctx.translate(-this.width/2, -this.height/2);
  //set the font size relative to the scaling of the molecule
  ctx.font = 'bold ' + 2.5*scale + 'px sans-serif';
  //center the atom numbers over the x,y coordinates of the atoms
  ctx.textAlign = 'center';
  ctx.textBaseline = 'middle';
  //iterate through all atoms
  for (var i = 0, ii=caffeine.atoms.length; i<ii; i++) {
   var atom = caffeine.atoms[i];
   if(showatomnumbers){
      //draw a white circle behind the atom number
      ctx.fillStyle = "white";
      ctx.beginPath();
      ctx.arc(atom.x, atom.y, scale*1.8, 0, Math.PI*2, true);
      ctx.fill();
      //draw the atom number in the color
      ctx.fillStyle = ChemDoodle.ELEMENT[atom.label].jmolColor;
      //the "+scale*0.7" is to adjust the text on the y-center of the atoms.
      ctx.fillText(i + 1, atom.x, atom.y);
   }
  } 
  //restore the ctx settings so that application next time will start from scratch
  ctx.restore();
 }
 tutorial3_testmol.repaint();

December 26, 2011

Chemdoodle web components tricks #2: showing atom numbers

Due to the release of ChemDoodle web components version 4.5 showing atom numbers have been much simplified, and this post has been superseded by a new script in a new post Chemdoodle web components tricks #4: showing atom numbers with the altLabel property

My second try on using the chemdoodle web compontents. Very nice GPL licensed javascript library for chemical 2D (and 3D if you have a webGL supporting browser).

One thing missing from the API is atom numbers (i.e. numbering from 1,2...10 for all atoms). Here I'll show how to display atom numbers using a viewercanvas.

Standard canvas with a caffeine molecule and colored atom labels (from my last post on scaling):

And here, the same canvas with atom numbers shown:

The javascript code for the standard scaled canvas:

 
  var tutorial2_testmol = new ChemDoodle.ViewerCanvas('tutorial2_testmol', 300, 300);
  tutorial2_testmol.specs.atoms_useJMOLColors = true;
  var caffeineMolFile = 'Molecule Name\n  CHEMDOOD08070920033D 0   0.00000     0.00000     0\n[Insert Comment Here]\n 14 15  0  0  0  0  0  0  0  0  1 V2000\n   -0.3318    2.0000    0.0000   O 0  0  0  1  0  0  0  0  0  0  0  0\n   -0.3318    1.0000    0.0000   C 0  0  0  1  0  0  0  0  0  0  0  0\n   -1.1980    0.5000    0.0000   N 0  0  0  1  0  0  0  0  0  0  0  0\n    0.5342    0.5000    0.0000   C 0  0  0  1  0  0  0  0  0  0  0  0\n   -1.1980   -0.5000    0.0000   C 0  0  0  1  0  0  0  0  0  0  0  0\n   -2.0640    1.0000    0.0000   C 0  0  0  4  0  0  0  0  0  0  0  0\n    1.4804    0.8047    0.0000   N 0  0  0  1  0  0  0  0  0  0  0  0\n    0.5342   -0.5000    0.0000   C 0  0  0  1  0  0  0  0  0  0  0  0\n   -2.0640   -1.0000    0.0000   O 0  0  0  1  0  0  0  0  0  0  0  0\n   -0.3318   -1.0000    0.0000   N 0  0  0  1  0  0  0  0  0  0  0  0\n    2.0640   -0.0000    0.0000   C 0  0  0  2  0  0  0  0  0  0  0  0\n    1.7910    1.7553    0.0000   C 0  0  0  4  0  0  0  0  0  0  0  0\n    1.4804   -0.8047    0.0000   N 0  0  0  1  0  0  0  0  0  0  0  0\n   -0.3318   -2.0000    0.0000   C 0  0  0  4  0  0  0  0  0  0  0  0\n  1  2  2  0  0  0  0\n  3  2  1  0  0  0  0\n  4  2  1  0  0  0  0\n  3  5  1  0  0  0  0\n  3  6  1  0  0  0  0\n  7  4  1  0  0  0  0\n  4  8  2  0  0  0  0\n  9  5  2  0  0  0  0\n 10  5  1  0  0  0  0\n 10  8  1  0  0  0  0\n  7 11  1  0  0  0  0\n  7 12  1  0  0  0  0\n 13  8  1  0  0  0  0\n 13 11  2  0  0  0  0\n 10 14  1  0  0  0  0\nM  END\n> \n07-08-2009\n';
  var caffeine = ChemDoodle.readMOL(caffeineMolFile);
  // get the dimension of the molecule
  var size = caffeine.getDimension();
  // find the scale by taking the minimum of the canvas/size ratios
  var scale = Math.min(tutorial2_testmol.width/size.x, tutorial2_testmol.height/size.y);
  // load the molecule first (this function automatically sets scale, so we need to change specs after)
  tutorial2_testmol.loadMolecule(caffeine);
  // change the specs.scale value to the scale calculated, shrinking it slightly so that text is not cut off
  tutorial2_testmol.specs.scale = scale*.9;
  // repaint the canvas
  tutorial2_testmol.repaint(); 

and the javascript code for the canvas with the atom numbers:

 
  var tutorial2_testmol2 = new ChemDoodle.ViewerCanvas('tutorial2_testmol2', 300, 300);
  tutorial2_testmol2.specs.atoms_useJMOLColors = true;
  var caffeineMolFile = 'Molecule Name\n  CHEMDOOD08070920033D 0   0.00000     0.00000     0\n[Insert Comment Here]\n 14 15  0  0  0  0  0  0  0  0  1 V2000\n   -0.3318    2.0000    0.0000   O 0  0  0  1  0  0  0  0  0  0  0  0\n   -0.3318    1.0000    0.0000   C 0  0  0  1  0  0  0  0  0  0  0  0\n   -1.1980    0.5000    0.0000   N 0  0  0  1  0  0  0  0  0  0  0  0\n    0.5342    0.5000    0.0000   C 0  0  0  1  0  0  0  0  0  0  0  0\n   -1.1980   -0.5000    0.0000   C 0  0  0  1  0  0  0  0  0  0  0  0\n   -2.0640    1.0000    0.0000   C 0  0  0  4  0  0  0  0  0  0  0  0\n    1.4804    0.8047    0.0000   N 0  0  0  1  0  0  0  0  0  0  0  0\n    0.5342   -0.5000    0.0000   C 0  0  0  1  0  0  0  0  0  0  0  0\n   -2.0640   -1.0000    0.0000   O 0  0  0  1  0  0  0  0  0  0  0  0\n   -0.3318   -1.0000    0.0000   N 0  0  0  1  0  0  0  0  0  0  0  0\n    2.0640   -0.0000    0.0000   C 0  0  0  2  0  0  0  0  0  0  0  0\n    1.7910    1.7553    0.0000   C 0  0  0  4  0  0  0  0  0  0  0  0\n    1.4804   -0.8047    0.0000   N 0  0  0  1  0  0  0  0  0  0  0  0\n   -0.3318   -2.0000    0.0000   C 0  0  0  4  0  0  0  0  0  0  0  0\n  1  2  2  0  0  0  0\n  3  2  1  0  0  0  0\n  4  2  1  0  0  0  0\n  3  5  1  0  0  0  0\n  3  6  1  0  0  0  0\n  7  4  1  0  0  0  0\n  4  8  2  0  0  0  0\n  9  5  2  0  0  0  0\n 10  5  1  0  0  0  0\n 10  8  1  0  0  0  0\n  7 11  1  0  0  0  0\n  7 12  1  0  0  0  0\n 13  8  1  0  0  0  0\n 13 11  2  0  0  0  0\n 10 14  1  0  0  0  0\nM  END\n> \n07-08-2009\n';
  var caffeine = ChemDoodle.readMOL(caffeineMolFile);
  // get the dimension of the molecule
  var size = caffeine.getDimension();
  // find the scale by taking the minimum of the canvas/size ratios
  var scale = Math.min(tutorial2_testmol2.width/size.x, tutorial2_testmol2.height/size.y);
  // load the molecule first (this function automatically sets scale, so we need to change specs after)
  tutorial2_testmol2.loadMolecule(caffeine);
  // change the specs.scale value to the scale calculated, shrinking it slightly so that text is not cut off
  tutorial2_testmol2.specs.scale = scale*.9;
  tutorial2_testmol2.repaint();
  var molcenter = caffeine.getCenter();
  
  tutorial2_testmol2.drawChildExtras = function(ctx){
    ctx.translate(this.width/2, this.height/2);
    ctx.rotate(this.specs.rotateAngle);
    ctx.scale(this.specs.scale, this.specs.scale);
    ctx.translate(-this.width/2, -this.height/2);
    //draw atom numbers
    ctx.font = 'bold ' + 2.5*scale + 'px sans-serif';
    ctx.textAlign = 'center';
    ctx.textBaseline = 'middle';
 //iterate through all atoms
    for (var i = 0, ii=caffeine.atoms.length; i<ii; i++) {
      var atom = caffeine.atoms[i];
      //draw a white circle behind the atom number
      ctx.fillStyle = "white";
      ctx.strokeStyle = "white";
      ctx.beginPath();
      ctx.arc(atom['x'],atom['y'],scale*1.8,0,Math.PI*2,true);
      ctx.fill();
      ctx.stroke();
      //draw the atom number
   thisatomcolor = ChemDoodle.ELEMENT[atom['label']].jmolColor;
      ctx.fillStyle = thisatomcolor;
      ctx.fillText(i + 1, atom['x'], atom['y']);
    }
  }

  tutorial2_testmol2.repaint();

The trickiest part about this is that while the first four lines in the drawChildExtras function (as suggested by Kevin in the comments and the fourth line is needed to simplify the atom locations), do scale the drawing on the canvas correctly, they do not scale or translate the molecule. That's why I had to subtract the x and y coordinates of the center of the molecule to get the positions right.

December 23, 2011

Chemdoodle web components tricks #1: scaling molecules

I have played around a bit with chemdoodle web compontents. Very nice GPL licensed javascript library for chemical 2D (and 3D if you have a webGL supporting browser).

However, the standard API still is missing some functionality. But you can fix most things by using standard javascript (at least as long as you're only using a standard viewer canvas, without zoom functionality).

Here I will show how you can scale a molecule to the size of the canvas. The standard viewercanvas will not scale molecules to the size of the canvas. So I added this to my code.

Standard canvas with a caffeine molecule and colored atom labels:

And here, the same canvas with the molecule scaled by the size of the canvas:

The javascript code for the standard canvas:

 
  var testmol = new ChemDoodle.ViewerCanvas('testmol', 300, 300);
  testmol.specs.atoms_useJMOLColors = true;
  var caffeineMolFile = 'Molecule Name\n  CHEMDOOD08070920033D 0   0.00000     0.00000     0\n[Insert Comment Here]\n 14 15  0  0  0  0  0  0  0  0  1 V2000\n   -0.3318    2.0000    0.0000   O 0  0  0  1  0  0  0  0  0  0  0  0\n   -0.3318    1.0000    0.0000   C 0  0  0  1  0  0  0  0  0  0  0  0\n   -1.1980    0.5000    0.0000   N 0  0  0  1  0  0  0  0  0  0  0  0\n    0.5342    0.5000    0.0000   C 0  0  0  1  0  0  0  0  0  0  0  0\n   -1.1980   -0.5000    0.0000   C 0  0  0  1  0  0  0  0  0  0  0  0\n   -2.0640    1.0000    0.0000   C 0  0  0  4  0  0  0  0  0  0  0  0\n    1.4804    0.8047    0.0000   N 0  0  0  1  0  0  0  0  0  0  0  0\n    0.5342   -0.5000    0.0000   C 0  0  0  1  0  0  0  0  0  0  0  0\n   -2.0640   -1.0000    0.0000   O 0  0  0  1  0  0  0  0  0  0  0  0\n   -0.3318   -1.0000    0.0000   N 0  0  0  1  0  0  0  0  0  0  0  0\n    2.0640   -0.0000    0.0000   C 0  0  0  2  0  0  0  0  0  0  0  0\n    1.7910    1.7553    0.0000   C 0  0  0  4  0  0  0  0  0  0  0  0\n    1.4804   -0.8047    0.0000   N 0  0  0  1  0  0  0  0  0  0  0  0\n   -0.3318   -2.0000    0.0000   C 0  0  0  4  0  0  0  0  0  0  0  0\n  1  2  2  0  0  0  0\n  3  2  1  0  0  0  0\n  4  2  1  0  0  0  0\n  3  5  1  0  0  0  0\n  3  6  1  0  0  0  0\n  7  4  1  0  0  0  0\n  4  8  2  0  0  0  0\n  9  5  2  0  0  0  0\n 10  5  1  0  0  0  0\n 10  8  1  0  0  0  0\n  7 11  1  0  0  0  0\n  7 12  1  0  0  0  0\n 13  8  1  0  0  0  0\n 13 11  2  0  0  0  0\n 10 14  1  0  0  0  0\nM  END\n> \n07-08-2009\n';
  var caffeine = ChemDoodle.readMOL(caffeineMolFile);
  testmol.loadMolecule(caffeine);

and the javascript code for the canvas with the scaled molecule:

 
  var testmol2 = new ChemDoodle.ViewerCanvas('testmol2', 300, 300);
  testmol2.specs.atoms_useJMOLColors = true;
  var caffeineMolFile = 'Molecule Name\n  CHEMDOOD08070920033D 0   0.00000     0.00000     0\n[Insert Comment Here]\n 14 15  0  0  0  0  0  0  0  0  1 V2000\n   -0.3318    2.0000    0.0000   O 0  0  0  1  0  0  0  0  0  0  0  0\n   -0.3318    1.0000    0.0000   C 0  0  0  1  0  0  0  0  0  0  0  0\n   -1.1980    0.5000    0.0000   N 0  0  0  1  0  0  0  0  0  0  0  0\n    0.5342    0.5000    0.0000   C 0  0  0  1  0  0  0  0  0  0  0  0\n   -1.1980   -0.5000    0.0000   C 0  0  0  1  0  0  0  0  0  0  0  0\n   -2.0640    1.0000    0.0000   C 0  0  0  4  0  0  0  0  0  0  0  0\n    1.4804    0.8047    0.0000   N 0  0  0  1  0  0  0  0  0  0  0  0\n    0.5342   -0.5000    0.0000   C 0  0  0  1  0  0  0  0  0  0  0  0\n   -2.0640   -1.0000    0.0000   O 0  0  0  1  0  0  0  0  0  0  0  0\n   -0.3318   -1.0000    0.0000   N 0  0  0  1  0  0  0  0  0  0  0  0\n    2.0640   -0.0000    0.0000   C 0  0  0  2  0  0  0  0  0  0  0  0\n    1.7910    1.7553    0.0000   C 0  0  0  4  0  0  0  0  0  0  0  0\n    1.4804   -0.8047    0.0000   N 0  0  0  1  0  0  0  0  0  0  0  0\n   -0.3318   -2.0000    0.0000   C 0  0  0  4  0  0  0  0  0  0  0  0\n  1  2  2  0  0  0  0\n  3  2  1  0  0  0  0\n  4  2  1  0  0  0  0\n  3  5  1  0  0  0  0\n  3  6  1  0  0  0  0\n  7  4  1  0  0  0  0\n  4  8  2  0  0  0  0\n  9  5  2  0  0  0  0\n 10  5  1  0  0  0  0\n 10  8  1  0  0  0  0\n  7 11  1  0  0  0  0\n  7 12  1  0  0  0  0\n 13  8  1  0  0  0  0\n 13 11  2  0  0  0  0\n 10 14  1  0  0  0  0\nM  END\n> \n07-08-2009\n';
  var caffeine = ChemDoodle.readMOL(caffeineMolFile);
  // get the dimension of the molecule
  var size = caffeine.getDimension();
  // find the scale by taking the minimum of the canvas/size ratios
  var scale = Math.min(testmol2.width/size.x, testmol2.height/size.y);
  // load the molecule first (this function automatically sets scale, so we need to change specs after)
  testmol2.loadMolecule(caffeine);
  // change the specs.scale value to the scale calculated, shrinking it slightly so that text is not cut off
  testmol2.specs.scale = scale*.9;
  // repaint the canvas
  testmol2.repaint(); 

The code for the scaled figure updated after the comment below from Kevin Theisen. So now the code is simpler, but requires a repaint() command instead. Thanks to Kevin for this.

December 4, 2011

Circular fingerprints in CDK

I've been working on some code that computes circular fingerprints for all atoms in a molecule with the help of the CDK java library.
The code implements the circular fingerprints from a paper by Xing et al. (J. Chem. Inf. Comput. Sci. 2003, 43, 870-879).

If anyone who has the skills and time would turn this into a proper CDK function which could be implemented into the CDK library that would probably be helpful for other people. Here follows my code with assumes that it's applied in java class where this is an AtomContainer.

 
// Calculates the atom environment descriptor of all atoms
public void calculateAtomEnvironmentDescriptor() throws CDKException{
        //maxLevel is how many bonds from the atom we will count atom types
 int maxLevel = 5;
 IAtomType[] atomTypes = null;
 int natoms = this.getAtomCount();
 //do atom type matching
 IAtomTypeMatcher atm = SybylAtomTypeMatcher.getInstance(NoNotificationChemObjectBuilder.getInstance());
 InputStream ins = this.getClass().getClassLoader().getResourceAsStream("org/openscience/cdk/dict/data/sybyl-atom-types.owl");
 AtomTypeFactory factory = AtomTypeFactory.getInstance(ins,"owl",NoNotificationChemObjectBuilder.getInstance());
 atomTypes = factory.getAllAtomTypes();
 //map atomtypes to the atomIndex integer array
 TreeMap map = new TreeMap();
    for (int i = 0; i < atomTypes.length; i++) { 
        map.put(atomTypes[i].getAtomTypeName(),new Integer(i));
    }
    int[] atomIndex = new int[natoms]; //array of atom type integers
    for (int i = 0; i < natoms; i++) {
        try {
            IAtomType a = atm.findMatchingAtomType(this,this.getAtom(i));
            if ( a != null) {
                Object mappedType = map.get(a.getAtomTypeName());
                if (mappedType != null) 
                    atomIndex[i] = ((Integer) mappedType).intValue();
                else {
                    //System.out.println(a.getAtomTypeName() + " not found in " + map);
                    atomIndex[i] = -1;
                }    
            } else //atom type not found 
             atomIndex[i] = -1;
        } catch (Exception x) {
            x.printStackTrace();
            throw new CDKException(x.getMessage() + "\ninitConnectionMatrix");
        }                
    }
    //compute bond distances between all atoms
    int[][] aMatrix = PathTools.computeFloydAPSP(AdjacencyMatrix.getMatrix(this));
    //assign values to the results arrays for all atoms
 int L = (atomTypes.length +1) ;
 int [][] result = new int[natoms][L*(maxLevel)+2]; //create result array
 for (int i = 0; i < natoms; i++) {
  //for every atom, iterate through its connections to all other atoms
  for (int j=0; j < natoms; j++) {
   if (aMatrix[i][j] == 0) 
    result[i][1] = atomIndex[j]; //atom j is atom i
         else if (aMatrix[i][j] > 0 && aMatrix[i][j] <= maxLevel){ //j is not atom i and bonds less or equal to maxlevel
          if (atomIndex[j] >= 0) //atom type defined in factory
              result[i][L*(aMatrix[i][j]-1)+atomIndex[j]+2]++;
          else if (atomIndex[j] == -1) //-1, unknown  type
           result[i][L*(aMatrix[i][j]-1)+(L-1)+2]++;
         }
  }
  //checksum for easy comparison        
  for (int j = 1; j < result[i].length; j++) result[i][0] += result[i][j];
 }
 if (debug == 1){
        //print out the names of atom types for reference use
     System.out.print("Sum\tAtomType\t");
  for (int j=0; j < factory.getSize(); j++) {
   System.out.print(j);
   System.out.print(".");
      System.out.print(atomTypes[j].getAtomTypeName());
      System.out.print("\t");
     }       
     System.out.println("");
     for (int i = 0; i < natoms; i++)
      System.out.println(Arrays.toString(result[i]));
 }
}
Big thanks to Nina Jeliazkova who gave me a preliminary version of the code from her old Ambit code. I have rewritten her code so now it's bug free and faster.

December 3, 2011

Why SMARTCyp works - challenges in reactivity models applied to P450s

During the last couple of years I've been developing the SMARTCyp methodology (www.farma.ku.dk/smartcyp) for prediction of site-of-metabolism in drug metabolism mediated by the cytochromes P450 enzyme family. Here I will describe why the SMARTCyp approach works, and why many other models applied to the same problem often fails when they try to include reactivity.


  1. SMARTCyp is simple
    The application of the reactivity model is through simple fragment matching towards a library of SMARTS strings. Each SMARTS string describes a number of reference calculation and assigns the average energy of these reference calculations to any atom that matches.
    The advantage of this is that we get rid of the problem of generating many different 3D structures to make sure the conformation is correct, and we get rid of the problem of possible hydrogen bond donors/acceptors messing up reactivity calculations by causing strange structures (this can happen when the oxygen radical attacks a site close to a donor/acceptor).
  2. SMARTCyp reactivities are very accurate
    For the so far 200+ fragments we have computed the reactivities are as good as can possibly be done today (DFT with B3LYP, large basis sets and transition states verified by frequency calculations). As long as the fragment of interest in a drug molecule is relatively similar, the predicted reactivity will be quite good (>1 kcal/mol error is unusual).
    This is because our calculations use a full heme model, and compute the transition state for the actual reaction mechanism, no simplified model or other assumptions are made. Other reactivity models use intermediates and also do not use a full heme model. Such models are relatively accurate for hydrogen abstraction reactions (that is hydroxylation of aliphatic carbon atoms or dealkylation reactions), but for any other type of reactions the models are usually pretty bad.
  3. P450 enzymes are flexible
    While it has been shown that CYP3A4 is highly flexible and reactivity is the only major determinant of site-of-metabolism in this enzyme, it's my belief that most other CYPs also are more flexible than most people believe. The crystal structures of several of the human CYPs do not have active sites large enough to explain all the substrates known to be metabolized by them. Hence, the approximation that reactivity is very important is very good for most CYPs (however, it's not sufficient for CYP2D6 and CYP2C9). Unpublished data show that SMARTCyp can be applied to all isoforms except 2D6 and 2C9 with high accuracy. The two exceptions stem from the fact that these are the only two isoforms with charged amino acids in the active site, leading to large contributions from binding. But even these binding contributions can be described by a simple 2D model.
    We recently showed that even for CYP2D6 we can build a good model based on SMARTCyp reactivities using only three descriptors (implemented in SMARTCyp 2.0). While some may claim such a model is too simple, the simplicity is actually an advantage. Simple models are usually robust, the dependence on data sets can be smaller if you have fewer descriptors that you're trying to fit. And hence we have less noise in the data. When compared to ensemble docking supplemented with reactivities from SMARTCyp 1.5, the new SMARTCyp 2D6 model is as accurate,but with a fraction of the computational cost, and much much easier to apply for a medicinal chemist.