[ensembl-dev] GenomicAlignBlock restriction

[Benoît Ballester]

Hi Compara/Ensembl, 

I'd like to have confirmations about some of the GenomicAlignBlock API logic. 

Given a human slice (300 to 400bp), I want to fetch a GAB to see if I could find an alignment across species. I want to know if that slice is aligneable in a couple of species. 


I use $gab_adaptor->fetch_all_by_MethodLinkSpeciesSet_Slice to fetch the GenomicAlignBlock. 

But this returns the entire GenomicAlign for the given block ! Right ? 


    #-- Get GAB                                                                                                                                     
    my $slice = $slice_adaptor->fetch_by_region( 'chromosome', $seq_region_name, $seq_region_start, $seq_region_end );                  
                                                                                       
    my $gabs = $gab_adaptor->fetch_all_by_MethodLinkSpeciesSet_Slice($mlss, $slice);                                                       
                                                                                                                                                                                                                                                                                                
    if (scalar @$gabs > 0){ #-- There is a GAB                                                                                                     
                                                                                                                                                    
        foreach my $gab (@$gabs) {                                                                                                                 
          #-- do something                                                                                                                                          
	  # but this will return the entire GenomicAlign
	}
     }


So, correct me if I am wrong, I can restrict the GAB using the method restrict_between_reference_positions giving the start/end coordinates  of my initial slice, and the boolean $skip_empty_GenomicAligns ? 


Once I have the GAB I restrict it : 

         my $gab2 = $gab->restrict_between_reference_positions($rf->seq_region_start, $rf->seq_region_end, undef ,1); 


Could you please confirm that what I am doing is correct. 
I want to do that from a human slice, and find if there is an alignment for a couple of species. I want to skip empty alignments, or long gaps. If this the correct way ? 

Thanks for the replies, 


Ben


--
Benoît Ballester, PhD
Inserm U1090, TAGC
Campus de Luminy
13288 Marseille Cedex 9
France
+33 4 91 82 87 39