Algorithme d'alignement global.

On se donne deux chaines de caractères à aligner le mieux possible. On autorise l'insertion de gaps. Par exemple pour aligner ACTTGTGA et CTGTGAAA on pourrait utiliser:
ACTTGTG--A
-CT-GTGAAA
Si on se donne un score de +2 pour une correspondance exacte et -1 sinon, on obtiendrait ici un score total de:
-1+2+2-1+2+2+2-1-1+2

Notations: si $S$ est une chaines de caractères, on note $\vert S\vert$ la longueur de $S$ et $S[i]$ le $i$-ème caractère de la chaine $S$.

Hypothèse: on suppose que le score de l'alignement de deux gaps "-" est négatif. Ceci assure qu'on n'insère un nombre fini de "-" pour obtenir le score maximum, le probleme de maximisation est alors une recherche de maximum dans un nombre fini de possibilités et admet donc une solution (pas forcément unique).

Si $S$ et $T$ deux chaines de longueur $n$ et $m$, on note $V(n,m)$ le score maximum d'alignement possible pour $S$ et $T$. On va calculer $V(n,m)$ par un procédé de récurrence à deux indices qui ressemble dans son concept au calcul des coefficients du triangle de Pascal.

On va donc calculer tous les $V(i,j)$ pour $i\leq n$ et $j \leq m$. Supposons qu'on souhaite calculer $V(i,j)$ avec $i$ et $j$ non nuls. On regarde quelle peut être la dernière paire de caractères alignés. Vu l'hypothèse faire sur la fonction score, il y a 3 possibilités:

• on aligne $S[i]$ avec $T[j]$, si c'est le cas:
  
  $$V(i,j)=V(i-1,j-1)+\sigma (S[i],T[j])$$
  
  

• on aligne $S[i]$ avec "-", si c'est le cas:
  
  $$V(i,j)=V(i-1,j)+\sigma (S[i],-)$$
  
  

• on aligne "-" avec $T[j]$, si c'est le cas:
  
  $$V(i,j)=V(i,j-1)+\sigma (-,S[j])$$
  
  

On en déduit donc que pour $i>0$ et $j>0$:

$$V(i,j)=\mbox{max} ( V(i-1,j-1)+\sigma (S[i],T[j]), V(i-1,j)+\sigma (S[i],-),V(i,j-1)+\sigma (-,T[j]) )$$ (1)

Il suffit ensuite d'initialiser $V(i,0)$ et $V(0,j)$ pour $i\leq n$ et $j \leq m$ ce qui est relativement simple puisque $V(i,0)$ ou $V(j,0)$ est l'alignement de $i$ lettres avec "-":

$$V(i,0)=\sum _{k=1}^i \sigma (S[k],-), \quad V(0,j)=\sum _{k=1}^j \sigma (-,T[k])$$

Par exemple si on se donne S=ACTG et T=AACC, on obtient (avec la fonction score valant 2 et -1):

	j	0	1	2	3	4
i			A	A	C	C
0		0	-1	-2	-3	-4
1	A	-1	2	1	0	-1
2	C	-2	1	1	3	2
3	T	-3	0	0	2	2
4	G	-4	-1	-1	1	1

On remplit la base ($i=0$ et $j=0$) puis par exemple ligne par ligne en appliquant la formule (1).

L'alignement optimal a donc pour score 1.

En C++, on utilisera un tableau d'entiers à deux indices, qu'on peut créer de la manière suivante:
vector<int> ligne(m+1,0);
vector< vector<int> > V(n+1,ligne);
On initialise ensuite V[0][j] et V[j][0] (faire une boucle à chaque fois). On exécute alors une double boucle (l'une imbriquée dans l'autre), par exemple on fait varier i de 1 à $n$ puis à i fixé, j de 1 à $m$ et on applique la formule de récurrence (1)

Amélioration: on remarque que la formule (1) ne fait intervenir que la ligne précédente ($i-1$) et la ligne actuelle ($i$), on peut donc se contenter de stocker seulement les 2 lignes Vprecedent[m+1] et Vactuel[m+1]

On peut ensuite chercher à retracer les alignements de score 1. Pour cela on part de $V(n,m)$ et on retrace dans (1) quel(s) terme(s) ont donné lieu au max. Cela revient à regarder dans les cases au-dessus, à gauche, ou en diagonale à gauche au-dessus quel(s) est (sont) la (les) case(s) telles que la différence entre les 2 scores correspond bien au score de l'alignement des 2 lettres (ou gap) qu'on rajoute. Par exemple ici, on part de $(4,4)$, on peut remonter en $(3,3)$ ou $(3,4)$ puis on se retouve en $(2,3)$ dans les deux cas. Ensuite on ne peut passer qu'en $(1,2)$ même si $(2,2)$ a la même valeur que $(1,2)$ car le score de $(2,3)$ est plus grand que les précédents il y a donc eu alignement de 2 lettres. Depuis $(1,2)$ on a à nouveau 2 possibilités $(1,1)$ et $(0,1)$ et on termine en $(0,0)$ ce qui donne 4 possibilités d'alignement de score 1:
A-CTG -ACTG A-CTG -ACTG
AACC- AACC- AAC-C AAC-C

Pour programmer la recherche des chemins, on note dans un tableau au fur et à mesure du calcul de V[i][j] les directions qui ont donné lieu au max de la formule. On peut pour cela utiliser un type structuré, par exemple:

struct element {
  int valeur;
  bool en_haut;
  bool a_gauche;
  bool diag;
};

vector<element> ligne(m+1);
vector< vector<element> > V(n+1,ligne);

Si le score V[i][j] est obtenu en utilisant la diagonale, on écrira par exemple:
V[i][j].diag=true;
On crée ensuite une fonction récursive calcule_chemin dont le rôle est d'ajouter à une liste d'alignements les alignements dont on connait la fin, plus précisément qui se terminent à la position (pos_a,pos_b) par les chaines de caractères fin_a et fin_b:

void calcule_chemin(const string & a,const string & b, 
		    const vector< vector<element> > & v,
		    string fin_a, string fin_b, int pos_a, int pos_b,
		    vector<string> & alignements){
  ...
}

Dans calcule_chemin, on teste d'abord si pos_a ou pos_b est nul. Par exemple si pos_a est nul, l'alignement est obtenu en rajoutant pos_b carcatères - à fin_a et en rajoutant le début (les pos_b premiers caractères) de b à fin_b.

Si pos_a et pos_b sont non nuls, on teste les trois possibilités (par en haut, par la gauche, par la diagonale) et on appelle récursivement la fonction d'affichage d'alignement par exemple pour un déplacement vertical:

  if (v[pos_a][pos_b].en_haut)
    calcule_chemin(a,b,v,a[pos_a-1]+fin_a,'-'+fin_b,pos_a-1,pos_b,alignements);