Minimap2: query vs target or target vs query in read overlap model have very different results

Created on 8 Oct 2019  路  3Comments  路  Source: lh3/minimap2

I have a very big reads data from pacbio seq2, and I want to get overlap for all reads. So I split the huge read fasta file, and use minimap2 ava-pb to map each pair of sub files. But I found when I use A as query and B as target will get nothing , but use B as query and A as target will get a good results. Why ??

my code:
minimap2 -t 56 -x ava-pb 031.fa 000.fa > 031_vs_000.paf
this cmd get a good results file

[M::mm_idx_gen::73.558*1.94] collected minimizers
[M::mm_idx_gen::81.744*4.10] sorted minimizers
[M::main::81.745*4.10] loaded/built the index for 250000 target sequence(s)
[M::mm_mapopt_update::90.531*3.80] mid_occ = 357
[M::mm_idx_stat] kmer size: 19; skip: 5; is_hpc: 1; #seq: 250000
[M::mm_idx_stat::95.469*3.65] distinct minimizers: 242745064 (44.26% are singletons); average occurrences: 3.701; average spacing: 4.351
[M::worker_pipeline::153.218*15.64] mapped 32255 sequences
[M::worker_pipeline::187.490*18.34] mapped 31983 sequences
[M::worker_pipeline::220.958*20.28] mapped 33009 sequences
[M::worker_pipeline::252.800*21.56] mapped 33889 sequences
[M::worker_pipeline::286.707*22.49] mapped 33790 sequences
[M::worker_pipeline::319.235*23.42] mapped 34288 sequences
[M::worker_pipeline::345.419*22.96] mapped 35404 sequences
[M::worker_pipeline::355.369*22.34] mapped 15382 sequences
[M::main] Version: 2.17-r941
[M::main] CMD: minimap2 -t 56 -x ava-pb 031.fa 000.fa
[M::main] Real time: 357.013 sec; CPU: 7941.179 sec; Peak RSS: 21.388 GB

minimap2 -t 56 -x ava-pb 000.fa 031.fa > 000_vs_031.paf
this cmd get a empty output paf

[M::mm_idx_gen::70.207*1.95] collected minimizers
[M::mm_idx_gen::78.683*4.33] sorted minimizers
[M::main::78.683*4.33] loaded/built the index for 250000 target sequence(s)
[M::mm_mapopt_update::88.403*3.96] mid_occ = 372
[M::mm_idx_stat] kmer size: 19; skip: 5; is_hpc: 1; #seq: 250000
[M::mm_idx_stat::93.258*3.81] distinct minimizers: 233769422 (45.47% are singletons); average occurrences: 3.633; average spacing: 4.374
[M::worker_pipeline::101.071*5.86] mapped 31360 sequences
[M::worker_pipeline::105.144*7.67] mapped 31463 sequences
[M::worker_pipeline::109.309*9.48] mapped 32001 sequences
[M::worker_pipeline::113.724*10.76] mapped 32064 sequences
[M::worker_pipeline::118.359*12.07] mapped 31908 sequences
[M::worker_pipeline::122.208*13.44] mapped 32472 sequences
[M::worker_pipeline::126.029*14.71] mapped 32316 sequences
[M::worker_pipeline::128.879*15.59] mapped 26416 sequences
[M::main] Version: 2.17-r941
[M::main] CMD: minimap2 -t 56 -x ava-pb 000.fa 031.fa
[M::main] Real time: 130.109 sec; CPU: 2010.615 sec; Peak RSS: 17.669 GB

my subfasta file is very normal
000.fa

Nx      Size    Number
N90     8491    152518
N80     13211   117916
N70     17485   93628
N60     20827   74277
N50     23564   57542
N40     26399   42639
N30     29694   29362
N20     34052   17646
N10     40813   7596
Total Length    3714217773      250000
Maximum Length  110245
Minimum Length  50

031.fa

Nx      Size    Number
N90     9065    153364
N80     14104   119225
N70     18485   95239
N60     21673   75819
N50     24399   58842
N40     27290   43686
N30     30747   30171
N20     34974   18222
N10     41607   7899
Total Length    3908445356      250000
Maximum Length  132603
Minimum Length  50

Most helpful comment

You should set --dual=yes, if query fasta file is different from target fasta file.
In default parameters of -x ava-pb, --dual is set to no to reduce half of alignment work.

All 3 comments

You should set --dual=yes, if query fasta file is different from target fasta file.
In default parameters of -x ava-pb, --dual is set to no to reduce half of alignment work.

@wzboy1984 Thank you very much, I will try it : )

As @wzboy1984 said (Thank you!).

Was this page helpful?
0 / 5 - 0 ratings