Dada2: silva 138

[bernt-matthias]

Silva 138 was just released. Do you plan to put it centrally on https://benjjneb.github.io/dada2/training.html? Or should users generate the dada versions on their own?

[benjjneb]

We will, but it may not happen that quickly.

This is in part because we have been relying on initial processing done by the mothur team of the Silva releases into database files more refined for 16S sequence assignment. Our processing after that point is in formalized in an automated R script (dada2:::makeTaxonomyFasta_Silva) but that script starts with the mothur-formatted Silva files rather than the release files.

TBH we should probably think about pulling that whole workflow into the R script, but it's not there right now.

[wbazant]

I've written such a script:
https://github.com/wbazant/silva-for-dada2/blob/master/get_and_format_silva_reference_for_dada2.sh

which tries to replace mothur's preprocessing + dada2:::makeTaxonomyFasta_Silva.

Differently to mothur's file as described by http://blog.mothur.org/2018/01/10/SILVA-v132-reference-files/, it doesn't "screen the sequences for those that span the 27f and 1492r primer region, have 5 or fewer ambiguous base calls, and that are unique" and the starting point is SILVA's exports, not their ARB files.

It reproduces ad hoc things that were done elsewhere, and maybe adds one or two more. In general it skips entries and adjusts headers to make the two target fastas from SILVA's reference fasta.

There are some sample results in the same repository as the script I've linked, I've looked at them a bit and convinced myself the results are vaguely okay. It could be used as check whether this reference is as acceptable to use as one prepared by mothur.

@benjjneb if you think it's a good enough replacement, we could archive the script and its outputs on zenodo and link to it, or add it to this repository. If it makes it more palatable, I could rewrite it into a single Perl script instead of the bash+Perl format it is currently in.

Also if I don't understand about something that needs to be done to prepare the reference files, and that my script doesn't do, I really want to find out :)

[benjjneb]

Silva v138 dada2 training fastas (generated using our pre-existing workflow) are now available from Zenodo and from the training fastas page.

[bernt-matthias]

Thanks for the notice. I'm already integrating this into Galaxy https://github.com/galaxyproject/tools-iuc/pull/2941

[pdcountway]

Silva v138 dada2 training fastas (generated using our pre-existing workflow) are now available from Zenodo and from the training fastas page.

@benjjneb - I normally use the PR2 reference database for Euks, but wanted to try processing the same data with Silva v138. It appears that the assignTaxonomy reference file: silva_nr_v138_train_set.fa on Zenodo includes 14,887 Euks, but the corresponding addSpecies reference file: silva_species_assignment_v138.fa doesn't contain any. Is that intentional?

[wbazant]

@pdcountway yes, have a look at how the file is prepared,
https://github.com/benjjneb/dada2/blob/c1adc9403fc5b8a988fde60f84283566f476747c/R/taxonomy.R#L535

[benjjneb]

Is that intentional?

The assignSpecies exact matching method was really designed for use with short-read 16S sequencing data in bacteria and archaea. I don't have nearly as much experience with Euks and so I am hesitant to recommend exact matching as the way to go there -- I mean it probably won't be wrong, but it could be too conservative in Euks, I just don't know.

So, yes, it was intentionally limited to the domain in which I felt I knew what I was doing.