[gmx-developers] creating a contact matrix from .xpm file from g_hbond

[Andrew Jewett]

I wrote a simple program in python can reformat the hydrogen-bond
information stored in the .xpm files created by "g_hbond -hbm" in a
way that is a little bit easier to read (both for people and for
programs).

I have been testing it using the CVS version of g_hbond (version
3.3.1), and it works.

For every frame in the trajectory, this program generates a
"connectivity matrix" (or "contact matrix?").  The following example
shows the connectivity between the different amino acids in a
simulation consisting of two short peptides (each 11 residues long).
Because the bonds tend to follow the main diagonal, you can see that
the two strands are in a (mostly) parallel fashion (as opposed to
anti-parallel).

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

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.

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:

 This program can consume lots of memory, due to the order that data is
 saved to the xpm file, especially for long trajectories.  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.  You can run xpm2matrix
 multiple times and concatenate the results together.