0001 function metabModel = extractMetModel(model,metabNames,nLayers,allCompFlag,nRxnsMetThr)
0002
0003
0004
0005
0006
0007
0008
0009
0010
0011
0012
0013
0014
0015
0016
0017
0018
0019
0020
0021 if (nargin < 5)
0022 nRxnsMetThr = 100;
0023 end
0024
0025
0026 nRxnsMet = full(sum(model.S~=0,2));
0027 baseMetNames = parseMetNames(model.mets);
0028 for i = 1:length(baseMetNames)
0029 nRxnsMetComp(i) = sum(nRxnsMet(strcmp(baseMetNames,baseMetNames{i})));
0030 end
0031 nRxnsMetComp = nRxnsMetComp';
0032 selLowDegMet = nRxnsMetComp < nRxnsMetThr;
0033 model.S = model.S(selLowDegMet,:);
0034 model.mets = model.mets(selLowDegMet);
0035
0036 if (~iscell(metabNames))
0037 tmpMetName = metabNames;
0038 clear metabNames;
0039 metabNames{1} = tmpMetName;
0040 end
0041
0042 if (allCompFlag)
0043 allMetNames = parseMetNames(model.mets);
0044 metabNames = parseMetNames(metabNames);
0045 else
0046 allMetNames = model.mets;
0047 end
0048
0049 selMets = find(ismember(allMetNames,metabNames));
0050
0051 metS = model.S(selMets,:);
0052 [nMet,tmp] = size(metS);
0053 if (nMet > 1)
0054 selRxns = any(full(metS) ~= 0)';
0055 else
0056 selRxns = (full(metS) ~= 0)';
0057 end
0058 for i = 1:nLayers+1
0059 metS = model.S(selMets,:);
0060 [nMet,tmp] = size(metS);
0061 if (nMet > 1)
0062 selRxns = any(full(metS) ~= 0)';
0063 else
0064 selRxns = (full(metS) ~= 0)';
0065 end
0066
0067 if (isfield(model,'c'))
0068 selRxns = selRxns & ~ (model.c == 1);
0069 end
0070 nextLayerMets = find(any(model.S(:,selRxns) ~= 0,2));
0071 selMets = union(selMets,nextLayerMets);
0072 end
0073
0074 rxnNames = model.rxns(selRxns);
0075 metNames = model.mets(selMets);
0076
0077 metabModel = extractSubNetwork(model,rxnNames,metNames);