diff --git a/core/checkModelStruct.m b/core/checkModelStruct.m index 6cadafcf..3fb1e15f 100755 --- a/core/checkModelStruct.m +++ b/core/checkModelStruct.m @@ -124,6 +124,10 @@ function checkModelStruct(model,throwErrors,trimWarnings) EM='If "grRules" field exists, the model should also contain a "genes" field'; dispEM(EM,throwErrors); else + %Erroneous grRules that start/end with OR/AND + EM='The following reaction(s) have grRules that start or end with ''OR'' or ''AND'':'; + dispEM(EM,throwErrors,model.rxns(startsWith(model.grRules,{'or ','and '}) | endsWith(model.grRules,{' or',' and'})),trimWarnings); + %grRules that are not in genes field geneList = getGenesFromGrRules(model.grRules); geneList = setdiff(unique(geneList),model.genes); if ~isempty(geneList) diff --git a/core/getExchangeRxns.m b/core/getExchangeRxns.m index fa197b21..d058f689 100755 --- a/core/getExchangeRxns.m +++ b/core/getExchangeRxns.m @@ -26,8 +26,8 @@ % positive flux value would imply uptake, % but reaction bounds are not considered % 'out' reactions where the boundary metabolite -% is the substrate of the reaction, a -% positive flux value would imply uptake, +% is the product of the reaction, a +% negative flux value would imply uptake, % but reaction bounds are not considered. % % Output: diff --git a/core/setParam.m b/core/setParam.m index e24ad2ca..734853b7 100755 --- a/core/setParam.m +++ b/core/setParam.m @@ -39,6 +39,9 @@ else rxnList=convertCharArray(rxnList); end +if isempty(rxnList) + return; +end %Allow to set several parameters to the same value if numel(rxnList)~=numel(params) && numel(params)~=1 @@ -59,23 +62,16 @@ %as we do not want to throw errors if matches fail indexes=zeros(numel(rxnList),1); -for i=1:numel(rxnList) - index=find(strcmp(rxnList{i},model.rxns),1); - if ~isempty(index) - indexes(i)=index; - else - indexes(i)=-1; - EM=['Reaction ' rxnList{i} ' is not present in the reaction list']; - dispEM(EM,false); - end +[Lia,Locb] = ismember(rxnList,model.rxns); +indexes(Lia) = Locb; +if any(~Lia) + params(~Lia)=[]; + indexes(~Lia)=[]; + paramType(~Lia)=[]; + dispEM('Reactions not present in model, will be ignored:',false,rxnLise(~Lia)); end -%Remove the reactions that were not found -params(indexes==-1)=[]; -indexes(indexes==-1)=[]; -paramType(indexes==-1)=[]; %Change the parameters - if ~isempty(indexes) if contains(paramType,'obj') model.c=zeros(numel(model.c),1); % parameter is changed, not added diff --git a/doc/core/checkModelStruct.html b/doc/core/checkModelStruct.html index 2f029e1f..920a7c8b 100644 --- a/doc/core/checkModelStruct.html +++ b/doc/core/checkModelStruct.html @@ -183,331 +183,335 @@

SOURCE CODE ^'If "grRules" field exists, the model should also contain a "genes" field'; 0125 dispEM(EM,throwErrors); 0126 else -0127 geneList = getGenesFromGrRules(model.grRules); -0128 geneList = setdiff(unique(geneList),model.genes); -0129 if ~isempty(geneList) -0130 problemGrRules = model.rxns(contains(model.grRules,geneList)); -0131 problemGrRules = strjoin(problemGrRules(:),'; '); -0132 EM=['The reaction(s) "' problemGrRules '" contain the following genes in its "grRules" field, but these are not in the "genes" field:']; -0133 dispEM(EM,throwErrors,geneList); -0134 end -0135 end -0136 end -0137 if isfield(model,'rxnComps') -0138 if ~isnumeric(model.rxnComps) -0139 EM='The "rxnComps" field must be of type "double"'; -0140 dispEM(EM,throwErrors); -0141 end -0142 end -0143 if isfield(model,'inchis') -0144 if ~iscellstr(model.inchis) -0145 EM='The "inchis" field must be a cell array of strings'; -0146 dispEM(EM,throwErrors); -0147 end -0148 end -0149 if isfield(model,'metSmiles') -0150 if ~iscellstr(model.metSmiles) -0151 EM='The "metSmiles" field must be a cell array of strings'; -0152 dispEM(EM,throwErrors); -0153 end -0154 end -0155 if isfield(model,'metFormulas') -0156 if ~iscellstr(model.metFormulas) -0157 EM='The "metFormulas" field must be a cell array of strings'; -0158 dispEM(EM,throwErrors); -0159 end -0160 end -0161 if isfield(model,'metCharges') -0162 if ~isnumeric(model.metCharges) -0163 EM='The "metCharges" field must be a double'; -0164 dispEM(EM,throwErrors); -0165 end -0166 end -0167 if isfield(model,'metDeltaG') -0168 if ~isnumeric(model.metDeltaG) -0169 EM='The "metDeltaG" field must be a double'; -0170 dispEM(EM,throwErrors); -0171 end -0172 end -0173 if isfield(model,'subSystems') -0174 for i=1:numel(model.subSystems) -0175 if ~iscell(model.subSystems{i,1}) -0176 EM='The "subSystems" field must be a cell array'; -0177 dispEM(EM,throwErrors); -0178 end -0179 end -0180 end -0181 if isfield(model,'eccodes') -0182 if ~iscellstr(model.eccodes) -0183 EM='The "eccodes" field must be a cell array of strings'; -0184 dispEM(EM,throwErrors); -0185 end -0186 end -0187 if isfield(model,'unconstrained') -0188 if ~isnumeric(model.unconstrained) -0189 EM='The "unconstrained" field must be of type "double"'; -0190 dispEM(EM,throwErrors); -0191 end -0192 end -0193 if isfield(model,'rxnNotes') -0194 if ~iscellstr(model.rxnNotes) -0195 EM='The "rxnNotes" field must be a cell array of strings'; -0196 dispEM(EM,throwErrors); -0197 end -0198 end -0199 if isfield(model,'rxnReferences') -0200 if ~iscellstr(model.rxnReferences) -0201 EM='The "rxnReferences" field must be a cell array of strings'; -0202 dispEM(EM,throwErrors); -0203 end -0204 end -0205 if isfield(model,'rxnConfidenceScores') -0206 if ~isnumeric(model.rxnConfidenceScores) -0207 EM='The "rxnConfidenceScores" field must be a double'; -0208 dispEM(EM,throwErrors); -0209 end -0210 end -0211 if isfield(model,'rxnDeltaG') -0212 if ~isnumeric(model.rxnDeltaG) -0213 EM='The "rxnDeltaG" field must be a double'; -0214 dispEM(EM,throwErrors); -0215 end -0216 end -0217 -0218 %Empty strings -0219 if isempty(model.id) -0220 EM='The "id" field cannot be empty'; -0221 dispEM(EM,throwErrors); -0222 end -0223 if any(cellfun(@isempty,model.rxns)) -0224 EM='The model contains empty reaction IDs'; +0127 %Erroneous grRules that start/end with OR/AND +0128 EM='The following reaction(s) have grRules that start or end with ''OR'' or ''AND'':'; +0129 dispEM(EM,throwErrors,model.rxns(startsWith(model.grRules,{'or ','and '}) | endsWith(model.grRules,{' or',' and'})),trimWarnings); +0130 %grRules that are not in genes field +0131 geneList = getGenesFromGrRules(model.grRules); +0132 geneList = setdiff(unique(geneList),model.genes); +0133 if ~isempty(geneList) +0134 problemGrRules = model.rxns(contains(model.grRules,geneList)); +0135 problemGrRules = strjoin(problemGrRules(:),'; '); +0136 EM=['The reaction(s) "' problemGrRules '" contain the following genes in its "grRules" field, but these are not in the "genes" field:']; +0137 dispEM(EM,throwErrors,geneList); +0138 end +0139 end +0140 end +0141 if isfield(model,'rxnComps') +0142 if ~isnumeric(model.rxnComps) +0143 EM='The "rxnComps" field must be of type "double"'; +0144 dispEM(EM,throwErrors); +0145 end +0146 end +0147 if isfield(model,'inchis') +0148 if ~iscellstr(model.inchis) +0149 EM='The "inchis" field must be a cell array of strings'; +0150 dispEM(EM,throwErrors); +0151 end +0152 end +0153 if isfield(model,'metSmiles') +0154 if ~iscellstr(model.metSmiles) +0155 EM='The "metSmiles" field must be a cell array of strings'; +0156 dispEM(EM,throwErrors); +0157 end +0158 end +0159 if isfield(model,'metFormulas') +0160 if ~iscellstr(model.metFormulas) +0161 EM='The "metFormulas" field must be a cell array of strings'; +0162 dispEM(EM,throwErrors); +0163 end +0164 end +0165 if isfield(model,'metCharges') +0166 if ~isnumeric(model.metCharges) +0167 EM='The "metCharges" field must be a double'; +0168 dispEM(EM,throwErrors); +0169 end +0170 end +0171 if isfield(model,'metDeltaG') +0172 if ~isnumeric(model.metDeltaG) +0173 EM='The "metDeltaG" field must be a double'; +0174 dispEM(EM,throwErrors); +0175 end +0176 end +0177 if isfield(model,'subSystems') +0178 for i=1:numel(model.subSystems) +0179 if ~iscell(model.subSystems{i,1}) +0180 EM='The "subSystems" field must be a cell array'; +0181 dispEM(EM,throwErrors); +0182 end +0183 end +0184 end +0185 if isfield(model,'eccodes') +0186 if ~iscellstr(model.eccodes) +0187 EM='The "eccodes" field must be a cell array of strings'; +0188 dispEM(EM,throwErrors); +0189 end +0190 end +0191 if isfield(model,'unconstrained') +0192 if ~isnumeric(model.unconstrained) +0193 EM='The "unconstrained" field must be of type "double"'; +0194 dispEM(EM,throwErrors); +0195 end +0196 end +0197 if isfield(model,'rxnNotes') +0198 if ~iscellstr(model.rxnNotes) +0199 EM='The "rxnNotes" field must be a cell array of strings'; +0200 dispEM(EM,throwErrors); +0201 end +0202 end +0203 if isfield(model,'rxnReferences') +0204 if ~iscellstr(model.rxnReferences) +0205 EM='The "rxnReferences" field must be a cell array of strings'; +0206 dispEM(EM,throwErrors); +0207 end +0208 end +0209 if isfield(model,'rxnConfidenceScores') +0210 if ~isnumeric(model.rxnConfidenceScores) +0211 EM='The "rxnConfidenceScores" field must be a double'; +0212 dispEM(EM,throwErrors); +0213 end +0214 end +0215 if isfield(model,'rxnDeltaG') +0216 if ~isnumeric(model.rxnDeltaG) +0217 EM='The "rxnDeltaG" field must be a double'; +0218 dispEM(EM,throwErrors); +0219 end +0220 end +0221 +0222 %Empty strings +0223 if isempty(model.id) +0224 EM='The "id" field cannot be empty'; 0225 dispEM(EM,throwErrors); 0226 end -0227 if any(cellfun(@isempty,model.mets)) -0228 EM='The model contains empty metabolite IDs'; +0227 if any(cellfun(@isempty,model.rxns)) +0228 EM='The model contains empty reaction IDs'; 0229 dispEM(EM,throwErrors); 0230 end -0231 if any(cellfun(@isempty,model.comps)) -0232 EM='The model contains empty compartment IDs'; +0231 if any(cellfun(@isempty,model.mets)) +0232 EM='The model contains empty metabolite IDs'; 0233 dispEM(EM,throwErrors); 0234 end -0235 EM='The following metabolites have empty names:'; -0236 dispEM(EM,throwErrors,model.mets(cellfun(@isempty,model.metNames)),trimWarnings); -0237 -0238 if isfield(model,'genes') -0239 if any(cellfun(@isempty,model.genes)) -0240 EM='The model contains empty gene IDs'; -0241 dispEM(EM,throwErrors); -0242 end -0243 end -0244 -0245 %Validate format of ids -0246 fields = {'rxns';'mets';'comps';'genes'}; -0247 fieldNames = {'reaction';'metabolite';'compartment';'gene'}; -0248 fieldPrefix = {'R_';'M_';'C_';'G_'}; -0249 for i=1:numel(fields) -0250 try -0251 numIDs = ~startsWith(model.(fields{i}),regexpPattern('^[a-zA-Z_]')); -0252 catch -0253 numIDs = []; -0254 end -0255 if any(numIDs) -0256 EM = ['The following ' fieldNames{i} ' identifiers do not start '... -0257 'with a letter or _ (conflicting with SBML specifications). '... -0258 'This does not impact RAVEN functionality, but be aware that '... -0259 'exportModel will automatically add ' fieldPrefix{i} ... -0260 ' prefixes to all ' fieldNames{i} ' identifiers:']; -0261 dispEM(EM,false,{model.(fields{i}){numIDs}},trimWarnings); -0262 end -0263 end -0264 -0265 %Duplicates -0266 EM='The following reaction IDs are duplicates:'; -0267 dispEM(EM,throwErrors,model.rxns(duplicates(model.rxns)),trimWarnings); -0268 EM='The following metabolite IDs are duplicates:'; -0269 dispEM(EM,throwErrors,model.mets(duplicates(model.mets)),trimWarnings); -0270 EM='The following compartment IDs are duplicates:'; -0271 dispEM(EM,throwErrors,model.comps(duplicates(model.comps)),trimWarnings); -0272 if isfield(model,'genes') -0273 EM='The following genes are duplicates:'; -0274 dispEM(EM,throwErrors,model.genes(duplicates(model.genes)),trimWarnings); -0275 end -0276 metInComp=strcat(model.metNames,'[',model.comps(model.metComps),']'); -0277 EM='The following metabolites already exist in the same compartment:'; -0278 dispEM(EM,throwErrors,metInComp(duplicates(metInComp)),trimWarnings); -0279 -0280 %Elements never used (print only as warnings -0281 EM='The following reactions are empty (no involved metabolites):'; -0282 dispEM(EM,false,model.rxns(~any(model.S,1)),trimWarnings); -0283 EM='The following metabolites are never used in a reaction:'; -0284 dispEM(EM,false,model.mets(~any(model.S,2)),trimWarnings); -0285 if isfield(model,'genes') -0286 EM='The following genes are not associated to a reaction:'; -0287 dispEM(EM,false,model.genes(~any(model.rxnGeneMat,1)),trimWarnings); -0288 end -0289 I=true(numel(model.comps),1); -0290 I(model.metComps)=false; -0291 EM='The following compartments contain no metabolites:'; -0292 dispEM(EM,false,model.comps(I),trimWarnings); -0293 -0294 %Contradicting bounds -0295 EM='The following reactions have contradicting bounds (lower bound is higher than upper bound):'; -0296 dispEM(EM,throwErrors,model.rxns(model.lb>model.ub),trimWarnings); -0297 EM='The following reactions have lower and upper bounds that indicate reversibility, but are indicated as irreversible in model.rev:'; -0298 dispEM(EM,false,model.rxns(model.lb < 0 & model.ub > 0 & model.rev==0),trimWarnings); -0299 -0300 %Multiple or no objective functions not allowed in SBML L3V1 FBCv2 -0301 if numel(find(model.c))>1 -0302 EM='Multiple objective functions found. This might be intended, but results in FBCv2 non-compliant SBML file when exported'; -0303 dispEM(EM,false,model.rxns(find(model.c)),trimWarnings); -0304 elseif ~any(model.c) -0305 EM='No objective function found. This might be intended, but results in FBCv2 non-compliant SBML file when exported'; -0306 dispEM(EM,false); -0307 end -0308 -0309 %Mapping of compartments -0310 if isfield(model,'compOutside') -0311 EM='The following compartments are in "compOutside" but not in "comps":'; -0312 dispEM(EM,throwErrors,setdiff(model.compOutside,[{''};model.comps]),trimWarnings); -0313 end -0314 -0315 %Met names which start with number -0316 I=false(numel(model.metNames),1); -0317 for i=1:numel(model.metNames) -0318 index=strfind(model.metNames{i},' '); -0319 if any(index) -0320 if any(str2double(model.metNames{i}(1:index(1)-1))) -0321 I(i)=true; -0322 end -0323 end -0324 end -0325 EM='The following metabolite names begin with a number directly followed by space, which could potentially cause problems:'; -0326 dispEM(EM,false,model.metNames(I),trimWarnings); -0327 -0328 %Non-parseable composition -0329 if isfield(model,'metFormulas') -0330 [~, ~, exitFlag]=parseFormulas(model.metFormulas,true,false); -0331 EM='The composition for the following metabolites could not be parsed:'; -0332 dispEM(EM,false,model.mets(exitFlag==-1),trimWarnings); -0333 end -0334 -0335 %Check if there are metabolites with different names but the same MIRIAM -0336 %codes -0337 if isfield(model,'metMiriams') -0338 miriams=containers.Map(); -0339 for i=1:numel(model.mets) -0340 if ~isempty(model.metMiriams{i}) -0341 %Loop through and add for each miriam -0342 for j=1:numel(model.metMiriams{i}.name) -0343 %Get existing metabolite indexes -0344 current=strcat(model.metMiriams{i}.name{j},'/',model.metMiriams{i}.value{j}); -0345 if isKey(miriams,current) -0346 existing=miriams(current); -0347 else -0348 existing=[]; -0349 end -0350 miriams(current)=[existing;i]; -0351 end -0352 end -0353 end -0354 -0355 %Get all keys -0356 allMiriams=keys(miriams); -0357 -0358 hasMultiple=false(numel(allMiriams),1); -0359 for i=1:numel(allMiriams) -0360 if numel(miriams(allMiriams{i}))>1 -0361 %Check if they all have the same name -0362 if numel(unique(model.metNames(miriams(allMiriams{i}))))>1 -0363 if ~regexp(allMiriams{i},'^sbo\/SBO:') % SBO terms are expected to be multiple -0364 hasMultiple(i)=true; -0365 end -0366 end -0367 end -0368 end -0369 -0370 %Print output -0371 EM='The following MIRIAM strings are associated to more than one unique metabolite name:'; -0372 dispEM(EM,false,allMiriams(hasMultiple),trimWarnings); -0373 end -0374 -0375 %Check if there are metabolites with different names but the same InChI -0376 %codes -0377 if isfield(model,'inchis') -0378 inchis=containers.Map(); -0379 for i=1:numel(model.mets) -0380 if ~isempty(model.inchis{i}) -0381 %Get existing metabolite indexes -0382 if isKey(inchis,model.inchis{i}) -0383 existing=inchis(model.inchis{i}); -0384 else -0385 existing=[]; -0386 end -0387 inchis(model.inchis{i})=[existing;i]; -0388 end -0389 end -0390 -0391 %Get all keys -0392 allInchis=keys(inchis); -0393 -0394 hasMultiple=false(numel(allInchis),1); -0395 for i=1:numel(allInchis) -0396 if numel(inchis(allInchis{i}))>1 -0397 %Check if they all have the same name -0398 if numel(unique(model.metNames(inchis(allInchis{i}))))>1 -0399 hasMultiple(i)=true; -0400 end -0401 end -0402 end -0403 -0404 %Print output -0405 EM='The following InChI strings are associated to more than one unique metabolite name:'; -0406 dispEM(EM,false,allInchis(hasMultiple),trimWarnings); -0407 end -0408 -0409 % %Check if there are metabolites with different names but the same SMILES -0410 % if isfield(model,'metSmiles') -0411 % metSmiles=containers.Map(); -0412 % for i=1:numel(model.mets) -0413 % if ~isempty(model.metSmiles{i}) -0414 % %Get existing metabolite indexes -0415 % if isKey(metSmiles,model.metSmiles{i}) -0416 % existing=metSmiles(model.metSmiles{i}); -0417 % else -0418 % existing=[]; -0419 % end -0420 % metSmiles(model.metSmiles{i})=[existing;i]; -0421 % end -0422 % end -0423 % -0424 % %Get all keys -0425 % allmetSmiles=keys(metSmiles); -0426 % -0427 % hasMultiple=false(numel(metSmiles),1); -0428 % for i=1:numel(metSmiles) -0429 % if numel(metSmiles(metSmiles{i}))>1 -0430 % %Check if they all have the same name -0431 % if numel(unique(model.metNames(metSmiles(allmetSmiles{i}))))>1 -0432 % hasMultiple(i)=true; -0433 % end -0434 % end -0435 % end -0436 % -0437 % %Print output -0438 % EM='The following metSmiles strings are associated to more than one unique metabolite name:'; -0439 % dispEM(EM,false,allmetSmiles(hasMultiple),trimWarnings); -0440 % end -0441 end -0442 -0443 function I=duplicates(strings) -0444 I=false(numel(strings),1); -0445 [J, K]=unique(strings); -0446 if numel(J)~=numel(strings) -0447 L=1:numel(strings); -0448 L(K)=[]; -0449 I(L)=true; -0450 end -0451 end +0235 if any(cellfun(@isempty,model.comps)) +0236 EM='The model contains empty compartment IDs'; +0237 dispEM(EM,throwErrors); +0238 end +0239 EM='The following metabolites have empty names:'; +0240 dispEM(EM,throwErrors,model.mets(cellfun(@isempty,model.metNames)),trimWarnings); +0241 +0242 if isfield(model,'genes') +0243 if any(cellfun(@isempty,model.genes)) +0244 EM='The model contains empty gene IDs'; +0245 dispEM(EM,throwErrors); +0246 end +0247 end +0248 +0249 %Validate format of ids +0250 fields = {'rxns';'mets';'comps';'genes'}; +0251 fieldNames = {'reaction';'metabolite';'compartment';'gene'}; +0252 fieldPrefix = {'R_';'M_';'C_';'G_'}; +0253 for i=1:numel(fields) +0254 try +0255 numIDs = ~startsWith(model.(fields{i}),regexpPattern('^[a-zA-Z_]')); +0256 catch +0257 numIDs = []; +0258 end +0259 if any(numIDs) +0260 EM = ['The following ' fieldNames{i} ' identifiers do not start '... +0261 'with a letter or _ (conflicting with SBML specifications). '... +0262 'This does not impact RAVEN functionality, but be aware that '... +0263 'exportModel will automatically add ' fieldPrefix{i} ... +0264 ' prefixes to all ' fieldNames{i} ' identifiers:']; +0265 dispEM(EM,false,{model.(fields{i}){numIDs}},trimWarnings); +0266 end +0267 end +0268 +0269 %Duplicates +0270 EM='The following reaction IDs are duplicates:'; +0271 dispEM(EM,throwErrors,model.rxns(duplicates(model.rxns)),trimWarnings); +0272 EM='The following metabolite IDs are duplicates:'; +0273 dispEM(EM,throwErrors,model.mets(duplicates(model.mets)),trimWarnings); +0274 EM='The following compartment IDs are duplicates:'; +0275 dispEM(EM,throwErrors,model.comps(duplicates(model.comps)),trimWarnings); +0276 if isfield(model,'genes') +0277 EM='The following genes are duplicates:'; +0278 dispEM(EM,throwErrors,model.genes(duplicates(model.genes)),trimWarnings); +0279 end +0280 metInComp=strcat(model.metNames,'[',model.comps(model.metComps),']'); +0281 EM='The following metabolites already exist in the same compartment:'; +0282 dispEM(EM,throwErrors,metInComp(duplicates(metInComp)),trimWarnings); +0283 +0284 %Elements never used (print only as warnings +0285 EM='The following reactions are empty (no involved metabolites):'; +0286 dispEM(EM,false,model.rxns(~any(model.S,1)),trimWarnings); +0287 EM='The following metabolites are never used in a reaction:'; +0288 dispEM(EM,false,model.mets(~any(model.S,2)),trimWarnings); +0289 if isfield(model,'genes') +0290 EM='The following genes are not associated to a reaction:'; +0291 dispEM(EM,false,model.genes(~any(model.rxnGeneMat,1)),trimWarnings); +0292 end +0293 I=true(numel(model.comps),1); +0294 I(model.metComps)=false; +0295 EM='The following compartments contain no metabolites:'; +0296 dispEM(EM,false,model.comps(I),trimWarnings); +0297 +0298 %Contradicting bounds +0299 EM='The following reactions have contradicting bounds (lower bound is higher than upper bound):'; +0300 dispEM(EM,throwErrors,model.rxns(model.lb>model.ub),trimWarnings); +0301 EM='The following reactions have lower and upper bounds that indicate reversibility, but are indicated as irreversible in model.rev:'; +0302 dispEM(EM,false,model.rxns(model.lb < 0 & model.ub > 0 & model.rev==0),trimWarnings); +0303 +0304 %Multiple or no objective functions not allowed in SBML L3V1 FBCv2 +0305 if numel(find(model.c))>1 +0306 EM='Multiple objective functions found. This might be intended, but results in FBCv2 non-compliant SBML file when exported'; +0307 dispEM(EM,false,model.rxns(find(model.c)),trimWarnings); +0308 elseif ~any(model.c) +0309 EM='No objective function found. This might be intended, but results in FBCv2 non-compliant SBML file when exported'; +0310 dispEM(EM,false); +0311 end +0312 +0313 %Mapping of compartments +0314 if isfield(model,'compOutside') +0315 EM='The following compartments are in "compOutside" but not in "comps":'; +0316 dispEM(EM,throwErrors,setdiff(model.compOutside,[{''};model.comps]),trimWarnings); +0317 end +0318 +0319 %Met names which start with number +0320 I=false(numel(model.metNames),1); +0321 for i=1:numel(model.metNames) +0322 index=strfind(model.metNames{i},' '); +0323 if any(index) +0324 if any(str2double(model.metNames{i}(1:index(1)-1))) +0325 I(i)=true; +0326 end +0327 end +0328 end +0329 EM='The following metabolite names begin with a number directly followed by space, which could potentially cause problems:'; +0330 dispEM(EM,false,model.metNames(I),trimWarnings); +0331 +0332 %Non-parseable composition +0333 if isfield(model,'metFormulas') +0334 [~, ~, exitFlag]=parseFormulas(model.metFormulas,true,false); +0335 EM='The composition for the following metabolites could not be parsed:'; +0336 dispEM(EM,false,model.mets(exitFlag==-1),trimWarnings); +0337 end +0338 +0339 %Check if there are metabolites with different names but the same MIRIAM +0340 %codes +0341 if isfield(model,'metMiriams') +0342 miriams=containers.Map(); +0343 for i=1:numel(model.mets) +0344 if ~isempty(model.metMiriams{i}) +0345 %Loop through and add for each miriam +0346 for j=1:numel(model.metMiriams{i}.name) +0347 %Get existing metabolite indexes +0348 current=strcat(model.metMiriams{i}.name{j},'/',model.metMiriams{i}.value{j}); +0349 if isKey(miriams,current) +0350 existing=miriams(current); +0351 else +0352 existing=[]; +0353 end +0354 miriams(current)=[existing;i]; +0355 end +0356 end +0357 end +0358 +0359 %Get all keys +0360 allMiriams=keys(miriams); +0361 +0362 hasMultiple=false(numel(allMiriams),1); +0363 for i=1:numel(allMiriams) +0364 if numel(miriams(allMiriams{i}))>1 +0365 %Check if they all have the same name +0366 if numel(unique(model.metNames(miriams(allMiriams{i}))))>1 +0367 if ~regexp(allMiriams{i},'^sbo\/SBO:') % SBO terms are expected to be multiple +0368 hasMultiple(i)=true; +0369 end +0370 end +0371 end +0372 end +0373 +0374 %Print output +0375 EM='The following MIRIAM strings are associated to more than one unique metabolite name:'; +0376 dispEM(EM,false,allMiriams(hasMultiple),trimWarnings); +0377 end +0378 +0379 %Check if there are metabolites with different names but the same InChI +0380 %codes +0381 if isfield(model,'inchis') +0382 inchis=containers.Map(); +0383 for i=1:numel(model.mets) +0384 if ~isempty(model.inchis{i}) +0385 %Get existing metabolite indexes +0386 if isKey(inchis,model.inchis{i}) +0387 existing=inchis(model.inchis{i}); +0388 else +0389 existing=[]; +0390 end +0391 inchis(model.inchis{i})=[existing;i]; +0392 end +0393 end +0394 +0395 %Get all keys +0396 allInchis=keys(inchis); +0397 +0398 hasMultiple=false(numel(allInchis),1); +0399 for i=1:numel(allInchis) +0400 if numel(inchis(allInchis{i}))>1 +0401 %Check if they all have the same name +0402 if numel(unique(model.metNames(inchis(allInchis{i}))))>1 +0403 hasMultiple(i)=true; +0404 end +0405 end +0406 end +0407 +0408 %Print output +0409 EM='The following InChI strings are associated to more than one unique metabolite name:'; +0410 dispEM(EM,false,allInchis(hasMultiple),trimWarnings); +0411 end +0412 +0413 % %Check if there are metabolites with different names but the same SMILES +0414 % if isfield(model,'metSmiles') +0415 % metSmiles=containers.Map(); +0416 % for i=1:numel(model.mets) +0417 % if ~isempty(model.metSmiles{i}) +0418 % %Get existing metabolite indexes +0419 % if isKey(metSmiles,model.metSmiles{i}) +0420 % existing=metSmiles(model.metSmiles{i}); +0421 % else +0422 % existing=[]; +0423 % end +0424 % metSmiles(model.metSmiles{i})=[existing;i]; +0425 % end +0426 % end +0427 % +0428 % %Get all keys +0429 % allmetSmiles=keys(metSmiles); +0430 % +0431 % hasMultiple=false(numel(metSmiles),1); +0432 % for i=1:numel(metSmiles) +0433 % if numel(metSmiles(metSmiles{i}))>1 +0434 % %Check if they all have the same name +0435 % if numel(unique(model.metNames(metSmiles(allmetSmiles{i}))))>1 +0436 % hasMultiple(i)=true; +0437 % end +0438 % end +0439 % end +0440 % +0441 % %Print output +0442 % EM='The following metSmiles strings are associated to more than one unique metabolite name:'; +0443 % dispEM(EM,false,allmetSmiles(hasMultiple),trimWarnings); +0444 % end +0445 end +0446 +0447 function I=duplicates(strings) +0448 I=false(numel(strings),1); +0449 [J, K]=unique(strings); +0450 if numel(J)~=numel(strings) +0451 L=1:numel(strings); +0452 L(K)=[]; +0453 I(L)=true; +0454 end +0455 end
Generated by m2html © 2005
\ No newline at end of file diff --git a/doc/core/getExchangeRxns.html b/doc/core/getExchangeRxns.html index 78d4e7a3..09e722ab 100644 --- a/doc/core/getExchangeRxns.html +++ b/doc/core/getExchangeRxns.html @@ -54,8 +54,8 @@

DESCRIPTION ^SOURCE CODE ^% positive flux value would imply uptake, 0027 % but reaction bounds are not considered 0028 % 'out' reactions where the boundary metabolite -0029 % is the substrate of the reaction, a -0030 % positive flux value would imply uptake, +0029 % is the product of the reaction, a +0030 % negative flux value would imply uptake, 0031 % but reaction bounds are not considered. 0032 % 0033 % Output: diff --git a/doc/core/setParam.html b/doc/core/setParam.html index 73b63b7e..96f0a107 100644 --- a/doc/core/setParam.html +++ b/doc/core/setParam.html @@ -111,96 +111,92 @@

SOURCE CODE ^else 0040 rxnList=convertCharArray(rxnList); 0041 end -0042 -0043 %Allow to set several parameters to the same value -0044 if numel(rxnList)~=numel(params) && numel(params)~=1 -0045 EM='The number of parameter values and the number of reactions must be the same'; -0046 dispEM(EM); -0047 end -0048 -0049 if length(rxnList)>1 -0050 if length(paramType)==1 -0051 paramType(1:length(rxnList))=paramType; -0052 end -0053 if length(params)==1 -0054 params(1:length(rxnList))=params; +0042 if isempty(rxnList) +0043 return; +0044 end +0045 +0046 %Allow to set several parameters to the same value +0047 if numel(rxnList)~=numel(params) && numel(params)~=1 +0048 EM='The number of parameter values and the number of reactions must be the same'; +0049 dispEM(EM); +0050 end +0051 +0052 if length(rxnList)>1 +0053 if length(paramType)==1 +0054 paramType(1:length(rxnList))=paramType; 0055 end -0056 end -0057 -0058 %Find the indexes for the reactions in rxnList. Do not use getIndexes -0059 %as we do not want to throw errors if matches fail -0060 indexes=zeros(numel(rxnList),1); -0061 -0062 for i=1:numel(rxnList) -0063 index=find(strcmp(rxnList{i},model.rxns),1); -0064 if ~isempty(index) -0065 indexes(i)=index; -0066 else -0067 indexes(i)=-1; -0068 EM=['Reaction ' rxnList{i} ' is not present in the reaction list']; -0069 dispEM(EM,false); -0070 end -0071 end -0072 -0073 %Remove the reactions that were not found -0074 params(indexes==-1)=[]; -0075 indexes(indexes==-1)=[]; -0076 paramType(indexes==-1)=[]; -0077 %Change the parameters -0078 -0079 if ~isempty(indexes) -0080 if contains(paramType,'obj') -0081 model.c=zeros(numel(model.c),1); % parameter is changed, not added -0082 end -0083 for j=1:length(paramType) -0084 if strcmpi(paramType{j},'eq') +0056 if length(params)==1 +0057 params(1:length(rxnList))=params; +0058 end +0059 end +0060 +0061 %Find the indexes for the reactions in rxnList. Do not use getIndexes +0062 %as we do not want to throw errors if matches fail +0063 indexes=zeros(numel(rxnList),1); +0064 +0065 [Lia,Locb] = ismember(rxnList,model.rxns); +0066 indexes(Lia) = Locb; +0067 if any(~Lia) +0068 params(~Lia)=[]; +0069 indexes(~Lia)=[]; +0070 paramType(~Lia)=[]; +0071 dispEM('Reactions not present in model, will be ignored:',false,rxnLise(~Lia)); +0072 end +0073 +0074 %Change the parameters +0075 if ~isempty(indexes) +0076 if contains(paramType,'obj') +0077 model.c=zeros(numel(model.c),1); % parameter is changed, not added +0078 end +0079 for j=1:length(paramType) +0080 if strcmpi(paramType{j},'eq') +0081 model.lb(indexes(j))=params(j); +0082 model.ub(indexes(j))=params(j); +0083 end +0084 if strcmpi(paramType{j},'lb') 0085 model.lb(indexes(j))=params(j); -0086 model.ub(indexes(j))=params(j); -0087 end -0088 if strcmpi(paramType{j},'lb') -0089 model.lb(indexes(j))=params(j); -0090 end -0091 if strcmpi(paramType{j},'ub') -0092 model.ub(indexes(j))=params(j); -0093 end -0094 if strcmpi(paramType{j},'obj') -0095 model.c(indexes(j))=params(j); +0086 end +0087 if strcmpi(paramType{j},'ub') +0088 model.ub(indexes(j))=params(j); +0089 end +0090 if strcmpi(paramType{j},'obj') +0091 model.c(indexes(j))=params(j); +0092 end +0093 if strcmpi(paramType{j},'rev') +0094 %Non-zero values are interpreted as reversible +0095 model.rev(indexes(j))=params(j)~=0; 0096 end -0097 if strcmpi(paramType{j},'rev') -0098 %Non-zero values are interpreted as reversible -0099 model.rev(indexes(j))=params(j)~=0; -0100 end -0101 if strcmpi(paramType{j},'var') -0102 if params(j) < 0 -0103 model.lb(indexes(j)) = params(j) * (1+var/200); -0104 model.ub(indexes(j)) = params(j) * (1-var/200); -0105 else -0106 model.lb(indexes(j)) = params(j) * (1-var/200); -0107 model.ub(indexes(j)) = params(j) * (1+var/200); -0108 end -0109 end -0110 if strcmpi(paramType{j},'unc') -0111 if isfield(model.annotation,'defaultLB') -0112 lb = model.annotation.defaultLB; -0113 else -0114 lb = -1000; -0115 end -0116 if isfield(model.annotation,'defaultUB') -0117 ub = model.annotation.defaultUB; -0118 else -0119 ub = 1000; -0120 end -0121 model.lb(indexes(j)) = lb; -0122 model.ub(indexes(j)) = ub; -0123 end -0124 end -0125 end -0126 if any(ismember(paramType,{'lb','ub','unc'})) -0127 invalidBound = model.lb(indexes) > model.ub(indexes); -0128 if any(invalidBound) -0129 error(['Invalid set of bounds for reaction(s): ', strjoin(model.rxns(indexes(invalidBound)),', '), '.']) -0130 end -0131 end +0097 if strcmpi(paramType{j},'var') +0098 if params(j) < 0 +0099 model.lb(indexes(j)) = params(j) * (1+var/200); +0100 model.ub(indexes(j)) = params(j) * (1-var/200); +0101 else +0102 model.lb(indexes(j)) = params(j) * (1-var/200); +0103 model.ub(indexes(j)) = params(j) * (1+var/200); +0104 end +0105 end +0106 if strcmpi(paramType{j},'unc') +0107 if isfield(model.annotation,'defaultLB') +0108 lb = model.annotation.defaultLB; +0109 else +0110 lb = -1000; +0111 end +0112 if isfield(model.annotation,'defaultUB') +0113 ub = model.annotation.defaultUB; +0114 else +0115 ub = 1000; +0116 end +0117 model.lb(indexes(j)) = lb; +0118 model.ub(indexes(j)) = ub; +0119 end +0120 end +0121 end +0122 if any(ismember(paramType,{'lb','ub','unc'})) +0123 invalidBound = model.lb(indexes) > model.ub(indexes); +0124 if any(invalidBound) +0125 error(['Invalid set of bounds for reaction(s): ', strjoin(model.rxns(indexes(invalidBound)),', '), '.']) +0126 end +0127 end
Generated by m2html © 2005
\ No newline at end of file diff --git a/doc/io/exportToExcelFormat.html b/doc/io/exportToExcelFormat.html index bbb64211..70b7328f 100644 --- a/doc/io/exportToExcelFormat.html +++ b/doc/io/exportToExcelFormat.html @@ -41,9 +41,6 @@

DESCRIPTION ^SOURCE CODE ^% sorted alphabetically by their identifiers (optional, 0014 % default false) 0015 % -0016 % No checks are made regarding the correctness of the model. Use -0017 % checkModelStruct to identify problems in the model structure. -0018 % -0019 % Usage: exportToExcelFormat(model, fileName, sortIds) -0020 -0021 if nargin<2 || isempty(fileName) -0022 [fileName, pathName] = uiputfile('*.xlsx', 'Select file for model export',[model.id '.xlsx']); -0023 if fileName == 0 -0024 error('You should provide a file location') -0025 else -0026 fileName = fullfile(pathName,fileName); -0027 end -0028 end -0029 fileName=char(fileName); -0030 if nargin<3 -0031 sortIds=false; +0016 % Usage: exportToExcelFormat(model, fileName, sortIds) +0017 +0018 if nargin<2 || isempty(fileName) +0019 [fileName, pathName] = uiputfile('*.xlsx', 'Select file for model export',[model.id '.xlsx']); +0020 if fileName == 0 +0021 error('You should provide a file location') +0022 else +0023 fileName = fullfile(pathName,fileName); +0024 end +0025 end +0026 fileName=char(fileName); +0027 if nargin<3 +0028 sortIds=false; +0029 end +0030 if sortIds==true +0031 model=sortIdentifiers(model); 0032 end -0033 if sortIds==true -0034 model=sortIdentifiers(model); -0035 end -0036 -0037 addList = matlab.addons.installedAddons; -0038 if any(strcmpi(addList.Name,'Text Analytics Toolbox')) -0039 error(['exportToExcelFormat is incompatible with MATLAB Text Analytics Toolbox. ' ... -0040 'Further instructions => https://github.com/SysBioChalmers/RAVEN/issues/55#issuecomment-1514369299']) -0041 end -0042 -0043 [~, A, B]=fileparts(fileName); -0044 -0045 %If a path was used call on exportToTabDelimited instead -0046 if ~any(A) || ~any(B) -0047 exportToTabDelimited(model,fileName); -0048 return; -0049 end -0050 -0051 if ~strcmpi(B,'.xlsx') -0052 EM='As of RAVEN version 1.9, only export to xlsx format is supported'; -0053 dispEM(EM); -0054 end -0055 -0056 import java.io.File; -0057 import java.io.FileOutputStream; -0058 import java.io.IOException; -0059 -0060 %Remove the output file if it already exists -0061 if isfile(fileName) -0062 delete(fileName); -0063 end -0064 -0065 %Load an empty workbook -0066 wb=loadWorkbook(fileName,true); -0067 -0068 %Construct equations -0069 model.equations=constructEquations(model,model.rxns,true); -0070 -0071 %Check if it should print genes -0072 if isfield(model,'grRules') -0073 rules=model.grRules; -0074 else -0075 rules=[]; -0076 end -0077 -0078 %Check if the model has default upper/lower bounds. This determines if -0079 %those values should be printed or not -0080 hasDefaultLB=false; -0081 hasDefaultUB=false; -0082 if isfield(model,'annotation') -0083 if isfield(model.annotation,'defaultLB') -0084 hasDefaultLB=true; -0085 end -0086 if isfield(model.annotation,'defaultUB') -0087 hasDefaultUB=true; -0088 end -0089 end -0090 -0091 %Add the RXNS sheet -0092 -0093 %Create the header row -0094 headers={'#';'ID';'NAME';'EQUATION';'EC-NUMBER';'GENE ASSOCIATION';'LOWER BOUND';'UPPER BOUND';'OBJECTIVE';'COMPARTMENT';'MIRIAM';'SUBSYSTEM';'REPLACEMENT ID';'NOTE';'REFERENCE';'CONFIDENCE SCORE'}; -0095 -0096 %Add empty comments -0097 emptyColumn=cell(numel(model.rxns),1); -0098 rxnSheet=emptyColumn; -0099 -0100 %Add the model fields -0101 rxnSheet=[rxnSheet model.rxns]; -0102 -0103 if isfield(model,'rxnNames') -0104 rxnSheet=[rxnSheet model.rxnNames]; -0105 else -0106 rxnSheet=[rxnSheet emptyColumn]; -0107 end -0108 -0109 rxnSheet=[rxnSheet model.equations]; -0110 -0111 if isfield(model,'eccodes') -0112 rxnSheet=[rxnSheet model.eccodes]; -0113 else -0114 rxnSheet=[rxnSheet emptyColumn]; -0115 end -0116 -0117 if ~isempty(rules) -0118 rxnSheet=[rxnSheet rules]; -0119 else -0120 rxnSheet=[rxnSheet emptyColumn]; -0121 end -0122 -0123 lb=emptyColumn; -0124 ub=emptyColumn; -0125 objective=emptyColumn; -0126 rxnMiriams=emptyColumn; -0127 -0128 for i=1:numel(model.rxns) -0129 if isfield(model,'lb') -0130 if hasDefaultLB==true -0131 if model.rev(i)==1 -0132 %If reversible, print only if different than defaultLB -0133 if model.lb(i) ~= model.annotation.defaultLB -0134 lb{i}=model.lb(i); -0135 end -0136 else -0137 %If irreversible, print only for non-zero values -0138 if model.lb(i)~=0 -0139 lb{i}=model.lb(i); -0140 end -0141 end -0142 else -0143 lb{i}=model.lb(i); -0144 end -0145 end -0146 -0147 if isfield(model,'ub') -0148 if hasDefaultUB==true -0149 if model.ub(i) ~= model.annotation.defaultUB -0150 ub{i}=model.ub(i); -0151 end -0152 else -0153 ub{i}=model.ub(i); -0154 end -0155 end -0156 -0157 if isfield(model,'c') -0158 if model.c(i)~=0 -0159 objective{i}=model.c(i); -0160 end -0161 end -0162 -0163 if isfield(model,'rxnMiriams') -0164 if ~isempty(model.rxnMiriams{i}) -0165 toPrint=[]; -0166 for j=1:numel(model.rxnMiriams{i}.name) -0167 toPrint=[toPrint strtrim(model.rxnMiriams{i}.name{j}) '/' strtrim(model.rxnMiriams{i}.value{j}) ';']; -0168 end -0169 rxnMiriams{i}=toPrint(1:end-1); -0170 end -0171 end -0172 end -0173 -0174 rxnSheet=[rxnSheet lb]; -0175 rxnSheet=[rxnSheet ub]; -0176 rxnSheet=[rxnSheet objective]; -0177 -0178 if isfield(model,'rxnComps') -0179 rxnSheet=[rxnSheet model.comps(model.rxnComps)]; -0180 else -0181 rxnSheet=[rxnSheet emptyColumn]; -0182 end -0183 -0184 rxnSheet=[rxnSheet rxnMiriams]; -0185 -0186 subsystems=''; -0187 if isfield(model,'subSystems') -0188 for i=1:numel(model.subSystems) -0189 if ~isempty(model.subSystems{i,1}) -0190 subsystems{i,1}=strjoin(model.subSystems{i,1},';'); -0191 else -0192 subsystems{i,1}=''; -0193 end -0194 end -0195 rxnSheet=[rxnSheet subsystems]; -0196 else -0197 rxnSheet=[rxnSheet emptyColumn]; -0198 end -0199 -0200 %For REPLACEMENT ID which isn't in the model -0201 rxnSheet=[rxnSheet emptyColumn]; -0202 -0203 if isfield(model,'rxnNotes') -0204 rxnSheet=[rxnSheet model.rxnNotes]; -0205 else -0206 rxnSheet=[rxnSheet emptyColumn]; -0207 end -0208 -0209 if isfield(model,'rxnReferences') -0210 rxnSheet=[rxnSheet model.rxnReferences]; -0211 else -0212 rxnSheet=[rxnSheet emptyColumn]; -0213 end -0214 -0215 if isfield(model,'rxnConfidenceScores') -0216 rxnSheet=[rxnSheet num2cell(model.rxnConfidenceScores)]; -0217 else -0218 rxnSheet=[rxnSheet emptyColumn]; -0219 end -0220 -0221 wb=writeSheet(wb,'RXNS',0,headers,[],rxnSheet); -0222 -0223 headers={'#';'ID';'NAME';'UNCONSTRAINED';'MIRIAM';'COMPOSITION';'InChI';'COMPARTMENT';'REPLACEMENT ID';'CHARGE'}; -0224 -0225 metSheet=cell(numel(model.mets),numel(headers)); -0226 -0227 for i=1:numel(model.mets) -0228 metSheet{i,2}=[model.metNames{i} '[' model.comps{model.metComps(i)} ']']; -0229 -0230 if isfield(model,'metNames') -0231 metSheet(i,3)=model.metNames(i); -0232 end -0233 -0234 if isfield(model,'unconstrained') -0235 if model.unconstrained(i)~=0 -0236 metSheet{i,4}=true; -0237 end -0238 end -0239 -0240 if isfield(model,'metMiriams') -0241 if ~isempty(model.metMiriams{i}) -0242 toPrint=[]; -0243 for j=1:numel(model.metMiriams{i}.name) -0244 toPrint=[toPrint strtrim(model.metMiriams{i}.name{j}) '/' strtrim(model.metMiriams{i}.value{j}) ';']; -0245 end -0246 metSheet{i,5}=toPrint(1:end-1); -0247 end -0248 end -0249 -0250 % Making sure that only these metFormulas are exported, which don't -0251 % have InChI strings -0252 if isfield(model,'metFormulas') -0253 if isfield(model,'inchis') -0254 if isempty(model.inchis{i}) -0255 metSheet(i,6)=model.metFormulas(i); -0256 end -0257 else -0258 metSheet(i,6)=model.metFormulas(i); -0259 end -0260 end -0261 -0262 if isfield(model,'inchis') -0263 metSheet(i,7)=model.inchis(i); -0264 end -0265 -0266 if isfield(model,'metComps') -0267 metSheet(i,8)=model.comps(model.metComps(i)); -0268 end -0269 -0270 metSheet(i,9)=model.mets(i); -0271 -0272 if isfield(model,'metCharges') -0273 metSheet{i,10}=model.metCharges(i); -0274 end -0275 end -0276 -0277 wb=writeSheet(wb,'METS',1,headers,[],metSheet); -0278 -0279 %Add the COMPS sheet -0280 -0281 %Create the header row -0282 headers={'#';'ABBREVIATION';'NAME';'INSIDE';'MIRIAM'}; -0283 -0284 compSheet=cell(numel(model.comps),numel(headers)); -0285 -0286 for i=1:numel(model.comps) -0287 compSheet(i,2)=model.comps(i); -0288 -0289 if isfield(model,'compNames') -0290 compSheet(i,3)=model.compNames(i); -0291 end -0292 -0293 if isfield(model,'compOutside') -0294 compSheet(i,4)=model.compOutside(i); -0295 end -0296 -0297 if isfield(model,'compMiriams') -0298 if ~isempty(model.compMiriams{i}) -0299 toPrint=[]; -0300 for j=1:numel(model.compMiriams{i}.name) -0301 toPrint=[toPrint strtrim(model.compMiriams{i}.name{j}) '/' strtrim(model.compMiriams{i}.value{j}) ';']; -0302 end -0303 compSheet{i,5}=toPrint(1:end-1); -0304 end -0305 end -0306 end -0307 -0308 wb=writeSheet(wb,'COMPS',2,headers,[],compSheet); -0309 -0310 %Add the GENES sheet -0311 if isfield(model,'genes') -0312 %Create the header row -0313 headers={'#';'NAME';'MIRIAM';'SHORT NAME';'COMPARTMENT'}; -0314 -0315 geneSheet=cell(numel(model.genes),numel(headers)); -0316 -0317 for i=1:numel(model.genes) -0318 geneSheet(i,2)=model.genes(i); -0319 -0320 if isfield(model,'geneMiriams') -0321 if ~isempty(model.geneMiriams{i}) -0322 toPrint=[]; -0323 for j=1:numel(model.geneMiriams{i}.name) -0324 toPrint=[toPrint strtrim(model.geneMiriams{i}.name{j}) '/' strtrim(model.geneMiriams{i}.value{j}) ';']; -0325 end -0326 geneSheet{i,3}=toPrint(1:end-1); -0327 end -0328 end -0329 if isfield(model,'geneShortNames') -0330 geneSheet(i,4)=model.geneShortNames(i); -0331 end -0332 if isfield(model,'geneComps') -0333 geneSheet(i,5)=model.comps(model.geneComps(i)); -0334 end -0335 end -0336 -0337 wb=writeSheet(wb,'GENES',3,headers,[],geneSheet); -0338 end -0339 -0340 %Add the MODEL sheet -0341 -0342 %Create the header row -0343 headers={'#';'ID';'NAME';'TAXONOMY';'DEFAULT LOWER';'DEFAULT UPPER';'CONTACT GIVEN NAME';'CONTACT FAMILY NAME';'CONTACT EMAIL';'ORGANIZATION';'NOTES'}; -0344 -0345 modelSheet=cell(1,numel(headers)); -0346 -0347 if ~isfield(model,'annotation') -0348 model.annotation = []; -0349 end -0350 -0351 if isfield(model,'id') -0352 modelSheet{1,2}=model.id; -0353 else -0354 modelSheet{1,2}='blankID'; -0355 end -0356 if isfield(model,'name') -0357 modelSheet{1,3}=model.name; -0358 else -0359 modelSheet{1,3}='blankName'; -0360 end -0361 if isfield(model.annotation,'taxonomy') -0362 modelSheet{1,4}=model.annotation.taxonomy; -0363 end -0364 if isfield(model.annotation,'defaultLB') -0365 modelSheet{1,5}=model.annotation.defaultLB; -0366 end -0367 if isfield(model.annotation,'defaultUB') -0368 modelSheet{1,6}=model.annotation.defaultUB; -0369 end -0370 if isfield(model.annotation,'givenName') -0371 modelSheet{1,7}=model.annotation.givenName; -0372 end -0373 if isfield(model.annotation,'familyName') -0374 modelSheet{1,8}=model.annotation.familyName; -0375 end -0376 if isfield(model.annotation,'email') -0377 modelSheet{1,9}=model.annotation.email; -0378 end -0379 if isfield(model.annotation,'organization') -0380 modelSheet{1,10}=model.annotation.organization; -0381 end -0382 if isfield(model.annotation,'note') -0383 modelSheet{1,11}=model.annotation.note; -0384 end -0385 -0386 if isfield(model,'genes') -0387 wb=writeSheet(wb,'MODEL',4,headers,[],modelSheet); -0388 else -0389 wb=writeSheet(wb,'MODEL',3,headers,[],modelSheet); -0390 end -0391 -0392 %Open the output stream -0393 out = FileOutputStream(fileName); -0394 wb.write(out); -0395 out.close(); -0396 end +0033 +0034 checkModelStruct(model); +0035 +0036 addList = matlab.addons.installedAddons; +0037 if any(strcmpi(addList.Name,'Text Analytics Toolbox')) +0038 error(['exportToExcelFormat is incompatible with MATLAB Text Analytics Toolbox. ' ... +0039 'Further instructions => https://github.com/SysBioChalmers/RAVEN/issues/55#issuecomment-1514369299']) +0040 end +0041 +0042 [~, A, B]=fileparts(fileName); +0043 +0044 %If a path was used call on exportToTabDelimited instead +0045 if ~any(A) || ~any(B) +0046 exportToTabDelimited(model,fileName); +0047 return; +0048 end +0049 +0050 if ~strcmpi(B,'.xlsx') +0051 EM='As of RAVEN version 1.9, only export to xlsx format is supported'; +0052 dispEM(EM); +0053 end +0054 +0055 import java.io.File; +0056 import java.io.FileOutputStream; +0057 import java.io.IOException; +0058 +0059 %Remove the output file if it already exists +0060 if isfile(fileName) +0061 delete(fileName); +0062 end +0063 +0064 %Load an empty workbook +0065 wb=loadWorkbook(fileName,true); +0066 +0067 %Construct equations +0068 model.equations=constructEquations(model,model.rxns,true); +0069 +0070 %Check if it should print genes +0071 if isfield(model,'grRules') +0072 rules=model.grRules; +0073 else +0074 rules=[]; +0075 end +0076 +0077 %Check if the model has default upper/lower bounds. This determines if +0078 %those values should be printed or not +0079 hasDefaultLB=false; +0080 hasDefaultUB=false; +0081 if isfield(model,'annotation') +0082 if isfield(model.annotation,'defaultLB') +0083 hasDefaultLB=true; +0084 end +0085 if isfield(model.annotation,'defaultUB') +0086 hasDefaultUB=true; +0087 end +0088 end +0089 +0090 %Add the RXNS sheet +0091 +0092 %Create the header row +0093 headers={'#';'ID';'NAME';'EQUATION';'EC-NUMBER';'GENE ASSOCIATION';'LOWER BOUND';'UPPER BOUND';'OBJECTIVE';'COMPARTMENT';'MIRIAM';'SUBSYSTEM';'REPLACEMENT ID';'NOTE';'REFERENCE';'CONFIDENCE SCORE'}; +0094 +0095 %Add empty comments +0096 emptyColumn=cell(numel(model.rxns),1); +0097 rxnSheet=emptyColumn; +0098 +0099 %Add the model fields +0100 rxnSheet=[rxnSheet model.rxns]; +0101 +0102 if isfield(model,'rxnNames') +0103 rxnSheet=[rxnSheet model.rxnNames]; +0104 else +0105 rxnSheet=[rxnSheet emptyColumn]; +0106 end +0107 +0108 rxnSheet=[rxnSheet model.equations]; +0109 +0110 if isfield(model,'eccodes') +0111 rxnSheet=[rxnSheet model.eccodes]; +0112 else +0113 rxnSheet=[rxnSheet emptyColumn]; +0114 end +0115 +0116 if ~isempty(rules) +0117 rxnSheet=[rxnSheet rules]; +0118 else +0119 rxnSheet=[rxnSheet emptyColumn]; +0120 end +0121 +0122 lb=emptyColumn; +0123 ub=emptyColumn; +0124 objective=emptyColumn; +0125 rxnMiriams=emptyColumn; +0126 +0127 for i=1:numel(model.rxns) +0128 if isfield(model,'lb') +0129 if hasDefaultLB==true +0130 if model.rev(i)==1 +0131 %If reversible, print only if different than defaultLB +0132 if model.lb(i) ~= model.annotation.defaultLB +0133 lb{i}=model.lb(i); +0134 end +0135 else +0136 %If irreversible, print only for non-zero values +0137 if model.lb(i)~=0 +0138 lb{i}=model.lb(i); +0139 end +0140 end +0141 else +0142 lb{i}=model.lb(i); +0143 end +0144 end +0145 +0146 if isfield(model,'ub') +0147 if hasDefaultUB==true +0148 if model.ub(i) ~= model.annotation.defaultUB +0149 ub{i}=model.ub(i); +0150 end +0151 else +0152 ub{i}=model.ub(i); +0153 end +0154 end +0155 +0156 if isfield(model,'c') +0157 if model.c(i)~=0 +0158 objective{i}=model.c(i); +0159 end +0160 end +0161 +0162 if isfield(model,'rxnMiriams') +0163 if ~isempty(model.rxnMiriams{i}) +0164 toPrint=[]; +0165 for j=1:numel(model.rxnMiriams{i}.name) +0166 toPrint=[toPrint strtrim(model.rxnMiriams{i}.name{j}) '/' strtrim(model.rxnMiriams{i}.value{j}) ';']; +0167 end +0168 rxnMiriams{i}=toPrint(1:end-1); +0169 end +0170 end +0171 end +0172 +0173 rxnSheet=[rxnSheet lb]; +0174 rxnSheet=[rxnSheet ub]; +0175 rxnSheet=[rxnSheet objective]; +0176 +0177 if isfield(model,'rxnComps') +0178 rxnSheet=[rxnSheet model.comps(model.rxnComps)]; +0179 else +0180 rxnSheet=[rxnSheet emptyColumn]; +0181 end +0182 +0183 rxnSheet=[rxnSheet rxnMiriams]; +0184 +0185 subsystems=''; +0186 if isfield(model,'subSystems') +0187 for i=1:numel(model.subSystems) +0188 if ~isempty(model.subSystems{i,1}) +0189 subsystems{i,1}=strjoin(model.subSystems{i,1},';'); +0190 else +0191 subsystems{i,1}=''; +0192 end +0193 end +0194 rxnSheet=[rxnSheet subsystems]; +0195 else +0196 rxnSheet=[rxnSheet emptyColumn]; +0197 end +0198 +0199 %For REPLACEMENT ID which isn't in the model +0200 rxnSheet=[rxnSheet emptyColumn]; +0201 +0202 if isfield(model,'rxnNotes') +0203 rxnSheet=[rxnSheet model.rxnNotes]; +0204 else +0205 rxnSheet=[rxnSheet emptyColumn]; +0206 end +0207 +0208 if isfield(model,'rxnReferences') +0209 rxnSheet=[rxnSheet model.rxnReferences]; +0210 else +0211 rxnSheet=[rxnSheet emptyColumn]; +0212 end +0213 +0214 if isfield(model,'rxnConfidenceScores') +0215 rxnSheet=[rxnSheet num2cell(model.rxnConfidenceScores)]; +0216 else +0217 rxnSheet=[rxnSheet emptyColumn]; +0218 end +0219 +0220 wb=writeSheet(wb,'RXNS',0,headers,[],rxnSheet); +0221 +0222 headers={'#';'ID';'NAME';'UNCONSTRAINED';'MIRIAM';'COMPOSITION';'InChI';'COMPARTMENT';'REPLACEMENT ID';'CHARGE'}; +0223 +0224 metSheet=cell(numel(model.mets),numel(headers)); +0225 +0226 for i=1:numel(model.mets) +0227 metSheet{i,2}=[model.metNames{i} '[' model.comps{model.metComps(i)} ']']; +0228 +0229 if isfield(model,'metNames') +0230 metSheet(i,3)=model.metNames(i); +0231 end +0232 +0233 if isfield(model,'unconstrained') +0234 if model.unconstrained(i)~=0 +0235 metSheet{i,4}=true; +0236 end +0237 end +0238 +0239 if isfield(model,'metMiriams') +0240 if ~isempty(model.metMiriams{i}) +0241 toPrint=[]; +0242 for j=1:numel(model.metMiriams{i}.name) +0243 toPrint=[toPrint strtrim(model.metMiriams{i}.name{j}) '/' strtrim(model.metMiriams{i}.value{j}) ';']; +0244 end +0245 metSheet{i,5}=toPrint(1:end-1); +0246 end +0247 end +0248 +0249 % Making sure that only these metFormulas are exported, which don't +0250 % have InChI strings +0251 if isfield(model,'metFormulas') +0252 if isfield(model,'inchis') +0253 if isempty(model.inchis{i}) +0254 metSheet(i,6)=model.metFormulas(i); +0255 end +0256 else +0257 metSheet(i,6)=model.metFormulas(i); +0258 end +0259 end +0260 +0261 if isfield(model,'inchis') +0262 metSheet(i,7)=model.inchis(i); +0263 end +0264 +0265 if isfield(model,'metComps') +0266 metSheet(i,8)=model.comps(model.metComps(i)); +0267 end +0268 +0269 metSheet(i,9)=model.mets(i); +0270 +0271 if isfield(model,'metCharges') +0272 metSheet{i,10}=model.metCharges(i); +0273 end +0274 end +0275 +0276 wb=writeSheet(wb,'METS',1,headers,[],metSheet); +0277 +0278 %Add the COMPS sheet +0279 +0280 %Create the header row +0281 headers={'#';'ABBREVIATION';'NAME';'INSIDE';'MIRIAM'}; +0282 +0283 compSheet=cell(numel(model.comps),numel(headers)); +0284 +0285 for i=1:numel(model.comps) +0286 compSheet(i,2)=model.comps(i); +0287 +0288 if isfield(model,'compNames') +0289 compSheet(i,3)=model.compNames(i); +0290 end +0291 +0292 if isfield(model,'compOutside') +0293 compSheet(i,4)=model.compOutside(i); +0294 end +0295 +0296 if isfield(model,'compMiriams') +0297 if ~isempty(model.compMiriams{i}) +0298 toPrint=[]; +0299 for j=1:numel(model.compMiriams{i}.name) +0300 toPrint=[toPrint strtrim(model.compMiriams{i}.name{j}) '/' strtrim(model.compMiriams{i}.value{j}) ';']; +0301 end +0302 compSheet{i,5}=toPrint(1:end-1); +0303 end +0304 end +0305 end +0306 +0307 wb=writeSheet(wb,'COMPS',2,headers,[],compSheet); +0308 +0309 %Add the GENES sheet +0310 if isfield(model,'genes') +0311 %Create the header row +0312 headers={'#';'NAME';'MIRIAM';'SHORT NAME';'COMPARTMENT'}; +0313 +0314 geneSheet=cell(numel(model.genes),numel(headers)); +0315 +0316 for i=1:numel(model.genes) +0317 geneSheet(i,2)=model.genes(i); +0318 +0319 if isfield(model,'geneMiriams') +0320 if ~isempty(model.geneMiriams{i}) +0321 toPrint=[]; +0322 for j=1:numel(model.geneMiriams{i}.name) +0323 toPrint=[toPrint strtrim(model.geneMiriams{i}.name{j}) '/' strtrim(model.geneMiriams{i}.value{j}) ';']; +0324 end +0325 geneSheet{i,3}=toPrint(1:end-1); +0326 end +0327 end +0328 if isfield(model,'geneShortNames') +0329 geneSheet(i,4)=model.geneShortNames(i); +0330 end +0331 if isfield(model,'geneComps') +0332 geneSheet(i,5)=model.comps(model.geneComps(i)); +0333 end +0334 end +0335 +0336 wb=writeSheet(wb,'GENES',3,headers,[],geneSheet); +0337 end +0338 +0339 %Add the MODEL sheet +0340 +0341 %Create the header row +0342 headers={'#';'ID';'NAME';'TAXONOMY';'DEFAULT LOWER';'DEFAULT UPPER';'CONTACT GIVEN NAME';'CONTACT FAMILY NAME';'CONTACT EMAIL';'ORGANIZATION';'NOTES'}; +0343 +0344 modelSheet=cell(1,numel(headers)); +0345 +0346 if ~isfield(model,'annotation') +0347 model.annotation = []; +0348 end +0349 +0350 if isfield(model,'id') +0351 modelSheet{1,2}=model.id; +0352 else +0353 modelSheet{1,2}='blankID'; +0354 end +0355 if isfield(model,'name') +0356 modelSheet{1,3}=model.name; +0357 else +0358 modelSheet{1,3}='blankName'; +0359 end +0360 if isfield(model.annotation,'taxonomy') +0361 modelSheet{1,4}=model.annotation.taxonomy; +0362 end +0363 if isfield(model.annotation,'defaultLB') +0364 modelSheet{1,5}=model.annotation.defaultLB; +0365 end +0366 if isfield(model.annotation,'defaultUB') +0367 modelSheet{1,6}=model.annotation.defaultUB; +0368 end +0369 if isfield(model.annotation,'givenName') +0370 modelSheet{1,7}=model.annotation.givenName; +0371 end +0372 if isfield(model.annotation,'familyName') +0373 modelSheet{1,8}=model.annotation.familyName; +0374 end +0375 if isfield(model.annotation,'email') +0376 modelSheet{1,9}=model.annotation.email; +0377 end +0378 if isfield(model.annotation,'organization') +0379 modelSheet{1,10}=model.annotation.organization; +0380 end +0381 if isfield(model.annotation,'note') +0382 modelSheet{1,11}=model.annotation.note; +0383 end +0384 +0385 if isfield(model,'genes') +0386 wb=writeSheet(wb,'MODEL',4,headers,[],modelSheet); +0387 else +0388 wb=writeSheet(wb,'MODEL',3,headers,[],modelSheet); +0389 end +0390 +0391 %Open the output stream +0392 out = FileOutputStream(fileName); +0393 wb.write(out); +0394 out.close(); +0395 end
Generated by m2html © 2005
\ No newline at end of file diff --git a/doc/io/readYAMLmodel.html b/doc/io/readYAMLmodel.html index 3a96f8c3..e8c684be 100644 --- a/doc/io/readYAMLmodel.html +++ b/doc/io/readYAMLmodel.html @@ -102,7 +102,7 @@

SOURCE CODE ^% leading spaces to avoid metaData to be concatenated. 0046 newLine=regexp(line_raw,'^ {6,}([\w\(\)].*)','tokenExtents'); 0047 brokenLine=find(~cellfun('isempty',newLine)); -0048 for i=1:numel(brokenLine) +0048 for i=flip(1:numel(brokenLine)) 0049 extraLine = char(line_raw(brokenLine(i))); 0050 extraLine = extraLine(newLine{brokenLine(i)}{1}(1):end); 0051 line_raw{brokenLine(i)-1} = strjoin({line_raw{brokenLine(i)-1},extraLine},' '); @@ -115,7 +115,7 @@

SOURCE CODE ^'.*:$',''); 0059 line_value = regexprep(line_value, '[^":]+: "?(.+)"?$','$1'); 0060 line_value = regexprep(line_value, '(")|(^ {4,}- )',''); -0061 +0061 line_value(strcmp(line_value,'''''')) = {''}; 0062 line_value(strcmp(line_value,line_raw)) = {''}; 0063 0064 modelFields = {'id',char();... @@ -197,7 +197,7 @@

SOURCE CODE ^for i=1:numel(line_key) 0146 tline_raw = line_raw{i}; diff --git a/doc/io/writeYAMLmodel.html b/doc/io/writeYAMLmodel.html index d2f385a9..8894eeb0 100644 --- a/doc/io/writeYAMLmodel.html +++ b/doc/io/writeYAMLmodel.html @@ -97,347 +97,350 @@

SOURCE CODE ^end 0040 -0041 %Sort identifiers alphabetically -0042 if sortIds == true -0043 model = sortIdentifiers(model); -0044 end -0045 -0046 %Simplify Miriam fields: -0047 if isfield(model,'metMiriams') -0048 [model.newMetMiriams,model.newMetMiriamNames] = extractMiriam(model.metMiriams); -0049 end -0050 if isfield(model,'rxnMiriams') -0051 [model.newRxnMiriams,model.newRxnMiriamNames] = extractMiriam(model.rxnMiriams); +0041 %Check that the model structure has no problems +0042 checkModelStruct(model); +0043 +0044 %Sort identifiers alphabetically +0045 if sortIds == true +0046 model = sortIdentifiers(model); +0047 end +0048 +0049 %Simplify Miriam fields: +0050 if isfield(model,'metMiriams') +0051 [model.newMetMiriams,model.newMetMiriamNames] = extractMiriam(model.metMiriams); 0052 end -0053 if isfield(model,'geneMiriams') -0054 [model.newGeneMiriams,model.newGeneMiriamNames] = extractMiriam(model.geneMiriams); +0053 if isfield(model,'rxnMiriams') +0054 [model.newRxnMiriams,model.newRxnMiriamNames] = extractMiriam(model.rxnMiriams); 0055 end -0056 if isfield(model,'compMiriams') -0057 [model.newCompMiriams,model.newCompMiriamNames] = extractMiriam(model.compMiriams); +0056 if isfield(model,'geneMiriams') +0057 [model.newGeneMiriams,model.newGeneMiriamNames] = extractMiriam(model.geneMiriams); 0058 end -0059 -0060 %Open file: -0061 fid = fopen(fileName,'wt'); -0062 if fid == -1 -0063 error(['Cannot write to ' fileName ', does the directory exist?']) -0064 end -0065 fprintf(fid,'---\n!!omap\n'); -0066 -0067 %Insert file header (metadata) -0068 writeMetadata(model,fid); +0059 if isfield(model,'compMiriams') +0060 [model.newCompMiriams,model.newCompMiriamNames] = extractMiriam(model.compMiriams); +0061 end +0062 +0063 %Open file: +0064 fid = fopen(fileName,'wt'); +0065 if fid == -1 +0066 error(['Cannot write to ' fileName ', does the directory exist?']) +0067 end +0068 fprintf(fid,'---\n!!omap\n'); 0069 -0070 %Metabolites: -0071 fprintf(fid,'- metabolites:\n'); -0072 for i = 1:length(model.mets) -0073 fprintf(fid,' - !!omap\n'); -0074 writeField(model, fid, 'mets', 'txt', i, ' - id', preserveQuotes) -0075 writeField(model, fid, 'metNames', 'txt', i, ' - name', preserveQuotes) -0076 writeField(model, fid, 'metComps', 'txt', i, ' - compartment', preserveQuotes) -0077 writeField(model, fid, 'metFormulas', 'txt', i, ' - formula', preserveQuotes) -0078 writeField(model, fid, 'metCharges', 'num', i, ' - charge', preserveQuotes) -0079 writeField(model, fid, 'inchis', 'txt', i, ' - inchis', preserveQuotes) -0080 writeField(model, fid, 'metSmiles', 'txt', i, ' - smiles', preserveQuotes) -0081 writeField(model, fid, 'metMiriams', 'txt', i, ' - annotation', preserveQuotes) -0082 writeField(model, fid, 'metDeltaG', 'num', i, ' - deltaG', preserveQuotes) -0083 writeField(model, fid, 'metNotes', 'txt', i, ' - notes', preserveQuotes) -0084 writeField(model, fid, 'metFrom', 'txt', i, ' - metFrom', preserveQuotes) -0085 end -0086 -0087 %Reactions: -0088 fprintf(fid,'- reactions:\n'); -0089 for i = 1:length(model.rxns) -0090 fprintf(fid,' - !!omap\n'); -0091 writeField(model, fid, 'rxns', 'txt', i, ' - id', preserveQuotes) -0092 writeField(model, fid, 'rxnNames', 'txt', i, ' - name', preserveQuotes) -0093 writeField(model, fid, 'S', 'txt', i, ' - metabolites', preserveQuotes) -0094 writeField(model, fid, 'lb', 'num', i, ' - lower_bound', preserveQuotes) -0095 writeField(model, fid, 'ub', 'num', i, ' - upper_bound', preserveQuotes) -0096 writeField(model, fid, 'grRules', 'txt', i, ' - gene_reaction_rule', preserveQuotes) -0097 writeField(model, fid, 'rxnFrom', 'txt', i, ' - rxnFrom', preserveQuotes) -0098 if model.c(i)~=0 -0099 writeField(model, fid, 'c', 'num', i, ' - objective_coefficient', preserveQuotes) -0100 end -0101 writeField(model, fid, 'eccodes', 'txt', i, ' - eccodes', preserveQuotes) -0102 writeField(model, fid, 'rxnReferences', 'txt', i, ' - references', preserveQuotes) -0103 writeField(model, fid, 'subSystems', 'txt', i, ' - subsystem', preserveQuotes) -0104 writeField(model, fid, 'rxnMiriams', 'txt', i, ' - annotation', preserveQuotes) -0105 writeField(model, fid, 'rxnDeltaG', 'num', i, ' - deltaG', preserveQuotes) -0106 writeField(model, fid, 'rxnConfidenceScores', 'num', i, ' - confidence_score', preserveQuotes) -0107 writeField(model, fid, 'rxnNotes', 'txt', i, ' - rxnNotes', preserveQuotes) -0108 end -0109 -0110 %Genes: -0111 fprintf(fid,'- genes:\n'); -0112 for i = 1:length(model.genes) -0113 fprintf(fid,' - !!omap\n'); -0114 writeField(model, fid, 'genes', 'txt', i, ' - id', preserveQuotes) -0115 writeField(model, fid, 'geneShortNames', 'txt', i, ' - name', preserveQuotes) -0116 writeField(model, fid, 'proteins', 'txt', i, ' - protein', preserveQuotes) -0117 writeField(model, fid, 'geneMiriams', 'txt', i, ' - annotation', preserveQuotes) -0118 end -0119 -0120 %Compartments: -0121 fprintf(fid,'- compartments: !!omap\n'); -0122 for i = 1:length(model.comps) -0123 writeField(model, fid, 'compNames', 'txt', i, ['- ' model.comps{i}], preserveQuotes) -0124 writeField(model, fid, 'compMiriams', 'txt', i, '- annotation', preserveQuotes) -0125 end -0126 -0127 -0128 %EC-model: -0129 if isfield(model,'ec') -0130 fprintf(fid,'- ec-rxns:\n'); -0131 for i = 1:length(model.ec.rxns) -0132 fprintf(fid,' - !!omap\n'); -0133 writeField(model.ec, fid, 'rxns', 'txt', i, '- id', preserveQuotes) -0134 writeField(model.ec, fid, 'kcat', 'num', i, '- kcat', preserveQuotes) -0135 writeField(model.ec, fid, 'source', 'txt', i, '- source', preserveQuotes) -0136 writeField(model.ec, fid, 'notes', 'txt', i, '- notes', preserveQuotes) -0137 writeField(model.ec, fid, 'eccodes', 'txt', i, '- eccodes', preserveQuotes) -0138 writeField(model.ec, fid, 'rxnEnzMat', 'txt', i, '- enzymes', preserveQuotes) -0139 end -0140 -0141 fprintf(fid,'- ec-enzymes:\n'); -0142 for i = 1:length(model.ec.genes) -0143 fprintf(fid,' - !!omap\n'); -0144 writeField(model.ec, fid, 'genes', 'txt', i, '- genes', preserveQuotes) -0145 writeField(model.ec, fid, 'enzymes', 'txt', i, '- enzymes', preserveQuotes) -0146 writeField(model.ec, fid, 'mw', 'num', i, '- mw', preserveQuotes) -0147 writeField(model.ec, fid, 'sequence', 'txt', i, '- sequence', preserveQuotes) -0148 writeField(model.ec, fid, 'concs', 'num', i, '- concs', preserveQuotes) -0149 end -0150 end -0151 -0152 %Close file: -0153 fclose(fid); +0070 %Insert file header (metadata) +0071 writeMetadata(model,fid); +0072 +0073 %Metabolites: +0074 fprintf(fid,'- metabolites:\n'); +0075 for i = 1:length(model.mets) +0076 fprintf(fid,' - !!omap\n'); +0077 writeField(model, fid, 'mets', 'txt', i, ' - id', preserveQuotes) +0078 writeField(model, fid, 'metNames', 'txt', i, ' - name', preserveQuotes) +0079 writeField(model, fid, 'metComps', 'txt', i, ' - compartment', preserveQuotes) +0080 writeField(model, fid, 'metFormulas', 'txt', i, ' - formula', preserveQuotes) +0081 writeField(model, fid, 'metCharges', 'num', i, ' - charge', preserveQuotes) +0082 writeField(model, fid, 'inchis', 'txt', i, ' - inchis', preserveQuotes) +0083 writeField(model, fid, 'metSmiles', 'txt', i, ' - smiles', preserveQuotes) +0084 writeField(model, fid, 'metMiriams', 'txt', i, ' - annotation', preserveQuotes) +0085 writeField(model, fid, 'metDeltaG', 'num', i, ' - deltaG', preserveQuotes) +0086 writeField(model, fid, 'metNotes', 'txt', i, ' - notes', preserveQuotes) +0087 writeField(model, fid, 'metFrom', 'txt', i, ' - metFrom', preserveQuotes) +0088 end +0089 +0090 %Reactions: +0091 fprintf(fid,'- reactions:\n'); +0092 for i = 1:length(model.rxns) +0093 fprintf(fid,' - !!omap\n'); +0094 writeField(model, fid, 'rxns', 'txt', i, ' - id', preserveQuotes) +0095 writeField(model, fid, 'rxnNames', 'txt', i, ' - name', preserveQuotes) +0096 writeField(model, fid, 'S', 'txt', i, ' - metabolites', preserveQuotes) +0097 writeField(model, fid, 'lb', 'num', i, ' - lower_bound', preserveQuotes) +0098 writeField(model, fid, 'ub', 'num', i, ' - upper_bound', preserveQuotes) +0099 writeField(model, fid, 'grRules', 'txt', i, ' - gene_reaction_rule', preserveQuotes) +0100 writeField(model, fid, 'rxnFrom', 'txt', i, ' - rxnFrom', preserveQuotes) +0101 if model.c(i)~=0 +0102 writeField(model, fid, 'c', 'num', i, ' - objective_coefficient', preserveQuotes) +0103 end +0104 writeField(model, fid, 'eccodes', 'txt', i, ' - eccodes', preserveQuotes) +0105 writeField(model, fid, 'rxnReferences', 'txt', i, ' - references', preserveQuotes) +0106 writeField(model, fid, 'subSystems', 'txt', i, ' - subsystem', preserveQuotes) +0107 writeField(model, fid, 'rxnMiriams', 'txt', i, ' - annotation', preserveQuotes) +0108 writeField(model, fid, 'rxnDeltaG', 'num', i, ' - deltaG', preserveQuotes) +0109 writeField(model, fid, 'rxnConfidenceScores', 'num', i, ' - confidence_score', preserveQuotes) +0110 writeField(model, fid, 'rxnNotes', 'txt', i, ' - rxnNotes', preserveQuotes) +0111 end +0112 +0113 %Genes: +0114 fprintf(fid,'- genes:\n'); +0115 for i = 1:length(model.genes) +0116 fprintf(fid,' - !!omap\n'); +0117 writeField(model, fid, 'genes', 'txt', i, ' - id', preserveQuotes) +0118 writeField(model, fid, 'geneShortNames', 'txt', i, ' - name', preserveQuotes) +0119 writeField(model, fid, 'proteins', 'txt', i, ' - protein', preserveQuotes) +0120 writeField(model, fid, 'geneMiriams', 'txt', i, ' - annotation', preserveQuotes) +0121 end +0122 +0123 %Compartments: +0124 fprintf(fid,'- compartments: !!omap\n'); +0125 for i = 1:length(model.comps) +0126 writeField(model, fid, 'compNames', 'txt', i, ['- ' model.comps{i}], preserveQuotes) +0127 writeField(model, fid, 'compMiriams', 'txt', i, '- annotation', preserveQuotes) +0128 end +0129 +0130 +0131 %EC-model: +0132 if isfield(model,'ec') +0133 fprintf(fid,'- ec-rxns:\n'); +0134 for i = 1:length(model.ec.rxns) +0135 fprintf(fid,' - !!omap\n'); +0136 writeField(model.ec, fid, 'rxns', 'txt', i, '- id', preserveQuotes) +0137 writeField(model.ec, fid, 'kcat', 'num', i, '- kcat', preserveQuotes) +0138 writeField(model.ec, fid, 'source', 'txt', i, '- source', preserveQuotes) +0139 writeField(model.ec, fid, 'notes', 'txt', i, '- notes', preserveQuotes) +0140 writeField(model.ec, fid, 'eccodes', 'txt', i, '- eccodes', preserveQuotes) +0141 writeField(model.ec, fid, 'rxnEnzMat', 'txt', i, '- enzymes', preserveQuotes) +0142 end +0143 +0144 fprintf(fid,'- ec-enzymes:\n'); +0145 for i = 1:length(model.ec.genes) +0146 fprintf(fid,' - !!omap\n'); +0147 writeField(model.ec, fid, 'genes', 'txt', i, '- genes', preserveQuotes) +0148 writeField(model.ec, fid, 'enzymes', 'txt', i, '- enzymes', preserveQuotes) +0149 writeField(model.ec, fid, 'mw', 'num', i, '- mw', preserveQuotes) +0150 writeField(model.ec, fid, 'sequence', 'txt', i, '- sequence', preserveQuotes) +0151 writeField(model.ec, fid, 'concs', 'num', i, '- concs', preserveQuotes) +0152 end +0153 end 0154 -0155 end -0156 -0157 function writeField(model,fid,fieldName,type,pos,name,preserveQuotes) -0158 %Writes a new line in the yaml file if the field exists and the field is -0159 %not empty at the correspoinding position. It's recursive for some fields -0160 %(metMiriams, rxnMiriams, and S) -0161 -0162 if isfield(model,fieldName) -0163 if strcmp(fieldName,'metComps') -0164 %metComps: write full name -0165 fieldName = 'comps'; -0166 pos = model.metComps(pos); -0167 end -0168 -0169 field = model.(fieldName); -0170 -0171 if strcmp(fieldName,'metMiriams') -0172 if ~isempty(model.metMiriams{pos}) -0173 fprintf(fid,' %s: !!omap\n',name); -0174 for i=1:size(model.newMetMiriams,2) -0175 %'i' represents the different miriam names, e.g. -0176 %kegg.compound or chebi -0177 if ~isempty(model.newMetMiriams{pos,i}) -0178 %As during the following writeField call the value of -0179 %'i' would be lost, it is temporarily concatenated to -0180 %'name' parameter, which will be edited later -0181 writeField(model, fid, 'newMetMiriams', 'txt', pos, [' - ' model.newMetMiriamNames{i} '_' sprintf('%d',i)], preserveQuotes) -0182 end -0183 end -0184 end -0185 -0186 elseif strcmp(fieldName,'rxnMiriams') -0187 if ~isempty(model.rxnMiriams{pos}) -0188 fprintf(fid,' %s: !!omap\n',name); -0189 for i=1:size(model.newRxnMiriams,2) -0190 if ~isempty(model.newRxnMiriams{pos,i}) -0191 writeField(model, fid, 'newRxnMiriams', 'txt', pos, [' - ' model.newRxnMiriamNames{i} '_' sprintf('%d',i)], preserveQuotes) -0192 end -0193 end -0194 end -0195 -0196 elseif strcmp(fieldName,'geneMiriams') -0197 if ~isempty(model.geneMiriams{pos}) -0198 fprintf(fid,' %s: !!omap\n',name); -0199 for i=1:size(model.newGeneMiriams,2) -0200 if ~isempty(model.newGeneMiriams{pos,i}) -0201 writeField(model, fid, 'newGeneMiriams', 'txt', pos, [' - ' model.newGeneMiriamNames{i} '_' sprintf('%d',i)], preserveQuotes) -0202 end -0203 end -0204 end -0205 -0206 elseif strcmp(fieldName,'compMiriams') -0207 if ~isempty(model.compMiriams{pos}) -0208 fprintf(fid,' %s: !!omap\n',name); -0209 for i=1:size(model.newCompMiriams,2) -0210 if ~isempty(model.newCompMiriams{pos,i}) -0211 writeField(model, fid, 'newCompMiriams', 'txt', pos, [' - ' model.newCompMiriamNames{i} '_' sprintf('%d',i)], preserveQuotes) -0212 end -0213 end -0214 end -0215 -0216 elseif strcmp(fieldName,'S') -0217 %S: create header & write each metabolite in a new line -0218 fprintf(fid,' %s: !!omap\n',name); -0219 if sum(field(:,pos) ~= 0) > 0 -0220 model.mets = model.mets(field(:,pos) ~= 0); -0221 model.coeffs = field(field(:,pos) ~= 0,pos); -0222 %Sort metabolites: -0223 [model.mets,order] = sort(model.mets); -0224 model.coeffs = model.coeffs(order); -0225 for i = 1:length(model.mets) -0226 writeField(model, fid, 'coeffs', 'num', i, [' - ' model.mets{i}], preserveQuotes) -0227 end -0228 end -0229 -0230 elseif strcmp(fieldName,'rxnEnzMat') -0231 %S: create header & write each enzyme in a new line -0232 fprintf(fid,' %s: !!omap\n',name); -0233 if sum(field(pos,:) ~= 0) > 0 -0234 model.enzymes = model.enzymes(field(pos,:) ~= 0); -0235 model.coeffs = field(pos,field(pos,:) ~= 0); -0236 %Sort metabolites: -0237 [model.enzymes,order] = sort(model.enzymes); -0238 model.coeffs = model.coeffs(order); -0239 for i = 1:length(model.enzymes) -0240 writeField(model, fid, 'coeffs', 'num', i, [' - ' model.enzymes{i}], preserveQuotes) -0241 end -0242 end -0243 -0244 elseif sum(strcmp({'subSystems','newMetMiriams','newRxnMiriams','newGeneMiriams','newCompMiriams','eccodes'},fieldName)) > 0 -0245 %eccodes/rxnNotes: if 1 write in 1 line, if more create header and list -0246 if strcmp(fieldName,'subSystems') -0247 list = field{pos}; %subSystems already comes in a cell array -0248 if isempty(list) -0249 return -0250 end -0251 elseif strcmp(fieldName,'newMetMiriams') -0252 index = str2double(regexprep(name,'^.+_','')); -0253 name = regexprep(name,'_\d+$',''); -0254 list = strsplit(model.newMetMiriams{pos,index},'; '); -0255 elseif strcmp(fieldName,'newRxnMiriams') -0256 index = str2double(regexprep(name,'^.+_','')); -0257 name = regexprep(name,'_\d+$',''); -0258 list = strsplit(model.newRxnMiriams{pos,index},'; '); -0259 elseif strcmp(fieldName,'newGeneMiriams') -0260 index = str2double(regexprep(name,'^.+_','')); -0261 name = regexprep(name,'_\d+$',''); -0262 list = strsplit(model.newGeneMiriams{pos,index},'; '); -0263 elseif strcmp(fieldName,'newCompMiriams') -0264 index = str2double(regexprep(name,'^.+_','')); -0265 name = regexprep(name,'_\d+$',''); -0266 list = strsplit(model.newCompMiriams{pos,index},'; '); -0267 elseif ~isempty(field{pos}) -0268 list = strrep(field{pos},' ',''); -0269 list = strsplit(list,';'); -0270 else -0271 return % empty, needs no line in file -0272 end -0273 list=strip(list); -0274 -0275 if length(list) == 1 && ~strcmp(list{1},'') && ~strcmp(fieldName,'subSystems') -0276 if preserveQuotes -0277 list = ['"' list{1} '"']; -0278 end -0279 if iscell(list) -0280 list=list{1}; +0155 %Close file: +0156 fclose(fid); +0157 +0158 end +0159 +0160 function writeField(model,fid,fieldName,type,pos,name,preserveQuotes) +0161 %Writes a new line in the yaml file if the field exists and the field is +0162 %not empty at the correspoinding position. It's recursive for some fields +0163 %(metMiriams, rxnMiriams, and S) +0164 +0165 if isfield(model,fieldName) +0166 if strcmp(fieldName,'metComps') +0167 %metComps: write full name +0168 fieldName = 'comps'; +0169 pos = model.metComps(pos); +0170 end +0171 +0172 field = model.(fieldName); +0173 +0174 if strcmp(fieldName,'metMiriams') +0175 if ~isempty(model.metMiriams{pos}) +0176 fprintf(fid,' %s: !!omap\n',name); +0177 for i=1:size(model.newMetMiriams,2) +0178 %'i' represents the different miriam names, e.g. +0179 %kegg.compound or chebi +0180 if ~isempty(model.newMetMiriams{pos,i}) +0181 %As during the following writeField call the value of +0182 %'i' would be lost, it is temporarily concatenated to +0183 %'name' parameter, which will be edited later +0184 writeField(model, fid, 'newMetMiriams', 'txt', pos, [' - ' model.newMetMiriamNames{i} '_' sprintf('%d',i)], preserveQuotes) +0185 end +0186 end +0187 end +0188 +0189 elseif strcmp(fieldName,'rxnMiriams') +0190 if ~isempty(model.rxnMiriams{pos}) +0191 fprintf(fid,' %s: !!omap\n',name); +0192 for i=1:size(model.newRxnMiriams,2) +0193 if ~isempty(model.newRxnMiriams{pos,i}) +0194 writeField(model, fid, 'newRxnMiriams', 'txt', pos, [' - ' model.newRxnMiriamNames{i} '_' sprintf('%d',i)], preserveQuotes) +0195 end +0196 end +0197 end +0198 +0199 elseif strcmp(fieldName,'geneMiriams') +0200 if ~isempty(model.geneMiriams{pos}) +0201 fprintf(fid,' %s: !!omap\n',name); +0202 for i=1:size(model.newGeneMiriams,2) +0203 if ~isempty(model.newGeneMiriams{pos,i}) +0204 writeField(model, fid, 'newGeneMiriams', 'txt', pos, [' - ' model.newGeneMiriamNames{i} '_' sprintf('%d',i)], preserveQuotes) +0205 end +0206 end +0207 end +0208 +0209 elseif strcmp(fieldName,'compMiriams') +0210 if ~isempty(model.compMiriams{pos}) +0211 fprintf(fid,' %s: !!omap\n',name); +0212 for i=1:size(model.newCompMiriams,2) +0213 if ~isempty(model.newCompMiriams{pos,i}) +0214 writeField(model, fid, 'newCompMiriams', 'txt', pos, [' - ' model.newCompMiriamNames{i} '_' sprintf('%d',i)], preserveQuotes) +0215 end +0216 end +0217 end +0218 +0219 elseif strcmp(fieldName,'S') +0220 %S: create header & write each metabolite in a new line +0221 fprintf(fid,' %s: !!omap\n',name); +0222 if sum(field(:,pos) ~= 0) > 0 +0223 model.mets = model.mets(field(:,pos) ~= 0); +0224 model.coeffs = field(field(:,pos) ~= 0,pos); +0225 %Sort metabolites: +0226 [model.mets,order] = sort(model.mets); +0227 model.coeffs = model.coeffs(order); +0228 for i = 1:length(model.mets) +0229 writeField(model, fid, 'coeffs', 'num', i, [' - ' model.mets{i}], preserveQuotes) +0230 end +0231 end +0232 +0233 elseif strcmp(fieldName,'rxnEnzMat') +0234 %S: create header & write each enzyme in a new line +0235 fprintf(fid,' %s: !!omap\n',name); +0236 if sum(field(pos,:) ~= 0) > 0 +0237 model.enzymes = model.enzymes(field(pos,:) ~= 0); +0238 model.coeffs = field(pos,field(pos,:) ~= 0); +0239 %Sort metabolites: +0240 [model.enzymes,order] = sort(model.enzymes); +0241 model.coeffs = model.coeffs(order); +0242 for i = 1:length(model.enzymes) +0243 writeField(model, fid, 'coeffs', 'num', i, [' - ' model.enzymes{i}], preserveQuotes) +0244 end +0245 end +0246 +0247 elseif sum(strcmp({'subSystems','newMetMiriams','newRxnMiriams','newGeneMiriams','newCompMiriams','eccodes'},fieldName)) > 0 +0248 %eccodes/rxnNotes: if 1 write in 1 line, if more create header and list +0249 if strcmp(fieldName,'subSystems') +0250 list = field{pos}; %subSystems already comes in a cell array +0251 if isempty(list) +0252 return +0253 end +0254 elseif strcmp(fieldName,'newMetMiriams') +0255 index = str2double(regexprep(name,'^.+_','')); +0256 name = regexprep(name,'_\d+$',''); +0257 list = strsplit(model.newMetMiriams{pos,index},'; '); +0258 elseif strcmp(fieldName,'newRxnMiriams') +0259 index = str2double(regexprep(name,'^.+_','')); +0260 name = regexprep(name,'_\d+$',''); +0261 list = strsplit(model.newRxnMiriams{pos,index},'; '); +0262 elseif strcmp(fieldName,'newGeneMiriams') +0263 index = str2double(regexprep(name,'^.+_','')); +0264 name = regexprep(name,'_\d+$',''); +0265 list = strsplit(model.newGeneMiriams{pos,index},'; '); +0266 elseif strcmp(fieldName,'newCompMiriams') +0267 index = str2double(regexprep(name,'^.+_','')); +0268 name = regexprep(name,'_\d+$',''); +0269 list = strsplit(model.newCompMiriams{pos,index},'; '); +0270 elseif ~isempty(field{pos}) +0271 list = strrep(field{pos},' ',''); +0272 list = strsplit(list,';'); +0273 else +0274 return % empty, needs no line in file +0275 end +0276 list=strip(list); +0277 +0278 if length(list) == 1 && ~strcmp(list{1},'') && ~strcmp(fieldName,'subSystems') +0279 if preserveQuotes +0280 list = ['"' list{1} '"']; 0281 end -0282 fprintf(fid,' %s: %s\n',name,list); -0283 elseif ischar(list) && strcmp(fieldName,'subSystems') -0284 if preserveQuotes -0285 list = ['"' list '"']; -0286 end -0287 fprintf(fid,' %s: %s\n',name,list); -0288 elseif length(list) > 1 || strcmp(fieldName,'subSystems') -0289 if preserveQuotes -0290 for j=1:numel(list) -0291 list{j} = ['"' list{j} '"']; -0292 end -0293 end -0294 fprintf(fid,' %s:\n',name); -0295 for i = 1:length(list) -0296 fprintf(fid,'%s - %s\n',regexprep(name,'(^\s*).*','$1'),list{i}); -0297 end -0298 end -0299 -0300 elseif sum(pos) > 0 -0301 %All other fields: -0302 if strcmp(type,'txt') -0303 value = field{pos}; -0304 if preserveQuotes && ~isempty(value) -0305 value = ['"',value,'"']; -0306 end -0307 elseif strcmp(type,'num') -0308 if isnan(field(pos)) -0309 value = []; -0310 else -0311 value = sprintf('%.15g',full(field(pos))); -0312 end -0313 end -0314 if ~isempty(value) -0315 fprintf(fid,' %s: %s\n',name,value); +0282 if iscell(list) +0283 list=list{1}; +0284 end +0285 fprintf(fid,' %s: %s\n',name,list); +0286 elseif ischar(list) && strcmp(fieldName,'subSystems') +0287 if preserveQuotes +0288 list = ['"' list '"']; +0289 end +0290 fprintf(fid,' %s: %s\n',name,list); +0291 elseif length(list) > 1 || strcmp(fieldName,'subSystems') +0292 if preserveQuotes +0293 for j=1:numel(list) +0294 list{j} = ['"' list{j} '"']; +0295 end +0296 end +0297 fprintf(fid,' %s:\n',name); +0298 for i = 1:length(list) +0299 fprintf(fid,'%s - %s\n',regexprep(name,'(^\s*).*','$1'),list{i}); +0300 end +0301 end +0302 +0303 elseif sum(pos) > 0 +0304 %All other fields: +0305 if strcmp(type,'txt') +0306 value = field{pos}; +0307 if preserveQuotes && ~isempty(value) +0308 value = ['"',value,'"']; +0309 end +0310 elseif strcmp(type,'num') +0311 if isnan(field(pos)) +0312 value = []; +0313 else +0314 value = sprintf('%.15g',full(field(pos))); +0315 end 0316 end -0317 end -0318 end -0319 end -0320 -0321 function writeMetadata(model,fid) -0322 % Writes model metadata to the yaml file. This information will eventually -0323 % be extracted entirely from the model, but for now, many of the entries -0324 % are hard-coded defaults for HumanGEM. -0325 -0326 fprintf(fid, '- metaData:\n'); -0327 if isfield(model,'id') -0328 fprintf(fid, ' id: "%s"\n', model.id); -0329 else -0330 fprintf(fid, ' id: "blankID"\n'); -0331 end -0332 if isfield(model,'name') -0333 fprintf(fid, ' name: "%s"\n',model.name); -0334 else -0335 fprintf(fid, ' name: "blankName"\n'); -0336 end -0337 if isfield(model,'version') -0338 fprintf(fid, ' version: "%s"\n',model.version); +0317 if ~isempty(value) +0318 fprintf(fid,' %s: %s\n',name,value); +0319 end +0320 end +0321 end +0322 end +0323 +0324 function writeMetadata(model,fid) +0325 % Writes model metadata to the yaml file. This information will eventually +0326 % be extracted entirely from the model, but for now, many of the entries +0327 % are hard-coded defaults for HumanGEM. +0328 +0329 fprintf(fid, '- metaData:\n'); +0330 if isfield(model,'id') +0331 fprintf(fid, ' id: "%s"\n', model.id); +0332 else +0333 fprintf(fid, ' id: "blankID"\n'); +0334 end +0335 if isfield(model,'name') +0336 fprintf(fid, ' name: "%s"\n',model.name); +0337 else +0338 fprintf(fid, ' name: "blankName"\n'); 0339 end -0340 fprintf(fid, ' date: "%s"\n',datestr(now,29)); % 29=YYYY-MM-DD -0341 if isfield(model,'annotation') -0342 if isfield(model.annotation,'defaultLB') -0343 fprintf(fid, ' defaultLB: "%g"\n', model.annotation.defaultLB); -0344 end -0345 if isfield(model.annotation,'defaultUB') -0346 fprintf(fid, ' defaultUB: "%g"\n', model.annotation.defaultUB); +0340 if isfield(model,'version') +0341 fprintf(fid, ' version: "%s"\n',model.version); +0342 end +0343 fprintf(fid, ' date: "%s"\n',datestr(now,29)); % 29=YYYY-MM-DD +0344 if isfield(model,'annotation') +0345 if isfield(model.annotation,'defaultLB') +0346 fprintf(fid, ' defaultLB: "%g"\n', model.annotation.defaultLB); 0347 end -0348 if isfield(model.annotation,'givenName') -0349 fprintf(fid, ' givenName: "%s"\n', model.annotation.givenName); +0348 if isfield(model.annotation,'defaultUB') +0349 fprintf(fid, ' defaultUB: "%g"\n', model.annotation.defaultUB); 0350 end -0351 if isfield(model.annotation,'familyName') -0352 fprintf(fid, ' familyName: "%s"\n', model.annotation.familyName); +0351 if isfield(model.annotation,'givenName') +0352 fprintf(fid, ' givenName: "%s"\n', model.annotation.givenName); 0353 end -0354 if isfield(model.annotation,'authors') -0355 fprintf(fid, ' authors: "%s"\n', model.annotation.authors); +0354 if isfield(model.annotation,'familyName') +0355 fprintf(fid, ' familyName: "%s"\n', model.annotation.familyName); 0356 end -0357 if isfield(model.annotation,'email') -0358 fprintf(fid, ' email: "%s"\n', model.annotation.email); +0357 if isfield(model.annotation,'authors') +0358 fprintf(fid, ' authors: "%s"\n', model.annotation.authors); 0359 end -0360 if isfield(model.annotation,'organization') -0361 fprintf(fid, ' organization: "%s"\n',model.annotation.organization); +0360 if isfield(model.annotation,'email') +0361 fprintf(fid, ' email: "%s"\n', model.annotation.email); 0362 end -0363 if isfield(model.annotation,'taxonomy') -0364 fprintf(fid, ' taxonomy: "%s"\n', model.annotation.taxonomy); +0363 if isfield(model.annotation,'organization') +0364 fprintf(fid, ' organization: "%s"\n',model.annotation.organization); 0365 end -0366 if isfield(model.annotation,'note') -0367 fprintf(fid, ' note: "%s"\n', model.annotation.note); +0366 if isfield(model.annotation,'taxonomy') +0367 fprintf(fid, ' taxonomy: "%s"\n', model.annotation.taxonomy); 0368 end -0369 if isfield(model.annotation,'sourceUrl') -0370 fprintf(fid, ' sourceUrl: "%s"\n', model.annotation.sourceUrl); +0369 if isfield(model.annotation,'note') +0370 fprintf(fid, ' note: "%s"\n', model.annotation.note); 0371 end -0372 end -0373 if isfield(model,'ec') -0374 if model.ec.geckoLight -0375 geckoLight = 'true'; -0376 else -0377 geckoLight = 'false'; -0378 end -0379 fprintf(fid,' geckoLight: "%s"\n',geckoLight); -0380 end -0381 end +0372 if isfield(model.annotation,'sourceUrl') +0373 fprintf(fid, ' sourceUrl: "%s"\n', model.annotation.sourceUrl); +0374 end +0375 end +0376 if isfield(model,'ec') +0377 if model.ec.geckoLight +0378 geckoLight = 'true'; +0379 else +0380 geckoLight = 'false'; +0381 end +0382 fprintf(fid,' geckoLight: "%s"\n',geckoLight); +0383 end +0384 end
Generated by m2html © 2005
\ No newline at end of file diff --git a/doc/solver/optimizeProb.html b/doc/solver/optimizeProb.html index ebc32e00..6a0138bd 100644 --- a/doc/solver/optimizeProb.html +++ b/doc/solver/optimizeProb.html @@ -144,191 +144,194 @@

SOURCE CODE ^structUpdate(solverparams,params); -0094 -0095 % Restructering problem according to gurobi format -0096 if isfield(prob, 'csense') -0097 prob.sense = renameparams(prob.csense, {'L','G','E'}, {'<','>','='}); -0098 prob = rmfield(prob, {'csense'}); -0099 end -0100 if isfield(prob, 'osense') -0101 osense = prob.osense; -0102 prob.modelsense = renameparams(num2str(prob.osense), {'1','-1'}, {'min','max'}); -0103 prob = rmfield(prob, {'osense'}); -0104 end -0105 [prob.obj, prob.rhs, prob.vtype] = deal(prob.c, prob.b, prob.vartype); -0106 prob = rmfield(prob, {'c','b','vartype'}); -0107 -0108 resG = gurobi(prob,solverparams); -0109 -0110 try -0111 % Name output fields the same as COBRA does -0112 res.full = resG.x; -0113 res.obj = resG.objval; -0114 res.origStat = resG.status; -0115 if isfield(resG,{'pi','rc'}) -0116 res.dual = -resG.pi*osense; -0117 res.rcost = -resG.rc*osense; -0118 end -0119 if milp && strcmp(resG.status, 'TIME_LIMIT') -0120 % If res has the objval field, it succeeded, regardless of -0121 % time_limit status -0122 resG.status = 'OPTIMAL'; -0123 end -0124 switch resG.status -0125 case 'OPTIMAL' -0126 res.stat = 1; -0127 case 'UNBOUNDED' -0128 res.stat = 2; -0129 otherwise -0130 res.stat = 0; -0131 end -0132 if ~milp -0133 res.vbasis = resG.vbasis; -0134 res.cbasis = resG.cbasis; -0135 else -0136 res.mipgap = resG.mipgap; -0137 end -0138 catch -0139 res.stat = 0; -0140 res.origStat = resG.status; % useful information to have -0141 end -0142 %% Use GLPK using RAVEN-provided binary -0143 case 'glpk' -0144 solverparams.scale = 1; % Auto scaling -0145 %solverparams.tmlim = defaultparams.timeLimit; -0146 solverparams.tolbnd = defaultparams.feasTol; -0147 solverparams.toldj = defaultparams.optTol; -0148 solverparams.tolint = defaultparams.intTol; -0149 solverparams.tolobj = defaultparams.objTol; -0150 solverparams.msglev = 0; % Level of verbosity -0151 solverparams = structUpdate(solverparams,params); -0152 -0153 prob.csense = renameparams(prob.csense, {'L','G','E'}, {'U','L','S'}); -0154 -0155 if milp -0156 solverparams.tmlim = solverparams.tmlim*10; -0157 solverparams.msglev = 1; % Level of verbosity -0158 disp('Issues have been observed when using GLPK for MILP solving. Be advised to carefully observe the results, or us another solver.') -0159 end -0160 -0161 % Ensure that RAVEN glpk binary is used, return to original -0162 % directory afterwards -0163 [ravenDir,currDir]=findRAVENroot(); -0164 cd(fullfile(ravenDir,'software','GLPKmex')) -0165 [xopt, fmin, errnum, extra] = glpk(prob.c, prob.A, prob.b, prob.lb, prob.ub, prob.csense, prob.vartype, prob.osense, solverparams); -0166 cd(currDir) -0167 -0168 switch errnum % 1 = undefined; 2 = feasible; 3 = infeasible; 4 = no feasible solution; 5 = optimal; 6 = no unbounded solution -0169 case 5 -0170 res.stat = 1; % Optimal -0171 case 2 -0172 res.stat = 2; % Feasible, but not optimal -0173 otherwise -0174 res.stat = 0; -0175 end -0176 res.origStat = errnum; -0177 res.full = xopt; -0178 res.obj = fmin; -0179 res.dual = -extra.lambda*prob.osense; -0180 res.rcost = -extra.redcosts*prob.osense; -0181 %% Use scip -0182 case {'soplex','scip'} % Old 'soplex' option also allowed -0183 [xopt,fval,exitflag] = scip([], prob.c, prob.A,-prob.b, prob.b, prob.lb, prob.ub, prob.vartype); -0184 -0185 % [x,fval,exitflag,stats] = scip(H, f, A, rl, ru, lb, ub, xtype, sos, qc, nl, x0, opts) -0186 % -0187 % Input arguments*: -0188 % H - quadratic objective matrix (sparse, optional [NOT TRIL / TRIU]) -0189 % f - linear objective vector -0190 % A - linear constraint matrix (sparse) -0191 % rl - linear constraint lhs -0192 % ru - linear constraint rhs -0193 % lb - decision variable lower bounds -0194 % ub - decision variable upper bounds -0195 % xtype - string of variable integrality ('c' continuous, 'i' integer, 'b' binary) -0196 % sos - SOS structure with fields type, index and weight (see below) -0197 % qc - Quadratic Constraints structure with fields Q, l, qrl and qru (see below) -0198 % nl - Nonlinear Objective and Constraints structure (see below) -0199 % x0 - primal solution -0200 % opts - solver options (see below) -0201 % -0202 % Return arguments: -0203 % x - solution vector -0204 % fval - objective value at the solution -0205 % exitflag - exit status (see below) -0206 % stats - statistics structure -0207 % -0208 % Option Fields (all optional, see also optiset for a list): -0209 % solverOpts - specific SCIP options (list of pairs of parameter names and values) -0210 % maxiter - maximum LP solver iterations -0211 % maxnodes - maximum nodes to explore -0212 % maxtime - maximum execution time [s] -0213 % tolrfun - primal feasibility tolerance -0214 % display - solver display level [0-5] -0215 % probfile - write problem to given file -0216 % presolvedfile - write presolved problem to file -0217 % -0218 % Return Status: -0219 % 0 - Unknown -0220 % 1 - User Interrupted -0221 % 2 - Node Limit Reached -0222 % 3 - Total Node Limit Reached -0223 % 4 - Stall Node Limit Reached -0224 % 5 - Time Limit Reached -0225 % 6 - Memory Limit Reached -0226 % 7 - Gap Limit Reached -0227 % 8 - Solution Limit Reached -0228 % 9 - Solution Improvement Limit Reached -0229 % 10 - Restart Limit Reached -0230 % 11 - Problem Solved to Optimality -0231 % 12 - Problem is Infeasible -0232 % 13 - Problem is Unbounded -0233 % 14 - Problem is Either Infeasible or Unbounded -0234 -0235 res.origStat = exitflag; -0236 res.full = xopt; -0237 res.obj = fval; -0238 -0239 switch exitflag -0240 case 11 -0241 res.stat = 1; -0242 case [5, 6, 7, 8, 9, 10, 13] -0243 res.stat = 2; -0244 otherwise -0245 res.stat = 0; -0246 end -0247 otherwise -0248 error('RAVEN solver not defined or unknown. Try using setRavenSolver(''solver'').'); -0249 end -0250 if res.stat>0 -0251 res.full=res.full(1:size(prob.a,2)); +0093 if ~isempty(getCurrentTask) % If run in parallel, then one thread per gurobi +0094 solverparams.Threads=1; +0095 end +0096 solverparams = structUpdate(solverparams,params); +0097 +0098 % Restructering problem according to gurobi format +0099 if isfield(prob, 'csense') +0100 prob.sense = renameparams(prob.csense, {'L','G','E'}, {'<','>','='}); +0101 prob = rmfield(prob, {'csense'}); +0102 end +0103 if isfield(prob, 'osense') +0104 osense = prob.osense; +0105 prob.modelsense = renameparams(num2str(prob.osense), {'1','-1'}, {'min','max'}); +0106 prob = rmfield(prob, {'osense'}); +0107 end +0108 [prob.obj, prob.rhs, prob.vtype] = deal(prob.c, prob.b, prob.vartype); +0109 prob = rmfield(prob, {'c','b','vartype'}); +0110 +0111 resG = gurobi(prob,solverparams); +0112 +0113 try +0114 % Name output fields the same as COBRA does +0115 res.full = resG.x; +0116 res.obj = resG.objval; +0117 res.origStat = resG.status; +0118 if isfield(resG,{'pi','rc'}) +0119 res.dual = -resG.pi*osense; +0120 res.rcost = -resG.rc*osense; +0121 end +0122 if milp && strcmp(resG.status, 'TIME_LIMIT') +0123 % If res has the objval field, it succeeded, regardless of +0124 % time_limit status +0125 resG.status = 'OPTIMAL'; +0126 end +0127 switch resG.status +0128 case 'OPTIMAL' +0129 res.stat = 1; +0130 case 'UNBOUNDED' +0131 res.stat = 2; +0132 otherwise +0133 res.stat = 0; +0134 end +0135 if ~milp +0136 res.vbasis = resG.vbasis; +0137 res.cbasis = resG.cbasis; +0138 else +0139 res.mipgap = resG.mipgap; +0140 end +0141 catch +0142 res.stat = 0; +0143 res.origStat = resG.status; % useful information to have +0144 end +0145 %% Use GLPK using RAVEN-provided binary +0146 case 'glpk' +0147 solverparams.scale = 1; % Auto scaling +0148 %solverparams.tmlim = defaultparams.timeLimit; +0149 solverparams.tolbnd = defaultparams.feasTol; +0150 solverparams.toldj = defaultparams.optTol; +0151 solverparams.tolint = defaultparams.intTol; +0152 solverparams.tolobj = defaultparams.objTol; +0153 solverparams.msglev = 0; % Level of verbosity +0154 solverparams = structUpdate(solverparams,params); +0155 +0156 prob.csense = renameparams(prob.csense, {'L','G','E'}, {'U','L','S'}); +0157 +0158 if milp +0159 solverparams.tmlim = solverparams.tmlim*10; +0160 solverparams.msglev = 1; % Level of verbosity +0161 disp('Issues have been observed when using GLPK for MILP solving. Be advised to carefully observe the results, or us another solver.') +0162 end +0163 +0164 % Ensure that RAVEN glpk binary is used, return to original +0165 % directory afterwards +0166 [ravenDir,currDir]=findRAVENroot(); +0167 cd(fullfile(ravenDir,'software','GLPKmex')) +0168 [xopt, fmin, errnum, extra] = glpk(prob.c, prob.A, prob.b, prob.lb, prob.ub, prob.csense, prob.vartype, prob.osense, solverparams); +0169 cd(currDir) +0170 +0171 switch errnum % 1 = undefined; 2 = feasible; 3 = infeasible; 4 = no feasible solution; 5 = optimal; 6 = no unbounded solution +0172 case 5 +0173 res.stat = 1; % Optimal +0174 case 2 +0175 res.stat = 2; % Feasible, but not optimal +0176 otherwise +0177 res.stat = 0; +0178 end +0179 res.origStat = errnum; +0180 res.full = xopt; +0181 res.obj = fmin; +0182 res.dual = -extra.lambda*prob.osense; +0183 res.rcost = -extra.redcosts*prob.osense; +0184 %% Use scip +0185 case {'soplex','scip'} % Old 'soplex' option also allowed +0186 [xopt,fval,exitflag] = scip([], prob.c, prob.A,-prob.b, prob.b, prob.lb, prob.ub, prob.vartype); +0187 +0188 % [x,fval,exitflag,stats] = scip(H, f, A, rl, ru, lb, ub, xtype, sos, qc, nl, x0, opts) +0189 % +0190 % Input arguments*: +0191 % H - quadratic objective matrix (sparse, optional [NOT TRIL / TRIU]) +0192 % f - linear objective vector +0193 % A - linear constraint matrix (sparse) +0194 % rl - linear constraint lhs +0195 % ru - linear constraint rhs +0196 % lb - decision variable lower bounds +0197 % ub - decision variable upper bounds +0198 % xtype - string of variable integrality ('c' continuous, 'i' integer, 'b' binary) +0199 % sos - SOS structure with fields type, index and weight (see below) +0200 % qc - Quadratic Constraints structure with fields Q, l, qrl and qru (see below) +0201 % nl - Nonlinear Objective and Constraints structure (see below) +0202 % x0 - primal solution +0203 % opts - solver options (see below) +0204 % +0205 % Return arguments: +0206 % x - solution vector +0207 % fval - objective value at the solution +0208 % exitflag - exit status (see below) +0209 % stats - statistics structure +0210 % +0211 % Option Fields (all optional, see also optiset for a list): +0212 % solverOpts - specific SCIP options (list of pairs of parameter names and values) +0213 % maxiter - maximum LP solver iterations +0214 % maxnodes - maximum nodes to explore +0215 % maxtime - maximum execution time [s] +0216 % tolrfun - primal feasibility tolerance +0217 % display - solver display level [0-5] +0218 % probfile - write problem to given file +0219 % presolvedfile - write presolved problem to file +0220 % +0221 % Return Status: +0222 % 0 - Unknown +0223 % 1 - User Interrupted +0224 % 2 - Node Limit Reached +0225 % 3 - Total Node Limit Reached +0226 % 4 - Stall Node Limit Reached +0227 % 5 - Time Limit Reached +0228 % 6 - Memory Limit Reached +0229 % 7 - Gap Limit Reached +0230 % 8 - Solution Limit Reached +0231 % 9 - Solution Improvement Limit Reached +0232 % 10 - Restart Limit Reached +0233 % 11 - Problem Solved to Optimality +0234 % 12 - Problem is Infeasible +0235 % 13 - Problem is Unbounded +0236 % 14 - Problem is Either Infeasible or Unbounded +0237 +0238 res.origStat = exitflag; +0239 res.full = xopt; +0240 res.obj = fval; +0241 +0242 switch exitflag +0243 case 11 +0244 res.stat = 1; +0245 case [5, 6, 7, 8, 9, 10, 13] +0246 res.stat = 2; +0247 otherwise +0248 res.stat = 0; +0249 end +0250 otherwise +0251 error('RAVEN solver not defined or unknown. Try using setRavenSolver(''solver'').'); 0252 end -0253 end -0254 -0255 function s_merged=structUpdate(s_old,s_new) -0256 %Remove overlapping fields from first struct; -0257 %Obtain all unique names of remaining fields; -0258 %Merge both structs -0259 s_merged = rmfield(s_old, intersect(fieldnames(s_old), fieldnames(s_new))); -0260 names = [fieldnames(s_merged); fieldnames(s_new)]; -0261 s_merged = cell2struct([struct2cell(s_merged); struct2cell(s_new)], names, 1); -0262 end -0263 -0264 function paramlist = renameparams(paramlist,old,new) -0265 if ~iscell(paramlist) -0266 wasNoCell = true; -0267 paramlist={paramlist}; -0268 else -0269 wasNoCell = false; -0270 end -0271 for i=1:numel(old) -0272 paramlist = regexprep(paramlist,old{i},new{i}); +0253 if res.stat>0 +0254 res.full=res.full(1:size(prob.a,2)); +0255 end +0256 end +0257 +0258 function s_merged=structUpdate(s_old,s_new) +0259 %Remove overlapping fields from first struct; +0260 %Obtain all unique names of remaining fields; +0261 %Merge both structs +0262 s_merged = rmfield(s_old, intersect(fieldnames(s_old), fieldnames(s_new))); +0263 names = [fieldnames(s_merged); fieldnames(s_new)]; +0264 s_merged = cell2struct([struct2cell(s_merged); struct2cell(s_new)], names, 1); +0265 end +0266 +0267 function paramlist = renameparams(paramlist,old,new) +0268 if ~iscell(paramlist) +0269 wasNoCell = true; +0270 paramlist={paramlist}; +0271 else +0272 wasNoCell = false; 0273 end -0274 if wasNoCell -0275 paramlist=paramlist{1}; +0274 for i=1:numel(old) +0275 paramlist = regexprep(paramlist,old{i},new{i}); 0276 end -0277 end +0277 if wasNoCell +0278 paramlist=paramlist{1}; +0279 end +0280 end
Generated by m2html © 2005
\ No newline at end of file diff --git a/doc/solver/solveLP.html b/doc/solver/solveLP.html index 8c462e72..9cec93be 100644 --- a/doc/solver/solveLP.html +++ b/doc/solver/solveLP.html @@ -151,8 +151,8 @@

SOURCE CODE ^'Invalid defintion of objective function in model.c.') 0071 end 0072 %Check for valid S-matrix -0073 invalidS = ~isfinite(model.S); -0074 if any(any(invalidS)) +0073 if ~allfinite(model.S) +0074 invalidS = ~isfinite(model.S); 0075 error(['Invalid coefficients defined for reaction(s): ', strjoin(model.rxns(any(invalidS)),', '), '.']) 0076 end 0077 diff --git a/io/exportToExcelFormat.m b/io/exportToExcelFormat.m index a4c92f8f..7274c387 100755 --- a/io/exportToExcelFormat.m +++ b/io/exportToExcelFormat.m @@ -13,9 +13,6 @@ function exportToExcelFormat(model,fileName,sortIds) % sorted alphabetically by their identifiers (optional, % default false) % -% No checks are made regarding the correctness of the model. Use -% checkModelStruct to identify problems in the model structure. -% % Usage: exportToExcelFormat(model, fileName, sortIds) if nargin<2 || isempty(fileName) @@ -34,6 +31,8 @@ function exportToExcelFormat(model,fileName,sortIds) model=sortIdentifiers(model); end +checkModelStruct(model); + addList = matlab.addons.installedAddons; if any(strcmpi(addList.Name,'Text Analytics Toolbox')) error(['exportToExcelFormat is incompatible with MATLAB Text Analytics Toolbox. ' ... diff --git a/io/readYAMLmodel.m b/io/readYAMLmodel.m index ac3f256a..bc93d78c 100755 --- a/io/readYAMLmodel.m +++ b/io/readYAMLmodel.m @@ -45,7 +45,7 @@ % leading spaces to avoid metaData to be concatenated. newLine=regexp(line_raw,'^ {6,}([\w\(\)].*)','tokenExtents'); brokenLine=find(~cellfun('isempty',newLine)); -for i=1:numel(brokenLine) +for i=flip(1:numel(brokenLine)) extraLine = char(line_raw(brokenLine(i))); extraLine = extraLine(newLine{brokenLine(i)}{1}(1):end); line_raw{brokenLine(i)-1} = strjoin({line_raw{brokenLine(i)-1},extraLine},' '); @@ -58,7 +58,7 @@ line_value = regexprep(line_raw, '.*:$',''); line_value = regexprep(line_value, '[^":]+: "?(.+)"?$','$1'); line_value = regexprep(line_value, '(")|(^ {4,}- )',''); - +line_value(strcmp(line_value,'''''')) = {''}; line_value(strcmp(line_value,line_raw)) = {''}; modelFields = {'id',char();... @@ -140,7 +140,7 @@ geneMiriams=cell(100000,3); genMirNo=1; subSystems=cell(100000,2); subSysNo=1; eccodes=cell(100000,2); ecCodeNo=1; -equations=cell(100000,3); equatiNo=1; +equations=cell(100000,3); equatiNo=1; for i=1:numel(line_key) tline_raw = line_raw{i}; diff --git a/io/writeYAMLmodel.m b/io/writeYAMLmodel.m index 2b2bd612..2042c5e3 100755 --- a/io/writeYAMLmodel.m +++ b/io/writeYAMLmodel.m @@ -38,6 +38,9 @@ function writeYAMLmodel(model,fileName,preserveQuotes,sortIds) model = ravenCobraWrapper(model); end +%Check that the model structure has no problems +checkModelStruct(model); + %Sort identifiers alphabetically if sortIds == true model = sortIdentifiers(model); diff --git a/solver/optimizeProb.m b/solver/optimizeProb.m index cf33d92d..3aa3bce2 100755 --- a/solver/optimizeProb.m +++ b/solver/optimizeProb.m @@ -90,6 +90,9 @@ solverparams.FeasibilityTol = defaultparams.feasTol; solverparams.OptimalityTol = defaultparams.optTol; solverparams.Presolve = 2; + if ~isempty(getCurrentTask) % If run in parallel, then one thread per gurobi + solverparams.Threads=1; + end solverparams = structUpdate(solverparams,params); % Restructering problem according to gurobi format diff --git a/solver/solveLP.m b/solver/solveLP.m index ddd123d7..2fba3805 100755 --- a/solver/solveLP.m +++ b/solver/solveLP.m @@ -70,8 +70,8 @@ error('Invalid defintion of objective function in model.c.') end %Check for valid S-matrix -invalidS = ~isfinite(model.S); -if any(any(invalidS)) +if ~allfinite(model.S) + invalidS = ~isfinite(model.S); error(['Invalid coefficients defined for reaction(s): ', strjoin(model.rxns(any(invalidS)),', '), '.']) end