0001 function model = readCbModel(fileName,defaultBound,fileType,modelDescription,compSymbolList,compNameList)
0002
0003
0004
0005
0006
0007
0008
0009
0010
0011
0012
0013
0014
0015
0016
0017
0018
0019
0020
0021
0022
0023
0024
0025
0026
0027
0028
0029
0030
0031
0032
0033
0034
0035
0036
0037
0038
0039
0040
0041
0042
0043
0044
0045
0046
0047
0048
0049
0050
0051
0052
0053
0054
0055
0056
0057
0058
0059
0060
0061
0062
0063
0064
0065
0066
0067
0068
0069
0070
0071
0072
0073
0074
0075
0076
0077
0078 if (nargin < 2)
0079 defaultBound = 1000;
0080 else
0081 if (isempty(defaultBound))
0082 defaultBound = 1000;
0083 end
0084 end
0085 if (nargin < 3)
0086 fileType = 'SBML';
0087 if (isempty(fileType))
0088 fileType = 'SBML';
0089 end
0090 end
0091
0092
0093
0094
0095 if (nargin < 1)
0096 [fileNameFull,filePath] = uigetfile({'*.xml';'*.sto'});
0097 if (fileNameFull)
0098 [t1,t2,t3,t4,tokens] = regexp(fileNameFull,'(.*)\.(\w*)');
0099 noPathName = tokens{1}{1};
0100 fileTypeTmp = tokens{1}{2};
0101 fileName = [filePath noPathName];
0102 else
0103 model = 0;
0104 return;
0105 end
0106 switch fileTypeTmp
0107 case 'xml'
0108 fileType = 'SBML';
0109 case 'sto'
0110 fileType = 'SimPheny';
0111 otherwise
0112 error('Cannot process this file type');
0113 end
0114 end
0115
0116 if (nargin < 4)
0117 if (exist('filePath'))
0118 modelDescription = noPathName;
0119 else
0120 modelDescription = fileName;
0121 end
0122 end
0123
0124 if (nargin < 5)
0125 compSymbolList = {};
0126 compNameList = {};
0127 end
0128
0129 switch fileType
0130 case 'SBML',
0131 if isempty(regexp(fileName,'\.xml$', 'once'))
0132 model = readSBMLCbModel([fileName '.xml'],defaultBound,compSymbolList,compNameList);
0133 else
0134 model = readSBMLCbModel(fileName,defaultBound,compSymbolList,compNameList);
0135 end
0136 case 'SimPheny',
0137 model = readSimPhenyCbModel(fileName,defaultBound,compSymbolList,compNameList);
0138 case 'SimPhenyPlus',
0139 model = readSimPhenyCbModel(fileName,defaultBound,compSymbolList,compNameList);
0140 model = readSimPhenyGprCmpd(fileName,model);
0141 case 'SimPhenyText',
0142 model = readSimPhenyCbModel(fileName,defaultBound,compSymbolList,compNameList);
0143 model = readSimPhenyGprText([fileName '_gpra.txt'],model);
0144 otherwise,
0145 error('Unknown file type');
0146 end
0147
0148
0149 model = checkReversibility(model);
0150
0151
0152 checkCobraModelUnique(model);
0153
0154 model.b = zeros(length(model.mets),1);
0155
0156 model.description = modelDescription;
0157
0158
0159
0160
0161 function model = checkReversibility(model)
0162
0163 selRev = (model.lb < 0 & model.ub > 0);
0164 model.rev(selRev) = 1;
0165
0166
0167 function model = readSBMLCbModel(fileName,defaultBound,compSymbolList,compNameList)
0168
0169 if ~(exist(fileName,'file'))
0170 error(['Input file ' fileName ' not found']);
0171 end
0172
0173 if isempty(compSymbolList)
0174 compSymbolList = {'c','m','v','x','e','t','g','r','n','p'};
0175 compNameList = {'Cytosol','Mitochondria','Vacuole','Peroxisome','Extra-organism','Pool','Golgi Apparatus','Endoplasmic Reticulum','Nucleus','Periplasm'};
0176 end
0177
0178
0179 modelSBML = TranslateSBML(fileName);
0180
0181
0182 model = convertSBMLToCobra(modelSBML,defaultBound,compSymbolList,compNameList);
0183
0184
0185 function model = readSimPhenyCbModel(baseName,defaultBound,compSymbolList,compNameList)
0186
0187
0188
0189
0190
0191
0192
0193
0194
0195
0196
0197
0198
0199
0200
0201
0202
0203 if (nargin < 2)
0204 defaultBound = 1000;
0205 end
0206
0207 if ~(exist([baseName '.met'],'file') & exist([baseName '.rxn'],'file') & exist([baseName '.sto'],'file'))
0208 error('One or more input files not found');
0209 end
0210
0211 if isempty(compSymbolList)
0212 compSymbolList = {'c','m','v','x','e','t','g','r','n','p'};
0213 compNameList = {'Cytosol','Mitochondria','Vacuole','Peroxisome','Extra-organism','Pool','Golgi Apparatus','Endoplasmic Reticulum','Nucleus','Periplasm'};
0214 end
0215
0216
0217 fid = fopen([baseName '.met']);
0218 cnt = 0;
0219 while 1
0220 tline = fgetl(fid);
0221 if ~ischar(tline), break, end
0222 if (~isempty(regexp(tline,'^\d', 'once')))
0223 cnt = cnt + 1;
0224 fields = splitString(tline,'\t');
0225 mets{cnt} = fields{2};
0226
0227
0228
0229
0230 comp{cnt,1} = fields{4};
0231
0232 compSymb = 0;
0233 TF = strcmp('comp{cnt}', compNameList);
0234 for n = 1:length(compNameList)
0235 if TF(n)
0236 compSymb = compSymbolList{n};
0237 end
0238 end
0239 if (isempty(compSymb))
0240 compSymb = comp{cnt};
0241 end
0242 if (~isempty(regexp(mets{cnt},'\(', 'once')))
0243 mets{cnt} = strrep(mets{cnt},'(','[');
0244 mets{cnt} = strrep(mets{cnt},')',']');
0245 else
0246 mets{cnt} = [mets{cnt} '[' compSymb ']'];
0247 end
0248 metNames{cnt} = fields{3};
0249 end
0250 end
0251 fclose(fid);
0252
0253 mets = columnVector(mets);
0254 metNames = columnVector(metNames);
0255
0256
0257 fid = fopen([baseName '.rxn']);
0258 cnt = 0;
0259 startRxns = false;
0260 while 1
0261 tline = fgetl(fid);
0262 if ~ischar(tline), break, end
0263 if (regexp(tline,'^REACTION'))
0264 startRxns = true;
0265 end
0266 if (startRxns & ~isempty(regexp(tline,'^\d', 'once')))
0267 cnt = cnt + 1;
0268 fields = splitString(tline,'\t');
0269 rxns{cnt} = fields{2};
0270 rxnNames{cnt} = fields{3};
0271 revStr{cnt} = fields{4};
0272 lb(cnt) = str2num(fields{5});
0273 ub(cnt) = str2num(fields{6});
0274 c(cnt) = str2num(fields{7});
0275 end
0276 end
0277 fclose(fid);
0278
0279 revStr = columnVector(revStr);
0280 rev = strcmp(revStr,'Reversible');
0281 rxns = columnVector(rxns);
0282 rxnNames = columnVector(rxnNames);
0283 lb = columnVector(lb);
0284 ub = columnVector(ub);
0285 c = columnVector(c);
0286 lb(lb < -defaultBound) = -defaultBound;
0287 ub(ub > defaultBound) = defaultBound;
0288
0289
0290 fid = fopen([baseName '.sto']);
0291 fid2 = fopen('load_simpheny.tmp','w');
0292 while 1
0293 tline = fgetl(fid);
0294 if ~ischar(tline), break, end
0295
0296
0297 if (~isempty(regexp(tline,'^[-0123456789.]', 'once')))
0298 fprintf(fid2,[tline '\n']);
0299 else
0300
0301
0302 end
0303 end
0304 fclose(fid);
0305 fclose(fid2);
0306 S = load('load_simpheny.tmp');
0307
0308
0309
0310
0311
0312
0313
0314
0315
0316
0317
0318
0319
0320
0321
0322
0323
0324
0325 model.mets = mets;
0326
0327 model.metComps = comp;
0328
0329 model.metNames = metNames;
0330 model.rxns = removeDeletedTags(rxns);
0331 model.rxnNames = rxnNames;
0332 model.rev = rev;
0333 model.lb = lb;
0334 model.ub = ub;
0335 model.c = c;
0336 model.S = sparse(S);
0337
0338
0339 delete('load_simpheny.tmp');
0340
0341
0342 function list = removeDeletedTags(list)
0343
0344
0345
0346
0347
0348
0349 for i = 1:length(list)
0350 item = list{i};
0351 ind = strfind(item,' [deleted');
0352 if (~isempty(ind))
0353 list{i} = item(1:(ind-1));
0354 end
0355 end
0356
0357
0358
0359 function model = readSimPhenyGprCmpd(baseName,model)
0360
0361 [rxnInfo,rxns,allGenes] = readSimPhenyGPR([baseName '_gpr.txt']);
0362
0363 nRxns = length(model.rxns);
0364
0365
0366 rxnGeneMat = sparse(nRxns,length(allGenes));
0367 h = waitbar(0,'Constructing GPR mapping ...');
0368 for i = 1:nRxns
0369 rxnID = find(ismember(rxns,model.rxns{i}));
0370 if (~isempty(rxnID))
0371 if mod(i,10) == 0
0372 waitbar(i/nRxns,h);
0373 end
0374 [tmp,geneInd] = ismember(rxnInfo(rxnID).genes,allGenes);
0375 rxnGeneMat(i,geneInd) = 1;
0376 rules{i} = rxnInfo(rxnID).rule;
0377 grRules{i} = rxnInfo(rxnID).gra;
0378 grRules{i} = regexprep(grRules{i},'\s{2,}',' ');
0379 grRules{i} = regexprep(grRules{i},'( ','(');
0380 grRules{i} = regexprep(grRules{i},' )',')');
0381 subSystems{i} = rxnInfo(rxnID).subSystem;
0382 for j = 1:length(geneInd)
0383
0384 rules{i} = strrep(rules{i},['x(' num2str(j) ')'],['x(' num2str(geneInd(j)) '_TMP_)']);
0385 end
0386 rules{i} = strrep(rules{i},'_TMP_','');
0387 else
0388 rules{i} = '';
0389 grRules{i} = '';
0390 subSystems{i} = '';
0391 end
0392 end
0393 if ( regexp( version, 'R20') )
0394 close(h);
0395 end
0396
0397
0398 [metInfo,mets] = readSimPhenyCMPD([baseName '_cmpd.txt']);
0399
0400 baseMets = parseMetNames(model.mets);
0401 nMets = length(model.mets);
0402 h = waitbar(0,'Constructing metabolite lists ...');
0403 for i = 1:nMets
0404 if mod(i,10) == 0
0405 waitbar(i/nMets,h);
0406 end
0407 metID = find(ismember(mets,baseMets{i}));
0408 if (~isempty(metID))
0409 metFormulas{i} = metInfo(metID).formula;
0410 else
0411 metFormulas{i} = '';
0412 end
0413 end
0414 if ( regexp( version, 'R20') )
0415 close(h);
0416 end
0417
0418 model.rxnGeneMat = rxnGeneMat;
0419 model.rules = columnVector(rules);
0420 model.grRules = columnVector(grRules);
0421 model.genes = columnVector(allGenes);
0422 model.metFormulas = columnVector(metFormulas);
0423 model.subSystems = columnVector(subSystems);
0424
0425