lecture code viewer

downloads

Code

#!/usr/local/bin/perl

$\ = "\n";

my $header;
my %bp;

open(FH,"/home/martink/work/ucsc/hg17/fasta/chr22.fa");
while($line = <FH>) {
  chomp $line;
  # skip the header
  next if $line =~ /^>/;
  # we're in the sequence
  for $bp (split("",$line)) {
    $bp{$bp}++;
  }
}

%bpstats = ();
for $bp (keys %bp) {
  print qq/$bp $bp{$bp}/;
  if($bp =~ /[atcg]/) {
      $bpstats{repeat} += $bp{$bp};
  } elsif ($bp =~ /n/i) {
      $bpstats{padding} += $bp{$bp};
  } elsif ($bp =~ /[GC]/) {
      $bpstats{gc} += $bp{$bp};
  } elsif ($bp =~ /[AT]/) {
      $bpstats{at} += $bp{$bp};
  } else {
      $bpstats{unk} += $bp{$bp};
  }
}

for $statistic (sort {$bpstats{$b} <=> $bpstats{$a}} keys %bpstats) {
    print qq($statistic $bpstats{$statistic});
}