Answers.Body: <p>Install (pip) biopython if u don't have it.</p>
  <pre><code>from Bio import SeqIO
  with open('/path/ids.txt', 'r') as fin:
      to_keep = set([])
      for line in fin:
          to_keep.add(line.strip().split(' ')[0])
  with open('/path/output.fa', 'w') as fout:
      for record in SeqIO.parse('/path/sequences.fa','fasta'):
          if record.id in to_keep:
              SeqIO.write(record, fout, 'fasta')
  </code></pre>
  
Answers.ContentLicense: CC BY-SA 4.0
Answers.CreationDate: 2021-02-04T20:55:08.470Z
Answers.LastActivityDate: 2021-02-04T20:55:08.470Z
Answers.UserDisplayName: Liam McIntyre
Answers.Body: <p>I tend to do this using the <a href="https://bioinformatics.stackexchange.com/a/2651/298">FastaToTbl and TblToFasta scripts I have posted before</a>. <code>FastaToTbl</code> will convert fastq sequences to single-line, tab-separated recors which makes them far easier to parse. Using these two scripts, you can do:</p>
  <pre><code>$ FastaToTbl file.fa | grep -vwf idList | TblToFasta 
  &gt;m64041_200717_231916/77/1941_50622 id=3115 
  TGAGCCGCACCTCCTCAATCTGCAGCAGTTGAGGCCACTACCCTTCTGCTCAATGGTTCC
  TGCAGACTTTATCATAGTCACTCACACTTGTCCATGAG
  &gt;m64041_200717_231916/162/0_6615 id=857 
  CAGACTTTATCATAGTCACTCACACTTGTCCATGAGGCAGACTTTATCATAGTCACTCAC
  ACTTG
  </code></pre>
  <p>The <code>-v</code> option will tell <code>grep</code> to return non-matching lines, and the <code>-f</code> option tells it to read the patterns to search for from a file. In this case, we give it the file with the IDs to skip (<code>idList</code>). The <code>-w</code> ensures we only find &quot;whole word&quot; matches, so that something like <code>m64041_200717_231916/77/1941_50622 id=31151111</code> will not count as a match for <code>m64041_200717_231916/77/1941_50622 id=31151111</code>.</p>
  <hr />
  <p>Another approach would be to just use <code>awk</code> directly:</p>
  <pre><code>$ awk '(NR==FNR){s[&quot;&gt;&quot;$0]++} (/^&gt;/){ s[$0] ? a=0 : a=1}a' idList file.fa 
  &gt;m64041_200717_231916/77/1941_50622 id=3115
  TGAGCCGCACCTCCTCAATCTGCAGCAGTTGAGGCCACTACCCTTCTGCTCAATGGTTCCTGCAGACTTTATCATAGTCACTCACACTTGTCCATGAG
  &gt;m64041_200717_231916/162/0_6615 id=857
  CAGACTTTATCATAGTCACTCACACTTGTCCATGAGGCAGACTTTATCATAGTCACTCACACTTG
  </code></pre>
  <p>Here is the same <code>awk</code> script slightly expanded to make it easier to understand:</p>
  <pre><code> awk '{
          if (NR==FNR){idsToSkip[&quot;&gt;&quot;<span class="math-container">$0]++} 
          else if (/^&gt;/){ 
              if($</span>0 in idsToSkip){
                  wantThisId = 0
              }
              else{
                  wantThisId = 1
              }
          }
          if(wantThisId == 1){
              print
          }
        }' idList file.fa 
  </code></pre>
  
Answers.ContentLicense: CC BY-SA 4.0
Answers.CreationDate: 2021-02-05T13:55:38.527Z
Answers.LastActivityDate: 2021-02-05T13:55:38.527Z
Answers.UserDisplayName: terdon
Body: <p>I have a fasta file that contains sequence reads and sequence id file that needed to be removed from the fasta file. I have done this earlier, but since id contains a space my piece of code is not working and I don't know how to change it since I'm still learning. fasta file looks like this.</p>
  <pre><code>&gt;m64041_200717_231916/74/12723_24868 id=30
  CCGCACCTCCTCAATCTGCAGCAGTTGAGGCCACTACCCTTCTGCTCAATGGTTCCTGCAGACTTTATCATAGTCACTCACACTTGTCCA
  &gt;m64041_200717_231916/77/1941_50622 id=3115
  TGAGCCGCACCTCCTCAATCTGCAGCAGTTGAGGCCACTACCCTTCTGCTCAATGGTTCCTGCAGACTTTATCATAGTCACTCACACTTGTCCATGAG
  &gt;m64041_200717_231916/105/20691_65844 id=488
  AGGCCACTACCCTTCTGCTCAATGGTTCCTGCAGACTTTATCATAGTCACTCACACTTGTCCAAGACTTTAT
  &gt;m64041_200717_231916/108/17414_67048 id=4956
  TGAGGCCACTACCCTTCTGCTCAATGGTTCCTGCAGACTTTATCATAGTCACTCACACTTGTCCATGAGGCAGACTTTATCATAGTCACTCACACTTG
  &gt;m64041_200717_231916/162/0_6615 id=857
  CAGACTTTATCATAGTCACTCACACTTGTCCATGAGGCAGACTTTATCATAGTCACTCACACTTG
  </code></pre>
  <p>IDList file is like this.</p>
  <pre><code>m64041_200717_231916/74/12723_24868 id=30
  m64041_200717_231916/108/17414_67048 id=4956
  m64041_200717_231916/105/20691_65844 id=488
  </code></pre>
  <p>I need to remove these sequences from that fasta file. I have a code I used earlier. But it is not working since id contains a space (I think). Previous code is like this.</p>
  <pre><code># use as follows
  # ./filter_seq.bash &lt;fastafile&gt; &lt;list_of_ids_you_want_to_remove&gt;
  # ./filter_seq.bash all.fasta bacterial_IDlist
  
  grep &quot;&gt;&quot; <span class="math-container">$1 | sed 's/&gt;//g' | sort | uniq &gt; all_del
  sort $</span>2 | uniq &gt; query_del
  comm -23 all_del query_del | sort | uniq &gt; filtered_del
  perl -ne 'if(/^&gt;(\S+)/){<span class="math-container">$c=$</span>i{<span class="math-container">$1}}$</span>c?print:chomp;<span class="math-container">$i{$_}=1 if @ARGV' filtered_del $</span>1 &gt; filtered.fasta
  rm *_del
  </code></pre>
  <p>I would be grateful if someone can help me.</p>
  
Comments.ContentLicense: CC BY-SA 4.0
Comments.CreationDate: 2021-02-04T18:49:32.340Z
Comments.Text: I would highly recommend getting rid of the spaces first using something like `sed 's, ,_,g' -I in.fa`. Spaces in any kind of genetic data file are bad news and will cause downstream errors. It will make further analysis much easier.
Comments.UserDisplayName: user438383
Comments.ContentLicense: CC BY-SA 4.0
Comments.CreationDate: 2021-02-05T14:02:28.323Z
Comments.Text: @user438383 why do you say that? Spaces in fasta headers are incredibly common, allowed by the standard and should not cause any issues at all. I can't think of any tool that would have trouble parsing a fasta header like that and would even say that if a tool cannot do this, it is not fit for purpose. As an example, have a look at pretty much any fasta sequence from refseq. Here's human tp53: https://www.ncbi.nlm.nih.gov/nuccore/NG_017013.2?report=fasta.
Comments.UserDisplayName: terdon
Tags: <bash><linux><filtering>
Title: How to remove sequences from a fasta file using a sequence ID list which contains a space within the id?
ThreadText: # How to remove sequences from a fasta file using a sequence ID list which contains a space within the id?
  
  (Asked by: username_0 on 2021-02-04)
  
  I have a fasta file that contains sequence reads and sequence id file that needed to be removed from the fasta file. I have done this earlier, but since id contains a space my piece of code is not working and I don't know how to change it since I'm still learning. fasta file looks like this.
  
  ```
  >m64041_200717_231916/74/12723_24868 id=30
  CCGCACCTCCTCAATCTGCAGCAGTTGAGGCCACTACCCTTCTGCTCAATGGTTCCTGCAGACTTTATCATAGTCACTCACACTTGTCCA
  >m64041_200717_231916/77/1941_50622 id=3115
  TGAGCCGCACCTCCTCAATCTGCAGCAGTTGAGGCCACTACCCTTCTGCTCAATGGTTCCTGCAGACTTTATCATAGTCACTCACACTTGTCCATGAG
  >m64041_200717_231916/105/20691_65844 id=488
  AGGCCACTACCCTTCTGCTCAATGGTTCCTGCAGACTTTATCATAGTCACTCACACTTGTCCAAGACTTTAT
  >m64041_200717_231916/108/17414_67048 id=4956
  TGAGGCCACTACCCTTCTGCTCAATGGTTCCTGCAGACTTTATCATAGTCACTCACACTTGTCCATGAGGCAGACTTTATCATAGTCACTCACACTTG
  >m64041_200717_231916/162/0_6615 id=857
  CAGACTTTATCATAGTCACTCACACTTGTCCATGAGGCAGACTTTATCATAGTCACTCACACTTG
  
  ```
  
  IDList file is like this.
  
  ```
  m64041_200717_231916/74/12723_24868 id=30
  m64041_200717_231916/108/17414_67048 id=4956
  m64041_200717_231916/105/20691_65844 id=488
  
  ```
  
  I need to remove these sequences from that fasta file. I have a code I used earlier. But it is not working since id contains a space (I think). Previous code is like this.
  
  ```
  # use as follows
  # ./filter_seq.bash <fastafile> <list_of_ids_you_want_to_remove>
  # ./filter_seq.bash all.fasta bacterial_IDlist
  
  grep ">" $1 | sed 's/>//g' | sort | uniq > all_del
  sort $2 | uniq > query_del
  comm -23 all_del query_del | sort | uniq > filtered_del
  perl -ne 'if(/^>(\S+)/){$c=$i{$1}}$c?print:chomp;$i{$_}=1 if @ARGV' filtered_del $1 > filtered.fasta
  rm *_del
  
  ```
  
  I would be grateful if someone can help me.
  
  **Tags:** bash, linux, filtering
  
  > **Comments:**
  > **username_1**: I would highly recommend getting rid of the spaces first using something like `sed 's, ,_,g' -I in.fa`. Spaces in any kind of genetic data file are bad news and will cause downstream errors. It will make further analysis much easier.  
  > **username_2**: @username_1 why do you say that? Spaces in fasta headers are incredibly common, allowed by the standard and should not cause any issues at all. I can't think of any tool that would have trouble parsing a fasta header like that and would even say that if a tool cannot do this, it is not fit for purpose. As an example, have a look at pretty much any fasta sequence from refseq. Here's human tp53: https://www.ncbi.nlm.nih.gov/nuccore/NG_017013.2?report=fasta.  
  
  ---
  
  ## 2 answers
  
  ### Answer by username_3 on 2021-02-04. Score: 1
  
  Install (pip) biopython if u don't have it.
  
  ```
  from Bio import SeqIO
  with open('/path/ids.txt', 'r') as fin:
      to_keep = set([])
      for line in fin:
          to_keep.add(line.strip().split(' ')[0])
  with open('/path/output.fa', 'w') as fout:
      for record in SeqIO.parse('/path/sequences.fa','fasta'):
          if record.id in to_keep:
              SeqIO.write(record, fout, 'fasta')
  
  ```
  
  
  ---
  
  ### Answer by username_2 on 2021-02-05. Score: 0
  
  I tend to do this using the [FastaToTbl and TblToFasta scripts I have posted before](https://bioinformatics.stackexchange.com/a/2651/298). `FastaToTbl` will convert fastq sequences to single-line, tab-separated recors which makes them far easier to parse. Using these two scripts, you can do:
  
  ```
  $ FastaToTbl file.fa | grep -vwf idList | TblToFasta 
  >m64041_200717_231916/77/1941_50622 id=3115 
  TGAGCCGCACCTCCTCAATCTGCAGCAGTTGAGGCCACTACCCTTCTGCTCAATGGTTCC
  TGCAGACTTTATCATAGTCACTCACACTTGTCCATGAG
  >m64041_200717_231916/162/0_6615 id=857 
  CAGACTTTATCATAGTCACTCACACTTGTCCATGAGGCAGACTTTATCATAGTCACTCAC
  ACTTG
  
  ```
  
  The `-v` option will tell `grep` to return non-matching lines, and the `-f` option tells it to read the patterns to search for from a file. In this case, we give it the file with the IDs to skip (`idList`). The `-w` ensures we only find "whole word" matches, so that something like `m64041_200717_231916/77/1941_50622 id=31151111` will not count as a match for `m64041_200717_231916/77/1941_50622 id=31151111`.
  
  ---
  
  Another approach would be to just use `awk` directly:
  
  ```
  $ awk '(NR==FNR){s[">"$0]++} (/^>/){ s[$0] ? a=0 : a=1}a' idList file.fa 
  >m64041_200717_231916/77/1941_50622 id=3115
  TGAGCCGCACCTCCTCAATCTGCAGCAGTTGAGGCCACTACCCTTCTGCTCAATGGTTCCTGCAGACTTTATCATAGTCACTCACACTTGTCCATGAG
  >m64041_200717_231916/162/0_6615 id=857
  CAGACTTTATCATAGTCACTCACACTTGTCCATGAGGCAGACTTTATCATAGTCACTCACACTTG
  
  ```
  
  Here is the same `awk` script slightly expanded to make it easier to understand:
  
  ```
   awk '{
          if (NR==FNR){idsToSkip[">"$0]++} 
          else if (/^>/){ 
              if($0 in idsToSkip){
                  wantThisId = 0
              }
              else{
                  wantThisId = 1
              }
          }
          if(wantThisId == 1){
              print
          }
        }' idList file.fa 
  
  ```
  
  
  ---
  
  