Fluxes only sampling at zero

[jacobmartin2022]

I'm having an issue that only comes up with certain reactions - it seems that when a reaction flux is bounded to be only in a single direction (upper bound <= 0 or lower bound >= 0), then that reaction ONLY returns a flux of zero. I believe that the system isn't overconstrained because the COBRA toolbox gpSampler function gives me the results I'd expect using the same input model struct (but is very slow). I've attached the model I'm using, along with the sample structs generated from gpSampleStruct = gpSampler(gapsplit_model, 100); and gapsplitSampleStruct = gapsplit(gapsplit_model, 100)'

Here's my .mat file with my model struct input along with those outputs: 2020_10_27_gapsplit_testing.zip

Thanks, Jake

[pauljensen]

Hi Jake,

Which reactions in the model are giving you problems? I don't see any single-sided reactions with LB or UB at zero; there are some fixed to zero and some with > 0:

>> [gapsplit_model.ub'; gapsplit_model.lb']

ans =

  Columns 1 through 30

   88.7391   88.7391   88.7391   88.7391   88.7391   88.7391   88.7391   88.7391   -1.0000   88.7391   88.7391   88.7391   88.7391   88.7391   88.7391   88.7391   88.7391   88.7391   88.7391   88.7391   88.7391   88.7391   88.7391   88.7391   88.7391   88.7391   88.7391   88.7391   88.7391   88.7391
  -88.7391  -88.7391  -88.7391  -88.7391  -88.7391  -88.7391  -88.7391    1.0000  -88.7391  -88.7391    1.0000  -88.7391  -88.7391  -88.7391  -88.7391  -88.7391  -88.7391  -88.7391    1.0000  -88.7391  -88.7391  -88.7391  -88.7391  -88.7391  -88.7391  -88.7391  -88.7391  -88.7391  -88.7391  -88.7391

  Columns 31 through 60

   88.7391   88.7391   88.7391   88.7391   88.7391   88.7391   88.7391   88.7391   88.7391   88.7391   88.7391   88.7391   88.7391   88.7391   88.7391   88.7391   88.7391   88.7391   88.7391   88.7391   88.7391   88.7391   88.7391   88.7391   88.7391   88.7391   88.7391   88.7391   88.7391   88.7391
  -88.7391    1.0000  -88.7391  -88.7391  -88.7391  -88.7391  -88.7391  -88.7391  -88.7391  -88.7391  -88.7391  -88.7391  -88.7391  -88.7391  -88.7391  -88.7391  -88.7391    1.0000  -88.7391  -88.7391  -88.7391  -88.7391  -88.7391  -88.7391  -88.7391  -88.7391  -88.7391  -88.7391    1.0000  -88.7391

  Columns 61 through 90

   88.7391   88.7391   88.7391   88.7391   88.7391   88.7391   88.7391   88.7391         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0
  -88.7391    1.0000    1.0000  -88.7391  -88.7391  -88.7391  -88.7391  -88.7391         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0

  Columns 91 through 120

         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0
         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0

  Columns 121 through 150

         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0
         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0

  Columns 151 through 180

         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0   88.7391
         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0  -88.7391

  Columns 181 through 195

         0         0         0         0   88.7391         0         0         0         0   88.7391         0         0         0   88.7391   88.7391
         0         0         0         0  -88.7391         0         0         0         0    1.0000         0         0         0  -88.7391  -88.7391

[jacobmartin2022]

Rxns 8 & 9 for example - I forgot that I'd tried to set the bounds to strictly positive or negative by making lb = 1 for forward-only reactions and ub = -1 for reverse only instead of zero, but I had the same behavior of only getting sampled fluxes at 0 when I set ub/lb at zero. So for example, for Rxn #8, my bounds are [1, 88.7], but gapsplit seems to only sample at a flux of zero (same as if I gave bounds of [0, 88.7], though the Cobra gpSample tool seems to be able to sample positive flux values (so I don't think I overconstrained the system somehow).

Let me know if that makes sense or if I can clarify anything, and thanks for the help!

-Jake

On Mon, Nov 9, 2020 at 2:38 PM Paul Jensen notifications@github.com wrote:

Hi Jake,

Which reactions in the model are giving you problems? I don't see any single-sided reactions with LB or UB at zero; there are some fixed to zero and some with > 0:

[gapsplit_model.ub'; gapsplit_model.lb']

ans =

Columns 1 through 30

88.7391 88.7391 88.7391 88.7391 88.7391 88.7391 88.7391 88.7391 -1.0000 88.7391 88.7391 88.7391 88.7391 88.7391 88.7391 88.7391 88.7391 88.7391 88.7391 88.7391 88.7391 88.7391 88.7391 88.7391 88.7391 88.7391 88.7391 88.7391 88.7391 88.7391 -88.7391 -88.7391 -88.7391 -88.7391 -88.7391 -88.7391 -88.7391 1.0000 -88.7391 -88.7391 1.0000 -88.7391 -88.7391 -88.7391 -88.7391 -88.7391 -88.7391 -88.7391 1.0000 -88.7391 -88.7391 -88.7391 -88.7391 -88.7391 -88.7391 -88.7391 -88.7391 -88.7391 -88.7391 -88.7391

Columns 31 through 60

88.7391 88.7391 88.7391 88.7391 88.7391 88.7391 88.7391 88.7391 88.7391 88.7391 88.7391 88.7391 88.7391 88.7391 88.7391 88.7391 88.7391 88.7391 88.7391 88.7391 88.7391 88.7391 88.7391 88.7391 88.7391 88.7391 88.7391 88.7391 88.7391 88.7391 -88.7391 1.0000 -88.7391 -88.7391 -88.7391 -88.7391 -88.7391 -88.7391 -88.7391 -88.7391 -88.7391 -88.7391 -88.7391 -88.7391 -88.7391 -88.7391 -88.7391 1.0000 -88.7391 -88.7391 -88.7391 -88.7391 -88.7391 -88.7391 -88.7391 -88.7391 -88.7391 -88.7391 1.0000 -88.7391

Columns 61 through 90

88.7391 88.7391 88.7391 88.7391 88.7391 88.7391 88.7391 88.7391 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -88.7391 1.0000 1.0000 -88.7391 -88.7391 -88.7391 -88.7391 -88.7391 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0

Columns 91 through 120

     0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0
     0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0

Columns 121 through 150

     0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0
     0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0

Columns 151 through 180

     0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0   88.7391
     0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0         0  -88.7391

Columns 181 through 195

     0         0         0         0   88.7391         0         0         0         0   88.7391         0         0         0   88.7391   88.7391
     0         0         0         0  -88.7391         0         0         0         0    1.0000         0         0         0  -88.7391  -88.7391

— You are receiving this because you authored the thread. Reply to this email directly, view it on GitHub https://github.com/jensenlab/gapsplit/issues/3#issuecomment-724264660, or unsubscribe https://github.com/notifications/unsubscribe-auth/AIRQZPKIELVRCDJSEXOVVYTSPBHNLANCNFSM4TBRHPFQ .

[pauljensen]

I'm not reproducing the bug:

>> s = gapsplit(gapsplit_model,100);
Calculating feasible ranges for variables. (0.67 seconds)
Targeting 1 primary and 3 secondary variables.
Sampling LP model with 70/195 unblocked variables ( 35.90%).

Samples   Coverage   MinGap   Median   MaxGap     Elapsed     Remaining   Infeasible
 10/100     42.88%   0.2472   0.5000   1.0000        0.05          0.44            1
 20/100     57.29%   0.1662   0.3281   1.0000        0.09          0.36            3
 30/100     63.97%   0.1080   0.2500   1.0000        0.13          0.30            4
 40/100     70.74%   0.0879   0.2225   1.0000        0.17          0.26            6
 50/100     74.11%   0.0855   0.1903   1.0000        0.22          0.22           11
 60/100     78.05%   0.0809   0.1655   1.0000        0.27          0.18           15
 70/100     82.13%   0.0733   0.1413   0.5005        0.31          0.13           17
 80/100     83.91%   0.0668   0.1252   0.5005        0.35          0.09           17
 90/100     85.23%   0.0623   0.1203   0.5005        0.40          0.04           19
100/100     86.45%   0.0511   0.1113   0.5005        0.46          0.00           21
>> s.samples(:,[8 9])

ans =

    1.0000  -34.3098
   88.7391  -67.6195
   88.7391  -67.6195
   88.7391  -50.9647
   88.7391  -30.7505
   88.7391  -67.6195
   44.8695  -67.6195
   22.9577  -67.6195
    1.0000  -15.8911
   88.7391  -50.9647
    1.0000  -17.7058
   88.7391  -67.6195
    1.0000  -67.6195
    1.0000  -67.6195
   12.9340  -40.8576
   66.8043  -67.6195
   15.8724   -9.3529
   88.7391  -34.7174
   88.7391  -67.6195
    1.0000  -67.6195
   29.8291  -67.6195
    6.9767  -67.6195
    1.0000  -67.6195
   88.7391  -67.6195
   88.7391  -51.1685
   88.7391  -18.2663
    1.0000  -67.6195
   60.4657  -67.6195
   88.7391  -59.2921
    1.0000  -59.2921
   88.7391  -10.0991
    1.0000  -67.6195
    1.0000  -67.6195
    1.0000  -62.3531
   22.0778  -67.6195
   14.3335  -24.5084
   88.7391  -67.6195
   49.4807  -67.6195
    1.0000  -67.6195
   88.7391  -67.6195
   83.1827  -45.9111
    1.0000  -67.6195
    1.0000  -67.6195
   88.7391   -8.4456
    1.0000  -67.6195
    1.0000  -67.6195
   88.7391  -67.6195
    1.0000  -42.6372
   88.7391  -67.6195
    1.0000  -35.8166
   88.7391  -67.6195
    1.0000  -67.6195
   88.7391  -67.6195
   88.7391  -67.6195
    1.0000  -67.6195
    1.0000  -67.6195
   36.3798  -50.9293
   73.9010  -26.0078
   53.2647  -67.6195
    1.0000  -67.6195
   88.7391  -67.6195
   88.7391   -5.1765
    1.0000  -67.6195
   88.7391  -67.6195
    1.0000  -67.6195
   88.7391  -67.6195
   71.1373  -67.6195
    1.0000  -20.1997
   88.7391  -67.6195
   88.7391  -67.6195
   88.7391  -67.6195
    1.0000  -67.6195
   79.9218  -55.2303
   88.7391  -67.6195
    1.0000  -67.6195
    1.0000  -67.6195
    1.0000  -67.6195
   40.6653  -12.9951
   16.9902  -65.0513
   88.7391  -67.6195
   88.7391  -13.5845
   55.9654  -38.3371
   78.0758  -67.6195
   88.7391  -67.6195
   88.7391  -67.6195
   26.3934  -67.6195
   88.7391  -67.6195
    1.0000  -67.6195
    1.0000  -67.6195
    1.0000  -67.6195
   53.4776  -67.6195
    1.0000  -67.6195
   88.7391  -67.6195
    1.0000  -67.6195
   88.7391  -48.4202
   88.7391  -67.6195
    1.0000  -59.3940
   88.7391  -67.6195
    1.0000  -54.9298
    1.0000  -67.6195

The strange thing is the FVA bounds in the structure returned by gapsplit. For the results you provided, these reactions (8 & 9) are somehow fixed at zero:

>> [gapsplitSampleStruct.maxval'; gapsplitSampleStruct.minval']

ans =

  Columns 1 through 14

   88.7391   59.1594   88.7391   88.7391   88.7391   88.7391  -13.2940         0         0   88.7391         0   88.7391   88.7391   88.7391
  -88.7391  -88.7391  -88.7391  -88.7391  -88.7391  -88.7391  -13.2940         0         0  -88.7391         0  -88.7391  -88.7391  -88.7391

  Columns 15 through 28

   88.7391   88.7391   88.7391   88.7391         0   88.7391   88.7391   88.7391   88.7391   88.7391   88.7391   88.7391  -23.9805   88.7391
  -88.7391  -88.7391  -88.7391  -88.7391         0  -88.7391  -88.7391  -88.7391  -88.7391  -88.7391  -88.7391  -81.9952  -23.9805  -88.7391

  Columns 29 through 42

   88.7391   88.7391   88.7391         0   88.7391   88.7391   88.7391   88.7391   88.7391   88.7391   88.7391   88.7391   88.7391   88.7391
  -88.7391  -88.7391  -88.7391         0  -88.7391  -88.7391  -88.7391  -88.7391  -88.7391  -88.7391  -88.7391  -88.7391  -88.7391  -88.7391

  Columns 43 through 56

   88.7391   88.7391   88.7391   88.7391   88.7391         0   81.9952   88.7391   88.7391   88.7391   88.7391   88.7391    8.2384   88.7391
  -88.7391  -88.7391  -88.7391  -88.7391  -88.7391         0  -88.7391  -88.7391  -44.5000  -88.7391  -88.7391  -88.7391  -79.5007  -88.7391

  Columns 57 through 70

   88.7391   88.7391         0   88.7391   88.7391         0         0   88.7391   88.7391   88.7391   59.1594   88.7391         0         0
  -88.7391  -88.7391         0  -88.7391  -88.7391         0         0  -88.7391  -88.7391  -88.7391  -88.7391  -88.7391         0         0

  Columns 71 through 84

         0         0         0         0         0         0         0         0         0         0         0         0         0         0
         0         0         0         0         0         0         0         0         0         0         0         0         0         0

  Columns 85 through 98

         0         0         0         0         0         0         0         0         0         0         0         0         0         0
         0         0         0         0         0         0         0         0         0         0         0         0         0         0

  Columns 99 through 112

         0         0         0         0         0         0         0         0         0         0         0         0         0         0
         0         0         0         0         0         0         0         0         0         0         0         0         0         0

  Columns 113 through 126

         0         0         0         0         0         0         0         0         0         0         0         0         0         0
         0         0         0         0         0         0         0         0         0         0         0         0         0         0

  Columns 127 through 140

         0         0         0         0         0         0         0         0         0         0         0         0         0         0
         0         0         0         0         0         0         0         0         0         0         0         0         0         0

  Columns 141 through 154

         0         0         0         0         0         0         0         0         0         0         0         0         0         0
         0         0         0         0         0         0         0         0         0         0         0         0         0         0

  Columns 155 through 168

         0         0         0         0         0         0         0         0         0         0         0         0         0         0
         0         0         0         0         0         0         0         0         0         0         0         0         0         0

  Columns 169 through 182

         0         0         0         0         0         0         0         0         0         0         0   88.7391         0         0
         0         0         0         0         0         0         0         0         0         0         0  -88.7391         0         0

  Columns 183 through 195

         0         0   88.7391         0         0         0         0         0         0         0         0   88.7391    0.8498
         0         0  -88.7391         0         0         0         0         0         0         0         0  -88.7391    0.8498

However, this isn't the case when I run it again:

>> [s.maxval'; s.minval']

ans =

  Columns 1 through 14

   88.7391   59.1594   88.7391   88.7391   88.7391   88.7391  -13.2940   88.7391   -1.0000   88.7391   67.6195   88.7391   88.7391   88.7391
  -88.7391  -88.7391  -88.7391  -88.7391  -88.7391  -88.7391  -13.2940    1.0000  -67.6195  -88.7391    1.0000  -88.7391  -88.7391  -88.7391

  Columns 15 through 28

   88.7391   88.7391   88.7391   88.7391   88.7391   88.7391   88.7391   88.7391   88.7391   88.7391   88.7391   88.7391  -23.9805   88.7391
  -88.7391  -88.7391  -88.7391  -88.7391    1.0000  -88.7391  -88.7391  -88.7391  -88.7391  -88.7391  -88.7391  -81.9952  -23.9805  -88.7391

  Columns 29 through 42

   88.7391   88.7391   88.7391   88.7391   88.7391   88.7391   88.7391   88.7391   88.7391   88.7391   88.7391   88.7391   88.7391   88.7391
  -88.7391  -88.7391  -88.7391    1.0000  -88.7391  -88.7391  -88.7391  -88.7391  -88.7391  -88.7391  -88.7391  -88.7391  -88.7391  -88.7391

  Columns 43 through 56

   88.7391   88.7391   88.7391   88.7391   88.7391   88.7391   81.9952   88.7391   88.7391   88.7391   88.7391   88.7391    8.2384   88.7391
  -88.7391  -88.7391  -88.7391  -88.7391  -88.7391    1.0000  -88.7391  -88.7391  -44.5000  -88.7391  -88.7391  -88.7391  -79.5007  -88.7391

  Columns 57 through 70

   88.7391   88.7391   88.7391   88.7391   88.7391   88.7391   88.7391   88.7391   88.7391   88.7391   59.1594   88.7391         0         0
  -88.7391  -88.7391    1.0000  -88.7391  -88.7391    1.0000    1.0000  -88.7391  -88.7391  -88.7391  -88.7391  -88.7391         0         0

  Columns 71 through 84

         0         0         0         0         0         0         0         0         0         0         0         0         0         0
         0         0         0         0         0         0         0         0         0         0         0         0         0         0

  Columns 85 through 98

         0         0         0         0         0         0         0         0         0         0         0         0         0         0
         0         0         0         0         0         0         0         0         0         0         0         0         0         0

  Columns 99 through 112

         0         0         0         0         0         0         0         0         0         0         0         0         0         0
         0         0         0         0         0         0         0         0         0         0         0         0         0         0

  Columns 113 through 126

         0         0         0         0         0         0         0         0         0         0         0         0         0         0
         0         0         0         0         0         0         0         0         0         0         0         0         0         0

  Columns 127 through 140

         0         0         0         0         0         0         0         0         0         0         0         0         0         0
         0         0         0         0         0         0         0         0         0         0         0         0         0         0

  Columns 141 through 154

         0         0         0         0         0         0         0         0         0         0         0         0         0         0
         0         0         0         0         0         0         0         0         0         0         0         0         0         0

  Columns 155 through 168

         0         0         0         0         0         0         0         0         0         0         0         0         0         0
         0         0         0         0         0         0         0         0         0         0         0         0         0         0

  Columns 169 through 182

         0         0         0         0         0         0         0         0         0         0         0   88.7391         0         0
         0         0         0         0         0         0         0         0         0         0         0  -88.7391         0         0

  Columns 183 through 195

         0         0   88.7391         0         0         0         0   88.7391         0         0         0   88.7391    0.8498
         0         0  -88.7391         0         0         0         0    1.0000         0         0         0  -88.7391    0.8498

These fields are calculated by a call to fluxVariability from the Cobra Toolbox; for some reason FVA thought those reactions are fixed. What solver are you using?

[jacobmartin2022]

I have

>> changeCobraSolver(); Defined solvers are: CBT_LP_SOLVER: glpk CBT_MILP_SOLVER: glpk CBT_QP_SOLVER: pdco

[pauljensen]

I switched to those solvers and it runs fine. I'm using Matlab 2018a, which still supported them. What version of Matlab and Cobra Toolbox are you using?

[jacobmartin2022]

I'm on Matlab 2020a, and not sure which version of Cobra but I only downloaded it on this machine a few months ago.

[pauljensen]

If you run fluxVariability on the model, are the min/max fluxes for reactions 8 & 9 zero?

[jacobmartin2022]

Yes, and it also gives me

Warning: The model is infeasible > In findBlockedIrrRxns (line 50) In fluxVariability (line 428)

[pauljensen]

Can you try fluxVariability(gapsplit_model, 0.0)?

[jacobmartin2022]

That’s actually what I did - without passing the ‘0’ argument, I got an error in the deal function assigning outputs.

On Nov 9, 2020, at 5:17 PM, Paul Jensen notifications@github.com wrote:



Can you try fluxVariability(gapsplit_model, 0.0)?

— You are receiving this because you authored the thread. Reply to this email directly, view it on GitHub https://github.com/jensenlab/gapsplit/issues/3#issuecomment-724338939, or unsubscribe https://github.com/notifications/unsubscribe-auth/AIRQZPNZJZ7PI3FYHDPGA7TSPB2AFANCNFSM4TBRHPFQ .

[pauljensen]

Is this the same model that was in the .mat file you originally posted? FVA works fine for me.

[jacobmartin2022]

It is, that's odd. When I run fluxVariability(model), I get

Error using deal (line 37)
The number of outputs should match the number of inputs.

Error in fluxVariability (line 119)
    advind, threads, heuristics, useMtFVA] = deal(funParams{:});

I'll try to update Cobra and see if that makes a difference

[pauljensen]

Ok. The second input to fluxVariability must be required in your version of Cobra.

I think you've found your problem. For whatever reason, fluxVariability is giving the incorrect feasible bounds for your model. You might want to bring this up to the Cobra team. Gapsplit uses these bounds when placing samples, so FVA bounds of [0,0] would put all of the samples to zero.

As a workaround, if you know what the FVA bounds should be, you can pass them to gapsplit with the 'minval' and 'maxval' parameters. If these are given as inputs, gapsplit will use them instead of calling fluxVariability.

[jacobmartin2022]

Sounds good - thanks for all the help!