Facing problem with BLASTx using BioPython in Google Colab

Question

I'm trying to BLASTX a sequence using Biopython in Google Colab. Here is the code I am using. I'm not able to understand why I am getting an empty list [] for alignments.

result_handle = NCBIWWW.qblast("blastx", "nr", test_seq.seq, expect=1e-05)
blast_record = NCBIXML.read(result_handle)
alignments = []
for alignment in tqdm(blast_record.alignments):
  for hsp in alignment.hsps:
      alignments.append({
          "title": alignment.title,
          "length": alignment.length,
          "e_value": hsp.expect,
          "score": hsp.score,
          "query": hsp.query,
          "match": hsp.match,
          "subject": hsp.sbjct
      })
len(alignments) # its coming 0 because its empty

As a note, here is the test sequence I am using, which is giving proper results in BLASTX on NCBI server in NCBI NR database. Any help would be very helpful.

>test_seq
GCAGCAGCAGCAGCAGCAGCGCGGCGCGCGACTGCCTGGCCATGAATGCTTTCTCCAAGCTGCTGCGCGT
GCTGACCACGCCCTCCATGCTGCACAACCTGTGCAACAGCGACCAAGGCGTGCACGTGCTGATCCGGCTC
ATCAAGCTGAGCTCGCCGCCAGAGACGGCTTTCATTGCTGACCAGGTCCTGCACGACTTGGCTGACATGG
TGGTCAGCAAGGGCGGCGCTCTGGTGGTCTGCACCCTGGTCGACCGCCTGGCCGTCGCGGCCAGCTGCGA
CAACCCGGGGGTGGTGGAGCGGGCCACGGGCGCGGTTGTGCGCATGTGCAGCCACCTGTGCCAGCACCAT
CGCGCGCTGCAGCAGTGTGCGCGGCACCAGTACCTGGGTCCCAAGGTCGTTGAGCTGCTGCAGTGCGGCC
TGCCCGCTGTGGATGCCTTCAAGATGGGGTTTGAGCTGGTCCAGAACGCGCCCGCGCTGGCGGCCACCCG
CCACGGCAGCCGCACGCTGCAGCTGCTCACGGAGATTCGGGCCGCCGATGCGCGGTCGCGAGAGATGCTG
CGAGACGTGGTCTGCCACGCGTCCAAGGCGCTGCAGGGCCACTTTGGCCAGCTGGCGGCCAAGGGGGAGA
CGGGAGGGCCCGCCCTGGTCCGCGCCCTGCTGGACTACCTGCTGGCGGAGCGGGAGGACGACTGGGTGGT
GCAGATCACCTGGGAGCTGCTGGCGGCGGGCAGTGAGCTGCGGCGGGACCGGGAGGCACATGACTTGCTG
GTGCACAGCCTGGCCCTGCCGCTGTGGCGCCCCACCGAGGCGCAGCAGCACTGGGCGGCGCTTGCGGCCG
CTGGCGTGGATGTGGCAGCCCTCAAGCCCAAGGTGGATGCAGCCAAGCAGCACCAGCCGCGGCTTCCGGA
GCGGGTGCTGGCATCGCTGTCCGAGGTGCAGCCTTCGTCACCTGCAGAGGTCAGCGCCAAGCAGCGCACG
AGCTGTCTTCCCCGCTCCGCATGCTCAGGGCCGGGCTGCTCCTGCCTGCCCGGAGTGGTGAGCGGGATCC
TGACGTGCCTCCCTGCCCCTCGCTTGTCCAATGATGTTCTATATTCCATCGTCGCAGGCTGTGCGGCGCT
TCCCTTACTATCTGCGGGATGGAAGCCAGGCGGCTCGCAACTCTGTGCTCAGAAAGCAGCAGGGCGCAGG
GCAGCCAGCCGTGCCCCCGCCGCCCCCGCCGCGGGAGCAGCAGCAGCCCAAGGCGCCCCCTGGCTTTGGG
GCGGCAGGAGCCGCTGCAGCGCCAGCGGCTCCTGCGGCTGCCAAGGTGCCCCCACCCGGCTTTGCGGCAG
CGCCGTCGGCGGCGCCAGCAGCAGCGCCAGCGGCCATCACTTTCCCAAATGGCGCAAGCGTGCCCATGCC
GCAGCCCTATGGCGCCAGCAGCAGCAGCGCTGCGCCTGCAGCCAGCAAGGCGGTTCCGCCGCCAGCCCCA
CAGCAGCCAGCTGTCGTAGCCGCTGCTGCGATGGTGCCCCCGCCCTCGCCCCTCAAGCCTGCAGCCCCTC
CGCAGCAGCCAGTGGCACTGCCACCGCCTGCCCCCATGCAGCCAGCAGCAGCGCCCCTGGTGCCTGCACC
TGTGCCGCAGCCCGCGCCTGTGCAGCCACCAGTGGCAGCGCCCCGGCCCCCGGCACCTGTGCCGCAGCCT
GCTGTGCAGCCACCTGTGGCGGTGCGGCCACCAGTGCCACCACCTGCGCCTGTGGCACCCGGCATGAATG
CAGCAGCCCGGCCATTTGTGCAGGCAGTGCCTGCCCCTGCGGCGCAACCCCTGCCCCCGCACCTGCAGCC
CAACCTGCAGCAGCAGCAACAGCCGCAGCAGCCCGCGCACCTTCAAGCGCTGCCGCCGCACCTGCGGCAG
CCGGTGGCAGCGGCAGCAGGCAGCAGTGCCGGGCAGGCGGGCACGGCAGCACCCATCACGCCTGTGCCGG
CACAGCACCCTGCGGGCGGCTTCACCGCGCCGCCTGCACCGCAGCAGCTGCTGCAGCAGCTGCTTGGCAG
CGCTGCCCCGCGGCCCGCGGTCCCGGCAGCGGCGCCAGCAGCGGTGCGGCCGCTTGCTCCGCCGCAGGCA
GCAGCAGCTCAGGGCAGCGCAGCTGCCTTTGATGGCTACGGCAGCCTCACAGAAATGCTGTCTGAGCTAA
CTGCCGCTGCGGCGCCACAGCCGGTGGCACCTGCACGCCCAGCCTACCGCCCGGCCTACCGCCCGCCTGG
AGCAGCTGCTTACCGCCCGCCCACCGCCGCGGCGCCGCCGCCACCGCCCGCGCCAGCGGCCGCACCTGCG
CCTGTTTCTGTGCCTGTGGGCATTCCAGCCCCCTCGCCTGCAACCCGGTTCCCCATGCCCGGCGCGCCAC
GGCCCGCACCGGTGCCAGCAGCTGCGGCGCCGCCGCCGCCGGCAGCGCCCGCCCGCCCGGCGGTGCCAGC
GGTGCCGGTGGCATCGGCGACCAACTGGACCTGCAGGGTGTGCACCTACCGGCATGAGGGAGCGGAGGCT
CAGTTCCTGTCCTGTGCCATGTGCGGCTCCCCACAGTCCTAGGCAGCCTGCTGTGGCAGTCCCTAGGCTG
GCTGCTGTGGCAGTCCTGGGCAGGCCGCTTGGTGGCTAGCTGGCTTCGCTGGCTCGGCAGCACTTCAACG
CTCCCACAACCTTCCATTGTTCCTTGCACTGCAACCCTTGCCCCCTGATCGCATCGCTCGCCCTTGCCTC
TTGTCATTCATGCTGCACCAACCCAGCCCCTGTGCATGTAACACCAGAGATGACCCCCCTGTTGTATTGG
GCTCGCTTGCTTTGAGGCACGGGGGTTGTTTCACAGGAGGAATGAACAGGGGGGCACATGGGACGGCTGG
GAAGAAGCAGCACAAAGAGCAACAGGAACAGAATCATCAAATGAGCATGTGGTGGTACACTTGACTTCGG
GGGATGAAGAAGAAGAGGCAGCCAGGGGTGCTGGGGCGGAGGGGGGTGTCGATCCGAGGCATCATGCCGC
CAGCGGTAAGGGTGCAAGCAGTGGGGGTGCGTATTTAGACAGCAGGGCATTCGTGGGGCAAGGAAGAGAT
GGAGGAGTGTGGCGGGGCTGGCAGGGGGTTCAATGCAGCCTTCCTGCCGCACACTTTGCAGACATGCATT
TTCCAGGGAAGAGACGAGAGGAAGAGCAGTGCGACGGGGCGGGACAGCGAAGGGGCGGGCGCACGGGGGC
GCAGTCCCGCCCATGGCGCAGCCGCCACGGCCGGAACAAACATTTTGTTTTTGATGGTACACAAACAAAC
AAACAGGAGCGATGTTAAGAGTTAGGGAGCGGAGCAATGCAGAGGTGGGTGGGCAGGTAGCAAGGCAGGG
AGCGGGTATATGCCGGGGGCGGGCCACAGTGCAGTGGCGGCATGGCGGCAGCGCCACGCTGCCGCCGCCC
ACACCGACGGCTGGGCGCCGCCACTTCAGTCGGCGAGCGTCGCGCGCGCCAGCTGCGCTTGCACAGCCGC
CGTGCCGCCGCCCTGCTCGCCAGATGCAGTGCGCACCGGCTGGGCAGCGGAAGCAGACTCCTCAGACTGC
TCCAGAGAGCTGGCAGAGTCCACCTTGCCCTTGCTGCGGCTGCCGTTGCGGCTCAGCATGTGCAGCGGCC
TGCGGCGAAGCAAGGGCGGTGCGGCGGGGGCCGCAGCCAGGGTCGGCATCCAAGCCGCAGCCGGCAGCAG
CGTGCGGTACAGAGGGCGGCTGGCCTCGGCCTGCCGGCGCTGCCAAGCATGGGTGCTATGAGGAAGCAGG
CCTGCAGAGTTGCAGGCCGGCAGGCAGGGCAGCGCCAATCAGGGCGGGCACCAGCAGTGCCAGGCTCGCC
AGCAACAAGCTGCAGCCCTTTGCCACACTCACCGCGTGTACGGCGGCTGCGACTGCTGCGACTGCTGCTG
CGCAAACGGCGTGCCGCTGTAGTGGCCGCCGCCGCCGCCGCCAGCGGGGTGCTGCGGCGGGAATGCCGGC
GATGCCATGCCGGTGCTTGTGGGGTAGCCGTAGGCAGCTGCGCCGGGCACATAGCCGCCCATCTGGCCTC
CCATCATCATGTACGGGCTGGGGGCCATGCCCCCCATGCCAGCCATTGGGCTCATGCCACCGTACTGCTG
GGCCATGCCGGCATACTGCATGCCAAAGCCGGGGCCGCCCGGCACCATGCCTGCATGCGGGTAGCCGAAC
ATGCCAAAGCCGGCGGCCGCGTAGGGGCTGGGCGGGCTGCCAAACAGCGCGCCGGGTGAGCCCACCACGG
CGCCGGGCGGGCCGCCAAAGGGCGGCAGGCCCATGCCGCCGCGGCCGCCGGCGCCGCTGGTGCTG

EDIT: By the way, to test the code, I tried running Biopython in Jupyter locally. When I tried with BLASTp with same code for a different sequence, I got this error:

URLError: <urlopen error [Errno 11001] getaddrinfo failed>

I'm thinking there is some problem with the scraper implemented in Biopython. When attempting BLASTX, the code keeps on running for an unreasonable amount of time (>10 min. I had to stop the Kernel). Let me know if it makes sense.

Thanks

Answer 1

I've been able to replicate the issue you've had (or at least produce a similar issue) with the following code:

#!/usr/bin/env python

## initial setup code from 'www_blast' biopython example
## (Brad Chapman / Biopython)

# standard library
from io import StringIO

# biopython
from Bio import SeqIO
from Bio.Blast import NCBIWWW
from Bio.Blast import NCBIXML

from tqdm import tqdm # progress bar

# first get the sequence we want to parse from a FASTA file
test_seq = next(SeqIO.parse("test_seq.fasta", "fasta"))

# this fetch query takes a *long* time
result_handle = NCBIWWW.qblast("blastn", "nr",
                               test_seq.format("fasta"),
                               expect=1e-05)

# save the results for later, in case we want to look at it
with open("m_cold_blast.out", "w") as save_file:
    blast_results = result_handle.read()
    save_file.write(blast_results)

blast_record = NCBIXML.read(result_handle)

alignments = []
for alignment in tqdm(blast_record.alignments):
  for hsp in alignment.hsps:
      alignments.append({
          "title": alignment.title,
          "length": alignment.length,
          "e_value": hsp.expect,
          "score": hsp.score,
          "query": hsp.query,
          "match": hsp.match,
          "subject": hsp.sbjct
      })
len(alignments) # result: 0

Looking at the message in the saved XML output - I'd recommend doing similar with saving and looking at intermediate data when debugging issues - I see the following:

<Iteration_message>[blastsrv4.REAL]: Error: CPU usage limit was exceeded, resulting in SIGXCPU (24).</Iteration_message>

So my guess is that your query string is too long.

After further investigation, I've come up with the following code, which is able to properly process a shorter string (280bp):

#!/usr/bin/env python

## initial setup code from 'www_blast' biopython example
## (Brad Chapman / Biopython)

# standard library
from io import StringIO

# biopython
from Bio import SeqIO
from Bio.Blast import NCBIWWW
from Bio.Blast import NCBIXML

from tqdm import tqdm # progress bar

# first get the sequence we want to parse from a FASTA file
test_seq = next(SeqIO.parse("test_seq_short.fasta", "fasta"))

# this takes a long time
result_handle = NCBIWWW.qblast("blastx", "nr",
                               test_seq.format("fasta"),
                               expect=1e-05)

# save the results for later, in case we want to look at it
with open("m_cold_blast.out", "w") as save_file:
    blast_results = result_handle.read()
    save_file.write(blast_results)

string_result_handle = StringIO(blast_results)
blast_record = NCBIXML.read(string_result_handle)

alignments = []
for alignment in tqdm(blast_record.alignments):
  for hsp in alignment.hsps:
      alignments.append({
          "title": alignment.title,
          "length": alignment.length,
          "e_value": hsp.expect,
          "score": hsp.score,
          "query": hsp.query,
          "match": hsp.match,
          "subject": hsp.sbjct
      })

>>> len(alignments)
4
>>> alignments
[{'title': 'gb|PRW60132.1| 60S ribosomal L15 [Chlorella sorokiniana]', 'length': 1087, 'e_value': 2.9921e-51, 'score': 466.0, 'query': 'SSSSSSAARDCLAMNAFSKLLRVLTTPSMLHNLCNSDQGVHVLIRLIKLSSPPETAFIADQVLHDLADMVVSKGGALVVCTLVDRLAVAASC', 'match': 'SSSSSSAARDCLAMNAFSKLLRVLTTPSMLHNLCNSDQGVHVLIRLIKLSSPPETAFIADQVLHDLADMVVSKGGALVVCTLVDRLAVAASC', 'subject': 'SSSSSSAARDCLAMNAFSKLLRVLTTPSMLHNLCNSDQGVHVLIRLIKLSSPPETAFIADQVLHDLADMVVSKGGALVVCTLVDRLAVAASC'}, {'title': 'gb|KAI7846133.1| hypothetical protein COHA_000394 [Chlorella ohadii]', 'length': 864, 'e_value': 1.83868e-41, 'score': 394.0, 'query': 'AARDCLAMNAFSKLLRVLTTPSMLHNLCNSDQGVHVLIRLIKLSSPPETAFIADQVLHDLADMVVSKGGALVVCTLVDRLAVAASC', 'match': 'A  DCLAM +FSKLLRVLT+ SMLHNLC SDQGVHVL+RLIKLSSPPETAFIADQV+HDLADMV SKGGALVVCTLVDRLAVAASC', 'subject': 'AGHDCLAMTSFSKLLRVLTSASMLHNLCMSDQGVHVLVRLIKLSSPPETAFIADQVMHDLADMVTSKGGALVVCTLVDRLAVAASC'}, {'title': 'ref|XP_005845358.1| hypothetical protein CHLNCDRAFT_137175 [Chlorella variabilis] >gb|EFN53256.1| hypothetical protein CHLNCDRAFT_137175 [Chlorella variabilis]', 'length': 1059, 'e_value': 2.1613e-11, 'score': 171.0, 'query': 'CLAMNAFSKLLRVLTTPSMLHNLCNSDQGVHVLIRLIKLSSPPETAFIADQVLHDLADMVVSKGGALVVCTLVDRLAVAAS', 'match': 'C A+ A+S+LL+ LT P+M H L  SD+GV+VL++L+KL SP ET +IA+Q    + +   SK G+L V  LVDRLA  A+', 'subject': 'CTALRAYSELLKALT-PAM-HQLATSDRGVYVLVKLVKLGSPMETVYIAEQATPMVVEWSGSKDGSLFVTQLVDRLAQVAT'}, {'title': 'gb|KAI3436075.1| hypothetical protein D9Q98_002133 [Chlorella vulgaris]', 'length': 1091, 'e_value': 4.23077e-06, 'score': 132.0, 'query': 'CLAMNAFSKLLRVLTTPSMLHNLCNSDQGVHVLIRLIKLSSPPETAFIADQVLHDLADMVVSKGGALVVCTLVDRLAVAAS', 'match': 'C AM A+S+LL+ L     +HNL   D+GVH L++L+   S  E  ++A+Q    + D + SK G L    L+ RLA +A+', 'subject': 'CTAMMAYSELLKALAL--AMHNLSMRDKGVHALVKLVDTGSTLEAVYLAEQATPLMLDWMGSKDGILFAVKLIHRLAQSAA'}]

Comparing with the code you have provided, the only difference I can see is that my working result text is wrapped inside a StringIO Structure.