shield_lock
JinDaGe
← Back to skills
extension
Category: Development & EngineeringNo API key required

bio-systems-biology-metabolic-reconstruction

Build genome-scale metabolic models from genome sequences using CarveMe and gapseq for automated reconstruction. Generate draft models ready for curation and analysis. Use when creating metabolic models for organisms without existing models.

personAuthor: jakexiaohubgithub
open_in_newRepositorydownloadDownload package

Version Compatibility

Reference examples tested with: CarveMe 1.6+, gapseq 1.2+, COBRApy 0.29+, DIAMOND 2.1+, Python 3.10+

Before using code patterns, verify installed versions match. If versions differ:

  • Python: pip show <package> then help(module.function) to check signatures
  • CLI: <tool> --version then <tool> --help to confirm flags

If code throws ImportError, AttributeError, or TypeError, introspect the installed package and adapt the example to match the actual API rather than retrying.

Note: CarveMe needs an LP solver (academic CPLEX/Gurobi; SCIP is a slow open-source fallback) and a DIAMOND install, and its universal model is BiGG-derived so the BiGG-model release matters. gapseq is cloned from GitHub (not pip-installable), emits ModelSEED-namespace models, and its reference DB version matters. Model predictions are only comparable within the same tool, DB, and namespace.

Metabolic Reconstruction

"Build a metabolic model from my organism's genome" -> Map the annotated proteome/genome to reactions in a reference database, assemble a draft network with a biomass reaction, and gap-fill so it can grow on a chosen medium.

  • CLI: carve genome.faa -o model.xml (CarveMe, top-down); gapseq doall genome.fna (gapseq, bottom-up)

The governing principle: a draft is a hypothesis, and "it grows" is guaranteed by construction

Automated reconstruction produces a DRAFT, not a finished model. The single most misleading signal is growth: CarveMe and gapseq GAP-FILL specifically to force biomass production on a chosen medium, so a draft that grows proves nothing biological - it was made to grow. Consequences:

  • Gap-filled reactions are the least-evidenced part of the model (added to close a hole, not because homology supports them), and the gap-fill medium determines what gets added. Gap-filling on the wrong medium bakes in the wrong reactions. Flag gap-filled reactions as low-confidence and record the medium.
  • Draft quality is bounded by the annotation and the reference database. CarveMe can only ever include reactions in the BiGG universe (biased toward well-studied organisms); gapseq's homology thresholds and pathway logic set its floor. Peripheral/novel metabolism is systematically underrepresented.
  • The draft is the START of curation, not the end. Different tools produce markedly different models from the same genome, and successive models of the same organism (iJR904 -> iJO1366 -> iML1515) give different predictions. Reconstruction feeds model-curation, never bypasses it.

Decision: CarveMe vs gapseq vs ModelSEED

| Goal | Tool | Why / trade-off | |------|------|-----------------| | Fast draft(s) for well-studied bacteria; batch/community | CarveMe (carve) | top-down MILP carving of a curated BiGG universe; minutes; universe is simulation-ready but BiGG-centric; universal biomass; weak transporters | | Non-model/environmental clade; carbon-source & fermentation phenotypes | gapseq | bottom-up homology + pathway-completeness; slower, more transparent; better SCFA/carbon-use recovery; ModelSEED namespace complicates merging | | Fully-automated web pipeline (RAST annotation) | ModelSEED/KBase | template-based; convenient; template biomass and aggressive gap-fill can force implausible reactions | | Eukaryotes / fungi / actinomycetes | RAVEN (MATLAB) | KEGG/MetaCyc-based, template or de novo; MATLAB license; the eukaryote-capable option |

Do NOT treat "CarveMe and gapseq do the same thing, pick the faster one" as true: different philosophies, namespaces (BiGG vs ModelSEED), and failure modes. The choice is scientific. No single tool dominates - which is why consensus/ensemble reconstruction exists.

CarveMe (top-down)

pip install carveme          # also needs DIAMOND and an LP solver (CPLEX/Gurobi; SCIP fallback)

# Draft from a PROTEIN FASTA (default input). Raw/GenBank genomes are NOT accepted.
carve genome.faa -o model.xml

# Gram type and universe are VALUES of -u/--universe, NOT --grampos/--gramneg flags.
carve genome.faa -o model.xml -u grampos    # {bacteria (default), grampos, gramneg, archaea, cyanobacteria}

# Gap-fill to force growth on a medium (opt-in; records what was added for that medium).
carve genome.faa -o model.xml --gapfill M9
carve genome.faa -o model.xml -u gramneg --gapfill M9,LB   # multiple media

# Nucleotide input instead of protein, or download by accession:
carve genome.fna --dna -o model.xml

Community reconstruction uses a SEPARATE merge_community command (not carve); see systems-biology/community-metabolic-modeling.

gapseq (bottom-up, pathway-evidence)

git clone https://github.com/jotech/gapseq && cd gapseq && ./gapseq test   # cloned, not pip; check deps

# One-shot: find + find-transport + draft + fill
./gapseq doall genome.fna

# Or the explicit steps (note find-transport is its OWN subcommand, not `find -t`):
./gapseq find -p all genome.fna          # -> genome-all-Reactions.tbl, genome-all-Pathways.tbl
./gapseq find-transport genome.fna       # -> genome-Transporter.tbl  (singular)
./gapseq draft -r genome-all-Reactions.tbl -t genome-Transporter.tbl \
               -p genome-all-Pathways.tbl -c genome.fna   # -> genome-draft.RDS, genome-rxnWeights.RDS
./gapseq fill -m genome-draft.RDS -n dat/media/M9.csv \
              -c genome-rxnWeights.RDS -g genome-rxnXgenes.RDS   # -> genome.xml / genome.RDS

Load and Sanity-Check the Draft

Goal: Read the draft, confirm it grows on the gap-fill medium, and inventory the parts most likely to be wrong.

Approach: Load the SBML into COBRApy, report network size and gene coverage, test growth, and count orphan (gene-less) reactions and exchanges - the draft's soft spots before curation.

import cobra

model = cobra.io.read_sbml_model('model.xml')
print(f'reactions={len(model.reactions)} metabolites={len(model.metabolites)} genes={len(model.genes)}')
print(f'grows on gap-fill medium: {model.slim_optimize() > 1e-3}')   # true by construction if gap-filled
orphans = [r for r in model.reactions if not r.genes]   # no GPR: gap-filled, spontaneous, or transport
print(f'orphan (gene-less) reactions: {len(orphans)}  exchanges: {len(model.exchanges)}')
# Typical bacterial draft: ~1000-2500 reactions. Far outside that range flags an annotation problem.

Namespaces (the silent killer of model comparison)

# Reaction/metabolite IDs come from the tool's reference DB: CarveMe = BiGG, gapseq/ModelSEED =
# ModelSEED (seed.*), RAVEN = KEGG/MetaCyc. Two models in different namespaces cannot be merged or
# compared directly. Reconcile through MetaNetX/MNXref (MNXM* metabolites, MNXR* reactions) BEFORE
# any cross-tool merge or community build. This BiGG-vs-ModelSEED split is exactly why community
# modeling of CarveMe + gapseq outputs breaks without reconciliation.

Common Errors

| Symptom | Cause | Fix | |---------|-------|-----| | carve errors on a genome file | GenBank/nucleotide passed where protein FASTA expected | supply a protein FASTA, or add --dna for nucleotide | | --grampos/--gramneg not recognized | those are -u/--universe VALUES, not flags | carve ... -u grampos | | Draft cannot grow at all | no gap-filling requested, or wrong medium | add --gapfill <medium>; confirm the medium supplies biomass precursors | | Draft grows on everything / implausibly | gap-fill forced reactions for the chosen medium | flag gap-filled reactions low-confidence; re-gap-fill on the correct medium; curate | | Two models will not merge / IDs mismatch | different namespaces (BiGG vs ModelSEED) | reconcile via MetaNetX/MNXref before merging | | gapseq find -t fails | transport is the find-transport subcommand | use ./gapseq find-transport genome.fna | | Very few genes / tiny network | poor annotation or wrong input file | check the proteome/annotation; verify gene IDs |

Related Skills

  • systems-biology/model-curation - Curate, gap-fill deliberately, and validate the draft (the required next step)
  • systems-biology/flux-balance-analysis - Predict growth/flux once the model is trustworthy
  • systems-biology/community-metabolic-modeling - Combine reconstructions into a community model
  • genome-annotation/prokaryotic-annotation - Produce the annotated protein FASTA CarveMe/gapseq consume
  • database-access/ncbi-datasets-cli - Fetch genome/proteome inputs

References

  • Machado D, Andrejev S, Tramontano M, Patil KR. 2018. Fast automated reconstruction of genome-scale metabolic models for microbial species and communities. Nucleic Acids Res 46(15):7542-7553. (CarveMe)
  • Zimmermann J, Kaleta C, Waschina S. 2021. gapseq: informed prediction of bacterial metabolic pathways and reconstruction of accurate metabolic models. Genome Biol 22(1):81.
  • Henry CS, DeJongh M, Best AA, et al. 2010. High-throughput generation, optimization and analysis of genome-scale metabolic models. Nat Biotechnol 28(9):977-982. (ModelSEED)
  • Wang H, Marcisauskas S, Sanchez BJ, et al. 2018. RAVEN 2.0: a versatile toolbox for metabolic network reconstruction. PLoS Comput Biol 14(10):e1006541.
  • Thiele I, Palsson BO. 2010. A protocol for generating a high-quality genome-scale metabolic reconstruction. Nat Protoc 5(1):93-121.
  • Mendoza SN, Olivier BG, Molenaar D, Teusink B. 2019. A systematic assessment of current genome-scale metabolic reconstruction tools. Genome Biol 20(1):158.
  • Moretti S, Tran VDT, Mehl F, et al. 2021. MetaNetX/MNXref: unified namespace for metabolites and biochemical reactions. Nucleic Acids Res 49(D1):D570-D574.
  • Feist AM, Palsson BO. 2010. The biomass objective function. Curr Opin Microbiol 13(3):344-349.
  • Monk JM, Lloyd CJ, Brunk E, et al. 2017. iML1515, a knowledgebase that computes Escherichia coli traits. Nat Biotechnol 35(10):904-908.