--seq-mode for estimating MAGs contig abundance

[Maiya19724]

Hi, Thank you for the development of Metabuli! I have all MAG contigs from a single sample and want to estimate their relative abundance with Metabuli. Can I do this, and—if so—which --seq-mode should I use for contig input? Thanks!

[jaebeom-kim]

Sure you can! please try --seq-mode 3 :)

[Maiya19724]

Thanks a lot for your help!

[yazhinia]

Why not --seq-mode 1 ?

[jaebeom-kim]

Why not --seq-mode 1 ?

--seq-mode 1 will be also fine. It would generates some differences at subspecies ranks with default settings.
Contigs are similarly long as long-read sequences, and --seq-mode 3 affects some parameters considering the read length.
So I recommended --seq-mode 3.

[Maiya19724]

Hi @jaebeom-kim,

1. Missing SAR11 / SAR86 in my marine MAG test

   • Workflow: I ran metabuli@ed609df with default parameters (--seq-mode 1) on a marine MAG reference set, using your pre-built GTDB r214 database.
   • The output report contains no high-abundance marine taxa SAR11 or SAR86.
     However, when I annotate the same MAGs with GTDB, several of them are classified as SAR11 or SAR86.
     In the Metabuli output I cannot even find o__SAR86. ( You mentioned earlier that the two modes only cause possible subspecies-level differences. )
   • Question:
     • Is this a problem with my ocean reference MAGs set?
     • Or are SAR11 / SAR86 missing from the Metabuli GTDB r214 database?
     • Or could --seq-mode 1 have filtered them out, even though you said mode 1 vs 3 should only affect subspecies?

2. Interpretation of the “percent” column
   Am I correct that the percentages in the Metabuli report represent DNA relative abundance (the proportion of total DNA), not cellular relative abundance?
   If I need cell-level abundance, should I normalise by genome length (e.g., using lengths from GTDB r214 genome)?

Thank you very much for your time and support!

[jaebeom-kim]

Thank you!

For question 1:
Nothing is wrong with your data. SAR11 / SAR86 are not included in the Metabuli GTDB r214 database.
It is because only high-quality genomes have been used to build GTDB database.
The quality criteria for r214 was CheckM completeness > 90% and contamination < 5%.

From r220 database, which is now available as "gtdb+virus+human.tar.gz" in https://metabuli.steineggerlab.workers.dev/, the completeness and contamination values of CheckM"2" are used. This change resulted 10 genomes of SAR86 to be included in r220 Metabuli database, but not SAR11.

HOWEVER, you can add genomes of your interest to a pre-built database using updateDB module. So I recommend you downloading "gtdb+virus+human.tar.gz" and add genomes to it.
Please follow this part in the README and feel free to ask any helps. + I recommend using the latest Metabuli downloaded from here

For question 2:
Yes, the percentages are DNA relative abundance. They have to be normalized to get cell-level abundance.

[Maiya19724]

Thank you for the quick and clear response!

Question 1
So the Metabuli GTDB r214 database is not a 1-to-1 mirror of the official GTDB r214 release, correct?

• Can I list / inspect exactly which genomes are included in the current Metabuli build, then use updateDB to add the missing ones?
• If I add lower-quality genomes (e.g. SAR11 with lower completeness) to fill the gaps, could that negatively affect classification accuracy for other clades?

Question 2
In many workflows, a Kraken2 report is post-processed by Bracken to obtain cell-level abundance estimates. Metabuli produces a Kraken-style report, but the database may not directly interchangeable with a standard Kraken database.

• If I build a Bracken database myself (based on the GTDB r214), the taxon set may differ from your pre-built Metabuli DB. Would that discrepancy bias the Bracken output?
• Do you recommend an alternative method to convert Metabuli’s DNA-level abundances to cell-level abundances, ideally without losing the taxonomic resolution Metabuli provides?

Thanks again for all your work—Metabuli is a fantastic tool!

[jaebeom-kim]

Hi!

For question 1:

• Correct. Metabuli GTDB r214 database is not a 1-to-1 mirror of the official GTDB r214 release.
• In the case of GTDB r214 database, you can use database-report command to see it which genomes are included.
• I think adding genomes of high contamination rate can negatively affect classification accuracy. When it comes to completeness, 10 reads mapped to a complete genome and 10 reads mapped to 50% complete genome should be counted (normalized) differently.

For question 2:

• If you build Braken DB using superset of the Metabuli DB, I guess it would be fine. Please let me know if you manage to use Braken with Metabuli results. I have been curious for this.

Thank you so much. I would like to hear more from you.

[gilmahu]

This has been my other challenge with downstream use of metabulit results. @Maiya19724 and @jaebeom-kim has any of you manage to do a bracken type of relative abundance done after metabuli? Or in other words, how do you recommend we best use metabuli data downstream since it give good classification rates but so far it stops there.

[Maiya19724]

In my understanding, Metabuli outputs Kraken-style report. To get a relative-abundance table with Bracken, you first need a Bracken-formatted database that matches the Metabuli DB. Building that Bracken DB requires an underlying Kraken DB, so the order is: rebuild a Kraken DB from the same dataset used by Metabuli, convert it to a Bracken DB, and then run Bracken on Metabuli’s output to generate the abundance table.

Alternatively, you can download pre-built Kraken / Bracken databases from https://benlangmead.github.io/aws-indexes/k2. save the time of rebuilding everything locally, but keep in mind the genome set in those public indexes might not identical to the Metabuli database you ran, so abundance estimates may shift—sometimes noticeably—when the underlying reference differs.
I haven’t had time to validate these works, so you might want to run a small test first.
@gilmahu