too few hits

[fmaumusINRA]

Dear Soeding lab members,
I was trying to use mmseqs for tblastn-like analysis.
With a given aa query db and nt targetdb, I have results of blastall and blast+ giving approx 11k and 14k hits (evalue 0.001), respectively.
Aiming at speeding this process, I wanted to use easy-search.
At first, I was getting only about 500-700 hits with default settings and varying sensitivity.
I then tried out the "--max-seqs" option
--max-seqs INT Maximum results per query sequence allowed to pass the prefilter (affects sensitivity) [0.000]
It seems that by default the max seq number is set to 300. When I switched this to --max-seqs 100000, I am now getting 2.5K hits.
mmseqs easy-search prt.fa genome.fa out /tmp/ -s 7 --max-seqs 100000
That's much better but still far from blast.
Would you have any suggestion to address this discrepancy?

Thanks a lot for your help !

[milot-mirdita]

We don't (/can't) generate similar k-mers for nucleotide searches, so the sensitivity parameter doesn't really affect anything in this case. You could try reducing the k-mer size a bit, that might result in more hits getting passed through the prefilter.

[fmaumusINRA]

Thanks a lot for your response. I tried reducing with "-k 5" which improved a bit and then got an error when trying "-k 4".
But isn't it doing a BLASTP against translated ORFs anyhow?

[milot-mirdita]

Sorry you are right, I misread the initial thread. Maybe try reducing the --min-length of ORFs extracted from the genome the default 30aa (90 nucl) is quite long, that might be one issue. Could you share the two FASTA files, then I can take a look.

[fmaumusINRA]

My query is a reverse transcriptase (RT) domain. For instance

Athila41
LDAGVIYPISDSTWVSPVHCVPKKGGMTVVKNEKDELIPTRTITGHRMCIDYRKLNAASRKDHFPLPFIDQMLERLANHPYYCFLDGYSGFFQIPIHPNDQEKTTFTCPYGTFAYKRMPFGLCNAPATFQRCMTSIFSDLIEEMVEVFMDDFSVYGPSFSSCLLNLGRVLTRCEETNLVLNWEKCHFMVKEGIVLDHKISEKGIEVDKGKVEVMMQLQPPKTVKDIRSFLGHAGFYRRFIKDFA
It is about 240 aa so I wonder if --min-length will help.

My target is a plant genome. The genome is too big to share (1Gb). It is from here https://www.ncbi.nlm.nih.gov/Traces/wgs/?val=APLD01.
Thanks a lot for your attempts and suggestions!!!

[milot-mirdita]

--min-length would affect the genome in this case.

[milot-mirdita]

I get over 37k hits, but I also split the genome db creation from the search:

mmseqs createdb APLD01.1.fsa_nt.gz APLD01.2.fsa_nt.gz ntdb
mmseqs easy-search asd.fa ntdb res.m8 tmp -s 7 --max-seqs 30000000

--min-length is not changing much though.

[fmaumusINRA]

OK, maybe I over simplified the question because I was actually counting merged hits with a minimum size of 520 nt (empirical, a bit less than RT size).
I indeed get 37061 hits with mmseqs, 34268 remain after merging overlaping hits, and 1099 merged hits pass the size filter.
With blast+, I get 8159 hits, 6262 merged, and 2947 merged hits pass the size filter.
So yes, I get more raw hits with mmseqs but they are short :(
Sorry for misleading at first !

[milot-mirdita]

The ORF extraction step might be a bit too naive to deal well with a plant. I guess BLAST is able to extract longer fragments for some reason, but I don't know how.

[fmaumusINRA]

Thanks a lot for your help, see you soon with another question, I am exploring still