I'd like to create some databases with Viral Refseq Genomic.fna https://ftp.ncbi.nlm.nih.gov/refseq/release/viral/. But how do I construct the samplesheet? The `viral.1.1.genomic.fna.gz`contains 18000 genomes and they all have different taxids. Do I have to create 18000 separate fasta files (one for each virus) and add a line to the samplesheet for each one. Also, it will be very very time consuming to query

nucl_gb.accession2taxid.gz

for each accession number to get the taxid. Maybe I'm not understanding something. Thank you for your help

Hi Paul, Yes generally we would expect one entry per genome because some of the tools require this. That said, I think if the tool you select accepts a single FASTA file anyway, then we maybe able to 'hack' it by having a single row and a taxid set to 1 or something. What tool do you want to make a database for?

The question seems a lot about preparing the right inputs. Do we have some guidance on that already?

Yes, that's why I'm asking. I think conceptually a single FASTA as a row in the input sample may work if the tool can accept that. Ganon for example though requires one taxid per file, it doesn't do multiple genomes in one FASTA (at least as currently implemented). Also personally, dependibg in the context I would rarely recommend dumping the whole of refseq I to a single database. You can get massively skewed hits due to unequal representarion of speciies

no need to hack the pipeline just for my use case. I'm just trying to get back into the analysis of viral metagenomes. Creating databases (DIAMOND, BLAST, Bowtie, etc...) for these analyses can be time consuming and I thought this pipeline might help automate things.

Its not a hack per se! It might be a good thing to document for other users too, if it works :)

Kraken2 is on my list to try

Ok I think kraken2 and diamond may work with a single row in the input sheet

what taxid do I put in the samplesheet? I'm confused what the purpose of the taxid is in the samplesheet

In this case just try putting

1

I'm still not clear if that is what you really want Paul? If you had an easy way to split that fasta into separate files and create the sheet of taxid to genome, wouldn't that be preferable?

Ok maybe I misunderstood... I th9ught Paul explicitly didn't want that because it takes too long, or?

While it might be not too difficult to split the fasta into separate files, managing 18,000 files is going to be a pain. Then, it will be much harder to determine the taxid for each genome. Grepping accession numbers in the ncbi taxonomy nucl file takes a ridiculous amount of time to get the taxid.

Also splitting into the files is also not trivial because some genome shave multiple chromosomes or other weird extra elements

See, here is what I do for BLAST: makeblastdb -dbtype nucl -in viral.1.1.genomic.fna -parse_seqids -out viral.1.1.genomic Simple and quick. I just thought this pipeline would be similar and just give me a set of databases with all the tools so that I do not have to determine all the comamnd lines for each tool

You can also use NCBI download or entrez get taxonomy IDs from accessions via APIs but also the second can be just as slow as grepping

OK, got it. If the individual fasta sequence headers it might be possible to add the splitting to the pipeline. Look forward to hearing how far you get 🙂

Hrm that's true... I guess we could add support for that... Even if I disagree with the concept 😅😬

And I can imagine that many fasta are just a nightmare when it comes to header formatting. Definitely a can of worms.

It absolutely it is