sbioparamestim

Perform parameter estimation

sbioparamestim will be removed in a future release. Use sbiofit instead.

Syntax

[k, result]= sbioparamestim(modelObj, tspan, xtarget, observed_array, estimated_array)
[___]= sbioparamestim(___, observed_array, estimated_array, k0)
[___]= sbioparamestim(___, observed_array, estimated_array, k0, method)

Arguments

kVector of estimated parameter values. For all optimization methods except 'fminsearch', the parameters are constrained to be greater than or equal to 0.
resultStructure with fields that provide information about the progress of optimization.
modelObj

SimBiology® model object.

tspann-by-1 vector representing the time span of the target data xtarget.
xtarget

n-by-m matrix, where n is the number of time samples and m is the number of states to match during the simulation. The number of rows in xtarget must equal the number of rows in tspan.

observed_array

Either of the following:

  • Array of objects (species, compartment, or nonconstant parameter) in modelObj, whose values should be matched during the estimation process

  • Cell array of object names (species, compartment, or nonconstant parameter) in modelObj, whose values should be matched during the estimation process

    Note:   If duplicate names exist for any species or parameters, ensure there are no ambiguities by specifying either an array of objects or a cell array of qualified names, such as compartmentName.speciesName or reactionName.parameterName. For example, for a species named sp1 that is in a compartment named comp2, the qualified name is comp2.sp1.

The length of observed_array must equal the number of columns in xtarget. sbioparamestim assumes that the order of elements in observed_array is the same as the order of columns in xtarget.

estimated_array

Either of the following:

  • Array of objects (compartment, species, or parameter) in modelObj whose initial values should be estimated

  • Cell array of names of objects (compartment, species, or parameter) in modelObj whose initial values should be estimated

    Note:   If duplicate names exist for any compartments, species, or parameters, ensure there are no ambiguities by specifying either an array of objects or a cell array of qualified names, such as compartmentName.speciesName or reactionName.parameterName. For example, for a parameter named param1 scoped to a reaction named reaction1, the qualified name is reaction1.param1.

k0Numeric vector containing the initial values of compartments, species, or parameters to be estimated. The length of k0 must equal that of estimated_array. If you do not specify k0, or specify an empty vector for k0, then sbioparamestim takes initial values for compartments, species, or parameters from modelObj, or, if there are active variants, sbioparamestim uses any initial values specified in the active variants. For details about variants, see Variant object.
method

Optimization algorithm to use during the estimation process, specified by either of the following:

  • String specifying one of the following functions:

    • 'fminsearch'

    • 'lsqcurvefit'

    • 'lsqnonlin'

    • 'fmincon'

    • 'patternsearch'

    • 'patternsearch_hybrid'

    • 'ga'

    • 'ga_hybrid'

    • 'particleswarm'

    • 'particleswarm_hybrid'

    For descriptions of how sbioparamestim uses the previous functions, see the Function Descriptions table.

  • Two-element cell array, with the first element being one of the previous functions, and the second element being an options structure or object. Use an appropriate options structure or object for each method listed next.

    MethodOptions Structure or Object
    'fminsearch'optimset

    'fmincon'

    'lsqcurvefit'

    'lsqnonlin'

    'particleswarm'

    'particleswarm_hybrid'

    optimoptions
    'patternsearch'

    'patternsearch_hybrid'

    psoptimset
    'ga'

    'ga_hybrid'

    gaoptimset

      Tip   Use a two-element cell array to provide your own options structure for the optimization algorithm.

    If you have Parallel Computing Toolbox™, you can enable parallel computing for faster data fitting by:

    1. Opening a MATLAB® worker pool:

      parpool
    2. Setting the name-value pair argument 'UseParallel' to true in an options structure or object.

Function Descriptions

FunctionDescription
fminsearch

sbioparamestim uses the default options structure associated with fminsearch, except for:

Display = 'off'
TolFun = 1e-6* (Initial value of objective function)

    Note:   'fminsearch' is an unconstrained optimization method, which can result in negative values for parameters.

lsqcurvefit

Requires Optimization Toolbox™.

sbioparamestim uses the default options structure associated with lsqcurvefit, except for:

Display = 'off'
FinDiffRelStep = value of the SolverOptions.RelativeTolerance property of the configuration set associated with modelObj, with a minimum of eps^(1/3)
TolFun = 1e-6* (Initial value of objective function)
TypicalX = 1e-6* (Initial values of components to be estimated)
lsqnonlin

Requires Optimization Toolbox.

sbioparamestim uses the default options structure associated with lsqnonlin, except for:

Display = 'off'
FinDiffRelStep = value of the SolverOptions.RelativeTolerance property of the configuration set associated with modelObj, with a minimum of eps^(1/3)
TolFun = 1e-6* (Initial value of objective function)
TypicalX = 1e-6* (Initial values of components to be estimated)
fmincon

Requires Optimization Toolbox.

sbioparamestim uses the default options structure associated with fmincon, except for:

Algorithm = 'interior-point'
Display = 'off'
FinDiffRelStep = value of the SolverOptions.RelativeTolerance property of the configuration set associated with modelObj, with a minimum of eps^(1/3)
TolFun = 1e-6* (Initial value of objective function)
TypicalX = 1e-6* (Initial values of components to be estimated)
patternsearch

Requires Global Optimization Toolbox.

sbioparamestim uses the default options structure associated with patternsearch, except for:

Display = 'off'
TolFun = 1e-6* (Initial value of objective function)
TolMesh = 1.0e-3
Cache = 'on'
MeshAccel = 'on'
patternsearch_hybrid

Requires Global Optimization Toolbox.

sbioparamestim calls the patternsearch function with the additional option SearchMethod = {@searchlhs,10,15}. This option adds an additional search step that uses Latin hypercube sampling.

The sbioparamestim function uses the default options structure associated with patternsearch, except for:

Display = 'off'
TolFun = 1e-6* (Initial value of objective function)
TolMesh = 1.0e-3
Cache = 'on'
MeshAccel = 'on'
SearchMethod = {@searchlhs,10,15}
ga

Requires Global Optimization Toolbox.

sbioparamestim uses the default options structure associated with ga, except for:

Display = 'off'
TolFun = 1e-6* (Initial value of objective function)
PopulationSize = 10
Generations = 30
MutationFcn = @mutationadaptfeasible
ga_hybrid

Requires Global Optimization Toolbox.

sbioparamestim calls the ga function with the additional option HybridFcn = {@fmincon, fminopt}, where fminopt is the same set of default options sbioparamestim uses for fmincon. This option causes an additional gradient-based minimization after the genetic algorithm step ends.

The sbioparamestim function uses the default options structure associated with ga, except for:

Display = 'off'
TolFun = 1e-6* (Initial value of objective function)
PopulationSize = 10
Generations = 30
MutationFcn = @mutationadaptfeasible
HybridFcn = {@fmincon, structure of name/value pairs for fmincon}
particleswarm

Requires Global Optimization Toolbox.

sbioparamestim uses the following default options for particleswarm, except for:
Display          = 'off';
TolFun           = 1e-6*[Initial objective function value]
SwarmSize        = 10;
MaxIter          = 30;

particleswarm_hybrid

Requires Global Optimization Toolbox.

sbioparamestim calls the particleswarm function with the additional option HybridFcn = {@objFcn, options}. The objective function, objFcn, is one of these supported functions: patternsearch, fminsearch, fminunc, or fmincon. options is a structure of options for these functions and their values.

Display          = 'off';
TolFun           = 1e-6*[Initial objective function value]
SwarmSize        = 10;
MaxIter          = 30;
HybridFcn        = {@fmincon, [Fmincon Options, described above]}

    Note:   sbioparamestim does not support setting the Vectorized option to 'on' in algorithms that support this option.

Description

[k, result]= sbioparamestim(modelObj, tspan, xtarget, observed_array, estimated_array) estimates the initial values of compartments, species, and parameters of modelObj, a SimBiology model object, specified in estimated_array, so as to match the values of species and nonconstant parameters given by observed_array with the target state, xtarget, whose time variation is given by the time span tspan. If you have Optimization Toolbox installed, sbioparamestim uses the lsqnonlin function as the default method for the parameter estimation. If you do not have Optimization Toolbox installed, sbioparamestim uses the MATLAB function fminsearch as the default method for the parameter estimation.

[___]= sbioparamestim(___, observed_array, estimated_array, k0) specifies the initial values of compartments, species, and parameters listed in estimated_array.

[___]= sbioparamestim(___, observed_array, estimated_array, k0, method) specifies the optimization method to use.

Examples

Given a model and some target data, estimate all of its parameters without explicitly specifying any initial values:

  1. Load a model from the project, gprotein_norules.sbproj. The project contains two models, one for the wild-type strain (stored in variable m1), and one for the mutant strain (stored in variable m2). Load the G protein model for the wild-type strain.

    sbioloadproject gprotein_norules m1;
  2. Store the target data in a variable:

    Gt = 10000;
    tspan  = [0 10 30 60 110 210 300 450 600]';
    Ga_frac = [0 0.35 0.4 0.36 0.39 0.33 0.24 0.17 0.2]';
    xtarget = Ga_frac * Gt;
  3. Store all model parameters in an array:

    p_array = sbioselect(m1,'Type','parameter');
  4. Store the species that should match target:

    Ga = sbioselect(m1,'Type','species','Name','Ga');
    % In this example only one species is selected.
    % To match more than one targeted species data
    % replace with selected species array.
  5. Estimate the parameters:

    [k, result] = sbioparamestim(m1, tspan, xtarget, Ga, p_array)
    k =
    
        0.0100
        0.0000
        0.0004
        4.0000
        0.0040
        1.0000
        0.0000
        0.1100
    
    
    result = 
    
              fval: 1.4193e+06
          residual: [9x1 double]
          exitflag: 2
        iterations: 2
         funccount: 27
         algorithm: 'trust-region-reflective'
           message: [1x413 char]

Estimate parameters specified in p_array for species Ga using different algorithms. This example uses data from the first example.

[k1,r1] = sbioparamestim(m1,tspan,xtarget,Ga,p_array, ...
          {},'fmincon');
[k2,r2] = sbioparamestim(m1,tspan,xtarget,Ga,p_array, ...
          {},'patternsearch');
[k3,r3] = sbioparamestim(m1,tspan,xtarget,Ga,p_array, ...
          {},'ga');
[k4,r4] = sbioparamestim(m1,tspan,xtarget,Ga,p_array, ...
	  {},'particleswarm');

Estimate parameters specified in p_array for species Ga, and change default optimization options to use user-specified options. This example uses data from the first example.

myopt1 = optimset('Display','iter');
[k1,r1] = sbioparamestim(m1,tspan,xtarget, ...
          Ga,p_array,{},{'fmincon',myopt1});

myopt2 = psoptimset('TolMesh',1.0e-4); 
[k2,r2] = sbioparamestim(m1,tspan,xtarget, ...
          Ga,p_array,{},{'patternsearch',myopt2});

myopt3 = gaoptimset('PopulationSize',25, 'Generations', 10);
[k3,r3] = sbioparamestim(m1,tspan,xtarget, ...
          Ga,p_array,{},{'ga',myopt3});

myopt4 = optimoptions('particleswarm','Display','iter');
[k4,r4] = sbioparamestim(m1,tspan,xtarget,Ga,p_array,{},{'particleswarm',myopt4});

More About

expand all

Algorithms

sbioparamestim estimates parameters by attempting to minimize the discrepancy between simulation results and the data to fit. The minimization uses one of these optimization algorithms: fminsearch (from MATLAB); lsqcurvefit, lsqnonlinfit, or fmincon (from Optimization Toolbox); or patternsearch or ga (from Global Optimization Toolbox). All optimization methods require an objective function as an input. This objective function takes as input a vector of parameter values and returns an estimate of the discrepancy between simulation and data. When using lsqcurvefit or lsqnonlinfit as the optimization method, this objective function returns a vector of the residuals. For other optimization methods, the objective function returns the 2-norm of the residuals.

References

[1] Yi, T-M., Kitano, H., and Simon, M.I. (2003) A quantitative characterization of the yeast heterotrimeric G protein cycle. PNAS 100, 10764–10769.

Was this topic helpful?