Re: [bedtools-discuss] bamToBed question: output bed is empty

[Ivan Gregoretti]

Hi Lina,

I am not sure that an empty file is what you should expect but what I
think is that you should run samtools twice.

The first time, use the flag -f 64 to select the first in the pair.
The second time use -f 128 to get the second member of the pair.

Finally, run bamToBed on each file separately.

For information about flags see
http://picard.sourceforge.net/explain-flags.html

I hope this helps.

Ivan



Ivan Gregoretti, PhD

On Fri, Sep 14, 2012 at 11:28 AM, Linda C <linda...@gmail.com> wrote:
> Hi,
>
> I ran bowtie 2 to mapped pair-end fastq files to the genome. The output is
> sam file. I then ran samtools to convert it to bam file. Finally, I ran
> bedtools to convert bam to bed but I got an empty bed file. Here are all
> cmd lines I used:
>
> bowtie2 -x SacCer2 -1 s_1_1.fastq -2 s_1_2.fastq -S s_1.sam
> samtools view -bS s_1.sam > s_1.bam
> bamToBed -i s_1.bam > s_1.bed
>
>
> I'm wondering if it is because the file is pair-end so I have to use
> different parameter setting to convert bam file into bed. Thank you.
>
> Linda
>

[Aaron Quinlan]

Could you provide a peak at the first 10 lines of the BAM file?

samtools view s_1.bam | head

- Aaron

quinlanlab.org

[Ivan Gregoretti]

That is a bit short and won't let us see past the references.

Try

samtools view s_1.bam | head -30

(Aaron, sorry to meddle but I just could not resist.)

Ivan


Ivan Gregoretti, PhD


On Fri, Sep 14, 2012 at 11:56 AM, Linda C <linda...@gmail.com> wrote:
> I got this:
>
> [bam_header_read] invalid BAM binary header (this is not a BAM file).
> [main_samview] fail to read the header from "s_1.bam".
>
> And this is head of my sam file:
>
> @HD VN:1.0 SO:unsorted
> @SQ SN:chrI LN:230218
> @SQ SN:chrII LN:813184
> @SQ SN:chrIII LN:316620
> @SQ SN:chrIV LN:1531933
> @SQ SN:chrV LN:576874
> @SQ SN:chrVI LN:270161
> @SQ SN:chrVII LN:1090940
> @SQ SN:chrVIII LN:562643
> @SQ SN:chrIX LN:439888
> Any hints?
> Thank you,
> Kajia

[Ivan Gregoretti]

Actually, another diagnosis that could help is this

samtools flagstat s_1.bam

Ivan


Ivan Gregoretti, PhD

[Aaron Quinlan]

Not at all, thanks the point of a mailing list - community help.

[Ivan Gregoretti]

Hi Kajia,

There you have it. Your BAM file suffers some pathology.

Things to try

1) re run samtools view -bS s_1.sam > s_1.bam. You may have mistyped a
character at the command line.

2) If 1 produces the same result, try Picard instead of samttols.

3) Try visualising the BAM file using IGV. Note that you'll need to
sort and index the file. You can do that too with picard.


If the problems persist, I recommend emailing the samtools email list.
This is not really a BEDTools bug, in my opinion. Perhaps it will help
realise samtools developers that a user could use something a little
bit more informative than "this is not a BAM file".

By the way, since the .bam file is actually not a BAM file, before you
write to the samtools list, try this at the command line

file s_1.bam

That should tell what type of file that is.

Ivan


Ivan Gregoretti, PhD


On Fri, Sep 14, 2012 at 12:31 PM, Linda C <linda...@gmail.com> wrote:
> Hi Ivan
>
> It gave me the error msg:

>
> [bam_header_read] invalid BAM binary header (this is not a BAM file).
> [main_samview] fail to read the header from "s_1.bam".
>

> Should I convert sam to bam for the two files separately as you suggested?
> My goal is to upload the bgr files to UCSC genome browser. So I need to do:
> sam --> bam --> bed --> bgr (bedgraph). So if I deal the two files
> separately, evetually, I should get two tracks on the genome browser for
> each sample (as they are pair-end)?
>
> I'm totally new to the pair-end sequencing. Thank you for your input!