[gmx-developers] reposting the code for hbond analysis program (xpm2matrix)

[Andrew Jewett]

   I'm disappointed to learn that the source code for the h-bond
analysis script that I emailed you last month (xpm2matrix) was
reformatted and will not compile.  This time I'm sending it as an
attachment.  Please forgive me for the double-post.
A general description of this program is included below:

---------- Forwarded message ----------
Date: Aug 16, 2007 11:16 PM
Subject: creating a contact matrix from .xpm file from g_hbond


   I wrote a simple program in python to help interpret hydrogen bond
information from gromacs trajectories.  It generates a "connectivity matrix"
for every frame in the gromacs trajectory, indicating exactly which
groups of atoms are h-bonded together.  Actually, all it does is reformat
the output of g_hbond -hbm to make it easier to read.  The following
example shows the connectivity matrix between the various amino acids
in a simulation consisting of two short peptides (each 11 residues long).

0 0 0 1 0 1 0 0 0 0 0
0 0 1 0 0 0 0 0 0 0 0
0 1 0 0 1 0 0 0 0 0 0
0 0 0 0 0 1 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 1 0 1 0 0
0 0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 1 0 0
0 0 0 0 0 0 0 0 0 0 2

   Why is this useful?
In my own work, I haved used this matrix to quantify exactly how "parallel"
or "anti-parallel" is the interaction between two peptides.
In this example, you can see that the position of the bonds is correlated.
The bonds tend to follow the main diagonal, but the
upper-right or lower left regions of the matrix are filled with zeros.
>From you can see that the two strands are in a (mostly) parallel fashion
(as opposed to anti-parallel).

There are many other ways this matrix could be used, and I will send a couple
programs that do this.

All of this code has been tested using the CVS version of g_hbond (version
3.3.1).  You will need to have python installed (preferably in
/usr/bin) to run the script.
Windows users can change the name of the file to add a .py extension
and run it in through python interpreter.

One final note:  If Erik decides to reverse the order of the bonds in
the "hbond.ndx" file generated by g_hbond, then simply uncomment lines
149 and 150 of the xpm2matrix code included in this attachment.

------------------------------------------------------------------------------------

General description: (I will send a more detailed documentation file soon)

    This program reads .xpm files generated by "g_hbond -hbm ...",
 and, for each frame in the trajectory, prints out a 2-dimensional table
 indicating which pairs of atoms (or residues, or groups of atoms)
 are connected together (a connectivity matrix). This table may be easier
 to interpret than the original .xpm file.  As an input it requires a pair
 of files ("groups1.ndx" and "groups2.ndx" in this example) and each file
 contains a list of groups.  (Each group typically represents the atoms
 belonging to each amino acid or nucleotide in a molecule.)  Since there
 are two lists of groups, you can look at the connectivity matrix between
 two different molecules.  The the two lists of groups do not have to
 correspond to residues from different molecules.  They can be identical.
     The atoms in each group are not required to be donors or acceptors
 (other atoms will be ignored).  It is possible to tell xpm2matrix
 to only consider the donor atoms from one group, and the acceptor atoms
 from another group.   For example, one of list of groups
(say "groups2.ndx") might contain only a single group consisting of all
of the solvent atoms.  In that case, xpm2matrix would simply report a row
of numbers for each frame in the trajectory, each number indicating the
number of solvent h-bonds connected to each residue/group.

Output:

    For each frame in the trajectory used to generate the xpm file, this
 program generates a connectivity matrix (delimited by spaces and newlines):

c11 c12 c13 ... c1n
c21 c22 c23 ... c2n
c31 c32 c33 ... c3n
 :   :   :       :
cm1 cm2 cm3 ... cmn

 where cij = the number of hydrogen bonds connecting group i (from the
             first list of groups) with group j (from the second list)

Input:

 hbmap.xpm    An .xpm file created by "g_hbond -hbm" that you want to analyze
              (See documentation for g_hbond.)

Arguments:

 hbond.ndx    The ndx file created by "g_hbond -hbn ..."
 groups1.ndx     The number of hydrogen bonds is counted between every pair
 groups2.ndx     of groups in these two index files: groups1.ndx & groups2.ndx.
 DA or AD     The 4th argument DA / AD is optional.

              By using the "DA" and "AD" arguments, the user can limit
              the atoms that the program considers from groups1.ndx
              and groups2.ndx.  The DA command line argument causes
              xpm2matrix to ignore all atoms from the first list of
              groups which are not potential donor atoms, and all atoms
              from the second list which are not potential acceptors.
              (The reverse is true for the AD command.)

Performance Issues:

 If the program is performing very slowly (running out of memory),
 you can get around the problem by running g_hbond multiple times
 using the "-b" and "-e" flags, in order to generate multiple smaller
 .xpm files.  Then you can run xpm2matrix multiple times and concatenate
 the results together.
-------------- next part --------------
A non-text attachment was scrubbed...
Name: xpm2matrix
Type: application/octet-stream
Size: 13574 bytes
Desc: not available
URL: <http://maillist.sys.kth.se/pipermail/gromacs.org_gmx-developers/attachments/20070907/e2fea4b4/attachment.obj>