GROMACS hydrogen-bond analysis question

GROMACS version: 2023.3
GROMACS modification: No

Hi GROMACS community,

I would like to use the GROMACS hydrogen-bond analysis tools to calculate the number of hydrogen bonds involving water in my cellulose trajectory. Since my system is a cellulose polymer, I want to distinguish the hydrogen bonds contributed by water to the inner and outer cellulose layers. In the figure below, the red chains represent the inner group, while the black chains represent the outer group.

I have tested two different approaches:

  1. Calculate hydrogen bonds for each individual chain and then classify them: In this approach, I create an index group for each cellulose chain separately. I then calculate the hydrogen bonds involving water for each chain and sum the results according to the inner/outer classification shown above.

  2. Create inner and outer groups first, then calculate the hydrogen bonds: In this approach, I create two index groups at the beginning: one containing all inner chains and another containing all outer chains. I then use these two groups directly to calculate the hydrogen bonds involving water.

In my understanding, these two approaches should give the same results, because the inner and outer groups contain exactly the same chains in both cases. However, the final results show some differences. For example, for residues 1–41 and 42–81, I obtained the results shown below. The red is Method 1, while the black is Method 2:

I am not sure which method is correct, or what might be causing the difference between the two approaches. I would appreciate any suggestions or explanations from anyone familiar with this type of analysis.

Thank you!

The good news is that both of these approaches calculate the same total number of hydrogen bonds, so the primary calculations for both methods are consistent.

Based on these results, this discrepancy may be caused by mislabeling one of the index groups in either method. For residues 1-41, the difference between the two approaches is exactly 117 hydrogen bonds. Using the Method 1 totals, the average number of hydrogen bonds for an outer chain is ~120 (1434 / 12). This suggests that one of the outer chains was accidentally grouped into the inner layer.

Since this difference is much smaller for the second group of residues (42-82), it’s possible that the index selection accidentally matched a partial chain or assigned the residues to an incorrect chain.

I would recommend double-checking the index files. It might be helpful to compare the raw atom indices generated using both methods to see if there is any mismatch.

1 Like