Materials i Mètodes

L' objectiu del nostre treball és identificar les selenoproteïnes DI, SelJ i SelQ en el genoma de diferents protistes. El protocol del treball s'ha estructurat en diferents apartats que es corresponen amb el programari utilitzat.

Automatització

Per agilitzar tot el procés d'anàlisi i obtenció de dades, hem creat diversos programes. Un anomenat Blast.bash permet fer automàticament tots el BLAST d'un conjunt de querys contra una llista predefinida de genomes. Un altre anomenat AMRAgenomeannotation.bash incorpora tot el procés d'identificació de hits, predicció d'exons, t-coffee i BLASTp contra el conjunt no redundant de totes les proteïnes disponibles en NCBI, minimitzant l'interacció de l'usuari als punts crítics i accelerant notòriament el procés d'obtenció de dades.

Per a més informació, descàrregar els programes i tutorials clica aquí.

Seqüències i genomes

Obtenció de seqüències o querys

Les seqüències de SelJ i SelQ les hem obtingut de l'arxiu TARBALL, mentre que DI les hem obtingut a la base de dades SelenoDB.
SelJ és molt similar a una familia de proteïnes anomenada CristallinJ1, que trobem, per exemple, a L.loa a la base de dades de NCBI, i que també hem analitzat.
SelQ de T.gondii de l'arxiu TARBALL, no ens va donar cap hit. Vam buscar homólegs en altres organismes mitjançant NCBI BLAST i vam trobar un bon hit amb e-value 9e-14 en el cromosoma VIII de Neospora caninum Liverpool. Vam descarregar el genoma i amb aquest protocol vam aïllar una possible query SelQ de N.caninum. Aquí teniu el .gff i el t-coffee, així com la seqüència: 
>gi|325117666|emb|FR823390.1|:subseq(1685381,20000) Neospora caninum Liverpool complete genome, chromosome VIII:[translate(1)]
MKDFKKTEEFKELEKEAQDREKGVQPEHRRTWRFPGNSPLNPHFAPTYRPNVNDRYQIRRDRGGUC

DI és una familia de Iodothyronine deiodinases formada per tres membres en Homo sapiens: DI_1, DI_2 i DI_3, que són els que hem utilitzat en aquesta cerca. D'altre banda, al llarg del procés d'annotació de genomes hem vist que aquestes seqüències eren molt llunyanes als nostres organismes i també n'hem agafat dues més properes: DI de D.discoideum i DI1 de D.rerio (zebrafish) de la base de dades de NCBI. [7]

Obtenció de genomes

Els genomes de protistes que hem utilitzat els trobem al directori: 

/cursos/BI/genomes/protists/2012
D'altre banda hem descarregat el genoma del cromosoma VIII de Neospora caninum Liverpool de NCBI, per tal d'obtenir una nova query SelQ. Per a poder fer un BLAST d'aquest genoma hem hagut de fer una base de dades gràcies a la següent comanda: 
$ formatdb -i ./N.caninum_genome.fa -p F -n Neosdb.fa

tBLASTn

En primer lloc hem fet una cerca mitjançant tBLASTn de cadascuna de les selenoproteïnes al genoma de tots els protistes. tBLASTn és un programa que produeix alineaments locals de seqüències. En concret, permet comparar una seqüència proteica (query) amb una base de dades de nucelòtids.

Per tal d'automatitzar el procés hem emprat l'script blast.bash. Només cal donar permís d'execució amb chmod u+x i executar el programa a la terminal amb ./blast.bash. Més informació aquí. 

#!/bin/bash
export PATH=/cursos/BI/bin/ncbiblast/bin:$PATH 
cp /cursos/BI/bin/ncbiblast/.ncbirc ~/ 

for query_fa in ./querys/*; do {
query=`basename $query_fa`
nomproteina=${query%.fa}
llistaquerys=$llistaquerys" $nomproteina"
} done

mkdir ./blasts
for nomproteina in $llistaquerys; do {
for genome in A.laibachii_Nc14 A.rara C.fasciculata D.discoideum_AX4 D.fasciculatum F.cylindrus G.niphandrodes I.multifiliis_strain_G5 L.donovani_BPK282A1 L.tarentolae P.capsici P.polycephalum S.arctica T.congolense; do {
echo "fent $nomproteina en $genome"
blastall -p tblastn -i ./querys/$nomproteina.fa -d /cursos/BI/genomes/protists/2012/$genome/genome.fa -o blasts/$nomproteina.$genome.txt
} done
echo fet
} done

Un cop hem obtingut tots els blast, els hem mirat manualment i hem considerat oportu continuar la búsqueda de les querys en aquells genomes on obteniem un hit amb un e-score de 1e-4 o inferior. Considerem aquest llindar com a significatiu, ja que experimentalment el mínim score amb el que hem predit una selenoproteina era de 1e-4.

Regió genòmica del hit

Per tal de poder extreure primer la seqüència que conté la regió, i després extreure la regió a partir de la seqüència, necessitem l'ajuda d'aquests tres programes que acompanyen a exonerate. $year fa referència a l'any (e.g. 2012) mentre que $organismo fa referència al organisme que estem mirant (e.g. A.rara)

Fastaindex

Com hem automatitzat el procés, i alguns genomes al estar distribuits en forma de contigs ens donaven error al realitzar el fastasubseq, hem aplicat fastaindex a tots els genomes, amb la següent comanda: 
 export PATH=/cursos/BI/bin/exonerate/i386/bin:$PATH
fastaindex /cursos/BI/genomes/protists/$year/$organismo/genome.fa ./$organismo.index

Fastafetch

Seguidament, hem utilitzat el programa Fastafetch per tal d'extreure el contig que conté el hit d'interès amb la següent comanda, sent $hit la regió o contig on hem trobat el hit amb el BLAST: 
fastafetch /cursos/BI/genomes/protists/$year/$organismo/genome.fa $organismo.index $hit > ./$query/$organismo/database/$organismo.$query.$hit.db

Fastasubseq

Un cop hem extret el contig, hem utilitzat el programa Fastasubseq per a extreure el fragment que hem estimat que conté el nostre gen. Hem afegit unes 10kb a cada costat del hit, per assegurar-nos que agafàvem el gen d’interès. $inicio fa referència a 10kb abans d'on hem trobat el hit i $longitud serà 20000 (10kb per cada costat): 
fastasubseq ./$query/$organismo/database/$organismo.$query.$hit.db $inicio $longitud > ./$query/$organismo/fastasubseq/$query.$organismo.genome.fa

Predicció d'exons

Exonerate i Genewise són dos programes per a la predicció de gens. Presenten algorismes diferents i, per tant, pot ser que on un falli al fer l'anotació l'altre ho anoti correctament.

Exonerate

 Hem generat les anotacions dels gens que estem buscant amb exonerate, alineant el fragment de DNA que hem extret amb Fastasubseq amb la seqüència de DNA de la proteïna inicial. En aquest punt hem de modificar la nostra query i canviar la U de la selenocisteina per una X, i eliminar qualsevol altre símbol com per exemple una @, ja que d'altre banda el exonerate mostrarà un error. També és important que fem servir la següent comanda amb l'opció --exhaustive yes per tal de fer una cerca més acurada. 
export PATH=/cursos/BI/bin/exonerate/i386/bin:$PATH
exonerate -m p2g --exhaustive yes --showtargetgff -q ./query/$query -t ./$query/$organismo/fastasubseq/$query.$organismo.genome.fa > ./$query/$organismo/exonerate/$query.$organismo.exonerate.gff
egrep -w exon ./$query/$organismo/exonerate/$query.$organismo.exonerate.gff  > $query.$organismo.exo.gff
Amb aquesta comanda guardem el resultat del exonerate en .gff, alhora que ens deixa un altre arxiu preparat per al processament i obtenció de la seqüència petídica.

FastaseqfromGFF

Posteriorment, hem utilitzat el programa fastaseqfromGFF.pl amb el fitxer GFF obtingut amb l’exonerate per tal d’extreure la seqüència de cDNA 
export PATH=/cursos/BI/bin:$PATH
chmod u+x ./fastaseqfromGFF.pl
./fastaseqfromGFF.pl ./$query/$organismo/fastasubseq/$query.$organismo.genome.fa $query.$organismo.exo.gff > ./$query/$organismo/cDNA/$query.$organismo.cDNA

Fastatranslate

A continuació hem traduït el cDNA a proteïna amb el programa Fastatranslate,que acompanya a l'exonerate. 
fastatranslate ./$query/$organismo/cDNA/$query.$organismo.cDNA > ./$query/$organismo/proteina/$query.$organismo.fa
Els arxius que ens ha generat el Fastatranslate contenen 6 seqüències fasta, corresponents a les 6 possibles pautes de lectura (3 directes i 3 reverse). Per a continuar, hem de seleccionar una de les seqüències fasta i copiar-la en un arxiu nou. Una de les formes per escollir la pauta de lectura adequada és fixar-nos en el .gff que ha creat l'exonerate i agafar la mateixa pauta que es mostra en el .gff. Un cop hem seleccionat la pauta és interessant modificar els possibles * per X de cara a facilitar la interpretació de l'alineament amb t-coffee.

Si l'exonerate no produeix l'anotació sempre podem revisar les comandes o utilitzar genewise. A vegades, aconsegueix fer una anotació incomplerta, una de les causes podria ser la mida del fragment de DNA del Fastasubseq, però no esperem que sigui això un problema ja que hem agafat 20kb. Una altra forma d'intentar solucionar-ho seria manipular el .gff i modificar manualment la posició dels exons predits, agafant unes quantes pb abans, potser d'aquesta manera es soluciona el problema tot i que a nosaltres no ens ha funcionat.

Genewise

Amb el programa genewise generem una nova anotació del gen. La seva funció és la mateixa que l’exonerate, però com que no fa servir el mateix algoritme els resultats poden variar. Farem servir la següent comanda, tenint en compte que per anotacions en la cadena reversa hem d'afegir l'instrucció -trev: 
export PATH=/cursos/BI/bin:$PATH
export WISECONFIGDIR=/cursos/BI/soft/wise-2.2.0/wisecfg
genewise -pep -pretty -cdna -gff ./query/$query ./$query/$organismo/fastasubseq/$query.$organismo.genome.fa > ./$query/$organismo/genewise/$query.$organismo.genewise.gff
Genewise ens genera de forma automàtica la seqüència peptídica. Per a comparar aquesta seqüència amb la query amb t-coffee haurem de copiar-la a un nou arxiu manualment, o podem fer servir el programa GFFtoFA.pl que hem creat amb l'ajuda de Robert Castelo, per a fer-ho automàticament: 
chmod u+x ./GFFtoFA.pl
./GFFtoFA.pl < ./$query/$organismo/genewise/$query.$organismo.genewise.gff > ./$query/$organismo/genewise_prot/$query.$organismo.genewise.fa

t-coffee

Una vegada hem obtingut les seqüències de les proteïnes hem utilitzat el t-coffee per tal d'alinear la proteïna obtinguda amb la query inicial. D'aquesta manera podem veure si els residus d'aminoàcids de la proteïna inicial i la predita es corresponen, i fins i tot inferir homologia: 
t_coffee ./query/$query ./$query/$organismo/proteina/$query.$organismo.fa

BLASTp contra el conjunt no redundant de totes les proteïnes disponibles a NCBI.

Podem comprovar si la proteina que hem obtingut té homòlegs en altres espècies mitjançant una cerca BLASTP contra el conjunt no-redundant de totes les proteines disponibles a NCBI utilitzant el software netblast de la següent forma: 
export PATH=/cursos/BI/bin/netblast/bin:$PATH
blastcl3 -p blastp -i  ./$query/$organismo/proteina/$query.$organismo.fa -d nr > ./$query/$organismo/ncbi/$query.$organismo.ncbi

Aquesta eina és molt potent ja que si trobem selenoproteïnes homòlogues a la nostra seqüència predita en altres espècies és molt indicatiu de que veritablement és una selenoproteina, sempre i quan tingui l'aminoàcid selenocisteïna.

Alineaments múltiples i visualització

Un cop hem obtingut totes les seqüències, de cara a presentar certs resultats a la discussió hem fet uns alineaments múltiples mitjançant t-coffee. Per tal de visualitzar millor els alineaments realitzats hem utilitzat el programa Jalview.

SECISearch

 La cerca d’elements SECIS la fem a través de la plana del software SECISearch. Només hem cercat SECIS en els organismes en els que hem trobat selenoproteïnes, amb l'opció Pattern: Loose (canonical and non-canonical). [8]

Cerca de maquinària de síntesi de selenoproteïnes

 Un cop analitzats els genomes dels organismes, hem cercat si aquests organismes presenten la màquinaria de síntesis de selenoproteïnes per tal de corroborar si en efecte els organismes estudiats tenen o no selenoproteïnes. Aquesta màquinaria és necessària per sintetitzar les selenoproteïnes en els genomes que en contenen, pero es desconeix si poden tenir altres funcions, de forma que encara que l'organisme no tingui selenoproteïnes pot ser que algun membre de la maquinaria estigui conservat. S'han cercat les següents proteïnes amb diferents querys en les següents bases de dades:

SelenoDB

  1. eEFSec: C.elegans, D.melanogaster, H.sapiens.
  2. sbp2: D.melanogaster, H.sapiens, M.musculus.

NCBI

  1. secS: C.brenneri, D.discoideum_AX4, H.sapiens, Micromonas, M.musculus.
  2. pstk: Bacteria, T.cruzi, H.sapiens.
  3. secp43: H.sapiens, M.musculus.
  4. sps2: A.gambiae, D.melanogaster, H.sapiens, M.musculus.

Amb aquestes querys hem utilitzat el protocol esmentat previament i el programa AMRAgenomeannotation.bash, per tal de veure si la maquinària de síntesi de selenoproteïnes es present als organismes 2012, de la mateixa manera que van cercar la presència o no de selenoproteïnes.

tRNAscan-SE

Per a completar la cerca maquinària de síntesi de selenoproteïnes, hem buscat també si hi ha tRNAs de selenocisteïna en els organismes estudiats. Aquesta cerca s'ha basat en el programa tRNAscan-SE 1.21, un programa que escaneja el genomes i prediu els tRNAs. L'output d'aquest programa proporciona informació de tots els tRNAs presents en el genoma d'un organisme i entre d'altres ens monstra el seu anticodó. Sabent que el codó de la selenocisteïna és TGA, buscarem el seu anticodó TCA en cadascun dels fitxers que ens ha proporcionat el programa. La presència d'aquest anticodó indica l'existència del tRNA per a selenocisteïna (SeCys), si es tracta d'un codó STOP no genera un anticodó tRNA ja que no és transcriu.

En primer lloc hem descarregat el programa i l'hem instalat. Després amb la següent comanda hem executat el programa i hem scannejat els genomes:tRNAscan-SE [-options] genome.fa, i hem mirat si hi ha o no anticodó TCA. Per tal d'accelerar el procés d'obtenció de de dades hem creat el següent script .bash, més informació aquí: 

#!/bin/bash
mkdir ~/tRNAscan
               for genome in A.laibachii_Nc14 A.rara C.fasciculata D.discoideum_AX4 D.fasciculatum F.cylindrus G.niphandrodes I.multifiliis_strain_G5 L.donovani_BPK282A1 L.tarentolae P.capsici P.polycephalum S.arctica T.congolense; do {
echo "fent $genome"                     
./tRNAscan-SE /cursos/BI/genomes/protists/2012/$genome/genome.fa > ~/tRNAscan/$genome.out
                                      } done
 for genome2 in T.cruzi T.parva P.sojae E.histolytica M.brevicollis; do {
echo "fent $genome2"                     
./tRNAscan-SE /cursos/BI/genomes/protists/2010/$genome2/genome.fa > ~/tRNAscan/$genome2.out 
                                      } done
for genome3 in D.purpureum L.mexicana; do {
echo "fent $genome3"                     
./tRNAscan-SE /cursos/BI/genomes/protists/2009/$genome3/genome.fa > ~/tRNAscan/$genome3.out 
                                      } done
for genome4 in P.tetraurelia; do {
echo "fent $genome4"                     
./tRNAscan-SE /cursos/BI/genomes/protists/2008/$genome4/genome.fa > ~/tRNAscan/$genome4.out 
                                      } done
echo fet

Selenoprofiles

 Per a verificar si el nostre anàlisi ha estat correcte hem emprat el programa Selenoprofiles. Aquest programa l'ha fet Marco Mariotti del CRG, podeu visitar la plana web http://big.crg.cat/services/selenoprofiles per a més informació. El programa fa psi-tblastn, exonerate, genewise i diferents filtrats per a tal d'annotar i predir selenoproteïnes als genomes. L'instal·lació del programa al UNIX de l'aula d'informàtica pot ser difícil, ja que necessita altres programes previs i més disc dur del disponible per a la base de dades de NCBI. Intentarem ara explicar els passos bàsics per a instal·lar-ho als ordinadors de l'aula d'informàtica.

En primer lloc descarreguem el programa i el descomprimim. Per a la seva instal·lació necesitem a més descarregar el software MAFFT, el descomprimim i l'instal·lem seguint les instruccions que trobem a la web del software (Recordeu que no som usuaris root i haurem de modificar l'arxiu Makefile cambiant PREFIX = /usr/local per PREFIX = /home/your_home/somewhere i BINDIR = $(PREFIX)/bin per BINDIR = /home/your_home/bin). D'altra banda també és necessari descarregar la base de dades de NCBI, aquest arxiu ocupa uns ~3Gb i trigarà unes 3-4 hores. Recomanem que poseu l'arxiu a un disc dur que tingui una capacitat per a uns 20-30Gb. Després de descarregar l'arxiu, l'heu de descomprimir i l'arxiu final ocuparà uns 9,6Gb, d'altra banda a vegades es necessari aplicar un -formatdb a l'arxiu descomprimit i axiò genera 10Gb extres.

Un cop tenim la base de dades de NCBI descomprimida i hem instalat el MAFFT podem passar a instalar el selenoprofiles. Anem la carpeta de selenoprofiles_2_installation, abans d'instal·lar el programa hem d'assegurar-nos que estiguin disponibles tots els programes necessaris per a la seva execució al bash terminal. Per això hem de posar al terminal les següents comandes: 

        $ export PATH=/cursos/BI/bin/ncbiblast/bin:$PATH             # pel NCBI Blast
        $ cp /cursos/BI/bin/ncbiblast/.ncbirc ~/                     # pel NCBI Blast 
        $ export PATH=/cursos/BI/bin/exonerate/i386/bin:$PATH        # per l'exonerate
        $ export PATH=/cursos/BI/bin:$PATH                           # pel GeneWise, el fastaseqfromGFF.pl i el t_coffee
        $ export WISECONFIGDIR=/cursos/BI/soft/wise-2.2.0/wisecfg    # pel GeneWise
        $PATH=PATH:~/your_home/bin                                   # pel MAFFT

Un cop hem aplicat aquestes comandes, instal·lem Selenoprofiles amb la següent comanda: 

$ ./install_selenoprofiles.py -no_nr -no_go

Un cop instal·lat, a la carpeta selenoprofiles_2_installation s'haura creat el fitxer Selenoprofiles. Per a executar-lo emprem la següent comanda (fem l'exemple de D.discoideum_AX4): 

./Selenoprofiles results_folder /cursos/BI/genomes/protists/2012/D.discoideum_AX4/genome.fa -species "Dictyostelium_discoideum_AX4"

Abans però, ens pot interessar definir quines selenoproteïnes volem que busqui, en aquest punt modificarem l'arxiu selenoprofiles_2.config que es troba a la mateixa carpeta selenoprofiles_2_installation, fent els següents canvis. 

DE
profile=machinery,eukaryotic
A
profile=grupbioinfo
...
Al mateix arxiu hem d'escribir la linea:
 
families_set.grupbioinfo = SelQ,SelJ,DI                     #aqui definiu quines son les vostres selenoproteines

Després d'executar Selenoprofiles, el programa ens crearà diverses carpetes per a cada organisme amb els arxius de la predicció que ha fet. Dins la carpeta output trobarem els resultats de l'anàlisi. Per a mirar els nostres resultats feu clic a discussió.