From 54c0c1fad6e238ac6d7c150e05c806a7c70a4019 Mon Sep 17 00:00:00 2001 From: Finlay Clark Date: Wed, 2 Sep 2026 13:22:08 +0100 Subject: [PATCH 1/9] Update update small param notebook --- environment.yaml | 6 +- .../small_molecule_parameterisation.ipynb | 266 ++++++++++-------- 2 files changed, 148 insertions(+), 124 deletions(-) diff --git a/environment.yaml b/environment.yaml index b48b799..d10e1b1 100644 --- a/environment.yaml +++ b/environment.yaml @@ -3,8 +3,10 @@ channels: - conda-forge dependencies: - jupyterlab - - openff-toolkit-examples>=0.17.0 # New version needed for NAGL compat. - - openff-interchange>=0.4.7 # New version needed for NAGL compat. + - openff-toolkit-examples>=0.19.0 # Need at least 0.17.0 for NAGL compat. + - openff-interchange>=0.5.4 # Need at least 0.4.7 for NAGL compat. + - packmol # Needed to pack water box only + - ambertools # Needed for AM1-BCC charge assignment comparison only - MDAnalysis - prolif - py3dmol # For 3D visualization of of interactions in prolif diff --git a/notebooks/small_molecule_parameterisation.ipynb b/notebooks/small_molecule_parameterisation.ipynb index e3182d2..cc08160 100644 --- a/notebooks/small_molecule_parameterisation.ipynb +++ b/notebooks/small_molecule_parameterisation.ipynb @@ -18,9 +18,10 @@ "\n", "| Action | Software|\n", "|--|--|\n", + "| [Go from smiles to simulation in few lines of code](#showcase) | OpenFF Toolkit, OpenFF Interchange, OpenMM\n", "| [Load and inspect a force field](#loading_ff) | OpenFF Toolkit\n", "| [Create a representation of your chemical system](#topology) | OpenFF Toolkit\n", - "| [Parameterise your system and run a quick simulation](#interchange) | OpenFF Interchange, OpenMM\n", + "| [Parameterise your system and run a quick simulation in water](#interchange) | OpenFF Interchange, OpenMM\n", "| [Rapidly assign partial charges with a graph neural network model](#gnn_charges) | OpenFF Toolkit, OpenFF NAGL Models\n", "| [Review what you've learnt](#summary) | \n", "| [Check out other OpenFF tutorials](#further_materials) | \n", @@ -40,6 +41,91 @@ "\n" ] }, + { + "cell_type": "markdown", + "id": "81bb0b4f", + "metadata": {}, + "source": [ + "\n", + "## 0. You can go from SMILES to simulation in a few lines of code" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "cef30154", + "metadata": {}, + "outputs": [], + "source": [ + "# Go from SMILES -> simulation input with OpenFF\n", + "from openff.toolkit import ForceField, Molecule, Topology\n", + "\n", + "molecule = Molecule.from_smiles(\"CC(=O)Nc1ccc(cc1)O\")\n", + "molecule.generate_conformers(n_conformers=1)\n", + "topology = Topology.from_molecules([molecule])\n", + "\n", + "force_field = ForceField(\"openff-2.3.0.offxml\")\n", + "interchange = force_field.create_interchange(topology)\n", + "interchange.minimize()\n", + "\n", + "openmm_system = interchange.to_openmm_system()\n", + "openmm_topology = interchange.to_openmm_topology()\n", + "openmm_positions = interchange.positions.to_openmm()" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "11ec8a77", + "metadata": {}, + "outputs": [], + "source": [ + "# Run the simulation with OpenMM\n", + "import openmm\n", + "\n", + "temperature = 298.15 * openmm.unit.kelvin\n", + "friction_coefficient = 1.0 / openmm.unit.picosecond\n", + "step_size = 2.0 * openmm.unit.femtosecond\n", + "\n", + "simulation = openmm.app.Simulation(\n", + " openmm_topology,\n", + " openmm_system,\n", + " openmm.LangevinIntegrator(temperature, friction_coefficient, step_size),\n", + ")\n", + "simulation.context.setPositions(openmm_positions)\n", + "simulation.context.setVelocitiesToTemperature(simulation.integrator.getTemperature())\n", + "\n", + "simulation.reporters.append(\n", + " openmm.app.DCDReporter(file=\"trajectory_showcase.dcd\", reportInterval=100)\n", + ")\n", + "simulation.step(10000)" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "a2983975", + "metadata": {}, + "outputs": [], + "source": [ + "# Load the trajectory with MDAnalysis and visualise with nglview\n", + "import MDAnalysis as mda\n", + "import nglview\n", + "\n", + "u = mda.Universe(openmm_topology, \"trajectory_showcase.dcd\")\n", + "\n", + "view = nglview.show_mdanalysis(u)\n", + "view" + ] + }, + { + "cell_type": "markdown", + "id": "e110064b", + "metadata": {}, + "source": [ + "That's it! You've run a vacuum simulation for paracetamol. Below and in the next notebook, we'll go into more detail on each of the steps in the OpenFF cell and show how you can set up more complex systems, but this is mainly for your understanding and you rarely need much more code than shown above." + ] + }, { "cell_type": "markdown", "id": "32c4b5c9", @@ -48,7 +134,11 @@ "\n", "## 1. Force fields are specified in `.offxml` files and can be loaded with the `ForceField` class\n", "\n", - "OpenFF's force fields use the The SMIRKS Native Open Force Field (SMIRNOFF) [specification](https://openforcefield.github.io/standards/standards/smirnoff/) and are conventionally encoded in `.offxml` files. The spec fully describes the contents of a SMIRNOFF force field, how parameters should be applied, and several other important usage details. You could implement a SMIRNOFF engine in your own code, but conveniently the OpenFF Toolkit already provides this and a handful of utilities. Let's load up the latest OpenFF small molecule force field, OpenFF 2.2.1, and inspect its contents! This force field shares the code name \"Sage\" with all other force fields with the same major version number (2.x.x)." + "Let's dive into the details of what went on above. Here's a summary of how data flows through a workflow utilising OpenFF tools -- the OpenFF toolkit allows you to create `Molecule` and `ForceField` objects, which get combined into an `Interchange` object, which contains all the information needed to start a simulation. From there, you can create input for the simulation engine of your choice:\n", + "\n", + "\"Description\n", + "\n", + "Let's start with the `.offxml` force field file. OpenFF's force fields use the The SMIRKS Native Open Force Field (SMIRNOFF) [specification](https://openforcefield.github.io/standards/standards/smirnoff/). The spec fully describes the contents of a SMIRNOFF force field, how parameters should be applied, and several other important usage details. You could implement a SMIRNOFF engine in your own code, but conveniently the OpenFF Toolkit already provides this and a handful of utilities. Let's load up the latest OpenFF small molecule force field, OpenFF 2.3,0, and inspect its contents! This force field shares the code name \"Sage\" with all other force fields with the same major version number (2.x.x)." ] }, { @@ -67,7 +157,7 @@ "source": [ "from openff.toolkit import ForceField\n", "\n", - "sage = ForceField(\"openff-2.2.1.offxml\")\n", + "sage = ForceField(\"openff-2.3.0.offxml\")\n", "sage" ] }, @@ -76,7 +166,7 @@ "id": "52c754d3", "metadata": {}, "source": [ - "If you'd like to see the raw file on disk that's being parsed, [here's the file on GitHub](https://github.com/openforcefield/openff-forcefields/blob/main/openforcefields/offxml/openff-2.2.1.offxml)." + "If you'd like to see the raw file on disk that's being parsed, [here's the file on GitHub](https://github.com/openforcefield/openff-forcefields/blob/main/openforcefields/offxml/openff-2.3.0.offxml)." ] }, { @@ -170,7 +260,7 @@ "\n", "The toolkit uses these SMIRKS patterns and direct chemical perception to assign parameters to particular atoms (or bonds, angles, etc.).\n", "\n", - "We'll use OpenFF 2.2.1 for the remainder of this tutorial, but you can learn more about this and other SMIRNOFF force fields below:\n", + "We'll use OpenFF 2.3.0 for the remainder of this tutorial. This is OpenFF's latest small molecule force field and is a leading open-source small molecule force field which [performs comparably to other open-source force fields](https://doi-org.libproxy.ncl.ac.uk/10.1021/acs.jctc.3c00039). You can learn more about this and other SMIRNOFF force fields below:\n", "
\n", " Click here to learn about available and planned SMIRNOFF force fields\n", "\n", @@ -192,7 +282,7 @@ "\n", "The Sage line of force fields (`openff-2.y.z.offxml`) continued the process of fitting to more (and more diverse) QM datasets, but also included a re-fit of the Lennard-Jones parameters. Small molecule geometries and energies [improved, in general,](https://openforcefield.org/community/news/general/sage2.0.0-release/) significantly over Parsley. These improvements notably transferred to protein-ligand binding free energies despite Sage not being specifically fit to them. For more, see the [associated paper](https://pubs.acs.org/doi/10.1021/acs.jctc.3c00039).\n", "\n", - "[Subsequent releases](https://github.com/openforcefield/openff-forcefields/releases) used different fitting procedures and tweaks to parameter typing to improve performance and address issues with several specific chemistries. The latest release, **Sage 2.2.1 (`openff-2.2.1.offxml`) is the recommended force field for small molecule studies.**\n", + "[Subsequent releases](https://github.com/openforcefield/openff-forcefields/releases) used different fitting procedures and tweaks to parameter typing to improve performance and address issues with several specific chemistries. Notably, Sage 2.3.0 includes fast graph neural network charge assignment with Ash-GC, which is discussed later in this notebook. This charge model is trained to reproduce AM1-BCC charges without the typical $O(N^3)$ scaling, making it suitable for large (>> 100 atoms) molecules). The latest release, **Sage 2.3.0 (`openff-2.3.0.offxml`) is the recommended force field for small molecule studies.**\n", "\n", "## Ports\n", "\n", @@ -212,6 +302,13 @@ "\n", "OpenFF, in collaboration with Dave Cerutti of the Amber community, created a port of [ff14SB](https://pubs.acs.org/doi/10.1021/acs.jctc.5b00255), a popular Amber protein force field. There are some small numerical differences with how improper torsions are evaluated, but all other terms reproduce a canonical Amber source to high accuracy. **This is the only protein force field currently in SMIRONOFF (`.offxml`) format** and therefore the current recommendation for use with proteins. Primarily for technical reasons, porting other Amber force fields is not planned.\n", "\n", + "## Rosemary Alpha\n", + "\n", + "A future line of force fields from OpenFF (code name \"Rosemary\", starting with `openff-3.0.0.offxml`) is intended to handle small molecules and biopolymers in a _self-consistent_ manner. This is exciting as it will streamline simulations of proteins with non-cannonical amino acids! See the workshop [Simulating Post-Translationally Modified Proteins with the OpenFF Rosemary Alpha](https://github.com/openforcefield/2026-virtual-workshops/blob/main/ptm/ptm-workshop.ipynb). The first release will handle proteins, but future versions may cover nucleic acids. The performance, depending on the metrics used, is hoped to be comparable with existing Amber-family protein force fields. \n", + "\n", + "An pre-release version of Rosemary is available for testing as [`openff_no_water-3.0.0-alpha0.offxml`](https://github.com/openforcefield/openff-forcefields/blob/main/openforcefields/offxml/openff_no_water-3.0.0-alpha0.offxml). If you use it, see [the release notes](https://github.com/openforcefield/openff-forcefields/releases/tag/2025.10.1). \n", + "\n", + "\n", "## Non-main-line force fields\n", "\n", "### SMIRNOFF plugins\n", @@ -225,15 +322,7 @@ "## From OpenFF\n", "\n", "### Rosemary\n", - "A future line of force fields from OpenFF (code name \"Rosemary\", starting with `openff-3.0.0.offxml`) is intended to handle small molecules and biopolymers in a _self-consistent_ manner. The first release is expected to handle proteins, but future versions may cover nucleic acids. The performance, depending on the metrics used, is hoped to be comparable with existing Amber-family protein force fields.\n", - "\n", - "There is no specific release date planned for Rosemary, but it may be available in 2026 (a beta release candidate may also be publically available prior to the full release).\n", - "\n", - "### Graph net charge assignment\n", - "\n", - "TODO: UPDAATE AND MENTION 2.3\n", - "\n", - "The Sage 2.3.0 release is expected imminently and will include graph-convolutional neutral network (GCNN)-based charge assignment using [NAGL](https://github.com/openforcefield/openff-nagl) by default. The charge model is trained to reproduce AM1-BCC charges without the typical $O(N^3)$ scaling, making it suitable for large (>> 100 atoms) molecules). The [second release candidate](https://github.com/openforcefield/openff-forcefields/blob/main/openforcefields/offxml/openff-2.3.0-rc2.offxml) (which may or may not become the final version) is already available for you to try!\n", + "There is no specific release date planned for the first full version of Rosemary, but it may be available in late 2026.\n", "\n", "### Virtual sites\n", "\n", @@ -253,7 +342,7 @@ "\n", "## 2. The `Topology` class represents a chemical system containing one or more `Molecule`s\n", "\n", - "Now we've loaded our desired force field (OpenFF 2.2.1), we need to specify the chemical system we want to assign force field parameters to (\"parameterise\"). Our system will be represented by a `Topology`, which we will build from one or more `Molecule`s. \n", + "Now we've loaded our desired force field (OpenFF 2.3.0), we need to specify the chemical system we want to assign force field parameters to (\"parameterise\"). Our system will be represented by a `Topology`, which we will build from one or more `Molecule`s. \n", "\n", "As a simple example, let's build a `Topology` containing an small molecule with some features which illustrate how parameters are applied according to SMIRKS matches. We'll use the crotonate anion, but you could draw any molecule you like and convert it to a SMILES string using tools like ChemDraw and [MolView](https://molview.org/)." ] @@ -416,7 +505,7 @@ "metadata": {}, "source": [ "
\n", - " ✏️ Exercise: Build a Topology containing an MCL-1 ligand. Create the Molecule from an SDF file (take a look at the docstring of Molecule to see how this can be done). Also, see Molecule cookbook for all the ways to make a Molecule. The crystallographic MCL-1 ligand from PDB ID 6o6f is provided at ../structures/606f_ligand.sdf. Note that you don't need to use get_data_file_path as we already know the path.\n", + " ✏️ Exercise: Build a Topology containing an MCL-1 ligand. Create the Molecule from an SDF file (take a look at the docstring of Molecule to see how this can be done, noting that you can just pass the sdf path given below and don't need the get_data_file_path function). Also, see Molecule cookbook for all the ways to make a Molecule. The crystallographic MCL-1 ligand from PDB ID 6o6f is provided at ../structures/606f_ligand.sdf. Note that you don't need to use get_data_file_path as we already know the path.\n", "
\n" ] }, @@ -439,7 +528,7 @@ "id": "0c9e16d5", "metadata": {}, "source": [ - "We will cover creating a topology for a protein-ligand complex this afternoon." + "We will cover creating a topology for a protein-ligand complex in the next notebook." ] }, { @@ -452,17 +541,7 @@ "\n", "Now we've specified our force field and our chemical system using classes from the OpenFF Tools package (`ForceField`, `Molecule`, and `Topology`), and we want to apply our force field to our chemical topologies (parameterisation).\n", "\n", - "To do this, we'll use the `Interchange` class from the OpenFF Interchange package, which stores a fully-parameterised molecular system and provides methods to write out simulation-ready input files for a number of software packages. They key objective of Interchange is to provide an intermediate inspectable state after parameterisation and before conversion to an engine-specific format. For most users, an `Interchange` forms the bridge between the OpenFF ecosystem and their simulation software of choice. The current focus is applying SMIRNOFF force fields to chemical topologies and exporting the result to engines preferred by our users. In order of stability, OpenMM, GROMACS, Amber, and LAMMPS are supported. Future development may include support for CHARMM and other engines.\n", - "\n", - "Below is a summary of how data flows through a workflow utilising OpenFF tools, including where Interchange sits in the flow." - ] - }, - { - "cell_type": "markdown", - "id": "ddde1c2a", - "metadata": {}, - "source": [ - "\"Description" + "To do this, we'll use the `Interchange` class from the OpenFF Interchange package, which stores a fully-parameterised molecular system and provides methods to write out simulation-ready input files for a number of software packages. They key objective of Interchange is to provide an intermediate inspectable state after parameterisation and before conversion to an engine-specific format. For most users, an `Interchange` forms the bridge between the OpenFF ecosystem and their simulation software of choice. The current focus is applying SMIRNOFF force fields to chemical topologies and exporting the result to engines preferred by our users. In order of stability, OpenMM, GROMACS, Amber, and LAMMPS are supported. Future development may include support for CHARMM and other engines." ] }, { @@ -705,7 +784,7 @@ "id": "c1603b89", "metadata": {}, "source": [ - "We can see that the C=C bond (indices (1,2)) is associated with a potential key with the SMIRKS pattern `[#6X3:1]=[#6X3:2]` (specifying any two carbons bonded to 3 atoms connected by a double bond). Note that the (1,0) C-C bond is matched by the SMRIKS `[#6X3:1]-[#6X3:2]`, which specifies the atoms in the same way, showing that the parameters have been assigned by directly using information about the bond. This contrasts to traditional atom typing approaches, where information about the bond would be implicitly encoded in the atom types used to assign the parameters. Another example of this \"direct chemical perception\" is the assignment of the carboxylate carbon-oxygen bond parameters, which only match (triply-connected carbon) - (singly-connnected oxygen) bonds when the carbon is bonded to another singly-connected oxygen.\n", + "We can see that the C=C bond (indices (1,2)) is associated with a potential key with the SMIRKS pattern `[#6X3:1]=[#6X3:2]` (specifying any two carbons each bonded to 3 atoms and connected by a double bond). Note that the (1,0) C-C bond is matched by the SMRIKS `[#6X3:1]-[#6X3:2]`, which specifies the atoms in the same way, showing that the parameters have been assigned by directly using information about the bond. This contrasts to traditional atom typing approaches, where information about the bond would be implicitly encoded in the atom types used to assign the parameters. Another example of this \"direct chemical perception\" is the assignment of the carboxylate carbon-oxygen bond parameters, which only match (triply-connected carbon) - (singly-connnected oxygen) bonds when the carbon is bonded to another singly-connected oxygen.\n", "\n", "To see the actual parmeters specified for this bond, we can look up the `Potential` objects using the `PotentialKey`s." ] @@ -907,7 +986,7 @@ "id": "212385e0", "metadata": {}, "source": [ - "Here, we'll export to OpenMM and run a short simulation directly from the noteboook. We can create an OpenMM `Simulation` object from the `Interchange` and run for a specified wall clock time using `runForClockTime` (the simluation time will depend on how quickly it runs on your machine). We keep the volume ($V$), number of particles ($N$), and average temperature ($T$) (using the LangevinMiddleIntegrator) constant and the simulation corresponds to the $NVT$ ensemble." + "Here, we'll export to OpenMM and run a short simulation directly from the noteboook. We can create an OpenMM `Simulation` object from the `Interchange` and run for a specified wall clock time using `runForClockTime` (the simluation time will depend on how quickly it runs on your machine). We keep the volume ($V$), number of particles ($N$), and average temperature ($T$) (using the LangevinIntegrator) constant and the simulation corresponds to the $NVT$ ensemble." ] }, { @@ -927,7 +1006,7 @@ "import openmm\n", "import openmm.unit\n", "from openff.interchange import Interchange\n", - "import mdtraj\n", + "import MDAnalysis as mda\n", "import nglview\n", "\n", "\n", @@ -937,7 +1016,7 @@ " trajectory_name: str = \"small_mol_solvated.dcd\",\n", "):\n", " simulation = interchange.to_openmm_simulation(\n", - " integrator=openmm.LangevinMiddleIntegrator(\n", + " integrator=openmm.LangevinIntegrator(\n", " 300 * openmm.unit.kelvin,\n", " 1 / openmm.unit.picosecond,\n", " 0.002 * openmm.unit.picoseconds,\n", @@ -955,12 +1034,10 @@ " topology: Topology, filename: str = \"small_mol_solvated.dcd\"\n", ") -> nglview.NGLWidget:\n", " \"\"\"Visualise a trajectory using nglview.\"\"\"\n", - " traj = mdtraj.load(\n", - " filename,\n", - " top=mdtraj.Topology.from_openmm(topology.to_openmm()),\n", - " )\n", "\n", - " view = nglview.show_mdtraj(traj)\n", + " u = mda.Universe(topology.to_openmm(), filename)\n", + "\n", + " view = nglview.show_mdanalysis(u)\n", " view.add_representation(\"licorice\", selection=\"water\")\n", "\n", " return view\n", @@ -988,23 +1065,19 @@ "\n", "## 4. Graph Neural Networks Allow Fast Assignment of Partial Charges\n", "\n", - "You might notice that [Sage](https://github.com/openforcefield/openff-forcefields/blob/main/openforcefields/offxml/openff-2.2.1.offxml) doesn't contain tabulated charges for most atomic environments in the way it does for all other terms in the force field. Instead, it specifies:\n", + "You might notice that Sage force fields don't contain tabulated charges for most atomic environments in the way they do for all other terms in the force field.For example, [Sage 2.2.1](https://github.com/openforcefield/openff-forcefields/blob/main/openforcefields/offxml/openff-2.2.1.offxml) instead specifies:\n", "```\n", - "\n", + "\n", "```\n", - "which means that partial charges will be calculated using the common AM1-BCC method. Charges from a semi-empirical quantum chemistry calculation (Austin Model 1) are corrected (bond charge correction) to approximate charges obtained by fitting to the electrostatic potential at the HF/6-31G* level (see [Jakalian et al.](https://onlinelibrary.wiley.com/doi/10.1002/(SICI)1096-987X(20000130)21:2%3C132::AID-JCC5%3E3.0.CO;2-P)). Unfortunately, AM1-BCC scales 𝒪(N2) in the number of atoms N, making it prohibitively slow for large molecules and biopolymers.\n", + "which means that partial charges will be calculated using the common AM1-BCC method. Charges from a semi-empirical quantum chemistry calculation (Austin Model 1) are corrected (bond charge correction) to approximate charges obtained by fitting to the electrostatic potential at the HF/6-31G* level (see [Jakalian et al.](https://onlinelibrary.wiley.com/doi/10.1002/(SICI)1096-987X(20000130)21:2%3C132::AID-JCC5%3E3.0.CO;2-P)). Unfortunately, parameterisation with AM1-BCC using OpenEye or AmberTools scales 𝒪(N2) in the number of atoms N, making it prohibitively slow for large molecules and biopolymers.\n", "\n", - "Methods which assign partial charges using graph neural networks offer rapid assignment with better scaling. They also offer the possibility of going beyond traditionally affordable QM levels of theory by training to quickly reproduce charges from expensive calculations. For example, [EspalomaCharge](https://pubs.acs.org/doi/full/10.1021/acs.jpca.4c01287) is fit to AM1-BCC charges and offers 𝒪(N2) scaling, while [Adams et al.](https://chemrxiv.org/engage/chemrxiv/article-details/6839c94c3ba0887c33d2cd8e) trained models to reproduce atoms-in-molecules charges and electrostatic potentials obtained at a high level of theory. Here, we'll use OpenFF's [AshGC](https://zenodo.org/records/15770227/files/AshGC_methods_2025-06-30.pdf?download=1) model, which is trained to reproduce AM1-BCC charges." - ] - }, - { - "cell_type": "markdown", - "id": "ede41aef", - "metadata": {}, - "source": [ - "
\n", - " ⚠️ OpenFF 2.2.1 has not been explicitly trained and validated with AshGC charges. However, the 2.3.0 release will be, and is expected imminently. AshGC charges will be used as default and will be specified in the .offxml file, so there will be no need to call Molecule.assign_partial_charges as shown below.\n", - "
" + "Methods which assign partial charges using graph neural networks offer rapid assignment with better scaling. They also offer the possibility of going beyond traditionally affordable QM levels of theory by training to quickly reproduce charges from expensive calculations. For example, OpenFF's [AshGC](https://doi.org/10.1021/acs.jctc.6c00169) model is fit to AM1-BCC charges and offers 𝒪(N) scaling, while [Adams et al.](https://doi.org/10.1021/acs.jctc.5c01520) trained models to reproduce atoms-in-molecules charges and electrostatic potentials obtained at a high level of theory. AshGC is used by [Sage 2.3.0](https://github.com/openforcefield/openff-forcefields/blob/main/openforcefields/offxml/openff-2.3.0.offxml) -- if you inspect the file, you'll see:\n", + "```\n", + "\n", + "```\n", + "where \"openff-gnn-am1bcc-1.0.0.pt\" is the AshGC model.\n", + "\n", + "The GNN charge model is the main difference between Sage 2.2.1 and 2.3.0. Here, we'll compare the parameterisation speed and charges obtained using each force field." ] }, { @@ -1021,42 +1094,17 @@ }, "outputs": [], "source": [ - "from openff.toolkit import Molecule, ForceField\n", - "from openff.toolkit.utils.nagl_wrapper import NAGLToolkitWrapper\n", + "from openff.toolkit import Molecule\n", "\n", - "# Disable RDKit warnings to avoid misleading NAGL warnings\n", - "# (see https://github.com/openforcefield/openff-nagl/issues/198)\n", - "from rdkit import RDLogger\n", - "RDLogger.DisableLog('rdApp.*') \n", - "\n", - "# OpenFF NAGL store models as PyTorch files.\n", - "ASH_GC_MODEL = \"openff-gnn-am1bcc-0.1.0-rc.3.pt\"\n", "molecule = Molecule(\"../structures/6o6f_ligand.sdf\")" ] }, - { - "cell_type": "code", - "execution_count": null, - "id": "e0646ac8", - "metadata": { - "execution": { - "iopub.execute_input": "2025-10-09T15:21:01.492129Z", - "iopub.status.busy": "2025-10-09T15:21:01.491979Z", - "iopub.status.idle": "2025-10-09T15:21:01.498775Z", - "shell.execute_reply": "2025-10-09T15:21:01.498442Z" - } - }, - "outputs": [], - "source": [ - "molecule.assign_partial_charges?" - ] - }, { "cell_type": "markdown", "id": "7fee80b9", "metadata": {}, "source": [ - "First, let's assign charges the traditional way with AM1-BCC and check how long this takes..." + "First, let's parameterise with Sage 2.2.1, which uses the traditional AM1-BCC model, and check how long this takes..." ] }, { @@ -1074,10 +1122,8 @@ "outputs": [], "source": [ "%%time\n", - "molecule_am1bcc = Molecule(molecule)\n", - "molecule_am1bcc.assign_partial_charges(\n", - " partial_charge_method=\"am1bcc\",\n", - ")" + "sage221 = ForceField(\"openff-2.2.1.offxml\")\n", + "interchange_sage221 = Interchange.from_smirnoff(force_field=sage221, topology=molecule.to_topology())" ] }, { @@ -1085,7 +1131,9 @@ "id": "d751d584", "metadata": {}, "source": [ - "Now, let's try AshGC" + "Note that repeating these cells will show much faster assignment as partial charges are cached for a given molecule and charge method.\n", + "\n", + "Now, let's try Sage 2.3.0, which uses AshGC charges:" ] }, { @@ -1103,54 +1151,28 @@ "outputs": [], "source": [ "%%time\n", - "molecule_ashgc = Molecule(molecule)\n", - "molecule_ashgc.assign_partial_charges(\n", - " partial_charge_method=ASH_GC_MODEL,\n", - " toolkit_registry=NAGLToolkitWrapper(),\n", - ")" + "sage230 = ForceField(\"openff-2.3.0.offxml\")\n", + "interchange_sage230 = Interchange.from_smirnoff(force_field=sage230, topology=molecule.to_topology())" ] }, { "cell_type": "markdown", - "id": "39515e35", + "id": "e18fe6a2", "metadata": {}, "source": [ - "Finally, let's create an `Interchange` with our AshGC charges, making sure to specify `charge_from_molecules` so that we don't replace them with `AM1BCC` charges:" + "
\n", + " ✏️ Exercise: Compare the charges obtained with AM1-BCC and AshGC by inspecting the electrostatics collection of each Interchange object (see Section 3). How big are these differences on average? What is the largest difference? Which atom are these on? The np.max function may be useful.\n", + "
" ] }, { "cell_type": "code", "execution_count": null, - "id": "a81cf535", - "metadata": { - "execution": { - "iopub.execute_input": "2025-10-09T15:21:30.018680Z", - "iopub.status.busy": "2025-10-09T15:21:30.018523Z", - "iopub.status.idle": "2025-10-09T15:21:30.070775Z", - "shell.execute_reply": "2025-10-09T15:21:30.070474Z" - } - }, - "outputs": [], - "source": [ - "# normally when we call `ForceField.create_interchange` or `ForceField.create_openmm_system`, the toolkit will call\n", - "# AMBERTools or OEChem to assign partial charges, since that's what's in the force field file. A future OpenFF release\n", - "# which uses NAGL for charge assignment will encode this instruction in the force field file itself, but until that we\n", - "# can use the `charge_from_molecules` argument to tell it to use the charges that we just assigned# for more, see:\n", - "# https://docs.openforcefield.org/projects/toolkit/en/stable/api/generated/openff.toolkit.typing.engines.smirnoff.ForceField.html#openff.toolkit.typing.engines.smirnoff.ForceField.create_openmm_system\n", - "interchange = sage.create_interchange(\n", - " molecule_ashgc.to_topology(),\n", - " charge_from_molecules=[molecule_ashgc],\n", - ")" - ] - }, - { - "cell_type": "markdown", - "id": "e18fe6a2", + "id": "9791759e", "metadata": {}, + "outputs": [], "source": [ - "
\n", - " ✏️ Exercise: Compare the charges obtained with AM1-BCC and AshGC by looking at the Molecule.partial_charges attribute. How big are these differences on average? What is the largest difference? Which atom are these on? The np.max function may be useful.\n", - "
" + "interchange_sage230[\"Electrostatics\"].key_map" ] }, { @@ -1204,8 +1226,8 @@ "
\n", " ✏️ Extra Exercises: Based on the above tutorials, can you:\n", "
    \n", - "
  • Generate several conformers for one of your MCL-1 ligands and compute their relative energies using OpenFF 2.2.1?
    • \n", - "
  • Modify OpenFF 2.2.1 to change some of the parameters applied to one of your MCL-1 ligands? Minimise the ligand with this new force field and see how your changes influence the conformation.
    • \n", + "
  • Generate several conformers for one of your MCL-1 ligands and compute their relative energies using OpenFF 2.3.0?
    • \n", + "
  • Modify OpenFF 2.3.0 to change some of the parameters applied to one of your MCL-1 ligands? Minimise the ligand with this new force field and see how your changes influence the conformation.
    • \n", "
  • Analyse which parameters are shared and which are only applied to one or few molecules for a set of MCL-1 ligands?
    • \n", "
\n", "
" @@ -1228,7 +1250,7 @@ "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", - "version": "3.12.11" + "version": "3.13.15" }, "widgets": { "application/vnd.jupyter.widget-state+json": { From 896654dfdf0d15e7eabbdfe96567d1f2918cd024 Mon Sep 17 00:00:00 2001 From: Finlay Clark Date: Wed, 2 Sep 2026 13:52:27 +0100 Subject: [PATCH 2/9] Update notebooks with solutions and second notebook --- ...gand_complex_parameterisation_and_md.ipynb | 28 +- .../small_molecule_parameterisation.ipynb | 8 +- ...gand_complex_parameterisation_and_md.ipynb | 22 +- .../small_molecule_parameterisation.ipynb | 314 +++++++++--------- 4 files changed, 183 insertions(+), 189 deletions(-) diff --git a/notebooks/protein_ligand_complex_parameterisation_and_md.ipynb b/notebooks/protein_ligand_complex_parameterisation_and_md.ipynb index 084e5d5..48b5711 100644 --- a/notebooks/protein_ligand_complex_parameterisation_and_md.ipynb +++ b/notebooks/protein_ligand_complex_parameterisation_and_md.ipynb @@ -135,7 +135,7 @@ "\n", "## 2. OpenFF Toolkit Allows Us to Assemble the Topology\n", "\n", - "Conceptually, this step involves putting together the positions of all of the components of the system. We'll create a [`Topology`] to keep track of the contents of our system. As discussed in this morning's session, `Topology` represents a collection of molecules; it doesn't have any association with any force field parameters.\n", + "Conceptually, this step involves putting together the positions of all of the components of the system. We'll create a [`Topology`] to keep track of the contents of our system. As discussed in the previous notebook, `Topology` represents a collection of molecules; it doesn't have any association with any force field parameters.\n", "\n", "[`Topology`]: https://docs.openforcefield.org/projects/toolkit/en/stable/api/generated/openff.toolkit.topology.Topology.html" ] @@ -298,7 +298,7 @@ { "cell_type": "code", "execution_count": null, - "id": "147d6ec7", + "id": "5c71fa50", "metadata": { "tags": [ "placeholder" @@ -316,14 +316,14 @@ "\n", "## 3. We Can Assemble a Combined `ForceField` and use this to Parameterise the Whole System\n", "\n", - "Now that we've prepared our coordinates, we should choose the force field. For now, we don't have any single SMIRNOFF force field that can handle both proteins and small molecules; the Rosemary 3.0.0 force field will support this, but it's not yet ready. As an alternative, we'll combine the AMBER-compatible [Sage] small molecule force field with the SMIRNOFF port of AMBER ff14SB. Note that Sage also includes the TIP3P water model, which is appropriate for AMBER ff14SB too.\n", + "Now that we've prepared our coordinates, we should choose the force field. For now, we don't have any single SMIRNOFF force field that can handle both proteins and small molecules. The Rosemary line of force fields (starting with `openff-3.0.0.offxml`) is intended to do exactly this. A pre-release version, [`openff_no_water-3.0.0-alpha0.offxml`](https://github.com/openforcefield/openff-forcefields/blob/main/openforcefields/offxml/openff_no_water-3.0.0-alpha0.offxml), is already available for testing (if you use it, see [the release notes](https://github.com/openforcefield/openff-forcefields/releases/tag/2025.10.1)). There is no specific release date planned for the first full version, but it may be available later in 2026. As an alternative, we'll combine the AMBER-compatible [Sage] small molecule force field with the SMIRNOFF port of AMBER ff14SB. Note that Sage also includes the TIP3P water model, which is appropriate for AMBER ff14SB too.\n", "\n", "When we combine multiple SMIRNOFF force fields into one, we provide them in an order from general to specific. Sage includes parameters that could be applied to a protein, but they're general across all molecules; ff14SB's parameters are specific to proteins. Since the Toolkit always applies the last parameters that match a moiety, this order makes sure the right parameters get assigned.\n", "\n", "[Sage]: https://openforcefield.org/force-fields/force-fields/#sage\n", "\n", "
\n", - "⚠️ Warning: If your small molecule has an amino acid substructure in it, the specific patterns in the ff14SB force field will override the general ones from openff-2.2.1.offxml. This is the SMIRNOFF format being applied correctly, but some users may find this surprising, especially since terminal caps like ACE and NME are relatively small substructures and will sometimes appear in ligands.\n", + "⚠️ Warning: If your small molecule has an amino acid substructure in it, the specific patterns in the ff14SB force field will override the general ones from openff-2.3.0.offxml. This is the SMIRNOFF format being applied correctly, but some users may find this surprising, especially since terminal caps like ACE and NME are relatively small substructures and will sometimes appear in ligands.\n", "
\n" ] }, @@ -343,16 +343,16 @@ "from openff.toolkit import ForceField\n", "\n", "# Assemble the combined force field\n", - "sage_ff14sb = ForceField(\"openff-2.2.1.offxml\", \"ff14sb_off_impropers_0.0.3.offxml\")" + "sage_ff14sb = ForceField(\"openff-2.3.0.offxml\", \"ff14sb_off_impropers_0.0.3.offxml\")" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ - "We now have a `Topology`, which stores the chemical information of the system, and a `ForceField`, which maps chemistry to force field parameters. To parametrize the system, we combine these two objects into an [`Interchange`], as discussed in this morning's session.\n", + "We now have a `Topology`, which stores the chemical information of the system, and a `ForceField`, which maps chemistry to force field parameters. To parametrize the system, we combine these two objects into an [`Interchange`], as discussed in the previous notebook.\n", "\n", - "An `Interchange` represents a completely parameterised molecular mechanics system. Partial charges are computed here according to the instructions in the force field, and this is where virtual sites required by the force field will be introduced. This all happens behind the scenes; all we have to do is combine an abstract chemical description with a force field. This makes it easy to change water model or force field, as the chemistry being modelled is completely independent of the model itself.\n", + "An `Interchange` represents a completely parameterised molecular mechanics system. Partial charges are computed here according to the instructions in the force field: Sage 2.3.0 assigns the ligand's charges with the AshGC graph neural network model (see section 4 of the previous notebook), while ff14SB supplies library charges for the protein, water, and ions. This is also where virtual sites required by the force field will be introduced. This all happens behind the scenes; all we have to do is combine an abstract chemical description with a force field. This makes it easy to change water model or force field, as the chemistry being modelled is completely independent of the model itself.\n", "\n", "[`Interchange`]: https://docs.openforcefield.org/projects/interchange/en/stable/_autosummary/openff.interchange.components.interchange.Interchange.html" ] @@ -377,14 +377,14 @@ "cell_type": "markdown", "metadata": {}, "source": [ - "*(This should take about a minute, largely because of the complexity of the AMBER protein force field port. In the future, this should be faster.)*" + "*(This should take well under a minute, most of it spent on the chemical perception required by the AMBER protein force field port. It used to be considerably slower: assigning the ligand's AM1-BCC charges was the bottleneck, and AshGC has removed it.)*" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ - "While that runs, let's recap. We've constructed a `Topology` out of a number of `Molecule` objects, each of which represents a particular chemical independent of any model details. The `Topology` then represents an entire chemical system, which in theory could be modelled in any number of ways. Our `Topology` also includes atom positions and box vectors, but if we thought that was too concrete for our use case we could leave them out and add them after parameterisation.\n", + "Before we simulate, let's recap. We've constructed a `Topology` out of a number of `Molecule` objects, each of which represents a particular chemical independent of any model details. The `Topology` then represents an entire chemical system, which in theory could be modelled in any number of ways. Our `Topology` also includes atom positions and box vectors, but if we thought that was too concrete for our use case we could leave them out and add them after parameterisation.\n", "\n", "Separately, we've constructed a `ForceField` by combining a general SMIRNOFF force field with a protein-specific SMIRNOFF force field. A SMIRNOFF force field is a bunch of rules for applying force field parameters to chemicals via SMARTS patterns. The force field includes everything needed to compute an energy: parameters, charges, functional forms, non-bonded methods and cutoffs, virtual sites, and so on.\n", "\n", @@ -462,8 +462,8 @@ "FRICTION_COEFFICIENT = 1 / openmm.unit.picosecond\n", "TIMESTEP = 0.002 * openmm.unit.picoseconds\n", "\n", - "# Construct and configure a LangevinMiddleIntegrator at 300 K with an appropriate friction constant and time-step\n", - "integrator = openmm.LangevinMiddleIntegrator(\n", + "# Construct and configure a LangevinIntegrator at 300 K with an appropriate friction constant and time-step\n", + "integrator = openmm.LangevinIntegrator(\n", " TEMPERATURE,\n", " FRICTION_COEFFICIENT,\n", " TIMESTEP,\n", @@ -711,7 +711,7 @@ { "cell_type": "code", "execution_count": null, - "id": "c1992869", + "id": "4a9349fa", "metadata": { "tags": [ "placeholder" @@ -939,7 +939,7 @@ { "cell_type": "code", "execution_count": null, - "id": "7e6c7e5e", + "id": "145abf20", "metadata": { "tags": [ "placeholder" @@ -955,7 +955,7 @@ "metadata": {}, "source": [ "
\n", - " ✏️ Exercise: Repeat this entire notebook using a ligand from docked to MCL-1 during this morning's session. (Hint: You'll need to convert the pdbqt files to sdf files using obabel, adding protons as appropriate for pH 7. This will look something like obabel docked_ligand.pdbqt -opdb | obabel -ipdb -osdf -p 7.0 -O docked_ligand.sdf. Make sure to use the docked coordinates! An example docked pdbqt file is provided at ../structures/docked_ligand.pdbqt) Is the binding pose stable? Are similar interactions formed by the docked ligand and the crystallographic ligand? Which do you think is likely to bind more strongly? What would be required to answer these questions robustly?\n", + " ✏️ Exercise: Repeat this entire notebook using a ligand docked to MCL-1. (Hint: You'll need to convert the pdbqt files to sdf files using obabel, adding protons as appropriate for pH 7. This will look something like obabel docked_ligand.pdbqt -opdb | obabel -ipdb -osdf -p 7.0 -O docked_ligand.sdf. Make sure to use the docked coordinates! An example docked pdbqt file is provided at ../structures/docked_ligand.pdbqt) Is the binding pose stable? Are similar interactions formed by the docked ligand and the crystallographic ligand? Which do you think is likely to bind more strongly? What would be required to answer these questions robustly?\n", "
" ] }, diff --git a/notebooks/small_molecule_parameterisation.ipynb b/notebooks/small_molecule_parameterisation.ipynb index cc08160..257b044 100644 --- a/notebooks/small_molecule_parameterisation.ipynb +++ b/notebooks/small_molecule_parameterisation.ipynb @@ -512,7 +512,7 @@ { "cell_type": "code", "execution_count": null, - "id": "1e1444fc", + "id": "8fcccecc", "metadata": { "tags": [ "placeholder" @@ -829,7 +829,7 @@ { "cell_type": "code", "execution_count": null, - "id": "e794a5c2", + "id": "df6bf585", "metadata": { "tags": [ "placeholder" @@ -853,7 +853,7 @@ { "cell_type": "code", "execution_count": null, - "id": "86804051", + "id": "4d49b1fd", "metadata": { "tags": [ "placeholder" @@ -1178,7 +1178,7 @@ { "cell_type": "code", "execution_count": null, - "id": "a7c3e001", + "id": "c6a0b74d", "metadata": { "tags": [ "placeholder" diff --git a/notebooks_with_solutions/protein_ligand_complex_parameterisation_and_md.ipynb b/notebooks_with_solutions/protein_ligand_complex_parameterisation_and_md.ipynb index 42107d1..894b967 100644 --- a/notebooks_with_solutions/protein_ligand_complex_parameterisation_and_md.ipynb +++ b/notebooks_with_solutions/protein_ligand_complex_parameterisation_and_md.ipynb @@ -135,7 +135,7 @@ "\n", "## 2. OpenFF Toolkit Allows Us to Assemble the Topology\n", "\n", - "Conceptually, this step involves putting together the positions of all of the components of the system. We'll create a [`Topology`] to keep track of the contents of our system. As discussed in this morning's session, `Topology` represents a collection of molecules; it doesn't have any association with any force field parameters.\n", + "Conceptually, this step involves putting together the positions of all of the components of the system. We'll create a [`Topology`] to keep track of the contents of our system. As discussed in the previous notebook, `Topology` represents a collection of molecules; it doesn't have any association with any force field parameters.\n", "\n", "[`Topology`]: https://docs.openforcefield.org/projects/toolkit/en/stable/api/generated/openff.toolkit.topology.Topology.html" ] @@ -322,14 +322,14 @@ "\n", "## 3. We Can Assemble a Combined `ForceField` and use this to Parameterise the Whole System\n", "\n", - "Now that we've prepared our coordinates, we should choose the force field. For now, we don't have any single SMIRNOFF force field that can handle both proteins and small molecules; the Rosemary 3.0.0 force field will support this, but it's not yet ready. As an alternative, we'll combine the AMBER-compatible [Sage] small molecule force field with the SMIRNOFF port of AMBER ff14SB. Note that Sage also includes the TIP3P water model, which is appropriate for AMBER ff14SB too.\n", + "Now that we've prepared our coordinates, we should choose the force field. For now, we don't have any single SMIRNOFF force field that can handle both proteins and small molecules. The Rosemary line of force fields (starting with `openff-3.0.0.offxml`) is intended to do exactly this. A pre-release version, [`openff_no_water-3.0.0-alpha0.offxml`](https://github.com/openforcefield/openff-forcefields/blob/main/openforcefields/offxml/openff_no_water-3.0.0-alpha0.offxml), is already available for testing (if you use it, see [the release notes](https://github.com/openforcefield/openff-forcefields/releases/tag/2025.10.1)). There is no specific release date planned for the first full version, but it may be available later in 2026. As an alternative, we'll combine the AMBER-compatible [Sage] small molecule force field with the SMIRNOFF port of AMBER ff14SB. Note that Sage also includes the TIP3P water model, which is appropriate for AMBER ff14SB too.\n", "\n", "When we combine multiple SMIRNOFF force fields into one, we provide them in an order from general to specific. Sage includes parameters that could be applied to a protein, but they're general across all molecules; ff14SB's parameters are specific to proteins. Since the Toolkit always applies the last parameters that match a moiety, this order makes sure the right parameters get assigned.\n", "\n", "[Sage]: https://openforcefield.org/force-fields/force-fields/#sage\n", "\n", "
\n", - "⚠️ Warning: If your small molecule has an amino acid substructure in it, the specific patterns in the ff14SB force field will override the general ones from openff-2.2.1.offxml. This is the SMIRNOFF format being applied correctly, but some users may find this surprising, especially since terminal caps like ACE and NME are relatively small substructures and will sometimes appear in ligands.\n", + "⚠️ Warning: If your small molecule has an amino acid substructure in it, the specific patterns in the ff14SB force field will override the general ones from openff-2.3.0.offxml. This is the SMIRNOFF format being applied correctly, but some users may find this surprising, especially since terminal caps like ACE and NME are relatively small substructures and will sometimes appear in ligands.\n", "
\n" ] }, @@ -349,16 +349,16 @@ "from openff.toolkit import ForceField\n", "\n", "# Assemble the combined force field\n", - "sage_ff14sb = ForceField(\"openff-2.2.1.offxml\", \"ff14sb_off_impropers_0.0.3.offxml\")" + "sage_ff14sb = ForceField(\"openff-2.3.0.offxml\", \"ff14sb_off_impropers_0.0.3.offxml\")" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ - "We now have a `Topology`, which stores the chemical information of the system, and a `ForceField`, which maps chemistry to force field parameters. To parametrize the system, we combine these two objects into an [`Interchange`], as discussed in this morning's session.\n", + "We now have a `Topology`, which stores the chemical information of the system, and a `ForceField`, which maps chemistry to force field parameters. To parametrize the system, we combine these two objects into an [`Interchange`], as discussed in the previous notebook.\n", "\n", - "An `Interchange` represents a completely parameterised molecular mechanics system. Partial charges are computed here according to the instructions in the force field, and this is where virtual sites required by the force field will be introduced. This all happens behind the scenes; all we have to do is combine an abstract chemical description with a force field. This makes it easy to change water model or force field, as the chemistry being modelled is completely independent of the model itself.\n", + "An `Interchange` represents a completely parameterised molecular mechanics system. Partial charges are computed here according to the instructions in the force field: Sage 2.3.0 assigns the ligand's charges with the AshGC graph neural network model (see section 4 of the previous notebook), while ff14SB supplies library charges for the protein, water, and ions. This is also where virtual sites required by the force field will be introduced. This all happens behind the scenes; all we have to do is combine an abstract chemical description with a force field. This makes it easy to change water model or force field, as the chemistry being modelled is completely independent of the model itself.\n", "\n", "[`Interchange`]: https://docs.openforcefield.org/projects/interchange/en/stable/_autosummary/openff.interchange.components.interchange.Interchange.html" ] @@ -383,14 +383,14 @@ "cell_type": "markdown", "metadata": {}, "source": [ - "*(This should take about a minute, largely because of the complexity of the AMBER protein force field port. In the future, this should be faster.)*" + "*(This should take well under a minute, most of it spent on the chemical perception required by the AMBER protein force field port. It used to be considerably slower: assigning the ligand's AM1-BCC charges was the bottleneck, and AshGC has removed it.)*" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ - "While that runs, let's recap. We've constructed a `Topology` out of a number of `Molecule` objects, each of which represents a particular chemical independent of any model details. The `Topology` then represents an entire chemical system, which in theory could be modelled in any number of ways. Our `Topology` also includes atom positions and box vectors, but if we thought that was too concrete for our use case we could leave them out and add them after parameterisation.\n", + "Before we simulate, let's recap. We've constructed a `Topology` out of a number of `Molecule` objects, each of which represents a particular chemical independent of any model details. The `Topology` then represents an entire chemical system, which in theory could be modelled in any number of ways. Our `Topology` also includes atom positions and box vectors, but if we thought that was too concrete for our use case we could leave them out and add them after parameterisation.\n", "\n", "Separately, we've constructed a `ForceField` by combining a general SMIRNOFF force field with a protein-specific SMIRNOFF force field. A SMIRNOFF force field is a bunch of rules for applying force field parameters to chemicals via SMARTS patterns. The force field includes everything needed to compute an energy: parameters, charges, functional forms, non-bonded methods and cutoffs, virtual sites, and so on.\n", "\n", @@ -468,8 +468,8 @@ "FRICTION_COEFFICIENT = 1 / openmm.unit.picosecond\n", "TIMESTEP = 0.002 * openmm.unit.picoseconds\n", "\n", - "# Construct and configure a LangevinMiddleIntegrator at 300 K with an appropriate friction constant and time-step\n", - "integrator = openmm.LangevinMiddleIntegrator(\n", + "# Construct and configure a LangevinIntegrator at 300 K with an appropriate friction constant and time-step\n", + "integrator = openmm.LangevinIntegrator(\n", " TEMPERATURE,\n", " FRICTION_COEFFICIENT,\n", " TIMESTEP,\n", @@ -998,7 +998,7 @@ "metadata": {}, "source": [ "
\n", - " ✏️ Exercise: Repeat this entire notebook using a ligand from docked to MCL-1 during this morning's session. (Hint: You'll need to convert the pdbqt files to sdf files using obabel, adding protons as appropriate for pH 7. This will look something like obabel docked_ligand.pdbqt -opdb | obabel -ipdb -osdf -p 7.0 -O docked_ligand.sdf. Make sure to use the docked coordinates! An example docked pdbqt file is provided at ../structures/docked_ligand.pdbqt) Is the binding pose stable? Are similar interactions formed by the docked ligand and the crystallographic ligand? Which do you think is likely to bind more strongly? What would be required to answer these questions robustly?\n", + " ✏️ Exercise: Repeat this entire notebook using a ligand docked to MCL-1. (Hint: You'll need to convert the pdbqt files to sdf files using obabel, adding protons as appropriate for pH 7. This will look something like obabel docked_ligand.pdbqt -opdb | obabel -ipdb -osdf -p 7.0 -O docked_ligand.sdf. Make sure to use the docked coordinates! An example docked pdbqt file is provided at ../structures/docked_ligand.pdbqt) Is the binding pose stable? Are similar interactions formed by the docked ligand and the crystallographic ligand? Which do you think is likely to bind more strongly? What would be required to answer these questions robustly?\n", "
" ] }, diff --git a/notebooks_with_solutions/small_molecule_parameterisation.ipynb b/notebooks_with_solutions/small_molecule_parameterisation.ipynb index 90a77cd..ca6cad3 100644 --- a/notebooks_with_solutions/small_molecule_parameterisation.ipynb +++ b/notebooks_with_solutions/small_molecule_parameterisation.ipynb @@ -18,9 +18,10 @@ "\n", "| Action | Software|\n", "|--|--|\n", + "| [Go from smiles to simulation in few lines of code](#showcase) | OpenFF Toolkit, OpenFF Interchange, OpenMM\n", "| [Load and inspect a force field](#loading_ff) | OpenFF Toolkit\n", "| [Create a representation of your chemical system](#topology) | OpenFF Toolkit\n", - "| [Parameterise your system and run a quick simulation](#interchange) | OpenFF Interchange, OpenMM\n", + "| [Parameterise your system and run a quick simulation in water](#interchange) | OpenFF Interchange, OpenMM\n", "| [Rapidly assign partial charges with a graph neural network model](#gnn_charges) | OpenFF Toolkit, OpenFF NAGL Models\n", "| [Review what you've learnt](#summary) | \n", "| [Check out other OpenFF tutorials](#further_materials) | \n", @@ -40,6 +41,91 @@ "\n" ] }, + { + "cell_type": "markdown", + "id": "81bb0b4f", + "metadata": {}, + "source": [ + "\n", + "## 0. You can go from SMILES to simulation in a few lines of code" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "cef30154", + "metadata": {}, + "outputs": [], + "source": [ + "# Go from SMILES -> simulation input with OpenFF\n", + "from openff.toolkit import ForceField, Molecule, Topology\n", + "\n", + "molecule = Molecule.from_smiles(\"CC(=O)Nc1ccc(cc1)O\")\n", + "molecule.generate_conformers(n_conformers=1)\n", + "topology = Topology.from_molecules([molecule])\n", + "\n", + "force_field = ForceField(\"openff-2.3.0.offxml\")\n", + "interchange = force_field.create_interchange(topology)\n", + "interchange.minimize()\n", + "\n", + "openmm_system = interchange.to_openmm_system()\n", + "openmm_topology = interchange.to_openmm_topology()\n", + "openmm_positions = interchange.positions.to_openmm()" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "11ec8a77", + "metadata": {}, + "outputs": [], + "source": [ + "# Run the simulation with OpenMM\n", + "import openmm\n", + "\n", + "temperature = 298.15 * openmm.unit.kelvin\n", + "friction_coefficient = 1.0 / openmm.unit.picosecond\n", + "step_size = 2.0 * openmm.unit.femtosecond\n", + "\n", + "simulation = openmm.app.Simulation(\n", + " openmm_topology,\n", + " openmm_system,\n", + " openmm.LangevinIntegrator(temperature, friction_coefficient, step_size),\n", + ")\n", + "simulation.context.setPositions(openmm_positions)\n", + "simulation.context.setVelocitiesToTemperature(simulation.integrator.getTemperature())\n", + "\n", + "simulation.reporters.append(\n", + " openmm.app.DCDReporter(file=\"trajectory_showcase.dcd\", reportInterval=100)\n", + ")\n", + "simulation.step(10000)" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "a2983975", + "metadata": {}, + "outputs": [], + "source": [ + "# Load the trajectory with MDAnalysis and visualise with nglview\n", + "import MDAnalysis as mda\n", + "import nglview\n", + "\n", + "u = mda.Universe(openmm_topology, \"trajectory_showcase.dcd\")\n", + "\n", + "view = nglview.show_mdanalysis(u)\n", + "view" + ] + }, + { + "cell_type": "markdown", + "id": "e110064b", + "metadata": {}, + "source": [ + "That's it! You've run a vacuum simulation for paracetamol. Below and in the next notebook, we'll go into more detail on each of the steps in the OpenFF cell and show how you can set up more complex systems, but this is mainly for your understanding and you rarely need much more code than shown above." + ] + }, { "cell_type": "markdown", "id": "32c4b5c9", @@ -48,7 +134,11 @@ "\n", "## 1. Force fields are specified in `.offxml` files and can be loaded with the `ForceField` class\n", "\n", - "OpenFF's force fields use the The SMIRKS Native Open Force Field (SMIRNOFF) [specification](https://openforcefield.github.io/standards/standards/smirnoff/) and are conventionally encoded in `.offxml` files. The spec fully describes the contents of a SMIRNOFF force field, how parameters should be applied, and several other important usage details. You could implement a SMIRNOFF engine in your own code, but conveniently the OpenFF Toolkit already provides this and a handful of utilities. Let's load up the latest OpenFF small molecule force field, OpenFF 2.2.1, and inspect its contents! This force field shares the code name \"Sage\" with all other force fields with the same major version number (2.x.x)." + "Let's dive into the details of what went on above. Here's a summary of how data flows through a workflow utilising OpenFF tools -- the OpenFF toolkit allows you to create `Molecule` and `ForceField` objects, which get combined into an `Interchange` object, which contains all the information needed to start a simulation. From there, you can create input for the simulation engine of your choice:\n", + "\n", + "\"Description\n", + "\n", + "Let's start with the `.offxml` force field file. OpenFF's force fields use the The SMIRKS Native Open Force Field (SMIRNOFF) [specification](https://openforcefield.github.io/standards/standards/smirnoff/). The spec fully describes the contents of a SMIRNOFF force field, how parameters should be applied, and several other important usage details. You could implement a SMIRNOFF engine in your own code, but conveniently the OpenFF Toolkit already provides this and a handful of utilities. Let's load up the latest OpenFF small molecule force field, OpenFF 2.3,0, and inspect its contents! This force field shares the code name \"Sage\" with all other force fields with the same major version number (2.x.x)." ] }, { @@ -67,7 +157,7 @@ "source": [ "from openff.toolkit import ForceField\n", "\n", - "sage = ForceField(\"openff-2.2.1.offxml\")\n", + "sage = ForceField(\"openff-2.3.0.offxml\")\n", "sage" ] }, @@ -76,7 +166,7 @@ "id": "52c754d3", "metadata": {}, "source": [ - "If you'd like to see the raw file on disk that's being parsed, [here's the file on GitHub](https://github.com/openforcefield/openff-forcefields/blob/main/openforcefields/offxml/openff-2.2.1.offxml)." + "If you'd like to see the raw file on disk that's being parsed, [here's the file on GitHub](https://github.com/openforcefield/openff-forcefields/blob/main/openforcefields/offxml/openff-2.3.0.offxml)." ] }, { @@ -170,7 +260,7 @@ "\n", "The toolkit uses these SMIRKS patterns and direct chemical perception to assign parameters to particular atoms (or bonds, angles, etc.).\n", "\n", - "We'll use OpenFF 2.2.1 for the remainder of this tutorial, but you can learn more about this and other SMIRNOFF force fields below:\n", + "We'll use OpenFF 2.3.0 for the remainder of this tutorial. This is OpenFF's latest small molecule force field and is a leading open-source small molecule force field which [performs comparably to other open-source force fields](https://doi-org.libproxy.ncl.ac.uk/10.1021/acs.jctc.3c00039). You can learn more about this and other SMIRNOFF force fields below:\n", "
\n", " Click here to learn about available and planned SMIRNOFF force fields\n", "\n", @@ -192,7 +282,7 @@ "\n", "The Sage line of force fields (`openff-2.y.z.offxml`) continued the process of fitting to more (and more diverse) QM datasets, but also included a re-fit of the Lennard-Jones parameters. Small molecule geometries and energies [improved, in general,](https://openforcefield.org/community/news/general/sage2.0.0-release/) significantly over Parsley. These improvements notably transferred to protein-ligand binding free energies despite Sage not being specifically fit to them. For more, see the [associated paper](https://pubs.acs.org/doi/10.1021/acs.jctc.3c00039).\n", "\n", - "[Subsequent releases](https://github.com/openforcefield/openff-forcefields/releases) used different fitting procedures and tweaks to parameter typing to improve performance and address issues with several specific chemistries. The latest release, **Sage 2.2.1 (`openff-2.2.1.offxml`) is the recommended force field for small molecule studies.**\n", + "[Subsequent releases](https://github.com/openforcefield/openff-forcefields/releases) used different fitting procedures and tweaks to parameter typing to improve performance and address issues with several specific chemistries. Notably, Sage 2.3.0 includes fast graph neural network charge assignment with Ash-GC, which is discussed later in this notebook. This charge model is trained to reproduce AM1-BCC charges without the typical $O(N^3)$ scaling, making it suitable for large (>> 100 atoms) molecules). The latest release, **Sage 2.3.0 (`openff-2.3.0.offxml`) is the recommended force field for small molecule studies.**\n", "\n", "## Ports\n", "\n", @@ -212,6 +302,13 @@ "\n", "OpenFF, in collaboration with Dave Cerutti of the Amber community, created a port of [ff14SB](https://pubs.acs.org/doi/10.1021/acs.jctc.5b00255), a popular Amber protein force field. There are some small numerical differences with how improper torsions are evaluated, but all other terms reproduce a canonical Amber source to high accuracy. **This is the only protein force field currently in SMIRONOFF (`.offxml`) format** and therefore the current recommendation for use with proteins. Primarily for technical reasons, porting other Amber force fields is not planned.\n", "\n", + "## Rosemary Alpha\n", + "\n", + "A future line of force fields from OpenFF (code name \"Rosemary\", starting with `openff-3.0.0.offxml`) is intended to handle small molecules and biopolymers in a _self-consistent_ manner. This is exciting as it will streamline simulations of proteins with non-cannonical amino acids! See the workshop [Simulating Post-Translationally Modified Proteins with the OpenFF Rosemary Alpha](https://github.com/openforcefield/2026-virtual-workshops/blob/main/ptm/ptm-workshop.ipynb). The first release will handle proteins, but future versions may cover nucleic acids. The performance, depending on the metrics used, is hoped to be comparable with existing Amber-family protein force fields. \n", + "\n", + "An pre-release version of Rosemary is available for testing as [`openff_no_water-3.0.0-alpha0.offxml`](https://github.com/openforcefield/openff-forcefields/blob/main/openforcefields/offxml/openff_no_water-3.0.0-alpha0.offxml). If you use it, see [the release notes](https://github.com/openforcefield/openff-forcefields/releases/tag/2025.10.1). \n", + "\n", + "\n", "## Non-main-line force fields\n", "\n", "### SMIRNOFF plugins\n", @@ -225,15 +322,7 @@ "## From OpenFF\n", "\n", "### Rosemary\n", - "A future line of force fields from OpenFF (code name \"Rosemary\", starting with `openff-3.0.0.offxml`) is intended to handle small molecules and biopolymers in a _self-consistent_ manner. The first release is expected to handle proteins, but future versions may cover nucleic acids. The performance, depending on the metrics used, is hoped to be comparable with existing Amber-family protein force fields.\n", - "\n", - "There is no specific release date planned for Rosemary, but it may be available in 2026 (a beta release candidate may also be publically available prior to the full release).\n", - "\n", - "### Graph net charge assignment\n", - "\n", - "TODO: UPDAATE AND MENTION 2.3\n", - "\n", - "The Sage 2.3.0 release is expected imminently and will include graph-convolutional neutral network (GCNN)-based charge assignment using [NAGL](https://github.com/openforcefield/openff-nagl) by default. The charge model is trained to reproduce AM1-BCC charges without the typical $O(N^3)$ scaling, making it suitable for large (>> 100 atoms) molecules). The [second release candidate](https://github.com/openforcefield/openff-forcefields/blob/main/openforcefields/offxml/openff-2.3.0-rc2.offxml) (which may or may not become the final version) is already available for you to try!\n", + "There is no specific release date planned for the first full version of Rosemary, but it may be available in late 2026.\n", "\n", "### Virtual sites\n", "\n", @@ -253,7 +342,7 @@ "\n", "## 2. The `Topology` class represents a chemical system containing one or more `Molecule`s\n", "\n", - "Now we've loaded our desired force field (OpenFF 2.2.1), we need to specify the chemical system we want to assign force field parameters to (\"parameterise\"). Our system will be represented by a `Topology`, which we will build from one or more `Molecule`s. \n", + "Now we've loaded our desired force field (OpenFF 2.3.0), we need to specify the chemical system we want to assign force field parameters to (\"parameterise\"). Our system will be represented by a `Topology`, which we will build from one or more `Molecule`s. \n", "\n", "As a simple example, let's build a `Topology` containing an small molecule with some features which illustrate how parameters are applied according to SMIRKS matches. We'll use the crotonate anion, but you could draw any molecule you like and convert it to a SMILES string using tools like ChemDraw and [MolView](https://molview.org/)." ] @@ -416,7 +505,7 @@ "metadata": {}, "source": [ "
\n", - " ✏️ Exercise: Build a Topology containing an MCL-1 ligand. Create the Molecule from an SDF file (take a look at the docstring of Molecule to see how this can be done). Also, see Molecule cookbook for all the ways to make a Molecule. The crystallographic MCL-1 ligand from PDB ID 6o6f is provided at ../structures/606f_ligand.sdf. Note that you don't need to use get_data_file_path as we already know the path.\n", + " ✏️ Exercise: Build a Topology containing an MCL-1 ligand. Create the Molecule from an SDF file (take a look at the docstring of Molecule to see how this can be done, noting that you can just pass the sdf path given below and don't need the get_data_file_path function). Also, see Molecule cookbook for all the ways to make a Molecule. The crystallographic MCL-1 ligand from PDB ID 6o6f is provided at ../structures/606f_ligand.sdf. Note that you don't need to use get_data_file_path as we already know the path.\n", "
\n" ] }, @@ -425,12 +514,6 @@ "execution_count": null, "id": "4673ba20", "metadata": { - "execution": { - "iopub.execute_input": "2025-10-09T15:20:46.525200Z", - "iopub.status.busy": "2025-10-09T15:20:46.525111Z", - "iopub.status.idle": "2025-10-09T15:20:46.541355Z", - "shell.execute_reply": "2025-10-09T15:20:46.540944Z" - }, "tags": [ "solution" ] @@ -447,7 +530,7 @@ "id": "0c9e16d5", "metadata": {}, "source": [ - "We will cover creating a topology for a protein-ligand complex this afternoon." + "We will cover creating a topology for a protein-ligand complex in the next notebook." ] }, { @@ -460,17 +543,7 @@ "\n", "Now we've specified our force field and our chemical system using classes from the OpenFF Tools package (`ForceField`, `Molecule`, and `Topology`), and we want to apply our force field to our chemical topologies (parameterisation).\n", "\n", - "To do this, we'll use the `Interchange` class from the OpenFF Interchange package, which stores a fully-parameterised molecular system and provides methods to write out simulation-ready input files for a number of software packages. They key objective of Interchange is to provide an intermediate inspectable state after parameterisation and before conversion to an engine-specific format. For most users, an `Interchange` forms the bridge between the OpenFF ecosystem and their simulation software of choice. The current focus is applying SMIRNOFF force fields to chemical topologies and exporting the result to engines preferred by our users. In order of stability, OpenMM, GROMACS, Amber, and LAMMPS are supported. Future development may include support for CHARMM and other engines.\n", - "\n", - "Below is a summary of how data flows through a workflow utilising OpenFF tools, including where Interchange sits in the flow." - ] - }, - { - "cell_type": "markdown", - "id": "ddde1c2a", - "metadata": {}, - "source": [ - "\"Description" + "To do this, we'll use the `Interchange` class from the OpenFF Interchange package, which stores a fully-parameterised molecular system and provides methods to write out simulation-ready input files for a number of software packages. They key objective of Interchange is to provide an intermediate inspectable state after parameterisation and before conversion to an engine-specific format. For most users, an `Interchange` forms the bridge between the OpenFF ecosystem and their simulation software of choice. The current focus is applying SMIRNOFF force fields to chemical topologies and exporting the result to engines preferred by our users. In order of stability, OpenMM, GROMACS, Amber, and LAMMPS are supported. Future development may include support for CHARMM and other engines." ] }, { @@ -713,7 +786,7 @@ "id": "c1603b89", "metadata": {}, "source": [ - "We can see that the C=C bond (indices (1,2)) is associated with a potential key with the SMIRKS pattern `[#6X3:1]=[#6X3:2]` (specifying any two carbons bonded to 3 atoms connected by a double bond). Note that the (1,0) C-C bond is matched by the SMRIKS `[#6X3:1]-[#6X3:2]`, which specifies the atoms in the same way, showing that the parameters have been assigned by directly using information about the bond. This contrasts to traditional atom typing approaches, where information about the bond would be implicitly encoded in the atom types used to assign the parameters. Another example of this \"direct chemical perception\" is the assignment of the carboxylate carbon-oxygen bond parameters, which only match (triply-connected carbon) - (singly-connnected oxygen) bonds when the carbon is bonded to another singly-connected oxygen.\n", + "We can see that the C=C bond (indices (1,2)) is associated with a potential key with the SMIRKS pattern `[#6X3:1]=[#6X3:2]` (specifying any two carbons each bonded to 3 atoms and connected by a double bond). Note that the (1,0) C-C bond is matched by the SMRIKS `[#6X3:1]-[#6X3:2]`, which specifies the atoms in the same way, showing that the parameters have been assigned by directly using information about the bond. This contrasts to traditional atom typing approaches, where information about the bond would be implicitly encoded in the atom types used to assign the parameters. Another example of this \"direct chemical perception\" is the assignment of the carboxylate carbon-oxygen bond parameters, which only match (triply-connected carbon) - (singly-connnected oxygen) bonds when the carbon is bonded to another singly-connected oxygen.\n", "\n", "To see the actual parmeters specified for this bond, we can look up the `Potential` objects using the `PotentialKey`s." ] @@ -760,12 +833,6 @@ "execution_count": null, "id": "d11c4ced", "metadata": { - "execution": { - "iopub.execute_input": "2025-10-09T15:20:47.518198Z", - "iopub.status.busy": "2025-10-09T15:20:47.518102Z", - "iopub.status.idle": "2025-10-09T15:20:47.519617Z", - "shell.execute_reply": "2025-10-09T15:20:47.519362Z" - }, "tags": [ "solution" ] @@ -794,12 +861,6 @@ "execution_count": null, "id": "2fc49bda", "metadata": { - "execution": { - "iopub.execute_input": "2025-10-09T15:20:47.520438Z", - "iopub.status.busy": "2025-10-09T15:20:47.520341Z", - "iopub.status.idle": "2025-10-09T15:20:47.521771Z", - "shell.execute_reply": "2025-10-09T15:20:47.521555Z" - }, "tags": [ "solution" ] @@ -946,7 +1007,7 @@ "id": "212385e0", "metadata": {}, "source": [ - "Here, we'll export to OpenMM and run a short simulation directly from the noteboook. We can create an OpenMM `Simulation` object from the `Interchange` and run for a specified wall clock time using `runForClockTime` (the simluation time will depend on how quickly it runs on your machine). We keep the volume ($V$), number of particles ($N$), and average temperature ($T$) (using the LangevinMiddleIntegrator) constant and the simulation corresponds to the $NVT$ ensemble." + "Here, we'll export to OpenMM and run a short simulation directly from the noteboook. We can create an OpenMM `Simulation` object from the `Interchange` and run for a specified wall clock time using `runForClockTime` (the simluation time will depend on how quickly it runs on your machine). We keep the volume ($V$), number of particles ($N$), and average temperature ($T$) (using the LangevinIntegrator) constant and the simulation corresponds to the $NVT$ ensemble." ] }, { @@ -966,7 +1027,7 @@ "import openmm\n", "import openmm.unit\n", "from openff.interchange import Interchange\n", - "import mdtraj\n", + "import MDAnalysis as mda\n", "import nglview\n", "\n", "\n", @@ -976,7 +1037,7 @@ " trajectory_name: str = \"small_mol_solvated.dcd\",\n", "):\n", " simulation = interchange.to_openmm_simulation(\n", - " integrator=openmm.LangevinMiddleIntegrator(\n", + " integrator=openmm.LangevinIntegrator(\n", " 300 * openmm.unit.kelvin,\n", " 1 / openmm.unit.picosecond,\n", " 0.002 * openmm.unit.picoseconds,\n", @@ -994,12 +1055,10 @@ " topology: Topology, filename: str = \"small_mol_solvated.dcd\"\n", ") -> nglview.NGLWidget:\n", " \"\"\"Visualise a trajectory using nglview.\"\"\"\n", - " traj = mdtraj.load(\n", - " filename,\n", - " top=mdtraj.Topology.from_openmm(topology.to_openmm()),\n", - " )\n", "\n", - " view = nglview.show_mdtraj(traj)\n", + " u = mda.Universe(topology.to_openmm(), filename)\n", + "\n", + " view = nglview.show_mdanalysis(u)\n", " view.add_representation(\"licorice\", selection=\"water\")\n", "\n", " return view\n", @@ -1027,23 +1086,19 @@ "\n", "## 4. Graph Neural Networks Allow Fast Assignment of Partial Charges\n", "\n", - "You might notice that [Sage](https://github.com/openforcefield/openff-forcefields/blob/main/openforcefields/offxml/openff-2.2.1.offxml) doesn't contain tabulated charges for most atomic environments in the way it does for all other terms in the force field. Instead, it specifies:\n", + "You might notice that Sage force fields don't contain tabulated charges for most atomic environments in the way they do for all other terms in the force field.For example, [Sage 2.2.1](https://github.com/openforcefield/openff-forcefields/blob/main/openforcefields/offxml/openff-2.2.1.offxml) instead specifies:\n", "```\n", - "\n", + "\n", "```\n", - "which means that partial charges will be calculated using the common AM1-BCC method. Charges from a semi-empirical quantum chemistry calculation (Austin Model 1) are corrected (bond charge correction) to approximate charges obtained by fitting to the electrostatic potential at the HF/6-31G* level (see [Jakalian et al.](https://onlinelibrary.wiley.com/doi/10.1002/(SICI)1096-987X(20000130)21:2%3C132::AID-JCC5%3E3.0.CO;2-P)). Unfortunately, AM1-BCC scales 𝒪(N2) in the number of atoms N, making it prohibitively slow for large molecules and biopolymers.\n", + "which means that partial charges will be calculated using the common AM1-BCC method. Charges from a semi-empirical quantum chemistry calculation (Austin Model 1) are corrected (bond charge correction) to approximate charges obtained by fitting to the electrostatic potential at the HF/6-31G* level (see [Jakalian et al.](https://onlinelibrary.wiley.com/doi/10.1002/(SICI)1096-987X(20000130)21:2%3C132::AID-JCC5%3E3.0.CO;2-P)). Unfortunately, parameterisation with AM1-BCC using OpenEye or AmberTools scales 𝒪(N2) in the number of atoms N, making it prohibitively slow for large molecules and biopolymers.\n", "\n", - "Methods which assign partial charges using graph neural networks offer rapid assignment with better scaling. They also offer the possibility of going beyond traditionally affordable QM levels of theory by training to quickly reproduce charges from expensive calculations. For example, [EspalomaCharge](https://pubs.acs.org/doi/full/10.1021/acs.jpca.4c01287) is fit to AM1-BCC charges and offers 𝒪(N2) scaling, while [Adams et al.](https://chemrxiv.org/engage/chemrxiv/article-details/6839c94c3ba0887c33d2cd8e) trained models to reproduce atoms-in-molecules charges and electrostatic potentials obtained at a high level of theory. Here, we'll use OpenFF's [AshGC](https://zenodo.org/records/15770227/files/AshGC_methods_2025-06-30.pdf?download=1) model, which is trained to reproduce AM1-BCC charges." - ] - }, - { - "cell_type": "markdown", - "id": "ede41aef", - "metadata": {}, - "source": [ - "
\n", - " ⚠️ OpenFF 2.2.1 has not been explicitly trained and validated with AshGC charges. However, the 2.3.0 release will be, and is expected imminently. AshGC charges will be used as default and will be specified in the .offxml file, so there will be no need to call Molecule.assign_partial_charges as shown below.\n", - "
" + "Methods which assign partial charges using graph neural networks offer rapid assignment with better scaling. They also offer the possibility of going beyond traditionally affordable QM levels of theory by training to quickly reproduce charges from expensive calculations. For example, OpenFF's [AshGC](https://doi.org/10.1021/acs.jctc.6c00169) model is fit to AM1-BCC charges and offers 𝒪(N) scaling, while [Adams et al.](https://doi.org/10.1021/acs.jctc.5c01520) trained models to reproduce atoms-in-molecules charges and electrostatic potentials obtained at a high level of theory. AshGC is used by [Sage 2.3.0](https://github.com/openforcefield/openff-forcefields/blob/main/openforcefields/offxml/openff-2.3.0.offxml) -- if you inspect the file, you'll see:\n", + "```\n", + "\n", + "```\n", + "where \"openff-gnn-am1bcc-1.0.0.pt\" is the AshGC model.\n", + "\n", + "The GNN charge model is the main difference between Sage 2.2.1 and 2.3.0. Here, we'll compare the parameterisation speed and charges obtained using each force field." ] }, { @@ -1060,42 +1115,17 @@ }, "outputs": [], "source": [ - "from openff.toolkit import Molecule, ForceField\n", - "from openff.toolkit.utils.nagl_wrapper import NAGLToolkitWrapper\n", + "from openff.toolkit import Molecule\n", "\n", - "# Disable RDKit warnings to avoid misleading NAGL warnings\n", - "# (see https://github.com/openforcefield/openff-nagl/issues/198)\n", - "from rdkit import RDLogger\n", - "RDLogger.DisableLog('rdApp.*') \n", - "\n", - "# OpenFF NAGL store models as PyTorch files.\n", - "ASH_GC_MODEL = \"openff-gnn-am1bcc-0.1.0-rc.3.pt\"\n", "molecule = Molecule(\"../structures/6o6f_ligand.sdf\")" ] }, - { - "cell_type": "code", - "execution_count": null, - "id": "e0646ac8", - "metadata": { - "execution": { - "iopub.execute_input": "2025-10-09T15:21:01.492129Z", - "iopub.status.busy": "2025-10-09T15:21:01.491979Z", - "iopub.status.idle": "2025-10-09T15:21:01.498775Z", - "shell.execute_reply": "2025-10-09T15:21:01.498442Z" - } - }, - "outputs": [], - "source": [ - "molecule.assign_partial_charges?" - ] - }, { "cell_type": "markdown", "id": "7fee80b9", "metadata": {}, "source": [ - "First, let's assign charges the traditional way with AM1-BCC and check how long this takes..." + "First, let's parameterise with Sage 2.2.1, which uses the traditional AM1-BCC model, and check how long this takes..." ] }, { @@ -1113,10 +1143,8 @@ "outputs": [], "source": [ "%%time\n", - "molecule_am1bcc = Molecule(molecule)\n", - "molecule_am1bcc.assign_partial_charges(\n", - " partial_charge_method=\"am1bcc\",\n", - ")" + "sage221 = ForceField(\"openff-2.2.1.offxml\")\n", + "interchange_sage221 = Interchange.from_smirnoff(force_field=sage221, topology=molecule.to_topology())" ] }, { @@ -1124,7 +1152,9 @@ "id": "d751d584", "metadata": {}, "source": [ - "Now, let's try AshGC" + "Note that repeating these cells will show much faster assignment as partial charges are cached for a given molecule and charge method.\n", + "\n", + "Now, let's try Sage 2.3.0, which uses AshGC charges:" ] }, { @@ -1142,54 +1172,28 @@ "outputs": [], "source": [ "%%time\n", - "molecule_ashgc = Molecule(molecule)\n", - "molecule_ashgc.assign_partial_charges(\n", - " partial_charge_method=ASH_GC_MODEL,\n", - " toolkit_registry=NAGLToolkitWrapper(),\n", - ")" + "sage230 = ForceField(\"openff-2.3.0.offxml\")\n", + "interchange_sage230 = Interchange.from_smirnoff(force_field=sage230, topology=molecule.to_topology())" ] }, { "cell_type": "markdown", - "id": "39515e35", + "id": "e18fe6a2", "metadata": {}, "source": [ - "Finally, let's create an `Interchange` with our AshGC charges, making sure to specify `charge_from_molecules` so that we don't replace them with `AM1BCC` charges:" + "
\n", + " ✏️ Exercise: Compare the charges obtained with AM1-BCC and AshGC by inspecting the electrostatics collection of each Interchange object (see Section 3). How big are these differences on average? What is the largest difference? Which atom are these on? The np.max function may be useful.\n", + "
" ] }, { "cell_type": "code", "execution_count": null, - "id": "a81cf535", - "metadata": { - "execution": { - "iopub.execute_input": "2025-10-09T15:21:30.018680Z", - "iopub.status.busy": "2025-10-09T15:21:30.018523Z", - "iopub.status.idle": "2025-10-09T15:21:30.070775Z", - "shell.execute_reply": "2025-10-09T15:21:30.070474Z" - } - }, - "outputs": [], - "source": [ - "# normally when we call `ForceField.create_interchange` or `ForceField.create_openmm_system`, the toolkit will call\n", - "# AMBERTools or OEChem to assign partial charges, since that's what's in the force field file. A future OpenFF release\n", - "# which uses NAGL for charge assignment will encode this instruction in the force field file itself, but until that we\n", - "# can use the `charge_from_molecules` argument to tell it to use the charges that we just assigned# for more, see:\n", - "# https://docs.openforcefield.org/projects/toolkit/en/stable/api/generated/openff.toolkit.typing.engines.smirnoff.ForceField.html#openff.toolkit.typing.engines.smirnoff.ForceField.create_openmm_system\n", - "interchange = sage.create_interchange(\n", - " molecule_ashgc.to_topology(),\n", - " charge_from_molecules=[molecule_ashgc],\n", - ")" - ] - }, - { - "cell_type": "markdown", - "id": "e18fe6a2", + "id": "9791759e", "metadata": {}, + "outputs": [], "source": [ - "
\n", - " ✏️ Exercise: Compare the charges obtained with AM1-BCC and AshGC by looking at the Molecule.partial_charges attribute. How big are these differences on average? What is the largest difference? Which atom are these on? The np.max function may be useful.\n", - "
" + "interchange_sage230[\"Electrostatics\"].key_map" ] }, { @@ -1197,36 +1201,26 @@ "execution_count": null, "id": "8add64c6", "metadata": { - "execution": { - "iopub.execute_input": "2025-10-09T15:21:30.072193Z", - "iopub.status.busy": "2025-10-09T15:21:30.072052Z", - "iopub.status.idle": "2025-10-09T15:21:30.074497Z", - "shell.execute_reply": "2025-10-09T15:21:30.074218Z" - }, "tags": [ "solution" ] }, "outputs": [], "source": [ - "# Compare charges assigned with AM1-BCC and AshGC...\n", + "# Compare charges assigned with AM1-BCC (Sage 2.2.1) and AshGC (Sage 2.3.0)...\n", "import numpy as np\n", - "print(f\"AM1 BCC charges: {molecule_am1bcc.partial_charges}\")\n", - "print(f\"AshGC charges: {molecule_ashgc.partial_charges}\")\n", - "\n", - "differences = molecule_am1bcc.partial_charges - molecule_ashgc.partial_charges\n", - "differences_by_atom_index = {idx: diff.magnitude for idx, diff in enumerate(differences)}\n", - "print(f\"Differences by atom index: {differences_by_atom_index}\")\n", "\n", - "max_difference = np.max(np.abs(differences.magnitude))\n", - "print(f\"Max difference: {max_difference} e\")\n", + "charges_am1bcc = np.array([c.m for c in interchange_sage221[\"Electrostatics\"].charges.values()])\n", + "charges_ashgc = np.array([c.m for c in interchange_sage230[\"Electrostatics\"].charges.values()])\n", "\n", - "mean_difference = np.mean(np.abs(differences.magnitude))\n", - "print(f\"Mean absolute difference: {mean_difference} e\")\n", + "differences = charges_am1bcc - charges_ashgc\n", + "print(f\"Mean absolute difference: {np.mean(np.abs(differences)):.4f} e\")\n", "\n", - "# Get the atom index with the largest difference\n", - "atom_index = np.argmax(np.abs(differences.magnitude))\n", - "print(f\"Largest absolute difference is for atom index {atom_index}, which is a {molecule_ashgc.atoms[atom_index].symbol} atom\")" + "max_index = int(np.argmax(np.abs(differences)))\n", + "print(\n", + " f\"Largest absolute difference is {np.abs(differences[max_index]):.4f} e, \"\n", + " f\"for atom index {max_index}, which is a {molecule.atoms[max_index].symbol} atom\"\n", + ")" ] }, { @@ -1266,8 +1260,8 @@ "
\n", " ✏️ Extra Exercises: Based on the above tutorials, can you:\n", "
    \n", - "
  • Generate several conformers for one of your MCL-1 ligands and compute their relative energies using OpenFF 2.2.1?
    • \n", - "
  • Modify OpenFF 2.2.1 to change some of the parameters applied to one of your MCL-1 ligands? Minimise the ligand with this new force field and see how your changes influence the conformation.
    • \n", + "
  • Generate several conformers for one of your MCL-1 ligands and compute their relative energies using OpenFF 2.3.0?
    • \n", + "
  • Modify OpenFF 2.3.0 to change some of the parameters applied to one of your MCL-1 ligands? Minimise the ligand with this new force field and see how your changes influence the conformation.
    • \n", "
  • Analyse which parameters are shared and which are only applied to one or few molecules for a set of MCL-1 ligands?
    • \n", "
\n", "
" @@ -1290,7 +1284,7 @@ "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", - "version": "3.12.11" + "version": "3.13.15" }, "widgets": { "application/vnd.jupyter.widget-state+json": { From 0d64f476f1faa085241f664c046d8415ecd756bb Mon Sep 17 00:00:00 2001 From: Finlay Clark Date: Wed, 2 Sep 2026 14:22:33 +0100 Subject: [PATCH 3/9] Make minor tweaks to the notebooks --- ...gand_complex_parameterisation_and_md.ipynb | 22 ++++++++++++++----- .../small_molecule_parameterisation.ipynb | 20 +++++------------ ...gand_complex_parameterisation_and_md.ipynb | 16 ++++++++++++-- .../small_molecule_parameterisation.ipynb | 12 +--------- 4 files changed, 37 insertions(+), 33 deletions(-) diff --git a/notebooks/protein_ligand_complex_parameterisation_and_md.ipynb b/notebooks/protein_ligand_complex_parameterisation_and_md.ipynb index 48b5711..6d3d11b 100644 --- a/notebooks/protein_ligand_complex_parameterisation_and_md.ipynb +++ b/notebooks/protein_ligand_complex_parameterisation_and_md.ipynb @@ -298,7 +298,7 @@ { "cell_type": "code", "execution_count": null, - "id": "5c71fa50", + "id": "82b9912f", "metadata": { "tags": [ "placeholder" @@ -421,6 +421,16 @@ "interchange.to_gromacs(prefix=\"complex\")" ] }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "# Check the new gromacs output\n", + "!ls" + ] + }, { "cell_type": "markdown", "metadata": { @@ -711,7 +721,7 @@ { "cell_type": "code", "execution_count": null, - "id": "4a9349fa", + "id": "ed983a3e", "metadata": { "tags": [ "placeholder" @@ -939,7 +949,7 @@ { "cell_type": "code", "execution_count": null, - "id": "145abf20", + "id": "d99ede30", "metadata": { "tags": [ "placeholder" @@ -988,19 +998,21 @@ "* Using OpenMM, we never had to leave Python to set up the simulation.\n", "* With Interchange, using OpenMM, GROMACS, Amber or LAMMPS is simple!\n", "* MDAnalysis and ProLIF allows us to perform varied analyses of our trajectories.\n", + "* The Rosemary force field is likely coming soon, and will allow easy set of simulations with post-translationally modified proteins. Check out [this workshop](https://github.com/openforcefield/2026-virtual-workshops/blob/main/ptm/ptm-workshop.ipynb)!\n", "\n", "\n", "## 7. There's Lots More to OpenFF!\n", "\n", "A variety of example notebooks for OpenFF software are provided [here](https://docs.openforcefield.org/en/latest/examples.html). A few which are particularly relevant are:\n", "\n", + "- [Simulating Post-Translationally Modified Proteins with the OpenFF Rosemary Alpha](https://github.com/openforcefield/2026-virtual-workshops/blob/main/ptm/ptm-workshop.ipynb)\n", "- [Host-guest systems](https://docs.openforcefield.org/en/latest/examples/openforcefield/openff-interchange/host-guest/host_guest.html)\n", "- [Protein-ligand-water systems with Interchange](https://docs.openforcefield.org/en/latest/examples/openforcefield/openff-interchange/protein_ligand/protein_ligand.html). This has a lot of overlap with the current notebook, but there are several extra details not covered here.\n", "\n", "\n", "## 8. Beyond OpenFF\n", "\n", - "You can parameterise your complex and run molecular dynamics -- so what's next? If you're interested in quantiatively assessing the binding affinity of your ligand for your target, then [alchemical (and path-based) free energy calculations are the gold-standard method](https://www.nature.com/articles/s42004-023-01019-9). [Open Free Energy](https://openfree.energy/) is another [Open Molecular Software Foundation](https://omsf.io/) initiative, which develops open-source tools for binding free energy calculations. Head to their [tutorials](https://docs.openfree.energy/en/latest/tutorials/index.html) to learn more! However, these calculations are computationally demanding. If you're interested in a relatively fast (but relatively inaccurate) ranking of the binding affinities of a set of ligands, methods such as MM/GBSA may be appropriate." + "You can parameterise your complex and run molecular dynamics -- so what's next? If you're interested in quantiatively assessing the binding affinity of your ligand for your target, then [alchemical (and path-based) free energy calculations are the gold-standard method](https://livecomsjournal.org/index.php/livecoms/article/view/v2i1e18378). [Open Free Energy](https://openfree.energy/) is another [Open Molecular Software Foundation](https://omsf.io/) initiative, which develops open-source tools for binding free energy calculations. We'll learn about OpenFE tomorrow! However, these calculations are computationally demanding. If you're interested in a relatively fast (but relatively inaccurate) ranking of the binding affinities of a set of ligands, methods such as MM/GBSA may be appropriate. Affinity prediction methods based on deep learning such such as Boltz-2 are appealingly fast, [but perform poorly on systems dissimilar to those they are trained on, and often show inflated performance on benchmarks due to data leakage](https://www.biorxiv.org/content/10.64898/2026.06.29.735309v1.abstract)." ] } ], @@ -1021,7 +1033,7 @@ "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", - "version": "3.12.11" + "version": "3.12.14" }, "widgets": { "application/vnd.jupyter.widget-state+json": { diff --git a/notebooks/small_molecule_parameterisation.ipynb b/notebooks/small_molecule_parameterisation.ipynb index 257b044..9baa0f1 100644 --- a/notebooks/small_molecule_parameterisation.ipynb +++ b/notebooks/small_molecule_parameterisation.ipynb @@ -512,7 +512,7 @@ { "cell_type": "code", "execution_count": null, - "id": "8fcccecc", + "id": "f2e4d2fd", "metadata": { "tags": [ "placeholder" @@ -829,7 +829,7 @@ { "cell_type": "code", "execution_count": null, - "id": "df6bf585", + "id": "82f8bb89", "metadata": { "tags": [ "placeholder" @@ -853,7 +853,7 @@ { "cell_type": "code", "execution_count": null, - "id": "4d49b1fd", + "id": "38f8788d", "metadata": { "tags": [ "placeholder" @@ -1168,17 +1168,7 @@ { "cell_type": "code", "execution_count": null, - "id": "9791759e", - "metadata": {}, - "outputs": [], - "source": [ - "interchange_sage230[\"Electrostatics\"].key_map" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "id": "c6a0b74d", + "id": "25e7a681", "metadata": { "tags": [ "placeholder" @@ -1250,7 +1240,7 @@ "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", - "version": "3.13.15" + "version": "3.12.14" }, "widgets": { "application/vnd.jupyter.widget-state+json": { diff --git a/notebooks_with_solutions/protein_ligand_complex_parameterisation_and_md.ipynb b/notebooks_with_solutions/protein_ligand_complex_parameterisation_and_md.ipynb index 894b967..3586394 100644 --- a/notebooks_with_solutions/protein_ligand_complex_parameterisation_and_md.ipynb +++ b/notebooks_with_solutions/protein_ligand_complex_parameterisation_and_md.ipynb @@ -427,6 +427,16 @@ "interchange.to_gromacs(prefix=\"complex\")" ] }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "# Check the new gromacs output\n", + "!ls" + ] + }, { "cell_type": "markdown", "metadata": { @@ -1031,19 +1041,21 @@ "* Using OpenMM, we never had to leave Python to set up the simulation.\n", "* With Interchange, using OpenMM, GROMACS, Amber or LAMMPS is simple!\n", "* MDAnalysis and ProLIF allows us to perform varied analyses of our trajectories.\n", + "* The Rosemary force field is likely coming soon, and will allow easy set of simulations with post-translationally modified proteins. Check out [this workshop](https://github.com/openforcefield/2026-virtual-workshops/blob/main/ptm/ptm-workshop.ipynb)!\n", "\n", "\n", "## 7. There's Lots More to OpenFF!\n", "\n", "A variety of example notebooks for OpenFF software are provided [here](https://docs.openforcefield.org/en/latest/examples.html). A few which are particularly relevant are:\n", "\n", + "- [Simulating Post-Translationally Modified Proteins with the OpenFF Rosemary Alpha](https://github.com/openforcefield/2026-virtual-workshops/blob/main/ptm/ptm-workshop.ipynb)\n", "- [Host-guest systems](https://docs.openforcefield.org/en/latest/examples/openforcefield/openff-interchange/host-guest/host_guest.html)\n", "- [Protein-ligand-water systems with Interchange](https://docs.openforcefield.org/en/latest/examples/openforcefield/openff-interchange/protein_ligand/protein_ligand.html). This has a lot of overlap with the current notebook, but there are several extra details not covered here.\n", "\n", "\n", "## 8. Beyond OpenFF\n", "\n", - "You can parameterise your complex and run molecular dynamics -- so what's next? If you're interested in quantiatively assessing the binding affinity of your ligand for your target, then [alchemical (and path-based) free energy calculations are the gold-standard method](https://www.nature.com/articles/s42004-023-01019-9). [Open Free Energy](https://openfree.energy/) is another [Open Molecular Software Foundation](https://omsf.io/) initiative, which develops open-source tools for binding free energy calculations. Head to their [tutorials](https://docs.openfree.energy/en/latest/tutorials/index.html) to learn more! However, these calculations are computationally demanding. If you're interested in a relatively fast (but relatively inaccurate) ranking of the binding affinities of a set of ligands, methods such as MM/GBSA may be appropriate." + "You can parameterise your complex and run molecular dynamics -- so what's next? If you're interested in quantiatively assessing the binding affinity of your ligand for your target, then [alchemical (and path-based) free energy calculations are the gold-standard method](https://livecomsjournal.org/index.php/livecoms/article/view/v2i1e18378). [Open Free Energy](https://openfree.energy/) is another [Open Molecular Software Foundation](https://omsf.io/) initiative, which develops open-source tools for binding free energy calculations. We'll learn about OpenFE tomorrow! However, these calculations are computationally demanding. If you're interested in a relatively fast (but relatively inaccurate) ranking of the binding affinities of a set of ligands, methods such as MM/GBSA may be appropriate. Affinity prediction methods based on deep learning such such as Boltz-2 are appealingly fast, [but perform poorly on systems dissimilar to those they are trained on, and often show inflated performance on benchmarks due to data leakage](https://www.biorxiv.org/content/10.64898/2026.06.29.735309v1.abstract)." ] } ], @@ -1064,7 +1076,7 @@ "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", - "version": "3.12.11" + "version": "3.12.14" }, "widgets": { "application/vnd.jupyter.widget-state+json": { diff --git a/notebooks_with_solutions/small_molecule_parameterisation.ipynb b/notebooks_with_solutions/small_molecule_parameterisation.ipynb index ca6cad3..5a75e02 100644 --- a/notebooks_with_solutions/small_molecule_parameterisation.ipynb +++ b/notebooks_with_solutions/small_molecule_parameterisation.ipynb @@ -1186,16 +1186,6 @@ "" ] }, - { - "cell_type": "code", - "execution_count": null, - "id": "9791759e", - "metadata": {}, - "outputs": [], - "source": [ - "interchange_sage230[\"Electrostatics\"].key_map" - ] - }, { "cell_type": "code", "execution_count": null, @@ -1284,7 +1274,7 @@ "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", - "version": "3.13.15" + "version": "3.12.14" }, "widgets": { "application/vnd.jupyter.widget-state+json": { From 2d6bae7283c2ae259feca23ab50280d146efd54b Mon Sep 17 00:00:00 2001 From: Finlay Clark Date: Wed, 2 Sep 2026 14:22:45 +0100 Subject: [PATCH 4/9] Update README for 2026 --- README.md | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/README.md b/README.md index 84396c1..864db43 100644 --- a/README.md +++ b/README.md @@ -1,8 +1,8 @@ -# 2025 CCPBioSim Training Week materials +# 2026 CCPBioSim Training Week materials -[![CI](https://github.com/openforcefield/ccpbiosim-2025/actions/workflows/ci.yaml/badge.svg)](https://github.com/openforcefield/ccpbiosim-2025/actions/workflows/ci.yaml) +[![CI](https://github.com/openforcefield/ccpbiosim-2025/actions/workflows/ci.yaml/badge.svg)](https://github.com/openforcefield/ccpbiosim-workshop/actions/workflows/ci.yaml) -These tutorials were delivered at the 2025 CCPBioSim training week, but are suitable for self-guided learning. +These tutorials were delivered at the 2026 CCPBioSim training week, but are suitable for self-guided learning. Presenters: @@ -40,6 +40,6 @@ $ mamba activate openff-env ## Acknowledgements -Most of the material for the notebook [Parameterising small molecules with OpenFF](notebooks/small_molecule_parameterisation.ipynb) was adapted from the [2023 CCPBioSim Workshop Open Force Field Sessions](https://github.com/openforcefield/ccpbiosim-2023?) created by Matt Thompson and Jeff Wagner. +Most of the material for the notebook [Parameterising small molecules with OpenFF](notebooks/small_molecule_parameterisation.ipynb) was adapted from the [2023 CCPBioSim Workshop Open Force Field Sessions](https://github.com/openforcefield/ccpbiosim-2023?) created by Matt Thompson and Jeff Wagner, as well as the [Simulating Post-Translationally Modified Proteins with the OpenFF Rosemary Alpha Workshop](https://github.com/openforcefield/2026-virtual-workshops/blob/main/ptm/ptm-workshop.ipynb) from Ashley Mitchell. Most of the material for the notebook [Parameterisation, molecular dynamics, and basic trajectory analysis for a protein-ligand complex](notebooks/protein_ligand_complex_parameterisation_and_md.ipynb) was adapted from the OpenFF [toolkit showcase](https://docs.openforcefield.org/en/latest/examples/openforcefield/openff-toolkit/toolkit_showcase/toolkit_showcase.html) and the [ProLIF Ligand-protein MD tutorial](https://prolif.readthedocs.io/en/latest/notebooks/md-ligand-protein.html#ligand-protein-md). From cc6d0782e1576761a674cf95c390d476cf111bdf Mon Sep 17 00:00:00 2001 From: Finlay Clark Date: Wed, 2 Sep 2026 14:49:12 +0100 Subject: [PATCH 5/9] Fix typos and trivial errors Also add more waters to the box to get better density. --- .gitignore | 13 ++++++ README.md | 2 +- ...gand_complex_parameterisation_and_md.ipynb | 6 +-- .../small_molecule_parameterisation.ipynb | 29 +++++++------ ...gand_complex_parameterisation_and_md.ipynb | 6 +-- .../small_molecule_parameterisation.ipynb | 43 ++++++++++--------- 6 files changed, 59 insertions(+), 40 deletions(-) diff --git a/.gitignore b/.gitignore index b7faf40..c80a5f0 100644 --- a/.gitignore +++ b/.gitignore @@ -205,3 +205,16 @@ cython_debug/ marimo/_static/ marimo/_lsp/ __marimo__/ + +# Files generated by running the workshop notebooks +notebooks/*.dcd +notebooks_with_solutions/*.dcd +!notebooks/trajectory_gpu.dcd +!notebooks_with_solutions/trajectory_gpu.dcd +notebooks*/ligand.prmtop +notebooks*/ligand.inpcrd +notebooks*/ligand_pointenergy.in +notebooks*/complex.gro +notebooks*/complex.top +notebooks*/complex_pointenergy.mdp +notebooks*/topology.json diff --git a/README.md b/README.md index 864db43..b728388 100644 --- a/README.md +++ b/README.md @@ -1,6 +1,6 @@ # 2026 CCPBioSim Training Week materials -[![CI](https://github.com/openforcefield/ccpbiosim-2025/actions/workflows/ci.yaml/badge.svg)](https://github.com/openforcefield/ccpbiosim-workshop/actions/workflows/ci.yaml) +[![CI](https://github.com/openforcefield/ccpbiosim-workshop/actions/workflows/ci.yaml/badge.svg)](https://github.com/openforcefield/ccpbiosim-workshop/actions/workflows/ci.yaml) These tutorials were delivered at the 2026 CCPBioSim training week, but are suitable for self-guided learning. diff --git a/notebooks/protein_ligand_complex_parameterisation_and_md.ipynb b/notebooks/protein_ligand_complex_parameterisation_and_md.ipynb index 6d3d11b..a2ff4d2 100644 --- a/notebooks/protein_ligand_complex_parameterisation_and_md.ipynb +++ b/notebooks/protein_ligand_complex_parameterisation_and_md.ipynb @@ -449,7 +449,7 @@ "\n", "### 4.1 Configure and run the simulation\n", "\n", - "Here, we'll use a Langevin thermostat at 300 Kelvin and a 2 fs time step. We'll write the structure to disk every 10 steps. In contrast to the previous notebook, we'll add a MonteCarloBarostat to fix the pressure, while allowing the volume to fluctuate. Our simulation corresponds to the $NPT$ ensemble." + "Here, we'll use a Langevin thermostat at 300 Kelvin and a 2 fs time step. We'll write the structure to disk every 50 steps. In contrast to the previous notebook, we'll add a MonteCarloBarostat to fix the pressure, while allowing the volume to fluctuate. Our simulation corresponds to the $NPT$ ensemble." ] }, { @@ -998,7 +998,7 @@ "* Using OpenMM, we never had to leave Python to set up the simulation.\n", "* With Interchange, using OpenMM, GROMACS, Amber or LAMMPS is simple!\n", "* MDAnalysis and ProLIF allows us to perform varied analyses of our trajectories.\n", - "* The Rosemary force field is likely coming soon, and will allow easy set of simulations with post-translationally modified proteins. Check out [this workshop](https://github.com/openforcefield/2026-virtual-workshops/blob/main/ptm/ptm-workshop.ipynb)!\n", + "* The Rosemary force field is likely coming soon, and will allow easy set-up of simulations with post-translationally modified proteins. Check out [this workshop](https://github.com/openforcefield/2026-virtual-workshops/blob/main/ptm/ptm-workshop.ipynb)!\n", "\n", "\n", "## 7. There's Lots More to OpenFF!\n", @@ -1012,7 +1012,7 @@ "\n", "## 8. Beyond OpenFF\n", "\n", - "You can parameterise your complex and run molecular dynamics -- so what's next? If you're interested in quantiatively assessing the binding affinity of your ligand for your target, then [alchemical (and path-based) free energy calculations are the gold-standard method](https://livecomsjournal.org/index.php/livecoms/article/view/v2i1e18378). [Open Free Energy](https://openfree.energy/) is another [Open Molecular Software Foundation](https://omsf.io/) initiative, which develops open-source tools for binding free energy calculations. We'll learn about OpenFE tomorrow! However, these calculations are computationally demanding. If you're interested in a relatively fast (but relatively inaccurate) ranking of the binding affinities of a set of ligands, methods such as MM/GBSA may be appropriate. Affinity prediction methods based on deep learning such such as Boltz-2 are appealingly fast, [but perform poorly on systems dissimilar to those they are trained on, and often show inflated performance on benchmarks due to data leakage](https://www.biorxiv.org/content/10.64898/2026.06.29.735309v1.abstract)." + "You can parameterise your complex and run molecular dynamics -- so what's next? If you're interested in quantitatively assessing the binding affinity of your ligand for your target, then [alchemical (and path-based) free energy calculations are the gold-standard method](https://livecomsjournal.org/index.php/livecoms/article/view/v2i1e18378). [Open Free Energy](https://openfree.energy/) is another [Open Molecular Software Foundation](https://omsf.io/) initiative, which develops open-source tools for binding free energy calculations. However, these calculations are computationally demanding. If you're interested in a relatively fast (but relatively inaccurate) ranking of the binding affinities of a set of ligands, methods such as MM/GBSA may be appropriate. Affinity prediction methods based on deep learning such as Boltz-2 are appealingly fast, [but perform poorly on systems dissimilar to those they are trained on, and often show inflated performance on benchmarks due to data leakage](https://www.biorxiv.org/content/10.64898/2026.06.29.735309v1.abstract)." ] } ], diff --git a/notebooks/small_molecule_parameterisation.ipynb b/notebooks/small_molecule_parameterisation.ipynb index 9baa0f1..e9e4f45 100644 --- a/notebooks/small_molecule_parameterisation.ipynb +++ b/notebooks/small_molecule_parameterisation.ipynb @@ -18,7 +18,7 @@ "\n", "| Action | Software|\n", "|--|--|\n", - "| [Go from smiles to simulation in few lines of code](#showcase) | OpenFF Toolkit, OpenFF Interchange, OpenMM\n", + "| [Go from SMILES to simulation in a few lines of code](#showcase) | OpenFF Toolkit, OpenFF Interchange, OpenMM\n", "| [Load and inspect a force field](#loading_ff) | OpenFF Toolkit\n", "| [Create a representation of your chemical system](#topology) | OpenFF Toolkit\n", "| [Parameterise your system and run a quick simulation in water](#interchange) | OpenFF Interchange, OpenMM\n", @@ -82,6 +82,7 @@ "source": [ "# Run the simulation with OpenMM\n", "import openmm\n", + "import openmm.app\n", "\n", "temperature = 298.15 * openmm.unit.kelvin\n", "friction_coefficient = 1.0 / openmm.unit.picosecond\n", @@ -138,7 +139,7 @@ "\n", "\"Description\n", "\n", - "Let's start with the `.offxml` force field file. OpenFF's force fields use the The SMIRKS Native Open Force Field (SMIRNOFF) [specification](https://openforcefield.github.io/standards/standards/smirnoff/). The spec fully describes the contents of a SMIRNOFF force field, how parameters should be applied, and several other important usage details. You could implement a SMIRNOFF engine in your own code, but conveniently the OpenFF Toolkit already provides this and a handful of utilities. Let's load up the latest OpenFF small molecule force field, OpenFF 2.3,0, and inspect its contents! This force field shares the code name \"Sage\" with all other force fields with the same major version number (2.x.x)." + "Let's start with the `.offxml` force field file. OpenFF's force fields use the SMIRKS Native Open Force Field (SMIRNOFF) [specification](https://openforcefield.github.io/standards/standards/smirnoff/). The spec fully describes the contents of a SMIRNOFF force field, how parameters should be applied, and several other important usage details. You could implement a SMIRNOFF engine in your own code, but conveniently the OpenFF Toolkit already provides this and a handful of utilities. Let's load up the latest OpenFF small molecule force field, OpenFF 2.3.0, and inspect its contents! This force field shares the code name \"Sage\" with all other force fields with the same major version number (2.x.x)." ] }, { @@ -260,7 +261,7 @@ "\n", "The toolkit uses these SMIRKS patterns and direct chemical perception to assign parameters to particular atoms (or bonds, angles, etc.).\n", "\n", - "We'll use OpenFF 2.3.0 for the remainder of this tutorial. This is OpenFF's latest small molecule force field and is a leading open-source small molecule force field which [performs comparably to other open-source force fields](https://doi-org.libproxy.ncl.ac.uk/10.1021/acs.jctc.3c00039). You can learn more about this and other SMIRNOFF force fields below:\n", + "We'll use OpenFF 2.3.0 for the remainder of this tutorial. This is OpenFF's latest small molecule force field and is a leading open-source small molecule force field which [performs comparably to other open-source force fields](https://doi.org/10.1021/acs.jctc.3c00039). You can learn more about this and other SMIRNOFF force fields below:\n", "
\n", " Click here to learn about available and planned SMIRNOFF force fields\n", "\n", @@ -282,7 +283,7 @@ "\n", "The Sage line of force fields (`openff-2.y.z.offxml`) continued the process of fitting to more (and more diverse) QM datasets, but also included a re-fit of the Lennard-Jones parameters. Small molecule geometries and energies [improved, in general,](https://openforcefield.org/community/news/general/sage2.0.0-release/) significantly over Parsley. These improvements notably transferred to protein-ligand binding free energies despite Sage not being specifically fit to them. For more, see the [associated paper](https://pubs.acs.org/doi/10.1021/acs.jctc.3c00039).\n", "\n", - "[Subsequent releases](https://github.com/openforcefield/openff-forcefields/releases) used different fitting procedures and tweaks to parameter typing to improve performance and address issues with several specific chemistries. Notably, Sage 2.3.0 includes fast graph neural network charge assignment with Ash-GC, which is discussed later in this notebook. This charge model is trained to reproduce AM1-BCC charges without the typical $O(N^3)$ scaling, making it suitable for large (>> 100 atoms) molecules). The latest release, **Sage 2.3.0 (`openff-2.3.0.offxml`) is the recommended force field for small molecule studies.**\n", + "[Subsequent releases](https://github.com/openforcefield/openff-forcefields/releases) used different fitting procedures and tweaks to parameter typing to improve performance and address issues with several specific chemistries. Notably, Sage 2.3.0 includes fast graph neural network charge assignment with AshGC, which is discussed later in this notebook. This charge model is trained to reproduce AM1-BCC charges, but scales 𝒪(N) rather than the 𝒪(N2) of common AM1-BCC implementations, making it suitable for large (>> 100 atoms) molecules. The latest release, **Sage 2.3.0 (`openff-2.3.0.offxml`) is the recommended force field for small molecule studies.**\n", "\n", "## Ports\n", "\n", @@ -300,13 +301,13 @@ "\n", "## ff14SB\n", "\n", - "OpenFF, in collaboration with Dave Cerutti of the Amber community, created a port of [ff14SB](https://pubs.acs.org/doi/10.1021/acs.jctc.5b00255), a popular Amber protein force field. There are some small numerical differences with how improper torsions are evaluated, but all other terms reproduce a canonical Amber source to high accuracy. **This is the only protein force field currently in SMIRONOFF (`.offxml`) format** and therefore the current recommendation for use with proteins. Primarily for technical reasons, porting other Amber force fields is not planned.\n", + "OpenFF, in collaboration with Dave Cerutti of the Amber community, created a port of [ff14SB](https://pubs.acs.org/doi/10.1021/acs.jctc.5b00255), a popular Amber protein force field. There are some small numerical differences with how improper torsions are evaluated, but all other terms reproduce a canonical Amber source to high accuracy. **This is the only protein force field currently in SMIRNOFF (`.offxml`) format** and therefore the current recommendation for use with proteins. Primarily for technical reasons, porting other Amber force fields is not planned.\n", "\n", "## Rosemary Alpha\n", "\n", - "A future line of force fields from OpenFF (code name \"Rosemary\", starting with `openff-3.0.0.offxml`) is intended to handle small molecules and biopolymers in a _self-consistent_ manner. This is exciting as it will streamline simulations of proteins with non-cannonical amino acids! See the workshop [Simulating Post-Translationally Modified Proteins with the OpenFF Rosemary Alpha](https://github.com/openforcefield/2026-virtual-workshops/blob/main/ptm/ptm-workshop.ipynb). The first release will handle proteins, but future versions may cover nucleic acids. The performance, depending on the metrics used, is hoped to be comparable with existing Amber-family protein force fields. \n", + "A future line of force fields from OpenFF (code name \"Rosemary\", starting with `openff-3.0.0.offxml`) is intended to handle small molecules and biopolymers in a _self-consistent_ manner. This is exciting as it will streamline simulations of proteins with non-canonical amino acids! See the workshop [Simulating Post-Translationally Modified Proteins with the OpenFF Rosemary Alpha](https://github.com/openforcefield/2026-virtual-workshops/blob/main/ptm/ptm-workshop.ipynb). The first release will handle proteins, but future versions may cover nucleic acids. The performance, depending on the metrics used, is hoped to be comparable with existing Amber-family protein force fields. \n", "\n", - "An pre-release version of Rosemary is available for testing as [`openff_no_water-3.0.0-alpha0.offxml`](https://github.com/openforcefield/openff-forcefields/blob/main/openforcefields/offxml/openff_no_water-3.0.0-alpha0.offxml). If you use it, see [the release notes](https://github.com/openforcefield/openff-forcefields/releases/tag/2025.10.1). \n", + "A pre-release version of Rosemary is available for testing as [`openff_no_water-3.0.0-alpha0.offxml`](https://github.com/openforcefield/openff-forcefields/blob/main/openforcefields/offxml/openff_no_water-3.0.0-alpha0.offxml). If you use it, see [the release notes](https://github.com/openforcefield/openff-forcefields/releases/tag/2025.10.1). \n", "\n", "\n", "## Non-main-line force fields\n", @@ -326,7 +327,7 @@ "\n", "### Virtual sites\n", "\n", - "Another release from OpenFF may include some virtual site parameters with off-center charges. No release date is planned, but the most of the supporting infrastructure is currently in place and some early studies have shown promise for better representing electrostatics of chemistries such as halogens and aromatic nitrogens.\n", + "Another release from OpenFF may include some virtual site parameters with off-center charges. No release date is planned, but most of the supporting infrastructure is currently in place and some early studies have shown promise for better representing electrostatics of chemistries such as halogens and aromatic nitrogens.\n", "\n", "## From you!\n", "\n", @@ -505,7 +506,7 @@ "metadata": {}, "source": [ "
\n", - " ✏️ Exercise: Build a Topology containing an MCL-1 ligand. Create the Molecule from an SDF file (take a look at the docstring of Molecule to see how this can be done, noting that you can just pass the sdf path given below and don't need the get_data_file_path function). Also, see Molecule cookbook for all the ways to make a Molecule. The crystallographic MCL-1 ligand from PDB ID 6o6f is provided at ../structures/606f_ligand.sdf. Note that you don't need to use get_data_file_path as we already know the path.\n", + " ✏️ Exercise: Build a Topology containing an MCL-1 ligand. Create the Molecule from an SDF file (take a look at the docstring of Molecule to see how this can be done, noting that you can just pass the sdf path given below and don't need the get_data_file_path function). Also, see Molecule cookbook for all the ways to make a Molecule. The crystallographic MCL-1 ligand from PDB ID 6o6f is provided at ../structures/6o6f_ligand.sdf. Note that you don't need to use get_data_file_path as we already know the path.\n", "
\n" ] }, @@ -567,6 +568,7 @@ "outputs": [], "source": [ "molecule = Molecule.from_smiles(\"C/C=C/C(=O)[O-]\")\n", + "molecule.generate_conformers(n_conformers=1)\n", "topology = molecule.to_topology()" ] }, @@ -929,7 +931,7 @@ "# the volume.\n", "topology = pack_box(\n", " molecules=[solute, water],\n", - " number_of_copies=[1, 1000],\n", + " number_of_copies=[1, 1400],\n", " box_vectors=3.5 * UNIT_CUBE * unit.nanometer,\n", ")\n", "\n", @@ -1004,6 +1006,7 @@ "outputs": [], "source": [ "import openmm\n", + "import openmm.app\n", "import openmm.unit\n", "from openff.interchange import Interchange\n", "import MDAnalysis as mda\n", @@ -1065,9 +1068,9 @@ "\n", "## 4. Graph Neural Networks Allow Fast Assignment of Partial Charges\n", "\n", - "You might notice that Sage force fields don't contain tabulated charges for most atomic environments in the way they do for all other terms in the force field.For example, [Sage 2.2.1](https://github.com/openforcefield/openff-forcefields/blob/main/openforcefields/offxml/openff-2.2.1.offxml) instead specifies:\n", + "You might notice that Sage force fields don't contain tabulated charges for most atomic environments in the way they do for all other terms in the force field. For example, [Sage 2.2.1](https://github.com/openforcefield/openff-forcefields/blob/main/openforcefields/offxml/openff-2.2.1.offxml) instead specifies:\n", "```\n", - "\n", + "\n", "```\n", "which means that partial charges will be calculated using the common AM1-BCC method. Charges from a semi-empirical quantum chemistry calculation (Austin Model 1) are corrected (bond charge correction) to approximate charges obtained by fitting to the electrostatic potential at the HF/6-31G* level (see [Jakalian et al.](https://onlinelibrary.wiley.com/doi/10.1002/(SICI)1096-987X(20000130)21:2%3C132::AID-JCC5%3E3.0.CO;2-P)). Unfortunately, parameterisation with AM1-BCC using OpenEye or AmberTools scales 𝒪(N2) in the number of atoms N, making it prohibitively slow for large molecules and biopolymers.\n", "\n", @@ -1240,7 +1243,7 @@ "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", - "version": "3.12.14" + "version": "3.13.15" }, "widgets": { "application/vnd.jupyter.widget-state+json": { diff --git a/notebooks_with_solutions/protein_ligand_complex_parameterisation_and_md.ipynb b/notebooks_with_solutions/protein_ligand_complex_parameterisation_and_md.ipynb index 3586394..2999559 100644 --- a/notebooks_with_solutions/protein_ligand_complex_parameterisation_and_md.ipynb +++ b/notebooks_with_solutions/protein_ligand_complex_parameterisation_and_md.ipynb @@ -455,7 +455,7 @@ "\n", "### 4.1 Configure and run the simulation\n", "\n", - "Here, we'll use a Langevin thermostat at 300 Kelvin and a 2 fs time step. We'll write the structure to disk every 10 steps. In contrast to the previous notebook, we'll add a MonteCarloBarostat to fix the pressure, while allowing the volume to fluctuate. Our simulation corresponds to the $NPT$ ensemble." + "Here, we'll use a Langevin thermostat at 300 Kelvin and a 2 fs time step. We'll write the structure to disk every 50 steps. In contrast to the previous notebook, we'll add a MonteCarloBarostat to fix the pressure, while allowing the volume to fluctuate. Our simulation corresponds to the $NPT$ ensemble." ] }, { @@ -1041,7 +1041,7 @@ "* Using OpenMM, we never had to leave Python to set up the simulation.\n", "* With Interchange, using OpenMM, GROMACS, Amber or LAMMPS is simple!\n", "* MDAnalysis and ProLIF allows us to perform varied analyses of our trajectories.\n", - "* The Rosemary force field is likely coming soon, and will allow easy set of simulations with post-translationally modified proteins. Check out [this workshop](https://github.com/openforcefield/2026-virtual-workshops/blob/main/ptm/ptm-workshop.ipynb)!\n", + "* The Rosemary force field is likely coming soon, and will allow easy set-up of simulations with post-translationally modified proteins. Check out [this workshop](https://github.com/openforcefield/2026-virtual-workshops/blob/main/ptm/ptm-workshop.ipynb)!\n", "\n", "\n", "## 7. There's Lots More to OpenFF!\n", @@ -1055,7 +1055,7 @@ "\n", "## 8. Beyond OpenFF\n", "\n", - "You can parameterise your complex and run molecular dynamics -- so what's next? If you're interested in quantiatively assessing the binding affinity of your ligand for your target, then [alchemical (and path-based) free energy calculations are the gold-standard method](https://livecomsjournal.org/index.php/livecoms/article/view/v2i1e18378). [Open Free Energy](https://openfree.energy/) is another [Open Molecular Software Foundation](https://omsf.io/) initiative, which develops open-source tools for binding free energy calculations. We'll learn about OpenFE tomorrow! However, these calculations are computationally demanding. If you're interested in a relatively fast (but relatively inaccurate) ranking of the binding affinities of a set of ligands, methods such as MM/GBSA may be appropriate. Affinity prediction methods based on deep learning such such as Boltz-2 are appealingly fast, [but perform poorly on systems dissimilar to those they are trained on, and often show inflated performance on benchmarks due to data leakage](https://www.biorxiv.org/content/10.64898/2026.06.29.735309v1.abstract)." + "You can parameterise your complex and run molecular dynamics -- so what's next? If you're interested in quantitatively assessing the binding affinity of your ligand for your target, then [alchemical (and path-based) free energy calculations are the gold-standard method](https://livecomsjournal.org/index.php/livecoms/article/view/v2i1e18378). [Open Free Energy](https://openfree.energy/) is another [Open Molecular Software Foundation](https://omsf.io/) initiative, which develops open-source tools for binding free energy calculations. However, these calculations are computationally demanding. If you're interested in a relatively fast (but relatively inaccurate) ranking of the binding affinities of a set of ligands, methods such as MM/GBSA may be appropriate. Affinity prediction methods based on deep learning such as Boltz-2 are appealingly fast, [but perform poorly on systems dissimilar to those they are trained on, and often show inflated performance on benchmarks due to data leakage](https://www.biorxiv.org/content/10.64898/2026.06.29.735309v1.abstract)." ] } ], diff --git a/notebooks_with_solutions/small_molecule_parameterisation.ipynb b/notebooks_with_solutions/small_molecule_parameterisation.ipynb index 5a75e02..8005dff 100644 --- a/notebooks_with_solutions/small_molecule_parameterisation.ipynb +++ b/notebooks_with_solutions/small_molecule_parameterisation.ipynb @@ -18,7 +18,7 @@ "\n", "| Action | Software|\n", "|--|--|\n", - "| [Go from smiles to simulation in few lines of code](#showcase) | OpenFF Toolkit, OpenFF Interchange, OpenMM\n", + "| [Go from SMILES to simulation in a few lines of code](#showcase) | OpenFF Toolkit, OpenFF Interchange, OpenMM\n", "| [Load and inspect a force field](#loading_ff) | OpenFF Toolkit\n", "| [Create a representation of your chemical system](#topology) | OpenFF Toolkit\n", "| [Parameterise your system and run a quick simulation in water](#interchange) | OpenFF Interchange, OpenMM\n", @@ -82,6 +82,7 @@ "source": [ "# Run the simulation with OpenMM\n", "import openmm\n", + "import openmm.app\n", "\n", "temperature = 298.15 * openmm.unit.kelvin\n", "friction_coefficient = 1.0 / openmm.unit.picosecond\n", @@ -138,7 +139,7 @@ "\n", "\"Description\n", "\n", - "Let's start with the `.offxml` force field file. OpenFF's force fields use the The SMIRKS Native Open Force Field (SMIRNOFF) [specification](https://openforcefield.github.io/standards/standards/smirnoff/). The spec fully describes the contents of a SMIRNOFF force field, how parameters should be applied, and several other important usage details. You could implement a SMIRNOFF engine in your own code, but conveniently the OpenFF Toolkit already provides this and a handful of utilities. Let's load up the latest OpenFF small molecule force field, OpenFF 2.3,0, and inspect its contents! This force field shares the code name \"Sage\" with all other force fields with the same major version number (2.x.x)." + "Let's start with the `.offxml` force field file. OpenFF's force fields use the SMIRKS Native Open Force Field (SMIRNOFF) [specification](https://openforcefield.github.io/standards/standards/smirnoff/). The spec fully describes the contents of a SMIRNOFF force field, how parameters should be applied, and several other important usage details. You could implement a SMIRNOFF engine in your own code, but conveniently the OpenFF Toolkit already provides this and a handful of utilities. Let's load up the latest OpenFF small molecule force field, OpenFF 2.3.0, and inspect its contents! This force field shares the code name \"Sage\" with all other force fields with the same major version number (2.x.x)." ] }, { @@ -260,7 +261,7 @@ "\n", "The toolkit uses these SMIRKS patterns and direct chemical perception to assign parameters to particular atoms (or bonds, angles, etc.).\n", "\n", - "We'll use OpenFF 2.3.0 for the remainder of this tutorial. This is OpenFF's latest small molecule force field and is a leading open-source small molecule force field which [performs comparably to other open-source force fields](https://doi-org.libproxy.ncl.ac.uk/10.1021/acs.jctc.3c00039). You can learn more about this and other SMIRNOFF force fields below:\n", + "We'll use OpenFF 2.3.0 for the remainder of this tutorial. This is OpenFF's latest small molecule force field and is a leading open-source small molecule force field which [performs comparably to other open-source force fields](https://doi.org/10.1021/acs.jctc.3c00039). You can learn more about this and other SMIRNOFF force fields below:\n", "
\n", " Click here to learn about available and planned SMIRNOFF force fields\n", "\n", @@ -282,7 +283,7 @@ "\n", "The Sage line of force fields (`openff-2.y.z.offxml`) continued the process of fitting to more (and more diverse) QM datasets, but also included a re-fit of the Lennard-Jones parameters. Small molecule geometries and energies [improved, in general,](https://openforcefield.org/community/news/general/sage2.0.0-release/) significantly over Parsley. These improvements notably transferred to protein-ligand binding free energies despite Sage not being specifically fit to them. For more, see the [associated paper](https://pubs.acs.org/doi/10.1021/acs.jctc.3c00039).\n", "\n", - "[Subsequent releases](https://github.com/openforcefield/openff-forcefields/releases) used different fitting procedures and tweaks to parameter typing to improve performance and address issues with several specific chemistries. Notably, Sage 2.3.0 includes fast graph neural network charge assignment with Ash-GC, which is discussed later in this notebook. This charge model is trained to reproduce AM1-BCC charges without the typical $O(N^3)$ scaling, making it suitable for large (>> 100 atoms) molecules). The latest release, **Sage 2.3.0 (`openff-2.3.0.offxml`) is the recommended force field for small molecule studies.**\n", + "[Subsequent releases](https://github.com/openforcefield/openff-forcefields/releases) used different fitting procedures and tweaks to parameter typing to improve performance and address issues with several specific chemistries. Notably, Sage 2.3.0 includes fast graph neural network charge assignment with AshGC, which is discussed later in this notebook. This charge model is trained to reproduce AM1-BCC charges, but scales 𝒪(N) rather than the 𝒪(N2) of common AM1-BCC implementations, making it suitable for large (>> 100 atoms) molecules. The latest release, **Sage 2.3.0 (`openff-2.3.0.offxml`) is the recommended force field for small molecule studies.**\n", "\n", "## Ports\n", "\n", @@ -300,13 +301,13 @@ "\n", "## ff14SB\n", "\n", - "OpenFF, in collaboration with Dave Cerutti of the Amber community, created a port of [ff14SB](https://pubs.acs.org/doi/10.1021/acs.jctc.5b00255), a popular Amber protein force field. There are some small numerical differences with how improper torsions are evaluated, but all other terms reproduce a canonical Amber source to high accuracy. **This is the only protein force field currently in SMIRONOFF (`.offxml`) format** and therefore the current recommendation for use with proteins. Primarily for technical reasons, porting other Amber force fields is not planned.\n", + "OpenFF, in collaboration with Dave Cerutti of the Amber community, created a port of [ff14SB](https://pubs.acs.org/doi/10.1021/acs.jctc.5b00255), a popular Amber protein force field. There are some small numerical differences with how improper torsions are evaluated, but all other terms reproduce a canonical Amber source to high accuracy. **This is the only protein force field currently in SMIRNOFF (`.offxml`) format** and therefore the current recommendation for use with proteins. Primarily for technical reasons, porting other Amber force fields is not planned.\n", "\n", "## Rosemary Alpha\n", "\n", - "A future line of force fields from OpenFF (code name \"Rosemary\", starting with `openff-3.0.0.offxml`) is intended to handle small molecules and biopolymers in a _self-consistent_ manner. This is exciting as it will streamline simulations of proteins with non-cannonical amino acids! See the workshop [Simulating Post-Translationally Modified Proteins with the OpenFF Rosemary Alpha](https://github.com/openforcefield/2026-virtual-workshops/blob/main/ptm/ptm-workshop.ipynb). The first release will handle proteins, but future versions may cover nucleic acids. The performance, depending on the metrics used, is hoped to be comparable with existing Amber-family protein force fields. \n", + "A future line of force fields from OpenFF (code name \"Rosemary\", starting with `openff-3.0.0.offxml`) is intended to handle small molecules and biopolymers in a _self-consistent_ manner. This is exciting as it will streamline simulations of proteins with non-canonical amino acids! See the workshop [Simulating Post-Translationally Modified Proteins with the OpenFF Rosemary Alpha](https://github.com/openforcefield/2026-virtual-workshops/blob/main/ptm/ptm-workshop.ipynb). The first release will handle proteins, but future versions may cover nucleic acids. The performance, depending on the metrics used, is hoped to be comparable with existing Amber-family protein force fields. \n", "\n", - "An pre-release version of Rosemary is available for testing as [`openff_no_water-3.0.0-alpha0.offxml`](https://github.com/openforcefield/openff-forcefields/blob/main/openforcefields/offxml/openff_no_water-3.0.0-alpha0.offxml). If you use it, see [the release notes](https://github.com/openforcefield/openff-forcefields/releases/tag/2025.10.1). \n", + "A pre-release version of Rosemary is available for testing as [`openff_no_water-3.0.0-alpha0.offxml`](https://github.com/openforcefield/openff-forcefields/blob/main/openforcefields/offxml/openff_no_water-3.0.0-alpha0.offxml). If you use it, see [the release notes](https://github.com/openforcefield/openff-forcefields/releases/tag/2025.10.1). \n", "\n", "\n", "## Non-main-line force fields\n", @@ -326,7 +327,7 @@ "\n", "### Virtual sites\n", "\n", - "Another release from OpenFF may include some virtual site parameters with off-center charges. No release date is planned, but the most of the supporting infrastructure is currently in place and some early studies have shown promise for better representing electrostatics of chemistries such as halogens and aromatic nitrogens.\n", + "Another release from OpenFF may include some virtual site parameters with off-center charges. No release date is planned, but most of the supporting infrastructure is currently in place and some early studies have shown promise for better representing electrostatics of chemistries such as halogens and aromatic nitrogens.\n", "\n", "## From you!\n", "\n", @@ -505,7 +506,7 @@ "metadata": {}, "source": [ "
\n", - " ✏️ Exercise: Build a Topology containing an MCL-1 ligand. Create the Molecule from an SDF file (take a look at the docstring of Molecule to see how this can be done, noting that you can just pass the sdf path given below and don't need the get_data_file_path function). Also, see Molecule cookbook for all the ways to make a Molecule. The crystallographic MCL-1 ligand from PDB ID 6o6f is provided at ../structures/606f_ligand.sdf. Note that you don't need to use get_data_file_path as we already know the path.\n", + " ✏️ Exercise: Build a Topology containing an MCL-1 ligand. Create the Molecule from an SDF file (take a look at the docstring of Molecule to see how this can be done, noting that you can just pass the sdf path given below and don't need the get_data_file_path function). Also, see Molecule cookbook for all the ways to make a Molecule. The crystallographic MCL-1 ligand from PDB ID 6o6f is provided at ../structures/6o6f_ligand.sdf. Note that you don't need to use get_data_file_path as we already know the path.\n", "
\n" ] }, @@ -569,6 +570,7 @@ "outputs": [], "source": [ "molecule = Molecule.from_smiles(\"C/C=C/C(=O)[O-]\")\n", + "molecule.generate_conformers(n_conformers=1)\n", "topology = molecule.to_topology()" ] }, @@ -870,17 +872,17 @@ "bond1_indices, bond2_indices = (3,4), (3,5)\n", "smiles = {\"anion\": \"C/C=C/C(=O)[O-]\", \"neutral\": \"C/C=C/C(=O)O\"}\n", "\n", + "# Use fresh variable names so we don't overwrite `molecule`, `topology` and `interchange`\n", "for name, smiles in smiles.items():\n", " print(f\"\\n{name} ({smiles}):\")\n", - " molecule = Molecule.from_smiles(smiles)\n", - " topology = molecule.to_topology()\n", - " interchange = Interchange.from_smirnoff(\n", + " protonation_state = Molecule.from_smiles(smiles)\n", + " protonation_state_interchange = Interchange.from_smirnoff(\n", " force_field=sage,\n", - " topology=topology,\n", + " topology=protonation_state.to_topology(),\n", " )\n", - " collection = interchange.collections[\"Bonds\"]\n", - " bond1 = collection[bond1_indices]\n", - " bond2 = collection[bond2_indices]\n", + " bonds = protonation_state_interchange.collections[\"Bonds\"]\n", + " bond1 = bonds[bond1_indices]\n", + " bond2 = bonds[bond2_indices]\n", " print(f\" Bond {bond1_indices}: {bond1}\")\n", " print(f\" Bond {bond2_indices}: {bond2}\")" ] @@ -950,7 +952,7 @@ "# the volume.\n", "topology = pack_box(\n", " molecules=[solute, water],\n", - " number_of_copies=[1, 1000],\n", + " number_of_copies=[1, 1400],\n", " box_vectors=3.5 * UNIT_CUBE * unit.nanometer,\n", ")\n", "\n", @@ -1025,6 +1027,7 @@ "outputs": [], "source": [ "import openmm\n", + "import openmm.app\n", "import openmm.unit\n", "from openff.interchange import Interchange\n", "import MDAnalysis as mda\n", @@ -1086,9 +1089,9 @@ "\n", "## 4. Graph Neural Networks Allow Fast Assignment of Partial Charges\n", "\n", - "You might notice that Sage force fields don't contain tabulated charges for most atomic environments in the way they do for all other terms in the force field.For example, [Sage 2.2.1](https://github.com/openforcefield/openff-forcefields/blob/main/openforcefields/offxml/openff-2.2.1.offxml) instead specifies:\n", + "You might notice that Sage force fields don't contain tabulated charges for most atomic environments in the way they do for all other terms in the force field. For example, [Sage 2.2.1](https://github.com/openforcefield/openff-forcefields/blob/main/openforcefields/offxml/openff-2.2.1.offxml) instead specifies:\n", "```\n", - "\n", + "\n", "```\n", "which means that partial charges will be calculated using the common AM1-BCC method. Charges from a semi-empirical quantum chemistry calculation (Austin Model 1) are corrected (bond charge correction) to approximate charges obtained by fitting to the electrostatic potential at the HF/6-31G* level (see [Jakalian et al.](https://onlinelibrary.wiley.com/doi/10.1002/(SICI)1096-987X(20000130)21:2%3C132::AID-JCC5%3E3.0.CO;2-P)). Unfortunately, parameterisation with AM1-BCC using OpenEye or AmberTools scales 𝒪(N2) in the number of atoms N, making it prohibitively slow for large molecules and biopolymers.\n", "\n", @@ -1274,7 +1277,7 @@ "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", - "version": "3.12.14" + "version": "3.13.15" }, "widgets": { "application/vnd.jupyter.widget-state+json": { From ab12283d6a18862774c48c027700fbfdcbf8fdd1 Mon Sep 17 00:00:00 2001 From: "github-actions[bot]" Date: Tue, 8 Sep 2026 20:48:49 +0000 Subject: [PATCH 6/9] Add generated notebook.md and answerless notebook files [skip ci] --- ..._ligand_complex_parameterisation_and_md.md | 55 +-- .../small_molecule_parameterisation.md | 340 ++++++++++-------- ...gand_complex_parameterisation_and_md.ipynb | 6 +- .../small_molecule_parameterisation.ipynb | 8 +- 4 files changed, 232 insertions(+), 177 deletions(-) diff --git a/notebooks-rendered/protein_ligand_complex_parameterisation_and_md.md b/notebooks-rendered/protein_ligand_complex_parameterisation_and_md.md index 5a99167..4cf04e4 100644 --- a/notebooks-rendered/protein_ligand_complex_parameterisation_and_md.md +++ b/notebooks-rendered/protein_ligand_complex_parameterisation_and_md.md @@ -78,10 +78,6 @@ view - /opt/conda/envs/openff-env/lib/python3.12/site-packages/nglview/__init__.py:12: UserWarning: pkg_resources is deprecated as an API. See https://setuptools.pypa.io/en/latest/pkg_resources.html. The pkg_resources package is slated for removal as early as 2025-11-30. Refrain from using this package or pin to Setuptools<81. - import pkg_resources - - NGLWidget() @@ -94,7 +90,7 @@ view ## 2. OpenFF Toolkit Allows Us to Assemble the Topology -Conceptually, this step involves putting together the positions of all of the components of the system. We'll create a [`Topology`] to keep track of the contents of our system. As discussed in this morning's session, `Topology` represents a collection of molecules; it doesn't have any association with any force field parameters. +Conceptually, this step involves putting together the positions of all of the components of the system. We'll create a [`Topology`] to keep track of the contents of our system. As discussed in the previous notebook, `Topology` represents a collection of molecules; it doesn't have any association with any force field parameters. [`Topology`]: https://docs.openforcefield.org/projects/toolkit/en/stable/api/generated/openff.toolkit.topology.Topology.html @@ -220,14 +216,14 @@ view ## 3. We Can Assemble a Combined `ForceField` and use this to Parameterise the Whole System -Now that we've prepared our coordinates, we should choose the force field. For now, we don't have any single SMIRNOFF force field that can handle both proteins and small molecules; the Rosemary 3.0.0 force field will support this, but it's not yet ready. As an alternative, we'll combine the AMBER-compatible [Sage] small molecule force field with the SMIRNOFF port of AMBER ff14SB. Note that Sage also includes the TIP3P water model, which is appropriate for AMBER ff14SB too. +Now that we've prepared our coordinates, we should choose the force field. For now, we don't have any single SMIRNOFF force field that can handle both proteins and small molecules. The Rosemary line of force fields (starting with `openff-3.0.0.offxml`) is intended to do exactly this. A pre-release version, [`openff_no_water-3.0.0-alpha0.offxml`](https://github.com/openforcefield/openff-forcefields/blob/main/openforcefields/offxml/openff_no_water-3.0.0-alpha0.offxml), is already available for testing (if you use it, see [the release notes](https://github.com/openforcefield/openff-forcefields/releases/tag/2025.10.1)). There is no specific release date planned for the first full version, but it may be available later in 2026. As an alternative, we'll combine the AMBER-compatible [Sage] small molecule force field with the SMIRNOFF port of AMBER ff14SB. Note that Sage also includes the TIP3P water model, which is appropriate for AMBER ff14SB too. When we combine multiple SMIRNOFF force fields into one, we provide them in an order from general to specific. Sage includes parameters that could be applied to a protein, but they're general across all molecules; ff14SB's parameters are specific to proteins. Since the Toolkit always applies the last parameters that match a moiety, this order makes sure the right parameters get assigned. [Sage]: https://openforcefield.org/force-fields/force-fields/#sage
-⚠️ Warning: If your small molecule has an amino acid substructure in it, the specific patterns in the ff14SB force field will override the general ones from openff-2.2.1.offxml. This is the SMIRNOFF format being applied correctly, but some users may find this surprising, especially since terminal caps like ACE and NME are relatively small substructures and will sometimes appear in ligands. +⚠️ Warning: If your small molecule has an amino acid substructure in it, the specific patterns in the ff14SB force field will override the general ones from openff-2.3.0.offxml. This is the SMIRNOFF format being applied correctly, but some users may find this surprising, especially since terminal caps like ACE and NME are relatively small substructures and will sometimes appear in ligands.
@@ -236,12 +232,12 @@ When we combine multiple SMIRNOFF force fields into one, we provide them in an o from openff.toolkit import ForceField # Assemble the combined force field -sage_ff14sb = ForceField("openff-2.2.1.offxml", "ff14sb_off_impropers_0.0.3.offxml") +sage_ff14sb = ForceField("openff-2.3.0.offxml", "ff14sb_off_impropers_0.0.3.offxml") ``` -We now have a `Topology`, which stores the chemical information of the system, and a `ForceField`, which maps chemistry to force field parameters. To parametrize the system, we combine these two objects into an [`Interchange`], as discussed in this morning's session. +We now have a `Topology`, which stores the chemical information of the system, and a `ForceField`, which maps chemistry to force field parameters. To parametrize the system, we combine these two objects into an [`Interchange`], as discussed in the previous notebook. -An `Interchange` represents a completely parameterised molecular mechanics system. Partial charges are computed here according to the instructions in the force field, and this is where virtual sites required by the force field will be introduced. This all happens behind the scenes; all we have to do is combine an abstract chemical description with a force field. This makes it easy to change water model or force field, as the chemistry being modelled is completely independent of the model itself. +An `Interchange` represents a completely parameterised molecular mechanics system. Partial charges are computed here according to the instructions in the force field: Sage 2.3.0 assigns the ligand's charges with the AshGC graph neural network model (see section 4 of the previous notebook), while ff14SB supplies library charges for the protein, water, and ions. This is also where virtual sites required by the force field will be introduced. This all happens behind the scenes; all we have to do is combine an abstract chemical description with a force field. This makes it easy to change water model or force field, as the chemistry being modelled is completely independent of the model itself. [`Interchange`]: https://docs.openforcefield.org/projects/interchange/en/stable/_autosummary/openff.interchange.components.interchange.Interchange.html @@ -250,9 +246,9 @@ An `Interchange` represents a completely parameterised molecular mechanics syste interchange = sage_ff14sb.create_interchange(topology) ``` -*(This should take about a minute, largely because of the complexity of the AMBER protein force field port. In the future, this should be faster.)* +*(This should take well under a minute, most of it spent on the chemical perception required by the AMBER protein force field port. It used to be considerably slower: assigning the ligand's AM1-BCC charges was the bottleneck, and AshGC has removed it.)* -While that runs, let's recap. We've constructed a `Topology` out of a number of `Molecule` objects, each of which represents a particular chemical independent of any model details. The `Topology` then represents an entire chemical system, which in theory could be modelled in any number of ways. Our `Topology` also includes atom positions and box vectors, but if we thought that was too concrete for our use case we could leave them out and add them after parameterisation. +Before we simulate, let's recap. We've constructed a `Topology` out of a number of `Molecule` objects, each of which represents a particular chemical independent of any model details. The `Topology` then represents an entire chemical system, which in theory could be modelled in any number of ways. Our `Topology` also includes atom positions and box vectors, but if we thought that was too concrete for our use case we could leave them out and add them after parameterisation. Separately, we've constructed a `ForceField` by combining a general SMIRNOFF force field with a protein-specific SMIRNOFF force field. A SMIRNOFF force field is a bunch of rules for applying force field parameters to chemicals via SMARTS patterns. The force field includes everything needed to compute an energy: parameters, charges, functional forms, non-bonded methods and cutoffs, virtual sites, and so on. @@ -277,6 +273,21 @@ interchange.to_gromacs(prefix="complex") +```python +# Check the new gromacs output +!ls +``` + + complex.gro + complex.top + complex_pointenergy.mdp + protein_ligand_complex_parameterisation_and_md.ipynb + small_molecule_parameterisation.ipynb + topology.json + trajectory_gpu.dcd + + + All that remains is to tell OpenMM the details about how we want to integrate and record data for the simulation, and then to put everything together and run it! The steps are: 1. Configure and run the simulation @@ -288,7 +299,7 @@ All that remains is to tell OpenMM the details about how we want to integrate an ### 4.1 Configure and run the simulation -Here, we'll use a Langevin thermostat at 300 Kelvin and a 2 fs time step. We'll write the structure to disk every 10 steps. In contrast to the previous notebook, we'll add a MonteCarloBarostat to fix the pressure, while allowing the volume to fluctuate. Our simulation corresponds to the $NPT$ ensemble. +Here, we'll use a Langevin thermostat at 300 Kelvin and a 2 fs time step. We'll write the structure to disk every 50 steps. In contrast to the previous notebook, we'll add a MonteCarloBarostat to fix the pressure, while allowing the volume to fluctuate. Our simulation corresponds to the $NPT$ ensemble. ```python @@ -299,8 +310,8 @@ PRESSURE = 1 * openmm.unit.atmosphere FRICTION_COEFFICIENT = 1 / openmm.unit.picosecond TIMESTEP = 0.002 * openmm.unit.picoseconds -# Construct and configure a LangevinMiddleIntegrator at 300 K with an appropriate friction constant and time-step -integrator = openmm.LangevinMiddleIntegrator( +# Construct and configure a LangevinIntegrator at 300 K with an appropriate friction constant and time-step +integrator = openmm.LangevinIntegrator( TEMPERATURE, FRICTION_COEFFICIENT, TIMESTEP, @@ -348,10 +359,10 @@ describe_state( ) ``` - Original state has energy 14441301.95 kJ/mol with maximum force 1367209099.85 kJ/(mol nm) + Original state has energy 14441257.72 kJ/mol with maximum force 1367206715.93 kJ/(mol nm) - Minimized state has energy -434035.42 kJ/mol with maximum force 2426.16 kJ/(mol nm) + Minimized state has energy -435737.31 kJ/mol with maximum force 2458.94 kJ/(mol nm) ### 4.3 Run a short simulation @@ -468,7 +479,7 @@ ax.set_ylabel(r'RMSD (Å)') -![png](output_39_1.png) +![png](output_40_1.png) @@ -519,7 +530,7 @@ ax.set_ylabel(r'RMSD (Å)') -![png](output_41_1.png) +![png](output_42_1.png) @@ -627,7 +638,7 @@ df = fp.to_dataframe() ```
- ✏️ Exercise: Repeat this entire notebook using a ligand from docked to MCL-1 during this morning's session. (Hint: You'll need to convert the pdbqt files to sdf files using obabel, adding protons as appropriate for pH 7. This will look something like obabel docked_ligand.pdbqt -opdb | obabel -ipdb -osdf -p 7.0 -O docked_ligand.sdf. Make sure to use the docked coordinates! An example docked pdbqt file is provided at ../structures/docked_ligand.pdbqt) Is the binding pose stable? Are similar interactions formed by the docked ligand and the crystallographic ligand? Which do you think is likely to bind more strongly? What would be required to answer these questions robustly? + ✏️ Exercise: Repeat this entire notebook using a ligand docked to MCL-1. (Hint: You'll need to convert the pdbqt files to sdf files using obabel, adding protons as appropriate for pH 7. This will look something like obabel docked_ligand.pdbqt -opdb | obabel -ipdb -osdf -p 7.0 -O docked_ligand.sdf. Make sure to use the docked coordinates! An example docked pdbqt file is provided at ../structures/docked_ligand.pdbqt) Is the binding pose stable? Are similar interactions formed by the docked ligand and the crystallographic ligand? Which do you think is likely to bind more strongly? What would be required to answer these questions robustly?
@@ -644,16 +655,18 @@ df = fp.to_dataframe() * Using OpenMM, we never had to leave Python to set up the simulation. * With Interchange, using OpenMM, GROMACS, Amber or LAMMPS is simple! * MDAnalysis and ProLIF allows us to perform varied analyses of our trajectories. +* The Rosemary force field is likely coming soon, and will allow easy set-up of simulations with post-translationally modified proteins. Check out [this workshop](https://github.com/openforcefield/2026-virtual-workshops/blob/main/ptm/ptm-workshop.ipynb)! ## 7. There's Lots More to OpenFF! A variety of example notebooks for OpenFF software are provided [here](https://docs.openforcefield.org/en/latest/examples.html). A few which are particularly relevant are: +- [Simulating Post-Translationally Modified Proteins with the OpenFF Rosemary Alpha](https://github.com/openforcefield/2026-virtual-workshops/blob/main/ptm/ptm-workshop.ipynb) - [Host-guest systems](https://docs.openforcefield.org/en/latest/examples/openforcefield/openff-interchange/host-guest/host_guest.html) - [Protein-ligand-water systems with Interchange](https://docs.openforcefield.org/en/latest/examples/openforcefield/openff-interchange/protein_ligand/protein_ligand.html). This has a lot of overlap with the current notebook, but there are several extra details not covered here. ## 8. Beyond OpenFF -You can parameterise your complex and run molecular dynamics -- so what's next? If you're interested in quantiatively assessing the binding affinity of your ligand for your target, then [alchemical (and path-based) free energy calculations are the gold-standard method](https://www.nature.com/articles/s42004-023-01019-9). [Open Free Energy](https://openfree.energy/) is another [Open Molecular Software Foundation](https://omsf.io/) initiative, which develops open-source tools for binding free energy calculations. Head to their [tutorials](https://docs.openfree.energy/en/latest/tutorials/index.html) to learn more! However, these calculations are computationally demanding. If you're interested in a relatively fast (but relatively inaccurate) ranking of the binding affinities of a set of ligands, methods such as MM/GBSA may be appropriate. +You can parameterise your complex and run molecular dynamics -- so what's next? If you're interested in quantitatively assessing the binding affinity of your ligand for your target, then [alchemical (and path-based) free energy calculations are the gold-standard method](https://livecomsjournal.org/index.php/livecoms/article/view/v2i1e18378). [Open Free Energy](https://openfree.energy/) is another [Open Molecular Software Foundation](https://omsf.io/) initiative, which develops open-source tools for binding free energy calculations. However, these calculations are computationally demanding. If you're interested in a relatively fast (but relatively inaccurate) ranking of the binding affinities of a set of ligands, methods such as MM/GBSA may be appropriate. Affinity prediction methods based on deep learning such as Boltz-2 are appealingly fast, [but perform poorly on systems dissimilar to those they are trained on, and often show inflated performance on benchmarks due to data leakage](https://www.biorxiv.org/content/10.64898/2026.06.29.735309v1.abstract). diff --git a/notebooks-rendered/small_molecule_parameterisation.md b/notebooks-rendered/small_molecule_parameterisation.md index 7a01b94..9641901 100644 --- a/notebooks-rendered/small_molecule_parameterisation.md +++ b/notebooks-rendered/small_molecule_parameterisation.md @@ -11,9 +11,10 @@ This is the first of two jupyter notebooks on handling force fields using [Open | Action | Software| |--|--| +| [Go from SMILES to simulation in a few lines of code](#showcase) | OpenFF Toolkit, OpenFF Interchange, OpenMM | [Load and inspect a force field](#loading_ff) | OpenFF Toolkit | [Create a representation of your chemical system](#topology) | OpenFF Toolkit -| [Parameterise your system and run a quick simulation](#interchange) | OpenFF Interchange, OpenMM +| [Parameterise your system and run a quick simulation in water](#interchange) | OpenFF Interchange, OpenMM | [Rapidly assign partial charges with a graph neural network model](#gnn_charges) | OpenFF Toolkit, OpenFF NAGL Models | [Review what you've learnt](#summary) | | [Check out other OpenFF tutorials](#further_materials) | @@ -33,27 +34,102 @@ Most of this material was adapted from the [2023 CCPBioSim Workshop Open Force F + +## 0. You can go from SMILES to simulation in a few lines of code + + +```python +# Go from SMILES -> simulation input with OpenFF +from openff.toolkit import ForceField, Molecule, Topology + +molecule = Molecule.from_smiles("CC(=O)Nc1ccc(cc1)O") +molecule.generate_conformers(n_conformers=1) +topology = Topology.from_molecules([molecule]) + +force_field = ForceField("openff-2.3.0.offxml") +interchange = force_field.create_interchange(topology) +interchange.minimize() + +openmm_system = interchange.to_openmm_system() +openmm_topology = interchange.to_openmm_topology() +openmm_positions = interchange.positions.to_openmm() +``` + + +```python +# Run the simulation with OpenMM +import openmm +import openmm.app + +temperature = 298.15 * openmm.unit.kelvin +friction_coefficient = 1.0 / openmm.unit.picosecond +step_size = 2.0 * openmm.unit.femtosecond + +simulation = openmm.app.Simulation( + openmm_topology, + openmm_system, + openmm.LangevinIntegrator(temperature, friction_coefficient, step_size), +) +simulation.context.setPositions(openmm_positions) +simulation.context.setVelocitiesToTemperature(simulation.integrator.getTemperature()) + +simulation.reporters.append( + openmm.app.DCDReporter(file="trajectory_showcase.dcd", reportInterval=100) +) +simulation.step(10000) +``` + + +```python +# Load the trajectory with MDAnalysis and visualise with nglview +import MDAnalysis as mda +import nglview + +u = mda.Universe(openmm_topology, "trajectory_showcase.dcd") + +view = nglview.show_mdanalysis(u) +view +``` + + + + + + /opt/conda/envs/openff-env/lib/python3.12/site-packages/MDAnalysis/coordinates/DCD.py:171: DeprecationWarning: DCDReader currently makes independent timesteps by copying self.ts while other readers update self.ts inplace. This behavior will be changed in 3.0 to be the same as other readers. Read more at https://github.com/MDAnalysis/mdanalysis/issues/3889 to learn if this change in behavior might affect you. + warnings.warn("DCDReader currently makes independent timesteps" + + + + NGLWidget(max_frame=99) + + +That's it! You've run a vacuum simulation for paracetamol. Below and in the next notebook, we'll go into more detail on each of the steps in the OpenFF cell and show how you can set up more complex systems, but this is mainly for your understanding and you rarely need much more code than shown above. + ## 1. Force fields are specified in `.offxml` files and can be loaded with the `ForceField` class -OpenFF's force fields use the The SMIRKS Native Open Force Field (SMIRNOFF) [specification](https://openforcefield.github.io/standards/standards/smirnoff/) and are conventionally encoded in `.offxml` files. The spec fully describes the contents of a SMIRNOFF force field, how parameters should be applied, and several other important usage details. You could implement a SMIRNOFF engine in your own code, but conveniently the OpenFF Toolkit already provides this and a handful of utilities. Let's load up the latest OpenFF small molecule force field, OpenFF 2.2.1, and inspect its contents! This force field shares the code name "Sage" with all other force fields with the same major version number (2.x.x). +Let's dive into the details of what went on above. Here's a summary of how data flows through a workflow utilising OpenFF tools -- the OpenFF toolkit allows you to create `Molecule` and `ForceField` objects, which get combined into an `Interchange` object, which contains all the information needed to start a simulation. From there, you can create input for the simulation engine of your choice: + +Description of image + +Let's start with the `.offxml` force field file. OpenFF's force fields use the SMIRKS Native Open Force Field (SMIRNOFF) [specification](https://openforcefield.github.io/standards/standards/smirnoff/). The spec fully describes the contents of a SMIRNOFF force field, how parameters should be applied, and several other important usage details. You could implement a SMIRNOFF engine in your own code, but conveniently the OpenFF Toolkit already provides this and a handful of utilities. Let's load up the latest OpenFF small molecule force field, OpenFF 2.3.0, and inspect its contents! This force field shares the code name "Sage" with all other force fields with the same major version number (2.x.x). ```python from openff.toolkit import ForceField -sage = ForceField("openff-2.2.1.offxml") +sage = ForceField("openff-2.3.0.offxml") sage ``` - + -If you'd like to see the raw file on disk that's being parsed, [here's the file on GitHub](https://github.com/openforcefield/openff-forcefields/blob/main/openforcefields/offxml/openff-2.2.1.offxml). +If you'd like to see the raw file on disk that's being parsed, [here's the file on GitHub](https://github.com/openforcefield/openff-forcefields/blob/main/openforcefields/offxml/openff-2.3.0.offxml). Each section of a force field is stored in memory within `ParameterHandler` objects, which can be looked up with brackets (just like looking up values in a dictionary): @@ -65,13 +141,13 @@ vdw_handler = sage["vdW"] vdw_handler ``` - ['Constraints', 'Bonds', 'Angles', 'ProperTorsions', 'ImproperTorsions', 'vdW', 'Electrostatics', 'LibraryCharges', 'ToolkitAM1BCC'] + ['Constraints', 'Bonds', 'Angles', 'ProperTorsions', 'ImproperTorsions', 'vdW', 'Electrostatics', 'LibraryCharges', 'NAGLCharges'] - + @@ -88,7 +164,7 @@ print(f"vdw_handler parameters: {vdw_handler.parameters}") vdw_handler cutoff: 9.0 angstrom vdw_handler combining rules: Lorentz-Berthelot vdw_handler scale14: 0.5 - vdw_handler parameters: [, , , , , , , , , , , , , , , , , , , , , , , , , , , , , , , , , , , , , ] + vdw_handler parameters: [, , , , , , , , , , , , , , , , , , , , , , , , , , , , , , , , , , , , , ] From here you can inspect all the way down to individual parameters, which are stored in custom objects (in this case, `vdWType`). Let's look at the type with id `n16`, which looks like a generic carbon with four bonded neighbors: @@ -102,7 +178,7 @@ vdw_type - + @@ -110,7 +186,7 @@ Note that the type contains both the physical parameters (sigma and epsilon, for The toolkit uses these SMIRKS patterns and direct chemical perception to assign parameters to particular atoms (or bonds, angles, etc.). -We'll use OpenFF 2.2.1 for the remainder of this tutorial, but you can learn more about this and other SMIRNOFF force fields below: +We'll use OpenFF 2.3.0 for the remainder of this tutorial. This is OpenFF's latest small molecule force field and is a leading open-source small molecule force field which [performs comparably to other open-source force fields](https://doi.org/10.1021/acs.jctc.3c00039). You can learn more about this and other SMIRNOFF force fields below:
Click here to learn about available and planned SMIRNOFF force fields @@ -132,7 +208,7 @@ The Parsley line of force fields (`openff-1.y.z.offxml`) was OpenFF's [first ful The Sage line of force fields (`openff-2.y.z.offxml`) continued the process of fitting to more (and more diverse) QM datasets, but also included a re-fit of the Lennard-Jones parameters. Small molecule geometries and energies [improved, in general,](https://openforcefield.org/community/news/general/sage2.0.0-release/) significantly over Parsley. These improvements notably transferred to protein-ligand binding free energies despite Sage not being specifically fit to them. For more, see the [associated paper](https://pubs.acs.org/doi/10.1021/acs.jctc.3c00039). -[Subsequent releases](https://github.com/openforcefield/openff-forcefields/releases) used different fitting procedures and tweaks to parameter typing to improve performance and address issues with several specific chemistries. The latest release, **Sage 2.2.1 (`openff-2.2.1.offxml`) is the recommended force field for small molecule studies.** +[Subsequent releases](https://github.com/openforcefield/openff-forcefields/releases) used different fitting procedures and tweaks to parameter typing to improve performance and address issues with several specific chemistries. Notably, Sage 2.3.0 includes fast graph neural network charge assignment with AshGC, which is discussed later in this notebook. This charge model is trained to reproduce AM1-BCC charges, but scales 𝒪(N) rather than the 𝒪(N2) of common AM1-BCC implementations, making it suitable for large (>> 100 atoms) molecules. The latest release, **Sage 2.3.0 (`openff-2.3.0.offxml`) is the recommended force field for small molecule studies.** ## Ports @@ -150,7 +226,14 @@ Existing main-line OpenFF force fields are fit against TIP3P water, so use of ot ## ff14SB -OpenFF, in collaboration with Dave Cerutti of the Amber community, created a port of [ff14SB](https://pubs.acs.org/doi/10.1021/acs.jctc.5b00255), a popular Amber protein force field. There are some small numerical differences with how improper torsions are evaluated, but all other terms reproduce a canonical Amber source to high accuracy. **This is the only protein force field currently in SMIRONOFF (`.offxml`) format** and therefore the current recommendation for use with proteins. Primarily for technical reasons, porting other Amber force fields is not planned. +OpenFF, in collaboration with Dave Cerutti of the Amber community, created a port of [ff14SB](https://pubs.acs.org/doi/10.1021/acs.jctc.5b00255), a popular Amber protein force field. There are some small numerical differences with how improper torsions are evaluated, but all other terms reproduce a canonical Amber source to high accuracy. **This is the only protein force field currently in SMIRNOFF (`.offxml`) format** and therefore the current recommendation for use with proteins. Primarily for technical reasons, porting other Amber force fields is not planned. + +## Rosemary Alpha + +A future line of force fields from OpenFF (code name "Rosemary", starting with `openff-3.0.0.offxml`) is intended to handle small molecules and biopolymers in a _self-consistent_ manner. This is exciting as it will streamline simulations of proteins with non-canonical amino acids! See the workshop [Simulating Post-Translationally Modified Proteins with the OpenFF Rosemary Alpha](https://github.com/openforcefield/2026-virtual-workshops/blob/main/ptm/ptm-workshop.ipynb). The first release will handle proteins, but future versions may cover nucleic acids. The performance, depending on the metrics used, is hoped to be comparable with existing Amber-family protein force fields. + +A pre-release version of Rosemary is available for testing as [`openff_no_water-3.0.0-alpha0.offxml`](https://github.com/openforcefield/openff-forcefields/blob/main/openforcefields/offxml/openff_no_water-3.0.0-alpha0.offxml). If you use it, see [the release notes](https://github.com/openforcefield/openff-forcefields/releases/tag/2025.10.1). + ## Non-main-line force fields @@ -165,19 +248,11 @@ https://github.com/jthorton/de-forcefields ## From OpenFF ### Rosemary -A future line of force fields from OpenFF (code name "Rosemary", starting with `openff-3.0.0.offxml`) is intended to handle small molecules and biopolymers in a _self-consistent_ manner. The first release is expected to handle proteins, but future versions may cover nucleic acids. The performance, depending on the metrics used, is hoped to be comparable with existing Amber-family protein force fields. - -There is no specific release date planned for Rosemary, but it may be available in 2026 (a beta release candidate may also be publically available prior to the full release). - -### Graph net charge assignment - -TODO: UPDAATE AND MENTION 2.3 - -The Sage 2.3.0 release is expected imminently and will include graph-convolutional neutral network (GCNN)-based charge assignment using [NAGL](https://github.com/openforcefield/openff-nagl) by default. The charge model is trained to reproduce AM1-BCC charges without the typical $O(N^3)$ scaling, making it suitable for large (>> 100 atoms) molecules). The [second release candidate](https://github.com/openforcefield/openff-forcefields/blob/main/openforcefields/offxml/openff-2.3.0-rc2.offxml) (which may or may not become the final version) is already available for you to try! +There is no specific release date planned for the first full version of Rosemary, but it may be available in late 2026. ### Virtual sites -Another release from OpenFF may include some virtual site parameters with off-center charges. No release date is planned, but the most of the supporting infrastructure is currently in place and some early studies have shown promise for better representing electrostatics of chemistries such as halogens and aromatic nitrogens. +Another release from OpenFF may include some virtual site parameters with off-center charges. No release date is planned, but most of the supporting infrastructure is currently in place and some early studies have shown promise for better representing electrostatics of chemistries such as halogens and aromatic nitrogens. ## From you! @@ -187,7 +262,7 @@ Anybody can write a SMIRNOFF force field! This workshop doesn't have time to cov ## 2. The `Topology` class represents a chemical system containing one or more `Molecule`s -Now we've loaded our desired force field (OpenFF 2.2.1), we need to specify the chemical system we want to assign force field parameters to ("parameterise"). Our system will be represented by a `Topology`, which we will build from one or more `Molecule`s. +Now we've loaded our desired force field (OpenFF 2.3.0), we need to specify the chemical system we want to assign force field parameters to ("parameterise"). Our system will be represented by a `Topology`, which we will build from one or more `Molecule`s. As a simple example, let's build a `Topology` containing an small molecule with some features which illustrate how parameters are applied according to SMIRKS matches. We'll use the crotonate anion, but you could draw any molecule you like and convert it to a SMILES string using tools like ChemDraw and [MolView](https://molview.org/). @@ -201,15 +276,7 @@ molecule - - - /opt/conda/envs/openff-env/lib/python3.12/site-packages/nglview/__init__.py:12: UserWarning: pkg_resources is deprecated as an API. See https://setuptools.pypa.io/en/latest/pkg_resources.html. The pkg_resources package is slated for removal as early as 2025-11-30. Refrain from using this package or pin to Setuptools<81. - import pkg_resources - - - - -![svg](output_12_2.svg) +![svg](output_17_0.svg) @@ -282,7 +349,7 @@ topology_with_water.molecule(0), topology_with_water.molecule(1), topology_with_
- ✏️ Exercise: Build a Topology containing an MCL-1 ligand. Create the Molecule from an SDF file (take a look at the docstring of Molecule to see how this can be done). Also, see Molecule cookbook for all the ways to make a Molecule. The crystallographic MCL-1 ligand from PDB ID 6o6f is provided at ../structures/606f_ligand.sdf. Note that you don't need to use get_data_file_path as we already know the path. + ✏️ Exercise: Build a Topology containing an MCL-1 ligand. Create the Molecule from an SDF file (take a look at the docstring of Molecule to see how this can be done, noting that you can just pass the sdf path given below and don't need the get_data_file_path function). Also, see Molecule cookbook for all the ways to make a Molecule. The crystallographic MCL-1 ligand from PDB ID 6o6f is provided at ../structures/6o6f_ligand.sdf. Note that you don't need to use get_data_file_path as we already know the path.
@@ -293,7 +360,7 @@ mcl1_mol = Molecule("../structures/6o6f_ligand.sdf") top = mcl1_mol.to_topology() ``` -We will cover creating a topology for a protein-ligand complex this afternoon. +We will cover creating a topology for a protein-ligand complex in the next notebook. ## 3. `Interchange` objects contain fully parameterised systems with all the information needed to start a simulation @@ -302,15 +369,12 @@ Now we've specified our force field and our chemical system using classes from t To do this, we'll use the `Interchange` class from the OpenFF Interchange package, which stores a fully-parameterised molecular system and provides methods to write out simulation-ready input files for a number of software packages. They key objective of Interchange is to provide an intermediate inspectable state after parameterisation and before conversion to an engine-specific format. For most users, an `Interchange` forms the bridge between the OpenFF ecosystem and their simulation software of choice. The current focus is applying SMIRNOFF force fields to chemical topologies and exporting the result to engines preferred by our users. In order of stability, OpenMM, GROMACS, Amber, and LAMMPS are supported. Future development may include support for CHARMM and other engines. -Below is a summary of how data flows through a workflow utilising OpenFF tools, including where Interchange sits in the flow. - -Description of image - First, let's recreate our `molecule` and `topology` in case you overwrote them during the previous exercises: ```python molecule = Molecule.from_smiles("C/C=C/C(=O)[O-]") +molecule.generate_conformers(n_conformers=1) topology = molecule.to_topology() ``` @@ -440,7 +504,7 @@ SVG(mol_with_atom_index(molecule)) -![svg](output_40_0.svg) +![svg](output_44_0.svg) @@ -469,7 +533,7 @@ collection.key_map -We can see that the C=C bond (indices (1,2)) is associated with a potential key with the SMIRKS pattern `[#6X3:1]=[#6X3:2]` (specifying any two carbons bonded to 3 atoms connected by a double bond). Note that the (1,0) C-C bond is matched by the SMRIKS `[#6X3:1]-[#6X3:2]`, which specifies the atoms in the same way, showing that the parameters have been assigned by directly using information about the bond. This contrasts to traditional atom typing approaches, where information about the bond would be implicitly encoded in the atom types used to assign the parameters. Another example of this "direct chemical perception" is the assignment of the carboxylate carbon-oxygen bond parameters, which only match (triply-connected carbon) - (singly-connnected oxygen) bonds when the carbon is bonded to another singly-connected oxygen. +We can see that the C=C bond (indices (1,2)) is associated with a potential key with the SMIRKS pattern `[#6X3:1]=[#6X3:2]` (specifying any two carbons each bonded to 3 atoms and connected by a double bond). Note that the (1,0) C-C bond is matched by the SMRIKS `[#6X3:1]-[#6X3:2]`, which specifies the atoms in the same way, showing that the parameters have been assigned by directly using information about the bond. This contrasts to traditional atom typing approaches, where information about the bond would be implicitly encoded in the atom types used to assign the parameters. Another example of this "direct chemical perception" is the assignment of the carboxylate carbon-oxygen bond parameters, which only match (triply-connected carbon) - (singly-connnected oxygen) bonds when the carbon is bonded to another singly-connected oxygen. To see the actual parmeters specified for this bond, we can look up the `Potential` objects using the `PotentialKey`s. @@ -480,16 +544,16 @@ for topology_key, potential_key in collection.key_map.items(): print(f"{topology_key} -> {potential}") ``` - atom_indices=(0, 1) bond_order=None -> parameters={'k': , 'length': } map_key=None - atom_indices=(0, 6) bond_order=None -> parameters={'k': , 'length': } map_key=None - atom_indices=(0, 7) bond_order=None -> parameters={'k': , 'length': } map_key=None - atom_indices=(0, 8) bond_order=None -> parameters={'k': , 'length': } map_key=None - atom_indices=(1, 2) bond_order=None -> parameters={'k': , 'length': } map_key=None - atom_indices=(1, 9) bond_order=None -> parameters={'k': , 'length': } map_key=None - atom_indices=(2, 3) bond_order=None -> parameters={'k': , 'length': } map_key=None - atom_indices=(2, 10) bond_order=None -> parameters={'k': , 'length': } map_key=None - atom_indices=(3, 4) bond_order=None -> parameters={'k': , 'length': } map_key=None - atom_indices=(3, 5) bond_order=None -> parameters={'k': , 'length': } map_key=None + atom_indices=(0, 1) bond_order=None -> parameters={'k': , 'length': } map_key=None + atom_indices=(0, 6) bond_order=None -> parameters={'k': , 'length': } map_key=None + atom_indices=(0, 7) bond_order=None -> parameters={'k': , 'length': } map_key=None + atom_indices=(0, 8) bond_order=None -> parameters={'k': , 'length': } map_key=None + atom_indices=(1, 2) bond_order=None -> parameters={'k': , 'length': } map_key=None + atom_indices=(1, 9) bond_order=None -> parameters={'k': , 'length': } map_key=None + atom_indices=(2, 3) bond_order=None -> parameters={'k': , 'length': } map_key=None + atom_indices=(2, 10) bond_order=None -> parameters={'k': , 'length': } map_key=None + atom_indices=(3, 4) bond_order=None -> parameters={'k': , 'length': } map_key=None + atom_indices=(3, 5) bond_order=None -> parameters={'k': , 'length': } map_key=None So our C=C bond (indices (1,2)) has a force constant of 904 kcal mol-1 Å-2 and an equilibrium bond length of 1.37 Å. Note that the [`ForceField.label_molecules`](https://docs.openforcefield.org/projects/toolkit/en/stable/api/generated/openff.toolkit.typing.engines.smirnoff.ForceField.html#openff.toolkit.typing.engines.smirnoff.ForceField.label_molecules) method is also useful for checking which parameters will be applied to your molecule. @@ -531,31 +595,31 @@ collection.key_map bond1_indices, bond2_indices = (3,4), (3,5) smiles = {"anion": "C/C=C/C(=O)[O-]", "neutral": "C/C=C/C(=O)O"} +# Use fresh variable names so we don't overwrite `molecule`, `topology` and `interchange` for name, smiles in smiles.items(): print(f"\n{name} ({smiles}):") - molecule = Molecule.from_smiles(smiles) - topology = molecule.to_topology() - interchange = Interchange.from_smirnoff( + protonation_state = Molecule.from_smiles(smiles) + protonation_state_interchange = Interchange.from_smirnoff( force_field=sage, - topology=topology, + topology=protonation_state.to_topology(), ) - collection = interchange.collections["Bonds"] - bond1 = collection[bond1_indices] - bond2 = collection[bond2_indices] + bonds = protonation_state_interchange.collections["Bonds"] + bond1 = bonds[bond1_indices] + bond2 = bonds[bond2_indices] print(f" Bond {bond1_indices}: {bond1}") print(f" Bond {bond2_indices}: {bond2}") ``` anion (C/C=C/C(=O)[O-]): - Bond (3, 4): parameters={'k': , 'length': } map_key=None - Bond (3, 5): parameters={'k': , 'length': } map_key=None + Bond (3, 4): parameters={'k': , 'length': } map_key=None + Bond (3, 5): parameters={'k': , 'length': } map_key=None neutral (C/C=C/C(=O)O): - Bond (3, 4): parameters={'k': , 'length': } map_key=None - Bond (3, 5): parameters={'k': , 'length': } map_key=None + Bond (3, 4): parameters={'k': , 'length': } map_key=None + Bond (3, 5): parameters={'k': , 'length': } map_key=None Finally, `interchange.box` and `interchange.velocities` are `None`, although `interchange.positions` is populated because we passed a topology with a molecule that had a defined conformer, so `from_smirnoff` set atomic positions from this information: @@ -568,7 +632,19 @@ interchange.positions, interchange.box, interchange.velocities - (None, None, None) + (, + None, + None) @@ -594,7 +670,7 @@ for atom in water.atoms: # the volume. topology = pack_box( molecules=[solute, water], - number_of_copies=[1, 1000], + number_of_copies=[1, 1400], box_vectors=3.5 * UNIT_CUBE * unit.nanometer, ) @@ -606,11 +682,11 @@ interchange.topology.n_atoms, interchange.box, interchange.positions.shape - (3012, + (4212, , - (3012, 3)) + (4212, 3)) @@ -623,7 +699,7 @@ interchange.to_amber(prefix="ligand") /opt/conda/envs/openff-env/lib/python3.12/site-packages/openff/interchange/components/mdconfig.py:502: UserWarning: Ambiguous failure while processing constraints. Constraining h-bonds as a stopgap. warnings.warn( - /opt/conda/envs/openff-env/lib/python3.12/site-packages/openff/interchange/components/mdconfig.py:430: 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. + /opt/conda/envs/openff-env/lib/python3.12/site-packages/openff/interchange/components/mdconfig.py:432: 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( @@ -638,17 +714,18 @@ interchange.to_amber(prefix="ligand") complex_pointenergy.mdp topology.json ligand.inpcrd trajectory.dcd ligand.prmtop trajectory_gpu.dcd - ligand_pointenergy.in + ligand_pointenergy.in trajectory_showcase.dcd -Here, we'll export to OpenMM and run a short simulation directly from the noteboook. We can create an OpenMM `Simulation` object from the `Interchange` and run for a specified wall clock time using `runForClockTime` (the simluation time will depend on how quickly it runs on your machine). We keep the volume ($V$), number of particles ($N$), and average temperature ($T$) (using the LangevinMiddleIntegrator) constant and the simulation corresponds to the $NVT$ ensemble. +Here, we'll export to OpenMM and run a short simulation directly from the noteboook. We can create an OpenMM `Simulation` object from the `Interchange` and run for a specified wall clock time using `runForClockTime` (the simluation time will depend on how quickly it runs on your machine). We keep the volume ($V$), number of particles ($N$), and average temperature ($T$) (using the LangevinIntegrator) constant and the simulation corresponds to the $NVT$ ensemble. ```python import openmm +import openmm.app import openmm.unit from openff.interchange import Interchange -import mdtraj +import MDAnalysis as mda import nglview @@ -658,7 +735,7 @@ def run_openmm( trajectory_name: str = "small_mol_solvated.dcd", ): simulation = interchange.to_openmm_simulation( - integrator=openmm.LangevinMiddleIntegrator( + integrator=openmm.LangevinIntegrator( 300 * openmm.unit.kelvin, 1 / openmm.unit.picosecond, 0.002 * openmm.unit.picoseconds, @@ -676,12 +753,10 @@ def visualise_traj( topology: Topology, filename: str = "small_mol_solvated.dcd" ) -> nglview.NGLWidget: """Visualise a trajectory using nglview.""" - traj = mdtraj.load( - filename, - top=mdtraj.Topology.from_openmm(topology.to_openmm()), - ) - view = nglview.show_mdtraj(traj) + u = mda.Universe(topology.to_openmm(), filename) + + view = nglview.show_mdanalysis(u) view.add_representation("licorice", selection="water") return view @@ -691,8 +766,12 @@ run_openmm(interchange) visualise_traj(interchange.topology) ``` + /opt/conda/envs/openff-env/lib/python3.12/site-packages/MDAnalysis/coordinates/DCD.py:171: DeprecationWarning: DCDReader currently makes independent timesteps by copying self.ts while other readers update self.ts inplace. This behavior will be changed in 3.0 to be the same as other readers. Read more at https://github.com/MDAnalysis/mdanalysis/issues/3889 to learn if this change in behavior might affect you. + warnings.warn("DCDReader currently makes independent timesteps" + + - NGLWidget(max_frame=49) + NGLWidget(max_frame=25)
@@ -702,116 +781,79 @@ visualise_traj(interchange.topology) ## 4. Graph Neural Networks Allow Fast Assignment of Partial Charges -You might notice that [Sage](https://github.com/openforcefield/openff-forcefields/blob/main/openforcefields/offxml/openff-2.2.1.offxml) doesn't contain tabulated charges for most atomic environments in the way it does for all other terms in the force field. Instead, it specifies: +You might notice that Sage force fields don't contain tabulated charges for most atomic environments in the way they do for all other terms in the force field. For example, [Sage 2.2.1](https://github.com/openforcefield/openff-forcefields/blob/main/openforcefields/offxml/openff-2.2.1.offxml) instead specifies: ``` ``` -which means that partial charges will be calculated using the common AM1-BCC method. Charges from a semi-empirical quantum chemistry calculation (Austin Model 1) are corrected (bond charge correction) to approximate charges obtained by fitting to the electrostatic potential at the HF/6-31G* level (see [Jakalian et al.](https://onlinelibrary.wiley.com/doi/10.1002/(SICI)1096-987X(20000130)21:2%3C132::AID-JCC5%3E3.0.CO;2-P)). Unfortunately, AM1-BCC scales 𝒪(N2) in the number of atoms N, making it prohibitively slow for large molecules and biopolymers. +which means that partial charges will be calculated using the common AM1-BCC method. Charges from a semi-empirical quantum chemistry calculation (Austin Model 1) are corrected (bond charge correction) to approximate charges obtained by fitting to the electrostatic potential at the HF/6-31G* level (see [Jakalian et al.](https://onlinelibrary.wiley.com/doi/10.1002/(SICI)1096-987X(20000130)21:2%3C132::AID-JCC5%3E3.0.CO;2-P)). Unfortunately, parameterisation with AM1-BCC using OpenEye or AmberTools scales 𝒪(N2) in the number of atoms N, making it prohibitively slow for large molecules and biopolymers. -Methods which assign partial charges using graph neural networks offer rapid assignment with better scaling. They also offer the possibility of going beyond traditionally affordable QM levels of theory by training to quickly reproduce charges from expensive calculations. For example, [EspalomaCharge](https://pubs.acs.org/doi/full/10.1021/acs.jpca.4c01287) is fit to AM1-BCC charges and offers 𝒪(N2) scaling, while [Adams et al.](https://chemrxiv.org/engage/chemrxiv/article-details/6839c94c3ba0887c33d2cd8e) trained models to reproduce atoms-in-molecules charges and electrostatic potentials obtained at a high level of theory. Here, we'll use OpenFF's [AshGC](https://zenodo.org/records/15770227/files/AshGC_methods_2025-06-30.pdf?download=1) model, which is trained to reproduce AM1-BCC charges. +Methods which assign partial charges using graph neural networks offer rapid assignment with better scaling. They also offer the possibility of going beyond traditionally affordable QM levels of theory by training to quickly reproduce charges from expensive calculations. For example, OpenFF's [AshGC](https://doi.org/10.1021/acs.jctc.6c00169) model is fit to AM1-BCC charges and offers 𝒪(N) scaling, while [Adams et al.](https://doi.org/10.1021/acs.jctc.5c01520) trained models to reproduce atoms-in-molecules charges and electrostatic potentials obtained at a high level of theory. AshGC is used by [Sage 2.3.0](https://github.com/openforcefield/openff-forcefields/blob/main/openforcefields/offxml/openff-2.3.0.offxml) -- if you inspect the file, you'll see: +``` + +``` +where "openff-gnn-am1bcc-1.0.0.pt" is the AshGC model. -
- ⚠️ OpenFF 2.2.1 has not been explicitly trained and validated with AshGC charges. However, the 2.3.0 release will be, and is expected imminently. AshGC charges will be used as default and will be specified in the .offxml file, so there will be no need to call Molecule.assign_partial_charges as shown below. -
+The GNN charge model is the main difference between Sage 2.2.1 and 2.3.0. Here, we'll compare the parameterisation speed and charges obtained using each force field. ```python -from openff.toolkit import Molecule, ForceField -from openff.toolkit.utils.nagl_wrapper import NAGLToolkitWrapper - -# Disable RDKit warnings to avoid misleading NAGL warnings -# (see https://github.com/openforcefield/openff-nagl/issues/198) -from rdkit import RDLogger -RDLogger.DisableLog('rdApp.*') +from openff.toolkit import Molecule -# OpenFF NAGL store models as PyTorch files. -ASH_GC_MODEL = "openff-gnn-am1bcc-0.1.0-rc.3.pt" molecule = Molecule("../structures/6o6f_ligand.sdf") ``` - -```python -molecule.assign_partial_charges? -``` - -First, let's assign charges the traditional way with AM1-BCC and check how long this takes... +First, let's parameterise with Sage 2.2.1, which uses the traditional AM1-BCC model, and check how long this takes... ```python %%time -molecule_am1bcc = Molecule(molecule) -molecule_am1bcc.assign_partial_charges( - partial_charge_method="am1bcc", -) +sage221 = ForceField("openff-2.2.1.offxml") +interchange_sage221 = Interchange.from_smirnoff(force_field=sage221, topology=molecule.to_topology()) ``` - CPU times: user 64.1 ms, sys: 7.24 ms, total: 71.3 ms - Wall time: 22.6 s + CPU times: user 318 ms, sys: 1.66 ms, total: 319 ms + Wall time: 30.8 s -Now, let's try AshGC +Note that repeating these cells will show much faster assignment as partial charges are cached for a given molecule and charge method. + +Now, let's try Sage 2.3.0, which uses AshGC charges: ```python %%time -molecule_ashgc = Molecule(molecule) -molecule_ashgc.assign_partial_charges( - partial_charge_method=ASH_GC_MODEL, - toolkit_registry=NAGLToolkitWrapper(), -) +sage230 = ForceField("openff-2.3.0.offxml") +interchange_sage230 = Interchange.from_smirnoff(force_field=sage230, topology=molecule.to_topology()) ``` - CPU times: user 1.13 s, sys: 32.9 ms, total: 1.17 s - Wall time: 1.1 s - - -Finally, let's create an `Interchange` with our AshGC charges, making sure to specify `charge_from_molecules` so that we don't replace them with `AM1BCC` charges: - + CPU times: user 1.84 s, sys: 70.3 ms, total: 1.91 s + Wall time: 1.56 s -```python -# normally when we call `ForceField.create_interchange` or `ForceField.create_openmm_system`, the toolkit will call -# AMBERTools or OEChem to assign partial charges, since that's what's in the force field file. A future OpenFF release -# which uses NAGL for charge assignment will encode this instruction in the force field file itself, but until that we -# can use the `charge_from_molecules` argument to tell it to use the charges that we just assigned# for more, see: -# https://docs.openforcefield.org/projects/toolkit/en/stable/api/generated/openff.toolkit.typing.engines.smirnoff.ForceField.html#openff.toolkit.typing.engines.smirnoff.ForceField.create_openmm_system -interchange = sage.create_interchange( - molecule_ashgc.to_topology(), - charge_from_molecules=[molecule_ashgc], -) -```
- ✏️ Exercise: Compare the charges obtained with AM1-BCC and AshGC by looking at the Molecule.partial_charges attribute. How big are these differences on average? What is the largest difference? Which atom are these on? The np.max function may be useful. + ✏️ Exercise: Compare the charges obtained with AM1-BCC and AshGC by inspecting the electrostatics collection of each Interchange object (see Section 3). How big are these differences on average? What is the largest difference? Which atom are these on? The np.max function may be useful.
```python -# Compare charges assigned with AM1-BCC and AshGC... +# Compare charges assigned with AM1-BCC (Sage 2.2.1) and AshGC (Sage 2.3.0)... import numpy as np -print(f"AM1 BCC charges: {molecule_am1bcc.partial_charges}") -print(f"AshGC charges: {molecule_ashgc.partial_charges}") -differences = molecule_am1bcc.partial_charges - molecule_ashgc.partial_charges -differences_by_atom_index = {idx: diff.magnitude for idx, diff in enumerate(differences)} -print(f"Differences by atom index: {differences_by_atom_index}") +charges_am1bcc = np.array([c.m for c in interchange_sage221["Electrostatics"].charges.values()]) +charges_ashgc = np.array([c.m for c in interchange_sage230["Electrostatics"].charges.values()]) -max_difference = np.max(np.abs(differences.magnitude)) -print(f"Max difference: {max_difference} e") +differences = charges_am1bcc - charges_ashgc +print(f"Mean absolute difference: {np.mean(np.abs(differences)):.4f} e") -mean_difference = np.mean(np.abs(differences.magnitude)) -print(f"Mean absolute difference: {mean_difference} e") - -# Get the atom index with the largest difference -atom_index = np.argmax(np.abs(differences.magnitude)) -print(f"Largest absolute difference is for atom index {atom_index}, which is a {molecule_ashgc.atoms[atom_index].symbol} atom") +max_index = int(np.argmax(np.abs(differences))) +print( + f"Largest absolute difference is {np.abs(differences[max_index]):.4f} e, " + f"for atom index {max_index}, which is a {molecule.atoms[max_index].symbol} atom" +) ``` - AM1 BCC charges: [-0.12740000000000004 -0.8283 -0.8283 -0.33690000000000003 -0.633 0.14839999999999998 0.19979999999999998 0.20379999999999998 -0.09940000000000004 -0.13500000000000004 -0.15200000000000002 -0.09040000000000004 -0.09040000000000004 -0.11100000000000004 -0.07100000000000004 -0.03410000000000003 -0.12800000000000003 -0.13800000000000004 -0.07840000000000004 0.9082 0.012399999999999965 0.04109999999999996 0.059599999999999966 -0.11570000000000004 -0.11560000000000004 -0.09030000000000005 -0.053300000000000035 -0.03870000000000003 0.04569999999999996 0.04569999999999996 0.056699999999999966 0.056699999999999966 0.03769999999999996 0.03769999999999996 0.039199999999999964 0.039199999999999964 0.14199999999999996 0.12199999999999996 0.049699999999999966 0.049699999999999966 0.049699999999999966 0.049699999999999966 0.15599999999999997 0.18099999999999997 0.047699999999999965 0.047699999999999965 0.13699999999999998 0.16099999999999998 0.04469999999999996 0.04469999999999996 0.08069999999999995] elementary_charge - AshGC charges: [-0.10291757708524957 -0.8212745068046976 -0.8212745068046976 -0.3278780457947184 -0.656814920661204 0.13626032793784842 0.22055995190406547 0.2206000658299993 -0.09947414970134988 -0.14171894168590798 -0.17288624024128213 -0.09598572826122537 -0.09598569845890298 -0.0818095085594584 -0.1024471399757792 -0.04326643323635354 -0.12058025872444406 -0.08199206268524423 -0.07742916321491494 0.9052312495734762 0.027429900186903337 0.0331591638352941 0.0954736592795919 -0.10021315043901696 -0.1429999051067759 -0.07631889259552255 -0.09621359681820169 -0.03022590180968537 0.044725027921445226 0.044725027921445226 0.03946309262777076 0.03946309262777076 0.03401340270305381 0.03401340270305381 0.04491722309852347 0.04491722309852347 0.14514142000938163 0.11216316843295798 0.047582616153008794 0.047582616153008794 0.047582616153008794 0.047582616153008794 0.15584553504253135 0.13253731751704917 0.056406603249556875 0.056406603249556875 0.14515182101989493 0.1600255640771459 0.05733717704082236 0.05733717704082236 0.05607166612411246] elementary_charge - Differences by atom index: {0: np.float64(-0.024482422914750473), 1: np.float64(-0.007025493195302435), 2: np.float64(-0.007025493195302435), 3: np.float64(-0.009021954205281624), 4: np.float64(0.023814920661203942), 5: np.float64(0.012139672062151552), 6: np.float64(-0.020759951904065488), 7: np.float64(-0.016800065829999322), 8: np.float64(7.41497013498349e-05), 9: np.float64(0.006718941685907948), 10: np.float64(0.020886240241282106), 11: np.float64(0.005585728261225331), 12: np.float64(0.0055856984589029435), 13: np.float64(-0.029190491440541644), 14: np.float64(0.031447139975779165), 15: np.float64(0.009166433236353508), 16: np.float64(-0.007419741275555974), 17: np.float64(-0.05600793731475581), 18: np.float64(-0.0009708367850850969), 19: np.float64(0.0029687504265237807), 20: np.float64(-0.015029900186903372), 21: np.float64(0.00794083616470586), 22: np.float64(-0.03587365927959193), 23: np.float64(-0.015486849560983076), 24: np.float64(0.027399905106775868), 25: np.float64(-0.0139811074044775), 26: np.float64(0.04291359681820165), 27: np.float64(-0.008474098190314663), 28: np.float64(0.0009749720785547367), 29: np.float64(0.0009749720785547367), 30: np.float64(0.01723690737222921), 31: np.float64(0.01723690737222921), 32: np.float64(0.003686597296946155), 33: np.float64(0.003686597296946155), 34: np.float64(-0.005717223098523509), 35: np.float64(-0.005717223098523509), 36: np.float64(-0.003141420009381668), 37: np.float64(0.009836831567041973), 38: np.float64(0.002117383846991172), 39: np.float64(0.002117383846991172), 40: np.float64(0.002117383846991172), 41: np.float64(0.002117383846991172), 42: np.float64(0.0001544649574686252), 43: np.float64(0.0484626824829508), 44: np.float64(-0.00870660324955691), 45: np.float64(-0.00870660324955691), 46: np.float64(-0.00815182101989495), 47: np.float64(0.0009744359228540667), 48: np.float64(-0.0126371770408224), 49: np.float64(-0.0126371770408224), 50: np.float64(0.02462833387588749)} - Max difference: 0.05600793731475581 e - Mean absolute difference: 0.013057460803529108 e - Largest absolute difference is for atom index 17, which is a C atom + Mean absolute difference: 0.0131 e + Largest absolute difference is 0.0560 e, for atom index 17, which is a C atom @@ -834,8 +876,8 @@ A variety of example notebooks for OpenFF software are provided [here](https://d
✏️ Extra Exercises: Based on the above tutorials, can you:
    -
  • Generate several conformers for one of your MCL-1 ligands and compute their relative energies using OpenFF 2.2.1?
    • -
  • Modify OpenFF 2.2.1 to change some of the parameters applied to one of your MCL-1 ligands? Minimise the ligand with this new force field and see how your changes influence the conformation.
    • +
  • Generate several conformers for one of your MCL-1 ligands and compute their relative energies using OpenFF 2.3.0?
    • +
  • Modify OpenFF 2.3.0 to change some of the parameters applied to one of your MCL-1 ligands? Minimise the ligand with this new force field and see how your changes influence the conformation.
  • Analyse which parameters are shared and which are only applied to one or few molecules for a set of MCL-1 ligands?
diff --git a/notebooks/protein_ligand_complex_parameterisation_and_md.ipynb b/notebooks/protein_ligand_complex_parameterisation_and_md.ipynb index a2ff4d2..0491f8a 100644 --- a/notebooks/protein_ligand_complex_parameterisation_and_md.ipynb +++ b/notebooks/protein_ligand_complex_parameterisation_and_md.ipynb @@ -298,7 +298,7 @@ { "cell_type": "code", "execution_count": null, - "id": "82b9912f", + "id": "cddfa873", "metadata": { "tags": [ "placeholder" @@ -721,7 +721,7 @@ { "cell_type": "code", "execution_count": null, - "id": "ed983a3e", + "id": "0b516cda", "metadata": { "tags": [ "placeholder" @@ -949,7 +949,7 @@ { "cell_type": "code", "execution_count": null, - "id": "d99ede30", + "id": "7b22239d", "metadata": { "tags": [ "placeholder" diff --git a/notebooks/small_molecule_parameterisation.ipynb b/notebooks/small_molecule_parameterisation.ipynb index e9e4f45..7686b8e 100644 --- a/notebooks/small_molecule_parameterisation.ipynb +++ b/notebooks/small_molecule_parameterisation.ipynb @@ -513,7 +513,7 @@ { "cell_type": "code", "execution_count": null, - "id": "f2e4d2fd", + "id": "f9639715", "metadata": { "tags": [ "placeholder" @@ -831,7 +831,7 @@ { "cell_type": "code", "execution_count": null, - "id": "82f8bb89", + "id": "5807541e", "metadata": { "tags": [ "placeholder" @@ -855,7 +855,7 @@ { "cell_type": "code", "execution_count": null, - "id": "38f8788d", + "id": "0c76d7cb", "metadata": { "tags": [ "placeholder" @@ -1171,7 +1171,7 @@ { "cell_type": "code", "execution_count": null, - "id": "25e7a681", + "id": "17b0325a", "metadata": { "tags": [ "placeholder" From 5ccacec283629398d0b3c332507c460d34d9c88b Mon Sep 17 00:00:00 2001 From: James Gebbie-Rayet Date: Thu, 17 Sep 2026 20:46:02 +0100 Subject: [PATCH 7/9] Delete .github/dependabot.yml --- .github/dependabot.yml | 16 ---------------- 1 file changed, 16 deletions(-) delete mode 100644 .github/dependabot.yml diff --git a/.github/dependabot.yml b/.github/dependabot.yml deleted file mode 100644 index ef2e5af..0000000 --- a/.github/dependabot.yml +++ /dev/null @@ -1,16 +0,0 @@ -# To get started with Dependabot version updates, you'll need to specify which -# package ecosystems to update and where the package manifests are located. -# Please see the documentation for all configuration options: -# https://docs.github.com/code-security/dependabot/dependabot-version-updates/configuration-options-for-the-dependabot.yml-file - -version: 2 -updates: - - package-ecosystem: "github-actions" - directory: "/" - schedule: - # Check for updates to GitHub Actions every day - interval: "daily" - time: "09:00" - timezone: "UTC" - assignees: - - "jimboid" From c495e9b6097232f28ee68e3201dffa833e1edebf Mon Sep 17 00:00:00 2001 From: James Gebbie-Rayet Date: Thu, 17 Sep 2026 20:46:28 +0100 Subject: [PATCH 8/9] Add Renovate configuration for dependency management --- .github/renovate.json | 48 +++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 48 insertions(+) create mode 100644 .github/renovate.json diff --git a/.github/renovate.json b/.github/renovate.json new file mode 100644 index 0000000..0b3e5d8 --- /dev/null +++ b/.github/renovate.json @@ -0,0 +1,48 @@ +{ + "$schema": "https://docs.renovatebot.com/renovate-schema.json", + "extends": [ + "config:best-practices", + "helpers:pinGitHubActionDigestsToSemver" + ], + "timezone": "Europe/London", + "schedule": ["after 5am and before 5pm every weekday"], + "dependencyDashboard": true, + "assignees": ["jimboid"], + "minimumReleaseAge": "3 days", + "packageRules": [ + { + "matchManagers": ["github-actions"], + "groupName": "GitHub Actions", + "addLabels": ["ci", "dependencies"] + } + ], + "customManagers": [ + { + "customType": "regex", + "managerFilePatterns": [ + "/(^|/)Dockerfile[^/]*$/" + ], + "matchStringsStrategy": "recursive", + "matchStrings": [ + "pip3? install[\\s\\S]*?(?(?:[\\w.-]+==[\\w.-]+[\\s\\\\]*)+)", + "(?[\\w.-]+)==(?[\\w.-]+)" + ], + "datasourceTemplate": "pypi", + "versioningTemplate": "pep440" + }, + { + "customType": "regex", + "managerFilePatterns": [ + "/(^|/)Dockerfile[^/]*$/" + ], + "matchStringsStrategy": "recursive", + "matchStrings": [ + "(?:mamba|conda) install[\\s\\S]*?(?(?:(?:[\\w-]+::)?[\\w.-]+=[\\w.-]+[\\s\\\\]*)+)", + "(?:(?[\\w-]+)::)?(?[\\w.-]+)=(?[\\w.-]+)" + ], + "datasourceTemplate": "conda", + "registryUrlTemplate": "https://api.anaconda.org/package/{{#if channel}}{{channel}}{{else}}conda-forge{{/if}}/", + "versioningTemplate": "pep440" + } + ] +} From 241a49b4c14c50a82ee205a15a7625256f14416d Mon Sep 17 00:00:00 2001 From: James Gebbie-Rayet Date: Thu, 17 Sep 2026 20:56:18 +0100 Subject: [PATCH 9/9] Update dependencies in Dockerfile and comment out script --- docker/Dockerfile | 16 ++++++++++++---- 1 file changed, 12 insertions(+), 4 deletions(-) diff --git a/docker/Dockerfile b/docker/Dockerfile index a567703..7e2bb6e 100644 --- a/docker/Dockerfile +++ b/docker/Dockerfile @@ -9,10 +9,18 @@ USER $NB_USER WORKDIR $HOME # Install workshop deps -RUN conda install -c conda-forge "openff-toolkit-examples>=0.17.0" "openff-interchange>=0.4.7" MDAnalysis prolif py3dmol openbabel -y - -COPY docker/fix-nglview.sh /tmp/fix-nglview.sh -RUN bash /tmp/fix-nglview.sh +RUN conda install -y -c conda-forge \ + conda-forge::openff-toolkit-examples=0.19.0 \ + conda-forge::openff-interchange=0.5.4 \ + conda-forge::packmol=21.2.1 \ + conda-forge::ambertools=26.0 \ + conda-forge::MDAnalysis=2.10.0 \ + conda-forge::prolif=2.2.1 \ + conda-forge::py3dmol=2.5.4 \ + conda-forge::openbabel=3.2.1 + +#COPY docker/fix-nglview.sh /tmp/fix-nglview.sh +#RUN bash /tmp/fix-nglview.sh # Get workshop files and move them to jovyan directory. COPY --chown=1000:100 . .