6
$\begingroup$

I'm trying to get the FastA files for some accessions (like NC_001416.1). I did not managed to find an FTP server or direct link to these files (I want to get it from command line with wget, not from a web browser). But I found an "equivalent" file in ftp://ftp.ncbi.nlm.nih.gov/sra/refseq/NC_001416.1

This file is an SRA file that should be processed with the SRA-toolkit. In particular, fastq-dump seems to be the tools of choice.

I ran

fastq-dump --fasta NC_001416.1

and I got a NC_001416.1.fasta file. But the file contains the sequence I'm looking for in chunks of 5000 bases. The IDs in the output file are:

>NC_001416.1.1 length=5000
>NC_001416.1.2 length=5000
>NC_001416.1.3 length=5000
>NC_001416.1.4 length=5000
>NC_001416.1.5 length=5000
>NC_001416.1.6 length=5000
>NC_001416.1.7 length=5000
>NC_001416.1.8 length=5000
>NC_001416.1.9 length=5000
>NC_001416.1.10 length=3502

I need that sequence under a single ID:

>NC_001416.1 length=48502

I cannot find the proper parameters to supply to fastq-dump to achieve the desired result.

Does anyone know how could I get the FastA file I'm expecting? How to process properly the SRA file or where I could download the FastA file directly.

$\endgroup$
  • 1
    $\begingroup$ Out of curiosity, why do you write "FastA"? I've never seen it written that way before, only "FASTA" or "fasta". Why the final capitalized A? $\endgroup$ – terdon Sep 19 '18 at 14:34
  • 1
    $\begingroup$ No reason, really. It's a personal way of life XD FASTA comes from Fast-All, because it originally worked with all alphabets. There was also FAST-P for proteins and FAST-N for nucleotides. So FASTA is a shorthand for Fast-All. And I just used CamelCase to name it. FastA. Same for FastQ. First full word capitalized, second word reduced to a single letter and capitalized too. $\endgroup$ – Poshi Sep 20 '18 at 9:34
  • $\begingroup$ @Poshi as much as I follow your logic, I find really bad not to use the "standard" .fasta. using your own capitalization generates confusion. $\endgroup$ – Kamil S Jaron Sep 21 '18 at 9:19
  • 1
    $\begingroup$ In file names I use the standard .fasta, .fa or .fas. In natural language, I use what it feels better for me. I'm sorry if it confuses you. To me, this capitalization is clearer. If we have to go purists, the proper way to refer the format is FASTA, all capitals (check en.wikipedia.org/wiki/FASTA_format or zhanglab.ccmb.med.umich.edu/FASTA). But that notation hurt my eyes. It's like it's shouting to me. $\endgroup$ – Poshi Sep 21 '18 at 9:43
  • $\begingroup$ typing it as FastA is kind of fun. It is a nice reminder that the file is pronounced like "fast-eh?" instead of "fast-uh". $\endgroup$ – conchoecia Sep 21 '18 at 22:18
5
$\begingroup$

Finally, I found an alternative to the SRA translation: a link that works! For those of you interested in knowing how to download FastA files from NCBI using an accession number, try the following link:

https://www.ncbi.nlm.nih.gov/search/api/sequence/${accession}/?report=fasta

Using wget to download the accession used as example:

wget -O NC_001416.1.fasta "https://www.ncbi.nlm.nih.gov/search/api/sequence/NC_001416.1/?report=fasta"
$\endgroup$
  • $\begingroup$ This works, but is it documented anywhere, and how to get a gzip? $\endgroup$ – Ciro Santilli 新疆改造中心法轮功六四事件 Dec 10 '18 at 20:38
  • 1
    $\begingroup$ I didn't found that in the documentation, but playing with the web. If you need a gzip, my advise is to download first and gzip later (or do it at the same time with pipes, but the gzip is on your side either way). $\endgroup$ – Poshi Dec 11 '18 at 8:08
  • $\begingroup$ OK. gzipping locally is problematic because the high cost operation is the network transfer. $\endgroup$ – Ciro Santilli 新疆改造中心法轮功六四事件 Dec 11 '18 at 8:13
  • $\begingroup$ I'm sorry, I didn't found any other option. Try to find information about the report parameter. Maybe there is something like gzfasta. But I cannot point you to any documentation :-? $\endgroup$ – Poshi Dec 11 '18 at 8:15
2
$\begingroup$

Here is a Biopython solution via Entrez utilies, simply copy the code into a python script and execute it:

from Bio import Entrez, SeqIO

Entrez.email = 'me@email.com'

def get_fasta_from_ids(ids_list):
    handle = Entrez.efetch(db='nucleotide', id=','.join(ids_list), rettype='fasta')
    for record in SeqIO.parse(handle, 'fasta'):
        yield record

SeqIO.write(get_fasta_from_ids(['NC_001416.1']), 'NC_001416.1.fa', 'fasta')

You can add more ids to the input list to retrieve multiple records at once.

$\endgroup$
  • $\begingroup$ Thanks for you code snippet, but this is also a solution that needs extra software to be installed (BioPython and this snippet itself), while I'm looking for something that can be run out of the box in any system with no extra software if possible. $\endgroup$ – Poshi Sep 21 '18 at 8:39
  • $\begingroup$ Is the Entrez.email actually doing anything? Does the Entrez package need it to be set in order to run? $\endgroup$ – terdon Sep 21 '18 at 8:55
  • $\begingroup$ @terdon Yes Entrez.email is required (I do use my actual email address although I've never received an email about this). Entrez.api_key can be also be added to speed up queries. BTW, I agree that using the E-utilities (as in your solution and mine) is the most flexible approach, whether you do it via python or something else is mostly a matter of taste $\endgroup$ – Chris_Rands Sep 21 '18 at 8:59
1
$\begingroup$

You can either do this directly as you show in your answer or, for a more sophisticated and flexible approach, use NCBI's edirect tool:

esearch -db nucleotide -query 'NC_001416.1' | efetch -format fasta > NC_001416.1.fa
$\endgroup$
  • $\begingroup$ That's a quite elegant solution, but it involves installing new software, and I'm trying to generate my scripts with as less extra software as possible. Althought this is a more robust solution than having a link that can change over time. $\endgroup$ – Poshi Sep 21 '18 at 8:37
  • $\begingroup$ @Poshi fair enough, I'm glad you found a solution. You might want to accept your answer to indicate that it is the one you preferred. $\endgroup$ – terdon Sep 21 '18 at 8:56
0
$\begingroup$

You can add the --origfmt option, and see what comes out. Either you will have what you want, or perhaps they uploaded the sequences that way.

See: https://ncbi.github.io/sra-tools/fastq-dump.html

$\endgroup$
  • $\begingroup$ With that parameter I have the exact same output :-( I think I already tried everything that made sense in my head. $\endgroup$ – Poshi Sep 19 '18 at 13:11
  • $\begingroup$ Did you also try to set the --fasta to 0, so that there is no line wrapping? $\endgroup$ – Macspider Sep 19 '18 at 13:31
  • $\begingroup$ Yes, same effect, oly difference is that instead of a line of ID and several lines of sequence, I have on line of ID and on line of sequence. But still 10 sequences. $\endgroup$ – Poshi Sep 19 '18 at 13:36

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.