tests/blast_test.py

   1 from itertools import *
   2 import re
   3 import unittest
   4 import glob
   5 import os
   6 from testlib import testutil, SkipTest, PygrTestProgram
   7 from pygr import worldbase
   8 from pygr import sequence, cnestedlist, seqdb, blast, logger, parse_blast
   9 from pygr.nlmsa_utils import CoordsGroupStart,CoordsGroupEnd
  10 from pygr import translationDB
  11
  12 def check_results(results, correct, formatter, delta=0.01,
  13                   reformatCorrect=False, reformatResults=True):
  14     if reformatResults:
  15         results = reformat_results(results, formatter)
  16
  17     if reformatCorrect: # reformat these data too
  18         correct = reformat_results(correct, formatter)
  19     else:
  20         correct.sort()
  21
  22     # this is to help troubleshooting the mismatches if there are any
  23     mismatch = [ (a, b) for a, b in zip(correct, results) if
  24                  testutil.approximate_cmp([a], [b], delta)]
  25     if mismatch:
  26         logger.warn('blast mismatches found')
  27         for m in mismatch:
  28             logger.warn('%s != %s' % m)
  29
  30     # this is the actual test
  31     assert testutil.approximate_cmp(correct, results, delta) == 0
  32
  33 def check_results_relaxed_blastp(results, correct, formatter, delta=0.01,
  34                                  reformatCorrect=False, allowedLengthDiff=0,
  35                                  identityMin=0.6, reformatResults=True):
  36     if reformatResults:
  37         results = reformat_results(results, formatter)
  38
  39     if reformatCorrect: # reformat these data too
  40         correct = reformat_results(correct, formatter)
  41     else:
  42         correct.sort()
  43
  44     # Length of output
  45     assert abs(len(results) - len(correct)) <= allowedLengthDiff
  46
  47     # Format check
  48     key_re = re.compile('^[A-Z]{3}[A-Z0-9]?_[A-Z]{2,5}$')
  49     for result in results:
  50         assert key_re.search(result[0])
  51         assert key_re.search(result[1])
  52         assert (0. < result[2] and result[2] <= 1.)
  53
  54     # High-identity comparison
  55     results_high = []
  56     correct_high = []
  57     for result in results:
  58         if result[2] > identityMin:
  59             results_high.append(result)
  60     for result in correct:
  61         if result[2] > identityMin:
  62             correct_high.append(result)
  63     assert testutil.approximate_cmp(correct_high, results_high, delta) == 0
  64
  65 def check_results_relaxed_blastx(results, correct, formatter, delta=0.01,
  66                                  reformatCorrect=False, allowedLengthDiff=0,
  67                                  identityMin=0.6):
  68     results = reformat_results(results, formatter)
  69
  70     if reformatCorrect: # reformat these data too
  71         correct = reformat_results(correct, formatter)
  72     else:
  73         correct.sort()
  74
  75     # Length of output
  76     assert abs(len(results) - len(correct)) <= allowedLengthDiff
  77
  78     # Format check
  79     for result in results:
  80         assert 3 * result[0] == result[2]
  81         assert (0. < result[3] and result[3] <= 1.)
  82
  83     # High-identity comparison
  84     results_high = []
  85     correct_high = []
  86     for result in results:
  87         if result[3] > identityMin:
  88             results_high.append(result)
  89     for result in correct:
  90         if result[3] > identityMin:
  91             correct_high.append(result)
  92     assert testutil.approximate_cmp(correct_high, results_high, delta) == 0
  93
  94 def reformat_results(results, formatter):
  95     reffed = []
  96     for result in results:
  97         for t in result.edges(mergeMost=True):
  98             reffed.append(formatter(t))
  99     reffed.sort()
 100     return reffed
 101
 102 def pair_identity_tuple(t):
 103     'standard formatter for blast matches'
 104     return (t[0].id, t[1].id, t[2].pIdentity())
 105
 106 class BlastBase(unittest.TestCase):
 107     def setUp(self):
 108         hbb1_mouse = testutil.datafile('hbb1_mouse.fa')
 109         hbb1_mouse_rc = testutil.datafile('hbb1_mouse_rc.fa')
 110         sp_hbb1 = testutil.datafile('sp_hbb1')
 111         gapping = testutil.datafile('gapping.fa')
 112
 113         self.dna = seqdb.SequenceFileDB(hbb1_mouse)
 114         self.dna_rc = seqdb.SequenceFileDB(hbb1_mouse_rc)
 115         self.prot = seqdb.SequenceFileDB(sp_hbb1)
 116         self.gapping = seqdb.SequenceFileDB(gapping)
 117
 118     def tearDown(self):
 119         'do the RIGHT thing... close resources that have been opened!'
 120         self.dna.close()
 121         self.dna_rc.close()
 122         self.prot.close()
 123         self.gapping.close()
 124
 125
 126 _multiblast_results = None
 127
 128 class Blast_Test(BlastBase):
 129     """
 130     Test basic BLAST stuff (using blastp).
 131     """
 132     def test_blastp(self):
 133         "Testing blastp"
 134         if not testutil.blast_enabled():
 135             raise SkipTest, "no BLAST installed"
 136
 137         blastmap = blast.BlastMapping(self.prot, verbose=False)
 138         results = blastmap[self.prot['HBB1_XENLA']]
 139
 140         check_results_relaxed_blastp([results], blastp_correct_results,
 141                                      pair_identity_tuple,
 142                                      allowedLengthDiff=2)
 143
 144     def test_repr(self):
 145         blastmap = blast.BlastMapping(self.prot, verbose=False)
 146         assert '<BlastMapping' in repr(blastmap)
 147
 148     def test_no_query(self):
 149         blastmap = blast.BlastMapping(self.dna, verbose=False)
 150         try:
 151             blastmap()
 152             assert 0, "should fail before this"
 153         except ValueError:
 154             pass
 155
 156     def test_both_seq_and_db(self):
 157         "Testing db arg present"
 158         blastmap = blast.BlastMapping(self.prot, verbose=False)
 159         seq = self.prot['HBB1_XENLA']
 160
 161         try:
 162             blastmap(seq=seq, queryDB=self.prot)
 163             assert 0, "should fail before this"
 164         except ValueError:
 165             pass
 166
 167     def test_multiblast(self):
 168         "testing multi sequence blast"
 169         results = self.get_multiblast_results()
 170         check_results_relaxed_blastp(results, correct_multiblast_results,
 171                                      None, reformatResults=False,
 172                                      allowedLengthDiff=10)
 173
 174     def get_multiblast_results(self):
 175         """return saved results or generate them if needed;
 176         results are saved so we only do this time-consuming operation once"""
 177         global _multiblast_results
 178
 179         if not testutil.blast_enabled():
 180             raise SkipTest, "no BLAST installed"
 181
 182         if not _multiblast_results:
 183             logger.info("running expensive multiblast")
 184             blastmap = blast.BlastMapping(self.prot, verbose=False)
 185             al = cnestedlist.NLMSA('blasthits', 'memory', pairwiseMode=True,
 186                                    bidirectional=False)
 187
 188             blastmap(al=al, queryDB=self.prot) # all vs all
 189
 190             al.build() # construct the alignment indexes
 191             results = [al[seq] for seq in self.prot.values()]
 192             _multiblast_results = reformat_results(results,
 193                                                    pair_identity_tuple)
 194
 195         return _multiblast_results
 196
 197     def test_multiblast_single(self):
 198         "Test multi-sequence BLAST results, for BLASTs run one by one."
 199         if not testutil.blast_enabled():
 200             raise SkipTest, "no BLAST installed"
 201
 202         blastmap = blast.BlastMapping(self.prot, verbose=False)
 203         al = cnestedlist.NLMSA('blasthits', 'memory', pairwiseMode=True,
 204                                bidirectional=False)
 205
 206         for seq in self.prot.values():
 207             blastmap(seq, al) # all vs all, one by one
 208
 209         al.build() # construct the alignment indexes
 210         results = [al[seq] for seq in self.prot.values()]
 211         results_multi = self.get_multiblast_results()
 212         # Strict check must work here even on live BLAST output
 213         check_results(results, results_multi, pair_identity_tuple)
 214
 215     def test_multiblast_long(self):
 216         "testing multi sequence blast with long db"
 217         if not testutil.blast_enabled():
 218             raise SkipTest, "no BLAST installed"
 219
 220         longerFile = testutil.datafile('sp_all_hbb')
 221
 222         sp_all_hbb = seqdb.SequenceFileDB(longerFile)
 223         blastmap = blast.BlastMapping(self.prot, verbose=False)
 224         al = cnestedlist.NLMSA('blasthits', 'memory', pairwiseMode=True,
 225                                bidirectional=False)
 226         blastmap(None, al, queryDB=sp_all_hbb) # all vs all
 227         al.build() # construct the alignment indexes
 228
 229     def test_maskEnd(self):
 230         """
 231         This tests against a minor bug in cnestedlist where maskEnd
 232         is used to clip the end to the mask region.
 233         """
 234         if not testutil.blast_enabled():
 235             raise SkipTest, "no BLAST installed"
 236
 237         db = self.gapping
 238         blastmap = blast.BlastMapping(db)
 239         ungapped = db['ungapped']
 240         gapped = db['gapped']
 241         results = blastmap[gapped]
 242
 243         results[ungapped]
 244
 245     def test_no_bidirectional(self):
 246         if not testutil.blast_enabled():
 247             raise SkipTest, "no BLAST installed"
 248
 249         db = self.gapping
 250         gapped = db['gapped']
 251         ungapped = db['ungapped']
 252
 253         blastmap = blast.BlastMapping(db)
 254         al = blastmap(queryDB=db)
 255         slice = al[gapped]
 256
 257         found_once = False
 258         for src, dest, edge in al[gapped].edges():
 259             if src == gapped[0:40] and dest == ungapped[0:40]:
 260                 assert not found_once, \
 261                        "BLAST results should not be bidirectional"
 262                 found_once = True
 263
 264         assert found_once, "should have found this match exactly once!"
 265
 266     def test_formatdb_fail(self):
 267         db = self.gapping
 268         try:
 269             blastmap = blast.BlastMapping(db, filepath='foobarbaz.fa',
 270                                           blastReady=True,
 271                                           showFormatdbMessages=False)
 272             assert 0, "should not reach this point"
 273         except IOError:                 # should fail with 'cannot build'
 274             pass
 275
 276         remnants = glob.glob('foobarbaz.fa.n??')
 277         for filename in remnants:
 278             os.unlink(filename)
 279
 280     def test_seq_without_db(self):
 281         "Check that sequences without associated DBs work as query strings"
 282         blastmap = blast.BlastMapping(self.prot, verbose=False)
 283
 284         seq = self.prot['HBB1_XENLA']
 285         seq_no_db = sequence.Sequence(str(seq), 'HBB1_XENLA_no_db')
 286         slice = blastmap(seq=seq_no_db)[seq_no_db]
 287         assert len(slice)
 288
 289 class Blastx_Test(BlastBase):
 290     def test_blastx(self):
 291         "Testing blastx"
 292         if not testutil.blast_enabled():
 293             raise SkipTest, "no BLAST installed"
 294
 295         blastmap = blast.BlastxMapping(self.prot, verbose=False)
 296
 297         correct = [(143, 143, 429, 0.53146853146853146),
 298                    (143, 145, 429, 0.28275862068965518),
 299                    (143, 145, 429, 0.28965517241379313),
 300                    (143, 145, 429, 0.29655172413793102),
 301                    (143, 145, 429, 0.30344827586206896),
 302                    (144, 144, 432, 0.4513888888888889),
 303                    (144, 144, 432, 0.4513888888888889),
 304                    (145, 145, 435, 0.45517241379310347),
 305                    (145, 145, 435, 0.51034482758620692),
 306                    (146, 142, 438, 0.35616438356164382),
 307                    (146, 146, 438, 0.4589041095890411),
 308                    (146, 146, 438, 0.46575342465753422),
 309                    (146, 146, 438, 0.4726027397260274),
 310                    (146, 146, 438, 0.4726027397260274),
 311                    (146, 146, 438, 0.4863013698630137),
 312                    (146, 146, 438, 0.59589041095890416),
 313                    (146, 146, 438, 0.62328767123287676),
 314                    (146, 146, 438, 0.66438356164383561),
 315                    (146, 146, 438, 0.74657534246575341),
 316                    (146, 146, 438, 0.91095890410958902),
 317                    (146, 146, 438, 0.97945205479452058)]
 318
 319         results = blastmap[self.dna['gi|171854975|dbj|AB364477.1|']]
 320         check_results_relaxed_blastx(results, correct,
 321                       lambda t:(len(t[0]), len(t[1]), len(t[0].sequence),
 322                                 t[2].pIdentity()), allowedLengthDiff=2)
 323
 324     def test_blastx_rc(self):
 325         "Testing blastx with negative frames"
 326         if not testutil.blast_enabled():
 327             raise SkipTest, "no BLAST installed"
 328
 329         blastmap = blast.BlastxMapping(self.prot, verbose=False)
 330
 331         correct = [(143, 143, 429, 0.53146853146853146),
 332                    (143, 145, 429, 0.28275862068965518),
 333                    (143, 145, 429, 0.28965517241379313),
 334                    (143, 145, 429, 0.29655172413793102),
 335                    (143, 145, 429, 0.30344827586206896),
 336                    (144, 144, 432, 0.4513888888888889),
 337                    (144, 144, 432, 0.4513888888888889),
 338                    (145, 145, 435, 0.45517241379310347),
 339                    (145, 145, 435, 0.51034482758620692),
 340                    (146, 142, 438, 0.35616438356164382),
 341                    (146, 146, 438, 0.4589041095890411),
 342                    (146, 146, 438, 0.46575342465753422),
 343                    (146, 146, 438, 0.4726027397260274),
 344                    (146, 146, 438, 0.4726027397260274),
 345                    (146, 146, 438, 0.4863013698630137),
 346                    (146, 146, 438, 0.59589041095890416),
 347                    (146, 146, 438, 0.62328767123287676),
 348                    (146, 146, 438, 0.66438356164383561),
 349                    (146, 146, 438, 0.74657534246575341),
 350                    (146, 146, 438, 0.91095890410958902),
 351                    (146, 146, 438, 0.97945205479452058)]
 352
 353         results = blastmap[self.dna_rc['hbb1_mouse_RC']]
 354         check_results_relaxed_blastx(results, correct,
 355                       lambda t:(len(t[0]), len(t[1]), len(t[0].sequence),
 356                                 t[2].pIdentity()), allowedLengthDiff=2)
 357
 358         results = blastmap[self.dna_rc['hbb1_mouse_RC_2']]
 359         check_results_relaxed_blastx(results, correct,
 360                       lambda t:(len(t[0]), len(t[1]), len(t[0].sequence),
 361                                 t[2].pIdentity()), allowedLengthDiff=2)
 362
 363         results = blastmap[self.dna_rc['hbb1_mouse_RC_3']]
 364         check_results_relaxed_blastx(results, correct,
 365                       lambda t:(len(t[0]), len(t[1]), len(t[0].sequence),
 366                                 t[2].pIdentity()), allowedLengthDiff=2)
 367
 368     def test_repr(self):
 369         blastmap = blast.BlastxMapping(self.prot, verbose=False)
 370         assert '<BlastxMapping' in repr(blastmap)
 371
 372     def test_blastx_no_blastp(self):
 373         if not testutil.blast_enabled():
 374             raise SkipTest, "no BLAST installed"
 375
 376         blastmap = blast.BlastxMapping(self.prot, verbose=False)
 377
 378         try:
 379             results = blastmap(self.prot['HBB1_MOUSE'])
 380             raise AssertionError('failed to trap blastp in BlastxMapping')
 381         except ValueError:
 382             pass
 383
 384     def test_no_query(self):
 385         blastmap = blast.BlastxMapping(self.prot, verbose=False)
 386         try:
 387             blastmap()
 388             assert 0, "should fail before this"
 389         except ValueError:
 390             pass
 391
 392     def test_both_seq_and_db(self):
 393         "Testing blastp"
 394         blastmap = blast.BlastxMapping(self.prot, verbose=False)
 395         seq = self.prot['HBB1_XENLA']
 396
 397         try:
 398             blastmap(seq=seq, queryDB=self.prot)
 399             assert 0, "should fail before this"
 400         except ValueError:
 401             pass
 402
 403     def test_translation_db_in_results_of_db_search(self):
 404         """
 405         Test that the NLMSA in a BlastxMapping properly picks up the
 406         translationDB from the query sequence dict.
 407         """
 408         blastmap = blast.BlastxMapping(self.prot, verbose=False)
 409         results = blastmap(queryDB=self.dna)
 410
 411         tdb = translationDB.get_translation_db(self.dna)
 412         assert tdb.annodb in results.seqDict.dicts
 413
 414     def test_translation_db_in_results_of_seq_search(self):
 415         """
 416         Test that the NLMSA in a BlastxMapping properly picks up the
 417         translationDB from a single input sequence.
 418         """
 419         blastmap = blast.BlastxMapping(self.prot, verbose=False)
 420
 421         query_seq = self.dna['gi|171854975|dbj|AB364477.1|']
 422         results = blastmap(seq=query_seq)
 423
 424         tdb = translationDB.get_translation_db(self.dna)
 425         assert tdb.annodb in results.seqDict.dicts
 426
 427     def test_translated_seqs_in_results(self):
 428         """
 429         Only NLMSASlices for the query sequence should show up in
 430         BlastxMapping.__getitem__, right?
 431         """
 432         blastmap = blast.BlastxMapping(self.prot, verbose=False)
 433
 434         query_seq = self.dna['gi|171854975|dbj|AB364477.1|']
 435         results = blastmap[query_seq]
 436
 437         tdb = translationDB.get_translation_db(self.dna)
 438         annodb = tdb.annodb
 439
 440         for slice in results:
 441             assert slice.seq.id in annodb, '%s not in annodb!' % slice.seq.id
 442
 443     def test_non_consumable_results(self):
 444         blastmap = blast.BlastxMapping(self.prot, verbose=False)
 445
 446         query_seq = self.dna['gi|171854975|dbj|AB364477.1|']
 447         results = blastmap[query_seq]
 448
 449         x = list(results)
 450         y = list(results)
 451
 452         assert len(x), x
 453         assert x == y, "BlastxMapping.__getitem__ should return list"
 454
 455 class Tblastn_Test(BlastBase):
 456     def test_tblastn(self):
 457         "tblastn test"
 458         if not testutil.blast_enabled():
 459             raise SkipTest, "no BLAST installed"
 460
 461         blastmap = blast.BlastMapping(self.dna, verbose=False)
 462         correct = [(144, 144, 432, 0.451)]
 463
 464         result = blastmap[self.prot['HBB1_XENLA']]
 465         check_results_relaxed_blastx([result], correct,
 466                       lambda t:(len(t[1]), len(t[0]), len(t[1].sequence),
 467                                 t[2].pIdentity()))
 468
 469     def test_tblastn_no_blastx(self):
 470         if not testutil.blast_enabled():
 471             raise SkipTest, "no BLAST installed"
 472
 473         blastmap = blast.BlastMapping(self.prot)
 474         try:
 475             results = blastmap[self.dna['gi|171854975|dbj|AB364477.1|']]
 476             raise AssertionError('failed to trap blastx in BlastMapping')
 477         except ValueError:
 478             pass
 479
 480     def test_megablast(self):
 481         '''test megablast'''
 482         if not testutil.blast_enabled():
 483             raise SkipTest, "no BLAST installed"
 484
 485         blastmap = blast.MegablastMapping(self.dna, verbose=False)
 486         # must use copy of sequence to get "self matches" from NLMSA...
 487         query = seqdb.Sequence(str(self.dna['gi|171854975|dbj|AB364477.1|']),
 488                                'foo')
 489         try:
 490             result = blastmap[query]
 491         except OSError: # silently ignore missing RepeatMasker, megablast
 492             return
 493         found = [(len(t[0]), len(t[1])) for t in result.edges()]
 494         assert found == [(444, 444)]
 495
 496     def test_megablast_repr(self):
 497         blastmap = blast.MegablastMapping(self.dna, verbose=False)
 498         assert '<MegablastMapping' in repr(blastmap)
 499
 500     def test_bad_subject(self):
 501         "Test bad subjects"
 502
 503         correctCoords = ((12, 63, 99508, 99661),
 504                          (65, 96, 99661, 99754),
 505                          (96, 108, 99778, 99814),
 506                          (108, 181, 99826, 100045))
 507
 508         fp = file(testutil.datafile('bad_tblastn.txt'))
 509         try:
 510             p = parse_blast.BlastHitParser()
 511             it = iter(correctCoords)
 512             for ival in p.parse_file(fp):
 513                 if not isinstance(ival, (CoordsGroupStart, CoordsGroupEnd)):
 514                     assert (ival.src_start,ival.src_end,
 515                             ival.dest_start,ival.dest_end) \
 516                         == it.next()
 517         finally:
 518             fp.close()
 519
 520 class BlastParsers_Test(BlastBase):
 521     def test_blastp_parser(self):
 522         "Testing blastp parser"
 523         blastp_output = open(testutil.datafile('blastp_output.txt'), 'r')
 524
 525         seq_dict = { 'HBB1_XENLA' : self.prot['HBB1_XENLA'] }
 526         prot_index = blast.BlastIDIndex(self.prot)
 527         try:
 528             alignment = blast.read_blast_alignment(blastp_output, seq_dict,
 529                                                    prot_index)
 530             results = alignment[self.prot['HBB1_XENLA']]
 531         finally:
 532             blastp_output.close()
 533
 534         check_results([results], blastp_correct_results, pair_identity_tuple)
 535
 536     def test_multiblast_parser(self):
 537         "Testing multiblast parser"
 538         multiblast_output = open(testutil.datafile('multiblast_output.txt'),
 539                                  'r')
 540
 541         try:
 542             al = cnestedlist.NLMSA('blasthits', 'memory', pairwiseMode=True,
 543                                    bidirectional=False)
 544             al = blast.read_blast_alignment(multiblast_output, self.prot,
 545                                             blast.BlastIDIndex(self.prot), al)
 546         finally:
 547             multiblast_output.close()
 548         al.build()
 549         results = [al[seq] for seq in self.prot.values()]
 550
 551         check_results(results, correct_multiblast_results,
 552                       pair_identity_tuple)
 553
 554     def test_multiblast_parser_long(self):
 555         "Testing multiblast parser with long input"
 556         longerFile = testutil.datafile('sp_all_hbb')
 557         sp_all_hbb = seqdb.SequenceFileDB(longerFile)
 558
 559         filename = testutil.datafile('multiblast_long_output.txt')
 560         multiblast_output = open(filename, 'r')
 561         try:
 562             al = cnestedlist.NLMSA('blasthits', 'memory', pairwiseMode=True,
 563                                    bidirectional=False)
 564             al = blast.read_blast_alignment(multiblast_output, sp_all_hbb,
 565                                             self.prot, al)
 566         finally:
 567             multiblast_output.close()
 568         al.build()
 569
 570         results = []
 571         for seq in sp_all_hbb.values():
 572             try:
 573                 results.append(al[seq])
 574             except KeyError:
 575                 pass
 576         correctfile = file(testutil.datafile('multiblast_long_correct.txt'),
 577                            'r')
 578         try:
 579             correct = []
 580             for line in correctfile:
 581                 t = line.split()
 582                 correct.append((t[0], t[1], float(t[2])))
 583         finally:
 584             correctfile.close()
 585         check_results(results, correct, pair_identity_tuple)
 586
 587     def test_blastx_parser(self):
 588         "Testing blastx parser"
 589         blastx_output = open(testutil.datafile('blastx_output.txt'), 'r')
 590         seq_dict =  { 'gi|171854975|dbj|AB364477.1|' :
 591                       self.dna['gi|171854975|dbj|AB364477.1|'] }
 592         try:
 593             results = blast.read_blast_alignment(blastx_output,
 594                                                  seq_dict,
 595                                                  blast.BlastIDIndex(self.prot),
 596                                                  translateSrc=True)
 597         finally:
 598             blastx_output.close()
 599         correct = [(143, 143, 429, 0.53146853146853146),
 600                    (143, 145, 429, 0.28275862068965518),
 601                    (143, 145, 429, 0.28965517241379313),
 602                    (143, 145, 429, 0.29655172413793102),
 603                    (143, 145, 429, 0.30344827586206896),
 604                    (144, 144, 432, 0.4513888888888889),
 605                    (144, 144, 432, 0.4513888888888889),
 606                    (145, 145, 435, 0.45517241379310347),
 607                    (145, 145, 435, 0.51034482758620692),
 608                    (146, 142, 438, 0.35616438356164382),
 609                    (146, 146, 438, 0.4589041095890411),
 610                    (146, 146, 438, 0.46575342465753422),
 611                    (146, 146, 438, 0.4726027397260274),
 612                    (146, 146, 438, 0.4726027397260274),
 613                    (146, 146, 438, 0.4863013698630137),
 614                    (146, 146, 438, 0.59589041095890416),
 615                    (146, 146, 438, 0.62328767123287676),
 616                    (146, 146, 438, 0.66438356164383561),
 617                    (146, 146, 438, 0.74657534246575341),
 618                    (146, 146, 438, 0.91095890410958902),
 619                    (146, 146, 438, 0.97945205479452058)]
 620
 621         check_results([results], correct,
 622                       lambda t:(len(t[0]), len(t[1]), len(t[0].sequence),
 623                                 t[2].pIdentity()))
 624
 625     def test_tblastn_parser(self):
 626         "Testing tblastn parser"
 627         seq_dict = { 'HBB1_XENLA' : self.prot['HBB1_XENLA'] }
 628         dna_db = blast.BlastIDIndex(self.dna)
 629         tblastn_output = open(testutil.datafile('tblastn_output.txt'), 'r')
 630         try:
 631             al = blast.read_blast_alignment(tblastn_output, seq_dict,
 632                                             dna_db, translateDest=True)
 633             result = al[self.prot['HBB1_XENLA']]
 634         finally:
 635             tblastn_output.close()
 636         src, dest, edge = iter(result.edges()).next()
 637
 638         self.assertEqual(str(src),
 639             'LTAHDRQLINSTWGKLCAKTIGQEALGRLLWTYPWTQRYFSSFGNLNSADAVFHNEAVAAHGEK'
 640             'VVTSIGEAIKHMDDIKGYYAQLSKYHSETLHVDPLNFKRFGGCLSIALARHFHEEYTPELHAAY'
 641             'EHLFDAIADALGKGYH')
 642         self.assertEqual(str(dest),
 643             'LTDAEKAAVSGLWGKVNSDEVGGEALGRLLVVYPWTQRYFDSFGDLSSASAIMGNAKVKAHGKK'
 644             'VITAFNEGLNHLDSLKGTFASLSELHCDKLHVDPENFRLLGNMIVIVLGHHLGKDFTPAAQAAF'
 645             'QKVMAGVATALAHKYH')
 646         self.assertEqual(str(dest.sequence),
 647             'CTGACTGATGCTGAGAAGGCTGCTGTCTCTGGCCTGTGGGGAAAGGTGAACTCCGATGAAGTTG'
 648             'GTGGTGAGGCCCTGGGCAGGCTGCTGGTTGTCTACCCTTGGACCCAGAGGTACTTTGATAGCTT'
 649             'TGGAGACCTATCCTCTGCCTCTGCTATCATGGGTAATGCCAAAGTGAAGGCCCATGGCAAGAAA'
 650             'GTGATAACTGCCTTTAACGAGGGCCTGAATCACTTGGACAGCCTCAAGGGCACCTTTGCCAGCC'
 651             'TCAGTGAGCTCCACTGTGACAAGCTCCATGTGGATCCTGAGAACTTCAGGCTCCTGGGCAATAT'
 652             'GATCGTGATTGTGCTGGGCCACCACCTGGGCAAGGATTTCACCCCCGCTGCACAGGCTGCCTTC'
 653             'CAGAAGGTGATGGCTGGAGTGGCCACTGCCCTGGCTCACAAGTACCAC')
 654
 655         self.assertAlmostEqual(edge.pIdentity(), 0.451, 3)
 656
 657
 658 # not used currently
 659 def all_vs_all_blast_save():
 660     """
 661     Creates the blast files used during testing.
 662     Must be called before running the tests
 663     """
 664
 665     tempdir = testutil.TempDir('blast-test')
 666     testutil.change_pygrdatapath(tempdir.path)
 667
 668     sp_hbb1 = testutil.datafile('sp_hbb1')
 669     all_vs_all = testutil.tempdatafile('all_vs_all')
 670
 671     sp = seqdb.BlastDB(sp_hbb1)
 672     msa = cnestedlist.NLMSA(all_vs_all ,mode='w', pairwiseMode=True, bidirectional=False)
 673
 674     # get strong homologs, save alignment in msa for every sequence
 675     reader = islice(sp.iteritems(), None)
 676     for id, s in reader:
 677         sp.blast(s, msa, expmax=1e-10, verbose=False)
 678
 679     # done constructing the alignment, so build the alignment db indexes
 680     msa.build(saveSeqDict=True)
 681
 682     db = msa.seqDict.dicts.keys()[0]
 683     working, result = {}, {}
 684     for k in db.values():
 685         edges = msa[k].edges(minAlignSize=12, pIdentityMin=0.5)
 686         for t in edges:
 687             assert len(t[0]) >= 12
 688         tmpdict = dict(map(lambda x:(x, None), [(str(t[0]), str(t[1]), t[2].pIdentity(trapOverflow=False)) for t in edges]))
 689         result[repr(k)] = tmpdict.keys()
 690         result[repr(k)].sort()
 691
 692     # save it into worldbase
 693     data = testutil.TestData()
 694     data.__doc__ = 'sp_allvall'
 695     data.result = result
 696     worldbase.Bio.Blast = data
 697     worldbase.commit()
 698
 699     #return msa
 700
 701 ###
 702
 703 blastp_correct_results = \
 704         [('HBB1_XENLA', 'HBB0_PAGBO', 0.44055944055944057),
 705          ('HBB1_XENLA', 'HBB1_ANAMI', 0.45323741007194246),
 706          ('HBB1_XENLA', 'HBB1_CYGMA', 0.46715328467153283),
 707          ('HBB1_XENLA', 'HBB1_IGUIG', 0.48951048951048953),
 708          ('HBB1_XENLA', 'HBB1_MOUSE', 0.44444444444444442),
 709          ('HBB1_XENLA', 'HBB1_ONCMY', 0.39436619718309857),
 710          ('HBB1_XENLA', 'HBB1_PAGBO', 0.44055944055944057),
 711          ('HBB1_XENLA', 'HBB1_RAT', 0.45833333333333331),
 712          ('HBB1_XENLA', 'HBB1_SPHPU', 0.4825174825174825),
 713          ('HBB1_XENLA', 'HBB1_TAPTE', 0.47222222222222221),
 714          ('HBB1_XENLA', 'HBB1_TORMA', 0.33793103448275863),
 715          ('HBB1_XENLA', 'HBB1_TRICR', 0.49305555555555558),
 716          ('HBB1_XENLA', 'HBB1_UROHA', 0.3776223776223776),
 717          ('HBB1_XENLA', 'HBB1_VAREX', 0.5174825174825175),
 718          ('HBB1_XENLA', 'HBB1_XENBO', 0.96551724137931039),
 719          ('HBB1_XENLA', 'HBB1_XENLA', 1.0),
 720          ('HBB1_XENLA', 'HBB1_XENTR', 0.75),
 721          ('HBB1_XENLA', 'MYG_DIDMA', 0.32075471698113206),
 722          ('HBB1_XENLA', 'MYG_ELEMA', 0.26415094339622641),
 723          ('HBB1_XENLA', 'MYG_ERIEU', 0.32075471698113206),
 724          ('HBB1_XENLA', 'MYG_ESCGI', 0.28282828282828282),
 725          ('HBB1_XENLA', 'MYG_GALCR', 0.32075471698113206)]
 726
 727
 728
 729
 730
 731 correct_multiblast_results = \
 732         [('HBB0_PAGBO', 'HBB0_PAGBO', 1.0),
 733          ('HBB0_PAGBO', 'HBB1_ANAMI', 0.66896551724137931),
 734          ('HBB0_PAGBO', 'HBB1_CYGMA', 0.68493150684931503),
 735          ('HBB0_PAGBO', 'HBB1_IGUIG', 0.4863013698630137),
 736          ('HBB0_PAGBO', 'HBB1_MOUSE', 0.45205479452054792),
 737          ('HBB0_PAGBO', 'HBB1_ONCMY', 0.55172413793103448),
 738          ('HBB0_PAGBO', 'HBB1_PAGBO', 0.69178082191780821),
 739          ('HBB0_PAGBO', 'HBB1_RAT', 0.4589041095890411),
 740          ('HBB0_PAGBO', 'HBB1_SPHPU', 0.4589041095890411),
 741          ('HBB0_PAGBO', 'HBB1_TAPTE', 0.4863013698630137),
 742          ('HBB0_PAGBO', 'HBB1_TORMA', 0.31506849315068491),
 743          ('HBB0_PAGBO', 'HBB1_TRICR', 0.4375),
 744          ('HBB0_PAGBO', 'HBB1_UROHA', 0.4041095890410959),
 745          ('HBB0_PAGBO', 'HBB1_VAREX', 0.49315068493150682),
 746          ('HBB0_PAGBO', 'HBB1_XENBO', 0.43356643356643354),
 747          ('HBB0_PAGBO', 'HBB1_XENLA', 0.44055944055944057),
 748          ('HBB0_PAGBO', 'HBB1_XENTR', 0.4726027397260274),
 749          ('HBB0_PAGBO', 'MYG_DIDMA', 0.22222222222222221),
 750          ('HBB0_PAGBO', 'MYG_ELEMA', 0.20833333333333334),
 751          ('HBB0_PAGBO', 'MYG_ERIEU', 0.21527777777777779),
 752          ('HBB0_PAGBO', 'MYG_ESCGI', 0.25),
 753          ('HBB0_PAGBO', 'MYG_GALCR', 0.24305555555555555),
 754          ('HBB1_ANAMI', 'HBB0_PAGBO', 0.66896551724137931),
 755          ('HBB1_ANAMI', 'HBB1_ANAMI', 1.0),
 756          ('HBB1_ANAMI', 'HBB1_CYGMA', 0.75862068965517238),
 757          ('HBB1_ANAMI', 'HBB1_IGUIG', 0.47586206896551725),
 758          ('HBB1_ANAMI', 'HBB1_MOUSE', 0.45517241379310347),
 759          ('HBB1_ANAMI', 'HBB1_ONCMY', 0.59310344827586203),
 760          ('HBB1_ANAMI', 'HBB1_PAGBO', 0.75862068965517238),
 761          ('HBB1_ANAMI', 'HBB1_RAT', 0.48965517241379308),
 762          ('HBB1_ANAMI', 'HBB1_SPHPU', 0.46206896551724136),
 763          ('HBB1_ANAMI', 'HBB1_TAPTE', 0.48965517241379308),
 764          ('HBB1_ANAMI', 'HBB1_TORMA', 0.32413793103448274),
 765          ('HBB1_ANAMI', 'HBB1_TRICR', 0.41258741258741261),
 766          ('HBB1_ANAMI', 'HBB1_UROHA', 0.38620689655172413),
 767          ('HBB1_ANAMI', 'HBB1_VAREX', 0.48275862068965519),
 768          ('HBB1_ANAMI', 'HBB1_XENBO', 0.4460431654676259),
 769          ('HBB1_ANAMI', 'HBB1_XENLA', 0.45323741007194246),
 770          ('HBB1_ANAMI', 'HBB1_XENTR', 0.4689655172413793),
 771          ('HBB1_CYGMA', 'HBB0_PAGBO', 0.68493150684931503),
 772          ('HBB1_CYGMA', 'HBB1_ANAMI', 0.75862068965517238),
 773          ('HBB1_CYGMA', 'HBB1_CYGMA', 1.0),
 774          ('HBB1_CYGMA', 'HBB1_IGUIG', 0.5),
 775          ('HBB1_CYGMA', 'HBB1_MOUSE', 0.47945205479452052),
 776          ('HBB1_CYGMA', 'HBB1_ONCMY', 0.53103448275862064),
 777          ('HBB1_CYGMA', 'HBB1_PAGBO', 0.86986301369863017),
 778          ('HBB1_CYGMA', 'HBB1_RAT', 0.50684931506849318),
 779          ('HBB1_CYGMA', 'HBB1_SPHPU', 0.47945205479452052),
 780          ('HBB1_CYGMA', 'HBB1_TAPTE', 0.4726027397260274),
 781          ('HBB1_CYGMA', 'HBB1_TORMA', 0.33561643835616439),
 782          ('HBB1_CYGMA', 'HBB1_TRICR', 0.4375),
 783          ('HBB1_CYGMA', 'HBB1_UROHA', 0.36986301369863012),
 784          ('HBB1_CYGMA', 'HBB1_VAREX', 0.4863013698630137),
 785          ('HBB1_CYGMA', 'HBB1_XENBO', 0.45985401459854014),
 786          ('HBB1_CYGMA', 'HBB1_XENLA', 0.46715328467153283),
 787          ('HBB1_CYGMA', 'HBB1_XENTR', 0.47945205479452052),
 788          ('HBB1_CYGMA', 'MYG_ESCGI', 0.2361111111111111),
 789          ('HBB1_IGUIG', 'HBB0_PAGBO', 0.4863013698630137),
 790          ('HBB1_IGUIG', 'HBB1_ANAMI', 0.47586206896551725),
 791          ('HBB1_IGUIG', 'HBB1_CYGMA', 0.5),
 792          ('HBB1_IGUIG', 'HBB1_IGUIG', 1.0),
 793          ('HBB1_IGUIG', 'HBB1_MOUSE', 0.63013698630136983),
 794          ('HBB1_IGUIG', 'HBB1_ONCMY', 0.51034482758620692),
 795          ('HBB1_IGUIG', 'HBB1_PAGBO', 0.4863013698630137),
 796          ('HBB1_IGUIG', 'HBB1_RAT', 0.61643835616438358),
 797          ('HBB1_IGUIG', 'HBB1_SPHPU', 0.71232876712328763),
 798          ('HBB1_IGUIG', 'HBB1_TAPTE', 0.64383561643835618),
 799          ('HBB1_IGUIG', 'HBB1_TORMA', 0.36301369863013699),
 800          ('HBB1_IGUIG', 'HBB1_TRICR', 0.47916666666666669),
 801          ('HBB1_IGUIG', 'HBB1_UROHA', 0.64383561643835618),
 802          ('HBB1_IGUIG', 'HBB1_VAREX', 0.77397260273972601),
 803          ('HBB1_IGUIG', 'HBB1_XENBO', 0.4825174825174825),
 804          ('HBB1_IGUIG', 'HBB1_XENLA', 0.48951048951048953),
 805          ('HBB1_IGUIG', 'HBB1_XENTR', 0.49315068493150682),
 806          ('HBB1_IGUIG', 'MYG_DIDMA', 0.25179856115107913),
 807          ('HBB1_IGUIG', 'MYG_ERIEU', 0.28368794326241137),
 808          ('HBB1_IGUIG', 'MYG_ESCGI', 0.27659574468085107),
 809          ('HBB1_IGUIG', 'MYG_GALCR', 0.28368794326241137),
 810          ('HBB1_MOUSE', 'HBB0_PAGBO', 0.45205479452054792),
 811          ('HBB1_MOUSE', 'HBB1_ANAMI', 0.45517241379310347),
 812          ('HBB1_MOUSE', 'HBB1_CYGMA', 0.47945205479452052),
 813          ('HBB1_MOUSE', 'HBB1_IGUIG', 0.63013698630136983),
 814          ('HBB1_MOUSE', 'HBB1_MOUSE', 1.0),
 815          ('HBB1_MOUSE', 'HBB1_ONCMY', 0.50344827586206897),
 816          ('HBB1_MOUSE', 'HBB1_PAGBO', 0.4726027397260274),
 817          ('HBB1_MOUSE', 'HBB1_RAT', 0.9178082191780822),
 818          ('HBB1_MOUSE', 'HBB1_SPHPU', 0.65753424657534243),
 819          ('HBB1_MOUSE', 'HBB1_TAPTE', 0.76027397260273977),
 820          ('HBB1_MOUSE', 'HBB1_TORMA', 0.35616438356164382),
 821          ('HBB1_MOUSE', 'HBB1_TRICR', 0.52083333333333337),
 822          ('HBB1_MOUSE', 'HBB1_UROHA', 0.47945205479452052),
 823          ('HBB1_MOUSE', 'HBB1_VAREX', 0.6095890410958904),
 824          ('HBB1_MOUSE', 'HBB1_XENBO', 0.44444444444444442),
 825          ('HBB1_MOUSE', 'HBB1_XENLA', 0.44444444444444442),
 826          ('HBB1_MOUSE', 'HBB1_XENTR', 0.4589041095890411),
 827          ('HBB1_MOUSE', 'MYG_DIDMA', 0.29655172413793102),
 828          ('HBB1_MOUSE', 'MYG_ELEMA', 0.27586206896551724),
 829          ('HBB1_MOUSE', 'MYG_ERIEU', 0.30344827586206896),
 830          ('HBB1_MOUSE', 'MYG_ESCGI', 0.28965517241379313),
 831          ('HBB1_MOUSE', 'MYG_GALCR', 0.28275862068965518),
 832          ('HBB1_ONCMY', 'HBB0_PAGBO', 0.55172413793103448),
 833          ('HBB1_ONCMY', 'HBB1_ANAMI', 0.59310344827586203),
 834          ('HBB1_ONCMY', 'HBB1_CYGMA', 0.53103448275862064),
 835          ('HBB1_ONCMY', 'HBB1_IGUIG', 0.51034482758620692),
 836          ('HBB1_ONCMY', 'HBB1_MOUSE', 0.50344827586206897),
 837          ('HBB1_ONCMY', 'HBB1_ONCMY', 1.0),
 838          ('HBB1_ONCMY', 'HBB1_PAGBO', 0.56551724137931036),
 839          ('HBB1_ONCMY', 'HBB1_RAT', 0.50344827586206897),
 840          ('HBB1_ONCMY', 'HBB1_SPHPU', 0.46206896551724136),
 841          ('HBB1_ONCMY', 'HBB1_TAPTE', 0.50344827586206897),
 842          ('HBB1_ONCMY', 'HBB1_TORMA', 0.33793103448275863),
 843          ('HBB1_ONCMY', 'HBB1_TRICR', 0.41258741258741261),
 844          ('HBB1_ONCMY', 'HBB1_UROHA', 0.44827586206896552),
 845          ('HBB1_ONCMY', 'HBB1_VAREX', 0.48965517241379308),
 846          ('HBB1_ONCMY', 'HBB1_XENBO', 0.40140845070422537),
 847          ('HBB1_ONCMY', 'HBB1_XENLA', 0.39436619718309857),
 848          ('HBB1_ONCMY', 'HBB1_XENTR', 0.39310344827586208),
 849          ('HBB1_ONCMY', 'MYG_DIDMA', 0.25694444444444442),
 850          ('HBB1_ONCMY', 'MYG_ERIEU', 0.2361111111111111),
 851          ('HBB1_ONCMY', 'MYG_ESCGI', 0.25),
 852          ('HBB1_ONCMY', 'MYG_GALCR', 0.24305555555555555),
 853          ('HBB1_PAGBO', 'HBB0_PAGBO', 0.69178082191780821),
 854          ('HBB1_PAGBO', 'HBB1_ANAMI', 0.75862068965517238),
 855          ('HBB1_PAGBO', 'HBB1_CYGMA', 0.86986301369863017),
 856          ('HBB1_PAGBO', 'HBB1_IGUIG', 0.4863013698630137),
 857          ('HBB1_PAGBO', 'HBB1_MOUSE', 0.4726027397260274),
 858          ('HBB1_PAGBO', 'HBB1_ONCMY', 0.56551724137931036),
 859          ('HBB1_PAGBO', 'HBB1_PAGBO', 1.0),
 860          ('HBB1_PAGBO', 'HBB1_RAT', 0.4863013698630137),
 861          ('HBB1_PAGBO', 'HBB1_SPHPU', 0.4726027397260274),
 862          ('HBB1_PAGBO', 'HBB1_TAPTE', 0.46575342465753422),
 863          ('HBB1_PAGBO', 'HBB1_TORMA', 0.34931506849315069),
 864          ('HBB1_PAGBO', 'HBB1_TRICR', 0.4375),
 865          ('HBB1_PAGBO', 'HBB1_UROHA', 0.35616438356164382),
 866          ('HBB1_PAGBO', 'HBB1_VAREX', 0.4726027397260274),
 867          ('HBB1_PAGBO', 'HBB1_XENBO', 0.43356643356643354),
 868          ('HBB1_PAGBO', 'HBB1_XENLA', 0.44055944055944057),
 869          ('HBB1_PAGBO', 'HBB1_XENTR', 0.47945205479452052),
 870          ('HBB1_RAT', 'HBB0_PAGBO', 0.4589041095890411),
 871          ('HBB1_RAT', 'HBB1_ANAMI', 0.48965517241379308),
 872          ('HBB1_RAT', 'HBB1_CYGMA', 0.50684931506849318),
 873          ('HBB1_RAT', 'HBB1_IGUIG', 0.61643835616438358),
 874          ('HBB1_RAT', 'HBB1_MOUSE', 0.9178082191780822),
 875          ('HBB1_RAT', 'HBB1_ONCMY', 0.50344827586206897),
 876          ('HBB1_RAT', 'HBB1_PAGBO', 0.4863013698630137),
 877          ('HBB1_RAT', 'HBB1_RAT', 1.0),
 878          ('HBB1_RAT', 'HBB1_SPHPU', 0.66438356164383561),
 879          ('HBB1_RAT', 'HBB1_TAPTE', 0.76712328767123283),
 880          ('HBB1_RAT', 'HBB1_TORMA', 0.35616438356164382),
 881          ('HBB1_RAT', 'HBB1_TRICR', 0.52777777777777779),
 882          ('HBB1_RAT', 'HBB1_UROHA', 0.5),
 883          ('HBB1_RAT', 'HBB1_VAREX', 0.62328767123287676),
 884          ('HBB1_RAT', 'HBB1_XENBO', 0.45833333333333331),
 885          ('HBB1_RAT', 'HBB1_XENLA', 0.45833333333333331),
 886          ('HBB1_RAT', 'HBB1_XENTR', 0.45205479452054792),
 887          ('HBB1_RAT', 'MYG_DIDMA', 0.29655172413793102),
 888          ('HBB1_RAT', 'MYG_ELEMA', 0.28275862068965518),
 889          ('HBB1_RAT', 'MYG_ERIEU', 0.29655172413793102),
 890          ('HBB1_RAT', 'MYG_ESCGI', 0.28275862068965518),
 891          ('HBB1_RAT', 'MYG_GALCR', 0.27586206896551724),
 892          ('HBB1_SPHPU', 'HBB0_PAGBO', 0.4589041095890411),
 893          ('HBB1_SPHPU', 'HBB1_ANAMI', 0.46206896551724136),
 894          ('HBB1_SPHPU', 'HBB1_CYGMA', 0.47945205479452052),
 895          ('HBB1_SPHPU', 'HBB1_IGUIG', 0.71232876712328763),
 896          ('HBB1_SPHPU', 'HBB1_MOUSE', 0.65753424657534243),
 897          ('HBB1_SPHPU', 'HBB1_ONCMY', 0.46206896551724136),
 898          ('HBB1_SPHPU', 'HBB1_PAGBO', 0.4726027397260274),
 899          ('HBB1_SPHPU', 'HBB1_RAT', 0.66438356164383561),
 900          ('HBB1_SPHPU', 'HBB1_SPHPU', 1.0),
 901          ('HBB1_SPHPU', 'HBB1_TAPTE', 0.63698630136986301),
 902          ('HBB1_SPHPU', 'HBB1_TORMA', 0.38356164383561642),
 903          ('HBB1_SPHPU', 'HBB1_TRICR', 0.47916666666666669),
 904          ('HBB1_SPHPU', 'HBB1_UROHA', 0.54109589041095896),
 905          ('HBB1_SPHPU', 'HBB1_VAREX', 0.69178082191780821),
 906          ('HBB1_SPHPU', 'HBB1_XENBO', 0.48951048951048953),
 907          ('HBB1_SPHPU', 'HBB1_XENLA', 0.4825174825174825),
 908          ('HBB1_SPHPU', 'HBB1_XENTR', 0.4726027397260274),
 909          ('HBB1_TAPTE', 'HBB0_PAGBO', 0.4863013698630137),
 910          ('HBB1_TAPTE', 'HBB1_ANAMI', 0.48965517241379308),
 911          ('HBB1_TAPTE', 'HBB1_CYGMA', 0.4726027397260274),
 912          ('HBB1_TAPTE', 'HBB1_IGUIG', 0.64383561643835618),
 913          ('HBB1_TAPTE', 'HBB1_MOUSE', 0.76027397260273977),
 914          ('HBB1_TAPTE', 'HBB1_ONCMY', 0.50344827586206897),
 915          ('HBB1_TAPTE', 'HBB1_PAGBO', 0.46575342465753422),
 916          ('HBB1_TAPTE', 'HBB1_RAT', 0.76712328767123283),
 917          ('HBB1_TAPTE', 'HBB1_SPHPU', 0.63698630136986301),
 918          ('HBB1_TAPTE', 'HBB1_TAPTE', 1.0),
 919          ('HBB1_TAPTE', 'HBB1_TORMA', 0.34931506849315069),
 920          ('HBB1_TAPTE', 'HBB1_TRICR', 0.4861111111111111),
 921          ('HBB1_TAPTE', 'HBB1_UROHA', 0.51369863013698636),
 922          ('HBB1_TAPTE', 'HBB1_VAREX', 0.62328767123287676),
 923          ('HBB1_TAPTE', 'HBB1_XENBO', 0.4861111111111111),
 924          ('HBB1_TAPTE', 'HBB1_XENLA', 0.47222222222222221),
 925          ('HBB1_TAPTE', 'HBB1_XENTR', 0.45205479452054792),
 926          ('HBB1_TAPTE', 'MYG_DIDMA', 0.26277372262773724),
 927          ('HBB1_TAPTE', 'MYG_ERIEU', 0.27007299270072993),
 928          ('HBB1_TAPTE', 'MYG_ESCGI', 0.30344827586206896),
 929          ('HBB1_TAPTE', 'MYG_GALCR', 0.27007299270072993),
 930          ('HBB1_TORMA', 'HBB0_PAGBO', 0.31506849315068491),
 931          ('HBB1_TORMA', 'HBB1_ANAMI', 0.32413793103448274),
 932          ('HBB1_TORMA', 'HBB1_CYGMA', 0.33561643835616439),
 933          ('HBB1_TORMA', 'HBB1_IGUIG', 0.36301369863013699),
 934          ('HBB1_TORMA', 'HBB1_MOUSE', 0.35616438356164382),
 935          ('HBB1_TORMA', 'HBB1_ONCMY', 0.33793103448275863),
 936          ('HBB1_TORMA', 'HBB1_PAGBO', 0.34931506849315069),
 937          ('HBB1_TORMA', 'HBB1_RAT', 0.35616438356164382),
 938          ('HBB1_TORMA', 'HBB1_SPHPU', 0.38356164383561642),
 939          ('HBB1_TORMA', 'HBB1_TAPTE', 0.34931506849315069),
 940          ('HBB1_TORMA', 'HBB1_TORMA', 1.0),
 941          ('HBB1_TORMA', 'HBB1_TRICR', 0.31724137931034485),
 942          ('HBB1_TORMA', 'HBB1_UROHA', 0.29452054794520549),
 943          ('HBB1_TORMA', 'HBB1_VAREX', 0.35616438356164382),
 944          ('HBB1_TORMA', 'HBB1_XENBO', 0.34482758620689657),
 945          ('HBB1_TORMA', 'HBB1_XENLA', 0.33793103448275863),
 946          ('HBB1_TORMA', 'HBB1_XENTR', 0.33561643835616439),
 947          ('HBB1_TORMA', 'MYG_ESCGI', 0.25675675675675674),
 948          ('HBB1_TRICR', 'HBB0_PAGBO', 0.4375),
 949          ('HBB1_TRICR', 'HBB1_ANAMI', 0.41258741258741261),
 950          ('HBB1_TRICR', 'HBB1_CYGMA', 0.4375),
 951          ('HBB1_TRICR', 'HBB1_IGUIG', 0.47916666666666669),
 952          ('HBB1_TRICR', 'HBB1_MOUSE', 0.52083333333333337),
 953          ('HBB1_TRICR', 'HBB1_ONCMY', 0.41258741258741261),
 954          ('HBB1_TRICR', 'HBB1_PAGBO', 0.4375),
 955          ('HBB1_TRICR', 'HBB1_RAT', 0.52777777777777779),
 956          ('HBB1_TRICR', 'HBB1_SPHPU', 0.47916666666666669),
 957          ('HBB1_TRICR', 'HBB1_TAPTE', 0.4861111111111111),
 958          ('HBB1_TRICR', 'HBB1_TORMA', 0.31724137931034485),
 959          ('HBB1_TRICR', 'HBB1_TRICR', 1.0),
 960          ('HBB1_TRICR', 'HBB1_UROHA', 0.3611111111111111),
 961          ('HBB1_TRICR', 'HBB1_VAREX', 0.4513888888888889),
 962          ('HBB1_TRICR', 'HBB1_XENBO', 0.4861111111111111),
 963          ('HBB1_TRICR', 'HBB1_XENLA', 0.49305555555555558),
 964          ('HBB1_TRICR', 'HBB1_XENTR', 0.49305555555555558),
 965          ('HBB1_UROHA', 'HBB0_PAGBO', 0.4041095890410959),
 966          ('HBB1_UROHA', 'HBB1_ANAMI', 0.42857142857142855),
 967          ('HBB1_UROHA', 'HBB1_CYGMA', 0.36986301369863012),
 968          ('HBB1_UROHA', 'HBB1_IGUIG', 0.64383561643835618),
 969          ('HBB1_UROHA', 'HBB1_MOUSE', 0.51666666666666672),
 970          ('HBB1_UROHA', 'HBB1_ONCMY', 0.50420168067226889),
 971          ('HBB1_UROHA', 'HBB1_PAGBO', 0.38333333333333336),
 972          ('HBB1_UROHA', 'HBB1_RAT', 0.54166666666666663),
 973          ('HBB1_UROHA', 'HBB1_SPHPU', 0.54109589041095896),
 974          ('HBB1_UROHA', 'HBB1_TAPTE', 0.55833333333333335),
 975          ('HBB1_UROHA', 'HBB1_TORMA', 0.31034482758620691),
 976          ('HBB1_UROHA', 'HBB1_TRICR', 0.39316239316239315),
 977          ('HBB1_UROHA', 'HBB1_UROHA', 1.0),
 978          ('HBB1_UROHA', 'HBB1_VAREX', 0.59589041095890416),
 979          ('HBB1_UROHA', 'HBB1_XENBO', 0.42608695652173911),
 980          ('HBB1_UROHA', 'HBB1_XENLA', 0.41739130434782606),
 981          ('HBB1_UROHA', 'HBB1_XENTR', 0.40000000000000002),
 982          ('HBB1_UROHA', 'MYG_ERIEU', 0.27927927927927926),
 983          ('HBB1_VAREX', 'HBB0_PAGBO', 0.49315068493150682),
 984          ('HBB1_VAREX', 'HBB1_ANAMI', 0.48275862068965519),
 985          ('HBB1_VAREX', 'HBB1_CYGMA', 0.4863013698630137),
 986          ('HBB1_VAREX', 'HBB1_IGUIG', 0.77397260273972601),
 987          ('HBB1_VAREX', 'HBB1_MOUSE', 0.6095890410958904),
 988          ('HBB1_VAREX', 'HBB1_ONCMY', 0.48965517241379308),
 989          ('HBB1_VAREX', 'HBB1_PAGBO', 0.4726027397260274),
 990          ('HBB1_VAREX', 'HBB1_RAT', 0.62328767123287676),
 991          ('HBB1_VAREX', 'HBB1_SPHPU', 0.69178082191780821),
 992          ('HBB1_VAREX', 'HBB1_TAPTE', 0.62328767123287676),
 993          ('HBB1_VAREX', 'HBB1_TORMA', 0.35616438356164382),
 994          ('HBB1_VAREX', 'HBB1_TRICR', 0.4513888888888889),
 995          ('HBB1_VAREX', 'HBB1_UROHA', 0.59589041095890416),
 996          ('HBB1_VAREX', 'HBB1_VAREX', 1.0),
 997          ('HBB1_VAREX', 'HBB1_XENBO', 0.51048951048951052),
 998          ('HBB1_VAREX', 'HBB1_XENLA', 0.5174825174825175),
 999          ('HBB1_VAREX', 'HBB1_XENTR', 0.4726027397260274),
1000          ('HBB1_VAREX', 'MYG_DIDMA', 0.25531914893617019),
1001          ('HBB1_VAREX', 'MYG_ERIEU', 0.25531914893617019),
1002          ('HBB1_VAREX', 'MYG_ESCGI', 0.24822695035460993),
1003          ('HBB1_VAREX', 'MYG_GALCR', 0.24822695035460993),
1004          ('HBB1_XENBO', 'HBB0_PAGBO', 0.43356643356643354),
1005          ('HBB1_XENBO', 'HBB1_ANAMI', 0.4460431654676259),
1006          ('HBB1_XENBO', 'HBB1_CYGMA', 0.45985401459854014),
1007          ('HBB1_XENBO', 'HBB1_IGUIG', 0.4825174825174825),
1008          ('HBB1_XENBO', 'HBB1_MOUSE', 0.44444444444444442),
1009          ('HBB1_XENBO', 'HBB1_ONCMY', 0.40140845070422537),
1010          ('HBB1_XENBO', 'HBB1_PAGBO', 0.43356643356643354),
1011          ('HBB1_XENBO', 'HBB1_RAT', 0.45833333333333331),
1012          ('HBB1_XENBO', 'HBB1_SPHPU', 0.48951048951048953),
1013          ('HBB1_XENBO', 'HBB1_TAPTE', 0.4861111111111111),
1014          ('HBB1_XENBO', 'HBB1_TORMA', 0.34482758620689657),
1015          ('HBB1_XENBO', 'HBB1_TRICR', 0.4861111111111111),
1016          ('HBB1_XENBO', 'HBB1_UROHA', 0.38461538461538464),
1017          ('HBB1_XENBO', 'HBB1_VAREX', 0.51048951048951052),
1018          ('HBB1_XENBO', 'HBB1_XENBO', 1.0),
1019          ('HBB1_XENBO', 'HBB1_XENLA', 0.96551724137931039),
1020          ('HBB1_XENBO', 'HBB1_XENTR', 0.76388888888888884),
1021          ('HBB1_XENBO', 'MYG_DIDMA', 0.32075471698113206),
1022          ('HBB1_XENBO', 'MYG_ELEMA', 0.27358490566037735),
1023          ('HBB1_XENBO', 'MYG_ERIEU', 0.32075471698113206),
1024          ('HBB1_XENBO', 'MYG_GALCR', 0.32075471698113206),
1025          ('HBB1_XENLA', 'HBB0_PAGBO', 0.44055944055944057),
1026          ('HBB1_XENLA', 'HBB1_ANAMI', 0.45323741007194246),
1027          ('HBB1_XENLA', 'HBB1_CYGMA', 0.46715328467153283),
1028          ('HBB1_XENLA', 'HBB1_IGUIG', 0.48951048951048953),
1029          ('HBB1_XENLA', 'HBB1_MOUSE', 0.44444444444444442),
1030          ('HBB1_XENLA', 'HBB1_ONCMY', 0.39436619718309857),
1031          ('HBB1_XENLA', 'HBB1_PAGBO', 0.44055944055944057),
1032          ('HBB1_XENLA', 'HBB1_RAT', 0.45833333333333331),
1033          ('HBB1_XENLA', 'HBB1_SPHPU', 0.4825174825174825),
1034          ('HBB1_XENLA', 'HBB1_TAPTE', 0.47222222222222221),
1035          ('HBB1_XENLA', 'HBB1_TORMA', 0.33793103448275863),
1036          ('HBB1_XENLA', 'HBB1_TRICR', 0.49305555555555558),
1037          ('HBB1_XENLA', 'HBB1_UROHA', 0.3776223776223776),
1038          ('HBB1_XENLA', 'HBB1_VAREX', 0.5174825174825175),
1039          ('HBB1_XENLA', 'HBB1_XENBO', 0.96551724137931039),
1040          ('HBB1_XENLA', 'HBB1_XENLA', 1.0),
1041          ('HBB1_XENLA', 'HBB1_XENTR', 0.75),
1042          ('HBB1_XENLA', 'MYG_DIDMA', 0.32075471698113206),
1043          ('HBB1_XENLA', 'MYG_ELEMA', 0.26415094339622641),
1044          ('HBB1_XENLA', 'MYG_ERIEU', 0.32075471698113206),
1045          ('HBB1_XENLA', 'MYG_ESCGI', 0.28282828282828282),
1046          ('HBB1_XENLA', 'MYG_GALCR', 0.32075471698113206),
1047          ('HBB1_XENTR', 'HBB0_PAGBO', 0.4726027397260274),
1048          ('HBB1_XENTR', 'HBB1_ANAMI', 0.4689655172413793),
1049          ('HBB1_XENTR', 'HBB1_CYGMA', 0.47945205479452052),
1050          ('HBB1_XENTR', 'HBB1_IGUIG', 0.49315068493150682),
1051          ('HBB1_XENTR', 'HBB1_MOUSE', 0.4589041095890411),
1052          ('HBB1_XENTR', 'HBB1_ONCMY', 0.39310344827586208),
1053          ('HBB1_XENTR', 'HBB1_PAGBO', 0.47945205479452052),
1054          ('HBB1_XENTR', 'HBB1_RAT', 0.45205479452054792),
1055          ('HBB1_XENTR', 'HBB1_SPHPU', 0.4726027397260274),
1056          ('HBB1_XENTR', 'HBB1_TAPTE', 0.45205479452054792),
1057          ('HBB1_XENTR', 'HBB1_TORMA', 0.33561643835616439),
1058          ('HBB1_XENTR', 'HBB1_TRICR', 0.49305555555555558),
1059          ('HBB1_XENTR', 'HBB1_UROHA', 0.35616438356164382),
1060          ('HBB1_XENTR', 'HBB1_VAREX', 0.4726027397260274),
1061          ('HBB1_XENTR', 'HBB1_XENBO', 0.76388888888888884),
1062          ('HBB1_XENTR', 'HBB1_XENLA', 0.75),
1063          ('HBB1_XENTR', 'HBB1_XENTR', 1.0),
1064          ('HBB1_XENTR', 'MYG_DIDMA', 0.2857142857142857),
1065          ('HBB1_XENTR', 'MYG_ERIEU', 0.27067669172932329),
1066          ('HBB1_XENTR', 'MYG_ESCGI', 0.27272727272727271),
1067          ('HBB1_XENTR', 'MYG_GALCR', 0.2781954887218045),
1068          ('MYG_DIDMA', 'HBB0_PAGBO', 0.22222222222222221),
1069          ('MYG_DIDMA', 'HBB1_IGUIG', 0.25179856115107913),
1070          ('MYG_DIDMA', 'HBB1_MOUSE', 0.29655172413793102),
1071          ('MYG_DIDMA', 'HBB1_ONCMY', 0.25694444444444442),
1072          ('MYG_DIDMA', 'HBB1_RAT', 0.29655172413793102),
1073          ('MYG_DIDMA', 'HBB1_TAPTE', 0.26277372262773724),
1074          ('MYG_DIDMA', 'HBB1_VAREX', 0.25531914893617019),
1075          ('MYG_DIDMA', 'HBB1_XENBO', 0.32075471698113206),
1076          ('MYG_DIDMA', 'HBB1_XENLA', 0.32075471698113206),
1077          ('MYG_DIDMA', 'HBB1_XENTR', 0.2857142857142857),
1078          ('MYG_DIDMA', 'MYG_DIDMA', 1.0),
1079          ('MYG_DIDMA', 'MYG_ELEMA', 0.81045751633986929),
1080          ('MYG_DIDMA', 'MYG_ERIEU', 0.87581699346405228),
1081          ('MYG_DIDMA', 'MYG_ESCGI', 0.83552631578947367),
1082          ('MYG_DIDMA', 'MYG_GALCR', 0.83006535947712423),
1083          ('MYG_ELEMA', 'HBB0_PAGBO', 0.20833333333333334),
1084          ('MYG_ELEMA', 'HBB1_MOUSE', 0.27586206896551724),
1085          ('MYG_ELEMA', 'HBB1_RAT', 0.28275862068965518),
1086          ('MYG_ELEMA', 'HBB1_XENBO', 0.27358490566037735),
1087          ('MYG_ELEMA', 'HBB1_XENLA', 0.26415094339622641),
1088          ('MYG_ELEMA', 'MYG_DIDMA', 0.81045751633986929),
1089          ('MYG_ELEMA', 'MYG_ELEMA', 1.0),
1090          ('MYG_ELEMA', 'MYG_ERIEU', 0.82352941176470584),
1091          ('MYG_ELEMA', 'MYG_ESCGI', 0.83552631578947367),
1092          ('MYG_ELEMA', 'MYG_GALCR', 0.84313725490196079),
1093          ('MYG_ERIEU', 'HBB0_PAGBO', 0.21527777777777779),
1094          ('MYG_ERIEU', 'HBB1_IGUIG', 0.28368794326241137),
1095          ('MYG_ERIEU', 'HBB1_MOUSE', 0.30344827586206896),
1096          ('MYG_ERIEU', 'HBB1_ONCMY', 0.2361111111111111),
1097          ('MYG_ERIEU', 'HBB1_RAT', 0.29655172413793102),
1098          ('MYG_ERIEU', 'HBB1_TAPTE', 0.27007299270072993),
1099          ('MYG_ERIEU', 'HBB1_UROHA', 0.27927927927927926),
1100          ('MYG_ERIEU', 'HBB1_VAREX', 0.25531914893617019),
1101          ('MYG_ERIEU', 'HBB1_XENBO', 0.32075471698113206),
1102          ('MYG_ERIEU', 'HBB1_XENLA', 0.32075471698113206),
1103          ('MYG_ERIEU', 'HBB1_XENTR', 0.27067669172932329),
1104          ('MYG_ERIEU', 'MYG_DIDMA', 0.87581699346405228),
1105          ('MYG_ERIEU', 'MYG_ELEMA', 0.82352941176470584),
1106          ('MYG_ERIEU', 'MYG_ERIEU', 1.0),
1107          ('MYG_ERIEU', 'MYG_ESCGI', 0.83552631578947367),
1108          ('MYG_ERIEU', 'MYG_GALCR', 0.85620915032679734),
1109          ('MYG_ESCGI', 'HBB0_PAGBO', 0.25),
1110          ('MYG_ESCGI', 'HBB1_CYGMA', 0.2361111111111111),
1111          ('MYG_ESCGI', 'HBB1_IGUIG', 0.27659574468085107),
1112          ('MYG_ESCGI', 'HBB1_MOUSE', 0.28965517241379313),
1113          ('MYG_ESCGI', 'HBB1_ONCMY', 0.25),
1114          ('MYG_ESCGI', 'HBB1_RAT', 0.28275862068965518),
1115          ('MYG_ESCGI', 'HBB1_TAPTE', 0.3611111111111111),
1116          ('MYG_ESCGI', 'HBB1_TORMA', 0.25675675675675674),
1117          ('MYG_ESCGI', 'HBB1_VAREX', 0.24822695035460993),
1118          ('MYG_ESCGI', 'HBB1_XENLA', 0.28282828282828282),
1119          ('MYG_ESCGI', 'HBB1_XENTR', 0.27272727272727271),
1120          ('MYG_ESCGI', 'MYG_DIDMA', 0.83552631578947367),
1121          ('MYG_ESCGI', 'MYG_ELEMA', 0.83552631578947367),
1122          ('MYG_ESCGI', 'MYG_ERIEU', 0.83552631578947367),
1123          ('MYG_ESCGI', 'MYG_ESCGI', 1.0),
1124          ('MYG_ESCGI', 'MYG_GALCR', 0.84210526315789469),
1125          ('MYG_GALCR', 'HBB0_PAGBO', 0.24305555555555555),
1126          ('MYG_GALCR', 'HBB1_IGUIG', 0.28368794326241137),
1127          ('MYG_GALCR', 'HBB1_MOUSE', 0.28275862068965518),
1128          ('MYG_GALCR', 'HBB1_ONCMY', 0.24305555555555555),
1129          ('MYG_GALCR', 'HBB1_RAT', 0.27586206896551724),
1130          ('MYG_GALCR', 'HBB1_TAPTE', 0.27007299270072993),
1131          ('MYG_GALCR', 'HBB1_VAREX', 0.24822695035460993),
1132          ('MYG_GALCR', 'HBB1_XENBO', 0.32075471698113206),
1133          ('MYG_GALCR', 'HBB1_XENLA', 0.32075471698113206),
1134          ('MYG_GALCR', 'HBB1_XENTR', 0.2781954887218045),
1135          ('MYG_GALCR', 'MYG_DIDMA', 0.83006535947712423),
1136          ('MYG_GALCR', 'MYG_ELEMA', 0.84313725490196079),
1137          ('MYG_GALCR', 'MYG_ERIEU', 0.85620915032679734),
1138          ('MYG_GALCR', 'MYG_ESCGI', 0.84210526315789469),
1139          ('MYG_GALCR', 'MYG_GALCR', 1.0),
1140          ('PRCA_ANASP', 'PRCA_ANASP', 1.0),
1141          ('PRCA_ANASP', 'PRCA_ANAVA', 0.97222222222222221),
1142          ('PRCA_ANAVA', 'PRCA_ANASP', 0.97222222222222221),
1143          ('PRCA_ANAVA', 'PRCA_ANAVA', 1.0)]
1144
1145 ###
1146
1147 if __name__ == '__main__':
1148     PygrTestProgram(verbosity=2)
1149