convertToIrrev

 convertToIrrev
   Converts a model to irreversible form

 Input:
   model         a model structure
   rxns          cell array with the reactions so split (if reversible)
                 (optional, default model.rxns)

 Output:
   irrevModel    a model structure where reversible reactions have
                 been split into one forward and one reverse reaction
   matchRev      matching forward reaction to its backward reaction
   rev2irrev     forward and backward reactions for reversible reactions
   irrev2rev     matching all reactions back to original model

   The reverse reactions are saved as 'rxnID_REV'. A warning is shown if
   some reaction identifiers already end with '_REV'.

 Usage: [irrevModel,matchRev,rev2irrev,irrev2rev]=convertToIrrev(model,rxns)

0001 function [irrevModel,matchRev,rev2irrev,irrev2rev]=convertToIrrev(model,rxns)
0002 % convertToIrrev
0003 %   Converts a model to irreversible form
0004 %
0005 % Input:
0006 %   model         a model structure
0007 %   rxns          cell array with the reactions so split (if reversible)
0008 %                 (optional, default model.rxns)
0009 %
0010 % Output:
0011 %   irrevModel    a model structure where reversible reactions have
0012 %                 been split into one forward and one reverse reaction
0013 %   matchRev      matching forward reaction to its backward reaction
0014 %   rev2irrev     forward and backward reactions for reversible reactions
0015 %   irrev2rev     matching all reactions back to original model
0016 %
0017 %   The reverse reactions are saved as 'rxnID_REV'. A warning is shown if
0018 %   some reaction identifiers already end with '_REV'.
0019 %
0020 % Usage: [irrevModel,matchRev,rev2irrev,irrev2rev]=convertToIrrev(model,rxns)
0021 
0022 if nargin<2
0023     I=true(numel(model.rxns),1);
0024 else
0025     rxns=convertCharArray(rxns);
0026     I=getIndexes(model,rxns,'rxns',true);
0027 end
0028 
0029 irrevModel=model;
0030 
0031 revIndexesBool=model.rev~=0 & I;
0032 revIndexes=find(revIndexesBool);
0033 numOrigRxns=numel(model.rxns);
0034 numRevRxns=numel(revIndexes);
0035 
0036 if any(revIndexesBool)
0037     irrevModel.S=[model.S,model.S(:,revIndexes)*-1];
0038     irrevModel.rev(revIndexes)=0;
0039     irrevModel.rev=[irrevModel.rev;zeros(numRevRxns,1)];
0040     
0041     %Get the limits for all normal and reversible rxns
0042     ubNormal=irrevModel.ub;
0043     ubNormal(revIndexes(ubNormal(revIndexes)<0))=0;
0044     lbNormal=irrevModel.lb;
0045     lbNormal(revIndexes(lbNormal(revIndexes)<0))=0;
0046     ubRev=irrevModel.lb(revIndexes)*-1;
0047     ubRev(ubRev<0)=0;
0048     lbRev=irrevModel.ub(revIndexes)*-1;
0049     lbRev(lbRev<0)=0;
0050     irrevModel.ub=[ubNormal;ubRev];
0051     irrevModel.lb=[lbNormal;lbRev];
0052     
0053     %The objective coefficents should be zero for the backwards reversible
0054     %reactions unless they were negative in the original. In that case they
0055     %should be positive for the backwards reversible and deleted from the
0056     %original
0057     irrevC=zeros(numRevRxns,1);
0058     
0059     if any(irrevModel.c(revIndexes)<0)
0060         originalC=irrevModel.c(revIndexes);
0061         irrevC(irrevModel.c(revIndexes)<0)=originalC(originalC<0)*-1;
0062         irrevModel.c(irrevModel.c<0 & revIndexesBool)=0;
0063     end
0064     irrevModel.c=[irrevModel.c;irrevC];
0065     
0066     irrevModel.rxns=[irrevModel.rxns;strcat(irrevModel.rxns(revIndexes),'_REV')];
0067     irrevModel.rxnNames=[irrevModel.rxnNames;strcat(irrevModel.rxnNames(revIndexes),' (reversible)')];
0068     
0069     if isfield(irrevModel,'grRules')
0070         irrevModel.grRules=[irrevModel.grRules;irrevModel.grRules(revIndexes,:)];
0071     end
0072     if isfield(irrevModel,'rxnMiriams')
0073         irrevModel.rxnMiriams=[irrevModel.rxnMiriams;irrevModel.rxnMiriams(revIndexes,:)];
0074     end
0075     if isfield(irrevModel,'rxnGeneMat')
0076         irrevModel.rxnGeneMat=[irrevModel.rxnGeneMat;irrevModel.rxnGeneMat(revIndexes,:)];
0077     end
0078     if isfield(irrevModel,'subSystems')
0079         irrevModel.subSystems=[irrevModel.subSystems;irrevModel.subSystems(revIndexes)];
0080     end
0081     if isfield(irrevModel,'eccodes')
0082         irrevModel.eccodes=[irrevModel.eccodes;irrevModel.eccodes(revIndexes)];
0083     end
0084     if isfield(irrevModel,'rxnComps')
0085         irrevModel.rxnComps=[irrevModel.rxnComps;irrevModel.rxnComps(revIndexes)];
0086     end
0087     if isfield(irrevModel,'rxnFrom')
0088         irrevModel.rxnFrom=[irrevModel.rxnFrom;irrevModel.rxnFrom(revIndexes)];
0089     end
0090     if isfield(irrevModel,'rxnScores')
0091         irrevModel.rxnScores=[irrevModel.rxnScores;irrevModel.rxnScores(revIndexes)];
0092     end
0093     if isfield(irrevModel,'rxnNotes')
0094         irrevModel.rxnNotes=[irrevModel.rxnNotes;irrevModel.rxnNotes(revIndexes)];
0095     end
0096     if isfield(irrevModel,'rxnConfidenceScores')
0097         irrevModel.rxnConfidenceScores=[irrevModel.rxnConfidenceScores;irrevModel.rxnConfidenceScores(revIndexes)];
0098     end
0099     if isfield(irrevModel,'rxnDeltaG')
0100         irrevModel.rxnDeltaG=[irrevModel.rxnDeltaG;-irrevModel.rxnDeltaG(revIndexes)]; % Invert dG for reversed rxns
0101     end
0102     if isfield(irrevModel,'rxnReferences')
0103         irrevModel.rxnReferences=[irrevModel.rxnReferences;irrevModel.rxnReferences(revIndexes)];
0104     end
0105 end
0106     % Additional output
0107     if nargout>1
0108         irrev2rev = [transpose(1:numOrigRxns);revIndexes];
0109         rev2irrev = num2cell(transpose(1:numOrigRxns));
0110         newIdxs   = [revIndexes transpose(1:numRevRxns)];
0111         for i=1:numRevRxns
0112             rev2irrev{revIndexes(i)} = newIdxs(i,:);
0113         end
0114         matchRev  = zeros(numel(irrev2rev),1);
0115         matchRev(revIndexes) = (1:numRevRxns)+numOrigRxns;
0116     end
0117 end