view release on metacpan or  search on metacpan

Bio/Search/Tiling/MapTiling.pm  view on Meta::CPAN

    return $self->{"_algorithm"};
}

=head1 "PRIVATE" METHODS

=head2 Calculators

See L<Bio::Search::Tiling::MapTileUtils> for lower level
calculation methods.

=head2 _calc_coverage_map

 Title   : _calc_coverage_map
 Usage   : $tiling->_calc_coverage_map($type)
 Function: Calculates the coverage map for the object's associated
           hit from the perspective of the desired $type (see Args:) 
           and sets the coverage_map() property
 Returns : True on success
 Args    : optional scalar $type: one of 'hit'|'subject'|'query'
           default is 'query'
 Note    : The "coverage map" is an array with the following format:
           ( [ $component_interval => [ @containing_hsps ] ], ... ),
           where $component_interval is a closed interval (see 
           DESCRIPTION) of the form [$a0, $a1] with $a0 <= $a1, and
           @containing_hsps is an array of all HspI objects in the hit 
           which completely contain the $component_interval.
           The set of $component_interval's is a disjoint decomposition
           of the minimum set of minimal intervals that completely
           cover the hit's HSPs (from the perspective of the $type)
 Note    : This calculates the map for all strand/frame contexts available
           in the hit

=cut

sub _calc_coverage_map {
    my $self = shift;
    my ($type) = @_;
    $self->_check_type_arg(\$type);

    # obtain the [start, end] intervals for all hsps in the hit (relative
    # to the type)
    unless ($self->{'hsps'}) {
	$self->warn("No HSPs for this hit");
	return;
    }

    my (@map, @hsps, %filters, @intervals);
    

    # conversion here?
    my $c = $self->mapping($type);
    
    # create the possible maps 
    for my $context ($self->contexts($type)) {
	@map = ();
	@hsps = ($self->hsps)[$self->contexts($type, $context)];
	@intervals = get_intervals_from_hsps( $type, @hsps );
	# the "frame"
	my $f = ($intervals[0]->[0] - 1) % $c;

	# convert interval endpoints...
	for (@intervals) {
	    $$_[0] = ($$_[0] - $f + $c - 1)/$c;
	    $$_[1]  = ($$_[1] - $f)/$c;
	}
	
	# determine the minimal set of disjoint intervals that cover the
	# set of hsp intervals
	my @dj_set = interval_tiling(\@intervals);

	# decompose each disjoint interval into another set of disjoint 
	# intervals, each of which is completely contained within the
	# original hsp intervals with which it overlaps
	my $i=0;
	my @decomp;
	for my $dj_elt (@dj_set) {
	    my ($covering, $indices) = @$dj_elt;
	    my @covering_hsps = @hsps[@$indices];
	    my @coverers = @intervals[@$indices];
	    @decomp = decompose_interval( \@coverers );
	    for (@decomp) {
		my ($component, $container_indices) = @{$_};
		push @map, [ $component, 
			     [@covering_hsps[@$container_indices]] ];
	    }
	    1;
	}
    
	# unconvert the components:
#####
	foreach (@map) {
	    $$_[0][0] = $c*$$_[0][0] - $c + 1 + $f;
	    $$_[0][1] = $c*$$_[0][1] + $f;
	}
	foreach (@dj_set) {
	    $$_[0][0] = $c*$$_[0][0] - $c + 1 + $f;
	    $$_[0][1] = $c*$$_[0][1] + $f;
	}	    

	# sort the map on the interval left-ends
	@map = sort { $a->[0][0]<=>$b->[0][0] } @map;
	$self->{"coverage_map_${type}_${context}"} = [@map];
	# set the _contig_intersection attribute here (side effect)
	$self->{"_contig_intersection_${type}_${context}"} = [map { $$_[0] } @map];
    }

    return 1; # success
}

=head2 _calc_stats

 Title   : _calc_stats
 Usage   : $tiling->_calc_stats($type, $action, $context)
 Function: Calculates [estimated] tiling statistics (identities, conserved sites
           length) and sets the public accessors
 Returns : True on success
 Args    : scalar $type: one of 'hit', 'subject', 'query'
           default is 'query'
           optional scalar $action: requests calculation method
            currently one of 'exact', 'est', 'fast', 'max'
           option scalar $context: strand/frame context string