opencobra / cobrapy

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

Constraint upper bound overwritten when copying model #730

Open coltonlloyd opened 6 years ago

coltonlloyd commented 6 years ago

Problem description

Using the model.copy() with the gurobi interface causes the upper bound of constraints to be overwritten.

Code Sample

import cobra 

model = cobra.Model()
model.add_metabolites(cobra.Metabolite('a'))
model.metabolites.a.constraint.lb = -1
model.metabolites.a.constraint.ub = 2

print('----using glpk interface (w/ copied model)-----')
model.solver = 'glpk'
model2 = model.copy()
print('Constraint lower bound: ', model2.metabolites.a.constraint.lb)
print('Constraint upper bound: ', model2.metabolites.a.constraint.ub)
print('')

model.solver = 'gurobi'
print('')
print('----using gurobi interface (w/ original model)-----')
print('Constraint lower bound: ', model.metabolites.a.constraint.lb)
print('Constraint upper bound: ', model.metabolites.a.constraint.ub)

model3 = model.copy()
print('')
print('----using gurobi interface (w/ copied model)-----')
print('Constraint lower bound: ', model3.metabolites.a.constraint.lb)
print('Constraint upper bound: ', model3.metabolites.a.constraint.ub)

Actual Output

----using glpk interface (w/ copied model)----- Constraint lower bound: -1.0 Constraint upper bound: 2.0

----using gurobi interface (w/ original model)----- Constraint lower bound: -1 Constraint upper bound: 2

----using gurobi interface (w/ copied model)----- Constraint lower bound: -1.0 Constraint upper bound: -1.0

Expected Output

----using glpk interface (w/ copied model)----- Constraint lower bound: -1.0 Constraint upper bound: 2.0

----using gurobi interface (w/ original model)----- Constraint lower bound: -1 Constraint upper bound: 2

----using gurobi interface (w/ copied model)----- Constraint lower bound: -1.0 Constraint upper bound: 2.0

Dependency Information

Please paste the output of python -c "from cobra import show_versions; show_versions()" in the details block below.

System Information ================== OS Windows OS-release 10 Python 3.6.5 Package Versions ================ cobra 0.12.0 depinfo 1.3.0 future 0.16.0 numpy 1.14.2 optlang 1.4.2 pandas 0.22.0 pip 9.0.1 ruamel.yaml 0.14.12 setuptools 39.1.0 six 1.11.0 swiglpk 1.4.4 tabulate 0.8.2 wheel 0.31.0
akaviaLab commented 2 years ago

This is still valid (and valid for GLPK as well), and since the tests use copy a lot, it can alter results of running optimization of tests.

I've attached an example script file and output that shows how the constraints and fluxes are different when knocking out a reaction. Just change it to your directory and run it. Scratch41_Archive.zip

akaviaLab commented 2 years ago

Maybe this is a problem with reading from SBML? Are the constraint bounds supposed to be read as int or float? @matthiaskoenig Either way, the copy shouldn't change something from int to float.

cdiener commented 2 years ago

For Gurobi this is a bug in when reading back the constraints from LP format (slack variables don't get applied to bounds). So this can be fixed upstream. It should be noted that this only affects the display of the bound not the actual enforcement of the bound which is correct. In particular when optimizing model3 you one see the constraint taking on arbitrary values between -1 and 2. So not immensely serious but worth fixing for sure.

cdiener commented 2 years ago

@akaviaLab Note that flux solutions don't have to be unique. The only thing guaranteed with the optimizations you are running is (1) that the optimum should the same in original and copy, and (2) that no bound should be violated within the solver tolerance. The fluxes do not have to be equal. As far as I can tell everything is fine with GLPK except for a type switch which is a bit unlucky. We could cast them to float explicitly if that is helpful.

akaviaLab commented 2 years ago

@cdiener - thanks for explaining. Should the flux be float to begin with when reading from SBML? If not, how can we fix the type switch in GLPK?

cdiener commented 2 years ago

Definitely. We can enforce it in the base optlang interface and check the individual interfaces.

akaviaLab commented 2 years ago

Okay - can you deal with that please? I don't understand optlang that well.

cdiener commented 2 years ago

Will do.