Glide

This tutorial demonstrates common operations using the Glide API, including generating and reading grids; scoring poses or computing their energy; and performing individual steps of the docking funnel (conformational generation, docking, torsion and placement optimization, sampling, and post-docking minimization). It also provides examples of lower-level API calls for accessing information like rotatable bonds or hydrogen bonds.

See Sample scripts at the bottom of this document for some simplified but complete examples that use the Glide API.

Note that this document focuses on the Python wrappers for the C++ library, which is a relatively low-level way of doing Glide calculations. For a higher-level interface, see schrodinger.application.glide.runner.

As an alternative way of exploring the API, we ship a Jupyter notebook as part of the Schrödinger installation; it can be found under $SCHRODINGER/glide-v*/tutorial/notebooks/glide_api_example.ipynb (the exact location depends on the current Glide version number). This notebook goes through the docking funnel for a single ligand and creates some structure representations, tables, histograms, and other plots along the way. It can be useful for understanding the Glide methodology better, or for troubleshooting a ligand which is not producing a good pose.

Getting started

First, let’s import the module and get a Schrodinger Structure object to play with:

>>> from schrodinger import glide
>>> from schrodinger import structure
>>> from schrodinger.test import mmshare_data_file

>>> lig_file = mmshare_data_file('glide/factorXa_reference.maegz')
>>> lig_st = structure.Structure.read(lig_file)

Next, let’s get the basic Glide objects needed to score a pose: the Pose and the GridArchive. For now, we’ll get the GridArchive by reading a file:

# fromStructure() does topological analysis of the ligand: atom typing,
# rings, rotatable bonds, etc.
>>> lig = glide.Ligand.fromStructure(lig_st)

# A Pose is effectively a mutable Ligand that can be optimized or scored.
>>> pose = glide.Pose(lig)

# A GridArchive represents the contents of the grid.zip file and is needed
# for docking or scoring.
>>> grid_file = mmshare_data_file('glide/factorXa_grid.zip')
>>> grid = glide.GridArchive.readArchive(grid_file)

Licenses

Now let’s check out the licenses. This can be done with a context manager, either DockingLicense or GridgenLicense. For example:

with glide.DockingLicense():
    # do stuff requiring docking license

Note that the functions in the previous examples, such as readArchive() and fromStructure(), did not require a license; many lower-level functions in the Glide API don’t require a license, but major functionality such as scoring, confgen, docking, and gridgen do.

Best practice is to establish a single license context at a high level, such as in a main() function. This is both more concise and more efficient. That said, in the examples below, we’ll use multiple license contexts to emphasize which functions require a license.

Scoring

Next we can score the pose, using the default options:

>>> config = glide.Config()
>>> with glide.DockingLicense():
...     score_info = pose.computeScore(grid, config)
>>> score_info.gscore
-11.52...

We can inspect the scoring terms, which add up to the total score (see ScoreTerms):

>>> terms = score_info.terms
>>> terms.ecoul
-3.25...

>>> terms.lipo
-3.44...

Or we can convert the Pose back into a Structure, with Glide properties added:

>>> with glide.DockingLicense():
...     st = pose.getScoredStructure(grid, config)
>>> st.property['r_i_glide_gscore']
-11.52...

Now let’s go back and inspect some Ligand properties, such as rotatable bonds:

>>> len(lig.getRotatableBonds())
8

>>> rb0 = lig.getRotatableBonds()[0]
>>> rb0.atom1, rb0.atom2, rb0.type
(6, 35, <RotType.PERIPHERAL: 1>)

Now look at hydrogen bonds between pose and receptor. Note that this does not require a GridArchive, as it is also possible to create a Receptor from a Structure (see below under Grid generation for how to do that).

>>> recep = grid.getReceptor()
>>> hbonds = glide.get_hbonds(pose, recep)
>>> len(hbonds.hbonds)
5

>>> hb0 = hbonds.hbonds[0]
>>> hb0.lig_atom.getIndex(), hb0.recep_atom.getIndex(), hb0.type
(30, 1422, <HBondType.CHARGED_CHARGED: 2>)

Next, let’s compute the intermolecular energy and forces using the grid. Note that the returned forces are indexed at 0:

>>> with glide.DockingLicense():
...     ec, forces = glide.get_intermolecular_energy_and_forces(pose, grid)
>>> ec.coul
-32.70...

>>> ec.vdw
-41.12...

Constraints

Some of the functions below require a LigandConstraintData object as an argument. This contains all the information necessary to support all kinds of constraints, including core constraints, receptor-based constraints, and excluded volumes. To create this object based on the job configuration and grid:

>>> jcd = glide.JobConstraintData(grid, config)
>>> constraints = glide.LigandConstraintData(lig, jcd)

For jobs that don’t need any constraints, it is possible to create an default LigandConstraintData in a single step:

constraints = glide.LigandConstraintData()

Confgen

Doing a conformational search is as simple as creating an instance of a ConformerEnsemble from a Ligand.

For the rest of this example, we’ll use a new Config object with HTVS precision to keep the test quick. Config objects can be created from dicts, from filenames (a Glide input file) or from GlideJob objects.

>>> htvs_config = glide.Config({'PRECISION': 'HTVS'})
>>> with glide.DockingLicense():
...     confs = glide.ConformerEnsemble.runConfgen(lig, htvs_config, constraints)

The main purpose of a ConformerEnsemble is to be an input to docking, but we can also inspect or export individual conformers:

>>> conf = confs.getConformers()[0]
>>> print(conf.getTorsions())
(...)
>>> # (3.81, 0.79, 0.23, 3.66, 5.77, -0.278, 5.79, 0.05, -0.07)

conf_st = conf.getStructure()
conf_st.write('my_conf.maegz')

Docking

For docking, we use the Docker class, which requires a ConformerEnsemble as its main input.

Also, Docker, among others, requires a DataFiles object, which holds some of the global data, normally obtained from files in the data directory of the installation.

>>> datafiles = glide.DataFiles()
>>> with glide.DockingLicense():
...     docker = glide.Docker(grid, htvs_config, datafiles)
...     pose_ensemble = docker.dockConformerEnsemble(confs, constraints)

Here pose_ensemble is a RefinedPoseEnsemble. One thing we can do with it is extract the poses as Structure objects, for example to write them to a file:

with structure.StructureWriter('ref_poses.maegz') as writer:
    writer.extend(pose.getStructure() for pose in pose_ensemble.getPoses())

NOTE: these are the MAXKEEP refined poses, which have not been minimized yet.

Torsion and placement optimization

The next step after docking and refinement is a quick optimization of the ligand placement and the torsional angles of the rotatable bonds, which is implemented as a method of the RefinedPoseEnsemble. Typically, we only minimize a subset of the poses; let’s use 20 for this example (in a real Glide job, this number is obtained from the MAXREF input keyword).

Right after the energy minimization, we’ll compute Emodel, the scoring function that is used for picking the best poses for further processing. The computeEmodels() call also sorts the poses in the ensemble by Emodel.

>>> pose_ensemble.truncate(20)  # Keep only the best 20 poses.
>>> with glide.DockingLicense():
...      pose_ensemble.minimizePoses(grid, htvs_config, constraints)
...      pose_ensemble.computeEmodels(grid, htvs_config, datafiles, constraints)

Sampling

Next, a small number of the top-ranked poses are then subjected to a sampling procedure in which alternative local-minima core and rotamer-group torsion angles are examined to try to improve the energy score.

>>> with glide.DockingLicense():
...     pose_ensemble.resampleConformers(grid, htvs_config, constraints, datafiles)

Before going to the next stage, we need to turn some of the top RefinedPose objects into Pose objects. Let’s go with 5, which is the default value of the POSTDOCK_NPOSE keyword:

>>> poses = [glide.Pose(p) for p in pose_ensemble.getPoses()[:5]]

Post-docking Minimization

PDM is as simple as calling the minimize() method of a FullMinimizer object on each pose:

>>> minimizer = glide.FullMinimizer(grid, htvs_config, constraints)
>>> with glide.DockingLicense():
...     for pose in poses:
...         minimizer.minimize(pose)
<...>

Finally, as described previously in the Scoring section, we can convert our scored poses into Structure objects, and write them or examine them or use them as inputs to other APIs.

with structure.StructureWriter('poses.maegz') as writer:
    writer.extend(pose.getScoredStructure() for pose in poses)

Grid generation

Let’s create a Config object to define the location and size of the grid box:

>>> gridgen_config = glide.Config({
...     "GRID_CENTER":   [10.8, 24.3, 5.6],
...     "OUTERBOX":   25.0,
... })

Using the Config and a Structure, we can create a Receptor:

>>> recep_file = mmshare_data_file('glide/glide_test_rec.mae.gz')
>>> recep_st = structure.Structure.read(recep_file)
>>> recep = glide.Receptor.fromStructure(recep_st, gridgen_config)

Finally, run the grid generation, which returns a GridArchive. Note that this requires a gridgen license!:

>>> with glide.GridgenLicense():
...     new_grid = glide.GridArchive.runGridgen(recep, gridgen_config, datafiles)

And it’s done! The function calls are all silent by default, but it is possible to obtain human-readable information about the gridgen job:

>>> print(new_grid.getInfo())
Receptor atoms: 584
Sites: 47
Rough grid cells: (25, 25, 25)
Fine grid cells: (64, 64, 64)
Adaptive grid points: 41046
Force field: OPLS_2005
Partial charges: from force field
Charge/vdW scaling: from RECEP_VSCALE / RECEP_CCUT
Movable atoms: 0
Excluded volumes: 0
Peptide mode: false
Constraints: 0

Finally, let’s save the grid archive:

new_grid.writeZip('grid.zip')

But note that we don’t need to save it; if we want to use it for docking or scoring as part of the same process, we can use the GridArchive object as is, without having to involve the filesystem at all.

Sample scripts

These scripts are all found in the installation under $SCHRODINGER/internal/bin/samples.

1. gridgen_low_level.py: how to generate a grid from a Structure using the low-level API. The ligand structure is used to determine the size and position of the grid box.

2. docking_low_level.py: how to dock a Structure into a GridArchive using the low-level API.

3. workflow_example.py: toy example of the use of runner.run_docking_workflow(), using a StructureReader for input and a StructureWriter for output.

4. fragment_expansion_screening.py: perform fragment enumeration followed by docking and scoring of the enumerated molecules. Fragment enumeration is performed by the multiroute module, part of Pathfinder application. Docking and scoring is performed by Glide module.