The best way to convert BAM to fastq is to use htslib:

htscmd bamshuf -Ou input.bam tmp-prefix | htscmd bam2fq -s se.fq.gz - | gzip > pe.fq.gz

Firstly it is much faster and more lightweight and secondly it is not affected by the input alignment sorting.

I stumbled over this thread and try to use htscmd from the htslib. On an Ubuntu 12.04 it does is compiled (v 0.2 -1.2), but only htsfile, tabix and bgzip. Are "htscmd" additional scripts not created during compiling?

I try to also extract unpaired reads from a bam file created with the samtools -f4 to extract unmapped reads. Picard; picard-tools simply does not work with the FU= parameter to ADDITIONALLY extract unpaired reads, which is why I would like to try the htslib as suggested by 'lh3'. Help is very appreciated.

Philipp

Hi guys,

I am getting a similar error. Here is what I'm trying to do.

I have a lot of data but I want to focus on specific areas because the software I am using is crashing.

So I go and use samtools to isolate certain coordinates:

samtools view whatever_sorted.bam chr1:115240000-115270000 chr12:25300000-25500000 chr3:178800000-178990000 chr7:55080000-55280000 > fetched.sam

I get rid of the reads that are unpaired (that means they should have the second pair aligned somewhere else in the genome) with

awk '{if ($7 == "=") print $0}' fetched.sam > fetchedpaired.sam

From 1124 unpaired reads according to picard I went down to 90. so I still get the same

SAM validation error: ERROR: Found 90 unpaired mates

when I run Picard SamToFastq

Can it be due to the fact that some reads will aligned on the same chromosome but away from the first read so even if the 7th column is still "=" the reads are considered unpaired??? Or is there something else going on? How can I fix that?

I need to end up having separate files for Read1 and Read2.

Thank you in advance

I ran into the same issue. Here is a solution: you can use samtools to filter for properly paired reads (samtools view -hf 0x2) , this way Picard SamToFastq will not complain about unpaired mates. Here is an example to pipe samtools with Picard to get Read1 and Read2 files for chr21 only.

samtools view -hf 0x2 Sample.sorted.bam chr21 | \
java -jar picard.jar SamToFastq I=/dev/stdin \
FASTQ=Sample_chr21_r1.fq  SECOND_END_FASTQ=Sample_chr21_r2.fq \
UNPAIRED_FASTQ=Sample_chr21_up.fq

Note even though I included UNPAIRED_FASTQ, this file will always be empty because samtools filtered out read pairs where only one read map to chr21.