cobrapy icon indicating copy to clipboard operation
cobrapy copied to clipboard

FastCC results differ from respective MATLAB function

Open NantiaL opened this issue 3 years ago • 20 comments

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:

  • 2186 irreversible reactions
  • 1274 blocked reactions

Python:

  • 2192 irreversible reactions
  • blocked reactions differ after each run Used solver: glpk, python 3.8.5

NantiaL avatar Feb 21 '22 13:02 NantiaL

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

Midnighter avatar Feb 21 '22 16:02 Midnighter

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.

NantiaL avatar Feb 21 '22 16:02 NantiaL

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.

cdiener avatar Feb 21 '22 17:02 cdiener

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

synchon avatar Feb 23 '22 12:02 synchon

@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

babessell1 avatar Aug 03 '22 19:08 babessell1

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.

synchon avatar Aug 04 '22 20:08 synchon