Hi David!
Thanks for clarifying! Ill check this out -
Best,
mat
Am 1/4/11 6:10 PM, schrieb Weese, David:
Hi Mat,
the matePairId associates two reads to a mate-pair. Thus
there are at most 2 reads with the same matePairId. The
alignedReadStore stores matches of reads. Every read can have
multiple matches. Each match has a unique id. Two matches of
reads of the same mate-pair can form a mate-pair match. If so
they have the same pairMatchId which is unique for these 2
matches.
Example:
X1 X2 Y1 Y2 Y3
2 reads A and B have matches at positions Xi and Yj.
Theoretically 2*3 pairs of matches are possible. If only X1-Y1
and X2-Y3 form a mate pair match (due to their insert size)
the alignedReadStore could look as follows:
read, position, pairMatchId
A, X1, 23
A, X2, 42
B, Y1, 23
B, Y3, 42
The pairMatchIds can
be chosen arbitrarily (non-sequential) as long as
there are exactly two matches with the
same pairMatchId.
There is a
function compactPairMatchIds that renumerates
the pairMatchIds sequentially from 0 and returns the
number of pair matches so that you can store
additional information in an extra array using the
pairMatchId.
If you need two
matches of all pair matches you can call:
sortAlignedReads(me.alignedReadStore,
SortPairMatchId());
to resort the
alignedReadStore after the pairMatchId so that matches
of a pair match are adjacent.
Cheers,
David
PS: Does the error
still occur? If yes, please create a ticket with SAM
and source files.
--
David Weese weese@inf.fu-berlin.de
Institut für
Informatik Phone: +49 30 838
75246
Takustraße 9 Algorithmic
Bioinformatics
14195 Berlin Room
021
Am 04.01.2011 um 16:47 schrieb Mat:
Hi Manuel!
Thanks for your reply. I can send you some source later,
but i think its a more simple problem. I figured out
that the pairMatchId of the FragmentStore's
alignedReadStore is in general not the same as the
matePairId from the matePairStore by looking into the
source. What still is missing is: I would like to
know on which contig the two reads of a mate-pair
aligned. So far i got:
- First i iterate over the alignedReadStore to get a
readId
- i get the matePairId via the readStore
- then i get the mate's id via the matePairStore
- now i need a function to search for a readId in
alignedReadStore efficiently - is there something like
that?
The code looks now like this:
...
sortAlignedReads(fragStore.alignedReadStore,
SortReadId());
for (unsigned int n=0;n<nrAlignments;++n){
//get the current readId
TValue
currentReadId=fragStore.alignedReadStore[n].readId;
//get the correct mate-pair id
TValue pair_id =
fragStore.readStore[currentReadId].matePairId;
if(pair_id!=TReadStoreElement::INVALID_ID){
++mateCounter;
std::cout << "read-id: " <<
currentReadId << "n: " << n <<
std::endl;
TValue r1 =
fragStore.matePairStore[pair_id].readId[0];
TValue r2 =
fragStore.matePairStore[pair_id].readId[1];
TValue other = 0;
//now r1 or r2 should equal currentReadId
if(r1==currentReadId){
other = r2;
std::cout << "other is r2:"
<< other << std::endl;
//now get the alignId
TValue alignId;
...how to do that?
So how do i get the alignId(s) of a mate's readId?
Thanks!
best,
mat
Am 1/4/11 4:16 PM, schrieb Holtgrewe, Manuel:
It turns out that the both pair_id's never match and the loop breaks with:
Trying to access an element behind the last one!
<<
This assertion is in the operator[]/value(str, pos) function for the String class. Did you try to run the program within a debugger to see where the problem occurs in your program?
Please provide a complete, compiling problem and a minimal SAM file (should be 2 lines after the header, right?) that allow us to reproduce the problem.
Bests,
Manuel
_______________________________________________
seqan-dev mailing list
seqan-dev@lists.fu-berlin.de
https://lists.fu-berlin.de/listinfo/seqan-dev
<ATT00001..txt>
_______________________________________________
seqan-dev mailing list
seqan-dev@lists.fu-berlin.de
https://lists.fu-berlin.de/listinfo/seqan-dev
|
