[X,Y]=meshgrid(-3:.2:3,-3:.2:3);
Z=cos(X)+sin(Y)+3+randn(size(X))/10;
[row,col]=size(Z);
b=Z(:);A=[X(:).^3 X(:).^2 X(:) Y(:).^3 Y(:).^2 Y(:) ones(size(Y(:)))];
x=A\b;
Zhat=reshape(A*x,row,col);
subplot(2,2,1);mesh(X,Y,Z);title('original surface');
subplot(2,2,2);surf(X,Y,Zhat);shading interp; title('least mean squared method');
subplot(2,1,2);mesh(X,Y,Z-Zhat);title('error surface');