Back to index

python-biopython  1.60
Classes | Functions | Variables
Bio.pairwise2 Namespace Reference

Classes

class  align
class  identity_match
class  dictionary_match
class  affine_penalty

Functions

def _align
def _make_score_matrix_generic
def _make_score_matrix_fast
def _recover_alignments
def _find_start
def _find_global_start
def _find_local_start
def _clean_alignments
def _pad_until_equal
def _lpad_until_equal
def _pad
def _lpad
def rint
def calc_affine_penalty
def print_matrix
def format_alignment
def _test

Variables

int MAX_ALIGNMENTS = 1000
tuple align = align()
int _PRECISION = 1000

Function Documentation

def Bio.pairwise2._align (   sequenceA,
  sequenceB,
  match_fn,
  gap_A_fn,
  gap_B_fn,
  penalize_extend_when_opening,
  penalize_end_gaps,
  align_globally,
  gap_char,
  force_generic,
  score_only,
  one_alignment_only 
) [private]

Definition at line 311 of file pairwise2.py.

00311 
00312            one_alignment_only):
00313     if not sequenceA or not sequenceB:
00314         return []
00315 
00316     if (not force_generic) and isinstance(gap_A_fn, affine_penalty) \
00317     and isinstance(gap_B_fn, affine_penalty):
00318         open_A, extend_A = gap_A_fn.open, gap_A_fn.extend
00319         open_B, extend_B = gap_B_fn.open, gap_B_fn.extend
00320         x = _make_score_matrix_fast(
00321             sequenceA, sequenceB, match_fn, open_A, extend_A, open_B, extend_B,
00322             penalize_extend_when_opening, penalize_end_gaps, align_globally,
00323             score_only)
00324     else:
00325         x = _make_score_matrix_generic(
00326             sequenceA, sequenceB, match_fn, gap_A_fn, gap_B_fn,
00327             penalize_extend_when_opening, penalize_end_gaps, align_globally,
00328             score_only)
00329     score_matrix, trace_matrix = x
00330 
00331     #print "SCORE"; print_matrix(score_matrix)
00332     #print "TRACEBACK"; print_matrix(trace_matrix)
00333          
00334     # Look for the proper starting point.  Get a list of all possible
00335     # starting points.
00336     starts = _find_start(
00337         score_matrix, sequenceA, sequenceB,
00338         gap_A_fn, gap_B_fn, penalize_end_gaps, align_globally)
00339     # Find the highest score.
00340     best_score = max([x[0] for x in starts])
00341 
00342     # If they only want the score, then return it.
00343     if score_only:
00344         return best_score
00345     
00346     tolerance = 0  # XXX do anything with this?
00347     # Now find all the positions within some tolerance of the best
00348     # score.
00349     i = 0
00350     while i < len(starts):
00351         score, pos = starts[i]
00352         if rint(abs(score-best_score)) > rint(tolerance):
00353             del starts[i]
00354         else:
00355             i += 1
00356     
00357     # Recover the alignments and return them.
00358     x = _recover_alignments(
00359         sequenceA, sequenceB, starts, score_matrix, trace_matrix,
00360         align_globally, penalize_end_gaps, gap_char, one_alignment_only)
00361     return x

Here is the call graph for this function:

Here is the caller graph for this function:

def Bio.pairwise2._clean_alignments (   alignments) [private]

Definition at line 702 of file pairwise2.py.

00702 
00703 def _clean_alignments(alignments):
00704     # Take a list of alignments and return a cleaned version.  Remove
00705     # duplicates, make sure begin and end are set correctly, remove
00706     # empty alignments.
00707     unique_alignments = []
00708     for align in alignments :
00709         if align not in unique_alignments :
00710             unique_alignments.append(align)
00711     i = 0
00712     while i < len(unique_alignments):
00713         seqA, seqB, score, begin, end = unique_alignments[i]
00714         # Make sure end is set reasonably.
00715         if end is None:   # global alignment
00716             end = len(seqA)
00717         elif end < 0:
00718             end = end + len(seqA)
00719         # If there's no alignment here, get rid of it.
00720         if begin >= end:
00721             del unique_alignments[i]
00722             continue
00723         unique_alignments[i] = seqA, seqB, score, begin, end
00724         i += 1
00725     return unique_alignments

Here is the caller graph for this function:

def Bio.pairwise2._find_global_start (   sequenceA,
  sequenceB,
  score_matrix,
  gap_A_fn,
  gap_B_fn,
  penalize_end_gaps 
) [private]

Definition at line 672 of file pairwise2.py.

00672 
00673                        score_matrix, gap_A_fn, gap_B_fn, penalize_end_gaps):
00674     # The whole sequence should be aligned, so return the positions at
00675     # the end of either one of the sequences.
00676     nrows, ncols = len(score_matrix), len(score_matrix[0])
00677     positions = []
00678     # Search all rows in the last column.
00679     for row in range(nrows):
00680         # Find the score, penalizing end gaps if necessary.
00681         score = score_matrix[row][ncols-1]
00682         if penalize_end_gaps:
00683             score += gap_B_fn(ncols, nrows-row-1)
00684         positions.append((score, (row, ncols-1)))
00685     # Search all columns in the last row.
00686     for col in range(ncols-1):
00687         score = score_matrix[nrows-1][col]
00688         if penalize_end_gaps:
00689             score += gap_A_fn(nrows, ncols-col-1)
00690         positions.append((score, (nrows-1, col)))
00691     return positions

Here is the caller graph for this function:

def Bio.pairwise2._find_local_start (   score_matrix) [private]

Definition at line 692 of file pairwise2.py.

00692 
00693 def _find_local_start(score_matrix):
00694     # Return every position in the matrix.
00695     positions = []
00696     nrows, ncols = len(score_matrix), len(score_matrix[0])
00697     for row in range(nrows):
00698         for col in range(ncols):
00699             score = score_matrix[row][col]
00700             positions.append((score, (row, col)))
00701     return positions

Here is the caller graph for this function:

def Bio.pairwise2._find_start (   score_matrix,
  sequenceA,
  sequenceB,
  gap_A_fn,
  gap_B_fn,
  penalize_end_gaps,
  align_globally 
) [private]

Definition at line 657 of file pairwise2.py.

00657 
00658                 penalize_end_gaps, align_globally):
00659     # Return a list of (score, (row, col)) indicating every possible
00660     # place to start the tracebacks.
00661     if align_globally:
00662         if penalize_end_gaps:
00663             starts = _find_global_start(
00664                 sequenceA, sequenceB, score_matrix, gap_A_fn, gap_B_fn, 1)
00665         else:
00666             starts = _find_global_start(
00667                 sequenceA, sequenceB, score_matrix, None, None, 0)
00668     else:
00669         starts = _find_local_start(score_matrix)
00670     return starts

Here is the call graph for this function:

Here is the caller graph for this function:

def Bio.pairwise2._lpad (   s,
  char,
  n 
) [private]

Definition at line 749 of file pairwise2.py.

00749 
00750 def _lpad(s, char, n):
00751     # Prepend n chars to the beginning of s.
00752     return char*n + s

Here is the caller graph for this function:

def Bio.pairwise2._lpad_until_equal (   s1,
  s2,
  char 
) [private]

Definition at line 735 of file pairwise2.py.

00735 
00736 def _lpad_until_equal(s1, s2, char):
00737     # Add char to the beginning of s1 or s2 until they are equal
00738     # length.
00739     ls1, ls2 = len(s1), len(s2)
00740     if ls1 < ls2:
00741         s1 = _lpad(s1, char, ls2-ls1)
00742     elif ls2 < ls1:
00743         s2 = _lpad(s2, char, ls1-ls2)
00744     return s1, s2

Here is the call graph for this function:

Here is the caller graph for this function:

def Bio.pairwise2._make_score_matrix_fast (   sequenceA,
  sequenceB,
  match_fn,
  open_A,
  extend_A,
  open_B,
  extend_B,
  penalize_extend_when_opening,
  penalize_end_gaps,
  align_globally,
  score_only 
) [private]

Definition at line 446 of file pairwise2.py.

00446 
00447     align_globally, score_only):
00448     first_A_gap = calc_affine_penalty(1, open_A, extend_A,
00449                                       penalize_extend_when_opening)
00450     first_B_gap = calc_affine_penalty(1, open_B, extend_B,
00451                                       penalize_extend_when_opening)
00452 
00453     # Create the score and traceback matrices.  These should be in the
00454     # shape:
00455     # sequenceA (down) x sequenceB (across)
00456     lenA, lenB = len(sequenceA), len(sequenceB)
00457     score_matrix, trace_matrix = [], []
00458     for i in range(lenA):
00459         score_matrix.append([None] * lenB)
00460         trace_matrix.append([[None]] * lenB)
00461 
00462     # The top and left borders of the matrices are special cases
00463     # because there are no previously aligned characters.  To simplify
00464     # the main loop, handle these separately.
00465     for i in range(lenA):
00466         # Align the first residue in sequenceB to the ith residue in
00467         # sequence A.  This is like opening up i gaps at the beginning
00468         # of sequence B.
00469         score = match_fn(sequenceA[i], sequenceB[0])
00470         if penalize_end_gaps:
00471             score += calc_affine_penalty(
00472                 i, open_B, extend_B, penalize_extend_when_opening)
00473         score_matrix[i][0] = score
00474     for i in range(1, lenB):
00475         score = match_fn(sequenceA[0], sequenceB[i])
00476         if penalize_end_gaps:
00477             score += calc_affine_penalty(
00478                 i, open_A, extend_A, penalize_extend_when_opening)
00479         score_matrix[0][i] = score
00480 
00481     # In the generic algorithm, at each row and column in the score
00482     # matrix, we had to scan all previous rows and columns to see
00483     # whether opening a gap might yield a higher score.  Here, since
00484     # we know the penalties are affine, we can cache just the best
00485     # score in the previous rows and columns.  Instead of scanning
00486     # through all the previous rows and cols, we can just look at the
00487     # cache for the best one.  Whenever the row or col increments, the
00488     # best cached score just decreases by extending the gap longer.
00489 
00490     # The best score and indexes for each row (goes down all columns).
00491     # I don't need to store the last row because it's the end of the
00492     # sequence.
00493     row_cache_score, row_cache_index = [None]*(lenA-1), [None]*(lenA-1)
00494     # The best score and indexes for each column (goes across rows).
00495     col_cache_score, col_cache_index = [None]*(lenB-1), [None]*(lenB-1)
00496 
00497     for i in range(lenA-1):
00498         # Initialize each row to be the alignment of sequenceA[i] to
00499         # sequenceB[0], plus opening a gap in sequenceA.
00500         row_cache_score[i] = score_matrix[i][0] + first_A_gap
00501         row_cache_index[i] = [(i, 0)]
00502     for i in range(lenB-1):
00503         col_cache_score[i] = score_matrix[0][i] + first_B_gap
00504         col_cache_index[i] = [(0, i)]
00505         
00506     # Fill in the score_matrix.
00507     for row in range(1, lenA):
00508         for col in range(1, lenB):
00509             # Calculate the score that would occur by extending the
00510             # alignment without gaps.
00511             nogap_score = score_matrix[row-1][col-1]
00512             
00513             # Check the score that would occur if there were a gap in
00514             # sequence A.
00515             if col > 1:
00516                 row_score = row_cache_score[row-1]
00517             else:
00518                 row_score = nogap_score - 1   # Make sure it's not the best.
00519             # Check the score that would occur if there were a gap in
00520             # sequence B.  
00521             if row > 1:
00522                 col_score = col_cache_score[col-1]
00523             else:
00524                 col_score = nogap_score - 1
00525 
00526             best_score = max(nogap_score, row_score, col_score)
00527             best_score_rint = rint(best_score)
00528             best_index = []
00529             if best_score_rint == rint(nogap_score):
00530                 best_index.append((row-1, col-1))
00531             if best_score_rint == rint(row_score):
00532                 best_index.extend(row_cache_index[row-1])
00533             if best_score_rint == rint(col_score):
00534                 best_index.extend(col_cache_index[col-1])
00535 
00536             # Set the score and traceback matrices.
00537             score = best_score + match_fn(sequenceA[row], sequenceB[col])
00538             if not align_globally and score < 0:
00539                 score_matrix[row][col] = 0
00540             else:
00541                 score_matrix[row][col] = score
00542             trace_matrix[row][col] = best_index
00543 
00544             # Update the cached column scores.  The best score for
00545             # this can come from either extending the gap in the
00546             # previous cached score, or opening a new gap from the
00547             # most previously seen character.  Compare the two scores
00548             # and keep the best one.
00549             open_score = score_matrix[row-1][col-1] + first_B_gap
00550             extend_score = col_cache_score[col-1] + extend_B
00551             open_score_rint, extend_score_rint = \
00552                              rint(open_score), rint(extend_score)
00553             if open_score_rint > extend_score_rint:
00554                 col_cache_score[col-1] = open_score
00555                 col_cache_index[col-1] = [(row-1, col-1)]
00556             elif extend_score_rint > open_score_rint:
00557                 col_cache_score[col-1] = extend_score
00558             else:
00559                 col_cache_score[col-1] = open_score
00560                 if (row-1, col-1) not in col_cache_index[col-1]:
00561                     col_cache_index[col-1] = col_cache_index[col-1] + \
00562                                              [(row-1, col-1)]
00563 
00564             # Update the cached row scores.
00565             open_score = score_matrix[row-1][col-1] + first_A_gap
00566             extend_score = row_cache_score[row-1] + extend_A
00567             open_score_rint, extend_score_rint = \
00568                              rint(open_score), rint(extend_score)
00569             if open_score_rint > extend_score_rint:
00570                 row_cache_score[row-1] = open_score
00571                 row_cache_index[row-1] = [(row-1, col-1)]
00572             elif extend_score_rint > open_score_rint:
00573                 row_cache_score[row-1] = extend_score
00574             else:
00575                 row_cache_score[row-1] = open_score
00576                 if (row-1, col-1) not in row_cache_index[row-1]:
00577                     row_cache_index[row-1] = row_cache_index[row-1] + \
00578                                              [(row-1, col-1)]
00579                     
00580     return score_matrix, trace_matrix
    

Here is the call graph for this function:

Here is the caller graph for this function:

def Bio.pairwise2._make_score_matrix_generic (   sequenceA,
  sequenceB,
  match_fn,
  gap_A_fn,
  gap_B_fn,
  penalize_extend_when_opening,
  penalize_end_gaps,
  align_globally,
  score_only 
) [private]

Definition at line 365 of file pairwise2.py.

00365 
00366     score_only):
00367     # This is an implementation of the Needleman-Wunsch dynamic
00368     # programming algorithm for aligning sequences.
00369     
00370     # Create the score and traceback matrices.  These should be in the
00371     # shape:
00372     # sequenceA (down) x sequenceB (across)
00373     lenA, lenB = len(sequenceA), len(sequenceB)
00374     score_matrix, trace_matrix = [], []
00375     for i in range(lenA):
00376         score_matrix.append([None] * lenB)
00377         trace_matrix.append([[None]] * lenB)
00378 
00379     # The top and left borders of the matrices are special cases
00380     # because there are no previously aligned characters.  To simplify
00381     # the main loop, handle these separately.
00382     for i in range(lenA):
00383         # Align the first residue in sequenceB to the ith residue in
00384         # sequence A.  This is like opening up i gaps at the beginning
00385         # of sequence B.
00386         score = match_fn(sequenceA[i], sequenceB[0])
00387         if penalize_end_gaps:
00388             score += gap_B_fn(0, i)
00389         score_matrix[i][0] = score
00390     for i in range(1, lenB):
00391         score = match_fn(sequenceA[0], sequenceB[i])
00392         if penalize_end_gaps:
00393             score += gap_A_fn(0, i)
00394         score_matrix[0][i] = score
00395 
00396     # Fill in the score matrix.  Each position in the matrix
00397     # represents an alignment between a character from sequenceA to
00398     # one in sequence B.  As I iterate through the matrix, find the
00399     # alignment by choose the best of:
00400     #    1) extending a previous alignment without gaps
00401     #    2) adding a gap in sequenceA
00402     #    3) adding a gap in sequenceB
00403     for row in range(1, lenA):
00404         for col in range(1, lenB):
00405             # First, calculate the score that would occur by extending
00406             # the alignment without gaps.
00407             best_score = score_matrix[row-1][col-1]
00408             best_score_rint = rint(best_score)
00409             best_indexes = [(row-1, col-1)]
00410 
00411             # Try to find a better score by opening gaps in sequenceA.
00412             # Do this by checking alignments from each column in the
00413             # previous row.  Each column represents a different
00414             # character to align from, and thus a different length
00415             # gap.
00416             for i in range(0, col-1):
00417                 score = score_matrix[row-1][i] + gap_A_fn(i, col-1-i)
00418                 score_rint = rint(score)
00419                 if score_rint == best_score_rint:
00420                     best_score, best_score_rint = score, score_rint
00421                     best_indexes.append((row-1, i))
00422                 elif score_rint > best_score_rint:
00423                     best_score, best_score_rint = score, score_rint
00424                     best_indexes = [(row-1, i)]
00425             
00426             # Try to find a better score by opening gaps in sequenceB.
00427             for i in range(0, row-1):
00428                 score = score_matrix[i][col-1] + gap_B_fn(i, row-1-i)
00429                 score_rint = rint(score)
00430                 if score_rint == best_score_rint:
00431                     best_score, best_score_rint = score, score_rint
00432                     best_indexes.append((i, col-1))
00433                 elif score_rint > best_score_rint:
00434                     best_score, best_score_rint = score, score_rint
00435                     best_indexes = [(i, col-1)]
00436 
00437             score_matrix[row][col] = best_score + \
00438                                      match_fn(sequenceA[row], sequenceB[col])
00439             if not align_globally and score_matrix[row][col] < 0:
00440                 score_matrix[row][col] = 0
00441             trace_matrix[row][col] = best_indexes
00442     return score_matrix, trace_matrix

Here is the call graph for this function:

Here is the caller graph for this function:

def Bio.pairwise2._pad (   s,
  char,
  n 
) [private]

Definition at line 745 of file pairwise2.py.

00745 
00746 def _pad(s, char, n):
00747     # Append n chars to the end of s.
00748     return s + char*n

Here is the caller graph for this function:

def Bio.pairwise2._pad_until_equal (   s1,
  s2,
  char 
) [private]

Definition at line 726 of file pairwise2.py.

00726 
00727 def _pad_until_equal(s1, s2, char):
00728     # Add char to the end of s1 or s2 until they are equal length.
00729     ls1, ls2 = len(s1), len(s2)
00730     if ls1 < ls2:
00731         s1 = _pad(s1, char, ls2-ls1)
00732     elif ls2 < ls1:
00733         s2 = _pad(s2, char, ls1-ls2)
00734     return s1, s2

Here is the call graph for this function:

def Bio.pairwise2._recover_alignments (   sequenceA,
  sequenceB,
  starts,
  score_matrix,
  trace_matrix,
  align_globally,
  penalize_end_gaps,
  gap_char,
  one_alignment_only 
) [private]

Definition at line 583 of file pairwise2.py.

00583 
00584                         penalize_end_gaps, gap_char, one_alignment_only):
00585     # Recover the alignments by following the traceback matrix.  This
00586     # is a recursive procedure, but it's implemented here iteratively
00587     # with a stack.
00588     lenA, lenB = len(sequenceA), len(sequenceB)
00589     tracebacks = [] # list of (seq1, seq2, score, begin, end)
00590     in_process = [] # list of ([same as tracebacks], prev_pos, next_pos)
00591 
00592     # sequenceA and sequenceB may be sequences, including strings,
00593     # lists, or list-like objects.  In order to preserve the type of
00594     # the object, we need to use slices on the sequences instead of
00595     # indexes.  For example, sequenceA[row] may return a type that's
00596     # not compatible with sequenceA, e.g. if sequenceA is a list and
00597     # sequenceA[row] is a string.  Thus, avoid using indexes and use
00598     # slices, e.g. sequenceA[row:row+1].  Assume that client-defined
00599     # sequence classes preserve these semantics.
00600 
00601     # Initialize the in_process stack
00602     for score, (row, col) in starts:
00603         if align_globally:
00604             begin, end = None, None
00605         else:
00606             begin, end = None, -max(lenA-row, lenB-col)+1
00607             if not end:
00608                 end = None
00609         # Initialize the in_process list with empty sequences of the
00610         # same type as sequenceA.  To do this, take empty slices of
00611         # the sequences.
00612         in_process.append(
00613             (sequenceA[0:0], sequenceB[0:0], score, begin, end,
00614              (lenA, lenB), (row, col)))
00615         if one_alignment_only:
00616             break
00617     while in_process and len(tracebacks) < MAX_ALIGNMENTS:
00618         seqA, seqB, score, begin, end, prev_pos, next_pos = in_process.pop()
00619         prevA, prevB = prev_pos
00620         if next_pos is None:
00621             prevlen = len(seqA)
00622             # add the rest of the sequences
00623             seqA = sequenceA[:prevA] + seqA
00624             seqB = sequenceB[:prevB] + seqB
00625             # add the rest of the gaps
00626             seqA, seqB = _lpad_until_equal(seqA, seqB, gap_char)
00627             
00628             # Now make sure begin is set.
00629             if begin is None:
00630                 if align_globally:
00631                     begin = 0
00632                 else:
00633                     begin = len(seqA) - prevlen
00634             tracebacks.append((seqA, seqB, score, begin, end))
00635         else:
00636             nextA, nextB = next_pos
00637             nseqA, nseqB = prevA-nextA, prevB-nextB
00638             maxseq = max(nseqA, nseqB)
00639             ngapA, ngapB = maxseq-nseqA, maxseq-nseqB
00640             seqA = sequenceA[nextA:nextA+nseqA] + gap_char*ngapA + seqA
00641             seqB = sequenceB[nextB:nextB+nseqB] + gap_char*ngapB + seqB
00642             prev_pos = next_pos
00643             # local alignment stops early if score falls < 0
00644             if not align_globally and score_matrix[nextA][nextB] <= 0:
00645                 begin = max(prevA, prevB)
00646                 in_process.append(
00647                     (seqA, seqB, score, begin, end, prev_pos, None))
00648             else:
00649                 for next_pos in trace_matrix[nextA][nextB]:
00650                     in_process.append(
00651                         (seqA, seqB, score, begin, end, prev_pos, next_pos))
00652                     if one_alignment_only:
00653                         break
00654                     
00655     return _clean_alignments(tracebacks)

Here is the call graph for this function:

Here is the caller graph for this function:

def Bio.pairwise2._test ( ) [private]
Run the module's doctests (PRIVATE).

Definition at line 856 of file pairwise2.py.

00856 
00857 def _test():
00858     """Run the module's doctests (PRIVATE)."""
00859     print "Runing doctests..."
00860     import doctest
00861     doctest.testmod(optionflags=doctest.IGNORE_EXCEPTION_DETAIL)
00862     print "Done"

def Bio.pairwise2.calc_affine_penalty (   length,
  open,
  extend,
  penalize_extend_when_opening 
)

Definition at line 809 of file pairwise2.py.

00809 
00810 def calc_affine_penalty(length, open, extend, penalize_extend_when_opening):
00811     if length <= 0:
00812         return 0
00813     penalty = open + extend * length
00814     if not penalize_extend_when_opening:
00815         penalty -= extend
00816     return penalty

Here is the caller graph for this function:

def Bio.pairwise2.format_alignment (   align1,
  align2,
  score,
  begin,
  end 
)
format_alignment(align1, align2, score, begin, end) -> string

Format the alignment prettily into a string.

Definition at line 834 of file pairwise2.py.

00834 
00835 def format_alignment(align1, align2, score, begin, end):
00836     """format_alignment(align1, align2, score, begin, end) -> string
00837 
00838     Format the alignment prettily into a string.
00839 
00840     """
00841     s = []
00842     s.append("%s\n" % align1)
00843     s.append("%s%s\n" % (" "*begin, "|"*(end-begin)))
00844     s.append("%s\n" % align2)
00845     s.append("  Score=%g\n" % score)
00846     return ''.join(s)
00847 
00848 
00849 # Try and load C implementations of functions.  If I can't,
00850 # then just ignore and use the pure python implementations.
00851 try:
    from cpairwise2 import rint, _make_score_matrix_fast
def Bio.pairwise2.print_matrix (   matrix)
print_matrix(matrix)

Print out a matrix.  For debugging purposes.

Definition at line 817 of file pairwise2.py.

00817 
00818 def print_matrix(matrix):
00819     """print_matrix(matrix)
00820 
00821     Print out a matrix.  For debugging purposes.
00822 
00823     """
00824     # Transpose the matrix and get the length of the values in each column.
00825     matrixT = [[] for x in range(len(matrix[0]))]
00826     for i in range(len(matrix)):
00827         for j in range(len(matrix[i])):
00828             matrixT[j].append(len(str(matrix[i][j])))
00829     ndigits = map(max, matrixT)
00830     for i in range(len(matrix)):
00831         #Using string formatting trick to add leading spaces,
00832         print " ".join("%*s " % (ndigits[j], matrix[i][j]) \
00833                        for j in range(len(matrix[i])))

def Bio.pairwise2.rint (   x,
  precision = _PRECISION 
)

Definition at line 754 of file pairwise2.py.

00754 
00755 def rint(x, precision=_PRECISION):
00756     return int(x * precision + 0.5)

Here is the caller graph for this function:


Variable Documentation

Definition at line 753 of file pairwise2.py.

Definition at line 305 of file pairwise2.py.

Definition at line 152 of file pairwise2.py.