vg haplotypes errors about mc graph top-level chain being a loop

[ASLeonard]

Hi,

I was trying to create a haplotype index for a 77-sample mc graph (built with v3.1) using vg 1.70. The steps were

vg gbwt --gbz-input graph.gbz --r-index graph.ri
vg index--no-nested-distance --dist-name graph.dist graph.gbz

vg haplotypes --haplotype-output graph.hapl --distance-index graph.dist --r-index graph.ri graph.gbz

But the haplotypes command reports an issue with top-level chains being a loop. This error seems to have appeared around v1.68 in c73f635. Going back to vg 1.67 gives the error

error: [vg haplotypes] HaplotypePartitioner::partition_haplotypes(): there are 1020 top-level chains and 1957 weakly connected components; haplotype sampling cannot be used with this graph

I believe asking cactus to make the hapl index does use reference tips to upweight nodes, so I can also check that later if necessary, but otherwise this should almost the same procedure.

Best,
Alex

[jltsiren]

Haplotype sampling as quite strict requirements for the structure of the graph. There must be one top-level chain in the snarl decomposition for each weakly connected component in the graph. And each top-level chain must be linear rather than a loop. Older versions of vg did not check for the latter, and they could produce nonsensical haplotypes from loops. Newer versions check for the presence of loops first, before realizing that the number of top-level chains does not match the number of components.

Minigraph–Cactus is supposed to create graphs like that. At least the default/clipped graph is supposed to be compatible with vg haplotypes. I'm not sure about the full graph.

[cwatt]

I'm also getting this error using a minigraph-cactus graph that was clipped with the default 10kb parameter.

[jltsiren]

@cwatt Which error message are you getting? The one about a loop or the one with the mismatch between the number of top-level chains and weakly connected components?

In the latter case, what are the numbers in the error message? Is the number of top-level chains larger or smaller than the number of weakly connected components?

If the number of top-level chains is larger, this is probably a bug in Minigraph–Cactus. If it's smaller, I don't really know what could be causing it. Maybe if there are very short contigs (at most 1024 bp?) in the reference and no other sequence aligns to them, you could get graph components that consist of a single node. Such components might not show at all in the snarl decomposition, making the number of top-level chains smaller than the number of components. In that case, you may be able to solve the problem by removing such contigs from the reference and then rebuilding the graph.

[cwatt]

@jltsiren Apologies, should have included specifics. The loop error.

Here is the command used:
vg haplotypes -v 2 -t 16 -H graph.hapl graph.gbz
Here is the error:
error: [vg haplotypes] HaplotypePartitioner::partition_haplotypes(): top-level chain 49 is a loop; haplotype sampling cannot be used with this graph

I also used the --no-nested-distance option when creating the distance index and the graph was constructed using minigraph-cactus v2.9.9 using a 10kb clipping threshold. It's also worth mentioning this is a maize genome graph.

[jltsiren]

@cwatt It sounds like the issue is the graph itself. If maize genome is full of repetitive elements and transposons and the number of top-level chains / graph components is much higher than the number of chromosomes, it's plausible that you have at least one component that is essentially a large cycle.

I'm not too familiar with the tools used to debug graph construction issues. @glennhickey might be able to help.

[cwatt]

@jltsiren You're right that the number of top-level chains is much higher than the number of chromosomes. Maize genome graphs are very fragmented due to a high proportion of repetitive sequence.

@glennhickey Is there a way to detect and remove cyclic components in a graph? Could cycle removal be worth including as a QC step prior to downstream mapping for complex pangenomes like we encounter in the plant sciences?