Universita degli Studi di RomaROMA TRE

FACOLTA DI INGEGNERIA

Corso di Laurea Magistrale in

INGEGNERIA GESTIONALE E DELL’AUTOMAZIONE

INFERENZA DI APLOTIPI:

UN BRANCH AND BOUND BASATO

SU UNA FORMULAZIONE DI

SET COVERING

Tesi di Laurea in Ricerca Operativa

Candidato:Pietro Leonori

Relatrice:

ProfessoressaGaia Nicosia

Corelatrice:

DottoressaPaola Bertolazzi

Sessione di Maggio

Anno Accademico 2005/2006

Indice

Introduzione 3

1 Introduzione al problema biologico 6

1.1 La genetica classica . . . . . . . . . . . . . . . . . . . . 6

1.2 La genetica moderna . . . . . . . . . . . . . . . . . . . 9

1.3 Le leggi mendeliane applicate agli esseri umani . . . . . 11

1.4 La scoperta del DNA . . . . . . . . . . . . . . . . . . . 12

1.5 La differenziazione genetica . . . . . . . . . . . . . . . 15

1.6 Introduzione alla biologia computazionale . . . . . . . . 16

2 Bioinformatica e problemi di analisi degli aplotipi 19

2.1 Cromosomi, genotipi, aplotipi e SNP . . . . . . . . . . 19

2.2 Il problema di inferenza di aplotipi . . . . . . . . . . . 20

2.3 Il modello biologico: la parsimonia pura e la filogenia

perfetta . . . . . . . . . . . . . . . . . . . . . . . . . . 22

2.4 Identificazione dei tag SNP . . . . . . . . . . . . . . . . 23

2.5 Il Progetto di Ricerca HapMap ed il Progetto Genoma 24

3 Metodi per la soluzione del problema dell’inferenza di

aplotipi 26

3.1 I primi approcci al problema dell’Inferenza di Aplotipi 27

3.2 Gusfield e la formulazione di PLI . . . . . . . . . . . . 29

3.2.1 Dimensioni della formulazione TIP . . . . . . . 31

3.3 Algoritmi statistici . . . . . . . . . . . . . . . . . . . . 33

3.4 Algoritmi polinomiali . . . . . . . . . . . . . . . . . . . 35

1

4 La formulazione di Set Covering puro 37

4.1 Il passo di Fourier e la formulazione di Covering puro . 37

4.2 L’algoritmo di Branch and Cut . . . . . . . . . . . . . 40

4.3 L’algoritmo HI PARS . . . . . . . . . . . . . . . . . . . 43

4.4 Il preprocessamento . . . . . . . . . . . . . . . . . . . . 45

5 Dettagli implementativi e risultati computazionali 48

5.1 Dettagli implementativi . . . . . . . . . . . . . . . . . 49

5.2 Osservazioni sulla complessita delle istanze . . . . . . . 53

5.3 Definizione delle classi di istanze del test . . . . . . . . 55

5.4 Analisi dei risultati . . . . . . . . . . . . . . . . . . . . 58

Conclusioni 60

Ringraziamenti 62

2

Introduzione

La bioinformatica e la disciplina scientifica dedicata alla modellazio-

ne ed alla soluzione dei problemi biologici secondo modelli teorici e

tecniche informatiche. Essa costituisce un tentativo di descrivere dal

punto di vista numerico e statistico i fenomeni biologici, in particolare

problemi legati al DNA, ai cromosomi ed allo studio delle somiglianze

e delle differenze genetiche tra gli esseri umani.

Il DNA e una molecola costituita da due catene di elementi piu

piccole detti nucleotidi ed e l’elemento biologico che caratterizza uni-

vocamente un individuo. Il DNA e situato nel nucleo delle cellule di un

organismo ed e organizzato in cormosomi. I geni sono sottosequenze

delle molecole di DNA situate in posizioni note dei cromosomi e sono

gli elementi che definiscono i caratteri genetici e che costituiscono le

unita ereditarie degli organismi viventi.

Le moderne tecniche di sequenziamento permettono di ottenere

l’esatta sequenza dei nucleotidi di una moelcola di DNA e quindi di

studiare le caratteristiche genetiche di un individuo o di una popolazio-

ne. Per supportare gli studi del patrimonio genetico degli esseri umani,

sono nati negli ultimi decenni diversi progetti, seguendo il piu famoso

“Progetto Genoma Umano”, con l’intento di riconoscere e studiare le

funzioni dei circa 40.000 geni che caratterizzano la specie umana.

Nell’uomo, in quanto organismo diploide, i cromosomi sono presen-

ti in coppie omologhe dette aplotipi, una ereditata dal padre, l’altra

dalla madre. L’analisi degli apolotipi e di fondamentale importanza

per lo studio di problemi come il rigetto nel trapianto di organi, l’ana-

lisi dell’ereditarieta delle malattie genetiche e lo sviluppo delle reazioni

3

dei pazienti a particolari cure farmacologiche. I prossimi sviluppi dello

studio del genoma richiedono, quindi, che sia disponibile una mappa

degli aplotipi che permetta di seguire, all’interno di una popolazione,

la distribuzione di particolari geni e le loro linee ereditarie.

Leggere i valori degli aplotipi di un individuo e ancora molto co-

stoso e complicato, mentre e piu semplice leggere una informazione

sintetica dei due cromosomi omologhi, detta genotipo, in cui alcune

posizioni risultano ambigue. Per ovviare a questo problema all’inter-

no della bioinformatica sono stati studiati molti algoritmi, dal primo

algoritmo di Clark del ’90, per calcolare gli aplotipi di un individuo a

partire dal suo genotipo.

Il problema di calcolare gli aplotipi di una popolazione data e detto

di “Inferenza di Aplotipi”. Tale problema viene solitamente trattato

secondo il modello della parsimonia pura, secondo cui la soluzione mi-

gliore e quella con il minor numero possibile di aplotipi: se l’obbiettivo

dell’Inferenza di Aplotipi e, dato un insieme di genotipi che rappresen-

ta una popolazione, trovare un insieme di aplotipi che spieghi tutti i

genotipi dati, il modello biologico della parsimonia pura dice che la so-

luzione ottima del problema e quella per cui la cardinalita dell’insieme

di aplotipi trovati e minima.

Obiettivo di questa tesi e di studiare un modello di Programma-

zione Lineare Intera per il problema di Inferenza di Aplotipi con par-

simonia pura, che nasce dall’applicazione della procedura di Fourier-

Motzkin al modello di Gusfield. Le variabili del modello presentato

sono molto ridotte rispetto alle variabili del modello di Gusfield e i

suoi vincoli, seppur in numero maggiore, sono facilmente identificabi-

li. L’algoritmo implementato per testare le qualita della nuova mo-

dellazione si basa su una procedura di ricerca dei tagli, procedura che

permette, nel caso medio, di ridurre notevolmente il numero dei vincoli

considerati per ottenere la soluzione ottima.

Nel primo capitolo vengono rivisitate le principali tappe stori-

che dello studio della genetica, dai primi esperimenti di Mendel alle

applicazioni della biologia computazionale.

4

Il secondo capitolo introduce ai problemi di bioinformatica legati

allo studio dei genotipi, degli aplotipi e dei cromosomi.

Il terzo capitolo ripercorre le tappe principali dello studio del pro-

blema di Inferenza di Aplotipi, dai primi approcci basati sulla regola

di Clark, al modello di Gusfield, fino ai modelli polinomiali e statistici.

Nel quarto capitolo viene introdotto il modello di set covering puro

presentato in questa tesi per la soluzione del problema di Inferenza

di Aplotipi. Per la soluzione di tale modello viene implementato un

algoritmo di Branch and Bound.

Nel quinto capitolo l’algoritmo implementato per la soluzione del

problema viene testato su insiemi di dati simulati e vengono analizzati

i risultati ottenuti.

5

Capitolo 1

Introduzione al problema

biologico

In questo capitolo vengono riassunti i concetti principali di biologia

da cui nasce il problema che si vuole risolvere in questo lavoro di tesi

secondo la loro recente evoluzione storica. Dai primi studi genetici

di Mendel e dagli studi della genetica moderna sulla locazione fisica

dei geni e sulle loro funzioni, alla scoperta e del DNA e della sua

importanza nell’analisi della differenziazione genetica nell’uomo.

1.1 La genetica classica

Fin dall’inizio della storia dell’uomo, ci si e domandati come i caratteri

venissero trasmessi da una generazione alla successiva. Nonostante

che i figli assomiglino spesso piu ad un genitore piuttosto che all’altro,

la maggior parte della progenie sembra essere una mescolanza delle

caratteristiche di entrambi i genitori.

Nel 1865 un monaco agostiniano di nome Gregor Mendel scoprı che

i caratteri individuali sono determinati da fattori discreti, piu avan-

ti denominati geni, che vengono ereditati dai genitori. Mendel inizio

partendo da genitori con un noto passato genetico, in modo da avere

una base per paragonare gli schemi di ereditarieta nella progenie otte-

6

Figura 1.1: Caratteri della pianta di pisello dolce

nuta. Successivamente analizzo come le caratteristiche dei genitori si

tramandavano alle generazioni successive.

Gli studi di Mendel si basavano sull’analisi delle piante di pisello

dolce, per le quali e facile distinguere alcune caratteristiche traman-

date geneticamente e che hanno un sistema riproduttivo che facilita

gli studi evolutivi. Piuttosto che considerarle come un’unica entita,

Mendel focalizzo la sua attenzione su sette singoli caratteri che poteva

prontamente distinguere. Scoprı che ciascun carattere presentava due

forme alternative, per esempio, il colore dei semi poteva essere verde

o giallo. Analizzando i risultati di diversi incroci, Mendel concluse

che ogni forma alternativa di un carattere e determinata da forme

alternative di un gene.

Per seguire le tracce dell’ereditarieta dei geni da genitori a figli

Mendel doveva innanzitutto essere sicuro di quali geni fossero forniti

da ciascun genitore. Dal momento che le piante di pisello sono natu-

ralmente auto fecondanti, le specie di linea pura erano prontamente

disponibili. Ogni varieta conteneva solo una delle possibili forme del

gene, che furono chiamate successivamente alleli, che determinava un

carattere. Linee pure di piante con semi gialli producevano solo pro-

genie con i semi gialli (con l’allele “semi gialli”), mentre linee pure

di piante con semi verdi producevano unicamente progenie con semi

verdi. Dai risultati dei successivi esperimenti, Mendel dedusse che le

7

linee pure delle piante dovevano avere due copie dello stesso gene per

ciascun carattere.

In generale, la progenie sembra essere una mescolanza delle carat-

teristiche dei genitori, ma non era vero per i caratteri delle piante che

Mendel aveva scelto di studiare: quando venivano incrociate le linee

pure, non venivano prodotte progenie con caratteri mescolati. Per

esempio, ci si potrebbe aspettare che un incrocio tra una linea pura di

piante dai semi verdi e una linea pura dai semi gialli avrebbe prodotto

una progenie con semi di un colore intermedio tra il verde e il giallo,

mentre questo incrocio produceva progenie di un solo colore, il giallo.

Non si riscontravano ne sfumature intermedie di colore, ne tracce di

verde.

Dall’analisi della progenie ottenuta incrociando piante della gene-

razione ibrida fu chiaro che per ogni caratteristica osservata si potesse-

ro distinguere un carattere dominante, che si manifestava nella prima

generazione ibrida, e un carattere recessivo. Il fatto che il carattere

recessivo in questa progenie fosse riscontrabile statisticamente su una

pianta su quattro era un chiaro segno di come venisse tramandato il

patrimonio genetico: Mendel ipotizzo per primo che in ogni indivi-

duo sono presenti due copie di ogni gene, ma che le cellule sessuali

(spermatozoi e ovuli) ne contengono solo una e che l’intero patrimo-

nio genetico viene ricostituito al momento della fecondazione. Queste

sue ipotesi furono il primo passo verso il riconoscimento del fatto che

nel patrimonio genetico i geni sono raggruppati in cromosomi e che

negli organismi diploidi, come l’uomo e le piante di pisello dolce, sono

presenti due copie di ogni cromosoma.

Mendel pubblico le sue ricerche, Experiments in Plant Hybridisa-

tion, nel 1865 e ne invio copie ad eminenti scienziati di numerosi paesi.

Tuttavia le sue conoscenze astratte dei geni non furono apprezzate dai

naturalisti del tempo, a cui era stato insegnato principalmente ad os-

servare e catalogare gli esseri viventi. Il lavoro di Mendel fu trascurato

fino al 1900, quando diversi scienziati europei confermarono separata-

mente i suoi risultati. Poco prima della riscoperta del lavoro di Mendel,

8

furono portati avanti attenti studi del comportamento dei cromosomi

durante la formazione delle cellule sessuali (meiosi).

La meiosi e il meccanismo di divisione mediante il quale una cellula

eucariota con corredo cromosomico diploide (in cui sono presenti due

copie di ogni cromosoma) da origine a cellule con corredo cromosomico

aploide (in cui e presente una sola copia di ogni cromosoma). La

diversita genetica viene mantenuta dalla riproduzione sessuata, che

comprende cicli cellulari con ricombinazione dell’informazione genetica

proveniente da cellule di due organismi differenti (padre e madre).

Nella meiosi la ricombinazione del materiale genetico provenien-

te dai due genitori avviene durante il processo di crossing-over (o

semplicemente ricombinazione), durante il quale coppie di cromosomi

scambiano porzioni omologhe di materiale genetico ([15]).

1.2 La genetica moderna

Thomas Hunt Morgan i suoi studenti della Columbia University apri-

rono l’era della genetica moderna dimostrando le basi fisiche dell’e-

reditarieta, grazie a studi che iniziarono nei primi anni del ventesimo

secolo e che culminarono nel 1933 con l’assegnazione del Premio Nobel

per la medicina, per essere stati tra i pionieri dell’embriologia e della

morfologia speriementale.

Invece delle piante di pisello coltivate da Mendel, il team della Co-

lumbia University studio l’ereditarieta nella comunissima mosca della

frutta, per la quale, pero, fu piu complesso trovare caratteri facilmen-

te identificabili che potessero essere studiati singolarmente. Alla fine,

scoprirono un’unica mosca maschio dagli occhi bianchi, che si distin-

gueva dai suoi pari, normalmente dagli occhi rossi. Un incrocio tra il

maschio mutante e la femmina dagli occhi rossi dava origine solamente

a progenie con gli occhi rossi e i mutanti con gli occhi bianchi ricom-

parivano nella generazione successiva, secondo il classico schema di un

carattere recessivo. Tuttavia, il carattere degli occhi bianchi veniva

9

riscontrato esclusivamente nei maschi della seconda generazione. Cosı

giunsero alla conclusione che gli occhi bianchi sono un carattere reces-

sivo legato al sesso e che il gene del colore degli occhi doveva essere fisi-

camente collocato nel cromosoma X, cromosoma che contraddistingue

gli individui di sesso maschile.

Man mano che Morgan e i suoi collaboratori identificavano nelle

mosche della frutta sempre piu caratteri ereditari, notarono che spesso

le mosche presentavano particolari combinazioni di caratteri. Questo li

indusse a pensare che determinati geni fossero connessi e che venissero

ereditati insieme come un unico blocco. Essi identificarono quattro di

questi blocchi, o “gruppi di connessione”, cioe nello stesso numero dei

cromosomi appaiati osservati al microscopio.

Il gruppo di Morgan sfrutto il fenomeno della concatenazione per

costruire una mappa dei cromosomi delle mosche della frutta. Scopri-

rono che geni connessi a volte si separano nella meiosi, nel processo

di ricombinazione e ne dedussero che la frequenza con cui due geni si

separano fornisce una misura della distanza relativa tra di essi su un

cromosoma: geni distanti si ricombinano frequentemente, geni vicini

si ricombinano raramente e sono strettamente legati.

Oltre alla posizione fisica dei geni ed alla loro trasmissione e ri-

combinazione nei processi riproduttivi, l’attenzione in quegli anni era

concentrata anche sul fatto che la chiave della teoria dell’evoluzio-

ne, attraverso la selezione naturale, e la variazione dei nuovi caratteri

che insorgono spontaneamente e rendono un organismo piu competi-

tivo nella lotta per la sopravvivenza. Queste variazioni che insorgono

mutando il patrimonio genetico di un individuo con meccanismi indi-

pendenti dai processi ereditari sono comunemente chiamate mutazioni.

Guidati dalla pubblicazione de “L’origine delle specie” di Darwin, nel

1859 erano stati allestiti dei laboratori sul campo dove gli scienziati

potevano studiare le caratteristiche uniche degli organismi che si erano

evoluti per occupare differenti ecosistemi. Comunque le osservazioni

sul campo non avrebbero potuto spiegare l’origine delle variazioni o

come nuovi caratteri vengono ereditati.

10

La nuova sotto-disciplina dell’evoluzione sperimentale si sviluppo

a cavallo del ventesimo secolo con lo scopo di ricreare l’evoluzione in

esperimenti controllati con piante coltivate e animali. Divenne in breve

chiaro che le mutazioni genetiche sono la sorgente delle variazioni e

che la genetica mendeliana offriva un metodo statistico per analizzare

l’ereditarieta delle nuove mutazioni. All’inizio del 1920 gli evoluzionisti

sperimentali erano diventati la prima generazione di genetisti.

1.3 Le leggi mendeliane applicate agli es-

seri umani

Nonostante che le leggi di Mendel fossero testate inizialmente sulle

piante di pisello e sulle mosche della frutta, fu da subito evidente che

potevano essere applicate a tutti gli esseri viventi. Come le mutazioni

avevano fornito le chiavi per capire la genetica delle mosche della frut-

ta, cosı gli alberi genealogici di famiglie affette da malattie offrirono i

primi casi di studio dell’ereditarieta mendeliana negli uomini.

L’ereditarieta di tipo recessivo fu descritta per prima nelle malattie

come l’alcaptonuria (1902) e l’albinismo (1903). Fra le prime malat-

tie di tipo dominante scoperte ci furono la brachidattilia (1905), la

cataratta congenita (1906) e il morbo di Huntington (1913). La di-

strofia muscolare di Duchenne (1913), il daltonismo (1914) e l’emofilia

(1916) furono le prime malattie per le quali si scprı che erano legate

al sesso. Il semplice concetto dell’ereditarieta del colore degli occhi -

marrone dominante, blu recessivo - fu pubblicato nel 1907 (tuttavia

attualmente gli scienziati credono che vi siano implicati vari geni).

Ci fu un ovvio interesse verso l’applicazione delle leggi di Men-

del all’agricoltura e le sue idee furono anche abbracciate dal movi-

mento dell’eugenetica, che si poneva l’obiettivo di migliorare la specie

umana attraverso unioni controllate. Gli eugenetisti incoraggiavano

i matrimoni tra persone dal buon bagaglio genetico e scoraggiavano

la riproduzione nelle coppie geneticamente inadatte, ma adoperarono

11

erroneamente i semplici schemi dominante/recessivo per spiegare com-

portamenti complessi e disturbi mentali, per i quali adesso sappiamo

essere coinvolti molti geni, e fallirono anche nel giustificare gli effetti

ambientali sullo sviluppo umano. Negli Stati Uniti restrittive legisla-

zioni eugenetiche rifletterono pregiudizi politici e sociali, piuttosto che

risultati genetici. La descrizione eugenetica della vita umana fu in

seguito screditata dalle orribili conseguenze della ricerca nazista della

purezza della razza.

1.4 La scoperta del DNA

Il DNA, acido deossiribonucleico, e stato riconosciuto come uno dei

principali componenti chimici del nucleo negli anni in cui venivano

pubblicati i lavori di Mendel e Darwin. Tuttavia ancora nei primi

anni del 1900, le proteine erano considerate le uniche molecole capaci

di trasportare grandi quantita di informazione ereditaria di genera-

zione in generazione. Sebbene fosse gia noto che il DNA e molecola

molto lunga, sembrava che le sue quattro componenti chimiche fossero

assemblate secondo un modello fisso, come un polimero sintetico, e

ancora non erano state trovate le sue specifiche funzioni cellulari. Le

proteine, d’altra parte, erano gia state riconosciute come importanti

enzimi e componenti strutturali delle cellule viventi. Si sapeva inoltre

che le proteine sono polimeri di numerosi aminoacidi e si credeva che

l’alfabeto dei 20 aminoacidi delle proteine potesse spiegare meglio le

caratterizzazioni univoche degli individui.

Nel 1902 Archibald Garrod descrisse la malattia ereditaria dell’al-

captonuria come un errore congenito del metabolismo. Egli propose

che una mutazione ad un gene potesse causare un difetto specifico

nella via biochimica per l’eliminazione dei liquidi residui. Il fenoti-

po della malattia - urine scure - e una conseguenza di questo difetto.

Questa ipotesi venne rigorosamente dimostrata nel 1941 da George

Beadle e da Edward Tatum, usando come organismo la Neurospora,

12

la semplice muffa del pane. Inizialmente scoprirono che le muffe espo-

ste a radiazione perdono la capacita di produrre nutrimenti essenziali,

il che rallenta o addirittura blocca la loro crescita, e successivamen-

te scoprirono che la crescita poteva essere ripristinata fornendo alla

muffa mutata degli additivi. Considerarono che una mutazione de-

ve essere responsabile della disattivazione dell’enzima (proteina) ne-

cessario per sintetizzare il nutrimento e che quindi un gene doveva

essere responsabile dell’informazione necessaria alla formazione della

proteina.

Negli anni venti fu dimostrato sperimentalmente che un ceppo in-

nocuo di batteri poteva diventare infettivo quando veniva mescolato

ad un ceppo virulento che era stato ucciso. Il batterio morto appa-

rentemente attivava dei processi chimici che trasformavano il batterio

innocuo in infettivo. Questo principio, detto principio trasformante,

sembrava essere legato ad un gene. Una equipe di scienziati guidata da

Oswald Avery al Rockfeller Institute, seguı rigorosamente quest’espe-

rimento negli anni ’40. Trovarono che un estratto puro del principio

trasformante non veniva colpito dal trattamento con enzimi che dige-

riscono proteine, ma era eliminato da un enzima che digeriva DNA.

Eppure molti scienziati furono lenti nell’accettare che la molecola che

trasporta il patrimonio genetico e il DNA e non le proteine.

L’avvento dei microscopi elettronici permise di osservare come agi-

sce un virus batterico: il virus attacca un batterio ospite iniettandogli

i propri geni attraverso la sua coda. Nel 1952 Alfred Hershey mostro

che il DNA da solo e responsabile della riproduzione di nuovi virus

all’interno di una cellula infettata. Questo forniva un’inconfutabile

conferma ai precedenti esperimenti di Avery secondo i quali un gene

e fatto di DNA e faceva emergere che i virus, cosı come i batteri, pos-

sono essere usati come modelli per studiare principi universali della

genetica.

Il modello che viene attualmente riconosciuto come primo modello

accurato della struttura del DNA fu proposto da James Watson e

Francis Crick. Nel loro articolo “La struttura molecolare degli acidi

13

nucleici: una struttura per l’acido deossiribonucleico”, pubblicato nel

1953 sul giornale scientifico Nature e basato sulle immagini ottenuti

da Rosalind Franklin con la tecnica della diffrazione a raggi X, per la

prima volta viene descritta la struttura a doppia elica del DNA. Dal

loro lavoro fu chiaro che il DNA e una lunga catena di monomeri, detti

nucleotidi, costituiti da un gruppo fosfato, da uno zucchero pentosio

(desossiribosio) e da una base azotata. I nomi dei nucleotidi del DNA

sono solitamente abbreviati da una lettera che indica da quale base

azotata sono costituiti: A per adenina, C per citosina, G per guanina

e T per timina. Nel DNA la quantita di guanina e uguale a quella

della citosina e la quantita di adenina e uguale a quella di timina,

questo indusse a pensare che nella catena di DNA le basi azotate

fossero sempre in coppie A:T e C:G e Watson e Crick proposero nel

loro modello che tali basi costituissero i ‘pioli’ della doppia elica e che

questa struttura ad informazione duplicata servisse per facilitare la

duplicazione delle catene di DNA.

A seguito di queste ultime due scoperte venne enunciato nel 1958,

per poi essere successivamente riformulato e pubblicato in una pubbli-

cazione su Nature nel 1970, il dogma centrale della biologia molecola-

re, secondo cui il flusso dell’informazione genetica e monodirezionale

e parte dal DNA per arrivare alle proteine. In questo processo pos-

sono essere identificati tre fasi principali: l’informazione genetica vie-

ne conservata nel DNA, trascritta nell’RNA (acido ribonucleico) che

viene successivamente tradotto in proteine, la forma operativa dell’in-

formazione contenuta nel genoma. Il flusso di informazioni e canaliz-

zato attraverso un codice per cui tre nucleotidi (un codone) codifica-

no un aminoacido. Gli amminoacidi presenti negli organismi viventi

sono numerosissimi ma solo venti di essi interessano nel processo di

trasmissione del codice genetico.

Se Mendel aveva descritto il gene come una porzione definita di

eredita che influenza caratteri visibili, studi successivi all’enunciato

del dogma centrale portarono ad una riformulazione del concetto di

gene come sottosequenza discreta di DNA che codifica una proteina

14

Figura 1.2: La doppia elica del DNA

e che comincia con un codone di partenza e finisce con un codone

di chiusura. L’analisi dei geni fece un grande passo avanti con la

scoperta dei metodi per determinare l’esatta sequenza di nucleotidi

che compone uno specifico gene, detti metodi di sequenziamento del

DNA.

1.5 La differenziazione genetica

Analizzando le sequenze di DNA di diversi individui appartenenti alla

stessa specie si puo osservare che non e possibile trovare due sequenze

identiche, ma che comunque le differenze che si possono riscontrare

sono circa di un nucleotide ogni 1000.

Tali differenze sono dovute alle mutazioni che da sempre guidano

e caratterizzano l’evoluzione delle specie e possono influire positiva-

mente o negativamente sull’evoluzione della famiglia di un individuo,

ma anche essere ininfluenti. Le mutazioni possono essere dovute all’e-

sposizione a raggi X, ad altre radiazioni ionizzanti (come per esempio

i raggi ultravioletti, componenti della luce solare), o possono esse-

re dovute a sostanze mutagene o anche ad alcune fasi imperfette del

processo stesso di replicazione del DNA.

Fino agli anni ’50 la maggior parte dei biologi pensavano che i geni

fossero unita stabili. La nozione che il DNA potesse essere danneggia-

to e poi riparato fu introdotta da ricercatori che stavano cercando di

spiegare il comportamento anomalo di alcuni microbi che stavano stu-

diando: colture che erano evidentemente state uccise dall’esposizione

15

ai raggi ultravioletti, si riprendevano dopo essere state messe vicino

ad una finestra e dei mutanti facevano curiosamente la loro comparsa

molto tempo dopo l’esposizione ad agenti mutageni. Studi condotti

su organismi di natura diversa, dai batteri agli esseri umani, hanno

rivelato un numero impressionante di enzimi in grado di riparare dan-

ni causati da agenti mutageni ambientali o da errori nella replicazione

del DNA. Senza questi enzimi, il danno al DNA causerebbe livelli di

mutazione intollerabili.

Da un lato, il fatto che disturbi causati da enzimi di riparazione

difettosi accorciano la durata della vita mette in luce il ruolo centra-

le della riparazione del DNA per la sopravvivenza di un organismo.

Dall’altro, il fatto che il DNA non sia sempre riparato o corretto e pu-

re centrale per la sopravvivenza della specie, perche queste anomalie

hanno come risultato delle mutazioni e senza mutazioni l’evoluzione

sarebbe impossibile.

Gli studi sul patrimonio genetico e sulle mutazioni suggeriscono

l’idea che se si conoscessero le funzioni e la posizione di tutti i geni

nel DNA si potrebbe cercare di indurre delle mutazioni puntuali per

prevenire o curare le malattie genetiche o si potrebbero studiare le

reazioni di un individuo ad una particolare cura farmacologica oppure,

piu semplicemente, si potrebbero approfondire le nostre conoscenze sul

funzionamento del corpo umano.

1.6 Introduzione alla biologia computa-

zionale

Il codice genetico e l’elemento biologico che distingue univocamente

un individuo da ogni altro essere della sua specie e da ogni altra forma

di vita. Nonostante cio individui appartenenti alla stessa specie ne

condividono in genere circa il 99,9%: le differenze che rendono unico

l’individuo risiedono in variazioni minime del DNA.

Per la ricerca delle basi biologiche di tali differenze e delle loro

16

implicazioni nel funzionamento del corpo umano sono nati negli anni

novanta diversi progetti che si ponevano l’obiettivo di supportare lo

studio della diagnostica e della cura delle malattie genetiche.

La dimensione di questi progetti e delle loro aspirazioni sono indice

del sempre maggior interesse che viene riposto nella Biologia Compu-

tazionale che, secondo il National Institute of Healt degli Stati Uniti,

e la “scienza che studia lo sviluppo e l’applicazione di metodi teorici di

analisi dei dati, modelli matematici e tecniche di simulazione compu-

tazionale per lo studio di sistemi biologici, comportamentali e sociali”.

Essa e anche definita come campo interdisciplinare che applica tecni-

che informatiche, di matematica applicata e di statistica per creare

modelli finalizzati alla soluzione di problemi ispirati dalla biologia. I

campi principali di applicazione della biologia computazionale sono:

• la bioinformatica, che applica algoritmi e tecniche statistiche

ad insiemi di dati provenienti da analisi biologiche, dati che

tipicamente consistono in sequenze di DNA, RNA o proteine.

• la genomica computazionale, che utilizza tecnologie di microarray

sul DNA per analizzare i geni espressi nelle tipologie di cellula

umana

• la modellazione molecolare, che studia modelli teorici e tecniche

computazionali per modellare o simulare il comportamento di

molecole di diverse dimensioni, da molecole di pochi atomi alle

grandi molecole presenti nei sistemi biologici

• la biologia dei sistemi, che vuole modellare reti di interazioni

biologiche su grande scala, spesso con il supporto delle equazioni

differenziali

• la predizione della struttura delle proteine e la genomica strut-

turale, che vogliono produrre modelli strutturali accurati per

strutture proteiniche tridimensionali che non sono state risolte

sperimentalmente

17

• la biochimica computazinale e la biofisica computazionale, che

fanno ampio uso della modellazione strutturale e di tecniche di

simulazione per mettere in luce le funzioni cinetiche e termodi-

namiche delle proteine.

L’inferenza di aplotipi, oggetto di questo lavoro di tesi, e un proble-

ma di biologia computazionale che vuole approfondire le conoscenze

sul genoma umano per supportare progetti futuri sulla diagnostica del-

le malattie genetiche e sullo studio dell’evoluzione e dell’ereditarieta

dei geni.

18

Capitolo 2

Bioinformatica e problemi di

analisi degli aplotipi

In questo capitolo vengono introdotti alcuni problemi relativi alla gene-

tica ed alla biologia computazionale, tra cui il problema dell’inferenza

di aplotipi, la cui soluzione e oggetto principale di questo lavoro di

tesi.

2.1 Cromosomi, genotipi, aplotipi e SNP

I cromosomi sono corpuscoli presenti nei nuclei delle cellule eucariote

costituiti da un filamento di DNA e da proteine e il loro numero va-

ria di specie in specie. Quasi tutte le cellule delle specie viventi sono

diploidi e quindi caratterizzate dalla presenza di coppie omologhe di

cromosomi, 23 coppie per la specie umana, di cui una copia e ereditata

dal padre, l’altra dalla madre. Nell’uomo le uniche cellule aploidi, con

una sola copia di ognuno dei 23 cromosomi, sono le cellule germina-

li, che hanno il compito di trasmettere i caratteri genetici ai figli e

che ricompongono l’intero patrimonio genetico del nuovo individuo al

momento della fecondazione.

Il genotipo e il profilo genetico di un individuo, ovvero la totalia dei

geni presenti nel suo genoma, e contribuisce con i fattori ambientali

19

a determinare il suo fenotipo, l’insieme delle caratteristiche rilevabili

dall’esterno. La combinazioni di varianti alleliche lungo una singola

copia di un cromosoma e, invece, chiamata aplotipo e, per semplicita,

possiamo considerarla come la parte del patrimonio genetico ereditata

da uno solo dei due genitori. In questo caso per allele si intende il

valore che ha un nucleotide in una particolare sequenza.

Uno SNP (Single Nucleotide Polymorphism, polimorfismo di un

singolo nucleotide) e un singolo sito del DNA dove esattamente due

(su quattro) differenti nucleotidi vengono riscontrati nella grande mag-

gioranza della popolazione e rappresenta la piu piccola unita di dif-

ferenziazione genetica misurabile tra due individui. Un aplotipo puo

rappresentare le sequenza completa del DNA, oppure, come avviene

piu comunemente, i valori degli SNP in una data regione del DNA.

Il problema della ricerca della sequenza completa degli SNP per

una delle due copie di un dato cromosoma (aplotipo) in un genoma

diploide e detta inferenza di aplotipi o phasing e i due aplotipi di un

individuo sono detti fase del relativo genotipo.

2.2 Il problema di inferenza di aplotipi

Nella biochimica sono noti diversi metodi di sequenziamento di DNA

che permettono di ottenere il genotipo di un individuo, mentre otte-

nere i due aplotipi e ancora molto costoso. Dato che nei problemi di

diagnosi medica, filogenia e sviluppo di nuovi farmaci e importante

avere non tanto il genotipo quanto la sua fase, dagli inizi degli anni

20

’90 diversi scienziati, nell’ambito della disciplina scientifica della bio-

informatica, affrontano problemi legati allo studio della fase del profilo

genetico non con metodi biochimici, ma con strumenti computazionali.

Il problema di inferenza per un singolo individuo (SIH, Single In-

dividual Haplotyping) puo essere enunciato nel seguente modo:

Problema di inferenza per un singolo individuo: da-

to un insieme di frammenti del DNA di un individuo, tro-

vare una coppia di aplotipi che sia il piu possibile consi-

stente con tali frammenti rimuovendo gli errori legati al

sequenziamento, alle ripetizioni e ad altre cause.

Questo problema e strettamente legato alle tecniche di sequenzia-

mento attualmente in uso, che richiedono la ricostruzione della sequen-

za di DNA da un insieme di piccoli frammenti, e risente del fatto che

le tecniche per ottenere tali frammenti spesso inducono errori e mu-

tazioni che alterano le catene di nucleotidi. Il problema SIH puo non

avere soluzione anche quando non ci siano errori se l’informazione e

insufficiente. Per esempio, se un insieme di frammenti non condivi-

de SNP con gli altri frammenti possiamo solo ricostruire parzialmente

gli aplotipi, ma non possiamo sapere come unirli in soli due aplotipi.

Rizzi et al. in [28] dimostrano che il problema e polinomiale quan-

do ogni frammento copre un insieme di SNP consecutivi (frammento

senza interruzioni) mentre in generale e NP-difficile.

L’approccio piu comune al calcolo della fase degli aplotipi e basato

sullo studio di una popolazione, piu che di un singolo. Il problema

dell’inferenza di aplotipi (HI, Haplotype Inference), che verra trattato

piu approfonditamente nel prossimo capitolo, vuole calcolare gli aplo-

tipi da un insieme di individui, o popolazione, rappresentato come

insieme di genotipi e puo essere enunciato nel seguente modo:

Inferenza di Aplotipi Dato un un insieme G = {g1, ..., gm}di m genotipi di n SNP, vogliamo trovare 2 ·m aplotipi che

risolvano i genotipi dati.

21

La soluzione di questo problema, da questa semplice formulazione,

puo essere trovata sia prendendo casulmente due aplotipi per ogni ge-

notipo, sia andando a studiare singolarmente i 2n possibili aplotipi per

trovare quali sia giusto scegliere. Per questo e importante completare

il problema con un modello biologico.

2.3 Il modello biologico: la parsimonia

pura e la filogenia perfetta

Il modello biologico piu comunemente usato e la parsimonia pura (PP),

secondo cui la soluzione migliore e quella con cardinalita minima. La

parsimonia pura segue il principio del Rasoio di Ockham, secondo

cui la soluzione piu semplice in natura e solitamente la piu giusta,

avallato nel nostro caso dal fatto che in natura il numero degli aplotipi

riscontrati in una popolazione e molto inferiore al numero dei possibili

aplotipi.

Con questo modello biologico, il problema puo essere cosı enuncia-

to:

Inferenza di Aplotipi con Parsimonia Pura (HIPP):

dato un insieme G = {g1, ..., gm} di m genotipi di n SNP,

vogliamo trovare un insieme di aplotipi H tale che ogni

genotipo in G sia spiegato almeno da una coppia di aplotipi

in H e vogliamo che H sia di cardinalita minima.

Nel prossimo capitolo verranno studiati vari modelli per la solu-

zione del problema di HIPP, in particolare il modello HIc sviluppato

a partire dal modello di Gusfield ed implementato per la prima volta

per questo lavoro di tesi.

La teoria coalescente vuole invece definire un modello biologico

seguendo l’evoluzione della catena di DNA. Una coalescenza e un pro-

cesso che fornisce una storia di come un insieme di aplotipi si e evoluta

in una popolazione. Supponendo che tutti i geni o gli alleli di una data

22

popolazione sono ereditati da un unico antenato, noto come l’antenato

comune piu recente, e dato il modello dei siti infiniti, secondo cui le

mutazioni avvengono sempre su siti distinti, possiamo rappresentare

la storia evolutiva degli aplotipi come un albero filogenetico. Senza

andare oltre gli interessi del problema di HI, un albero filogenetico e

un albero che ha per nodi degli aplotipi, teoricamente a diversi livelli

evolutivi, e per archi delle mutazioni, dove archi diversi rappresentano

mutazioni su loci diversi.

Il problema dell’Inferenza di Aplotipi con Filogenia Perfetta vuole

che una soluzione ottenuta per il problema di HI possa essere disposta

su di un albero filogenetico, il che implica che non sia prevista ricombi-

nazione nell’evoluzione degli individui. Questo modello biologico viene

affrontato gia in [2] da Gusfield e Orzack.

2.4 Identificazione dei tag SNP

Per ovviare ai problemi dovuti alla grande quantita di dati che e ne-

cessario processare nel risolvere problemi legati al DNA, molti scien-

ziati si sono dedicati alla ricerca di sottoinsiemi minimi di SNP che

permettano comunque di rappresentare le associazioni tra malattie e

geni. Questo processo viene chiamato identificazione dei tag SNP e in

generale gli SNP selezionati vengono indicati come tag SNP e quelli

non selezionati come tagged SNP.

Sia S = {s1, ..., sn} un insieme di n SNP in una data regione e

D = {h1, ..., hm} un insieme di m aplotipi consistenti con gli n SNP.

Supponiamo che il numero massimo di tag SNP sia k e che la funzione

f(T ′, D) valuti quanto il sottoinsieme T ′ ⊂ S rappresenta bene D,

allora il problema puo essere enunciato come:

Identificazione dei tag SNP: dato un insieme di SNP S,

un insieme di aplotipi A, un numero massimo di tag SNP

k e una funzione f(T ′, D), vogliamo calcolare un insieme

23

di tag SNP

T = argmax{T ′ t.c. T ′⊂S&|T ′|≤k}f(T ′, D).

Phil Hyoun Lee, in [22], divide gli algoritmi per la soluzione di

questo problema in quattro categorie basate sull’approccio secondo

cui misurano l’informazione degli aplotipi:

• differenza tra aplotipi

• associazione tra coppie di SNP

• predizione basata sui tag SNP

• associazione al fenotipo.

2.5 Il Progetto di Ricerca HapMap ed il

Progetto Genoma

Il progetto di biologia computazionale sicuramente piu noto e piu am-

bizioso e il Progetto Genoma Umano che nasce nell’Ottobre del 1990

come progetto di collaborazione internazionale con lo scopo di creare

una mappa del genoma umano, ovvero di descrivere la struttura, la

posizione e la funzione dei circa 40.000 geni che caratterizzano la spe-

cie umana. Lo studio del genoma implica il sequenziamento del DNA,

cioe l’identificazione della sequenza dei 3 miliardi di coppie di basi

azotate che ne compongono la molecola. Obiettivo finale del progetto

Genoma Umano e la comprensione della funzione dei geni e di quali

malattie possono derivare dalla loro alterazione.

Con ambizioni simili nasce nell’Ottobre del 2002 il progetto inter-

nazionale HapMap, per coordinare una collaborazione internazionale

tra Giappone, Gran Bretagna, Canada, Cina, Nigeria e Stati Uniti

mirata ad identificare e catalogare somiglianze e differenze nella gene-

tica degli esseri umani. A differenza del progetto genoma, il progetto

24

HaplMap mira a creare una mappa non dei genotipi, quanto degli

aplotipi delle popolazioni analizzate. Il suo obiettivo principale e di

mettere a disposizione di chiunque abbia interesse ad analizzarli dati

grazie ai quali sia possibile trovare i geni che influiscono sulla salute,

sulle malattie e sulla reazione degli individui ai medicinali ed ai fattori

ambientali.

25

Capitolo 3

Metodi per la soluzione del

problema dell’inferenza di

aplotipi

In questo capitolo vengono presentati i principali modelli ed algoritmi

per la soluzione dell’inferenza di aplotipi. Dalla regola di Clark, al

modello TIP di Gusfiel, agli algoritmi statistici e polinomiali.

Prima di entrare nel dettaglio del problema faro chiarezza sulla

terminologia che verra da qui in poi utilizzata, ricordando che il pro-

blema trattato e relativo all’essere umano, organismo di natura diploi-

de, nel cui patrimonio genetico sono quindi presenti due copie di ogni

cromosoma.

Quasi tutti gli SNP hanno solo due alleli (su quattro), di cui uno e

l’allele comune in natura, l’altro e la mutazione. Possiamo quindi rap-

presentare un aplotipo di n SNP come un vettore binario h ∈ {0, 1}n

dove h(k) = 0 se il k-esimo SNP di h presenta l’allele comune, h(k) = 1

se presenta l’allele mutato. Un genotipo (di n SNP), la cui fase sia

data dalla coppia di aplotipi (h′, h′′), puo essere rappresentato con un

26

vettore ternario g ∈ {0, 1, 2}n dove

g(k) =

0 se h′(k) = h′′(k) = 0

1 se h′(k) = h′′(k) = 1

2 se h′(k) 6= h′′(k).

Se h′(k) = h′′(k) si dice che il sito k e omozigote, altrimenti e

eterozigote.

Definiamo l’operatore di somma ⊕ : {0, 1} → {0, 1, 2} nel seguente

modo:0⊕ 0 = 0

1⊕ 1 = 1

0⊕ 1 = 2

1⊕ 0 = 2

Questo operatore viene utilizzato anche per indicare che due aplo-

tipi h′ e h′′ (entrambi di dimensione n) generano un genotipo g (an-

ch’esso di dimensione n): h′ ⊕ h′′ = g, se

g(p) = h′(p)⊕ h′′(p) ∀p ∈ {1, ..., n}.

3.1 I primi approcci al problema dell’In-

ferenza di Aplotipi

I primi approcci alla soluzione del problema dell’inferenza di aplotipi

(HI) si basavano sull’applicazione della regola di Clark. Dato un’insie-

me di genotipi ed un insieme iniziale di aplotipi, viene calcolato, con

un processo iterativo, un nuovo aplotipo ad ogni passo attraverso una

regola data, fino a che non divenga impossibile applicarla ulteriormen-

te. Clark presento il primo algoritmo basato su questo principio in

[14] nel ’90, sulla base della seguente regola:

Regola di Clark (i-esimo passo) : dato l’insieme Hi degli

aplotipi della soluzione al passo i-esimo e l’insieme Gi dei

27

genotipi non ancora risolti, determiniamo un genotipo g ∈Gi e un aplotipo h ∈ Hi per il quale esista h′ ∈ {0, 1}n tale

che h⊕ h′ = g.

Se esistono un genotipo g ed un aplotipo h che soddisfano

queste caratteristiche, allora Gi+1 = Gi\g e Hi+1 = Hi∪h′.

Questa regola viene applicata fino a che

1. Gi = ∅

oppure

2. non divenga impossibile trovare g e h che soddisfano le caratte-

ristiche richieste.

Secondo il modello biologico della massima parsimonia, l’algoritmo

di Clark ricerca per ogni passo una soluzione localmente ottima (algo-

ritmo di tipo greedy), cercando una soluzione (locale) che incrementi

il meno possibile il numero degli aplotipi, senza vincoli relativi alla

soluzione globale. La mancanza di un modello di soluzione globale im-

plica che la soluzione trovata puo essere diversa, ad ogni applicazione

dell’algoritmo, a seconda della scelte effettuate ad ogni passo se ad

una iterazione si ha piu di un genotipo candidato per essere risolto o

piu di un aplotipo compatibile con il genotipo candidato.

La seconda condizione di interruzione implica che l’algoritmo po-

trebbe finire senza aver risolto tutti i genotipi. Inoltre bisogna consi-

derare che l’insieme iniziale degli aplotipi H0 potrebbe essere vuoto.

H0 viene solitamente inizializzato secondo due semplici osservazioni:

· se un genotipo non presenta alcun 2 e compatibile con un solo

aplotipo e tale aplotipo deve stare necesariamente nella soluzione

· se un genotipo presenta un solo 2 esiste una sola coppia di

aplotipi che lo possono aver generato e tale coppia deve stare

necessariamente nella soluzione.

28

Gli aplotipi che secondo queste osservazioni sono necessariamente nel-

la soluzione vengono messi in H0, ma, se non ci sono genotipi con al

piu un 2, H0 = ∅ e l’algoritmo non parte perche non e possibile ap-

plicare la regola di Clark. Per questo in lavori successivi a quello di

Clark, come anche in [2], si cerca di definire quale sia il modello ge-

netico che guida il suo algoritmo, distinguendo tra la minimizzazione

della dimensione dell’insieme H e la massimizzazione dell’insieme dei

genotipi risolti. Gusfield in [3] suggerisce che entrambi gli obiettivi

possono essere perseguiti applicando piu volte l’algoritmo e variando

il criterio di scelta di g e di h ad ogni passo. Anche in questo modo,

pero, non si possono dare garanzie sulla qualita della soluzione finale.

3.2 Gusfield e la formulazione di PLI

Gusfield in [2] e in [3] presenta un modello per la soluzione del pro-

blema di HI basato sulla Programmazione Lineare Intera (PLI). De-

scriviamo per prima la formulazione TIP, che Gusfield stesso definisce

concettuale perche sarebbe generalmente poco pratica da utilizzare.

Successivamente illustreremo delle semplici osservazioni che rendono

questa formulazione praticamente applicabile su insiemi di dati di reale

interesse biologico.

Sia gi un vettore ternario (gi ∈ {0, 1, 2}n) che rappresenta un ge-

notipo di n SNP e sia li il numero dei suoi siti eterozigoti (li = |{k ∈{1, ..., n} : gi(k) = 2}|, numero dei siti di gi con valore 2). Esi-

stono esattamente 2li aplotipi compatibili con gi e 2li−1 coppie pij

(j ∈ {1, ..., 2li−1}) di aplotipi che potrebbero aver generato gi. Consi-

deriamo la totalita di queste coppie e generiamo per ciascuna di esse

una variabile yij. Tali variabili avranno il seguente significato:

yij =

1 se nella soluzione si risolve il genotipo gi con la coppia pij

0 altrimenti

Per ogni coppia pij = (hu, hv) si devono poi considerare separata-

29

mente gli aplotipi hu ed hv che la compongono, ricordando che il fatto

che pij sia una coppia di gi vuol dire che hu ⊕ hv = gi. Generiamo per

ogni aplotipo hk una variabile xk con il seguente significato:

xk =

1 se nella soluzione si prende l’aplotipo hk

0 altrimenti

Per quanto riguarda i vincoli, vogliamo che ogni genotipo gi sia

risolto da almeno una coppia di aplotipi pij e che, se viene scelta una

coppia, nella soluzione siano presenti i due aplotipi che la compongono.

Sara allora necessario introdurre due tipologie di vincoli:

• vincoli di copertura, per assicurare che per ogni genotipo sia presa

almeno una coppia

∑j yij ≥ 1 ∀gi ∈ G (3.1)

• vincoli di attivazione per forzare nella soluzione gli aplotipi di

ogni coppia presa xu ≥ yij

xv ≥ yij

per ogni coppia pij = (hu, hv) di gi,∀gi ∈ G.

(3.2)

Variabili e vincoli vanno creati allo stesso modo per ogni genotipo

gi ∈ G tenendo presente che se un aplotipo hk e compatibile con piu

di un genotipo, comunque sara rappresentato da una sola variabile xk.

Dobbiamo infine imporre che tutte le variabili siano intere e binarie: x ∈ {0, 1} ∀x ∈ X

y ∈ {0, 1} ∀y ∈ Y(3.3)

(dove Y e l’insieme di tutte le variabili di tipo yij).

Per quanto rigurda la funzione obiettivo, vogliamo minimizzare

il numero degli aplotipi nella soluzione, quindi possiamo esplicitarla

30

come:

min∑

xk∈X

xk (3.4)

dove X e l’insieme di tutte le variabili di tipo x.

L’insieme dei vincoli (3.1), (3.2) e (3.3), con la funzione obiettivo

vista, definiscono il problema secondo la formulazione TIP. Cerchiamo

ora di raffinare il problema per evitare vincoli e variabili che non in-

fluiscono sul risultato del problema. La prima idea e che se gli aplotipi

hu e hv di una coppia pij non sono parte di altre coppie, non sono

quindi compatibili con altri genotipi, possiamo escludere dal problema

di PLI le variabili xu, xv e yij perche ogni altra coppia del genotipo gi

sara non peggiore rispetto alla funzione obiettivo. Se in questo modo

eliminiamo tutte le coppie per un genotipo gi, allora vuol dire che gi

e non compatibile con i genotipi di G \ gi, i.e., esiste almeno un sito

k ∈ {1, ..., n} dove gi(k) = 0 e gj(k) = 1 (o gi(k) = 1 e gj(k) = 0) per

ogni gj ∈ {G \ gi}, allora la scelta della coppia che risolve gi non in-

fluisce sul valore della soluzione. Gusfield indico questa formulazione

ridotta con RTIP.

Riassumendo, la formulazione di PLI di Gusfield per il problema

di Inferenza di Aplotipi con Parsimonia Pura e:

min∑

x∈X x

∑j yij ≥ 1 ∀gi ∈ G

xu ≥ yij

xv ≥ yij

∀pij = (hu, hv) t.c. hu ⊕ hv = gi,∀gi ∈ G

x ∈ {0, 1} ∀x ∈ X

y ∈ {0, 1} ∀y ∈ Y

(3.5)

3.2.1 Dimensioni della formulazione TIP

Vogliamo ora esaminare la formulazione TIP di Gusfield per avere una

stima delle sue dimensioni.

Se li e il numero di siti eterozigoti del genotipo gi, allora gli aplotipi

31

compatibili con gi sono 2li e compongono 2li−1 coppie. Il numero delle

variabili per ogni genotipo e O(2li) per quanto riguarda gli aplotipi e

O(2li−1) per quanto riguarda le coppie. Se il numero dei genotipi e

m e il numero dei siti di ogni genotipo e n, il numero degli aplotipi

sara nh ≈ O(2n) e il numero delle coppie np ≈ O(2n−1). Sappiamo che

complessivamente il numero di aplotipi non puo essere maggiore di 2n,

che e il numero di vettori binari di dimensione n. Per quanto riguarda

i vincoli, ne avremo m di copertura (3.1) e 2 ·np (O(2n)) di attivazione

(3.2). I vincoli di interezza (3.3) infine sono nh + np (O(2n)).

Le dimensioni esponenziali del modello lo rendono in molti casi

inapplicabile anche riducendolo alla formulazione RTIP.

Come abbiamo visto il parametro fondamentale che rende un pro-

blema risolvibile o meno non e tanto il numero di genotipi ne il numero

di SNP, quanto il numero massimo di siti eterozigoti rilevati in una

istanza del problema. Per caratterizzare computazionalmente il pro-

blema in funzione di questo parametro introduciamo una notazione

per le classi di istanze: un problema appartiene alla classe delle istan-

ze (k, l)-limitate se la relativa matrice dei genotipi di input ha al piu

k ‘2’ per ogni riga e l ‘2’ per ogni colonna. Un * al posto di k indica

istanze limitate solo per quanto riguarda le colonne e un * al posto di

l indica istanze limitate solo per quanto riguarda le righe.

E stato dimostrato in [19] che il problema HIPP appartiene alla

classe dei problemi APX-difficili. D’altra parte in [12] viene dimostrato

come la classe dei problemi (3,3)-limitati sia APX-difficile, mentre le

classi (2,*) e (*,1) sono risolvibili in tempo polinomiale come anche la

classe (*,2) sotto particolari proprieta di compatibilita tra gli aplotipi.

Per quanto riguarda questo lavoro di tesi prenderemo in esame

problemi di dimensioni molto maggiori, quindi ci basta sapere che

il problema e APX-difficile, il che vuol dire che, volendo risolverlo

all’ottimo, non esistono algoritmi polinomiali per risolverlo (a meno di

dimostrare che P = NP) su una macchina di touring deterministica.

32

3.3 Algoritmi statistici

Oltre ai modelli di PLI, un’altro approccio molto utilizzato per la solu-

zione del problema HIPP e quello del modello statistico. Gli elementi

principali di un modello per generare dati relativi agli aplotipi e ai

genotipi sono:

• G = (g1, ..., gm) e l’insieme dei genotipi osservati

• H = (H1, ..., Hm) sono le corrispondenti coppie di aplotipi inco-

gnite, dove Hi = (hi1, hi2) per i = 1, ...,m

• ΘH = (θ1, ..., θv) indica il valore delle v frequenze non note degli

aplotipi (dove v e la cardinalita dell’insieme di tutti gli aplotipi)

• l’equilibrio di Hardy-Weinberg, che dice che la probabilita di

osservare un genotipo e pari al prodotto delle probabilita di

osservare gli aplotipi che lo costituiscono:

pr(g) =∑

(h,h′):h⊕h′=g

pr(h) · pr(h′),

dove pr(h) e la probabilita di trovare l’aplotipo h.

L’introduzione dei modelli statistici e dovuta principalmente al-

la presenza di alcuni difetti presenti nella formulazione presentata in

[16]. Infatti l’algoritmo illustrato in [16] richiede un insieme iniziale

di aplotipi risolti e la sua soluzione dipende fortemente dall’ordine in

cui vengono risolti gli aplotipi. L’approccio illustrato in [3] diminuisce

l’importanza di tale dipendenza, rendendo la procedura piu affidabile.

L’idea principale dei modelli statistici e che gli aplotipi hanno una di-

stribuzione di probabilita non nota nella popolazione in esame e che

i genotipi osservati di ogni individuo sono semplicemente combinazio-

ni di due aplotipi presi casualmente dalla popolazione. L’obiettivo

dell’approcico statistico al problema di HI e quindi di stimare le fre-

quenze degli aplotipi in modo che sia possibile calcolare facilmente gli

33

aplotipi di ogni individuo basandosi sulle distribuzioni di probabilita e

su assunzioni biologiche (come per esempio l’assunzione di casualita).

Mostriamo ora, come esempi della formulazione di questo problema,

la formulazione di Maximum Likelihood inference, presentata in [17] e

la formulazione di Bayesian Inference presentata in [18].

Frequencies Haplotype Inference problem (FHI):

Input: un insieme G di genotipi.

Output: l’insieme delle frequenze degli aplotipi {h1, ..., hv}(dove v e il numero di tutti i possibili aplotipi) che massi-

mizza la funzione di somiglianza (likelihood) per l’insieme

di genotipi osservato.

Bayesian Frequencies Haplotype Inference problem

(BFHI):

Input: un insieme G di genotipi e una distribuzione a

priori delle frequenze dei genotipi.

Output: la distribuzione a posteriori delle frequenze degli

aplotipi dato G.

Il poblema FHI viene presentato in [17], dove un algoritmo di

Expectation-Maximization (EM) viene proposto per stimare per gli

aplotipi le frequenze che massimizzano la funzione di somiglianza del-

l’insieme dei genotipi. Sia ΘHt l’insieme delle frequenze degli aplotipi

e Gt l’insieme delle probabilita di tutti i genotipi al passo t. L’algorit-

mo EM assegna un valore iniziale ΘH0 (un insieme iniziale possibile di

frequenze e quello che corrisponde all’assunzione che tutti i possibili

aplotipi siano equiprobabili). Basandosi su ΘH0 puo essere facilmen-

te calcolato il valore atteso di un genotipo di G e puo quindi essere

facilmente calcolato il valore di G1. I valori attesi dei genotipi in G1

vengono usati successivamente per stimare le frequenze gli aplotipi, il

che porta a calcolare ΘH1. L’algoritmo consiste nell’iterare i due passi

fino alla convergenza (finche la differenza tra ΘHt+1 ed ΘHt non sia piu

piccola di un valore predefinito). Ad ogni iterazione il valore di ΘHt

34

viene migliorato massimizzando la funzione di somiglianza dell’insie-

me G. Si possono prendere valori iniziali differenti per aumentare la

possibilita di ottenere un ottimo globale.

Il problema BFHI prende in consideraizone due tipi di probabilita:

a priori e a posteriori. La probabilita a priori puo essere interpretata

come cio che si conosce gia della variabile. La probabilita a posteriori

e la probabilita che viene calcolata in funzione dell’insieme di genotipi

G dato, delle probabilita a priori e di una funzione di somiglianza

che tiene conto del teorema di Bayes. Viene utilizzata con l’obiettivo

di attribuire incertezza piu che casualita ad una quantita incerta. Si

applica il teorema di Bayes, moltiplicando le probabilita a priori per la

funzione di somiglianza e normalizzando, per ottenere la probabilita

a posteriori, che e la probabilita condizionata dell’incertezza dei dati

forniti.

3.4 Algoritmi polinomiali

Il problema dell’inferenza di aplotipi con parsimonia pura e stato anche

affrontato cercando una formulazione alternativa quella di Gusfield per

cercare di ridurre non tanto la complessita dell’algoritmo di soluzione,

tenendo sempre presente che e un problema della classe di complessita

APX-difficile, quanto nelle dimensioni del problema da risolvere.

Brown e Harrower in [13] presentano un modello di rappresenta-

zione del problema polinomiale per numero di variabili e vincoli nel

numero dei genotipi e di SNP dell’input.

Ogni genotipo gi e spiegato dagli aplotipi h2i−1 e h2i e le variabili

yi,k rappresentano il valore dell’aplotipo i nella posizione k. I vincoli

y2i−1,k ⊕ y2i,k = gi(k) ∀1 ≤ i ≤ m, ∀k (3.6)

assicurano che ogni genotipo sia spiegato propriamente. Per contare

gli aplotipi unici che spiegano i genotipi di input, vengono usate delle

variabili binarie per indicare le differenze tra coppie di aplotipi. Per

35

ogni coppia di aplotipi (i, j), 1 ≤ i < j ≤ 2m, una variabile di,j e

uguale ad 1 quando hi 6= hj. Se hi 6= hj c’e almeno una posizione k

per la quale hi(k) 6= hj(k) = 0 oppura hi(k) 6= hj(k) = 1 e allora le

seguenti equazioni forzeranno di,j a 1:

di,j ≥ yi,k − yj,k ∀ 1 ≤ i < j ≤ 2m, ∀kdi,j ≥ yj,k − yi,k ∀ 1 ≤ i < j ≤ 2m, ∀k.

(3.7)

Per ogni aplotipo hi si ha una variabile xi che e pari a 1 se hi e

unico in (h1, ..., hi). Questa condizione e garantita dall’equazione:

xi ≥ 2− i +∑i−1

j=1 dj,i ∀1 ≤ i ≤ 2m. (3.8)

La funzione obiettivo min∑2m

i=1 xi porta il problema a minimizzare

il numero di aplotipi unici, quindi anche il numero totale di aplotipi

ed i vincoli di interezza

xi, yi,k, di,j ∈ {0, 1} ∀i, j e k (3.9)

completano il problema di PLI.

Il problema cosı strutturato ha, rispetto all’approccio illustrato da

Gusfield, il grande merito di essere comunque trattabile, mentre le

formulazioni TIP ed RTIP portano spesso a problemi con un numero

eccessivo di variabili e di vincoli.

36

Capitolo 4

La formulazione di Set

Covering puro

4.1 Il passo di Fourier e la formulazione

di Covering puro

Obiettivo centrale di questa tesi e di dimostrare le proprieta di una

formulazione del problema HIPP proposta da A. Godi in [1] che riporta

il modello di PLI di Gusfiel illustrato nel paragrafo 3.2 ad un problema

di Set Covering puro.

Il modello di Gusfield prevede l’utilizzo di due tipi di variabili: va-

riabili x per gli aplotipi e variabili y per le coppie, ma nella funzione

obiettivo sono presenti solo le variabili x. Questa struttura suggerisce

di applicare la procedura di Fourier-Motzkin per proiettare il poliedro

associato al modello in uno spazio di dimensioni minori. La proce-

dura di Fourier-Motzkin fu inventata da Fourier (1827) e riscoperta

successivamente da Dines (1928) e Motzkin (1936).

Diamo un’idea generale di come funziona l’algoritmo di Fourier-

Motzkin. Sia I un insieme di equazioni e x una variabile. Sia Ci una

disequazione in I dove x presenta un coefficiente negativo −ci e sia

Cj una disequazione in I dove x presenta un coefficiente positivo cj,

37

allora

ciCj − cjCi

e una disequazione che non contiene x. L’insieme di tutte le equazioni

di questo tipo e di tutte le disequazioni in I che non contengono la

variabile x e il risultato di un passo di Fourier che elimina x da I.

Se x non ha nessun coefficiente positivo o non ha nessun coefficiente

negativo nelle disequazioni di I, allora il passo di Fourier per eliminare

x da I cancella semplicemente tutte le disequazioni che contengono x.

L’algoritmo di Fourier-Motzkin consiste nell’applicare iterativamente il

passo di Fourier fino a che non si presenta una delle seguenti condizioni:

- l’insieme I contiene una disequazione non ammissibile

- tutte le variabili sono state eliminate e I non contiene contrad-

dizioni.

Se ad un certo passo dell’algoritmo nell’insieme si trova una di-

sequazione non ammissibile, allora il problema iniziale non ammette

soluzioni. Nel secondo caso il problema iniziale ammette soluzioni. La

complessita asintotica dell’algoritmo per un insieme di m disequazioni

e n variabili e O(m2n).

Applichiamo ora la procedura di Fourier-Motzkin alla Formula-

zione TIP di Gusfield descritta in (3.5) con l’intento di eliminare le

variabili y. Prendiamo per esempio il genotipo g = 20021 e i relativi

aplotipi h1 = 00001, h2 = 10011, h3 = 00011 e h4 = 10001 che forma-

no le coppie p1 = (h1, h2) e p2 = (h3, h4). Abbiamo quattro variabili

di aplotipo xi, i ∈ {1, ..., 4}, e due variabili di coppia y1 per p1 e y2

per p2. I relativi vincoli sono:

y1 + y2 ≥ 1 (4.1)

x1 − y1 ≥ 0 (x1 ≥ y1) (4.2)

x2 − y1 ≥ 0 (x2 ≥ y1) (4.3)

38

x3 − y2 ≥ 0 (x3 ≥ y2) (4.4)

x4 − y2 ≥ 0 (x4 ≥ y2). (4.5)

Sommando la (4.1) alla (4.2) e alla (4.4) otteniamo

x1 + x3 ≥ 1.

Ripetendo l’operazione per le quattro combinazioni possibili otte-

niamo l’insieme di vincoli:

x1 + x3 ≥ 1

x1 + x4 ≥ 1

x2 + x3 ≥ 1

x2 + x4 ≥ 1

Ripetendo la procedura per ogni genotipo, otteniamo la formula-

zione che chiameremo di Covering (HIc):

Formulazione HIc

min∑

∀xk∈X

xi

∑(hu,hv):hu⊕hv=g

∑hk∈{hu,hv} xk ≥ 1 ∀g ∈ G

xk ∈ {0, 1} ∀xk ∈ X(4.6)

Dato un genotipo, i vincoli a lui associati sono del tipo:

date tutte le coppie, si prende un aplotipo per ogni cop-

pia, la somma delle variabili di tali aplotipi deve essere

maggiore o uguale a 1.

Perche la formulazione sia completa devo prendere per ogni geno-

tipo tutti i vincoli possibili di questo tipo. In questa formulazione

39

ho eliminato le variabili di tipo y, ma a costo di aumentare in modo

esponenziale il numero dei vincoli.

Per sapere quanti saranno i vincoli per ogni genotipo cerchiamo di

rappresentarli nel seguente modo: per un genotipo gi che ha li ‘2’ ho

2li aplotipi e 2li−1 coppie. Vogliamo sapere quante sono sono le diverse

combinazioni del tipo visto precedentemente. Possiamo rappresentare

un vincolo con un vettore binario di 2li−1 elementi dove il j-esimo ele-

mento mi dice se nel vincolo ho preso il primo elemento della coppia

pi,j (0) o il secondo (1). E facile vedere come ogni vettore binario di

questo tipo rappresenta uno di questi vincoli, che sono quindi 22li−1.

Se gia la formulazione TIP portava il problema ad essere ingestibile

con un piccolo aumento delle dimensioni, questa formulazione sara an-

cora meno gestibile nella sua forma completa. Per questo nei prossimi

paragrafi si cerca un modo per estrapolare il minor numero possibile

di vincoli necessari per avere una soluzione ammissibile.

4.2 L’algoritmo di Branch and Cut

La Programmazione Lineare Intera (PLI) e un problema noto e molto

utilizzato in Ricerca Operativa ed appartiene alla classe dei proble-

mi NP-difficili. Nella teoria della complessita, l’inseme NP (e i suoi

sottoinsiemi, come l’insieme dei problemi NP-difficili) e l’insieme dei

problemi di decisione che possono essere risolti in tempo polinomia-

le su una macchina di Touring non deterministica e possono essere

verificati in tempo polinomiale da una macchina di Touring determi-

nistica. Dato che i calcolatori su cui lavoriamo sono implementazioni

di macchine di Touring deterministiche non ci e possibile risolvere un

problema NP in tempo polinomiale. In realta determinare se se sia

possibile risolvere problemi NP in tempo polinomiale (determinare se

P = NP) e una delle grandi questioni ancora aperte dell’informatica

teorica.

Un problema di PLI e composto principalmente da:

40

- una funzione obiettivo del tipo:

min∑

xi∈X

cixi (o max∑

xi∈X

cixi)

- vincoli che definiscono il problema

∑j

∑xi∈X

aijxi ≥ bj

- vincoli di bound per le variabili

lbi ≤ xi ≤ ubi ∀xi ∈ X

- e vincoli di interezza, non necessariamente per tutte le variabili

xi intera ∀xi ∈ XI ⊆ X.

Per le caratteristiche del problema in esame sara sufficiente d’ora

in avanti considerare esclusivamente problemi di minimo con variabili

binarie.

L’algoritmo piu utilizzato per risolvere questo tipo di problemi e

il Branch and Bound. L’algoritmo di Branch and Bound si basa sul-

l’idea che un problema di Programmazione Lineare (non intera) puo

essere risolto in tempo polinomiale e quindi si puo cercare di risolvere

un problema di PLI non considerando i vincoli di interezza (rilassan-

do i vincoli di interezza del problema) con un algoritmo polinomiale

e vedere se la soluzione ottenuta e o no intera. Se la soluzione e in-

tera abbiamo una soluzione ottima perche una soluzione ottima del

rilassamento di un problema non puo essere peggiore di una soluzione

ottima del problema non rilassato. Se la soluzione e non intera sara

comunnque, per quanto detto prima, un limite inferiore oltre il quale

il problema non rilassato non potra andare e avremo che almeno una

variabile x ha valore non intero. Possiamo allora copiare il problema in

41

due sottoproblemi ed imporre da una parte il vincolo x = 0 e dall’altra

x = 1 (trattando esclusivamente variabili binarie), risolvere i relativi

rilassamenti e confrontare i risultati. Il processo di divisione di un

problema in due sottoproblemi secondo la procedura sopra indicata e

detto di branching e, se iterato su piu variabili, porta a sviluppare un

albero binario, detto albero di Branch and Bound.

Tale albero ha per ogni ramo come vincolo che viene imposto al

relativo sottoalbero una variabile che viene forzata a 0 o ad 1 e per ogni

foglia il problema originale a cui sono stati aggiunti vincoli che forzano

una o piu variabili. Confrontando i risultati ottenuti per le foglie si

puo cercare di orientare lo sviluppo dell’albero al fine di ottenere un

risultato intero ottimo per il problema originale.

La complessita asintotica dell’algortimo di Branch and Bound e

esponenziale dal momento che l’albero che viene creato ha nel caso

peggiore (quando viene costruito interamente) 2v foglie (per variabili

binarie) dove v e il numero delle variabili del problema.

Una variante dell’algoritmo, chiamata algoritmo di Branch and

Cut, prevede la possibilita ad ogni nodo dell’albero di introdurre dei

tagli. I tagli sono vincoli del problema che rendono non ammissibile la

soluzione corrente non intera, ma non influiscono sulle soluzioni intere.

Introducendo un taglio si e sicuri che non verra piu presa la soluzione

non intera attuale e che restano ammissibili tutte le soluzioni intere

del problema originario.

La formulazione TIP di Gusfield porta ad un problema di PLI

con un numero esponenziale di variabili e di vincoli. Risolvere tale

problema con l’algoritmo di Branch and Bound vuol dire risolvere con

un algoritmo polinomiale un insieme di problemi di PL, che nel caso

peggiore sono in numero esponenziale, sapendo che le dimensioni dei

problemi sono sempre esponenziali per numero di vincoli e di variabili.

Rispetto alla formulazione TIP, nella formulazione HIc viene ridot-

to il numero delle variabili in quanto vengono eliminate tutte le va-

riabili di coppia, ma aumentano esponenzialmente i vincoli da O(2n)

a O(22n). Quest’aumento sarebbe chiaramente insostenibile per una

42

formulazione completa.

L’algoritmo che e stato utilizzato in questo lavoro di tesi mira a in-

trodurre nel problema il minor numero possibile di vincoli, in modo da

non compromettere la ricerca dell’ottimo e di contenere le dimensioni

del problema nel caso medio.

4.3 L’algoritmo HI PARS

Partendo dalla formulazione HIc, dal problema completo Pc vogliamo

calcolare un sottoinsieme il piu piccolo possibile di vincoli che ci porti

ad una soluzione ammissibile. Una soluzione trovata in questo modo

e sicuramente ottima perche la soluzione di un rilassamento di un

problema di PLI e non peggiore di una soluzione del problema non

rilassato.

All’i-esimo nodo dell’albero di Branch and Bound avremo un pro-

blema Pi con le stesse variabili del problema Pc, ma con un sottoin-

sieme dei suoi vincoli. Ci serve ora una procedura che ci indichi se

la soluzione x∗ ottima per Pi e o no ammissibile per Pc e restituisca

eventualmente uno o piu vincoli di Pc che taglino x∗.

Per fare cio introduciamo il seguente teorema:

Teorema AG: dato un problema Pc costruito secondo la

formulazione HIc e una soluzione intera x∗, x∗ e ammissibile

per Pc se per ogni genotipo g del problema vale la seguente

disequazione:

∑(hu,hv):hu⊕hv=g

min{x∗u, x∗v} ≥ 1. (4.7)

E facile intuire che se per un genotipo g ho preso almeno una

coppia (hu, hv) allora tutti i vincoli relativi a g sono soddisfatti perche

in tutti e presente o xu o xv. Se per g non ho preso nessuna coppia

esiste almeno un vincolo di Pc che non e soddisfatto e che taglia x∗.

43

Facciamo un semplice esempio. Sia g un genotipo con tre 2, siano

{h1, h2}, {h3, h4}, {h5, h6} e {h7, h8} le sue coppie di aplotipi. Sup-

poniamo di non prendere la prima coppia e sia, senza perdita di ge-

neralita, x1 = 1 e x2 = 0. Supponiamo, seguendo lo stesso ragio-

namento che sia x3 = x5 = x7 = 1 e x4 = x6 = x8 = 0. Allora

x∗ = (1, 0, 1, 0, 1, 0, 1, 0) e

x2 + x4 + x6 + x8 ≥ 1

e un vincolo di Pc ed e un taglio per x∗.

Chiameremo Procedura di AG la procedura grazie alla quale

viene verificato, dati un problema Pc e una soluzione x∗ e seguendo

il teorema AG, se esistono in Pc tagli alla soluzione x∗ e che even-

tualmente restituisce un insieme CAG di tali tagli, al piu uno per ogni

genotipo.

La procedura AG puo essere scritta in questo modo:

• sia CAG un insieme vuoto di vincoli

• per ogni genotipo g:

– se∑

(hu,hv):hu⊕hv=g min{x∗u, x∗v} ≥ 1 non ci sono tagli per x∗

relativamente a g (ho gia preso almeno una coppia per g)

– altrimenti viene aggiunto a CAG il vincolo

∑(hu,hv):hu⊕hv=g

xargmin{x∗u,x∗v} ≥ 1.

Diamo ora le linee guida dell’algoritmo di Branch and Cut che e

stato implementato ai fini di questa tesi per risolvere il problema HIc

utilizzando la procedura AG per la ricerca dei tagli.

• Il problema iniziale P0 viene inizializzato con vincoli che im-

pongono di prendere gli aplotipi dell’insieme iniziale H0 visto in

44

3.1:

xk ≥ 1 ∀hk ∈ H0 (4.8)

e viene introdotto nell’insieme dei problemi aperti.

• finche non viene trovata una soluzione ottima

– viene scelto un problema Pi dall’insieme dei problemi aperti

– vengono rilassati i vincoli di interezza e viene calcolata la

soluzione ottima (non necessariamente intera) x∗

– finche ci sono nuovi tagli con la procedura AG, vengono

aggiunti i nuovi tagli e si risolve nuovamente il problema

– quando non ci sono piu nuovi tagli se la soluzione e intera

si ha un ottimo locale

– se la soluzione e non intera si generano due nuovi problemi,

applicando la procedura di branching su una delle variabili

frazionarie nella soluzione ottenuta, e tali problemi vengono

aggiunti all’insieme dei problemi aperti.

4.4 Il preprocessamento

Date le dimensioni del problema, e molto importante verificare a priori

se e come sia possibile ridurre il numero dei vincoli e delle variabili.

Consideriamo l’insieme dei genotipi G. Sappiamo che i genotipi

che hanno al piu un 2 sono univocamente determinati e possono essere

considerati a parte. Chiamiamo Gud l’insieme di tali genotipi. Sia

Gnc l’insieme dei genotipi non compatibili con gli altri genotipi di G

(Gnc = {gi ∈ G : ∃p ∈ {1, ..., n} t.c. ∀gj ∈ G \ gi, gi(p) + gj(p) = 1}).Come gia detto nel paragrafo 3.2, un genotipo non compatibile con

nessun altro genotipo di G puo essere processato indipendentemente

dagli altri perche la scelta degli aplotipi che lo spiegano non inficia la

soluzione degli altri genotipi, dal momento che il criterio secondo cui

vengono scelti gli aplotipi e la parsimonia pura.

45

Sia Ha l’iniseme degli aplotipi che spiegano i genotipi di Gud, pos-

siamo dire che se un genotipo puo essere spiegato da due aplotipi

contenuti in Ha, allora e inutile generare nuovi vincoli per quel geno-

tipo. Sia Gs l’insieme di tali aplotipi, Gs = {gi ∈ G : ∃{h′, h′′} ∈Ha t.c. h′ ⊕ h′′ = gi}.

Dall’insieme G abbiamo ora calcolato facilmente un insieme di

aplotipi obbligati Ha e il sottoinsieme di genotipi rimanenti Gr =

G \ {Gud ∪Gnc ∪Gs} sui quali applicheremo l’algoritmo di Branch &

Cut.

Seguendo la distinzione tra la formulazione TIP di Gusfield e la

formulaizone RTIP, considereremo per ognuno dei genotipi di Gr non

tutte le possibili coppie, ma solo quelle coppie per le quali almeno uno

dei due aplotipi e compatibile con almeno un altro genotipo di Gr o

di Gud.

Per dare un idea del possibile risparmio di tempo e spazio che

questa procedura puo apportare consideriamo i genotipi

g1 = 222121

g2 = 222220

g3 = 111202.

Dato che g2 ha cinque 2, e compatibile con 25 aplotipi (24 coppie).

Nella formulazione HIc il numero complessivo di vincoli relativi a g2

e 224= 216 = 65.536. Ma:

g1 e g2 non hanno aplotipi in comune (non sono compatibili) dato che

g1(6) = 1 e g2(6) = 0.

g1 e g3 hanno in comune solo l’aplotipo h1 = 111101

g2 e g3 hanno in comune gli aplotipi h2 = 111000 e h3 = 111100

Complessivamente prendo quindi in esame una coppia di aplotipi

per g1, due coppie per g2 e tre per g3 senza compromettere l’ottimalita

della formulazione.

46

Cio non deve illudere che ogni istanza sottoposta al programma

sia facilmente risolvibile in tempo breve. Ma, come sara illustrato nel

prossimo capitolo, il piu delle volte l’algoritmo arriva ad una soluzione

ottima in tempi molto brevi. Resta il fatto che se l’ultimo valore di g1

fosse stato uno 0 avremmo avuto 24 aplotipi in comune tra g1 e g2 e

di conseguenza 24 coppie l’uno, quindi 224= 65.536 vincoli ciascuno.

Avremmo piu di cinquanta variabili per piu di centomila vincoli nella

formulazione completa per un problema molto simile al precedente e

che comunque resta di dimensioni molto ridotte per numero di SNP.

47

Capitolo 5

Dettagli implementativi e

risultati computazionali

In questo capitolo vengono presentati i principali dettagli implemen-

tativi del programma HI_PARS, implementato in C per risolvere il pro-

blema dell’Inferenza di Aplotipi secondo il modello HIc presentato nel

capitolo precedente e con l’algoritmo di Branch and Cut. Successi-

vamente vengono presentati i risultati ottenuti dall’ambiente di test

sviluppato per verificare l’efficienza del programma.

Per gestire i nodi dell’albero di Branch and Cut, risolvere i rilas-

samenti lineari dei problemi di PLI ed aggiungere i tagli al problema

e stato utilizzato il programma CPLEX 8.0 dell’azienda di software di

ottimizzazione ILOG, messo a disposizione dall’Istituto di Analisi dei

Sistemi ed Informatica (IASI) del Consiglio Nazionale delle Ricerche,

presso cui e stato svolto il tirocinio da cui nasce questo lavoro di tesi.

Il software e stato eseguito su un computer dello IASI, in parti-

colare HYDRA03, terza “testa” del cluster Hydra. HYDRA03, come

gli altri nodi del cluster, ha un processore P4 a 3.6 GHz e 2 GB di

memoria RAM.

48

5.1 Dettagli implementativi

L’ambiente di test sviluppato per questo lavoro di tesi si divide in:

• generazione degli aplotipi con il programma ms di R.R.Hudson

• generazione delle istanze di genotipi

• inizializzazione delle strutture dati su CPLEX

• soluzione del problema secondo l’algoritmo presentato.

Il programma ms, presentato in [21] di R. R. Hudson, permette di

generare campioni di dati come stringhe binarie secondo diversi mo-

delli neutri ed e stato implementato con l’obiettivo di permettere di

investigare le proprieta statistiche di tali campioni, per dare stime di

test statistici e in generale per aiutare l’interpretazione di insiemi di

dati polimorfici. In molte pubblicazioni sulla biologia computaziona-

le, come per esempio [3], [1] e [20], viene utilizzato ms per genera-

re gli insiemi di aplotipi di partenza. Tali aplotipi sono stati presi

successivamente a due a due casualmente per generare i genotipi.

Nel nostro problema gli aplotipi sono rappresentati come array bi-

nari e, per contenere le dimensioni delle strutture dati che li rappresen-

tano, vengono memorizzati in un albero binario. Partendo dalla radice

dell’albero tutti gli aplotipi che iniziano con uno 0 appartengono al sot-

toalbero sinistro, tutti gli aplotipi che iniziano con un 1 appartengono

al sottoalbero destro. Iterando per tutti gli SNP, per ogni aplotipo dif-

ferente ho una foglia differente, che viene indicizzata con un numero

progressivo. Tali indici rappresenteranno le variabili del modello im-

plementato. Nel file bintree.h ho definito lo struct bin_tree per la

rappresentazione dell’albero binario e le funzioni che mi permettono di

inserire e gestire facilmente nuovi aplotipi. In particolare la funzione

int insHaplo(int* h, int dim, bin_tree* tree)

aggiunge all’albero binario un nuovo aplotipo (array binario) e resti-

tuisce l’intero che lo identifica univocamente.

49

I genotipi sono rappresentati con un apposito struct, lo struct geno,

in cui sono memorizzati il valore del genotipo (come stringa ternaria)

ed i relativi aplotipi, rappresentati con gli interi che li identificano

nell’albero binario. Nel file geno.h sono definiti lo struct geno ed

alcune funzioni utili per gestire tale struct. In particolare la funzione

commonH(geno *gi, geno *gj, bin_tree *tree)

calcola, dati due genotipi gi e gj, l’insieme Hcomm degli aplotipi che

hanno in comune e aggiunge a gi e a gj gli aplotipi di Hcomm con la

funzione

addH(geno *gi, int **hcom, int hsize, bin_tree *tree).

Quest’ultima, per ogni aplotipo hu ∈ Hcomm, se hu non e ancora stato

aggiunto agl genotipo gi, calcola l’aplotipo hv complementare di hu

rispetto a gi (hv tale che hu ⊕ hv = gi) ed aggiunge entrambi alla lista

dei genotipi di gi.

La funzione commonH viene utilizzata nell’inizializzazione dei dati

per implementare il preprocessamento indicato da Gusfield ed illustra-

to nel paragrafo 3.2.

Facciamo un esempio dell’applicazione di queste funzioni. Dati

g1 = 022

g2 = 220

la funzione commonH calcola che i due genotipi hanno in comune gli

aplotipi

h1 = 000

h2 = 010.

Se h1 non e stato gia agigunto a g1, addH calcola h3 = 011 (h1⊕h3 = g1)

e aggiunge entrambi a g1. Ripetendo l’operazione per h2 e per g2 avro

gli aplotipi

h4 = 001 (h2 ⊕ h4 = g1)

50

h5 = 110 (h1 ⊕ h5 = g2)

h6 = 100 (h2 ⊕ h6 = g2).

Figura 5.1: Albero binario degli aplotipi

Nel file toCplex.h ho definito lo struct problem, che contiene tut-

ti i valori necessari per caricare il problema CPLEX e alcune funzioni

necessarie per la creazione, aggiornamento e gestione di un problema.

La gestione dei problemi in CPLEX e abbastanza complessa e le funzio-

ni di toCplex.h offrono un’interfaccia abbastanza semplice, studiata

appositamente per il problema in questione.

Per vedere con quali parametri CPLEX reppresenta un problema di

PLI (consideriamo un problema di minimo, con vincoli sempre ‘≥ 1’ e

variabili binarie), consideriamo la notazione comune:

min cT x

Ax ≥ 1

x ∈ {0, 1}(5.1)

Il numero delle righe e delle colonne sono dati dai numcols e

numrows. Il vettore objective rappresenta il vettore c, il vettore

rhs rappresenta la parte destra dei vincoli e i questo caso e un vettore

di tutti 1, il vettore sense rappresenta il senso di ogni vincolo, qui

sempre ‘≥’ (‘G’).

51

Gli array matbeg, matval, matcnt e matval rappresentano la ma-

trice A per colonne e rappresentano solo gli elementi diversi da 0.

matcnt indica quanti valori sono diversi da 0 per ogni colonna, in

matval sono registrati i valori di tutti gli elementi non zero presi per

colonne di A, matbeg indica il primo elemento di ogni colonna sul

vettore matval e matind indica gli indici di riga degli elementi in

matval.

Per esempio, la matrice

1 0 0 5

0 2 4 0

0 3 0 6

e rappresentata dai vettori

matval {1, 2, 3, 4, 5, 6}matind {0, 1, 2, 1, 0, 2}matbeg {0, 1, 3, 4}matcnt {1, 2, 1, 2}.

Il problema viene inizializzato con dei vincoli che impongano nella

soluzione gli aplotipi obbligati dell’insieme H0 visto nel paragrafo 3.1

e poi risolto con il programma HI_PARS, che implementa l’algoritmo

descritto nel paragrafo 4.3. All’i-esimo nodo dell’albero di branching,

dato il problema Pi, viene calcolata la soluzione del rilassamento li-

neare x∗i e viene applicata la procedura presentata in 4.3 per verificare

se esitano tagli alla soluzione data. Sia Ci l’insieme dei tagli alla so-

luzione x∗i . Se Ci 6= ∅ aggiungo a Pi i tagli di Ci e risolvo nuovamente

il rilassamento lineare. Altrimenti non ho nuovi tagli da aggiungere

e se la soluzione e intera ho trovato un ottimo locale e posso chiude-

re il nodo, altrimenti scelgo una variabile non intera sulla quale fare

branching. Posso chiudere un nodo anche quando il valore della fun-

zione obiettivo per la soluzione del rilassamento lineare e minore de

miglior ottimo locale finora trovato, quindi le soluzioni che troverei

sviluppando tale nodo sarebbero sicuramente subottime.

52

Quando ho chiuso tutti i nodi dell’albero sono sicuro che la solu-

zione migliore finora trovata e la soluzione ottima al problema.

5.2 Osservazioni sulla complessita delle

istanze

Prima di presentare le istanze del problema che sono state prese in

analisi, le modalita di generazione di tali istanze e i risultati che pos-

sono essere estrapolati dalle soluzioni da esse ottenute, sono necessa-

rie alcune osservazioni per la comprensione dell’ambiente di testing

sviluppato.

Partendo dalla formulazione TIP presentata da Gusfiled in [3], che

e stata utilizzata come base di partenza anche in questo lavoro di tesi,

la dimensione di esponenzialita del problema, come visto nel capitolo

precedente, si manifesta sotto due principali aspetti:

• la formulazione TIP porta ad un problema di PLI, la cui com-

plessita e nota essere esponenziale

• i problemi di PLI costruiti secondo tale formulazione sono a lo-

ro volta di dimensioni esponenziali per numero di variabili e di

vincoli.

La dimensione fondamentale del problema e data dal numero di siti

polimorfi di ogni genotipo (siti in cui il genotipo presenta un 2) che

e non nota a priori e limitata superiormente dal numero n di SNP

di ogni genotipo. Secondo la formulazione TIP, per ogni genotipo

vengono generati O(2l) aplotipi, dove l e il numero di 2 del genotipo,

e O(2l) vincoli.

Consideriamo quindi di difficile soluzione un insieme di genotipi

che contenga almeno un genotipo con un numero elevato di 2. Per

dare un’idea della dimensione della difficolta, restando sempre nella

formulazione TIP, per un singolo genotipo con un numero di 2 pari a

53

16 vengono generate circa 64 mila variabili aplotipo (e 32 mila varia-

bili di coppia) per circa 64 mila vincoli. Se il numero delle variabili al

considerare gli altri genotipi puo non crescere sensibilmente se geno-

tipi diversi hanno in comune piu aplotipi, stante il fatto che il limite

superiore del numero degli aplotipi e 2n, i vincoli complessivi saranno

la somma dei vincoli di ogni genotipo.

Riprendendo i concetti esplicati nelle introduzioni al problema bio-

logico ed al problema computazionale ricordiamo che l’Inferenza di

Aplotipi mira a estrapolare da un insieme di individui e secondo un

particolare modello biologico l’influenza che hanno avuto i genitori

di ogni individuo nella costituzione del suo patrimonio genetico. Dal

punto di vista computazionale, da un insieme di genotipi vogliamo

calcolare uno degli insieme di aplotipi che possa averli generati con

maggior probabilita (non e detto che ne esita uno solo), secondo un

modello biologico dato, nel nostro caso la massima parsimonia.

Per un individuo non prendiamo in esame l’intero patrimonio ge-

netico, dal momento che eccede per dimensioni dalle nostre capacita

computazionali, ma un sottoinsieme limitato e al caso pratico mirato

dei suoi SNP. Per ogni aplotipo non sara possibile, per quanto con-

cerne questo lavoro di tesi, distinguere il sesso dell’individuo da cui

abbiamo preso il genotipo, ne tantomeno possiamo distinguere a quale

dei due genitori, padre e madre, appartengano gli aplotipi che lo hanno

generato.

Dato che il numero degli SNP considerato e molto piccolo rispetto

al numero totale dei siti dell’intero patrimonio genetico, possiamo dire

che un genotipo ed un aplotipo in questo caso non corrispondano ad

un unico individuo, ma ad un insieme di individui, una famiglia, che

condividono un percorso evolutivo. Questo giustifica il fatto che tre

diversi aplotipi possano generare tre diversi genotipi, senza dare adi-

to ad ipotesi di moralita e pratiche procreative inusuali, e che in un

stesso insieme possano essere presenti due o piu genotipi uguali, non

necessariamente appartenenti allo stesso individuo.

In un problema reale la dimensione della difficolta e funzione del

54

numero degli SNP presi in esame, del numero degli individui e di

quanto siano diversi gli aplotipi che hanno generato ogni genotipo

(quindi di quanti siti eterozigoti abbia ogni genotipo). Ricordando

l’esempio nel paragrafo 4.4 possiamo dire che, per quanto visto in

questo lavoro di tesi, il fatto che una istanza del problema, insieme di

individui presi casualmente da una popolazione piu o meno eterogenea,

sia o no risolvibile per classi di problemi con piu di 20 SNP e puramente

casulae: per controllarne la difficolta si dovrebbe dare una dimensione

della diversita tra aplotipi ed interferire sul loro accoppiamento.

5.3 Definizione delle classi di istanze del

test

Il problema di Inferenza di Aplotipi puo essere visto nel seguente modo:

• esiste in una popolazione un insieme di aplotipi

• tali aplotipi si ricombinano generando dei genotipi

• nostro intento e di riscostruire dai genotipi l’insieme iniziale degli

aplotipi.

Indichiamo come soluzione l’insieme originario di aplotipi, istanza

del problema l’insieme dei genotipi generati da una soluzione e come

risultato l’insieme degli aplotipi calcolati con HI_PARS da una istanza.

L’accuratezza di un risultato e pari alla percentuale dei suoi aplotipi

che sono presenti anche nella relativa soluzione.

La generazione degli aplotipi e stata effettuata utilizzando il pro-

gramma ms di R. R. Hudson presentato in [21]. ms permette di generare

campioni di dati come stringhe binarie secondo diversi modelli neutri

ed e stato implementato con l’obiettivo di permettere di investigare le

proprieta statistiche di tali campioni, per dare stime di test statistici

e in generale per aiutare l’interpretazione di insiemi di dati polimor-

fici. In molte pubblicazioni sulla biologia computazionale, come per

55

esempio [3], [1] e [20], viene utilizzato ms per generare gli insiemi di

aplotipi di partenza.

Lasciando invariati il numero dei genotipi e degli SNP per ogni set

di dati (50 genotipi per 30 SNP) sono state generate classi di istanze

in funzione di:

• cardinalita dell’insieme degli aplotipi di partenza: 10, 30, 50

• fattore di ricombinazione; 0, 4, 16, 40

La cardinalita degli aplotipi di partenza e stata fatta variare con-

siderando che, a livello teorico e a meno di ripetizioni (che comunque

come abbiamo visto possono essere presenti), il legame tra il numero

di aplotipi e il numero di genotipi di una istanza e il seguente:

• dato un insieme di m genotipi, il massimo numero degli aplotipi

che possono aver generato tali genotipi e 2m e corrisponde al

caso in cui siano stati presi tutti aplotipi distinti

• dato un insieme di aplotipi di cardinalita a, il numero massi-

mo di genotipi diversi che possiamo generare e a(a + 1)/2 (se

consideriamo possibile ottenere genotipi senza nessun 2, quindi

ottenuti accoppiamento due copie dello stesso aplotipo).

Il numero minimo degli aplotipi distinti necessari per generare 50

genotipi distinti e 10 (con 9 si hanno massimo 45 genotipi) e se il

numero di aplotipi utilizzati inizialmente e minimo, il modello del-

la parsimonia pura ci porta presumibilmente a ricostruire l’insieme

originario.

Se consideriamo l’altro estremo, quindi m genotipi ottenuti da 2m

aplotipi diversi, e possibile che esista un insieme di aplotipi di cardi-

nalita minore a 2m che risolve l’istanza e a cui portera il modello della

parsimonia pura.

Per fare un esempio consideriamo i seguenti aplotipi:

h1 = 0000, h2 = 1100, h3 = 0010, h4 = 1000

56

e generiamo i genotipi

g1 = h1 ⊕ h2 = 2200, g2 = h3 ⊕ h4 = 2020.

E facile vedere che la soluzione

g1 = h1 ⊕ h2, g2 = h1 ⊕ h5

dove h5 = 1010, e migliore dal punto di vista della parsimonia, ma

peggiore per quanto riguarda l’accuratezza e che la soluzione

g1 = 0100⊕ 1000, g2 = 0000⊕ 1010

e ammissibile per quanto riguarda il problema computazionale ed ha

la stessa cardinalita della soluzione che stiamo cercando pur essendone

diversa.

Possiamo dire che all’aumentare del rapporto tra il numero degli

aplotipi ed il numero dei genotipi, aumenta la possibilita di trova-

re piu soluzioni ottime o addirittura soluzioni di cardinalita mino-

re. Per questo ci aspettiamo che l’accuratezza dei risultati sia minore

all’aumentare del numero degli aplotipi generatori.

Per quanto riguarda la ricombinazione, ms permette di generare

campioni secondo il modello dei siti finiti e ricombinazione uniforme.

In questo modello il numero dei siti tra i quali puo essere presente

ricombinazione e finito e specificato dall’utente. Nonostante il numero

finito di siti tra i quali puo avvenire la ricombinazione, il processo di

mutazione segue ancora il modello dei siti infiniti, secondo cui non

possono avvenire due mutazioni sullo stesso sito. Lo switch −r per-

mette di specificare il parametro r che indica il rapporto di cross-over.

Come anche in [3] r viene posto uguale a 0, 4, 16 e 40. Non si sa quale

sia il valore di r che corrisponde al reale livello di ricombinazione in

natura, ma il range usato qui dovrebbe essere sufficientemente ampio

da simulare livelli realistici di ricombinazione.

57

5.4 Analisi dei risultati

Nell’appendice A sono riassunti i risultati ottenuti per le classi di istan-

za illustrate nel capitolo precedente. Osservando tali risultati possiamo

vedere come le aspettative teoriche siano state confermate dai risultati

sperimentali.

Per i campioni di dati ottenuti a partire da 10 aplotipi l’accuratezza

e sempre del 100%, perche, come abbiamo visto nel capitolo preceden-

te, il numero di aplotipi e il minore possibile per poter generare 50

genotipi diversi. Possiamo notare, dalla cardinalita degli insiemi ot-

tenuti, che la casualita della generazione dei genotipi puo portare ad

ottenere risultati di cardinalita minore delle soluzioni iniziali, quindi

a non considerare uno o piu aplotipi. L’accuratezza, assecondando le

ipotesi teoriche precedenti, diminuisce all’aumentare del numero degli

aplotipi di partenza: quanto piu aumenta il numero degli aplotipi, fis-

sato il numero dei genotipi, piu aumenta la possibilita di avere diverse

soluzioni ugualmente ottime per la parsimonia pura o addirittura di ot-

tenere risultati parzialmente differenti dalla soluzione cercata, seppur

di cardinalita minore.

L’aumentare della ricombinazione contribuisce a diminuire la per-

centuale di accuratezza perche una maggiore ricombinazione comporta

una maggiore differenziazione degli aplotipi, quindi un aumento della

difficolta delle singole istanze.

La colonna dei ‘2’ indica per ogni istanza il numero l di siti eterozi-

goti del genotipo con piu siti eterozigoti. Ricordiamo che, per quanto

riguarda la modellazione HIc implementata, il numero delle variabili e

O(2l) e il numero dei vincoli e O(22l−1). Osservando le dimensioni dei

problemi trattati possiamo osservare il buon funzionamento nel caso

medio del preprocessamento per quanto riguarda le variabili e della

procedura di ricerca dei tagli per quanto riguarda i vincoli. Comun-

que non e possibile formalizzare tale vantaggio perche dipende dalla

singola istanza e non ha fondamenti teorici rigorosi.

Prima conseguenza delle dimensioni ridotte dei problemi consi-

58

derati e il tempo di processamento che si mantiene, per la grande

maggioranza delle istanze, inferiore ai cinque secondi.

59

Conclusioni

Lo sviluppo di una mappa degli aplotipi e un problema attuale, ancora

aperto e fondamentale per lo sviluppo delle prossime fasi dello stu-

dio del genoma umano per applicazioni in campo medico diagnostico,

chirurgico e filogenetico. Il problema dell’inferenza di aplotipi viene

studiato come problema di biologia computazionale per ottenere infor-

mazioni relative agli aplotipi secondo modelli teorici e implementazioni

basate su tecniche informatiche.

Obiettivo di questa tesi e lo studio delle capacita applicative del

modello HIc presentato per la soluzione del problema di inferenza di

aplotipi con parsimonia pura.

Dai risultati ottenuti per campioni di dati appositamente generati,

abbiamo visto come il modello, seppur di dimenzioni molto maggiori

rispetto al comunemente usato modello TIP di Gusfield, con la funzio-

ne di generazione dei tagli implementata, permette di ottenere risul-

tati limitando i tempi di soluzione e il numero dei vincoli considerati,

quindi le dimensioni dei problemi.

I limiti riscontrati sono di poco inferiori a quelli osservati nei lavori

precedenti sviluppati a partire dalla formulazione di PLI presentata

da Gusfield, in quanto non varia l’ordine di grandezza delle variabili

che cresce esponenzialmente all’aumentare della differenziazione dei

genotipi.

Le diverse classi di istanza generate e risolte hanno messo in evi-

denza come il modello della parsimonia pura sia efficace quanto piu sia

grande il rapporto tra la cardinalita dell’insieme dei genotipi e la car-

dinalita dell’insieme iniziale degli aplotipi e che, quindi, per ottenere

60

risultati con una buona accuratezza e necessario tarare le dimensioni

degli insiemi dei genotipi da risolvere secondo studi statistici riguardo

la differenziazione degli aplotipi nella popolazione in esame.

Per migliorare le prestazioni dell’algoritmo implementato, sviluppi

futuri nello studio del problema secondo il modello dato dovrebbero

studiare una finzione di pricing che permetta di ridurre il numero delle

variabili generate secondo uno studio dei costi ridotti legati all’introdu-

zione di tali variabili nel problema di PLI. Puo inoltre essere raffinata

la fase di preprocessamento per cercare di ridurre a priori variabili che

non influiranno sulla soluzione.

Maggiori informazioni relativamente all’ereditarieta e quindi alla

natura degli aplotipi possono essere considerate prendendo in esame

campioni di dati relativi non a singoli individui presi casualmente al-

l’interno di una popolazione, quanto studiando piccoli nuclei familiari

collegati da un passato evolutivo comune.

61

Ringraziamenti

Dedico questa tesi ai miei genitori, entrambi, per il materiale geneti-

co che mi hanno tramandato e per il supporto morale ed economico

durante i lunghi anni di ingegneria.

Alla professoressa Nicosia che mi ha permesso di svolgere il tirocinio

allo IASI, alla Dottoressa Bertolazzi, che ha saputo darmi la forza per

motivare i miei risultati e ad Alessandra, che probabilmente ha fatto

la scelta migliore.

Allo IASI, dal Direttore ai centri di gestione delle informazioni:

Mauro per l’analogico, i Maestri per il digitale.

A chi ogni giorno mi regala un emozione.

Agli amici, che e sempre meglio non metterli, che c’e sempre uno

che ti scordi, uno che non vorresti mettere, uno che dovresti, uno che

non puoi, uno che non lo sai, magari lo saprai, ma per ora meglio

far finta di niente. Amici, quelli abusati dalle canzonette leggere, che

quando hai bisogno ci sono e quando no non importa, che fa piacere

comunque, anche se ci pensi ogni tanto e quando potrai ricomincia

tutto. Qelli che non ci pensi, che neanche te ne rendi conto, che la

vita senza sarebbe peggio. Certo, si sopravvive, ma peggio. Quelli

nuovi, vecchi. Chi si rende conto se non va. No? Ometto.

Diciamo a chi vuole.

62

Bibliografia

[1] A. Godi (2005), An Application of Integer Linear Programming to

Haplotyping Inference by Parsimony Problem, Tesi di Dottorato

di Ricerca presso l’Universitaa degli Studi di Roma ROMA TRE.

[2] Dan Gusfield and Steven H. Orzack (2005), Haplotype Inference,

In CRC Handbook on Bioin-formatics, S. Aluru Editor.

[3] Gusfield D. (2003), Haplotype Inference by Pure Parsimo-

ny, UC Davis Computer Science Technical Report CSE-2003-2,

University of California, Davis.

[4] V. Bafna, D. Gusfield, S. Hannenhalli and S. Yooseph (2004), A

Note on Efficient Computation of Haplotypes via Perfect Phylo-

geny, Journal of Computational Biology, Oct 2004, Vol. 11, No. 5

: 858 -866.

[5] V. Bafna, D. Gusfield, G. Lancia and S. Yooseph (2003), Ha-

plotyping as Perfect Phylogeny: a Direct Approach, Journal of

Computational Biology, Jun 2003, Vol. 10, No. 3-4 : 323 -340.

[6] Lynce et al. (2006), Assessing the efficacy of haplotype inferen-

ce by pure parsimony on biological data, INESC-ID Technical

Report, December 2006, Lisboa.

[7] L. Wang and Y, Xu (2003), Haplotype Inference by Maximum

Parsimony, Bioinformatics , vol. 19, no. 14, pp. 1773-1780.

63

[8] Lancia (2004), Integer Programming Models for Computational

Biology Problems, Journal of Computer Science and Technology,

v.19 n.1, p.60-77.

[9] Lancia, Counting Variables in Gusfield’s IP Model, manoscritto.

[10] Lancia, Pinotti and Rizzi, Haplotyping Populations: Complexity

and Approximations, Technical Report DIT-02-080, Informatica

e Telecomunicazioni, University of Trento.

[11] Lancia and Rizzi, A Polynomial Case of the Parsimony Haplo-

typing Problem (2006), Operation Research Letters 34(3), pp.

289-295, Elsevier.

[12] L. van Iersel, J. Keijsper, S. Kelk and L. Stougie (2007), Sho-

relines of islands of tractability: Algorithms for parsimony and

perfect phylogeny haplotyping problems. Accepted (April 2007)

and to appear in IEEE/ACM TCBB (IEEE/ACM Transactions

on Computational Biology and Bioinformatics).

[13] D. Brown and I. Harrower (2006), Integer programming approa-

ches to haplotype inference by pure parsimony, IEEE/ACM Tran-

sactions on Computational Biology and Bioinformatics, 3(2):141-

154, April-June 2006.

[14] P. Bonizzoni, G. Della Vedova, R. Dondi and J. Li (2003), The Ha-

plotyping Problem: An Overview of Computational Models and

Solutions, Journal of Computer Science and Technology, vol.18

no.6 P.675-688.

[15] Contributori di Wikipedia, ’Meiosi’, Wikipedia, L’enciclopedia

libera, 18 aprile 2007, 10:04 UTC,

mhttp://it.wikipedia.org/w/index.php?title=Meiosi&oldid=8232529

[in data 18 aprile 2007]

64

[16] Andrew G. Clark (1990), Inference of Haplotypes from PCR-

amplified Samples of Diploid Populations, Mol. Biol. Evol, 7,

111-122.

[17] L.Excoffier and M. Slatkin, Maximum-likelihood estimation

of molecular haplotype frequencies in a diploid population,

Molecular Biology and Evolution, 12(5):921-927,1995.

[18] M. Stephens, N.J. Smith and P. Donnely, A new Statistical

Method for Haplotype Reconstruction from Population Data,

American Journal of Human Genetics, 68:978-989, 2001.

[19] Vineet Bafna, Sorin Istrail, Giuseppe Lancia and Romeo Riz-

zi (2005), Polynomial and APX-hard cases of the individual

haplotyping problem, Theoretical Computer Science, 335(1),

109-125.

[20] Bertolazzi, Godi, Labbe and Tininini (2007), Solving haplotyping

inference parsimony problem using a new basic polynomial for-

mulation, in via di pubblicazione su Computers and Mathematics

with Applications.

[21] R. R. Hudson, Generating samples under a wright-fisher neutral

model og genetic variation, Bioinformatics, 18(2):337, Febbraio

2002.

[22] Phil Hyoun Lee, Computational Haplotype Analysis: An

overview of computational methods in genetic variation study.

[23] J. Marchini, D. Cutler, N. Patterson, M. Stephens, E. Eskin, E.

Halperin, S. Lin, Z.S. Qin, H.M. Munro, G.R. Abecasis, P. Don-

nelly, and International HapMap Consortium (2006) A Compa-

rison of Phasing Algorithms for Trios and Unrelated Individuals.

Amercan Journal of Human Genetics, 78:437-450.

[24] Jonathan Marchini, Phasing Algorithm Benchmark Datasets,

http://www.stats.ox.ac.uk/ marchini/phaseoff.html, 2007.

65

[25] Stephens M, Donnelly P (2003) A comparison of Bayesian me-

thods for haplotype reconstruction from population genotype data.

Am J Hum Genet 73:1162-1169.

[26] Lin S, Cutler DJ, Zwick ME, Chakravarti A (2002) Haploty-

pe inference in random population samples, Am J Hum Genet

71:1129-1137

[27] Romeo Rizzi, Vineet Bafna, Sorin Istrail and Giuseppe Lancia

(2002), Practical Algorithms and Fixed-Parameter Tractability

for the Single Individual SNP Haplotyping problem, 2nd Work-

shop on Algorithms in Bioinformatics (WABI), Lecture Notes in

Computer Science, 2452:29-43, Springer.

[28] G. Lancia, V. Bafna, S. Istrail, R. Lippert and R. Scharz (2001),

SNPs problems, complexity and algorithms, in Proceedings of An-

nual European Symposium on Algorithms (ESA), volume 2161 of

Lecture Notes in Computer Science, pp 182.193, Springer.

66