--- tags: Random title: Command-line blast example --- # Command-line blast example [toc] ## Making conda env ```bash conda create -n blast -c conda-forge -c bioconda -c defaults blast conda activate blast ``` ## Making example files ```bash printf ">NR_024570.1 Escherichia coli strain U 5/41 16S ribosomal RNA, partial sequence AGTTTGATCATGGCTCAGATTGAACGCTGGCGGCAGGCCTAACACATGCAAGTCGAACGGTAACAGGAAGCAGCTTGCTGCTTTGCTGACGAGTGGCGGACGGGTGAGTAATGTCTGGGAAACTGCCTGATGGAGGGGGATAACTACTGGAAACGGTAGCTAATACCGCATAACGTCGCAAGCACAAAGAGGGGGACCTTAGGGCCTCTTGCCATCGGATGTGCCCAGATGGGATTAGCTAGTAGGTGGGGTAACGGCTCACCTAGGCGACGATCCCTAGCTGGTCTGAGAGGATGA >NR_169460.1 Pseudomonas aylmerensis strain S1E40 16S ribosomal RNA, partial sequence CTCAGATTGAACGCTGGCGGCAGGCCTAACACATGCAAGTCGAGCGGTAGAGAGAAGCTTGCTTCTCTTGAGAGCGGCGGACGGGTGAGTAATGCCTAGGAATCTGCCTGGTAGTGGGGGATAACGTTCGGAAACGGACGCTAATACCGCATACGTCCTACGGGAGAAAGCAGGGGACCTTCGGGCCTTGCGCTATCAGATGAGCCTAGGTCGGATTAGCTAGTTGGTGGGGTAATGGCTCACCAAGGCGACGATCCGTAACTGGTCTGAGAGGATGATCAGTCACACTGGAACTGA " > refs.fa printf ">Q1 AGTTTGATCATGGCTCAGATTGAACGCTGGCGGCAGGCCTAACACATGCAAGTCGAACGGTAACAGGAAGCAGCTTGCGGGGGTGCTGACGAGTGGCGGACGGGTGAGTAATGTCTGGGAAACTGCCTGATGGAGGGGGATAACTACTGGAAACGGTAGCTAATACCGCATAACGTCGCAAGCACAAAGAGGGGGACCTTAGGGCCTCTTGCCCCCCCATGTGCCCAGATGGGATTAGCTAGTAGGTGGGGTAACGGCTCACCTAGGCGACGATCCCTAGCTGGTCTGAGAGGATGA >Q2 CTCAGATTGAACGCTGGCGGCAGGCCTAACACATGCAAGTCGAGCGGTAGAGAGAAGCTTGCTTCTCTTGAGAGCGGCGGACCCCCGAGTAATGCCTAGGAATCTGCCTGGTAGTGGGGGATAACGTTCGGAAACGGACGCTAATACCGCATACGTCCTACGGGAGAAAGCAGGGGACCTTCGGGCCTTGCGCTATCAGATGAGGGGGGGTCGGATTAGCTAGTTGGTGGGGTAATGGCTCACCAAGGCGACGATCCGTAACTGGTCTGAGAGGATGATCAGTCACACTGGAACTGA >Q3_exact AGTTTGATCATGGCTCAGATTGAACGCTGGCGGCAGGCCTAACACATGCAAGTCGAACGGTAACAGGAAGCAGCTTGCTGCTTTGCTGACGAGTGGCGGACGGGTGAGTAATGTCTGGGAAACTGCCTGATGGAGGGGGATAACTACTGGAAACGGTAGCTAATACCGCATAACGTCGCAAGCACAAAGAGGGGGACCTTAGGGCCTCTTGCCATCGGATGTGCCCAGATGGGATTAGCTAGTAGGTGGGGTAACGGCTCACCTAGGCGACGATCCCTAGCTGGTCTGAGAGGATGA >Q4_exact CTCAGATTGAACGCTGGCGGCAGGCCTAACACATGCAAGTCGAGCGGTAGAGAGAAGCTTGCTTCTCTTGAGAGCGGCGGACGGGTGAGTAATGCCTAGGAATCTGCCTGGTAGTGGGGGATAACGTTCGGAAACGGACGCTAATACCGCATACGTCCTACGGGAGAAAGCAGGGGACCTTCGGGCCTTGCGCTATCAGATGAGCCTAGGTCGGATTAGCTAGTTGGTGGGGTAATGGCTCACCAAGGCGACGATCCGTAACTGGTCTGAGAGGATGATCAGTCACACTGGAACTGA " > queries.fa ``` ## Blasting ### With few references With a small number of references we want to blast against, it's fine to just provide the references as a fasta file like so: ```bash blastn -query queries.fa -subject refs.fa > blast-results.txt ``` Typically we want things in table form, here's one way to run it for that and sort by bitscore: ```bash blastn -query queries.fa -subject refs.fa -outfmt "6 qseqid qlen sseqid slen length pident evalue bitscore" -max_hsps 1 -max_target_seqs 1 | sort -nrk 8 > blast-results.tmp # adding a header cat <(printf "qseqid\tqlen\tsseqid\tslen\tlength\tpident\tevalue\tbitscore\n") blast-results.tmp > blast-results.tsv rm blast-results.tmp column -ts $'\t' blast-results.tsv ``` ```bash # qseqid qlen sseqid slen length pident evalue bitscore # Q4_exact 297 NR_169460.1 297 297 100.000 5.96e-161 549 # Q3_exact 297 NR_024570.1 297 297 100.000 5.96e-161 549 # Q2 297 NR_169460.1 297 297 97.306 1.31e-147 505 # Q1 297 NR_024570.1 297 297 97.306 1.31e-147 505 ``` ### With many references When we have many references, it is more efficient to create a blast database like so: ```bash makeblastdb -dbtype nucl -in refs.fa -out refs-blast-db ``` Then we could run the blast like so: ```bash blastn -query queries.fa -db refs-blast-db -outfmt "6 qseqid qlen sseqid slen length pident evalue bitscore" -max_hsps 1 -max_target_seqs 1 | sort -nrk 8 > blast-results2.tmp # adding a header cat <(printf "qseqid\tqlen\tsseqid\tslen\tlength\tpident\tevalue\tbitscore\n") blast-results2.tmp > blast-results2.tsv rm blast-results2.tmp column -ts $'\t' blast-results2.tsv ``` ```bash # qseqid qlen sseqid slen length pident evalue bitscore # Q4_exact 297 NR_169460.1 297 297 100.000 5.96e-161 549 # Q3_exact 297 NR_024570.1 297 297 100.000 5.96e-161 549 # Q2 297 NR_169460.1 297 297 97.306 1.31e-147 505 # Q1 297 NR_024570.1 297 297 97.306 1.31e-147 505 ```