Home > RAVEN > haveFlux.m

haveFlux

PURPOSE ^

haveFlux

SYNOPSIS ^

function I=haveFlux(model,cutOff,rxns)

DESCRIPTION ^

 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

CROSS-REFERENCE INFORMATION ^

This function calls: This function is called by:

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'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)));

Generated on Mon 06-Jan-2014 14:58:12 by m2html © 2005