constrainFluxData
Constrains fluxes to the data that is provided in the fluxData
structure, which itself is read by loadFluxData from data/fluxData.tsv.
Input:
model an ecModel in GECKO 3 format (with ecModel.ec structure)
fluxData structure with flux data
conds sampling condition
Ptot total protein (g/gDCW)
grRate growth rate (1/h)
exchFluxes exchange fluxes (mmol/gDCWh)
exchMets exchanged metabolites, matching exchFluxes
exchRxnIDs exchange reaction IDs, matching exchMets
condition either index number or name of the sample condition in
fluxData.conds (Optional, default = 1)
maxMinGrowth 'max' if the provided growth rate should be set as
maximum growth rate (= upper bound), or 'min' if it
should be set as minimum growth rate (= lower bound).
The latter option is suitable if minimization of
prot_pool_exchange is used as objective function. (Opt,
default = 'max')
looseStrictFlux how strictly constrained the exchange fluxes should be,
optional, default = 'loose'
'loose' if the exchange fluxes should be constraint
only by the "outer bounds". If exchFluxes(i)
> 0, LB = 0 and UB = exchFluxes(i). If
exchFluxes(i) < 0, LB = exchFluxes(i) and
UB = 0
0-100 LB and UB constraints are set with a specified
percentage of variance around exchFluxes. If 10
is specified, LB = exchFluxes*0.95 and UB =
exchFluxes*1.05. This allows for 10% variance
around the exchFluxes values, but strictly
forces a flux through the exchRxns.
modelAdapter a loaded model adapter (Optional, will otherwise use
the default model adapter)
Output:
model an ecModel where fluxes are constraint
Note: If a provided constraint is either -1000 or 1000, then the function
will update the reaction lower and upper bound to either allow uptake or
excretion, irrespective of what option is given as the looseStrictFlux
parameter.
Usage:
model = constrainFluxData(model, fluxData, condition, maxMinGrowth, looseStrictFlux, modelAdapter)0001 function model = constrainFluxData(model, fluxData, condition, maxMinGrowth, looseStrictFlux, modelAdapter) 0002 % constrainFluxData 0003 % Constrains fluxes to the data that is provided in the fluxData 0004 % structure, which itself is read by loadFluxData from data/fluxData.tsv. 0005 % 0006 % Input: 0007 % model an ecModel in GECKO 3 format (with ecModel.ec structure) 0008 % fluxData structure with flux data 0009 % conds sampling condition 0010 % Ptot total protein (g/gDCW) 0011 % grRate growth rate (1/h) 0012 % exchFluxes exchange fluxes (mmol/gDCWh) 0013 % exchMets exchanged metabolites, matching exchFluxes 0014 % exchRxnIDs exchange reaction IDs, matching exchMets 0015 % condition either index number or name of the sample condition in 0016 % fluxData.conds (Optional, default = 1) 0017 % maxMinGrowth 'max' if the provided growth rate should be set as 0018 % maximum growth rate (= upper bound), or 'min' if it 0019 % should be set as minimum growth rate (= lower bound). 0020 % The latter option is suitable if minimization of 0021 % prot_pool_exchange is used as objective function. (Opt, 0022 % default = 'max') 0023 % looseStrictFlux how strictly constrained the exchange fluxes should be, 0024 % optional, default = 'loose' 0025 % 'loose' if the exchange fluxes should be constraint 0026 % only by the "outer bounds". If exchFluxes(i) 0027 % > 0, LB = 0 and UB = exchFluxes(i). If 0028 % exchFluxes(i) < 0, LB = exchFluxes(i) and 0029 % UB = 0 0030 % 0-100 LB and UB constraints are set with a specified 0031 % percentage of variance around exchFluxes. If 10 0032 % is specified, LB = exchFluxes*0.95 and UB = 0033 % exchFluxes*1.05. This allows for 10% variance 0034 % around the exchFluxes values, but strictly 0035 % forces a flux through the exchRxns. 0036 % modelAdapter a loaded model adapter (Optional, will otherwise use 0037 % the default model adapter) 0038 % 0039 % Output: 0040 % model an ecModel where fluxes are constraint 0041 % 0042 % Note: If a provided constraint is either -1000 or 1000, then the function 0043 % will update the reaction lower and upper bound to either allow uptake or 0044 % excretion, irrespective of what option is given as the looseStrictFlux 0045 % parameter. 0046 % 0047 % Usage: 0048 % model = constrainFluxData(model, fluxData, condition, maxMinGrowth, looseStrictFlux, modelAdapter) 0049 0050 if nargin < 6 || isempty(modelAdapter) 0051 modelAdapter = ModelAdapterManager.getDefault(); 0052 if isempty(modelAdapter) 0053 error('Either send in a modelAdapter or set the default model adapter in the ModelAdapterManager.') 0054 end 0055 end 0056 params = modelAdapter.getParameters(); 0057 0058 if nargin < 5 || isempty(looseStrictFlux) 0059 looseStrictFlux = 'loose'; 0060 end 0061 0062 if nargin < 4 || isempty(maxMinGrowth) 0063 maxMinGrowth = 'max'; 0064 end 0065 0066 if nargin < 2 || isempty(fluxData) 0067 fluxData = loadFluxData(fullfile(params.path,'data','fluxData.tsv'),modelAdapter); 0068 end 0069 0070 if nargin < 3 || isempty(condition) 0071 condition = 1; 0072 elseif ~isnumeric(condition) 0073 idx = find(strcmp(fluxData.conds,condition)); 0074 if isempty(condition) 0075 error(['Condition ' condition ' cannot be found in fluxData']) 0076 else 0077 condition = idx; 0078 end 0079 end 0080 0081 % Select the exchange fluxes of specified condition 0082 fluxData.exchFluxes = fluxData.exchFluxes(condition,:); 0083 0084 %Set original c-source to zero 0085 model = setParam(model,'eq',params.c_source,0); 0086 %Set growth 0087 switch maxMinGrowth 0088 case 'max' 0089 model = setParam(model,'lb',params.bioRxn,0); 0090 model = setParam(model,'ub',params.bioRxn,fluxData.grRate(condition)); 0091 case 'min' 0092 model = setParam(model,'lb',params.bioRxn,fluxData.grRate(condition)); 0093 model = setParam(model,'ub',params.bioRxn,1000); 0094 end 0095 0096 negFlux = lt(fluxData.exchFluxes,0); % less than 0 0097 ub = fluxData.exchFluxes(~negFlux); 0098 posFlux = fluxData.exchRxnIDs(~negFlux); 0099 lb = fluxData.exchFluxes(negFlux); 0100 negFlux = fluxData.exchRxnIDs(negFlux); 0101 0102 switch looseStrictFlux 0103 case 'loose' 0104 model = setParam(model,'lb',negFlux,lb); 0105 model = setParam(model,'ub',negFlux,0); 0106 model = setParam(model,'lb',posFlux,0); 0107 model = setParam(model,'ub',posFlux,ub); 0108 otherwise 0109 model = setParam(model,'var',fluxData.exchRxnIDs,fluxData.exchFluxes,looseStrictFlux); 0110 end 0111 extremeFlux = find(abs(fluxData.exchFluxes)==1000); 0112 if any(extremeFlux) 0113 exchFlux = fluxData.exchFluxes(extremeFlux); 0114 if any(exchFlux==-1000) 0115 model = setParam(model,'lb',fluxData.exchRxnIDs(extremeFlux(exchFlux==-1000)),-1000); 0116 model = setParam(model,'ub',fluxData.exchRxnIDs(extremeFlux(exchFlux==-1000)),0); 0117 end 0118 if any(exchFlux==1000) 0119 model = setParam(model,'lb',fluxData.exchRxnIDs(extremeFlux(exchFlux==1000)),0); 0120 model = setParam(model,'ub',fluxData.exchRxnIDs(extremeFlux(exchFlux==1000)),1000); 0121 end 0122 end 0123 end