tests/tabix_test.py

   1 #!/usr/bin/env python
   2 '''unit testing code for pysam.
   3
   4 Execute in the :file:`tests` directory as it requires the Makefile
   5 and data files located there.
   6 '''
   7
   8 import sys, os, shutil, gzip
   9 import pysam
  10 import unittest
  11 import itertools
  12 import subprocess
  13
  14 def checkBinaryEqual( filename1, filename2 ):
  15     '''return true if the two files are binary equal.'''
  16     if os.path.getsize( filename1 ) !=  os.path.getsize( filename2 ):
  17         return False
  18
  19     infile1 = open(filename1, "rb")
  20     infile2 = open(filename2, "rb")
  21
  22     def chariter( infile ):
  23         while 1:
  24             c = infile.read(1)
  25             if c == "": break
  26             yield c
  27
  28     found = False
  29     for c1,c2 in itertools.izip( chariter( infile1), chariter( infile2) ):
  30         if c1 != c2: break
  31     else:
  32         found = True
  33
  34     infile1.close()
  35     infile2.close()
  36     return found
  37
  38 class TestIndexing(unittest.TestCase):
  39     filename = "example.gtf.gz"
  40     filename_idx = "example.gtf.gz.tbi"
  41
  42     def setUp( self ):
  43
  44         self.tmpfilename = "tmp_%i.gtf.gz" % id(self)
  45         shutil.copyfile( self.filename, self.tmpfilename )
  46
  47     def testIndexPreset( self ):
  48         '''test indexing via preset.'''
  49
  50         pysam.tabix_index( self.tmpfilename, preset = "gff" )
  51         checkBinaryEqual( self.tmpfilename + ".tbi", self.filename_idx )
  52
  53     def tearDown( self ):
  54         os.unlink( self.tmpfilename )
  55         os.unlink( self.tmpfilename + ".tbi" )
  56
  57 class TestCompression(unittest.TestCase):
  58     filename = "example.gtf.gz"
  59     filename_idx = "example.gtf.gz.tbi"
  60
  61     def setUp( self ):
  62
  63         self.tmpfilename = "tmp_%i.gtf" % id(self)
  64         infile = gzip.open( self.filename, "r")
  65         outfile = open( self.tmpfilename, "w" )
  66         outfile.write( "".join(infile.readlines()) )
  67         outfile.close()
  68         infile.close()
  69
  70     def testIndexPreset( self ):
  71         '''test indexing via preset.'''
  72
  73         pysam.tabix_index( self.tmpfilename, preset = "gff" )
  74         checkBinaryEqual( self.tmpfilename + ".gz", self.filename )
  75         checkBinaryEqual( self.tmpfilename + ".gz.tbi", self.filename_idx )
  76
  77     def tearDown( self ):
  78         os.unlink( self.tmpfilename + ".gz" )
  79         os.unlink( self.tmpfilename + ".gz.tbi" )
  80
  81 class TestIteration( unittest.TestCase ):
  82
  83     filename = "example.gtf.gz"
  84
  85     def setUp( self ):
  86
  87         self.tabix = pysam.Tabixfile( self.filename )
  88         lines = gzip.open(self.filename).readlines()
  89         # creates index of contig, start, end, adds content without newline.
  90         self.compare = [
  91             (x[0][0], int(x[0][3]), int(x[0][4]), x[1])
  92             for x in [ (y.split("\t"), y[:-1]) for y in lines ] ]
  93
  94     def getSubset( self, contig = None, start = None, end = None):
  95
  96         if contig == None:
  97             # all lines
  98             subset = [ x[3] for x in self.compare ]
  99         else:
 100             if start != None and end == None:
 101                 # until end of contig
 102                 subset = [ x[3] for x in self.compare if x[0] == contig and x[2] > start ]
 103             elif start == None and end != None:
 104                 # from start of contig
 105                 subset = [ x[3] for x in self.compare if x[0] == contig and x[1] <= end ]
 106             elif start == None and end == None:
 107                 subset = [ x[3] for x in self.compare if x[0] == contig ]
 108             else:
 109                 # all within interval
 110                 subset = [ x[3] for x in self.compare if x[0] == contig and \
 111                                min( x[2], end) - max(x[1], start) > 0 ]
 112
 113         return subset
 114
 115     def checkPairwise( self, result, ref ):
 116
 117         result.sort()
 118         ref.sort()
 119
 120         a = set(result)
 121         b = set(ref)
 122
 123         self.assertEqual( len(result), len(ref),
 124                           "unexpected number of results: %i, expected %i, differences are %s: %s" \
 125                               % (len(result), len(ref),
 126                                  a.difference(b),
 127                                  b.difference(a) ))
 128
 129         for x, d in enumerate( zip( result, ref )):
 130
 131             self.assertEqual( d[0], d[1],
 132                               "unexpected results in pair %i: '%s', expected '%s'" % \
 133                                   (x,
 134                                    d[0],
 135                                    d[1]) )
 136
 137
 138     def testAll( self ):
 139         result = list(self.tabix.fetch())
 140         ref = self.getSubset( )
 141         self.checkPairwise( result, ref )
 142
 143     def testPerContig( self ):
 144         for contig in ("chr1", "chr2", "chr1", "chr2" ):
 145             result = list(self.tabix.fetch( contig ))
 146             ref = self.getSubset( contig )
 147             self.checkPairwise( result, ref )
 148
 149     def testPerContigToEnd( self ):
 150
 151         end = None
 152         for contig in ("chr1", "chr2", "chr1", "chr2" ):
 153             for start in range( 0, 200000, 1000):
 154                 result = list(self.tabix.fetch( contig, start, end ))
 155                 ref = self.getSubset( contig, start, end )
 156                 self.checkPairwise( result, ref )
 157
 158     def testPerContigFromStart( self ):
 159
 160         start = None
 161         for contig in ("chr1", "chr2", "chr1", "chr2" ):
 162             for end in range( 0, 200000, 1000):
 163                 result = list(self.tabix.fetch( contig, start, end ))
 164                 ref = self.getSubset( contig, start, end )
 165                 self.checkPairwise( result, ref )
 166
 167     def testPerContig( self ):
 168
 169         start, end  = None, None
 170         for contig in ("chr1", "chr2", "chr1", "chr2" ):
 171             result = list(self.tabix.fetch( contig, start, end ))
 172             ref = self.getSubset( contig, start, end )
 173             self.checkPairwise( result, ref )
 174
 175     def testPerInterval( self ):
 176
 177         start, end  = None, None
 178         for contig in ("chr1", "chr2", "chr1", "chr2" ):
 179             for start in range( 0, 200000, 2000):
 180                 for end in range( start, start + 2000, 500):
 181                     result = list(self.tabix.fetch( contig, start, end ))
 182                     ref = self.getSubset( contig, start, end )
 183                     self.checkPairwise( result, ref )
 184
 185
 186     def testInvalidIntervals( self ):
 187
 188         self.assertRaises( ValueError, self.tabix.fetch, "chr1", 0, -10)
 189         self.assertRaises( ValueError, self.tabix.fetch, "chr1", -10, 200)
 190         self.assertRaises( ValueError, self.tabix.fetch, "chr1", 200, 0)
 191         self.assertRaises( ValueError, self.tabix.fetch, "chr1", -10, -20)
 192         self.assertRaises( ValueError, self.tabix.fetch, "chrUn" )
 193
 194     def testGetContigs( self ):
 195         self.assertEqual( sorted(self.tabix.contigs), ["chr1", "chr2"] )
 196         # check that contigs is read-only
 197         self.assertRaises( AttributeError, setattr, self.tabix, "contigs", ["chr1", "chr2"] )
 198
 199 class TestParser( unittest.TestCase ):
 200
 201     filename = "example.gtf.gz"
 202
 203     def setUp( self ):
 204
 205         self.tabix = pysam.Tabixfile( self.filename )
 206         self.compare = [ x[:-1].split("\t") for x in gzip.open( self.filename, "r") ]
 207
 208     def testGTF( self ):
 209
 210         for x, r in enumerate(self.tabix.fetch( parser = pysam.asGTF() )):
 211             self.assertEqual( "\t".join( self.compare[x]), str(r) )
 212
 213     def testTuple( self ):
 214
 215         for x, r in enumerate(self.tabix.fetch( parser = pysam.asTuple() )):
 216             self.assertEqual( self.compare[x], list(r) )
 217
 218             self.assertEqual( len(self.compare[x]), len(r) )
 219             for c in range(0,len(r)):
 220                 self.assertEqual( self.compare[x][c], r[c] )
 221
 222 if __name__ == "__main__":
 223     unittest.main()
 224
 225