Skip to content

Fix multi chain RMSD analysis#64

Merged
hannahbaumann merged 66 commits intomainfrom
fix_rmsd_multichain
Feb 16, 2026
Merged

Fix multi chain RMSD analysis#64
hannahbaumann merged 66 commits intomainfrom
fix_rmsd_multichain

Conversation

@hannahbaumann
Copy link
Contributor

@hannahbaumann hannahbaumann commented Dec 18, 2025

Fixes #30

@hannahbaumann hannahbaumann changed the title Fix multi chain RMSD analysis [WIP] Fix multi chain RMSD analysis Dec 18, 2025
@hannahbaumann hannahbaumann linked an issue Dec 19, 2025 that may be closed by this pull request
@hannahbaumann hannahbaumann self-assigned this Dec 19, 2025
@hannahbaumann
Copy link
Contributor Author

Check if this works for triclinic boxes. From Irfan: There might also be code in MDAnalysis that does the math also working for triclinic boxes.

@hannahbaumann hannahbaumann marked this pull request as ready for review January 16, 2026 10:30
@IAlibay IAlibay self-assigned this Feb 2, 2026
@hannahbaumann
Copy link
Contributor Author

pre-commit.ci autofix

@IAlibay IAlibay mentioned this pull request Feb 13, 2026
Copy link
Member

@IAlibay IAlibay left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Couple of small things - main question is about the fractional to cartesian return, should that be on box or box.T?


# this only works for orthogonal boxes
ag.positions -= np.rint(vec / box) * box
frac = np.linalg.solve(box.T, vec) # fractional coordinates
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

This seems to make sense.

However I'm not 100% sure about the shift calculation below, do you have a reference for this / where it was inspired from?

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

This was inspired from the wrap_positions function here:
https://gitlab.com/ase/ase/-/blob/master/ase/geometry/geometry.py

assert_allclose(
output["protein_RMSD"][0][:6],
[0, 30.620948, 31.158894, 1.045068, 30.735975, 30.999849],
[0, 1.089747, 1.006143, 1.045068, 1.476353, 1.332893],
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Were these hand calculated or just a regression test for how the numbers changed?

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

This is just a regression test for how the numbers changed.

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Ok - not here but maybe when you move things to separate classes, it might be good to have at least one test where you hand verify the accuracy. I find that this can be important in sanity checking things.

# this only works for orthogonal boxes
ag.positions -= np.rint(vec / box) * box
frac = np.linalg.solve(box.T, vec) # fractional coordinates
shift = np.dot(np.round(frac), box) # nearest image, then compute shift
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Looking at https://www.homepages.ucl.ac.uk/~rmhajc0/frorth.pdf, should this not be applied to box.T here? (equation 17)

hannahbaumann and others added 7 commits February 16, 2026 13:01
Co-authored-by: Irfan Alibay <IAlibay@users.noreply.github.com>
Co-authored-by: Irfan Alibay <IAlibay@users.noreply.github.com>
Co-authored-by: Irfan Alibay <IAlibay@users.noreply.github.com>
Co-authored-by: Irfan Alibay <IAlibay@users.noreply.github.com>
Co-authored-by: Irfan Alibay <IAlibay@users.noreply.github.com>
@hannahbaumann
Copy link
Contributor Author

I realized that the test for the ClosestImageShift class was on a frame that didn't need shifting, so I updated the test to first shift the ligand arbitrarily into a different image and then to see whether we get it back to the same COM distance as before.

Copy link
Member

@IAlibay IAlibay left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

This looks good to me! Let's see how it goes in practice :)

@hannahbaumann hannahbaumann merged commit d67cfc2 into main Feb 16, 2026
7 checks passed
@hannahbaumann hannahbaumann deleted the fix_rmsd_multichain branch February 16, 2026 15:57
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

None yet

Projects

None yet

Development

Successfully merging this pull request may close these issues.

Analysis does not work properly for proteins with multiple chains

3 participants