Package Bio :: Module pairwise2
[hide private]
[frames] | no frames]

Module pairwise2

source code

Pairwise sequence alignment using a dynamic programming algorithm.

This provides functions to get global and local alignments between two sequences. A global alignment finds the best concordance between all characters in two sequences. A local alignment finds just the subsequences that align the best. Local alignments must have a positive score to be reported and they will not be extended for 'zero counting' matches. This means a local alignment will always start and end with a positive counting match.

When doing alignments, you can specify the match score and gap penalties. The match score indicates the compatibility between an alignment of two characters in the sequences. Highly compatible characters should be given positive scores, and incompatible ones should be given negative scores or 0. The gap penalties should be negative.

The names of the alignment functions in this module follow the convention <alignment type>XX where <alignment type> is either "global" or "local" and XX is a 2 character code indicating the parameters it takes. The first character indicates the parameters for matches (and mismatches), and the second indicates the parameters for gap penalties.

The match parameters are:

CODE  DESCRIPTION
x     No parameters. Identical characters have score of 1, otherwise 0.
m     A match score is the score of identical chars, otherwise mismatch
      score.
d     A dictionary returns the score of any pair of characters.
c     A callback function returns scores.

The gap penalty parameters are:

CODE  DESCRIPTION
x     No gap penalties.
s     Same open and extend gap penalties for both sequences.
d     The sequences have different open and extend gap penalties.
c     A callback function returns the gap penalties.

All the different alignment functions are contained in an object align. For example:

>>> from Bio import pairwise2
>>> alignments = pairwise2.align.globalxx("ACCGT", "ACG")

will return a list of the alignments between the two strings. For a nice printout, use the format_alignment method of the module:

>>> from Bio.pairwise2 import format_alignment
>>> print(format_alignment(*alignments[0]))
ACCGT
| ||
A-CG-
  Score=3
<BLANKLINE>

All alignment functions have the following arguments:

The other parameters of the alignment function depend on the function called. Some examples:

To see a description of the parameters for a function, please look at the docstring for the function via the help function, e.g. type help(pairwise2.align.localds) at the Python prompt.

Classes [hide private]
  identity_match
Create a match function for use in an alignment.
  dictionary_match
Create a match function for use in an alignment.
  affine_penalty
Create a gap function for use in an alignment.
Functions [hide private]
 
_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)
Return optimal alignments between two sequences (PRIVATE).
source code
 
_make_score_matrix_generic(sequenceA, sequenceB, match_fn, gap_A_fn, gap_B_fn, penalize_end_gaps, align_globally, score_only)
Generate a score and traceback matrix (PRIVATE).
source code
 
_recover_alignments(sequenceA, sequenceB, starts, best_score, score_matrix, trace_matrix, align_globally, gap_char, one_alignment_only, gap_A_fn, gap_B_fn, reverse=False)
Do the backtracing and return a list of alignments (PRIVATE).
source code
 
_find_start(score_matrix, best_score, align_globally)
Return a list of starting points (score, (row, col)) (PRIVATE).
source code
 
_reverse_matrices(score_matrix, trace_matrix)
Reverse score and trace matrices (PRIVATE).
source code
 
_clean_alignments(alignments)
Take a list of alignments and return a cleaned version (PRIVATE).
source code
 
_finish_backtrace(sequenceA, sequenceB, ali_seqA, ali_seqB, row, col, gap_char)
Add remaining sequences and fill with gaps if necessary (PRIVATE).
source code
 
_find_gap_open(sequenceA, sequenceB, ali_seqA, ali_seqB, end, row, col, col_gap, gap_char, score_matrix, trace_matrix, in_process, gap_fn, target, index, direction, best_score, align_globally)
Find the starting point(s) of the extended gap (PRIVATE).
source code
 
calc_affine_penalty(length, open, extend, penalize_extend_when_opening)
Calculate a penality score for the gap function.
source code
 
print_matrix(matrix)
Print out a matrix for debugging purposes.
source code
 
format_alignment(align1, align2, score, begin, end)
Format the alignment prettily into a string.
source code
 
_python_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)
Generate a score and traceback matrix according to Gotoh (PRIVATE).
source code
 
_python_rint(x, precision=1000)
Print number with declared precision.
source code
Variables [hide private]
  MAX_ALIGNMENTS = 1000
  align = align()
  _PRECISION = 1000
  __package__ = 'Bio'
Function Details [hide private]

_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)

source code 

Return optimal alignments between two sequences (PRIVATE).

This method either returns a list of optimal alignments (with the same score) or just the optimal score.

_make_score_matrix_generic(sequenceA, sequenceB, match_fn, gap_A_fn, gap_B_fn, penalize_end_gaps, align_globally, score_only)

source code 

Generate a score and traceback matrix (PRIVATE).

This implementation according to Needleman-Wunsch allows the usage of general gap functions and is rather slow. It is automatically called if you define your own gap functions. You can force the usage of this method with force_generic=True.

_recover_alignments(sequenceA, sequenceB, starts, best_score, score_matrix, trace_matrix, align_globally, gap_char, one_alignment_only, gap_A_fn, gap_B_fn, reverse=False)

source code 

Do the backtracing and return a list of alignments (PRIVATE).

Recover the alignments by following the traceback matrix. This is a recursive procedure, but it's implemented here iteratively with a stack.

sequenceA and sequenceB may be sequences, including strings, lists, or list-like objects. In order to preserve the type of the object, we need to use slices on the sequences instead of indexes. For example, sequenceA[row] may return a type that's not compatible with sequenceA, e.g. if sequenceA is a list and sequenceA[row] is a string. Thus, avoid using indexes and use slices, e.g. sequenceA[row:row+1]. Assume that client-defined sequence classes preserve these semantics.

_find_start(score_matrix, best_score, align_globally)

source code 

Return a list of starting points (score, (row, col)) (PRIVATE).

Indicating every possible place to start the tracebacks.

_clean_alignments(alignments)

source code 

Take a list of alignments and return a cleaned version (PRIVATE).

Remove duplicates, make sure begin and end are set correctly, remove empty alignments.

format_alignment(align1, align2, score, begin, end)

source code 

Format the alignment prettily into a string.

IMPORTANT: Gap symbol must be "-" (or ['-'] for lists)!

Since Biopython 1.71 identical matches are shown with a pipe character, mismatches as a dot, and gaps as a space.

Prior releases just used the pipe character to indicate the aligned region (matches, mismatches and gaps).

Also, in local alignments, if the alignment does not include the whole sequences, now only the aligned part is shown, together with the start positions of the aligned subsequences. The start positions are 1-based; so start position n is the n-th base/amino acid in the un-aligned sequence.

NOTE: This is different to the alignment's begin/end values, which give the Python indices (0-based) of the bases/amino acids in the aligned sequences.

_python_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)

source code 

Generate a score and traceback matrix according to Gotoh (PRIVATE).

This is an implementation of the Needleman-Wunsch dynamic programming algorithm as modified by Gotoh, implementing affine gap penalties. In short, we have three matrices, holding scores for alignments ending in (1) a match/mismatch, (2) a gap in sequence A, and (3) a gap in sequence B, respectively. However, we can combine them in one matrix, which holds the best scores, and store only those values from the other matrices that are actually used for the next step of calculation. The traceback matrix holds the positions for backtracing the alignment.