1
$\begingroup$

I'm leveraging minimap2 to map select genes from short-read microbial fastq metagenomics zip files.

minimap2 --heap-sort=yes --splice-flank=no metagenome.fastq.gz fastaDB.fas 

What I noticed is the low tolerance of mismatch when I need to maximise mismatch because microbes have loads of genetic diversity. Thus the collection of references sequences are never really a reference, in fact the "microbe reference genome" is a misnomer, because in situ diversity of microbes is huge.

Question How do you relax the mapping stringency criteria for minimap2? If you are familiar with Blast parameters this would help alot here. I know blast but I've never parameterised it.

The options are this are:

Alignment options

-A Matching score [2]

-B INT Mismatching penalty [4]

-O INT1[,INT2] Gap open penalty [4,24].

-s INT Minimal peak DP alignment score to output [40].

--score-N INT Score of a mismatch involving ambiguous bases 1

Anything else that is relevant.

Note, -f is a famous option, not sure it helps here.

I'm going to code up a hyperparameterisation unit-test, but understanding the parameters first would really, really help.

Alternative approaches welcome as well.


Plugging minimap2. I code high performance (hadoop-style, back-tracking, graph, traverse tall buildings in a stride etc ...). The coder that wrote minimap has knowledge beyond, i.e. very fast run times directly from gzip file. I assume they are indexing then chunking the raw zip file then unzipping each 'chunk' JIT (just in time).


Thanks @MaximilianPress great answer. I'll fill in the backstory.

I'm trapping total genetic diversity across Genbank and the metagenomes. My protocol to do both:

  • Blasting out Genbank and then
  • piping the output across metagenomes via in this case minimap2.

Hence speed is key. So I ain't replacing Blast with minimap2 but leveraging minimap2 via the output of blast.

However on reflection of @MaximilianPress's discussion the answer is to generate a gold standard of diversity to assess alternative approaches.

$\endgroup$

2 Answers 2

2
$\begingroup$

First of all, use short read option

The short read preset (-x sr) for minimap2 is heavily recommended for doing such mapping:

minimap2 -x sr

Heng Li on minimap2 for related genomes

Heng Li (minimap2 author) has some suggestions on a somewhat related github issue (though it is now 4 years old):

  1. Minimap2 is not really intended for this purpose. For short reads he recommends stampy, but it is clearly quite old (python2.7) and does not seem to have much recent developer history.
  2. Failing that, he suggests an option set such as follows:

If you want to give minimap2 a try anyway, you can use the following setting: minimap2 -axsr -A2 -B4 -O4,24 -E2,1 If minimap2 is not sensitive enough, reduce -k and -w at the cost of mapping speed. Reducing -w should have relatively lower impact on speed, but reducing -k should be more sensitive.

I would suggest engaging with Heng, he is quite responsive in general and may have more suggestions.

Other aligners for metagenomics

I would suggest that you look for more possibilities for alignment, or possibly not using alignment in the first place. However, this is very sensitive to your ultimate application, which I'm not sure that I understand from your post.

I was pleased to see that others have reported good results using minimap2 for metagenomics, which I was skeptical about, having been out of the field for a bit. Usually in metagenomics there is an expected tradeoff between sensitivity and specificity for the reasons that you state (this paper is old but it uses BLAST, so it is a sort of best case for sensitivity on unrelated genomes at the level of divergence that they test).

There are of course new shotgun metagenomic aligners coming out all the time, here is one recent example that I know nothing about but which may be interesting as a comparison.

If you are really looking for a rich parameter set that allows low-stringency mappings, I would suggest going back to BLAST, though it will be much slower. minimap2 and related algorithms are to some extent dependent on the ability to find short seed matches, so there will necessarily be less than perfect sensitivity in your application.

$\endgroup$
1
  • 1
    $\begingroup$ Many thanks. Cool answer. I've addressed the issues you've below the final "-------" $\endgroup$
    – M__
    Commented Oct 2, 2023 at 3:25
1
$\begingroup$

Just to follow-up on the performance of minimap2: the solution @MaximillianPress described works. Thousands of hits are obtained for comparatively obscure genes (some will be duplicates). Parameterisation gives a 2-log10 increase (astonishing).

The issue is that its:

  • relatively slow - precisely as described by @MP;
  • very, very sensitive to parameterisation, just a single integer changes the output by a couple of thousand hits.

I can see why there was concern about using it - such notable sensitive is an issue - and alternatively its clear why other investigators have had great success with it - hit the "sweet-spot" and its a bonanza. For my purposes it's fine - it's not a quantitation of complete genetic diversity, as long as I capture a reasonable proportion of diversity that is fine. Clearly, they'll be times the "sweet-spot" is missed but thats still okay (famous last words :-) ).

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.