amplab / snap

Scalable Nucleotide Alignment Program -- a fast and accurate read aligner for high-throughput sequencing data
https://www.microsoft.com/en-us/research/project/snap/
Apache License 2.0
288 stars 66 forks source link

Read alignment duplication in output #48

Closed cboursnell closed 4 years ago

cboursnell commented 9 years ago

Snap was run with this command:

snap paired index left.fq right.fq -o output.bam -s 0 1000 -H 300000 -h 2000 -d 30 -t 12 -b -M -D 5 -om 5 -omax 10

This is one such example of the alignment duplication I am getting:

FCC2HFRACXX:7:2313:8305:100035#TGNCCAAT 99  Sb02g028080.1   1553    70  86M2I4M2I6M =   1608    155 CAGTAACGCTATTCGTCGCAGTGTTTTGTTTTTAAAAATGTGTAAATAACTGTGCACTGTATATATGTATGTATGTATGTATGTATGTATATATGCATAT    [^^cc`QaWQ[ceXb`aeabe^eb`^f]Yaeehd_aafhh_c^ae[Sa`edddhhhfdd_cddhddd^`^V^bdcdacc`aacZaacZ]aaaacaac`_a    RG:Z:FASTQ  PL:Z:Illumina   PU:Z:pu LB:Z:lb SM:Z:sm PG:Z:SNAP   NM:i:4
FCC2HFRACXX:7:2313:8305:100035#TGNCCAAT 147 Sb02g028080.1   1608    70  31M2I4M2I59M1I1M    =   1553    -155    ACTGTATATATGTATGTATGTATGTATGTATGTATATATGCATATTAAATCAACAGATAGGAGCCATTCCGTAGCCTCTAAGTGGAATCGTTAGGTACGT    Tb`b]]b_Zb]UZZZZcbabb^degdegbb_S]d_b\\S__Z\S\M[bbXaXIa^XXXXea_ggcV`a_Ycfd\gbhgfgd_^d_fdgeac`caca^^_Z    RG:Z:FASTQ  PL:Z:Illumina   PU:Z:pu LB:Z:lb SM:Z:sm PG:Z:SNAP   NM:i:5
FCC2HFRACXX:7:2313:8305:100035#TGNCCAAT 355 Sb02g028080.1   1553    0   86M2I4M2I6M =   1608    155 CAGTAACGCTATTCGTCGCAGTGTTTTGTTTTTAAAAATGTGTAAATAACTGTGCACTGTATATATGTATGTATGTATGTATGTATGTATATATGCATAT    [^^cc`QaWQ[ceXb`aeabe^eb`^f]Yaeehd_aafhh_c^ae[Sa`edddhhhfdd_cddhddd^`^V^bdcdacc`aacZaacZ]aaaacaac`_a    RG:Z:FASTQ  PL:Z:Illumina   PU:Z:pu LB:Z:lb SM:Z:sm PG:Z:SNAP   NM:i:4
FCC2HFRACXX:7:2313:8305:100035#TGNCCAAT 403 Sb02g028080.1   1608    0   31M2I4M2I59M1I1M    =   1553    -155    ACTGTATATATGTATGTATGTATGTATGTATGTATATATGCATATTAAATCAACAGATAGGAGCCATTCCGTAGCCTCTAAGTGGAATCGTTAGGTACGT    Tb`b]]b_Zb]UZZZZcbabb^degdegbb_S]d_b\\S__Z\S\M[bbXaXIa^XXXXea_ggcV`a_Ycfd\gbhgfgd_^d_fdgeac`caca^^_Z    RG:Z:FASTQ  PL:Z:Illumina   PU:Z:pu LB:Z:lb SM:Z:sm PG:Z:SNAP   NM:i:5
FCC2HFRACXX:7:2313:8305:100035#TGNCCAAT 355 Sb02g028080.1   1553    0   86M2I4M2I6M =   1608    155 CAGTAACGCTATTCGTCGCAGTGTTTTGTTTTTAAAAATGTGTAAATAACTGTGCACTGTATATATGTATGTATGTATGTATGTATGTATATATGCATAT    [^^cc`QaWQ[ceXb`aeabe^eb`^f]Yaeehd_aafhh_c^ae[Sa`edddhhhfdd_cddhddd^`^V^bdcdacc`aacZaacZ]aaaacaac`_a    RG:Z:FASTQ  PL:Z:Illumina   PU:Z:pu LB:Z:lb SM:Z:sm PG:Z:SNAP   NM:i:4
FCC2HFRACXX:7:2313:8305:100035#TGNCCAAT 403 Sb02g028080.1   1608    0   31M2I4M2I59M1I1M    =   1553    -155    ACTGTATATATGTATGTATGTATGTATGTATGTATATATGCATATTAAATCAACAGATAGGAGCCATTCCGTAGCCTCTAAGTGGAATCGTTAGGTACGT    Tb`b]]b_Zb]UZZZZcbabb^degdegbb_S]d_b\\S__Z\S\M[bbXaXIa^XXXXea_ggcV`a_Ycfd\gbhgfgd_^d_fdgeac`caca^^_Z    RG:Z:FASTQ  PL:Z:Illumina   PU:Z:pu LB:Z:lb SM:Z:sm PG:Z:SNAP   NM:i:5
FCC2HFRACXX:7:2313:8305:100035#TGNCCAAT 355 Sb02g028080.1   1553    0   86M2I4M2I6M =   1608    155 CAGTAACGCTATTCGTCGCAGTGTTTTGTTTTTAAAAATGTGTAAATAACTGTGCACTGTATATATGTATGTATGTATGTATGTATGTATATATGCATAT    [^^cc`QaWQ[ceXb`aeabe^eb`^f]Yaeehd_aafhh_c^ae[Sa`edddhhhfdd_cddhddd^`^V^bdcdacc`aacZaacZ]aaaacaac`_a    RG:Z:FASTQ  PL:Z:Illumina   PU:Z:pu LB:Z:lb SM:Z:sm PG:Z:SNAP   NM:i:4
FCC2HFRACXX:7:2313:8305:100035#TGNCCAAT 403 Sb02g028080.1   1608    0   31M2I4M2I59M1I1M    =   1553    -155    ACTGTATATATGTATGTATGTATGTATGTATGTATATATGCATATTAAATCAACAGATAGGAGCCATTCCGTAGCCTCTAAGTGGAATCGTTAGGTACGT    Tb`b]]b_Zb]UZZZZcbabb^degdegbb_S]d_b\\S__Z\S\M[bbXaXIa^XXXXea_ggcV`a_Ycfd\gbhgfgd_^d_fdgeac`caca^^_Z    RG:Z:FASTQ  PL:Z:Illumina   PU:Z:pu LB:Z:lb SM:Z:sm PG:Z:SNAP   NM:i:5

The same pair of reads has aligned to this transcript 4 times in the same location.

This error doesn't occur when the -om 5 -omax 10 options are omitted

bolosky commented 9 years ago

Can you tell me where I can get the fasta for the index you used? I can probably repro it pretty quickly from there.

It’s clearly related to the multiple alignment finding code, since all but one are marked as secondary alignments (and, as you observed, it doesn’t happen if you don’t ask for multiple alignments).

--Bill

From: Chris Boursnell [mailto:notifications@github.com] Sent: Friday, February 6, 2015 4:18 AM To: amplab/snap Subject: [snap] Read alignment duplication in output (#48)

Snap was run with this command:

snap paired index left.fq right.fq -o output.bam -s 0 1000 -H 300000 -h 2000 -d 30 -t 12 -b -M -D 5 -om 5 -omax 10

This is one such example of the alignment duplication I am getting:

FCC2HFRACXX:7:2313:8305:100035#TGNCCAAT 99 Sb02g028080.1 1553 70 86M2I4M2I6M = 1608 155 CAGTAACGCTATTCGTCGCAGTGTTTTGTTTTTAAAAATGTGTAAATAACTGTGCACTGTATATATGTATGTATGTATGTATGTATGTATATATGCATAT [^^ccQaWQ[ceXbaeabe^eb^f]Yaeehd_aafhh_c^ae[Saedddhhhfdd_cddhddd^^V^bdcdaccaacZaacZ]aaaacaac`_a RG:Z:FASTQ PL:Z:Illumina PU:Z:pu LB:Z:lb SM:Z:sm PG:Z:SNAP NM:i:4

FCC2HFRACXX:7:2313:8305:100035#TGNCCAAT 147 Sb02g028080.1 1608 70 31M2I4M2I59M1I1M = 1553 -155 ACTGTATATATGTATGTATGTATGTATGTATGTATATATGCATATTAAATCAACAGATAGGAGCCATTCCGTAGCCTCTAAGTGGAATCGTTAGGTACGT Tbb]]b_Zb]UZZZZcbabb^degdegbb_S]d_b\\S__Z\S\M[bbXaXIa^XXXXea_ggcVaYcfd\gbhgfgd^d_fdgeac`caca^^_Z RG:Z:FASTQ PL:Z:Illumina PU:Z:pu LB:Z:lb SM:Z:sm PG:Z:SNAP NM:i:5

FCC2HFRACXX:7:2313:8305:100035#TGNCCAAT 355 Sb02g028080.1 1553 0 86M2I4M2I6M = 1608 155 CAGTAACGCTATTCGTCGCAGTGTTTTGTTTTTAAAAATGTGTAAATAACTGTGCACTGTATATATGTATGTATGTATGTATGTATGTATATATGCATAT [^^ccQaWQ[ceXbaeabe^eb^f]Yaeehd_aafhh_c^ae[Saedddhhhfdd_cddhddd^^V^bdcdaccaacZaacZ]aaaacaac`_a RG:Z:FASTQ PL:Z:Illumina PU:Z:pu LB:Z:lb SM:Z:sm PG:Z:SNAP NM:i:4

FCC2HFRACXX:7:2313:8305:100035#TGNCCAAT 403 Sb02g028080.1 1608 0 31M2I4M2I59M1I1M = 1553 -155 ACTGTATATATGTATGTATGTATGTATGTATGTATATATGCATATTAAATCAACAGATAGGAGCCATTCCGTAGCCTCTAAGTGGAATCGTTAGGTACGT Tbb]]b_Zb]UZZZZcbabb^degdegbb_S]d_b\\S__Z\S\M[bbXaXIa^XXXXea_ggcVaYcfd\gbhgfgd^d_fdgeac`caca^^_Z RG:Z:FASTQ PL:Z:Illumina PU:Z:pu LB:Z:lb SM:Z:sm PG:Z:SNAP NM:i:5

FCC2HFRACXX:7:2313:8305:100035#TGNCCAAT 355 Sb02g028080.1 1553 0 86M2I4M2I6M = 1608 155 CAGTAACGCTATTCGTCGCAGTGTTTTGTTTTTAAAAATGTGTAAATAACTGTGCACTGTATATATGTATGTATGTATGTATGTATGTATATATGCATAT [^^ccQaWQ[ceXbaeabe^eb^f]Yaeehd_aafhh_c^ae[Saedddhhhfdd_cddhddd^^V^bdcdaccaacZaacZ]aaaacaac`_a RG:Z:FASTQ PL:Z:Illumina PU:Z:pu LB:Z:lb SM:Z:sm PG:Z:SNAP NM:i:4

FCC2HFRACXX:7:2313:8305:100035#TGNCCAAT 403 Sb02g028080.1 1608 0 31M2I4M2I59M1I1M = 1553 -155 ACTGTATATATGTATGTATGTATGTATGTATGTATATATGCATATTAAATCAACAGATAGGAGCCATTCCGTAGCCTCTAAGTGGAATCGTTAGGTACGT Tbb]]b_Zb]UZZZZcbabb^degdegbb_S]d_b\\S__Z\S\M[bbXaXIa^XXXXea_ggcVaYcfd\gbhgfgd^d_fdgeac`caca^^_Z RG:Z:FASTQ PL:Z:Illumina PU:Z:pu LB:Z:lb SM:Z:sm PG:Z:SNAP NM:i:5

FCC2HFRACXX:7:2313:8305:100035#TGNCCAAT 355 Sb02g028080.1 1553 0 86M2I4M2I6M = 1608 155 CAGTAACGCTATTCGTCGCAGTGTTTTGTTTTTAAAAATGTGTAAATAACTGTGCACTGTATATATGTATGTATGTATGTATGTATGTATATATGCATAT [^^ccQaWQ[ceXbaeabe^eb^f]Yaeehd_aafhh_c^ae[Saedddhhhfdd_cddhddd^^V^bdcdaccaacZaacZ]aaaacaac`_a RG:Z:FASTQ PL:Z:Illumina PU:Z:pu LB:Z:lb SM:Z:sm PG:Z:SNAP NM:i:4

FCC2HFRACXX:7:2313:8305:100035#TGNCCAAT 403 Sb02g028080.1 1608 0 31M2I4M2I59M1I1M = 1553 -155 ACTGTATATATGTATGTATGTATGTATGTATGTATATATGCATATTAAATCAACAGATAGGAGCCATTCCGTAGCCTCTAAGTGGAATCGTTAGGTACGT Tbb]]b_Zb]UZZZZcbabb^degdegbb_S]d_b\\S__Z\S\M[bbXaXIa^XXXXea_ggcVaYcfd\gbhgfgd^d_fdgeac`caca^^_Z RG:Z:FASTQ PL:Z:Illumina PU:Z:pu LB:Z:lb SM:Z:sm PG:Z:SNAP NM:i:5

The same pair of reads has aligned to this transcript 4 times in the same location.

— Reply to this email directly or view it on GitHubhttps://github.com/amplab/snap/issues/48.

blahah commented 9 years ago

@cboursnell is this still happening?