Toolbar

Summary

Bio::EnsEMBL::Compara::Production::GenomicAlignBlock::Mavid

Package variables

Included modules

Inherit

Bio::EnsEMBL::Hive::Process

Synopsis

Description

Methods

Methods description

    Title   :   fetch_input
    Usage   :   $self->fetch_input
    Function:   Fetches input data for repeatmasker from the database
    Returns :   none
    Args    :   none

Methods code

sub constraints_file {

  my $self = shift;
  $self->{'_constraints_file'} = shift if(@_);
  return $self->{'_constraints_file'};

}

sub dumpConstraints {

  my $self = shift;
  my $file = $self->fasta_file($WORKDIR . "/contraints_file");
  
}

1;

}

sub dumpFasta {

  my $self = shift;

#  $self->{'comparaDBA'}->dbc->disconnect_when_inactive(1);

  $WORKDIR = $self->worker_temp_directory;

  my $sra = $self->{'comparaDBA'}->get_SyntenyRegionAdaptor;

  my $sr = $sra->fetch_by_dbID($self->synteny_region_id);
  
  my $file = $self->fasta_file($WORKDIR . "/seqs.fa");
  my $masked_file = $WORKDIR . "/seqs.fa.masked";

  open F, ">$file" || throw("Couldn't open $file");
  open MF, ">$masked_file" || throw("Couldn't open $masked_file");
  foreach my $dfr (@{$sr->children}) {
    # WARNING this is a hack. It won't work at all on self comparisons!!!
    # !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
    $self->{'_dnafrag_regions'}{$dfr->dnafrag_id} = $dfr;
    $dfr->retain;
    my $slice = $dfr->slice;
    print F ">" . $dfr->dnafrag_id . "\n";
    print MF ">" . $dfr->dnafrag_id . "\n";
    my $seq = $slice->seq;
    $seq =~ s/(.{80})/$1\n/g;
    chomp $seq;
    print F $seq,"\n";
    $seq = $slice->get_repeatmasked_seq->seq;
    $seq =~ s/(.{80})/$1\n/g;
    chomp $seq;
    print MF $seq,"\n";
  }
  close F;
  close MF;

#  $self->{'comparaDBA'}->dbc->disconnect_when_inactive(0);
  $sr->release;
 
  if ($self->tree_file) {
    $self->dumpTreeFile;
  }
  return 1;

}

sub dumpTreeFile {

  my $self = shift;
  my $tree_file = $self->tree_file;
  open F, $tree_file || throw("Can not open $tree_file");
  my $newick = "";
  while (<F>) {
    chomp;
    if (/^\s*(.*)\s*$/) {
      $newick .= $1;
    }
  }
  close F;
  my $tree = Bio::EnsEMBL::Compara::DBSQL::ProteinTreeAdaptor->parse_newick_into_tree($newick);
  $self->update_node_names($tree);
  $tree_file = $WORKDIR . "/tree_file";
  open F, ">$tree_file";
  print F $tree->newick_simple_format,"\n";
  close F;
  $tree->release;
  $self->tree_file($tree_file);
  return $tree_file;

}

sub fasta_file {

  my $self = shift;
  $self->{'_fasta_file'} = shift if(@_);
  return $self->{'_fasta_file'};

}

sub fetch_input {

  my( $self) = @_;

  $self->{'comparaDBA'} = Bio::EnsEMBL::Compara::Production::DBSQL::DBAdaptor->new(-DBCONN=>$self->db->dbc);

  $self->get_params($self->parameters);
  $self->get_params($self->input_id);


  return 1;

}

sub get_params {

  my $self         = shift;
  my $param_string = shift;

  return unless($param_string);
  print("parsing parameter string : ",$param_string,"\n");

  my $params = eval($param_string);
  return unless($params);

  if(defined($params->{'mavid_constraints'})) {
    $self->mavid_constraints($params->{'mavid_constraints'});
  }
  if(defined($params->{'synteny_region_id'})) {
    $self->synteny_region_id($params->{'synteny_region_id'});
  }
  if(defined($params->{'method_link_species_set_id'})) {
    $self->method_link_species_set_id($params->{'method_link_species_set_id'});
  }
  if(defined($params->{'tree_file'})) {
    $self->tree_file($params->{'tree_file'});
  }
  return 1;

}

sub mavid_constraints {

  my $self = shift;
  $self->{'_mavid_constraints'} = shift if(@_);
  return $self->{'_mavid_constraints'};

}

sub method_link_species_set_id {

  my $self = shift;
  $self->{'_method_link_species_set_id'} = shift if(@_);
  return $self->{'_method_link_species_set_id'};
}

##########################################
#
# internal methods
#
##########################################

}

sub run {

  my $self = shift;
  $self->dumpFasta;
  if ($self->mavid_constraints) {
    $self->dumpConstraints;
  }
  print "fasta_file: ",$self->fasta_file,"\n";
  my $runnable = new Bio::EnsEMBL::Analysis::Runnable::Mavid
    (-workdir => $WORKDIR,
     -fasta_file => $self->fasta_file,
     -tree_file => $self->tree_file,
     -constraints_file => $self->constraints_file,
     -analysis => $self->analysis);
  $self->{'_runnable'} = $runnable;
  $runnable->run_analysis;
#  $self->global_cleanup;

}

sub synteny_region_id {

  my $self = shift;
  $self->{'_synteny_region_id'} = shift if(@_);
  return $self->{'_synteny_region_id'};

}

sub tree_file {

  my $self = shift;
  $self->{'_tree_file'} = shift if(@_);
  return $self->{'_tree_file'};

}

sub update_node_names {

  my $self = shift;
  my $tree = shift;
  my %gdb_id2dfr;
  foreach my $dfr (values %{$self->{'_dnafrag_regions'}}) {
    $gdb_id2dfr{$dfr->dnafrag->genome_db->dbID} = $dfr->dnafrag_id;
  }

  foreach my $leaf (@{$tree->get_all_leaves}) {
    if (defined $gdb_id2dfr{$leaf->name}) {
      $leaf->name($gdb_id2dfr{$leaf->name});
    } else {
#      $leaf->parent->remove_child($leaf);
      $leaf->parent->retain->remove_child($leaf);
    }
  }
  $tree = $tree->minimize_tree;
  if ($tree->get_child_count == 1) {
    my $child = $tree->children->[0];
    $child->parent->merge_children($child);
    $child->disavow_parent;
  }
  return $tree;

}

sub write_output {

  my ($self) = @_;

  my $mlssa = $self->{'comparaDBA'}->get_MethodLinkSpeciesSetAdaptor;
  my $mlss = $mlssa->fetch_by_dbID($self->method_link_species_set_id);
  my $gaba = $self->{'comparaDBA'}->get_GenomicAlignBlockAdaptor;
  foreach my $gab (@{$self->{'_runnable'}->output}) {
    foreach my $ga (@{$gab->genomic_align_array}) {
#      $ga->method_link_species_set_id($self->method_link_species_set_id);
      $ga->method_link_species_set($mlss);
      my $dfr = $self->{'_dnafrag_regions'}{$ga->dnafrag_id};
      $ga->dnafrag_id($dfr->dnafrag_id);
      $ga->dnafrag($dfr->dnafrag);
      $ga->dnafrag_start($dfr->dnafrag_start);
      $ga->dnafrag_end($dfr->dnafrag_end);
      $ga->dnafrag_strand($dfr->dnafrag_strand);
      $ga->level_id(1);
      $dfr->release;
      unless (defined $gab->length) {
        $gab->length(length($ga->aligned_sequence));
      }
    }
#    $gab->method_link_species_set_id($self->method_link_species_set_id);
    $gab->method_link_species_set($mlss);
    $gaba->store($gab);
  }
  return 1;
}

##########################################
#
# getter/setter methods
# 
##########################################

#sub input_dir {
#  my $self = shift;
#  $self->{'_input_dir'} = shift if(@_);
#  return $self->{'_input_dir'};
#

}

General documentation

Abel Ureta-Vidal <abel@ebi.ac.uk>

The rest of the documentation details each of the object methods.
Internal methods are usually preceded with a _