Automatically Getting The Ncbi Taxonomy Id From The Genbank Identifier

Question

[Alf]

The question is whether, given a (long) list of Genbank identifiers, is possible to get the ncbi taxonomy identifier for each one. I know it may seem very easy, but I have not found any web service which makes this, and I wouldn't like to do this manually.

Answer 1

[Pierre Lindenbaum]

NCBI efetch can use an accession number instead of a gi. and the XML/Fasta returned by efetch contains the taxonomy-ID:

$ curl "https://eutils.ncbi.nlm.nih.gov/entrez/eutils/efetch.fcgi?db=nuccore&id=HQ844023&rettype=fasta&retmode=xml" 

 http://www.ncbi.nlm.nih.gov/dtd/NCBI_TSeq.dtd">
 <TSeqSet>
<TSeq>
  <TSeq_seqtype value="nucleotide"/>
  <TSeq_gi>341832806</TSeq_gi>
  <TSeq_accver>HQ844023.1</TSeq_accver>
  <TSeq_taxid>1056490</TSeq_taxid>
  <TSeq_orgname>Rotavirus A HC91xUK reassortant (UKg9KC-1)</TSeq_orgname>
  <TSeq_defline>Rotavirus A HC91xUK reassortant (UKg9KC-1) NSP3 protein gene, complete cds</TSeq_defline>
  <TSeq_length>942</TSeq_length>
  <TSeq_sequence>(...)</TSeq_sequence>
</TSeq>

</TSeqSet>

So you can use efetch and egrep to get the taxonomic-id for each of your ACNs.

for ACC in A00002 X53307 BB145968 CAA42669 V00181  AH002406  HQ844023
do
   echo -n -e "$ACC\t"
   curl -s "https://eutils.ncbi.nlm.nih.gov/entrez/eutils/efetch.fcgi?db=nuccore&id=${ACC}&rettype=fasta&retmode=xml" |\
   grep TSeq_taxid |\
   cut -d '>' -f 2 |\
   cut -d '<' -f 1 |\
   tr -d "\n"
   echo
 done


A00002    9913
X53307    1423
BB145968    10090
CAA42669    9913
V00181    7154
AH002406    538120
HQ844023    1056490

[Alf]

Man, the script kick asses. I only was looking for some orientation, but the script is _exactly_ what I want. Thanks a lot!

[Steve Moss]

I modified your script slightly to retrieve these data for some accession IDs that I failed to recover when running blastdbcmd -entry_batch against my local blast database install, even though I could get information on them when using Entrez (my local data was the latest release from the FTP). I ran the following, where failed_accessions.txt was a list of 662 accession IDs, one per line:

cat failed_accessions.txt | xargs ./get_failed_acc.sh > failed_acctax.txt

The get_failed_acc.sh script was a modification of yours as follows:

for ACC in $@
do
    echo -n "$ACC\t"
    curl -s "http://eutils.ncbi.nlm.nih.gov/entrez/eutils/efetch.fcgi?db=nuccore&id=${ACC}&rettype=fasta&retmode=xml" |\
    grep TSeq_taxid |\
    cut -d '>' -f 2 |\
    cut -d '<' -f 1 |\
    tr -d "\n"
echo
sleep 0.25
done

I added the sleep in so I didn't hammer the server too much!

[toddknutson]

NCBI is now using https:// instead of http://

Make sure you include the "s" in your link!

[wintergb]

Hi thanks for the script! I am new to linux and I got a question. I wanna use my own accession numbers to replace "A00002 X53307 BB145968 CAA42669 V00181 AH002406 HQ844023". But my accession numbers are in a .CSV file. There are several hundreds of them. How can I copy them to the for loop? Or is it possible to read the .csv file directly? Thanks a lot!

[steven.vanmaldeghem]

Just do tr ',' ' ' < filename.csv | xargs ./get_failed_acc.sh. The tr command will replace all comma's by spaces an convert them to arguments.

[Naren]

This may work only if the url goes https not http, as ncbi has turned https only in 2017.

[Pierre Lindenbaum]

I've fixed this, thanks.

[vishwaas1704]

How to use this script for multiple file as input, I want to extract the locus_tag by using start and stop genomic position

For example ACC=accession.txt START=start.txt STOP=stop.txt

curl -s "https://eutils.ncbi.nlm.nih.gov/entrez/eutils/efetch.fcgi?db=nuccore&id=${ACC}&seq_start=$(START)&$(STOP)&rettype=gb" |\ grep locus_tag |\

I could only manage with - for ACC in $(cat accession.txt), and could not make a nested for loop to take other variables from input files.

If the loop works, i am getting duplicate/triplicates hits like, without accession number. So, I could not use the retrieved data :( curl -s "https://eutils.ncbi.nlm.nih.gov/entrez/eutils/efetch.fcgi?db=nuccore&id=AE006468.2&seq_start=8718&seq_stop=9319&rettype=gb" | grep "/locus_tag "

Output is
                 /locus_tag="STM0008"
                 /locus_tag="STM0008"
                 /locus_tag="STM0008"

What to do ?

Answer 2

[Yogesh Pandit]

The BioIDMapper package in R might be helpful

 library(BioIDMapper)
 data(glist)

 # 1=GI number, 31=NCBI Taxid
 myMap <- bio.convert(glist, 1, 31)

Quick Guide for the BioIDMapper package

[marcelavillegasp]

Hi!!

I didn't know about this package... !

Can you explain me please, how should I load my ID_list file? I have been tried with the function data("myfile") and read.table("myfile") but does not work...

Thank you so much :)

Answer 3

[Neilfws]

There's a discussion of this problem on the Bioperl mailing list.

The best solution is to use ELink from EUtils. If you have a GenBank GI number such as 341926284, you construct a query linking nucleotide and taxonomy like this:

http://eutils.ncbi.nlm.nih.gov/entrez/eutils/elink.fcgi?dbfrom=nucleotide&db=taxonomy&id=341926284

Result: taxid 9606 (Homo sapiens).

 
http://www.ncbi.nlm.nih.gov/entrez/query/DTD/eLink_101123.dtd"> 
    <eLinkResult>

        <LinkSet> 
            <DbFrom>nuccore</DbFrom> 
            <IdList> 
                <Id>341926284</Id> 
        </IdList> 
        <LinkSetDb> 
            <DbTo>taxonomy</DbTo> 
            <LinkName>nuccore_taxonomy</LinkName> 
            <Link> 
                <Id>9606</Id> 
            </Link> 
        </LinkSetDb> 
    </LinkSet> 
</eLinkResult>

A second solution is to use the file gi2taxonid which can be downloaded from the NCBI FTP site (I forget where exactly; see the Bioperl thread or search Google for it).

Of course if by GenBank identifier, you mean something other than GI, then you'll need a separate initial query to find the GI for whatever identifier you have.

[Naga]

Here is the link to download the gi2taxid table ftp://ftp.ncbi.nih.gov/pub/taxonomy/

Answer 4

[-_-]

First, you need to map accession numbers (GI is deprecated) to tax ids based on nucl_*accession2taxid.gz files from here, ftp://ftp.ncbi.nih.gov/pub/taxonomy/accession2taxid/. Then based on a tax id, you can trace its whole lineage.

The whole NCBI taxonomy database is not that big. I have written some code to convert NCBI taxdump into lineages identified by tax ids, https://github.com/zyxue/ncbitax2lin. You may find it useful.

Answer 5

[Will]

Just for completeness ...

If you have a LARGE number of IDS and don't want to be limited by the EUtils pipe there is a set of files in the ftp directory: ftp://ftp.ncbi.nih.gov/pub/taxonomy/ that map TaxIDs to various other identifiers. This would only be helpful if you had too many IDS to submit to EUtils.

[annamaria.asplund]

Where do I find accession numbers?