Home > RAVEN > contractModel.m

contractModel

PURPOSE ^

contractModel

SYNOPSIS ^

function [reducedModel removedRxns]=contractModel(model)

DESCRIPTION ^

 contractModel
   Contracts a model by grouping all identical reactions. Similar to the
   deleteDuplicates part in simplifyModel but more care is taken here
   when it comes to gene associations

    model           a model structure

    reducedModel    a model structure with grouped reactions
   removedRxns     cell array with the duplicate reactions

   NOTE: This code might not work for advanced grRules strings
         that involve nested expressions of 'and' and 'or'.

   Usage: [reducedModel removedRxns]=contractModel(model)

   Rasmus Agren, 2011-06-13

CROSS-REFERENCE INFORMATION ^

This function calls: This function is called by:

SOURCE CODE ^

0001 function [reducedModel removedRxns]=contractModel(model)
0002 % contractModel
0003 %   Contracts a model by grouping all identical reactions. Similar to the
0004 %   deleteDuplicates part in simplifyModel but more care is taken here
0005 %   when it comes to gene associations
0006 %
0007 %    model           a model structure
0008 %
0009 %    reducedModel    a model structure with grouped reactions
0010 %   removedRxns     cell array with the duplicate reactions
0011 %
0012 %   NOTE: This code might not work for advanced grRules strings
0013 %         that involve nested expressions of 'and' and 'or'.
0014 %
0015 %   Usage: [reducedModel removedRxns]=contractModel(model)
0016 %
0017 %   Rasmus Agren, 2011-06-13
0018 %
0019 
0020 %First sort the model so that reversible reactions are in the same
0021 %direction
0022 modelS=sortModel(model);
0023 
0024 %Get a list of duplicate reactions
0025 x=[modelS.S; model.rev']';
0026 [B,I,J] = unique(x,'rows','first');
0027 
0028 duplicateRxns=setdiff(1:numel(model.rxns),I);
0029 mergeTo=I(J(duplicateRxns));
0030 
0031 %Now add all the info from this one. Print a warning if they have different
0032 %bounds or objective function coefficients. Uses the widest bounds and largest
0033 %magnitude of objective coefficient
0034 for i=1:numel(duplicateRxns)
0035     if model.lb(duplicateRxns(i))<model.lb(mergeTo(i))
0036        fprintf(['WARNING: Duplicate reaction ' model.rxns{duplicateRxns(i)} ' has wider lower bound. Uses the most negative/smallest lower bound\n']); 
0037        model.lb(mergeTo(i))=model.lb(duplicateRxns(i));
0038     end
0039     if model.ub(duplicateRxns(i))>model.ub(mergeTo(i))
0040        fprintf(['WARNING: Duplicate reaction ' model.rxns{duplicateRxns(i)} ' has wider upper bound. Uses the most positive/largest upper bound\n']); 
0041        model.ub(mergeTo(i))=model.ub(duplicateRxns(i));
0042     end
0043     if abs(model.c(duplicateRxns(i)))>abs(model.c(mergeTo(i)))
0044        fprintf(['WARNING: Duplicate reaction ' model.rxns{duplicateRxns(i)} ' has a larger objective function coefficient. Uses the largest coefficient\n']); 
0045        model.c(mergeTo(i))=model.c(duplicateRxns(i));
0046     end
0047     
0048     %Genes are added as 'or'
0049     commonGenes=find(model.rxnGeneMat(duplicateRxns(i),:) & model.rxnGeneMat(mergeTo(i),:));
0050     newGenes=model.rxnGeneMat(duplicateRxns(i),:);
0051     newGenes(commonGenes)=0;
0052     if isfield(model,'rxnGeneMat')
0053        model.rxnGeneMat(mergeTo(i),:)=model.rxnGeneMat(mergeTo(i),:)+newGenes;
0054     end
0055     if isfield(model,'grRules')
0056         if any(model.grRules{duplicateRxns(i)})
0057            if any(model.grRules{mergeTo(i)})
0058                %Split both grStrings on ' or ' and then put together union
0059                %with ' or '
0060                rules1=regexp(model.grRules{mergeTo(i)},' or ','split');
0061                rules2=regexp(model.grRules{duplicateRxns(i)},' or ','split');
0062                allRules=union(rules1,rules2);
0063                
0064                %Probably not the nicest way to do this
0065                model.grRules{mergeTo(i)}=allRules{1};
0066                for j=2:numel(allRules)
0067                    model.grRules{mergeTo(i)}=[model.grRules{mergeTo(i)} ' or ' allRules{j}];
0068                end
0069            else
0070                model.grRules{mergeTo(i)}=model.grRules{duplicateRxns(i)};
0071            end
0072         end
0073     end
0074     if isfield(model,'eccodes')
0075         if any(model.eccodes{duplicateRxns(i)})
0076            if any(model.eccodes{mergeTo(i)})
0077                %Split on ';' and put together the union with ';'
0078                codes1=regexp(model.eccodes{mergeTo(i)},';','split');
0079                codes2=regexp(model.eccodes{duplicateRxns(i)},';','split');
0080                codes=union(codes1,codes2);
0081                model.eccodes{mergeTo(i)}=codes{1};
0082                for j=2:numel(codes)
0083                   model.eccodes{mergeTo(i)}=[model.eccodes{mergeTo(i)} ';' codes{j}]; 
0084                end
0085            else
0086                model.eccodes{mergeTo(i)}=model.eccodes{duplicateRxns(i)};
0087            end
0088         end
0089     end
0090 end
0091 
0092 %Delete the duplicate reactions
0093 reducedModel=removeRxns(model,duplicateRxns);
0094 removedRxns=model.rxns(duplicateRxns);
0095 end

Generated on Tue 23-Apr-2013 15:18:37 by m2html © 2005