diff --git a/core/haveFlux.m b/core/haveFlux.m index c60df45f..01215959 100755 --- a/core/haveFlux.m +++ b/core/haveFlux.m @@ -1,23 +1,25 @@ function I=haveFlux(model,cutOff,rxns) % haveFlux -% Checks which reactions can carry a (positive or negative) flux. -% Is used as a faster version of getAllowedBounds if it is only interesting +% Checks which reactions can carry a (positive or negative) flux. Is used +% as a faster version of getAllowedBounds if it is only interesting % whether the reactions can carry a flux or not % +% Input: % model a model structure % cutOff the flux value that a reaction has to carry to be % identified as positive (optional, default 10^-8) % rxns either a cell array of IDs, a logical vector with the -% same number of elements as metabolites in the model, -% of a vector of indexes (optional, default model.rxns) +% same number of elements as metabolites in the model, or a +% vector of indexes (optional, default model.rxns) % -% I logical array with true if the corresponding -% reaction can carry a flux +% Output: +% I logical array with true if the corresponding reaction can +% carry a flux % -% NOTE: If a model has +/- Inf bounds then those are replaced with an -% arbitary large value of +/- 10000 prior to solving +% If a model has +/- Inf bounds then those are replaced with an arbitary +% large value of +/- 10000 prior to solving % -% Usage: I=haveFlux(model,cutOff, rxns) +% Usage: I = haveFlux(model, cutOff, rxns) if nargin<2 cutOff=10^-6; @@ -58,10 +60,7 @@ %Loop through and maximize then minimize each rxn if it does not already %have a flux Z=zeros(numel(smallModel.c),1); -tryMin = find(J == false); -isRev = smallModel.rev - -parfor +hsSolOut=[]; for i=[1 -1] for j=1:numel(J) if J(j)==false @@ -69,7 +68,7 @@ if i==1 || smallModel.rev(mixIndexes(j))~=0 smallModel.c=Z; smallModel.c(mixIndexes(j))=i; - sol=solveLP(smallModel,0); + [sol, hsSolOut]=solveLP(smallModel,0,[],hsSolOut); if any(sol.x) J(abs(sol.x(mixIndexes))>cutOff)=true; end diff --git a/doc/core/haveFlux.html b/doc/core/haveFlux.html index fc497959..4fdbd8b6 100644 --- a/doc/core/haveFlux.html +++ b/doc/core/haveFlux.html @@ -28,24 +28,26 @@

SYNOPSIS ^DESCRIPTION ^

 haveFlux
-   Checks which reactions can carry a (positive or negative) flux.
-   Is used as a faster version of getAllowedBounds if it is only interesting
+   Checks which reactions can carry a (positive or negative) flux. Is used
+   as a faster version of getAllowedBounds if it is only interesting
    whether the reactions can carry a flux or not
 
+ Input:
    model       a model structure
    cutOff      the flux value that a reaction has to carry to be
                identified as positive (optional, default 10^-8)
    rxns        either a cell array of IDs, a logical vector with the
-               same number of elements as metabolites in the model,
-               of a vector of indexes (optional, default model.rxns)
+               same number of elements as metabolites in the model, or a
+               vector of indexes (optional, default model.rxns)
 
-   I           logical array with true if the corresponding
-               reaction can carry a flux
+ Output:
+   I           logical array with true if the corresponding reaction can
+               carry a flux
 
-   NOTE: If a model has +/- Inf bounds then those are replaced with an
-   arbitary large value of +/- 10000 prior to solving
+ If a model has +/- Inf bounds then those are replaced with an arbitary
+ large value of +/- 10000 prior to solving
 
- Usage: I=haveFlux(model,cutOff, rxns)
+ Usage: I = haveFlux(model, cutOff, rxns)

CROSS-REFERENCE INFORMATION ^

@@ -62,86 +64,85 @@

CROSS-REFERENCE INFORMATION ^
 <h2><a name=SOURCE CODE ^

0001 function I=haveFlux(model,cutOff,rxns)
 0002 % haveFlux
-0003 %   Checks which reactions can carry a (positive or negative) flux.
-0004 %   Is used as a faster version of getAllowedBounds if it is only interesting
+0003 %   Checks which reactions can carry a (positive or negative) flux. Is used
+0004 %   as a faster version of getAllowedBounds if it is only interesting
 0005 %   whether the reactions can carry a flux or not
 0006 %
-0007 %   model       a model structure
-0008 %   cutOff      the flux value that a reaction has to carry to be
-0009 %               identified as positive (optional, default 10^-8)
-0010 %   rxns        either a cell array of IDs, a logical vector with the
-0011 %               same number of elements as metabolites in the model,
-0012 %               of a vector of indexes (optional, default model.rxns)
-0013 %
-0014 %   I           logical array with true if the corresponding
-0015 %               reaction can carry a flux
-0016 %
-0017 %   NOTE: If a model has +/- Inf bounds then those are replaced with an
-0018 %   arbitary large value of +/- 10000 prior to solving
-0019 %
-0020 % Usage: I=haveFlux(model,cutOff, rxns)
-0021 
-0022 if nargin<2
-0023     cutOff=10^-6;
-0024 end
-0025 if isempty(cutOff)
-0026     cutOff=10^-6;
-0027 end
-0028 if nargin<3
-0029     rxns=model.rxns;
-0030 elseif ~islogical(rxns) && ~isnumeric(rxns)
-0031     rxns=convertCharArray(rxns);
-0032 end
-0033 
-0034 %This is since we're maximizing for the sum of fluxes, which isn't possible
-0035 %when there are infinite bounds
-0036 model.lb(model.lb==-inf)=-10000;
-0037 model.ub(model.ub==inf)=10000;
-0038 
-0039 %Get the reaction IDs. A bit of an awkward way, but fine.
-0040 indexes=getIndexes(model,rxns,'rxns');
-0041 rxns=model.rxns(indexes);
-0042 
-0043 %First remove all dead-end reactions
-0044 smallModel=simplifyModel(model,false,false,true,true);
-0045 
-0046 %Get the indexes of the reactions in the connected model
-0047 indexes=getIndexes(smallModel,intersect(smallModel.rxns,rxns),'rxns');
-0048 J=false(numel(indexes),1);
-0049 mixIndexes=indexes(randperm(numel(indexes)));
-0050 
-0051 %Maximize for all fluxes first in order to get fewer rxns to test
-0052 smallModel.c=ones(numel(smallModel.c),1);
-0053 sol=solveLP(smallModel);
-0054 if ~isempty(sol.x)
-0055     J(abs(sol.x(mixIndexes))>cutOff)=true;
-0056 end
-0057 
-0058 %Loop through and maximize then minimize each rxn if it does not already
-0059 %have a flux
-0060 Z=zeros(numel(smallModel.c),1);
-0061 tryMin = find(J == false);
-0062 isRev = smallModel.rev
-0063 
-0064 parfor 
-0065 for i=[1 -1]
-0066     for j=1:numel(J)
-0067         if J(j)==false
-0068             %Only minimize for reversible reactions
-0069             if i==1 || smallModel.rev(mixIndexes(j))~=0
-0070                 smallModel.c=Z;
-0071                 smallModel.c(mixIndexes(j))=i;
-0072                 sol=solveLP(smallModel,0);
-0073                 if any(sol.x)
-0074                     J(abs(sol.x(mixIndexes))>cutOff)=true;
-0075                 end
-0076             end
-0077         end
-0078     end
-0079 end
-0080 
-0081 I=ismember(rxns,smallModel.rxns(mixIndexes(J)));
-0082 end
+0007 % Input: +0008 % model a model structure +0009 % cutOff the flux value that a reaction has to carry to be +0010 % identified as positive (optional, default 10^-8) +0011 % rxns either a cell array of IDs, a logical vector with the +0012 % same number of elements as metabolites in the model, or a +0013 % vector of indexes (optional, default model.rxns) +0014 % +0015 % Output: +0016 % I logical array with true if the corresponding reaction can +0017 % carry a flux +0018 % +0019 % If a model has +/- Inf bounds then those are replaced with an arbitary +0020 % large value of +/- 10000 prior to solving +0021 % +0022 % Usage: I = haveFlux(model, cutOff, rxns) +0023 +0024 if nargin<2 +0025 cutOff=10^-6; +0026 end +0027 if isempty(cutOff) +0028 cutOff=10^-6; +0029 end +0030 if nargin<3 +0031 rxns=model.rxns; +0032 elseif ~islogical(rxns) && ~isnumeric(rxns) +0033 rxns=convertCharArray(rxns); +0034 end +0035 +0036 %This is since we're maximizing for the sum of fluxes, which isn't possible +0037 %when there are infinite bounds +0038 model.lb(model.lb==-inf)=-10000; +0039 model.ub(model.ub==inf)=10000; +0040 +0041 %Get the reaction IDs. A bit of an awkward way, but fine. +0042 indexes=getIndexes(model,rxns,'rxns'); +0043 rxns=model.rxns(indexes); +0044 +0045 %First remove all dead-end reactions +0046 smallModel=simplifyModel(model,false,false,true,true); +0047 +0048 %Get the indexes of the reactions in the connected model +0049 indexes=getIndexes(smallModel,intersect(smallModel.rxns,rxns),'rxns'); +0050 J=false(numel(indexes),1); +0051 mixIndexes=indexes(randperm(numel(indexes))); +0052 +0053 %Maximize for all fluxes first in order to get fewer rxns to test +0054 smallModel.c=ones(numel(smallModel.c),1); +0055 sol=solveLP(smallModel); +0056 if ~isempty(sol.x) +0057 J(abs(sol.x(mixIndexes))>cutOff)=true; +0058 end +0059 +0060 %Loop through and maximize then minimize each rxn if it does not already +0061 %have a flux +0062 Z=zeros(numel(smallModel.c),1); +0063 hsSolOut=[]; +0064 for i=[1 -1] +0065 for j=1:numel(J) +0066 if J(j)==false +0067 %Only minimize for reversible reactions +0068 if i==1 || smallModel.rev(mixIndexes(j))~=0 +0069 smallModel.c=Z; +0070 smallModel.c(mixIndexes(j))=i; +0071 [sol, hsSolOut]=solveLP(smallModel,0,[],hsSolOut); +0072 if any(sol.x) +0073 J(abs(sol.x(mixIndexes))>cutOff)=true; +0074 end +0075 end +0076 end +0077 end +0078 end +0079 +0080 I=ismember(rxns,smallModel.rxns(mixIndexes(J))); +0081 end
Generated by m2html © 2005
\ No newline at end of file