Blast
Introducció
Section titled “Introducció”Basic Local Alignment Search Tool (BLAST) és l’eina principal de l’NCBI per comparar una seqüència d’ADN o proteïna amb altres seqüències de diverses bases de dades.
La cerca BLAST és una de les maneres fonamentals d’aprendre sobre una proteïna o un gen: la cerca revela quines seqüències relacionades estan presents en el mateix organisme i en altres organismes.
Blast+
Section titled “Blast+”Els programes BLAST hostetjats al lloc web de NCBI són extremadament populars.
Com a alternativa, també pots utilitzar BLAST+, un conjunt d’executables de línia d’ordres per a ús local. Aquesta alternativa és molt interessant, ja que així pots configurar bases de dades personalitzades, fer cerques massives (és a dir, utilitzar un nombre molt gran de consultes), implementar estratègies de cerca complexes mitjançant scripts personalitzats o utilitzar el teu clúster d’ordinadors per millorar el rendiment.
En aquest enllaç tens el manual: BLAST Command Line Applications User Manual.
Entorno de trabajo
Section titled “Entorno de trabajo”Per treballar amb blast utilitzarem la imatge ncbi\blast.
Descarrega l’script blast.sh:
$ wget https://gitlab.com/xtec/bio-blast/-/raw/main/blast.sh$ chmod u+x blast.shExecuta l’script - el primer cop tardará més temps en executar-se perquè ha de descarregar la imatge:
$ ./blast.sh
...Status: Downloaded newer image for ncbi/blast:latestroot@8e54f06bb10e:/blast#L’script també crea unes carpetes compartides entre el contenidor i la teva màquina:
$ ls blast...Per sortir del contenidor executa exit:
$ exitCrea una base de dades personalitzada
Section titled “Crea una base de dades personalitzada”BLAST et permet crear bases de dades personalitzades
Per realitzar la consulta utilitzarem una base de dades amb 7 seqüències de proteïnes del tauró nodrissa.
Primer hem de descarregar els fitxers fasta de la base de dades protein mitjançant Entrez:
$ ./blast.sh$ efetch -db protein -format fasta -id Q90523,P80049,P83981,P83982,P83983,P83977,P83984,P83985,P27950 > fasta/nurse-shark-proteins.fastaPots veure que les seqüències s’han guardat en el fitxer nurse-shack-proteins.fasta:
$ cat fasta/nurse-shark-proteins.fasta | head -n 1
>sp|P27950.1|NDK_GINCI RecName: Full=Nucleoside diphosphate kinase; Short=NDK; Short=NDP kinaseA continuació crea una base de dades BLAST amb el nom nurse-shark-proteins:
$ makeblastdb -in fasta/nurse-shark-proteins.fasta -dbtype prot -parse_seqids -out blastdb_custom/nurse-shark-proteins -title "Nurse shark proteins" -taxid 7801 -blastdb_version 5
Building a new DB, current time: 03/18/2024 17:08:17New DB name: /blast/blastdb_custom/nurse-shark-proteinsNew DB title: Nurse shark proteinsSequence type: ProteinKeep MBits: TMaximum file size: 3000000000BAdding sequences from FASTA; added 7 sequences in 0.00194597 seconds.Com que els fitxers fasta són pocs i petits, la base de dades es construeix en un moment.
Pots veure els fitxers que formen la base de dades:
$ ls blastdb_custom/
nurse-shark-proteins.pdb nurse-shark-proteins.pin nurse-shark-proteins.pog nurse-shark-proteins.pot nurse-shark-proteins.ptfnurse-shark-proteins.phr nurse-shark-proteins.pjs nurse-shark-proteins.pos nurse-shark-proteins.psq nurse-shark-proteins.ptoPots verificar la base de dades BLAST que acabes de crear executar l’ordre següent per mostrar les accessions, la longitud de la seqüència i el nom comú de les seqüències a la base de dades.
$ blastdbcmd -entry all -db nurse-shark-proteins -outfmt "%a %l %T"
P27950.1 151 7801P80049.1 132 7801P83977.1 95 7801P83981.1 53 7801P83984.1 190 7801P83985.1 195 7801Q90523.1 106 7801Consultar la base de dades
Section titled “Consultar la base de dades”Per realitzar una consulta, primer hem d’obtenir la seqüència de consulta.
En el nostre exemple serà la proteína P01349 del tauró tigre de sorra que es va seqüenciar el 1986:
$ efetch -db protein -format fasta -id P01349 > queries/P01349.fastaMira el contingut del fitxer fasta que acabes de baixar:
$ cat queries/P01349.fasta
>sp|P01349.2|RELX_CARTA RecName: Full=Relaxin; Contains: RecName: Full=Relaxin B chain; Contains: RecName: Full=Relaxin A chainQLCGRGFIRAIIFACGGSRWATSPAMSIKCCIYGCTKKDISVLCA continuació pots executar la consulta blastp:
$ blastp -query queries/P01349.fasta -db nurse-shark-proteins
...>P80049.1 RecName: Full=Fatty acid-binding protein, liver; AltName: Full=Liver-typefatty acid-binding protein; Short=L-FABPLength=132
Score = 14.2 bits (25), Expect = 0.96, Method: Compositional matrix adjust. Identities = 3/9 (33%), Positives = 6/9 (67%), Gaps = 0/9 (0%)
Query 2 LCGRGFIRA 10 +C R ++RSbjct 123 VCTREYVRE 131...El resultat d’aquest consulta és la identificació de la seqüència de proteïnes P80049.1 com una coincidència amb una puntuació de 14,2 i un valor E de 0,96.
Per un anàlisi posterior pots utilitzar l’opció -out per guardar la sortida en un fitxer.
$ blastp -query queries/P01349.fasta -db nurse-shark-proteins -out results/blastp.outProtein Data Bank
Section titled “Protein Data Bank”Les bases de dades personals són útils quan estàs treballant amb seqüències que no són públiques, però en la majoria dels casos es treballa amb bases de dades públiques.
Primer has de descarregar la base de dades d’aminoàcids Protein Data Bank (pdbaa)
El fitxer update_blastdb.pl` és un script de Perl que proporciona l’NCBI que ja està incorporat i configurat a les variables d’entorn.
L’has d’executar dins de la carpeta blastdb!
Fem servir l’opció -source gcp perquè descarregui les dades desde Google Cloud enlloc del servidor FTP del NCBI:
$ cd blastdb$ update_blastdb.pl -source gcp pdbaaConnected to GCP$ cd ..Pots veure els fitxers que composen la base de dades s’han descarregat a la carpeta blastdb:
$ ls blastdb
pdbaa.pdb pdbaa.phr pdbaa.pin pdbaa.pog pdbaa.pos pdbaa.pot pdbaa.ppd pdbaa.ppi pdbaa.psq pdbaa.ptf pdbaa.ptoA continuació podem consultar la seqüència P01349.fasta contra la base de dades PDB tal com es mostra a continuació:
$ blastp -query queries/P01349.fasta -db pdbaa...Sequences producing significant alignments: (Bits) Value
2FHW_B Chain B, Relaxin 3 (Prorelaxin H3) (Insulin-like peptide I... 38.1 1e-052H8B_B Chain B, Insulin-like 3 [synthetic construct] 34.7 3e-042K6U_B Chain B, Insulin-like 3 B chain [synthetic construct] 30.8 0.010...En aquest cas hi ha moltes més coincidències i amb un score més alt.
Bases de dades
Section titled “Bases de dades”Tal com has vist abans update_blastdb.pl és un script en perl que et permet gestionar les bases de dades.
Executa update_blastdb.pl sense arguments perquè et mostri el menú d’ajuda, i poder veure totes les bases de dades disponibles
$ update_blastdb.pl --showall18S_fungal_sequencesBetacoronavirus...1.- Descarrega i descomprimeix la base de dades swissprot i realitza una consulta BLAST a la proteïna P01929 que guardi els resultats en un fitxer a la carpeta results.
update_blastdb.pl --source gcp --decompress swissprot{% sol %}
efetch -db protein -format fasta -id P01929 > queries/P01929.fastablastp -query queries/P01929.fasta -db swissprot -out results/P01929_blast_results.txtAlguns resultats obtinguts:
Query= sp|P01929.1|HBA_LEOFU RecName: Full=Hemoglobin subunit alpha;AltName: Full=Alpha-globin; AltName: Full=Hemoglobin alpha chain;Contains: RecName: Full=Hemopressin
Length=141 Score ESequences producing significant alignments: (Bits) Value
P01929.1 RecName: Full=Hemoglobin subunit alpha; AltName: Full=Al... 282 1e-98P63107.2 RecName: Full=Hemoglobin subunit alpha; AltName: Full=Al... 277 1e-96P67817.2 RecName: Full=Hemoglobin subunit alpha; AltName: Full=Al... 276 2e-96P67818.1 RecName: Full=Hemoglobin subunit alpha; AltName: Full=Al... 276 3e-96P01923.1 RecName: Full=Hemoglobin subunit alpha; AltName: Full=Al... 275 7e-96{% endsol %}
Algoritme
Section titled “Algoritme”Els algorismes de semblança global optimitzen l’alineament global de dues seqüències, i són els més adequats per trobar coincidències que consisteixen en trams llargs de poca similitud.
En canvi, els algorismes de similitud local com BLAST identifiquen alineaments relativament curts, i són útils per a la cerca en bases de dades perquè moltes seqüències de consulta (“query”) tenen dominis, llocs actius o altres motius que tenen regions locals de semblança amb altres proteïnes ( però no globals),
Normalment, les bases de dades també tenen fragments d’ADN i seqüències de proteïnes que es poden alinear localment a una consulta.
Llista
Section titled “Llista”En la primera fase, l’algorisme BLASTP compila una llista de “paraules” de longitud fixa w que es deriven de la seqüència de consulta.
Per a les cerques de proteïnes, la mida de la paraula normalment té un valor predeterminat de 3 aminoàcids, encara que es pot modificar. Com que hi ha 20 aminoàcids, hi ha 20^3 = 8000 paraules possibles. Això és així perquè es tracta d’un cas de variacions amb repetició (és important l’ordre de cada aminoàcid i es poden repetir).
Per exemple si la seqüència “query” és la “human beta globin” (NP_000509.1) la seqüència incloura aquests aminoàcdis …VTALWGKVNVD… i algunes de les paraules derivades de la seqüència de consulta seran VTA TAL ALW ``LWG WGK GKV KVN VNV NVD.
Per tant, si una seqüència té 100 aminoàcids tindrem 98 paraules de consulta de mida 3.
A continuació, per cada paraula de consulta es genera una llista de paraules objectiu (“target”) de la mateixa mida amb la seva puntuació de similitud derivada d’una matriu com pot ser BLOSUM62.
Aquesta llista està ordenada pel valor de similitud entre la paraula de consulta i la paraula target.
Per exemple, per a la paraula de consulta LWG es pot generar aquesta llista:
LWG 4+11+6=21 LWA 4+11+0=15IWG 2+11+6=19 LYG 4+ 2+6=12MWG 2+11+6=19 --------------VWG 1+11+6=18 LFG 4+ 1+6=11FWG 0+11+6=17 FWS 0+11+0=11AWG 0+11+6=17 AWS -1+11+0=10LWS 4+11+0=15 CWS -1+11+0=10LWN 4+11+0=15 IWC 2+11-3=10No totes les paraules es fan servir, sinó aquelles que estiguin per sobre d’un llindar establert T (“threshold”). En el nostre exemple hem establert el llindar a 12.
Si reduïm el valor de T augmentem la sensibilitat de la cerca (tindrem més “hits” perquè tindrem més paraules target) a canvi de que la cerca trigui més temps en executar-se.
Scan y extensiones
Section titled “Scan y extensiones”A continuació, per cada paraula de consulta es busca a l’índex de la base de dades totes les paraules target associades a la paraula de consulta amb una puntuació superior al llindar T. Si una paraula target és de mida 3, i suposem que tots els aminoàcids tenen la mateixa freqüència d’aparició (que sabem que no és cert), una paraula target tindrà 1 “hit” cada 8000 seqüències de la base de dades. Si una seqüència te un “hit”, es comença a fer una alieament en ambdues direccions tal com es mostra en aquest exemple amb l’alpha globin:
LTPEEKSAVTA[LWG]KV--NVDEVGGEALGRLLVVYPWTQRFFESFGDLSTPDAVMGNPKV HBBL+P +K+ V A[ WG]KV + E G EAL R+ + +P T+ +F F D G+ +VLSPADKTNVKA[AWG]KVGAHAGEYGAEALERMFLSFPTTKTYFPHF------DLSHGSAQV HBATal com pots veure a l’exemple, primer l’extensió es fa sense buits i després amb buits si la puntuació( “score”) ha arribat a un llindar, per obtenir parells de segments d’alta puntuació (HSP).
Si una puntuació HSP supera una puntuació de tall S particular, aquest alineament local passa a formar part del resultat.
Si durant l’extensió la puntuació cau per sota d’un límit l’alineament es descarta la seqüència.
Traceback
Section titled “Traceback”Un cop es tenen totes les seqüències amb un puntuació HSP superior al valor de tall S, s’assignen les ubicacions d’insercions, supressions i desajustos de les seqüències seleccionades. També s’apliquen unes estadístiques que descriurem més endavant.
BLASTN
Section titled “BLASTN”Per a BLASTN, la primera fase és lleugerament diferent perquè l’algoritme exigeix coincidències exactes de paraules (no hi ha matrius de similitut en l’alineament de seqüències d’ADN).
La mida predeterminada de la paraula és 11 – es pot ajustar-la a valors de 7 o 15. Reduir la longitud de la paraula té el mateix efecte que reduir la puntuació llindar en consultes blastp.
Threshold
Section titled “Threshold”A continuació es demostrar l’efecte de diferents nivells de llindar en una cerca BLASTP canviant el paràmetre T del seu valor predeterminat (11) a un rang d’altres valors en una consulta de betaglobina humana:
{% image “blastp.png” %}
Python
Section titled “Python”Amb Python pots interactuar directament amb les eines BLAST+ tal com s’explica a {% link “/python/docker/” %}.
A continuació tens un exemple; que funcionarà fora del contenidor si has seguit totes les indicacions prèvies:
import osfrom python_on_whales import docker
CWD = os.getcwd()
def blast(command): output = docker.run( image="ncbi/blast", remove=True, volumes=[(f"{CWD}/blast/blastdb", "/blast/blastdb"), (f"{CWD}/blast/blastdb_custom", "/blast/blastdb_custom"), (f"{CWD}/blast/fasta", "/blast/fasta"), (f"{CWD}/blast/queries", "/blast/queries")], command=command) return output
output = blast(["efetch", "-db", "protein", "-format", "fasta", "-id", "NP_000509"])
print(output)Pots veure que la funció blast et retorna el fitxer fasta de l’hemoglobina beta:
$ python3 blastdemo.py>NP_000509.1 hemoglobin subunit beta [Homo sapiens]MVHLTPEEKSAVTALWGKVNVDEVGGEALGRLLVVYPWTQRFFESFGDLSTPDAVMGNPKVKAHGKKVLGAFSDGLAHLDNLKGTFATLSELHCDKLHVDPENFRLLGNVLVCVLAHHFGKEFTPPVQAAYQKVVAGVANALAHKYHAl realitzar una consulta amb blastp pots utilitzar l’argument -outfmt per especificar un format que puguis utilitzar directament amb Python.
A continuació es mostren tots els formats disponibles en JSON (amb la comanda blastp -help pots veure els formats disponibles):
$ blastp --help... 12 = Seqalign (JSON), 13 = Multiple-file BLAST JSON, 15 = Single-file BLAST JSON,...Per analitzar un resultat, pots indicar que la sortida sigui en BLAST JSON tal com es mostra a continuació:
record = blast(["blastp", "-query", "/blast/queries/P01349.fasta", "-db", "pdbaa", "-outfmt", "15"])record = json.loads(record)hits = record["BlastOutput2"][0]["report"]["results"]["search"]["hits"]print(json.dumps(hits[0], indent=4))Pots veure com tens accés a tota la informació al convertir l’string JSON en un diccionari:`
{ "num": 1, "description": [ { "id": "pdb|2FHW|B", "accession": "2FHW_B", "title": "Chain B, Relaxin 3 (Prorelaxin H3) (Insulin-like peptide INSL7) (Insulin-like peptide 7) [synthetic construct]",...Cerca BLAST
Section titled “Cerca BLAST”A continuació es descriu els diferents passos per fer una consulta.
Especificar la seqüència d’interès
Section titled “Especificar la seqüència d’interès”Una cerca BLAST comença amb la selecció d’una seqüència d’ADN o proteïna a analitzar.
La seqüència es pot passar directament – per exemple exemple en el format FASTA – o mitjançant el número d’accés – per exemple, un número d’identificació RefSeq o GenBank (GI).
La cerca BLAST també et permet seleccionar un subconjunt d’una seqüència de consulta sencera, com ara una regió o un domini d’interès.
Seleccionar el programa
Section titled “Seleccionar el programa”La família de programes NCBI BLAST inclou cinc programes principals:
-
El programa BLASTP (protein) compara una seqüència de consulta d’aminoàcids amb una base de dades de seqüències de proteïnes. Tingues en compte que per a aquest tipus de cerca hi ha paràmetres opcionals que són especialment rellevants per a les cerques de proteïnes, com ara l’elecció de diverses matrius de puntuació PAM i BLOSUM.
-
El programa BLASTN (nucleotid) s’utilitza per comparar una seqüència de consulta de nucleòtids amb una base de dades de seqüències de nucleòtids.
Ja saps que qualsevol seqüència d’ADN es pot transcriure i traduir en sis marcs de lectura potencials (tres a la cadena superior i tres a la cadena inferior). Pert tant, donada una seqüència d’ADN, tenim 6 seqüències proteïnes potencials.
Els altres tres programes utilitzen aquest fet per realitzar alineaments proteïna-proteïna.
-
El programa BLASTX compara una seqüència de consulta de nucleòtids, traduïda a 6 proteïnes potencials, amb una base de dades de seqüències de proteïnes.
Per tant, si tens una seqüència d’ADN i vols saber si codifica alguna proteïna pots fer una cerca
blastx. -
El programa TBLASTN compara una seqüència de consulta de proteïnes, traduïda a 6 seqüencies d’ADN, amb una base de dades de seqüències de nucleòtids.
Per tant, si tens una base de dades d’ADN i vols saber si quines seqüències de la base de dades codifiquen una proteïna, pots fer una cerca tblastn.
-
El programa TBLASTX compara una seqüència de consulta de nucleòtids, traduïda a 6 proteïnes, que a la seva vegada es tradueix cada proteïna a 6 seqüencies d’ADN, amb una base de seqüències de nucleòtids.
Per exemple, si tens una seqüència d’ADN sense coincidències de base de dades evidents i vols saber si codifica una proteïna amb coincidències de base de dades distants i estadísticament significatives en una base de dades d’etiquetes de seqüències expressades, una cerca tblastx és més sensible que blastn i, per tant, útil per revelar gens que codifiquen proteïnes homòlogues a la teva consulta.
Selecció d’una base de dades
Section titled “Selecció d’una base de dades”Les bases de dades disponibles per a la cerca BLAST es mostren a cada pàgina BLAST.
-
Per a les cerques de bases de dades de proteïnes (BLASTP i BLASTX), l’opció predeterminada és la base de dades no redundant (nr). Consisteix en els registres de proteïnes combinats de GenBank, el Protein Data Bank (PDB), SwissProt, PIR i PRF. Una altra opció és cercar només proteïnes Refeq.
-
Per a les cerques de bases de dades d’ADN (BLASTN, TBLASTN, TBLASTX) l’opció predeterminada és cercar a la base de dades de nucleòtids nr/nt. Això inclou seqüències de nucleòtids de GenBank, EMBL, DDBJ, PDB i RefSeq. No obstant això, la base de dades nr no té registres de l’EST, el lloc etiquetat per seqüències (STS), la seqüència de genoma sencer (WGS), la seqüència d’enquesta del genoma (GSS), l’assemblatge del transcriptoma (TSA), les patents o la seqüència genòmica d’alt rendiment. bases de dades (HTGS).
Altres opcions que s’utilitzen habitualment inclouen la base de dades genòmica més transcripció humana (o ratolí) o la base de dades EST.
Les bases de dades nr es deriven fusionant diverses bases de dades principals de proteïnes o ADN. Aquestes bases de dades sovint contenen seqüències idèntiques. Generalment, la base de dades nr conserva només una d’aquestes seqüències, juntament amb diversos números d’accés. (Encara que dues seqüències de la base de dades nr semblen ser idèntiques, almenys haurien de tenir alguna diferència subtil.) Les bases de dades nr són sovint els llocs preferits per cercar la majoria de seqüències disponibles.
Selecció de paràmetres de cerca opcionals
Section titled “Selecció de paràmetres de cerca opcionals”TODO
Revisar taules:
https://www.ncbi.nlm.nih.gov/books/NBK279684/
NCBI-BLAST API
Section titled “NCBI-BLAST API”La API de URL común NCBI-BLAST te permite realizar búsquedas de forma remota. Soporta lo mismo comandos en el servidor web NCBI y en la instalación de un proveedor de nube. Esto te permite ejecutar búsquedas en el servidor web del NCBI o en un proveedor de nube sin tener que modificar el código.
La API REST de l’NCBI et permet realitzar consultes BLAST+, però amb un límit de 1 consulta cada 10 segons.
- Necessites fer cerques BLAST de gran volum?
- Tens dades de seqüència propietaries per cercar i no pots utilitzar el lloc web de NCBI BLAST?
- Voleu evitar la feina de muntar i mantenir BLAST autònom?
NCBI ofereix ElasticBLAST que pot executar les vostres cerques al núvol.
Una altra opció és que tu mateix construeixis un aplicació web que permeti fer les consultes en un servidor.
Es tracta que no només provis els exemples proporcionats en aquestes activitats sinó que usis el teu nou servidor local per altres BLAST que consideris interessants.
Ja sigui cercant exemples de BLAST ja realitzats de fonts fiables o bé experimentant amb gens i organismes que et resultin interessants.
Algunes fonts que pots provar si vols:
Millores
Pots crear una interfície web.