merenlab / anvio

An analysis and visualization platform for 'omics data
http://merenlab.org/software/anvio
GNU General Public License v3.0
413 stars 142 forks source link

[DISCUSSION] Clarifying bitscore heuristic and KO profile check in anvi-run-kegg-kofams #2259

Closed jma1991 closed 1 month ago

jma1991 commented 2 months ago

Hello,

Could you please provide more details about the bitscore heuristic used in the anvi-run-kegg-kofams module? I am currently developing a stand-alone script and need to implement this heuristic. However, I'm confused because in my kofam_scan outputs, no single gene seems to have multiple hits to the same KO identifier. Given this observation, why does the module check if all hits correspond to unique KO profiles?

Website documentation:

For every gene without a KOfam annotation, we examine all the hits with an e-value below X and a bitscore above Y percent of the threshold. If those hits are all to a unique KOfam profile, then we annotate the gene call with that KO.

Function documentation:

To find these annotations, we adopt the following heuristic: For every gene without a KOfam annotation, we examine all the hits with an e-value below X and a bitscore above Y% of the threshold. If those hits are all to a unique KO profile, then we annotate the gene call with that KO.

Thanks! James

ivagljiva commented 2 months ago

It can happen sometimes that a single gene sequence has multiple matches to the same HMM model, especially if that gene sequence has multiple similar protein domains. Here is a random example I took from a test database (as reported in the HMMER output table):

992                  -          K23537               -            8.3e-84  276.7   0.1   2.9e-44  146.1   0.0   2.1   1   1   1   2   2   2   2 -
992                  -          K23547               -            9.4e-63  207.0   0.3   7.7e-34  111.6   0.1   3.0   2   1   0   2   2   2   2 -

The first column is the gene caller ID (992), and it has two separate hits for the same KOfam (K23537). Neither of them are great hits (the bit score threshold for this KO is 578.67, compared to 276.7 and 207.0 coming from the 6th column in the table). Since this particular KO is an ATP-binding protein, I would venture a guess that the multiple hits are related to its ATP-binding domain(s).

The reason you don't see this in the kofamscan output is because kofamscan simply doesn't report multiple hits to the same KO for a given gene. And I would say that you will rarely (if ever) see multiple annotations for the same KO for the same gene call in anvi'o, either, since this situation only tends to happen for weaker matches that are often below the bit score threshold, and if they pass the annotation heuristic in this case, we only add one annotation to the unique KO so it never results in multiple hits being added to the database.

And to more directly answer your question, we look for unique hits as an additional confidence check since these are weak matches and we want to be careful to avoid introducing garbage. If the matching KOs are not all unique, then we know this is a case of 'non-specific' matching and we shouldn't annotate :)

I hope that makes sense.

jma1991 commented 1 month ago

Hi @ivagljiva

Thank you so much for the detailed response. I believe that has cleared up my confusion.

Best, James