<div dir="ltr"><div dir="ltr">Thank you, David.<div><br></div><div>Since I emailed the mailing list, I may have found a solution independently: please see <a href="https://www.biostars.org/p/361700/#362231">https://www.biostars.org/p/361700/#362231</a>.</div><div><br></div><div>But if I find your solution to be better in some way(s), I'll use it instead.</div><div>Thank you for taking the time to explain your solution in detail, I truly appreciate it.<br></div><div><br></div><div>Sincerely,</div><div>Anand</div><div><br></div></div></div><br><div class="gmail_quote"><div dir="ltr" class="gmail_attr">On Tue, Feb 5, 2019 at 7:01 AM <<a href="mailto:emboss-request@mailman.open-bio.org">emboss-request@mailman.open-bio.org</a>> wrote:<br></div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex">Send EMBOSS mailing list submissions to<br>
        <a href="mailto:emboss@mailman.open-bio.org" target="_blank">emboss@mailman.open-bio.org</a><br>
<br>
To subscribe or unsubscribe via the World Wide Web, visit<br>
        <a href="http://mailman.open-bio.org/mailman/listinfo/emboss" rel="noreferrer" target="_blank">http://mailman.open-bio.org/mailman/listinfo/emboss</a><br>
or, via email, send a message with subject or body 'help' to<br>
        <a href="mailto:emboss-request@mailman.open-bio.org" target="_blank">emboss-request@mailman.open-bio.org</a><br>
<br>
You can reach the person managing the list at<br>
        <a href="mailto:emboss-owner@mailman.open-bio.org" target="_blank">emboss-owner@mailman.open-bio.org</a><br>
<br>
When replying, please edit your Subject line so it is more specific<br>
than "Re: Contents of EMBOSS digest..."<br>
<br>
<br>
Today's Topics:<br>
<br>
   1. Re: run parameters for inexact matching with EMBOSS needle<br>
      (David Mathog)<br>
   2. Re: run parameters for inexact matching with EMBOSS needle<br>
      (David Mathog)<br>
<br>
<br>
----------------------------------------------------------------------<br>
<br>
Message: 1<br>
Date: Mon, 04 Feb 2019 09:57:06 -0800<br>
From: David Mathog <<a href="mailto:mathog@caltech.edu" target="_blank">mathog@caltech.edu</a>><br>
To: <a href="mailto:emboss@mailman.open-bio.org" target="_blank">emboss@mailman.open-bio.org</a><br>
Subject: Re: [EMBOSS] run parameters for inexact matching with EMBOSS<br>
        needle<br>
Message-ID: <<a href="mailto:1af376ce7957ba022462adfcb19765f2@saf.bio.caltech.edu" target="_blank">1af376ce7957ba022462adfcb19765f2@saf.bio.caltech.edu</a>><br>
Content-Type: text/plain; charset=UTF-8; format=flowed<br>
<br>
On 03-Feb-2019 04:18, Anandkumar Surendrarao wrote:<br>
<br>
>    - I want to search with FASTA formatted DNA sequence queries that <br>
> are <=<br>
>    100nt in length,<br>
>    - against FASTA formatted genomic DNA sequence,<br>
>    - for global end-to-end matches  (using Needleman Wunsch type global<br>
>    search, not Smith Waterman type local matching)<br>
>    - but allowing >= 80% identity and >= 80% of query length coverage<br>
<br>
Not sure about needle, except that I expect it would be very slow.<br>
<br>
However, I made a modified version of Lastal that added command line <br>
parameters for identity/query/target cutoff for almost exactly the same <br>
application.  It is available here in the "modified programs" directory.<br>
<br>
   <a href="https://sourceforge.net/projects/lemonade-assemble/" rel="noreferrer" target="_blank">https://sourceforge.net/projects/lemonade-assemble/</a><br>
<br>
There is a Centos 6 64 bit binary there, for any other platform you will <br>
need to recompile. Lastal does not have NW exactly, but it does have a <br>
-T parameter which does something which in many cases is similar.  -T 0 <br>
is a local alignment, and -T 1 "extends to the end" (starting from the <br>
local alignment).  My use case involved Illumina reads, and what would <br>
happen is if there was a mismatch at 3 or 4 bp from the end of that read <br>
in the alignment against the target (in this case mRNA transcripts) -T 0 <br>
would truncate the alignment at the position before the mismatch, so the <br>
read would not be seen as a full length match. With -T 1 it would pick <br>
up the extra matching bases.  Anyway, let's say you have<br>
<br>
   reads.fasta   # bunch of 100bp reads<br>
   genome.fasta  # self explanatory<br>
<br>
CPUS=40<br>
lastdb -P $CPUS -w 10 genome_10 genome.fasta<br>
lastal -P $CPUS -T1 -I80 -Y80 -m250 -O50 -8 1 -9 1 \<br>
   -f BlastTab genome_10 reads.fasta > results.txt 2>errors.log<br>
<br>
<br>
The lastdb command makes a database but uses only every 10th position.  <br>
This will make the database smaller and the search much faster.  These <br>
positions are only used as seeds, so if the expected match is very good <br>
-w10 gives the same results as -w1 (default).  If the match is expected <br>
to be poor use -w 1 and expect much longer run times.  Set CPUS to match <br>
your system.<br>
<br>
The lastal command looks for 80% identity (-I80), 80% coverage on query <br>
(-Y80), allows up to 250 overlapping matches (-m250, see below), allows <br>
offset overlap alignment as<br>
long as the overlap is at least 50bp (-O50), requires overlap alignments <br>
to go right<br>
to the ends (-8 1 -9 1), and -f BlastTab gives a blast tabular format <br>
output file (except when it is on the other strand the query positions <br>
are flipped rather than the target).  Offset overlap is like this:<br>
<br>
     -------------------------<br>
               |||||||||||||||  <- -O50 == this must be >= 50bp<br>
               ------------------------------<br>
<br>
This only concerns alignments which overlap the ends.   The -m parameter <br>
is a standard lastal command line option:<br>
<br>
   -m: maximum initial matches per query position (10)<br>
<br>
What it means is that if there are a lot of alignments in exactly the <br>
same place by default one only retrieves a few of them.  Increasing the <br>
value to 250 means many more will be retrieved.  Conceivably if there is <br>
super high coverage you might want to use an even larger number.  Best <br>
to first try a small subset to determine the appropriate -m value for <br>
your case.<br>
<br>
I don't know if it would be faster if the query/database were the other <br>
way around.<br>
<br>
Regards,<br>
<br>
David Mathog<br>
<a href="mailto:mathog@caltech.edu" target="_blank">mathog@caltech.edu</a><br>
Manager, Sequence Analysis Facility, Biology Division, Caltech<br>
<br>
<br>
------------------------------<br>
<br>
Message: 2<br>
Date: Mon, 04 Feb 2019 10:32:20 -0800<br>
From: David Mathog <<a href="mailto:mathog@caltech.edu" target="_blank">mathog@caltech.edu</a>><br>
To: <a href="mailto:emboss@mailman.open-bio.org" target="_blank">emboss@mailman.open-bio.org</a><br>
Subject: Re: [EMBOSS] run parameters for inexact matching with EMBOSS<br>
        needle<br>
Message-ID: <<a href="mailto:f5553bd8c0b5781ade9bafca9bc40ed2@saf.bio.caltech.edu" target="_blank">f5553bd8c0b5781ade9bafca9bc40ed2@saf.bio.caltech.edu</a>><br>
Content-Type: text/plain; charset=UTF-8; format=flowed<br>
<br>
On 04-Feb-2019 09:57, David Mathog wrote:<br>
> lastal -P $CPUS -T1 -I80 -Y80 -m250 -O50 -8 1 -9 1 \<br>
>   -f BlastTab genome_10 reads.fasta > results.txt 2>errors.log<br>
<br>
I should add that for my case it was -I93 -J100.  The reads were only <br>
76bp and<br>
it was common to have a few bp mismatches on the ends, hence 93% <br>
identity.  -J instead<br>
of -Y because it was transcripts searching a reads database rather than <br>
the other way around.  80% is pretty low if the short sequences should <br>
be present 1:1 in the assembled genome, even for a highly polymorphic <br>
genome where only one haplotype might be present.<br>
<br>
Regards,<br>
<br>
David Mathog<br>
<a href="mailto:mathog@caltech.edu" target="_blank">mathog@caltech.edu</a><br>
Manager, Sequence Analysis Facility, Biology Division, Caltech<br>
<br>
<br>
------------------------------<br>
<br>
Subject: Digest Footer<br>
<br>
_______________________________________________<br>
EMBOSS mailing list<br>
<a href="mailto:EMBOSS@mailman.open-bio.org" target="_blank">EMBOSS@mailman.open-bio.org</a><br>
<a href="http://mailman.open-bio.org/mailman/listinfo/emboss" rel="noreferrer" target="_blank">http://mailman.open-bio.org/mailman/listinfo/emboss</a><br>
<br>
------------------------------<br>
<br>
End of EMBOSS Digest, Vol 121, Issue 2<br>
**************************************<br>
</blockquote></div>