new vg combine, adding new options

[xinzilan]

Changelog Entry

To be copied to the draft changelog by merger:

• vg combine gained options for stitching chunked assemblies and graphs: -f/--connect-fragments (now the default mode), -s/--shared-nodes, -m/--merge, and -u/--fusion.
• vg concat is removed (use updated vg combine instead)

Description

This PR adds four new options to vg combine to support our chunked-assembly stitching workflow.

• -f/--connect-fragments (default mode): Treats paths that share the same sample/locus/haplotype/phase block but differ only by subrange — e.g. CHM13#0#chr22 and CHM13#0#chr22[17475777] — as fragments of one path. It sorts the fragments by start offset, concatenates them in order, adds the missing inter-fragment edges, and rewrites the path name to a single [start-end] subrange covering the union. This runs by default whenever none of -c/-p/-m is given.

• -s/--shared-nodes: Modifies the default (-f) mode so node IDs are not renumbered when merging. Node IDs appearing in both inputs must carry the same sequence; one copy is kept. This is intended for vg chunk output that shares boundary nodes across chunks, where the graph (rather than the assembly) was chunked. Fragments are merged into one path, trimming duplicated boundary steps. Valid only with the default mode; cannot be combined with -c/-p/-m.

• -m/--merge: Overlap-aware stitching for chunked-assembly graphs that share REFERENCE-sense paths (e.g. CHM13#0#chr22). Inputs are bucketed by REFERENCE identity (sample/locus/haplotype), so chunks of different chromosomes (e.g. CHM13#0#chr21 and CHM13#0#chr22) can be passed together — each reference is stitched independently and kept as a separate path in the output. Within a reference, it sorts inputs by REFERENCE start offset, validates that each chunk boundary is a single chain node shared by every path (and equal to the REFERENCE end), trims each right chunk's leading node by the REFERENCE-offset overlap, connects the left end to the trimmed right start, then merges fragments like -f.

• -u/--fusion (requires -m): Instead of adding an edge between the left end node and the right chunk's trimmed start node, fuses them into a single node carrying the concatenated sequence.

[faithokamoto]

Are we ready to remove vg concat, then? #4932 (comment)

[adamnovak]

Yeah I think this will replace vg concat.

[adamnovak]

I found a few problems with this.

1. There aren't any tests. There need to be tests for different edge cases of path handling (paths that overlap between input graphs, paths with gaps between them, paths that use phase blocks and don't really want to be connected, paths that are in pieces in their original input graphs) so that we can ensure that the code actually handles paths in the way we expect and produces the appropriate errors when graphs shouldn't be able to be combined.
2. Because there aren't any tests, I think the code is able to do the wrong things with some of those situations. I think it's able to produce paths where their subranges don't actually reflect the total length of all their visits, for example. It's also claiming to be able to change the sense of paths from HAPLOTYPE to REFERENCE.
3. The comments often substitute a summary of the code for an explanation of the abstraction that the code is supposed to provide to its caller to work with. This includes the user-facing documentation for the options.
4. In general the different modes don't share code where it seems like they should, and don't explain why they can't. It looks like -f mode was implemented to combine graphs and join up paths, and then -m mode was implemented to combine fragments along a reference-sense backbone and join up paths, and besides the path-joining function they share basically nothing, with similar problems solved by main force with completely different code. I would expect to see one main implementation of loading and combining the graphs, which is able to do different things at different points by e.g. checking a condition and deciding whether or not to enforce constraints around a single REFERENCE set of paths. I'm not sure the reference mode can't just be a pass, rather than its own nearly-complete implementation.
5. We're approaching the level of complexity where we no longer want all the code to live in the subcommand file. There probably should be a widget class in src or some new algorithms in src/algorithms to hold the guts of this.

[adamnovak]

Just passing argc and argv to all of these (and implicitly optind) was probably a bad idea for cat_proto_graphs(), and it gets worse and worse the more we use it.

I think we did it so we could use the argc/argv version of get_input_file(), but since we also have versions of that function and its friends that work OK with string filenames and "-" for stdin, what we should really do is something like:

std::vector<std::string> filenames_to_combine;
while (optind < argc) {
        filenames_to_combine.push_back(get_input_file_name(optind, argc, argv));
}

And then we would pass a vector of filenames to all the combine functions and they can use:

for (const std::string& filename : filenames_to_combine) {
    get_input_file(filename, [&](std::istream& stream) {
        // Use the file
    });
}

Or something vaguely similar with numbered indexing if they need to handle the first file specially.

[adamnovak]

This should probably check if the loaded graph happens to implement PathPositionHandleGraph and if so use get_path_length().

PathPositionHandleGraph* pos_graph = dynamic_cast<PathPositionHandleGraph*>(graph);
if (pos_graph != nullptr) {
    return pos_graph->get_path_length(path);
}

(We also might consider eventually promoting that method up to PathHandleGraph since I think some graphs track that for you but don't do full efficient position lookups.)

[adamnovak]

This looks like a summary where an explanation is needed. Some of the other, new functions also seem to have this problem.

The right place to set out how a function works, with an overview of the steps it is going to take, is inside the function body, because how the function does what it does is its own business. The point of a function is to use abstraction to keep all that business inside it as much as possible, so that whoever is using the function doesn't need to take up any of their 7 +/- 2 working memory slots thinking about it, and can thus keep larger problems in their head than they otherwise could.

I usually do that interleaved with the code for each step:

// Load the first graph as the destination accumulator
...
for (...) {
    // For each subsequent graph
    ...
    // Renumber its node IDs out of the way
   ...
}
...

That way it's easy to make sure the summary changed when the code changes. But you could also put an outline of the plan at the top of the function.

If you're putting a comment right above a function to introduce it, that comment should explain the why of the function, rather than the how. (It also ought to be a doc comment with /// or /** ... */ delimiters, so documentation tools will pick it up.) It should say things like under what circumstances one may or may not call this function, and what the function is intended to accomplish, and what the caller can rely on about the result. Something like:

/**
 * Combine graphs while connecting paths.
 *
 * Takes argc and argv, and expects flags to have already been parsed,
 * so optind must be pointing to the filename of the first graph on the command line.
 *
 * Everything will end up in a single graph, which comes out via standard output.
 * Graph nodes will end up renumbered so node IDs do not collide.
 * Paths with identical names will end up concatenated together, with the required edges in the graph.
 * Returns 0 on success.
 */

If instead you put a summary of the function's code, anyone who might need to use the function is left to try to read the summary, figure out why the particular things it describes being done are being done and what their consequences will be, and then figure out if those consequences will guarantee that the things they need are going to be true when the function returns. Also, they could just read the code, which has a good chance of being clearer about what the function does in confusing circumstances than the summary can be.

AI code generators in particular love to substitute summaries for explanations, because they're fundamentally translation systems at heart and find it simple to translate code to English and challenging to use a theory of mind to produce text about why the code is there in the first place. If you want to try and use one for this, you might have to wrestle with it or adjust the project's BOTS.md to encourage it.

[adamnovak]

There should be a comment explaining what the 5 fields mean here.

[adamnovak]

Isn't this second-guessing the distinction that the source graph is making between subrange and phase block? If the source graph has genuine unpositioned phase blocks 1, 2, 3, etc., won't we be mistaking those for start offsets here?

[adamnovak]

This again looks like summary rather than explanation.

[adamnovak]

If this is about what the flag is meant for it shouldn't be a parenthetical in the middle attached to a sentence about something else.

[adamnovak]

This is too much like a summary of the code. Instead, we want to talk here about the important differences between the modes from the point of view of someone who has graph pieces to put together and needs to figure out which mode is right for them. They are going to assume that nodes get trimmed by the correct amount to accomplish the task; they need to know what task each of the flags represents.

[adamnovak]

It would be better if the logic described things in terms of the function's own arguments, and not the CLI flags which ought to live on a completely different layer of the stack of abstractions.

[adamnovak]

Didn't we define another class earlier that holds all the path metadata? Maybe there should be a shared struct-of-path-metadata?

[faithokamoto]

I've taken care of removing vg concat, so that now the two changes will merge at the same time. Shouldn't affect any of the vg combine updates.

[faithokamoto]

I did a quick check through the old issues, and this may be relevant to #2688?

[adamnovak]

Yeah, this should fix #2688.

[faithokamoto]

Ah, this line is going to need to be replaced by vg combine. I defer the decision of how to do that to the people who understand what is happening here.

vg/test/t/06_vg_index.t

Line 149 in d9ea751

vg concat x.vg y.vg z.vg >q.vg

[faithokamoto]

#2677 may also be relevant, if this PR means we could do the workflow Adam said was impossible here: #2677 (comment)