Bioinformatics & Computational Biology: Chapter 3 — Pairwise alignment

3.1 Three words that are not interchangeable #

3.1.1 Identity #

Identity is the fraction of aligned positions with the same letter. The denominator is part of the claim. Count gapped columns and every gap lowers the percentage; ignore them and gaps vanish from the number; divide by the shorter sequence and an unaligned tail is hidden. The same columns can be reported as 75% and as 100%.

Three matches in four columns is not three hundred in four hundred, though both can be written “75%”. Quote the counts and the length.

3.1.2 Similarity #

Similarity is not a standard percentage. It usually means identical letters, or conservative substitutes under a named matrix or chemical group. Sometimes the word is used for the raw score. Report instead the score under a named matrix and gap penalty, the identity count with its denominator, and the fraction of each sequence covered.

“Sixty percent similar” is none of those. It is not a probability and not a claim about structure. A biased stretch can agree while the rest of the two proteins do not.

3.1.3 Homology #

Homology means common ancestry. It is not a percentage, so “forty percent homologous” is a misuse. The defensible sentence is that the sequences are 40% identical over a stated alignment, more than a stated chance model expects, and homology is therefore inferred. The inference can be wrong. The percentage still only describes the alignment.

A shared domain does not make two multi-domain proteins homologous end to end. Analogy is similarity without common ancestry. Sequences rarely prove analogy; failing to support homology is a different result. High identity is evidence, not the definition, and low identity does not prove that two proteins are unrelated. Orthology and paralogy are Chapter 6.

3.2 What an alignment is #

The study of What an alignment is provides fundamental principles governing cellular architecture, molecular interactions, and physiological regulation.

🔬 Key Illustrated Scientific Figure

Full-Screen Zoom

Click to enlarge on this page

3.2.1 Matches, mismatches, and gaps #

An alignment puts in one column the letters proposed to come from one ancestral letter. A match is the same letter. A mismatch is a different letter, read as a substitution if the homology is accepted. A gap is a letter against nothing: an insertion in one sequence or a deletion in the other. The alignment does not say which, so “indel” is the honest word.

A gap is part of the hypothesis. Move it and the identity, any motif, a tree, and a homology model all change. Below, a match scores +1, a mismatch −1, and each gapped residue −1. Protein searches replace +1/−1 with a log-odds matrix and usually replace the linear gap with an affine penalty. The recurrence is the same idea. The optimal alignment is not.

3.2.2 Global, local, and the cases in between #

Global alignment (Needleman–Wunsch) scores every residue of both sequences, end gaps included. Use it when the sequences are homologous along roughly their whole length.

Local alignment (Smith–Waterman) scores the best segment pair and ignores the rest. Use it for a shared domain, a short motif in a long protein, or a read against a chromosome. A global alignment of a 50-residue query to a 2,000-residue protein spends its score on gaps. Database search is a local problem.

Semi-global alignment charges internal gaps and not end gaps. That fits a read inside a reference, or a shotgun overlap. Only the boundary conditions change.

3.3 Dynamic programming #

3.3.1 Filling a matrix from smaller problems #

Let F(i, j) be the best score of the first i letters of x with the first j letters of y. That alignment ends by pairing the last letters, or by gapping the last letter of x, or by gapping the last letter of y:

F(i, j) = max{ F(i−1, j−1) + s(xᵢ, yⱼ), F(i−1, j) + gap, F(i, j−1) + gap }.

Here s is +1 or −1 and the gap penalty is −1. A prefix against an empty string is only gaps: F(i, 0) = i × gap, F(0, j) = j × gap, and F(0, 0) = 0. Each cell uses only cells above or to the left, so one row-by-row pass is enough, and F(n, m) is the best global score. A pointer records which choice won. From the corner back to the origin, a diagonal writes both letters, up writes a letter of x against a gap, and left writes a gap against a letter of y. A tie means both alignments are optimal. Time and memory scale with n × m, which is why a database search does not fill this matrix against every subject.

score(a, b) = +1 if a equals b, otherwise −1
gap = −1
F(0, 0) = 0
for i in 1..n: F(i, 0) = i * gap
for j in 1..m: F(0, j) = j * gap
for i in 1..n:
    for j in 1..m:
        F(i, j) = max(
            F(i−1, j−1) + score(x_i, y_j),
            F(i−1, j) + gap,
            F(i, j−1) + gap)
traceback from F(n, m) to F(0, 0)

Smith–Waterman takes the maximum of those three candidates and zero, keeps the borders at zero, starts traceback at the highest cell anywhere, and stops at a zero.

3.3.2 A worked global alignment #

Take x = `ATGC` and y = `AGC`, match +1, mismatch −1, gap −1.

🔬 Key Illustrated Scientific Figure

Full-Screen Zoom

Click to enlarge on this page

Check the G of x against the G of y. Diagonal is F(T, A) + 1 = 0 + 1 = 1. Up is 0 − 1 = −1. Left is −1 − 1 = −2. The cell is 1, from the diagonal. Check A of x against G of y: diagonal −2, up −3, left 1 − 1 = 0, so the cell is 0 from a gap.

The corner is +2. Traceback is uniquely diagonal, diagonal, up, diagonal. In reverse the alignment is

A T G C
A - G C

Score +1 − 1 + 1 + 1 = +2. Three identical pairs sit in four columns: identity is 3/4 if gapped columns count and 3/3 if they do not. Say which. The path asserts an indel at T. It is the only optimal path under this scoring.

3.3.3 The same strings, aligned locally #

Borders are zero, and no cell is allowed to stay negative.

The best cell is +2. Traceback steps C with C, then G with G, then meets a zero and stops.

G C
G C

The A/A cell, score +1, is a separate local alignment, cut off by the zeros on the T row. A program that reports only the best local hit will omit it. Under +1/−1/−1 the two optima share a score and are not the same alignment: the global path gains +1 for A/A and pays −1 for the gap. With a gap penalty of −2, that global alignment scores +1, worse than local `GC` at +2. Local alignment may drop the awkward residue. Global alignment may not.

3.4 Affine gaps #

3.4.1 Opening and extension #

A linear penalty charges the same for every gapped residue, so one insertion of length ten costs as much as ten insertions of length one. A loop, an intron, or a slipped-strand run is usually one event. That event is what should be expensive.

An affine penalty, in the convention used here, costs open + extend × k for a gap of length k. A large opening and a small extension make one long gap cheaper than many short ones. Classical protein BLAST pairs BLOSUM62 with open 11 and extend 1, so the first gapped residue costs 12 and each further residue costs 1. Those are search settings, not chemical constants, and they change the λ and K of Chapter 4. E-values at different gap costs cannot be compared.

Gotoh’s method keeps three matrices, for paths ending on a substitution, a gap in x, or a gap in y. Entering a gap pays the opening; staying pays the extension. An intron does not literally cost “11 plus its length”, but the split stops one indel being scored as many single-residue gaps. A profile HMM makes the same split: a transition into insert or delete, versus a self-loop.

3.5 Key terms, summary, and self-check #

Key terms. Identity; similarity; homology; analogy; match; mismatch; indel; global alignment; local alignment; semi-global alignment; dynamic programming; Needleman–Wunsch; Smith–Waterman; traceback; affine gap.

Summary. Identity is a fraction on an alignment and needs a stated denominator. Similarity should be replaced by a score, an identity count, and a coverage. Homology is inferred common ancestry and is not a percentage. Global alignment scores both strings in full; local alignment scores the best segment pair. Dynamic programming fills each cell from three smaller cells. On `ATGC` against `AGC` at +1/−1/−1, the unique global alignment is `ATGC` over `A-GC` at score +2, and the optimal local alignment is `GC` over `GC`. Affine gaps charge an opening once so that one long indel costs less than many short ones.

Self-check.

1. How can `ATGC` over `A-GC` be reported both as 75% identity and as 100% identity?

2. Why is “the domains are 30% homologous” not a usable sentence?

3. In the global matrix, which predecessor supplies the +1 for G against G, and what are the other two candidates?

4. Why does the local traceback omit the initial A, even though A matched A?

5. Which gap model makes twelve separate one-residue gaps cheaper than one gap of length twelve, and which model reverses that?