<p>You can do this using <a href="https://bioinf.shenwei.me/seqkit/" rel="nofollow noreferrer"><code>seqkit</code></a> as follows:</p>
<pre><code>seqkit grep -r -n -p '.*Pseudomonas.*' temp.fa
</code></pre>
<p>To explain a little, <code>seqkit grep</code> will allow you to search FASTA/Q files by sequence name or sequence itself. In this instance:</p>
<ul>
<li><code>-r</code> tells that the pattern is a regular expression</li>
<li><code>-n</code> to match by full name instead of just id</li>
<li><code>-p</code> to specify the regular expression pattern to search</li>
</ul>
<p>There are a bunch of other tools that come with <code>seqkit</code> that can prove indispensable for basic sequence analyses.</p>

CC BY-SA 4.0
2021-02-19T13:08:43.647Z
Hi @vkkodali, thank you for your time posting it. I did not know about seqkit (living and learning). I'm testing it out and looks like a great tool. Thanks again.
Fabricio
CC BY-SA 4.0
2021-02-19T13:46:59.780Z
Hi again @vkkodali, another question (sorry, I'm a novice in the world of bioinformatics). Is there a way to retrieve the whole sequence header or ID using seqkit? I filtered the sequences that belong to Pseudomonas and the fasta file contains 38K entries of Pseudomonas. I would like to extract the sequence ID to count how many species of Pseudomonas are in this file (38K entries).

illustration:

>AB001446.1.1538 Pseudomonas amygdali pv. mori
AACTGAAGAGTT...

>AB001440.1.1538 Pseudomonas coronafaciens pv. atropurpurea
AACTGAAGAGTTT...
> ...

Thank you for your time and attention.
Fabricio
CC BY-SA 4.0
2021-02-19T17:01:50.373Z
So you just want the `AB001446.1.1538` portion of the header? So many ways to do that... Easiest would be something like `grep '^>' in.fa | cut -f1 -d ' ' | sed 's/>//'`. Another way is to do `seqkit faidx in.fa` and that creates a `in.fa.fai` which is essentially a tab-delimited file with the first column containing the identifiers.
vkkodali
CC BY-SA 4.0
2021-02-19T20:22:37.973Z
Hi, can be the whole header or portion of it. The portion is the genus and species name (e.g., Pseudomonas putid). I would like to extract the species name from the filtered Psedumonas file (e.g., Pseudomonas amygdali pv. mori, Pseudomonas coronafaciens pv. atropurpurea etc.). I would like to know how many different species of Pseudomonas are there. Thanks again for your time. Fabricio
Fabricio
CC BY-SA 4.0
2021-02-19T20:58:59.913Z
This can be done using a combination of `grep` and `sed` but I'd be making assumptions about how the organism names in the headers for all your 9k sequences are formatted. Perhaps you should describe the entire issue with an example of what you are starting with and what you would like to see in the end and I can update my answer?
vkkodali
CC BY-SA 4.0
2021-02-22T16:37:39.873Z
Hi @vkkodali, due to the limitation of characters in the comment section here. I wrote an answer which actually is a follow-up (explanation) for the same issue. Thank you for your time. Fabricio.
Fabricio
CC BY-SA 4.0
2021-02-18T02:01:19.337Z
2021-02-18T02:01:19.337Z
vkkodali
<p>Using Biopython:</p>
<pre><code>from Bio import SeqIO

for record in SeqIO.parse(&quot;input.fa&quot;, &quot;fasta&quot;):
    if record.description.split()[-2:] == [&quot;Pseudomonas&quot;, &quot;aeruginosa&quot;]:
        print(record.format(&quot;fasta&quot;))
</code></pre>

CC BY-SA 4.0
2021-02-17T19:47:16.303Z
Hi Chris, right after the post I realize that the issue was not properly described. So, this database has ~9k entries of Pseudomonas species, i.e., hundreds of P. species. The command above asks for the name of each Pseudomonas entry, correct ("Pseudomonas", "X")? Is there a way to retrieve the Pseudomonas entries without specifying all Pseudomonas names (I do not have a list name)? Thanks.
Fabricio
CC BY-SA 4.0
2021-02-17T19:50:49.113Z
@Fabricio Hi Fabricio, do you understand the above code? If so you can easily modify the `if` line to search for entries with "Pseudomonas"
Chris_Rands
CC BY-SA 4.0
2021-02-17T21:19:49.307Z
@Fabricio Add this at the top of the code, after the other `import` statement: `import re` . Change the `if` expression to this: `if re.findall(r'pseudomonas', record.description, re.I):` , in order to search in the fasta header for the word `pseudomonas`, in case-insensitive manner.
Timur Shtatland
CC BY-SA 4.0
2021-02-17T19:35:43.563Z
2021-02-17T19:35:43.563Z
Chris_Rands
<p>&quot;describing the entire issue&quot;.</p>
<p>I have a database (silva 16S formatted by DADA2 - fasta format) that contains ~300k entries for bacteria (16S rRNA gene). e.g.,</p>
<blockquote>
<p>AB001441.1.1538 Pseudomonas syringae pv. broussonetiae
AACTGAAGAGTTTGATCATGGCTCAGATTGAACGCTGGCGGCAGG...
KR054096.1.1224 Marinococcus halophilus
GAAGGTAAAGGCTTCCCCAAGGGCGACGATGGGTAGCCGACTTGA...
AB681405.1.1478 Paenibacillus barengoltzii
GACGAACGCTGGCGGCGTGCCTAATACATGCAAGTCGAGCGGAGT...
...</p>
</blockquote>
<p>I would like to retrieve sequences for some genus, e.g., Pseudomonas. After running the seqkit command (mentioned by you above) I ended up with ~9 of Pseudomonas entries in a fasta file (example below).</p>
<blockquote>
<p>AB001441.1.1538 Pseudomonas syringae pv. broussonetiae
AACTGAAGAGTTTGATCATGGCTCAGATTGAACGCTGGCGGCAGGCCTAACA...
AB001442.1.1538 Pseudomonas amygdali pv. eriobotryae
AACTGAAGAGTTTGATCATGGCTCAGATTGAACGCTGGCGGCAGGCCTAACA...
...</p>
</blockquote>
<p>This pool of 9k sequences of Pseudomonas has several species of Pseudomonas (e.g., P. putida, P. syringe Et Cetera). Since some species of Pseudomonas have multiple sequences in the database (e.g., P. putida has 1.8k entries) I would like to know how many different species of Pseudomonas are there (in the 9k Pseudomonas entries), by names (i.e., in the filtered Psedumonas fasta file). At the end of the day, I would like to say that in this database I have 9k entries of Pseudomonas which represents 'X' number of Pseudomonas species.</p>
<p>What I'm thinking... Retrieve the full ID (AB001442.1.1538 Pseudomonas amygdali pv. eriobotryaeall) or if possible the portion of the ID (headers) that contains the species names (e.g. Psudomonas amygdali pv. eriobotryaeall) and with this, I can list the species names and count them. Does this explanation help to understand what I'm trying to do?</p>
<p>Let me know if you need more info.</p>
<p>Thank you for your time.
Fabricio</p>

CC BY-SA 4.0
2021-02-22T16:34:06.940Z
2021-02-22T16:34:06.940Z
Fabricio
<p>First of all, sorry if this question has been posted previously (I could not found a solution accessing the previous Q&amp;A).</p>
<p>I have a fasta file as follow:</p>
<pre><code>&gt;LNIV02000036.519060.520603 Pseudomonas aeruginosa
GAACTGAAGAGTTTGATCATGGCTCAGATTGAACGCTGGCGGCAGGCCTAA...

&gt;LDLH01000045.361.1876 Pseudomonas aeruginosa
AGAGTTTGATCATGGCTCAGATTGAACGCTGGCGGCAGGCCTAACACATG...

&gt;LNMQ01000087.5825.7394 Enterococcus faecium
TTTTTATGAGAGTTTGATCCTGGCTCAGGACGAACGCTGGCGGCGTGC...

...
</code></pre>
<p>I would like to filter sequences or extract entries for a specific species, let's say Pseudomonas (I have others in mind...). I tried a few things but I think the main problem is my skills associated with bioinformatics. I would like to create a file content Pseudomonas entries (~9k entries in the current database). Any suggestion or tip?</p>
<p>Thank you for your time.</p>
<p>Best,<br />
Fabricio</p>

<fasta><database>
Filtering Sequences (entries) by headers ID from a Fasta file database
# Filtering Sequences (entries) by headers ID from a Fasta file database

(Asked by: username_0 on 2021-02-17)

First of all, sorry if this question has been posted previously (I could not found a solution accessing the previous Q&A).
I have a fasta file as follow:
>LNIV02000036.519060.520603 Pseudomonas aeruginosa
GAACTGAAGAGTTTGATCATGGCTCAGATTGAACGCTGGCGGCAGGCCTAA...

>LDLH01000045.361.1876 Pseudomonas aeruginosa
AGAGTTTGATCATGGCTCAGATTGAACGCTGGCGGCAGGCCTAACACATG...

>LNMQ01000087.5825.7394 Enterococcus faecium
TTTTTATGAGAGTTTGATCCTGGCTCAGGACGAACGCTGGCGGCGTGC...

...

I would like to filter sequences or extract entries for a specific species, let's say Pseudomonas (I have others in mind...). I tried a few things but I think the main problem is my skills associated with bioinformatics. I would like to create a file content Pseudomonas entries (~9k entries in the current database). Any suggestion or tip?
Thank you for your time.
Best,
Fabricio

**Tags:** fasta, database

---

## 3 answers

### Answer by username_1 on 2021-02-18. Score: 4

You can do this using seqkit as follows:
seqkit grep -r -n -p '.*Pseudomonas.*' temp.fa

To explain a little, seqkit grep will allow you to search FASTA/Q files by sequence name or sequence itself. In this instance:

-r tells that the pattern is a regular expression
-n to match by full name instead of just id
-p to specify the regular expression pattern to search

There are a bunch of other tools that come with seqkit that can prove indispensable for basic sequence analyses.

> **Comments:**
> **username_0**: Hi @vkkodali, thank you for your time posting it. I did not know about seqkit (living and learning). I'm testing it out and looks like a great tool. Thanks again.  
> **username_0**: Hi again @vkkodali, another question (sorry, I'm a novice in the world of bioinformatics). Is there a way to retrieve the whole sequence header or ID using seqkit? I filtered the sequences that belong to Pseudomonas and the fasta file contains 38K entries of Pseudomonas. I would like to extract the sequence ID to count how many species of Pseudomonas are in this file (38K entries).  illustration:  >AB001446.1.1538 Pseudomonas amygdali pv. mori AACTGAAGAGTT...  >AB001440.1.1538 Pseudomonas coronafaciens pv. atropurpurea AACTGAAGAGTTT... > ...  Thank you for your time and attention.  
> **username_1**: So you just want the `AB001446.1.1538` portion of the header? So many ways to do that... Easiest would be something like `grep '^>' in.fa | cut -f1 -d ' ' | sed 's/>//'`. Another way is to do `seqkit faidx in.fa` and that creates a `in.fa.fai` which is essentially a tab-delimited file with the first column containing the identifiers.  
> **username_0**: Hi, can be the whole header or portion of it. The portion is the genus and species name (e.g., Pseudomonas putid). I would like to extract the species name from the filtered Psedumonas file (e.g., Pseudomonas amygdali pv. mori, Pseudomonas coronafaciens pv. atropurpurea etc.). I would like to know how many different species of Pseudomonas are there. Thanks again for your time. Fabricio  
> **username_1**: This can be done using a combination of `grep` and `sed` but I'd be making assumptions about how the organism names in the headers for all your 9k sequences are formatted. Perhaps you should describe the entire issue with an example of what you are starting with and what you would like to see in the end and I can update my answer?  
> **username_0**: Hi @vkkodali, due to the limitation of characters in the comment section here. I wrote an answer which actually is a follow-up (explanation) for the same issue. Thank you for your time. Fabricio.  

---

### Answer by username_2 on 2021-02-17. Score: 2

Using Biopython:
from Bio import SeqIO

for record in SeqIO.parse("input.fa", "fasta"):
    if record.description.split()[-2:] == ["Pseudomonas", "aeruginosa"]:
        print(record.format("fasta"))

> **Comments:**
> **username_0**: Hi Chris, right after the post I realize that the issue was not properly described. So, this database has ~9k entries of Pseudomonas species, i.e., hundreds of P. species. The command above asks for the name of each Pseudomonas entry, correct ("Pseudomonas", "X")? Is there a way to retrieve the Pseudomonas entries without specifying all Pseudomonas names (I do not have a list name)? Thanks.  
> **username_2**: @Fabricio Hi Fabricio, do you understand the above code? If so you can easily modify the `if` line to search for entries with "Pseudomonas"  
> **username_3**: @Fabricio Add this at the top of the code, after the other `import` statement: `import re` . Change the `if` expression to this: `if re.findall(r'pseudomonas', record.description, re.I):` , in order to search in the fasta header for the word `pseudomonas`, in case-insensitive manner.  

---

### Answer by username_0 on 2021-02-22. Score: 0

"describing the entire issue".

I have a database (silva 16S formatted by DADA2 - fasta format) that contains ~300k entries for bacteria (16S rRNA gene). e.g.,

> AB001441.1.1538 Pseudomonas syringae pv. broussonetiae
> AACTGAAGAGTTTGATCATGGCTCAGATTGAACGCTGGCGGCAGG...
> KR054096.1.1224 Marinococcus halophilus
> GAAGGTAAAGGCTTCCCCAAGGGCGACGATGGGTAGCCGACTTGA...
> AB681405.1.1478 Paenibacillus barengoltzii
> GACGAACGCTGGCGGCGTGCCTAATACATGCAAGTCGAGCGGAGT...
> ...

I would like to retrieve sequences for some genus, e.g., Pseudomonas. After running the seqkit command (mentioned by you above) I ended up with ~9 of Pseudomonas entries in a fasta file (example below).

> AB001441.1.1538 Pseudomonas syringae pv. broussonetiae
> AACTGAAGAGTTTGATCATGGCTCAGATTGAACGCTGGCGGCAGGCCTAACA...
> AB001442.1.1538 Pseudomonas amygdali pv. eriobotryae
> AACTGAAGAGTTTGATCATGGCTCAGATTGAACGCTGGCGGCAGGCCTAACA...
> ...

This pool of 9k sequences of Pseudomonas has several species of Pseudomonas (e.g., P. putida, P. syringe Et Cetera). Since some species of Pseudomonas have multiple sequences in the database (e.g., P. putida has 1.8k entries) I would like to know how many different species of Pseudomonas are there (in the 9k Pseudomonas entries), by names (i.e., in the filtered Psedumonas fasta file). At the end of the day, I would like to say that in this database I have 9k entries of Pseudomonas which represents 'X' number of Pseudomonas species.

What I'm thinking... Retrieve the full ID (AB001442.1.1538 Pseudomonas amygdali pv. eriobotryaeall) or if possible the portion of the ID (headers) that contains the species names (e.g. Psudomonas amygdali pv. eriobotryaeall) and with this, I can list the species names and count them. Does this explanation help to understand what I'm trying to do?

Let me know if you need more info.

Thank you for your time.
Fabricio


---

