Breve introducción a nuestro trabajo
En este trabajo buscamos homólogos de una familia de selenoproteínas en genomas recién secuenciados de una serie de protistas. Para llevarlo a cabo, hemos creado un protocolo, usando un programario y una serie de bases de datos para poder obtener y tratar la información necesaria. Finalmente conseguimos crear un programa para automatizar este proceso, pero, pese a eso, explicaremos paso a paso el procedimiento.
En la búsqueda manual de selenoproteínas en los protistas del 2011 hemos usado diversas queries. Para DSBA cogimos la query que nos proporcionaron en la web de bioinformática que era de Ciona intestinalis y para SelI cogimos la query de SelenoDB. Para SelI sólo había 3 queries: humana, de ratón y un homólogo con cisteína en drosófila. La búsqueda manual la hicimos con la query humana y después con la de ratón. Como estas queries eran de organismos demasiado lejanos filogenéticamente, hicimos un BLASTp en la base de datos del NCBI y buscamos queries de organismos más cercanos. Para la selenoproteína I también cogimos como query la selenoproteína de Phytophora infestans anotada por compañeros de años anteriores en sus respectivos trabajos y para dsbA cogimos las queries de Arthroderma gypseum y Psychromonas ingrahamii
Las queries elegidas finalmente se pueden encontrar en el apartado de introducción de esta página.
PROGRAMA AUTOMÁTICO PARA BUSCAR SELENOPROTEÍNAS
Para facilitar la búsqueda de selenoproteínas en protistas hemos creado un programa automatizado que nos permite buscar la presencia de selenoproteínas usando cualquier query de partida y en cualquier genoma. Este programa está dividido en diversas partes. El programa principal contiene el proceso automatizado hasta el Tcoffee. Posteriormente hicimos una serie de programas opcionales: uno para analizar proteínas no redundantes con el programa BLASTP.sh, otro para buscar elementos SECIS llamado SECIS.sh y un último programa llamado multiTCOFFEE.sh que permite hacer alineamientos múltiples para todas las proteínas encontradas usando cada una de las queries con la query de partida. Este último programa permite interpretar de un modo más visual los resultados ya que con él podemos ver las regiones en las que las secuencias proteicas de los diversos genomas alinean con la query.
La descripción de los programas está incluida en los archivos .sh en forma de comentarios indicando de forma detallada qué hace cada paso.
Para poder poner en marcha el programa es necesario seguir estos pasos:
- Crear una carpeta con el nombre de Semillas que contenga dentro las queries a analizar en formato FASTA (.fa)
- Tener en el mismo directorio en el que hemos guardado el programa principal los programas Evalues.pl y fastaseqfromGFF.pl.
BÚSQUEDA MANUAL DE SELENOPROTEÍNAS
Obtención del material genómico
Los genomas de los protistas secuenciados en los últimos años han sido proporcionados por los profesores de la asignatura. Están organizados en carpetas en función del año en el que han sido secuenciados. Nosotros hemos analizado manualmente los genomas que se encontraban en la carpeta de 2011 accediendo mediante la siguiente ruta:
/cursos/BI/genomes/protists/2011/Nombreprotista/genome.faJunto a esto, procederemos a exportar los diversos programas que necesitaremos para nuestro trabajo, mediante los comandos
$ export PATH=/cursos/BI/bin/ncbiblast/bin:$PATH # para el NCBI Blast $ cp /cursos/BI/bin/ncbiblast/.ncbirc ~/ # para el NCBI Blast $ export PATH=/cursos/BI/bin/exonerate/i386/bin:$PATH # para el exonerate $ export WISECONFIGDIR=/cursos/BI/soft/wise-2.2.0/wisecfg # para el GeneWise $ export PATH=/cursos/BI/bin:$PATH # para el GeneWise, el t_coffee y el fastaseqfromGFF.pl.
Basic Local Alignment Search Tool
Hacemos los tblastn de nuestra query contra todos los genomas presentes en la carpeta de protistas de 2011 usando el comando siguiente:
$ blastall -p tblastn -i query.fa -d genoma_protista.fa -o queryVSgenoma.tblastn
Para hacerlo más rápido hicimos un mini-programa para automatizar el proceso:
#!/bin/bash
export PATH=/cursos/BI/bin/ncbiblast/bin:$PATH
cp /cursos/BI/bin/ncbiblast/.ncbirc ~/
for genome in A_taiwanensis B_hominis C_merolae C_muris C_owczarzaki C_parvum E_siliculosus H_arabidopsidis N_gruberi P_pallidum P_tricornutum P_ultimum S_parasitica T_trahens ; do {
blastdb=`grep $genome /cursos/BI/genomes/protists/genomes_list_and_info.tab | cut -f 4`
blastall -p tblastn -i dsba.fa -d $blastdb -o dsbavs$genome.tblastn.fa
blastall -p tblastn -i selI.fa -d $blastdb -o selIvs$genome.tblastn
} done
Así, este programa va automáticamente a la carpeta que contiene los genomas de nuestros organismos y corta la columna 4, que es la que contiene el nombre de los protistas a analizar. Los resultados los guardamos en un archivo que llamamos selenoproteinaVSprotista.tblastn
También hicimos otro tBLASTn usando la opción -m9 y -F F:
-m 9 muestra las posiciones que ocupan los distintos alineamientos del BLAST. Esta opción nos facilita el análisis de los resultados obtenidos ya que el output es en forma de tabla.
-F F permite considerar las low complexity regions de nuestros genomas. El BLAST, por defecto, las sustituye por Xs.
En este caso el comando usado es:
$ blastall -p tblastn -i query.fa -d genoma_protista.fa -o queryVSgenoma.tblastn.fa -F F -m 9
También hicimos un programa que automatizase este paso:
#!/bin/bash
export PATH=/cursos/BI/bin/ncbiblast/bin:$PATH # pel NCBI Blast
cp /cursos/BI/bin/ncbiblast/.ncbirc ~/ # pel NCBI Blast
for genome in A_taiwanensis B_hominis C_merolae C_muris C_owczarzaki C_parvum E_siliculosus H_arabidopsidis N_gruberi P_pallidum P_tricornutum P_ultimum S_parasitica T_trahens ; do {
blastdb=`grep $genome /cursos/BI/genomes/protists/genomes_list_and_info.tab | cut -f 4`
blastall -p tblastn -i dsba.fa -d $blastdb -o dsbavs$genome.tblastn -F F -m9
blastall -p tblastn -i selI.fa -d $blastdb -o selIvs$genome.tblastn -F F -m9
} done
Una vez obtenidos los resultados de los BLAST hemos seguido con el proceso de análisis utilizando tan sólo aquellos hits que nos han dado resultados significativos, es decir, aquellos que tienen un E-value menor a 0'001. En los BLAST de Dsba veíamos como la selenocisteína de la query alineaba con una cisteína en algunos casos o con una selenocisteína. Sin embargo, en la SelI todos los BLAST finalizaban antes de la región que contenía la selenocisteína así que los siguientes pasos los hicimos para todos los hits significativos que se encontraban en contigs diferentes.
Extracción de la región genómica que potencialmente contiene el gen
Para los diversos protistas, cogemos el genoma de la carpeta del 2011 y hacemos los fastaindex usando el siguiente comando:
$ fastaindex /cursos/BI/genomes/protists/2011/Nombre.protista/genome.fa Protista.index ejemplo: $ fastaindex /cursos/BI/genomes/protists/2011/T.trahens/genome.fa T_trahens.index Este paso permite ordenar los genoma de los protistas.
Después de tener los fastaindex hacemos los fastafetch. Este paso nos permite extraer la región genómica de interés. Para hacer los fastafetch volvemos a usar la carpeta que contiene los genomas de los protistas, los fasta index que acabamos de hacer y los identificadores obtenidos en el tblastn que tenían un buen E-value. El comando usado es el siguiente:
$ fastafetch /cursos/BI/genomes/protists/2011/Nombre.protista/genome.fa Protista.index 'identificador' > Protista.fastafetch.fa ejemplo: $ fastafetch /cursos/BI/genomes/protists/2011/T.trahens/genome.fa T_trahens.index 'gi|298517488|gb|GL349486.1|' > T_trahens.fastafetch.fa
Después extraemos la región del genoma que contiene la secuencia de interés y para ello hacemos un fastasubseq en el que indicamos el inicio y la longitud de la región que queremos extraer. Para establecer el nucleótido de inicio cogimos el más pequeño de las columnas de s.start o s.end del tBLASTn y le restamos unos 2000 Nts. Para ello usamos el siguiente comando:
$ fastasubseq Protista.fastafetch.fa inicio longitud > Protista.fastasubseq.fa ejemplo: $ fastasubseq T_trahens.fastafetch.fa 21500 5000 > T_trahens.fastasubseq.fa
Anotación con exonerate o genewise
Exonerate y Genewise permiten obtener el modelo de un gen. Ambos tienen una función similar, aunque muestran ligeras diferencias en el resultado. Nosotros para este trabajo hemos recurrido especialmente a Exonerate, que se encarga de alinear la secuencia nucleotídica de la query con nuestra secuencia y predice su estructura exónica, sin mostrar sensibilidad por la direccionalidad del alineamiento (si las cadenas usadas son Forward o Reverse). Es importante recordar que para poder hacer el Exonerate primero debemos sustituir la U de la query por una X.
$ exonerate -m p2g --showtargetgff -q query.fa -t Protista.fastasubseq.fa > Protista.alineamientoquery.gff ejemplo: $ exonerate -m p2g --showtargetgff -q selI.fa -t T_trahens.fastasubseq.fa > T_trahens.alineamientoquery.gff
Genewise nos ofrece una alternativa a Exonerate. Compara una secuencia proteica con una nucleotídica de ADN, haciendo innecesario el uso de fastaseqfromGFF y fastatranslate, al obtener el cDNA y la proteína directamente. Sin embargo, a diferencia del Exonerate, requiere que se indique si las cadenas usadas son Forward o Reverse ya que el sentido de lectura cambia totalmente los tripletes, y, por tanto, los aminoácidos codificados. Para evitar tener que indicar el sentido de la cadena podemos usar la opcíon -both.
$ genewise -pep -pretty -cdna -gff -both selI.fa Protista.fastasubseq.fa > Protista.genewise.gff ejemplo: $ genewise -pep -pretty -cdna -gff -both selI.fa C.owc2.fastasubseq.fa > C.owc2.genewise.gff
Una vez obtenidos los resultados de estas alineaciones, en formato GFF, obtendremos la secuencia de cDNA correspondiente en formato fasta empleando un programa perl llamado fastaseqfromGFF.pl.
$ fastaseqfromGFF.pl Protista.fastasubseq.fa Protista.alineamientoquery.gff > Protista.cdna.fa ejemplo: $ fastaseqfromGFF.pl T_trahens.fastasubseq.fa T_trahens.alineamientoquery.gff > T_trahens.cdna.fa
Traducción de la secuencia de cDNA a proteínas
Ahora nos interesará obtener la secuencia proteica que se transcribe a partir del cDNA que hemos obtenido. Para ello usaremos el comando del fastatranslate que traduce el cDNA obtenido a proteínas usando los 6 posibles ORFs:
$ fastatranslate Protista.cdna.fa > translate.Protista.fa ejemplo: $ fastatranslate C.parvum.cdna.fa > translate.C.parvum.fa
Con los resultados obtenidos tras la introducción de éstos comandos, pasaremos a analizar el fichero, que contiene las 6 posibles proteinas, para determinar cu&aactue;l és la mejor. Guardaremos la mejor proteína en un nuevo fichero emacs que nombramos como Protista.proteina.aa.fa en el que tan solo copiaremos la proteína traducida que no contenga asteriscos ya que corresponde a la mejor secuencia traducida de los 6 ORFs posibles o, en el caso de no haber ninguna secuencia sin asteriscos, aquella que tenga menos.
$ emacs Protista.proteina.aa.fa ejemplo: $ C.parvum.proteina.aa.fa
Alineamiento global de la proteína resultante con la secuencia proteica original mediante el Tcoffee
Para hacer el alineamiento de la query con la proteína resultante usamos la secuencia traducida a proteínas que no contiene asteriscos. Usaremos el programa T-Coffee (Tree-based Consistency Objective Function for alignment Evaluation), programa que nos permite alinear progresivamente todas las parejas mientras cada paso de la alineación múltiple progresiva se va realizando, obteniendo como resultado un alineamiento simple y veloz pero con poca tendencia a cometer errores.
Para usar el Tcoffee con nuestros datos usamos el siguiente comando:$ t_coffee query.fa Protista.proteina.aa.fa > Protista.tcoffee.txt ejemplo: $ t_coffee selI.fa T_trahens.proteina.aa.fa > T_trahens.tcoffee.txt
BLAST recíproco en el NCBI
Para comprobar si con la proteína que hemos obtenido encontramos la selenoproteína de la query de partida hacemos un BLASTp en el NCBI usando como query el archivo Protista.proteina.aa.fa. De esta manera, además, podremos ver si las proteínas tienen homólogos en las bases de datos del NCBI.
Búsqueda de elementos SECIS
Hemos realizado la búsqueda de elementos SECIS mediante cortes de genoma con fastasubseq, de gran tamaño para intentar coger el SECIS dentro, y usando fastarevcomp si el alineamiento era en una cadena reverse. Posteriormente, usaremos el programa SECISearch para encontrar las secuencias SECIS en nuestros cortes de genoma, y crear una imagen de las mismas. Nos quedaremos con aquel resultado con un mayor cove score.