package GDxBase::GeneModel::Collection; =head1 Name GDxBase::GeneModel::Collection =head1 Description Last Revision: 25-Jul-06 Author: James Allen Email: james.allen@cimr.cam.ac.uk =cut use warnings; use strict; use Carp; use Config::IniFiles; use GDxBase::Globals; use Data::Dumper; use Bio::Graphics::Panel; use Class::AutoClass; our @ISA = "Class::AutoClass"; our (@AUTO_ATTRIBUTES, %DEFAULTS); @AUTO_ATTRIBUTES = qw( common_name build chromosome gene_id gene_name gene_start gene_stop collection_start collection_stop single_panel ignore_build show_gene_ref gene_models features image_dir http_image_dir config_file ); %DEFAULTS = ( gene_models => [], features => {} ); Class::AutoClass::declare(__PACKAGE__); ################################################################################ =head1 Methods =cut =head2 new Title: new Usage: GDxBase::GeneModel::Collection->new() Function: Returns: Args: common_name - mandatory !optionalise! build - mandatory !optionalise! chromosome - mandatory !optionalise! gene_id - mandatory gene_name - mandatory !optionalise! gene_start - mandatory !optionalise! gene_stop - mandatory !optionalise! collection_start - optional, same as gene_start by default collection_stop - optional, same as gene_stop by default show_gene_ref - optional, undef by default show_summary - optional, true by default =cut sub _init_self { my $self = shift; #confess "A taxonomy id (e.g. 9606) must be supplied" unless $self->tax_id; # Check for other required fields too. make sure dirs have trailing slash $self->collection_start($self->gene_start) unless $self->collection_start; $self->collection_stop($self->gene_stop) unless $self->collection_stop; } ################################################################################ =head2 add_gene_model Title: add_gene_model Usage: GDxBase::GeneModel::Collection->add_gene_model($gene_model) Function: Returns: Args: =cut sub add_gene_model { my $self = shift; my %args = @_; confess "'gene_model' argument required" unless exists $args{"gene_model"}; confess "'db_config_section' argument required" unless exists $args{"db_config_section"}; my $model_prefix = $args{"model_prefix"} || 'GM'; my $gene_model = $args{"gene_model"}; my $global_config = $self->_ini_file_section("GLOBAL"); my $db_config = $self->_ini_file_section($args{"db_config_section"}); my $source_config = $self->_ini_file_section($model_prefix."_".$gene_model->source); my $core_config = $self->_ini_file_section($model_prefix."_GFF"); my %config = (%$global_config, %$db_config, %$core_config, %$source_config); #%config = (%$db_config, %$source_config); $config{host} = $GDxBase::Globals::dbhandler->get_host($args{"db_config_section"}); my ($features, $model_start, $model_stop) = $gene_model->get_features($self, \%config); if (defined $model_start && defined $model_stop) { $self->collection_start($model_start) if $model_start < $self->collection_start; $self->collection_stop($model_stop) if $model_stop > $self->collection_stop; } if ($features && (scalar keys %{$features})) { push @{$self->gene_models}, [$gene_model, $features]; ${$self->features}{$gene_model->source} = $features; return (scalar keys %{$features}); } else { return(0); } } ################################################################################ =head2 render_gene_models Title: render_gene_models Usage: GDxBase::GeneModel::Collection->render_gene_models Function: Returns: Args: =cut sub render_gene_models { my $self = shift; my %args = @_; my $model_prefix = $args{"model_prefix"} || 'GM'; my (@images, @maps); # If the 'show_gene_ref' hasn't been explicitly set as a parameter, then # set it depending on whether the panel bounds are beyond the gene bounds. if (not defined($self->show_gene_ref)) { if (($self->collection_start != $self->gene_start) or ($self->collection_stop != $self->gene_stop)) { $self->show_gene_ref(1); } else { $self->show_gene_ref(0); } } # Get config data here for panel my $panel_config = $self->_ini_file_section($model_prefix."_PANEL_SETTINGS" ); if ($self->single_panel) { my $panel = $self->_render_panel($panel_config); $self->_render_scale($panel); if ($self->show_gene_ref) { my $gene_ref_config = $self->_ini_file_section($model_prefix."_GENE_REF" ); $self->_render_gene_ref($panel, $gene_ref_config); } my @gene_model_order; if (exists $args{"gene_model_order"}) { @gene_model_order = @{$args{"gene_model_order"}}; } else { @gene_model_order = sort @{keys %{$self->features}}; } foreach my $source (@gene_model_order) { if (exists ${$self->features}{$source}) { my $features = ${$self->features}{$source}; my $gene_model_config = $self->_ini_file_section( $model_prefix."_".$source ); $self->_render_tracks($panel, $features, $gene_model_config); } } my ($img, $map) = $self->_create_map($panel); push @images, $img; my $title = "Gene Models"; $title .= " - Build ".$self->build unless $self->ignore_build; my $map_hash; $$map_hash{"title"} = $title; $$map_hash{"location"} = $map; push @maps, $map_hash; } else { # With separate panels it's sometimes nice to combine related # gene models - we manage this by getting a list of models that # will appear underneath other models, on the same panel. A little # sleight of hand and hey presto... my $gene_model_merge; $gene_model_merge = $args{"gene_model_merge"} if exists $args{"gene_model_merge"}; my %to_be_merged; if ($gene_model_merge) { foreach my $source_list (values %{$gene_model_merge}) { foreach my $source (@{$source_list}) { $to_be_merged{$source} = 1; } } } my @gene_model_order; if (exists $args{"gene_model_order"}) { @gene_model_order = @{$args{"gene_model_order"}}; } else { @gene_model_order = sort @{keys %{$self->features}}; } foreach my $source (@gene_model_order) { if (exists ${$self->features}{$source}) { next if exists $to_be_merged{$source}; foreach my $fref (sort keys %{${$self->features}{$source}}){ my $panel = $self->_render_panel($panel_config); my $label; if ($self->chromosome eq $fref) { my $chr = $self->chromosome; $chr =~ s/chr//; $label = $self->common_name." Chromosome ".$chr; $label .= " Build ".$self->build unless $self->ignore_build; } else { $label = $fref; } $label .= ": ".$self->_insert_commas($self->collection_start); $label .= " - ".$self->_insert_commas($self->collection_stop); $label .= " (size ".$self->_insert_commas($self->collection_stop - $self->collection_start + 1).")"; $self->_render_scale($panel, $label); if ($self->show_gene_ref) { my $gene_ref_config = $self->_ini_file_section( $model_prefix."_GENE_REF" ); my ($ref_name, $ref_start, $ref_stop, $link); if ($self->chromosome eq $fref) { $ref_name = $self->gene_name." Gene"; $link = $self->gene_id; $ref_start = $self->gene_start; $ref_stop = $self->gene_stop; } else { $ref_name = $fref; $link = $fref; $ref_start = $self->collection_start; $ref_stop = $ref_start; foreach my $feature (@{${$self->features}{$source}{$fref}}) { $ref_stop = $feature->end if $feature->end > $ref_stop; } } my $gene_ref = Bio::SeqFeature::Generic->new( -display_name => $ref_name, -source_tag => $$gene_ref_config{"source"}, -start => $ref_start, -end => $ref_stop, -strand => 0, -tag => { "link" => $link }); $self->_render_gene_ref($panel, $gene_ref_config, $gene_ref); } my $gene_model_config = $self->_ini_file_section( $model_prefix."_".$source ); $self->_render_tracks($panel, ${$self->features}{$source}{$fref}, $gene_model_config); if ($gene_model_merge) { if (exists $$gene_model_merge{$source}) { foreach my $merge_source (@{$$gene_model_merge{$source}}) { if (exists ${$self->features}{$merge_source}) { $gene_model_config = $self->_ini_file_section( $model_prefix."_".$merge_source ); $self->_render_tracks($panel, ${$self->features}{$merge_source}{$fref}, $gene_model_config); } } } } my $title; if ($self->chromosome eq $fref) { $title = $source." Gene Model"; $title .= " - Build ".$self->build unless $self->ignore_build; } else { $title = $fref; } my ($img, $map) = $self->_create_map($panel,$self->build); push @images, $img; my $map_hash; $$map_hash{"title"} = $title; $$map_hash{"location"} = $map; push @maps, $map_hash; } } } } return (\@images, \@maps); } ################################################################################ =head2 _render_panel Title: _render_panel Usage: GeneModel->_render_panel($panel_config) Function: Create an empty Panel object Returns: A Bio::Graphics::Panel object Args: panel_config - mandatory, a hash of config params =cut sub _render_panel { my $self = shift; my ($panel_config) = @_; my %options; $options{-start} = $self->collection_start; $options{-stop} = $self->collection_stop; $options{-width} = $$panel_config{"width"}; $options{-pad_left} = $$panel_config{"pad_left"}; $options{-pad_right} = $$panel_config{"pad_right"}; $options{-key_style} = $$panel_config{"key_style"}; $options{-bgcolor} = $$panel_config{"bgcolor"}; $options{-grid} = 1 if $$panel_config{"gridcolor"} ne ""; $options{-gridcolor} = $$panel_config{"gridcolor"}; my $panel = Bio::Graphics::Panel->new(%options); return ($panel); } ################################################################################ =head2 _render_scale Title: _render_scale Usage: GeneModel->_render_scale($panel) Function: Add a scale track to a Panel object Returns: A Bio::Graphics::Panel object Args: panel - mandatory, a Bio::Graphics::Panel object =cut sub _render_scale { my $self = shift; my ($panel, $label) = @_; my $scale = Bio::SeqFeature::Generic->new( -start => $self->collection_start, -end => $self->collection_stop, -tag => {mouseover_text => "Scale"} ); $panel->unshift_track([$scale], -glyph => "arrow", -tick => 1, -height => 5, -fgcolor => "black", -double => 1, -font2color => "blue", -label => $label); } ################################################################################ =head2 _render_gene_ref Title: _render_gene_ref Usage: GeneModel->_render_gene_ref($panel) Function: Add a gene track to a Panel object Returns: A Bio::Graphics::Panel object Args: panel - mandatory, a Bio::Graphics::Panel object Notes: A track showing the gene is useful for reference if the transcripts extend beyond the start and/or stop points of the gene. Display is optional, controlled by the 'show_gene_ref' flag - if it is not supplied as a parameter, then the gene will be displayed if the bounds of the panel have extended beyond the bounds of the gene, for any of the gene models. =cut sub _render_gene_ref { my $self = shift; my ($panel, $gene_ref_config, $gene_ref) = @_; my %options; if ($$gene_ref_config{"key"}) { $options{-key} = $$gene_ref_config{"key"}; } else { my $key = $gene_ref->display_name.": "; $key .= $self->_insert_commas($gene_ref->start); $key .= " - ".$self->_insert_commas($gene_ref->end)." (size "; $key .= $self->_insert_commas($gene_ref->end - $gene_ref->start + 1).")"; $options{-key} = $key; } $options{-glyph} = $$gene_ref_config{"glyph"}; $options{-bgcolor} = $$gene_ref_config{"bgcolor"}; $options{-fgcolor} = $$gene_ref_config{"fgcolor"}; $panel->add_track($gene_ref, %options); } ################################################################################ =head2 _render_tracks Title: _render_tracks Usage: GeneModel->_render_tracks($panel, $features) Function: Add tracks representing features to a Panel Returns: A Bio::Graphics::Panel object Args: panel - mandatory, a Bio::Graphics::Panel object features - mandatory, a hash of feature arrays Notes: The appearance of the transcripts can be controlled in three ways: 1) Defaults for glyph style and colour are used in the absence of other information. 2) One or more of the defaults can be overridden when creating the GeneModel object. 3) Within the 'get_features' method of a subclass, a feature can be tagged with a style or colour; e.g. $feature->add_tag_value("bgcolor", "red"). See the _render_tracks code for details of what tag names are appropriate. These tags override the defaults and any constructor parameters. =cut sub _render_tracks { my $self = shift; my ($panel, $features, $track_config) = @_; # Panels can have multiple tracks and the ordering should be consistent - we # order alphabetically, with numbers coming after letters. There's scope for # allowing user-defined sorts in the future, but for now this is it. sub feature_sort { if ($a =~ /^[0-9]+/ or $b =~ /^[0-9]+/) { $b cmp $a; } else { $a cmp $b; } } # We allow each feature set to override any of the formatting options, # by the use of tags with the same name as the relevant option (except # for 'description' as that clashes with a pre-existing tag name). If an # option doesn't exist in the config file (ie doesn't exist in the # $track_config hash) then the default value will be used. $panel->add_track( $features, -key => $$track_config{"key"}, -glyph => sub { my $f=shift; if ($f->has_tag("glyph")) { join("", $f->get_tag_values("glyph")); } else { $$track_config{"glyph"}; } }, -bgcolor => sub { my $f=shift; if ($f->has_tag("bgcolor")) { join("", $f->get_tag_values("bgcolor")); } else { $$track_config{"bgcolor"}; } }, -fgcolor => sub { my $f=shift; if ($f->has_tag("fgcolor")) { join("", $f->get_tag_values("fgcolor")); } else { $$track_config{"fgcolor"}; } }, -font2color => sub { my $f=shift; if ($f->has_tag("fgcolor")) { join("", $f->get_tag_values("fgcolor")); } else { $$track_config{"fgcolor"}; } }, -utr_color => sub { my $f=shift; if ($f->has_tag("utr_color")) { join("", $f->get_tag_values("utr_color")); } else { $$track_config{"utr_color"}; } }, -label => sub { my $f=shift; return $f->display_name; }, -description => sub { my $f=shift; if ($f->has_tag("desc")) { join("", $f->get_tag_values("desc")); } }, -height => sub { my $f=shift; if ($f->has_tag("height")) { join("", $f->get_tag_values("height")); } else { 8; } }, -box_subparts => $$track_config{"box_subparts"}, -connector => $$track_config{"connector"}, -bump => $$track_config{"bump"}, -all_callbacks => $$track_config{"all_callbacks"}, -sort_order => "high_score|name", -implied_utrs => 0 ); } ################################################################################ =head2 create_map Title: _create_map Usage: GeneModel->_create_map($panel) Function: Generate an image, then an image map, from a Panel object Returns: A title for the panel, and the location of the map file Args: panel - mandatory, a Bio::Graphics::Panel object =cut sub _create_map { my $self = shift; my ($panel,$build) = @_; my ($boxes, $width, $height, $map_name, $img_dir, $http_img_dir, $map_file, $map); my ($img_file) = $self->_create_image($panel); $boxes = $panel->boxes; ($width, $height) = $panel->gd->getBounds; $map_name = "Map_for_$img_file"; $img_dir = $self->image_dir; $http_img_dir = $self->http_image_dir; $map_file = "$img_file"."_map"; $map = qq(\n); foreach my $box (@$boxes) { my $coordinates = "$box->[1], $box->[2], $box->[3], $box->[4]"; $map .= qq([0]->has_tag("mouseover_text")) { $mouseover_text .= join("", $box->[0]->get_tag_values("mouseover_text")); } elsif ($box->[0]->display_name) { $mouseover_text .= $box->[0]->display_name; } $mouseover_text .= " ("; $mouseover_text .= $self->_insert_commas($box->[0]->start); $mouseover_text .= " - "; $mouseover_text .= $self->_insert_commas($box->[0]->end); $mouseover_text .= ") "; $map .= qq(title="$mouseover_text" ); $map .= qq(alt="$mouseover_text" ); my $href = $self->_map_link($box->[0],$build); if ($href) { $map .= qq(href="$href" ); } $map .= qq(target="new" />\n); } $map .= qq(\n); $map .= qq(\n); open(OUT, ">$img_dir$map_file") || confess "Cannot open $img_dir$map_file"; print OUT $map; close(OUT); `chmod a+w $img_dir$map_file`; return ("$img_dir$img_file", "$img_dir$map_file"); } ################################################################################ =head2 _create_image Title: _create_image Usage: GeneModel->_create_image($panel) Function: Save a Panel object as an image file in a given directory Returns: Image file location Args: panel - mandatory, a Bio::Graphics::Panel object =cut sub _create_image { my $self = shift; my ($panel) = @_; my ($img_dir, $img_file); $img_dir = $self->image_dir; $img_file = "gene_model_" . time . int(rand 999999) . ".png"; $img_file =~ s/\s/_/g; open(OUT, ">$img_dir$img_file") || confess "Cannot open $img_dir$img_file"; print OUT $panel->png; close(OUT); `chmod a+w $img_dir$img_file`; return ($img_file); } ################################################################################ =head2 _map_link Title: _map_link Usage: GeneModel->_map_link($cfg, $feature) Function: Create a link that corresponds to a glyph in an image Returns: String containing a hyperlink Args: cfg - mandatory, a Config::IniFiles object feature - mandatory, a Bio::SeqFeature::Generic or a Bio::DB::GFF::Feature object =cut sub _map_link { my $self = shift; my ($feature,$build) = @_; my ($links, $href); if ($feature->source_tag && $feature->display_name) { my $common_name = lc($self->common_name); my $source = uc($feature->source_tag); $links = $self->_ini_file_section("LINKS $source"); if (keys %$links) { if (exists $$links{$common_name}) { $href = $$links{$common_name}; if($build=~ m/NCBI36/){ if($href=~ m/genome.ucsc.edu/ ){ $href=~ s/hg19/hg18/; } if($href=~ m/ensembl.org/){ $href=~ s/www.ensembl.org/may2009.archive.ensembl.org/; } } } else { $href = $$links{"link_out"}; } if (exists $$links{"link_out_type"}) { if ($$links{"link_out_type"} eq "name") { $href .= $feature->display_name; } elsif ($$links{"link_out_type"} eq "position") { my $chr = $self->chromosome; $chr =~ s/chr//; $href .= "Chr".$chr.":".$self->collection_start."-".$self->collection_stop; } elsif ($$links{"link_out_type"} eq "other") { my $link = eval {join("", $feature->get_tag_values("link"))}; $href .= $link if $link; } } else { $href .= $feature->display_name; } } } # Make sure that links are xhtml compliant. $href =~ s/&/&/g if $href; return $href; } ################################################################################ =head2 _insert_commas Title: _insert_commas Usage: GeneModel->_insert_commas($number) Function: Format a number with commas to make it more readable Returns: String containing a comma-ified number Args: number - mandatory =cut sub _insert_commas { my $self = shift; my $number = shift; 1 while $number =~ s/^(.*\d)(\d{3})/$1,$2/; return $number; } ################################################################################ =head2 ini_file_section Title: ini_file_section Usage: Function: Returns: Args: =cut sub _ini_file_section { my $self = shift; my ($section) = @_; confess "'section' argument required" unless $section; my %config = (); my $cfg = $self->config_file; my @parameters = $cfg->Parameters($section); foreach my $parameter (@parameters) { my $value = $cfg->val($section, $parameter); if ($value =~ /\n/) { my @value = split("\n", $value); $config{$parameter} = \@value; } else { $config{$parameter} = $value; } } return(\%config); } ################################################################################ 1;