Re: [Seqan-dev] about extendSeed of Seqan

["Kehr, Birte" <Birte.Kehr@fu-berlin.de>]

Hi Beifang,

 

right now, we do not plan to implement traceback for extendSeed, though we might in the future.

 

If you do not fear to modify code, the most efficient way to get the *number* of matching positions is without traceback at all, i.e. also without global alignment:

While extending the seed, parts of a DP matrix is filled with scores. In addition to that score, you could store the number of matching positions for each matrix entry. In the code this would mean to replace occurrences of TScoreValue by Pair<TScoreValue, unsigned> (for the antidiagonals and tmp variables), and to increase the unsigned value if the max for one matrix entry is a match.

 

 

The current seeds-module uses end positions that are not consistent with the rest of SeqAn. For the rest of SeqAn, the position behind the last position is the end position, i.e. end - start = length. Also, diagonals are counted from left/bottom to right/top and not the other way around. Unfortunately, the seeds-module does it differently. That is why you have to take the negative of the diagonals.

 

This is also one of the reasons why the seeds2-module is already there (fixing these inconsistencies). As soon as all functionality is implemented and tested in the seeds2-module we will replace the seeds-module by it. The implementation of seed extension in the seeds2-module is already stable.

 

-Birte

 

 

From: Beifang Niu [mailto:neilniu.cn@gmail.com]
Sent: Donnerstag, 27. Oktober 2011 16:23
To: Kehr, Birte
Subject: Re: [Seqan-dev] about extendSeed of Seqan

 

Hi Birte,

 

Do you have plan to add traceback implement in extendSeed function to get match numbers or seed identity? maybe, it will consume less CPU time than the process extendSeed+globalAlignment.

 

I am confusing that why the minus dignonals "-leftDiagonal(seed) - 2, -rightDiagonal(seed) + 2"  were used in tutorial example. Can I use the same style in my real seed extension?

 

I used the simple score "Score<TScore, Simple> scoreMatrix(1, -2, -1, -6);" and I think    this kind of score is not helpful to get match numbers based on the return score from extendSeed.

 

thanks,

Beifang.

 

 

On Thu, Oct 27, 2011 at 3:41 PM, Kehr, Birte <Birte.Kehr@fu-berlin.de> wrote:

Yes, in the case of gapped X-drop extension you have to do globalAlignment. The function extendSeed does not do the traceback and does not determine the number of matching positions.
But it does compute maximal and minimal diagonals, such that you can band the global alignment (see the Tutorial example).

What kind of score are you using? Would the score of the extensions help you?

-Birte

________________________________________
From: Beifang Niu [neilniu.cn@gmail.com]

Sent: Thursday, October 27, 2011 11:15 PM
To: Kehr, Birte
Subject: Re: [Seqan-dev] about extendSeed of Seqan

Hi Birte,


Unfortunately, I have to use gapped X-drop extension.
Do I have to do globalAlignment to get the matched number?  I just need the matched number after gapped seed extension and there will be increase in computation time if I have to do globalAlignment for getting matched number.
Any ideas?

thanks,
Beifang.

On Thu, Oct 27, 2011 at 2:08 PM, Kehr, Birte <Birte.Kehr@fu-berlin.de<mailto:Birte.Kehr@fu-berlin.de>> wrote:
Hi Beifang,

the function extendSeed does not return the number of matched sequence positions.

I assume you have used ungapped X-drop extension? Then you can count matching positions by simply iterating over the infixes:

typedef typename Infix<TSeq>::Type TInfix;
TInfix infix1 = infix(seq1, leftPosition(seed, 0), rightPosition(seed, 0)+1);
TInfix infix2 = infix(seq2, leftPosition(seed, 1), rightPosition(seed, 1)+1);

unsigned count = 0;
for(int i = 0; i < length(seed); ++i)
{
  if (value(infix1, i) == value(infix2, i))
      ++count;
}

-Birte

________________________________________

From: Beifang Niu [neilniu.cn@gmail.com<mailto:neilniu.cn@gmail.com>]

Sent: Thursday, October 27, 2011 9:19 PM
To: Kehr, Birte
Subject: Re: seqan-dev Digest, Vol 25, Issue 6

Hi Birte,

Thank you for your prompt response but I didn't receive your reply from seqan development mail list.
 I did see the example for seed extension in the SeqAn-Tutorial. Now, I have other questions for you.
 How can I get the actual aligned bases number between two extended seeds after running extendSeeds?
 for example, sequences: ACGTAGTTT  and ACGTGGTTT , there is one seed GTTT, after the extension of left , I got the extension seeds:     ACGTAGTTT   and ACGTGGTTT ( there is only one mismatch) , the actual aligned bases number is  8.
 I want to get this number but I don;t know how to get it only running extendSeeds.



thank you,
Beifang.

On Tue, Oct 25, 2011 at 3:00 AM, <seqan-dev-request@lists.fu-berlin.de<mailto:seqan-dev-request@lists.fu-berlin.de><mailto:seqan-dev-request@lists.fu-berlin.de<mailto:seqan-dev-request@lists.fu-berlin.de>>> wrote:
Send seqan-dev mailing list submissions to

     seqan-dev@lists.fu-berlin.de<mailto:seqan-dev@lists.fu-berlin.de><mailto:seqan-dev@lists.fu-berlin.de<mailto:seqan-dev@lists.fu-berlin.de>>

To subscribe or unsubscribe via the World Wide Web, visit
     https://lists.fu-berlin.de/listinfo/seqan-dev
or, via email, send a message with subject or body 'help' to

     seqan-dev-request@lists.fu-berlin.de<mailto:seqan-dev-request@lists.fu-berlin.de><mailto:seqan-dev-request@lists.fu-berlin.de<mailto:seqan-dev-request@lists.fu-berlin.de>>

You can reach the person managing the list at

     seqan-dev-owner@lists.fu-berlin.de<mailto:seqan-dev-owner@lists.fu-berlin.de><mailto:seqan-dev-owner@lists.fu-berlin.de<mailto:seqan-dev-owner@lists.fu-berlin.de>>

When replying, please edit your Subject line so it is more specific
than "Re: Contents of seqan-dev digest..."


Today's Topics:

 1. about extendSeed of Seqan (Beifang Niu)
 2. Re: about extendSeed of Seqan (Kehr, Birte)


----------------------------------------------------------------------

Message: 1
Date: Mon, 24 Oct 2011 14:50:33 -0700

From: Beifang Niu <neilniu.cn@gmail.com<mailto:neilniu.cn@gmail.com><mailto:neilniu.cn@gmail.com<mailto:neilniu.cn@gmail.com>>>

Subject: [Seqan-dev] about extendSeed of Seqan

To: seqan-dev@lists.fu-berlin.de<mailto:seqan-dev@lists.fu-berlin.de><mailto:seqan-dev@lists.fu-berlin.de<mailto:seqan-dev@lists.fu-berlin.de>>
Message-ID:
     <CABnPkb9P5nAvhQD_4in6mo+Vtim6uhjWEFgzGMhxfb5XWuLQ2g@mail.gmail.com<mailto:CABnPkb9P5nAvhQD_4in6mo%2BVtim6uhjWEFgzGMhxfb5XWuLQ2g@mail.gmail.com><mailto:CABnPkb9P5nAvhQD_4in6mo%2BVtim6uhjWEFgzGMhxfb5XWuLQ2g@mail.gmail.com<mailto:CABnPkb9P5nAvhQD_4in6mo%252BVtim6uhjWEFgzGMhxfb5XWuLQ2g@mail.gmail.com>>>

Content-Type: text/plain; charset="iso-8859-1"

Hi,

I am trying to use Seqan library to do the MEM (max exact match) extension.
Firstly, I get the MEM of the two genome sequences using MUMMER3 and then I
want to use extendSeed of Seqan to do extension of MEMs.
Is there any examples for extendSeed function of seeds class?

thanks,
Beifang.
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <https://lists.fu-berlin.de/pipermail/seqan-dev/attachments/20111024/bca37918/attachment.htm>

------------------------------

Message: 2
Date: Tue, 25 Oct 2011 00:59:24 +0200

From: "Kehr, Birte" <Birte.Kehr@fu-berlin.de<mailto:Birte.Kehr@fu-berlin.de><mailto:Birte.Kehr@fu-berlin.de<mailto:Birte.Kehr@fu-berlin.de>>>

Subject: Re: [Seqan-dev] about extendSeed of Seqan

To: SeqAn Development <seqan-dev@lists.fu-berlin.de<mailto:seqan-dev@lists.fu-berlin.de><mailto:seqan-dev@lists.fu-berlin.de<mailto:seqan-dev@lists.fu-berlin.de>>>
Message-ID:
     <DAD226CB6878494EABEFD5215AA102015A8DE1CD68@exchange6.fu-berlin.de">DAD226CB6878494EABEFD5215AA102015A8DE1CD68@exchange6.fu-berlin.de<mailto:DAD226CB6878494EABEFD5215AA102015A8DE1CD68@exchange6.fu-berlin.de">DAD226CB6878494EABEFD5215AA102015A8DE1CD68@exchange6.fu-berlin.de><mailto:DAD226CB6878494EABEFD5215AA102015A8DE1CD68@exchange6.fu-berlin.de">DAD226CB6878494EABEFD5215AA102015A8DE1CD68@exchange6.fu-berlin.de<mailto:DAD226CB6878494EABEFD5215AA102015A8DE1CD68@exchange6.fu-berlin.de">DAD226CB6878494EABEFD5215AA102015A8DE1CD68@exchange6.fu-berlin.de>>>

Content-Type: text/plain; charset="us-ascii"

Hi Beifang,

you can find an example for seed extension in the SeqAn-Tutorial at
http://trac.mi.fu-berlin.de/seqan/wiki/Tutorial/Seed-and-Extend#SeedExtensionAndBandedAlignment.

You might also want to consider to use the seeds2 module instead of the seeds module since we plan to replace the seeds module by the seeds2 module. Unfortunately, there is no example on how to use the seeds2 module, yet.

-Birte

From: Beifang Niu [mailto:neilniu.cn@gmail.com<mailto:neilniu.cn@gmail.com><mailto:neilniu.cn@gmail.com<mailto:neilniu.cn@gmail.com>>]

Sent: Montag, 24. Oktober 2011 14:51

To: seqan-dev@lists.fu-berlin.de<mailto:seqan-dev@lists.fu-berlin.de><mailto:seqan-dev@lists.fu-berlin.de<mailto:seqan-dev@lists.fu-berlin.de>>

Subject: [Seqan-dev] about extendSeed of Seqan

Hi,

I am trying to use Seqan library to do the MEM (max exact match) extension.
Firstly, I get the MEM of the two genome sequences using MUMMER3 and then I want to use extendSeed of Seqan to do extension of MEMs.
Is there any examples for extendSeed function of seeds class?

thanks,
Beifang.
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <https://lists.fu-berlin.de/pipermail/seqan-dev/attachments/20111025/a7b787b0/attachment.htm>

------------------------------

_______________________________________________
seqan-dev mailing list

seqan-dev@lists.fu-berlin.de<mailto:seqan-dev@lists.fu-berlin.de><mailto:seqan-dev@lists.fu-berlin.de<mailto:seqan-dev@lists.fu-berlin.de>>

https://lists.fu-berlin.de/listinfo/seqan-dev


End of seqan-dev Digest, Vol 25, Issue 6
****************************************