**You are now following this question**

- You will see updates in your followed content feed.
- You may receive emails, depending on your communication preferences.

# Curve Fitting of large Data Measurement?

20 views (last 30 days)

Show older comments

I have a measurement data which I have uploded, and I am trying to utilize a function to fit the measurement with a damping sinusoidal order of the measurement dataset. But without success. Can anyone please point me in the right direction? The curve fitting app in matlab is unable to make it. I tried a damping function in it, it didn't work. Also the lsqcurvefit nothing is working. I tried a sine function with not success too.

Thanks in advance.

##### 2 Comments

Ahmed raafat
on 23 Jan 2022

I have download your data , and plot the second and third col

the repeating sequence leads to failing the curve fitting like this (Untitled.png)

you need to add another factor

Hossam Amin
on 23 Jan 2022

Maybe I forgot to mention. The data is composed of speed and torque measurements.

You would need to add the x-axis at time using linspace command.

### Accepted Answer

Star Strider
on 24 Jan 2022

Try this —

T1 = readtable('https://www.mathworks.com/matlabcentral/answers/uploaded_files/871405/Measurement%20Data.xlsx', 'VariableNamingRule','preserve');

format longE

x = T1.Time;

yc(:,1) = T1.('Torque[Nm]');

yc(:,1) = smoothdata(yc(:,1), 'sgolay', 250); % Remove Some Noise To Make The Fit Easier

yc(:,2) = T1.('RPM[1/min]');

for k = 1:size(yc,2)

yk = yc(:,k);

y = yk;

% y = detrend(y);

ym = mean(y); % Estimate offset

y = y - ym;

yu = max(y);

yl = min(y);

yr = (yu-yl); % Range of ‘y’

yz = y-yu+(yr/2);

zci = @(v) find(diff(sign(v))); % Returns Approximate Zero-Crossing Indices Of Argument Vector

zt = x(zci(y));

per = 2*mean(diff(zt)); % Estimate period

fit = @(b,x) b(1) .* exp(b(2).*x) .* (sin(2*pi*x./b(3) + 2*pi/b(4))) + b(5); % Objective Function to fit

fcn = @(b) norm(fit(b,x) - yk); % Least-Squares cost function

[s,nmrs] = fminsearch(fcn, [yr; -10; per; -1; ym]) % Minimise Least-Squares

xp = linspace(min(x),max(x), 500);

figure

plot(x,yk,':b', 'LineWidth',1.5)

hold on

plot(xp,fit(s,xp), '--r')

hold off

grid

xlabel('Time')

ylabel('Amplitude')

legend('Original Data', 'Fitted Curve')

text(0.1*max(xlim),0.7*min(ylim), sprintf('$y = %.3f\\cdot e^{%.3E\\cdot x}\\cdot sin(2\\pi\\cdot x\\cdot %.3E%+.3f) %+.4f$', [s(1:2); 1./s(3:4); s(5)]), 'Interpreter','latex')

end

Exiting: Maximum number of function evaluations has been exceeded
- increase MaxFunEvals option.
Current function value: 7.253750

s = 5×1

1.0e+00 *
2.054070381848720e-01
-1.211851657300356e-04
3.381604644061252e+03
-1.273801702648572e+00
-6.517492492161443e-03

nmrs =

7.253749849224849e+00

Exiting: Maximum number of function evaluations has been exceeded
- increase MaxFunEvals option.
Current function value: 6954.008039

s = 5×1

1.0e+00 *
5.565560752001775e+02
-3.434517388063775e-04
3.862185655766382e+03
-6.747881574969156e-01
3.384702768264112e+02

nmrs =

6.954008038834750e+03

Because it uses the least-square approach, the fit is dominated by the highest peaks, and fits them preferentially. Experiment with tweaking ‘fit’ to get different results. (This uses a slightly-modified version of my original code to process and fit the data.)

.

##### 25 Comments

Hossam Amin
on 24 Jan 2022

Edited: Hossam Amin
on 24 Jan 2022

Thanks for the code.

Yeah I noticed that you used a smoothing function to filter the signal, and then a for loop that goes through the torque and speed measurement one after the other.

You didn't change the 'fit' function though, which begs to the question, what if I had a different measurement like the the one I newly uploaded, do I need to change the cost function (fminsearch) or the filter or the 'fit', it's also damping but the behavior is slighly different.

I also tried to chane the x-axis to represent it in seconds using the linspace instead of milisecond that exist now, it the it failed to plot correctly.

Thanks for the help. I appreciate it.

Star Strider
on 24 Jan 2022

As always, my pleasure!

I changed the code slightly and created a separate variable ‘s0’ as the initial parameter estimates. Changing the second element to allowed a decent fit.

The smoothing was only to do the initial approximation, since the period (frequency) calculations are based on the signal zero-crossings, and excessive noise creates problems in that step. The actual fitting is done on the original signal, and the original signal is plotted with the fitted function.

Changing the ‘Time’ vector from milliseconds to seconds is straightforward, and produces a decent fit with the same initial parameter estimates and correspondingly different fitted parameters. I also retreved the variable names to use in various places in the code, including titling the plots.

T1 = readtable('https://www.mathworks.com/matlabcentral/answers/uploaded_files/871810/Measurement%20Data%20_2.xlsx', 'VariableNamingRule','preserve');

vn = T1.Properties.VariableNames; % Variable Names

ttlstr = vn([2 3]); % Variable Names Used In The Code

x = T1.Time * 1E-3; % Change 'Time' From Milliseconds To Seconds

yc(:,1) = T1.(ttlstr{1});

yc(:,1) = smoothdata(yc(:,1), 'sgolay', 250); % Remove Some Noise To Make The Fit Easier

yc(:,2) = T1.(ttlstr{2});

for k = 1:size(yc,2)

yk = yc(:,k);

y = yk;

% y = detrend(y);

ym = mean(y); % Estimate offset

y = y - ym;

yu = max(y);

yl = min(y);

yr = (yu-yl); % Range of ‘y’

yz = y-yu+(yr/2);

zci = @(v) find(diff(sign(v))); % Returns Approximate Zero-Crossing Indices Of Argument Vector

zt = x(zci(y));

per = 2*mean(diff(zt)); % Estimate period

fit = @(b,x) b(1) .* exp(b(2).*x) .* (sin(2*pi*x./b(3) + 2*pi/b(4))) + b(5); % Objective Function to fit

fcn = @(b) norm(fit(b,x) - yk); % Least-Squares cost function

s0 = [yr; -1E-2; per; -1; ym]; % Initial Parameter Estimates

[s,nmrs] = fminsearch(fcn, s0) % Minimise Least-Squares

xp = linspace(min(x),max(x), 500);

figure

plot(x,yk,':b', 'LineWidth',1.5)

hold on

plot(xp,fit(s,xp), '--r')

hold off

grid

xlabel('Time (s)')

ylabel('Amplitude')

title(ttlstr{k})

legend('Original Data', 'Fitted Curve')

text(0.1*max(xlim),0.7*min(ylim), sprintf('$y = %.3f\\cdot e^{%.3f\\cdot x}\\cdot sin(2\\pi\\cdot x\\cdot %.3f%+.3f) %+.4f$', [s(1:2); 1./s(3:4); s(5)]), 'Interpreter','latex')

end

s = 5×1

0.1789
-0.1048
3.1084
-1.2119
-0.0062

nmrs = 6.6842

s = 5×1

436.6429
-0.3349
3.5279
-2.0160
274.6755

nmrs = 6.5964e+03

This version should be reasonably robust to any other data provided to it.

.

Hossam Amin
on 24 Jan 2022

That's great yes.

One last inquiry please, if I needed to determine the height "ampliude" of each oscillation peak and the distance (x-displacement) between each two peaks, in which they are all visible on the plot, how can I do so using findpeaks or any other function. I tried using 'findpeaks' on the torque measurement but it doesn't do the job perfectly.

On the speed measurement, I was able to get the peaks of the first 4. But not the distances between them.

Thanks in advance.

The code I used. Currently it is set to read the torque measurement only.

z=importdata('Measurement Data_Torque.xlsx');

ind = 1:15000;

Time = z.data(ind,1);

y_meas = z.data(ind,2);

[Ypk,Xpk,Wpk,Ppk] = findpeaks(y_meas(ind));

plot(Time,y_meas,Time(Xpk),Ypk,'dr');

n_peaks=numel(Xpk);

%Vector length of interest

t_new = Time(Xpk);

y_new = Ypk;

%---Calculate Logarithmic Decrement, undamped natural frequency, damped natural frequency, damping ratio

Log_Dec = zeros(length(n_peaks));

for nn = 1:n_peaks-1 %

Log_Dec(nn)= log(y_new(nn)/y_new(nn+1)); % computed with n = 1 periods apart

end

Mean_dec = mean(Log_Dec); %Calculate Average Logarithmic Decrement

%Damping

damp_ratio_logdec = 1/sqrt(1+((2*pi/(Mean_dec))^2)); %Assesses Damping Constant

Star Strider
on 24 Jan 2022

Try this —

T1 = readtable('https://www.mathworks.com/matlabcentral/answers/uploaded_files/871810/Measurement%20Data%20_2.xlsx', 'VariableNamingRule','preserve');

vn = T1.Properties.VariableNames; % Variable Names

ttlstr = vn([2 3]); % Variable Names Used In The Code

x = T1.Time * 1E-3; % Change 'Time' From Milliseconds To Seconds

yc(:,1) = T1.(ttlstr{1});

yc(:,1) = smoothdata(yc(:,1), 'sgolay', 250); % Remove Some Noise To Make The Fit Easier

yc(:,2) = T1.(ttlstr{2});

for k = 1:size(yc,2)

yk = yc(:,k);

y = yk;

% y = detrend(y);

ym = mean(y); % Estimate offset

y = y - ym;

yu = max(y);

yl = min(y);

yr = (yu-yl); % Range of ‘y’

yz = y-yu+(yr/2);

zci = @(v) find(diff(sign(v))); % Returns Approximate Zero-Crossing Indices Of Argument Vector

zt = x(zci(y));

per = 2*mean(diff(zt)); % Estimate period

fit = @(b,x) b(1) .* exp(b(2).*x) .* (sin(2*pi*x./b(3) + 2*pi/b(4))) + b(5); % Objective Function to fit

fcn = @(b) norm(fit(b,x) - yk); % Least-Squares cost function

s0 = [yr; -1E-2; per; -1; ym]; % Initial Parameter Estimates

[s,nmrs] = fminsearch(fcn, s0) % Minimise Least-Squares

xp = linspace(min(x),max(x), 500);

figure

plot(x,yk,':b', 'LineWidth',1.5)

hold on

plot(xp,fit(s,xp), '--r')

hold off

grid

xlabel('Time (s)')

ylabel('Amplitude')

title(ttlstr{k})

legend('Original Data', 'Fitted Curve')

text(0.1*max(xlim),0.7*min(ylim), sprintf('$y = %.3f\\cdot e^{%.3f\\cdot x}\\cdot sin(2\\pi\\cdot x\\cdot %.3f%+.3f) %+.4f$', [s(1:2); 1./s(3:4); s(5)]), 'Interpreter','latex')

end

s = 5×1

0.1789
-0.1048
3.1084
-1.2119
-0.0062

nmrs = 6.6842

s = 5×1

436.6429
-0.3349
3.5279
-2.0160
274.6755

nmrs = 6.5964e+03

dt = mean(diff(x));

for k = 1:size(yc,2)

yk = yc(:,k);

[pks{k},locs{k}] = findpeaks(yk, 'MinPeakProminence',0.1*max(yk-mean(yk)), 'MinPeakDistance',750);

ipkdist{k} = diff(locs{k});

ipktime{k} = ipkdist{k} * dt;

dipktime = mean(ipktime{k})

figure

plot(x,yk,':b', 'LineWidth',1.5)

hold on

plot(x(locs{k}),yk(locs{k}), '^r', 'MarkerFaceColor','r')

hold off

grid

xlabel('Time (s)')

ylabel('Amplitude')

title(ttlstr{k})

ipktc = compose('\\Delta t = %.3f s\n\\downarrow', ipktime{k});

text(x(locs{k}(1:numel(locs{k})-1))+dipktime/2, pks{k}(1:numel(locs{k})-1), ipktc, 'Horiz','left','Vert','middle')

end

dipktime = 3.0478

dipktime = 3.0300

The time difference calculations may be robust to other data sets, hwoever since I have only used it with this one, I cannot be certain of that. The code might require some tweaking to work with the others.

.

Hossam Amin
on 25 Jan 2022

Thank you. I appreciate it.

Star Strider
on 25 Jan 2022

As always, my pleasure!

Hossam Amin
on 25 Jan 2022

I was wondering if the fit could be better though? Even for the very first one, because I had a closer analysis onto it and there is phase shift near the end and the amplitude is not nearly fit.

The problem with my data set is that the decay rate and the logarithmic decrement is not constant. (as shown in my code I put earlier). The damping fit function you have used tries to have an approximate by keeping the decay rate and lograithmic decrement constant.

To fix that, is it possible to use two sinus function with the exponential ? I haven't found any similar example on that though. I would be glad if you have the time to help me on this regard. Appreciate it.

Star Strider
on 25 Jan 2022

There are enough data so that the number of parameters estimated from it would not be a problem, however getting a nonlinear regression optimiisation routine to converge on it could be. (For that, one of the Global Optimization Toolbox functions. specifically ga, would likely be required.)

Extending it could go something like this —

fit = @(b,x) b(1) .* exp(b(2).*x) .* (sin(2*pi*x./b(3) + 2*pi/b(4))) + exp(b(5).*x) .* (sin(2*pi*x./b(6) + 2*pi/b(7))) + b(8);

fcn = @(b) norm(fit(b,x) - yk);

with ‘fcn’ (the fitness function or objective function) being the same with the Global Optimization Toolbox as it is here.

Note — There are now 8 parameters, so this is likely beyond the ability of fminsearch (that can only estimate about 7 parameters reliably), and even lsqcurvefit could find it challenging.

.

Hossam Amin
on 26 Jan 2022

OK nice I will try working around with all that information to get the best fit I want.

Lastly, sorry to keep asking so many questions, it's a thing I have to work on since it's not my expertise (data analysis).. :D

How can come up with a code that shows the fit decay rate like this image and also prints the equation of the fit on the plot.

Thanks in advance.

Star Strider
on 26 Jan 2022

My original code prints the equation of the fit on the plot, as does this version.

Adding a plot of the exponential decay envelope simply requires a few changes to the original code —

T1 = readtable('https://www.mathworks.com/matlabcentral/answers/uploaded_files/871810/Measurement%20Data%20_2.xlsx', 'VariableNamingRule','preserve');

vn = T1.Properties.VariableNames; % Variable Names

ttlstr = vn([2 3]); % Variable Names Used In The Code

x = T1.Time * 1E-3; % Change 'Time' From Milliseconds To Seconds

yc(:,1) = T1.(ttlstr{1});

yc(:,1) = smoothdata(yc(:,1), 'sgolay', 250); % Remove Some Noise To Make The Fit Easier

yc(:,2) = T1.(ttlstr{2});

for k = 1:size(yc,2)

yk = yc(:,k);

y = yk;

% y = detrend(y);

ym = mean(y); % Estimate offset

y = y - ym;

yu = max(y);

yl = min(y);

yr = (yu-yl); % Range of ‘y’

yz = y-yu+(yr/2);

zci = @(v) find(diff(sign(v))); % Returns Approximate Zero-Crossing Indices Of Argument Vector

zt = x(zci(y));

per = 2*mean(diff(zt)); % Estimate period

expenv = @(b,x) [-1 1]*b(1) .* exp(b(2).*x) + b(5); % Exponential Envelope Decay Function

fcn = @(b) norm(fit(b,x) - yk); % Least-Squares cost function

s0 = [yr; -1E-2; per; -1; ym]; % Initial Parameter Estimates

[s,nmrs] = fminsearch(fcn, s0) % Minimise Least-Squares

xp = linspace(min(x),max(x), 500);

figure

plot(x,yk,':b', 'LineWidth',1)

hold on

plot(xp,fit(s,xp), '--r')

plot(x, expenv(s,x), '-', 'Color',[1 0.5 0], 'LineWidth',2)

hold off

grid

xlabel('Time (s)')

ylabel('Amplitude')

title(ttlstr{k})

legend('Original Data', 'Fitted Curve', 'Exponential Envelope', 'Location','best')

text(0.1*max(xlim),0.7*min(ylim), sprintf('$y = %.3f\\cdot e^{%.3f\\cdot x}\\cdot sin(2\\pi\\cdot x\\cdot %.3f%+.3f) %+.4f$', [s(1:2); 1./s(3:4); s(5)]), 'Interpreter','latex')

end

s = 5×1

0.1789
-0.1048
3.1084
-1.2119
-0.0062

nmrs = 6.6842

s = 5×1

436.6429
-0.3349
3.5279
-2.0160
274.6755

nmrs = 6.5964e+03

Make appropriate changes to get the desired result.

.

Hossam Amin
on 26 Jan 2022

Hossam Amin
on 26 Jan 2022

... because what I will try to do next is try to estimate the error between both decay curves.

Star Strider
on 26 Jan 2022

Use the envelope function to determine it (findpeaks is also an option, however the noise would make that more difficult), and then use an appropriate nonlinear parameter estimation method to estimate its parameters.

Use a slightly edited version of ‘expenv’ for that:

expenv = @(b,x) b(1) .* exp(b(2).*x) + b(3); % Exponential Envelope Decay Function

Then use it similarly to the way I originally used ‘fit’ to estimate its parameters.

.

Hossam Amin
on 26 Jan 2022

Edited: Hossam Amin
on 26 Jan 2022

All right. Sounds good.

Is there a rule of thumb in chossing the "Initial Parameter Estimates" s0 = [yr; -1E-2; per; -1; ym]?

Because it requires a lot of trial and error for the different data measurements I have, what's the idea behind trying to choose correct values intuitively to get the fit faster?

Star Strider
on 27 Jan 2022

The ‘s0’ vector elements were chosen to attempt to provide a bit of accuracy with respect to the previous calculations.

Beginning with that many random numbers will rarely result in a decent fit to the data, since there are simply too many values they could take, and too many local minima the function could fall into. For the extended version of ‘fit’ the ‘s0’ vector should be similar for the additional elements, however not the same. I would consider multiplying the initial frequency estimate by 2 (equivalent to dividing the period by 2) and perhaps using random values for the rest of the additional estimates, since the frequency of the second sin function will drive the fit of the associated parameters.

The preferred method however would be to use one of the Global Optimization Toolbox functions since they will intelligently attempt to fit the function with a large number of different parameters, and (given enough time and perhaps a few re-starts) will almost always produce an acceptable solution in a well-posed problem.

Hossam Amin
on 27 Jan 2022

I see, well I don't have license for the Global Optimization Toolbox and the trial is only 30 Days. I guess this is sufficient for my need for now.

I will now subtract the fit curve from the signal in order to calculate the error.

And in order to do that, the best way is name the fit(s0,xp) with certain variable let's say Y_fit. I changed the xp to have the same length as the original data, so xp now = linspace(min(x),max(x), length(yk)). And when I did so, it didn't affect the fit, so we are good. And so now:

Sig_diff = yk - Y_fit and I believe it checks out.

Any insights or maybe a better method you propose ?

Star Strider
on 27 Jan 2022

A simple genetic algorithm is straightforward (although not trivial) to code, and although it will lack the sophistication of the ga function, it should be good enough to get the correct result.

If your research involves a significant amount of optimisations, especially difficult optimisations, the Global Optimization Toolbox is worth getting. If the objective is to create control system models using the input-output relationships of the system, the System Identification Toolbox is definitely worth getting.

There are several ways to calculate the goodness-of-fit. I am not certain where we are just now, however if the design is to compare those two vectors, I would use the fitlm function, with one being the independent and one being the dependent variable. That will produce all the statistics necessary to evaluate how well they match, and additional calls to the ‘mdl’ result will produce additional information, including parameter confidence intervals, and confidence intervals on the regression itself. An alternative is the regress function that will also provide statistics on the fit, however fitlm will likely be better for this problem. (The fitnlm function is also an option, however that requires developing a nonlinear model for the regression.)

These functions are better than simply getting the mean of the differences or even the mean of the squared differences, and doing only simple statistics on the result.

Hossam Amin
on 27 Jan 2022

Edited: Hossam Amin
on 27 Jan 2022

Great. Thank you for all your help and information. :D

I think I am gonna go with the fitlm or fitnlm for the time being, as it's not the focus of my research to dive in complex optimization algorithms. But it's worth the time and knowledge in the future.

Star Strider
on 28 Jan 2022

As always, my pleasure!

Just to be clear, I was considering fitnlm both with respect to the model regression (use the fminsearch results as the initial parameter estimates in order to get the statistics on the parameters and the fit) as well as comparing the envelope vectors, noting that it will then be necessary to provide a model of some sort to model the envelope functions against each other (most likely a simple exponential or something similar).

I would only consider fitlm to compare the two envelope vectors. They are likely similar enough that a linear hypothesis would work for testing them.

Hossam Amin
on 28 Jan 2022

Edited: Hossam Amin
on 28 Jan 2022

Aha, I see. You mean I can use fitnlm instead of fminsearch?.. I mean I looked into it and nlinfit and fitnlm are very similar to one another, I used [beta,R,J,CovB,MSE,ErrorModelInfo] = nlinfit(X,Y,modelfun,beta0) instead of the fminsearch ... I set the following:

- initial conditions of s0 as beta0,
- modelfun is the 'fit' which the exponential sinus damping function
- x and y are time and measurement signal
- the output beta is the variables of damping; b(1), b(2) ,,etc

The curve fit slightly changed but it is still a very good fit. Plus now I have the different other outputs in regards to the Covariance and MSE, etc.

Isn't that what you meant by considering fitnlm/nlinfit with respect to the model regression?

And that fitlm can be used for the envelope using a Model which has an exponential only function. However, I don't think I need to compare the two envelopes. By intuition I believe it is just the mirror of one another.

Star Strider
on 29 Jan 2022

‘Isn't that what you meant by considering fitnlm/nlinfit with respect to the model regression?’

Yes! My original code uses fminsearch because everyone has it, and it works for this problem. The ‘extended’ version of the ‘fit’ function would not be able to use fminsearch because it has more parameters than fminsearch can comfortably estimate.

‘I need to compare the two envelopes. By intuition I believe it is just the mirror of one another.’

I thought the objective of that part of the research was to compare one part of the envelope (for example the positive one) as estimated by the regression with the (positive) signal envelope provided by the envelope function. That’s the reason I suggested a regression function that can produce a number of relevant statistics.

.

Hossam Amin
on 29 Jan 2022

Right! .. so just to be clear when you say extended version of 'fit' function you mean the fitnlm/nlinfit and not the ModelFun which in your code is called 'fit' but has no relevance to the in-build matlab 'fit' function or the extended ones which requires the Statistics and Machine Learning Toolbox.

Secondly, on the other part, that might also be interesting to look at. But can you please show me an extended code example for it in regards to the envelope. Thank you.

Star Strider
on 29 Jan 2022

I was in the middle of providing a response to this when Windows 11 crashed. Again. Taking everything I was working on with it. This also occurred with Windows 10, and I was hoping that it would be solve in Windows 11, but of course not. The instability remains. I’m an Instrument-Rated Private Piliot and am truly grateful that Micro$oft doesn’t produce airplanes!

Comparing the envelopes —

Nrows = size(T1,1)

Nrows = 31000

vn = T1.Properties.VariableNames; % Variable Names

ttlstr = vn([2 3]); % Variable Names Used In The Code

x = T1.Time * 1E-3; % Change 'Time' From Milliseconds To Seconds

yc(:,1) = T1.(ttlstr{1});

yc(:,1) = smoothdata(yc(:,1), 'sgolay', 250); % Remove Some Noise To Make The Fit Easier

yc(:,2) = T1.(ttlstr{2});

for k = 1:size(yc,2)

yk = yc(:,k);

[ehi, elo] = envelope(yk, Nrows*0.2, 'peak');

y = yk;

% y = detrend(y);

ym = mean(y); % Estimate offset

y = y - ym;

yu = max(y);

yl = min(y);

yr = (yu-yl); % Range of ‘y’

yz = y-yu+(yr/2);

zci = @(v) find(diff(sign(v))); % Returns Approximate Zero-Crossing Indices Of Argument Vector

zt = x(zci(y));

per = 2*mean(diff(zt)); % Estimate period

expenv = @(b,x) [-1 1]*b(1) .* exp(b(2).*x) + b(5); % Exponential Envelope Decay Function

fcn = @(b) norm(fit(b,x) - yk); % Least-Squares cost function

s0 = [yr; -1E-2; per; -1; ym]; % Initial Parameter Estimates

[s,nmrs] = fminsearch(fcn, s0) % Minimise Least-Squares

xp = linspace(min(x),max(x), 500);

figure

hp{1} = plot(x,yk,':b', 'LineWidth',1, 'DisplayName','Original Data');

hold on

hp{2} = plot(xp,fit(s,xp), '--r', 'Linewidth',1.5, 'DisplayName','Fitted Curve');

hp{3} = plot(x, expenv(s,x), '-', 'Color',[1 0.5 0], 'LineWidth',2, 'DisplayName','Exponential Regression Envelope');

hp{4} = plot(x, [ehi elo], '-g', 'LineWidth',1.5, 'DisplayName','Data Envelope');

hold off

grid

xlabel('Time (s)')

ylabel('Amplitude')

title(ttlstr{k})

legend([hp{1};hp{2};hp{3}(1);hp{4}(1)], 'Location','best')

% legend('Original Data', 'Fitted Curve', 'Exponential Envelope', 'Location','best')

expenvm = expenv(s,x);

envr = @(b,x) b(1).*exp(b(2).*x) + b(3);

B0 = [100; 0.01; mean(expenvm(:,2))];

envmdl{k} = fitnlm(ehi, expenvm(:,2), envr, B0);

[ernew,erci] = predict(envmdl{k}, ehi);

envmdl{k}

figure

plot(ehi, expenvm(:,2), ':b', 'LineWidth',2)

hold on

plot(ehi, [ernew erci], '-r', 'LineWidth',1.25)

hold off

grid

xlabel('Data Envelope')

ylabel('Exponential Regression Envelope')

title(ttlstr{k})

legend('Data', 'Regression', '95% Regression Confidence Intervals', 'Location','best')

end

s = 5×1

0.1789
-0.1048
3.1084
-1.2119
-0.0062

nmrs = 6.6842

ans =

Nonlinear regression model:
y ~ b1*exp(b2*x) + b3
Estimated Coefficients:

**Estimate****SE****tStat****pValue****_________****________****________****_______****b1**50.022 110.92 0.45099 0.652**b2**0.0085013 0.018827 0.45155 0.65159**b3**-50.011 110.92 -0.45088 0.65208 Number of observations: 31000, Error degrees of freedom: 30997 Root Mean Squared Error: 0.00673 R-Squared: 0.979, Adjusted R-Squared 0.979 F-statistic vs. constant model: 7.09e+05, p-value = 0s = 5×1

436.6429
-0.3349
3.5279
-2.0160
274.6755

nmrs = 6.5964e+03

ans =

Nonlinear regression model:
y ~ b1*exp(b2*x) + b3
Estimated Coefficients:

**Estimate****SE****tStat****pValue****________****__________****______****______****b1**1.2361 0.0039735 311.09 0**b2**0.010609 5.8863e-06 1802.3 0**b3**249.22 0.05228 4766.9 0 Number of observations: 31000, Error degrees of freedom: 30997 Root Mean Squared Error: 3.09 R-Squared: 0.999, Adjusted R-Squared 0.999 F-statistic vs. constant model: 1.2e+07, p-value = 0This illustrates one approach. Use whatever works in the context of your project.

At least this time Windows 11 managed to complete this without crashing!

.

Hossam Amin
on 30 Jan 2022

Oh! Sorry to hear that. Yeah I guess they should have no buisness with airplanes LOL.

OK so know I got you all right. Because it is true that the extended ModelFun or 'fit' which include two sinus functions till variable B(8) didn't work with the fminsearch, so I would need to use nlinfit or fitnlm with it "I still need to try it though".

And thanks for sharing the code and for your time.

Star Strider
on 30 Jan 2022

As always, my pleasure!

### More Answers (0)

### See Also

### Categories

### Community Treasure Hunt

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

Start Hunting!**An Error Occurred**

Unable to complete the action because of changes made to the page. Reload the page to see its updated state.

Select a Web Site

Choose a web site to get translated content where available and see local events and offers. Based on your location, we recommend that you select: .

You can also select a web site from the following list:

## How to Get Best Site Performance

Select the China site (in Chinese or English) for best site performance. Other MathWorks country sites are not optimized for visits from your location.

### Americas

- América Latina (Español)
- Canada (English)
- United States (English)

### Europe

- Belgium (English)
- Denmark (English)
- Deutschland (Deutsch)
- España (Español)
- Finland (English)
- France (Français)
- Ireland (English)
- Italia (Italiano)
- Luxembourg (English)

- Netherlands (English)
- Norway (English)
- Österreich (Deutsch)
- Portugal (English)
- Sweden (English)
- Switzerland
- United Kingdom (English)

### Asia Pacific

- Australia (English)
- India (English)
- New Zealand (English)
- 中国
- 日本Japanese (日本語)
- 한국Korean (한국어)