t/SearchIO/fasta.t

   1 # -*-Perl-*- Test Harness script for Bioperl
   2 # $Id: SearchIO_fasta.t 14995 2008-11-16 06:20:00Z cjfields $
   3
   4 use strict;
   5
   6 BEGIN {
   7         use lib '.';
   8     use Bio::Root::Test;
   9
  10     test_begin(-tests => 301);
  11
  12         use_ok('Bio::SearchIO');
  13 }
  14
  15 my ($searchio, $result, $hit, $hsp);
  16
  17 # simple FASTA
  18 $searchio = Bio::SearchIO->new(-format => 'fasta',
  19                                  -file   => test_input_file('HUMBETGLOA.FASTA'));
  20 $result = $searchio->next_result;
  21 like($result->database_name, qr/dros_clones.2.5/);
  22 is($result->database_letters, 112936249);
  23 is($result->database_entries, 657);
  24 is($result->algorithm, 'FASTN');
  25 is($result->algorithm_version, '3.3t08');
  26 is($result->query_name, "HUMBETGLOA");
  27 is($result->query_description, "Human haplotype C4 beta-globin gene, complete cds.");
  28 is($result->query_length, 3002);
  29 is($result->get_parameter('gapopen'), -16);
  30 is($result->get_parameter('gapext'), -4);
  31 is($result->get_parameter('ktup'), 6);
  32
  33 is($result->get_statistic('lambda'), 0.0823);
  34 is($result->get_statistic('dbletters'), 112936249);
  35 is($result->get_statistic('dbentries'), 657);
  36
  37 my @valid = ( [ 'BACR21I23', 73982, 'BACR21I23', '0.017', 44.2],
  38            [ 'BACR40P19', 73982, 'BACR40P19', '0.017', 44.2],
  39            [ 'BACR30L17', 32481, 'BACR30L17', '0.018', 44.1]);
  40 my $count = 0;
  41
  42 while( my $hit = $result->next_hit ) {
  43     my $d = shift @valid;
  44     is($hit->name, shift @$d);
  45     is($hit->length, shift @$d);
  46     is($hit->accession, shift @$d);
  47     float_is($hit->significance, shift @$d);
  48     is($hit->raw_score, shift @$d );
  49     is($hit->rank, $count + 1);
  50     if( $count == 0 ) {
  51         my $hsps_left = 1;
  52         while( my $hsp = $hit->next_hsp ) {
  53             is($hsp->query->start, 31);
  54             is($hsp->query->end, 289);
  55             is($hsp->query->strand, -1);
  56             is($hsp->hit->end, 65167);
  57             is($hsp->hit->start, 64902);
  58             is($hsp->hit->strand, 1);
  59             is($hsp->length('total'), 267);
  60             float_is($hsp->evalue, 0.017);
  61             is($hsp->score, 134.5);
  62             is($hsp->bits,44.2);
  63             is(sprintf("%.2f",$hsp->percent_identity), '57.30');
  64             is(sprintf("%.4f",$hsp->frac_identical('query')), 0.5907);
  65             is(sprintf("%.4f",$hsp->frac_identical('hit')), 0.5752);
  66                         # these are really UNGAPPED values not CONSERVED
  67                         # otherwise ident and conserved would be identical for
  68                         # nucleotide alignments
  69                         is(sprintf("%.4f",$hsp->frac_conserved('total')), '0.5918');
  70                         is(sprintf("%.4f",$hsp->frac_conserved('query')), '0.6100');
  71                         is(sprintf("%.4f",$hsp->frac_conserved('hit')), '0.5940');
  72             is($hsp->query->frame(), 0);
  73             is($hsp->hit->frame(), 0);
  74             is($hsp->gaps('query'), 8);
  75             is($hsp->gaps('hit'),1);
  76             is($hsp->gaps('total'), 9);
  77             is($hsp->query_string, 'GATTAAAACCTTCTGGTAAGAAAAGAAAAAATATATATATATATATATGTGTATATGTACACACATACATATACATATATATGCATTCATTTGTTGTTGTTTTTCTTAATTTGCTCATGCATGCTA----ATAAATTATGTCTAAAAATAGAAT---AAATACAAATCAATGTGCTCTGTGCATTA-GTTACTTATTAGGTTTTGGGAAACAAGAGGTAAAAAACTAGAGACCTCTTAATGCAGTCAAAAATACAAATAAATAAAAAGTCACTTACAACCCAAAGTGTGACTATCAATGGGGTAATCAGTGGTGTCAAATAGGAGGT');
  78             is($hsp->hit_string, 'GATGTCCTTGGTGGATTATGGTGTTAGGGTATATATATATATATATATATATATATATATATATATATATATATATATATATATATATATATATATATATATATATATATAATATAATACAAAATATAATACAAAATATAATACAAAATATAATACAAAATATAATACAAAATATAATACAAAATATAATACAAAATATAATATAAAATATAATATAAAATATAATATAAAATAAAATATAAAATAAAATATAAAATAAAATATAAAATAAAATATAAAATAAAATAT-AATATAAAATATAAAATAAAATATAATATAAAATATAATATAAAATATAATATAAAATATAATATAAAATA');
  79             is($hsp->homology_string, '                              :::::::::::::::::: : ::::: :: : : ::: ::::: ::::::::  ::  :: : :   : : : : :  ::    : :: ::   ::    : ::: :::     :::::: :::   ::::: ::  :::  :    :    : ::   :::  : ::   : :   : : :: :   :: : : :: : :       ::  : : ::: ::: ::  ::::: ::: : :  :: ::   ::: : : : ::: ::   '.' 'x60);
  80             is(join(' ', $hsp->seq_inds('query', 'nomatch',1)), '33 37 39 41 43 47-49 52 55 56 58 60 64 70 71 74 78 82 84 86 87 90-96 98 100 103 105 107 110-112 114 117 119 121-123 125 127-129 132 134 135 139 140 143 145-148 150-153 155 156 160 161 164 173 180-184 188 192 196-198 201 204 206-209 212 213 215 217 219 221 223-225 227 229 232 233 236 237 246 252 256 258 260 263 269 271');
  81             is(join(' ', $hsp->seq_inds('query', 'mismatch',1)), '33 37 41 43 47-49 52 55 56 58 60 64 70 71 74 78 82 84 86 87 90-96 98 100 103 105 107 110-112 114 117 119 121-123 125 127-129 132 134 135 139 140 143 145-148 150-153 155 156 160 161 164 173 180-184 188 192 196-198 201 204 206-209 212 213 215 217 219 221 223-225 227 229 232 233 236 237 246 252 256 258 260 263 269 271');
  82             is(join(' ', $hsp->seq_inds('query', 'conserved',1)), '31 32 34-36 38 40 42 44-46 50 51 53 54 57 59 61-63 65-69 72 73 75-77 79-81 83 85 88 89 97 99 101 102 104 106 108 109 113 115 116 118 120 124 126 130 131 133 136-138 141 142 144 149 154 157-159 162 163 165-172 174-179 185-187 189-191 193-195 199 200 202 203 205 210 211 214 216 218 220 222 226 228 230 231 234 235 238-245 247-251 253-255 257 259 261 262 264-268 270 272-289');
  83             is(join(' ', $hsp->seq_inds('hit', 'nomatch',1)), '64920 64922 64928 64931 64933 64935 64939 64945 64954 64955 64958 64959 64962 64964 64966-64968 64970 64972 64974 64976 64978 64979 64982-64985 64987 64990 64993-64995 64998-65001 65003 65007 65011-65015 65022 65026-65028 65034 65037 65038 65042 65043 65045-65048 65050-65053 65055 65058-65060 65064 65065 65067 65070-65072 65074 65076-65078 65080 65082 65085 65087-65089 65092 65094 65096 65099 65101 65103-65109 65112 65113 65115 65117 65121 65125 65128 65129 65135 65139 65141 65143 65144 65147 65150-65152 65156 65158 65161 65165');
  84             is(join(' ', $hsp->seq_inds('hit', 'mismatch',1)), '64920 64922 64928 64931 64933 64935 64939 64945 64954 64955 64958 64959 64962 64964 64966-64968 64970 64972 64974 64976 64978 64979 64982-64985 64987 64990 64993-64995 65003 65007 65011-65015 65022 65034 65037 65038 65042 65043 65045-65048 65050-65053 65055 65059 65060 65064 65065 65067 65070-65072 65074 65076-65078 65080 65082 65085 65087-65089 65092 65094 65096 65099 65101 65103-65109 65112 65113 65115 65117 65121 65125 65128 65129 65135 65139 65141 65143 65144 65147 65150-65152 65156 65158 65161 65165');
  85             is(join(' ', $hsp->seq_inds('hit', 'conserved',1)), '64902-64919 64921 64923-64927 64929 64930 64932 64934 64936-64938 64940-64944 64946-64953 64956 64957 64960 64961 64963 64965 64969 64971 64973 64975 64977 64980 64981 64986 64988 64989 64991 64992 64996 64997 65002 65004-65006 65008-65010 65016-65021 65023-65025 65029-65033 65035 65036 65039-65041 65044 65049 65054 65056 65057 65061-65063 65066 65068 65069 65073 65075 65079 65081 65083 65084 65086 65090 65091 65093 65095 65097 65098 65100 65102 65110 65111 65114 65116 65118-65120 65122-65124 65126 65127 65130-65134 65136-65138 65140 65142 65145 65146 65148 65149 65153-65155 65157 65159 65160 65162-65164 65166 65167');
  86             is(join(' ', $hsp->seq_inds('query', 'gap',1)), '141 170 194');
  87             is(join(' ', $hsp->seq_inds('hit', 'frameshift')), '');
  88             is(join(' ', $hsp->seq_inds('query', 'frameshift')), '');
  89             is($hsp->ambiguous_seq_inds, '');
  90             # note: the reason this is not the same percent id above
  91             # is we are calculating average percent id
  92             is(sprintf("%.2f",$hsp->get_aln->percentage_identity()), '59.30');
  93             $hsps_left--;
  94         }
  95         is($hsps_left, 0);
  96     }
  97     last if( $count++ > @valid );
  98 }
  99 is(@valid, 0);
 100
 101 $searchio = Bio::SearchIO->new(-format => 'fasta',
 102                                  -file   => test_input_file('cysprot1.FASTA'));
 103 $result = $searchio->next_result;
 104 like($result->database_name, qr/ecoli.aa/);
 105 is($result->database_letters, 1358987);
 106 is($result->database_entries, 4289);
 107 is($result->algorithm, 'FASTP');
 108 is($result->algorithm_version, '3.3t08');
 109 is($result->query_name, 'CYS1_DICDI');
 110 is($result->query_length, 343);
 111 is($result->get_parameter('gapopen'), -12);
 112 is($result->get_parameter('gapext'), -2);
 113 is($result->get_parameter('ktup'), 2);
 114
 115 is($result->get_statistic('lambda'), 0.1456);
 116 is($result->get_statistic('dbletters'), 1358987);
 117 is($result->get_statistic('dbentries'), 4289);
 118
 119
 120 @valid = ( [ 'gi|1787478|gb|AAC74309.1|', 512, 'AAC74309', 1787478, 1.2, 29.2],
 121            [ 'gi|1790635|gb|AAC77148.1|', 251, 'AAC77148', 1790635, 2.1, 27.4],
 122            [ 'gi|1786590|gb|AAC73494.1|', 94, 'AAC73494', 1786590, 2.1, 25.9]);
 123 $count = 0;
 124
 125 while( my $hit = $result->next_hit ) {
 126     my $d = shift @valid;
 127     is($hit->name, shift @$d);
 128     is($hit->length, shift @$d);
 129     is($hit->accession, shift @$d);
 130     is($hit->ncbi_gi, shift @$d);
 131     float_is($hit->significance, shift @$d);
 132     is($hit->raw_score, shift @$d );
 133
 134     if( $count == 0 ) {
 135         my $hsps_left = 1;
 136         while( my $hsp = $hit->next_hsp ) {
 137             is($hsp->query->start, 125);
 138             is($hsp->query->end, 305);
 139             is($hsp->query->strand, 0);
 140             is($hsp->hit->start, 2);
 141             is($hsp->hit->end, 181);
 142             is($hsp->hit->strand, 0);
 143             is($hsp->length('total'), 188);
 144             float_is($hsp->evalue, 1.2);
 145             is($hsp->score, 109.2);
 146             is($hsp->bits,29.2);
 147             is(sprintf("%.2f",$hsp->percent_identity), 23.94);
 148             is(sprintf("%.4f",$hsp->frac_identical('query')), 0.2486);
 149             is(sprintf("%.4f",$hsp->frac_identical('hit')), '0.2500');
 150             is(sprintf("%.4f",$hsp->frac_conserved('query')), '0.2707');
 151             is(sprintf("%.4f",$hsp->frac_conserved('hit')), '0.2722');
 152             # there is slight rounding different here so file says 26.012%
 153             # but with the rounding this ends up as 0.2606
 154             is(sprintf("%.4f",$hsp->frac_conserved('total')), '0.2606');
 155             is($hsp->query->frame(), 0);
 156             is($hsp->hit->frame(), 0);
 157             is($hsp->gaps('query'), 7);
 158             is($hsp->gaps('hit'), 8);
 159             is($hsp->gaps, 15);
 160             is($hsp->query_string, 'NKEAIFTDDLPVADYLDDEFINSIPTAFDWRTRGAVTPVKNQGQCGSCWSFSTT-GNV----EGQHFISQNKLVSLSEQNLVDCDHECME-YEGEEACDEGCNGGLQPNAYNYIIKNGGIQTESSYPYTAETGTQCNFNSANIGAKISNFTMIPKNETVMAGYIVSTGP-LAIAADAVEWQFYIGGVFDIPCNPNSLDHGILIVGYSAKNTIFRKNMPYWIVKNSWGADWGEQGYIYLRRGKNTCGVSNFVSTSII');
 161             is($hsp->hit_string, (' 'x29).'MKIRSQVGMVLNLDKCIGCHTCSVTCKNVWTSREGVEYAWFNNVETKPGQGF-PTDWENQEKYKGGWI--RKINGKLQPRMGNRAMLLGKIFANPHLPGIDDYYEPFDFDYQNLHTAPEG----SKSQPIARPRSLITGERMAKIEKGPNWEDDLGGEFDKLAKDKNFDN-IQKAMYSQFENTFMMYLPRLCEHCLNPACVATCPSGAIYKREEDGIVLIDQDKCRGWRMCITGCPYKKIYFNWKSGKSEKCIFCYPRIEAGQPTVCSETC');
 162             is($hsp->homology_string, '                              . :. :  : :  .: .: . :.:  ::    :: ..   :.. .   :..   : : .: :.:     .  :: :::   :  .  : : ..   :   .     .:.  :. .   .     :.. .     . ::  .:    . .:.  .:: ::   . ...:. :  . ::  .. :   .:                      '.' 'x60);
 163             # note: the reason this is not the same percent id above
 164             # is we are calculating average percent id
 165             is(sprintf("%.2f",$hsp->get_aln->percentage_identity()), 26.01);
 166             is(join(' ', $hsp->seq_inds('query', 'nomatch',1)), '126 129 131 132 134 136 137 140 143 145 149 154 157-159 163 165-167 171-173 175 177 183-187 189 190 193 197-199 201 202 204 205 207 209 212-214 216-218 220-224 228 229 232 234-236 238-242 246 248-252 254 257 260-263 265 269 270 274 277-279 281 287 289 290 292 295 296 299 301-303');
 167             is(join(' ', $hsp->seq_inds('query', 'mismatch',1)), '126 129 131 132 134 136 137 140 143 145 149 154 157-159 163 165-167 172 173 175 177 183-185 189 190 193 197-199 201 202 204 205 207 209 212-214 216-218 220-224 228 229 232 234-236 242 246 248-252 254 257 260-263 265 269 270 274 277-279 281 289 290 292 295 296 299 301-303');
 168             is(join(' ', $hsp->seq_inds('query', 'conserved',1)), '125 127 128 130 133 135 138 139 141 142 144 146-148 150-153 155 156 160-162 164 168-170 174 176 178-182 188 191 192 194-196 200 203 206 208 210 211 215 219 225-227 230 231 233 237 243-245 247 253 255 256 258 259 264 266-268 271-273 275 276 280 282-286 288 291 293 294 297 298 300 304 305');
 169             is(join(' ', $hsp->seq_inds('hit', 'nomatch',1)), '3 6 8 9 11 13 14 17 20 22 26 27 30-33 36 39-41 45 47-49 53 54 56 58 61 65-67 69 70 73 77-79 81 82 84 85 87 89 92-94 96-98 100-104 108 109 112 114-116 118 122 124-128 130 133 134 137-140 142 146 147 151 154-156 158 165 166 168 171 172 175 177-179');
 170             is(join(' ', $hsp->seq_inds('hit', 'mismatch',1)), '3 6 8 9 11 13 14 17 20 22 27 36 39-41 45 47-49 53 54 56 58 65-67 69 70 73 77-79 81 82 84 85 87 89 92-94 96-98 100-104 108 109 112 114-116 118 122 124-128 130 133 137-140 142 146 147 151 154-156 158 165 166 168 171 172 175 177-179');
 171             is(join(' ', $hsp->seq_inds('hit', 'conserved',1)), '2 4 5 7 10 12 15 16 18 19 21 23-25 28 29 34 35 37 38 42-44 46 50-52 55 57 59 60 62-64 68 71 72 74-76 80 83 86 88 90 91 95 99 105-107 110 111 113 117 119-121 123 129 131 132 135 136 141 143-145 148-150 152 153 157 159-164 167 169 170 173 174 176 180 181');
 172             is(join(' ', $hsp->seq_inds('query', 'gap',1)), '148 151 179 257');
 173             is(join(' ', $hsp->seq_inds('hit', 'frameshift')), '');
 174             is(join(' ', $hsp->seq_inds('query', 'frameshift')), '');
 175             $hsps_left--;
 176         }
 177         is($hsps_left, 0);
 178     }
 179     last if( $count++ > @valid );
 180 }
 181 is(@valid, 0);
 182
 183 is($result->hits, 8);
 184 $searchio = Bio::SearchIO->new(-format => 'fasta',
 185                                  -file   => test_input_file('cysprot_vs_gadfly.FASTA'));
 186 $result = $searchio->next_result;
 187 like($result->database_name, qr/gadflypep2/);
 188 is($result->database_letters, 7177762);
 189 is($result->database_entries, 14334);
 190 is($result->algorithm, 'FASTP');
 191 is($result->algorithm_version, '3.3t08');
 192 is($result->query_name, 'cysprot.fa');
 193 is($result->query_length, 2385);
 194 is($result->get_parameter('gapopen'), -12);
 195 is($result->get_parameter('gapext'), -2);
 196 is($result->get_parameter('ktup'), 2);
 197 is($result->get_parameter('matrix'), 'BL50');
 198
 199 is($result->get_statistic('lambda'), 0.1397);
 200 is($result->get_statistic('dbletters'), 7177762 );
 201 is($result->get_statistic('dbentries'), 14334);
 202
 203
 204 @valid = ( [ 'Cp1|FBgn0013770|pp-CT20780|FBan0006692', 341,
 205              'FBan0006692', '3.1e-59', 227.8],
 206            [ 'CG11459|FBgn0037396|pp-CT28891|FBan0011459', 336,
 207              'FBan0011459', '6.4e-41',  166.9],
 208            [ 'CG4847|FBgn0034229|pp-CT15577|FBan0004847', 390,
 209              'FBan0004847',  '2.5e-40', 165.2]);
 210 $count = 0;
 211
 212 while( my $hit = $result->next_hit ) {
 213     my $d = shift @valid;
 214
 215     is($hit->name, shift @$d);
 216     is($hit->length, shift @$d);
 217     is($hit->accession, shift @$d);
 218     float_is($hit->significance, shift @$d);
 219     is($hit->raw_score, shift @$d );
 220
 221     if( $count == 0 ) {
 222         my $hsps_left = 1;
 223         while( my $hsp = $hit->next_hsp ) {
 224             is($hsp->query->start, 1373);
 225             is($hsp->query->end, 1706);
 226             is($hsp->query->strand, 0);
 227             is($hsp->hit->start, 5);
 228             is($hsp->hit->end, 341);
 229             is($hsp->hit->strand, 0);
 230             is($hsp->length('total'), 345);
 231             float_is($hsp->evalue, 3.1e-59);
 232             is($hsp->score, 1170.6);
 233             is($hsp->bits,227.8);
 234             is(sprintf("%.2f",$hsp->percent_identity), 53.04);
 235             is(sprintf("%.4f",$hsp->frac_identical('query')), 0.5479);
 236             is(sprintf("%.4f",$hsp->frac_identical('hit')), '0.5430');
 237             is($hsp->query->frame(), 0);
 238             is($hsp->hit->frame(), 0);
 239             is($hsp->gaps('query'), 11);
 240             is($hsp->gaps('hit'), 8);
 241             is($hsp->gaps, 19);
 242             is($hsp->hit_string, (' 'x26).'MRTAVLLPLLAL----LAVAQA-VSFADVVMEEWHTFKLEHRKNYQDETEERFRLKIFNENKHKIAKHNQRFAEGKVSFKLAVNKYADLLHHEFRQLMNGFNYTLHKQLRAADESFKGVTFISPAHVTLPKSVDWRTKGAVTAVKDQGHCGSCWAFSSTGALEGQHFRKSGVLVSLSEQNLVDCSTKYGNNGCNGGLMDNAFRYIKDNGGIDTEKSYPYEAIDDSCHFNKGTVGATDRGFTDIPQGDEKKMAEAVATVGPVSVAIDASHESFQFYSEGVYNEPQCDAQNLDHGVLVVGFGTDESGED---YWLVKNSWGTTWGDKGFIKMLRNKENQCGIASASSYPLV');
 243             is($hsp->query_string, 'SNWGNNGYFLIERGKNMCGLAACASYPIPQVMNPTLILAAFCLGIASATLTFDHSLEAQWTKWKAMHNRLY-GMNEEGWRRAVWEKNMKMIELHNQEYREGKHSFTMAMNAFGDMTSEEFRQVMNGFQ---NRKPR------KGKVFQEPLFYEAPRSVDWREKGYVTPVKNQGQCGSCWAFSATGALEGQMFRKTGRLISLSEQNLVDCSGPQGNEGCNGGLMDYAFQYVQDNGGLDSEESYPYEATEESCKYNPKYSVANDTGFVDIPK-QEKALMKAVATVGPISVAIDAGHESFLFYKEGIYFEPDCSSEDMDHGVLVVGYGFESTESDNNKYWLVKNSWGEEWGMGGYVKMAKDRRNHCGIASAASYPTVMTPLLLLAVLCLGTALATPKFDQTFNAQWHQWKSTHRRLYGTNEE');
 244             # note: the reason this is not the same percent id above
 245             # is we are calculating average percent id
 246             is(sprintf("%.2f",$hsp->get_aln->percentage_identity()), 56.13);
 247             is(join(' ', $hsp->seq_inds('query', 'nomatch',1)), '1375 1377 1379 1381-1384 1391 1395-1397 1399 1400 1403 1404 1407 1408 1410 1412 1414 1415 1419 1422 1423 1429 1431 1433 1434 1440 1444 1447 1452 1457 1458 1473 1477 1480 1481 1483-1487 1495 1498 1501 1524 1530 1544-1546 1558 1580 1588-1592 1596 1607 1609 1630 1638 1658 1660 1662 1663 1665-1667 1677 1678 1681 1682 1688 1705');
 248             is(join(' ', $hsp->seq_inds('query', 'mismatch',1)), '1375 1377 1379 1395-1397 1399 1400 1403 1404 1407 1408 1410 1412 1414 1415 1419 1422 1423 1429 1431 1433 1434 1440 1444 1447 1452 1457 1458 1473 1477 1480 1481 1483-1487 1495 1498 1501 1524 1530 1544-1546 1558 1580 1588-1592 1596 1607 1609 1630 1638 1658 1660 1662 1663 1677 1678 1681 1682 1688 1705');
 249             is(join(' ', $hsp->seq_inds('query', 'conserved',1)), '1373 1374 1376 1378 1380 1385-1390 1392-1394 1398 1401 1402 1405 1406 1409 1411 1413 1416-1418 1420 1421 1424-1428 1430 1432 1435-1439 1441-1443 1445 1446 1448-1451 1453-1456 1459-1472 1474-1476 1478 1479 1482 1488-1494 1496 1497 1499 1500 1502-1523 1525-1529 1531-1543 1547-1557 1559-1579 1581-1587 1593-1595 1597-1606 1608 1610-1629 1631-1637 1639-1657 1659 1661 1664 1668-1676 1679 1680 1683-1687 1689-1704 1706');
 250             is(join(' ', $hsp->seq_inds('hit', 'nomatch',1)), '7 9 11 22-24 26 27 30 31 34 35 37 39 41-43 47 50 51 57 59 61 62 68 72 75 80 85 86 98-100 104 106-111 114 117 118 120-124 132 135 138 161 167 181-183 195 217 225-229 233 241 245 247 268 276 296 298 300 301 312 313 316 317 323 340');
 251             is(join(' ', $hsp->seq_inds('hit', 'mismatch',1)), '7 9 11 22-24 26 27 30 31 34 35 37 39 42 43 47 50 51 57 59 61 62 68 72 75 80 85 86 104 114 117 118 120-124 132 135 138 161 167 181-183 195 217 225-229 233 245 247 268 276 296 298 300 301 312 313 316 317 323 340');
 252             is(join(' ', $hsp->seq_inds('hit', 'conserved',1)), '5 6 8 10 12-21 25 28 29 32 33 36 38 40 44-46 48 49 52-56 58 60 63-67 69-71 73 74 76-79 81-84 87-97 101-103 105 112 113 115 116 119 125-131 133 134 136 137 139-160 162-166 168-180 184-194 196-216 218-224 230-232 234-240 242-244 246 248-267 269-275 277-295 297 299 302-311 314 315 318-322 324-339 341');
 253             is(join(' ', $hsp->seq_inds('query', 'gap',1)), '1413 1469 1474 1603');
 254             is(join(' ', $hsp->seq_inds('hit', 'frameshift')), '');
 255             is(join(' ', $hsp->seq_inds('query', 'frameshift')), '');
 256             $hsps_left--;
 257         }
 258         is($hsps_left, 0);
 259     }
 260     last if( $count++ > @valid );
 261 }
 262 is(@valid, 0);
 263 is($result->hits, 21);
 264
 265 # test on TFASTXY
 266 $searchio = Bio::SearchIO->new(-format => 'fasta',
 267                               -file   => test_input_file('5X_1895.FASTXY'));
 268 $result = $searchio->next_result;
 269 like($result->database_name, qr/yeast_nrpep.fasta/);
 270 is($result->database_letters, 4215311);
 271 is($result->database_entries, 9190);
 272 is($result->algorithm, 'FASTY');
 273 is($result->algorithm_version, '3.4t07');
 274 is($result->query_name, '5X_1895.fa');
 275 is($result->query_length, 7972);
 276 is($result->get_parameter('gapopen'), -14);
 277 is($result->get_parameter('gapext'), -2);
 278 is($result->get_parameter('ktup'), 2);
 279 is($result->get_parameter('matrix'), 'BL50');
 280
 281 is($result->get_statistic('lambda'), 0.1711);
 282 is($result->get_statistic('dbletters'), 4215311);
 283 is($result->get_statistic('dbentries'), 9190);
 284
 285
 286 @valid = ( [ 'NR_SC:SW-YNN2_YEAST', 1056, 'NR_SC:SW-YNN2_YEAST','1.6e-154', '547.0'],
 287            [ 'NR_SC:SW-MPCP_YEAST', 311, 'NR_SC:SW-MPCP_YEAST', '1.3e-25', 117.1],
 288            [ 'NR_SC:SW-YEO3_YEAST', 300, 'NR_SC:SW-YEO3_YEAST', '5.7e-05', 48.5]);
 289 $count = 0;
 290
 291 while( my $hit = $result->next_hit ) {
 292     my $d = shift @valid;
 293
 294     is($hit->name, shift @$d);
 295     is($hit->length, shift @$d);
 296     is($hit->accession, shift @$d);
 297     float_is($hit->significance, shift @$d);
 298     is($hit->raw_score, shift @$d );
 299
 300     if ( $count == 1 ) {
 301         my $hsps_left = 1;
 302         while( my $hsp = $hit->next_hsp ) {
 303             is($hsp->query->start, 6947);
 304             is($hsp->query->end, 7938);
 305             is($hsp->query->strand, -1);
 306             is($hsp->hit->start, 19);
 307             is($hsp->hit->end, 287);
 308             is($hsp->hit->strand, 0);
 309             is($hsp->length('total'), 330);
 310
 311             float_is($hsp->evalue, 1.3e-25);
 312             is($hsp->score, 563.4);
 313             is($hsp->bits,'117.1');
 314             is(sprintf("%.2f",$hsp->percent_identity), 54.24);
 315             is(sprintf("%.4f",$hsp->frac_identical('query')), '0.5413');
 316             is(sprintf("%.4f",$hsp->frac_identical('hit')), '0.6654');
 317             is($hsp->query->frame(), 0);
 318             is($hsp->hit->frame(), 0);
 319             is($hsp->gaps('query'), 0);
 320             is($hsp->gaps('hit'), 61);
 321             is($hsp->gaps, 61);
 322             is($hsp->query_string, 'RFALAGALGCAVTHGALTPVDVVKTRIQLEPEVYNRVGRFFNSS*GF*EL*GVVLMSQT\KGMVASFRQIIAKEGAGALLTGFGPTAVGYAIQGAFKFGG*VMMSLQITA*SRANLPISYEFWKKKAIDLVGVDKARENRQAIYLGASAIAEFFADIALCPLEATRIRLVSQPSFANGLSGGFLRILREEGPAAFYAGFGPILFKQVPYTMAKFAV*VDRTA*\YQTFG*YYRSYEVAVEKILKATGKSKDSLTGGQLTGLNLTSGLIAGLAAAVISQPADTLLSKINKTKGAPGQSTTSRLVQMAGQLGVSGLFTGMTTRLVMIGTLTAGQ');
 323             is($hsp->hit_string, 'KFALAGAIGCGSTHSSMVPIDVVKTRIQLEPTVYN-------------------------KGMVGSFKQIIAGEGAGALLTGFGPTLLGYSIQGAFKFGG-------------------YEVFKKFFIDNLGYDTASRYKNSVYMGSAAMAEFLADIALCPLEATRIRLVSQPQFANGLVGGFSRILKEEGIGSFYSGFTPILFKQIPYNIAKFLVFERASEF-YYGFAG------------------PKEKLSSTSTTLLNLLSGLTAGLAAAIVSQPADTLLSKVNKTKKAPGQSTVGLLAQLAKQLGFFGSFAGLPTRLVMVGTLTSLQ');
 324             is($hsp->homology_string, '.::::::.::. ::....:.::::::::::: :::                         ::::.::.:::: ::::::::::::: .::.:::::::::                   :: .::  :: .: : : . ....:.:..:.:::.:::::::::::::::::::.::::: ::: :::.::: ..::.:: ::::::.::..::: :    .   :  :.                    :..:.. . : ::: ::: ::::::..::::::::::.:::: ::::::.. :.:.: :::  : :.:. :::::.::::. :');
 325             # note: the reason this is not the same percent id above
 326             # is we are calculating average percent id
 327             is(sprintf("%.2f",$hsp->get_aln->overall_percentage_identity()),
 328                '53.92');
 329             is(sprintf("%.2f",$hsp->get_aln->average_percentage_identity()),
 330                '67.29');
 331             is(join(' ', $hsp->seq_inds('query', 'nomatch',1)), '6950-6952 6986-6988 7001-7003 7007-7012 7022-7024 7040-7042 7067-7069 7139-7141 7151-7153 7163-7165 7169-7171 7175-7177 7196-7255 7262-7267 7272-7277 7281-7292 7296-7298 7341-7343 7365-7367 7389-7391 7401-7403 7524-7526 7530-7532 7536-7538 7542-7544 7551-7553 7560-7565 7575-7577 7584-7640 7680-7682 7722-7724 7762-7833 7843-7845 7903-7905');
 332             is(join(' ', $hsp->seq_inds('query', 'mismatch',1)), '6950-6952 6986-6988 7001-7003 7007-7012 7022-7024 7040-7042 7067-7069 7139-7141 7151-7153 7163-7165 7169-7171 7175-7177 7196-7198 7253-7255 7262-7267 7272-7277 7281-7292 7296-7298 7341-7343 7365-7367 7389-7391 7401-7403 7524-7526 7530-7532 7536-7538 7542-7544 7551-7553 7560-7565 7575-7577 7680-7682 7722-7724 7843-7845 7903-7905');
 333             is(join(' ', $hsp->seq_inds('query', 'conserved',1)), '6947-6949 6953-6985 6989-7000 7004-7006 7013-7021 7025-7039 7043-7066 7070-7138 7142-7150 7154-7162 7166-7168 7172-7174 7178-7195 7256-7261 7269-7272 7278-7280 7293-7295 7299-7340 7344-7364 7368-7388 7392-7400 7404-7523 7527-7529 7533-7535 7539-7541 7545-7550 7554-7559 7566-7574 7578-7583 7641-7679 7683-7721 7725-7757 7759-7762 7834-7842 7846-7902 7906-7938');
 334             is(join(' ', $hsp->seq_inds('query', 'gap',1)), '');
 335             is(join(' ', $hsp->seq_inds('hit', 'nomatch',1)), '30 50 66 80 96 100 101 104 107 109 111 113 154 158 166 174 189 191-194 196 197 199 200 203 204 211 213 215 219 223 247 256 262 266 267 269 274 286');
 336             is(join(' ', $hsp->seq_inds('hit', 'mismatch',1)), '30 50 66 80 96 100 101 104 107 109 111 113 154 158 166 174 189 191-194 196 197 199 200 203 204 211 213 215 219 223 247 256 262 266 267 269 274 286');
 337             is(join(' ', $hsp->seq_inds('hit', 'conserved',1)), '19-29 31-49 51-65 67-79 81-95 97-99 102 103 105 106 108 110 112 114-153 155-157 159-165 167-173 175-188 190 195 198 201 202 205-210 212 214 216-218 220-222 224-246 248-255 257-261 263-265 268 270-273 275-285 287');
 338             is(join(' ', $hsp->seq_inds('hit', 'gap',1)), '53 93 203');
 339             is(join(' ', $hsp->seq_inds('query', 'frameshift', 1)), '7271-7273 7761-7763');
 340             is(join(' ', $hsp->seq_inds('hit', 'frameshift')), '');
 341             is($hsp->ambiguous_seq_inds, 'query');
 342             $hsps_left--;
 343         }
 344         is($hsps_left, 0);
 345     }
 346     last if( $count++ > @valid );
 347 }
 348 is(@valid, 0);
 349 is($result->hits, 58);
 350
 351 # test FASTA v35.04, params encoding changed
 352 # test on TFASTXY
 353 $searchio = Bio::SearchIO->new(-format => 'fasta',
 354                               -file   => test_input_file('BOSS_DROME.FASTP_v35_04'));
 355 $result = $searchio->next_result;
 356 like($result->database_name, qr/wormpep190/, 'TFASTXY');
 357 is($result->database_letters, 10449259);
 358 is($result->database_entries, 23771);
 359 is($result->algorithm, 'FASTP');
 360 is($result->algorithm_version, '35.04');
 361 is($result->query_name, 'BOSS_DROME');
 362 is($result->query_length, 896);
 363 is($result->get_parameter('gapopen'), -10);
 364 is($result->get_parameter('gapext'), -2);
 365 is($result->get_parameter('ktup'), 2);
 366 is($result->get_parameter('matrix'), 'BL50');
 367
 368 # wordsize is the same as ktup, not opt width, as we used to parse
 369 # is($result->get_parameter('wordsize'), 16);
 370
 371 # this is the range of the scoring matrix, not a filter (which is meant
 372 # to capture whether xS seg filtering used)
 373 # is($result->get_parameter('filter'), '15:-5');
 374
 375 is($result->get_statistic('lambda'), 0.122629);
 376 is($result->get_statistic('dbletters'), 10449259);
 377 is($result->get_statistic('dbentries'), 23771);
 378 is($result->get_statistic('effectivespace'),23771);
 379
 380 # bug 2473 - fasta3.4 parsing with -U option
 381
 382 $searchio = Bio::SearchIO->new(-format => 'fasta',
 383                                                           -file   => test_input_file('bug2473.fasta'));
 384
 385 while(my $res = $searchio->next_result) {
 386     is($res->query_name, 'total:39860_L:12096_-3:12346_0:617_+3:14801');
 387     is($res->query_description, '');
 388     is($res->query_length, 22);
 389     is($res->algorithm, 'FASTN');
 390 }
 391
 392 # bug 2537
 393
 394 $searchio = Bio::SearchIO->new(
 395    -format => 'fasta',
 396    -file   => test_input_file('in.fasta'),
 397 );
 398
 399 while ( my $result = $searchio->next_result() ) {
 400    while( my $hit = $result->next_hit ) {
 401       # process the Bio::Search::Hit::HitI object
 402       while( my $hsp = $hit->next_hsp ) {
 403          # process the Bio::Search::HSP::HSPI object
 404          is($hsp->num_identical, 30, "num_identical()");
 405          is($hsp->num_conserved, 30, "num_conserved()");
 406       }
 407    }
 408 }
 409
 410
 411 # bug 2937 and version 35.04
 412 $searchio = Bio::SearchIO->new(
 413    -format => 'fasta',
 414    -file   => test_input_file('bug2937.fasta'),
 415 );
 416 ok($searchio, 'bug 2937 and FASTA version 3.5');
 417
 418 while ( my $result = $searchio->next_result() ) {
 419    is ($result->algorithm_version, '35.04', 'algorithm version');
 420    is($result->query_name, 'ILTV-miR1', 'query name');
 421    is($result->query_description, '', 'query description');
 422    is($result->query_length, 70, 'query length');
 423    is($result->algorithm, 'FASTN', 'algorithm');
 424
 425    while( my $hit = $result->next_hit ) {
 426       # process the Bio::Search::Hit::HitI object
 427       while( my $hsp = $hit->next_hsp ) {
 428          # process the Bio::Search::HSP::HSPI object
 429          is($hsp->num_identical, 26, "num_identical()");
 430          is($hsp->num_conserved, 26, "num_conserved()");
 431          is($hsp->strand('hit'), '1', 'hsp->strand(hit)');
 432                  is($hsp->hit->strand, '1', 'hsp->hit->strand');
 433                  is($hsp->strand('query'), '-1', 'hsp->strand(query)');
 434                  is($hsp->query->strand, '-1', 'hsp->query->strand');
 435       }
 436    }
 437 }
 438
 439 $searchio = Bio::SearchIO->new( -format => 'fasta',
 440                                 -file => test_input_file('issue255_ssearch.fasta'),);
 441 ok($searchio, 'issue 255');
 442 while (my $result = $searchio->next_result() ) {
 443   while (my $hit = $result->next_hit ) {
 444     is ($hit->next_hsp->score, '424.1', 'issue 255: first hsp score');
 445   }
 446 }