Answers.Body: <p><strong>UPDATE:</strong> while procrastinating, I worked up the below throwaway code into a more generally applicable (but still bare-bones) <a href="https://github.com/maximilianpress/teaching_code/blob/master/exact_match.py" rel="nofollow noreferrer">script</a>.</p>
  <hr />
  <p><strong>ORIGINAL:</strong></p>
  <p>I feel like you might want to use something other than BLAST to do this.</p>
  <p>BLAST is not really designed to solve this problem, I think it's just easy to use. Also, your N-containing sequence (16nt) is shorter than the version without Ns (21 nt), which can't be correct if you care about &quot;exact&quot; matches.</p>
  <p>As a side note, the idea of exact matches is a bit tricky with Ns involved. The whole point of ambiguity codes is that things won't be exact.</p>
  <p>And all of this is on top of the short-string problem- most alignment programs won't even return results that don't meet some minimum alignment score, which is exacerbated a bit by the Ns. For how BLAST handles ambiguity codes, you can see section 7.2 <a href="https://www.ncbi.nlm.nih.gov/blast/BLAST_guide.pdf" rel="nofollow noreferrer">here</a>.</p>
  <p>Why not use regular expressions instead to search your sequence against the transcriptome?</p>
  <pre><code># python
  import re
  from Bio import SeqIO
  
  fwd = r'ATTG[ACGT]{4}GCAAACCA'
  rev = r'TGGTTTGC[ACGT]{4}CAAT'  # reverse complement
  
  # reference file
  seqfile = SeqIO.parse(&quot;your_reference_file.fa&quot;, &quot;fasta&quot;)  # a transcriptome reference fasta
  
  for seq in seqfile:
      # maybe use re.match() to get coords but easier to show it with findall
      a = re.findall(fwd, str(seq.seq))
      b = re.findall(rev, str(seq.seq))
      if len(a) &gt; 0:
          print(seq.id, a)
      if len(b) &gt; 0:
          print(seq.id, b)
  
  </code></pre>
  <p>When I put this in a file and run it against a test case with <code>ATTGAAAAGCAAACCA</code>, it works:</p>
  <pre><code>$ python exact_match.py 
  the_contig ['ATTGAAAAGCAAACCA']
  </code></pre>
  <p>Alternately, you could just blast the 16 non-N versions of that sequence and combine them, as you've shown.</p>
  
Answers.ContentLicense: CC BY-SA 4.0
Answers.CreationDate: 2020-09-26T10:15:53.720Z
Answers.LastActivityDate: 2020-09-26T21:00:54.343Z
Answers.LastEditDate: 2020-09-26T21:00:54.343Z
Answers.UserDisplayName: Maximilian Press
Body: <p>I'm trying to Blast for exact matches to the sequence: 'ATTGNNNNGCAAACCA' in the human transcriptome using NCBI Blast on its 'refseq_rna' database. However, when I do a basic query I get &quot;No significant similarity found.&quot;</p>
  <p>Note the N's in the middle of the sequence, if I get rid of those N's, yielding the sequence 'AGCGGATTGCAAAGCAAACCA', I get a match in the 3' UTR of the human MeCP2 gene (which is correct). I don't understand why adding the N's here make it not work.</p>
  <p>I looked at the help section and I believe this may be due to shorter sequences being less &quot;statistically significant&quot;. Statistical significance of output is completely irrelevant to me, so I set 'Expect Value' to 100000000000000. However, I still get the same result.</p>
  <p>I also tried to modify the HTML of the page to allow me to submit a word size of 4 to the form. This went through, but I still got &quot;No significant similarity found&quot;.</p>
  <p>Could someone help me out with this? I feel like searching for a short sequence shouldn't be this hard.</p>
  
Tags: <blast>
Title: NCBI blast for exact match of a short sequence
ThreadText: # NCBI blast for exact match of a short sequence
  
  (Asked by: username_0 on 2020-09-25)
  
  I'm trying to Blast for exact matches to the sequence: 'ATTGNNNNGCAAACCA' in the human transcriptome using NCBI Blast on its 'refseq\_rna' database. However, when I do a basic query I get "No significant similarity found."
  
  Note the N's in the middle of the sequence, if I get rid of those N's, yielding the sequence 'AGCGGATTGCAAAGCAAACCA', I get a match in the 3' UTR of the human MeCP2 gene (which is correct). I don't understand why adding the N's here make it not work.
  
  I looked at the help section and I believe this may be due to shorter sequences being less "statistically significant". Statistical significance of output is completely irrelevant to me, so I set 'Expect Value' to 100000000000000. However, I still get the same result.
  
  I also tried to modify the HTML of the page to allow me to submit a word size of 4 to the form. This went through, but I still got "No significant similarity found".
  
  Could someone help me out with this? I feel like searching for a short sequence shouldn't be this hard.
  
  **Tags:** blast
  
  
  ---
  
  ## 1 answers
  
  ### Answer by username_1 on 2020-09-26. Score: 2 | Accepted answer
  
  **UPDATE:** while procrastinating, I worked up the below throwaway code into a more generally applicable (but still bare-bones) [script](https://github.com/username_1/teaching_code/blob/master/exact_match.py).
  
  ---
  
  **ORIGINAL:**
  
  I feel like you might want to use something other than BLAST to do this.
  
  BLAST is not really designed to solve this problem, I think it's just easy to use. Also, your N-containing sequence (16nt) is shorter than the version without Ns (21 nt), which can't be correct if you care about "exact" matches.
  
  As a side note, the idea of exact matches is a bit tricky with Ns involved. The whole point of ambiguity codes is that things won't be exact.
  
  And all of this is on top of the short-string problem- most alignment programs won't even return results that don't meet some minimum alignment score, which is exacerbated a bit by the Ns. For how BLAST handles ambiguity codes, you can see section 7.2 [here](https://www.ncbi.nlm.nih.gov/blast/BLAST_guide.pdf).
  
  Why not use regular expressions instead to search your sequence against the transcriptome?
  
  ```
  # python
  import re
  from Bio import SeqIO
  
  fwd = r'ATTG[ACGT]{4}GCAAACCA'
  rev = r'TGGTTTGC[ACGT]{4}CAAT'  # reverse complement
  
  # reference file
  seqfile = SeqIO.parse("your_reference_file.fa", "fasta")  # a transcriptome reference fasta
  
  for seq in seqfile:
      # maybe use re.match() to get coords but easier to show it with findall
      a = re.findall(fwd, str(seq.seq))
      b = re.findall(rev, str(seq.seq))
      if len(a) > 0:
          print(seq.id, a)
      if len(b) > 0:
          print(seq.id, b)
  
  ```
  
  When I put this in a file and run it against a test case with `ATTGAAAAGCAAACCA`, it works:
  
  ```
  $ python exact_match.py 
  the_contig ['ATTGAAAAGCAAACCA']
  
  ```
  
  Alternately, you could just blast the 16 non-N versions of that sequence and combine them, as you've shown.
  
  
  ---
  
  