Queries

El objetivo es encontrar si los genomas de estos protistas contienen o no las Selenoproteínas T, V y W.Para esto, usaremos las secuencias de aminoácidos de estas Selenoproteínas en otras especies, y las alinearemos con los genomas de nuestros protistas.Las secuencias de otras especies usadas como modelo para el alineamiento las llamaremos “Query”.

El reino de los protistas es tan variado que a menudo, la distancia filogenética entre especies de protistas es tan grande como la distancia entre protistas y vertebrados. Es por esto que las Queries que usaremos provienen tanto de otros protistas como de vertebrados e invertebrados.

Antes que nada, cambiamos las U de estas secuencias query por una X.

Las Queries están almacenadas en ficheros en formato FASTA.

Sel T Sel V Sel W
Homo sapiens Homo sapiens Homo sapiens (SelW-1)
Cryptosporidium hominis (1) Danio rerio (SelW-1)
Aureococcus anophagefferens Artemia franciscana (SelW-2)
Leishmania baziliensis
Chlamydomonas reinhardtii
Drosophila melanogaster

Para asegurar que las queries que se van a utilizar son selenoproteínas T u homólogas de selT con cisteína, se realiza un alineamiento múltiple mediante el programa T_coffee. Se comprueba que efectivamente todas son SelT u homólogas, ya que se alinean todas las X (representan Selenocisteína) y C entre ellas (Archivo).

TBLASTN para alineamiento de secuencias

TBLASTN sirve para comparar nuestra proteína (query) vs una secuencia de nucleótidos (el genoma de los protistas). Por tanto lo usaremos para localizar la región donde se encuentra el posible gen de la selenoproteína.

Necesitaremos el software NCBI BLAST, y lo usaremos desde el shell del OS Unix con las siguientes comandas:

export PATH=/cursos/BI/bin/ncbiblast/bin:$PATH

cp /cursos/BI/bin/ncbiblast/.ncbirc ~/

blastall  -p  tblastn  -i  fichero_query.fa  -d  nombre_base_datos_BLAST  -o  fichero_salida –F F

  • Fichero_query.fa: nombre del fichero donde tenemos almacenada la query.
  • Nombre_base_datos_BLAST: en nuestro caso, la base de datos la encontramos en el siguiente enlace: /cursos/BI/genomes/protists/genomes_list_and_info.tab
  • Fichero_salida: nombre del fichero donde queremos que se nos almacene el resultado de la búsqueda en BLAST.

La comanda –F F sirve para autorizar las “Low Complexity regions” en nuestro alineamiento.

Existe una herramienta llamada “-m 9” que sirve para visualizar la información del documento BLAST de forma esquemàtica y tabulada. La comanda para el uso de la herramienta “-m 9” es:

blastall  -p  tblastn  -i  fitxer_query.fa  -d  nom_base_dades_BLAST  -o  fitxer_sortida –F F –m 9

Así pues, obtenemos un documento BLAST que nos alinea cada query con el genoma de uno de los protistas.

Al inicio de cada documento BLAST tenemos la lista de todos los Scaffolds que se han encontrado con esa query.

Se tiene que mirar el nombre del Scaffold y no el nombre del cromosoma, porque puede ser que el genoma de estos protistas no esté bien ensamblado.

Comprobamos en cada protista si hay algún Scaffold determinado que se ha alineado con todas o casi todas las Queries.

Si alguno de los Scaffolds lo encontramos en varias Queries: Miramos que en todos los alineamientos con las Queries tengan la misma pauta de lectura (Frame). Si es así, luego comprobamos si hay posiciones solapadas entre los fragmentos alineados. Entonces seleccionamos un fragmento de genoma de una longitud suficiente para que englobe el alineamiento de todas las Queries.

También seleccionaremos los Scaffolds que sólo se encuentren en una de las Queries pero que tengan un E-value bueno (criterio).

Nosotros buscamos un homólogo de la Selenoproteína query en el genoma del protista, y por tanto, nuestro hit tiene que contener la U alineada con:

  • U.
  • C (cisteína).
  • * (stop codon).

Tenemos que añadir unas 10.000Kb (o más) por lado, en cada Scaffold seleccionado. De esta forma nos aseguramos que abarcamos el gen entero.

Exonerate: Extracción de la región genómica

Ahora usaremos un programa que forma parte de un software llamado Exonerate.

Primero hay que conseguir que la base de datos del genoma de protista esté indexada según los Scaffolds. Para esto introducimos usamos la herramienta Fastaindex con las comandas siguientes:

formatdb -i ruta_fichero_base_datos.fa -p F -n nuevo_nombre_base_datos.fa

La comanda - p F significa que desactivamos la función de péptidos, y por tanto, sólo tenemos activada la función de nucleótidos.

fastaindex ruta_fichero_base_datos.fa nuevo_nombre_base_datos.index

Después usamos la función Fastafetch para escoger el Scaffold que nos interesa:

fastafetch ruta_fichero_base_datos.fa nuevo_nombre_base_datos.index nº_scaffold_deseado > fichero_nº_scaffold_deseado.fa

De esta forma conseguimos que el fichero de la base de datos contenga una sola secuencia.

Ahora podemos usar la herramienta Fastasubseq, que nos permitirá extraer el fragmento de genoma que queremos y almacenarlo en un fichero de formato FASTA.

Introducimos las siguientes comandas:

export PATH=/cursos/BI/bin/exonerate/i386/bin:$PATH

fastasubseq  nombre_base_datos_BLAST  start length  >  fichero_corte.fa

  • nombre_base_datos_BLAST: es el nombre de la base de datos en donde hemos hecho la búsqueda anterior.
  • Start: posición de inicio del corte.
  • Length: longitud del corte.
  • fichero_corte.fa: nombre del fichero en el que queremos almacenar el resultado del corte.

Ahora usaremos los programas Exonerate y fastaseqfromGFF.pl para obtener la sección correcta de genoma que contiene nuestro corte sin intrones. Las comandas son las siguientes:

export PATH=/cursos/BI/bin:$PATH

Realizaremos la acción en dos pasos. En el primer paso obtenemos un fichero con la anotación en GFF:

exonerate -m p2g --showtargetgff -q fichero_query.fa -t fichero_corte.fa –exhaustive yes| egrep –w exon > fichero_resultado_exonerate.gff

  • - m p2g: significa “model protein vs genome”.
  • q: documento query.
  • t: documento Target.
  • -- exhaustive yes: para que la búsqueda sea más exhaustiva.
  • | egrep –w exon: hacemos un pipe para almacenar sólo aquellas líneas del fichero que contengan la palabra “exon”. De esta forma, los intrones se descartan.
  • > fichero_resultado_exonerate.gff: con un redireccionamiento de salida, creamos un fichero donde almacenamos el resultado de toda la operación.

En el segundo paso, mediante el programa fastaseqfromGFF.pl y un redireccionamiento de salida (>) extraemos la secuencia en un fichero FASTA que llamaremos cdna.fa:

fastaseqfromGFF.pl fichero_corte.fa fichero_resultado_exonerate.gff > cdna.fa

El software Exonerate también incorpora la herramienta fastatranslate que sirve para traducir una secuencia de DNA a proteína en las 6 pautas de lectura posibles.

Nuestro fichero contiene la proteína sin intrones, de forma que el primer aminoácido del fichero es el aminoácido por el que empieza nuestra proteína extraída. Así, tenemos que coger sólo la primera de las 6 pautas de lectura. Lo conseguimos con la siguiente comanda:

fastatranslate -F 1 cdna.fa > prot.fa

  • prot.fa: nombre del fichero donde almacenamos la secuencia de aminoácidos de nuestra proteína extraída.
  • - F 1: indica que cogemos sólo la primera pauta de lectura.

El reino de los protistas es tan variado que a menudo, la distancia filogenética entre especies de protistas es tan grande como la distancia entre protistas y vertebrados. Es por esto que las Queries que usaremos provienen tanto de otros protistas como de vertebrados e invertebrados.

TCOFFEE: Alineamiento global queries vs proteína obtenida

Se alinean cada una de las proteínas obtenidas con todas las queries usadas para obtenerlas. Para esto usamos el programa TCOFFEE.

t_coffee fichero_query.fa prot.fa

Miramos si hay conservación o no entre las secuencias.

BLASTP: Búsqueda de proteínas homólogas

Comprobamos si la proteína obtenida tiene homólogos en otras especies mediante una búsqueda BLASTp contra el conjunto No-redundante de todas las proteínas disponibles en TCOFFEE:

Con BLASTp comparamos esta proteína con la base de datos de GenBank de nuestra especie target de protista. Esto también podemos usarlo en caso de que sólo hayamos encontrado un fragmento de la proteína: BLASTp nos permitirá encontrar la proteína completa.

  • Si sale en esta base de datos, entonces ya podemos anotar bien esta proteína en el genoma target.
  • Si no sale en esta base de datos, significa que tendremos que buscar la proteína en el GenBank de otro protista filogenéticamente parecido al nuestro. Así podremos ver si tiene proteínas homólogas, y también podremos anotarla.

Búsqueda de elementos SeCis

Usaremos el programa de la web SeCiSearch. Este programa es muy susceptible de proporcionar falsos positivos. Por esto es necesario acotar bien la región donde los vamos a buscar.

Para seleccionar dicha región, tenemos que partir desde la última posición de la proteína y avanzar en sentido 3' (la región UTR). En total seleccionamos una longitud de 5000Kb.

fastasubseq fichero_scaffold start length (5000) > fichero_salida_secis.fa

Si nuestra proteína se encuentra en la hebra reverse tendremos que aplicar la siguiente comanda en el fichero obtenido:

fastarevcomp fichero_salida_secis.fa > fichero_salida_secis_reverse.fa

En la web del programa SeCis Sarch introducimos el fichero, y comprobamos la existencia de SeCis con los diferentes Patterns que nos proporciona la web. Si nos aparece más de un SeCis, entonces tendremos que seleccionar el que tenga Cove Score más alto.

Automatización del proceso

Con tal de automatizar todos los pasos, creamos además dos programas distintos:

  • Blastsgrupdy.pl
  • El programa Blastsgrupdy.pl es un programa que realiza todos los blasts de nuestras familias de selenoproteínas a la vez. Este programa está escrito en formato shell, es decir, todas las comandas que se encuentran dentro del fichero serán ejecutados por la terminal. El programa necesita la carpeta TotesSeleno donde se encuentran todas las queries. El siguente paso es comparar las diferentes queries y la información de la base de datos (/cursos/BI/ genomes/protists_List_and_info.tab) obteniendo 4 documentos de BLAST, uno para cada selenoproteína (blastallT.txt, blastallV.txt, blastallW1.txt y blastallW2.txt).

  • new.pl
  • El programa new.pl es un programa confeccionado para reducir la tarea de realizar repetidas veces las mismas acciones. Este programa necesita 5 documentos en shell para realizar su función: exonerate1.pl, scaf.pl, noscaf.pl, start1.pl y start2.pl. La función de este programa es realizar el fastasubindex, fastasubfech, fastasubseq, exonerate, genwise, fastaseqfromGFF.pl, fastatranslate y t-coffee del exonerate y del genwise.

    El programa empieza con una serie de preguntas para tener toda la información necesaria para realizar todas sus funciones. Una vez terminadas las preguntas, el programa de perl new.pl es el encargado de activar los diferentes componentes para hacer la tarea.

    Al final se muestran los diferentes documentos dentro de una carpeta que tiene el nombre del protista estudiado que se encuentra dentro de la carpeta de la família de selenoproteína a la que pertenece.

Ambos programas nos dieron la misma información que la obtenida manualmente y por tanto podemos afirmar que funcionan.