Home > RAVEN > followChanged.m

followChanged

PURPOSE ^

followChanged

SYNOPSIS ^

function followChanged(model,fluxesA,fluxesB, cutOffChange, cutOffFlux, cutOffDiff, metaboliteList)

DESCRIPTION ^

 followChanged
    Prints fluxes and reactions for each of the reactions that results in
   different fluxes compared to the reference case.

   model           a model structure
   fluxesA         flux vector for the test case
   fluxesB         flux vector for the reference test
   cutOffChange    reactions where the fluxes differ by less than
                   this many percent won't be printed (opt, default 10^-8)
   cutOffFlux      reactions where the absolute value of both fluxes 
                   are below this value won't be printed (opt, 
                   default 10^-8)
   cutOffDiff      reactions where the fluxes differ by less than
                   cutOffDiff won't be printed (opt, default 10^-8)
   metaboliteList  cell array of metabolite names. Only reactions
                   involving any of these metabolites will be 
                   printed (opt)

   Usage: followChanged(model,fluxesA,fluxesB, cutOffChange, cutOffFlux,
           cutOffDiff, metaboliteList)

   Rasmus Agren, 2013-08-07

CROSS-REFERENCE INFORMATION ^

This function calls: This function is called by:

SOURCE CODE ^

0001 function followChanged(model,fluxesA,fluxesB, cutOffChange, cutOffFlux, cutOffDiff, metaboliteList)
0002 % followChanged
0003 %    Prints fluxes and reactions for each of the reactions that results in
0004 %   different fluxes compared to the reference case.
0005 %
0006 %   model           a model structure
0007 %   fluxesA         flux vector for the test case
0008 %   fluxesB         flux vector for the reference test
0009 %   cutOffChange    reactions where the fluxes differ by less than
0010 %                   this many percent won't be printed (opt, default 10^-8)
0011 %   cutOffFlux      reactions where the absolute value of both fluxes
0012 %                   are below this value won't be printed (opt,
0013 %                   default 10^-8)
0014 %   cutOffDiff      reactions where the fluxes differ by less than
0015 %                   cutOffDiff won't be printed (opt, default 10^-8)
0016 %   metaboliteList  cell array of metabolite names. Only reactions
0017 %                   involving any of these metabolites will be
0018 %                   printed (opt)
0019 %
0020 %   Usage: followChanged(model,fluxesA,fluxesB, cutOffChange, cutOffFlux,
0021 %           cutOffDiff, metaboliteList)
0022 %
0023 %   Rasmus Agren, 2013-08-07
0024 %
0025 
0026 %Checks if a cut off flux has been set
0027 if nargin<4
0028     cutOffChange=10^-8;
0029 end
0030 if nargin<5
0031     cutOffFlux=10^-8;
0032 end
0033 if nargin<6
0034     cutOffDiff=10^-8;
0035 end
0036 if nargin<7
0037     metaboliteList=[];
0038 end
0039 
0040 %If a metabolite list is to be used, then find all the reactions involving
0041 %any of those metabolites
0042 %Finds the metabolites
0043 if nargin>6
0044     reactionIndexes=[];
0045     for i=1:length(metaboliteList)
0046         metaboliteIndex=find(strcmpi(metaboliteList(i),model.metNames)); %Should use id maybe, setting
0047         if length(metaboliteIndex)>0
0048             [crap b]=find(model.S(metaboliteIndex,:));
0049             reactionIndexes=[reactionIndexes; b(:)];   
0050         else
0051             fprintf('Could not find any reactions with the metabolite %s\n\n',char(metaboliteList(i)))
0052         end
0053     end
0054     reactionIndexes=unique(reactionIndexes);
0055 else
0056     reactionIndexes=[1:length(fluxesA)]';
0057 end
0058 
0059 %Finds the reactions where either flux is at or above the cutOffFlux value
0060 in1=find(abs(fluxesA(reactionIndexes))>=cutOffFlux)';
0061 in2=find(abs(fluxesB(reactionIndexes))>=cutOffFlux)';
0062 ineither=reactionIndexes(unique([in1 in2]));
0063 
0064 %Keep only those solutions where the difference is larger than or equal to
0065 %cutOffDiff
0066 ineither=ineither(find(abs(fluxesA(ineither)-fluxesB(ineither))>=cutOffDiff));
0067 
0068 %Finds the reactions where the fluxes differ more than cutOffChange percent
0069 %First check those fluxes that are non-zero in solution1.x
0070 nonZeroFluxes=ineither(find(fluxesA(ineither)));
0071 quota=1+cutOffChange/100;
0072 larger=nonZeroFluxes(find((fluxesB(nonZeroFluxes)./fluxesA(nonZeroFluxes))>=(quota)))';
0073 smaller=nonZeroFluxes(find((fluxesB(nonZeroFluxes)./fluxesA(nonZeroFluxes))<(1/quota)))';
0074 fluxIndexes=[larger smaller];
0075 
0076 %Then add those where solution1 has a zero flux
0077 zeroFluxes=ineither(find(fluxesA(ineither)==0));
0078 fluxIndexes=unique([fluxIndexes zeroFluxes(find(abs(fluxesB(zeroFluxes))>=cutOffFlux))']);
0079 
0080 formulas=constructEquations(model,model.rxns(fluxIndexes));
0081 
0082 if nargin>4
0083     if nargin>5
0084         fprintf('These reactions have flux values that differ by more than %s percent, absolute values above %s, and a total difference above %s (%s reactions)\n\n',num2str(cutOffChange),num2str(cutOffFlux),num2str(cutOffDiff),num2str(length(formulas)));
0085     else
0086         fprintf('These reactions have flux values that differ by more than %s percent and absolute values above %s (%s reactions)\n\n',num2str(cutOffChange),num2str(cutOffFlux),num2str(length(formulas)));
0087     end
0088 else
0089     fprintf('These reactions have flux values that differ by more than %s percent (%s reactions)\n\n',num2str(cutOffChange),num2str(length(formulas)));
0090 end
0091 
0092 metaboliteNames=[];
0093 for i=1:length(metaboliteList)
0094    metaboliteNames=[metaboliteNames char(metaboliteList(i)) ' ']; 
0095 end
0096 
0097 if length(metaboliteNames)>0
0098     fprintf('Only prints reactions involving one or more of the following metabolites:\n%s\n\n',metaboliteNames)
0099 end
0100 
0101 for i=1:length(formulas)
0102     fluxText=['Flux: ' num2str(fluxesA(fluxIndexes(i))) ' Reference flux: ' num2str(fluxesB(fluxIndexes(i))) ' Difference: ' num2str(fluxesA(fluxIndexes(i))-fluxesB(fluxIndexes(i)))];
0103         fprintf('%s: %s\n\t%s\n\t%s\n\n', char(model.rxns(fluxIndexes(i))), char(formulas(i)),...
0104         char(model.rxnNames(fluxIndexes(i))),fluxText);
0105 end

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