lsqnonlin initial conditions (transcendental equation)

5 visualizzazioni (ultimi 30 giorni)
Hey all,
I am trying to fit a transcendental equation to some data points but I am coming across the issue where lsqnonlin ends up giving solutions ("possible minima") that heavily depend upon the choice of initial condition. I am wondering if this is because my choice of initial condition simply is incorrect or if there is something else going on and how I could better probe what exactly is happening. I have attached the data_points below and my code is:
a1 = 20;
a2 = 10^-4;
x_data = importdata("x_data.mat");
y_data = importdata("y_data.mat");
options = optimoptions(@lsqnonlin,'Algorithm','levenberg-marquardt','StepTolerance',1e-12,'FunctionTolerance',1e-12,'Display','iter');
a0 = [10^-3,1,5,2]; %initial guess%
m = @(a)sqrt(a2^2 + (x_data.*a(2)).^2 + a(3).*(1./a1).^a(4));
factor1 = @(a)sqrt((x_data+1i.*(m(a)./a(2)))./(x_data-1i.*(m(a)./a(2))));
b = @(a)(1/2) - 1i.*(y_data./a(2));
factor2 = @(a)double(gamma(sym(complex(1-b(a))))./gamma(sym(complex(b(a)))));
factor3 = @(a)(2./(a(1).*sqrt(x_data.^2+(m(a).^2)./a(2).^2))).^(2.*1i.*y_data./a(2));
func = @(a)factor1(a).*factor2(a).*factor3(a) - 1i;
[k,resnorm,residual,exitflag,output] = lsqnonlin(func,a0,[],[],options)
First-order Norm of Iteration Func-count Resnorm optimality Lambda step 0 5 30.8735 4.5e+05 0.01 1 13 25.0777 2.5e+05 10 0.0800468 2 18 22.909 1.49e+09 1 0.0593238 3 25 22.8683 2.07e+04 100 0.0010927 4 30 22.531 1.21e+05 10 0.00652635 5 35 22.5103 3.04e+08 1 0.0514446 6 45 22.3993 7.89e+04 100000 0.000507316 7 50 22.3525 3.14e+04 10000 0.000157084 8 55 22.3327 1.71e+04 1000 0.0005333 9 60 22.326 7.43e+03 100 0.00057971 10 65 22.3236 3.02e+08 10 0.00486204 11 75 22.3235 1.49e+03 1e+06 5.98625e-05 12 81 22.3235 583 1e+07 3.31727e-07 13 86 22.3235 224 1e+06 3.21627e-07 14 91 22.3235 117 100000 3.1936e-06 15 96 22.3233 327 10000 2.66618e-05 16 101 22.3227 1.21e+03 1000 9.88313e-05 17 106 22.3224 6.84e+08 100 0.000484278 18 115 22.3223 3.02e+08 1e+06 6.04123e-05 19 126 22.3223 79.5 1e+12 1.54083e-08 20 129 7.94303e-13 Local minimum possible. lsqnonlin stopped because the relative size of the current step is less than the value of the step size tolerance.
k =
0.0010 + 0.0000i 1.0263 + 0.0021i 4.9869 - 0.0016i 2.1974 + 0.0242i
resnorm = 22.3223
residual =
Columns 1 through 10 0.0866 + 0.0926i -0.2591 + 0.0924i 0.2009 - 2.1180i -0.1175 + 0.1314i -0.1175 + 0.1314i -0.9598 - 1.5939i 0.4355 + 0.0215i 0.4355 + 0.0215i 0.4276 - 0.0045i 0.4276 - 0.0045i Columns 11 through 20 0.1274 + 0.0405i 0.1274 + 0.0405i 0.0000 + 0.0000i -0.1159 - 0.0532i -0.3642 - 0.1519i -0.3642 - 0.1519i -0.3532 - 0.1716i -0.3532 - 0.1716i 0.7534 - 1.4662i 0.7534 - 1.4662i Columns 21 through 27 0.0908 - 0.1256i -0.1557 - 1.8665i -0.1557 - 1.8665i 0.2056 - 0.1334i 0.2056 - 0.1334i -0.0721 - 0.0905i -0.3719 - 0.1255i
exitflag = 4
output = struct with fields:
iterations: 20 funcCount: 129 stepsize: 7.9430e-13 cgiterations: [] firstorderopt: 79.5088 algorithm: 'levenberg-marquardt' message: 'Local minimum possible.↵lsqnonlin stopped because the relative size of the current step is less than↵the value of the step size tolerance.↵↵<stopping criteria details>↵↵Optimization stopped because the relative norm of the current step, 1.432370e-13,↵is less than options.StepTolerance = 1.000000e-12.' bestfeasible: [] constrviolation: []
Thanks.
  2 Commenti
Matt J
Matt J il 20 Mar 2023
Modificato: Matt J il 20 Mar 2023
It's pretty amazing that the optimization runs to completion at all, given that your objective is complex-valued. Generally speaking, that isn't allowed, but for exceptions discussed here:
bil
bil il 20 Mar 2023
Hi, Could it be the case that my function is still analytical? Regardless, do you have any suggestion regarding how else I could fit the parameters of the objective function to the data set?

Accedi per commentare.

Risposta accettata

Matt J
Matt J il 20 Mar 2023
You can split the model function into real and imaginary parts,
func=@(a) [real(func(a)); imag(func(a))];
[k,resnorm,residual,exitflag,output] = lsqnonlin(func,a0,[],[],options)
  4 Commenti
bil
bil il 20 Mar 2023
Modificato: bil il 20 Mar 2023
Could you elaborate on what you mean by well-arranged graphic. Do you mean like the explicit function itself?
Edit: Here is the function anyway
where
and
where are defined as originally in my code and the fitted parameters are and correspond to the x_data and y_data respectively.
Matt J
Matt J il 20 Mar 2023
Modificato: Matt J il 20 Mar 2023
When you say you are getting highly different solutions for different initial points, are you getting very different resnorms in each case as well? If not, the problem is simply ill-posed: you have non-unique solutions.
If you are getting very different resnorms, then you are probably getting stuck at local minima, in which case you could contemplate doing a 4D parameter search over the different k(i) combinations for a more accurate initial guess k0. This shouldn't be very computationally demanding, since the operations in func() are very well-vectorized. Can you write down lb,ub bounds on the parameters? If so, you can download ndgridVecs and create grid search data quite easily, e.g. as below,
ranges=arrayfun(@(l,u)linspace(l,u,30), lb,ub,'uni',0); %30-point discretization
[~,~,K1,K2,K3,K4]=ndgridVecs(1,1,ranges{:}); %grid search coordinates
values = vecnorm(func(K1,K2,K3,K4)); %rewrite func to support 4-argument input syntax
[~,i]=min(values,[],'all','linear');
k0=[K1(i), K2(i), K3(i), K4(i)]; %initial point

Accedi per commentare.

Più risposte (1)

John D'Errico
John D'Errico il 20 Mar 2023
This does not make sense.
x_data = importdata("x_data.mat");
x_data
x_data = 1×27
1.0e-42 * 0.7175 0.7175 0.7175 0.7175 0.7175 0.7175 0.7175 0.7175 0.7175 0.7175 0.7175 0.7175 0.7175 0.7175 0.7175 0.7175 0.7175 0.7175 0.7175 0.7175 0.7175 0.7175 0.7175 0.7175 0.7175 0.7175 0.7175
Your x_data vector is identically constant values.
unique(x_data)
ans = 7.1746e-43
Why do you expect any model to be fit here?
  1 Commento
bil
bil il 20 Mar 2023
Modificato: bil il 20 Mar 2023
Hi,
The transcendental equation should have multiple solutions for the same value of x. The simplest example is how can have multiple (actually infinitely many) solutions x for a single value of .

Accedi per commentare.

Categorie

Scopri di più su Mathematics in Help Center e File Exchange

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by