opencobra / cobrapy

COBRApy is a package for constraint-based modeling of metabolic networks.
http://opencobra.github.io/cobrapy/
GNU General Public License v2.0
465 stars 218 forks source link

FastCC results differ from respective MATLAB function #1154

Open NantiaL opened 2 years ago

NantiaL commented 2 years ago

Hello, I would like to raise a question regarding the FastCC implementation of Cobrapy. After multiple runs, FastCC outputs a different number of blocked reactions using the Recon1 model. It might be reasonable because optimization never delivers unique solutions. But, the respective MATLAB function from the COBRAToolbox always returns the same results after several runs. Since both implementations are based on the same algorithm, I guess both should give the same results. So, I think this might be something worths checking.

I noticed that the number of irreversible reactions detected in the FastCC implementation deviates from the number of irreversible reactions computed by the MATLAB function. Reactions that are irreversible in the backward direction are not reported as irreversible in the MATLAB command window.

For Recon1: MATLAB:

Python:

Midnighter commented 2 years ago

Thanks for the report @NantiaL. Can you please post the exact code that you used to generate these results?

NantiaL commented 2 years ago

Hey, here it is:

MATLAB:

is_active = fastcc(model, 1e-4);
inactiveRxns = setdiff(model.rxns, model.rxns(is_active));

Python:

from cobra import *
model = io.read_sbml_model('RECON1.xml')
model.solver= 'glpk'
fastcc(model)

In Python I also got different results over multiple runs after setting the zero_cutoff and/or the flux_threshold to 1e-4.

cdiener commented 2 years ago

From a cursory look it seems like there are some bugs in the FastCC coefficient setting. For instance a reaction with a forward flux of 1 and reverse flux of 1 (net flux 0) is considered active on our formulation.

synchon commented 2 years ago

Apologies for the delayed response. I'll take a look. Thanks for the report!

babessell1 commented 2 years ago

@synchon

Hello, I was wondering if you ever found any leads on what is causing this issue. I tried to diagnose myself but it got quite overwhelming for me to understand. This issue seems pretty important for research since results should be reproducible.

Thank you

synchon commented 2 years ago

Hi @babessell1 ! The problem is most likely in the problem formulation as @cdiener pointed out earlier. I have been quite busy lately and I'm not sure when I can take tackle it. Apologies for that as I had self-assigned it. I'll remove it and let someone else take a stab at it.