Omega.py

From DISI
Revision as of 22:42, 26 April 2013 by Oliv Eidam (talk | contribs) (Created page with "==Changing settings for Omega (OpenEye)== You may be interested to change the default settings for Omega (defined in omega.py) when generating conformations for docking. Maby...")
(diff) ← Older revision | Latest revision (diff) | Newer revision → (diff)
Jump to navigation Jump to search

Changing settings for Omega (OpenEye)

You may be interested to change the default settings for Omega (defined in omega.py) when generating conformations for docking. Mabye you want to increase the number of conformations sampled per ligand, or the force field used.

The easiest way to change the settings for Omega is to copy the default omega.py file into the directory where you want to generate ligand conformations using dbgen.csh.

Type:

which omega.py

Copy the output of the which command above into your current working directory:

cp /raid3/software/dockenv/etc/omega.py .

Rename that file to omega.parm:

mv omega.py omega.parm

Edit omega.parm in your favorite editor.

Run dbgen.csh with altered settings for omega:

dbgen.csh my_ligs.smi

Increase ligand conformation sampling or change the type of force field

An example of omega.py is shown at the bottom of this page. Settings you may want to consider changing are:


    omega.SetEnergyWindow(12.5)

    omega.SetMaxConfs(600)

    omega.SetRMSThreshold(0.80)
    omega.SetBuildForceField('mmff94s_NoEstat') 
    omega.SetSearchForceField('mmff94s_NoEstat')
  • you could increase the energy window: e.g. omega.SetEnergyWindow(25.0)

=> There is no benchmarking of this setting that I am aware of...

  • you could increase the number of output conformations: e.g. omega.SetMaxConfs(10000)

=> Molecules in ZINC have a cut a 600 conformations, which is certainly on the low side for flexible ligands

  • you could lower the RMS threshold between the confs: e.g. omega.SetRMSThreshold(0.40)

=> This results in finer sampling of conformational space

  • you could change the force field to included electrostatics: e.g. omega.SetBuildForceField('mmff94s') and omega.SetSearchForceField('mmff94s')

=> Internal benchmarkings on DUDE (by RGC) have shown that electrostatics improves enrichment for many targets.


Example of omega.py

#!/usr/bin/env python
import os
import sys
from openeye.oechem import *
from openeye.oeomega import *

def read_parm(parm):
    if os.path.exists(parm):
        f = open(parm, 'r')
        for line in f:
            line = line.strip()
            if line and line[0] != '#':
                eval('omega.%s' % line)
        f.close()

def write_molecule(outfile, mol):
    ofs = oemolostream(outfile)
    OEWriteMolecule(ofs, mol)
    ofs.close()

def set_defaults(omega):
    # File Parameters
    omega.SetCommentEnergy(True)
    omega.SetIncludeInput(False)
    omega.SetRotorOffset(True)
    omega.SetSDEnergy(False)
    omega.SetWarts(True)
    # 3D Construction Parameters
    omega.SetBuildForceField('mmff94s_NoEstat')
    omega.SetCanonOrder(True)
    omega.SetFixDeleteH(True)
    omega.SetDielectric(1.0)
    omega.SetExponent(1.0)
    omega.SetFixRMS(0.15)
    omega.SetFromCT(False)
    omega.SetFixMaxMatch(1)
    omega.SetFixUniqueMatch(True)
    # Structure Enumeration Parameters
    omega.SetEnumNitrogen(False)
    omega.SetEnumRing(False)
    # Torsion Driving Parameters
    omega.SetEnergyWindow(12.5)
    omega.SetMaxConfGen(100000)
    omega.SetMaxConfs(600)
    omega.SetMaxPoolSize(10000)
    omega.SetMaxRotors(-1)
    omega.SetMaxSearchTime(120.0)
    omega.SetRangeIncrement(5)
    omega.SetRMSThreshold(0.80)
    omega.SetSearchForceField('mmff94s_NoEstat')
    omega.SetTorsionDrive(True)

# Parse arguments here
parm_prefix = sys.argv[1]
infile = sys.argv[2]

parmfile = parm_prefix + ".parm"
noringparmfile = parm_prefix + "_noring.parm"

omega = OEOmega()
set_defaults(omega)

mol = OEMol()
inroot, inext = os.path.splitext(infile)
ifs = oemolistream(infile)
OEReadMolecule(ifs, mol)
ifs.close()
count, ringlist = OEDetermineRingSystems(mol)
rcf = open('ring_count', 'w')
rcf.write('%d\n' % count)
rcf.close()
if count == 0:
    omega.SetMaxConfs(30)
    read_parm(noringparmfile)
    outfile = "%s.%d.mol2" % (inroot, count) 
    if omega(mol):
        write_molecule(outfile, mol)
else:
    read_parm(parmfile)
    pred = OEPartPredAtom(ringlist)
    for i in xrange(1, count+1):
        pred.SelectPart(i)
        outfile = "%s.%d.mol2" % (inroot, i) 
        molcopy = OEMol(mol)
        if omega(molcopy, pred):
            write_molecule(outfile, molcopy)