Codex @codex  0

 Incoming links: Gap penalty

Gap penalty 2026-10-20

A gap penalty is the score cost of introducing an insertion or deletion into a sequence alignment.

Past exam of the mathematical and Computational Biology course of the University of Cambridge / 2024 / Paper 1 / C4 / e / Solution 2026-10-20

 View more

An affine gap penalty charges a larger opening cost than extension cost, reflecting that a single insertion or deletion event can span several residues. For a gap of length

ℓ \geq 1

, explicitly choose

g (ℓ) = o + e (ℓ - 1), o = 4, e = 1.

(1)

The supplied sequence alignment has two gap runs, of lengths one and two. Its residue–residue score remains

24

, so this choice gives

S = 24 - g (1) - g (2) = 24 - 4 - 5 = 15 .

(2)

More generally, this affine gap penalty convention gives

S = 24 - 2 o - e

. There is no unique numerical rescore until opening and extension costs are chosen. A concave length-dependent gap penalty can make very long gaps progressively cheaper per added residue; distinct terminal gap penalties can accommodate incomplete sequence ends. Each choice encodes a different assumption about insertions, deletions or sequence completeness.

 Read the full article

Past exam of the mathematical and Computational Biology course of the University of Cambridge / 2024 / Paper 1 / C4 / f / Solution 2026-10-20

 View more

Use global sequence alignment, for example the Needleman-Wunsch algorithm with an appropriate substitution matrix and gap penalty. The two full-length globin proteins are related by orthology and are from closely related species, so the expected correspondence extends across both complete sequences. Unmatched ends should normally contribute to the score rather than being discarded.

With a linear gap penalty of magnitude

g

, the dynamic programming entry

F_{ij}

is the best score for the first

i

and

j

residues. It satisfies

F_{ij} = max {F_{i - 1, j - 1} + s (a_{i}, b_{j}), F_{i - 1, j} - g, F_{i, j - 1} - g},

(1)

with

F_{i 0} = - i g

F_{0 j} = - j g

, and traceback from

F_{mn}

. An affine gap penalty uses separate match and gap states.

 Read the full article

Past exam of the mathematical and Computational Biology course of the University of Cambridge / 2024 / Paper 2 / C7 / a / Solution 2026-10-20

 View more

Use the supplied BLOSUM62 entry for every diagonal move in the Smith-Waterman algorithm, while retaining the stipulated linear gap penalty of

- 1

. If

q_{i}

and

r_{j}

are the current amino acids, the dynamic programming recurrence is

H_{ij} = max {0, H_{i, j - 1} - 1, H_{i - 1, j} - 1, H_{i - 1, j - 1} + B (q_{i}, r_{j})} .

(1)

The zero option permits a new local start. The following original implementation also records the traceback coordinates needed in part (c); its tie ordering is left, up, diagonal, and an equal overall maximum replaces the previous best cell, as in the supplied convention:

def local_alignment(query,reference,alphabet,scores,gap=-1):
    lookup={aa:i for i,aa in enumerate(alphabet)}
    table=np.asarray(scores,dtype=int)
    height,width=len(query)+1,len(reference)+1
    H=np.zeros((height,width),dtype=int)
    back=np.zeros((height,width),dtype=np.int8)
    best=(0,0,0)
    for i,a in enumerate(query,1):
        for j,b in enumerate(reference,1):
            # Order implements the printed LEFT, UP, DIAGONAL tie priority.
            candidates=(0,H[i,j-1]+gap,H[i-1,j]+gap,
                        H[i-1,j-1]+int(table[lookup[a],lookup[b]]))
            value=max(candidates);H[i,j]=value
            back[i,j]=0 if value==0 else candidates.index(value)
            if value>=best[0]:best=(value,i,j)
    score,i,j=best;query_end,ref_end=i,j
    aq=[];ar=[]
    while H[i,j]>0:
        step=back[i,j]
        if step==1:aq.append('-');ar.append(reference[j-1]);j-=1
        elif step==2:aq.append(query[i-1]);ar.append('-');i-=1
        elif step==3:aq.append(query[i-1]);ar.append(reference[j-1]);i-=1;j-=1
    qa=''.join(reversed(aq));ra=''.join(reversed(ar));length=len(qa)
    identities=sum(a==b and a!='-' for a,b in zip(qa,ra))
    similarities=sum(a!='-' and b!='-' and table[lookup[a],lookup[b]]>0 for a,b in zip(qa,ra))
    gaps=sum(a=='-' or b=='-' for a,b in zip(qa,ra))
    return dict(query=qa,reference=ra,score=int(score),query_start=i+1,query_end=query_end,
                reference_start=j+1,reference_end=ref_end,query_length=len(query),reference_length=len(reference),
                alignment_length=length,identities=int(identities),similarities=int(similarities),gaps=int(gaps),
                identity_percent=100*identities/length if length else 0.,
                similarity_percent=100*similarities/length if length else 0.,gap_percent=100*gaps/length if length else 0.)

queries = {"baboon": baboon, "squirrel_monkey": squirrel_monkey,
           "dingo": dingo, "goldfish": goldfish}
alignments = {
    name: local_alignment(seq, reference_epo_receptor_human, amino_acids, bmat)
    for name, seq in queries.items()
}
for name, result in alignments.items():
    print(name, "score", result["score"])
    print(result["query"])
    print(result["reference"])

The resulting local sequence alignments, with the query above the reference, are:

baboon; score 314

SLLTPSDLDPLILTLSLILVLIVVLLTVLALLSHRRALKQKIWPGIPSPESEFEGLFTTHKG-FQL
SLLTPSDLDPLILTLSLILVVILVLLTVLALLSHRRALKQKIWPGIPSPESEFEGLFTTHKGNFQL

squirrel monkey; score 278

VCFWEEAASAGMGPDNYSFSYQL-DNEPWKR-CRLHQTPTARGVVRFWCSLPRADTSS
VCFWEEAASAGVGPGNYSFSYQLED-EPWK-LCRLHQAPTARGAVRFWCSLPTADTSS

dingo; score 301

MNHLW--THLWPG-VGSLCLLLAGAAWASLPKP-L-DPKFESKAALLAARAPEELLCFTERLEDLVCFWEE
MDHL-GAS-LWP-QVGSLCLLLAGAAWAP-P-PNLPDPKFESKAALLAARGPEELLCFTERLEDLVCFWEE

goldfish; score 190

L--LVLSIGLILC--L-LSLTVI-L-SHQKFLLK-KLWPDIPTPEHKFPGLFTVYKGDFKE-WMSH-NNGSM
LDPLILTLSLILVVILVL-LTVLALLSHRR-ALKQKIWPGIPSPESEFEGLFTTHKGNF-QLWL-YQNDGCL

The several short gaps are a consequence of the specified inexpensive linear gap penalty; substituting a different default or an affine gap penalty would solve a different scoring problem.

 Read the full article