-
Notifications
You must be signed in to change notification settings - Fork 35
SepTop Protocol #1057
New issue
Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.
By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.
Already on GitHub? Sign in to your account
SepTop Protocol #1057
Changes from all commits
b24ea6b
ee11781
06b04a2
fd81c70
f26f886
69e6d8c
4f59bc5
101f49a
98059a0
080526f
0bd8b7d
40c3772
7faae6e
3fa964f
46465fb
c54c113
9e711fd
d06d7ca
7f10674
74c88ab
956f3f9
0516aad
9166110
e22c703
13769c6
5515687
c000832
60db5a6
9ed4ee6
6ac3004
77bfba0
efc351a
f031006
4f1f7de
cc2aa27
6587df9
4047bd7
bd8c2d1
feb6a97
9cd84d9
fa1e497
99a0465
95bc3c3
02f02f8
7f36893
86ad7b5
bb880a7
99a8386
aa4f79d
5bc8202
4061098
6c5e40a
a8a7008
bb5b2ea
763f154
8155483
7ff06c6
8b3b533
406ee07
b1c2e9d
6ed8181
a5b78a9
0ca65e6
dd29739
5f9a1fc
3e8ee0b
896a5f7
d5c167a
a885191
8e6d340
02243de
f44a35c
a3ca799
cceaa4f
cb6dba7
23f1896
40a4db3
406d89c
af84399
c171fbc
88ea11e
729864b
9c2837f
dcac3a8
47128df
578ca4b
7d17558
dac4491
551caec
c7d8733
e820815
369888a
7d91feb
07d55f2
15e9183
c445c11
777e79e
ed60b46
dbc3afd
11731fd
12ef50f
94d68f2
bda1ce4
e73b680
affb2bb
610b556
88840d2
0715458
aa0fdcc
78354fa
90ee5d5
b288874
c5c710a
9b84f9d
0735f21
aeea9f4
b7a3a2e
8ab0a32
30a2f3e
8c311e9
1d70d60
21088fd
8318698
e08f9e8
67c6cc5
d413072
06e89f4
47c00b0
5f7716a
3852ea3
6a2cdb7
90aac5c
1e0e137
3097669
9c1c8e5
c629e96
6ffc25e
2ae1ebb
b476858
87b855b
09dbd4c
1a1fe22
094e81a
2d6bbc1
ad4d12e
b2f80d1
577ac3e
c2ee83b
c3b089b
b1f0b40
24ec9c0
d604e67
f8c0b42
60c3fa2
e1ca7c7
55c9824
7bec29a
c1550eb
a22fd33
ea3baa1
2befe9b
7093825
c7648cf
e93a6fc
29fe940
97d7387
f54be85
97dda95
9d71f6d
93c197d
3fc4eed
1b1234b
c2a4d1e
3f79a88
a709fbc
b0f1e4f
5926240
bb26a61
ac358e4
b549394
25c4989
cca8e09
fff7331
6b9242b
a151bf0
5c10c2c
89e8035
1da79c2
d9e0e18
41fdeca
04da3bb
1f0c084
373512e
09b4207
62b2f5e
a528910
6483176
820e2ca
e1de8a6
feab5c2
6b86d54
c41cc57
9c706d9
3e2e4a7
014cc4f
d9f4eb0
32184b2
274fa25
97c5dc5
5af185a
1579150
6091880
b2cdfa9
8592684
30c3c17
8a62d40
5e26f21
b61a5f4
770725d
4c8d3f3
52a3f77
7425e06
6936be6
f7711c7
a609f78
b39157a
e43f23f
8fe93a4
33d4f25
28eaf50
881126f
26b0091
1b29d4e
686631b
10c1981
141eb97
dcd429d
c3d3536
c53821a
2ea6d37
1edbb43
e82f5e7
825ddc8
dc31215
d1128c3
f28ff8e
c805a97
a82ef53
fac9aee
dffae05
File filter
Filter by extension
Conversations
Jump to
Diff view
Diff view
There are no files selected for viewing
| Original file line number | Diff line number | Diff line change |
|---|---|---|
| @@ -0,0 +1,131 @@ | ||
| Separated Topologies Protocol | ||
| ============================= | ||
|
|
||
| Overview | ||
| -------- | ||
|
|
||
| The :class:`SepTopProtocol <.SepTopProtocol>` [1]_, [2]_ calculates the difference in binding free energy between two ligands. | ||
| This protocol essentially performs two absolute binding free energy calculations simultaneously in opposite directions, | ||
| by (alchemically) inserting one ligand into the binding site, while removing the other ligand at the same time. | ||
| In contrast to the :ref:`RelativeHybridTopologyProtocol <userguide_relative_hybrid_topology_protocol>`, the two ligand topologies are | ||
| completely separate (meaning there is no common core), making atom mapping unnecessary and allowing transformations between chemically diverse ligands. | ||
|
|
||
| The relative binding free energy is calculated through a thermodynamic cycle by transforming one ligand into the other ligand | ||
| both in the solvent and in the binding site. | ||
|
|
||
| Restraints are required to keep the weakly | ||
| coupled and fully decoupled ligand in the binding site region and thereby reduce the phase | ||
| space that needs to be sampled. In the :class:`SepTopProtocol <.SepTopProtocol>` | ||
| we apply orientational, or Boresch-style, restraints, as described below. | ||
|
|
||
| In this cycle, the interactions of one molecule are turned off while simultaneously turning on interactions of the other molecule both in the solvent and complex phases. | ||
| The relative binding free energy is then obtained via summation of free energy differences along the thermodynamic cycle. | ||
|
|
||
| .. figure:: img/septop_cycle.png | ||
| :scale: 50% | ||
|
|
||
| Thermodynamic cycle for the SepTop free energy protocol. | ||
|
|
||
| Scientific Details | ||
| ------------------ | ||
|
|
||
| Orientational restraints | ||
| ~~~~~~~~~~~~~~~~~~~~~~~~ | ||
|
|
||
| Orientational, or Boresch-style, restraints are automaticallly (unless manually specified) applied between three protein and three ligand atoms using one bond, | ||
| two angle, and three dihedral restraints. Reference atoms are picked based on different criteria, such as the root mean squared | ||
| fluctuation of the atoms in a short MD simulation, the secondary structure of the protein, and the distance between atoms, based on heuristics from Baumann et al. [2]_. | ||
|
|
||
| Partial annihilation scheme | ||
| ~~~~~~~~~~~~~~~~~~~~~~~~~~~ | ||
|
|
||
| In the :class:`SepTopProtocol <.SepTopProtocol>` the coulombic interactions of the molecules are fully turned off (annihilated) in the respective non-interacting end states. | ||
| The Lennard-Jones interactions are instead decoupled, meaning the intermolecular interactions are turned off, keeping the intramolecular Lennard-Jones interactions. | ||
|
|
||
| The lambda schedule | ||
| ~~~~~~~~~~~~~~~~~~~ | ||
|
|
||
| Molecular interactions are modified along an alchemical path using a discrete set of lambda windows. | ||
| For the transformation of ligand A to ligand B in the binding site, the following steps are carried out, starting with ligand A being fully interacting in the binding site while ligand B is decoupled. | ||
|
|
||
| 1. Insert the non-interacting dummy ligand B into the binding site and restrain it using orientational restraints. The contribution of the restraints is calculated analytically. | ||
|
Member
There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. Separate issue raised: the description of the lambda schedule is non symmetric, i.e. number 5 describes transferring the ligand A to solvent, but none of these refer to transferring ligand B from the solvent (I think it's what is meant by "insert the non-interacting dummy", but the language probably needs to be improved). |
||
| 2. Turn on the van der Waals (vdW) interactions of ligand B while also turning on orientational restraints on ligand A. | ||
| 3. Turn on the electrostatic interactions of ligand B while at the same time turning off the electrostatics of ligand A. | ||
| 4. Turn off vdW interactions of ligand A while simultaneously releasing restraints on ligand B. | ||
| 5. Release the restraints of the now dummy ligand A analytically and transfer the ligand into the solvent. | ||
|
|
||
| The lambda schedule in the solvent phase is similar to the one in the complex, except that a single harmonic distance restraint is | ||
| applied between the respective central atom in the two ligands to keep the ligands apart while doing the alchemical transformation. | ||
| A soft-core potential from Beutler et al. [3]_ is applied to the Lennard-Jones potential to avoid instablilites in intermediate lambda windows. | ||
| The lambda schedule is defined in the ``lambda_settings`` objects ``lambda_elec_A``, ``lambda_elec_B``, ``lambda_vdw_A``, ``lambda_vdw_B``, | ||
| ``lambda_restraints_A``, and ``lambda_restraints_B``. | ||
|
|
||
| Simulation overview | ||
| ~~~~~~~~~~~~~~~~~~~ | ||
|
|
||
| The :class:`.ProtocolDAG` of the :class:`SepTopProtocol <.SepTopProtocol>` contains :class:`.ProtocolUnit`\ s from both the complex and solvent transformations. | ||
| This means that both legs of the thermodynamic cycle are constructed and run sequentially in the same :class:`.ProtocolDAG`. This is different from the :class:`.RelativeHybridTopologyProtocol` where the :class:`.ProtocolDAG` only runs a single leg of a thermodynamic cycle. | ||
| If multiple ``protocol_repeats`` are run (default: ``protocol_repeats=3``), the :class:`.ProtocolDAG` contains multiple :class:`.ProtocolUnit`\ s of both complex and solvent transformations. | ||
|
Member
There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. Unclear from this text if that means that the same restraint will be used for all repeats. |
||
|
|
||
| Simulation steps | ||
| """""""""""""""" | ||
|
|
||
| Each :class:`.ProtocolUnit` (whether complex or solvent) carries out the following steps: | ||
|
|
||
| 1. Parameterize the system using `OpenMMForceFields <https://github.com/openmm/openmmforcefields>`_ and `Open Force Field <https://github.com/openforcefield/openff-forcefields>`_. | ||
| 2. Equilibrate the fully interacting system using a short MD simulation using the same approach as the :class:`.PlainMDProtocol` (in the solvent leg this will include rounds of NVT and NPT equilibration). | ||
| 3. Add restraints to the system: Orientational restraints in the complex, a single harmonic distance restraint in the solvent leg. | ||
| 4. Create an alchemical system. | ||
| 5. Minimize the alchemical system. | ||
| 6. Equilibrate and production simulate the alchemical system using the chosen multistate sampling method (under NPT conditions). | ||
| 7. Analyze results for the transformation. | ||
|
|
||
|
|
||
| .. note:: Three different types of multistate sampling (i.e. replica swapping between lambda states) methods can be chosen; HREX, SAMS, and independent (no lambda swaps attempted). | ||
|
Member
There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. Note: do we have tests for SAMS and Indendent? Last I heard SAMS was broken - we might need to fix this elsewhere.
Contributor
Author
There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. No tests for those currently. Will raise an issue for that. |
||
| By default the HREX approach is selected, this can be altered using ``solvent_simulation_settings.sampler_method`` or ``complex_simulation_settings.sampler_method`` (default: ``repex``). | ||
|
|
||
|
|
||
| Simulation details | ||
| """""""""""""""""" | ||
|
|
||
| Here are some details of how the simulation is carried out which are not detailed in the :class:`SepTopProtocol <.SepTopProtocol>`: | ||
|
|
||
| * The protocol applies a `LangevinMiddleIntegrator <https://openmmtools.readthedocs.io/en/latest/api/generated/openmmtools.mcmc.LangevinDynamicsMove.html>`_ which uses Langevin dynamics, with the LFMiddle discretization [4]_. | ||
| * A `Monte Carlo barostat <https://docs.openmm.org/latest/api-python/generated/openmm.openmm.MonteCarloBarostat.html>`_ is used in the NPT ensemble to maintain constant pressure. | ||
|
|
||
| Getting the free energy estimate | ||
| ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ | ||
|
|
||
| The free energy differences are obtained from simulation data using the `MBAR estimator <https://www.alchemistry.org/wiki/Multistate_Bennett_Acceptance_Ratio>`_ (multistate Bennett acceptance ratio estimator) as implemented in the `PyMBAR package <https://pymbar.readthedocs.io/en/master/mbar.html>`_. | ||
| Both the MBAR estimates of the two legs of the thermodynamic cycle, and the overall relative binding free energy (of the entire cycle) are obtained, | ||
| which is different compared to the results in the :class:`.RelativeHybridTopologyProtocol` where results from two legs of the thermodynamic cycle are obtained separately. | ||
|
|
||
| In addition to the estimates of the free energy changes and their uncertainty, the protocol also returns some metrics to help assess convergence of the results, these are detailed in the :ref:`multistate analysis section <multistate_analysis>`. | ||
|
|
||
| See Also | ||
| -------- | ||
|
|
||
| **Tutorials** | ||
|
|
||
| * :any:`Separated Topologies Free Energies tutorial <../../tutorials/septop_tutorial>` | ||
|
|
||
| **Cookbooks** | ||
|
|
||
| :ref:`Cookbooks <cookbooks>` | ||
|
|
||
| **API Documentation** | ||
|
|
||
| * :ref:`OpenMM Protocol Settings <openmm protocol settings api>` | ||
|
|
||
| References | ||
| ---------- | ||
|
|
||
| * `pymbar <https://pymbar.readthedocs.io/en/stable/>`_ | ||
| * `yank <http://getyank.org/latest/>`_ | ||
| * `OpenMMTools <https://openmmtools.readthedocs.io/en/stable/>`_ | ||
| * `OpenMM <https://openmm.org/>`_ | ||
|
|
||
| .. [1] Separated topologies--a method for relative binding free energy calculations using orientational restraints, G. Rocklin, D. Mobley, K. Dill; Chem Phys, 2013; 138(8):085104. doi: 10.1063/1.4792251. | ||
| .. [2] Broadening the Scope of Binding Free Energy Calculations Using a Separated Topologies Approach, H. Baumann, E. Dybeck, C. McClendon, F. Pickard IV, V. Gapsys, L. Pérez-Benito, D. Hahn, G. Tresadern, A. Mathiowetz, D. Mobley, J. Chem. Theory Comput., 2023, 19, 15, 5058–5076 | ||
| .. [3] Avoiding singularities and numerical instabilities in free energy calculations based on molecular simulations, T.C. Beutler, A.E. Mark, R.C. van Schaik, P.R. Greber, and W.F. van Gunsteren, Chem. Phys. Lett., 222 529–539 (1994) | ||
| .. [4] Unified Efficient Thermostat Scheme for the Canonical Ensemble with Holonomic or Isokinetic Constraints via Molecular Dynamics, Zhijun Zhang, Xinzijian Liu, Kangyu Yan, Mark E. Tuckerman, and Jian Liu, J. Phys. Chem. A 2019, 123, 28, 6056-6079 | ||
| Original file line number | Diff line number | Diff line change |
|---|---|---|
| @@ -0,0 +1,83 @@ | ||
| OpenMM Separated Topologies Protocol | ||
| ==================================== | ||
|
|
||
| .. _septop protocol api: | ||
|
|
||
| This section provides details about the OpenMM Separated Topologies Protocol | ||
| implemented in OpenFE. | ||
|
|
||
| Protocol API specification | ||
| -------------------------- | ||
|
|
||
| .. module:: openfe.protocols.openmm_septop.equil_septop_method | ||
hannahbaumann marked this conversation as resolved.
Show resolved
Hide resolved
|
||
|
|
||
| .. autosummary:: | ||
| :nosignatures: | ||
| :toctree: generated/ | ||
|
|
||
| SepTopProtocol | ||
| SepTopComplexSetupUnit | ||
| SepTopComplexRunUnit | ||
| SepTopSolventSetupUnit | ||
| SepTopSolventRunUnit | ||
| SepTopProtocolResult | ||
|
|
||
| Protocol Settings | ||
| ----------------- | ||
|
|
||
| Below are the settings which can be tweaked in the protocol. The default settings (accessed using :meth:`SepTopProtocol.default_settings`) will automatically populate settings which we have found to be useful for running a Separated Topologies free energy calculation. There will however be some cases (such as when calculating difficult to converge systems) where you will need to tweak some of the following settings. | ||
|
|
||
|
|
||
| .. module:: openfe.protocols.openmm_septop.equil_septop_settings | ||
|
|
||
| .. autopydantic_model:: SepTopSettings | ||
| :model-show-json: False | ||
| :model-show-field-summary: False | ||
| :model-show-config-member: False | ||
| :model-show-config-summary: False | ||
| :model-show-validator-members: False | ||
| :model-show-validator-summary: False | ||
| :field-list-validators: False | ||
| :inherited-members: SettingsBaseModel | ||
| :exclude-members: get_defaults | ||
| :member-order: bysource | ||
|
|
||
|
|
||
| Protocol Specific Settings Classes | ||
| ---------------------------------- | ||
|
|
||
| Below are Settings classes which are unique to the `SepTopProtocol`. | ||
|
|
||
|
|
||
| .. autopydantic_model:: AlchemicalSettings | ||
| :model-show-json: False | ||
| :model-show-field-summary: False | ||
| :model-show-config-member: False | ||
| :model-show-config-summary: False | ||
| :model-show-validator-members: False | ||
| :model-show-validator-summary: False | ||
| :field-list-validators: False | ||
| :inherited-members: SettingsBaseModel | ||
| :member-order: bysource | ||
|
|
||
| .. autopydantic_model:: LambdaSettings | ||
| :model-show-json: False | ||
| :model-show-field-summary: False | ||
| :model-show-config-member: False | ||
| :model-show-config-summary: False | ||
| :model-show-validator-members: False | ||
| :model-show-validator-summary: False | ||
| :field-list-validators: False | ||
| :inherited-members: SettingsBaseModel | ||
| :member-order: bysource | ||
|
|
||
| .. autopydantic_model:: SepTopEquilOutputSettings | ||
| :model-show-json: False | ||
| :model-show-field-summary: False | ||
| :model-show-config-member: False | ||
| :model-show-config-summary: False | ||
| :model-show-validator-members: False | ||
| :model-show-validator-summary: False | ||
| :field-list-validators: False | ||
| :inherited-members: SettingsBaseModel | ||
| :member-order: bysource | ||
| Original file line number | Diff line number | Diff line change |
|---|---|---|
| @@ -0,0 +1,6 @@ | ||
| { | ||
| "path": "../ExampleNotebooks/openmm_septop/septop_tutorial.ipynb", | ||
| "extra-media": [ | ||
| "../ExampleNotebooks/openmm_septop/septop_cycle.png" | ||
| ] | ||
| } |
| Original file line number | Diff line number | Diff line change |
|---|---|---|
| @@ -0,0 +1,23 @@ | ||
| **Added:** | ||
|
|
||
| * Added a new RBFE protocol based on Separated Topologies. | ||
|
|
||
| **Changed:** | ||
|
|
||
| * <news item> | ||
|
|
||
| **Deprecated:** | ||
|
|
||
| * <news item> | ||
|
|
||
| **Removed:** | ||
|
|
||
| * <news item> | ||
|
|
||
| **Fixed:** | ||
|
|
||
| * <news item> | ||
|
|
||
| **Security:** | ||
|
|
||
| * <news item> |
| Original file line number | Diff line number | Diff line change |
|---|---|---|
| @@ -0,0 +1,29 @@ | ||
| # This code is part of OpenFE and is licensed under the MIT license. | ||
| # For details, see https://github.com/OpenFreeEnergy/openfe | ||
| """ | ||
| Run SepTop free energy calculations using OpenMM and OpenMMTools. | ||
|
|
||
| """ | ||
|
|
||
| from .equil_septop_settings import ( | ||
| SepTopSettings, | ||
| ) | ||
|
|
||
| from .equil_septop_method import ( | ||
| SepTopProtocol, | ||
| SepTopProtocolResult, | ||
| SepTopComplexSetupUnit, | ||
| SepTopComplexRunUnit, | ||
| SepTopSolventSetupUnit, | ||
| SepTopSolventRunUnit, | ||
| ) | ||
|
|
||
| __all__ = [ | ||
| "SepTopProtocol", | ||
| "SepTopSettings", | ||
| "SepTopProtocolResult", | ||
| "SepTopComplexSetupUnit", | ||
| "SepTopSolventSetupUnit", | ||
| "SepTopSolventRunUnit", | ||
| "SepTopComplexRunUnit", | ||
| ] |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
TODO (issued raised): update the docs here (and in ABFEs) to reflect the new picking, i.e. it's getting closer to a mix between Hannah's appproach and MDRestraintsGenerator.