Residue rotation and UA translation axes for terminal residues

[ioanaapapa]

Summary

Terminal residue axes are no longer dependent on atom indices. The edge atom (i.e. heavy atom bonded to neighbour residue) is chosen alongside the average position of heavy atoms bonded to the edge atom and average position of all other heavy in the residue as the three points used for determining the axes. The centre of rotation (at res level)/translation (at UA level) and axes are then computed from the three points similar to the case of non-terminal residues.

Changes

Terminal residue treatment was changed

• Residue rotation and UA translation axes for terminal residues no longer use atom indices to emulate the backbone for non-terminal residues.
• The edge atom (E1 - red), defined as the heavy atom bonded to a neighbour residue is identified. All heavy atoms bonded to the edge atom (purple) are identified and their average position (C) is used as the third point for obtaining the centre of rotation. All other heavy atoms in the residue (green) are identified and used as the equivalent to the second edge(E2).

• The centre of rotation (O - blue) is the point where the perpendicular from the C (purple) meets the E1 (red)-E2 (green) vector.
• x-axis is along the centre of rotation (O) - first edge(E1) vector
• y axis is along the centre of rotation (O) - average position of bonded atoms (C) vector
• z-axis is perpendicular to the two

Edge cases have been introduced

• If there are no other heavy atoms in a terminal residue, other than the edge atom and those bonded to the edge atom, the centre of rotation is the edge atom and axes are defined as: x-axis is along the edge-bonded vector, y axis is an arbitrary perpendicular axis and z-axis is perpendicular to the two. (This uses the get_custom_axes function and is the same as case2 in get_bonded_axes for UAs)
• If the two edge atoms in a non-terminal residue are bonded (there is no backbone atom between them), the centre of rotation is located at their average position and axes are defined as: x-axis is along the centre-of-rotation-first edge vector, y-axis is arbitrary perpendicular and z-axis is perpendicular to the two.

Tests for terminal residues and edge cases

• Unit tests for terminal residues have been updated to reflect changes.
• New tests for the edge cases covered have been introduced

Documentation

• Vibrational entropy in theory section of the documentation was modified to reflect changes.

get_residue_axes and get_UA_axes structure change

• New functions get_terminal_axes and get_non_terminal_axes are introduced and called in get_residue_axes and get_UA_axes to improve code flow.

Impact

• This closes issue [Feature]: Axes for terminal residues #399.
• Terminal residue treatment is no longer dependent on atom indexing, which is arbitrary. The same result would be obtained using topologies of the same molecule with different indexes.
• Residue rovibrational results are now slightly higher compared to total entropy results obtained using the index-based backbone. No signficant difference observed for UA transvibrational results. (only treatment of terminal residues differs in table below)

[harryswift01]

Thank you for sending these changes through! I think this is a really good improvement, and it’s great to see the terminal residue handling becoming more robust and independent of atom indexing. I’ve just left a few comments around some small technical and maintainability points. If you’re able to have a look through them and let me know what you think, that would be great. Overall, I think this is looking really good and is very close to being ready to go in!

[harryswift01]

I wonder if it might be worth splitting the different cases into their own helper functions and calling them from the if statements? This function is getting quite large now, and there’s a fair bit of conditional logic to follow. I think separating the cases out could make the overall flow a bit clearer and easier to maintain as the logic grows.

[ioanaapapa]

This is a good idea as it would make the code easier to follow.

[ioanaapapa]

@harryswift01 Thank you, Harry! I will make the changes.