# "Solve" for optimization problem does not store linear parameters

3 views (last 30 days)
Samuele Bolotta on 27 Mar 2021
Answered: Matt J on 27 Mar 2021
I am fitting a function to some data I simulated. I managed to get intelligent constraints that help the fit quite a bit, even with a lot of noise.
This is the function; DF1 and DF2 are constant, c(1) and c(2) are linear, while lam(1), lam(2), lam(3) and lam(4) are nonlinear. I am following the procedure explained here (https://it.mathworks.com/help/optim/ug/nonlinear-data-fitting-problem-based-example.html#NonlinearDataFittingProblemBasedExample-4) to split linear and nonlinear parameters.
% Create a function that computes the value of the response at times t when the parameters are c and lam
temp1 = 1 - exp(-t / lam(1));
temp2 = exp(-t / lam(2));
temp3 = 1 - exp(-t / lam(3));
temp4 = exp(-t / lam(4));
diffun = (temp1 .* temp2) * DF1 * c(1) + (temp3 .* temp4) * DF2 * c(2);
This is the code that I came up with, but it's not working because it does not store correctly the two values for c.
To generate the data:
function [EPSC, IPSC, CPSC, t] = generate_current(G_max_chl, G_max_glu, EGlu, EChl, Vm, tau_rise_In, tau_decay_In, tau_rise_Ex, tau_decay_Ex,tmax)
dt = 0.1; % time step duration (ms)
t = 0:dt:tmax-dt;
% Compute compound current
IPSC = ((G_max_chl) .* ((1 - exp(-t / tau_rise_In)) .* exp(-t / tau_decay_In)) * (Vm - EChl));
EPSC = ((G_max_glu) .* ((1 - exp(-t / tau_rise_Ex)) .* exp(-t / tau_decay_Ex)) * (Vm - EGlu));
CPSC = IPSC + EPSC;
end
To fit the function:
% Simulated data
[EPSC,IPSC,CPSC,t] = generate_current(80,15,0,-70,-30,0.44,15,0.73,3,120);
ydata = awgn(CPSC,25,'measured'); % Add white noise
ydata = ydata';
t = t';
% Values
Vm = -30; DF1 = Vm +70; DF2 = Vm;
% Initial values for fitting
gmc = 40; gmg = 20; tde = 0.2; tdi = 8; tre = 1.56; tri = 3;
% Objective function
c = optimvar('c',2); % Linear parameters
lam = optimvar('lam',4); % Nonlinear parameters
% Bounds
c.LowerBound = [0, 0];
c.UpperBound = [200, 200];
lam.LowerBound = [0.16,7.4,1.1,2.6];
lam.UpperBound = [0.29,8.4,2.3,3.3];
x0.c = [gmc,gmg]; % Starting values
x0.lam = [tri,tdi,tre,tde]; % Starting values
% Create a function that computes the value of the response at times t when the parameters are c and lam
temp1 = 1 - exp(-t / lam(1));
temp2 = exp(-t / lam(2));
temp3 = 1 - exp(-t / lam(3));
temp4 = exp(-t / lam(4));
diffun = (temp1 .* temp2) * DF1 * c(1) + (temp3 .* temp4) * DF2 * c(2);
%Solve the problem
%To do so, first convert the fitvector function to an optimization expression using fcn2optimexpr.
F2 = fcn2optimexpr(@(x) fitvector(x,t,ydata),lam,'OutputSize',[length(t),1]);
% Create a new optimization problem with objective as the sum of squared differences between the converted fitvector function and the data y
ssqprob = optimproblem('Objective',sum((F2 - ydata).^2));
[sol,fval,exitflag,output] = solve(ssqprob,x0)
% Plot
resp = evaluate(diffun,sol);
hold on
plot(t,resp)
hold off
Fitvector function:
function yEst = fitvector(lam,xdata,ydata)
temp1 = 1 - exp(-xdata / lam(1));
temp2 = exp(-xdata / lam(2));
temp3 = 1 - exp(-xdata / lam(3));
temp4 = exp(-xdata / lam(4));
A(:,1) = temp1 .* temp2 * 40; % DF1 is 40
A(:,2) = temp3 .* temp4 * -30; % DF2 is -30
c = A\ydata; % solve A*c = y for linear parameters c
yEst = A*c; % return the estimated response based on c
Not sure what I am doing wrong, because I get this error:
Error using optim.problemdef.OptimizationExpression/evaluate
Second argument must contain values for variable 'c'.
Error in SplittFit (line 44)
resp = evaluate(diffun,sol);

Matt J on 27 Mar 2021
Have you looked at what sol contains? I'm willing to bet it contains values for lam, but not c. Your objective F2 is only a function of lam, after all, is it not?