transcript / samsa2

SAMSA pipeline, version 2.0. An open-source metatranscriptomics pipeline for analyzing microbiome data, built around DIAMOND and customizable reference databases.
GNU General Public License v3.0
54 stars 36 forks source link

counter scripts using blast format 6 with salltitles #64

Open bartns opened 3 years ago

bartns commented 3 years ago

From issue: Header parsing using blast format 6 with salltitles

Diamond can add the fasta header in the results table using format 6 with salltitles. Function and organism counts can therefor be obtained using the diamond output instead of parsing the reference database fasta file. This greatly improves speed and storage space needed.

We run diamond like this: diamond blastx --query input.fastq.gz --salltitles --block-size 12 --db ncbi-bact-refseq_28-01-2020_proteins.dmnd --index-chunks 1 --matrix BLOSUM62 --max-target-seqs 1 --outfmt 100 --out output_refseq --threads 20 convert to tsv: diamond view --outfmt 6 qseqid sseqid pident length mismatch gapopen qstart qend sstart send evalue bitscore stitle --out output.daa.tsv --daa output_refseq.daa

Another important change to the refseq counter is that organisms with just "sp." are added to other unknown species with the same genus

Also, "MULTISPECIES" entries are not counted in the error counter anymore.

Notes. This pull request only contains changes to the "diamond counter" scripts. So this will break the pipeline scripts. We have used are our own preprocessing pipelines up to diamond so only these needed modifications for our usage.

Another thing to do is to merge the first step of the line counting and the "database" part.

We did not compare the output results with the original version!

Hope this works/helps for you! Credits can go to @jjkoehorst and me.