Home > RAVEN > getEssentialRxns.m

getEssentialRxns

PURPOSE ^

getEssentialRxns

SYNOPSIS ^

function [essentialRxns, essentialRxnsIndexes]=getEssentialRxns(model,ignoreRxns)

DESCRIPTION ^

 getEssentialRxns
   Calculate the essential reactions for a model to be solvable

   model                   a model structure
   ignoreRxns              cell array of reaction IDs which should not be
                           checked (opt, default {})

   essentialRxns           cell array with the IDs of the essential reactions
   essentialRxnsIndexes    vector with the indexes of the essential reactions

   Essential reactions are those which, when constrained to 0, result in an
   infeasible problem.

   Usage: [essentialRxns, essentialRxnsIndexes]=getEssentialRxns(model,ignoreRxns)

   Rasmus Agren, 2013-08-01

CROSS-REFERENCE INFORMATION ^

This function calls: This function is called by:

SOURCE CODE ^

0001 function [essentialRxns, essentialRxnsIndexes]=getEssentialRxns(model,ignoreRxns)
0002 % getEssentialRxns
0003 %   Calculate the essential reactions for a model to be solvable
0004 %
0005 %   model                   a model structure
0006 %   ignoreRxns              cell array of reaction IDs which should not be
0007 %                           checked (opt, default {})
0008 %
0009 %   essentialRxns           cell array with the IDs of the essential reactions
0010 %   essentialRxnsIndexes    vector with the indexes of the essential reactions
0011 %
0012 %   Essential reactions are those which, when constrained to 0, result in an
0013 %   infeasible problem.
0014 %
0015 %   Usage: [essentialRxns, essentialRxnsIndexes]=getEssentialRxns(model,ignoreRxns)
0016 %
0017 %   Rasmus Agren, 2013-08-01
0018 %
0019 
0020 if nargin<2
0021     ignoreRxns={};
0022 end
0023 
0024 %Too make sure that it doesn't try to optimize for something
0025 model.c=zeros(numel(model.rxns),1);
0026 
0027 %First check that the problem is solvable
0028 [sol hsSolOut]=solveLP(model,1);
0029 
0030 if sol.stat==-1 || isempty(sol.x)
0031     dispEM('No feasible solution to the full model');
0032 end
0033 
0034 %Check which reactions have flux. Only those can be essential. This
0035 %is not the smallest list of reactions, but it's a fast way
0036 rxnsToCheck=setdiff(model.rxns(abs(sol.x)>10^-8),ignoreRxns);
0037 nToCheck=numel(rxnsToCheck);
0038 minimize=true;
0039 while 1
0040     if minimize==true
0041         sol=solveLP(setParam(model,'obj',rxnsToCheck,-1));
0042     else
0043         sol=solveLP(setParam(model,'obj',rxnsToCheck,1));
0044     end
0045     rxnsToCheck=intersect(rxnsToCheck,model.rxns(abs(sol.x)>10^-8));
0046     if numel(rxnsToCheck)>=nToCheck
0047         if minimize==true
0048             minimize=false;
0049         else
0050             break;
0051         end
0052     else
0053         nToCheck=numel(rxnsToCheck);
0054     end
0055 end
0056 
0057 essentialRxns={};
0058 for i=1:numel(rxnsToCheck)
0059    sol=solveLP(setParam(model,'eq',rxnsToCheck(i),0),0,[],hsSolOut);
0060    if sol.stat==-1 || isempty(sol.x)
0061        essentialRxns=[essentialRxns;rxnsToCheck(i)];
0062    end
0063 end
0064 
0065 [crap essentialRxnsIndexes]=ismember(essentialRxns,model.rxns);

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