GROMACS version:2021.6
GROMACS modification: No
i’ m using
echo -e "8\n16\n" | gmx hbond -f rest_3pc_MD${i}_topolmdtest_50dt.xtc -s topolmdtest.tpr -num 50ns_off_rest_3pc_MD${i}_topolmdtest_50dt_gmx_hbond_sc_wat_num.xvg -hbn 50ns_off_rest_3pc_MD${i}_topolmdtest_50dt_gmx_hbond_sc_wat.ndx -hbm 50ns_off_rest_3pc_MD${i}_topolmdtest_50dt_gmx_hbond_sc_wat_matrix.xpm -b 50000
to calculate the hydrogen bonds between protein and water, where 8 and 16 represent SideChain and SOL respectively.
I want to know the average number of hydrogen bonds that each residue’s side chain forms with the water molecule per frame.
First of all, i write a gmx_hbond_ndx_count_dat.py
gmx_hbond_ndx_count_dat_py.txt (2.3 KB)
file to calculate the number of each type of hbonds according to the file 50ns_off_rest_3pc_MD${i}_topolmdtest_50dt_gmx_hbond_sc_wat.ndx and file 50ns_off_rest_3pc_MD${i}_topolmdtest_50dt_gmx_hbond_sc_wat_matrix.xpm. The result file 50ns_off_rest_3pc_MD1_topolmdtest_50dt_gmx_hbond_sc_wat.dat is as follows:
30 31 6839 1
30 31 6883 1
30 31 6887 1
30 31 6907 1
30 31 6979 1
30 31 7031 1
30 31 7047 1
30 31 7131 1
30 31 7171 2
30 31 7179 2
30 31 7247 1
30 31 7275 3
30 31 7331 2
30 31 7335 1
30 31 7355 1
30 31 7387 1
...
224095 224096 187 1
224095 224096 188 1
224095 224096 440 2
224095 224096 441 1
224095 224096 704 3
224095 224096 919 2
224099 224100 73 1
224099 224100 187 4
224099 224100 188 1
224099 224100 540 1
224099 224100 969 1
224099 224100 983 2
first three columns are the last group in .ndx file, which are donors, hydrogen and receptors. the last column is the frequency of hydrogen bonds.
Then I use gmx_hbond_dist_on_res.py
gmx_hbond_dist_on_res_py.txt (2.1 KB)
to compare the atoms forming hydrogen bonds with the input1.txt
input1.txt (44.3 KB)
file and add up the total number of hydrogen bonds formed by each residue (side chain) in the 50ns_off_rest_3pc_MD1_topolmdtest_50dt_gmx_hbond_sc_wat.dat file. The result file 50ns_off_rest_3pc_MD1_topolmdtest_50dt_gmx_hbond_dist_on_res_sc_wat.dat is as follows:
01ACE 0 0.0
02VAL 0 0.0
03SER 13380 1.6722909636295462
04GLY 0 0.0
05SER 12725 1.5904261967254094
06PRO 0 0.0
07SER 13127 1.6406699162604674
08GLY 0 0.0
09PRO 0 0.0
10ARG 28430 3.5533058367704036
11GLY 0 0.0
12ASP 40320 5.039370078740157
13ARG 28307 3.537932758405199
14SER 11530 1.4410698662667167
15GLU 44833 5.6034245719285085
16GLY 0 0.0
17PRO 0 0.0
...
61VAL 0 0.0
62GLN 25551 3.1934758155230596
63GLN 23571 2.9460067491563553
64GLU 39128 4.8903887014123235
65GLU 43764 5.469816272965879
66GLU 45135 5.641169853768279
67VAL 0 0.0
68NME 0 0.0
The second column of the above results shows the total number of hydrogen bonds formed by each residue, and the third column shows the average number of hydrogen bonds formed over the number of frames.
But there’s a problem with the results. The side chain of Glu residues contains a carboxyl group, which theoretically can form up to five hydrogen bonds per frame, but most of the calculated results of Glu form more than five hydrogen bonds. For example, 15Glu forms an average of 5.6; 65Glu forms an average of 5.5.
I want to know why this happened.I’d be grateful if anyone could help with that.