Fit a model to data using lsqnonlin()

We are trying to fit a model :

=> f(x; alpha, beta) = alpha*exp(-x / beta)

with beta > 0 on some data (x, z).

The goal is therefore to compute the values of coefficients alpha and beta corresponding to the minimum squared error.

Here is the data generator used to generate (x, z) :

randn('seed',3); % always the same source of randomness

theta_1 = 5; % ground truth
theta_2 = 1; % ground truth

x = 0:0.3:19;
y = exp(-x/theta_1)-0.8*exp(-x/theta_2);
N = length(x);
noise = 0.03; % noise level (can be changed)
z = y + noise*randn(1, N);

figure(1), clf, plot(x,z,'o','linewidth',2); grid on;
save 'data0.mat'

We have used lsqnonlin

% Optimization of alpha and beta coefficients
fun = @(v)v(1) * exp(-x / v(2)) - z;
 
lb = [-Inf, 0.0001];
ub = [Inf, Inf];
x0 = [0, 0.0001];
[ab, f_val] = lsqnonlin(fun, x0, lb, ub);
 

But it seems like the coefficients we are finding are not the right one.

Would you have an idea about what's wrong ?

Thanks a lot.

Answer 1

[Matt J]

Your initial guess is way off. With x0=[1,1] I get something very reasonable looking.

[Quentin Pradelle]

Well, that was efficient. Thanks a lot.

Also I'm generally confused about how to determine initial guesses, do you have any advice regarding this ?

[Matt J]

There is no general method for choosing the initial guess, but it shouldn't be random. In this case, you know that the thetas are on the order of 1 or 5, not .0001. A value of 0.0001 is, in fact, so far away that it causes the exp() calculation to underflow for the values x = 0:0.3:19 that you are considering,

>> exp(-19/.0001)

ans =

     0

>> isequal(0,ans)

ans =

  logical

   1