converting sam to bam format samtools
1
0
Entering edit mode
7.6 years ago
Hans ▴ 130

Hello

This is my first time runing bwa and samtools. I am not so similar with the terminology.

I have generated sam file with the command:

bwa sampe ref.fa aln_sa1.sai aln_sa2.sai read1.fq read2.fq > aln-pe.sam

The read files are actually files that were processed by the servicing lab and have only unique reads. they look like this:

>@SBS123:70:B0563ABXX:3:1101:1246:2041 1:N:0:GCCAAT
NTAAGGCTACCTAAAGGAGCTCATGCANCATTANCTATGGTANTTGTCATCTTCACTTTGTACAGTACTGCAGCATTGGGTCAGCTCTTCTTTTTCCATGT

The sam file head looks like this:

@SQ     SN:UCW_Tu-k45_contig_21575;tu-k61_contig_17060  LN:1805
@SQ     SN:UCW_Tu-k25_contig_8865;tu-k31_contig_8795    LN:3898

When I run:

samtools view -bS aln.sam > aln.bam

I get this message:

[sam_header_line_parse] expected '@XY', got [@SBS123:70:B0563ABXX:3:1101:1206:2056    83    UCW_Tu-k41_contig_1687    3010    60    101M    =    2977    -134    CTGATGACCCATGTGTGATGAACAGCTACCTCTTTGTTCAAGAGCTAGTGATGGAGCTAGAAATGTCNACTTGTAGTGTAAGTTGGCAAAATGACTGCTTN    *    XT:A:U    NM:i:5    SM:i:25    AM:i:25    X0:i:1    X1:i:0    XM:i:5    XO:i:0    XG:i:0    MD:Z:1C8T56C2A29T0]
Hint: The header tags must be tab-separated.
[samopen] no @SQ lines in the header.
[sam_read1] missing header? Abort!

Any suggestions are welcome

Thank you

samtools sam bam • 3.5k views
ADD COMMENT
0
Entering edit mode

can you please show us the output of 'head -n 20 aln-pe.sam'

ADD REPLY
0
Entering edit mode

Thank you for the quick response. The next 20 lines will look the same as the first two. I forgot to say this is Rseq analysis against 88000 transcriptome contigs.

ADD REPLY
0
Entering edit mode

Hello

The problem was solved by removing the @ signs from the reads in the aln.sam file. I did it with the command

sed -i 's/@SBS/SBS/g' aln.sam

Thank you for your help

ADD REPLY
3
Entering edit mode
7.6 years ago

It looks like your aligner didn't remove the @ signs from in front of the read names! That'll cause the problem you're seeing. Can you show the first read in the fastq file?

ADD COMMENT
0
Entering edit mode

Thank you for the quick response. I do not have the original fastq files. I only got the processed fasta files. I can try remove the @ with sed. Thank you

ADD REPLY

Login before adding your answer.

Traffic: 701 users visited in the last hour
Help About
FAQ
Access RSS
API
Stats

Use of this site constitutes acceptance of our User Agreement and Privacy Policy.

Powered by the version 2.3.6