CaviTracer: seed channel searches per chamber, report chamber links, and name outputs after what they hold - #2270
Conversation
…er_radius, and min_volume defaults to 10 A^3. A channel narrower than the probe that traced it is not something that probe can actually pass, so reporting it by default invited the reading that it was open; the old default of 0.0 filtered nothing and left every such entry in the list. min_volume=10, roughly the van der Waals volume of a water molecule, drops the degenerate stubs the same way.
…to avoid dedup's cost explosion when a cavity is searched from more than one seed: on 1tqn with 26 seeds the pathfinding took 238s against 6s single-seed. The opening centres and radii are grown on append rather than rebuilt per candidate, so a small kept set pays nothing for the change. Measured 238s -> 31s on that load, with bit-for-bit identical output on all three regression cases (linb_h, linb_xray, porin).
…epest tetrahedron - findDeepestTetrahedra gives a cavity one seed, the far end of its geodesic depth field. That is proper site for a compact pocket but problematic for a large or branched ones, especially when a small probe connects the real pockets into as single cavity, so the deepest point can land in one of them only. Every channel is then at least as long as the cavity is deep, and the mouths of the other lobes go unreported, because a route from that far seed is absorbed at a nearer mouth before it arrives. findChambers carves the cleared void a second time with seed_radius (default max(1.4, inner_radius)) and labels the connected pieces that survive; each is a chamber and seeds one search. then setStartingTetrahedraFromChambers picks each chamber's widest tetrahedron among those no shallower than min_depth; The pass is skipped when start_point is given or when seed_radius does not exceed inner_radius, where the carve can prune nothing and would silently move the seed from the deepest to the widest tetrahedron rather than being a no-op.
Per-object files now carry the start point they were traced from, as out_sp0_chl3.pqr, so everything belonging to one void globs together, and a run told only a directory drops the run stem altogether: channels.pqr, cavities.pqr and pores.pqr beside sp0_chl3.pqr, since a placeholder stem repeated on every file names nothing. Pores are written as pores rather than under the channel writer's defaults, which had them colliding with the channels once the stem was gone, and a cavity file states its volume, depth and tetrahedron count in a REMARK, its atoms being fixed-size markers rather than a measurement. selectChannelBySelection reads the new names back through one shared parser, so it finds the residue and parameter rows of any object kind and reports what it found, with file_prefix for the runs whose file names no longer carry one.
The get*ParametersMultipleFrames wrappers wrote every frame under a single name, so the frames either overwrote one another or piled up under colliding keys, and only one frame's parameters could be read back. Each frame now writes its own file, named _frame<i> or _model<i> after the new trajectory argument, matching what the residue wrappers already do so that the parameters and the lining of one frame stay findable under the same prefix.
min_chamber_tetrahedra becomes seed_volume (default 30 A^3), measuring a chamber on the same scale a cavity is measured on so that the two can share a threshold, and max_seeds now keeps the largest chambers rather than the widest. The chamber carve stops at min_depth, so a chamber is a buried site rather than a lobe holding its own mouth, and a cavity where no chamber qualifies is seeded at its widest deep tetrahedron instead of its deepest one, which no longer caps every channel of that cavity at a narrow seed. Cavities are ordered largest first, so cavity 0 is the biggest void wherever it is named. min_volume defaults to None and min_tetrahedra/max_tetrahedra move to the advanced options: what gets reported is deliberately separate from where a search may start.
A deep chamber often has no opening of its own and reaches the solvent through a shallower one, and its route out was dropped as a duplicate of that chamber's shorter channel, leaving the deep site unreported. Such access is now reported as a link: the route from the deep chamber, cut where it joins the shallower one, returned in details['links'], written beside the channels, and measured by getLinkParameters and getLinkResidueNames with their multi-frame wrappers. Every channel and link carries the start point it was traced from, in the file name, in the REMARK and as an attribute, so the objects of one void can be analysed together. The log gains a table of the search sites with their volume, depth and object counts, including the sites that report nothing because every way out of them is narrower than bottleneck.
The bottleneck docstring is reworded, saying the same thing more plainly. Ångström is spelled A with a ring in the messages that used a bare A, and the surface cavity line now reports how many cavities were found alongside how many survived min_depth and what is done with them.
|
@karolamik13 pls have a look. This resolves the 1tqn trouble observed in https://github.com/prody/ProDy-website/blob/master/tutorials/cavitracer_tutorial/singlePDB.rst. We now have roughly the same # of channels but also a few linkers (new channel types, aside from pores and channels) connecting the cavities but not going to the surface via a duplicate path. Importantly, the bottlenecks are now physically sensible, and paths starts from appropriate locations: |
karolamik13
left a comment
karolamik13
left a comment
There was a problem hiding this comment.
Looks good. I tested the code on different levels.
2 participants
Summary
Channel detection previously started every search at a single point — the deepest tetrahedron of a cavity. That is the right site for a compact pocket, but not for a large or branched one, and especially not when a narrow probe fuses several real pockets into one cavity. This PR replaces the single seed with one seed per chamber, adds a way to report how those chambers connect, and makes the defaults, output names and logs consistent with the multi-seed picture.
Changes implemented here:
1) Channels are seeded per chamber, not at the single deepest point. With one seed, every channel of a cavity was at least as long as the cavity is deep, and the mouths of the other lobes went unreported — a route from the far seed is absorbed at a nearer mouth before it arrives. The cleared void is now carved a second time at seed_radius, and each surviving connected piece is a chamber that seeds its own search. Chambers qualify by volume and depth (seed_volume, min_depth) rather than by tetrahedron count, so a chamber is measured on the same scale as a cavity, and max_seeds bounds the run on structures whose chambers are least trustworthy. Cavities are ordered largest first, so cavity 0 is the biggest void wherever it sits.
2) Chamber links are reported. A deep chamber often has no opening of its own and reaches solvent only through a shallower one; its route out then duplicates the shallower chamber's shorter channel and was dropped, leaving the deep site unreported. That access is now reported separately as a link — the route from the deep chamber, cut where it joins the shallower one — so the whole way out reads as link(deep → shallow) + channel(shallow → surface), and the link's bottleneck is the neck that actually governs passage. Links are returned in details['links'], written beside the channels, and measured by new getLinkParameters / getLinkResidueNames with their multi-frame wrappers.
3) Defaults now describe what a probe can actually pass. bottleneck defaults to inner_radius instead of 0: a channel narrower than the probe that traced it is not something that probe can pass, and reporting it by default invited the opposite reading. min_volume and the tetrahedron-count filters stay purely about what is reported, deliberately decoupled from where a search may start, so asking for only the large channels never removes the seeds that find them.
4) Output files are named after what they hold, and by their start point (start site). Per-object files carry the site they were traced from (out_sp0_chl3.pqr), so everything belonging to one void globs together; a run given only a directory drops the run stem entirely (channels.pqr, links.pqr, cavities.pqr, pores.pqr). Pores no longer collide with channels once the stem is gone, cavity files state their volume, depth and tetrahedron count in a REMARK, and selectChannelBySelection reads the new names back through one shared parser for any object kind.
5)Multi-frame parameter wrappers write one file per frame. They previously wrote every frame under a single name, so frames overwrote one another and only one frame's parameters could be read back. Each frame now writes
_frame<i>_ / _model<i>,matching what the residue wrappers already do, so parameters and lining for the same frame stay findable under one prefix.Behaviour changes to be aware of