Algoritmi per la
concatenazione di
match per il mapping
del cDNA
Flavia Rustici
Miguel Gutierrez
Michele Orbini Michelucci
Sommario
1.
2.
3.
Introduzione
Terminologia e definizioni
Match Chaining Problem su 2 sequenze
1.
2.
3.
Algoritmo base
Range Maximum Query (RMQ)
Algoritmo MCCM
1.
2.
4.
Match Chaining Problem su k sequenze
1.
Terminologia e definizioni (2)
2.
Concatenzione di match non overlapping
1.
3.
RMQ a dimensioni multiple
Concatenzione di match overlapping
1.
2.
3.
5.
Limitazione della distanza tra i match
Overlap
Riduzione a grafo
Algoritmo basato su geometria
Risultati Sperimentali
Riferimenti
2
Introduzione

Genomica comparativa
 branca
della biologia che si occupa del
confronto di sequenze genomiche
Confronto tra genomi di organismi in stretta
relazione
 Confronto di regioni che conservano la sintenia

(regioni in cui geni ortologhi occorrono nello stesso
ordine)

Stabilire le regioni codificanti e non codificanti di
sequenze (struttura esoni-introni)
3
Introduzione (2)

cDNA



DNA complementare ottenuto dall’mRNA nella fase di
trascrizione inversa
Contiene solo esoni
cDNA mapping


Obiettivo: trovare il gene (e riconoscere la struttura esoniintroni) nel genoma da cui il cDNA è stato estratto
Permette ulteriori analisi degli elementi regolatori del gene
(promotori, fattori di trascrizione, binding sites, etc.)
4
Introduzione (3)

Approccio software (anchor-based):
1.
Calcolo dei frammenti (segmenti uguali in
sequenze diverse)
2.
Calcolo di una catena globale di score
massimo di frammenti (ancore)
3.
Allineamento delle regioni comprese tra le
ancore
5
Sommario
1.
Introduzione
2.
Terminologia e definizioni
3.
Match Chaining Problem su 2 sequenze
1.
2.
3.
Algoritmo base
Range Maximum Query (RMQ)
Algoritmo MCCM
1.
2.
4.
Match Chaining Problem su k sequenze
1.
Terminologia e definizioni (2)
2.
Concatenzione di match non overlapping
1.
3.
RMQ a dimensioni multiple
Concatenzione di match overlapping
1.
2.
3.
5.
Limitazione della distanza tra i match
Overlap
Riduzione a grafo
Algoritmo basato su geometria
Risultati Sperimentali
Riferimenti
6
Terminologia e definizioni

Siano P e T due sequenze tali che




P: pattern (cDNA)
T: testo (DNA)
Sia M=(m1,m2, …,mn)una lista di match
Definizioni

Un match (frammento) mi tra due sequenze P e T è una tripla
(pi,ti,li) tale che P[ pi  pi  li  1]  T [ti ti  li  1]
T (DNA) = … G G A A A T A G T T A …
mi
P (cDNA) = … A C C A T A G T C G C C …
Esempio match mi tra 2 sequenze
7
Terminologia e definizioni(2)

Due match mi e mj hanno overlap sul pattern quando
pi  p j  pi  li

Due match mi e mj hanno overlap sul testo quando
ti  t j  ti  li
T (DNA) = … G G A C A C A G T ….….. G T A C C C G C C
mi
pj
mj
P (cDNA) = … G G A C A C A G T A C C C G C C …
overlap sulppattern
pEsempio
i
i+li
8
Terminologia e definizioni(3)

Possono esistere casi in cui due match mi e mj hanno overlap sia su
pattern che su testo
text_overlap(mi, mj) = max { 0, ti + li - tj}=2
T (DNA) = … G G A C A C A G T G T G T A C C C G C C
mi
mj
P (cDNA) = … G G A C A C A G T G T A C C C G C C …
pattern_overlap(mi, mj) = max { 0, pi + li - pj}=4
overlap(mi, mj) = max {pattern_overlap(mi,mj), text_overlap(mi,mj)}=4
9
Sommario
2.
Introduzione
Terminologia e definizioni
3.
Match Chaining Problem su 2 sequenze
1.
1.
2.
3.
Algoritmo base
Range Maximum Query (RMQ)
Algoritmo MCCM
1.
2.
4.
Match Chaining Problem su k sequenze
1.
Terminologia e definizioni (2)
2.
Concatenzione di match non overlapping
1.
3.
RMQ a dimensioni multiple
Concatenzione di match overlapping
1.
2.
3.
5.
Limitazione della distanza tra i match
Overlap
Riduzione a grafo
Algoritmo basato su geometria
Risultati Sperimentali
Riferimenti
10
Match Chaining Problem
Definizione del problema:

Input

M = {m1, m2, … , mn} una lista di match ordinati in base
alla loro posizione sul pattern pi
i j

pi  p j
i  j  ti  t j

Output

Catena di match C  { mc1 , mc2 ,, mc f } :

CM

i,j: mci ,mc j  C, i  j  (pci  pc j )  (tci  tc j )

C deve massimizzare una funzione di score
11
Match Chaining Problem (2)

score(C)


misura della “bontà” di una catena di frammenti C
Esempi

score( C )  f , dove f è il numero di match in C
f

score( C )   lci , somma delle lunghezze dei match in C
i 1

Algoritmo di base per la concatenazione
1. For i from 1 to n
find prev(mi)
O( n  g( n ))
compute score(mi)
2. Select the match with the largest score() and the construct the actual chain by
tracing prev() entries.
12
Range Maximum Query (RMQ)
Problema la cui soluzione sarà utile per il calcolo del match
predecessore ( prev(m) )

Problema: Dato un array A di n interi trovare l’indice di
massimo valore nel subarray A[i … j] per ogni i e j.

Soluzioni:


Tempo O(n)
Array statico
Esempio:
i
j
A … 4 2 8 1 3 5 …
RMQ( A, i, j) = 8
13
RMQ Dinamico
Definizione del problema:

Input

insieme S di elementi, inizialmente vuoto

ogni elemento I S ha associati 2 valori



a(I)
b(I)
Output RMQ(S, p, q)

elemento I di S tale che
I  S’ & b(I)= max { b(I’) | I’  S’}
dove S’ = {I’ | I’  S & p ≤ a(I’) ≤ q }  S
14
RMQ Dinamico: implementazione

Necessaria una struttura dati che supporti le
seguenti operazioni:


Ricerca
Aggiornamento (inserzione/cancellazione)
AVL TREE
Caratteristiche generali:
 Albero binario di ricerca bilanciato in altezza
 Ricerca e Aggiornamento O(log2 n), con n = n° nodi
correnti
6
3
1
8
4
9
15
Nozioni sugli alberi

Un albero è definito da un insieme di nodi V collegati tra
loro da un insieme di archi E tale che:
 esiste un nodo speciale detto radice (root)
 tutti i nodi tranne la radice hanno un padre
 ogni nodo ha un insieme di nodi figli


ad eccezione di un sottoinsieme di nodi detti foglia (leaf) che
non hanno figli
la relazione padre-figlio è rappresentata da un arco
che collega il nodo padre al nodo figlio
16
Albero(V,E) rappresentazione grafica
a
b
e
f
c
g
d
i
h
m
n
l
o
p
q
radice
nodo interno
r
s
foglia

Ogni nodo è radice di un sottoalbero
17
Nozioni sugli alberi(2)

Un albero è detto binario quando ogni nodo ha al più due figli

Dato un nodo v indichiamo con
 sx(v): il figlio sinistro di v
 dx(v): il figlio destro di v
 sinistra(v): il sottoalbero radicato nel figlio sinistro di v
 destra(v): il sottoalbero radicato nel figlio destro di v
Esempio:
a
sx(a) = b
dx(a) = c
b
sinistra(a)
c
destra(a)
18
Nozioni sugli alberi(3)

Un albero binario è detto di ricerca quando
ogni suo nodo v rispetta le seguenti proprietà:
v'  sin istra( v ) v'  v
v' '  destra( v ) v  v' '
6
3
1
8
4
9
Albero binario di ricerca ordinato su valori interi
19
Nozioni sugli alberi(4)

Un albero binario è detto bilanciato in altezza quando,
in ogni nodo, le altezze del sottoalbero destro e sinistro
differiscono di al più 1
 l’altezza di un albero binario bilanciato in altezza è O(log2 n)
a
hdx = 1
c
b
hsx = 3
d
b
a
d
e
hsx - hdx = 2
f
Albero binario non bilanciato
f
e
c
Albero binario bilanciato
20
Nozioni sugli alberi(5)

Un albero AVL è un albero binario di ricerca
bilanciato in altezza
6
3
1
8
4
9
Albero AVL
L’albero AVL è la struttura dati che utilizzeremo per
risolvere il problema dell’RMQ dinamico

Ricerca e Aggiornamento O(log n), con n = n° nodi correnti
21
RMQ Dinamico: implementazione(2)

Ogni nodo v dell’AVL-tree rappresenta un oggetto IS,
ed ha associato:



a(v), chiave (viene usata per l’ordinamento dei nodi nell’AVL)
b(v), rappresenta lo score
MAX_POINTER, puntatore al nodo con il massimo valore di b
nel sottoalbero radicato in v
6
a(v)
3
6
3
1
4
8
1
3
b(v)
3
MAX_POINTER(v)
4
7
9
5
22
RMQ Dinamico: esempio

Aggiornamento: insert(11,8)

Fase 1: ricerca posizione
b(11) = 8
6
b(MAX_POINTER(6))
b(MAX_POINTER(8)) = 57
b(MAX_POINTER(9))
579 <non
8 aggiornare
ha figli: aggiungo il
nodo come figlio dx
MAX_POINTER(8)
MAX_POINTER(9)
MAX_POINTER(6)
3
1
4
11
3
8
8
1
4
7
11
3
8
9
11
5
8
11
8
23
RMQ Dinamico: esempio

Aggiornamento: insert(11,8)

Fase 2: rotazione
Dopo l’inserimento il nodo 8 è
sbilanciato, c’è bisogno di una
ROTAZIONE.
6
3
1
4
3
8
1
4
7
3
9
h( sx(8) ) - h( dx(8) ) = -2
5
11
8
24
RMQ Dinamico: esempio

Aggiornamento: insert(11,8)

Fase 2: rotazione
6
8 è diventato nodo foglia, si deve
aggiornare MAX_POINTER(8).
3
1
4
3
9
1
4
7
8
3
5
11
8
Rotazione e aggiornamento dei MAX_POINTER
aumentano la complessità dell’inserzione?
25
RMQ Dinamico: aggiornamento

Complessità di un aggiornamento con rotazione

v1
v2
v2
Rotazione
v1
T1
T3
T2
T3
T1
T2
La rotazione richiede tempo O(1)
in quanto i puntatori da aggiornare
sono 3.
 Anche
l’aggiornamento
del
MAX_POINTER richiede tempo
O(1), in quanto riguarda solo v1 e
v2: basta confrontare i loro
MAX_POINTER con quelli delle
radici dei sottoalberi T1, T2 e T3.
La rotazione e l’aggiornamento dei MAX_POINTER non aumentano
l’ordine di complessità:

O(log n), dovuto a fase di ricerca della posizione in cui eseguire
l’inserzione/cancellazione.
26
RMQ Dinamico: algoritmo

rmq(S,p,d)

Fase 1:




Ricerca del nodo vleft tale che a(vleft) = min { a(I) | p ≤ a(I) , IS }
Ricerca del nodo vright tale che a(vright) = max { a(I) | a(I) ≤ q , IS }
Ricerca del nodo vlca, lowest common ancestor di vleft e vright
Fase 2:

Controllando i MAX_POINTER dei nodi lungo il cammino da vleft a
vright troviamo il nodo v con il massimo b tale che p ≤ a(v) ≤ q
27
RMQ Dinamico: esempio
p q

Ricerca: rmq(S, 4, 9)

Fase 1:



ricerca vleft Pleft = [ 6 , 3 , 4 ]
ricerca vright Pright= [ 6 , 9 ]
ricerca vlca (lowest common ancestor)
6<
> 94 → dx
sx
p=4
93 =
< 94 ,→
trovato
dx
vright
4 = 4 , trovato vleft
6
3
1
4
q=9
3
9
1
4
7
8
3
5
11
8
28
RMQ Dinamico: esempio

Fase2:

checkSX visita ricorsivamente tutti i nodi
nel cammino P1 = [ sx(vlca), …,vleft] per
trovare il massimo b nel range [p, a(vlca))
 checkDX è speculare per cammino
P2 = [ dx(vlca), …,vright] nel range (a(vlca), q]
checkSX(v,P1,m)::
begin
if (v == null) then return m;
if (isLeaf(v)) then return maxb(v,m);
if (dx(v)P1) then
Vlca
return checkSX(dx(v),P1,m);
fuori range
else //dx(v)P1
MAX_POINTER(v1) =
MAX_POINTER(v3)
 checkSX(dx(v1))
V1
fuori range
m1 = maxb(MAX_POINTER(dx(v)),v,m);
V2
return checkSX(sx(v),P1,m1);
end
Vright
V3
fuori range
Vleft
V4
O(log n)
max_tmp
V5
max_tmp = maxb { v3, MAX_POINTER(dx(v3)), max_tmp}
29
RMQ Dinamico: complessità
rmq(S, p, q)::
begin
Pleft = search_vleft(S, p);
Pright = search_vright(S, q);
Plca = compute_vlca(Pleft, Pright);
P1 = Pleft – Plca;
O(log n)
P2 = Pright – Plca;
O(log n)
maxleft = checkSX(vlca, P1, null);
maxright = checkDX(vlca, P2, null);
return maxb(maxleft, maxright, vlca);
O(1)
end
30
Sommario
1.
2.
3.
Introduzione
Terminologia e definizioni
Match Chaining Problem su 2 sequenze
2.
Algoritmo base
Range Maximum Query (RMQ)
3.
Algoritmo MCCM
1.
1.
2.
4.
Match Chaining Problem su k sequenze
1.
Terminologia e definizioni (2)
2.
Concatenzione di match non overlapping
1.
3.
RMQ a dimensioni multiple
Concatenzione di match overlapping
1.
2.
3.
5.
Limitazione della distanza tra i match
Overlap
Riduzione a grafo
Algoritmo basato su geometria
Risultati Sperimentali
Riferimenti
31
MCCM (Shibuya e Kurochkin 2003)

Input: M = {m1, m2, … , mn} una lista di match ordinati in base alla
loro posizione sul pattern pi

Output: catena di match C C={m
…,,mm
} :
{m
mc2c ,,
c1c, ,m
c cf}




CM
i,j: mci ,mc j  C, i  j
1
2
f
 (pci  pc j )  (tci  tc j )
C deve massimizzare una funzione di score
Caratteristiche:



Si applica a due sequenze genomiche (un genoma ed una
sequenza di cDNA)
Basato su Dynamic Range Maximum Query (RMQ)
Permette:



Mapping efficiente e accurato
Limitazione della distanza tra frammenti
Overlap di frammenti
32
MCCM: distanza massima tra match

È necessario limitare la distanza nel testo tra due match consecutivi
per evitare situazioni come questa
MAX_LEN
DNA
prev(mi)
mi
cDNA

Chiamiamo MAX_LEN la distanza massima ammissibile tra due
match nella catena C.
33
MCCM con MAX_LEN
f
score (C )   lci
i 1
1.
2.
R = struttura dati RMQ, inizialmente vuota;
per ogni i da 1 a nn
a.
per ogni match mj tale che j<i, pj+lj≤pi e mmjjR
R
insert(R, mm
tj+l
-1, -1,
score(m
j,j,
j));
insert(R,
tjj+l
score(m
j
j));
b.
c.
3.
4.
prev(mi)= rmq(R,
rmq(R, tit
- iMAX_LEN
-1, t-1,
- MAX_LEN
i-1);ti-1);
se (prev(mi) = null) allora score(mi)=li;
altrimenti score(mi)= score(prev(mi))+li;
m = <trova il match con il massimo score>;
C = <catena costruita tracciando i prev a
partire da m>
tempo
O(n log(n))
Non permette overlap!
34
MCCM con overlap

È utile permettere gli overlap perchè:


sono molto comuni nel mapping del cDNA, occorrendo
spesso alle estremità degli esoni del cDNA
aumentano la quantità di sequenza ricoperta dalla
catena
T (DNA) = … G G A C A C A G T G T G T A C C C G C C
mi
mj
P (cDNA) = … G G A C A C A G T G T A C C C G C C …
35
MCCM con overlap: idea

L’algoritmo considera 4 tipi di match candidati per prev(mi):
A.
B.
C.
D.
match che non si sovrappongono con mi
match mj che hanno pattern_overlap(mj,mi)>text_overlap(mj,mi)≥0
match mj che hanno pattern_overlap(mj,mi)=0 e text_overlap(mj,mi) >0
match mj che hanno text_overlap(mj,mi)>pattern_overlap(mj,mi)>0
S2
(cDNA)
mi
b2
b1
d
a
c
S1 (genoma)
36
MCCM con overlap: idea (2)
1.
Per ogni match mi :
a)
b)
2.
per ciascun insieme A, B, C e D si cerca il miglior
match candidato ad essere prev(mi)
prev(mi) sarà il “migliore” tra i 4 candidati
Si determina il match m con il massimo score e
si costruisce la catena a partire da m utilizzando
l’informazione prev()
f
f 1
i 1
i 1
score( C )   lci   overlap _ length( mci , mci1 )
37
MCCM con overlap

Ricerca del match candidato predecessore per mi di tipo A:
struttura dati RMQ dinamica ( R1 )
in R1 sono inseriti i match mj tali che




j<ie
pattern_overlap(mj, mi) = 0
rmq(R1, ti  MAX_LEN  1, ti  1)

MAX_LEN
: mi
: mj (j < i)  R1
T
P
Non ci sono match in R1 che abbiano overlap su pattern con mi
38
MCCM con overlap

Ricerca del match candidato predecessore per mi di tipo B:
struttura dati RMQ dinamica ( R2 )
in R2 sono inseriti i match mj tali che




ji e
pattern_overlap(mj, mi)  0
rmq(R2, ti  pi  MAX_LEN  1, ti  pi  1)

: mi
: mj (j < i)  R2
T
P
In R2 ci sono solo match che hanno overlap su pattern con mi
Con l’RMQ cerchiamo il match che massimizza lo score tra quelli
che hanno
pattern_overlap(mj, mi)  text_overlap(mj, mi)
39
MCCM con overlap

Ricerca del match candidato predecessore per mi di tipo C:

struttura dati coda ( Q ) implementata con albero binario di ricerca, che
supporta inserimenti, cancellazioni e ricerche (non c’è bisogno di rmq)
in Q sono inseriti i match mj tali che




ji
e
pattern_overlap(mj, mi) = 0
find_best(Q,ti)
: mi
: mj (j < i)  Q
T
P
Non ci sono match in Q che abbiano overlap su pattern con mi.
In Q andiamo a cercare l’elemento che massimizza lo score di mi
tra quelli che hanno
text_overlap(mj, mi)  0
41
MCCM con overlap

Ricerca del match candidato predecessore per mi di tipo D:

struttura dati RMQ dinamica ( R3 )
in R3 sono inseriti i match mj tali che



ji
e
pattern_overlap(mj, mi)  0
rmq(R3, ti  pi  1, ti  1)

: mi
: mj (j < i)  R3
T
P
In R3 ci sono solo match che hanno overlap su pattern con mi
Con l’RMQ cerchiamo il match che massimizza lo score tra quelli
che hanno
text_overlap(mj, mi)  pattern_overlap(mj, mi)
NB: non si garantisce che tj  ti
42
MCCM con overlap

Come si calcola il migliore tra i 4 match
candidati?

si calcola lo score che mi avrebbe usando
come predecessore di mi ognuno dei 4
candidati
 si sceglie quello che massimizza score(mi)
score( mi )  score( m )  li  overlap( m , mi )
m  prev( mi )
44
MCCM con overlap
R1, R2, R3 = strutture dati RMQ, inizialmente vuote;
Q = coda, inizialmente vuota;
per ogni i da 1 a n in ordine
1.
2.
3.
per ogni match mj tale che j<i, pj+lj≤pi e mjR1
a)




m(a) = rmq(R1,ti - MAX_LEN - 1, ti - 1);
m(b) = rmq(R2,ti - pi - MAX_LEN - 1, ti - pi - 1);
m(c) = find_best(Q,ti);
m(d) = rmq(R3,ti - pi, ti- 1);
se ((m(d)  null)AND(t(m(d))>ti )) allora m(d)=nil;
se ((m(a) = null)AND(m(b) = null) AND
(m(c) = null) AND(m(d) = null)) allora
b)
c)
d)
e)
f)
g)


h)
i)
j)
k)
4.
5.
insert(R1, mj, tj+lj-1, score(mj));
delete(R2, mj)
insert_queue((Q, mj, tj, tj+lj-1, score(mj)-tj-lj );
delete(R3, mj)
tempo
O(n log(n))
prev(mi)= nil;
score(mi)=li;
altrimenti prev(mi)= max_score(mi, m(a), m(b), m(c), m(d));
//max_score: restituisce il match che massimizza la funzione di score
per mi
score(mi) = score(prev(mi) + li – overlap(prev(mi), mi);
insert(R2, mi, ti-pi, score(mi)- pi - li);
insert(R3, mi, ti-pi, key(score(mi)- ti - lj, ti - pi);
m = <match con il massimo score>;
C = <catena costruita tracciando i prev a
partire da m>
45
Sommario
1.
2.
3.
Introduzione
Terminologia e definizioni
Match Chaining Problem su 2 sequenze
1.
2.
3.
Algoritmo base
Range Maximum Query (RMQ)
Algoritmo MCCM
1.
2.
4.
Match Chaining Problem su k sequenze
1.
Terminologia e definizioni (2)
2.
Concatenzione di match non overlapping
1.
3.
RMQ a dimensioni multiple
Concatenzione di match overlapping
1.
2.
3.
5.
Limitazione della distanza tra i match
Overlap
Riduzione a grafo
Algoritmo basato su geometria
Risultati Sperimentali
Riferimenti
47
Mapping su sequenze multiple


I genomi di specie in stretta relazione oppure i
genomi di ceppi differenti di una stessa specie,
condividono una grande quantità di sequenza
lo stesso avviene per una sequenza di cDNA
con il genoma da cui proviene
Estendiamo MCCM per mappare una sequenza
di cDNA su più sequenze genomiche per:


identificare geni comuni a più genomi
individuare le regioni di sintenia
48
Terminologia e definizioni

Definizioni
 Un match (frammento) m su k sequenze (k-1
sequenze genomiche e una di cDNA) può essere
rappresentato come un iper-rettangolo in uno spazio
k-dimensionale Rk, ossia una k-tupla
([p1 … q1], …, [pk … qk])
tale che S1[p1 … q1] = …= Sk[pk … qk]
S2
S3
q1
p1
S1
p2
q2
p3
q3
49
Terminologia e definizioni

Possiamo definire un match m anche dai suoi angoli estremi,
ossia da una coppia di k-tuple (beg(m),end(m)) tale che:



beg(m) = (beg(m).x1, … , beg(m).xk) = (p1, … , pk)
end(m) = (end(m).x1, … , end(m).xk) = (q1, … , qk)
Relazione << tra frammenti:

m’<<m (m’ precede m) sse
i  1..k

 1 ) beg( m' ).xi  beg( m ).xi

2 ) end ( m' ).xi  end ( m ).xi
m’ e m si dicono colineari se m’<<m oppure m<<m’

ossia, appaiono nello stesso ordine in tutte le sequenze
50
Problema più generale


Obiettivo: trovare una catena C di match non-overlaping
frammenti tale che l'ammontare di sequenza ricoperta dai
frammenti della catena sia massimizzato.
Input

M = {m1, m2, … , mn}, un insieme di n match su k sequenze
Output


trovare un sottoinsieme CC={f
, fc2, ,…
 {1m
m,c2fl,}, mc f } di match colineari tale
1
che:

elementi di C non si sovrappongono

la quantità di sequenza ricoperta dai frammenti è massimizzato.
f
score( C )   mci .length
i 1
51
Higher-dimensional Chaining


Gli algoritmi di approssimazione risolvono il
problema della concatenazione dei match per più
di 2 genomi

Problema: complessità quadratica nel numero n di match

Può essere un grande svantaggio per grandi valori di n
Alternativa: Algoritmo di Zhang

costruisce una catena ottima usando la divisione dello
spazio basandosi su kd-alberi
52
Terminologia e definizioni

KD-Albero: albero K-Dimensionale


Un Albero kd (k-dimensionale) è una struttura dati per il
partizionamento dello spazio, che organizza i punti su uno
spazio euclídeo a k dimensioni
Utilizza soltanto piani perpendicolari a uno degli assi del
sistema di coordinate
53
K-Dimensional Tree
Esempio di suddivisione di spazio per la creazione di kd-tree, con k = 3





La prima divisione (rosso) taglia
la cella principale (bianco) in 2
sottocelle
le sottocelle sono divise a loro
volta (verde) in due sottocelle
Queste 4 sottocelle sono
ancora divise (blu) a loro volta
in 2 sottocelle ciscuna
Poiché non ho più divisioni da
fare, le 8 celle finali sono
chiamate celle foglie
Le sfere gialle rappresentano i
nodi dell’albero
54
Higher-dimensional Chaining


Gli algoritmi di approssimazione risolvono il problema
della concatenazione dei match per più di 2 genomi

Problema: compessità quadratica nel numero n di match

Può essere un grande svantaggio per grandi valori di n
Alternativa: Algoritmo di Zhang


costruisce una catena ottima usando la divisione dello
spazio basandosi su kd-alberi
Problema: l’analisi rigorosa della complessità in tempo di
questo algoritmo è difficile perché la costruzione della
catena è integrata nella struttura dell’kd-tree
55
Higher-dimensional Chaining (2)


Alternativa: Algoritmo di Myers e Miller

Basato sul paradigma line-sweep

Ricerca ortogonale su Range Tree

Tempo O(n logk n)
Adesso presenteremo un algoritmo che migliora la
complessità dell’algoritmo di Myers e Miller
56
Higher-dimensional Chaining (3)

Ricerca ortogonale su range tree

Costruzione di range tree
57
Costruzione di un range tree

Dato un insieme S di n punti d-dimensionali, si
costruisce il suo range tree in questo modo:


Se d = 1: il range tree di S è un albero binario di ricerca di
altezza minima o un array ordinato sugli alementi di S
Se d > 1:



il range tree di S è un albero T binario di ricerca di altezza
minima con n foglie, la cui i-esima foglia più a sinistra contiene
l’elemento in S con l’i-esimo valore pù piccolo della coordinata x1
associamo ad ogni nodo interno v di T un sottoinsieme canonico
Cv S che contiene punti memorizzati nelle foglie del sottoalbero
radicato in v
Per ogni v,





lv è la più piccola coordinata x1 di ogni punto di Cv
hv è la più grande coordinata x1 di ogni punto di Cv
C*v  {(p.x 2 ..., p.x d )  R d -1 : (P.x 1 , p.x 2 ..., p.x d )  C v }
Ogni nodo interno v contiene lv, hv e un range tree a d-1
dimensioni costruito su C*v
Per ogni dimensione fissata d la struttura dati è costruita in
58
O(n logd-1 n) tempo e spazio
Higher-dimensional Chaining (3)

Ricerca ortogonale su range tree

Costruzione di range tree

Range Query
59
Range Query (RQ)

Dato un insieme S di punti in uno spazio a d dimensioni
Rd, una Range Query restituisce tutti i punti di S che si
trovano nel iper-rettangolo
R (p, q) = [l1... H1] × [L2... H2] × ... × [ld ... hd]

Se d = 1: la query può essere risolta in tempo O(log n),
dovuto alla ricerca binaria sull’albero

Se d > 1: si visita il range tree a partire dalla radice
60
Visita del range tree

per ogni nodo visitato v:

se v è una foglia:


restituiamo il punto corrispondente se si trova nell’iperrettangolo R(p, q)
se v è un nodo interno:


[l v . .. h v ]  [l1... h1 ]  {} non si fa niente
1 [l v . .. h v ]  [l1... h1 ] si ricerca ricorsivamente nel range tree di
dimensione d-1 memorizzato in v con l’iper-rettagolo
[l2... h2] × ... × [ld ... hd]


altrimenti, visitiamo ricorsivamente entrambi i figli di v
Tempo O (logd n + z), con z il numero di punti nell’iperrettangolo R(p, q)
61
Range Query (RQ)

Dato un insieme S di punti in uno spazio a d dimensioni
Rd, una Range Query restituisce tutti i punti di S che si
trovano nel iper-rettangolo
R (p, q) = [l1... H1] × [L2... H2] × ... × [ld ... hd]

Se d = 1: la query può essere risolta in tempo O(log n),
dovuto alla ricerca binaria sull’albero

Se d > 1: la query può essere risolta in tempo
O (logd n + z), con z il numero di punti nell’iper-rettangolo
R(p, q) visitando ricorsivamente il range tree a partire
dalla radice
62
Higher-dimensional Chaining (3)

Ricerca ortogonale su range tree

Costruzione di range tree

Range Query

Possiamo migliorare l’efficienza se come ultima
dimensione dell’albero se utiliziamo un albero di priorità
anziché un range tree, ottenendo un algoritmo di
concatenazione di complessità in tempo O(n logk-1 n)

Possiamo migliorare ancora l’algoritmo mediante la
tecnica di fractional cascading, che risparmia un fattore
logaritmico nel rispondere alle Range Query
63
Fractional Cascading


Vogliamo modificare la tecnica per RQ(0,q) per
rispondere efficientemente al Range Maximum
Query RMQ(0,q)
Consideriamo il caso d = 2

il range tree è un x-albero binario di ricerca di y-array

per ogni nodo v del x-albero indichiamo con v.L e v.R
rispettivamente i suoi figli sinistro e destro
l’y-array di v, indicato con Av, contiene tutti i punti in Cv
ordinati in modo crescente in base alla loro coordinata y
ogni elemento p  Av contiene due puntatori






Lptr punta all’elemento q1 più grande in Av.L tale che q1  p
Rptr punta all’elemento q2 più grande in Av.R tale che q2  p
se questi elementi non esistono i puntatori hanno valore
NULL
Esempio di questa struttura
64
Esempio
Range tree con Fractional Cascading per i punti nella figura in alto a destra
65
Esempio (2)




I nodi colorati sono quelli visitati dalla Range Query
I nodi tratteggiati sono i nodi canonici
I cerchi indicano i puntatori NULL
Gli elementi cololati dell’y-array dei nodi canonici sono i punti
all’interno del rettangolo di query
66
Fractional Cascading (2)


La ricerca di tutti i punti nel rettangolo R(0,(h1,h2))
avviene in 2 passi:
 Primo passo: la ricerca binaria sul y-array del
nodo radice del x-albero per localizzare il punto
ph2 più a destra, tale che ph2.y  [0... h2]
 Secondo
passo: visitiamo ricorsivamente
x-albero tenendo traccia dei puntatori per
localizzare la foglia più a destra ph1 tale che
ph1.x  [0... h1]
Durante la visita ricorsiva del x-albero,
identifichiamo un insieme di nodi che chiameremo
canonici, rispetto alla nostra query
67
Terminologia e definizioni

Un insieme di nodi canonici è il più piccolo
insieme P di nodi v1, ... , vl del x-albero tale che

l
l
P :  j 1 Av j   j 1 C v j contiene tutti i punti p  S tali che
p.x  [0...h1]
68
Fractional Cascading (3)


Se nessun p  P soddisfa p.y  [0... h1] analiziamo
i puntantori durante la visita del x-albero
Se incontriamo un nodo canonico vj ,l’elemento ej
puntato dell’ultimo puntatore visitato suddivide la
lista AVj nel modo seguente:





gli elementi a destra di ej non sono contenuti nell’iperrettangolo R (0, (h1, h2))
tutti gli altri elementi di AVj sono contenuti nell’iperrettangolo R (0, (h1, h2))
ej è detto split element.
Il numero di nodi canonici è O(log n)
Possiamo trovare tutti i nodi canonici e gli split
element del loro y-array in tempo O(log n)
69
Fractional Cascading (4)

Quindi, sia z il numero di punti nell’iperrettangolo R(0, q), utilizzando la tecnica di
fractional cascading possiamo realizzare i
Range Query in tempo:


se d = 2 O(log n + z)
se d > 2 O(logd-1 n + z)
70
RMQ su Range Tree

Possiamo rispondere al RMQ migliorando gli yarray nel fractional cascading implementandoli
con
una
coda
di
priorità
costruita
sull’ordinamento dei punti contenuti nell’ y-array


il rank space dei punti nel y-array è l’intervallo [1 …
N ], dove N è la dimensione dell’array
la coda di priorità supporta:

inserzioni, cancellazione, ricerca del predecessore e
successore in tempo O(log log N)
71
RMQ su Range Tree (2)

Per rispondere a RMQ(0,q):
1.
2.
si cercano tutti i nodi canonici
per ogni nodo canonico vj
1.
si cerca lo split element in AVj e sia rj il suo indice
all’interno del y-array


2.
3.
tutti i punti in AVj a sinistra dello split elemento (esso stesso
compreso) sono contenuti in R(0,q)
i restanti elementi di AVj non appartengono a R(0,q)
qj = predecessor(rj) è l’indice dell’elemento di massimo
score in AVj
l’elemento con il massimo score nel range(0,q)
sarà l’elemento che ha il massimo score tra i
massimi trovati su ogni nodo canonico
72
RMQ su Range Tree: complessità

Rispondere a una RMQ a due dimensioni
(d=2) costa in tempo O(log n log(log n))
 Il numero di nodi canonici è O(log n)


Ogni operazione sulla coda di priorità costa
O(log log n)
Se d ≥ 2
costa in tempo O(logd-1 n log(log n))
73
Teorema


Dati k>2 sequenze genomiche e n match si può
trovare una catena ottima globale in :

tempo O(n logk-2 n log (log n))

spazio O(n logk-2 n)
Dimostrazione (cenno)



Per ogni match eseguiamo una RMQ
Poichè i punti sono già ordinati sulla prima dimensione
l’RMQ necessita solo di d = k – 1 dimensioni
Per d ≥ 2

la ricerca dell’RMQ costa spazio O(n logd-1 n) e tempo
O(n logd-1 n log log n)

Pertanto, trovare la catena ottima globale costa in tempo
O(n logk-2 n log (log n)) e spazio O (n logk-2 n)
74
Sommario
1.
2.
3.
Introduzione
Terminologia e definizioni
Match Chaining Problem su 2 sequenze
1.
2.
3.
Algoritmo base
Range Maximum Query (RMQ)
Algoritmo MCCM
1.
2.
4.
Match Chaining Problem su k sequenze
1.
Terminologia e definizioni (2)
2.
Concatenzione di match non overlapping
1.
3.
RMQ a dimensioni multiple
Concatenzione di match overlapping
1.
2.
3.
5.
Limitazione della distanza tra i match
Overlap
Riduzione a grafo
Algoritmo basato su geometria
Risultati Sperimentali
Riferimenti
75
Terminologia e definizioni

Ad ogni match m è associato un valore di rarità r definito come
il numero massimo di volte che la sottostringa Si[pi … pi+l-1]
occorre nella stessa sequenza Si, per 1  i  k
S1
r=4
S2
S3

Dati due match m e m’ tali che
i  1..k beg( m' ).xi  beg( m ).xi
diremo che hanno overlap sse
end (m' ). xk  end (m). xk
&
end (m' ). xk  beg (m). xk
S1
m’
S2
m
S3
Non c’è ovelap tra m’ e m
76
Terminologia e definizioni

Per ogni coppia di frammenti m’ e m su k sequenze, dove la kesima è la sequenza di cDNA, la quantità di overlap sul
cDNA è calcolata come
end ( m' ).xk  beg( m ).xk  1
overlapk ( m' , m ) 
0


altrimenti
Relazione << tra frammenti:

m’<<m (m’ precede m) sse
i  1..k

se beg( m ).xk  end ( m' ).xk  end ( m ).xk
 1 ) beg( m' ).xi  beg( m ).xi

2 ) end ( m' ).xi  end ( m ).xi
m’ e m si dicono colineari se m’<<m oppure m<<m’

ossia, appaiono nello stesso ordine in tutte le sequenze
77
Match chaining problem su k sequenze


Input
 M = {m1, m2, … , mn} una lista di match su K ( 2) sequenze:

S1, S2, ... , Sk-1 sequenze di DNA

Sk sequenza di cDNA
Output

{m
mc2
Catena di match colineari C C
={m
…,m
, cmf cf} , cioè
c2 ,,
c1c,1 ,m
mc1
 mmc2c2
<<
mc f , tale che la funzione
m
<<…
c1 <<
cf
c f 1
i 1
i 1
score C    mi .length  overlapk (mi , mi 1 )
sia massimizzata
78
Soluzione con grafo

Il match chaining problem può essere ridotto al
problema della “ricerca del cammino di costo
massimo su un grafo orientato aciclico”

Un grafo G(V, E) è definito da una coppia di
insiemi



V, insieme dei nodi (o vertici)
E, insieme degli archi che collegano coppie di
nodi: E  (V, V)
Un grafo G(V, E) si dice pesato quando ad ogni
arco ei  E è associato un valore w(ei)
79
G(V,E) rappresentazione grafica
V = { a, b, c, d, e, f }
E = { (a, b), (b, c), (b, f), (c, f), (e, f),}
a
w((a, b))=5
5
w((b, c))=3
w((b, f))=8
w((c, f))=0
w((e, f))=2
d
8
b
3
f
2
e
0
c
non pesato
Grafo pesato
80
Nozioni sui grafi

Un grafo G(V, E) è orientato quando i suoi archi hanno un verso,
cioè sono percorribili in un solo senso
(vi, vj)  E, è un arco percorribile dal nodo vi al nodo vj
 (vi, vj)  (vj, vi)


Un cammino di lunghezza l in G(V, E) è dato da una sequenza di
vertici v0, v1, ... , vl tali che
(vi, vi+1)  E ,  i = 0, … , l – 1

Il costo di un cammino v0, v1, ... , vl su un grafo pesato G(V, E) è
dato dalla somma dei pesi degli archi che lo compongono
l 1
w(v0 , , vl )   w((vi , vi 1 ))
i 0
81
Nozioni sui grafi(2)

Un cammino è chiuso (ciclo) quando v0 = vl,
cioè quando inizia e termina sullo stesso nodo.

Un grafo aciclico è un grafo in cui non esistono
cammini chiusi
b
f
e
c
Grafo orientato ciclico
b
f
e
c
Grafo orientato aciclico
82
G(V, E) rappresentazione grafica
V = { a, b, c, d, e, f }
E = { (a, b), (b, c), (b, f), (c, f), (e, f),}
a
w((a, b))=5
5
w((b, c))=3
w((b, f))=8
w((c, f))=0
w((e, f))=2
d
8
b
3
f
0
c
2
e
esempio:
w(b,c,f,e)= 3 + 0 + 2 = 5
Grafo pesato orientato aciclico
83
Riduzione a grafo

Costruiamo un grafo G(V, E) tale che
 V rappresenta l’insieme dei match
 E rappresenta le relazioni di precedenza tra i match

ogni match in input è rappresentato come un nodo nel grafo G
 sia v’ il nodo che rappresenta m’
 sia v il nodo che rappresenta m



esiste in E un arco orientato (v’, v) solo se m’ << m
w((v’, v)) = m.length – overlapk(m’, m)
si aggiungono in V i nodi 0 e t tali che

il nodo 0



ha un arco uscente verso ogni altro nodo vi = mi del grafo
w((0, vi)) = mi.length
il nodo t


ha un arco entrante proveniente da ogni altro nodo vi del grafo
w((vi, t)) = 0
84
G(V, E) rappresentazione grafica
V = { m1, m2, m3, m4, m5, m6,}0,
0 }t }
),}(0,m1), … , (0,m6),
E = { (m2, m1), (m2, m6), (m2, m3), (m2, m5), (m3, m5), (m3, m6), (m6, m5),
)} (m1,t), … , (m6,t)}
3
m2
6
m4
0
0
5
5
0
7
0
8
m1
m5
2
m6
t
0
5
3
4
3
0
4
3
0
m3
“ricerca del cammino di costo massimo da
0 a t sul grafo pesato orientato aciclico”
85
Risoluzione del problema sul grafo
È
un problema classico della teoria
dei grafi risolto in tempo O(n2), con
programmazione dinamica

Inefficiente con un numero di frammenti
n molto alto
 La
riduzione costa comunque tempo
O(n2)
86
Sommario
1.
2.
3.
Introduzione
Terminologia e definizioni
Match Chaining Problem su 2 sequenze
1.
2.
3.
Algoritmo base
Range Maximum Query (RMQ)
Algoritmo MCCM
1.
2.
4.
Match Chaining Problem su k sequenze
1.
Terminologia e definizioni (2)
2.
Concatenzione di match non overlapping
1.
3.
5.
Limitazione della distanza tra i match
Overlap
RMQ a dimensioni multiple
Concatenzione di match overlapping
1.
Riduzione a grafo
2.
Algoritmo basato su geometria
3.
Risultati Sperimentali
Riferimenti
87
Soluzione geometrica

Usando una soluzione basata sulla geometria è
possibile ridurre la complessità in tempo

Concetto di regione (k-dimensionale):
S2
(cDNA)
C(m)
m

A(m) = Rect (0,beg(m))
B(m) = ([beg(m).x1 … end(m).x1], [0 … beg(m).x2])
C(m) = ([0 … end(m).x1], [beg(m).x2 … end(m).x2])

AB(m) = A(m)  B(m)


A(m)
B(m)
S1 (genoma)
88
Idea

ordinare i match in input m1, ..., mn e processarli in
base ad un criterio di ordinamento.

Per ogni match mi:
 cercare il miglior match predecessore tra quelli già
processati, distinguendo in base al fatto che abbia o
meno sovrapposizione sul pattern (cDNA)
 calcolare score(mi)
 marcare il match mi come già processato (attivo)
89
Considerazioni iniziali

Consideriamo

i match come coppie di punti a k dimensioni
(beg(m),end(m)) con



beg(m) = beg(m).x1, … , beg(m).xk
end(m) = end(m).x1, … , end(m).xk
D1, D2 sono strutture dati semi-dinamiche che
supportano RMQ a 2k dimensioni con attivazione
per manipolare i match (RMQD1,RMQD2)


RMQ con attivazione: nel calcolo dell’RMQ si considerano
solo i match attivi, ignorando quelli inattivi
activate(D, m, p): procedura che attiva il match m (già
inserito nella struttura dati D) con priorità p
90
Algoritmo
Ordiniamo i punti di inizio dei match in base alla loro
coordinata x1
Tutti i match sono caricati come inattivi sia in D1 che in
D2
Processiamo i match come a simulare una linea che
percorre i punti sulla coordinata x1 (line-sweep
paradigm)
1.
2.
3.

4.
Quando il match mi viene processato:

gli viene associato prev(mi)

viene calcolato score(mi)
 mi viene attivato in D1 e in D2
Prendiamo il match di massimo score e costruiamo la
catena utilizzando l’informazione sul predecessore
91
Ricerca del predecessore

Si sfrutta il concetto di regione per individuare i match
già attivi che sono in relazione << con il match corrente
mi, e tra questi distinguiamo quelli che hanno overlap
sul pattern da quelli che non lo hanno:

mj << mi sse



beg(mj)  A(mi) and
end(mj)  { A(mi)  B(mi)  C(mi) }
mj ha overlapk con mi sse

end(mj)  C(mi)
S2
(cDNA)
m1
Esempio:
 mj è un match non (ancora) attivo
 non è vero che m1 << mi
 m2 << mi e non ha overlapk
 m3 << mi ed ha overlapk
mi
C(mi )
mj
m3
m2
A(mi )
B(mi )
S1 (genoma)
92
Ricerca del predecessore(2)

La ricerca dei match candidati viene fatta utilizzando RMQ a 2k
dimensioni su due regioni disgiunte

le prime k dimensioni devono assicurare che beg(m’)  A(mi)


le ultime k dimensioni devono assicurare, distinguendo i due casi, che

end(m’)  {AB(mi)} oppure

end(m’)  C(mi)
m1



m2


= RMQD1(A(mi), {AB(mi)})
per la ricerca del match di massimo
score che non ha sovrapposizione sul
pattern (cDNA) con il match corrente mi
end(m1).xk< beg(mi).xk
= RMQD2(A(mi), C(mi))
per la ricerca del match di massimo
score che ha sovrapposizione sul
pattern (cDNA) con il match corrente mi
beg(mi).xk≤ end(m2).xk ≤ end(mi).xk
mi
S2
(cDNA)
C(mi )
m2
A(mi )
m1
B(mi )
S1 (genoma)
93
Ricerca del predecessore(3)
Tra i 2 match candidati si sceglie quello che massimizza lo score di
mi:




score1 = m1.score + mi.length
score2 = m2.score + mi.length – overlapk(m2, mi)
if (score1 ≥ score2) then



prev(mi)=m1
score(mi)=score1
else

prev(mi)=m2

score(mi)=score2
mi
S2
(cDNA)
C(mi )
m2
A(mi )
m1
B(mi )
S1 (genoma)
94
Attivazione di un match

Dopo aver stabilito il prev(mi) e score(mi), si deve attivare mi:

in D1 il match mi viene attivato con priorità score(mi)

in D2 il match mi viene attivato con priorità score(mi) – end(mi).xk
Perchè usiamo priorità diverse?


Quando si calcola RMQD1 si cercano frammenti senza overlap su cDNA con il
match corrente mi, perciò quello che massimizza score(mi) sarà quello col più alto
score


score(mi) = score(m1) + mi.length
Con RMQD2 si cerca il match migliore tra quelli che hanno overlap con il match
corrente mi, perciò non è detto che il match di score massimo in D2 sia quello che
massimizza lo score(mi)

score(mi) = score(m2) + mi.length – overlapk(m2, m)
= score(m2) + mi.length – ( end(m2).xk – beg(mi).xk )
= score(m2) – end(m2).xk + mi.length + beg(mi).xk
95
Ricerca della catena ottima

Esempio

m1


Input: M = {m1, m2, m3, m4, m5, m6, m7, m8}
5
m2
11
m3
20
m4
12
m5
9
m6
23
m7
16
m8
8
Output: C = {m1, m2, m4, m7, m6}
Questa parte dell’algoritmo ha costo O(n)
NB: l’ordine dei match in input non è quello con cui sono processati nella fase precedente
96
Complessità

La complessità dell’algoritmo dipende dalla struttura dati
utilizzata per implementare l’RMQ con attivazione in D1 e D2.

Range tree:


eseguire RMQ costa O(log2k-1 n log (log n)) tempo e O(log2k-1 n) spazio
L’algoritmo esegue 2 RMQ per ogni match, quindi calcola per O(n)
volte l’RMQ
La complessità totale è
O(n log2k-1 n log (log n)) tempo
e
O(n log2k-1 n) spazio
97
Abbassiamo la complessità

È possibile ridurre la complessità dell’algoritmo riducendo le
dimensioni dell’RMQ in D1 e D2

il primo range [0 … beg(mi).x1] può essere ignorato perché i match sono
processati in ordine crescente secondo la coordinata x1 del punto di
inizio


tutti i match già attivi rispettano questo range
il k-esimo range [0 … beg(mi).xk] può essere ignorato

in D1:


con il k-esimo range si vuole un match m1 tale che beg(m1).xk < beg(mi).xk, ma
questa condizione è ridondante in quanto il 2k-esimo range controlla che
end(m1).xk < beg(mi).xk
in D2 si possono verificare 2 casi:


beg(m2).xk < beg(mi).xk , OK
beg(m2).xk ≥ beg(mi).xk, in questo caso m2 deve essere scartato. Questo non
compromette l’ottimalità della catena poiché non esiste un altro match che termina
in C(mi) che possa dare ad mi uno score migliore di uno che temina in AB(mi)
La complessità scende a O(n log2k-3 n log (log n)) tempo
e O(n log2k-3 n) spazio
98
Casi speciali

L’algoritmo può essere specializzato per
1.
usare in input match di tipo multi_MUM, riducendo la
complessità
2.
limitare la quantità di overlap ammissibile tra due
match consecutivi nella catena di output, riducendo la
complessità
99
Casi speciali: multi-MUM

Maximal Multiple Unique Match: match con valore di rarità 1, cioè
per ogni sequenza Si la sottostringa Si[pi … pi + l - 1] occorre in Si
esattamente una sola volta
Conseguenze:

basta controllare le k dimensioni sulla
fine dei match in quanto
S2

(cDNA)

anche il primo range [0 … end(m).x1]
può essere ignorato, a causa del
paradigma line-sweep

S1
(genoma)

end(m’) < end(m)  beg(m’) < beg(m)
ogni m’ con end(m’).x1 > end(m).x1 avrà
beg(m’).x1 > beg(m).x1, quindi m’ non è
stato ancora attivato
D1 e D2 contengono punti a k - 1
dimensioni (end(m).x2,…,end(m).xk)
La complessità si riduce a O(n logk-2 n log (log n))
100
Casi speciali: limitazione
sull’overlap
Assunzioni:
 l: lunghezza minima di un match
 l -1: massimo overlap fra due match consecutivi nella catena
Per trovare il predecessore m’ del match m, imponiamo che
i  1k end (m' ).xi  beg (m).xi  l
ma poiché l è la lunghezza minima di un match, sicuramente avremo che
i  1k end (m' ).xi  beg (m).xi  l
m
S2
(cDNA)
 beg (m' )  A(m)
C(m)
l
A(m)
B(m)
101
l
S1 (genoma)
Casi speciali: limitazione sull’overlap(2)
i  1k end (m' ).xi  beg (m).xi  l
 beg (m' )  A(m)
S2
(cDNA)
v(m)
m
l
OK
OK
A (m)
OK
l : lunghezza minima di un match
l
S1 (genoma)
Conseguenze:
 basta controllare le k dimensioni
sulla fine dei match!
 se applichiamo il paradigma linesweep in base al punto virtuale
v(m) = (beg(m).x1+l,…,beg(m).xk+l)
si può eliminare anche la prima
dimensione
 D1 e D2 contengono punti a k - 1
dimensioni (end(m).x2,…,end(m).xk )
La complessità si riduce a O(n logk-2 n log (log n))
102
Sommario
1.
2.
3.
Introduzione
Terminologia e definizioni
Match Chaining Problem su 2 sequenze
1.
2.
3.
Algoritmo base
Range Maximum Query (RMQ)
Algoritmo MCCM
1.
2.
4.
Match Chaining Problem su k sequenze
1.
Terminologia e definizioni (2)
2.
Concatenzione di match non overlapping
1.
3.
RMQ a dimensioni multiple
Concatenzione di match overlapping
2.
Riduzione a grafo
Algoritmo basato su geometria
3.
Risultati Sperimentali
1.
5.
Limitazione della distanza tra i match
Overlap
Riferimenti
103
Risultati sperimentali

Per gli esperimenti sono stati utilizzati come
sequenze:
 il cDNA del D.melanogaster (23 Mbp)
 cromosoma 2R del
D.melanogaster (21
Mbp)
 cromosoma 2R del D.simulans (20 Mbp)

Con il programma ramaco si sono generati
frammenti di lunghezza minima 15.

Sono stati confrontati 3 algoritmi su differenti
valori di rarità massima dei match in input:
 soluzione basata su grafo (graph)
 soluzione geometrica (geom1)
 soluzione geometrica con limite di overlap
(geom2)
104
Riferimenti
1.
Match Chaining Algorithms for cDNA Mapping di
T.Shibuya e I.Kurochkin, proceedings of 3rd Workshop
on Algorithms in Bioinformatics (WABI), 462-475,
2003.
2.
Chaining Algorithms and applications to comparative
genomics di M.I.Abouelhoda e E.Ohlebusch, Journal
of Discrete Algorithms 3(2-4), 321-341, 2005.
3.
A chaining Algorithm for Mapping cDNA Sequences to
Multiple Genomic Sequences di M.Abouelhoda,
proceedings of 14th International Symposium on String
Processing and Information Retrieval (SPIRE), 1-13,
2007.
105
Scarica

Lucidi