FIX: Handle a non-linear FrameIterator in HydrogenBondAnalysis #5202
Conversation
|
Oops it seems like my auto-formatter went a bit wild - despite still passing Black. Will clean up. |
BradyAJohnston
force-pushed
the
fix-hbond-iterator
branch
from
695a668 to
84e63d2
Compare
Codecov Report✅ All modified and coverable lines are covered by tests. Additional details and impacted files@@ Coverage Diff @@
## develop #5202 +/- ##
========================================
Coverage 92.72% 92.72%
========================================
Files 180 180
Lines 22475 22475
Branches 3190 3190
========================================
+ Hits 20840 20841 +1
Misses 1177 1177
+ Partials 458 457 -1 ☔ View full report in Codecov by Sentry. 🚀 New features to boost your workflow:
|
orbeckst
left a comment
orbeckst
left a comment
There was a problem hiding this comment.
Looks good!
Optionally: Look at the performance, perhaps there's a faster way to do the lookup.
| count_lookup = dict(zip(indices, tmp_counts)) | ||
| return np.array([count_lookup.get(i, 0) for i in range(len(self.frames))]) |
There was a problem hiding this comment.
Looking up each frame looks slow. Perhaps there's some numpy magic (take???) ?
There was a problem hiding this comment.
The only really faster approach I could figure out would be this:
if self.start is None:
counts = np.zeros(len(self.frames), dtype=int)
positions = np.searchsorted(self.frames, indices)
counts[positions] = tmp_counts
return countsBut this assumes the self.frames to be sorted. Would this always be the case, given the FrameIterator could be a non-sorted sequence of frames?
There was a problem hiding this comment.
I am not sure if self.frames is sorted, possibly not when using run(frames=[2, 3, 0, 7, 6]). Maybe do a quick test?
Perhaps one could sort frames and rearrange counts in the same way and then un-sort everything again before returning?
There was a problem hiding this comment.
Okay I just revisited it and the whole function can actually be simplified which I just pushed.
def count_by_time(self):
"""Counts the number of hydrogen bonds per timestep.
Returns
-------
counts : numpy.ndarray
Contains the total number of hydrogen bonds found at each timestep.
Can be used along with :attr:`HydrogenBondAnalysis.times` to plot
the number of hydrogen bonds over time.
"""
hbond_frames = self.results.hbonds[:, 0].astype(int)
frame_unique, frame_counts = np.unique(hbond_frames, return_counts=True)
counts = np.zeros(max(self.frames) + 1, dtype=int)
counts[frame_unique] = frame_counts
return counts[self.frames]There was a problem hiding this comment.
Alternatively could ensure we only make a small-as-necessary results array:
def count_by_time(self):
"""Counts the number of hydrogen bonds per timestep.
Returns
-------
counts : numpy.ndarray
Contains the total number of hydrogen bonds found at each timestep.
Can be used along with :attr:`HydrogenBondAnalysis.times` to plot
the number of hydrogen bonds over time.
"""
hbond_frames = self.results.hbonds[:, 0].astype(int)
frame_unique, frame_counts = np.unique(hbond_frames, return_counts=True)
frame_min = min(self.frames)
frame_max = max(self.frames)
counts = np.zeros(frame_max - frame_min + 1, dtype=int)
counts[frame_unique - frame_min] = frame_counts
return counts[self.frames - frame_min]There was a problem hiding this comment.
ended up updating it for the second approach using min and max
There was a problem hiding this comment.
The last one looks like a good idea. It looks to me that you can do things like
plot(frames, counts)or
zip(frames, counts)because there's a one-to-one correspondence between the two arrays, right?
There was a problem hiding this comment.
Yep that is the case. I believe it was already mostly the case for a linear range but this ensures it is the same for all frame combinations.
There was a problem hiding this comment.
Might be worth looking to see if other analysis / results classes need a similar tweak
| """ | ||
| hbond_frames = self.results.hbonds[:, 0].astype(int) | ||
| frame_unique, frame_counts = np.unique(hbond_frames, return_counts=True) | ||
| frame_min, frame_max = self.frames.min(), self.frames.max() |
There was a problem hiding this comment.
good, you used np min/max — they are much faster than applying python min/max to arrays
2 participants
Fixes #5200
Changes made in this Pull Request:
count_by_time()to handle returning values on a non-linear FrameIteratorLLM / AI generated code disclosure
LLMs or other AI-powered tools (beyond simple IDE use cases) were used in this contribution: no
PR Checklist
package/CHANGELOGfile updated?package/AUTHORS? (If it is not, add it!)Developers Certificate of Origin
I certify that I can submit this code contribution as described in the Developer Certificate of Origin, under the MDAnalysis LICENSE.
📚 Documentation preview 📚: https://mdanalysis--5202.org.readthedocs.build/en/5202/