Quadratic interpolation with Vandermonde matrix - why the use of PINV?

Hi,

I am a beginner in Matlab and I have to understand a code.

I don't understand this part where they want to do a quadratic interpolation of unevenly spaced points into a grid

for i = 1:X
    for j = 1:Y
        closest_pts = [x,y,z, zeros(length(x),1)];
   % Find distance of each data point to the current grid point
            closest_pts(:,4)= sqrt((closest_pts(:,1)-pos(1)).^2 + (closest_pts(:,2)-pos(2)).^2);

   % Use surface fitting to find interpolated displacement at 'pos'
              Xf = closest_pts(1:num,1); Yf = closest_pts(1:num,2); Uf = closest_pts(1:num,3);
              P = zeros(num,6);
              P(:,1) = 1;
              P(:,2) = Xf;
              P(:,3) = Yf;
              P(:,4) = Xf.^2;
              P(:,5) = Yf.^2;
              P(:,6) = Xf.*Yf;
              %This is a standard method of finding the least squares
              %solution for a system of equation (see 'Linear Algebra With
              %Applications' by Otto Bretscher)
              a = pinv(P'*P)*P'*Uf;  
              interp = a(1) + a(2)*pos(1) + a(3)*pos(2) + a(4)*(pos(1)^2)+...
                  a(5)*(pos(2)^2) + a(6)*pos(1)*pos(2);

I do understand that interp has the form ax^2+bxy+cy^2+dx+ey+f=z because it is a quadratic interpolation, and I think it is thanks to Vandermonde matrix due to the construction of P...

but I don't understand : a = pinv(P'*P)*P'*Uf;

Maybe it is just gaps in my knowledge of mathematics...

Thank you for your help!!

[David Goodmanson]

Hello Aude, I don't really know why they are doing it this way either. You want to solve P*a = Uf for 'a' so you multiply both sides by P' to get

(P'*P)*a = P'*Uf.

Then you invert (P'*P) -- which is square -- using 'pinv' and take it to the other side. That's one way, which I believe is an older method. Simply using the backslash a = P\Uf also gives a least squares solution and lets Matlab decide the best way. (Which it does with the QR matrix decomposition, while pinv uses a different one, the SVD decomposition).

[Walter Roberson]

pinv() can work with singular matrices, which \ cannot.

[John D'Errico]

If this is what it tells you to do:

a = pinv(P'*P)*P'*Uf;

I find it amazing that a text that is in its 5th edition now has such a poor piece of advice. Sadly Amazon won't let me scroll down that far to see what it says to verify your claim that this is how one solves the least squares problem using matrices. And there is no way I'd spend any money at all to buy a book that has such poor advice, just to learn what it does say.

I can only hope that what is written in there in fact is something like: "While this is a method often used, it is also a BAD thing to do to a linear system of equations (for reasons listed below), and here is are some far better ways to solve the problem."

My advice is, IF it does tell you to use that solution, without explaining why it is the wrong thing to do, then return the book and get your money back.