com.actelion.research.chem.docking.DockingEngine Maven / Gradle / Ivy
Go to download
Show more of this group Show more artifacts with this name
Show all versions of openchemlib Show documentation
Show all versions of openchemlib Show documentation
Open Source Chemistry Library
package com.actelion.research.chem.docking;
import java.util.ArrayList;
import java.util.HashMap;
import java.util.HashSet;
import java.util.List;
import java.util.Map;
import java.util.Random;
import java.util.Set;
import java.util.TreeSet;
import org.openmolecules.chem.conf.gen.ConformerGenerator;
import com.actelion.research.calc.Matrix;
import com.actelion.research.chem.Canonizer;
import com.actelion.research.chem.Coordinates;
import com.actelion.research.chem.Molecule;
import com.actelion.research.chem.Molecule3D;
import com.actelion.research.chem.StereoMolecule;
import com.actelion.research.chem.conf.Conformer;
import com.actelion.research.chem.conf.ConformerSet;
import com.actelion.research.chem.conf.ConformerSetGenerator;
import com.actelion.research.chem.docking.scoring.AbstractScoringEngine;
import com.actelion.research.chem.docking.scoring.ChemPLP;
import com.actelion.research.chem.docking.scoring.IdoScore;
import com.actelion.research.chem.forcefield.mmff.ForceFieldMMFF94;
import com.actelion.research.chem.forcefield.mmff.PositionConstraint;
import com.actelion.research.chem.interactionstatistics.InteractionAtomTypeCalculator;
import com.actelion.research.chem.io.pdb.converter.MoleculeGrid;
import com.actelion.research.chem.optimization.OptimizerLBFGS;
import com.actelion.research.chem.phesa.DescriptorHandlerShapeOneConf;
import com.actelion.research.chem.phesa.MolecularVolume;
import com.actelion.research.chem.phesa.PheSAAlignment;
import com.actelion.research.chem.phesa.PheSAAlignment.PheSAResult;
import com.actelion.research.chem.phesa.PheSAMolecule;
public class DockingEngine {
public enum ScoringFunction {CHEMPLP,IDOSCORE;}
public enum StartPosition {PHESA, RANDOM;}
private static final int DEFAULT_NR_MC_STEPS = 50;
private static final int DEFAULT_START_POSITIONS = 15;
private static final double BOLTZMANN_FACTOR = 1.2; //as for AutoDock Vina
public static final double GRID_DIMENSION = 6.0;
public static final double GRID_RESOLUTION = 0.5;
public static final double MINI_CUTOFF = 100; //if energy higher after MC step, don't minimize
private Matrix rotation; //for initial prealignment to native ligand
private Coordinates origCOM;
private int mcSteps;
private Random random;
private AbstractScoringEngine engine;
private int startPositions;
private StereoMolecule nativeLigand;
private StartPosition startPosition;
public DockingEngine(StereoMolecule rec, StereoMolecule nativeLig, int mcSteps, int startPositions,
ScoringFunction scoringFunction, StartPosition startPosition) {
this.startPosition = startPosition;
nativeLigand = new Molecule3D(nativeLig);
nativeLigand.ensureHelperArrays(Molecule.cHelperCIP);
Molecule3D receptor = new Molecule3D(rec);
receptor.ensureHelperArrays(Molecule.cHelperCIP);
MolecularVolume molVol = new MolecularVolume(nativeLigand);
origCOM = new Coordinates(molVol.getCOM());
Conformer conf = new Conformer(nativeLigand);
rotation = PheSAAlignment.preProcess(conf, molVol);
this.startPositions = startPositions;
preprocess(receptor,nativeLigand);
MoleculeGrid grid = new MoleculeGrid(nativeLigand,DockingEngine.GRID_RESOLUTION,
new Coordinates(DockingEngine.GRID_DIMENSION,DockingEngine.GRID_DIMENSION,
DockingEngine.GRID_DIMENSION));
Set bindingSiteAtoms = new HashSet();
if(scoringFunction==ScoringFunction.CHEMPLP) {
DockingEngine.getBindingSiteAtoms(receptor, bindingSiteAtoms, grid, true);
engine = new ChemPLP(receptor,bindingSiteAtoms,grid);
}
else if(scoringFunction==ScoringFunction.IDOSCORE) {
DockingEngine.getBindingSiteAtoms(receptor, bindingSiteAtoms, grid, false);
engine = new IdoScore(receptor,bindingSiteAtoms, getReceptorAtomTypes(receptor),grid);
}
this.mcSteps = mcSteps;
this.random = new Random(LigandPose.SEED);
}
public DockingEngine(StereoMolecule receptor, StereoMolecule nativeLigand) {
this(receptor,nativeLigand,DEFAULT_NR_MC_STEPS,DEFAULT_START_POSITIONS,ScoringFunction.CHEMPLP, StartPosition.PHESA);
}
private double getStartingPositions(StereoMolecule mol, ConformerSet initialPos) {
ConformerSetGenerator confSetGen = new ConformerSetGenerator(100,ConformerGenerator.STRATEGY_LIKELY_RANDOM, false,
LigandPose.SEED);
ConformerSet confSet = confSetGen.generateConformerSet(mol);
double eMin = Double.MAX_VALUE;
Map ffOptions = new HashMap();
ffOptions.put("dielectric constant", 80.0);
if(startPosition==StartPosition.PHESA) { //make PheSA Prealignment
DescriptorHandlerShapeOneConf dhs = new DescriptorHandlerShapeOneConf();
PheSAMolecule refVol = dhs.createDescriptor(nativeLigand);
TreeSet alignments = new TreeSet((e1,e2) -> {
return Double.compare(e1.getSim(), e2.getSim());});
for(Conformer conformer : confSet) {
if(conformer!=null) {
PheSAMolecule fitVol = dhs.createDescriptor(conformer.toMolecule());
double sim = dhs.getSimilarity(refVol, fitVol);
PheSAResult result = new PheSAResult(dhs.getPreviousAlignment()[0],dhs.getPreviousAlignment()[1],
sim);
alignments.add(result);
}
}
Set sortedAlignments = alignments.descendingSet();
for(PheSAResult res : sortedAlignments) {
StereoMolecule alignedMol = res.getFitMol();
if(initialPos.size()>=startPositions)
break;
ForceFieldMMFF94 mmff = new ForceFieldMMFF94(alignedMol, ForceFieldMMFF94.MMFF94SPLUS, ffOptions);
PositionConstraint constraint = new PositionConstraint(alignedMol,50,0.2);
mmff.addEnergyTerm(constraint);
mmff.minimise();
Conformer ligConf = new Conformer(alignedMol);
initialPos.add(ligConf);
if(initialPos.size()>=startPositions)
break;
double e = mmff.getTotalEnergy();
if(e=startPositions)
break;
double e = mmff.getTotalEnergy();
if(e ffOptions = new HashMap();
ffOptions.put("dielectric constant", 80.0);
if(ForceFieldMMFF94.table(ForceFieldMMFF94.MMFF94SPLUS)==null)
ForceFieldMMFF94.initialize(ForceFieldMMFF94.MMFF94SPLUS);
ConformerSet startPoints = new ConformerSet();
double eMin = getStartingPositions(mol, startPoints);
for(Conformer ligConf : startPoints) {
for(double[] transform : PheSAAlignment.initialTransform(1)) {
Conformer newLigConf = new Conformer(ligConf);
PheSAAlignment.rotateMol(newLigConf, transform);
LigandPose pose = initiate(newLigConf,eMin);
double energy = mcSearch(pose);
if(energy bindingSiteAtoms, MoleculeGrid grid,
boolean includeHydrogens) {
int[] gridSize = grid.getGridSize();
int atoms = receptor.getAtoms();
if(includeHydrogens)
atoms = receptor.getAllAtoms();
for(int i=0;i0 && x0 && y0 && z ffOptions = new HashMap();
ffOptions.put("dielectric constant", 80.0);
if(ForceFieldMMFF94.table(ForceFieldMMFF94.MMFF94SPLUS)==null)
ForceFieldMMFF94.initialize(ForceFieldMMFF94.MMFF94SPLUS);
Molecule3D nativePose = new Molecule3D(nativeLigand);
new Canonizer(nativePose);
ForceFieldMMFF94 mmff = new ForceFieldMMFF94(nativePose, ForceFieldMMFF94.MMFF94SPLUS, ffOptions);
PositionConstraint constraint = new PositionConstraint(nativePose,50,0.2);
mmff.addEnergyTerm(constraint);
mmff.minimise();
ConformerSetGenerator confSetGen = new ConformerSetGenerator(100,ConformerGenerator.STRATEGY_LIKELY_RANDOM, false,
LigandPose.SEED);
ConformerSet confSet = confSetGen.generateConformerSet(nativePose);
double eMin = Double.MAX_VALUE;
ConformerSet initialPos = new ConformerSet();
for(Conformer conformer : confSet) {
if(conformer!=null) {
StereoMolecule conf = conformer.toMolecule();
conf.ensureHelperArrays(Molecule.cHelperParities);
mmff = new ForceFieldMMFF94(conf, ForceFieldMMFF94.MMFF94SPLUS, ffOptions);
constraint = new PositionConstraint(conf,50,0.2);
mmff.addEnergyTerm(constraint);
mmff.minimise();
Conformer ligConf = new Conformer(conf);
initialPos.add(ligConf);
if(initialPos.size()>=startPositions)
break;
double e = mmff.getTotalEnergy();
if(e © 2015 - 2025 Weber Informatics LLC | Privacy Policy