Commit 6112e049 authored by Dan T Andrews's avatar Dan T Andrews
Browse files

Added some functionality from StaticGoldenPathAdaptor. Added methods...

Added some functionality from StaticGoldenPathAdaptor.  Added methods fetch_Slice_by_contig, fetch_Slice_by_clone, fetch_Slice_by_gene, and three attendant methods get_chr_start_end_of_contig, fetch_Slice_by_chr_start_end and get_Gene_chr_bp.  Checked generally that these are returning what is expected sequence-wise.
parent 822e87fc
......@@ -44,6 +44,7 @@ methods. Internal methods are usually preceded with a _
package Bio::EnsEMBL::DBSQL::SliceAdaptor;
use Bio::EnsEMBL::Utils::Eprof qw(eprof_start eprof_end);
use vars qw(@ISA);
use strict;
......@@ -204,3 +205,273 @@ sub fetch_all_similarity_features_above_pid{
return @out;
}
=head2 get_chr_start_end_of_contig
Title : get_chr_start_end_of_contig
Usage :
Function: returns the chromosome name, absolute start and absolute end of the
specified contig
Returns : returns chr,start,end
Args : contig id
=cut
sub get_chr_start_end_of_contig {
my ($self,$contigid) = @_;
if( !defined $contigid ) {
$self->throw("Must have contig id to fetch Slice of contig");
}
my $type = $self->db->static_golden_path_type()
or $self->throw("No assembly type defined");
my $sth = $self->db->prepare("SELECT c.name,
a.chr_start,
a.chr_end,
a.chromosome_id
FROM assembly a, contig c
WHERE c.name = '$contigid'
AND c.contig_id = a.contig_id
AND a.type = '$type'"
);
$sth->execute();
my ($contig,$start,$end,$chr_name) = $sth->fetchrow_array;
if( !defined $contig ) {
$self->throw("Contig $contigid is not on the golden path of type $type");
}
return ($chr_name,$start,$end);
}
=head2 fetch_Slice_by_chr_start_end
Title : fetch_Slice_by_chr_start_end
Usage :
Function: create a Slice based on a segment of a chromosome and
start/end
Example :
Returns : A Slice
Args : chromosome, start, end (in Chromosome coordinates)
=cut
sub fetch_Slice_by_chr_start_end {
my ($self,$chr,$start,$end) = @_;
if( !defined $end ) { # Why defined? Is '0' a valid end?
$self->throw("must provide chr, start and end");
}
if( $start > $end ) {
$self->throw("start must be less than end: parameters $chr:$start:$end");
}
my $slice;
&eprof_start('Slice: staticcontig build');
my $type = $self->db->static_golden_path_type();
eval {
$slice = Bio::EnsEMBL::Slice->new(
-chr_name => $chr,
-chr_start => $start,
-chr_end => $end,
-assembly_type => $type,
-adaptor => $self->db->get_SliceAdaptor
);
} ;
if( $@ ) {
$self->throw("Unable to build a slice for $chr, $start,$end\n\nUnderlying exception $@\n");
}
&eprof_end('Slice: staticcontig build');
return $slice;
}
=head2 fetch_Slice_by_contig
Title : fetch_Slice_by_contig
Usage : $slice = $slice_adaptor->fetch_Slice_by_contig('AC000012.00001',1000);
Function: Creates a slice of the specified slice adaptor object. If a context size is given, the slice is extended by that number of basepairs on either side of the contig. Throws if the contig is not golden.
Returns : Slice object
Args : contig id, [context size in bp]
=cut
sub fetch_Slice_by_contig{
my ($self,$contigid,$size) = @_;
if( !defined $size ) {$size=0;}
my ($chr_name,$start,$end) = $self->get_chr_start_end_of_contig($contigid);
return $self->fetch_Slice_by_chr_start_end( $chr_name,
$start-$size,
$end+$size
);
}
=head2 fetch_Slice_by_clone
Title : fetch_Slice_by_clone
Usage : $slice = $slice_adaptor->fetch_Slice_by_clone('AC000012',1000);
Function: Creates a Slice of the specified object. If a context size is given, the Slice is extended by that number of basepairs on either side of the clone. Throws if the clone is not golden.
Returns : Slice object
Args : clone id, [context size in bp]
=cut
sub fetch_Slice_by_clone{
my ($self,$clone,$size) = @_;
if( !defined $clone ) {
$self->throw("Must have clone to fetch Slice of clone");
}
if( !defined $size ) {$size=0;}
my $type = $self->db->static_golden_path_type()
or $self->throw("No assembly type defined");
my $sth = $self->db->prepare("SELECT c.name,
a.chr_start,
a.chr_end,
a.chromosome_id
FROM assembly a,
contig c,
clone cl
WHERE c.clone_id = cl.clone_id
AND cl.name = '$clone'
AND c.contig_id = a.contig_id
AND a.type = '$type'
ORDER BY a.chr_start"
);
$sth->execute();
my ($contig,$start,$end,$chr_name);
my $counter;
my $first_start;
while ( my @row=$sth->fetchrow_array){
$counter++;
($contig,$start,$end,$chr_name)=@row;
if ($counter==1){$first_start=$start;}
}
if( !defined $contig ) {
$self->throw("Clone is not on the golden path. Cannot build Slice");
}
my $slice = $self->fetch_Slice_by_chr_start_end( $chr_name,
$first_start-$size,
$end+$size
);
$slice->adaptor->db($self->db);
return $slice;
}
=head2 get_Gene_chr_bp
Title : get_Gene_chr_bp
Usage :
Function:
Returns :
Args :
=cut
sub get_Gene_chr_bp {
my ($self,$geneid) = @_;
my $type = $self->db->static_golden_path_type()
or $self->throw("No assembly type defined");
my $sth = $self->db->prepare("SELECT
if(a.contig_ori=1,(e.contig_start-a.contig_start+a.chr_start),
(a.chr_start+a.contig_end-e.contig_end)),
if(a.contig_ori=1,(e.contig_end-a.contig_start+a.chr_start),
(a.chr_start+a.contig_end-e.contig_start)),
a.chromosome_id
FROM exon e,
transcript tr,
exon_transcript et,
assembly a,
gene_stable_id gsi
WHERE e.exon_id=et.exon_id
AND et.transcript_id =tr.transcript_id
AND a.contig_id=e.contig_id
AND a.type = '$type'
AND tr.gene_id = gsi.gene_id
AND gsi.stable_id = '$geneid';"
);
$sth->execute();
my ($start,$end,$chr);
my @start;
while ( my @row=$sth->fetchrow_array){
($start,$end,$chr)=@row;
push @start,$start;
push @start,$end;
}
my @start_sorted=sort { $a <=> $b } @start;
$start=shift @start_sorted;
$end=pop @start_sorted;
return ($chr,$start,$end);
}
=head2 fetch_Slice_by_gene
Title : fetch_Slice_by_gene
Usage : $slice = $slice_adaptor->fetch_Slice_by_gene('ENSG00000012123',1000);
Function: Creates a slice of the specified object. If a context size is given, the slice is extended by that number of basepairs on either side of the gene. Throws if the gene is not golden.
Returns : Slice object
Args : gene id, [context size in bp]
=cut
sub fetch_Slice_by_gene{
my ($self,$geneid,$size) = @_;
if( !defined $geneid ) {
$self->throw("Must have gene id to fetch Slice of gene");
}
if( !defined $size ) {$size=0;}
my ($chr_name,$start,$end) = $self->get_Gene_chr_bp($geneid);
if( !defined $start ) {
my $type = $self->adaptor->db->static_golden_path_type()
or $self->throw("No assembly type defined");
$self->throw("Gene is not on the golden path '$type'. Cannot build Slice.");
}
return $self->fetch_Slice_by_chr_start_end( $chr_name,
$start-$size,
$end+$size
);
}
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment