34 views (last 30 days)

Hello everyone,

I'm a bit stuck on how the curve fitting toolbox works. I attach the data with which I am working, where "mx_inelastic_breather_fitting" is the variable y and "time_inelastic_breather_fitting" is the variable x. The fitting I want to do has the functional form y=a*cos(2*pi*f*x)*exp(-x/tp). I know that both f and tp can be found more easily if the initial values f = 470.5e9 and tp = 10e-12 are used.

Any advice?

John D'Errico
on 16 Jul 2020 at 11:26

Edited: John D'Errico
on 16 Jul 2020 at 11:30

Without even looking at your data, just your question, you need to understand how any such optimization tool works. Pretend you are a blind person, tasked with finding the location of lowest altitude on the surface of the earth. Yes, I know, it happens to be deep under water in the Pacific ocean. I'll give you scuba gear.

ll you are allowed to do is to use a cane, to infer the local clope of the surface. Then you will walk downhill. I told you to use the scuba gear! Things will get wet.

But now, suppose you were started out in some unfortunate place. Perhaps the vicinity of the Dead Sea? Thus a local depression where as you walk downhill, while you may get wet, you will also get stuck in the wrong place?

The point is, you need to use intelligent starting values for any such solver. Start it in a bad (wildly wrong) place, and you should expect to get complete crapola as a result.

DISCLAIMER: No blind individuals werre harmed in this thought experiment.

Having said all that, now we can look at your data.

plot(x,y)

Sigh. Can I suggest you learn to work in better units? Numbers with hugely varying scales will always cause you numerical problems on a computer. So instead of working with numbers that are on the order of 10^-11 seconds, use picoseconds or at least nanoseconds. Since a picosecond is 10^-12 seconds, that seems perfect here.

xp = x/1e-12;

If nanoseconds would make you happier, then use them instead.

Next, consider shifting/translating the model. Since your data starts at roughly 18, then translate the model to reflect that. Otherwise, you will be using exponentials that start at 18. That causes you to have strange values appear in the coefficient a.

Sorry though. I won't use the app here, as I cannot as easily show and explain the results.

ft = fittype('a*cos(2*pi*f*(x-x0)).*exp(-(x-x0)/tp)','indep','x','problem','x0')

ft

ft =

General model:

ft(a,f,tp,x0,x) = a*cos(2*pi*f*(x-x0)).*exp(-(x-x0)/tp)

Here, when I plot the data, I can infer intelligent start points for the parameters.

plot(xp,y)

a is the value at x == x0. That should be on the order of 1.

f and tp are rate parameters. tp is an exponential decay. 10 should be reasonable.

The most sensitive one is f, which as you have written it, should be roughly the inverse of the period, so a frequency. We can think of f as the number of oscillations per picosecond. Since the period of oscillation seems to be roughly 2, 0.5 is an intelligent start point. To get a decent estimate of the period of oscillation, I'll use the intersections tool written by Doug Schwarz, as found on the file exchange.

xpint = intersections(xp,y,[0 40],[0 0])

xpint =

18.6739787975366

19.7766155889863

20.8497137648191

21.9085954944223

22.9711077650663

24.0312152815146

25.0899375680738

26.1465450240289

27.2047855420838

28.2608584039546

29.3180979049269

30.3759822488665

31.4383534867009

32.4978593885987

33.5605930233276

34.622453287097

35.6820730748939

36.7440919564783

37.8046527675794

diff(xpint(1:2:end))

ans =

2.17573496728251

2.12139400024718

2.11882980300753

2.11484797400998

2.11331236284309

2.12025558177397

2.12223953662676

2.12148005156622

2.12257969268551

The period seems to be vaguely constant, though perhaps the first period is off from the remainder.

1/median(diff(xpint(1:2:end)))

ans =

0.471388153206563

Time to fit a model.

mdl = fit(xp,y,ft,'start',[1,.5,10],'problem',18)

mdl =

General model:

mdl(x) = a*cos(2*pi*f*(x-x0))*exp(-(x-x0)/tp)

Coefficients (with 95% confidence bounds):

a = 0.9925 (0.9156, 1.07)

f = 0.4445 (0.4405, 0.4485)

tp = 3.367 (2.971, 3.764)

Problem parameters:

x0 = 18

That seems to fit. At least, when I plot it, things seem to be somewhat reaonable.

plot(mdl)

hold on

plot(xp,y)

Now getting reasonable. This model suffers from a flaw, in that it assumes a fixed period of oscillation in that decay. Instead, your oscillations seem to be also slowing down. You could repair that, using a time varying coefficient f. Perhaps f(xp) might take the form f0 - f1*xp.

ft = fittype('a*cos(2*pi*(f0 - f1*(x - x0)).*(x-x0)).*exp(-(x-x0)/tp)','indep','x','problem','x0')

ft =

General model:

ft(a,f0,f1,tp,x0,x) = a*cos(2*pi*(f0 - f1*(x - x0)).*(x-x0)).*exp(-(x-x0)/tp)

mdl = fit(xp,y,ft,'start',[1,.5,.01,10],'problem',18)

mdl =

General model:

mdl(x) = a*cos(2*pi*(f0 - f1*(x - x0)).*(x-x0)).*exp(-(x-x0)/tp)

Coefficients (with 95% confidence bounds):

a = 0.9983 (0.9463, 1.05)

f0 = 0.4179 (0.4129, 0.4229)

f1 = -0.005473 (-0.006309, -0.004638)

tp = 3.494 (3.218, 3.771)

Problem parameters:

x0 = 18

figure

plot(mdl)

hold on

plot(xp,y)

That was a bit of a hack in the model. It does better, but we could surely do better yet with some effort. The problem is your fit will be jerked around by that first period. Since the absolute numbers are the largest there, the fit will be most impacted by errors in that first period. That would suggest a simple weighted fit might be as good, where the weight on the latter part of the series is increased.

Opportunities for recent engineering grads.

Apply Today
## 2 Comments

## Direct link to this comment

https://uk.mathworks.com/matlabcentral/answers/565526-problems-using-the-curve-fitting-app#comment_938039

⋮## Direct link to this comment

https://uk.mathworks.com/matlabcentral/answers/565526-problems-using-the-curve-fitting-app#comment_938039

## Direct link to this comment

https://uk.mathworks.com/matlabcentral/answers/565526-problems-using-the-curve-fitting-app#comment_938414

⋮## Direct link to this comment

https://uk.mathworks.com/matlabcentral/answers/565526-problems-using-the-curve-fitting-app#comment_938414

Sign in to comment.