This is a test version of Biostars. For the public version, visit https://www.biostars.org.
sam to bam conversion problem

I use bwa mem to align reads. And got several log info at the head of the samfile. Like this:

[M::main_mem] read71286 sequences(10000200 bp)...
[M::mem_pestat]# candidate unique pairs for(FF, FR, RF, RR):(0,5009,0,0) 
[M::mem_pestat]skip orientation FF as there are not enough pairs 
[M::mem_pestat]analyzing insert size distribution for orientation FR...
[M::mem_pestat](25,50,75) percentile:(160,170,174)
[M::mem_pestat]low and high boundaries for computing mean and std. 
[M::mem_pestat]mean and std. dev:(165.99,9.13)
[M::mem_pestat]low and high boundaries for proper pairs:(118,216) 
[M::mem_pestat]skip orientation RF as there are not enough pairs 
[M::mem_pestat]skip orientation RR as there are not enough pairs 
[M::mem_pestat]processed71286 reads in6.210 CPU sec,6.219 real sec

Then I want to convert the sam file into bam, when I use samtools view to do that, It has an error :

# [samopen] no @SQ lines in the header.
# [sam_read1] missing header? Abort!

so I cut the first 11 lines in the head, and then run it again, and I got another error, 

# [samopen] SAM header is present: 25 sequences.
# Parse error at line 71433: sequence and quality are inconsistent

So there are mutiple this log info scattered in the sam file ,so I can not convert into bam.

How could this happen? is there other ways I can convert to bam file?

Thank you in advance.

bwa breakdancer
• 5,161 views • link
updated by mittu1602 mittu1602 • written by jinxinhao1988 jinxinhao1988

Please add the commands you used for bwa mem.

• 0 views • link
written by WouterDeCoster WouterDeCoster

I use the basic command of mem align.

bwa mem test.fasta 1.fq 2.fq > aln_pe.sam

• 0 views • link
written by jinxinhao1988 jinxinhao1988

Sure you didn't use nohup or '&>` or something similar? Because this looks like combining stderr and stdout together.

• 0 views • link
written by WouterDeCoster WouterDeCoster

I use nohup and & in the bwa mem running. but don't use in the samtools running cause I want to see the error info. Sorry I don't quite understand the meaning of stderr and stdout combination.

• 0 views • link
written by jinxinhao1988 jinxinhao1988

WouterDeCoster means standard log message displayed by BWA got concatenated with your standard .sam output file instead of adding to nohup file. Thus you were able to see BWA log message as first 11 lines of your .sam file.

I think you should rerun the BWA mem again without nohup or anything similar followed by samtools view.

And again if you get this error

[samopen] SAM header is present: 25 sequences.

Parse error at line 71433: sequence and quality are inconsistent

then there is surely something wrong in your reads file as error "sequence and quality are inconsistent" means number of total bases in your sequence is not same as number of quality representing those bases

• 0 views • link
written by Tm Tm

I use nohup and & in the bwa mem running.

Well there's your problem.

I use the basic command of mem align.

bwa mem test.fasta 1.fq 2.fq > aln_pe.sam

Next time, when someone asks you about the commands you used, give the commands you used and not a piece of it. You left out important information.

• 0 views • link
written by WouterDeCoster WouterDeCoster

2 answers

Try following the commands below with same parameters:\

bwa mem ucsc.hg19.fasta 1.fastq 2.fastq -M   > test.alnpe.sam
samtools view -Sb test.alnpe.sam >test_bam.alnpe.bam
samtools "sort" -@ 24 test_bam.alnpe.bam >test_bam.alnpe.sort.bam

Good luck!!

• 0 views • link
written by mittu1602 mittu1602

Simplified, without intermediate files:

bwa mem ref.fasta R1.fastq R2.fastq | samtools sort -o sorted_result.bam -
• 0 views • link
written by WouterDeCoster WouterDeCoster

sorry, but how this samtools command is transforming a SAM to a BAM file? You don't use view, just sort. I would appreacite the help, thanks!

• 0 views • link
written by gubrins gubrins

Samtools is smart and will figure out the output format needed based on the extension (.bam) used.

• 0 views • link
written by GenoMax GenoMax

last question, what is the - symbol at the end? If I do that, a lot of temporary files are created

• 0 views • link
written by gubrins gubrins

Temp files are for the sorting process. They will be deleted once the sorting is done. - indicates STDIN i.e. read from there. bwa mem writes to STDOUT, which becomes STDIN for the samtools.

• 0 views • link
written by GenoMax GenoMax

you have not mentioned about the parameters you are using for sam to bam conversion, however you can try using simple basic command:

samtools view -bS in.sam > out.bam

Or you can check out similar post addressing the same issue here

• 0 views • link
written by Tm Tm

The bwa align I use basic command line : bwa mem test.fasta 1.fq 2.fa > aln_pe.sam

and the conversion I use exact like the one you write in the comments.

• 0 views • link
written by jinxinhao1988 jinxinhao1988

one of your file is in fastq and the other in fasta? or just a typo error?

• 0 views • link
written by mittu1602 mittu1602

Sorry, that's a typo error. I am trying to run the bwa mem without nohup as toralmanvar said. I think maybe that's the problem.

• 0 views • link
written by jinxinhao1988 jinxinhao1988

Log in to answer this question.

More posts like this
  • Error of SAM file
    written by Pish • Pish
    I can not view the SAM file or convert SAM to BAM using samtools. As below; > user@computer:$ samtools view aln_orn5v2_4marker.sam | head > [W::sam_read1] …
  • BWA mem skip orientation
    written by williamsbrian5064 williamsbrian5064
    Hi, So I am not entirely sure if this is an issue or not. So to give a little background, I used bwa, picard, samtools, …
  • there are not enough pairs message on screen when running BWA MEM
    written by S AR S AR
    ``` [M::mem_pestat] low and high boundaries for proper pairs: (150, 353) [M::mem_pestat] skip orientation RF as there are not enough pairs [M::mem_pestat] skip orientation RR …
  • sam file error
    written by 2yuna710 • 2yuna710
    After sequences mapping, I got Sam format files. So, I wanted to sort the files by samtools. But I got some error messages as below. …
  • Error when sorting a SAM file with samtools sort, sam file created with bwa mem
    written by mpyle • mpyle
    Hi, I created a .sam file using bwa mem, but when I try to sort the file using samtools sort I get this error: [W::sam_read1] …
  • bwa -mem log messages
    written by aleka aleka
    Hello, I am using bwa mem to aligh paired ends and I get the following errors: [M::process] read 1057202 sequences (80000122 bp)... [M::mem_pestat] # candidate …
  • Skip orientation RF
    written by Folder40g Folder40g
    Hi, I'm using BWA to align some exome sequences: ``` /home/Programas/bwa-0.7.10/bwa index -p genome -a bwtsw /local/Referencias/Homo_sapiens/UCSC/hg19/Sequence/BWAIndex/version0.7.10-r789/genome.fa /home/Programas/bwa-0.7.10/bwa mem -t 10 /local/Referencias/Homo_sapiens/UCSC/hg19/Sequence/BWAIndex/version0.7.10-r789/genome /local/CLL/Raw_data/MONONUCLEARCells_R1.fastq.gz /local/CLL/Raw_data/MONONUCLEARCells_R2.fastq.gz > …
  • BWA MEM 0.7.12 Error : [gzread] <fd:4>: invalid distance code
    written by gsr9999 gsr9999
    Dear BioStars Leaders, I have been running BWA MEM 0.7.12 on a handful of paired-end sample fastq files and it worked out pretty well. I …
  • BWA MEM indexing
    written by david.4.ray • david.4.ray
    I'm a relative newbie and trying t work my way through using bwa mem. I've successfully aligned PE reads to a reference genome previously using …
  • Error while running BWA mem
    written by ravi.uhdnis ravi.uhdnis
    Hi, I am running BWA mem for aligning my PE reads against Human reference genome (GRCh38) although it is running but I am encountering with …