t/SearchIO/hmmer.t

   1 # -*-Perl-*- Test Harness script for Bioperl
   2
   3 use strict;
   4
   5 BEGIN {
   6     use lib '.';
   7     use Bio::Root::Test;
   8
   9     test_begin( -tests => 313 );
  10
  11     use_ok('Bio::SearchIO');
  12 }
  13
  14 my $searchio = Bio::SearchIO->new(
  15     -format => 'hmmer',
  16     -file   => test_input_file('hmmpfam.out')
  17 );
  18 my $result;
  19
  20 while ( $result = $searchio->next_result ) {
  21     is( ref($result),
  22         'Bio::Search::Result::HMMERResult',
  23         'Check for the correct result reference type'
  24     );
  25     is( $result->algorithm,         'HMMPFAM', 'Check algorithm' );
  26     is( $result->algorithm_version, '2.1.1',   'Check algorithm version' );
  27     is( $result->hmm_name,          'pfam',    'Check hmm_name' );
  28     is( $result->sequence_file,
  29         '/home/birney/src/wise2/example/road.pep',
  30         'Check sequence_file'
  31     );
  32     is( $result->query_name,        'roa1_drome', 'Check query_name' );
  33     is( $result->query_description, '',           'Check query_description' );
  34     is( $result->num_hits(),        2,            'Check num_hits' );
  35     my ( $hsp, $hit );
  36
  37     if ( $hit = $result->next_model ) {
  38         is( $hit->name,      'SEED',  'Check hit name' );
  39         is( $hit->raw_score, '146.1', 'Check hit raw_score' );
  40         float_is( $hit->significance, 6.3e-40, 'Check hit significance' );
  41         is( ref($hit), 'Bio::Search::Hit::HMMERHit',
  42             'Check for the correct hit reference type' );
  43         is( $hit->num_hsps, 1, 'Check num_hsps' );
  44
  45         if ( defined( $hsp = $hit->next_domain ) ) {
  46             is( $hsp->hit->start,   1,    'Check for hit hmmfrom value' );
  47             is( $hsp->hit->end,     77,   'Check for hit hmm to value' );
  48             is( $hsp->query->start, 33,   'Check for query alifrom value' );
  49             is( $hsp->query->end,   103,  'Check for query ali to value' );
  50             is( $hsp->score,        71.2, 'Check for hsp score' );
  51             float_is( $hsp->evalue, 2.2e-17, 'Check for hsp c-Evalue' );
  52             is( $hsp->query_string,
  53                 'LFIGGLDYRTTDENLKAHFEKWGNIVDVVVMKD-----PRTKRSRGFGFITYSHSSMIDEAQK--SRpHKIDGRVVEP',
  54                 'Check for query string'
  55             );
  56             is( $hsp->gaps('query'), 7, 'Check for number of gaps in query' );
  57             is( $hsp->hit_string,
  58                 'lfVgNLppdvteedLkdlFskfGpivsikivrDiiekpketgkskGfaFVeFeseedAekAlealnG-kelggrklrv',
  59                 'Check for hit string'
  60             );
  61             is( $hsp->homology_string,
  62                 'lf+g+L + +t+e Lk++F+k G iv++ +++D     + t++s+Gf+F+++  ++  + A +    +++++gr+++ ',
  63                 'Check for homology string'
  64             );
  65             is( length( $hsp->homology_string ),
  66                 length( $hsp->hit_string ),
  67                 'Check if homology string and hit string have an equal lenght'
  68             );
  69             is( length( $hsp->query_string ),
  70                 length( $hsp->homology_string ),
  71                 'Check if query string and homology string have an equal lenght'
  72             );
  73         }
  74     }
  75     if ( defined( $hit = $result->next_model ) ) {
  76         if ( defined( $hsp = $hit->next_domain ) ) {
  77             is( $hsp->hit->start,   1,    'Check for hit hmmfrom value' );
  78             is( $hsp->hit->end,     77,   'Check for hit hmm to value' );
  79             is( $hsp->query->start, 124,  'Check for query alifrom value' );
  80             is( $hsp->query->end,   194,  'Check for query ali to value' );
  81             is( $hsp->score,        75.5, 'Check for hsp score' );
  82             float_is( $hsp->evalue, 1.1e-18, 'Check for hsp c-Evalue' );
  83             is( $hsp->query_string,
  84                 'LFVGALKDDHDEQSIRDYFQHFGNIVDINIVID-----KETGKKRGFAFVEFDDYDPVDKVVL-QKQHQLNGKMVDV',
  85                 'Check for query string'
  86             );
  87             is( $hsp->gaps('query'), 6, 'Check for number of gaps in query' );
  88             is( $hsp->hit_string,
  89                 'lfVgNLppdvteedLkdlFskfGpivsikivrDiiekpketgkskGfaFVeFeseedAekAlealnGkelggrklrv',
  90                 'Check for hit string'
  91             );
  92             is( $hsp->homology_string,
  93                 'lfVg L  d +e+ ++d+F++fG iv+i+iv+D     ketgk +GfaFVeF++++ ++k +     ++l+g+ + v',
  94                 'Check for homology string'
  95             );
  96             is( length( $hsp->homology_string ),
  97                 length( $hsp->hit_string ),
  98                 'Check if homology string and hit string have an equal lenght'
  99             );
 100             is( length( $hsp->query_string ),
 101                 length( $hsp->homology_string ),
 102                 'Check if query string and homology string have an equal lenght'
 103             );
 104         }
 105         last;
 106     }
 107 }
 108 $searchio = Bio::SearchIO->new(
 109     -format => 'hmmer',
 110     -file   => test_input_file('hmmsearch.out')
 111 );
 112 while ( $result = $searchio->next_result ) {
 113     is( ref($result),
 114         'Bio::Search::Result::HMMERResult',
 115         'Check for the correct result reference type'
 116     );
 117     is( $result->algorithm,         'HMMSEARCH',  'Check algorithm' );
 118     is( $result->algorithm_version, '2.0',        'Check algorithm version' );
 119     is( $result->hmm_name,          'HMM [SEED]', 'Check hmm_name' );
 120     is( $result->sequence_file, 'HMM.dbtemp.29591', 'Check sequence_file' );
 121     is( $result->database_name, 'HMM.dbtemp.29591', 'Check database_name' );
 122     is( $result->query_name,    'SEED',             'Check query_name' );
 123     is( $result->query_description, '',   'Check query_description' );
 124     is( $result->num_hits(),        1215, 'Check num_hits' );
 125     my $hit = $result->next_model;
 126     is( $hit->name, 'Q91581', 'Check hit name' );
 127     is( $hit->description,
 128         'Q91581 POLYADENYLATION FACTOR 64 KDA SUBUN',
 129         'Check for hit description'
 130     );
 131     float_is( $hit->significance, 2e-31, 'Check hit significance' );
 132     is( $hit->raw_score, 119.7, 'Check hit raw_score' );
 133     my $hsp = $hit->next_domain;
 134     is( $hsp->score, 119.7, 'Check for hsp score' );
 135     float_is( $hsp->evalue, 2e-31, 'Check for hsp c-Evalue' );
 136     is( $hsp->query->start,    18,       'Check for query alifrom value' );
 137     is( $hsp->query->end,      89,       'Check for query ali to value' );
 138     is( $hsp->hit->start,      1,        'Check for hit hmmfrom value' );
 139     is( $hsp->hit->end,        77,       'Check for hit hmm to value' );
 140     is( $hsp->query->seq_id(), 'SEED',   'Check for query seq_id' );
 141     is( $hsp->hit->seq_id(),   'Q91581', 'Check for hit seq_id' );
 142 }
 143
 144 $searchio = Bio::SearchIO->new(
 145     -format => 'hmmer',
 146     -file   => test_input_file('L77119.hmmer')
 147 );
 148
 149 while ( $result = $searchio->next_result ) {
 150     is( ref($result),
 151         'Bio::Search::Result::HMMERResult',
 152         'Check for the correct result reference type'
 153     );
 154     is( $result->algorithm,         'HMMPFAM',    'Check algorithm' );
 155     is( $result->algorithm_version, '2.2g',       'Check algorithm version' );
 156     is( $result->hmm_name,          'Pfam',       'Check hmm_name' );
 157     is( $result->sequence_file,     'L77119.faa', 'Check sequence_file' );
 158     is( $result->query_name, 'gi|1522636|gb|AAC37060.1|',
 159         'Check query_name' );
 160     is( $result->query_description,
 161         'M. jannaschii predicted coding region MJECS02 [Methanococcus jannaschii]',
 162         'Check query_description'
 163     );
 164     is( $result->num_hits(), 1, 'Check num_hits' );
 165     my $hit = $result->next_hit;
 166     is( $hit->name, 'Methylase_M', 'Check hit name' );
 167     is( $hit->description,
 168         'Type I restriction modification system, M',
 169         'Check for hit description'
 170     );
 171     float_is( $hit->significance, 0.0022, 'Check hit significance' );
 172     is( $hit->raw_score, -105.2, 'Check hit raw_score' );
 173     my $hsp = $hit->next_hsp;
 174     is( $hsp->score, -105.2, 'Check for hsp score' );
 175     float_is( $hsp->evalue, 0.0022, 'Check for hsp evalue' );
 176     is( $hsp->query->start, 280, 'Check for query alifrom value' );
 177     is( $hsp->query->end,   481, 'Check for query ali to value' );
 178     is( $hsp->hit->start,   1,   'Check for hit hmmfrom value' );
 179     is( $hsp->hit->end,     279, 'Check for hit hmm to value' );
 180     is( $hsp->query->seq_id(),
 181         'gi|1522636|gb|AAC37060.1|', 'Check for query seq_id' );
 182     is( $hsp->hit->seq_id(), 'Methylase_M', 'Check for hit seq_id' );
 183     is( $hsp->hit_string,
 184         'lrnELentLWavADkLRGsmDaseYKdyVLGLlFlKYiSdkFlerrieieerktdtesepsldyakledqyeqlededlekedfyqkkGvFilPsqlFwdfikeaeknkldedigtdldkifseledqialgypaSeedfkGlfpdldfnsnkLgskaqarnetLtelidlfselelgtPmHNG-dfeelgikDlfGDaYEYLLgkFAeneGKsGGeFYTPqeVSkLiaeiLtigqpsegdfsIYDPAcGSGSLllqaskflgehdgkrnaisyYGQEsn',
 185         'Check for hiy string'
 186     );
 187     is( $hsp->query_string,
 188         'NTSELDKKKFAVLLMNR--------------LIFIKFLEDK------GIV---------PRDLLRRTYEDY---KKSNVLI-NYYDAY-L----KPLFYEVLNTPEDER--KENIRT-NPYYKDIPYL---N-G-------GLFRSNNV--PNELSFTIKDNEIIGEVINFLERYKFTLSTSEGsEEVELNP-DILGYVYEKLINILAEKGQKGLGAYYTPDEITSYIAKNT-IEPIVVE----------------RFKEIIK--NWKINDINF----ST',
 189         'Check for query string'
 190     );
 191     is( $hsp->homology_string,
 192         ' ++EL+++  av+   R              L+F K++ dk      +i+         p +   + +++y   ++   ++ ++y ++      + lF++++   e ++  ++++ + +    ++      + +       Glf ++++  ++ +s+   +ne ++e+i+ +++ +++     G++ +el   D++G +YE L+   Ae   K+ G +YTP e++  ia+ + i+  ++                  +++ ++    k+n+i +    s+',
 193         'Check for homology string'
 194     );
 195     is( join( ' ', $hsp->seq_inds( 'query', 'nomatch', 1 ) ),
 196         '280 288 289 293-295 300 304 311 313-315 317 324-326 332 335 337 344-346 348 355 358-361 364-366 372 379 383-385 389 396 400 404-408 412 416 417 422 426 429-431 434-436 439 441 446 450 451 455 459 460 463 464 468 471 472 478',
 197         'Check for nomatch indices in query'
 198     );
 199     is( join( ' ', $hsp->seq_inds( 'hit', 'nomatch', 1 ) ),
 200         '1 9 10 14-16 18-31 35 39 42-47 51-59 61 63-65 67 72-74 77-79 82 86 89-94 96 103-105 107 110 111 116 118 120-123 126-131 133 135-141 145 150 151 154 158-160 164 171 175 179-183 187 191-193 198 202 205-207 210-212 215 217 222 226 227 231 233 236 237 240-257 261 264-267 273 275-278',
 201         'Check for nomatch indices in hit'
 202     );
 203     is( join( ' ', $hsp->seq_inds( 'query', 'gap', 1 ) ),
 204         '296 306 309 321 328 334 335 350 356 366-368 376 417 456 463 470 479',
 205         'Check for gap indices in query'
 206     );
 207     is( join( ' ', $hsp->seq_inds( 'hit', 'gap', 1 ) ),
 208         '', 'Check for gap indices in hit' );
 209 }
 210
 211 $searchio = Bio::SearchIO->new(
 212     -format => 'hmmer',
 213     -file   => test_input_file('cysprot1b.hmmsearch')
 214 );
 215
 216 while ( $result = $searchio->next_result ) {
 217     is( ref($result),
 218         'Bio::Search::Result::HMMERResult',
 219         'Check for the correct result reference type'
 220     );
 221     is( $result->algorithm,         'HMMSEARCH', 'Check algorithm' );
 222     is( $result->algorithm_version, '2.2g',      'Check algorithm version' );
 223     is( $result->hmm_name,
 224         'Peptidase_C1.hmm [Peptidase_C1]',
 225         'Check hmm_name'
 226     );
 227     is( $result->database_name,   'cysprot1b.fa', 'Check database_name' );
 228     is( $result->sequence_file,   'cysprot1b.fa', 'Check sequence_file' );
 229     is( $result->query_name,      'Peptidase_C1', 'Check query_name' );
 230     is( $result->query_accession, 'PF00112',      'Check query_accession' );
 231     is( $result->query_description,
 232         'Papain family cysteine protease',
 233         'Check query_description'
 234     );
 235     is( $result->num_hits(), 4, 'Check num_hits' );
 236     my $hit = $result->next_hit;
 237     is( $hit->name,        'CATL_RAT', 'Check hit name' );
 238     is( $hit->description, '',         'Check for hit description' );
 239     float_is( $hit->significance, 2e-135, 'Check hit significance' );
 240     is( $hit->raw_score, 449.4, 'Check hit raw_score' );
 241     my $hsp = $hit->next_hsp;
 242     is( $hsp->score, 449.4, 'Check for hsp score' );
 243     float_is( $hsp->evalue, 2e-135, 'Check for hsp evalue' );
 244     is( $hsp->query->start, 1,   'Check for query alifrom value' );
 245     is( $hsp->query->end,   337, 'Check for query ali to value' );
 246     is( $hsp->hit->start,   114, 'Check for hit hmmfrom value' );
 247     is( $hsp->hit->end,     332, 'Check for hit hmm to value' );
 248     is( $hsp->query->seq_id(), 'Peptidase_C1', 'Check for query seq_id' );
 249     is( $hsp->hit->seq_id(),   'CATL_RAT',     'Check for hit seq_id' );
 250     is( $hsp->hit_string,
 251         'IPKTVDWRE-KG-CVTPVKNQG-QCGSCWAFSASGCLEGQMFLKT------GKLISLSEQNLVDCSH-DQGNQ------GCNG-GLMDFAFQYIKE-----NGGLDSEESY-----PYE----AKD-------------------GSCKYR-AEYAV-----ANDTGFVDIPQQ-----EKALMKAVATVGPISVAMDASHPS---LQFYSSG-------IYYEP---NCSSK---DLDHGVLVVGYGYEG-T------------------------------------DSNKDKYWLVKNSWGKEWGMDGYIKIAKDRN----NHCGLATAASYPI',
 252         'Check for hiy string'
 253     );
 254     is( $hsp->homology_string,
 255         '+P+++DWRe kg  VtpVK+QG qCGSCWAFSa g lEg+ ++kt      gkl+sLSEQ+LvDC++ d gn+      GCnG Glmd Af+Yik+     NgGl++E++Y     PY+    +kd                   g+Cky+  + ++     a+++g++d+p++     E+al+ka+a++GP+sVa+das+ s    q+Y+sG       +Y+++    C+++   +LdH+Vl+VGYG e+                                      ++++ +YW+VKNSWG++WG++GY++ia+++n    n+CG+a+ asypi',
 256         'Check for homology string'
 257     );
 258     is( $hsp->query_string,
 259         'lPesfDWReWkggaVtpVKdQGiqCGSCWAFSavgalEgryciktgtkawggklvsLSEQqLvDCdgedygnngesCGyGCnGGGlmdnAfeYikkeqIsnNgGlvtEsdYekgCkPYtdfPCgkdggndtyypCpgkaydpndTgtCkynckknskypktyakikgygdvpynvsTydEealqkalaknGPvsVaidasedskgDFqlYksGendvgyGvYkhtsageCggtpfteLdHAVliVGYGteneggtfdetssskksesgiqvssgsngssgSSgssgapiedkgkdYWIVKNSWGtdWGEnGYfriaRgknksgkneCGIaseasypi',
 260         'Check for query string'
 261     );
 262     $hit = $result->next_hit;
 263     is( $hit->name,        'CATL_HUMAN', 'Check hit name' );
 264     is( $hit->description, '',           'Check for hit description' );
 265     float_is( $hit->significance, 6.1e-134, 'Check hit significance' );
 266     is( $hit->raw_score, 444.5, 'Check hit raw_score' );
 267 }
 268
 269 # test for bug 2632 - CS lines should get ignored without breaking the parser
 270 $searchio = Bio::SearchIO->new(
 271     -format => 'hmmer',
 272     -file   => test_input_file('hmmpfam_cs.out')
 273 );
 274 $result = $searchio->next_result;
 275 my $hit = $result->next_hit;
 276 my $hsp = $hit->next_hsp;
 277 is( $hsp->seq_str,
 278     'CGV-GFIADVNNVANHKIVVQALEALTCMEHRGACSADRDSGDGAGITTAIPWNLFQKSLQNQNIKFEQnDSVGVGMLFLPAHKLKES--KLIIETVLKEENLEIIGWRLVPTVQEVLGKQAYLNKPHVEQVFCKSSNLSKDRLEQQLFLVRKKIEKYIGINGKDwaheFYICSLSCYTIVYKGMMRSAVLGQFYQDLYHSEYTSSFAIYHRRFSTNTMPKWPLAQPMR---------FVSHNGEINTLLGNLNWMQSREPLLQSKVWKDRIHELKPITNKDNSDSANLDAAVELLIASGRSPEEALMILVPEAFQNQPDFA-NNTEISDFYEYYSGLQEPWDGPALVVFTNGKV-IGATLDRNGL-RPARYVIT----KDNLVIVSSES',
 279     'Check for hsp seq_str'
 280 );
 281
 282 # Tests for hmmer3 output here
 283 $searchio = Bio::SearchIO->new(
 284     -format  => 'hmmer',
 285     -file    => test_input_file('hmmscan.out'),
 286     -verbose => 1
 287 );
 288 is( ref($searchio), 'Bio::SearchIO::hmmer3',
 289     'Check if correct searchio object is returned' );
 290 while ( $result = $searchio->next_result ) {
 291     is( ref($result),
 292         'Bio::Search::Result::HMMERResult',
 293         'Check for the correct result reference type'
 294     );
 295     is( $result->algorithm,         'HMMSCAN', 'Check algorithm' );
 296     is( $result->algorithm_version, '3.0',     'Check algorithm version' );
 297     is( $result->hmm_name,
 298         '/data/biodata/HMMerDB/Pfam.hmm',
 299         'Check hmm_name'
 300     );
 301     is( $result->sequence_file, 'BA000019.orf1.fasta',
 302         'Check sequence_file' );
 303     is( $result->query_name, 'BA000019.orf1', 'Check query_name' );
 304     is( $result->query_length, '198', 'Check query_length' );
 305     is( $result->query_description, '', 'Check query_description' );
 306     is( $result->num_hits(),        1,  'Check num_hits' );
 307     my ( $hsp, $hit );
 308
 309     if ( $hit = $result->next_model ) {
 310         is( ref($hit), 'Bio::Search::Hit::HMMERHit',
 311             'Check for the correct hit reference type' );
 312         is( $hit->name, 'Peripla_BP_2', 'Check hit name' );
 313         is( $hit->description,
 314             'Periplasmic binding protein',
 315             'Check for hit description'
 316         );
 317         is( $hit->raw_score, '105.2', 'Check hit raw_score' );
 318         float_is( $hit->significance, 6e-30, 'Check hit significance' );
 319         is( $hit->num_hsps, 1, 'Check num_hsps' );
 320
 321         if ( defined( $hsp = $hit->next_domain ) ) {
 322             is( ref($hsp), 'Bio::Search::HSP::HMMERHSP',
 323                 'Check for correct hsp reference type' );
 324             is( $hsp->hit->start,   59,  'Check for hit hmmfrom value' );
 325             is( $hsp->hit->end,     236, 'Check for hit hmm to value' );
 326             is( $hsp->query->start, 2,   'Check for query alifrom value' );
 327             is( $hsp->query->end,   173, 'Check for query ali to value' );
 328             is( $hsp->score, '105.0', 'Check for hsp score' );
 329             float_is( $hsp->evalue, 1.5e-33, 'Check for hsp c-Evalue' );
 330             is( $hsp->query_string,
 331                 'LKPDLIIGREYQ---KNIYNQLSNFAPTVLVDWGSF-TSFQDNFRYIAQVLNEEEQGKLVLQQYQKRIRDLQDRMGERlQKIEVSVIGFSGQSIKSLNR-DAVFNQVLDDAGIKRIsIQKNQQERYLEISIENLNKYDADVLFVINE---SKEQLYPDLKNPLWHHLRAVKKQQVYVVNQ',
 332                 'Check for query string'
 333             );
 334             is( $hsp->hit_string,
 335                 'lkPDlvivsafgalvseieellelgipvvavessstaeslleqirllgellgeedeaeelvaelesridavkaridsl-kpktvlvfgyadegikvvfgsgswvgdlldaaggeni-iaeakgseseeisaEqilaadpdviivsgrgedtktgveelkenplwaelpAvkngrvyllds',
 336                 'Check for hit string'
 337             );
 338             is( $hsp->homology_string,
 339                 'lkPDl+i+ +++   ++i+++l++ +p+v v+  s+  s+++ +r ++++l+ee++++ + +++++ri+++++r  +  ++ +v+v+g+++ +ik+++  +  ++++ld+ag++ i i++++++ + eis+E+++++d+dv++v       k+ +   ++nplw +l+Avk+++vy++++',
 340                 'Check for homology string'
 341             );
 342         }
 343     }
 344 }
 345
 346 $searchio = Bio::SearchIO->new(
 347     -format  => 'hmmer',
 348     -file    => test_input_file('hmmsearch3.out'),
 349     -verbose => 1
 350 );
 351 while ( $result = $searchio->next_result ) {
 352     is( ref($result),
 353         'Bio::Search::Result::HMMERResult',
 354         'Check for the correct result reference type'
 355     );
 356     is( $result->algorithm,         'HMMSEARCH', 'Check algorithm' );
 357     is( $result->algorithm_version, '3.0',       'Check algorithm version' );
 358     is( $result->hmm_name,          'Kv9.hmm',   'Check hmm_name' );
 359     is( $result->sequence_file,
 360         '/home/pboutet/Desktop/databases/nr_May26',
 361         'Check sequence_file'
 362     );
 363     is( $result->query_name,        'Kv9', 'Check query_name' );
 364     is( $result->query_length,      '481', 'Check query_length' );
 365     is( $result->query_description, '',    'Check query_description' );
 366     is( $result->num_hits(),        2,     'Check num_hits' );
 367
 368     while ( my $hit = $result->next_model ) {
 369     }
 370 }
 371
 372 $searchio = Bio::SearchIO->new(
 373     -format  => 'hmmer',
 374     -file    => test_input_file('hmmscan_multi_domain.out'),
 375     -verbose => 1
 376 );
 377
 378 my @multi_hits = (
 379     [   'PPC',
 380         'Bacterial pre-peptidase C-terminal domain',
 381         '111.0', 3.1e-32, 6,
 382         [   [ 4,  59, 117,  183,  0.5,  0.16 ],
 383             [ 12, 58, 347,  388,  -0.6, 0.36 ],
 384             [ 1,  69, 470,  549,  71.3, 1.3e-23 ],
 385             [ 15, 25, 582,  603,  -3.2, 2 ],
 386             [ 13, 36, 987,  1019, -1.1, 0.5 ],
 387             [ 1,  69, 1087, 1168, 54.4, 2.4e-18 ]
 388         ]
 389     ],
 390     [   'HemolysinCabind',
 391         'Hemolysin-type calcium-binding repeat (2 cop',
 392         '47.9', 4.7e-13, 3,
 393         [   [ 2, 13, 1214, 1225, 5.9,  0.0026 ],
 394             [ 1, 18, 1231, 1248, 10.8, 6.8e-5 ],
 395             [ 4, 18, 1243, 1257, 11.4, 4.3e-05 ]
 396         ]
 397     ]
 398 );
 399
 400 while ( $result = $searchio->next_result ) {
 401     is( ref($result),
 402         'Bio::Search::Result::HMMERResult',
 403         'Check for the correct result reference type'
 404     );
 405     is( $result->algorithm,         'HMMSCAN', 'Check algorithm' );
 406     is( $result->algorithm_version, '3.0',     'Check algorithm version' );
 407     is( $result->hmm_name,
 408         '/data/biodata/HMMerDB/Pfam-A.hmm',
 409         'Check hmm_name'
 410     );
 411     is( $result->sequence_file, 'BA000019.orf37.fasta',
 412         'Check sequence_file' );
 413     is( $result->query_name, 'BA000019.orf37', 'Check query_name' );
 414     is( $result->query_length, '1418', 'Check query_length' );
 415     is( $result->query_description, '', 'Check query_description' );
 416     is( $result->num_hits(),        2,  'Check num_hits' );
 417     my ( $hsp, $hit );
 418
 419     while ( $hit = $result->next_model ) {
 420         my @expected = @{ shift @multi_hits };
 421         is( ref($hit), 'Bio::Search::Hit::HMMERHit',
 422             'Check for the correct hit reference type' );
 423         is( $hit->name,        shift @expected, 'Check hit name' );
 424         is( $hit->description, shift @expected, 'Check for hit description' );
 425         is( $hit->raw_score,   shift @expected, 'Check hit raw_score' );
 426         float_is(
 427             $hit->significance,
 428             shift @expected,
 429             'Check hit significance'
 430         );
 431         is( $hit->num_hsps, shift @expected, 'Check num_hsps' );
 432         my @hsp_list = @{ shift @expected };
 433
 434         while ( defined( $hsp = $hit->next_domain ) ) {
 435             my @hsp_exp = @{ shift @hsp_list };
 436             is( ref($hsp), 'Bio::Search::HSP::HMMERHSP',
 437                 'Check for correct hsp reference type' );
 438             is( $hsp->hit->start,
 439                 shift @hsp_exp,
 440                 'Check for hit envfrom value'
 441             );
 442             is( $hsp->hit->end, shift @hsp_exp,
 443                 'Check for hit env to value' );
 444             is( $hsp->query->start,
 445                 shift @hsp_exp,
 446                 'Check for query hmmfrom value'
 447             );
 448             is( $hsp->query->end,
 449                 shift @hsp_exp,
 450                 'Check for query hmm to value'
 451             );
 452             is( $hsp->score, shift @hsp_exp, 'Check for hsp score' );
 453             float_is( $hsp->evalue, shift @hsp_exp,
 454                 'Check for hsp c-Evalue' );
 455         }
 456     }
 457 }
 458
 459 $searchio = Bio::SearchIO->new(
 460     -format  => 'hmmer',
 461     -file    => test_input_file('hmmscan_sec_struct.out'),
 462     -verbose => 1
 463 );
 464
 465 @multi_hits = (
 466     [   'HTH_AraC',
 467         'Bacterial regulatory helix-turn-helix proteins, Ara',
 468         '41.3', 6.7e-11, 2,
 469         [   [ 'siadiAeevgfSpsyfsrlFkkytGvt', 'SLMELSRQVGLNDCTLKRGFRLVFDTT' ],
 470             [   'nwsiadiAeevgf-SpsyfsrlFkkytGvtPsqyr',
 471                 'EINISQAARRVGFsSRSYFATAFRKKFGINPKEFL'
 472             ]
 473         ]
 474     ],
 475     [   'PKSI-KS_m3',
 476         '', '38.2', 3.8e-12, 2,
 477         [   [ 'GPSvtVDTACSSSLvA', 'GPSVTVDTLCSSSLVA' ],
 478             [ 'GPSvtVDTACSSSLv',  'GPNLVIDSACSSALV' ]
 479         ]
 480     ],
 481     [   'DUF746',
 482         'Domain of Unknown Function (DUF746)',
 483         '13.9', 0.023, 2,
 484         [   [   'rllIrlLsqplslaeaadqlgtdegiiak',
 485                 'EILIRNLENPPSLMELSRQVGLNDCTLKR'
 486             ],
 487             [ 'plslaeaadqlgtdeg', 'EINISQAARRVGFSSR' ]
 488         ]
 489     ]
 490 );
 491
 492 while ( $result = $searchio->next_result ) {
 493     is( ref($result),
 494         'Bio::Search::Result::HMMERResult',
 495         'Check for the correct result reference type'
 496     );
 497     is( $result->algorithm,         'HMMSCAN',    'Check algorithm' );
 498     is( $result->algorithm_version, '3.0',        'Check algorithm version' );
 499     is( $result->hmm_name,          'Pfam-A.hmm', 'Check hmm_name' );
 500     is( $result->sequence_file, 'BA000019.orf8.fasta',
 501         'Check sequence_file' );
 502     is( $result->query_name, 'BA000019.orf8', 'Check query_name' );
 503     is( $result->query_length, '348', 'Check query_length' );
 504     is( $result->query_description, '', 'Check query_description' );
 505     is( $result->num_hits(),        3,  'Check num_hits' );
 506     my ( $hsp, $hit );
 507
 508     while ( $hit = $result->next_model ) {
 509         my @expected = @{ shift @multi_hits };
 510         is( ref($hit), 'Bio::Search::Hit::HMMERHit',
 511             'Check for the correct hit reference type' );
 512         is( $hit->name,        shift @expected, 'Check hit name' );
 513         is( $hit->description, shift @expected, 'Check for hit description' );
 514         is( $hit->raw_score,   shift @expected, 'Check hit raw_score' );
 515         float_is(
 516             $hit->significance,
 517             shift @expected,
 518             'Check hit significance'
 519         );
 520         is( $hit->num_hsps, shift @expected, 'Check num_hsps' );
 521         my @hsp_list = @{ shift @expected };
 522
 523         while ( defined( $hsp = $hit->next_domain ) ) {
 524             my @hsp_exp = @{ shift @hsp_list };
 525             is( ref($hsp), 'Bio::Search::HSP::HMMERHSP',
 526                 'Check for correct hsp reference type' );
 527             is( $hsp->hit_string,   shift @hsp_exp, 'Check hit sequence' );
 528             is( $hsp->query_string, shift @hsp_exp, 'Check query sequence' );
 529         }
 530     }
 531 }
 532
 533 # Make sure that you can also directly call the hmmer2 and hmmer3 subclasses
 534 $searchio = Bio::SearchIO->new(
 535     -format => 'hmmer2',
 536     -file   => test_input_file('hmmpfam.out')
 537 );
 538 is( ref($searchio), 'Bio::SearchIO::hmmer2',
 539     'Check if loading hmmpfam output via the hmm2 parser directly works' );
 540 is( ref( $searchio->next_result ),
 541     'Bio::Search::Result::HMMERResult',
 542     'Check for the correct result reference type'
 543 );
 544
 545 $searchio = Bio::SearchIO->new(
 546     -format => 'hmmer2',
 547     -file   => test_input_file('hmmsearch.out')
 548 );
 549 is( ref($searchio), 'Bio::SearchIO::hmmer2',
 550     'Check if loading hmmsearch2 output via the hmm2 parser directly works' );
 551 is( ref( $searchio->next_result ),
 552     'Bio::Search::Result::HMMERResult',
 553     'Check for the correct result reference type'
 554 );
 555
 556 $searchio = Bio::SearchIO->new(
 557     -format => 'hmmer3',
 558     -file   => test_input_file('hmmscan.out')
 559 );
 560 is( ref($searchio), 'Bio::SearchIO::hmmer3',
 561     'Check if loading hmmscan output via the hmm3 parser directly works' );
 562 is( ref( $searchio->next_result ),
 563     'Bio::Search::Result::HMMERResult',
 564     'Check for the correct result reference type'
 565 );
 566
 567 $searchio = Bio::SearchIO->new(
 568     -format => 'hmmer',
 569     -file   => test_input_file('hmmsearch3.out')
 570 );
 571 is( ref($searchio), 'Bio::SearchIO::hmmer3',
 572     'Check if loading hmmsearch3 output via the hmm3 parser directly works' );
 573 is( ref( $searchio->next_result ),
 574     'Bio::Search::Result::HMMERResult',
 575     'Check for the correct result reference type'
 576 );
 577
 578 # Make sure that you can also specify the -version parameter directly
 579 $searchio = Bio::SearchIO->new(
 580     -format  => 'hmmer',
 581     -file    => test_input_file('hmmpfam.out'),
 582     -version => 2
 583 );
 584 is( ref($searchio), 'Bio::SearchIO::hmmer2',
 585     'Check if selecting the correct hmmpfam parser using -version works' );
 586 is( ref( $searchio->next_result ),
 587     'Bio::Search::Result::HMMERResult',
 588     'Check for the correct result reference type'
 589 );
 590
 591 $searchio = Bio::SearchIO->new(
 592     -format  => 'hmmer',
 593     -file    => test_input_file('hmmsearch.out'),
 594     -version => 2
 595 );
 596 is( ref($searchio), 'Bio::SearchIO::hmmer2',
 597     'Check if selecting the correct hmmsearch2 parser using -version works' );
 598 is( ref( $searchio->next_result ),
 599     'Bio::Search::Result::HMMERResult',
 600     'Check for the correct result reference type'
 601 );
 602
 603 $searchio = Bio::SearchIO->new(
 604     -format  => 'hmmer3',
 605     -file    => test_input_file('hmmscan.out'),
 606     -version => 3
 607 );
 608 is( ref($searchio), 'Bio::SearchIO::hmmer3',
 609     'Check if selecting the correct hmmscan parser using -version works' );
 610 is( ref( $searchio->next_result ),
 611     'Bio::Search::Result::HMMERResult',
 612     'Check for the correct result reference type'
 613 );
 614
 615 $searchio = Bio::SearchIO->new(
 616     -format  => 'hmmer',
 617     -file    => test_input_file('hmmsearch3.out'),
 618     -version => 3
 619 );
 620 is( ref($searchio), 'Bio::SearchIO::hmmer3',
 621     'Check if selecting the correct hmmsearch3 parser using -version works' );
 622 is( ref( $searchio->next_result ),
 623     'Bio::Search::Result::HMMERResult',
 624     'Check for the correct result reference type'
 625 );
 626
 627 my $pipestr = "cat " . test_input_file('hmmpfam.out') . " |";
 628 open( my $pipefh, $pipestr );
 629
 630 $searchio = Bio::SearchIO->new(
 631     -format => 'hmmer',
 632     -fh     => $pipefh
 633 );
 634 is( ref($searchio), 'Bio::SearchIO::hmmer2',
 635     'Check if reading from a pipe works' );
 636 $result = $searchio->next_result;
 637 is( ref($result),
 638     'Bio::Search::Result::HMMERResult',
 639     'Check for the correct result reference type'
 640 );
 641 is( $result->num_hits(), 2, 'Check num_hits' );
 642
 643 # bug 3376
 644 {
 645     my $in = Bio::SearchIO->new(
 646         -format => 'hmmer',
 647         -file   => test_input_file('pfamOutput-bug3376.out')
 648     );
 649     my $result = $in->next_result;
 650     my $hit    = $result->next_hit;
 651     my $hsp    = $hit->next_hsp;
 652     is( $hsp->hit_string,
 653         "svfqqqqssksttgstvtAiAiAigYRYRYRAvtWnsGsLssGvnDnDnDqqsdgLYtiYYsvtvpssslpsqtviHHHaHkasstkiiikiePr",
 654         "bug3376"
 655     );
 656 }
 657 # end bug 3376
 658
 659 # bug 3421 - making sure a full line of dashes in an HSP is parsed correctly
 660 {
 661     my $in = Bio::SearchIO->new(
 662         -format => 'hmmer',
 663         -file   => test_input_file('hmmpfam_HSPdashline.txt')
 664     );
 665     my $result = $in->next_result;
 666     my $hit    = $result->next_hit;
 667     my $hsp    = $hit->next_hsp;
 668     is( $hsp->length, "561",
 669         "bug3421 - Check if can correctly parse an HSP with line full of dashes"
 670     );
 671 }
 672 # end bug 3421
 673
 674 # bug 3302
 675 {
 676     my $in = Bio::SearchIO->new(
 677         -format => 'hmmer',
 678         -file   => test_input_file('hmmpfam_multiresult.out')
 679     );
 680     my $result = $in->next_result;
 681     $result = $in->next_result;
 682     my $hit = $result->next_hit;
 683     is( $hit->name, "IS66_ORF3.uniq", "bug3302 - Check if can parse multiresult hmmer" );
 684 }
 685 # end bug 3302
 686
 687 # HMMER 3.1 nhmmer output
 688 {
 689     my $in = Bio::SearchIO->new(
 690         -format  => 'hmmer',
 691         -version => 3,
 692         -file    => test_input_file('nhmmer-3.1.out')
 693     );
 694     my $result = $in->next_result;
 695     is( $result->algorithm_version, '3.1b1', 'Check nhmmer algorithm version' );
 696
 697     my $hit = $result->next_hit;
 698     is( $hit->name, "seq1", "Check nhmmer hit name" );
 699     is( $hit->description, "Description of seq1", "Check nhmmer hit description" );
 700     is( $hit->significance, "3.2e-48", "Check nhmmer hit significance" );
 701     is( $hit->score,        "148.2",   "Check nhmmer hit score" );
 702
 703     my $hsp = $hit->next_hsp;
 704     is( $hsp->score, "148.2",   "Check nhmmer hsp score" );
 705     is( $hsp->significance, "3.2e-48",   "Check nhmmer hsp score" );
 706     is( $hsp->start('hit'), "1",   "Check nhmmer hsp hit start" );
 707     is( $hsp->end('hit'), "151",   "Check nhmmer hsp hit end" );
 708     is( $hsp->start('query'), "258",   "Check nhmmer hsp query start" );
 709     is( $hsp->end('query'), "411",   "Check nhmmer hsp query end" );
 710 }
 711 # end HMMER 3.1 nhmmer output
 712