0001 function [finalState,finalInputs1States,finalInputs2States] = solveBooleanRegModel(model,initialState,inputs1States,inputs2States)
0002
0003
0004
0005
0006
0007
0008
0009
0010
0011
0012
0013
0014
0015
0016
0017
0018
0019
0020
0021
0022
0023
0024 finalInputs1States = [];
0025 [selExc,selUpt] = findExcRxns(model);
0026 for i = 1:length(model.regulatoryInputs1)
0027 met = model.regulatoryInputs1{i};
0028 fullS = full(model.S);
0029 rxnID = intersect(find(fullS(findMetIDs(model,met),:)),find(selExc));
0030 if model.lb(rxnID) < 0
0031 finalInputs1States(i,1) = true;
0032 else
0033 finalInputs1States(i,1) = false;
0034 end
0035 end
0036
0037
0038 drGenes = {};
0039 for i = 1:length(model.regulatoryGenes)
0040 if initialState(i) == false
0041 drGenes{length(drGenes)+1} = model.regulatoryGenes{i};
0042 end
0043 end
0044 drGenes = intersect(model.genes,drGenes);
0045 modelDR = deleteModelGenes(model,drGenes);
0046 fbasolDR = optimizeCbModel(modelDR,'max',true);
0047
0048
0049 finalInputs2States = [];
0050 if ~any(fbasolDR.x)
0051 finalInputs2States = false.*ones(length(model.regulatoryInputs2),1);
0052 else
0053 for i = 1:length(model.regulatoryInputs2)
0054 rxnFlux = fbasolDR.x(findRxnIDs(model,model.regulatoryInputs2{i}));
0055 if rxnFlux == 0
0056 finalInputs2States(i,1) = false;
0057 else
0058 finalInputs2States(i,1) = true;
0059 end
0060 end
0061 end
0062
0063
0064
0065 ruleList = parseRegulatoryRules(model);
0066
0067 geneStates = [];
0068 for i = 1:length(model.regulatoryRules)
0069 geneStates(i,1) = eval(ruleList{i});
0070 end
0071
0072 finalState = geneStates;
0073
0074
0075
0076 function ruleList = parseRegulatoryRules(model)
0077
0078 ruleList = cell(length(model.regulatoryRules),1);
0079
0080 for i = 1:length(model.regulatoryRules)
0081 fields = splitString(model.regulatoryRules{i},'[\s~|&()]');
0082 newFields = fields;
0083 for j = 1:length(fields)
0084 word = fields{j};
0085 if strcmp(word,'true')
0086 newFields{j} = 'true';
0087 elseif strcmp(word,'false')
0088 newFields{j} = 'false';
0089 else
0090 if any(strcmp(word,model.regulatoryGenes))
0091 index = find(strcmp(word,model.regulatoryGenes));
0092 newFields{j} = ['initialState(',num2str(index),')'];
0093 elseif any(strcmp(word,model.regulatoryInputs1))
0094 index = find(strcmp(word,model.regulatoryInputs1));
0095 newFields{j} = ['inputs1States(',num2str(index),')'];
0096 elseif any(strcmp(word,model.regulatoryInputs2))
0097 index = find(strcmp(word,model.regulatoryInputs2));
0098 newFields{j} = ['inputs2States(',num2str(index),')'];
0099 else
0100 newFields{j} = '';
0101 end
0102 end
0103 end
0104 newRule = model.regulatoryRules{i};
0105 for j = 1:length(fields)
0106 newRule = strrep(newRule,fields{j},newFields{j});
0107 end
0108 ruleList{i} = newRule;
0109 end
0110
0111