OD1 taxonomic calls depending upon number of sequences in FASTA

[smith_cc (Chad Smith)]

The version is qiime2-2018.6.

1. When I run sequence A01 by itself, feature_classifier returns a g__Herbaspirillum classification.
2. When I run the exact same sequence with sequence A02, feature_classifier returns p__OD classification.
3. feature_classifier on A02 always returns p_OD, even when A02 is run by itself.

the commands are:
qiime tools import --type FeatureData[Sequence] --input-path A01A02.fasta --output-path A01A02.qza
qiime feature-classifier classify-sklearn --i-classifier /databases/qiime2/gg-13-8-99-nb-classifier.qza --i-reads A01A02.qza --o-classification gg.qza
qiime tools export gg.qza --output-dir gg

When I run Qiime1 with RDP classifier with Greengenes 13_8 99 db, I get g__Herbaspirillum for both A01 and A02.

The FASTA:

plate_999-A01 plate_184-A01 Contig - A01
TCGACGGCAGCATGGGAGCTTGCTCCTGATGGCGAGTGGCGAACGGGTGAGTAATATATC
GGAACGTGCCCTAGAGTGGGGGATAACTAGTCGAAAGACTAGCTAATACCGCATACGATC
TACGGATGAAAGTGGGGGATCTCAAGACCTCATGCTCCTGGAGCGGCCGATATCTGATTA
GCTAGTTGGTGGGGTAAAAGCCTACCAAGGCAACGATCAGTAGCTGGTCTGAGAGGACGA
CCAGCCACACTGGGACTGAGACACGGCCCAGACTCCTACGGGAGGCAGCAGTGGGGAATT
TTGGACAATGGGGGCAACCCTGATCCAGCAATGCCGCGTGAGTGAAGAAGGCCTTCGGGT
TGTAAAGCTCTTTTGTCAGGGAAGAAACGGTTGTCTCTAATAATATTACTAATGACGGTA
CCTGAAGAATAAGCACCGGCTAACTACGTGCCAGCAGCCGCGGTAATACGTAGGGTGCAA
GCGTTAATCGGAATTACTGGGCGTAAAGCGTGCGCAGGCGGTTGTGTAAGTCAGATGTGA
AATCCCCGGGCTCAACCTGGGAATTGCATTTGAGACTGCACGGCTAGAGTGTGTCAGAGG
GGGGGTAGAATTCCACGTGTAGCAGTGAAATGCGTAGATATGTGGAGGAATACCGATGGC
GAAGGCAGCCCCCTGGGATAACACTGACGCTCATGCACGAAAGCGTGGGGAGCAAACAGG
ATTAGATACCCTGGTAGTCCACGCCCTAAACGATGTCTACTAGTTGTCGGGTCTTAATTG
CCTTGGTAACGCAGCTAACGCGTGAAGTAGACCGCCTGGGGAGTACGGTCGCAAGATTAA
AACTCAAAGGAATTGACGGGGACCCGCACAAGCGGTGGATGATGTGGATTAATTCGATGC
AACGCGAAAAACCTTACCTACCCTTGACATGGATGGAATCCCGAAGAGATTTGGGAGTGC
TCGAAAGAGAACCATCACACAGGTGCTGCATGGCTGTCGTCAGCTCGTGTCGTGAGATGT
TGGGTTAAGTCCCGCAACGAGCGCAACCCTTGTCATTAGTTGCTACGAAAGGGCACTCTA
ATGAGACTGCCGGTGACAAACCGGAGGAAGGTGGGGATGACGTCAAGTCCTCATGGCCCT
TATGGGTAGGGCTTCACACGTCATACAATGGTACATACAGAGGGCCGCCAACCCGCGAGG
GGGAGCTAATCCCAGAAAGTGTATCGTAGTCCGGATTGCAGTCTGCAACTCGACTGCATG
AAGTTGGAATCGCTAGTAATCGCGGATCAGCATGTCGCGGTGAATACGTTCCCGGGTCTT
GTACACACCGCCCGTCACACCATGGGAGCGGGTTTACCAGAAGTG
plate_999-A02 plate_184-A02 Contig - A02
CTGATGGCGAGTGGCGAACGGGTGAGTAATATATCGGAACGTGCCCTAGTAGTGGGGGAT
AACTAGTCGAAAGACTAGCTAATACCGCATACGATCTACGGATGAAAGCGGGGGATCTCA
CGACCTCATGCTATTGGAGCGGCCGATATCTGATTAGCTAGTTGGTGGGGTAAAAGCCTA
CCAAGGCGACGATCAGTAGCTGGTCTGAGAGGACGACCAGCCACACTGGGACTGAGACAC
GGCCCAGACTCCTACGGGAGGCAGCAGTGGGGAATTTTGGACAATGGGGGCAACCCTGAT
CCAGCAATGCCGCGTGAGTGAAGAAGGCCTTCGGGTTGTAAAGCTCTTTTGTCAGGGAAG
AAACGGTAGTATCTAATACATATTGGTAATGACGGTACCTGAAGAATAAGCACCGGCTAA
CTACGTGCCAGCAGCCGCGGTAATACGTAGGGTGCAAGCGTTAATCGGAATTACTGGGCG
TAAAGCGTGCGCAGGCGGTTTTGTAAGTCAGATGTGAAATCCCCGGGCTCAACCTGGGAA
CTGCATTTGAGACTGCCCGGCTAGAGTGTGTCAGAGGGGGGTAGAATTCCACGTGTAGCA
GTGAAATGCGTAGATATGTGGAGGAATACCGATGGCGAAGGCAGCCCCCTGGGACAACAC
TGACGCTCATGCACGAAAGCGTGGGGAGCAAACAGGATTAGATACCCTGGTAGTCCACGC
CCTAAACGATGTCAACTAGTTGTCGGGTCTTAATTGACTTGGTAACGCAGCTAACGCGTG
AAGTAGACCGCCTGGGGAGTACGGTCGCAAGATTAGAACTCAAAGGAATTGACGGGGACC
CGCACAAGCGGTGGATGATGTGGATTAATTCTTGCTGGCGAAAAACCTTACCTACCCTTG
ACATGGTGGAATCCCGAAGAGATAGTGAGTGCTCCCTTAGAACCGCACCAGGTGCTGCAT
GGTGTCGTCAGCTCGTGTCTGAATGTTGGGTTAATCCCGCAACGAGCGCAACCCTTG

Any input appreciated! I am running sets of 96 samples of full-length 16S rRNA sequences and getting tons of p_OD1 classifications that are suspect based on the above.Thanks!

[Nicholas_Bokulich]

Hi @smith_cc,
Neat problem, thanks for reporting!

The issue here seems to be that the classifier is getting confused about the orientation that your reads are in. There is a --p-read-orientation parameter to manually set the orientation relative to the reference sequences — but if this is not used, by default classify-sklearn will attempt to determine the orientation (based on the first 100 sequences I believe) and assume that all reads are in the same orientation.

Your reads appear to be in the same orientation as the reference sequences, so when I classify these with --p-read-orientation same, I get the following classifications:

Feature ID	Taxon	Confidence
A01	k__Bacteria; p__Proteobacteria; c__Betaproteobacteria; o__Burkholderiales; f__Oxalobacteraceae; g__Herbaspirillum; s__	0.9849395619546033
A02	k__Bacteria; p__Proteobacteria; c__Betaproteobacteria; o__Burkholderiales; f__Oxalobacteraceae; g__; s__	0.7882450394455649

Regardless of whether I am only inputting a01, a02, or both.

However, otherwise I replicate your classification issue when I do not manually set the read orientation parameter, or when I set it to reverse (big clue there).

So your problem should be solved by setting the read orientation parameter to same… unless if your reads are in fact in mixed orientations and this test subset does not reflect that.

The other question is why are your reads confusing the read orientation detector? Two guesses:

1. Your sequences are full 16S, and have the word “contig” in the header line, so I am guessing these are metagenome assemblies or something along those lines. If these contain any sequence information that falls outside of the 16S (or even outside the sequences present in the reference database, which are probably based on amplicon sequences rather than full genome sequences from which 16S has been excised), it could confuse the classifier. If this is the case, you could try trimming your 16S contigs using commonly used primers (e.g., 27f/1492r) and see if that helps.
2. Again, if the full set of query sequences is in mixed orientation, that will cause problems. We do not currently have a workaround for mixed-orientation reads (since this is somewhat atypical), other than attempting to comb them into the same orientation (outside of QIIME 2). Always very happy for contributions to QIIME 2 if you have a solution you would like to share

Let us know if that solves your problem! And please do tell us more about your data if you care to share, I am very curious how you are assembling these reads.