Skip to content

Commit

Permalink
More documentation
Browse files Browse the repository at this point in the history
  • Loading branch information
ym2877 committed Jul 6, 2023
1 parent 695ab80 commit dff5815
Show file tree
Hide file tree
Showing 5 changed files with 74 additions and 6 deletions.
17 changes: 15 additions & 2 deletions pymgpipe/diet.py
Original file line number Diff line number Diff line change
Expand Up @@ -16,7 +16,7 @@ def get_available_diets():
"""Returns all diets that come pre-packaged with pymgpipe"""
return [f.split('.txt')[0] for f in os.listdir(resource_filename("pymgpipe", "resources/diets/"))]

def get_adapted_diet(
def _get_adapted_diet(
diet, essential_metabolites=None, micronutrients=None, vaginal=False, threshold=0.8
):
if not isinstance(diet, pd.DataFrame):
Expand Down Expand Up @@ -403,6 +403,7 @@ def get_adapted_diet(

# Removes any diet
def remove_diet(model):
"""Removes existing diet from model"""
model = load_model(model)
print("Removing diet from model...")
for d in get_reactions(model, regex="Diet_EX_.*"):
Expand All @@ -411,6 +412,7 @@ def remove_diet(model):

# Finds diet current set in model
def get_diet(model):
"""Returns existing diet from model"""
model = load_model(model)
print("Fetching diet from model...")

Expand All @@ -431,6 +433,17 @@ def add_diet_to_model(
threshold=0.8,
check=True
):
"""Add pymgpipe-adapated diet to model as defined by original mgPipe paper (see README for more details)
Args:
model (optlang.interface.model): LP problem
diet (pandas.DataFrame | str): Path to diet or dataframe
essential_metabolites (list): Custom of essential metabolites (uses pre-defined list by default)
micronutrients (list): Custom list of micronutrients (uses pre-defined list by default)
vaginal (bool): Whether or not this is a vaginal diet
threshold (float): Value between 0 and 1 that defines how strict the diet constraints are (with 1 being the least strict)
check (bool): Check whether or not this diet is feasible (can take some time depending on size of model)
"""
model = load_model(model)

print("\nAttempting to add diet...")
Expand Down Expand Up @@ -496,7 +509,7 @@ def add_diet_to_model(
if micronutrients is not None:
print("Using custom set of micronutrients...")

d = get_adapted_diet(diet_df, essential_metabolites, micronutrients, vaginal, threshold)
d = _get_adapted_diet(diet_df, essential_metabolites, micronutrients, vaginal, threshold)

logger.info("Adding %s diet to model..." % diet)
added = []
Expand Down
11 changes: 10 additions & 1 deletion pymgpipe/io.py
Original file line number Diff line number Diff line change
Expand Up @@ -72,6 +72,7 @@ def _get_optlang_interface(solver):
# Loads cobra file and returns cobrapy model
# RETURNS- cobra model
def load_cobra_model(file, solver="gurobi"):
"""Loads and returns COBRA model"""
_, ext = path.splitext(file)
read_func = _read_funcs[ext]
try:
Expand All @@ -85,8 +86,11 @@ def load_cobra_model(file, solver="gurobi"):
return model


# Warning- COBRA won't save any additional custom constraints to the model (i.e. coupling constraints)
def write_cobra_model(model, file):
"""Writes COBRA model to file
Warning: COBRA won't save any custom modifications added to the model (i.e. coupling constraints, diet, etc). For this reason, we recommend working and using the corresponding optlang LP problems.
"""
_, ext = path.splitext(file)
write_func = _write_funcs[ext]
try:
Expand All @@ -99,6 +103,10 @@ def write_cobra_model(model, file):
# Loads either LP file or cobra file and returns optlang model representing underlying optimization problem
# RETURNS- optlang model
def load_model(path, solver="gurobi"):
"""Loads optlang.interface.Model from either an LP file (.lp or .mps) or any of the available COBRA file types (.xml, .mat, etc)
Returns: optlang.interface.Model
"""
if isinstance(path, cobra.Model):
path.solver.name = path.name
return path.solver
Expand Down Expand Up @@ -136,6 +144,7 @@ def load_model(path, solver="gurobi"):


def write_lp_problem(model, out_file=None, compress=True, force=True):
"""Writes optlang.interface.Model out to file (will compress by default)"""
out_file = "./" + model.name + ".xml" if out_file is None else out_file

if compress and not (out_file.endswith(".gz") or out_file.endswith(".7z")):
Expand Down
4 changes: 2 additions & 2 deletions pymgpipe/main.py
Original file line number Diff line number Diff line change
Expand Up @@ -18,7 +18,7 @@
from .io import load_cobra_model, write_lp_problem, write_cobra_model, suppress_stdout
from .utils import load_dataframe, remove_reverse_vars
from .coupling import add_coupling_constraints
from .metrics import compute_diversity_metrics
from .metrics import _compute_diversity_metrics
from .logger import logger

cobra_config = Configuration()
Expand Down Expand Up @@ -266,7 +266,7 @@ def _inner(
write_cobra_model(pymgpipe_model, model_out)

if compute_metrics:
metrics = compute_diversity_metrics(pymgpipe_model)
metrics = _compute_diversity_metrics(pymgpipe_model)
if metrics is None or len(metrics)==0:
logger.warning('Unable to compute diversity metrics for %s'%pymgpipe_model.name)

Expand Down
2 changes: 1 addition & 1 deletion pymgpipe/metrics.py
Original file line number Diff line number Diff line change
@@ -1,7 +1,7 @@
from .utils import get_abundances
import cobra

def compute_diversity_metrics(model):
def _compute_diversity_metrics(model):
assert isinstance(model, cobra.Model), '`model` needs to be COBRA model to compute diversity metrics.'
print('Computing metrics for %s...'%model.name)

Expand Down
46 changes: 46 additions & 0 deletions pymgpipe/utils.py
Original file line number Diff line number Diff line change
Expand Up @@ -32,6 +32,19 @@ def solve_model(
flux_threshold=None,
ex_only=True,
):
"""Solves optlang.interface.Model
Args:
model (optlang.interface.model): LP problem
regex (str): Regex string corresponding to list of reactions you want to return
reactions (list): List of reactions you want to return
ex_only (bool): Only return fluxes for exchange reactions (will be overridden by either `regex` or `reactions` param)
flux_threshold (float): Any flux below this threshold value will be set to 0
Notes:
`presolve` and `method` are both solver-specific parameters. Please refer to [optlang documentation](https://optlang.readthedocs.io/en/latest/) for more info.
"""

model = load_model(model, solver)
model.configuration.verbosity = verbosity
model.configuration.presolve = presolve
Expand Down Expand Up @@ -76,6 +89,16 @@ def _get_fluxes_from_model(model, reactions=None, regex=None, threshold=1e-5):


def get_reactions(model, reactions=None, regex=None, include_reverse=False):
"""Fetches reactions from model
Args:
model (optlang.interface.model): LP problem
regex (str): Regex string corresponding to list of reactions you want to return
reactions (list): List of reactions you want to return
include_reverse (bool): Whether or not you want to return reverse variables as well (if model contains reverse variables)
Returns: List of optlang.interface.Variable
"""
model = load_model(model)
r = []
if reactions is not None and len(reactions) > 0:
Expand Down Expand Up @@ -112,11 +135,34 @@ def get_reactions(model, reactions=None, regex=None, include_reverse=False):
return r

def get_net_reactions(model, reactions=None, regex=None):
"""Fetches NET reactions from model
Args:
model (optlang.interface.model): LP problem
regex (str): Regex string corresponding to list of reactions you want to return
reactions (list): List of reactions you want to return
Notes:
By default, COBRA adds a reverse variable for every forward variable within the model.
Therefore, to calculate the flux going through `reaction_A` for example, you must substract the forward and reverse fluxes like so- `reaction_A - reaction_A_reverse`
"""
rxns = get_reactions(model,reactions,regex,include_reverse=True)
net = {rxns[i].name:rxns[i]-rxns[i+1] for i in range(0,len(rxns),2)}
return net

def constrain_reactions(model, flux_map, threshold=0.0):
"""Constrains reactions within model to specific value (or range of values)
Args:
model (optlang.interface.model): LP problem
flux_map (dictionary): Dictionary defining flux value for each reaction you want to constrain
threshold (float): Fixed value defining flexibility threshold (not a percentage!)
Example:
Passing in the following dictionary {'EX_reactionA': 400, 'EX_reactionB': 100} with a threshold of 50 will set the bounds for these reactions like so-
`350 <= EX_reactionA <= 450`
`50 <= EX_reactionB <= 150`
"""
model = load_model(model)
if isinstance(flux_map, pd.Series):
flux_map = flux_map.to_dict()
Expand Down

1 comment on commit dff5815

@github-actions
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Coverage

Coverage Report
FileStmtsMissCoverMissing
pymgpipe
   coupling.py39685%41, 64, 67–68, 73, 88
   diet.py1031982%17, 23, 397–398, 451–467, 480–491, 496–499, 508, 510, 539
   fva.py1445959%32, 36, 88–92, 98, 107–113, 116–117, 120–121, 139–149, 165, 174–235, 244
   io.py1063765%17, 55, 61–66, 81–82, 99–100, 116, 119, 125–128, 133–141, 151, 155–156, 161–166, 170–176, 185–186
   main.py1484470%112, 115, 167, 210–212, 243–297, 304, 309–310, 323
   metrics.py242112%5–35
   modeling.py145795%36, 39, 54–57, 129, 131
   nmpc.py64592%95, 132–134, 136
   utils.py23510754%58–59, 62, 65, 80, 86, 106, 110, 125–128, 134, 149–151, 166–184, 188–198, 206, 219–220, 222–223, 245–246, 251–253, 291–294, 296–310, 316, 327–332, 335–339, 359–360, 371–413
pymgpipe/tests
   test_e2e.py90199%197
TOTAL130930677% 

Tests Skipped Failures Errors Time
23 0 💤 0 ❌ 0 🔥 1m 12s ⏱️

Please sign in to comment.