4.3.1. Hydrogen Bond Analysis — MDAnalysis.analysis.hydrogenbonds.hbond_analysis¶
- Author:
Paul Smith
- Year:
2019
- Copyright:
GNU Public License v3
New in version 1.0.0.
This module provides methods to find and analyse hydrogen bonds in a Universe.
The HydrogenBondAnalysis class is a new version of the original
MDAnalysis.analysis.hbonds.HydrogenBondAnalysis class from the module
MDAnalysis.analysis.hbonds.hbond_analysis, which itself was modeled after the VMD
HBONDS plugin.
Please cite [Smith2019] if you use this module in addition to the normal MDAnalysis citations.
4.3.1.1. Input¶
- Required:
universe : an MDAnalysis Universe object
- Options:
donors_sel [None] : Atom selection for donors. If None, then will be identified via the topology.
hydrogens_sel [None] : Atom selection for hydrogens. If None, then will be identified via charge and mass.
acceptors_sel [None] : Atom selection for acceptors. If None, then will be identified via charge.
d_h_cutoff (Å) [1.2] : Distance cutoff used for finding donor-hydrogen pairs
d_a_cutoff (Å) [3.0] : Distance cutoff for hydrogen bonds. This cutoff refers to the D-A distance.
d_h_a_angle_cutoff (degrees) [150] : D-H-A angle cutoff for hydrogen bonds.
update_selections [True] : If true, will update atom selections at each frame.
4.3.1.2. Output¶
frame : frame at which a hydrogen bond was found
donor id : atom id of the hydrogen bond donor atom
hydrogen id : atom id of the hydrogen bond hydrogen atom
acceptor id : atom id of the hydrogen bond acceptor atom
distance (Å): length of the hydrogen bond
angle (degrees): angle of the hydrogen bond
Hydrogen bond data are returned in a numpy.ndarray on a “one line, one observation” basis
and can be accessed via HydrogenBondAnalysis.results.hbonds:
results = [
[
<frame>,
<donor index (0-based)>,
<hydrogen index (0-based)>,
<acceptor index (0-based)>,
<distance>,
<angle>
],
...
]
4.3.1.3. Example use of HydrogenBondAnalysis¶
The simplest use case is to allow HydrogenBondAnalysis to guess the acceptor and hydrogen atoms, and to
identify donor-hydrogen pairs via the bonding information in the topology:
import MDAnalysis
from MDAnalysis.analysis.hydrogenbonds.hbond_analysis import HydrogenBondAnalysis as HBA
u = MDAnalysis.Universe(psf, trajectory)
hbonds = HBA(universe=u)
hbonds.run()
It is also possible to specify which hydrogens and acceptors to use in the analysis. For example, to find all hydrogen bonds in water:
import MDAnalysis
from MDAnalysis.analysis.hydrogenbonds.hbond_analysis import HydrogenBondAnalysis as HBA
u = MDAnalysis.Universe(psf, trajectory)
hbonds = HBA(universe=u, hydrogens_sel='resname TIP3 and name H1 H2', acceptors_sel='resname TIP3 and name OH2')
hbonds.run()
Alternatively, hydrogens_sel and acceptors_sel may be generated via the guess_hydrogens and
guess_acceptors. This selection strings may then be modified prior to calling run, or a subset of
the universe may be used to guess the atoms. For example, find hydrogens and acceptors belonging to a protein:
import MDAnalysis
from MDAnalysis.analysis.hydrogenbonds.hbond_analysis import HydrogenBondAnalysis as HBA
u = MDAnalysis.Universe(psf, trajectory)
hbonds = HBA(universe=u)
hbonds.hydrogens_sel = hbonds.guess_hydrogens("protein")
hbonds.acceptors_sel = hbonds.guess_acceptors("protein")
hbonds.run()
Slightly more complex selection strings are also possible. For example, to find hydrogen bonds involving a protein and any water molecules within 10 Å of the protein (which may be useful for subsequently finding the lifetime of protein-water hydrogen bonds or finding water-bridging hydrogen bond paths):
import MDAnalysis
from MDAnalysis.analysis.hydrogenbonds.hbond_analysis import HydrogenBondAnalysis as HBA
u = MDAnalysis.Universe(psf, trajectory)
hbonds = HBA(universe=u)
protein_hydrogens_sel = hbonds.guess_hydrogens("protein")
protein_acceptors_sel = hbonds.guess_acceptors("protein")
water_hydrogens_sel = "resname TIP3 and name H1 H2"
water_acceptors_sel = "resname TIP3 and name OH2"
hbonds.hydrogens_sel = f"({protein_hydrogens_sel}) or ({water_hydrogens_sel} and around 10 not resname TIP3})"
hbonds.acceptors_sel = f"({protein_acceptors_sel}) or ({water_acceptors_sel} and around 10 not resname TIP3})"
hbonds.run()
To calculate the hydrogen bonds between different groups, for example a
protein and water, one can use the between keyword. The
following will find protein-water hydrogen bonds but not protein-protein
or water-water hydrogen bonds:
import MDAnalysis
from MDAnalysis.analysis.hydrogenbonds.hbond_analysis import (
HydrogenBondAnalysis as HBA)
u = MDAnalysis.Universe(psf, trajectory)
hbonds = HBA(
universe=u,
between=['resname TIP3', 'protein']
)
protein_hydrogens_sel = hbonds.guess_hydrogens("protein")
protein_acceptors_sel = hbonds.guess_acceptors("protein")
water_hydrogens_sel = "resname TIP3 and name H1 H2"
water_acceptors_sel = "resname TIP3 and name OH2"
hbonds.hydrogens_sel = f"({protein_hydrogens_sel}) or ({water_hydrogens_sel}"
hbonds.acceptors_sel = f"({protein_acceptors_sel}) or ({water_acceptors_sel}"
hbonds.run()
It is further possible to compute hydrogen bonds between several groups with
with use of between. If in the above example,
between=[[‘resname TIP3’, ‘protein’], [‘protein’, ‘protein’]], all
protein-water and protein-protein hydrogen bonds will be found, but
no water-water hydrogen bonds.
In order to compute the hydrogen bond lifetime, after finding hydrogen bonds
one can use the lifetime function:
...
hbonds.run()
tau_timeseries, timeseries = hbonds.lifetime()
It is highly recommended that a topology with bond information is used to
generate the universe, e.g PSF, TPR, or PRMTOP files. This is the only
method by which it can be guaranteed that donor-hydrogen pairs are correctly
identified. However, if, for example, a PDB file is used instead, a
donors_sel may be provided along with a hydrogens_sel and the
donor-hydrogen pairs will be identified via a distance cutoff,
d_h_cutoff:
import MDAnalysis
from MDAnalysis.analysis.hydrogenbonds.hbond_analysis import (
HydrogenBondAnalysis as HBA)
u = MDAnalysis.Universe(pdb, trajectory)
hbonds = HBA(
universe=u,
donors_sel='resname TIP3 and name OH2',
hydrogens_sel='resname TIP3 and name H1 H2',
acceptors_sel='resname TIP3 and name OH2',
d_h_cutoff=1.2
)
hbonds.run()
4.3.1.4. The class and its methods¶
- class MDAnalysis.analysis.hydrogenbonds.hbond_analysis.HydrogenBondAnalysis(universe, donors_sel=None, hydrogens_sel=None, acceptors_sel=None, between=None, d_h_cutoff=1.2, d_a_cutoff=3.0, d_h_a_angle_cutoff=150, update_selections=True)[source]¶
Perform an analysis of hydrogen bonds in a Universe.
Set up atom selections and geometric criteria for finding hydrogen bonds in a Universe.
- Parameters:
universe (Universe) – MDAnalysis Universe object
donors_sel (str) – Selection string for the hydrogen bond donor atoms. If the universe topology contains bonding information, leave
donors_selas None so that donor-hydrogen pairs can be correctly identified.hydrogens_sel (str) – Selection string for the hydrogen bond hydrogen atoms. Leave as None to guess which hydrogens to use in the analysis using
guess_hydrogens. Ifhydrogens_selis left as None, also leavedonors_selas None so that donor-hydrogen pairs can be correctly identified.acceptors_sel (str) – Selection string for the hydrogen bond acceptor atoms. Leave as None to guess which atoms to use in the analysis using
guess_acceptorsbetween (List (optional),) – Specify two selection strings for non-updating atom groups between which hydrogen bonds will be calculated. For example, if the donor and acceptor selections include both protein and water, it is possible to find only protein-water hydrogen bonds - and not protein-protein or water-water - by specifying between=[“protein”, “SOL”]`. If a two-dimensional list is passed, hydrogen bonds between each pair will be found. For example, between=[[“protein”, “SOL”], [“protein”, “protein”]]` will calculate all protein-water and protein-protein hydrogen bonds but not water-water hydrogen bonds. If None, hydrogen bonds between all donors and acceptors will be calculated.
d_h_cutoff (float (optional)) – Distance cutoff used for finding donor-hydrogen pairs. Only used to find donor-hydrogen pairs if the universe topology does not contain bonding information
d_a_cutoff (float (optional)) – Distance cutoff for hydrogen bonds. This cutoff refers to the D-A distance.
d_h_a_angle_cutoff (float (optional)) – D-H-A angle cutoff for hydrogen bonds, in degrees.
update_selections (bool (optional)) – Whether or not to update the acceptor, donor and hydrogen lists at each frame.
Note
It is highly recommended that a universe topology with bond information is used, as this is the only way that guarantees the correct identification of donor-hydrogen pairs.
New in version 2.0.0: Added between keyword
- results.hbonds¶
A
numpy.ndarraywhich contains a list of all observed hydrogen bond interactions. See Output for more information.New in version 2.0.0.
- hbonds¶
Alias to the
results.hbondsattribute.Deprecated since version 2.0.0: Will be removed in MDAnalysis 3.0.0. Please use
results.hbondsinstead.[Smith2019]Paul Smith, Robert M. Ziolek, Elena Gazzarrini, Dylan M. Owen, and Christian D. Lorenz. On the interaction of hyaluronic acid with synovial fluid lipid membranes. Phys. Chem. Chem. Phys., 21:9845–9857, 2019. doi:10.1039/C9CP01532A.