problem sorting bed and bam files

Question

[gubrins]

Heys,

I am having problems when trying to calculate coverage levels of a bed file with genomeCoverageBed (v2.3), as my bam and bed files are sorted differently:

cut -f 1 combined_after_scaffolding.bed |uniq
    scaffold_1 scaffold_10
samtools view -H both_not_filtered_sorted.bam
    àHD VN:1.6 SO:queryname àSQ SN:scaffold_1 LN:233710445 àSQ SN:scaffold_2 LN:214025719

However, when I try to sort it with the genome.file like this:

bedtools sort -g genomefile.txt -i combined_after_scaffolding.bed > sorted.bed

My computer rans out of memory (250Gb). Any solutions to that?

Thanks in advance!

Answer 1

[ATpoint]

First of all genomeCoverageBed takes either a bam or bed file, not both. In any case, the bam must be sorted by position so use samtools sort -o sorted.bam in.bam and for sorting bed files one generally uses sort -k1,1 -k2,2n in.bed > sorted.bed. That is just the "normal" GNU sort. On mac one can get that by installing coreutils, e.g. via brew. bedtools is great for my things but not for sorting, don't use it for that. The mentioned sort is just GNU sort which is fast and efficient. See its manual for options to play with memory and threads to make it faster.

[gubrins]

thanks ATpoint for answering! Both bam and bed were sorted as you mentioned. However, when I try to run genomeCoverageBed with a bam file, I get this error message: Input error: Chromosome scaffold_17 found in non-sequential lines. This suggests that the input file is not sorted correctly.

With the bed file I manage to obtain a coverage file. However, when I run bedGraphToBigWig, I obtain this error: end (297) before start (343) line 344 of coverage.bed

Not sure why, in my coverage file I have endings that are previous to the beginnings, as for example this:

scaffold_1 343 297 233710445 1.2708e-06

Any idea of what is going on?