In protein simulation, it is sometimes necessary to analyze the angle between different spirals and use it to characterize the different states of the protein. This document discusses how to use GROMACS to calculate the angle between different spirals.

Use gmx bundle to calculate the angle between different spirals 1

If you use AMBER trajectory, you need to use the .parm7 file to convert the .nc trajectory to the trajectory format supported by GROMACS. GROMACS supports trajectory files in the following formats: xtc trr cpt trj gro g96 pdb tng.

Script to convert .nc track to .dcd format (.dcd format does not contain time information, but it is not very important, you can also add it manually at the end)

cpptraj strip.parm7 <<EOF trajin file.nc trajout file.dcd dcd goEOF

In addition, the .tpr file is also required to run the gmx bundle.

Create an index file that defines the axis vector

Use make_ndx (or g_select) to create an index file, and specify the centroid groups that define the two ends of the spiral axis. Run
make_ndx -f conf.gro

Assuming that the residue numbers of the two helices to be processed are 43–52 and 62–73, respectively, if you want to use CA to define the centroid group, you can enter it when prompted
r 43-47 | r 62-66 & a CAr 48-52 | r 69-73 & a CAq

The above two sets of choices respectively define the bottom and top of the two spirals, each using 5 residues. Note that the number of atoms in each group must be the same, and it must be an even number, otherwise you will get an error when you run the following command.

Calculate the tilt angle

Use gmx bundle to calculate the inclination angle of the two spiral axes relative to the average axis
g_bundle -f traj.trr -s topol.tpr -n index.ndx -na 2 -oa

Use -na 2 because there are only two axes, and use -oa to output the PDB file of the axis coordinates for viewing or verifying the calculation results.

When prompted to select a group, just specify the two groups defined in the previous step.

Analyze the output file

The calculation result of the tilt angle is stored in bun_tilt.xvg by default, but the angle between each spiral axis and the average spiral axis is saved, not the required angle between the two spiral axes.

Only when the two spiral axes are in the same plane, the angle between the two spiral axes is equal to the angle between the first spiral axis and the average axis plus the angle between the second spiral axis and the average axis. Angle, that is, the sum of the second and third columns of the file. At this time, you can simply use the awk script to calculate the included angle. An example is as follows:
awk'/^[^#@]/ {print 1,2+$3 }'bun_tilt.xvg> angle.dat

When the two spiral axes to be studied are in different planes, the result obtained by the above method is incorrect. In this case, the two coordinate axis directions saved in axes.pdb can be used to calculate. The sample awk script is as follows:

# Language: bashawk 'BEGIN{rad2deg=45/atan2(1,1)}/t=/ {t=$NF}/ATOM +1/{x1 = $6; y1 = $7; z1 = $8}/ATOM +3 /{x1 -= $6; y1 -= $7; z1 -= $8}/ATOM +4/{x2 = $6; y2 = $7; z2 = $8}/ATOM +6/(x2 -= $6; y2 -= $7 ; z2 -= $8)/TER/ {print t, Acos((x1*x2+y1*y2+z1*z2)/sqrt((x1^2+y1^2+z1^2)*(x2^2+ y2^2+z2^2)))*rad2deg)function Acos(x) {if(x==0) return 2*atan2(1,1) return atan2(sqrt(1/x^2-1), 1 ) }'axes.pdb