A dispetto di quanto detto nella sezione precedente esiste un modo
molto elegante che permette non solo di calcolare il costo di un allineamento
con complessità nello spazio lineare, ma anche di ricostruire un allineamento
con tale costo in tempo quadratico ().
Per poter ottenere questo risultato dobbiamo appliccare un approccio
basato sul dividi et impera (in inglese divide and conquer).
Prima di procedere dobbiamo notare che sebbene noi abbiamo utilizzato
la convenzione di cercare un allineamento ottimo con il massimo punteggio
a partire dall'inizio delle due sequenze fino al termine, è comunque
equivalente allineare le due sequenze a partire dal fondo. In pratica,
se consideriamo, che le nostre due sequenze partano da 0 e siano lunghe
rispettivamente e
(gli indici da
e
), allora
le relazioni di ricorrenza precedenti per il calcolo di un allineamento
con punteggio ottimale divengono
per l'inizializzazione:
Come esempio dell'equivalenza degli algoritmo, vediamo in
Figura , lo stesso costo e il medesimo all'allineamento
per le due sequenze 'ALS; e 'AS'.
Traducendo il tutto in codice Python, in Figura ,
vediamo l'analogo di linear_cost che però calcola il costo
di un allineamento tra due sequenze, occupando uno spazio lineare,
a partire dal termine e arrivando all'inizio (come effettuato manualmente
in Figura
).
![]() |
![]() |
Adesso che siamo convinti che sia equivalente allineare dall'inizio o dalla fine le due sequenze, possiamo passare a vedere come utilizzare un approccio dividi et impera.
Se noi consideriamo la linea che taglia
a metà la matrice dello score rispetto alla prima sequenza, abbiamo
che esiste un punto, in cui un allineamento con punteggio ottimo
attraversa quella riga. Ricordandoci che siamo in grado di calcolare
il costo di un allineamento sia partendo dall'inizio delle sequenze, sia
dalla fine, possiamo valutare qual'è il punto in cui il nostro
allineamento ottimo attraversa la frontiera, semplicemente
calcolando la posizione per cui il punteggio è massimo tra i
vettori che ci ritornano calcolando gli allineamenti parziali.
Più precisamante, se abbiamo che le due sequenze sono s1=s1[0:m]
e s1=s2[0:n] e la matrice dei punteggi è C
possiamo verificare che la colonna (jmax) in cui il nostro allineamento
attraversa la riga i-esima (C[i][jmax]),
la possiamo calcolare come
>>> s1,s2='ABCDX','ACLD' >>> B # backtrace matrix calulated as usual [['0', 'R', 'R', 'R', 'R'], ['C', 'D', 'R', 'R', 'R'], ['C', 'C', 'D', 'D', 'D'], ['C', 'C', 'D', 'R', 'R'], ['C', 'C', 'C', 'D', 'D'], ['C', 'C', 'C', 'D', 'C']] >>> # best path from the bottom-right corner: C D (R D) C D ... >>> def find_jmax(sTop,sBot): ... jmax=0 ... smax=sTop[0]+sBot[0] ... for j in range(1,len(sTop)): ... if sTop[j]+sBot[j]>= smax : ... smax=sTop[j]+sBot[j] ... jmax=j ... return jmax ... >>> len1=len(s1) >>> for i in range(len1,0,-1): # from the last to the first row ... (sTop,b)=linear_cost(s1[0:i],s2) # first half ... (sBot,br)=back_linear_cost(s1[i:],s2) # second half ... print i,s1[0:i],s1[i:],sTop,sBot ... print "backtrace",B[i][find_jmax(sTop,sBot)] ... 5 ABCDX [-5.0, -3.0, -1.0, -1.0, 0.0] [-4.0, -3.0, -2.0, -1.0, 0.0] backtrace C 4 ABCD X [-4.0, -2.0, 0.0, 0.0, 1.0] [-4.0, -3.0, -2.0, -1.0, -1.0] backtrace D 3 ABC DX [-3.0, -1.0, 1.0, 0.0, -1.0] [-3.0, -2.0, -1.0, 0.0, -2.0] backtrace R 2 AB CDX [-2.0, 0.0, 0.0, -1.0, -2.0] [-1.0, 0.0, -1.0, -1.0, -3.0] backtrace C 1 A BCDX [-1.0, 1.0, 0.0, -1.0, -2.0] [-1.0, -1.0, -2.0, -2.0, -4.0] backtrace D
L'aproccio è quindi quello di partire dal calcolo della colonna
() in cui l'allineamento attraversa la riga centrale (mid) della matrice,
dopo di che il problema è scomposto in due sottoproblemi più piccoli
che riguardano la valutazione del passaggio dell'allineamento ottimale
nel quarto superiore sinistro e nel quarto inferiore destro.
Questa procedura viene iterata dimezzando ogni volta l'area della
matrice che dobbiamo calcolare. La divisione in sottoproblemi
piccoli riduce drasticamente la complessità computazionale
nel tempo passando da una complessità
a
, come
mostrato in Figura
.
La tecnica del divide and conquer, può quindi essere schematizzata come
L'agoritmo che implementa questo approccio basato sul divide and conquer
è riportato in Figura .
È particolarmente elegante il fatto che pur utilizzando soltanto
dei vettori lunghi quanto la seconda sequenza si possa ottenere il costo
ed anche un allineamento ottimo con un tempo che rimane ancora
.
Per meglio capire il funzionamento di quest'approccio ricorsivo
in Figura riportiamo le varie fasi dell'algoritmo
utilizzando le due sequenze 'ABCDX' e 'ACLD'. Ciascuna chiamata è etichettata
da 0,T,B in dipendenza del fatto che siamo al livello iniziale,
si considera la parte sopra della sequenza (top), o quella finale (bottom).
Inotre la funzione ricorsiva, a seconda del caso che ci si trovi nel caso
base o in quello generale, ritornerà immediatamente il risultato richiesto,
oppure genererà a sua volta due chiamate ricorsive.
Per poter meglio visualizzare le chiamate annidate le riportiamo
indentate (Figura
), per cui per esempio
0.T.T.B identifica il fatto che a partire dall'inizio ('0')
abbiamo diviso il problema in sottoproblemi di dimensione
e che
la chiamata attuale (0.T.T.B) è stata ottenuta prendendo la:
metà sotto a destra ('B') della metà sopra a sinistra ('T')
della metà sopra a sinistra ('T') della matrice di score.
![]()
|
![]() |