Download Notebook View in GitHub Open in Google Colab
Mixing Sage and Amber force fields: BRD4 benchmark
This example applies SMIRNOFF-format parameters to BRD4 inhibitors from the living review on binding free energy benchmark systems by Mobley and Gilson. The BRD4 system comes from the accompanying GitHub repository.
# Retrieve protein and ligand files for BRD4 and a docked inhibitor from the benchmark systems GitHub repository
# https://github.com/MobleyLab/benchmarksets
import requests
repo_url = (
"https://raw.githubusercontent.com/MobleyLab/benchmarksets/master/input_files/"
)
sources = {
"receptor.pdb": repo_url + "BRD4/pdb/BRD4.pdb",
"ligand.pdb": repo_url + "BRD4/pdb/ligand-1.pdb",
"ligand.sdf": repo_url + "BRD4/sdf/ligand-1.sdf",
}
for filename, url in sources.items():
r = requests.get(url)
open(filename, "w").write(r.text)
from openff.toolkit import ForceField, Molecule, Quantity, Topology
Combining parametrization results from different force fields requires the use of Interchange.combine
, which is a newer feature and has not been tested for every conceivable combination of inputs. It may not work as expected for all systems. If extending this workflow, please validate results to avoid surprises in simulation results!
ligand_molecule = Molecule.from_file("ligand.sdf")
sage = ForceField("openff-2.2.0.offxml")
ligand_system = sage.create_interchange(topology=ligand_molecule.to_topology())
receptor_topology = Topology.from_pdb("receptor.pdb")
ff14sb = ForceField("ff14sb_off_impropers_0.0.4.offxml")
receptor_system = ff14sb.create_interchange(topology=receptor_topology)
complex_system = receptor_system.combine(ligand_system)
# Make pseudo-vacuum because GROMACS does not support proper vacuum simulations
complex_system.box = Quantity([5, 5, 5], "nanometer")
/home/runner/micromamba/envs/openff-docs-examples/lib/python3.11/site-packages/openff/interchange/operations/_combine.py:77: InterchangeCombinationWarning: Interchange object combination is complex and likely to produce strange results outside of a limited set of use cases it has been tested in. Any workflow using this method is not guaranteed to be suitable for production or stable between versions. Use with extreme caution and thoroughly validate results!
warnings.warn(
/home/runner/micromamba/envs/openff-docs-examples/lib/python3.11/site-packages/openff/interchange/operations/_combine.py:57: InterchangeCombinationWarning: Found electrostatics 1-4 scaling factors of 5/6 with slightly different rounding (0.833333 and 0.8333333333). This likely stems from OpenFF using more digits in rounding 1/1.2. The value of 0.8333333333 will be used, which may or may not introduce small errors.
warnings.warn(
Export to OpenMM
By default, Interchange.to_openmm()
returns an openmm.System
object. There is also Interchange.to_openmm_simulation()
which requires an openmm.Integrator
as an argument and returns a openmm.app.Simulation
object containing information about the topology, box vectors, and particle positions.
See the Interchange user docs for more optional features.
complex_system.to_openmm()
<openmm.openmm.System; proxy of <Swig Object of type 'OpenMM::System *' at 0x7ff72c334600> >
Export to Amber
Amber parameter/topology (.prmtop
), coordinate (.inpcrd
) and run control (.in
) files can each be written out separately with individual Interchange.to_x
methods. These are wrapped up into a convenient Interchange.to_amber
method which handles them all at once.
See the Interchange user docs for more.
complex_system.to_amber("complex")
/home/runner/micromamba/envs/openff-docs-examples/lib/python3.11/site-packages/openff/interchange/components/mdconfig.py:473: UserWarning: Ambiguous failure while processing constraints. Constraining h-bonds as a stopgap.
warnings.warn(
/home/runner/micromamba/envs/openff-docs-examples/lib/python3.11/site-packages/openff/interchange/components/mdconfig.py:401: SwitchingFunctionNotImplementedWarning: A switching distance 8.0 angstrom was specified by the force field, but Amber does not implement a switching function. Using a hard cut-off instead. Non-bonded interactions will be affected.
warnings.warn(
Export to GROMACS
Gromacs topology (.top
), coordinate (.gro
) and run control (.mdp
) files can each be written out separately with individual Interchange.to_x
methods. These are wrapped up into a convenient Interchange.to_gromacs
method which handles them all at once.
See the Interchange user docs for more.
complex_system.to_gromacs(prefix="complex")
/home/runner/micromamba/envs/openff-docs-examples/lib/python3.11/site-packages/openff/interchange/components/mdconfig.py:473: UserWarning: Ambiguous failure while processing constraints. Constraining h-bonds as a stopgap.
warnings.warn(