nwrec.cpp

unix,linux下编译。用于蛋白质

/***
Needleman-Wunch recursions

Notation: i,j are prefix lengths so are in
ranges i = [0,|A|] and j = [0,|B|].

Profile positions are in ranges [0,|A|-1]
and [0,|B|-1] so prefix length i corresponds
to position (i-1) in the profile, and similarly
for j.

Terminal gap scoring
--------------------
Terminal gaps are scored as with open [close]
penalties only at the left [right] terminal,
as follows:

      0  i
	  |  |
	A XXXXX...
	B ---XX...

          i |A|-1
          |  |
	A ...XXXXX
	B ...XX---

In these examples, open / close penalty at position
i is  included, but close / open penalty at |A|-1 /
0 is not included.

This is implemented by setting the open [close] 
penalty to zero in the first [last] position of
each profile.

Consider adding a column to a sub-alignment. After the
column is added, there are i letters from A and j letters
from B.

The column starts a left-terminal gap if:
	Delete with i=1, j=0 or
	Insert with i=0, j=1.

The column ends a left-terminal gap if:
	Match following Delete with j=1, or
	Match following Insert with i=1.

The column starts a right-terminal gap if:
	Delete following a Match and i=|A|, or
	Insert following a Match and j=|B|.

The column ends a right-terminal gap if:
	Match with i=|A|, j=|B| following Delete or Insert.
	
RECURSION RELATIONS
===================

         i-1
          |
DD	A ..X X
	B ..- -

MD	A ..X X
	B ..X -

D(i,j) = max
			D(i-1,j) + e
			M(i-1,j) + goA(i-1)
Valid for:
	i = [1,|A|-1]
	j = [1,|B|]

I(i,j) By symmetry with D(i,j).

       i-2
        | i-1
		| |
MM	A ..X X
	B ..X X

DM	A ..X X
	B ..- X

IM  A ..- X
	B ..X X
	    | |
		| j-1
	   j-2

M(i,j) = L(i-1,j-1) + max
			M(i-1,j-1)
			D(i-1,j-1) + gcA(i-2)
			I(i-1,j-1) + gcB(j-2)
Valid for:
	i = [2,|A|]
	j = [2,|B|]

Equivalently:

M(i+1,j+1) = L(i,j) + max
			M(i,j)
			D(i,j) + gcA(i-1)
			I(i,j) + gcB(j-1)

Valid for:
	i = [1,|A|-1]
	j = [1,|B|-1]

Boundary conditions
===================

A XXXX
B ----
	D(0,0) = -infinity

	D(i,0) = ie
		i = [1,|A|]

	D(0,j) = -infinity
		j = [0,|B|]

I(0,0), I(0,j) and I(i,0) by symmetry with D.

	M(0,0) = 0
	M(i,0) = -infinity, i > 0
	M(0,j) = -infinity, j > 0

A X
B -
	D(1,0) = e
	D(1,j) = -infinity, j = [1,|B|]
		(assuming no I-D allowed).

	D(0,1) = -infinity
	D(1,1) = -infinity
	D(i,1) = max.
***/