I found something that may help you in the following link. It seems that if you convert the trajectory into a gromacs format, you can use g_energy module with option energygrp.

http://gromacs.5086.x6.nabble.com/Calculate-a-single-hydrogen-bond-energy-in-the-polypeptide-molecule-td4402617.html

Equations 9, 10 an 11 in this paper http://arxiv.org/abs/0904.1365 show how we calculate hydrogen bond energy (in the paper 1 energy unit = 2.62 kJ/mol). This energy should be interpreted as the free energy of a hydrogen bond, i.e. taking solvent effects into account, not only the electrostatics of the interaction itself. This has been verified by the fact that implicit water simulations of small peptides give pretty realistic folding temperatures using this hydrogen bond term (and some other energy terms too of course, as described in the paper). If you're only interested in hydrogen bonds it's probably easiest to implement this calculation yourself, otherwise our simulation code is available here http://cbbp.thep.lu.se/activities/profasi/.

A slightly simpler, but probably less accurate approach is that employed by DSSP http://en.wikipedia.org/wiki/DSSP_%28protein%29 .

Thank you Simon, this seems to be what I am looking for. Thank you too, Hiqmet, I hadn't thought of GROMACS (I have no experience with it).

A follow-up question, on your method, Simon. The hydrogen bond term in your model describes backbone-backbone (bb-bb) and backbone-side chain (sc-bb) hydrogen bonds, giving the appropriate values for bond strength. What about side chain-side chain (sc-sc) hydrogen bonds? Should I use the sc-bb strength value for them? Or is the side chain potential term better for these?

I am glad to be of help!

There are two main reasons we didn't include sc-sc hydrogen bonds in the hydrogen bond potential. First, the charged side changes have a net charge and will thus have qualitatively different properties. Second, we have reason to beleive they are weaker than the back bone hydrogen bonds, mainly because the involved atoms have a higher degree of solvation, i.e. in a model with only implicit solvent these should be treated as weaker than bb-bb bonds (which is also why bb-sc is weaker than bb-bb). In the end we decided it was enough to have a general fuzzier term for sc-sc interactions, that models effects from both general Coulomb interactions and sc-sc h-bonds. As I said our goal was to get correct folding thermodynamics.

If you are only interested in counting hydrogen bonds in a given structure, it is probably sufficient to use the same energy scale for all h-bonds, and use a cut-off to decide if they are present (we used 1.03 e.u. as cutoff in that paper, which is similar to the DSSP cutoff in kcal/mol). If you are weighing their importance somehow, using the sc-bb value for sc-sc bonds, or something even weaker is probably more appropriate (extrapolating the difference 3.0 -> 2.3 linearly a value of 1.6 e.u for epsilon_{sc-sc} might be reasonale?).

Thanks again, Simon!

you should use Pisa on EMBL-EBI. This is the online program recommended by RCSB when depositing PDB files to determine the biological relevant oligomer. The program gives a break down of all contact energies between amino acids including hydrogen bonds, salt bridges and hyfrophobic interactions.

@Simon Mitternacht

You are using directed hydrogen bond term in your force field. What do you think about physical sense of this directionality? Are there some none crystallographic arguments?

@Alexey Niktin

I didn't take part in the development of the details of this H-bond term. But I know my colleagues did some experimenting on how strong the directionality should be, and that the measure of success was realistic beta-sheet (and presumably alpha-helix) geometries -- which I guess is a crystallographic argument to a large extent. Of course, making the energy minimum narrower will reduce the conformational entropy of low-energy states and thus affect the free energy landscape, which is what we are always optimizing for.

Quantum mechanical calculations might also support having directional H-bonds, but that is beyond my expertise.

@Simon Mitternacht

Thanks for the reply. Can we talk about this field in more detail?

I can not see whether is it possible to send personal messages in this system.

My e-mail on the page http://www.biomolecular-modeling.com/Contacts.html

Best wishes,

Alexei