+ fixed issue with negative strand in lookup table
authorFabio <fabio@congo.fml.local>
Wed, 15 Oct 2008 14:48:25 +0000 (16:48 +0200)
committerFabio <fabio@congo.fml.local>
Wed, 15 Oct 2008 14:48:25 +0000 (16:48 +0200)
README
doc/qpalma-manual.tex
qpalma/Lookup.py
setup.py
tests/test_qpalma.py

diff --git a/README b/README
index f704ae1..121c233 100644 (file)
--- a/README
+++ b/README
@@ -1,4 +1,4 @@
-Installation instructions for QPalma Version 1.0
+Installation instructions for QPalma Version 0.9
   
 The package requires version 2.5 or newer of Python, and is built from 
 source, so the header files and libraries for Python must be installed, 
@@ -25,11 +25,11 @@ Pythongrid and Genefinding are two utility packages from our lab. These can be o
 
 Once you have the above packages installed you can start installing QPalma.
 
-Assuming you downloaded QPalma-1.0.tar.gz type:
+Assuming you downloaded QPalma-0.9.tar.gz type:
 
-   tar -xzvf QPalma-1.0.tar.gz
+   tar -xzvf QPalma-0.9.tar.gz
 
-   cd QPalma-1.0
+   cd QPalma-0.9
 
    python setup.py build
 
@@ -56,4 +56,4 @@ Further information
 Contact
 -------
 
-   fabio@tuebingen.mpg.de
+   qpalma@tuebingen.mpg.de
index 806c3d3..b4d9585 100644 (file)
@@ -67,8 +67,8 @@ For training \QP you need one of the following optimization toolkits:
 \item Install the Pythongrid and the Genefinding tool packages
 \item Update your PYTHONPATH variable to point to the above packages
 \item Unpack the QPalma tarball via
-\item[$\rightarrow$] tar -xzvf QPalma-1.0.tar.gz
-\item Enter the QPalma-1.0 directory and type:
+\item[$\rightarrow$] tar -xzvf QPalma-0.9.tar.gz
+\item Enter the QPalma-0.9 directory and type:
 \item[$\rightarrow$] python setup.py build
 \end{enumerate}
 \noindent
index d9ea7be..baa6419 100644 (file)
@@ -52,24 +52,27 @@ class LookupTable:
       assert chromo in self.chromo_list
       assert start <= stop
 
-      if strand == '-':
-         total_size = self.seqInfo.getFragmentSize(chromo)
-         startP   = total_size - stop
-         stopP    = total_size - start
-         _stop     = stopP
-         _start    = startP
-
       chromo_idx = self.chromo_map[chromo]
 
+
       if strand == '+':
          currentSequence         = self.genomicSequences[chromo_idx][start:stop]
          currentAcceptorScores   = self.acceptorScoresPos[chromo_idx][start:stop]
          currentDonorScores      = self.donorScoresPos[chromo_idx][start:stop]
       elif strand == '-':
-         currentSequence         = self.genomicSequences[chromo_idx][_start:_stop]
+         currentSequence         = self.genomicSequences[chromo_idx][start:stop]
          currentSequence = reverse_complement(currentSequence)
-         currentAcceptorScores   = self.acceptorScoresNeg[chromo_idx][start:stop]
-         currentDonorScores      = self.donorScoresNeg[chromo_idx][start:stop]
+
+         total_size = self.seqInfo.getFragmentSize(chromo)
+         _start = total_size - stop
+         _stop = total_size - start
+         currentAcceptorScores   = self.acceptorScoresNeg[chromo_idx][_start:_stop]
+         currentDonorScores      = self.donorScoresNeg[chromo_idx][_start:_stop]
+
+         #gt_tuple_pos = [p for p,e in enumerate(currentSequence) if (p>0 and p<len(currentSequence)-1 and e=='g' ) and (currentSequence[p+1]=='t' or currentSequence[p+1]=='c')]
+         #ag_tuple_pos = [p for p,e in enumerate(currentSequence) if p>1 and currentSequence[p-1]=='a' and currentSequence[p]=='g']
+         #acc_pos = [p for p,e in enumerate(currentAcceptorScores) if e != -inf and p > 1]
+         #don_pos = [p for p,e in enumerate(currentDonorScores) if e != -inf and p > 0]
 
       return currentSequence, currentAcceptorScores, currentDonorScores
 
@@ -86,9 +89,7 @@ class LookupTable:
 
       print 'lt total_size %d' % self.seqInfo.getFragmentSize(chromo)
 
-      strand = '+'
-      genomicSeq,currentAcc,currentDon = self.seqInfo.get_seq_and_scores(chromo,strand,genomicSeq_start,genomicSeq_stop)
-      #print len(currentAcc),len(currentDon)
+      genomicSeq,currentAcc,currentDon = self.seqInfo.get_seq_and_scores(chromo,'+',genomicSeq_start,genomicSeq_stop)
       genomicSeq = genomicSeq.lower()
 
       chromo_idx = self.chromo_map[chromo]
@@ -97,8 +98,14 @@ class LookupTable:
       self.acceptorScoresPos[chromo_idx]   = currentAcc
       self.donorScoresPos[chromo_idx]      = currentDon
 
-      strand = '-'
-      genomicSeq,currentAcc,currentDon = self.seqInfo.get_seq_and_scores(chromo,strand,genomicSeq_start,genomicSeq_stop)
+      dummy,currentAcc,currentDon = self.seqInfo.get_seq_and_scores(chromo,'-',genomicSeq_start,genomicSeq_stop)
+
+      #genomicSeq = reverse_complement(self.genomicSequences[chromo_idx])
+      #gt_tuple_pos = [p for p,e in enumerate(genomicSeq) if (p>0 and p<len(genomicSeq)-1 and e=='g' ) and (genomicSeq[p+1]=='t' or genomicSeq[p+1]=='c')]
+      #ag_tuple_pos = [p for p,e in enumerate(genomicSeq) if p>1 and genomicSeq[p-1]=='a' and genomicSeq[p]=='g']
+      #acc_pos = [p for p,e in enumerate(currentAcc) if e != -inf and p > 1]
+      #don_pos = [p for p,e in enumerate(currentDon) if e != -inf and p > 0]
+      #pdb.set_trace()
 
       self.acceptorScoresNeg[chromo_idx]   = currentAcc
       self.donorScoresNeg[chromo_idx]      = currentDon
index a628a00..22e923d 100644 (file)
--- a/setup.py
+++ b/setup.py
@@ -149,7 +149,7 @@ class ExtBuilder(_build_py):
 setup (
    name = 'QPalma',
    description = 'A package for optimal alignment of short reads',
-   version = '1.0', 
+   version = '0.9', 
    long_description = '''
    A
    ''', 
index 41e30ba..e06c5d2 100644 (file)
@@ -335,9 +335,9 @@ def check_example(chromo,strand,b,e,seqInfo,lt1):
    dna,acc,don = seqInfo.get_seq_and_scores(1,strand,b,e)
 
    if lt1 != None:
-      _dna,_acc,_don= lt1.get_seq_and_scores(1,strand,b,e)
+      _dna,_acc,_don= lt1.get_seq_and_scores(chromo,strand,b,e)
    else:
-      _dna,_acc,_don = seqInfo.get_seq_and_scores(1,strand,b,e)
+      _dna,_acc,_don = seqInfo.get_seq_and_scores(chromo,strand,b,e)
 
    print 'Current interval: (%d,%d), current strand: %s'%(b,e,strand)
    print 'Results for dna,acc,don: %s %s %s'%(str(dna==_dna),str(acc==_acc),str(don==_don))
@@ -362,7 +362,8 @@ def simple_check(settings,seqInfo,lt1):
    chromo = 1
    total_size = seqInfo.getFragmentSize(chromo)
 
-   for strand in ['+','-']:
+   #for strand in ['+','-']:
+   for strand in ['-']:
       for (b,e) in [(206,874),(545,874),(999,1234),(1000,total_size-1000),(3,total_size-3),(0,total_size)]:
          check_example(chromo,strand,b,e,seqInfo,lt1)
 
@@ -392,7 +393,7 @@ def checks():
    lt1 = LookupTable(settings)
 
    print 'Checking with toy data...'
-   #simple_check(settings,seqInfo,lt1)
+   simple_check(settings,seqInfo,lt1)
 
    settings = {}
 
@@ -410,9 +411,9 @@ def checks():
    seqInfo = SeqSpliceInfo(accessWrapper,settings['allowed_fragments'])
 
    lt1 = None
-   #lt1 = LookupTable(settings)
+   lt1 = LookupTable(settings)
 
-   #print 'Checking with real data...'
+   print 'Checking with real data...'
    simple_check(settings,seqInfo,lt1)