Home > RAVEN > printModelStats.m

printModelStats

PURPOSE ^

printModelStats

SYNOPSIS ^

function printModelStats(model, printModelIssues, printDetails)

DESCRIPTION ^

 printModelStats
   prints some statistics about a model to the screen

   model               a model structure
   printModelIssues    true if information about unconnected
                       reactions/metabolites and elemental balancing
                       should be printed (opt, default false)
   printDetails        true if detailed information should be printed
                       about model issues. Only used if printModelIssues
                       is true (opt, default true)

   Usage: printModelStats(model,printModelIssues, printDetails)

   Rasmus Agren, 2013-08-01

CROSS-REFERENCE INFORMATION ^

This function calls: This function is called by:

SOURCE CODE ^

0001 function printModelStats(model, printModelIssues, printDetails)
0002 % printModelStats
0003 %   prints some statistics about a model to the screen
0004 %
0005 %   model               a model structure
0006 %   printModelIssues    true if information about unconnected
0007 %                       reactions/metabolites and elemental balancing
0008 %                       should be printed (opt, default false)
0009 %   printDetails        true if detailed information should be printed
0010 %                       about model issues. Only used if printModelIssues
0011 %                       is true (opt, default true)
0012 %
0013 %   Usage: printModelStats(model,printModelIssues, printDetails)
0014 %
0015 %   Rasmus Agren, 2013-08-01
0016 %
0017 
0018 if nargin<2
0019     printModelIssues=false;
0020 end
0021 if nargin<3
0022     printDetails=true;
0023 end
0024 
0025 fprintf(['Network statistics for ' model.id ': ' model.description '\n']);
0026 
0027 %Get which reactions are present in each compartment
0028 rxnComps=sparse(numel(model.rxns),numel(model.comps));
0029 
0030 %For each compartment, find the metabolites that are present in that
0031 %compartment and then the reactions they are involved in
0032 for i=1:numel(model.comps)
0033     [crap I]=find(model.S(model.metComps==i,:));
0034     rxnComps(I,i)=1;
0035 end
0036 
0037 if isfield(model,'eccodes')
0038     fprintf(['EC-numbers\t\t\t' num2str(numel(unique(model.eccodes))) '\n']);
0039 end
0040 
0041 %Print information about genes
0042 if isfield(model,'genes')
0043     fprintf(['Genes*\t\t\t\t' num2str(numel(model.genes)) '\n']);
0044     %Find the genes in each compartment
0045     for i=1:numel(model.comps)
0046         [crap I]=find(model.rxnGeneMat(rxnComps(:,i)==1,:));
0047         fprintf(['\t' model.compNames{i} '\t' num2str(numel(unique(I))) '\n']);
0048     end
0049 end
0050 
0051 %Print information about reactions
0052 fprintf(['\nReactions*\t\t\t' num2str(numel(model.rxns)) '\n']);
0053 for i=1:numel(model.comps)
0054     fprintf(['\t' model.compNames{i} '\t' num2str(sum(rxnComps(:,i))) '\n']);
0055 end
0056 
0057 %Removes the effect of compartments and removes duplicate reactions
0058 temp=model;
0059 temp.comps(:)={'s'}; %Set all compartments to be the same
0060 equ=constructEquations(sortModel(temp,true,true),temp.rxns,false);
0061 
0062 fprintf(['Unique reactions**\t' num2str(numel(unique(equ))) '\n']);
0063 
0064 %Print information about metabolites
0065 fprintf(['\nMetabolites\t\t\t' num2str(numel(model.mets)) '\n']);
0066 for i=1:numel(model.comps)
0067     fprintf(['\t' model.compNames{i} '\t' num2str(sum(model.metComps==i)) '\n']);
0068 end
0069 
0070 fprintf(['Unique metabolites\t' num2str(numel(unique(model.metNames))) '\n']);
0071 
0072 fprintf('\n* Genes and reactions are counted for each compartment if any of the corresponding metabolites are in that compartment. The sum may therefore not add up to the total number.\n');
0073 fprintf('** Unique reactions are defined as being biochemically unique (no compartmentalization)\n');
0074 
0075 %Also print some potential problems if there are any
0076 if printModelIssues==true
0077     fprintf(['\nShort model quality summary for ' model.id ': ' model.description '\n']);
0078     
0079     %Check that all the metabolites are being used
0080     involvedMat=model.S;
0081     involvedMat(involvedMat~=0)=1;
0082     usedMets=sum(involvedMat,2);
0083     notPresent=find(usedMets==0);
0084     if ~isempty(notPresent)
0085         errorText=['Non-used metabolites\t' num2str(numel(notPresent)) '\n'];
0086         if printDetails==true
0087             for i=1:numel(notPresent)
0088                 errorText=[errorText '\t(' model.mets{notPresent(i)} ') ' model.metNames{notPresent(i)} '\n'];
0089             end
0090             errorText=[errorText '\n'];
0091         end
0092         fprintf(errorText);
0093     end
0094     
0095     %Check if there are empty reactions
0096     usedRxns=sum(involvedMat,1);
0097     notUsed=find(usedRxns==0);
0098     if ~isempty(notUsed)
0099         errorText=['Empty reactions\t' num2str(numel(notUsed)) '\n'];
0100         if printDetails==true
0101             for i=1:numel(notUsed)
0102                 errorText=[errorText '\t' model.rxns{notUsed(i)} '\n'];
0103             end
0104             errorText=[errorText '\n'];
0105         end
0106         fprintf(errorText);
0107     end
0108     
0109     %Check if there are dead-end reactions/metabolites
0110     [reducedModel, deletedReactions, deletedMetabolites]=simplifyModel(model,true,false,false,true);
0111     
0112     if ~isempty(deletedReactions)
0113         errorText=['Dead-end reactions\t' num2str(numel(deletedReactions)) '\n'];
0114         if printDetails==true
0115             for i=1:numel(deletedReactions)
0116                 errorText=[errorText '\t' deletedReactions{i} '\n'];
0117             end
0118             errorText=[errorText '\n'];
0119         end
0120         fprintf(errorText);
0121     end
0122     
0123     %Ignore non-used metabolites
0124     deletedMetabolites=setdiff(deletedMetabolites,model.mets(notPresent));
0125     %Must map to indexes in order to print names
0126     deletedMetabolites=find(ismember(model.mets,deletedMetabolites));
0127     if ~isempty(deletedMetabolites)
0128         errorText=['Dead-end metabolites\t' num2str(numel(deletedMetabolites)) '\n'];
0129         if printDetails==true
0130             for i=1:numel(deletedMetabolites)
0131                 errorText=[errorText '\t(' model.mets{deletedMetabolites(i)} ') ' model.metNames{deletedMetabolites(i)} '\n'];
0132             end
0133             errorText=[errorText '\n'];
0134         end
0135         fprintf(errorText);
0136     end
0137     
0138     balanceStructure=getElementalBalance(model);
0139     
0140     notParsed=find(balanceStructure.balanceStatus<0);
0141     notBalanced=find(balanceStructure.balanceStatus==0);
0142     
0143     if ~isempty(notParsed)
0144         errorText=['Reactions which could not be elementally balanced\t' num2str(numel(notParsed)) '\n'];
0145         if printDetails==true
0146             for i=1:numel(notParsed)
0147                 errorText=[errorText '\t' model.rxns{notParsed(i)} '\n'];
0148             end
0149             errorText=[errorText '\n'];
0150         end
0151         fprintf(errorText);
0152     end
0153     if ~isempty(notBalanced)
0154         errorText=['Reactions which are elementally unbalanced\t' num2str(numel(notBalanced)) '\n'];
0155         if printDetails==true
0156             names=strcat(balanceStructure.elements.names,{', '});
0157             for i=1:numel(notBalanced)
0158                 badOnes=sprintf('%s', names{abs(balanceStructure.leftComp(notBalanced(i),:)-balanceStructure.rightComp(notBalanced(i),:))>10^-7});
0159                 errorText=[errorText '\t' model.rxns{notBalanced(i)} '\t' badOnes(1:end-2) '\n'];
0160             end
0161             errorText=[errorText '\n'];
0162         end
0163         fprintf(errorText);
0164     end   
0165 end
0166 end

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