Identify Hydrogen Bonds¶

Use the functionality of hbonds.py in order to identify hydrogen bonds in a molecular dynamics trajectory.

Attempts to describe a hydrogen bond in the context of molecular simulations has yielded a number of definitions with no true consensus cite{prada-gracia_quest_2013} especially since the geometry of hydrogen bonds has some dependence on the system being studied. Luzar and Chandler proposed the geometric criterion such that a hydrogen bond exists if the distance between the donor, D, and acceptor, A, atoms is less than 3.5 AA~and the angle formed by D–H$cdotcdotcdot$A is less than 30degree. cite{luzar_effect_1996} The definition of Luzar and Chandler is easily visualized for trajectories using the texttt{hbonds} representation of the Visual Molecular Dynamics (VMD) software package which allows us to directly check the validity of identified hydrogen bonds.

Some definitions:

donor (D): atom covalently bonded to hydrogen
acceptor (A) : atom which ‘accepts’ hydrogen bond

A diagram of an hbond:

D–H - - A

Criterion:

Distance between D and A below some distance
Angle between DHA less than some cut-off

Note: You will need to properly annotate any residue that you think might participate in a hydrogen bond. See Summary of Annotations.

Classes¶

class hbonds.System(traj, gro, begin=0, end=-1, skip=1, t=None)¶

__init__(traj, gro, begin=0, end=-1, skip=1, t=None)¶

Load in the trajectory

Parameters:	traj (str or NoneType) – GROMACS trajectory file (.xtc or .trr), can be None if t is provided gro (str) – GROMACS coordinate file (.gro) begin (int) – First frame to analyze end (int) – Last frame to analyze skip (int) – Only analyze every skip frames t (object) – mdtraj trajectory object. If this is provided, traj and gro will not be loaded

identify_hbonds(cut, angle)¶

Identify hydrogen bonds based on geometric criteria. If the angle formed by D-H–A and the distance between donor and acceptor are less than the cut-offs, then we consider a hydrogen bond to exist

Parameters:	cut (float) – cut-off distance (nm) angle (float) – cut-off angle (degrees)

plot_hbonds(show=True, save=True, savename='hbonds.png')¶

Plot the total number of hydrogen bonds per frame as a function of time

Parameters:	show (bool) – display plot save (bool) – save plot under savename savename (str) – name under which to save plot

set_eligible(res, atoms, acceptor_only=False, donors_only=False)¶

Define atoms that will be included in h-bond calculation. Make sure res is updated with proper annotations or no atoms will be selected!

Parameters:	res (str) – residue to include in calculation atoms (str or list) – atoms from residue to include in calculation. Use ‘all’ to include all potential hbonding atoms acceptor_only (bool) – restrict residue so it can only accept hydrogen bonds donor_only – restrict residue so it can only donate hydrogen bonds

Examples¶

# imports
from LLC_Membranes.analysis import hbonds

# Initialize system
sys = hbonds.System(traj, gro)

# Define residues involved in hbonds of interest
residues = ['HII', 'HOH'] # hydrogen bonds between LLC monomer and water

# Restrict system to hbonding atoms that are a part of HII and HOH
for i, r, in enumerate(residues):
    sys.set_eligible(r, 'all')  # restrict to all hbonding atoms of each residue

# Apply geometric criteria to identify hydrogen bonds
distance = 0.35  # distance cut-off
angle = 30
sys.identify_hbonds(distance, angle)

# Plot total hydrogen bonds as a function of time
sys.plot_hbonds()

Command Line Functionality¶

Use geometric criteria to identify hydrogen bonds

usage: hbonds.py [-h] [-t TRAJ] [-g GRO] [-b BEGIN] [-e END] [-sk SKIP]
                 [-r RESIDUES [RESIDUES ...]] [-a ATOMS [ATOMS ...]]
                 [-d DISTANCE] [-angle ANGLE_CUT] [-l LOAD]
                 [-acc ACCEPTORS [ACCEPTORS ...]]
                 [-donors DONORS [DONORS ...]] [-s SAVENAME] [-noshow]

Named Arguments¶

`-t, --traj`	Path to input file Default: “wiggle.trr”
`-g, --gro`	Coordinate file Default: “wiggle.gro”
`-b, --begin`	End frame Default: 0
`-e, --end`	Start frame Default: -1
`-sk, --skip`	Skip every skip frames Default: 1
`-r, --residues`	Residues to include in h-bond search Default: [‘HII’]
`-a, --atoms`	Atoms to include for each residue. Each list of atoms must be passed with a separate -a flag for each residue inargs.residues
`-d, --distance`	Maximum distance between acceptor and donor atoms Default: 0.35
`-angle, --angle_cut`
	Maximum DHA angle to be considered an H-bond Default: 30
`-l, --load`	Load pickled system object Default: False
`-acc, --acceptors`
	If you only want hbonds with acceptoratoms of a residue, use this flag. 1 for this condition, otherwise 0. Input as a list of 0 and 1 in the same order as args.residues Default: False
`-donors, --donors`
	If you only want hbonds with donoratoms of a residue, use this flag. 1 for this condition, otherwise 0. Input as a list of 0 and 1 in the same order as args.residues Default: False
`-s, --savename`	Name of pickled object to save aftercalcualtions are complete Default: “hbonds.pl”
`-noshow, --noshow`
	Do not display plots Default: False

Example Usage¶

# hydrogen bonds between atoms O3 and O4 of the monomer HII with all atoms of water
hbonds.py -t trajectory.xtc -g coordinates.gro -r HII HOH -a O3 O4 -a all