nonlinear fit with function a*exp((q1+​i*q2)*x)/s​qrt(x)

I am trying to do a nonlinear fit of the attached data, which is the blue decaying sine or cosine wave in the attached image.

I want to fit the data into the function of a*exp((q1+i*q2)*x)/sqrt(x), where a, q1 and q2 are real numbers.

I have tried to use lsqnonlin, but could not limit those three fitting variables to be real. And the fitted results do not make sense.

Any help will be much appreciated.

[Matt J]

We need to see the code you tried, so that we can recommend modifications.

[Alex Sha]

Hi, refer to the result below:

Root of Mean Square Error (RMSE): 3.84515212861012E-7
Sum of Squared Residual: 2.08471247979384E-11
Correlation Coef. (R): 0.894601522262208
R-Square: 0.800311883633859

Parameter          	Best Estimate
--------------------	-------------
a	3.5614255240146E-9
q1	-569495.123566796
q2	13936802.0739758

[Jiong Yang]

Hi Alex,

Thank you for sharing the results.

From a physical point of view, q1 and q2 should be both positive, and q2 should be smaller than q1 with a rough ratio of q1 ~ 10*q2. It may be a good idea to fit from the second peak, rather than the first peak.

Also, could you share the matlab code?

Best Regards

[Alex Sha]

Hi, if you want both of q1 and a2 to be positive, the results below may be you want:

Root of Mean Square Error (RMSE): 4.84330718082384E-7
Sum of Squared Residual: 3.30752504714258E-11
Correlation Coef. (R): 0.822474332455383
R-Square: 0.676464027547928

Parameter          	Best Estimate
--------------------	-------------
a	1.35682343774127E-9
q1	9.79296293081736E-11
q2	13826159.1347865

[Jiong Yang]

Thank you, Alex. Problem solved already.

Answer 1

[Bjorn Gustavsson]

This seems like you need to take a standard step from the physicsist's complex representation (which admittedly we often use a bit carelessly as a very convenient(!) shorthand). My guess is that you need to do something like this:

% your model-function: the real part of a modified damped oscillating
% exponential:
curve_fcn = @(pars,x) real((pars(1)+1i*pars(2))*exp((pars(3)+1i*pars(4))*x)./sqrt(x));
err_fcn = @(pars,x,y,fcn) sum((y-fcn(pars,x)).^2);
par0 = [-0.1378e-8, -0.0541e-8, 4.1643e+03, 2.3100e7]; % this might have to be adjusted to get a good enough start-guess
parBest = fminsearch(@(pars) err_fcn(pars,x,y,curve_fcn),par0);
subplot(2,1,1)
plot(x,y,x,curve_fcn(parBest,x))
subplot(2,1,2)
plot(x,y-curve_fcn(parBest,x))

You might have to iterate the optimization-step with the parameter estimate parBest as the new input for par0 to get a improved fits.

You might consider a properly weighted sum-of-square optimization instead of the straight-forward sum-of-squares, but in order to do that you'll need estimates of the uncertainties (standard-deviation) of each measurement point.

HTH

[Jiong Yang]

Hi Bjorn,

Thank you very much for the code. I made a minot adjustment by adding "options = optimset('MaxFunEvals', 1e3)", and I got fitting results of q1 = -9.6476e5 and q2 = 1.3983e7. Please find the attached matlab code file and plot. The fitting seems to be quite correct.

However, from a physical point of view, q1 and q2 should be both positive, and q2 should be smaller than q1 with a rough ratio of q1 ~ 10*q2.

Could you help? Thanks a lot in advance.

Best

[Bjorn Gustavsson]

You must have made a typo somewhere, if both q1 and q2 should be positive you'll have a growing exponential factor - while it seems your signal decays. After a couple of re-runs I got the best fitting parameters:

parBest = [0.1744e-8    0.2007e-8, -3.5048e+05, 1.3371e+07];

That fit might be close to the best to this data for this model-function. The double-peak at the end indicates that the data is more complicated than that simple function can produce. Exactly how or what is "a very difficult question". More data is required to do this experimentally, or a far more thorough modeling derivation has to be made.

The reason for using two components for the amplitude is that you have a small phase-shift of the oscillating part relative to x=0. To fit for that you need a model for that.

HTH

[Jiong Yang]

It really helps. Thank you very much!

Now I have a almost perfect fit for the first 4 peaks. I am not too worried about the double peak at the end, as the data is the inverse fourier transform of the fourier transform of the original data.