Patching not working?

[yuliamostovoy]

Hi, I'm trying to do something relatively simple and can't figure out why it isn't working. Any help would be greatly appreciated! So I have a high-quality local assembly (11Mb) from mouse consisting of three contigs generated from HiFi reads. I want to use the reference chromosome to scaffold the three contigs and fill in the small gaps (like a few kb each) between them.

I ran ragtag scaffold and it successfully turned the three contigs into one scaffold with two N-gaps, as expected, so that's great.

Then I ran ragtag patch and I keep getting back the same scaffold containing the three contigs with two N-gaps, instead of the N-gaps being filled in by the reference chromosome. I'm sure that this area of the reference does not have N-gaps, so it should be working, right? My command looks like this:

ragtag.py patch ragtag_output2/ragtag.scaffold.fasta ../ref/chr7_mm39.fa --aligner minimap2 -o ragtag_output3

Things I've tried include increasing the -i parameter by a lot, adding '--fill-only', setting -q 1, but none of that changed the outcome. I'm a step away from giving up and just doing this manually. Thanks for your help!

[malonge]

Hi there,

This is a beta version of patch, so there may very well be a bug. That said, not all gaps/joins are eligible for patching even when everything is working properly. For example, if the region is too repetitive, RagTag may not find a confident patch.

The way to check this is to rerun with the --debug option. Then let me know how the ragtag.patch.debug.*.paf files look. In particular, ragtag.patch.debug.merged.paf will contain filtered and merged alignments. If any subset of these alignments are useful for patching, they will be in ragtag.patch.debug.useful.paf. Knowing the contents of this file should lead us in the right direction.

[malonge]

Also, I recommend using Nucmer (default) instead of minimap2, unless you have a particular justification for using miniamp2. I find that it generally works better for patching. It may not make a difference for this particular problem, but it usually helps to make the patch boundaries more accurate.

[yuliamostovoy]

Thanks! I reran patch with Nucmer and used the debug option. The 'useful.paf' file is indeed empty. I'm attaching all the ragtag.patch.debug.*.paf files. Based on the UCSC browser, it doesn't look like these two gap regions are especially repetitive, but they do have what seems like a typical assortment of LINE and SINE elements. The HiFi coverage just dropped out in these regions, which made the contigs break up.
debug.paf.tar.gz

[malonge]

Great - I will take a look at the alignments and see what is going on. Could you also provide the coordinates for the gaps? Maybe the AGP file from the scaffolding step?

[yuliamostovoy]

Sure, it's attached. Thanks again!
ragtag.scaffold.agp.gz

[malonge]

Hi there,

Unfortunately I am having a hard time viewing some of the PAF files - it looks like some of them are binary. Could you try sending them again?

[yuliamostovoy]

Oh, weird, it looks like one of them got corrupted when I downloaded it from the server. Hopefully this works better.
debug.paf.tar.gz

[malonge]

I took a look and it is not clear to me why there are not more informative alignments after merging. Would you be willing to email me (malonge11@gmail.com) the assembly so I could run things on my end and debug? And if you could point me to where I could download the reference genome that would also be helpful.

Thanks

[yuliamostovoy]

Sure, will do. Thanks a lot!

[malonge]

Hi there,

So I found a few interesting things causing problems. There were two main issues:

1. the target is well-contained within the query

Currently, ragtag assumes that the query and target sequences are roughly homologous. Because the target sequence is a local assembly representing a portion of chromosome 7, this assumption was broken, and that had two consequences:

a. The unique anchor filtering (UAF) was broken. UAF finds the best match for each query position, kind of like delta-filter -q. Accordingly, things may be repetitive with respect to the target. So a bunch of repeats from all over the query (chr7) were non-uniquely mapping to the target (the local assembly), but they weren't removed because they were unique with respect to the query.

b. Because there is so much query sequence up and downstream of the target sequence (everything flanking the local assembly), this exceeded the -i parameter. So much overhanging sequence can indicate an error when patching homologous sequences.

2. there is an insertion in the query

I found a ~150k insertion in the local assembly that broke one of the alignments (its more like ~50kb of the reference was replaced with ~200kb of novel sequence, or at least that is how it appears to me ).

The solution to the first problem is to cut out the homologous portion of chr7 and use that as the query instead of using the entire chromosome. The solution to the second problem is to set -d 200000. This worked for me on my end.

With regard to (1), I will make a note to address this in the next version. At the least, I will make this caveat clear in the documentation. But I may also add options for this sort of scenario. I'd like to leave the issue open until then.

[yuliamostovoy]

Thank you for looking into this!

[malonge]

You're welcome. I want to make one correction upon reflection, more as a note to myself: (1b) is actually not correct. This was only a problem because of the insertion mentioned in (2).