haveFlux Checks which reactions can carry a (positive or negative) flux. Is used as a faster version of getAllowedBounds if it's only interesting whether the reactions can carry a flux or not model a model structure cutOff the flux value that a reaction has to carry to be identified as positive (opt, 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 (opt, default model.rxns) I logical array with true if the corresponding reaction can carry a flux Usage: I=haveFlux(model,cutOff, rxns) Rasmus Agren, 2013-04-19
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's 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 (opt, 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 (opt, default model.rxns) 0013 % 0014 % I logical array with true if the corresponding 0015 % reaction can carry a flux 0016 % 0017 % Usage: I=haveFlux(model,cutOff, rxns) 0018 % 0019 % Rasmus Agren, 2013-04-19 0020 % 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 end 0031 0032 %Get the reaction IDs. A bit of an awkward way, but fine. 0033 indexes=getIndexes(model,rxns,'rxns'); 0034 rxns=model.rxns(indexes); 0035 0036 %First remove all dead-end reactions 0037 smallModel=simplifyModel(model,false,false,true,true); 0038 0039 %Get the indexes of the reactions in the connected model 0040 indexes=getIndexes(smallModel,intersect(smallModel.rxns,rxns),'rxns'); 0041 J=false(numel(indexes),1); 0042 mixIndexes=indexes(randperm(numel(indexes))); 0043 0044 %Maximize for all fluxes first in order to get fewer rxns to test 0045 smallModel.c=ones(numel(smallModel.c),1); 0046 sol=solveLP(smallModel); 0047 J(abs(sol.x(mixIndexes))>cutOff)=true; 0048 0049 %Loop through and maximize then minimize each rxn if it doesn't already have a flux 0050 Z=zeros(numel(smallModel.c),1); 0051 hsSolOut=[]; 0052 for i=[1 -1] 0053 for j=1:numel(J) 0054 if J(j)==false 0055 %Only minimize for reversible reactions 0056 if i==1 || smallModel.rev(mixIndexes(j))~=0 0057 smallModel.c=Z; 0058 smallModel.c(mixIndexes(j))=i; 0059 [sol hsSolOut]=solveLP(smallModel,0,[],hsSolOut); 0060 if any(sol.x) 0061 J(abs(sol.x(mixIndexes))>cutOff)=true; 0062 end 0063 end 0064 end 0065 end 0066 end 0067 0068 I=ismember(rxns,smallModel.rxns(mixIndexes(J)));