This is probably really simple, but how do you go from a BED file of genomic features to a wiggle file which contains the number of features at a given position on the chromosome? So, every position in the wiggle file gives the total number of features that overlap with that position in the chromosome.
In additions to @David's solution, if your data is large, you can get directly to bigwig format by using
bedtools genomecov (genomeCoverageBed) and
then bedGraphToBigWig
By the way... is there also a possibility to distinguish between the strands? With my script it's possible, but absolutely not memory efficient for a whole genome, when working with NGS data... but with genomeCoverageBed it's not doable, is it?
ah ok thx, didn't know that. But anyways, I always use a script, since I also want to normalize for multiple hits... I work with short RNAseq and there it happens a lot. And that is not possible with such tools.
I think the simplest way is to write a small perl-script...
Copy this to a file called something like bed2wig.pl:
use strict;
use warnings;
my $file = $ARGV[0];
my %hash = ();
open(FILE, "<$file") || die "cannot open $file\n";
while(<FILE>){
chomp;
my ($chr, $start, $end, $id, $score, $strand) = split(/\s+/,$_);
for(my $i=$start; $i<=$end; $i++){
$hash{$chr}{$i}++;
}
}
close(FILE);
foreach my $chr (keys %hash){
print "variableStep chrom=$chr\n";
foreach my $i (sort {$a<=>$b} keys %{$hash{$chr}}){
print "$i\t$hash{$chr}{$i}\n";
}
}
Now call it like ./bed2wig.pl yourFile.bed > yourFile.wig
It's not very memory efficient for big files (because of the hash), but it should give you an idea of how to do it.
Log in to answer this question.
