CobraPy and MatLab Cobra Toolbox showing different number of "active" reactions

[Josh Loecker]

Hi all,

I am working with a model and the number of "active" genes determined by CobraPy (python) and Cobra Toolbox (matlab) are different. Say I have this code in python

model_path = "/Users/joshl/Downloads/my_model.mat"

model = cobra.io.read_matlab_model(model_path)

model.solver = "gurobi"

model.optimize()

active_reactions = 0

for reaction in model.reactions:

    if reaction.flux < -1e-6 or reaction.flux > 1e-6:

        active_reactions += 1

print(active_reactions)

I get a result of 355 "active" genes in this case. If I do the same thing in matlab (with the code below)

load('/Users/joshl/Downloads/my_model.mat');  % Provides a "model" variable

changeCobraSolver('gurobi');

optimizeCbModel(model);

sum(solution.v<-1e-6);

sum(solution.v>1e-6);

The final two lines show "570" and "573", respectively, for a total of 1,143 "active" reactions. 

Does anyone have an idea on why the number of active reactions is different in the two programs? I assumed since they are both using Gurobi as the solver, the difference would be negligible.

Python CobraPy Version: 0.26.3

MatLab Cobra Toolbox Version: 3.4

Gurobi Version: 10.0.2

[Thomas Pfau]

Hi,

How big is your model? If you got a very big model, you might well have a lot of cycles in the model, which can either be active or inactive in any solution which could explain the large discrepancy of active reactions.
Other than that, it's entirely possible that the two different toolboxes have different default settings for the solvers, which can easily lead to different solutions.

You could try things like parsimonious FBA to also minimize the reaction flux while maintaining the optimal value which should give you more similar solutions as to the number of active reactions, but even there this is not a guarantee.

Best

Thomas

--

---
You received this message because you are subscribed to the Google Groups "COBRA Toolbox" group.
To unsubscribe from this group and stop receiving emails from it, send an email to cobra-toolbo...@googlegroups.com.
To view this discussion on the web visit https://groups.google.com/d/msgid/cobra-toolbox/4c91e05c-af4e-4973-95f2-5271265fd175n%40googlegroups.com.

[Josh Loecker]

The following are some statistics about my model, I'm not sure if it is considered "large" or not

1,154 total genes

4,532 total reactions

Thank you for the advise with parsimonious FBA, I will try that to see if it helps show a comparable number of active reactions