Polinomial approximation of surface

I have a massive of [x y z] coordinates of points of some surface. I need to get the polynomial approximation of this surface. Standard Matlab function polyfit does not approach like others. I have found polyfitn (+ polyvaln) function via file exchange which applies to the similar task, but it does not satisfy me.

All the functions that I have found alloy to approximate surface with formula z = f(x,y), but I need to get f(x,y,z) = 0. Approximate polynomial will be like this:

A11*x^2 + A22*y^2 + A33*z^2 + A12*x*y + A13*x*z + A23*y*z + A1*x + A2*y + A3*z + A0 = 0;

in case of 2nd order. Can anyone help me?

[dpb]

If you have the Statistics Toolbox, you can use the linear regression tool

doc linearmodel.fit

If not you can write the explicit model coefficient matrix and use the backslash operator to solve the least squares equations

doc mldivide

Answer 1

[Roger Stafford]

You can try using matlab's 'svd' function. It will find a set of coefficients A11, A22, A33, ...,A0 for which the sum of the squares of the left sides of your quadratic expression for each of your (x,y,z) points is a minimum subject to the constraint that the sum of the squares of the coefficients be one. There are many ways of defining a "best" fit but this is one definition that has an easy solution.

Let X, Y, Z be column vectors for your given points. That is, the n-th elements of the three vectors are, respectively, the x, y, and z coordinates for the n-th point. Then form the matrix M and use 'svd' to do a singular value decomposition of M:

 M = [X.^2,Y.^2,Z.^2,X.*Y,X.*Z,Y.*Z,X,Y,Z,ones(size(X))];
 [~,S,V] = svd(M,0);
 A = V(:,10);

The vector A which is the rightmost column of V will contain the desired list of coefficients, A11 A22, etc. It is a normalized eigenvector and the sum of its squares will accordingly be one. The smallest singular value that occurs at the lower right side of the main diagonal of the 10-by-10 matrix S, S(10,10), can be considered a measure of the quality of fit.

[Alexandra Nicole O'Connor]

Would this method still work if the equation were not equal to zero? I have a data set of before and after coordinates modeled by the following equation: x1 = a1 + A11*X1 + A12X2 + B111X1X1 + B112X1X2 + B121X2X1 + B122X2X2 where X is undeformed and x is deformed and I need to find the corresponding coefficients.