SEQanswers

Go Back   SEQanswers > Bioinformatics > Bioinformatics



Similar Threads
Thread Thread Starter Forum Replies Last Post
BWA-MEM and X0/X1/etc flags mturchin20 Bioinformatics 5 07-15-2015 06:54 AM
bwa mem with paired/unpaired reads mchaisso Bioinformatics 3 09-09-2013 10:43 AM
BWA-MEM: output mapped reads larger than input reads Kennels Bioinformatics 9 09-05-2013 05:27 PM
why the properly paired results from bwa mem is an odd number Pengfei Liu Bioinformatics 2 08-24-2013 06:26 PM
bwa mem segfault; bwa bwasw breaks MarkDuplicates ElenaN Bioinformatics 2 06-30-2013 10:23 PM

Reply
 
Thread Tools
Old 09-13-2013, 04:16 AM   #1
Linnea
Member
 
Location: Uppsala, Sweden

Join Date: Mar 2010
Posts: 23
Exclamation bwa-mem paired reads problem

Hi,

I have aligned several re-sequenced individuals (Illumina HiSeq, paired-end with IS ~450bp) to a draft assembly, and when looking at the alignments in IGV, it seemed like almost all of them had their mate mapped on a different scaffold.

I extracted some read pairs from the bam file and got this weired result:
HWI-ST0866:1861D42ACXX:2:1313:18722:3523 pR1 scaffold37985 44312 60 100M scaffold54822 4963 0 AATAAATCTACCCCAGCCAGCCTAAAGCCACTCCCTCCTGTCCCATCACTTCATGCCTTTTCACAAAGTCCCTTCTCAGCCTTTTTGCAGCAGCAATCAG @BCFDDFFHHHHHJJIIIFIIJIIJGIGDIIIJJJJIJI<BGHIFJIJJIIJIIJJJJHIJJCE;CDHCC?CCEHFEFFFFE3>@CDDDDDB@DBDACCD NM:i:0 AS:i:100 XS:i:0
HWI-ST0866:1861D42ACXX:2:1313:18722:3523 prR2 scaffold37985 44635 60 100M scaffold39009 30358 0 TGAAGGGTTAAAAACACATGGAGCCAGCATGTTCTGATCACTGAGAAACAACACTGTGAATCAGCACAGCAGCAAATGGGAGTGACAGGCACTGGAGATG CDDBACEECDECCBFFHFHHHIJIJIGIJIJJJJJJIHHIJJJJJIJIHDDHGIHIGIIGDJJIJJIIIIGIFJJJJJJJJJJJJJJHHHHGFFFFFCCC NM:i:1 AS:i:95 XS:i:0

The first read in the pair is located on scaffold37985 and says it has its mate on scaffold54822, but when looking at the mate (second line) it's also located on scaffold37985, but reports to have a mate on scaffold39009!
To me this seems so wrong, why can't they just form a proper pair, when they obviously have perfect matches close to each other on the same scaffold??? Am I missing something here?

I ran bwa mem on 8 cores with otherwise default settings :
bwa mem -t 8 reference readfile1.fastq readfile2.fastq >out.sam

Any input on this would be much appreciated!
Linnéa
Linnea is offline   Reply With Quote
Old 09-16-2013, 10:11 AM   #2
ak352
Junior Member
 
Location: Europe

Join Date: Sep 2013
Posts: 7
Default

Hi Linnea,
Are you sure these are the only two mappings of this fragment? Could you post the results of

grep -w "HWI-ST0866:1861D42ACXX:2:1313:18722:3523" file.bam
__________________
Please answer my unsolved question - http://seqanswers.com/forums/showthread.php?t=33740
ak352 is offline   Reply With Quote
Old 09-16-2013, 10:01 PM   #3
Linnea
Member
 
Location: Uppsala, Sweden

Join Date: Mar 2010
Posts: 23
Default

Thanks for your reply, actually that was the result from grep.. so no other mappings.

However, I think I've just found the problem anyway.. When I checked if the fastq files were sorted properly, I used "head" and "less", and all the read names looked the same between the two files, but checking the line numbers for the pair above showed that they are in fact not on the same line. Turned out that while the beginning of the files were nicely sorted, the rest are not

So no problem with bwa mem, only with my data!
Now I only need to figure out how the h* that could happen...

Sorry for posting!
Linnea is offline   Reply With Quote
Reply

Thread Tools

Posting Rules
You may not post new threads
You may not post replies
You may not post attachments
You may not edit your posts

BB code is On
Smilies are On
[IMG] code is On
HTML code is Off




All times are GMT -8. The time now is 02:03 PM.


Powered by vBulletin® Version 3.8.9
Copyright ©2000 - 2020, vBulletin Solutions, Inc.
Single Sign On provided by vBSSO