search

opencobra / memote

memote – the genome-scale metabolic model test suite
https://memote.readthedocs.io/
Apache License 2.0
125 stars 26 forks source link

Inconsistency in finding blocked reactions #639

Open hvdinh16 opened 5 years ago

hvdinh16 commented 5 years ago

Hi,

I have run memote report snapshot for my model. For some reason, in Network Topology, test Universally Blocked Reactions, the number is different when I run with cobrapy manually:

cobra.flux_analysis_find_blocked_reactions(model, reaction_list=model.reactions, open_exchanges=True)

The number when I run with cobrapy is smaller than memote's number. Is there any difference between the two tests? Which number is more reliable?

Edit: I have rerun memote test the second time and it got me the number in cobrapy. I am not sure if there is any issue in my end. If there aren't, this test result is also different between memote runs.

Best, Hoang

Midnighter commented 5 years ago

Did you run each of those with the same versions? Did you use the same solvers?

Midnighter commented 5 years ago

Because internally memote actually calls find_blocked_reactions so unless solvers changed, or the cobra version changed, I'd be very surprised to see a difference.

hvdinh16 commented 5 years ago

Hi @Midnighter,

It is indeed due to the solver. I think the issue is associated with using glpk. Here are the outputs: Selection_084

Midnighter commented 5 years ago

This is currently not handled very elegantly in cobrapy and we are trying to address that with https://github.com/opencobra/cobrapy/pull/825. Basically, (depending on the solver) the feasibility and optimality thresholds are somewhere around 1E-06 or 1E-07 but the default value for find_blocked_reactions is 1E-09.

Can you do me a favour, please, and set the zero_cutoff=1E-06 and see if you then get the same results for all the solvers? If that is the case, it will be addressed in the next cobrapy release.

hvdinh16 commented 5 years ago

After setting the zero_cutoff to 1e-6, the results obtained are still the same. image

Midnighter commented 5 years ago

Okay, odd. Please take a look at my notebook blocked_reactions.ipynb.gz where I get also 677 blocked reactions as a result. I don't have the other solvers on this computer right now, so I can't test them. Big difference is also that I'm working on Python 3.6 and your syntax looks like Python 2.

hvdinh16 commented 5 years ago

I was able to reproduce my previous results and your results. The only difference is how the function is loaded. This works: from cobra.flux_analysis import find_blocked_reactions result = find_blocked_reactions()

This will cause the miscalculation: import cobra result = cobra.flux_analysis.find_blocked_reactions() blocked_reactions.ipynb.tar.gz

ChristianLieven commented 5 years ago

@hvdinh16 could you quickly let us know if this is still an issue? If it is could you open another issue on Cobrapy and reference it here, as it seems that this is where the error originates.