{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# Calculating hydrogen bonds: advanced selections\n", "\n", "We will find all intramolcular hydrogen bonds of a protein without passing any atom selections. We will also look at how to use more advanced atom selections than we saw in [Calculating hydrogen bonds: the basiscs](hbonds.ipynb).\n", "\n", "**Last updated:** June 23, 2021 with MDAnalysis 2.0.0-dev\n", "\n", "**Minimum version of MDAnalysis:** 2.0.0-dev0\n", "\n", "**Packages required:**\n", " \n", "* MDAnalysis (Michaud-Agrawal *et al.*, 2011, Gowers *et al.*, 2016)\n", "* MDAnalysisTests\n", "* [numpy](https://numpy.org)\n", "\n", "**See also**\n", "\n", "* [Calculating hydrogen bonds: the basics](hbonds.ipynb)\n", "* [Calculating hydrogen bond lifetimes](hbonds-lifetimes.ipynb)\n", "\n", "
\n", " \n", "**Note**\n", "\n", "Please cite [Smith et al. (2018)](http://dx.doi.org/10.1039/C9CP01532A) when using HydrogenBondAnaysis in published work.\n", "\n", "
" ] }, { "cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [], "source": [ "import numpy as np\n", "\n", "import MDAnalysis as mda\n", "from MDAnalysis.tests.datafiles import PSF, DCD\n", "from MDAnalysis.analysis.hydrogenbonds import HydrogenBondAnalysis\n" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Loading files" ] }, { "cell_type": "code", "execution_count": 2, "metadata": {}, "outputs": [], "source": [ "u = mda.Universe(PSF, DCD)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The test files we will be working with here feature adenylate kinase (AdK), a phosophotransferase enzyme. (Beckstein *et al.*, 2009)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Find all hydrogen bonds\n", "The simplest use case is to allow `HydrogenBondAnalysis` to guess the acceptor and hydrogen atoms, then to identify donor-hydrogen pairs via the bond information in the topology.\n", "\n", "Accetptor and hydrogen atoms are guessed based on the mass and partial charge of the atoms.\n" ] }, { "cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [], "source": [ "hbonds = HydrogenBondAnalysis(universe=u)" ] }, { "cell_type": "code", "execution_count": 4, "metadata": {}, "outputs": [ { "data": { "application/vnd.jupyter.widget-view+json": { "model_id": "b250ac8556114eda9bba2f2675130383", "version_major": 2, "version_minor": 0 }, "text/plain": [ " 0%| | 0/98 [00:00\n", " \n", "**Note**\n", " \n", "The Universe we are analysing has the water removed. The above example is for illustrative purposes only.\n", " \n", "" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Hydrogen bonds between specific groups\n", "\n", "To calculate the hydrogen bonds between different groups, for example protein and water, one can use the `between` keyword. Below we will find protein-water and protein-protein hydrogen bonds, but not water-water hydrogen bonds. We do this by passing atom selection for the groups we wish to find hydrogen bonds between." ] }, { "cell_type": "code", "execution_count": 10, "metadata": {}, "outputs": [], "source": [ "hbonds = HydrogenBondAnalysis(\n", " universe=u,\n", " between=[\n", " [\"protein\", \"resname TIP3\"], # for protein-water hbonds\n", " [\"protein\", \"protein\"] # for protein-protein hbonds\n", " ]\n", ")\n" ] }, { "cell_type": "code", "execution_count": 11, "metadata": {}, "outputs": [ { "data": { "application/vnd.jupyter.widget-view+json": { "model_id": "dc286e109cf14ea388e324ba70946368", "version_major": 2, "version_minor": 0 }, "text/plain": [ " 0%| | 0/98 [00:00\n", " \n", "**Note**\n", " \n", "The Universe we are analysing has the water removed. The above example is for illustrative purposes only.\n", " \n", "" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## References\n", "\n", "[1] Richard J. Gowers, Max Linke, Jonathan Barnoud, Tyler J. E. Reddy, Manuel N. Melo, Sean L. Seyler, Jan Domański, David L. Dotson, Sébastien Buchoux, Ian M. Kenney, and Oliver Beckstein.\n", "MDAnalysis: A Python Package for the Rapid Analysis of Molecular Dynamics Simulations.\n", "Proceedings of the 15th Python in Science Conference, pages 98–105, 2016.\n", "00152.\n", "URL: https://conference.scipy.org/proceedings/scipy2016/oliver_beckstein.html, doi:10.25080/Majora-629e541a-00e.\n", "\n", "[2] Naveen Michaud-Agrawal, Elizabeth J. Denning, Thomas B. Woolf, and Oliver Beckstein.\n", "MDAnalysis: A toolkit for the analysis of molecular dynamics simulations.\n", "Journal of Computational Chemistry, 32(10):2319–2327, July 2011.\n", "00778.\n", "URL: http://doi.wiley.com/10.1002/jcc.21787, doi:10.1002/jcc.21787.\n", "\n", "[3] Paul Smith, Robert M. Ziolek, Elena Gazzarrini, Dylan M. Owen, and Christian D. Lorenz.\n", "On the interaction of hyaluronic acid with synovial fluid lipid membranes.\n", "Phys. Chem. Chem. Phys., 21(19):9845-9857, 2018.\n", "URL: http://dx.doi.org/10.1039/C9CP01532A\n", "\n", "[4] Oliver Beckstein, Elizabeth J. Denning, Juan R. Perilla, and Thomas B. Woolf.\n", "Zipping and Unzipping of Adenylate Kinase: Atomistic Insights into the Ensemble of OpenClosed Transitions.\n", "Journal of Molecular Biology, 394(1):160–176, November 2009.\n", "00107.\n", "URL: https://linkinghub.elsevier.com/retrieve/pii/S0022283609011164, doi:10.1016/j.jmb.2009.09.009.\n", "\n" ] } ], "metadata": { "celltoolbar": "Raw Cell Format", "kernelspec": { "display_name": "Python [conda env:mda-user-guide]", "language": "python", "name": "conda-env-mda-user-guide-py" }, "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.8.5" }, "toc": { "base_numbering": 1, "nav_menu": {}, "number_sections": false, "sideBar": true, "skip_h1_title": false, "title_cell": "Table of Contents", "title_sidebar": "Contents", "toc_cell": false, "toc_position": {}, "toc_section_display": true, "toc_window_display": true }, "widgets": { "application/vnd.jupyter.widget-state+json": { "state": { "0234c4febb954879a30d246dea830afa": { "model_module": "@jupyter-widgets/controls", "model_module_version": "1.5.0", "model_name": "DescriptionStyleModel", "state": { "description_width": "" } }, "02c2a05dba9d4ea0bf1840387af43607": { "model_module": "@jupyter-widgets/controls", "model_module_version": "1.5.0", "model_name": "HTMLModel", "state": { "layout": "IPY_MODEL_5533b16b291a40789f736f3d53916ea2", "style": "IPY_MODEL_6ab80580d1314dccb4fad4c8f4964d2f", "value": " 98/98 [00:09<00:00, 9.31it/s]" } }, "050ae1e0b9db4d41ab07d70f7340e4ba": { "model_module": "@jupyter-widgets/controls", "model_module_version": "1.5.0", "model_name": "HTMLModel", "state": { "layout": "IPY_MODEL_ba52338f1cb540c5b2e3bea90a82135d", "style": "IPY_MODEL_be3afb50421e4d99aabea002929a8215", "value": "100%" } }, "05a42445d9a544dca7c7b9d5613ff433": { "model_module": "@jupyter-widgets/base", "model_module_version": "1.2.0", "model_name": "LayoutModel", "state": {} }, "1894730a30084ab6b591f7d7fa0e76d0": { "model_module": "@jupyter-widgets/controls", "model_module_version": "1.5.0", "model_name": "DescriptionStyleModel", "state": { "description_width": "" } }, "18f91f27ac184092a093f5f0226b7ce9": { "model_module": "@jupyter-widgets/base", "model_module_version": "1.2.0", "model_name": "LayoutModel", "state": {} }, "2fb5cee414374def9e06e658182fec15": { "model_module": "@jupyter-widgets/controls", "model_module_version": "1.5.0", "model_name": "ProgressStyleModel", "state": { "description_width": "" } }, "3e6ee95e18b4465582fd42a7f995901d": { "model_module": "@jupyter-widgets/base", "model_module_version": "1.2.0", "model_name": "LayoutModel", "state": {} }, "49331ea4ddde4282a74f1cb6c1ebe776": { "model_module": "@jupyter-widgets/controls", "model_module_version": "1.5.0", "model_name": "DescriptionStyleModel", "state": { "description_width": "" } }, "5533b16b291a40789f736f3d53916ea2": { "model_module": "@jupyter-widgets/base", "model_module_version": "1.2.0", "model_name": "LayoutModel", "state": {} }, "59b9ba3247a74156a18e50e38b317f42": { "model_module": "@jupyter-widgets/controls", "model_module_version": "1.5.0", "model_name": "HTMLModel", "state": { "layout": "IPY_MODEL_aa7b0762cb3643c29b40dd0b3abf7650", "style": "IPY_MODEL_1894730a30084ab6b591f7d7fa0e76d0", "value": "100%" } }, "5a55198aba0a45728d4188874af326e7": { "model_module": "@jupyter-widgets/controls", "model_module_version": "1.5.0", "model_name": "HTMLModel", "state": { "layout": "IPY_MODEL_b8e364e4ea6944cf821f8bf3045ef4c1", "style": "IPY_MODEL_bb54c783014b491fa46c88100189b4a7", "value": " 98/98 [00:03<00:00, 31.66it/s]" } }, "5ff94c37bd9c4f6a85940ba2ca2b8cc8": { "model_module": "@jupyter-widgets/base", "model_module_version": "1.2.0", "model_name": "LayoutModel", "state": {} }, "6132218a274d4a39b8f0eb456894c5e1": { "model_module": "@jupyter-widgets/controls", "model_module_version": "1.5.0", "model_name": "HBoxModel", "state": { "children": [ "IPY_MODEL_89db791279874d749dd4ca3340a3b5e2", "IPY_MODEL_e37df7b9be4a46a9b4e0996eee7231d7", "IPY_MODEL_6ddd96ffd6ee432aad8ae92986dfa40b" ], "layout": "IPY_MODEL_dfe10a9d20f744719cc9bc5dc47a03cf" } }, "682c51aa834e47c6b314a5a57019b83b": { "model_module": "@jupyter-widgets/base", "model_module_version": "1.2.0", "model_name": "LayoutModel", "state": {} }, "6ab80580d1314dccb4fad4c8f4964d2f": { "model_module": "@jupyter-widgets/controls", "model_module_version": "1.5.0", "model_name": "DescriptionStyleModel", "state": { "description_width": "" } }, "6ddd96ffd6ee432aad8ae92986dfa40b": { "model_module": "@jupyter-widgets/controls", "model_module_version": "1.5.0", "model_name": "HTMLModel", "state": { "layout": "IPY_MODEL_5ff94c37bd9c4f6a85940ba2ca2b8cc8", "style": "IPY_MODEL_49331ea4ddde4282a74f1cb6c1ebe776", "value": " 98/98 [00:09<00:00, 10.24it/s]" } }, "6fb8a5a9a0a84c3eb71cb03afb3c65f0": { "model_module": "@jupyter-widgets/base", "model_module_version": "1.2.0", "model_name": "LayoutModel", "state": {} }, "78722f3bfe9742b29450019f4031a8f6": { "model_module": "@jupyter-widgets/controls", "model_module_version": "1.5.0", "model_name": "FloatProgressModel", "state": { "bar_style": "success", "layout": "IPY_MODEL_6fb8a5a9a0a84c3eb71cb03afb3c65f0", "max": 98, "style": "IPY_MODEL_cc2b120548e3467fb923f411cb52b24e", "value": 98 } }, "80d76bac7f894582bf10091d3cee3bda": { "model_module": "@jupyter-widgets/controls", "model_module_version": "1.5.0", "model_name": "ProgressStyleModel", "state": { "description_width": "" } }, "89db791279874d749dd4ca3340a3b5e2": { "model_module": "@jupyter-widgets/controls", "model_module_version": "1.5.0", "model_name": "HTMLModel", "state": { "layout": "IPY_MODEL_ff76f3ec14104842abae9d4b63379007", "style": "IPY_MODEL_0234c4febb954879a30d246dea830afa", "value": "100%" } }, "a62ebbb093504dbc90318cde0396d094": { "model_module": "@jupyter-widgets/controls", "model_module_version": "1.5.0", "model_name": "DescriptionStyleModel", "state": { "description_width": "" } }, "aa7b0762cb3643c29b40dd0b3abf7650": { "model_module": "@jupyter-widgets/base", "model_module_version": "1.2.0", "model_name": "LayoutModel", "state": {} }, "b250ac8556114eda9bba2f2675130383": { "model_module": "@jupyter-widgets/controls", "model_module_version": "1.5.0", "model_name": "HBoxModel", "state": { "children": [ "IPY_MODEL_59b9ba3247a74156a18e50e38b317f42", "IPY_MODEL_78722f3bfe9742b29450019f4031a8f6", "IPY_MODEL_daa341a4dd454d31a9f20b94654836d4" ], "layout": "IPY_MODEL_d1f2b3d6c1f943f983d95dc076c10486" } }, "b8e364e4ea6944cf821f8bf3045ef4c1": { "model_module": "@jupyter-widgets/base", "model_module_version": "1.2.0", "model_name": "LayoutModel", "state": {} }, "ba52338f1cb540c5b2e3bea90a82135d": { "model_module": "@jupyter-widgets/base", "model_module_version": "1.2.0", "model_name": "LayoutModel", "state": {} }, "bb54c783014b491fa46c88100189b4a7": { "model_module": "@jupyter-widgets/controls", "model_module_version": "1.5.0", "model_name": "DescriptionStyleModel", "state": { "description_width": "" } }, "be3afb50421e4d99aabea002929a8215": { "model_module": "@jupyter-widgets/controls", "model_module_version": "1.5.0", "model_name": "DescriptionStyleModel", "state": { "description_width": "" } }, "c5d77db92ba44e0ab4beb82dd9c39493": { "model_module": "@jupyter-widgets/controls", "model_module_version": "1.5.0", "model_name": "FloatProgressModel", "state": { "bar_style": "success", "layout": "IPY_MODEL_682c51aa834e47c6b314a5a57019b83b", "max": 98, "style": "IPY_MODEL_80d76bac7f894582bf10091d3cee3bda", "value": 98 } }, "ca060a426ae7417e9af161559ffd6ecc": { "model_module": "@jupyter-widgets/controls", "model_module_version": "1.5.0", "model_name": "HBoxModel", "state": { "children": [ "IPY_MODEL_cc2c0bdcf8f54624be570a8ed6355795", "IPY_MODEL_e13fdfbb8d1444b3842f444e47fcf35f", "IPY_MODEL_5a55198aba0a45728d4188874af326e7" ], "layout": "IPY_MODEL_f0d2ed5157c642bd983bd287bac048fd" } }, "cc2b120548e3467fb923f411cb52b24e": { "model_module": "@jupyter-widgets/controls", "model_module_version": "1.5.0", "model_name": "ProgressStyleModel", "state": { "description_width": "" } }, "cc2c0bdcf8f54624be570a8ed6355795": { "model_module": "@jupyter-widgets/controls", "model_module_version": "1.5.0", "model_name": "HTMLModel", "state": { "layout": "IPY_MODEL_18f91f27ac184092a093f5f0226b7ce9", "style": "IPY_MODEL_a62ebbb093504dbc90318cde0396d094", "value": "100%" } }, "ceb88130dc8940afbaa53ad0c57f9645": { "model_module": "@jupyter-widgets/base", "model_module_version": "1.2.0", "model_name": "LayoutModel", "state": {} }, "d1f2b3d6c1f943f983d95dc076c10486": { "model_module": "@jupyter-widgets/base", "model_module_version": "1.2.0", "model_name": "LayoutModel", "state": {} }, "daa341a4dd454d31a9f20b94654836d4": { "model_module": "@jupyter-widgets/controls", "model_module_version": "1.5.0", "model_name": "HTMLModel", "state": { "layout": "IPY_MODEL_ceb88130dc8940afbaa53ad0c57f9645", "style": "IPY_MODEL_e15db1b4181c4b009979903348388ca1", "value": " 98/98 [00:09<00:00, 10.47it/s]" } }, "dc286e109cf14ea388e324ba70946368": { "model_module": "@jupyter-widgets/controls", "model_module_version": "1.5.0", "model_name": "HBoxModel", "state": { "children": [ "IPY_MODEL_050ae1e0b9db4d41ab07d70f7340e4ba", "IPY_MODEL_c5d77db92ba44e0ab4beb82dd9c39493", "IPY_MODEL_02c2a05dba9d4ea0bf1840387af43607" ], "layout": "IPY_MODEL_ddc08bfd042f451e9a86c69f19937920" } }, "ddc08bfd042f451e9a86c69f19937920": { "model_module": "@jupyter-widgets/base", "model_module_version": "1.2.0", "model_name": "LayoutModel", "state": {} }, "dfe10a9d20f744719cc9bc5dc47a03cf": { "model_module": "@jupyter-widgets/base", "model_module_version": "1.2.0", "model_name": "LayoutModel", "state": {} }, "e13fdfbb8d1444b3842f444e47fcf35f": { "model_module": "@jupyter-widgets/controls", "model_module_version": "1.5.0", "model_name": "FloatProgressModel", "state": { "bar_style": "success", "layout": "IPY_MODEL_05a42445d9a544dca7c7b9d5613ff433", "max": 98, "style": "IPY_MODEL_fb3c902982a241be823bc7c7fc10ad70", "value": 98 } }, "e15db1b4181c4b009979903348388ca1": { "model_module": "@jupyter-widgets/controls", "model_module_version": "1.5.0", "model_name": "DescriptionStyleModel", "state": { "description_width": "" } }, "e37df7b9be4a46a9b4e0996eee7231d7": { "model_module": "@jupyter-widgets/controls", "model_module_version": "1.5.0", "model_name": "FloatProgressModel", "state": { "bar_style": "success", "layout": "IPY_MODEL_3e6ee95e18b4465582fd42a7f995901d", "max": 98, "style": "IPY_MODEL_2fb5cee414374def9e06e658182fec15", "value": 98 } }, "f0d2ed5157c642bd983bd287bac048fd": { "model_module": "@jupyter-widgets/base", "model_module_version": "1.2.0", "model_name": "LayoutModel", "state": {} }, "fb3c902982a241be823bc7c7fc10ad70": { "model_module": "@jupyter-widgets/controls", "model_module_version": "1.5.0", "model_name": "ProgressStyleModel", "state": { "description_width": "" } }, "ff76f3ec14104842abae9d4b63379007": { "model_module": "@jupyter-widgets/base", "model_module_version": "1.2.0", "model_name": "LayoutModel", "state": {} } }, "version_major": 2, "version_minor": 0 } } }, "nbformat": 4, "nbformat_minor": 4 }