search

opencobra / cobrapy

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

model.summary() is inconsistent with and without 'with' context #900

Closed bdelepine closed 4 years ago

bdelepine commented 5 years ago

Problem description

Hello,

I want to use a context to safely change the C source, compute an FBA and then display the summary. I was surprised to see that the changes introduced in the context were ignored by the summary. I think this unexpected behavior because 1/ this is different than the behavior without using context, 2/ this leads to meaningless values if we additionally ask for FVA in the summary.

++

Code Sample

import cobra
import cobra.test

# With context
model = cobra.test.create_test_model("ecoli")
with model as m:
    m.reactions.EX_glc__D_e.bounds = (0, 1000)
    m.reactions.EX_ac_e.bounds = (-10, 1000)
    solution = m.optimize()
    ans1 = m.summary()
    ans2 = m.summary(solution=solution)
    ans3 = m.summary(solution=solution, fva=0.95)

# Without context
m = cobra.test.create_test_model("ecoli")
m.reactions.EX_glc__D_e.bounds = (0, 1000)
m.reactions.EX_ac_e.bounds = (-10, 1000)
solution = m.optimize()
ans11 = m.summary()
ans22 = m.summary(solution=solution)
ans33 = m.summary(solution=solution, fva=0.95)

#
print(ans1, ans11)
print(ans2, ans22)
print(ans3, ans33)

Actual Output

ans1 and ans11 are different:

IN_FLUXES IN_FLUXES OUT_FLUXES OUT_FLUXES            OBJECTIVES           OBJECTIVES
    ID       FLUX       ID        FLUX                   ID                  FLUX   
     o2_e    19.8      h2o_e      50.1     BIOMASS_Ec_iJO1366_core_53p95M   0.982   
    nh4_e    10.6      co2_e      19.7                                        nan   
 glc__D_e      10      fe3_e      9.03                                        nan   
    fe2_e    9.04                  nan                                        nan   
     pi_e   0.948                  nan                                        nan    
 IN_FLUXES IN_FLUXES OUT_FLUXES OUT_FLUXES            OBJECTIVES           OBJECTIVES
    ID       FLUX       ID        FLUX                   ID                  FLUX   
    ac_e       10      h2o_e      16.4     BIOMASS_Ec_iJO1366_core_53p95M   0.247   
    o2_e     9.33      co2_e      9.85                                        nan   
     h_e     7.73                  nan                                        nan   
   nh4_e     2.67                  nan                                        nan   
    pi_e    0.238                  nan                                        nan

ans2 and ans22 are equivalent:

IN_FLUXES IN_FLUXES OUT_FLUXES OUT_FLUXES            OBJECTIVES           OBJECTIVES
    ID       FLUX       ID        FLUX                   ID                  FLUX   
    ac_e       10      h2o_e      16.4     BIOMASS_Ec_iJO1366_core_53p95M   0.247   
    o2_e     9.33      co2_e      9.85                                        nan   
     h_e     7.73                  nan                                        nan   
   nh4_e     2.67                  nan                                        nan   
    pi_e    0.238                  nan                                        nan    
 IN_FLUXES IN_FLUXES OUT_FLUXES OUT_FLUXES            OBJECTIVES           OBJECTIVES
    ID       FLUX       ID        FLUX                   ID                  FLUX   
    ac_e       10      h2o_e      16.4     BIOMASS_Ec_iJO1366_core_53p95M   0.247   
    o2_e     9.33      co2_e      9.85                                        nan   
     h_e     7.73                  nan                                        nan   
   nh4_e     2.67                  nan                                        nan   
    pi_e    0.238                  nan                                        nan

ans3 and ans33 are different. ans3 is incoherent:

IN_FLUXES IN_FLUXES IN_FLUXES IN_FLUXES    OUT_FLUXES   OUT_FLUXES OUT_FLUXES OUT_FLUXES            OBJECTIVES           OBJECTIVES
    ID       FLUX    FLUX_MIN  FLUX_MAX        ID          FLUX     FLUX_MIN   FLUX_MAX                 ID                  FLUX   
     ac_e        10      -0      -1.82            h2o_e    16.4       40.4        549     BIOMASS_Ec_iJO1366_core_53p95M   0.247   
     o2_e      9.33    16.1        270            co2_e    9.85       15.8       21.7                                        nan   
      h_e      7.73   -14.5        994            for_e       0          0       5.92                                        nan   
    nh4_e      2.67    10.1       14.2           urea_e       0          0       2.06                                        nan   
     pi_e     0.238     0.9       1.49         glyclt_e       0          0       1.82                                        nan   
    so4_e    0.0623   0.235       1.34            gly_e       0          0       1.64                                        nan   
      k_e    0.0483   0.182      0.192          acald_e       0          0       1.28                                        nan   
    fe2_e   0.00397   -4.98      1e+03            pyr_e       0          0       1.28                                        nan   
    fe3_e  1.36e-20    4.99     -1e+03         ser__L_e       0          0       1.13                                        nan   
 glc__D_e        -0    9.51         10            h2s_e       0          0        1.1                                        nan   
[...]

 IN_FLUXES IN_FLUXES IN_FLUXES IN_FLUXES  OUT_FLUXES OUT_FLUXES OUT_FLUXES OUT_FLUXES            OBJECTIVES           OBJECTIVES
    ID       FLUX    FLUX_MIN  FLUX_MAX      ID        FLUX     FLUX_MIN   FLUX_MAX                 ID                  FLUX   
    ac_e         10     9.53       10         h2o_e    16.4       15.2        514     BIOMASS_Ec_iJO1366_core_53p95M   0.247   
    o2_e       9.33     8.91      258         co2_e    9.85       8.97       10.4                                        nan   
     h_e       7.73     6.43    1e+03         for_e       0          0       1.39                                        nan   
   nh4_e       2.67     2.54     3.62        urea_e       0          0       0.54                                        nan   
    pi_e      0.238    0.227     0.37      glyclt_e       0          0      0.475                                        nan   
   so4_e     0.0623   0.0592    0.356         gly_e       0          0      0.376                                        nan   
   fe2_e    0.00397    -1.25      994       acald_e       0          0       0.33                                        nan   
   fe3_e   1.36e-20     1.25     -994         pyr_e       0          0       0.33                                        nan   
[...]

Expected Output

I would have expected that the summary with and without context to be equivalent.

Dependency Information

System Information ================== OS Linux OS-release 3.10.0-957.27.2.el7.x86_64 Python 3.6.8 Package Versions ================ cobra 0.16.0 depinfo 1.3.0 future 0.16.0 numpy 1.15.4 optlang 1.4.3 pandas 0.23.4 pip 19.1.1 python-libsbml-experimental 5.18.0 ruamel.yaml 0.16.5 setuptools 41.0.1 six 1.11.0 swiglpk 1.4.4 wheel 0.33.4
Midnighter commented 5 years ago

Across two different solutions I would only expect the objective value and fairly closely the in-fluxes to be equal. However, as you demonstrate, this is not the case. I will try to reproduce these but it looks quite problematic.

cdiener commented 5 years ago

Another thing that pops out is that in many cases FLUX_MIN > FLUX_MAX...

gregmedlock commented 5 years ago

I recently noticed that the flux_min/flux_max returned by FVA had the same issue--the cobrapy version I'm using for this project is 0.14.2, so the problem might be quite old.

cdiener commented 5 years ago

@gregmedlock I can't reproduce this bug with FVA alone. Could you make a new bug report with a reproducible example and indicate which solver you are using? Thanks!

cdiener commented 5 years ago

The problem with here is that model.summary() does not actually return a string representation of the summary but a Summary object that does nothing initially. The summary gets calculated when ans1 is printed and not in the with block. You can easily see that when printing ans1-ans3 within the with block which will give you correct results.

In your case ans1 uses the model after reset when print is called thus giving the unconstrained model. ans2 is using the solution from within the with block which is correct, thus everything works. ans3 uses the solution from within the with block which is correct but runs FVA on the model thus the discrepancy.

I feel like this is very odd behavior. I would expect summaries to be immutable. They should be calculated on call and not on print which can create all kinds of problems. Preferably Summary objects should be a one-shot thing. They only get calculated on init and then just store the results without any chance of re-computation.

Midnighter commented 5 years ago

Ah, thank you for investigating that Christian :smiley: That is actually what I proposed to do with my additional changes to the summary PR. So I guess it's urgent that I finish up that work.

cdiener commented 4 years ago

Do we have planned fixes for this? This has been broken since several releases so might be a good idea to at least do a workaround.

Midnighter commented 4 years ago

Sorry, completely dropped off my radar. Please poke me again in mid May if I haven't submitted anything by then.

cdiener commented 4 years ago

Sure thing. Thank you!

he-hai commented 4 years ago

Sorry, completely dropped off my radar. Please poke me again in mid May if I haven't submitted anything by then.

Maybe a hint because it seems to work fine in the old version, please see below with version 0.15.4:

import cobra
import cobra.test

# With context
model = cobra.test.create_test_model("ecoli")
with model as m:
    m.reactions.EX_glc_e.bounds = (0, 1000)
    m.reactions.EX_ac_e.bounds = (-10, 1000)
    solution = m.optimize()
    ans1 = m.summary()
    ans2 = m.summary(solution=solution)
    ans3 = m.summary(solution=solution, fva=0.95)

# Without context
m = cobra.test.create_test_model("ecoli")
m.reactions.EX_glc_e.bounds = (0, 1000)
m.reactions.EX_ac_e.bounds = (-10, 1000)
solution = m.optimize()
ans11 = m.summary()
ans22 = m.summary(solution=solution)
ans33 = m.summary(solution=solution, fva=0.95)

#
print(ans1, ans11)
print(ans2, ans22)
print(ans3, ans33)

Output

ERROR:root:'M_sertrna[sec]_c' is not a valid SBML 'SId'.
IN FLUXES          OUT FLUXES    OBJECTIVES
-----------------  ------------  ----------------------
ac_e    10         h2o_e  16.4   Ec_biomass_i...  0.247
o2_e     9.33      co2_e   9.85
h_e      7.73
nh4_e    2.67
pi_e     0.238
so4_e    0.0623
k_e      0.0483
fe2_e    0.00397
mg2_e    0.00214
ca2_e    0.00129
cl_e     0.00129
cu2_e    0.000175
mn2_e    0.000171
zn2_e    8.43e-05
ni2_e    7.98e-05
mobd_e   3.19e-05
IN FLUXES          OUT FLUXES    OBJECTIVES
-----------------  ------------  ----------------------
ac_e    10         h2o_e  16.4   Ec_biomass_i...  0.247
o2_e     9.33      co2_e   9.85
h_e      7.73
nh4_e    2.67
pi_e     0.238
so4_e    0.0623
k_e      0.0483
fe2_e    0.00397
mg2_e    0.00214
ca2_e    0.00129
cl_e     0.00129
cu2_e    0.000175
mn2_e    0.000171
zn2_e    8.43e-05
ni2_e    7.98e-05
mobd_e   3.19e-05
IN FLUXES                                OUT FLUXES                             OBJECTIVES
---------------------------------------  -------------------------------------  ----------------------
id           Flux  Range                 id                 Flux  Range         Ec_biomass_i...  0.247
------  ---------  --------------------  ---------------  ------  ------------
ac_e    10         [9.53, 10]            h2o_e             16.4   [15.2, 514]
o2_e     9.33      [8.91, 258]           co2_e              9.85  [8.97, 10.4]
h_e      7.73      [6.43, 1e+03]         fe3_e              0     [-1.25, 994]
nh4_e    2.67      [2.54, 3.62]          for_e              0     [0, 1.39]
pi_e     0.238     [0.227, 0.37]         urea_e             0     [0, 0.54]
so4_e    0.0623    [0.0592, 0.356]       glyclt_e           0     [0, 0.475]
k_e      0.0483    [0.0458, 0.0483]      gly_e              0     [0, 0.376]
fe2_e    0.00397   [-1.25, 994]          acald_e            0     [0, 0.33]
mg2_e    0.00214   [0.00204, 0.00214]    pyr_e              0     [0, 0.33]
ca2_e    0.00129   [0.00122, 0.00129]    h2s_e              0     [0, 0.297]
cl_e     0.00129   [0.00122, 0.00129]    ala__L_e           0     [0, 0.276]
cu2_e    0.000175  [0.000167, 0.000175]  etoh_e             0     [0, 0.276]
mn2_e    0.000171  [0.000162, 0.000171]  glyc__R_e          0     [0, 0.276]
zn2_e    8.4e-05   [8e-05, 8.4e-05]      lac__D_e           0     [0, 0.276]
ni2_e    8e-05     [7.6e-05, 8e-05]      ser__L_e           0     [0, 0.27]
mobd_e   3.2e-05   [3e-05, 3.2e-05]      asp__L_e           0     [0, 0.247]
                                         mal__L_e           0     [0, 0.237]
                                         asn__L_e           0     [0, 0.228]
                                         ura_e              0     [0, 0.228]
                                         lac__L_e           0     [0, 0.223]
                                         succ_e             0     [0, 0.22]
                                         dha_e              0     [0, 0.216]
                                         glyald_e           0     [0, 0.216]
                                         alltn_e            0     [0, 0.191]
                                         glyc_e             0     [0, 0.191]
                                         akg_e              0     [0, 0.188]
                                         ala__D_e           0     [0, 0.185]
                                         hom__L_e           0     [0, 0.177]
                                         12ppd__R_e         0     [0, 0.177]
                                         12ppd__S_e         0     [0, 0.177]
                                         etha_e             0     [0, 0.175]
                                         4abut_e            0     [0, 0.172]
                                         thr__L_e           0     [0, 0.172]
                                         glu__L_e           0     [0, 0.17]
                                         3hpp_e             0     [0, 0.165]
                                         cit_e              0     [0, 0.163]
                                         acser_e            0     [0, 0.16]
                                         xan_e              0     [0, 0.154]
                                         glyc3p_e           0     [0, 0.143]
                                         thym_e             0     [0, 0.143]
                                         hxan_e             0     [0, 0.141]
                                         gua_e              0     [0, 0.138]
                                         val__L_e           0     [0, 0.138]
                                         pro__L_e           0     [0, 0.136]
                                         ade_e              0     [0, 0.135]
                                         cys__L_e           0     [0, 0.135]
                                         ptrc_e             0     [0, 0.129]
                                         orn_e              0     [0, 0.128]
                                         agm_e              0     [0, 0.116]
                                         arg__L_e           0     [0, 0.115]
                                         leu__L_e           0     [0, 0.11]
                                         his__L_e           0     [0, 0.104]
                                         15dap_e            0     [0, 0.104]
                                         lys__L_e           0     [0, 0.104]
                                         hxa_e              0     [0, 0.103]
                                         ile__L_e           0     [0, 0.0997]
                                         quin_e             0     [0, 0.0955]
                                         alaala_e           0     [0, 0.092]
                                         g3pe_e             0     [0, 0.0913]
                                         cgly_e             0     [0, 0.0897]
                                         uri_e              0     [0, 0.0826]
                                         g3pg_e             0     [0, 0.0802]
                                         cytd_e             0     [0, 0.0778]
                                         xtsn_e             0     [0, 0.0708]
                                         tyr__L_e           0     [0, 0.0689]
                                         ins_e              0     [0, 0.068]
                                         phe__L_e           0     [0, 0.0662]
                                         adn_e              0     [0, 0.0661]
                                         thymd_e            0     [0, 0.0661]
                                         indole_e           0     [0, 0.0656]
                                         gthrd_e            0     [0, 0.059]
                                         trp__L_e           0     [0, 0.0544]
                                         LalaDgluMdap_e     0     [0, 0.0427]
                                         LalaDgluMdap...    0     [0, 0.0367]
                                         5mtr_e             0     [0, 0.0354]
                                         spmd_e             0     [0, 0.0354]
                                         anhgm_e            0     [0, 0.033]
                                         enter_e            0     [0, 0.0216]
                                         feenter_e          0     [0, 0.0216]
                                         pheme_e            0     [0, 0.0185]
                                         kdo2lipid4_e       0     [0, 0.00703]
                                         lipa_e             0     [0, 0.00533]
                                         enlipa_e           0     [0, 0.0052]
                                         lipa_cold_e        0     [0, 0.00514]
                                         colipa_e           0     [0, 0.00343]
                                         acolipa_e          0     [0, 0.00333]
                                         colipap_e          0     [0, 0.00249]
                                         eca4colipa_e       0     [0, 0.0023]
ERROR:root:'M_sertrna[sec]_c' is not a valid SBML 'SId'.
IN FLUXES          OUT FLUXES    OBJECTIVES
-----------------  ------------  ----------------------
ac_e    10         h2o_e  16.4   Ec_biomass_i...  0.247
o2_e     9.33      co2_e   9.85
h_e      7.73
nh4_e    2.67
pi_e     0.238
so4_e    0.0623
k_e      0.0483
fe2_e    0.00397
mg2_e    0.00214
ca2_e    0.00129
cl_e     0.00129
cu2_e    0.000175
mn2_e    0.000171
zn2_e    8.43e-05
ni2_e    7.98e-05
mobd_e   3.19e-05
IN FLUXES          OUT FLUXES    OBJECTIVES
-----------------  ------------  ----------------------
ac_e    10         h2o_e  16.4   Ec_biomass_i...  0.247
o2_e     9.33      co2_e   9.85
h_e      7.73
nh4_e    2.67
pi_e     0.238
so4_e    0.0623
k_e      0.0483
fe2_e    0.00397
mg2_e    0.00214
ca2_e    0.00129
cl_e     0.00129
cu2_e    0.000175
mn2_e    0.000171
zn2_e    8.43e-05
ni2_e    7.98e-05
mobd_e   3.19e-05
IN FLUXES                                OUT FLUXES                             OBJECTIVES
---------------------------------------  -------------------------------------  ----------------------
id           Flux  Range                 id                 Flux  Range         Ec_biomass_i...  0.247
------  ---------  --------------------  ---------------  ------  ------------
ac_e    10         [9.53, 10]            h2o_e             16.4   [15.2, 514]
o2_e     9.33      [8.91, 258]           co2_e              9.85  [8.97, 10.4]
h_e      7.73      [6.43, 1e+03]         fe3_e              0     [-1.25, 994]
nh4_e    2.67      [2.54, 3.62]          for_e              0     [0, 1.39]
pi_e     0.238     [0.227, 0.37]         urea_e             0     [0, 0.54]
so4_e    0.0623    [0.0592, 0.356]       glyclt_e           0     [0, 0.475]
k_e      0.0483    [0.0458, 0.0483]      gly_e              0     [0, 0.376]
fe2_e    0.00397   [-1.25, 994]          acald_e            0     [0, 0.33]
mg2_e    0.00214   [0.00204, 0.00214]    pyr_e              0     [0, 0.33]
ca2_e    0.00129   [0.00122, 0.00129]    h2s_e              0     [0, 0.297]
cl_e     0.00129   [0.00122, 0.00129]    ala__L_e           0     [0, 0.276]
cu2_e    0.000175  [0.000167, 0.000175]  etoh_e             0     [0, 0.276]
mn2_e    0.000171  [0.000162, 0.000171]  glyc__R_e          0     [0, 0.276]
zn2_e    8.4e-05   [8e-05, 8.4e-05]      lac__D_e           0     [0, 0.276]
ni2_e    8e-05     [7.6e-05, 8e-05]      ser__L_e           0     [0, 0.27]
mobd_e   3.2e-05   [3e-05, 3.2e-05]      asp__L_e           0     [0, 0.247]
                                         mal__L_e           0     [0, 0.237]
                                         asn__L_e           0     [0, 0.228]
                                         ura_e              0     [0, 0.228]
                                         lac__L_e           0     [0, 0.223]
                                         succ_e             0     [0, 0.22]
                                         dha_e              0     [0, 0.216]
                                         glyald_e           0     [0, 0.216]
                                         alltn_e            0     [0, 0.191]
                                         glyc_e             0     [0, 0.191]
                                         akg_e              0     [0, 0.188]
                                         ala__D_e           0     [0, 0.185]
                                         hom__L_e           0     [0, 0.177]
                                         12ppd__R_e         0     [0, 0.177]
                                         12ppd__S_e         0     [0, 0.177]
                                         etha_e             0     [0, 0.175]
                                         4abut_e            0     [0, 0.172]
                                         thr__L_e           0     [0, 0.172]
                                         glu__L_e           0     [0, 0.17]
                                         3hpp_e             0     [0, 0.165]
                                         cit_e              0     [0, 0.163]
                                         acser_e            0     [0, 0.16]
                                         xan_e              0     [0, 0.154]
                                         glyc3p_e           0     [0, 0.143]
                                         thym_e             0     [0, 0.143]
                                         hxan_e             0     [0, 0.141]
                                         gua_e              0     [0, 0.138]
                                         val__L_e           0     [0, 0.138]
                                         pro__L_e           0     [0, 0.136]
                                         ade_e              0     [0, 0.135]
                                         cys__L_e           0     [0, 0.135]
                                         ptrc_e             0     [0, 0.129]
                                         orn_e              0     [0, 0.128]
                                         agm_e              0     [0, 0.116]
                                         arg__L_e           0     [0, 0.115]
                                         leu__L_e           0     [0, 0.11]
                                         his__L_e           0     [0, 0.104]
                                         15dap_e            0     [0, 0.104]
                                         lys__L_e           0     [0, 0.104]
                                         hxa_e              0     [0, 0.103]
                                         ile__L_e           0     [0, 0.0997]
                                         quin_e             0     [0, 0.0955]
                                         alaala_e           0     [0, 0.092]
                                         g3pe_e             0     [0, 0.0913]
                                         cgly_e             0     [0, 0.0897]
                                         uri_e              0     [0, 0.0826]
                                         g3pg_e             0     [0, 0.0802]
                                         cytd_e             0     [0, 0.0778]
                                         xtsn_e             0     [0, 0.0708]
                                         tyr__L_e           0     [0, 0.0689]
                                         ins_e              0     [0, 0.068]
                                         phe__L_e           0     [0, 0.0662]
                                         adn_e              0     [0, 0.0661]
                                         thymd_e            0     [0, 0.0661]
                                         indole_e           0     [0, 0.0656]
                                         gthrd_e            0     [0, 0.059]
                                         trp__L_e           0     [0, 0.0544]
                                         LalaDgluMdap_e     0     [0, 0.0427]
                                         LalaDgluMdap...    0     [0, 0.0367]
                                         5mtr_e             0     [0, 0.0354]
                                         spmd_e             0     [0, 0.0354]
                                         anhgm_e            0     [0, 0.033]
                                         enter_e            0     [0, 0.0216]
                                         feenter_e          0     [0, 0.0216]
                                         pheme_e            0     [0, 0.0185]
                                         kdo2lipid4_e       0     [0, 0.00703]
                                         lipa_e             0     [0, 0.00533]
                                         enlipa_e           0     [0, 0.0052]
                                         lipa_cold_e        0     [0, 0.00514]
                                         colipa_e           0     [0, 0.00343]
                                         acolipa_e          0     [0, 0.00333]
                                         colipap_e          0     [0, 0.00249]
                                         eca4colipa_e       0     [0, 0.0023]
None None
None None
None None
cdiener commented 4 years ago

@Midnighter quick reminder.

Midnighter commented 4 years ago

I started working on this this morning :wink: Must be telepathy.