Give parameters to model.optimize()

[KurtzMarina]

Hi, I am new to Cobrapy.
I used this command in the COBRA toolbox in MATLAB:
solution= optimizeCbModel(model,'max',1e-6);

Is there a simple way to do it in Cobrapy also?
When I use model.optimize(), I don't get the same results as in MATLAB.

[cdiener]

Sorry for the late reply. The API of cobrapy and the cobratoolbox are pretty different, so generally you can't really translate code from one to the other. The first two arguments are easy to translate and are the defualt for model.optimize() anyway. The last one I am not sure because it is not apprent from the docs if it is optimized in a single objective with the initial one or a sperate one. Assuming it is separate and is the normal L2 minimization (in that case I am not sure why to use such small weights since that makes the problem unstable), you would have to add that step by hand since we don't have that by default. For instance:

from cobra.io import load_model
from cobra.util import fix_objective_as_constraint

model = load_model("textbook")
fix_objective_as_constraint(model)
obj = sum(1e-6 * rxn.flux_expression**2 for rxn in model.reactions)
model.objective = obj
model.objective_direction = "minimize"
sol = model.optimize()

sol would now contain the equivalent to what optimizeCBModel call did above.

[KurtzMarina]

Thanks for your answer. But when I run your solution I get the following error message:

Traceback (most recent call last):
File "/home/kurtz/testdict.py", line 42, in
cell1.objective = obj
File "/home/kurtz/venv/lib/python3.9/site-packages/cobra/core/model.py", line 1312, in objective
value = self.problem.Objective(value, sloppy=False)
File "/home/kurtz/venv/lib/python3.9/site-packages/optlang/glpk_interface.py", line 310, in init
raise ValueError(
ValueError: GLPK only supports linear objectives. Maximize
1.0e-6*(1.0ACALD - 1.0ACALD_reverse_fda2b)*2 + 1.0e-6(1.0ACALDt - 1.0ACALDt_reverse_858fa)*2 + 1.0e-6(1.0ACKr - 1.0ACKr_reverse_b49c0)*2 + 1.0e-6(1.0ACONTa - 1.0ACONTa_reverse_cad6d)*2 + 1.0e-6(1.0ACONTb - 1.0ACONTb_reverse_e198a)*2 + 1.0e-6(1.0ACt2r - 1.0ACt2r_reverse_64e05)*2 + 1.0e-6(1.0ADK1 - 1.0ADK1_reverse_a6f90)*2 + 1.0e-6(1.0AKGDH - 1.0AKGDH_reverse_08bdc)*2 + 1.0e-6(1.0AKGt2r - 1.0AKGt2r_reverse_5d500)*2 + 1.0e-6(1.0ALCD2x - 1.0ALCD2x_reverse_5d107)*2 + 1.0e-6(1.0ATPS4r - 1.0ATPS4r_reverse_64306)*2 + 1.0e-6(1.0Biomass_Ecoli_core - 1.0Biomass_Ecoli_core_reverse_2cdba)*2 + 1.0e-6(1.0CO2t - 1.0CO2t_reverse_7c42f)*2 + 1.0e-6(1.0CS - 1.0CS_reverse_8d7e9)*2 + 1.0e-6(1.0CYTBD - 1.0CYTBD_reverse_611ba)*2 + 1.0e-6(1.0D_LACt2 - 1.0D_LACt2_reverse_d0dc8)*2 + 1.0e-6(1.0ENO - 1.0ENO_reverse_40eea)*2 + 1.0e-6(1.0ETOHt2r - 1.0ETOHt2r_reverse_b3d52)*2 + 1.0e-6(1.0EX_ac_e - 1.0EX_ac_e_reverse_0be96)*2 + 1.0e-6(1.0EX_acald_e - 1.0EX_acald_e_reverse_c096e)*2 + 1.0e-6(1.0EX_akg_e - 1.0EX_akg_e_reverse_70d85)*2 + 1.0e-6(1.0EX_co2_e - 1.0EX_co2_e_reverse_d0466)*2 + 1.0e-6(1.0EX_etoh_e - 1.0EX_etoh_e_reverse_cc64f)*2 + 1.0e-6(1.0EX_for_e - 1.0EX_for_e_reverse_23269)*2 + 1.0e-6(1.0EX_fru_e - 1.0EX_fru_e_reverse_c3828)*2 + 1.0e-6(1.0EX_fum_e - 1.0EX_fum_e_reverse_e3432)*2 + 1.0e-6(1.0EX_glc__D_e - 1.0EX_glc__D_e_reverse_af641)*2 + 1.0e-6(1.0EX_gln__L_e - 1.0EX_gln__L_e_reverse_6a1a1)*2 + 1.0e-6(1.0EX_glu__L_e - 1.0EX_glu__L_e_reverse_42f6c)*2 + 1.0e-6(1.0EX_h2o_e - 1.0EX_h2o_e_reverse_3ced4)*2 + 1.0e-6(1.0EX_h_e - 1.0EX_h_e_reverse_3e0c5)*2 + 1.0e-6(1.0EX_lac__D_e - 1.0EX_lac__D_e_reverse_f95b4)*2 + 1.0e-6(1.0EX_mal__L_e - 1.0EX_mal__L_e_reverse_af154)*2 + 1.0e-6(1.0EX_nh4_e - 1.0EX_nh4_e_reverse_f9cc6)*2 + 1.0e-6(1.0EX_o2_e - 1.0EX_o2_e_reverse_efa94)*2 + 1.0e-6(1.0EX_pi_e - 1.0EX_pi_e_reverse_1fb09)*2 + 1.0e-6(1.0EX_pyr_e - 1.0EX_pyr_e_reverse_1f6de)*2 + 1.0e-6(1.0EX_succ_e - 1.0EX_succ_e_reverse_a9039)*2 + 1.0e-6(1.0FBA - 1.0FBA_reverse_84806)*2 + 1.0e-6(1.0FBP - 1.0FBP_reverse_bf2c9)*2 + 1.0e-6(1.0FORt2 - 1.0FORt2_reverse_89839)*2 + 1.0e-6(1.0FORti - 1.0FORti_reverse_18c06)*2 + 1.0e-6(1.0FRD7 - 1.0FRD7_reverse_e0172)*2 + 1.0e-6(1.0FRUpts2 - 1.0FRUpts2_reverse_58d3a)*2 + 1.0e-6(1.0FUM - 1.0FUM_reverse_d3642)*2 + 1.0e-6(1.0FUMt2_2 - 1.0FUMt2_2_reverse_7b4b6)*2 + 1.0e-6(1.0G6PDH2r - 1.0G6PDH2r_reverse_19ddf)*2 + 1.0e-6(1.0GAPD - 1.0GAPD_reverse_459c1)*2 + 1.0e-6(1.0GLCpts - 1.0GLCpts_reverse_a52ae)*2 + 1.0e-6(1.0GLNS - 1.0GLNS_reverse_59581)*2 + 1.0e-6(1.0GLNabc - 1.0GLNabc_reverse_1d82a)*2 + 1.0e-6(1.0GLUDy - 1.0GLUDy_reverse_fa4e7)*2 + 1.0e-6(1.0GLUN - 1.0GLUN_reverse_4ccdb)*2 + 1.0e-6(1.0GLUSy - 1.0GLUSy_reverse_6a00f)*2 + 1.0e-6(1.0GLUt2r - 1.0GLUt2r_reverse_3e88a)*2 + 1.0e-6(1.0GND - 1.0GND_reverse_eec5c)*2 + 1.0e-6(1.0H2Ot - 1.0H2Ot_reverse_aa560)*2 + 1.0e-6(1.0ICDHyr - 1.0ICDHyr_reverse_7f84b)*2 + 1.0e-6(1.0ICL - 1.0ICL_reverse_2f27e)*2 + 1.0e-6(1.0LDH_D - 1.0LDH_D_reverse_f8507)*2 + 1.0e-6(1.0MALS - 1.0MALS_reverse_d7382)*2 + 1.0e-6(1.0MALt2_2 - 1.0MALt2_2_reverse_a635f)*2 + 1.0e-6(1.0MDH - 1.0MDH_reverse_ee52c)*2 + 1.0e-6(1.0ME1 - 1.0ME1_reverse_9736c)*2 + 1.0e-6(1.0ME2 - 1.0ME2_reverse_2b0a2)*2 + 1.0e-6(1.0NADH16 - 1.0NADH16_reverse_330a2)*2 + 1.0e-6(1.0NADTRHD - 1.0NADTRHD_reverse_49725)*2 + 1.0e-6(1.0NH4t - 1.0NH4t_reverse_551ee)*2 + 1.0e-6(1.0O2t - 1.0O2t_reverse_4d957)*2 + 1.0e-6(1.0PDH - 1.0PDH_reverse_ca160)*2 + 1.0e-6(1.0PFK - 1.0PFK_reverse_d24a6)*2 + 1.0e-6(1.0PFL - 1.0PFL_reverse_af9ec)*2 + 1.0e-6(1.0PGI - 1.0PGI_reverse_27efc)*2 + 1.0e-6(1.0PGK - 1.0PGK_reverse_02696)*2 + 1.0e-6(1.0PGL - 1.0PGL_reverse_2bb6b)*2 + 1.0e-6(1.0PGM - 1.0PGM_reverse_fc9af)*2 + 1.0e-6(1.0PIt2r - 1.0PIt2r_reverse_1cd61)*2 + 1.0e-6(1.0PPC - 1.0PPC_reverse_e854a)*2 + 1.0e-6(1.0PPCK - 1.0PPCK_reverse_2557d)*2 + 1.0e-6(1.0PPS - 1.0PPS_reverse_1c319)*2 + 1.0e-6(1.0PTAr - 1.0PTAr_reverse_fce15)*2 + 1.0e-6(1.0PYK - 1.0PYK_reverse_bc8ff)*2 + 1.0e-6(1.0PYRt2 - 1.0PYRt2_reverse_a6e24)*2 + 1.0e-6(1.0RPE - 1.0RPE_reverse_a1b04)*2 + 1.0e-6(1.0RPI - 1.0RPI_reverse_853a1)*2 + 1.0e-6(1.0SUCCt2_2 - 1.0SUCCt2_2_reverse_206ca)*2 + 1.0e-6(1.0SUCCt3 - 1.0SUCCt3_reverse_7d82a)*2 + 1.0e-6(1.0SUCDi - 1.0SUCDi_reverse_480f4)*2 + 1.0e-6(1.0SUCOAS - 1.0SUCOAS_reverse_22958)*2 + 1.0e-6(1.0TALA - 1.0TALA_reverse_adfda)*2 + 1.0e-6(1.0THD2 - 1.0THD2_reverse_f65dd)*2 + 1.0e-6(1.0TKT1 - 1.0TKT1_reverse_a1021)*2 + 1.0e-6(1.0TKT2 - 1.0TKT2_reverse_7ebc7)*2 + 1.0e-6(1.0TPI - 1.0TPI_reverse_c2c3b)*2 + 1.0e-6(1.0ATPM - 1.0ATPM_reverse_5b752)**2 is not linear.

In Matlab with the cobratoolbox there is no problem with GLPK as solver.
I found this in the cobratoolbox documentation:

The remaining options work only with a valid QP solver:

0 Minimises the squared Euclidean Norm of internal fluxes. Typically 1e-6 works well.
min 1/2v′∗v
s.t. Sv=dxdt
c^Tv=f
lb≤v≤ub

So it's not the normal L2 minimization and that's why there isn't a way to do it in COBRApy?

[cdiener]

You do need a QP solver exactly the same as in Matlab (hence the note in their docs). GLPK does not support quadratic problems. You probably have another QP solver installed in your Matlab setup.

So it's not the normal L2 minimization and that's why there isn't a way to do it in COBRApy?

It is the normal L2 minimization in cobrapy as well, that is what the code above does. But minimizing the L2 norm is a quadratic programming problem. For cobrapy you can use cplex, gurobi or the free osqp solver (pip install osqp).