view release on metacpan or  search on metacpan

lib/Bio/GeneDesign/Graph.pm  view on Meta::CPAN

=head1 AUTHOR

Sarah Richardson <SMRichardson@lbl.gov>.

=cut

package Bio::GeneDesign::Graph;
require Exporter;

use GD::Graph::lines;
use GD::Graph::colour qw(sorted_colour_list);
use GD::Image;

use strict;
use warnings;

our $VERSION = 5.56;

use base qw(Exporter);
our @EXPORT_OK = qw(
  _make_graph
  _dotplot
  $VERSION
);
our %EXPORT_TAGS =  (GD => \@EXPORT_OK);


=head1 Functions

=head2 _make_graph()

=cut

sub _make_graph
{
  my ($arrref, $window, $orgname, $codon_t, $rscu_t, $revcodon_t) = @_;

  my @sizes = sort {$b <=> $a} map {length($_->seq)} @$arrref;
  my $maxlen = $sizes[0];

  my $graph = GD::Graph::lines->new(1024, 768);
  my @colors = reverse sorted_colour_list(29);
  $graph->set(
    x_label           => 'Window Position (Codon Offset)',
    y_label           => 'Average Relative Synonymous Codon Usage Value',
    title             => "Sliding window of $window using $orgname RSCU values",
    y_max_value       => 1,
    y_min_value       => 0,
    tick_length       => 3,
    y_tick_number     => 1,
    x_label_position  => 0.5,
    y_label_skip      => 0.1,
    x_label_skip      => int($maxlen/50),
    markers           => [1],
    line_width        => 2,
    marker_size       => 2,
    dclrs             => \@colors,
  ) || croak $graph->error;

  my $data = [];
  my @legend;
  my $first = 0;
  my %AAfams = map {$_ => scalar(@{$revcodon_t->{$codon_t->{$_}}})}
               keys %$codon_t;
  my %perc_t = map {$_ => $rscu_t->{$_} / $AAfams{$_}}
               keys %$codon_t;

  foreach my $seqobj (@$arrref)
  {
    my ($x, $y)  = index_codon_percentages($seqobj->seq, $window, \%perc_t);
    push @$data, $x if ($first == 0);
    push @$data, $y;
    $first++;
    push @legend, $seqobj->id;
  }
  $graph->set_legend(@legend);
  my $format = $graph->export_format;
  return ($graph->plot($data)->$format(), $format);
}

=head2 dotplot()

#NO UNIT TESTS

=cut

sub _dotplot
{
  my ($seq1, $seq2, $winsize, $stringency, $outfile) = @_;
  my $Lseq1 = length($seq1);
  my $Lseq2 = length($seq2);

  my $BitMap = GD::Image->new($Lseq1, $Lseq2);

  my $white = $BitMap->colorAllocate(255,255,255);
  my $black = $BitMap->colorAllocate(0,0,0);

  $BitMap->transparent($white);

  for (my $i = 0; $i < $Lseq1 - $winsize; $i++)
  {
    for (my $j = 0; $j < $Lseq2 - $winsize; $j++)
    {
      my $match = 0;
      for (my $w = 0; $w < $winsize; $w++)
      {
        if (substr($seq1, $i + $w, 1) eq substr($seq2, $j + $w, 1))
        {
          $match++;
        }
      }
      if (100 * ($match / $winsize) >= $stringency)
      {
        $BitMap->setPixel($i, $j, $black);
      }
    }
  }
  return $BitMap->png;

  #open   (my $IMG, '>', $outfile) or croak $!;
  #binmode $IMG;
  #print   $IMG $BitMap->png;
  #close   $IMG;
  #return;
}

=head2 index_codon_percentages()

Generates two arrays for x and y values of a graph of codon percentage values.

  in: dna sequence (string),
      window size (integer),
      codon percentage table (hash reference)
  out: x values (array reference), y values (array reference)

=cut