Model fit using fminunc based on measured data

Hi all,

I have measurement data that I would like to use to fit a model. I want to use the command "fminunc" for the fit, as I will later have a model with several coefficients (9) and variables (5). However, in order to get my script to run at all, I have reduced it and you can find it below as a minimal example.

These coefficients (in my minimal example these are A0 and EA) are to be fitted in such a way that the function including these fitted coefficients represents the entire result of my experiments as best as possible (each coefficient is a constant for the entire experimental data set). In my minimal example, I had four test runs at different temperatures (x1), which resulted in different lifetimes (x2).

Code:

kB = 8.617 * 1e-5; % in eV/K

x1 = [233; 264; 295; 326]; % temperatures

x2 = [420000; 970000; 3800000; 10000000]; % lifetimes

% Function to solve

fun = @(A0,EA) ((A0 .* exp(EA/kB./x1)) - x2).^2;

x0 = [2.7 * 1e10, 0.220]; % initial values

% Solve

[p, fval] = fminunc(fun,x0);

However, I only get an error message. Can anyone help me find out where my error is?

Many thanks in advance!

BR, Seb

Answer 1

[Bruno Luong]

kB = 8.617 * 1e-5; % in eV/K
x1 = [233; 264; 295; 326]; % temperatures
x2 = [420000; 970000; 3800000; 10000000]; % lifetimes
% Function to solve
fun = @(x) sum((x1 .* exp(x(2)/kB./x1) - x2).^2); % Change here
x0 = [2.7 * 1e10, 0.220];   % initial values
% Solve
[p, fval] = fminunc(fun,x0);

Local minimum possible.

fminunc stopped because it cannot decrease the objective function
along the current search direction.

[SebK]

Thank you for your response. I changed x1 to x(1) in line 5 because this should be a coefficient and not the temperature vector.

Solving the function with the fitted coefficients EA and A1 leads leads only to "Infs" in the lifetime result vector, unfortunatelly.

For valdidation reasons I fitted the data with lsqcurvefit and with this command I get the model fit which you can find attached. (as mentioned before I will have several coefficients and variables later so lsqcurvefit is only possible in this minimal example)

Can you imagine where the problem is with fminunc?

BR, Markus

kB = 8.617 * 1e-5; % in eV/K
x1 = [233; 264; 295; 326];
x2 = [420000; 970000; 3800000; 10000000];

% Function to solve
fun = @(x) sum((x(1) .* exp(x(2)/kB./x1) - x2).^2);   % x(1) instead of x1
x0 = [2.5 * 1e10, 0.200];   % initial values

% Solve
[p, fval] = fminunc(fun,x0);

% Optimised parameter values
EA = p(1);
A1 = p(2);

temperature = x1(1) : 0.01 : x1(4);   % temperature vector for plot
    for i = 1:length(temperature)
        cycles(i) = A1 * exp(EA/ kB / temperature(i)); % cycles correspond to lifetime
    end
    
figure
plot(temperature, cycles)

[Bruno Luong]

I belive the order is

EA = p(2);
A1 = p(1);

That having said it is well known that the exponetial fit often stuck in local minima if you don't have good estimation of x0.

[Bruno Luong]

kB = 8.617 * 1e-5; % in eV/K
x1 = [233; 264; 295; 326];
x2 = [420000; 970000; 3800000; 10000000];
% Function to solve
fun = @(x) sum((x(1) .* exp(x(2)/kB./x1) - x2).^2);   % x(1) instead of x1
P=polyfit(1./(kB*x1),log(x2),1);
A1 = exp(P(2));
EA = P(1);
x0 = [A1,EA];   % initial values
% Solve
[p, fval] = fminunc(fun,x0);

Local minimum possible.

fminunc stopped because it cannot decrease the objective function
along the current search direction.

% Optimised parameter values
EA = p(2);
A1 = p(1);
temperature = x1(1) : 0.01 : x1(4);   % temperature vector for plot
for i = 1:length(temperature)
    cycles(i) = A1 * exp(EA/ kB / temperature(i)); % cycles correspond to lifetime
end

figure
plot(temperature, cycles,'b',x1,x2,'or')

[SebK]

The result looks good but fminunc is not needed any more, isn't it?

So with more coefficients and variables I can't solve this equation any more. Or do I misunderstand your solution?

[Torsten]

Yes, you misunderstand the solution.

"polyfit" is called for the linearized equation and gives good initial guesses for the parameters.

These initial guesses are then used in a second step to solve the nonlinear fitting problem with "fminunc".

[Bruno Luong]

I only fix your original code.

I won't comment on general case beside what I told you: fitting exponantial is well known to have local minima, so you have to work through to overcome this and not throw the function and first guess without care.

[SebK]

Now it works. Great job!!!. Thank you two ;-)