assigning taxonomy to OTUs with a custom made database

[Sachia]

Dear QIIME community,

I am trying to assign taxonomy using a database we made ourselves (picking out organisms, and sequencing each of them).

for this custom DB I have :
1) a .fas file with all the sequences and a unique ID to each
2) a .txt file containing the unique IDs and the corresponding taxonomic information. the .txt file is separated by ";".
 
However, when I try:

parallel_assign_taxonomy_uclust.py -i ~/../repset_OTUs.fasta -r ~/../JELLY_no_align.fas -t ~/../FASTA_TAX.txt -o assign_custom -O 4 --uclust_max_accepts 10 --similarity 0.99

it never finishes, and when i check what is running with "top", I see nothing related to the parallel_assign command appearing. the output folder is being generated with a couple of temporary files, but no log file.

I am suspecting the DB format.. but I can't find more points on how the construct it other than making the .txt file ";" separated and have the same IDs in both files.

I have attached the two DB files to this message.

Thank you for your time,

Sachia

[Colin Brislawn]

Hello Sachia,

I don't know all that much database construction (a qiime dev can comment more), but I have two suggestions about the command you ran:

--similarity 0.99 is way too high. This means that only hits over 99% will be reported! Are you expecting 10 hits over 99%? Try leaving it at the default of 90%. This will still report the best hits, but is less stringent. 

--retain_temp_files could be useful to help you see how far the script got. 

--poll_directly may prevent the script from hanging. Sometimes I find this works better for me, when I can't see the qiime script running in top.

I'll be building a database like this soon and am also interested in any advice the qiime devs have.

Colin

[Sachia]

Hi Colin,

how are you? :)

thanks for the suggestions I will try them first think tomorrow.

as for the 0.99 similarity, I did it because the sequences are 18S and when we went over the results using 0.97 (both during clustering and later taxonomy) alot of things wasn't adding up with the taxonomy. Also, the sequences in the custom DB, many are more than 0.98 similar.

Based on that would you still take the OTUs (with 0.99 clustering) and assign taxonomy with 0.9?

Thanks for replying (again) to my stream of questions!

Sachia

[Colin Brislawn]

Based on that would you still take the OTUs (with 0.99 clustering) and assign taxonomy with 0.9?

Absolutely! This is because OTU picking and taxonomy assignment use this threshold differently. 

During OTU picking, each read only goes into one OTU, so you want it to be a very high match. 97% or up is good choice. During taxonomy assignment, each read is mapped to multiple taxa, so having several good hits is good. Actually, having several bad hits (say at 93%, 91%, and 91%) is still a good thing, because this accurately reflects the lack of a perfect match in your taxonomy database. 

Colin

 

[Sachia]

So during taxonomy assignment the threshold of 0.9 is only to evaluate which assignments to list in the result file?

ah, sorry for the stupid question, just want to be sure I get what you mean.

Sachia

[Colin Brislawn]

Yes. 0.9 is the minimum to be listed as one of the three hits (or 10 hits in your example). If you take a look at your log file, you will see that many of your hits are closer to 99% or 98%.

Colin

[Sachia]

Update: it still won´t work.

I ran:

parallel_assign_taxonomy_uclust.py -i ~/../repset_OTUs.fasta -r ~/../JELLY_no_align.fas -t ~/../FASTA_TAX.txt -o assign_custom -O 4 --poll_directly --retain_temp_files

and it seems like its the same. never see commands related to assign taxonomy, to job seems to hang and I just have a bunch of temp files in the output directory.

[Mauro Tutino]

Hi Sachia, 

I have been having the same problem but with parallel_assign_taxonomy_rdp.py.

When I use "assign_taxonomy.py -m rdp" I do not have such a problem.

Have you tried the not parallelized one?

In my case, since RDP perfectly works with assign_taxonomy.py, I suppose it is a problem related to installation and how it runs parallel jobs.

Unfortunately I work on a server and I cannot look directly at the problem.

Hope that helps.

Mauro.

[Sachia]

Hi Mauro,

I did try assign_taxonomy.py but with standard settings - and I get an error about specific lines in the script ending on "valueError need more than 1 value to unpack".

However, I did not try RDP, and I will give that a try. Could I ask how your ref_db and tax files look? like the layout, should I avoid specific things like "," or empty space etc.

thank you for the input, I'll keep you updated :)

Sachia

[Mauro Tutino]

Here are the first two fasta from Greengenes ref:

>1111883

GCTGGCGGCGTGCCTAACACATGTAAGTCGAACGGGACTGGGGGCAACTCCAGTTCAGTGGCAGACGGGTGCGTAACACGTGAGCAACTTGTCCGACGGCGGGGGATAGCCGGCCCAACGGCCGGGTAATACCGCGTACGCTCGTTTAGGGACATCCCTGAATGAGGAAAGCCGTAAGGCACCGACGGAGAGGCTCGCGGCCTATCAGCTAGTTGGCGGGGTAACGGCCCACCAAGGCGACGACGGGTAGCTGGTCTGAGAGGATGGCCAGCCACATTGGGACTGAGACACGGCCCAGACTCCTACGGGAGGCAGCAGTGGGGAATCTTGCGCAATGGCCGCAAGGCTGACGCAGCGACGCCGCGTGTGGGATGACGGCCTTCGGGTTGTAAACCACTGTCGGGAGGAACGAATACTCGGCTAGTCCGAGGGTGACGGTACCTCCAAAGGAAGCACCGGCTAACTCCGTGCCAGCAGCCGCGGTAATACGGAGGGTGCGAGCGTTGTCCGGAATCACTGGGCGTAAAGGGCGCGTAGGTGGCCCGTTAAGTGGCTGGTGAAATCCCGGGGCTCAACTCCGGGGCTGCCGGTCAGACTGGCGAGCTAGAGCACGGTAGGGGCAGATGGAATTCCCGGTGTAGCGGTGGAATGCGTAGATATCGGGAAGAATACCAGTGGCGAAGGCGTTCTGCTGGGCCGTTGCTGACACTGAGGCGCGACAGCGTGGGGAGCAAACAGGATTAGATACCCTGGTAGTCCACGCCGTAAACGATGGACACTAGACGTCGGGGGGAGCGACCCTCCCGGTGTCGTCGCTAACGCAGTAAGTGTCCCGCCTGGGGAGTACGGCCGCAAGGCTGAAACTCAAAGGAATTGACGGGGGCCCGCACAAGCGGTGGAGCATGTGGTTTAATTCGAAGCAACGCGAAGAACCTTACCTGGGCTTGACATGCTGGTGCAAGCCGGTGGAAACATCGGCCCCTCTTCGGAGCGCCAGCACAGGTGCTGCATGGCTGTCGTCAGCTCGTGTCGTGAGATGTTGGGTTAAGTCCCGCAACGAGCGCAACTCTCGCTCCCAGTTGCCAGCGGTTCGGCCGGGGACTCTGGGGGGACTGCCGGCGTTAAGCCGGAGGAAGGTGGGGACGACGTCAAGTCATCATGGCCCTTACGTCCAGGGCGACACACGTGCTACAATGCCTGGTACAGCGCGTCGCGAACTCGCAAGAGGGAGCCAATCGCCAAAAGCCGGGCTAAGTTCGGATTGTCGTCTGCAACTCGACGGCATGAAGCCGGAATCGCTAGTAATCGCGGATCAGCCACGCCGCGGTGAATACGTTCCCGGGCCTTGTACACACCGCCCGTCACGCCATGGAAGCCGGAGGGACCCGAAACCGGTGGGCCAACCGCAAGGGGGCAGCCGTCTAAGGT

>1111882

AGAGTTTGATCATGGCTCAGGATGAACGCTAGCGGCAGGCCTAACACATGCAAGTCGAGGGGTAGAGGCTTTCGGGCCTTGAGACCGGCGCACGGGTGCGTAACGCGTATGCAATCTGCCTTGTACTAAGGGATAGCCCAGAGAAATTTGGATTAATACCTTATAGTATATAGATGTGGCATCACATTTCTATTAAAGATTTATCGGTACAAGATGAGCATGCGTCCCATTAGCTAGTTGGTATGGTAACGGCATACCAAGGCAATGATGGGTAGGGGTCCTGAGAGGGAGATCCCCCACACTGGTACTGAGACACGGACCAGACTCCTACGGGAGGCAGCAGTGAGGAATATTGGTCAATGGGCGCAAGCCTGAACCAGCCATGCCGCGTGCAGGATGACGGTCCTATGGATTGTAAACTGCTTTTGTACGGGAAGAAACACTCCTACGTGTAGGGGCTTGACGGTACCGTAAGAATAAGGATCGGCTAACTCCGTGCCAGCAGCCGCGGTAATACGGAGGATCCAAGCGTTATCCGGAATCATTGGGTTTAAAGGGTCCGTAGGCGGTTTTATAAGTCAGTGGTGAAATCCGGCAGCTCAACTGTCGAACTGCCATTGATACTGTAGAACTTGAATTACTGTGAAGTAACTAGAATATGTAGTGTAGCGGTGAAATGCTTAGATATTACATGGAATACCAATTGCGAAGGCAGGTTACTAACAGTATATTGACGCTGATGGACGAAAGCGTGGGGAGCGAACAGGATTAGATACCCTGGTAGTCCACGCCGTAAACGATGGATACTAGCTGTTTGGCAGCAATGCTGAGTGGCTAAGCGAAAGTGTTAAGTATCCCACCTGGGGAGTACGAACGCAAGTTTGAAACTCAAAGGAATTGACGGGGGCCCGCACAAGCGGTGGAGCATGTGGTTTAATTCGATGATACGCGAGGAACCTTACCAGGGCTTAAATGTAGAGTGACAGGACTGGAAACAGTTTTTTCTTCGGACACTTTACAAGGTGCTGCATGGTTGTCGTCAGCTCGTGCCGTGAGGTGTCAGGTTAAGTCCTATAACGAGCGCAACCCCTGTTGTTAGTTGCCAGCGAGTAATGTCGGGAACTCTAACAAGACTGCCGGTGCAAACCGTGAGGAAGGTGGGGATGACGTCAAATCATCACGGCCCTTACGTCCTGGGCTACACACGTGCTACAATGGCCGGTACAGAGAGCAGCCACCTCGCGAGGGGGAGCGAATCTATAAAGCCGGTCACAGTTCGGATTGGAGTCTGCAACCCGACTCCATGAAGCTGGAATCGCTAGTAATCGGATATCAGCCATGATCCGGTGAATACGTTCCCGGGCCTTGTACACACCGCCCGTCAAGCCATGGAAGCTGGGGGTACCTGAAGTCGGTGACCGCAAGGAGCTGCCTAGGGTAAAACTGGTAACTGGGGCTAAGTCGTACAAGGTAGCCGTA

And here are the corresponding lines in the taxa file for those OTUs. There is a TAB between the OTU_ID and the taxonomy. Between taxa levels there are only spaces, no TABs.

1111883 k__Bacteria; p__Gemmatimonadetes; c__Gemm-1; o__; f__; g__; s__

1111882 k__Bacteria; p__Bacteroidetes; c__Flavobacteriia; o__Flavobacteriales; f__Flavobacteriaceae; g__Flavobacterium; s__

I gave a quick look at your files and you only used semicolons in the taxa file instead of TABs and there is the Windows character ^M in your fasta file. I am not sure that is a problem for uclust but it is always better using a file without those characters.

Let me know if that makes any difference.

Mauro

[Sachia]

Hi guys,

sorry for the delay, I'm back with some results.

making sure the ID string and taxonomy string was tab separated and making sure there was a space between each taxonomic level, in the tax string seems to have worked (this is all in the taxonomy.txt file). Then i ran:

assign_taxonomy.py -i repset_OTUs.fasta -o assign_RDP -m rdp -t JELLY_TAX.txt -r JELLY_no_align.fas

I have attached the log file and it seems to have processed the data the correct way?

However, I am not sure I am using the right settings, so I a few more questions if you don't mind.

1) The length of the sequences in the reference DB, should i trim them to fit the specific region of the amplicons?

2) Does it matter with the reference DB, should I preferably use the aligned sequences? or is this solely for later using UniFrac?

Good luck with doing this Colin I hope it will run smoothly and thanks for the suggestions Mauro :)


Cheers,

Sachia

[Mauro Tutino]

Hi Sachia,

glad to hear it worked.

I usually use the default settings with RDP.

1) I do not think trimming the ref makes any difference for the assignment taxonomy step. It can be done when aligning the sequences. I did not try it yet but it is something I want to do.

2) RDP does not accept aligned sequences.

Just one question, is your log file called repset99_OTUs_tax_assignment.log because you clustered the OTUs at 99% of similarity?

Mauro

[Sachia]

Hi Mauro,

sorry for dropping out for a while.
Yes we clustered at 99% because its 18S and less than 250 bp length.
I know there is a discrepancy between the commands I write in here and then in the log, I just forgot to write the real name here.

Again, thank you so much for helping out.

Cheers,
Sachia