Home > RAVEN > fitParameters.m

fitParameters

PURPOSE ^

fitParameters

SYNOPSIS ^

function [parameters fitnessScore exitFlag newModel]=fitParameters(model,xRxns,xValues,rxnsToFit,valuesToFit,parameterPositions,fitToRatio,initialGuess,plotFitting)

DESCRIPTION ^

 fitParameters
   Fits parameters such as maintenance ATP by quadratic programming

   model                 a model structure
   xRxns                 cell array with the IDs of the reactions that will be fixed for each data point
   xValues               matrix with the corresponding values for each
                         xRxns (columns are reactions)
   rxnsToFit             cell array with the IDs of reactions that will be fitted to
   valuesToFit           matrix with the corresponding values for each 
                         rxnsToFit (columns are reactions)
   parameterPositions    stucture that determines where the parameters are in the
                         stoichiometric matrix. Contains the fields:
       position          cell array of vectors where each element contains
                         the positions in the S-matrix for that parameter
       isNegative        cell array of vectors where the elements are true
                         if that position should be the negative of the
                         fitted value (to differentiate between
                         production/consumption)
    fitToRatio            if the ratio of simulated to measured values should
                         be fitted instead of the absolute value. Used to prevent 
                         large fluxes from having too large impact (opt,
                         default true)
   initialGuess          initial guess of the parameters (opt)
   plotFitting           true if the resulting fitting should be plotted
                         (opt, default false)

   parameters            fitted parameters in the same order as in
                         parameterPositions
   fitnessScore          the correponding residual sum of squares
   newModel              updated model structure with the fitted parameters

   Usage: [parameters fitnessScore exitFlag newModel]=fitParameters(model,...
           xRxns,xValues,rxnsToFit,valuesToFit,parameterPositions,fitToRatio,...
           initialGuess,plotFitting)

   Rasmus Agren, 2010-12-16

CROSS-REFERENCE INFORMATION ^

This function calls: This function is called by:

SUBFUNCTIONS ^

SOURCE CODE ^

0001 function [parameters fitnessScore exitFlag newModel]=fitParameters(model,xRxns,xValues,rxnsToFit,valuesToFit,parameterPositions,fitToRatio,initialGuess,plotFitting)
0002 % fitParameters
0003 %   Fits parameters such as maintenance ATP by quadratic programming
0004 %
0005 %   model                 a model structure
0006 %   xRxns                 cell array with the IDs of the reactions that will be fixed for each data point
0007 %   xValues               matrix with the corresponding values for each
0008 %                         xRxns (columns are reactions)
0009 %   rxnsToFit             cell array with the IDs of reactions that will be fitted to
0010 %   valuesToFit           matrix with the corresponding values for each
0011 %                         rxnsToFit (columns are reactions)
0012 %   parameterPositions    stucture that determines where the parameters are in the
0013 %                         stoichiometric matrix. Contains the fields:
0014 %       position          cell array of vectors where each element contains
0015 %                         the positions in the S-matrix for that parameter
0016 %       isNegative        cell array of vectors where the elements are true
0017 %                         if that position should be the negative of the
0018 %                         fitted value (to differentiate between
0019 %                         production/consumption)
0020 %    fitToRatio            if the ratio of simulated to measured values should
0021 %                         be fitted instead of the absolute value. Used to prevent
0022 %                         large fluxes from having too large impact (opt,
0023 %                         default true)
0024 %   initialGuess          initial guess of the parameters (opt)
0025 %   plotFitting           true if the resulting fitting should be plotted
0026 %                         (opt, default false)
0027 %
0028 %   parameters            fitted parameters in the same order as in
0029 %                         parameterPositions
0030 %   fitnessScore          the correponding residual sum of squares
0031 %   newModel              updated model structure with the fitted parameters
0032 %
0033 %   Usage: [parameters fitnessScore exitFlag newModel]=fitParameters(model,...
0034 %           xRxns,xValues,rxnsToFit,valuesToFit,parameterPositions,fitToRatio,...
0035 %           initialGuess,plotFitting)
0036 %
0037 %   Rasmus Agren, 2010-12-16
0038 %
0039 
0040 if nargin<7
0041     fitToRatio=true;
0042 end
0043 if nargin<8
0044     initialGuess=ones(numel(parameterPositions.position),1);
0045 end
0046 if isempty(initialGuess)
0047     initialGuess=ones(numel(parameterPositions.position),1);
0048 end
0049 if nargin<9
0050     plotFitting=false;
0051 end
0052 
0053 %Find the indexes of reactions that will be fitted
0054 [I rxnsToFitIndexes]=ismember(rxnsToFit,model.rxns);
0055 
0056 if ~all(I)
0057     throw(MException('','Could not find all reactions in rxnsToFit'));
0058 end
0059 
0060 %Find the indexes of reactions that will be used for constraints.
0061 [I xRxnsIndexes]=ismember(xRxns,model.rxns);
0062 
0063 if ~all(I)
0064     throw(MException('','Could not find all reactions in xRxns'));
0065 end
0066 
0067 [parameters fitnessScore exitFlag]=fminsearch(@(parameters) getRSS(parameters,model,xRxnsIndexes,xValues,rxnsToFitIndexes,valuesToFit,parameterPositions,fitToRatio),initialGuess);
0068 
0069 parameters=abs(parameters);
0070 
0071 if plotFitting==true
0072     %Set the resulting parameters
0073     [rss resultingFluxes newModel]=getRSS(parameters,model,xRxnsIndexes,xValues,rxnsToFitIndexes,valuesToFit,parameterPositions,true);
0074     plot(xValues,valuesToFit,'o',xValues,resultingFluxes,'-*');
0075 end
0076 end
0077 
0078 function [rss resultingFluxes newModel]=getRSS(parameters,model,xRxnsIndexes,xValues,rxnsToFitIndexes,valuesToFit,parameterPositions,fitToRatio)
0079 parameters=abs(parameters);
0080 
0081 %Set the parameters at the positions specified in parameterPositions
0082 for i=1:numel(parameterPositions.position)
0083     %Set positive
0084     model.S(parameterPositions.position{i}(parameterPositions.isNegative{i}==false))=parameters(i);
0085     
0086     %Set negative
0087     model.S(parameterPositions.position{i}(parameterPositions.isNegative{i}==true))=parameters(i)*-1;
0088 end
0089 
0090 %Also return an updated model
0091 newModel=model;
0092 
0093 %Loop through each data point, set xRxns to xValues and calculate the sum
0094 %of squares for the rxnsToFit
0095 rss=0;
0096 resultingFluxes=[];
0097 for i=1:size(xValues,1)
0098     %Fix for more xRxns!
0099     model.lb(xRxnsIndexes)=xValues(i,:);
0100     model.ub(xRxnsIndexes)=xValues(i);
0101     
0102     sol=solveLP(model);
0103     
0104     %Calculate the rss
0105     if fitToRatio==false
0106         rs=sol.x(rxnsToFitIndexes)'-valuesToFit(i,:);
0107     else
0108         rs=sol.x(rxnsToFitIndexes)'./valuesToFit(i,:)-ones(1,size(valuesToFit,2));
0109     end
0110     rss=rss+rs*rs';
0111     resultingFluxes=[resultingFluxes sol.x(rxnsToFitIndexes)];
0112 end
0113 end

Generated on Tue 16-Jul-2013 21:50:02 by m2html © 2005