For MATLAB R2014b and later, specifying parameter bounds during parameter estimation is directly supported by using "sbiofit" (from the command line) and the Fit Task (from the SimBiology Desktop):
1) Load your data
load(fullfile(matlabroot,'examples','simbio','data15.mat'))
2) Create a SimBiology model or use an existing model
pkmd = PKModelDesign;
pkc1 = addCompartment(pkmd,'Central');
pkc1.DosingType = 'Bolus';
pkc1.EliminationType = 'linear-clearance';
pkc1.HasResponseVariable = true;
model = construct(pkmd);
configset = getconfigset(model);
configset.CompileOptions.UnitConversion = true;
3) Define dosing
dose = sbiodose('dose');
dose.TargetName = 'Drug_Central';
dose.StartTime = 0;
dose.Amount = 10;
dose.AmountUnits = 'milligram';
dose.TimeUnits = 'hour';
4) Specify the parameters to estimate
paramsToEstimate = {'log(Central)','log(Cl_Central)'};
5) Specify the boundaries of the parameters
estimatedParams = estimatedInfo(paramsToEstimate,'InitialValue',[1 1],'Bounds',[1 5;0.5 2]);
6) Estimate the parameters
fitConst = sbiofit(model,gData,responseMap,estimatedParams,dose);
You can find the complete example in the attached Fitting_with_Bounds_Example.m.
For MATLAB R2014a and older the above approach above cannot be used. For these releases use the following approach:
The basic steps are as follows:
1) Load your .sbproj file
sbioloadproject('Tumor_Growth.sbproj', 'm1')
2) Load your data
ds = dataset('xlsfile', 'Data.xlsx') ;
3) Create vectors containing the input arguments for FMINCON. Among them you can specify upper and lower boundaries for the parameters as well as the initial values
p0 = [0.1, 0.1, 0.01] ; % initial estimate
lb = [1e-3 1e-3 1e-3] ; % lower bound
ub = [10 10 10] ; % upper bound
A = [] ; %
b = [] ; %
Aeq = [] ; %
beq = [] ; %
nlc = [] ; %
4) Define the optimization options
opt1 = optimoptions(@fmincon , ...
'Display' , 'iter' , ... %Iterations should be displayed
'Algorithm' , 'interior-point' , ... %Optimization Algorithm that should be used
'FinDiffRelStep' , m1Exp.SimulationOptions.RelativeTolerance ,.... %Step size factor
'TolFun' , 1e-3 ,... %Termination tolerance on the functionvalue, a positive scalar
'TolX' , 1e-6) ; %Termination tolerance on x,a positive scalar
5) Create your regression function
a) Set values of parameters to be estimated to value in current iteration, p
b) Set dose, amount and time
c) Call SIMULATE to simulate the model for this parameter set
d) Remove NaNs from the output
e) Calculate the sum of the squared errors
function f = objectiveFcn(p, yObs, m1Exp, ds)
% p - parameter
% yObs - Observed response
% yPred - predicted response
% f - value of the objective function
%%Simulation
m1Exp.InitialValues = p ; % Set values of parameters to be estimated to value in current iteration, p
d1 = getdose(m1Exp) ; % get dose
d1.Amount = ds.Dose(~isnan(ds.Dose)) ; % set amount
d1.Time = ds.Time(~isnan(ds.Dose)) ; % set time
m1Exp.SimulationOptions.OutputTimes = ds.Time ; % Set sampling time
[~, yPred] = simulate(m1Exp, d1) ; % Simulate
idx = ~isnan(yObs) ; % Identify non-nan values in yObs
%%Objective function: sum of squared errors
f = sum((yPred(idx) - yObs(idx)).^2) ;
end
6) Call your fitting algorithm
fh = @(p) objectiveFcn(p, ds.Tumor, m1Exp, ds) ;
pFit_FMC = fmincon(fh, p0, A, b, Aeq, beq, lb, ub, nlc, opt1 ) ;
You can also use another fitting algorithm that supports boundaries, for example lsqnonlin.
