vg map issue using pacbio reads

[shilpagarg]

Hi, I am aligning pacbio reads using vg map.

vg map -x true_yeast_graph_X.xg -g true_yeast_graph_X.gcsa -f yeast_pacbio.fastq > yeast_pacbio.gam

and I get

terminate called after throwing an instance of 'std::runtime_error'
what(): Cannot simplify Mapping with no position; must update position when removing leading deletions
Aborted

[ekg]

This may be related to the previous problem that Adam reported, where we have a failure during alignment simplification.

However, it isn't possible to know because this report does not include enough data to reproduce the failure. Would you please upload the input graph, the failing alignment, and write the commands which you used for indexing?

As an aside, do we know that the pacbio alignments are correct? What about simulating some reads with the right error rates and lengths and aligning them with given to vg map --compare? The output is a table that describes how accurate the reads are.

[ekg]

This is vg map --compare:

% vg map --help 2>&1 | grep compare
    -w, --compare         consider GAM input (-G) as truth, table of name, overlap with truth, identity, score, mapqual

If vg map is aligning a GAM that was simulated, then it measures how much of the realignment matches the original ground truth.

Using this with vg sim is critical to understand if we are getting correct alignments. Also, running this interactively will help show if the alignments are taking a reasonable amount of time. It may be that we need to set different parameters for long reads. We may, as long reads haven't been tested on the current version of vg's aligner.

[shilpagarg]

Here is the data:
https://transfer.sh/KWuYL/true-yeast-graph-x.vg
https://transfer.sh/m5Tiv/true-yeast-graph-x.xg

https://transfer.sh/r9Ovz/yeast-pacbio.fastq

[shilpagarg]

here is the back-trace:

#0  0x00007ffff5b14035 in __GI_raise (sig=6) at ../nptl/sysdeps/unix/sysv/linux/raise.c:64
#1  0x00007ffff5b1779b in __GI_abort () at abort.c:91
#2  0x00007ffff662421d in __gnu_cxx::__verbose_terminate_handler() () at /lustre/scratch110/sanger/gt1/build/gcc/gcc-4.9.2-src/libstdc++-v3/libsupc++/vterminate.cc:95
#3  0x00007ffff6622286 in __cxxabiv1::__terminate(void (*)()) () at /lustre/scratch110/sanger/gt1/build/gcc/gcc-4.9.2-src/libstdc++-v3/libsupc++/eh_terminate.cc:47
#4  0x00007ffff6621379 in __cxa_call_terminate () at /lustre/scratch110/sanger/gt1/build/gcc/gcc-4.9.2-src/libstdc++-v3/libsupc++/eh_call.cc:54
#5  0x00007ffff6621adb in __gxx_personality_v0 () at /lustre/scratch110/sanger/gt1/build/gcc/gcc-4.9.2-src/libstdc++-v3/libsupc++/eh_personality.cc:670
#6  0x00007ffff5eabf53 in _Unwind_RaiseException_Phase2 () from /software/hpag/gcc/4.9.2/lib64/libgcc_s.so.1
#7  0x00007ffff5eac477 in _Unwind_Resume () from /software/hpag/gcc/4.9.2/lib64/libgcc_s.so.1
#8  0x000000000075bc5a in vg::Mapper::align_banded(vg::Alignment const&, int, int, int, int) () at /software/hpag/gcc/4.9.2/include/c++/4.9.2/bits/stl_list.h:387
#9  0x000000000075a996 in vg::Mapper::align_multi_internal(bool, vg::Alignment const&, int, int, int, int, double&, int, std::vector<vg::MaximalExactMatch, std::allocator<vg::MaximalExactMatch> >*) () at src/mapper.cpp:2597
#10 0x000000000075ae8a in vg::Mapper::align_multi(vg::Alignment const&, int, int, int, int) () at src/mapper.cpp:2564
#11 0x000000000049bac9 in std::_Function_handler<void (vg::Alignment&), main_map(int, char**)::{lambda(vg::Alignment&)#5}>::_M_invoke(std::_Any_data const&, vg::Alignment&) ()
    at src/main.cpp:6111
#12 0x000000000078dd4b in vg::fastq_unpaired_for_each_parallel(std::string&, std::function<void (vg::Alignment&)>) [clone ._omp_fn.1] ()
    at /software/hpag/gcc/4.9.2/include/c++/4.9.2/functional:2439
#13 0x00007ffff60bf806 in gomp_thread_start () at /lustre/scratch110/sanger/gt1/build/gcc/gcc-4.9.2-src/libgomp/team.c:117
#14 0x00007ffff7762e9a in start_thread (arg=0x7fffeb70f700) at pthread_create.c:308
#15 0x00007ffff5bd136d in clone () at ../sysdeps/unix/sysv/linux/x86_64/clone.S:112
#16 0x0000000000000000 in ?? ()

[ekg]

Awesome, thanks for this. I'll check it out as soon as possible, hopefully tomorrow.

…

On Mon, Feb 13, 2017, 7:48 PM shilpagarg ***@***.***> wrote: here is the back-trace: #0 0x00007ffff5b14035 in __GI_raise (sig=6) at ../nptl/sysdeps/unix/sysv/linux/raise.c:64 #1 <#1> 0x00007ffff5b1779b in __GI_abort () at abort.c:91 #2 <#2> 0x00007ffff662421d in __gnu_cxx::__verbose_terminate_handler() () at /lustre/scratch110/sanger/gt1/build/gcc/gcc-4.9.2-src/libstdc++-v3/libsupc++/vterminate.cc:95 #3 <#3> 0x00007ffff6622286 in __cxxabiv1::__terminate(void ( *)()) () at /lustre/scratch110/sanger/gt1/build/gcc/gcc-4.9.2-src/libstdc++-v3/libsupc++/eh_terminate.cc:47 #4 <#4> 0x00007ffff6621379 in __cxa_call_terminate () at /lustre/scratch110/sanger/gt1/build/gcc/gcc-4.9.2-src/libstdc++-v3/libsupc++/eh_call.cc:54 #5 <#5> 0x00007ffff6621adb in __gxx_personality_v0 () at /lustre/scratch110/sanger/gt1/build/gcc/gcc-4.9.2-src/libstdc++-v3/libsupc++/eh_personality.cc:670 #6 <#6> 0x00007ffff5eabf53 in _Unwind_RaiseException_Phase2 () from /software/hpag/gcc/4.9.2/lib64/libgcc_s.so.1 #7 <#7> 0x00007ffff5eac477 in _Unwind_Resume () from /software/hpag/gcc/4.9.2/lib64/libgcc_s.so.1 #8 <#8> 0x000000000075bc5a in vg::Mapper::align_banded(vg::Alignment const&, int, int, int, int) () at /software/hpag/gcc/4.9.2/include/c++/4.9.2/bits/stl_list.h:387 #9 <#9> 0x000000000075a996 in vg::Mapper::align_multi_internal(bool, vg::Alignment const&, int, int, int, int, double&, int, std::vector<vg::MaximalExactMatch, std::allocatorvg::MaximalExactMatch >*) () at src/mapper.cpp:2597 #10 <#10> 0x000000000075ae8a in vg::Mapper::align_multi(vg::Alignment const&, int, int, int, int) () at src/mapper.cpp:2564 #11 <#11> 0x000000000049bac9 in std::_Function_handler<void (vg::Alignment&), main_map(int, char**)::{lambda(vg::Alignment&)#5 <#5>}>::_M_invoke(std::_Any_data const&, vg::Alignment&) () at src/main.cpp:6111 #12 <#12> 0x000000000078dd4b in vg::fastq_unpaired_for_each_parallel(std::string&, std::function<void (vg::Alignment&)>) [clone ._omp_fn.1] () at /software/hpag/gcc/4.9.2/include/c++/4.9.2/functional:2439 #13 <#13> 0x00007ffff60bf806 in gomp_thread_start () at /lustre/scratch110/sanger/gt1/build/gcc/gcc-4.9.2-src/libgomp/team.c:117 #14 <#14> 0x00007ffff7762e9a in start_thread (arg=0x7fffeb70f700) at pthread_create.c:308 #15 <#15> 0x00007ffff5bd136d in clone () at ../sysdeps/unix/sysv/linux/x86_64/clone.S:112 #16 <#16> 0x0000000000000000 in ?? () — You are receiving this because you commented. Reply to this email directly, view it on GitHub <#662 (comment)>, or mute the thread <https://github.com/notifications/unsubscribe-auth/AAI4EQJwNpkqN4v6jjkriopalh1HtsRQks5rcKV9gaJpZM4L-jAN> .

[ekg]

@shilpagarg this may not be quite quite enough--- how do I generate the index pair from the vg? Is it just vg index -x true-yeast-graph.xg -g true-yeast-graph.gcsa -k 16 true-yeast-graph.vg?

[shilpagarg]

vg index -x true_yeast_graph_X.xg -g true_yeast_graph_X.gcsa -k 8 -X 2 true_yeast_graph_X.vg

…

On Wed, Feb 15, 2017 at 3:59 PM, Erik Garrison ***@***.***> wrote: @shilpagarg <https://github.com/shilpagarg> this may not be quite quite enough--- how do I generate the index pair from the vg? Is it just vg index -x true-yeast-graph.xg -g true-yeast-graph.gcsa -k 16 true-yeast-graph.vg? — You are receiving this because you were mentioned. Reply to this email directly, view it on GitHub <#662 (comment)>, or mute the thread <https://github.com/notifications/unsubscribe-auth/AHEq52gxVNdKBCwBfqswGbum_v04sR3qks5rcxK_gaJpZM4L-jAN> .

[shilpagarg]

Here is the gcsa of the graph:

https://transfer.sh/kWSSf/true-yeast-graph-x.gcsa