Examples

This page will contain several example scripts/functions using the SBaddon package. So far there is only one.

Simple Parameter Estimation Function

The following two functions should be placed in the same file. They perform a very basic parameter estimation for a given model and given referencedata. These functions should be seen as an example on how MEX simulation functions can be used and not as "state-of-the-art" in parameter estimation.

MEXmodel: The MEX simulation function to use
timevector: Time vector to use for simulations (same as used to obtain the reference data)
reference: Matrix containing reference data used to calculate the cost in the costfunction
parameters: Cell-array containing the names of the parameters that are to be optimized
parameterv: Vector with initial guesses for the parameters

function[optvalues] = opt(MEXmodel,timevector,reference,parameters,parameterv)
% Global variables needed to pass reference data, timevector,
% model, and parameter names to cost function
global MEXmodel_global reference_global parameters_global timevector_global MEXmodel_global = MEXmodel; timevector_global = timevector; reference_global = reference; parameters_global = parameters; % Set desired options for the optimization function
optionsSimplexSB = [];
% Start the optimization function (here the simplexSB function from the SBTOOLBOX is used)
optvalues = simplexSB(@costfunction,parameterv,optionsSimplexSB);
return

function [cost] = costfunction(x)
global MEXmodel_global reference_global parameters_global timevector_global
% Define a full parameter vector for the MEX simulation function
% based on the varied parameters
pv = createparametervectorSBAO(MEXmodel_global,parameters_global,x)
% Set options for the simulation
optionsMEXmodel = [];
% Simulate the model
try
    simdata = feval(MEXmodel_global,timevector_global,[],pv,optionsMEXmodel);
catch
    % If the simulation failed this is most often due to completely wrong
    % parameter values, chosen by the optimization function. To make it nicer
    % you could make a check if the parameter values in "x" are within given
    % bounds before starting the simulation.
    cost = inf;
    return
end
% Get the simulation data you want to use for optimization
% (needs to have the same format as the referencedata). In this example it is
% assumed that all state variables are measured and available in the
% referencedata in the same order as in the simulation data
simulation = simdata.statevalues;
% Use sum of squared errors to determine the cost
cost = norm(simulation-reference_global);
return