Genome Informatics Laboratory at KRIBB

#!/usr/bin/perl

use Bio::SeqIO;

$window = 5000; # 10 kb
$shift = 5000;   # 2 kb

$seqIn = Bio::SeqIO->new(-file => $ARGV[0], -format => 'fasta');
$seqObj = $seqIn->next_seq();
$seqStr = $seqObj->seq();
$length = $seqObj->length();
@whole_seq = split //, $seqStr;

$cum_GC_skew = 0;

for ($i = 1; $i <= $length; $i += $shift) {
    $end = $i + $shift - 1;
    $end = $length if $end > $length;
    $frag = $seqObj->trunc($i, $end)->seq();
    $G = ($frag =~ s/g/G/ig);
    $C = ($frag =~ s/c/C/ig);
    $GC_skew = ($G - $C) / ($G + $C);
    $cum_GC_skew += $GC_skew;
    $GC = ($G + $C) / $window;
    print join "\t", $i, $GC_skew, $cum_GC_skew, $GC . "\n";
}