Seminari 4
Introducció als algorismes (II)

Bioinformàtica - 1er trimestre curs 2012/2013 - UPF

En el següent enllaç trobareu els apunts en format PDF corresponents a la part d'introducció als algorismes.

Aquesta pràctica conclou la introducció del tema de la programació i el desenvolupament d'algorismes. Els conceptes més importants que heu d'entendre són:

Eficiència dels algorismes

Quan ens enfrontem amb la tasca de dissenyar un algorisme podem trobar més d'una manera d'escriure aquest algorisme per resoldre el mateix problema.


Dos algorismes differents per al mateix problema poden ser també differents en com de ràpid resolen el problema. Aleshores es diu que un algorisme és més eficient que l'altre.


Per exemple, considerem el problema de cercar una paraula a un diccionari. La forma més simple de cercar una paraula a un diccionari és llegir la primera paraula del diccionari i comprovar si és la paraula que estem cercant. Si ho és, haurem acabat. Si no ho és, llegirem la segona paraula. Si la segona paraula és la que estem cercant, haurem acabat. Si no ho és, llegirem la tercera, i així repetidament fins trobar la paraula que busquem.


Si, per exemple, busquem d'aquesta forma la paraula tronc, quan l'haurem trobat haurem llegit totes les paraules del diccionari entre la lletra A i la lletra T a mes a mes de les que comencen per T i són alfabèticament anteriors a tronc. Qualsevol diccionari de butxaca sol tenir unes 50 000 paraules, per tant per buscar la paraula patio podem haver llegit fàcil ment unes 40 000 paraules.


ara pensem com fem realment nosaltres la cerca de paraules al diccionari. La manera en què ho fem es similar al que es coneix com cerca binaria.


cerca binaria: considerem un rang de paraules, per exemple les que hi ha de la lletra A a la lletra H, que anomenarem finestra.


inicialment la finestra correspon al rang de paraules de la A a la Z, es a dir, el diccionari sencer. Tal i com l'algorisme va executant-se, la finestra s'anirà escurçant, de vegades per l'esquerra, de vegades per la dreta (o per dalt i per baix si us imagineu la finestra verticalment). Arribarem a un punt en el qual o be la finestra conté una única paraula que és la que estem buscant, o be la finestra és buida i vol dir que la paraula no és al diccionari. 

@dic = ("Genis", "Mar", "Pep", "Robert", "Roderic", "Sergi", "Xavi");
$par = "Pep";
$min = 0;
$max = 6;
$trobada = 0;

while ($min <= $max && !$trobada) {

  $mig = int( ($min+$max)/2 );

  if ($dic[$mig] lt $par) {
    $min = $mig + 1;
  }
  else {

    if ($dic[$mig] gt $par) {
      $max = $mig - 1;
    }
    else {
      $trobada = 1;
    }

  }

}

if ($trobada == 1) {
  print "paraula trobada\n";
}
else {
  print "paraula no trobada\n";
}

En l'exemple anterior on buscavem la paraula tronc al diccionari de butxaca, haviem de llegir unes 40 000 paraules d'un diccionari de 50 000. Quantes hauriem de llegir, fent el mateix, si el diccionari tingués el doble de paraules, es a dir, 100 000 ?


assumint que el tamany del diccionari creix proporcionalment per cada lletra de l'abecedari, hauriem de llegir 80 000 paraules per trobar tronc. En un cas com aquest un diu que l'algorisme té un ordre de creixement lineal, i ho anotarem ${\cal O}(N)$. També es sol dir que la complexitat de l'algorisme és d'ordre lineal. Això implica que el nombre d'operacions que l'algorisme ha de fer es proporcional al nombre d'elements que ha de tractar l'algorisme.


ara pensem quin esforç extra hauria de fer l'algorisme de cerca binaria si en lloc de utilitzar-lo amb un diccionari de 50 000 paraules l'utilitzéssim amb un de 100 000.


doncs, simplement, una volta extra de la composició iterativa (el while). La cerca binaria pertany a la classe d'algorismes que que tenen un ordre de creixement logarítmic, que l'anotarem com ${\cal O}(\log N)$.

\includegraphics{logx}

és important adonar-se'n que la noció de complexitat, o ordre de creixement, s'interpreta en termes del comportament d'un algorisme en general. Podem pensar fàcilment d'algun cas en el que la cerca simple al diccionari sigui més eficient que la cerca binària, i no per aixo l'algorisme de cerca binaria deixa de ser més "eficient" que l'algorisme de cerca simple.


per tal de tenir una idea més precisa de fins a quin punt l'ordre de creixement determina l'eficiència dels algorismes, a la seguent taula trobarem diversos ordres (funcions) de creixement amb el seu valor per $N=1 000 000$.

Funció de creixement Valor per $N=1 000 000$
1 1
$\log N$ 13.8
$\sqrt N$ 1,000
$N$ 1,000,000
$N\log N$ 13,815,510
$N^2$ 1,000,000,000,000
$N^3$ 1,000,000,000,000,000,000
$2^N$ un número amb 693,148 dígits

considerem la següent matriu de 4 x 4 elements, que en llenguatge Perl l'especificarem com a un vector de vectors. 


        @v=([11,12,13,14],
            [21,22,23,24],
            [31,32,33,34],
            [41,42,43,44]);


la forma en la que ens referim a un element en concret d'una matriu en Perl es mitjançant la notació $v[2][3] on, en aquest cas, el 2 es refereix a la tercera posició a la primera dimensió i el 3 es refereix a la quarta posició de la segona dimensió, es a dir, l'element 34.


penseu ara un algorisme que, donada la matriu @v ens imprimeix la suma dels elements en cadascuna de les quatres "files" corresponents a la primera dimensió. Es a dir:


\begin{displaymath}11 + 12 + 13 + 14 \end{displaymath}


\begin{displaymath}21 + 22 + 23 + 24 \end{displaymath}


\begin{displaymath}31 + 32 + 33 + 34 \end{displaymath}


\begin{displaymath}41 + 42 + 43 + 44 \end{displaymath}

l'algorisme doncs, podria ser el seguent: 

$i = 0;

while ($i < 4) {
  $s = 0;
  $j = 0;
  while ($j < 4) {
    $s = $s + $v[$i][$j]; 
    $j = $j + 1;
  }
  print $s,"\n";
  $i = $i + 1;
}

quina és la complexitat (ordre de creixement) d'aquest algorisme ?

(a)
${\cal O}(N)$ lineal ?
(b)
${\cal O}(N^2)$ quadràtica ?
(c)
${\cal O}(N^3)$ cúbica ?
(d)
${\cal O}(2^N)$ exponencial ?

teniu en compte que si la matriu en lloc de ser de dimensions $n\times n$ (quadrada), fos $n\times m$ amb $n\neq m$ (rectangular), aleshores escriuriem l'ordre de creixement com:


\begin{displaymath}{\cal O}(nm) \end{displaymath}

Correspondència entre seqüències de símbols

motivació: en seqüències de DNA, RNA o aminoàcids, un grau similaritat alt entre seqüències normalment implica també un grau alt de similaritat estructural i/o funcional.


donades dues seqüències de símbols (per exemple de DNA): 

                         t={AATGC}
                         p={AGGC}

anomenarem com a problema de la correspondència entre seqüències de símbols (en anglés string matching problem) al problema de trobar si t apareix dins de p, o p dins de t, o en quina mesura t i p són similars.


aquest problema té dues variants principals:

Correspondència exacta

considerem una seqüència 

t={ATGCATAATGCGTCA}

i una seqüència 

p={ATAA}

en aquest cas la seqüència p es troba dins de t a partir del cinqué símbol i llavors direm que la seqüència p apareix amb desplaçament $s$ o que apareix a partir de la posició $s+1$ on, en aquest cas, $s=4$. També anomenarem a p el patró a buscar dins de t.

\includegraphics[width=\textwidth]{exactm}

en general, suposem que t té $n$ símbols i p en té $m$. El problema de la correspondència exacta entre seqüències de símbols correspon al problema de trobar tots els possibles desplaçaments $s$ que fan que els símbols 


                  t[s+1]t[s+2]t[s+3]...t[s+m]


siguin exactament els de p, es a dir: 


                  t[s+1] == p[0]
                  t[s+2] == p[1]
                  t[s+3] == p[2]
                  ...
                  t[s+m] == p[m-1]


quants possibles desplaçaments $s$ com a màxim podem arribar a trobar ?


\begin{displaymath}n-m+1 \end{displaymath}


quin seria l'algorisme més senzill per trobar tots els desplaçaments en els quals la seqüència p apareix exactament dins de t? 

$s=0;

while ($s <= $n-$m) {
  $i = 0;

  while ($i < $m  &&  $t[$s+$i] eq $p[$i]) {
    $i = $i + 1;
  }

  if ($i == $m) {
    print "patro trobat amb desplacament ",$s,"\n";
  }

  $s = $s + 1;
}

quina és la complexitat de l'algorisme anterior ?


\begin{displaymath}{\cal O}((n-m+1)m) \end{displaymath}

aquest algorisme simple no és el més eficient, existeixen algorismes amb ordre de creixement ${\cal O}(n+m)$.

Correspondència inexacta

trobar una correspondència inexacta, o aproximada, de dues seqüències de símbols és de tant o més interés que intentar trobar una correspondència exacta que moltes vegades no existeix. 

AASRPRSGVPAQSDSDPCQNLAATPIPSRPPSSQSADARQGRWGP
      || | |      |  |  | |  || |        
      SGAPGQRGEPGPQGHAGAPGPPGPPGSD

idea fonamental: assignar una puntuació (score) a cadascuna de les possibles correspondències inexactes. Aquesta puntuació serà la suma de puntuacions individuals de cada parell de símbols que estan a la mateixa posició. Per exemple, suposem que puntuem amb 1 quan dos símbols a la mateixa posició són iguals, i amb 0, si no ho són: 

                          ATCGCA
                           ATGTA
                          000101 = 2

buscarem, aleshores, la correspondència que ens maximitza la puntuació. Respecte a l'exemple anterior, podriem trobar una correspondència inexacta millor si introduïm salts (gaps) dins de la seqüència més curta: 

                          ATCGCA
                          AT-GTA
                          110101 = 4

un altre exemple podria ser el següent: 

AASRPRSGVPAQSDSDPCQNLAATPIPSRPPSSQSCQKCRADARQGRWGP
      || | |      |  |  | |  || | 
      SGAPGQRGEPGPQGHAGAPGPPGPPGSDGSPARKG

AASRPRSGVPAQSDSDPCQNLAATPIPSRPPSSQSCQKCRADARQGRWGP
      || | |      |  |  | |  || |         || |
      SGAPGQRGEPGPQGHAGAPGPPGPPGSDG-----SPARKG

el problema de trobar una correspondència inexacta òptima entre dues seqüències de símbols, introduïnt salts, es coneix a la biologia com el problema de l'alineament de seqüències, o sequence alignment problem en anglès.


una forma de trobar aquesta correspondència inexacta òptima seria per força bruta:

per tenir una idea de la complexitat d'aquest algorisme, simplement penseu que per una seqüència de $N$ símbols podem inserir salts en aproximadament $2^N$ llocs differents dins la seqüència.


per cada parell de seqüències, hem de generar tots els possibles alineaments, que son $N$. Per cada alineament hem de calcular la seva puntuació que requerirà de l'ordre de $N$ operacions. Es a dir, per cada parell haurem d'executar un conjunt d'accions amb una complexitat d'ordre ${\cal O}(N^2)$, i tenim $2^N\times 2^N=2^{2N}$ parells. Globalment, la complexitat del problema d'alinear seqüències es, en principi:


\begin{displaymath}{\cal O}(2^{2N}N^2) \end{displaymath}

mes concretament, s'ha de tenir en compte que molts d'aquests alineaments no tenen sentit des d'un punt de vista biològic, amb lo qual la complexitat real del problema de l'alineament de seqüències (o correspondència inexacta de seqüències) és (Waterman, 1984):


\begin{displaymath}{\cal O}\left(\frac{2^{2N}}{4\sqrt{N\pi}}\right) \end{displaymath}

per exemple, per dues seqüències de longitut $N=1,000$ hauriem de fer de l'ordre de $10^{600}$ operacions!!


aquest problema té una solució algorísmica més eficient (ordre ${\cal O}(nm)$) mitjançant el que es coneix com programació dinàmica (dynamic programming).


considerem que alineament de l'esquerra és òptim. La puntuació (score) total d'aquest alineament és la puntuació de les primeres 4 posicions més la puntuació de la última posició, tal i com s'especifica al mig. 

        AATGC            AATG   C            AATGC
        |  ||            |  | + |            |  ||
        AG-GC            AG-G   C            A-GGC

idea fonamental: si l'alineament de totes les posicions és òptim, aleshores l'alineament de les primeres quatre posicions també ho és.


si això no fos veritat de forma que, per exemple, alineant G i T donés una puntuació més alta, aleshores l'alineament de la dreta tindria una puntuació global més alta que el de l'esquerra, contradint el supòsit inicial.

Aquest argument és equivalent a considerar, per exemple, les següents distàncies en Km. dels diferents trajectes que podem fer amb tren entre les ciutats de Lleida, Valls, Tarragona, San Vicenç i Barcelona:

\includegraphics{lleidabcn}

i adonar-nos que si el trajecte més curt entre Lleida i Barcelona, passa per Valls, forçosament el trajecte més curt entre Lleida i Sant Vicenç també ha de passar per Valls. Plantejeu-vos, per exemple, que si no fos així i el trajecte més curt entre Lleida i San Vicenç passés per Tarragona, que implicaria això respecte al trajecte entre Lleida i Barcelona que passa per Valls?


tornant al problema de l'alineament de seqüències, es tracta d'explotar el fet de que la puntuació d'un alineament és la suma de les alineacions individuals de cada parell de símbols: 

        AATGC            A   A   T   G   C
        |  ||            | +   +   + | + |
        AG-GC            A   G   -   G   C

si ens fixem be, fent l'alineament símbol per símbol, ens adonarem que a cada pas només tenim 3 possibilitats:

l'algorisme de programació dinàmica per la correspondència inexacta, o alineament, de seqüències consisteix aleshores en construir la matriu sencera d'una seqüència contra l'altra.


Comencarem per la fila 1 columna 1 i, d'esquerra a dreta, i de dalt a baix, a cada posició considerarem les tres possibles propagacions.


De les tres possibles propagacions calcularem quines puntuacions generen cadascuna, i escollirem la propagació de màxima puntuació.


un cop hem construït la matriu sencera, començarem per la posició de la última fila i la última columna, i resseguint el camí de màxima puntuació arribarem a la posició de la primera fila i primera columna, recuperant l'alineament òptim al temps que escrivim els alineaments individuals d'acord amb el tipus de propagació dut a terme (diagonal, desplaçament vertical o desplaçament horitzontal).


per les seqüències anteriors AATGC i AGGC, intenteu construïr la matriu per trobar l'alineament òptim ficant la seqüència AATGC a la dimensió horitzontal i la seqüència AGGC a la vertical. Considereu una penalització de salt nul.la $w(-)=0$, i una puntuació de 1 quan dos símbols són iguals i de 0 quan no ho són. La matriu en qüestió seria aquesta:

\includegraphics{exalin}

la reconstrucció de l'alineament òptim és fa de la següent forma: 

                                        C
                                        |        diagonal
                                        C

                                      - C
                                        |        desp. vertical
                                      G C

                                    G - C
                                        |        diagonal
                                    G G C

                                  T G - C
                                        |        desp. horitzontal
                                  - G G C

                                A T G - C
                                        |        desp. horitzontal
                                - - G G C

                              A A T G - C
                                        |        diagonal (inici)
                              A - - G G C

com podem apreciar, el cost de tot el procés és a la construcció de la matriu, que té una complexitat ${\cal O}(nm)$ donat que les puntuacions es van reutilitzant tal i com anem construïnt la matriu.