cobrapy versus cameo (Ice Cream)

[aabversteeg]

Hello,

I have been trying to perform the cobrapy's Ice Cream example (https://cobrapy.readthedocs.io/en/latest/milp.html) in cameo:

from cameo import Model, Metabolite, Reaction

cone = Reaction("cone")
popsicle = Reaction("popsicle")
cone.variable_kind = "integer"
popsicle.variable_kind = "integer"

# constrainted to a budget
budget = Metabolite("budget")
budget._constraint_sense = "L"
budget._bound = starting_budget

cone.add_metabolites({budget: cone_production_cost})
popsicle.add_metabolites({budget: popsicle_production_cost})

# objective coefficient is the profit to be made from each unit
cone.objective_coefficient = cone_selling_price - cone_production_cost
popsicle.objective_coefficient = popsicle_selling_price - popsicle_production_cost

m = Model("lerman_ice_cream_co")
m.add_reactions((cone, popsicle))

print('Objective (%s): %s' % (m.objective.direction, m.objective.expression))
print('x:', m.optimize().x_dict)
print('y:', m.optimize().y_dict)

For cameo this results in the following:

Objective (max): 4.0*cone - 4.0*cone_reverse_29d6f + 1.0*popsicle - 1.0*popsicle_reverse_a0c73
x: OrderedDict([('cone', 0.0), ('popsicle', 0.0)])
y: OrderedDict([('cone', 0.0), ('popsicle', -0.6666666666666665)])

And for cobrapy in (I was not sure how to display the expression of the objective function, there are however both maximized):

x: {'popsicle': 1.0, 'cone': 33.0}
y: None

First of the fluxes differ, secondly the y-value correspond to metabolites and reactions for cobrapy and cameo, respectively. Clearly, I'm missing something here. I would greatly appreciate any explanation on this matter. Also I believe it would be good to state information explaining these differences on the cameo vs cobrapy page (http://cameo.bio/development.html#cameo-vs-cobrapy).

Thank you in advance

[hredestig]

Thanks for your interest in cobrapy and cameo. You are raising a topic that has become something of a pet-peeve here. That ice cream example is actually not very well corresponding to models of metabolic flux since it assigns a variable to the stoichiometry of products/reactants in a reaction, but those are of course more commonly thought of a fixed, i.e. it takes six CO2 to create one C6H12O6 and that is fixed, not something that one needs to model, hence rarely any need for the attribute ._constraint_sense of a Metabolite (the little underscore also indicates that it is private and not really meant to be fully supported)..

For what that documentation is trying to convey, the possibility to easily formulate problems and solve these using dedicated optimisers, please have a look at optlang instead and the basic example in the readme there which shows how to solve problems like the ice cream example without contrived analogies to metabolism.

The difference cameo/cobrapy is that cameo never implemented this as it directly uses optlang and thereby (intentionally) doesn't support for ._constraint_sense.

In the near future, that part of cobrapy documentation will be removed as cobrapy also makes the switch to optlang and the functionality will thereby be dropped so please just ignore that chapter, but cheers for pointing it out and yes you are right we should update the docs in cameo.

[aabversteeg]

Thank you for the clarification as well as for referring me to the optlang example. It is indeed a better example of imposing constraints on the model. I do not consider the docs lacking much. Just this one point seemed missing, but after all was only due to my misunderstanding.