Getting a maximum from a fitted curve (smoothing spline)?

Context: vibration experiments, plotting frequency vs displacement

What I have: data for freqency (x) and displacement (y). Plotted they form roughly a -x^2 shape. I used the curve fitting tool, with smooth spline selected to interpolate my data. The code returned was as follows:

% Fit
[xData1, yData1] = prepareCurveData( Frequency_UD, Displacement_UD );

% Set up fittype and options.
ft = fittype( 'smoothingspline' );

% Fit model to data.
[fitresult{1}, gof(1)] = fit( xData1, yData1, ft );

%plot
plot( fitresult{1},'c', xData1, yData1, 'k.');

What I want: I want to determine the resonant frequency (the point where displacement is greatest), so basically the peak of my graph. Obviously I can calculate it or visually locate it but thats not the point of the question.

For further clarity here is a picture of the graph I am working on

Answer 1

[John D'Errico]

Where does the problem lie? Are you surprised that the spline does not pass through the maximum data point? So are you just asking how to locate the maximum point, or are you worried that the smoothing spline seems to have missed that max?

Any smoothing tool will tend to reduce the peaks, increase the valleys. This is especially true of sharp peaks. Why? Because that is a location where the function is most inconsistent with the rest of the curve. So unless you have a nice, smooth, well-behaved maximum, then the smooth will tend to reduce it, to draw it down. Because it looks like it might be noise.

There are bumps in your data, clearly noise. So why should the smoothing tool not be willing to disbelieve just a little bit of that extreme point?

If you NEED to find the maximum, AND you believe that that peak is what you are looking for, then you NEED to use an interpolating spline, NOT a smoothing spline.

Ok, perhaps your question is not why does the curve not pass smoothly through the peak. Perhaps you are just trying to find out how to locate that maximal value.

How do you find a minimal point on a curve? The simple answer is to just use fminbnd, at least in one dimension. Or you could as easily use fminsearch. It will work fine here, though not quite as efficiently. But you want to find a maximum. The classic answer in optimization, as to how to find a max, is just find a minimum of the negative of your function. Multiply it by -1. Then use a minimization tool. Either fminbnd or fminsearch will suffice. (You could also use my SLM toolbox, which does have a tool in there that will locate the maximum or minimum point on a spline. But you don't need that here, and there is no reason to go that far.)

[Jacques Mclean]

no I understand the spline smooths the curve so it may miss some values. i should have been more clear:

My data is experimental, so the highest value I have recorded is not the highest possible value, that would take hours of experimentation to find. I want to extrapolate (interpolate?) what the highest value is from my data.

From what I take from your answer I need to redo my code with the interpolate tool rather than smooth spline in cftool?

[John D'Errico]

I lack your data here, so I cannot show how I might do so on your curve. But I can always make up my own data. So, lets see:

x = linspace(-3,3,24);
y = exp(-x.^2) + randn(size(x))/20;
plot(x,y,'-o')
grid on

Now, you and I both know the real peak lies at the (x,y) pair (0,1). (SHHHHSH! Don't tell anyone. It is a secret.)

But once your data becomes corrupted with noise, locating that peak in the (x,y) plane becomes non-trivial. At best, we need to make some assumptions, then work with what we have.

So, we might decide to use an interpolating spline. Simplest is to use spline here. (DO NOT USE PCHIP. A pchip interpolant will always put any extremum at a data point. It was effectively designed to do that.)

pp = spline(x,y);
[ymax,xloc] = slmpar(pp,'maxfun')
ymax =
       1.0826
xloc =
     -0.10842

I was lazy there, using my own slmpar to find the maximal point on the curve.

https://www.mathworks.com/matlabcentral/fileexchange/24443-slm-shape-language-modeling

I could have used fminbnd.

[xloc2,ymin2] = fminbnd(@(x) -ppval(pp,x),min(x),max(x))
ymin2 = -ymin2
ymin2 =
       1.0826
       
xloc2
xloc2 =
     -0.10842

Regardless we can find a maximal point on the fitted curve.

fnplt(pp)
hold on
plot(x,y,'go')
grid on
plot(xloc,ymax,'r*')

As you see, the spline interpolant found a local max that was not the same as one of the data points. It was not too much better in this case.

[xData1, yData1] = prepareCurveData( x,y);
ft = fittype( 'smoothingspline' );
fitresult = fit( xData1, yData1, ft );
plot(fitresult)
hold on
grid on
plot(x,y,'go')

[ymax,xloc] = slmpar(coeffvalues(fitresult),'maxfun')
ymax =
       1.0352
xloc =
    -0.070431

plot(xloc,ymax,'b*')

So I've found the peak of the smoothed approximant. Again, fminbnd would have sufficed as well to locate that point.

slmpar is part of my SLM toolbox.

https://www.mathworks.com/matlabcentral/fileexchange/24443-slm-shape-language-modeling

Remember that locating the true peak of a noisy curve will be only as good as your data.

[Jacques Mclean]

this was very helpful, thank you so much!

(also love your toolbox :) )

[John D'Errico]

:)