Home > io > readYAMLmodel.m

readYAMLmodel

PURPOSE ^

readYAMLmodel

SYNOPSIS ^

function model=readYAMLmodel(fileName, verbose)

DESCRIPTION ^

 readYAMLmodel
   Reads a yaml file matching (roughly) the cobrapy yaml structure

   Input:
   fileName    a model file in yaml file format. A dialog window will open
               if no file name is specified.
   verbose     set as true to monitor progress (optional, default false)

   Output:
   model       a model structure

 Usage: model = readYAMLmodel(fileName, verbose)

CROSS-REFERENCE INFORMATION ^

This function calls: This function is called by:

SUBFUNCTIONS ^

SOURCE CODE ^

0001 function model=readYAMLmodel(fileName, verbose)
0002 % readYAMLmodel
0003 %   Reads a yaml file matching (roughly) the cobrapy yaml structure
0004 %
0005 %   Input:
0006 %   fileName    a model file in yaml file format. A dialog window will open
0007 %               if no file name is specified.
0008 %   verbose     set as true to monitor progress (optional, default false)
0009 %
0010 %   Output:
0011 %   model       a model structure
0012 %
0013 % Usage: model = readYAMLmodel(fileName, verbose)
0014 if nargin<1 || isempty(fileName)
0015     [fileName, pathName] = uigetfile({'*.yml;*.yaml'}, 'Please select the model file');
0016     if fileName == 0
0017         error('You should select a model file')
0018     else
0019         fileName = fullfile(pathName,fileName);
0020     end
0021 end
0022 if nargin < 2
0023     verbose = false;
0024 end
0025 
0026 if ~isfile(fileName)
0027     error('Yaml file %s cannot be found', string(fileName));
0028 end
0029 
0030 if verLessThan('matlab','9.9') %readlines introduced 2020b
0031     fid=fopen(fileName);
0032     line_raw=cell(1000000,1);
0033     i=1;
0034     while ~feof(fid)
0035         line_raw{i}=fgetl(fid);
0036         i=i+1;
0037     end
0038     line_raw(i:end)=[];
0039     line_raw=string(line_raw);
0040     fclose(fid);
0041 else
0042     line_raw=readlines(fileName);
0043 end
0044 % If entry is broken of multiple lines, concatenate. Assumes at least 6
0045 % 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)
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},' ');
0052 end
0053 line_raw(brokenLine)=[];
0054 
0055 line_key = regexprep(line_raw,'^ *-? ([^:]+)(:)($| .*)','$1');
0056 line_key = regexprep(line_key,'(.*!!omap)|(---)|( {4,}.*)','');
0057 
0058 line_value = regexprep(line_raw, '.*:$','');
0059 line_value = regexprep(line_value, '[^":]+: "?(.+)"?$','$1');
0060 line_value = regexprep(line_value, '(")|(^ {4,}- )','');
0061 
0062 line_value(strcmp(line_value,line_raw)) = {''};
0063 
0064 modelFields =   {'id',char();...
0065                'name',char();...
0066         'description',char();...
0067             'version',char();...
0068                'date',char();...
0069          'annotation',struct();...
0070                'rxns',{};...
0071            'rxnNames',{};...
0072                'mets',{};...
0073            'metNames',{};...
0074                   'S',sparse([]);...
0075                  'lb',{};... %Changed to double in the end.
0076                  'ub',{};... %Changed to double in the end.
0077                 'rev',{};... %Changed to double in the end.
0078                   'c',[];...
0079                   'b',cell(0,0);... %Changed to double in the end.
0080               'genes',cell(0,0);...
0081             'grRules',cell(0,0);...
0082          'rxnGeneMat',sparse([]);...
0083            'rxnComps',cell(0,0);... %Changed to double in the end.
0084          'subSystems',cell(0,0);...
0085             'eccodes',cell(0,0);...
0086          'rxnMiriams',cell(0,0);...
0087           'rxnDeltaG',{};... %Changed to double in the end.
0088            'rxnNotes',cell(0,0);...
0089       'rxnReferences',cell(0,0);...
0090 'rxnConfidenceScores',cell(0,0);...
0091            'metComps',cell(0,0);... %Changed to double in the end.
0092              'inchis',cell(0,0);...
0093           'metSmiles',cell(0,0);...
0094         'metFormulas',cell(0,0);...
0095          'metMiriams',cell(0,0);...
0096           'metDeltaG',{};... %Changed to double in the end.
0097          'metCharges',cell(0,0);... %Changed to double in the end.
0098            'metNotes',cell(0,0);...
0099               'comps',cell(0,0);...
0100           'compNames',cell(0,0);...
0101         'compOutside',cell(0,0);...
0102           'geneComps',cell(0,0);... %Changed to double in the end.
0103         'geneMiriams',cell(0,0);...
0104      'geneShortNames',cell(0,0);...
0105        'proteins',cell(0,0);...
0106       'unconstrained',cell(0,0);... %Changed to double in the end.
0107             'metFrom',cell(0,0);...
0108             'rxnFrom',cell(0,0)};
0109 for i=1:size(modelFields,1)
0110     model.(modelFields{i,1})=modelFields{i,2};
0111 end
0112 
0113 % If GECKO model
0114 if any(contains(line_key,'geckoLight'))
0115     isGECKO=true;
0116     ecFields = {'geckoLight', false;...
0117                       'rxns', {};...
0118                       'kcat', {};...
0119                     'source', cell(0,0);...
0120                      'notes', cell(0,0);...
0121                    'eccodes', cell(0,0);...
0122                      'genes', cell(0,0);...
0123                    'enzymes', cell(0,0);...
0124                         'mw', cell(0,0);...
0125                   'sequence', cell(0,0);...
0126                      'concs', cell(0,0);...
0127                  'rxnEnzMat', []};
0128     for i=1:size(ecFields,1)
0129         model.ec.(ecFields{i,1})=ecFields{i,2};
0130     end
0131     ecGecko=cell(25000,2);      ecGeckoNo=1;
0132     enzStoich=cell(100000,3);   enzStoichNo=1;
0133 else
0134     isGECKO=false;
0135 end
0136 
0137 section = 0;
0138 metMiriams=cell(100000,3);   metMirNo=1;
0139 rxnMiriams=cell(100000,3);   rxnMirNo=1;
0140 geneMiriams=cell(100000,3);  genMirNo=1;
0141 subSystems=cell(100000,2);   subSysNo=1;
0142 eccodes=cell(100000,2);      ecCodeNo=1;
0143 equations=cell(100000,3);   equatiNo=1;
0144 
0145 for i=1:numel(line_key)
0146     tline_raw = line_raw{i};
0147     tline_key = line_key{i};
0148     tline_value = line_value{i};
0149     % import different sections
0150     switch tline_raw
0151         case '- metaData:'
0152             section = 1;
0153             if verbose
0154                 fprintf('\t%d\n', section);
0155             end
0156             continue % Go to next line
0157         case '- metabolites:'
0158             section = 2;
0159             if verbose
0160                 fprintf('\t%d\n', section);
0161             end
0162             pos=0;
0163             continue
0164         case '- reactions:'
0165             section = 3;
0166             if verbose
0167                 fprintf('\t%d\n', section);
0168             end
0169             pos=0;
0170             continue
0171         case '- genes:'
0172             section = 4;
0173             if verbose
0174                 fprintf('\t%d\n', section);
0175             end
0176             pos=0;
0177             continue
0178         case '- compartments: !!omap'
0179             section = 5;
0180             if verbose
0181                 fprintf('\t%d\n', section);
0182             end
0183             pos=0;
0184             continue
0185         case '- ec-rxns:'
0186             section = 6;
0187             if verbose
0188                 fprintf('\t%d\n', section);
0189             end
0190             pos=0;
0191             continue
0192         case '- ec-enzymes:'
0193             section = 7;
0194             if verbose
0195                 fprintf('\t%d\n', section);
0196             end
0197             pos=0;
0198             continue
0199     end
0200 
0201     % skip over empty keys
0202     if isempty(tline_raw) || (isempty(tline_key) && contains(tline_raw,'!!omap'))
0203         continue;
0204     end
0205     
0206     % import metaData
0207     if section == 1
0208         switch tline_key
0209             case {'short_name','id'} %short_name used by human-GEM
0210                 model.id = tline_value;
0211             case 'name'
0212                 model.name = tline_value;                
0213             case 'full_name' %used by human-GEM
0214                 model.description = tline_value;
0215             case 'version'
0216                 model.version = tline_value;
0217             case 'date'
0218                 model.date = tline_value;                
0219             case 'taxonomy'
0220                 model.annotation.taxonomy = tline_value;
0221             case {'description','note'} %description used by human-GEM
0222                 model.annotation.note = tline_value;
0223             case 'github'
0224                 model.annotation.sourceUrl = tline_value;
0225             case 'sourceUrl'
0226                 model.annotation.sourceUrl = tline_value;                
0227             case 'givenName'
0228                 model.annotation.givenName = tline_value;
0229             case 'familyName'
0230                 model.annotation.familyName = tline_value;
0231             case 'authors'
0232                 model.annotation.authorList = tline_value;
0233             case 'email'
0234                 model.annotation.email = tline_value;
0235             case 'organization'
0236                 model.annotation.organization = tline_value;
0237             case 'geckoLight'
0238                 if strcmp(tline_value,'true')
0239                     model.ec.geckoLight = true;
0240                 end
0241         end; continue
0242     end
0243 
0244     % import metabolites:
0245     if section == 2
0246         switch tline_key
0247             case 'id'
0248                 pos = pos + 1;
0249                 model = readFieldValue(model, 'mets', tline_value,pos);
0250                 readList=''; miriamKey='';
0251             case 'name'
0252                 model = readFieldValue(model, 'metNames', tline_value, pos);
0253                 readList=''; miriamKey='';
0254             case 'compartment'
0255                 model = readFieldValue(model, 'metComps', tline_value, pos);
0256                 readList=''; miriamKey='';
0257             case 'formula'
0258                 model = readFieldValue(model, 'metFormulas', tline_value, pos);
0259                 readList=''; miriamKey='';
0260             case 'charge'
0261                 model = readFieldValue(model, 'metCharges', tline_value, pos);
0262                 readList=''; miriamKey='';
0263             case 'notes'
0264                 model = readFieldValue(model, 'metNotes', tline_value, pos);
0265                 readList=''; miriamKey='';
0266             case 'inchis'
0267                 model = readFieldValue(model, 'inchis', tline_value, pos);
0268                 readList=''; miriamKey='';
0269             case 'smiles'
0270                 model = readFieldValue(model, 'metSmiles', tline_value, pos);
0271                 readList=''; miriamKey='';    
0272             case 'deltaG'
0273                 model = readFieldValue(model, 'metDeltaG', tline_value, pos);
0274                 readList=''; miriamKey='';                                
0275             case 'metFrom'
0276                 model = readFieldValue(model, 'metFrom', tline_value, pos);
0277                 readList=''; miriamKey='';
0278             case 'annotation'
0279                 readList = 'annotation';
0280             otherwise
0281                 switch readList
0282                     case 'annotation'
0283                         [metMiriams, miriamKey, metMirNo] = gatherAnnotation(pos,metMiriams,tline_key,tline_value,miriamKey,metMirNo);
0284                     otherwise
0285                         error(['Unknown entry in yaml file: ' tline_raw])
0286                 end
0287         end; continue
0288     end
0289 
0290     % import reactions:
0291     if section == 3
0292         switch tline_key
0293             case 'id'
0294                 pos = pos + 1;
0295                 model = readFieldValue(model, 'rxns', tline_value,pos);
0296                 readList=''; miriamKey='';
0297             case 'name'
0298                 model = readFieldValue(model, 'rxnNames', tline_value, pos);
0299                 readList=''; miriamKey='';
0300             case 'lower_bound'
0301                 model.lb(pos,1) = {tline_value};
0302                 readList=''; miriamKey='';
0303             case 'upper_bound'
0304                 model.ub(pos,1) = {tline_value};
0305                 readList=''; miriamKey='';
0306             case 'rev'
0307                 model.rev(pos,1) = {tline_value};
0308                 readList=''; miriamKey='';
0309             case 'gene_reaction_rule'
0310                 model = readFieldValue(model, 'grRules', tline_value, pos);
0311                 readList=''; miriamKey='';
0312             case 'rxnNotes'
0313                 model = readFieldValue(model, 'rxnNotes', tline_value, pos);
0314                 readList=''; miriamKey='';
0315             case 'rxnFrom'
0316                 model = readFieldValue(model, 'rxnFrom', tline_value, pos);
0317                 readList=''; miriamKey='';
0318             case 'deltaG'
0319                 model = readFieldValue(model, 'rxnDeltaG', tline_value, pos);
0320                 readList=''; miriamKey='';                
0321             case 'objective_coefficient'
0322                 model.c(pos,1) = 1;
0323                 readList=''; miriamKey='';
0324             case 'references'
0325                 model = readFieldValue(model, 'rxnReferences', tline_value, pos);
0326                 readList=''; miriamKey='';
0327             case 'confidence_score'
0328                 model = readFieldValue(model, 'rxnConfidenceScores', tline_value, pos);
0329                 readList=''; miriamKey='';
0330             case 'eccodes'
0331                 if isempty(tline_value)
0332                     readList = 'eccodes';
0333                 else
0334                     eccodes(ecCodeNo,1:2)={pos,tline_value};
0335                     ecCodeNo=ecCodeNo+1;
0336                 end
0337             case 'subsystem'
0338                 if isempty(tline_value)
0339                     readList = 'subsystem';
0340                 else
0341                     subSystems(subSysNo,1:2)={pos,tline_value};
0342                     subSysNo=subSysNo+1;
0343                 end
0344             case 'metabolites'
0345                 readList = 'equation';
0346             case 'annotation'
0347                 readList = 'annotation';
0348                 
0349             otherwise
0350                 switch readList
0351                     case 'eccodes'
0352                         eccodes(ecCodeNo,1:2)={pos,regexprep(tline_value,'^ +- "?(.*)"?$','$1')};
0353                         ecCodeNo=ecCodeNo+1;
0354                     case 'subsystem'
0355                         subSystems(subSysNo,1:2)={pos,regexprep(tline_value,'^ +- "?(.*)"?$','$1')};
0356                         subSysNo=subSysNo+1;
0357                     case 'annotation'
0358                         [rxnMiriams, miriamKey,rxnMirNo] = gatherAnnotation(pos,rxnMiriams,tline_key,tline_value,miriamKey,rxnMirNo);
0359                     case 'equation'
0360                         coeff = sscanf(tline_value,'%f');
0361                         equations(equatiNo,1:3)={pos,tline_key,coeff};
0362                         equatiNo=equatiNo+1;
0363                     otherwise
0364                         error(['Unknown entry in yaml file: ' tline_raw])
0365                 end                
0366         end; continue
0367     end
0368 
0369     % import genes:
0370     if section == 4
0371         switch tline_key
0372             case 'id'
0373                 pos = pos + 1;
0374                 model = readFieldValue(model, 'genes', tline_value, pos);
0375                 readList = '';
0376                 miriamKey = '';
0377             case 'name'
0378                 model = readFieldValue(model, 'geneShortNames', tline_value, pos);
0379             case 'protein'
0380                 model = readFieldValue(model, 'proteins', tline_value, pos);
0381             case 'annotation'
0382                 readList = 'annotation';
0383             otherwise
0384                 switch readList
0385                     case 'annotation'
0386                         [geneMiriams, miriamKey,genMirNo] = gatherAnnotation(pos,geneMiriams,tline_key,tline_value,miriamKey,genMirNo);
0387                     otherwise
0388                         error(['Unknown entry in yaml file: ' tline_raw])
0389                 end
0390         end; continue
0391     end
0392 
0393     % import compartments:
0394     if section == 5
0395         model.comps(end+1,1) = {tline_key};
0396         model.compNames(end+1,1) = {tline_value};
0397     end
0398 
0399     % import ec reaction info
0400     if section == 6
0401         switch tline_key
0402             case 'id'
0403                 pos = pos + 1;
0404                 model.ec = readFieldValue(model.ec, 'rxns', tline_value, pos);
0405                 readList='';
0406             case 'kcat'
0407                 model.ec = readFieldValue(model.ec, 'kcat', tline_value, pos);
0408                 readList='';
0409             case 'source'
0410                 model.ec = readFieldValue(model.ec, 'source', tline_value, pos);
0411                 readList='';
0412             case 'notes'
0413                 model.ec = readFieldValue(model.ec, 'notes', tline_value, pos);
0414                 readList='';
0415             case 'eccodes'
0416                 if isempty(tline_value)
0417                     readList = 'eccodes';
0418                 else
0419                     ecGecko(ecGeckoNo,1:2)={pos,tline_value};
0420                     ecGeckoNo=ecGeckoNo+1;
0421                 end
0422             case 'enzymes'
0423                 readList = 'enzStoich';
0424             otherwise
0425                 switch readList
0426                     case 'eccodes'
0427                         ecGecko(ecGeckoNo,1:2)={pos,regexprep(tline_value,'^ +- "?(.*)"?$','$1')};
0428                         ecGeckoNo=ecGeckoNo+1;
0429                     case 'enzStoich'
0430                         coeff = sscanf(tline_value,'%f');
0431                         enzStoich(enzStoichNo,1:3)={pos,tline_key,coeff};
0432                         enzStoichNo=enzStoichNo+1;
0433                     otherwise
0434                         error(['Unknown entry in yaml file: ' tline_raw])
0435                 end
0436         end; continue
0437     end
0438 
0439     % import ec enzyme info
0440     if section == 7
0441         switch tline_key
0442             case 'genes'
0443                 pos = pos + 1;
0444                 model.ec = readFieldValue(model.ec, 'genes', tline_value, pos);
0445             case 'enzymes'
0446                 model.ec = readFieldValue(model.ec, 'enzymes', tline_value, pos);
0447             case 'mw'
0448                 model.ec = readFieldValue(model.ec, 'mw', tline_value, pos);
0449             case 'sequence'
0450                 model.ec = readFieldValue(model.ec, 'sequence', tline_value, pos);
0451             case 'concs'
0452                 model.ec = readFieldValue(model.ec, 'concs', tline_value, pos);
0453             otherwise
0454                 error(['Unknown entry in yaml file: ' tline_raw])
0455         end; continue
0456     end
0457 end
0458 
0459 %Parse annotations
0460 if ~isempty(metMiriams)
0461     locs = cell2mat(metMiriams(:,1));
0462     for i=unique(locs)'
0463         model.metMiriams{i,1}.name=metMiriams(locs==i,2);
0464         model.metMiriams{i,1}.value=metMiriams(locs==i,3);
0465     end
0466 end
0467 if ~isempty(rxnMiriams)
0468     locs = cell2mat(rxnMiriams(:,1));
0469     for i=unique(locs)'
0470         model.rxnMiriams{i,1}.name=rxnMiriams(locs==i,2);
0471         model.rxnMiriams{i,1}.value=rxnMiriams(locs==i,3);
0472     end
0473 end
0474 if ~isempty(geneMiriams)
0475     locs = cell2mat(geneMiriams(:,1));
0476     for i=unique(locs)'
0477         model.geneMiriams{i,1}.name=geneMiriams(locs==i,2);
0478         model.geneMiriams{i,1}.value=geneMiriams(locs==i,3);
0479     end
0480 end
0481 
0482 %Parse subSystems
0483 if ~isempty(subSystems)
0484     locs = cell2mat(subSystems(:,1));
0485     for i=unique(locs)'
0486         model.subSystems{i,1}=subSystems(locs==i,2);
0487     end
0488 end
0489 
0490 %Parse ec-codes
0491 if ~isempty(eccodes)
0492     locs = cell2mat(eccodes(:,1));
0493     for i=unique(locs)'
0494         eccodesCat=strjoin(eccodes(locs==i,2),';');
0495         model.eccodes{i,1}=eccodesCat;
0496     end
0497     emptyEc=cellfun('isempty',model.eccodes);
0498     model.eccodes(emptyEc)={''};
0499 end
0500 
0501 % follow-up data processing
0502 if verbose
0503     fprintf('\nimporting completed\nfollow-up processing...');
0504 end
0505 [~, model.metComps] = ismember(model.metComps, model.comps);
0506 [~, model.geneComps] = ismember(model.geneComps, model.comps);
0507 [~, model.rxnComps] = ismember(model.rxnComps, model.comps);
0508 
0509 % Fill S-matrix
0510 rxnIdx = cellfun('isempty', equations(:,1));
0511 equations(rxnIdx,:) = '';
0512 rxnIdx = cell2mat(equations(:,1));
0513 [~,metIdx] = ismember(equations(:,2),model.mets);
0514 coeffs = cell2mat(equations(:,3));
0515 model.S=sparse(max(metIdx),max(rxnIdx));
0516 linearIndices = sub2ind([max(metIdx), max(rxnIdx)],metIdx,rxnIdx);
0517 model.S(linearIndices) = coeffs;
0518 
0519 % Convert strings to numeric
0520 model.metCharges = str2double(model.metCharges);
0521 model.lb = str2double(model.lb);
0522 model.ub = str2double(model.ub);
0523 model.rxnConfidenceScores = str2double(model.rxnConfidenceScores);
0524 model.b = zeros(length(model.mets),1);
0525 model.metDeltaG = str2double(model.metDeltaG);
0526 model.rxnDeltaG = str2double(model.rxnDeltaG);
0527 
0528 % Fill some other fields
0529 model.annotation.defaultLB = min(model.lb);
0530 model.annotation.defaultUB = max(model.ub);
0531 if numel(model.lb)<numel(model.rxns) %No LB reported = min
0532     model.lb(end+1:numel(model.rxns)-numel(model.lb),1) = double(model.annotation.defaultLB);
0533 end
0534 if numel(model.ub)<numel(model.rxns) %No UB reported = max
0535     model.ub(end+1:numel(model.rxns)-numel(model.ub),1) = double(model.annotation.defaultUB);
0536 end
0537 if ~all(cellfun('isempty',model.rev))
0538     model.rev = str2double(model.rev);
0539 else
0540     model.rev = [];
0541 end
0542 if numel(model.rev)<numel(model.rxns) %No rev reported, assume from LB and UB
0543     model.rev(end+1:numel(model.rxns)-numel(model.rev),1) = double(model.lb<0 & model.ub>0);
0544 end
0545 
0546 % Remove empty fields, otherwise fill to correct length
0547 % Reactions
0548 for i={'rxnNames','grRules','eccodes','rxnNotes','rxnReferences',...
0549        'rxnFrom','subSystems','rxnMiriams'} % Empty strings
0550    model = emptyOrFill(model,i{1},{''},'rxns');
0551 end
0552 for i={'c'} % Zeros
0553    model = emptyOrFill(model,i{1},0,'rxns',true);
0554 end
0555 for i={'rxnConfidenceScores','rxnDeltaG'} % NaNs
0556    model = emptyOrFill(model,i{1},NaN,'rxns');
0557 end
0558 for i={'rxnComps'} % Ones, assume first compartment
0559    model = emptyOrFill(model,i{1},1,'rxns');
0560 end
0561 % Metabolites
0562 for i={'metNames','inchis','metFormulas','metMiriams','metFrom','metSmiles','metNotes'} % Empty strings
0563    model = emptyOrFill(model,i{1},{''},'mets');
0564 end
0565 for i={'metCharges','unconstrained'} % Zeros
0566    model = emptyOrFill(model,i{1},0,'mets');
0567 end
0568 for i={'metDeltaG'} % % NaNs
0569     model = emptyOrFill(model,i{1},NaN,'mets');
0570  end
0571 for i={'metComps'} % Ones, assume first compartment
0572    model = emptyOrFill(model,i{1},1,'mets');
0573 end
0574 % Genes
0575 for i={'geneMiriams','geneShortNames','proteins'} % Empty strings
0576    model = emptyOrFill(model,i{1},{''},'genes');
0577 end
0578 for i={'geneComps'} % Ones, assume first compartment
0579    model = emptyOrFill(model,i{1},1,'genes');
0580 end
0581 % Comps
0582 for i={'compNames'} % Empty strings
0583    model = emptyOrFill(model,i{1},{''},'comps');
0584 end
0585 for i={'compOutside'} % First comp
0586    model = emptyOrFill(model,i{1},model.comps{1},'comps');
0587 end
0588 % Single fields are kept, even if empty
0589 % for i={'description','name','version','date','annotation'}
0590 %     if isempty(model.(i{1}))
0591 %         model = rmfield(model,i{1});
0592 %     end
0593 % end
0594 
0595 % Make rxnGeneMat fields and map to the existing model.genes field
0596 [genes, rxnGeneMat] = getGenesFromGrRules(model.grRules);
0597 model.rxnGeneMat = sparse(numel(model.rxns),numel(model.genes));
0598 [~,geneOrder] = ismember(genes,model.genes);
0599 if any(geneOrder == 0)
0600     error(['The grRules includes the following gene(s), that are not in '...
0601            'the list of model genes: ', genes{~geneOrder}])
0602 end
0603 model.rxnGeneMat(:,geneOrder) = rxnGeneMat;
0604 
0605 % Finalize GECKO model
0606 if isGECKO
0607     % Fill in empty fields and empty entries
0608     for i={'kcat','source','notes','eccodes'} % Even keep empty
0609         model.ec = emptyOrFill(model.ec,i{1},{''},'rxns',true);
0610     end
0611     for i={'enzymes','mw','sequence'}
0612         model.ec = emptyOrFill(model.ec,i{1},{''},'genes',true);
0613     end
0614     model.ec = emptyOrFill(model.ec,'concs',{'NaN'},'genes',true);
0615     model.ec = emptyOrFill(model.ec,'kcat',{'0'},'genes',true);
0616     % Change string to double
0617     for i={'kcat','mw','concs'}
0618         if isfield(model.ec,i{1})
0619             model.ec.(i{1}) = str2double(model.ec.(i{1}));
0620         end
0621     end
0622     % Fill rxnEnzMat
0623     rxnIdx              = cellfun('isempty', enzStoich(:,1));
0624     enzStoich(rxnIdx,:) = '';
0625     rxnIdx              = cell2mat(enzStoich(:,1));
0626     [~,enzIdx]          = ismember(enzStoich(:,2),model.ec.enzymes);
0627     coeffs              = cell2mat(enzStoich(:,3));
0628     model.ec.rxnEnzMat  = zeros(numel(model.ec.rxns), numel(model.ec.genes));
0629     linearIndices       = sub2ind([numel(model.ec.rxns), numel(model.ec.genes)], rxnIdx, enzIdx);
0630     model.ec.rxnEnzMat(linearIndices) = coeffs;
0631     %Parse ec-codes
0632     if ~isempty(ecGecko)
0633         locs = cell2mat(ecGecko(:,1));
0634         for i=unique(locs)'
0635             ecGeckoCat=strjoin(ecGecko(locs==i,2),';');
0636             model.ec.eccodes{i,1}=ecGeckoCat;
0637         end
0638         emptyEc=cellfun('isempty',model.ec.eccodes);
0639         model.ec.eccodes(emptyEc)={''};
0640     end
0641 end
0642 
0643 if verbose
0644     fprintf(' Done!\n');
0645 end
0646 end
0647 
0648 function model = emptyOrFill(model,field,emptyEntry,type,keepEmpty)
0649 if nargin<5
0650     keepEmpty=false;
0651 end
0652 if isnumeric(emptyEntry)
0653     emptyCells=isempty(model.(field));
0654 else
0655     emptyCells=cellfun('isempty',model.(field));
0656 end
0657 if all(emptyCells) && ~keepEmpty
0658     model = rmfield(model, field);
0659 elseif numel(model.(field))<numel(model.(type))
0660     model.(field)(end+1:numel(model.(type)),1)=emptyEntry;
0661 end
0662 end
0663 
0664 function model = readFieldValue(model, fieldName, value, pos)
0665 if numel(model.(fieldName))<pos-1
0666     model.(fieldName)(end+1:pos,1) = {''};
0667 end
0668 model.(fieldName)(pos,1) = {value};
0669 end
0670 
0671 function [miriams,miriamKey,entryNumber] = gatherAnnotation(pos,miriams,key,value,miriamKey,entryNumber)
0672 if isempty(key)
0673     key=miriamKey;
0674 else
0675     miriamKey=key;
0676 end
0677 if ~isempty(value)
0678     miriams(entryNumber,1:3) = {pos, key, strip(value)};
0679     entryNumber = entryNumber + 1;
0680 end
0681 end

Generated by m2html © 2005