diff --git a/cookbook/create_alchemical_network.ipynb b/cookbook/create_alchemical_network.ipynb index edebab0..959bc0c 100644 --- a/cookbook/create_alchemical_network.ipynb +++ b/cookbook/create_alchemical_network.ipynb @@ -13,7 +13,12 @@ "id": "2b11a786-0fab-4990-a2aa-242067f89287", "metadata": {}, "source": [ - "The final setup step is to compile all the different bits of information into a description of a single simulation campaign. This description takes the form of an [Alchemical Network]. Similarly to `LigandNetwork`, the `AlchemicalNetwork` class is a graph of all the transformations in the campaign; however, in an `AlchemicalNetwork`, these transformations include all the information needed to perform the transformation. By contrast, a `LigandNetwork` includes only the ligands themselves.\n", + "The final setup step is to compile all the different bits of information into a description of a single simulation campaign. This description takes the form of an [Alchemical Network].\n", + "Similar to `LigandNetwork`, the `AlchemicalNetwork` class is a graph of all the transformations in the campaign.\n", + "\n", + "However, in an `AlchemicalNetwork`, these transformations include all the information needed to perform the transformation. By contrast, a `LigandNetwork` includes only the ligands themselves.\n", + "\n", + "See the [User Guide: Alchemical Networks](https://docs.openfree.energy/en/stable/guide/setup/alchemical_network_model.html) for more background information.\n", "\n", "[Alchemical Network]: https://docs.openfree.energy/en/stable/reference/api/generated/openfe.AlchemicalNetwork.html" ] @@ -28,19 +33,7 @@ "tags": [] }, "source": [ - "## Setup" - ] - }, - { - "cell_type": "code", - "execution_count": 1, - "id": "ba6d0d2f-03df-463f-8976-c15a2bdf65e5", - "metadata": {}, - "outputs": [], - "source": [ - "import openfe, rdkit.Chem\n", - "from openff.units import unit\n", - "from openfe.protocols.openmm_rfe import RelativeHybridTopologyProtocol" + "## Start with a `LigandNetwork`" ] }, { @@ -48,18 +41,25 @@ "id": "6bfee4c7-0aaf-49a6-b35a-bdbe99444bab", "metadata": {}, "source": [ - "This cookbook assumes you've already loaded a `LigandNetwork`. For more information, see [Generate a Ligand Network Automatically]:\n", + "This cookbook assumes you've already loaded a `LigandNetwork`.\n", "\n", - "[Generate a Ligand Network Automatically]: https://docs.openfree.energy/en/stable/cookbook/generate_ligand_network.html" + "To learn how to build a `LigandNetwork`, see [Planning a Ligand Network]:\n", + "\n", + "[Planning a Ligand Network]: https://docs.openfree.energy/en/stable/cookbook/generate_ligand_network.html" ] }, { "cell_type": "code", - "execution_count": 2, + "execution_count": 1, "id": "4512c4d6-720f-4008-af53-a0b56b2ee1e9", "metadata": {}, "outputs": [], "source": [ + "import openfe\n", + "import rdkit.Chem\n", + "from openff.units import unit\n", + "from openfe.protocols.openmm_rfe import RelativeHybridTopologyProtocol\n", + "\n", "ligand_network = openfe.ligand_network_planning.generate_minimal_spanning_network(\n", " ligands=[\n", " openfe.SmallMoleculeComponent(mol) \n", @@ -85,7 +85,7 @@ }, { "cell_type": "code", - "execution_count": 3, + "execution_count": 2, "id": "940e2e19-4e7d-4853-ab5e-54626b3384d8", "metadata": {}, "outputs": [], @@ -105,7 +105,7 @@ }, { "cell_type": "code", - "execution_count": 4, + "execution_count": 3, "id": "d3054cfc-5040-4882-b3a3-92aba7ab5d1d", "metadata": { "slideshow": { @@ -149,7 +149,7 @@ }, { "cell_type": "code", - "execution_count": 5, + "execution_count": 4, "id": "733dee64-d23a-4eb0-87cf-c03065b8d2a6", "metadata": {}, "outputs": [], @@ -181,7 +181,7 @@ }, { "cell_type": "code", - "execution_count": 6, + "execution_count": 5, "id": "c49f36b1-c3c1-425e-b1de-24da25ed01c3", "metadata": {}, "outputs": [], @@ -189,19 +189,26 @@ "# In an RBFE, each edge includes two \"legs\": \n", "# one for the ligand complexed to the protein, \n", "# and the other for the ligand free in solution\n", + "\n", + "\n", "legs = {\n", - " \"solvent\": {\n", - " # Specify the components common to all systems in this leg\n", - " 'solvent': solvent,\n", - " },\n", - " \"complex\": {\n", - " # Specify the components common to all systems in this leg\n", - " 'solvent': solvent,\n", - " 'protein': protein,\n", - " }\n", - "}\n", - "transformations = []\n", + " # Specify the components common to all systems in this leg\n", + " \"solvent\": {'solvent': solvent},\n", + " # Specify the components common to all systems in this leg\n", + " \"complex\": {'solvent': solvent, 'protein': protein}\n", + "}" + ] + }, + { + "cell_type": "code", + "execution_count": 6, + "id": "77d71229-5dcc-40cf-a9ae-0724a7327578", + "metadata": {}, + "outputs": [], + "source": [ + "# create a list of transformations\n", "\n", + "transformations = []\n", "for mapping in ligand_network.edges:\n", " for leg, common_components in legs.items():\n", " system_a = openfe.ChemicalSystem(\n", @@ -229,7 +236,7 @@ " name=f\"easy_rbfe_{system_a.name}_{system_b.name}\"\n", " )\n", " \n", - " transformations.append(transformation)" + " transformations.append(transformation)\n" ] }, { @@ -324,8 +331,16 @@ "\n", "for n, transformation in enumerate(alchemical_network.edges):\n", " transformation_name = transformation.name or transformation.key\n", - " transformation.dump(transformations_dir / f\"{n}_{transformation_name}.json\")" + " transformation.to_json(transformations_dir / f\"{n}_{transformation_name}.json\")" ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "28d7af7c-0e8e-4379-91c2-c9d291c3923b", + "metadata": {}, + "outputs": [], + "source": [] } ], "metadata": { @@ -344,7 +359,7 @@ "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", - "version": "3.11.6" + "version": "3.13.11" }, "widgets": { "application/vnd.jupyter.widget-state+json": { diff --git a/networks/Preparing AlchemicalNetworks.ipynb b/networks/Preparing AlchemicalNetworks.ipynb deleted file mode 100644 index b415cb8..0000000 --- a/networks/Preparing AlchemicalNetworks.ipynb +++ /dev/null @@ -1,610 +0,0 @@ -{ - "cells": [ - { - "cell_type": "markdown", - "id": "4ed54201-aa1f-4ba9-bdb7-1c7f5a2cbd05", - "metadata": {}, - "source": [ - "# Preparing `AlchemicalNetwork`s" - ] - }, - { - "cell_type": "markdown", - "id": "6da80e66-6e17-407f-9ba4-a25ded4777bc", - "metadata": {}, - "source": [ - "This notebook will illustrate how to build `AlchemicalNetwork` objects,\n", - "from a starting point of chemical models stored in sdf and pdb files.\n", - "\n", - "An `AlchemicalNetwork` is used to represent an entire network of calculations,\n", - "and is composed of many smaller objects:\n", - "\n", - "- An `AlchemicalNetwork` composed of \n", - " - each node a `ChemicalSystem`\n", - " - each containing many components, such as `SmallMoleculeComponent`, `ProteinComponent`\n", - " - internally each Component usually wraps an RDKit representation\n", - " - each directed edge a `Transformation`, containing\n", - " - two `ChemicalSystem`s, the 'A' and 'B' side\n", - " - zero or more `Mapping` objects relating these two sides\n", - " - a `Protocol` defining the computational method to be applied to other items" - ] - }, - { - "cell_type": "code", - "execution_count": 1, - "id": "27b4aa61-ca82-44fa-be84-c802f1083a29", - "metadata": {}, - "outputs": [], - "source": [ - "# suppress `numba` warnings, if present\n", - "from numba.core.errors import NumbaWarning\n", - "import warnings\n", - "\n", - "warnings.simplefilter('ignore', category=NumbaWarning)" - ] - }, - { - "cell_type": "code", - "execution_count": 2, - "id": "6bb3f4a4-2135-494a-8365-9ef229dd124d", - "metadata": {}, - "outputs": [], - "source": [ - "import openfe\n", - "from gufe import AlchemicalNetwork\n", - "from openff.units import unit\n", - "from rdkit import Chem" - ] - }, - { - "cell_type": "markdown", - "id": "5d3285c5-0b7e-461f-a761-dabd01182775", - "metadata": {}, - "source": [ - "## Define `ChemicalSystem`s for network nodes\n", - "\n", - "We'll start by defining the nodes for our network.\n", - "A `ChemicalSystem` is made of one or more `Component`s. These can be one of `ProteinComponent`, `SmallMoleculeComponent`, or `SolventComponent`, and potentially others as needed. This design allows for memory efficient representation of large networks with perhaps hundreds or thousands of nodes, but perhaps far fewer variants in proteins, ligands, etc." - ] - }, - { - "cell_type": "markdown", - "id": "69822065", - "metadata": {}, - "source": [ - "### Reading Ligands\n", - "\n", - "The ligands are concatenated in a single sdf file, we'll read these using RDKit.\n", - "\n", - "Each of the ligands have been pre-docked into the protein and aligned to their common scaffold. It is important to recognize that any processing required to prepare ligand and protein structures for alchemical free energy calculations should be done *before* the steps we are taking here." - ] - }, - { - "cell_type": "code", - "execution_count": 3, - "id": "16527065", - "metadata": {}, - "outputs": [ - { - "data": { - "text/plain": [ - "[SmallMoleculeComponent(name=lig_ejm_31),\n", - " SmallMoleculeComponent(name=lig_ejm_42),\n", - " SmallMoleculeComponent(name=lig_ejm_43),\n", - " SmallMoleculeComponent(name=lig_ejm_45),\n", - " SmallMoleculeComponent(name=lig_ejm_46),\n", - " SmallMoleculeComponent(name=lig_ejm_47),\n", - " SmallMoleculeComponent(name=lig_ejm_48),\n", - " SmallMoleculeComponent(name=lig_ejm_50),\n", - " SmallMoleculeComponent(name=lig_ejm_54),\n", - " SmallMoleculeComponent(name=lig_ejm_55),\n", - " SmallMoleculeComponent(name=lig_jmc_23),\n", - " SmallMoleculeComponent(name=lig_jmc_27),\n", - " SmallMoleculeComponent(name=lig_jmc_28)]" - ] - }, - "execution_count": 3, - "metadata": {}, - "output_type": "execute_result" - } - ], - "source": [ - "ligands = [\n", - " openfe.SmallMoleculeComponent(m) for m in Chem.SDMolSupplier('data/tyk2_ligands.sdf', removeHs=False)\n", - "]\n", - "ligands" - ] - }, - { - "cell_type": "markdown", - "id": "0bb2cb69", - "metadata": {}, - "source": [ - "### Reading the protein\n", - "\n", - "The protein is supplied as a PDB file, readable via the `ProteinComponent.from_pdb_file` class method." - ] - }, - { - "cell_type": "code", - "execution_count": 4, - "id": "a31c89a4", - "metadata": {}, - "outputs": [ - { - "data": { - "text/plain": [ - "ProteinComponent(name=tyk2)" - ] - }, - "execution_count": 4, - "metadata": {}, - "output_type": "execute_result" - } - ], - "source": [ - "protein = openfe.ProteinComponent.from_pdb_file('./data/tyk2_protein.pdb', name='tyk2')\n", - "\n", - "protein" - ] - }, - { - "cell_type": "markdown", - "id": "e653c208-1f69-4f08-a38e-55b7123292d9", - "metadata": {}, - "source": [ - "### Defining the solvent\n", - "\n", - "We'll also need at least one `SolventComponent` to encode our choice of solvent and counterions, with concentration.\n", - "The concentration is defined as having units supplied by `openff.units`, this package is used to avoid confusion.\n", - "\n", - "\n", - "The `SolventComponent` doesn't actually perform any actual solvation (packing water molecules, ions); that is performed just before simulation time during `Protocol` execution." - ] - }, - { - "cell_type": "code", - "execution_count": 5, - "id": "3c47166d-143b-41f0-9dce-7805e56d7191", - "metadata": {}, - "outputs": [ - { - "data": { - "text/plain": [ - "SolventComponent(name=O, Na+, Cl-)" - ] - }, - "execution_count": 5, - "metadata": {}, - "output_type": "execute_result" - } - ], - "source": [ - "solvent = openfe.SolventComponent(positive_ion='Na', \n", - " negative_ion='Cl',\n", - " neutralize=True, \n", - " ion_concentration=0.15*unit.molar)\n", - "solvent" - ] - }, - { - "cell_type": "markdown", - "id": "b0b11357-2179-496b-9bff-89303e1c9c33", - "metadata": {}, - "source": [ - "### Build the `ChemicalSystem`s\n", - "\n", - "We can now construct the `ChemicalSystem`s we want represented in our network. Since we are planning to perform relative binding free energy (RBFE) calculations, we'll define both *complex* and *solvent* variants for each ligand.\n", - "\n", - "This produces a dictionary mapping the ligand name to the `ChemicalSystem` that contains that ligand.\n", - "There are two dictionaries, for complexed and solvated ligands respectively." - ] - }, - { - "cell_type": "code", - "execution_count": 6, - "id": "10edc8b3-5f45-4e37-b1d5-bd508601e352", - "metadata": {}, - "outputs": [ - { - "data": { - "text/plain": [ - "{'lig_ejm_31': ChemicalSystem(name=lig_ejm_31_complex, components={'ligand': SmallMoleculeComponent(name=lig_ejm_31), 'solvent': SolventComponent(name=O, Na+, Cl-), 'protein': ProteinComponent(name=tyk2)}),\n", - " 'lig_ejm_42': ChemicalSystem(name=lig_ejm_42_complex, components={'ligand': SmallMoleculeComponent(name=lig_ejm_42), 'solvent': SolventComponent(name=O, Na+, Cl-), 'protein': ProteinComponent(name=tyk2)}),\n", - " 'lig_ejm_43': ChemicalSystem(name=lig_ejm_43_complex, components={'ligand': SmallMoleculeComponent(name=lig_ejm_43), 'solvent': SolventComponent(name=O, Na+, Cl-), 'protein': ProteinComponent(name=tyk2)}),\n", - " 'lig_ejm_45': ChemicalSystem(name=lig_ejm_45_complex, components={'ligand': SmallMoleculeComponent(name=lig_ejm_45), 'solvent': SolventComponent(name=O, Na+, Cl-), 'protein': ProteinComponent(name=tyk2)}),\n", - " 'lig_ejm_46': ChemicalSystem(name=lig_ejm_46_complex, components={'ligand': SmallMoleculeComponent(name=lig_ejm_46), 'solvent': SolventComponent(name=O, Na+, Cl-), 'protein': ProteinComponent(name=tyk2)}),\n", - " 'lig_ejm_47': ChemicalSystem(name=lig_ejm_47_complex, components={'ligand': SmallMoleculeComponent(name=lig_ejm_47), 'solvent': SolventComponent(name=O, Na+, Cl-), 'protein': ProteinComponent(name=tyk2)}),\n", - " 'lig_ejm_48': ChemicalSystem(name=lig_ejm_48_complex, components={'ligand': SmallMoleculeComponent(name=lig_ejm_48), 'solvent': SolventComponent(name=O, Na+, Cl-), 'protein': ProteinComponent(name=tyk2)}),\n", - " 'lig_ejm_50': ChemicalSystem(name=lig_ejm_50_complex, components={'ligand': SmallMoleculeComponent(name=lig_ejm_50), 'solvent': SolventComponent(name=O, Na+, Cl-), 'protein': ProteinComponent(name=tyk2)}),\n", - " 'lig_ejm_54': ChemicalSystem(name=lig_ejm_54_complex, components={'ligand': SmallMoleculeComponent(name=lig_ejm_54), 'solvent': SolventComponent(name=O, Na+, Cl-), 'protein': ProteinComponent(name=tyk2)}),\n", - " 'lig_ejm_55': ChemicalSystem(name=lig_ejm_55_complex, components={'ligand': SmallMoleculeComponent(name=lig_ejm_55), 'solvent': SolventComponent(name=O, Na+, Cl-), 'protein': ProteinComponent(name=tyk2)}),\n", - " 'lig_jmc_23': ChemicalSystem(name=lig_jmc_23_complex, components={'ligand': SmallMoleculeComponent(name=lig_jmc_23), 'solvent': SolventComponent(name=O, Na+, Cl-), 'protein': ProteinComponent(name=tyk2)}),\n", - " 'lig_jmc_27': ChemicalSystem(name=lig_jmc_27_complex, components={'ligand': SmallMoleculeComponent(name=lig_jmc_27), 'solvent': SolventComponent(name=O, Na+, Cl-), 'protein': ProteinComponent(name=tyk2)}),\n", - " 'lig_jmc_28': ChemicalSystem(name=lig_jmc_28_complex, components={'ligand': SmallMoleculeComponent(name=lig_jmc_28), 'solvent': SolventComponent(name=O, Na+, Cl-), 'protein': ProteinComponent(name=tyk2)})}" - ] - }, - "execution_count": 6, - "metadata": {}, - "output_type": "execute_result" - } - ], - "source": [ - "complexed = {l.name: openfe.ChemicalSystem(components={'ligand': l,\n", - " 'solvent': solvent, \n", - " 'protein': protein}, \n", - " name=f\"{l.name}_complex\") \n", - " for l in ligands}\n", - "complexed" - ] - }, - { - "cell_type": "code", - "execution_count": 7, - "id": "c8a0b60d-bfae-4f99-9f4d-656c982433b9", - "metadata": { - "scrolled": false - }, - "outputs": [ - { - "data": { - "text/plain": [ - "{'lig_ejm_31': ChemicalSystem(name=lig_ejm_31_water, components={'ligand': SmallMoleculeComponent(name=lig_ejm_31), 'solvent': SolventComponent(name=O, Na+, Cl-)}),\n", - " 'lig_ejm_42': ChemicalSystem(name=lig_ejm_42_water, components={'ligand': SmallMoleculeComponent(name=lig_ejm_42), 'solvent': SolventComponent(name=O, Na+, Cl-)}),\n", - " 'lig_ejm_43': ChemicalSystem(name=lig_ejm_43_water, components={'ligand': SmallMoleculeComponent(name=lig_ejm_43), 'solvent': SolventComponent(name=O, Na+, Cl-)}),\n", - " 'lig_ejm_45': ChemicalSystem(name=lig_ejm_45_water, components={'ligand': SmallMoleculeComponent(name=lig_ejm_45), 'solvent': SolventComponent(name=O, Na+, Cl-)}),\n", - " 'lig_ejm_46': ChemicalSystem(name=lig_ejm_46_water, components={'ligand': SmallMoleculeComponent(name=lig_ejm_46), 'solvent': SolventComponent(name=O, Na+, Cl-)}),\n", - " 'lig_ejm_47': ChemicalSystem(name=lig_ejm_47_water, components={'ligand': SmallMoleculeComponent(name=lig_ejm_47), 'solvent': SolventComponent(name=O, Na+, Cl-)}),\n", - " 'lig_ejm_48': ChemicalSystem(name=lig_ejm_48_water, components={'ligand': SmallMoleculeComponent(name=lig_ejm_48), 'solvent': SolventComponent(name=O, Na+, Cl-)}),\n", - " 'lig_ejm_50': ChemicalSystem(name=lig_ejm_50_water, components={'ligand': SmallMoleculeComponent(name=lig_ejm_50), 'solvent': SolventComponent(name=O, Na+, Cl-)}),\n", - " 'lig_ejm_54': ChemicalSystem(name=lig_ejm_54_water, components={'ligand': SmallMoleculeComponent(name=lig_ejm_54), 'solvent': SolventComponent(name=O, Na+, Cl-)}),\n", - " 'lig_ejm_55': ChemicalSystem(name=lig_ejm_55_water, components={'ligand': SmallMoleculeComponent(name=lig_ejm_55), 'solvent': SolventComponent(name=O, Na+, Cl-)}),\n", - " 'lig_jmc_23': ChemicalSystem(name=lig_jmc_23_water, components={'ligand': SmallMoleculeComponent(name=lig_jmc_23), 'solvent': SolventComponent(name=O, Na+, Cl-)}),\n", - " 'lig_jmc_27': ChemicalSystem(name=lig_jmc_27_water, components={'ligand': SmallMoleculeComponent(name=lig_jmc_27), 'solvent': SolventComponent(name=O, Na+, Cl-)}),\n", - " 'lig_jmc_28': ChemicalSystem(name=lig_jmc_28_water, components={'ligand': SmallMoleculeComponent(name=lig_jmc_28), 'solvent': SolventComponent(name=O, Na+, Cl-)})}" - ] - }, - "execution_count": 7, - "metadata": {}, - "output_type": "execute_result" - } - ], - "source": [ - "solvated = {l.name: openfe.ChemicalSystem(components={'ligand': l, \n", - " 'solvent': solvent}, \n", - " name=f\"{l.name}_water\") \n", - " for l in ligands}\n", - "solvated" - ] - }, - { - "cell_type": "markdown", - "id": "7bc3a65a-979c-4114-b948-f7315796709d", - "metadata": {}, - "source": [ - "## Define `Transformation`s between `ChemicalSystem`s as network edges" - ] - }, - { - "cell_type": "markdown", - "id": "0a92b4c1-fffb-410c-9c76-177fe3ea0a87", - "metadata": {}, - "source": [ - "A `Transformation` is a directed edge between two `ChemicalSystem`s. It includes a `Protocol` parameterized with `Settings`, and optionally a `ComponentMapping`. \n", - "\n", - "The `Protocol` defines the actual computational method used to evaluate the `Transformation` to yield estimates for the free energy difference between the `ChemicalSystem`s.\n", - "\n", - "The `ComponentMapping` defines the atom mapping(s) between corresponding `Component`s in the two `ChemicalSystem`s. This is often critical for relative binding free energy calculations, since the choice of mapping can heavily influence convergence of the resulting estimates." - ] - }, - { - "cell_type": "markdown", - "id": "ace9de53-2527-47f5-8aae-98d9adef2c1d", - "metadata": {}, - "source": [ - "### Define the `Protocol` used for `Transformation` evaluation" - ] - }, - { - "cell_type": "markdown", - "id": "0338be2a-b861-4534-bc78-49c2e00ca2ac", - "metadata": {}, - "source": [ - "For this example, we'll use the same `Protocol` for all our `Transformation`s, with identical `Settings` for each." - ] - }, - { - "cell_type": "code", - "execution_count": 8, - "id": "183b23a6-1839-4586-bbbe-60b2defd2f01", - "metadata": {}, - "outputs": [], - "source": [ - "from openfe.protocols import openmm_rfe\n" - ] - }, - { - "cell_type": "markdown", - "id": "9e686c7c-9ef7-4e97-b384-63334eeb894e", - "metadata": {}, - "source": [ - "Any given `Protocol` has a `default_settings()` method, which can be used to get the default settings that are specific to that `Protocol`:" - ] - }, - { - "cell_type": "code", - "execution_count": 9, - "id": "d7bac455-ddc1-4735-96d7-5c6681ab3cbc", - "metadata": { - "scrolled": true, - "tags": [] - }, - "outputs": [ - { - "data": { - "text/plain": [ - "{'forcefield_settings': {'constraints': 'hbonds',\n", - " 'rigid_water': True,\n", - " 'remove_com': False,\n", - " 'hydrogen_mass': 4.0,\n", - " 'forcefields': ['amber/ff14SB.xml',\n", - " 'amber/tip3p_standard.xml',\n", - " 'amber/tip3p_HFE_multivalent.xml',\n", - " 'amber/phosaa10.xml'],\n", - " 'small_molecule_forcefield': 'openff-2.0.0'},\n", - " 'thermo_settings': {'temperature': 298.15 ,\n", - " 'pressure': 0.9869232667160129 ,\n", - " 'ph': None,\n", - " 'redox_potential': None},\n", - " 'system_settings': {'nonbonded_method': 'PME',\n", - " 'nonbonded_cutoff': 1.0 },\n", - " 'solvation_settings': {'solvent_model': 'tip3p',\n", - " 'solvent_padding': 1.2 },\n", - " 'alchemical_settings': {'lambda_functions': 'default',\n", - " 'lambda_windows': 11,\n", - " 'unsampled_endstates': False,\n", - " 'use_dispersion_correction': False,\n", - " 'softcore_LJ_v2': True,\n", - " 'softcore_electrostatics': True,\n", - " 'softcore_alpha': 0.85,\n", - " 'softcore_electrostatics_alpha': 0.3,\n", - " 'softcore_sigma_Q': 1.0,\n", - " 'interpolate_old_and_new_14s': False,\n", - " 'flatten_torsions': False},\n", - " 'alchemical_sampler_settings': {'online_analysis_interval': None,\n", - " 'n_repeats': 3,\n", - " 'sampler_method': 'repex',\n", - " 'online_analysis_target_error': 0.1 ,\n", - " 'online_analysis_minimum_iterations': 500,\n", - " 'flatness_criteria': 'logZ-flatness',\n", - " 'gamma0': 1.0,\n", - " 'n_replicas': 11},\n", - " 'engine_settings': {'compute_platform': None},\n", - " 'integrator_settings': {'timestep': 4 ,\n", - " 'collision_rate': 1.0 ,\n", - " 'n_steps': 250 ,\n", - " 'reassign_velocities': False,\n", - " 'splitting': 'V R O R V',\n", - " 'n_restart_attempts': 20,\n", - " 'constraint_tolerance': 1e-06,\n", - " 'barostat_frequency': 25 },\n", - " 'simulation_settings': {'equilibration_length': 2.0 ,\n", - " 'production_length': 5.0 ,\n", - " 'forcefield_cache': None,\n", - " 'minimization_steps': 5000,\n", - " 'output_filename': 'simulation.nc',\n", - " 'output_indices': 'all',\n", - " 'checkpoint_interval': 250 ,\n", - " 'checkpoint_storage': 'checkpoint.nc'}}" - ] - }, - "execution_count": 9, - "metadata": {}, - "output_type": "execute_result" - } - ], - "source": [ - "protocol_settings = openmm_rfe.RelativeHybridTopologyProtocol.default_settings()\n", - "protocol_settings.dict()" - ] - }, - { - "cell_type": "markdown", - "id": "d7898feb-419a-45cc-8837-53cd583d130b", - "metadata": {}, - "source": [ - "These can be edited, e.g. with:" - ] - }, - { - "cell_type": "code", - "execution_count": 10, - "id": "998ef869-7086-4c24-9c53-2d4e9553b6b9", - "metadata": { - "tags": [] - }, - "outputs": [], - "source": [ - "protocol_settings.thermo_settings.temperature = 299 * unit.kelvin" - ] - }, - { - "cell_type": "markdown", - "id": "70e558b5-93e4-4c91-8ff0-bc5b99646a71", - "metadata": {}, - "source": [ - "We can now produce a parameterized `RelativeHybridTopologyProtocol` instance:" - ] - }, - { - "cell_type": "code", - "execution_count": 11, - "id": "bd1d0a0d-00c9-44a1-820c-186b4ee05334", - "metadata": {}, - "outputs": [], - "source": [ - "protocol = openmm_rfe.RelativeHybridTopologyProtocol(protocol_settings)" - ] - }, - { - "cell_type": "markdown", - "id": "57239762-43b9-4501-b3ad-9eebbcd64a20", - "metadata": {}, - "source": [ - "### Build the `Transformation`s" - ] - }, - { - "cell_type": "markdown", - "id": "8ccdec56-519e-4293-86c4-60b0255cfcad", - "metadata": {}, - "source": [ - "We can now construct the `Transformation`s we want represented in our network.\n", - "\n", - "We'll use the predefined connections from the `tyk2` system from above as the basis for our choices here, but you could use any network planner of your choice to generate connections and use those instead." - ] - }, - { - "cell_type": "code", - "execution_count": 12, - "id": "4a73c2d6", - "metadata": {}, - "outputs": [], - "source": [ - "from openfe.setup import LomapAtomMapper\n", - "\n", - "mapper = LomapAtomMapper(element_change=False)\n", - "\n", - "connections = [(\"lig_ejm_31\", \"lig_ejm_50\"),\n", - " (\"lig_ejm_46\", \"lig_jmc_23\"),\n", - " (\"lig_ejm_31\", \"lig_ejm_55\"),\n", - " (\"lig_ejm_31\", \"lig_ejm_48\"),\n", - " (\"lig_ejm_31\", \"lig_ejm_54\"),\n", - " (\"lig_ejm_31\", \"lig_ejm_47\"),\n", - " (\"lig_ejm_31\", \"lig_ejm_46\"),\n", - " (\"lig_ejm_46\", \"lig_jmc_27\"),\n", - " (\"lig_ejm_46\", \"lig_jmc_28\"),\n", - " (\"lig_ejm_42\", \"lig_ejm_43\"),\n", - " (\"lig_ejm_31\", \"lig_ejm_42\"),\n", - " (\"lig_ejm_45\", \"lig_ejm_55\"),]" - ] - }, - { - "cell_type": "markdown", - "id": "cfd3fe16-1666-4b78-8263-becff81dfef1", - "metadata": {}, - "source": [ - "Since we are planning to perform relative binding free energy (RBFE) calculations, we'll define both *complex* and *solvent* variants for each `Transformation`:" - ] - }, - { - "cell_type": "code", - "execution_count": 13, - "id": "f3379b4a", - "metadata": {}, - "outputs": [], - "source": [ - "complexed_transformations = []\n", - "solvated_transformations = []\n", - "\n", - "for (ligA_name, ligB_name) in connections:\n", - " ligA = complexed[ligA_name]['ligand']\n", - " ligB = complexed[ligB_name]['ligand']\n", - " \n", - " mapping = next(mapper.suggest_mappings(ligA, ligB))\n", - " \n", - " complexed_transformations.append(\n", - " openfe.Transformation(stateA=complexed[ligA_name], \n", - " stateB=complexed[ligB_name], \n", - " mapping={'ligand': mapping},\n", - " protocol=protocol) \n", - " )\n", - " solvated_transformations.append(\n", - " openfe.Transformation(stateA=solvated[ligA_name], \n", - " stateB=solvated[ligB_name], \n", - " mapping={'ligand': mapping},\n", - " protocol=protocol) \n", - " )" - ] - }, - { - "cell_type": "markdown", - "id": "36c79a08-d4aa-4388-a0af-be51906f7672", - "metadata": {}, - "source": [ - "## Create the `AlchemicalNetwork`" - ] - }, - { - "cell_type": "markdown", - "id": "6ebc1174-b348-4e44-97f7-444c44c8d6d0", - "metadata": {}, - "source": [ - "An `AlchemicalNetwork` is simply the combination of `ChemicalSystem`s (nodes) and `Transformation`s (directed edges) that we want to evaluate $\\Delta G$s for. This data structure functions as a declaration of what you want to compute.\n", - "\n", - "We'll finish here by creating an `AlchemicalNetwork` from the collection of objects we've built so far." - ] - }, - { - "cell_type": "code", - "execution_count": 14, - "id": "b45baada-efa4-46f7-a2cb-6e24a9309f1e", - "metadata": {}, - "outputs": [ - { - "data": { - "text/plain": [ - "" - ] - }, - "execution_count": 14, - "metadata": {}, - "output_type": "execute_result" - } - ], - "source": [ - "network = AlchemicalNetwork(edges=(solvated_transformations + complexed_transformations), \n", - " nodes=(list(solvated.values()) + list(complexed.values())),\n", - " name=\"tyk2_relative_benchmark\")\n", - "network" - ] - }, - { - "cell_type": "markdown", - "id": "262c5908-640d-45dc-9ca4-6d912428b8ac", - "metadata": {}, - "source": [ - "That's it! We simply toss in all `Transformation`s (edges) and `ChemicalSystem`s (nodes) we want included in this `AlchemicalNetwork`, and optionally give it a name that means something to us (it need not be unique, but can be used to query for network(s) later)." - ] - }, - { - "cell_type": "markdown", - "id": "43849b09-5031-4eb0-989c-ceca2110f736", - "metadata": {}, - "source": [ - "We could have chosen here to leave the `nodes` argument off, since every `ChemicalSystem` we included was already represented among the `edges`, but we show it here for completeness. In this way, it's possible to include `ChemicalSystem`s in the network that aren't connected via any `Transformation`s to others, though in practice there isn't much utility in this." - ] - } - ], - "metadata": { - "kernelspec": { - "display_name": "Python 3 (ipykernel)", - "language": "python", - "name": "python3" - }, - "language_info": { - "codemirror_mode": { - "name": "ipython", - "version": 3 - }, - "file_extension": ".py", - "mimetype": "text/x-python", - "name": "python", - "nbconvert_exporter": "python", - "pygments_lexer": "ipython3", - "version": "3.10.11" - } - }, - "nbformat": 4, - "nbformat_minor": 5 -}