systems-biology-metabolic-reconstruction

name	systems-biology-metabolic-reconstruction
description	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.
tool_type	cli
primary_tool	CarveMe

Version Compatibility

Reference examples tested with: COBRApy 0.29+, NCBI BLAST+ 2.15+

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.

Metabolic Reconstruction

"Build a metabolic model for my organism from its genome" -> Generate a genome-scale metabolic model by mapping protein sequences to a universal reaction database, then gap-filling to ensure growth capability on specified media.

• CLI: carve (CarveMe) for automated reconstruction from protein FASTA
• CLI: gapseq find + gapseq draft + gapseq fill for pathway-based reconstruction

CarveMe (Recommended)

# Install CarveMe
pip install carveme

# Basic reconstruction from protein FASTA
carve genome.faa -o model.xml

# Specify output format
carve genome.faa -o model.xml --format sbml
carve genome.faa -o model.json --format json

# Gap-fill for specific media
carve genome.faa -o model.xml --gapfill M9

# Available media: M9, LB, M9[glc], M9[glyc], etc.

CarveMe Options

# Use diamond instead of blastp (faster)
carve genome.faa -o model.xml --diamond

# Specify organism type
carve genome.faa -o model.xml --grampos  # Gram-positive
carve genome.faa -o model.xml --gramneg  # Gram-negative (default)

# Initialize from template model
carve genome.faa -o model.xml --init M9

# Verbose output for debugging
carve genome.faa -o model.xml -v

gapseq (Alternative)

# Install gapseq
git clone https://github.com/jotech/gapseq
cd gapseq
./gapseq check  # Check dependencies

# Full reconstruction workflow
./gapseq find -p all genome.fasta  # Find metabolic pathways
./gapseq find -t all genome.fasta  # Find transporters
./gapseq draft -r genome-all-Reactions.tbl \
               -t genome-Transporters.tbl \
               -p genome-all-Pathways.tbl \
               -c genome.fasta
./gapseq fill -m genome-draft.RDS -c genome.fasta -n M9

Python API for CarveMe

import subprocess

def reconstruct_model(fasta_path, output_path, media='M9', grampos=False):
    '''Run CarveMe reconstruction

    Args:
        fasta_path: Path to protein FASTA file
        output_path: Output model file path (.xml or .json)
        media: Gap-filling media (M9, LB, etc.)
        grampos: True for Gram-positive organisms

    Model size expectations:
    - Bacteria: 1000-2500 reactions typical
    - Fungi: 1500-3000 reactions
    - Archaea: 800-1500 reactions
    '''
    cmd = ['carve', fasta_path, '-o', output_path, '--gapfill', media]

    if grampos:
        cmd.append('--grampos')

    subprocess.run(cmd, check=True)
    return output_path

Load and Inspect Draft Model

import cobra

model = cobra.io.read_sbml_model('model.xml')

print(f'Reactions: {len(model.reactions)}')
print(f'Metabolites: {len(model.metabolites)}')
print(f'Genes: {len(model.genes)}')

# Check if model can grow
solution = model.optimize()
print(f'Growth rate: {solution.objective_value:.4f}')

# List exchange reactions (available nutrients)
for rxn in model.exchanges[:10]:
    print(f'{rxn.id}: {rxn.reaction}')

Quality Metrics

def assess_model_quality(model):
    '''Basic quality assessment for draft model

    Returns metrics to evaluate reconstruction quality.
    '''
    metrics = {
        'reactions': len(model.reactions),
        'metabolites': len(model.metabolites),
        'genes': len(model.genes),
        'gene_reaction_ratio': len(model.reactions) / max(1, len(model.genes))
    }

    # Count reaction types
    metrics['exchanges'] = len(model.exchanges)
    metrics['transport'] = len([r for r in model.reactions if 'transport' in r.name.lower()])

    # Test growth
    sol = model.optimize()
    metrics['can_grow'] = sol.status == 'optimal' and sol.objective_value > 0.001

    # Gene-reaction rules
    metrics['orphan_reactions'] = len([r for r in model.reactions if not r.genes])

    return metrics

Multiple Genome Reconstruction

import os
from pathlib import Path

def batch_reconstruction(fasta_dir, output_dir, media='M9'):
    '''Reconstruct models for multiple genomes

    Use for comparative genomics or community modeling.
    '''
    os.makedirs(output_dir, exist_ok=True)

    for fasta in Path(fasta_dir).glob('*.faa'):
        output = Path(output_dir) / f'{fasta.stem}.xml'
        reconstruct_model(str(fasta), str(output), media=media)
        print(f'Completed: {fasta.name}')

Community Model Construction

def merge_models(model_paths, community_name='community'):
    '''Create community model from individual organisms

    For microbiome FBA, need to create a shared compartment
    for metabolite exchange between organisms.
    '''
    import cobra

    models = [cobra.io.read_sbml_model(p) for p in model_paths]

    # Add species prefix to all components
    for i, model in enumerate(models):
        species_id = f'sp{i+1}'
        for rxn in model.reactions:
            rxn.id = f'{species_id}_{rxn.id}'
        for met in model.metabolites:
            met.id = f'{species_id}_{met.id}'
        for gene in model.genes:
            gene.id = f'{species_id}_{gene.id}'

    # Merge into community model
    community = models[0].copy()
    for model in models[1:]:
        community.merge(model)

    return community

Related Skills

• systems-biology/model-curation - Validate and curate draft models
• systems-biology/flux-balance-analysis - Analyze reconstructed models
• database-access/entrez-fetch - Download genome sequences