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

Arguments

k Vector of estimated parameter values. For all optimization methods except 'fminsearch', the parameters are constrained to be greater than or equal to zero.

result Structure with fields that provide information about the progress of optimization.

modelObj

SimBiology model object.

tspan n-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 species objects or nonconstant parameter objects in modelObj, whose values should be matched during the estimation process
Cell array of names of species objects or nonconstant parameter objects 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.

k0 Numeric 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'
For descriptions of how sbioparamestim uses the previous functions, see the following table.
Two-element cell array, with the first element being one of the previous functions, and the second element being an options structure such as returned by optimset, psoptimset, or gaoptimset.

Function Description

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 terminates.

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, struct of name/value pairs for fmincon}

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.

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.

Examples

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

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;
```

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;

Store all model parameters in an array:

p_array = sbioselect(m1,'Type','parameter');

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.

Estimate the parameters:

[k, result] = sbioparamestim(m1, tspan, xtarget, Ga, p_array)

k =

    0.1988
    0.0000
    0.0045
    6.2859
    0.0040
    0.9726
    0.0000
    0.1164

result = 

          fval: 8.7248e+005
      residual: [9x1 double]
      exitflag: 2
    iterations: 2
     funccount: 27
     algorithm: 'large-scale: trust-region reflective Newton'
       message: [1x77 char]

Estimate parameters specified in p_array and species specified in sp_array 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')

Estimate parameters specified in p_array, species specified in sp_array, 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, ...
          sp_array, p_array, {},{'fmincon', myopt1});

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

myopt3.PopulationSize = 50; 
myopt3.Generations = 20; 
[k3,r3] = sbioparamestim(m1, tspan, xtarget, ...
          sp_array, p_array, {},{'ga', myopt3});