BLAST
En este capitulo vamos a explorar el uso de blast en local, es decir, instalado en nuestras máquinas. Lo primero que teneis que hacer, si no lo habeis hecho ya, es seguir las instrucciones para instalar Blast en vuestra máquina, podeis encontrarlo aqui.
En este tema solo se van a dar unas pinceladas de como funciona Blast y sus parámetro. Para el resto tendreis que visitar la documentación y su manual https://www.ncbi.nlm.nih.gov/books/NBK279690/ (esta disponible para descargar en PDF).
Homologia de secuencias
En biología, dos genes se definen como homólogos si derivan de un ancestro común. Pero la homología no termina ahí, ya que dos secuencias pueden ser homólogas ortólogas o parálogas, tal y como podemos observar en la Figura1.
Ortólogos Secuencias similares en dos os mas organismos diferentes que provienen de un proceso de especiación. Se espera que tengan la misma función biológica.
parálogos Secuecuencias similares dentro de un mismo organismo que proceden de un proceso de duplicación génica. No tienen porque tener la misma función
Alineamiento de secuencias
Un alineamiento de secuencias no es mas que el reordenamiento de dos cadenas de caracteres para conseguir el mayor numero de matches posibles entre ambas. El mejor alineamiento posible es aquel que muestra un mayor nivel de correspondencia entre ambas secuencias, y el menor numero de diferencias o discrepancias. Dos secuencias similares tienen una alta probabilidad de tener una misma función, pero siempre se cumple (recordar parálogos). Términos importantes:
identidad: El mismo aminoácido o nucleotido en la misma posicion en dos secuencias.
similaridad: Numero de positivos entre ambas secuencias, es decir, sustituciones de aminoácidos o nucleótidos con las mismas propiedades bioquímica en dos secuencias.
>Importante similaridad no es igual que homología
Si las secuencias son mayores de 100 nucleótidos o 100 aminoácidos, puedes decir que dos proteinas son homólogas si al menos el 25 % de los aminoácidos son idénticos. En el caso del DNA se requiere un porcentage de identidad del 70 %.
Hasta la fecha miles de genomas de diferentes especies se han secuenciado y están disponibles para que trabajemos con ellas, revelando información valiosa sobre sus genes y su función correspondiente. Para asignar funciones a los genes de nuevos genomas recien secuenciados se utilizan herramientas que buscan y miden los niveles de similitud entre millones de secuencias, ahí es donde entran en juegos herramientas como BLAST y otros alineadores de secuencias.
BLAST
BLAST es el acrónimo de Basic Local Alignment Search Tool. es el algoritmo mas empleado por el NCBI. es un programa informáJco de alineamiento de secuencias de tipo local, ya sea de ADN, ARN o de proteínas. El programa es capaz de comparar una secuencia problema (query) contra una gran canJdad de secuencias que se encuentren en una base de datos (subject).
Tipos de Blast
Podemos usar todos los tipos de blast con blast local, es decir, blastx, blastp, blastn y tblastn al instalarnos blast+ de NCBI en local. Recordad que:
- blastn: nucleotide > nucleotide (nucleotide blast)
- blastp : protein > protein (protein blast)
- tblastn : protein > translated nucleotide
- blastx: translated nucleotide > protein
Syntaxis general de blast
Cada vez que utilicemos blast usaremos un comando diferente segun el tipo de blast de los listados arriba vayamos a usar, es decir (blastp, blastn, blastx o tblastn). Además, es obligatorio indicar el archivo fasta que contiene la secuencia query (-query), el archivo con la/s secuencia/s subjects (-subject) e indicar el output si queremos guardar el resultado en un archivo con la opción –out
~$ [t]blast[pnx] -query fasta -subject fasta -out result.out
–query: secuencia problema fasta
–subject: archivo fasta con la/s secuencia/s que enfrentaremos a nuestra secuencia problema
–out: nombre del archivo de salida
veamos un ejemplo con blastp. Para ello usaremos los archivos fasta PBP1_staphylococcus.fasta como query PBP1.fasta como subject. Estos archivos contienen proteinas PBPs que son proteinas de union a la penicilina de diferentes bacterias.
~$ blastp -query PBP1_staphylococcus.fasta -subject PBP1.fasta -out result.out
...
Query= ASD49649.1 penicillin binding protein 1 [Staphylococcus aureus]
Length=744
Score E
Sequences producing significant alignments: (Bits) Value
ANW82078.1 penicillin binding protein 1 [Staphylococcus aureus] 1522 0.0
ASD49626.1 penicillin binding protein 1 [Staphylococcus aureus] 1522 0.0
ASD49630.1 penicillin binding protein 1 [Staphylococcus aureus] 1522 0.0
ASD49632.1 penicillin binding protein 1 [Staphylococcus aureus] 1522 0.0
ASD49633.1 penicillin binding protein 1 [Staphylococcus aureus] 1522 0.0
ASD49643.1 penicillin binding protein 1 [Staphylococcus aureus] 1522 0.0
ASD49646.1 penicillin binding protein 1 [Staphylococcus aureus] 1522 0.0
ASD49647.1 penicillin binding protein 1 [Staphylococcus aureus] 1522 0.0
Se obtiene un archivo con una estructura de datos muy similar a la que vemos cuando usamos la version web, sin embargo, el problema es que estas estructura es muy dificil de parsear. Las opciones de salida de datos que permite Blast son extensas y aqui solo vamos a mostrar algunas, os recomiendo que mireis el manual, si qureis profundizar en el tema.
En la Tabla 1 teneis algunas de las opciones disponibles comunes a todos los tipos diferentes de Blast.
Tabla 1. Opciones generales de Blast
| option | type | default value | description and notes |
|---|---|---|---|
| db | string | none | BLAST database name. |
| query | string | stdin | Query file name. |
| query_loc | string | none | Location on the query sequence (Format: start-stop) |
| out | string | stdout | Output file name |
| evalue | real | 10 | Expect value (E) for saving hits |
| subject | string | none | File with subject sequence(s) to search. |
| show_gis | flag | N/A | Show NCBI GIs in report. |
| num_descriptions | integer | 500 | Show one-line descriptions for this number of database sequences. |
| num_alignments | integer | 250 | Show alignments for this number of database sequences. |
| html | flag | N/A | Produce HTML output |
| gilist | string | none | Restrict search of database to GI’s listed in this file. Local searches only. |
| negative_gilist | string | none | Restrict search of database to everything except the GI’s listed in this file. Local searches only. |
| entrez_query | string | none | Restrict search with the given Entrez query. Remote searches only. |
| best_hit_overhang | real | none | Best Hit algorithm overhang value (recommended value: 0.1) |
| best_hit_score_edge | real | none | Best Hit algorithm score edge value (recommended value: 0.1) |
| dbsize | integer | none | Effective size of the database |
| searchsp | integer | none | Effective length of the search space |
| import_search_strategy | string | none | Search strategy file to read. |
| export_search_strategy | string | none | Record search strategy to this file. |
| parse_deflines | flag | N/A | Parse query and subject bar delimited sequence identifiers (e.g., gi|129295). |
| num_threads | integer | 1 | Number of threads (CPUs) to use in blast search. |
| remote | flag | N/A | Execute search on NCBI servers? |
output format
Si utilizamos la opción -outfmt podremos usar diferentes estilos de salida de datos (18 en total), nosotros usaremos la numero 6 (tabular), veamos que aspecto tiene:
~$ blastp -query PBP1_staphylococcus.fasta -subject PBP1.fasta -outfmt "6"
ASD49649.1 ANW82078.1 100.000 744 0 0 1 744 1 744 0.0 1522
ASD49649.1 ASD49626.1 100.000 744 0 0 1 744 1 744 0.0 1522
ASD49649.1 ASD49630.1 100.000 744 0 0 1 744 1 744 0.0 1522
ASD49649.1 ASD49632.1 100.000 744 0 0 1 744 1 744 0.0 1522
ASD49649.1 ASD49633.1 100.000 744 0 0 1 744 1 744 0.0 1522
ASD49649.1 ASD49643.1 100.000 744 0 0 1 744 1 744 0.0 1522
ASD49649.1 ASD49646.1 100.000 744 0 0 1 744 1 744 0.0 1522
ASD49649.1 ASD49647.1 100.000 744 0 0 1 744 1 744 0.0 1522
ASD49649.1 ASD49649.1 100.000 744 0 0 1 744 1 744 0.0 1522
ASD49649.1 ASD49651.1 100.000 744 0 0 1 744 1 744 0.0 1522
...
Pero quizas a veces necesitemos menos información. Utilizando la opción -outfmt tambien podemos especificar que campos son los que queremos mostrar, en este ejemplo, mostraremos las siguientes opciones que son de las mas importantes en cualquier analisis de blast:
- qseqid: query seq ID
- sseqid: subject seq ID
- pident: Percentage of identical matches
- length: Alignment length
- mismatch: Number of mismatches
- gapopen: Number of gap openings
- qstart: Start of alignment in query
- qend: End of alignment in query
- sstart: Start of alignment in subject
- send: End of alignment in subject
- evalue: Expect value (E) for saving hits
- bitscore: Bit score
como podeis ver ahora obtenemos un formato en columna mucho mas apetecible de ser parseado. las opciones mostradas por default en la opción 6 son por orden:
opciones que son de las mas usadas en cualquier análisis de blast:
qseqid
sseqid
qcovs
pident
evalue
~$ blastp -query PBP1_staphylococcus.fasta -subject PBP1.fasta -outfmt "6 qseqid sseqid qcovs pident evalue" ASD49649.1 ANW82078.1 100 100.000 0.0 ASD49649.1 ASD49626.1 100 100.000 0.0 ASD49649.1 ASD49630.1 100 100.000 0.0 ASD49649.1 ASD49632.1 100 100.000 0.0 ASD49649.1 ASD49633.1 100 100.000 0.0 ASD49649.1 ASD49643.1 100 100.000 0.0 ASD49649.1 ASD49646.1 100 100.000 0.0 ASD49649.1 ASD49647.1 100 100.000 0.0 ASD49649.1 ASD49649.1 100 100.000 0.0 ASD49649.1 ASD49651.1 100 100.000 0.0 ASD49649.1 ANW82077.1 100 99.866 0.0 ASD49649.1 ASD49623.1 100 99.866 0.0 …
tambien podemos controlar el valor de cut off del e-value por defecto es 10, pero no deberíamos de considerar resultados que no estén por debajo de 1xe-5, a menos que sepamos muy bien lo que estamos haciendo.
Por ejemplo veamos que pasa si filtramos por e-value de forma restrictiva utilizando la opción -evalue:
~$ blastp -query PBP1_staphylococcus.fasta -subject PBP1.fasta -evalue "0.0000001" -outfmt "6 qseqid sseqid qcovs pident evalue"
vereis que al ser mas restrictivo a nivel de e-value hemos obtenido una cantidad menor de hits.
Bases de datos de Blast
En vez de analizar mediante Blast un archivo Fasta, podemos generar o descargar de NCBI bases de datos de BLAST con secuencias de ADN o proteínas. La ventaja de usar bases de datos, es que la velocidad de analisis y comparación de secuencias de BLAST aumentará notablemente, especialmente cuanto utilizemos un conjunto de secuencias subject del orden de cientos de miles o millones.
Crear una base de datos en blast es muy sencillo:
~$ makeblastdb –in mydb.fsa –dbtype "nucl/prot" –out mydb
–in: el multifasta con las secuencias que formarán la base de datos.
–dbtype: Indica el tipo de datos que van a formar la base de datos, proteinas (prot) o nucleótidos (nucl).
–out: el nombre que le pondremos a la base de datos
Hagamos un ejemplo con nuestras secuencias PBPs, vamos a crear una base de datos que se llamará:
~$ makeblastdb -in PBP1.fasta -dbtype 'prot' -out PBP1
Building a new DB, current time: 04/24/2020 09:50:15
New DB name: /Users/JGL/Dropbox/programacion_bioinformatica/Tema_5_blast/PBP1
New DB title: PBP1.fasta
Sequence type: Protein
Keep MBits: T
Maximum file size: 1000000000B
Adding sequences from FASTA; added 54 sequences in 0.00509191 seconds.
Ahora hechemos un vistazo a los archivos que se han creado:
-rw-r--r--@ 1 JGL staff 6.8K Apr 24 09:50 PBP1.phr
-rw-r--r--@ 1 JGL staff 504B Apr 24 09:50 PBP1.pin
-rw-r--r--@ 1 JGL staff 34K Apr 24 09:50 PBP1.psq
Lo mas interesante es guardar las bases de datos en un directorio, con el mismo nombre. En este caso por ejemplo PBP1/.
Como podemos observar las bases de datos de Blast contienen tres clases de archivos, .phr .pin y .psq.
Una vez creada nuestra base de datos, hacer un blast frente a la misma es muy fácil, solo tenemos que
~$ blastp -query PBP1_staphylococcus.fasta -db PBP1/PBP1 -outfmt "6 qseqid sseqid qcovs pident evalue"
Tambien podeis bajaros bases de datos de NCBI y utlizarlas para hacer análisis en local, de aqui podeis bajaros la que querais y explorarlas:
ftp://ftp.ncbi.nlm.nih.gov/blast/db/
Podeis copias el link address de cualquiera de ellas y descargadlas en un directorio con wget o curl. Utilizando el comando Blastdbcmd podréis hacer análisis de las bases datos que os descarguéis, ¡os animo a explorarlo!.
