#acl ChrisSeidel:read,write,delete,revert All:read = Solexa Eland Export format to BED file = Here is a perl script which converts from Eland export format to BED format. {{{ #!/usr/bin/perl # Program to convert eland export format to BED format # for running MACS # Chris Seidel, June 2009 # # Requires tab delim file of chromosome or contig names # (eland fa match files) in the format: # UCSC_chr_name chr_length eland_name # corrects for alignments that go off the ends of the chrs # negative bases are trimmed to 1, # bases > chr_length are set to chr_length # (I know the former exist, I don't know if the latter exist) # results are not sorted, but can be sorted in linux by: # sort -o infile.bed -k 1,1 -k 2,2n infile.bed # (sort in place, first column, then by second column numeric) die("usage: $0 chrmap.txt eland_export.txt") unless(scalar(@ARGV) == 2); # create output filename $outfile = $ARGV[1]; $outfile =~ s/\.txt$/\.bed/; open(FOUT, ">$outfile") || die("can't open output file: $outfile"); # get info on chromosomes open(cmap, $ARGV[0]) || die("no chromosome name mapping file!"); %chrmap = {}; while($line = ){ chomp($line); ($newval, $size, $oldval) = split(/\t/, $line); $chrmap{$oldval} = $newval; $chrsize{$oldval} = $size; } # open input file open(fp, $ARGV[1]) || die("can't open eland file"); $lines = 0; while(){ chop; ++$lines; @bits = split(/\t/); # skip reads that didn't pass filtering next if($bits[21] eq "N"); # get match name $seqname = $bits[10]; # skip No Matches or QC failures # next if($seqname =~ /NM|QC/); # skip repeat matches # next if($seqname =~ /\d+:\d+:\d+/); # we're only interested in sequences that match our chrs next unless(exists($chrmap{$seqname})); $seqlen = length($bits[8]); $start = $bits[12]; $end = $start + $seqlen - 1; $strand = $bits[13]; # parse match descriptor $n = ($bits[14] =~ tr/[ACGTN]/[ACGTN]/); # skip reads beyond a certain threshold next if($n > 2); $read_code = "U".$n; # correct for alignments off the chromosome ends if( $start <= 0 ){ print STDERR "start less than or equal to 0: ", $start, "\n"; print STDERR join("\t", @bits), "\n"; $start = 1; } if($end > $chrsize{$seqname}){ print STDERR "end greater than chr end $chrsize{$seqname}: $end, diff: ", $end - $chrsize{$seqname}, "\n"; print STDERR join("\t", @bits), "\n"; $end = $chrsize{$seqname}; } if($strand eq "F"){ $strand = "+"; $color = "0,0,255"; } else{ $strand = "-"; $color = "255,0,0"; } $score = 0; print FOUT join("\t", $chrmap{$seqname}, $start, $end, $read_code, $score, $strand, $start, $end, $color), "\n"; # give some feedback print STDERR "$lines processed\n" if(!($lines % 100000)); } close(FOUT); print STDERR "output file: $outfile\n"; }}} it takes a companion file of chromosome descriptions. This file is used by the program to map the chromosome names used for the alignment - which are often big ugly names to nice short chromosome names. It also has the length of each chromosome. It's a tab delimited text file that you create for your organism to reflect the alignment values. {{{ chr2L 23011544 Drosophila_melanogaster.BDGP5.4.54.dna.chromosome.2L.fa chr2LHet 368872 Drosophila_melanogaster.BDGP5.4.54.dna.chromosome.2LHet.fa chr2R 21146708 Drosophila_melanogaster.BDGP5.4.54.dna.chromosome.2R.fa chr2RHet 3288761 Drosophila_melanogaster.BDGP5.4.54.dna.chromosome.2RHet.fa chr3L 24543557 Drosophila_melanogaster.BDGP5.4.54.dna.chromosome.3L.fa chr3LHet 2555491 Drosophila_melanogaster.BDGP5.4.54.dna.chromosome.3LHet.fa chr3R 27905053 Drosophila_melanogaster.BDGP5.4.54.dna.chromosome.3R.fa chr3RHet 2517507 Drosophila_melanogaster.BDGP5.4.54.dna.chromosome.3RHet.fa chr4 1351857 Drosophila_melanogaster.BDGP5.4.54.dna.chromosome.4.fa chrX 22422827 Drosophila_melanogaster.BDGP5.4.54.dna.chromosome.X.fa chrXHet 204112 Drosophila_melanogaster.BDGP5.4.54.dna.chromosome.XHet.fa chrYHet 347038 Drosophila_melanogaster.BDGP5.4.54.dna.chromosome.YHet.fa }}}