Error in SolveCobraQP

[jacobmartin2022]

I'm trying to use gapsplit to indirectly sample metabolite concentrations by sampling "cumulative" fluxes and then applying those to initial concentrations. I'm then applying upper and lower bounds on b in Sv = b to make sure my final concentrations are within experimental ranges. Since it seems like the Cobra Toolbox (and gapsplit) only allow either upper bounds or lower bounds, what I've been doing was doubling the linear problem, so I have `[S, -1.S] v <= [b_ub, -1 . b_lb]` However, when I do this, it seems unable to find a solution, even though I've confirmed that there's a feasible solution with

[x_init,resnorm,residual,exitflag,output,lambda] = ...
            lsqlin(C, d, A, b, [], [], ...
            tot_flux_lb, tot_flux_ub, [], lsqlin_options);

that there is a solution.

I attached an example workspace where I run

flux_struct_test = gapsplit(gs_model, n_samples, ...
                'primary', 'max',...
                'secondaryFrac', gs_metab_secondary_frac,...
                'minval', gs_minval,...
                'maxval', gs_maxval);

I pulled gs_minval and gs_maxval from a single homemade FVA, but it seems the same if I don't pass those varargin to gapsplit.

If I only solve the lower-bounded problem (bottom half of my double-bounded problem), gapsplit seems able to sample, but very slowly.

Lastly, the gurobi solver seems to be having more issues than pdco for this problem. Using the gurobi solver, I get the exception "'Gurobi error 10003: Element 0 of a double array is Nan or Inf.'" at line 713 of SolveCobraQP. If I'm using pdco, gapsplit just seems unable to find a solution (keeps entering the statement at line 263 of gapsplit and restarts the loop).

[jacobmartin2022]

2020_12_14_gapsplit_conc_sampling_test.zip

[jacobmartin2022]

Here is a updated workspace file for the problem written with slack variables.

Using pdco for the QP solver and running

flux_struct = gapsplit(gs_model, 10, 'minval', FVA_minval, 'maxval', FVA_maxval);

gives

Warning: Matrix is close to singular or badly scaled. Results may be inaccurate. RCOND =  4.100597e-18. 
> In pdco (line 1058)
  In solveCobraQP (line 609)
  In gapsplit (line 253)

If I use gurobi for the QP solver, I get stuck in the catch statement at line 714 of solveCobraQP with the exception message "Gurobi error 10003: Element 1 of a double array is Nan or Inf."

Again, when I solve the problem using lsqlin -

[x_init,~,~,exitflag] = lsqlin(C, d, [], [], gs_model.S, gs_model.b, gs_model.lb, gs_model.ub, [], lsqlin_options);

I get a solution that satisfies the system, so I don't think it's over-constrained.

Let me know if you need any other info or would want me to try something else.

2021_01_06_gapsplit_conc_sampling_test.zip

[pauljensen]

Hi Jake,

I just pushed a possible fix. Can you give it a try?

Thanks, Paul

On Jan 6, 2021, at 5:18 PM, jacobmartin2022 notifications@github.com<mailto:notifications@github.com> wrote:

Here is a updated workspace file for the problem written with slack variables.

Using pdco for the QP solver and running

flux_struct = gapsplit(gs_model, 10, 'minval', FVA_minval, 'maxval', FVA_maxval);

gives

Warning: Matrix is close to singular or badly scaled. Results may be inaccurate. RCOND = 4.100597e-18.

In pdco (line 1058) In solveCobraQP (line 609) In gapsplit (line 253)

If I use gurobi for the QP solver, I get stuck in the catch statement at line 714 of solveCobraQP with the exception message "Gurobi error 10003: Element 1 of a double array is Nan or Inf."

Again, when I solve the problem using lsqlin -

[x_init,~,~,exitflag] = lsqlin(C, d, [], [], gs_model.S, gs_model.b, gs_model.lb, gs_model.ub, [], lsqlin_options);

I get a solution that satisfies the system, so I don't think it's over-constrained.

Let me know if you need any other info or would want me to try something else.

2021_01_06_gapsplit_conc_sampling_test.ziphttps://urldefense.com/v3/__https://github.com/jensenlab/gapsplit/files/5778659/2021_01_06_gapsplit_conc_sampling_test.zip__;!!DZ3fjg!vGrMhNXdpj2r4XMK4jEiWkensEfuhQh2zx0UahVeaPSMiJ8djkdZwd2I2wR0rm6j$

— You are receiving this because you are subscribed to this thread. Reply to this email directly, view it on GitHubhttps://urldefense.com/v3/__https://github.com/jensenlab/gapsplit/issues/4*issuecomment-755773920__;Iw!!DZ3fjg!vGrMhNXdpj2r4XMK4jEiWkensEfuhQh2zx0UahVeaPSMiJ8djkdZwd2I29XiJgif$, or unsubscribehttps://urldefense.com/v3/__https://github.com/notifications/unsubscribe-auth/AAYNH3YEFLNFQLWP3URS56TSYTVVJANCNFSM4VHRNVSQ__;!!DZ3fjg!vGrMhNXdpj2r4XMK4jEiWkensEfuhQh2zx0UahVeaPSMiJ8djkdZwd2I23xcyX_x$.

[jacobmartin2022]

Some of those bounds were very large because those are species (water for example) whose change I don't want to constrain. I changed it so now those rows in the S matrix are just zeros - would this be a better way to treat those cases?

I pulled the fix and with pcdo for the QP solver, I get the same message.

With gurobi, in silveCobraQP I now get a value of 'INF_OR_UNBD' for resultgurobi.status at line 713, and then again the same status at line 766 - the comment here says this means the problem is infeasible.

2021_01_07_gapsplit_conc_sampling_test.zip

I've attached the updated model struct with the zero rows for ignored species.