Seeding by tree - #626
Merged
Merged
Seeding by tree#626
Conversation
Each include_samples entry may now be a (strain, match_date) tuple in
addition to a bare strain ID. A non-None match_date matches the seed in
on that date instead of its actual date; it may be before the actual
date, in which case the seed node retains its actual date via a negative
("in the future") node time relative to the match-day time-zero. This
lets widely-diverged seeds (e.g. Omicron BA.1/BA.2/BA.3) be matched in on
the same early day so a rudimentary ancestry tree forms before other
tree building, instead of matching to derived samples with many
reversions.
Also detach future (negative-time) nodes when setting up match_tsinfer,
so that other samples cannot copy from a node that lies in the future
relative to the current time-zero.
A seed strain in include_samples that is not present in the dataset now raises a ValueError up-front in extend() instead of being silently warned-and-skipped, so typos and bad configs fail fast. The no-alignment corner case is handled the same as for any other sample (no seed-specific path). detach_future_nodes now fully detaches negative-time nodes: in addition to removing incident edges, it drops mutations over those nodes and clears their NODE_IS_SAMPLE flag so the matcher cannot treat them as ancestors, recomputing mutation parents afterwards. Node IDs and times are preserved. Rewrite the detach_future_nodes tests with hand-drawn topologies (including recombination) following the conventions in test_tree_ops.py.
include_samples changes from a list of strains (optionally with a
match_date) to a list of groups, each a list of strain IDs inserted into
the ARG together as a single local tree.
Rather than matching each seed individually, infer a tree over the
group's haplotypes, match the haplotype of the group's inferred ancestor
against the ARG with recombination disallowed, and attach the whole group
at that single placement. Individual seeds from a major saltation are so
far diverged from the contemporaneous ARG that matching them one at a
time gives either spurious recombination or a match to some derived
sample with a pile of reversions; the inferred ancestor is much closer.
A group is inserted on the minimum date over its members; later-dated
members keep their real dates and so get negative ("future") node times,
using the same mechanism the old match_date override relied on.
Seeds no longer go through the MatchDb. Their matches are constructed
directly by forcing them down the ancestor's path, which removes the
hmm_match.cost = 0.5 hack and the risk of a group being split apart by
add_matching_results' (path, immediate_reversions) grouping key.
Main changes:
- tree_ops.infer_binary_topology(ts) now returns the (pi, tau) oriented
forest instead of mutating a table collection, and returns None for
fewer than two samples. tree_ops.infer_binary(ts, topology=None)
accepts a pre-inferred topology, so a group's topology can be inferred
once against the reference and reused when its mutations are re-mapped
against the matched path.
- New inference helpers: node_haplotypes, reference_haplotype,
path_haplotype, flat_group_ts and force_match.
- normalise_include_samples is replaced by check_include_samples and
seed_group_dates.
- add_matching_results is split into the MatchDb query/grouping and
add_sample_groups, which add_seed_groups also uses.
- create_mask_table is now called unconditionally in _extend, since a day
can consist solely of seeds and the retrospective query must not run
against a stale (or absent) used_samples table.
Backwards compatibility with the old include_samples format is
deliberately not preserved.
create_mask_table records the samples already in the ARG so that the retrospective query skips them, but it was only called when there were samples to match. On a day with no new samples the query therefore ran against the previous day's table, and could add a second copy of a sample that had just been inserted. The mask table is now rebuilt unconditionally in _extend. Update the retro group tests, which were pinning the old behaviour: the 2020-02-15 run has no samples of its own, so ERR4204459 and ERR4206593 were being re-added despite already being samples in the 2020-02-13 base ARG. Assert that samples_strain has no duplicates so this can't regress silently.
Codecov Report✅ All modified and coverable lines are covered by tests. Additional details and impacted files@@ Coverage Diff @@
## main #626 +/- ##
==========================================
+ Coverage 88.13% 88.67% +0.54%
==========================================
Files 12 12
Lines 4163 4319 +156
Branches 596 619 +23
==========================================
+ Hits 3669 3830 +161
+ Misses 363 361 -2
+ Partials 131 128 -3
Flags with carried forward coverage won't be shown. Click here to find out more. ☔ View full report in Codecov by Harness. 🚀 New features to boost your workflow:
|
The option is a list of seed groups, not a list of samples, so name it for what it is. Pure rename: check_include_samples becomes check_seed_groups and the error messages follow, but the behaviour is unchanged. Since the CLI splats [extend_parameters] straight into extend(), this renames the TOML config key too. An old config now fails with a TypeError from extend(), which is the same way any other unknown parameter is reported.
A seed group was represented as a bare tuple of strain IDs, with its date and per-day working state carried alongside in separate maps and lists. Give it a class instead, analogous to Sample: strains and date are its specification, and samples, topology and root are the state filled in on the day it is inserted. check_seed_groups now returns SeedGroups, absorbing seed_group_dates, which is removed. The group hash moves to a sample_group_id function shared with SampleGroup, so the two provably agree on the group_id written to node metadata. SeedGroup.sample_group() builds the SampleGroup used to insert the group, which will let add_seed_groups construct it with all fields known rather than patching path and immediate_reversions in afterwards.
Seeds were preprocessed in the same loop as ordinary samples and then diverted into a dict, with a seed_date map deciding which strains the day's metadata query should suppress and which to inject. Every strain in a seed group is now simply excluded from the standard pipeline on every date, and add_seed_groups preprocesses its own samples. This is the same set of samples by construction: a seed never reached the samples list before either. _extend now mentions seeds twice: one set-difference and one call. The per-day state that was carried in parallel lists lives on the SeedGroup, so the pending list of (group, root_sample) tuples goes away. New helpers, both reusing what was there: - prepare_samples wraps preprocess with the metadata decoration and the no-alignment drop, and is used by both pipelines. - infer_seed_group_root is the group ancestor inference lifted out of add_seed_groups, which makes it testable on its own. Seed samples are no longer merged back into the day's samples, so they no longer appear in the daily samples_processed statistics. A forced match is not an HMM match, and since force_match replaced the old cost = 0.5 sentinel a seed matching its parent exactly was being counted as an exact match.
Construct SampleGroup by keyword in add_matching_results, so that topology no longer has to be the trailing field. Let tmp_dataset take a date per sample, so that behaviour depending on samples having different dates can be tested on the 30bp genome rather than through the 21-day fixture. Add unit tests for infer_seed_group_root, which is now separable enough to test directly.
reference.fasta.fai is generated on demand when the bundled reference is read, and vim swap files turn up next to whatever is being edited. Both were showing up as untracked and are easy to commit by accident.
The group's local tree was built by a round trip: infer a topology against the reference, carry it on the group, synthesise a per-sample hmm_match against the matched path so that match_path_ts could rebuild the same tree from HMM mutations, and reuse the stored topology to keep the two inferences consistent. None of that is needed. The stored hmm_match is informational, so every member of a seed group now reports the match that placed the group. Once the tree is no longer derived from hmm_match it can be built straight from the samples' haplotypes, and each inference is a single infer_binary call: - infer_binary_topology and infer_binary go back to their previous form; the only remaining tree_ops change on this branch is detach_future_nodes. - flat_group_ts takes samples rather than bare haplotypes, and optionally a group_id, in which case the sample nodes carry the metadata attach_tree needs. It now writes each site's ancestral state from the supplied ancestral haplotype rather than from the ARG reference, which is a no-op for the reference-rooted call and is what makes it usable for the other. - SampleGroup.topology becomes SampleGroup.flat_ts, the star tree to infer from. MatchDb groups leave it None and build it from their HMM matches. - infer_seed_group_root returns just the ancestor's haplotype. - force_match and path_haplotype are gone. The parent haplotype comes from node_haplotypes, since a seed group attaches to a single node. The group's tree is still inferred twice, once against the reference to find the ancestor to match and once against the matched node's haplotype to build the tree that is attached. That second rooting is not optional: attach_tree hangs the group's root children off the matched node, so a tree built against the reference silently gives the seeds the parent's allele wherever the parent differs from the reference and the group does not. TestSeedGroupAttachment pins this down; it fails if the tree is built against the reference. Also drop the unreachable "haplotype is None" handling: preprocess always assigns a haplotype, and Dataset.alignment and Dataset.metadata share one sample_id_map, so a strain that passes check_seed_groups has an alignment.
Member
Author
|
This seems to be working, but hasn't been tested over the long term, so holding off on merging. One minor TODO is to fix the debug code to be able to handle nodes with negative times. Currently the ArgInfo class crashes on inputs with future nodes. |
jeromekelleher
marked this pull request as ready for review
July 30, 2026 13:03
Seed group ancestors were matched with NO_RECOMBINATION_NUM_MISMATCHES, which drives rho down to its 1e-200 floor and so effectively forbids recombination. That was introduced back when each seed was matched individually, where far diverged seeds gave spurious recombinations. We now match an inferred ancestor, which is much closer to the contemporaneous ARG, so there is no longer a reason to place it under different HMM parameters to everything else. add_seed_groups takes num_mismatches, and _extend passes its own value down. A group whose ancestor matches a recombinant path is now supported rather than asserted against: - path_haplotype comes back, and the group's local tree is inferred against the mosaic of the matched path rather than the first segment's parent. Without it the group re-derives the other parent's alleles below the recombination point. The genotypes are still correct, so it is only the redundant mutations that give it away. - characterise_recombinants runs on the group roots after match_samples, as it does for the day's samples in _extend, and breakpoint_intervals is propagated to the members alongside the shared hmm_match. is_recombinant is derived from that match, so without this add_sample_to_tables would write an empty interval list for a recombinant seed. attach_tree already attaches across the whole path and flags the group's root node as a recombinant, so nothing downstream needed changing. NO_RECOMBINATION_NUM_MISMATCHES stays for rematch_recombinant, which runs a forced no-recombination match deliberately, to compare against the real one. TestSeedGroupRecombinantAttachment covers the new path with a two-member group whose ancestor is a clean two-segment match at num_mismatches=2. No seed placement in the existing tests changes.
Member
Author
|
The last commit adds support for recombinant paths from the seed group, but I've not checked it thoroughly. If we decide to keep this (i.e., if Delta and Omicron match in without recombination) then it'll need some more work. |
Seed group members dated after time-zero get negative ("in the future")
node times. detach_future_nodes only clears NODE_IS_SAMPLE on the local
copy used for matching, so the persisted tree sequence still has these
flagged as samples and ArgInfo raised ValueError from np.bincount. This
happened before the quick=True guard, so no ArgInfo could be built at all.
Count only samples at time >= 0 in num_samples_per_day, which is indexed
by days-ago and so has no bin for a future sample. Also take the maximum
over sample dates for latest_sample, rather than assuming that node ID
order matches date order.
path_haplotype had no direct coverage at all, and node_haplotypes and reference_haplotype were only used as helpers within other tests, where a wrong result would have agreed with a wrong ARG. Cover the IUPAC encoding that distinguishes node_haplotypes from extract_haplotypes, the duplicate node workaround, and the mosaic assembly in path_haplotype: half-open segment intervals, gaps left as MISSING, repeated parents, overlaps and zero-width segments.
A seed inserted on day D gets a node time of exactly zero, so detach_future_nodes never touches it. It is nonetheless not a copying target for that day's ordinary samples, because they are matched against base_ts, the previous day's ARG, and seeds are appended only afterwards. It becomes matchable on the next processed date, exactly as an ordinary sample added on a given day is only copied from on the following day. Cover the same-day case with two identical samples, one of them seeded, along with the next-day complement that shows matching does work once the seed is visible. Also match against a node at time exactly zero, which is the boundary the existing tests skipped over. Correct the extend docstring, which claimed future-dated members become visible once their real date is reached; it is the day after.
Member
Author
|
I've spend some more time on this, and I'm fairly happy it's working well. I'm going to merge - we can make any further fixes we need in follow-ups. |
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment
2 participants
Add this suggestion to a batch that can be applied as a single commit.This suggestion is invalid because no changes were made to the code.Suggestions cannot be applied while the pull request is closed.Suggestions cannot be applied while viewing a subset of changes.Only one suggestion per line can be applied in a batch.Add this suggestion to a batch that can be applied as a single commit.Applying suggestions on deleted lines is not supported.You must change the existing code in this line in order to create a valid suggestion.Outdated suggestions cannot be applied.This suggestion has been applied or marked resolved.Suggestions cannot be applied from pending reviews.Suggestions cannot be applied on multi-line comments.Suggestions cannot be applied while the pull request is queued to merge.Suggestion cannot be applied right now. Please check back later.
Alternative to #625 - stacked on for simplicity