Analisi di sequenze genomiche: stringhe e algoritmi
Una stringa è una successione ordinata di caratteri o simboli derivata da un insieme finito chiamato alfabeto. Il termine sequenza viene normalmente usato come sinonimo di stringa.
- Alfabeto del DNA: {A, C, G, T}
- Alfabeto delle proteine: {20 caratteri indicanti gli amminoacidi}
La lunghezza s di una stringa, indicata come |s| (modulo di s), corrisponde al numero di caratteri che contiene. Il carattere presente nella posizione i viene indicato con s[i]. Gli indici dei caratteri vanno da 1 a |s|. Es. |4| Una stringa di lunghezza nulla è detta stringa vuota.
Sottosequenza e sottostringa
Benché stringa e sequenza abbiano lo stesso significato, i termini sottosequenza e sottostringa indicano concetti diversi:
- Una sottosequenza di s è una sequenza che può essere ottenuta da s rimuovendo alcuni caratteri. Data s=ATATAT, TTT è una sottosequenza di s; TAAA non lo è (posso solo rimuovere, non spostare). Se la sequenza t è una sottosequenza di s, allora s è una supersequenza di t.
- Una sottostringa di s è una stringa formata da caratteri consecutivi di s, ordinati come in s. Tutte le sottostringhe sono sottosequenze, ma non vale il viceversa: data s=AGTACA, TAC è una sottostringa di s; TTGAC non lo è. Se la sequenza t è una sottostringa di s, allora s è una superstringa di t.
Un intervallo in una stringa s è un insieme di indici consecutivi [i..j] tale che 1≤i
Il concetto di intervallo è utile per distinguere le sottostringhe ripetute (es.: TT ricorre più volte in CTTTAGCATTAA). Se i=j+1, si ha la stringa nulla.
La concatenazione di due stringhe s e t viene indicata con st ed è formata aggiungendo tutti i caratteri di t alla fine di s, nello stesso ordine in cui appaiono in t. La lunghezza di st è data da |s|+|t|.
Definiamo prefisso di s una sottostringa del tipo s[1..j] con 0≤j<|s|, mentre con suffisso di s ci si riferisce a una stringa del tipo s[i..|s|] con 1≤i<|s|+1.
Algoritmi e complessità
Un algoritmo è una sequenza ben definita e limitata di azioni volte a risolvere un problema. Deve poter trasformare un input (non nullo) in un qualche tipo di risposta (output) o di azione, e deve fornire la risposta in un tempo finito.
Il primo passaggio nella definizione di un algoritmo è la trasformazione della formulazione del problema dal linguaggio naturale a un linguaggio formale utilizzabile da una persona o da una macchina per risolvere il problema dato. Ogni passaggio deve essere formalizzato in maniera precisa e ogni passaggio deve poter essere eseguito da un sistema automatico.
Analisi di un algoritmo
- Verificare che il risultato fornito è corretto.
- Determinare il suo tempo di esecuzione.
L'efficienza di un algoritmo in termini di tempo di esecuzione è determinata dal numero di operazioni necessarie per eseguirlo. Quante più volte è ripetuta una stessa operazione, tanto maggiore sarà il tempo richiesto. La stima del tempo di esecuzione viene determinata in funzione della dimensione del problema da risolvere; la dimensione del problema viene descritta dal termine di ordine maggiore (funzione di complessità nel tempo). L'ordine di un algoritmo viene descritto dal suo comportamento asintotico.
- Un algoritmo O(n2) sarà più veloce di un algoritmo O(n) per valori piccoli di n.
- Per valori di n sufficientemente grandi, l’algoritmo O(n) sarà più veloce.
Algoritmi i cui tempi computazionali sono O(n) sono detti lineari, mentre quelli O(n2) sono detti quadratici.
Date due stringhe di 5 e 6 lettere, ci sono 5×6=30 possibili appaiamenti tra le varie lettere. Il numero di operazioni necessarie per esplorare tutti gli appaiamenti è proporzionale alla lunghezza delle due stringhe. L’algoritmo è di ordine O(nm).
Per molti problemi della bioinformatica sono stati messi a punto algoritmi che richiedono un numero di operazioni O(n2). Altri problemi sono risolvibili solo con algoritmi più complessi, O(n3), che richiedono tempi di esecuzione maggiori e quindi sono utili solo per valori di n relativamente piccoli.
Algoritmi di ordine O(kn), in cui n è l’esponente di una costante, non hanno generalmente alcuna applicazione pratica richiedendo un tempo di esecuzione che cresce esponenzialmente al crescere delle dimensioni del problema (ad es.: al crescere della lunghezza delle sequenze).
Classi dei problemi
- Appartengono alla classe P tutti i problemi di decisione che possono essere risolti da algoritmi di tipo polinomiale (polynomial time algorithm). I problemi in P sono ben risolti poiché per essi esiste un algoritmo efficiente.
- Appartengono alla classe NP tutti i problemi di decisione per cui, essendo nota una soluzione, è possibile verificare in tempo polinomiale che la risposta è affermativa (non-deterministic polynomial time algorithm). Di fatto nessuno è stato in grado di trovare un algoritmo efficiente per risolvere tali problemi, ma qualora si potesse provare la sua esistenza, la soluzione sarebbe verificabile in un tempo polinomiale.
Metodi predittivi
Il programma genetico dispone dei dati e delle informazioni per:
- Codificare/costruire tutte le proteine.
- Attivare/disattivare i geni necessari per reagire opportunamente agli eventi (promotori).
- Codificare per tutti gli altri elementi (rRNA, tRNA, telomeri, etc).
Uno degli aspetti più interessanti dell’analisi di un genoma consiste nell’identificare le regioni codificanti proteine. Nei procarioti le regioni codificanti possono essere facilmente identificate per la presenza di lunghi tratti privi di codoni di terminazione. Negli eucarioti questo approccio è reso difficile dalla presenza di sequenze introniche intercalate all’interno delle regioni codificanti. Dei 64 codoni, 3 rappresentano segnali di terminazione (TAA, TAG, TGA).
La stringa è abbastanza lunga per poter dire che è codificante? Se tutto è casuale, in una sequenza casuale di A, T, C e G si può supporre di trovare un segnale di terminazione in media ogni 21,3 codoni.
Leggiamo la sequenza e la sua complementare partendo dalla prima base e procedendo 3 basi alla volta (123, 456, 789, ...).
- Si tracci una barretta corta ogni volta che si incontra un codone d’inizio (ATG).
- Si tracci una barretta lunga ogni volta che si incontra un codone di terminazione.
- Si ripeta per ciascuna delle altre due fasi di lettura (234, 567, ...) e (345, 678, ...).
Le regioni estese prima dei codoni di terminazione sono dette fasi di lettura aperte (open reading frames o ORF). La presenza di una ORF non è sempre sufficiente ad identificare con certezza un gene; quelli che codificano per proteine corte possono confondersi con tratti di sequenza casualmente privi di codoni di terminazione. Non è sempre chiaro quale sia il punto d’inizio della regione codificante. Infatti, non è detto che il primo ATG sia necessariamente il primo codone codificante o che uno stop-codon non venga utilizzato per codificare un amminoacido e non si può nemmeno escludere la possibilità che un gene abbia due siti di inizio alternativi.
Per essere certi dell’ipotetica funzione codificante di una ORF e della sua esatta estensione non è sufficiente fare un’analisi dei codoni di inizio e di terminazione. Va ricercato un modulo, ovvero si ricerca uno o più pattern caratterizzanti e motivi funzionali. Le sequenze codificanti non sono infatti sequenze casuali; nelle proteine alcuni amminoacidi sono più abbondanti di altri, la maggior parte degli amminoacidi è codificata da diversi codoni, perciò le regioni codificanti hanno quindi sequenze particolari.
Alcuni approcci si basano quindi sull’osservazione che alcune regioni di poche basi hanno una probabilità molto diversa di essere trovate in regioni codificanti piuttosto che in regioni non codificanti.
GeneMark (Hidden Markov Models)
Supponiamo di avere una serie di esempi di sequenze che abbiano una particolare proprietà (es.: N basi a monte e a valle di un sito di splicing). Le sequenze possono essere descritte con una matrice che ha 2xN+1 colonne (numero di posizioni) e 4 righe (una per ogni base). Ogni elemento rappresenta la frequenza con cui ogni base è osservata nella data posizione (matrice sito-specifica). Se il numero di esempi è sufficientemente alto, i valori di frequenza approssimano quelli di probabilità.
Una volta definita la matrice sito-specifica, la ricerca, ad esempio, di un sito di splicing avviene facendo slittare la sequenza in cui si cerca rispetto alla matrice e moltiplicando tra loro le probabilità che le basi della sequenza siano presenti nella rispettiva posizione. Un valore alto indicherà un’alta probabilità che la sequenza faccia parte dell’insieme delle sequenze intorno alla regione di splicing.
Prendiamo ad esempio la sequenza GGTCACAACGTTAGG:
- Considerando la prima G come corrispondente alla prima posizione della tabella: 0.05×0.15×0.05×0.30×0.65×0.40×0.12×0.20×0.10 = 7.02×10-8
- Considerando la seconda G come corrispondente alla prima posizione della tabella: 0.05×0.15×0.80×0.45×0.25×0.32×0.12×0.30×0.20 = 1.56×10-6
- Considerando invece la sequenza AACCCACTA si ha una probabilità di avere il motivo: 0.30×0.50×0.80×0.30×0.25×0.32×0.20×0.44×0.65 = 1.65×10-4
Ma come si fa a sapere se un valore di probabilità è sufficientemente alto? Possiamo confrontarlo con il valore che ci aspetteremmo per una sequenza casuale assumendo che la sua composizione attesa sia 25% di A, 25% di T, 25% di C e 25% di G (0.259 = 3.81×10-6) mi dà la frequenza attesa. Solo la sequenza AACCCACTA avrebbe una probabilità significativa di appartenere alla regione descritta dalla matrice.
Poiché le probabilità sono minori di 1, il loro prodotto tende a generare numeri molto piccoli, soprattutto per N elevati. Un modo per ovviare a tale problema è utilizzare il log-odd, ovvero calcolare il rapporto tra frequenza osservata e frequenza attesa per ciascuna posizione e poi estrarne il logaritmo. Numeri minori di 1 generano logaritmi negativi e viceversa. Inoltre, poiché il logaritmo di un prodotto è la somma dei logaritmi, la probabilità finale è generata dalla somma dei log-odd.
Prendiamo ad esempio la sequenza GGTCACAACGTTAGG:
- Considerando la prima G come corrispondente alla prima posizione della tabella: -2.32-0.74-2.32+0.26+1.38+0.68-1.06-0.32-1.32 = -5.76
- Considerando la seconda G come corrispondente alla prima posizione della tabella: -2.32-0.74+1.68+0.85+0.00+0.36-1.06+0.26-0.32 = -1.29
Entrambe generano score negativi perché nella maggior parte delle loro posizioni compaiono basi poco probabili nella regione descritta dalla matrice.
- Al contrario, AACCCACTA avrà punteggio: 0.26+1.00+1.68+0.26+0.00+0.36-0.32+0.82+1.38 = 5.44
Attenzione al numero 0 (frequenza nulla di una base in una data posizione) in quanto log20 = -∞!
Nelle matrici sito-specifiche si assume che la probabilità di una base in una data posizione sia indipendente dal tipo di base che la precede; tuttavia, ad esempio, un sito di splicing è identificato dalla presenza del dinucleotide GA, quindi la probabilità che una A sia all’inizio di un introne dipende da quale base la precede. I processi statistici che descrivono il caso in cui un evento dipende da ciò che lo precede si chiamano processi di Markov. Questi metodi predittivi sono in grado di cogliere in maniera abbastanza efficiente i motivi nel genoma.
Allineamento di sequenze di acidi nucleici e proteine
L’allineamento tra due sequenze consiste nella determinazione di una relazione tra i residui della prima sequenza con quelli della seconda in modo da rendere massimo il grado di similarità o analogamente rendere minimo il numero di differenze. Esso stabilisce una relazione biunivoca tra due sequenze (o parti di esse) in modo da minimizzare il numero di operazioni necessarie per...