extractMetModel

PURPOSE ^

extractMetModel Create a subnetwork model around one or more metabolites

SYNOPSIS ^

function metabModel = extractMetModel(model,metabNames,nLayers,allCompFlag,nRxnsMetThr)

DESCRIPTION ^

extractMetModel Create a subnetwork model around one or more metabolites

 metabModel =  extractMetModel(model,metabNames,nLayers,allCompFlag,nRxnsMetThr)

INPUTS
 model             COBRA model structure
 metabNames        Metabolites to build subnetwork model around 
 nLayers           
 allCompFlag       Use all metabolites regardless of compartment

OPTIONAL INPUT
 nRxnsMetThr       Ignore metabolites which appear in more than nRxnMetThr
                   (Default = 100)

OUTPUT
 metabModel        COBRA model around one or more metabolites

 Markus Herrgard 3/1/06

CROSS-REFERENCE INFORMATION ^

This function calls: This function is called by:

SOURCE CODE ^

0001 function metabModel = extractMetModel(model,metabNames,nLayers,allCompFlag,nRxnsMetThr)
0002 %extractMetModel Create a subnetwork model around one or more metabolites
0003 %
0004 % metabModel =  extractMetModel(model,metabNames,nLayers,allCompFlag,nRxnsMetThr)
0005 %
0006 %INPUTS
0007 % model             COBRA model structure
0008 % metabNames        Metabolites to build subnetwork model around
0009 % nLayers
0010 % allCompFlag       Use all metabolites regardless of compartment
0011 %
0012 %OPTIONAL INPUT
0013 % nRxnsMetThr       Ignore metabolites which appear in more than nRxnMetThr
0014 %                   (Default = 100)
0015 %
0016 %OUTPUT
0017 % metabModel        COBRA model around one or more metabolites
0018 %
0019 % Markus Herrgard 3/1/06
0020 
0021 if (nargin < 5)
0022     nRxnsMetThr = 100;
0023 end
0024 
0025 % Filter out high degree metabolites
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);

Generated on Thu 21-Jun-2012 15:39:23 by m2html © 2003