 |  |  |
Hem treballat amb dos tipus de fitxers, un que conté els SNPs i l'altre que conté les anotacions de tots els gens del genoma per H.Sapiens.
Aquests fitxers els hem extret dels següents ubicacions:
/disc8/genomes/H.sapiens/golden_path_200405hg17/database/refGene.txt /disc8/genomes/H.sapiens/golden_path_200405hg17/database/snp.txt
Els fitxers d'anotacions dels gens (refGene.txt) tenen el següent format:
| identificador del gen | cromosoma | strand | inici transcripció | final transcripció | inici traducció | final traducció | num. exons | inici exons | final exons |
| chr1 | NM_000016 | + | 75902302 | 75941204 | 75902493 | 75940469 | 12 | 75902302,75906106,75910349,... | 75902523,75906194,75910447,... |
Per tal de poder treballar fàcilment amb aquest fitxer i enregistrar-lo en el programa, hem fet les següents modificacions amb comandes de Unix:
$ sort refGene.txt | uniq > ref1.txt
Per tal d'eliminar els gens que estan repetits.
$ gawk '{print $2" "$1" "$3" "$4" "$5" "$6" "$7" "$8" "$9" "$10}' ref1.txt > genoma.txt
D'aquesta forma ordenem cadascuna de les columnes del fitxer per tal de tenir el cromosoma a la primera columna.
$ egrep "chr1 [^0-9]" genoma.txt > chr1.txt
El mateix per a cadascun dels cromosomes. Per chrX i chrY:
$ egrep "chrX" genoma.txt > chr23.txt i $ egrep "chrY" genoma.txt > chr24.txt
D'aquesta manera tenim els gens ordenats per cromosoma i un fitxer per a cada cromosoma.
Una vegada tenim els fitxers dels gens de cada cromosoma ens interessa que aquests estiguin ordenats per inici de transcripció per tal de poder dur a terme una subfunció de cerca binària que necessitem dins el programa. Amb la següent comanda aconseguim aquest propòsit:
$ gawk '{print $4" "$1" "$2" "$3" "$5" "$6" "$7" "$8" "$9" "$10}'chr1.txt > chr1a.txt
$ sort chr1a.txt > chr1ord.txt | rm chr1.txt
$ gawk '{print $2" "$3" "$4" "$1" "$5" "$6" "$7" "$8" "$9" "$10}'chr1ord.txt > chr1.txt
Ara ja tenim els fitxers d'anotacions modificats per tal d'enregistrar-los en el nostre programa.
El fitxer de SNPs és massa gran per enregistrar-lo en el disc local en el qual estem treballant perquè ocupa 920Mb. Per aquest motiu hem decidit no modificar-lo i cridar-lo al mateix temps que executem el programa des de la seva ubicació original.
Per tal d'estudiar la distribució dels SNPs en els gens, hem creat un programa que ens permet recórrer el fitxer dels SNPs i buscar la posició on cau el SNP dins el fitxer dels gens.
El primer pas que hem realitzar ha estat enregistrar els fitxers dels cromosomes dins una matriu que hem anomenat @genoma. Aquesta matriu és de 3 dimensions. La primera dimensió conté el cromosoma, la segona els gens i la tercera cadascuna de les propietats que ens interessen de cada gen: identificador, inci transcripció, final transcripció, inici traducció, final traducció, número exons, i un vector anomenat @ife que conté les posicions d'inici i final de cadascun dels exons.
A continuació, obrim el fitxer dels SNPs i l'anem llegint línia per línia. Amb cadascuna de les línies que llegim, buscar la posició del SNP dins la matriu del genoma a través d'una subfunció que hem anomenat buscagens, i que consisteix en un algorisme de cerca binària, reduïnt a cada volta la finestra on buscar la posició del SNP.
Al final ens dóna la posició dins la matriu del gen on cau el SNP.
Finalment hem buscat en quina part del gen cau el SNP (intró, lloc donador de splicing, lloc acceptor de splicing, CDS, 5'UTR o 3'UTR). Això ho fem mitjançant una sèrie de composicions alternatives que defineixen cadascuna de les parts del gen.
El nostre OUTPUT l'hem redireccionat a un fitxer anomenat resultats.out que ens mostra el % de SNPs que cauen en cada regió definida.
Amb l'objectiu d'estudiar la quantitat de SNPs que conté cada cromosoma, hem creat un segon programa que es basa en:
Obrir el fitxer de SNPS, llegir línia per línia i enregistrar el número d'SNPs que cauen en cada cromosoma. L'OUTPUT del programa es un fitxer "contachr.out" que conté el total de SNPS de cada cromosoma.
Estudi comparatiu
Diferències en la distribució dels SNPs en gens relacionats i no relacionats amb malalties
PART 1: OBTENCIÓ DELS FITXERS
Necessitem, a partit del fitxer d'annotacions dels gens (RefGene.txt), obtenir dos fitxers que separin els gens entre aquells que tenen una malaltia asociada i els que no estan relacionats amb malaltia.
Per tal d'obtenir aquests fitxers hem realitzat els següents passos:
- Entrem a la pàgina d'Ensembl, i dins d'aquí a l'enllaç EnsMart.
- Seleccionem Ensembl Genes i H.Sapiens i cliquem 'Next'.
- Com que volem recuperar tot el genoma, desactivem la opció 'Limit to' i cliquem 'Next'.
- Dins del quadre de 'Regions' seleccionem 'Chromosome Name','Start Position','End Position' i 'Strand'.
- Dins del quadre 'Gene' activem les opcions 'Ensembl Gene ID','RefSeq DNA ID','Disease OMIM ID' i 'Disease Description'.
- Pel que fa a 'output format', seleccionem 'Text,tab separated', i 'gzip'.
- Al fitxer de sortida li diem hsapensembl.tsv.
A partir d'aquest fitxer i del RefGene.txt, realitzem les següents comandes amb Unix per tal d'obtenir els dos fitxer amb els gens relacionats i no relacionats amb malalties:
Ordenem el RefGene.txt per la columna de l'identificador.
$ sort -k 1,1 refGene.txt > refGene.orderedbyid.txt
Recuperem els gens que tenen entrades de OMIM. Tallem 13 columnes, i només si n'obtenim més de 11 vol dir que el gen té les entrades de OMIM.
$ cut -f 1,2,3,6,7,8,9,13,14,15,17,18,19 hsapmmusensmbl.tsv | perl -nae 'if (scalar(@F)>11 && $F[3]=~/\ANM/){print}' | cut -f 4 | perl -ne '/\A(\w+)/; print "$1\n"' | sort | uniq | wc -l
Hi ha mé annotacions quan ho ajuntem amb el RefSeq degut als transcrits alternatius. Amb la següent opció forcem que les columnes de l'output estiguin separades per tabulacions.
$ cut -f 1,2,3,6,7,8,9,13,14,15,17,18,19 hsapmmusensmbl.tsv | perl -nae 'if (scalar(@F)>11 && $F[3]=~/\ANM/){print}' | cut -f 4 | perl -ne '/\A(\w+)/; print "$1\n"' | sort | uniq | join -t "`echo -e '\t'`" - refGene.orderedbyid.txt | wc -l
Finalment extraiem els gens amb annotacions OMIM i sense annotacions OMIM i ho redireccionem a dos fitxers.
$ cut -f 1,2,3,6,7,8,9,13,14,15,17,18,19 hsapmmusensmbl.tsv | perl -nae 'if (scalar(@F)>11 && $F[3]=~/\ANM/){print}' | cut -f 4,11 | perl -nae '/(\w+)\.\d+\s+([^\n]+)/;if ($2 ne ""){print "$1\t$2\n"}' | sort | uniq > refseqwithomimdnds.txt
$ cut -f 1,2,3,6,7,8,9,13,14,15,17,18,19 hsapmmusensmbl.tsv | perl -nae 'if (scalar(@F)<=11 && $F[3]=~/\ANM/){print}' | cut -f 4,11 | perl -nae '/(\w+)\.\d+\s+([^\n]+)/;if ($2 ne ""){print "$1\t$2\n"}' | sort | uniq > refseqwithoutomimdnds.txt
Finalment dividim els fitxers per cromosoma per tal de poder enregistrar-los després en una matriu:
Per a gens relacionats amb malalties.
$ egrep "chr1 [^0-9]" genoma.txt > rchr1.txt
El mateix per a cadascun dels cromosomes. Per chrX i chrY:
$ egrep "chrX" genoma.txt > rchr23.txt i $ egrep "chrY" genoma.txt > rchr24.txt
I per als gens no relacionats amb malalties.
$ egrep "chr1 [^0-9]" genoma.txt > nchr1.txt
El mateix per a cadascun dels cromosomes. Per chrX i chrY:
$ egrep "chrX" genoma.txt > nchr23.txt i $ egrep "chrY" genoma.txt > nchr24.txt
Ja tenim els fitxers llestos per enregistrar-los al nostre programa.
Utilitzarem el mateix programa que en el cas anterior, però el nom d'entrada dels fitxers per enregistrar la matriu serà diferent:
rchr$n.txt pels gens relacionats amb malalties.
nchr$n.txt pels gens no relacionats amb malalties.
| | |