--- name: bio-systems-biology-metabolic-reconstruction description: Builds draft genome-scale metabolic models from an annotated genome using CarveMe (top-down carving of a BiGG universal model) or gapseq (bottom-up pathway-evidence reconstruction), then loads and sanity-checks the draft in COBRApy. Use when creating a model for an organism without one, choosing between CarveMe and gapseq, gap-filling to a target medium, understanding why a draft that grows is still only a hypothesis, handling BiGG-vs-ModelSEED namespace mismatch, or preparing a draft for curation and community modeling. tool_type: cli primary_tool: CarveMe --- ## 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 ` then `help(module.function)` to check signatures - CLI: ` --version` then ` --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) ```bash 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) ```bash 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. ```python 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) ```python # 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 `; 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.