Home > RAVEN > guessComposition.m

guessComposition

PURPOSE ^

guessComposition

SYNOPSIS ^

function [model, guessedFor, couldNotGuess]=guessComposition(model, printResults)

DESCRIPTION ^

 guessComposition
   Attempts to guess the composition of metabolites without information
   about elemental composition

   model               a model structure
   printResults        true if the output should be printed (opt, default true)

   model               a model structure with information about elemental
                       composition added
   guessedFor          indexes for the metabolites for which a composition
                       could be guessed
   couldNotGuess       indexes for the metabolites for which no
                       composition could be assigned

   This function works in a rather straight forward manner:

   1. Get the metabolites which lack composition and participates in
   at least one reaction where all other metabolites have composition information
   2. Loop through them and calculate their composition based on the rest
   of the involved metabolites. If there are any inconsistencies, so that
   a given metabolite should have different composition in different
   equations, then throw an error
   3. Go to 1

   This simple approach requires that the rest of the metabolites have
   correct composition information, and that the involved reactions are
   correct. The function will exit with an error on any inconsistencies,
   which means that it could also be used as a way of checking the model
   for errors. Note that just because this exits sucessfully, the
   calculated compositions could still be wrong (in case that the existing
   compositions were wrong)

   Usage: [newModel, guessedFor, couldNotGuess]=guessComposition(model, printResults)

   Rasmus Agren, 2013-11-06

CROSS-REFERENCE INFORMATION ^

This function calls: This function is called by:

SUBFUNCTIONS ^

SOURCE CODE ^

0001 function [model, guessedFor, couldNotGuess]=guessComposition(model, printResults)
0002 % guessComposition
0003 %   Attempts to guess the composition of metabolites without information
0004 %   about elemental composition
0005 %
0006 %   model               a model structure
0007 %   printResults        true if the output should be printed (opt, default true)
0008 %
0009 %   model               a model structure with information about elemental
0010 %                       composition added
0011 %   guessedFor          indexes for the metabolites for which a composition
0012 %                       could be guessed
0013 %   couldNotGuess       indexes for the metabolites for which no
0014 %                       composition could be assigned
0015 %
0016 %   This function works in a rather straight forward manner:
0017 %
0018 %   1. Get the metabolites which lack composition and participates in
0019 %   at least one reaction where all other metabolites have composition information
0020 %   2. Loop through them and calculate their composition based on the rest
0021 %   of the involved metabolites. If there are any inconsistencies, so that
0022 %   a given metabolite should have different composition in different
0023 %   equations, then throw an error
0024 %   3. Go to 1
0025 %
0026 %   This simple approach requires that the rest of the metabolites have
0027 %   correct composition information, and that the involved reactions are
0028 %   correct. The function will exit with an error on any inconsistencies,
0029 %   which means that it could also be used as a way of checking the model
0030 %   for errors. Note that just because this exits sucessfully, the
0031 %   calculated compositions could still be wrong (in case that the existing
0032 %   compositions were wrong)
0033 %
0034 %   Usage: [newModel, guessedFor, couldNotGuess]=guessComposition(model, printResults)
0035 %
0036 %   Rasmus Agren, 2013-11-06
0037 %
0038 
0039 if nargin<2
0040     printResults=true;
0041 end
0042 
0043 %The metabolites for which there is no elemental composition
0044 originalMissing=unique(model.metNames(cellfun(@isempty,model.metFormulas)));
0045 guessedFor={};
0046 
0047 %This is to keep track of if new metabolite compositions were predicted
0048 predicted=true;
0049 while predicted==true
0050     predicted=false;
0051 
0052     %Get the unique names (composition should be independent of compartment)
0053     %for the metabolites without composition
0054     metNames=unique(model.metNames(cellfun(@isempty,model.metFormulas)));
0055 
0056     %Parse the formulas in the model
0057     [elements, useMat, exitFlag]=parseFormulas(model.metFormulas, true,false);
0058 
0059     for i=1:numel(metNames)
0060         %Get the metabolites with this name. Not so neat, but this is a
0061         %fast function anyways
0062         mets=find(ismember(model.metNames,metNames(i)));
0063         
0064         currentComp=[];
0065 
0066         %Loop through the metabolites
0067         %-1: Could not assign due to missing info. -2: Could not assign due to contradiction
0068         %1: Composition assigned
0069         metStatus=-1;
0070         for j=1:numel(mets)
0071             %Get the reactions that the metabolite participates in
0072             [crap, I]=find(model.S(mets(j),:));
0073             if any(I)
0074                 for k=1:numel(I)
0075                     %Loop through the reactions and check if all other mets in them
0076                     %have known composition
0077                     eqn=model.S(:,I(k));
0078                     eqn(mets(j))=0;
0079                     if all(exitFlag(eqn~=0)==1)
0080                         %This means that all other mets had composition. Calculate
0081                         %the resulting composition for the unknown one
0082                         comp=useMat'*eqn;
0083                         
0084                         %This can result in round off errors if there are
0085                         %stoichiometries with many decimals. Ignore values
0086                         %below 10^-12
0087                         comp(abs(comp)<10^-12)=0;
0088 
0089                         %Check if the composition consist of both negative and
0090                         %positive values. If so, throw an error
0091                         if all(comp<=0) || all(comp>=0)
0092                             comp=abs(comp);
0093                             if isempty(currentComp)
0094                                 currentComp=comp;
0095                             end
0096                             %Also to deal with round off errors
0097                             if all(abs(currentComp-comp)<10^-10)
0098                                 metStatus=1;
0099                             else
0100                                 metStatus=-2;
0101                                 break;
0102 
0103                                 %%Check if there is an inconcistency
0104                                 %if any(currentComp~=comp)
0105                                 %    dispEM(['Could not predict composition of ' model.metNames{mets(i)} ],false);
0106                                 %end
0107                             end
0108                         else
0109                             %Check if there is an inconcistency
0110                             %if any(currentComp~=comp)
0111                             %    dispEM(['Could not predict composition of ' model.metNames{loopThrough(i)} ],false);
0112                             %end
0113                             metStatus=-2;
0114                             break;
0115                         end
0116                     end
0117                 end
0118                 %If there was contradictions in one compartment, then abort for
0119                 %all compartments
0120                 if metStatus==-2
0121                     break;
0122                 end
0123             else
0124                 %The metabolite doesn't participate, no composition can be
0125                 %calculated
0126                 metStatus=-1;
0127             end
0128         end
0129         %Check status of the metabolite
0130         switch metStatus
0131             case -2
0132                 dispEM(['Could not predict composition for "' metNames{i} '" due to inconsistencies'],false);
0133             case 1
0134                 %Calculate and add the composition
0135                 str=getCompString(elements,comp);
0136                 model.metFormulas(mets)={str};
0137                 if printResults==true
0138                     fprintf(['Predicted composition for "' metNames{i} '" to be ' str '\n']);
0139                 end
0140                 
0141                 %Keep track
0142                 guessedFor=[guessedFor;metNames(i)];
0143                 
0144                 predicted=true; %To loop again
0145         end        
0146     end
0147 end
0148 
0149 couldNotGuess=setdiff(originalMissing,guessedFor);
0150 end
0151 
0152 %Helper function for getting the composition string
0153 function str=getCompString(elements,comp)
0154     str='';
0155     
0156     for i=1:numel(comp)
0157        if comp(i)~=0
0158           if comp(i)==1
0159              str=[str  elements.abbrevs{i}];
0160           else
0161              str=[str  elements.abbrevs{i} num2str(comp(i))]; 
0162           end
0163        end
0164     end
0165 end

Generated on Mon 06-Jan-2014 14:58:12 by m2html © 2005