Probalign is a sequence alignment tool that calculates a maximum expected accuracy alignment using partition function posterior probabilities.[1] Base pair probabilities are estimated using an estimate similar to Boltzmann distribution. The partition function is calculated using a dynamic programming approach.
The following describes the algorithm used by probalign to determine the base pair probabilities.[2]
To score an alignment of two sequences two things are needed:
- a similarity function (e.g. PAM, BLOSUM,...)
- affine gap penalty:
The score of an alignment a is defined as:
Now the boltzmann weighted score of an alignment a is:
Where is a scaling factor.
The probability of an alignment assuming boltzmann distribution is given by
Where is the partition function, i.e. the sum of the boltzmann weights of all alignments.
Dynamic programming
[edit]
Let denote the partition function of the prefixes and . Three different cases are considered:
- the partition function of all alignments of the two prefixes that end in a match.
- the partition function of all alignments of the two prefixes that end in an insertion .
- the partition function of all alignments of the two prefixes that end in a deletion .
Then we have:
The matrixes are initialized as follows:
The partition function for the alignments of two sequences and is given by , which can be recursively computed:
- analogously
Base pair probability
[edit]
Finally the probability that positions and form a base pair is given by:
are the respective values for the recalculated with inversed base pair strings.