Derivative of spline with respect to y-values

Hello,

Consider a (x-data,y-data) spline interpolation problem. I am concerned with the problem of calculating the derivative of a spline function with respect to the perturbations of one of the y_i. This question was answered by @Bruno Luong in the following link

https://nl.mathworks.com/matlabcentral/answers/1925070-natural-cubic-spline-interpolation-of-y-values-how-to-get-derivative-of-the-spline-wrt-the-y-values#answer_1188400

However, this is very counterintuitive to me when Bruno says that "spline is linear to y values". Surprisingly, he is right but still it is hard for me to believe he is right. The following test justifies his claim:

L=-2*pi;R=2*pi;x1=linspace(L,R,10);y1=(sin(x1));x=linspace(a,R,1000);
f=@(x)interp1(x1,y1,x,'spline');
h=0.01;
Y1=y1;Y1(5)=Y1(5)+h;  % I perturb y(5)
g=@(x)(interp1(x1,Y1,x,'spline')-interp1(x1,y1,x,'spline'))./h;   % g is the answer
plot(x,f(x),'-k');hold on;plot(x1,f(x1),'r*');hold on;plot(x,g(x),'-b');  

We know that 'h' should be small. However, since Bruno says that spline is linear to the perturbations of y even when I use a very big h, I still get the same result. For instance,

h=10000;
Y1=y1;Y1(5)=Y1(5)+h;  % I perturb y(5)
g=@(x)(interp1(x1,Y1,x,'spline')-interp1(x1,y1,x,'spline'))./h;   % g is the answer
plot(x,f(x),'-k');hold on;plot(x1,f(x1),'r*');hold on;plot(x,g(x),'-b'); 

And, this proves his claim over the linearity of splines to y values. Since, I was doubting about this I also tested a different function as below:

% Sensitivity of sin(Ax) to A 
A=3;a=-2*pi;R=2*pi;x1=linspace(a,R,1000);h=0.01;
g=@(x)(sin((A+h).*x)-sin(A.*x))./h;
plot(x1,g(x1),'-k');hold on;plot(x1,x1.*cos(A.*x1),'-r');

And, for this problem it matters if you choose a bigger h, meaning that derivative of sin(Ax) w.r.t. A is xcos(Ax), not a linear relation as we know.

Do you have any intuitive explanation for this?

Of course I ham happy that Bruno is right and I was wrong. But, I find it difficult to believe this. Why? if we were concerned with Lagrange interpolation, for instance, then we know that

L(x) = sum (y_ j l_ j(x))

where {l_ 0(x), l_ 1(x), ....} are Lagrange basis polynomials. In this Lagrange interpolation problem,, it is clear that "L (x) has a linear dependence to y_ j". But, a cubic spline is more complicated than this and it is hard to believe the claimed linear relation by Bruno.

Any explanation is greatly appreciated!

Thanks in advance,

Babak

Answer 1

[John D'Errico]

Bruno is correct. Anyway, is there a good reason why he would lie? Just because something is more complicated, does not make it nonlinear. Look at a spline as a two step process. Then understand each step.

The coefficients of a spline can be written as the solution of a LINEAR systems of equations. That is, the coefficients of the spline are a linear combination of the y-values at the given points.

Next, the value of a spline, evaluated at any given point, is a LINEAR combination of the coefficients of the spline. Remember that a spline is just a cubic polynomial at heart. Well, actually a lot of little segments of cubic polynomials. But still just a cubic polynomial at heart. Do you accept that a cubic polynomial is linear in the coefficients of that cubic?

Therefore, a spline is trivially linear in the y values of your data, since each step is truly a linear combination of those y-values.

Be careful though, since it is not true that ALL curves commonly called splines are linear in y. For example, a pchip interpolant, which many think is a spline, is NOT linear in y. Close, but not exactly so, but if you are careful in how you construct the y-values, you can easily generate a set which shows the curve is not linear in y. What causes a problem in the case of pchip are the monotonicity constraints.

By the way, your "proof" of linearity is not truly a proof. It only shows that in that specific case, it appears the result behaves linearly. However, as I suggested, I can give an example where I apply a pchip interpolant to some data, and it will also seem linear, but then on another set of carefully constructed data, it will obviously not be linear. What matters is understanding how the curve is constructed. When you do, then you will also understand how to construct a counter-example of linearity for pchip.

Finally, your question goes about the derivative of a spline with respect to y, misunderstands the problem.

"And, for this problem it matters if you choose a bigger h, meaning that derivative of sin(Ax) w.r.t. A is xcos(Ax), not a linear relation as we know."

The issue is, that is NOT the derivative with respect to y. When you differentiate with respect to x, then yes, things can become nonlinear. I think you don't understand your own question.

And you don't need to believe me either. :-)

[Mohammad Shojaei Arani]

Hi John,

Thanks for your informative response!

Your response was helpful. But, still I find some issues with your argument and I will mention why. Before that let me first tell you that 4 months ago I had around 20 email exchange with Carl de Boor who has spent his entire life on splines and similar stuff. I asked him the same question I asked you. He then asked me to study his book "A practical guide to splines". He wrote to me that in the case of Lagrange interpolation the problem has linear dependence to y_j but in the case of cubic splines we do not have such a linear dependence. He never told me that for cubic splines also the answer is the same. I just mentioned this to tell you that I have spent some time in the past and have thought about this problem. Aparently this 'seemingly simple' problem seems to be very counterintuitive so that even a man whose life is about splines was not aware of this lineaer property. However, the question I asked has a lot of important applications in practical science and engineering. We really need it for optimizations. If we use this nice linear property then we can easily calculate gradient vector and Hessian matrix for our optimization problem and get a faster and more accurate solution.

You mentioned that :

Next, the value of a spline, evaluated at any given point, is a LINEAR combination of the coefficients of the spline. Remember that a spline is just a cubic polynomial at heart. Well, actually a lot of little segments of cubic polynomials. But still just a cubic polynomial at heart. Do you accept that a cubic polynomial is linear in the coefficients of that cubic?

Yes, I accept what you say. But, do you agree with me that we can use your argument for parabolic splines simply because parabolic splines are just quadratic pieces glued to gether? So, a parabolic spline is quardatic at heart. Now, let me do the same test using a parabolic spline as below:

L=-2*pi;R=2*pi;x1=linspace(L,R,10);y1=sin(x1);x=linspace(a,R,1000);

f=@(x)interp1(x1,y1,x,'spline');

H=[0.01 0.1 0.5];A=[0 0 0];i=1;

for h=H

Y1=y1;Y1(5)=Y1(5)+h; % I perturb y(5)

g=@(x)(fnval(spapi(3,x1,Y1),x)-fnval(spapi(3,x1,y),x))./h; % parabolic spline

A(i)=plot(x,g(x),'-');hold on;B=plot(x1,g(x1),'r.');

i=i+1;

end

legend([A B],{'h = 0.01','h = 0.1','h = 0.5','knots'});

How do you explain this? But, when I use cubic splines the 3 curves in above are identical, justifying the linear property.

[Bruno Luong]

@Mohammad Shojaei Arani "How do you explain this?"

Your (bug fixed) code does not reproduce what I see:

L=-2*pi;R=2*pi;x1=linspace(L,R,10);y1=sin(x1);x=linspace(L,R,1000); % Bruno change L not a
% f=@(x)interp1(x1,y1,x,'spline'); NOT used
H=[0.01 0.1 0.5];A=[0 0 0];i=1;
for h=H
    Y1=y1;Y1(5)=Y1(5)+h;  % I perturb y(5)
    g=@(x)(fnval(spapi(3,x1,Y1),x)-fnval(spapi(3,x1,y1),x))./h;   %  Bruno change y1 not y
    A(i)=plot(x,g(x),'-');hold on;B=plot(x1,g(x1),'r.');
    i=i+1;
end
legend([A B],{'h = 0.01','h = 0.1','h = 0.5','knots'});

[John D'Errico]

@Mohammad Shojaei Arani

First, I have no idea about exactly what you were talking to de Boor. There are MANY different types of splines, and you could easily have been confused. So throwing de Boor's name around is not relevant, since you clearly lack sufficient expertise in the area of mathematics to convince me that you truly understood what he told you. If I told you I once spent an afternoon talking to him about splines, would that have any relevance here? Of course not!

It is a provable fact that a simple cubic interpolating spline is linearly dependent on the y values passed to it. PROVABLE. PERIOD. I gave an outline of that proof in my answer. That is rock solid fact.

If you diverge into other classes of splines,for example, PCHIP. As I said, pchip is not linear in the function values, although in many cases, it will be! It is only when you push a pchip interpolant to certain limits, where the monotonicity constraints become active that it will become no longer linear.

Similarly, a smoothing spline, depending on how the smoothing parameter is chosen, need not be linear in the y-values, although for a fixed smoothing parameter, I think it should be. (I would need to think for a moment to be sure, and I won't do that here. Again, there are many different variations on smoothing splines.)

My point is, there are SOME splines which will be linear in the y-values, and that is a fact, pure and simple. I have been very specific here. If we are talking about the spline generated, for example, by the function spline, which is effectively what you have been using in your "tests", then it is linear in y.

As for your test using a quadratic spline, Bruno has already shown that your test is simply not valid. If you will try to make an argument, then I would suggest you make sure you understand what you are doing.