<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>

CC BY-SA 4.0
2021-02-04T20:55:08.470Z
2021-02-04T20:55:08.470Z
Liam McIntyre
<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>

CC BY-SA 4.0
2021-02-05T13:55:38.527Z
2021-02-05T13:55:38.527Z
terdon
<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>

CC BY-SA 4.0
2021-02-04T18:49:32.340Z
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.
user438383
CC BY-SA 4.0
2021-02-05T14:02:28.323Z
@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.
terdon
<bash><linux><filtering>
How to remove sequences from a fasta file using a sequence ID list which contains a space within the id?
# 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  
# ./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**: @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.  

---

## 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. 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


---

