Skip to content

Commit

Permalink
Merge pull request #4 from vivarium-collective/reorg
Browse files Browse the repository at this point in the history
Reorg
  • Loading branch information
eagmon authored Nov 8, 2024
2 parents 4f84855 + 676b7e5 commit e055815
Show file tree
Hide file tree
Showing 12 changed files with 853 additions and 650 deletions.
5 changes: 5 additions & 0 deletions .gitignore
Original file line number Diff line number Diff line change
@@ -0,0 +1,5 @@
.python-version
.ipynb_checkpoints/
spatio_flux.egg-info/
__pycache__/
out/
404 changes: 233 additions & 171 deletions demo/particle_comets.ipynb

Large diffs are not rendered by default.

27 changes: 19 additions & 8 deletions spatio_flux/__init__.py
Original file line number Diff line number Diff line change
Expand Up @@ -4,9 +4,7 @@
"""

from process_bigraph import ProcessTypes

# make type system
core = ProcessTypes()
from spatio_flux.processes import register_processes


def apply_non_negative(schema, current, update, core):
Expand All @@ -17,13 +15,26 @@ def apply_non_negative(schema, current, update, core):
positive_float = {
'_type': 'positive_float',
'_inherit': 'float',
'_apply': apply_non_negative
}
core.register('positive_float', positive_float)
'_apply': apply_non_negative}


bounds_type = {
'lower': 'maybe[float]',
'upper': 'maybe[float]'
'upper': 'maybe[float]'}


particle_type = {
'id': 'string',
'position': 'tuple[float,float]',
'size': 'float',
'local': 'map[float]',
'exchange': 'map[float]', # {mol_id: delta_value}
}
core.register_process('bounds', bounds_type)


def register_types(core):
core.register('positive_float', positive_float)
core.register('bounds', bounds_type)
core.register('particle', particle_type)

return register_processes(core)
13 changes: 13 additions & 0 deletions spatio_flux/processes/__init__.py
Original file line number Diff line number Diff line change
@@ -0,0 +1,13 @@
from spatio_flux.processes.dfba import DynamicFBA
from spatio_flux.processes.diffusion_advection import DiffusionAdvection
from spatio_flux.processes.particles import Particles, MinimalParticle


def register_processes(core):
core.register_process('DynamicFBA', DynamicFBA)
core.register_process('DiffusionAdvection', DiffusionAdvection)
core.register_process('Particles', Particles)
core.register_process('MinimalParticle', MinimalParticle)

return core

30 changes: 30 additions & 0 deletions spatio_flux/processes/comets.py
Original file line number Diff line number Diff line change
Expand Up @@ -20,6 +20,36 @@
'initial_min_max': {'glucose': (0, 10), 'acetate': (0, 0), 'biomass': (0, 0.1)},
}

## TODO -- maybe we need to make specific composites
class COMETS(Composite):
"""
This needs to declare what types of processes are in the composite.
"""
config_schema = {
'n_bins': 'tuple',
}

def __init__(self, config, core=None):
# set up the document here
state = {
'dFBA': {
'config': {
'n_bins': config['n_bins'],
}
},
'diffusion': {
'config': {
'something_else': config['n_bins'],
}
}
}

super().__init__(config, core=core)

# TODO -- this could be in Process.
def get_default(self):
return self.core.default(self.config_schema)


def run_comets(
total_time=10.0,
Expand Down
144 changes: 20 additions & 124 deletions spatio_flux/processes/dfba.py
Original file line number Diff line number Diff line change
Expand Up @@ -10,16 +10,12 @@
from cobra.io import load_model
from process_bigraph import Process, Composite
from bigraph_viz import plot_bigraph
from spatio_flux import core # import the core from the processes package
from spatio_flux.viz.plot import plot_time_series, plot_species_distributions_to_gif

# Suppress warnings
warnings.filterwarnings("ignore", category=UserWarning, module="cobra.util.solver")
warnings.filterwarnings("ignore", category=FutureWarning, module="cobra.medium.boundary_types")

# TODO -- can set lower and upper bounds by config instead of hardcoding
MODEL_FOR_TESTING = load_model('textbook')


class DynamicFBA(Process):
"""
Expand All @@ -28,21 +24,16 @@ class DynamicFBA(Process):
Parameters:
- model: The metabolic model for the simulation.
- kinetic_params: Kinetic parameters (Km and Vmax) for each substrate.
- biomass_reaction: The identifier for the biomass reaction in the model.
- substrate_update_reactions: A dictionary mapping substrates to their update reactions.
- biomass_identifier: The identifier for biomass in the current state.
- bounds: A dictionary of bounds for any reactions in the model.
TODO -- check units
"""

config_schema = {
'model_file': 'string',
'model': 'Any',
'kinetic_params': 'map[tuple[float,float]]',
'biomass_reaction': {
'_type': 'string',
'_default': 'Biomass_Ecoli_core'
},
'substrate_update_reactions': 'map[string]',
'biomass_identifier': 'string',
'bounds': 'map[bounds]',
Expand All @@ -51,9 +42,7 @@ class DynamicFBA(Process):
def __init__(self, config, core):
super().__init__(config, core)

if self.config['model_file'] == 'TESTING':
self.model = MODEL_FOR_TESTING
elif not 'xml' in self.config['model_file']:
if not 'xml' in self.config['model_file']:
# use the textbook model if no model file is provided
self.model = load_model(self.config['model_file'])
elif isinstance(self.config['model_file'], str):
Expand All @@ -78,9 +67,15 @@ def outputs(self):
'substrates': 'map[positive_float]'
}

# def interface(self):
# return {
# 'inputs': {'substrates': 'map[positive_float]'},
# 'outputs': {'substrates': 'map[positive_float]'},
# }

# TODO -- can we just put the inputs/outputs directly in the function?
def update(self, state, interval):
substrates_input = state['substrates']
def update(self, inputs, interval):
substrates_input = inputs['substrates']

for substrate, reaction_id in self.config['substrate_update_reactions'].items():
Km, Vmax = self.config['kinetic_params'][substrate]
Expand All @@ -93,13 +88,13 @@ def update(self, state, interval):
solution = self.model.optimize()
if solution.status == 'optimal':
current_biomass = substrates_input[self.config['biomass_identifier']]
biomass_growth_rate = solution.fluxes[self.config['biomass_reaction']]
biomass_growth_rate = solution.objective_value
substrate_update[self.config['biomass_identifier']] = biomass_growth_rate * current_biomass * interval

for substrate, reaction_id in self.config['substrate_update_reactions'].items():
flux = solution.fluxes[reaction_id] * current_biomass * interval
old_concentration = substrates_input[substrate]
new_concentration = max(old_concentration + flux, 0) # keep above 0
new_concentration = max(old_concentration + flux, 0) # keep above 0 -- TODO this should not happen
substrate_update[substrate] = new_concentration - old_concentration
# TODO -- assert not negative?
else:
Expand All @@ -113,15 +108,10 @@ def update(self, state, interval):
}


# register the process
core.register_process('DynamicFBA', DynamicFBA)


# Helper functions to get specs and states
def dfba_config(
model_file='textbook',
kinetic_params=None,
biomass_reaction='Biomass_Ecoli_core',
substrate_update_reactions=None,
biomass_identifier='biomass',
bounds=None
Expand All @@ -141,14 +131,19 @@ def dfba_config(
return {
'model_file': model_file,
'kinetic_params': kinetic_params,
'biomass_reaction': biomass_reaction,
'substrate_update_reactions': substrate_update_reactions,
'biomass_identifier': biomass_identifier,
'bounds': bounds
}


def get_single_dfba_spec(mol_ids=None, path=None, i=None, j=None):
def get_single_dfba_spec(
model_file='textbook',
mol_ids=None,
path=None,
i=None,
j=None,
):
"""
Constructs a configuration dictionary for a dynamic FBA process with optional path indices.
Expand Down Expand Up @@ -183,7 +178,7 @@ def build_path(mol_id):
return {
'_type': 'process',
'address': 'local:DynamicFBA',
'config': dfba_config(),
'config': dfba_config(model_file=model_file),
'inputs': {
'substrates': {mol_id: build_path(mol_id) for mol_id in mol_ids}
},
Expand Down Expand Up @@ -233,102 +228,3 @@ def get_spatial_dfba_state(
}


def run_dfba_single(
total_time=60,
mol_ids=None,
):
single_dfba_config = {
'dfba': get_single_dfba_spec(path=['fields']),
'fields': {
'glucose': 10,
'acetate': 0,
'biomass': 0.1
}
}

# make the simulation
sim = Composite({
'state': single_dfba_config,
'emitter': {'mode': 'all'}
}, core=core)

# save the document
sim.save(filename='single_dfba.json', outdir='out')

# simulate
print('Simulating...')
sim.update({}, total_time)

# gather results
dfba_results = sim.gather_results()

print('Plotting results...')
# plot timeseries
plot_time_series(
dfba_results,
# coordinates=[(0, 0), (1, 1), (2, 2)],
out_dir='out',
filename='dfba_single_timeseries.png',
)


def run_dfba_spatial(
total_time=60,
n_bins=(3, 3), # TODO -- why can't do (5, 10)??
mol_ids=None,
):
if mol_ids is None:
mol_ids = ['glucose', 'acetate', 'biomass']
composite_state = get_spatial_dfba_state(
n_bins=n_bins,
mol_ids=mol_ids,
)

# make the composite
print('Making the composite...')
sim = Composite({
'state': composite_state,
'emitter': {'mode': 'all'}
}, core=core)

# save the document
sim.save(filename='spatial_dfba.json', outdir='out')

# # save a viz figure of the initial state
# plot_bigraph(
# state=sim.state,
# schema=sim.composition,
# core=core,
# out_dir='out',
# filename='dfba_spatial_viz'
# )

# simulate
print('Simulating...')
sim.update({}, total_time)

# gather results
dfba_results = sim.gather_results()

print('Plotting results...')
# plot timeseries
plot_time_series(
dfba_results,
coordinates=[(0, 0), (1, 1), (2, 2)],
out_dir='out',
filename='dfba_timeseries.png',
)

# make video
plot_species_distributions_to_gif(
dfba_results,
out_dir='out',
filename='dfba_results.gif',
title='',
skip_frames=1
)


if __name__ == '__main__':
run_dfba_single()
run_dfba_spatial(n_bins=(8,8))
Loading

0 comments on commit e055815

Please sign in to comment.