1

我正在尝试创建一个拟合平面,以便在给定 X 和 Y 的情况下预测 Z 值。我使用了以下代码(3D 数据的最佳拟合平面),但预测的 Z 值非常不正确,坦率地说,不要没什么意义。平面填充有已知的 X、Y、Z 值。最后一行然后估计每个原始 X 和 Y 的 Z 值,这与实际 Z 值完全不同。估计的 Z 值包括在下面,而实际的 Z 值(见代码)一般在 30 到 40 之间。

估计值(距离):

ZHat = 129.6104 \ -45.8558 \ 157.4270 \ 79.9675 \ 185.5349 \ 56.1384 \ -6.7733 \ 29.1776 \ -4.7795 \ 59.1381 \ -0.7739 \ 95.5122 \ 38.3086 \ 2.0756 \ -58.7012 \ 76.5445 \ -111.0525 \ -44.5676 \ 29.8922 \ 38.1766

代码:

X = [34.5,32.8,35.4,33.4,33.2,36.1,35.8,35.3,37.6,34.1,31.8,36,34.8,33.7,36,33,33.8,30.6,34.6,29.3];
Y = [33.5,35.9,33.2,34.1,32.5,34.8,35.7,35.1,35.9,34.5,35.1,34.2,34.9,35.3,36.5,34.1,37,35.6,35,34.2];
Z = [41,39,34,35,30,46,41,43,31,42,31,39,24,32,35,26,29,34,39,34];
% First, remove NaNs. All the NaN (if any) will be in the same places,
% so we need do only one test. Testing each of X,Y,Z leads to confusion
% and leads to the expectation that the NaNs are NOT in the same places
% with some potential.
k = isfinite(X);
if all(k(:))
  % in this case there were no nans, so reshape the arrays into vectors
  % While X(k) would convert an array into a column vector anyway in
  % this case, it seems far more sane to do the reshape explicitly,
  % rather than let X(k) do it implicitly. Again, a source of ambiguity
  % removed.
  X = X(:);
  Y = Y(:);
  Z = Z(:);
else
  X = X(k);
  Y = Y(k);
  Z = Z(k);
end

% Combine X,Y,Z into one array
XYZ = [X,Y,Z];

% column means. This will allow us to do the fit properly with no
% constant term needed. It is necessary to do it this way, as
% otherwise the fit would not be properly scale independent
cm = mean(XYZ,1);

% subtract off the column means
XYZ0 = bsxfun(@minus,XYZ,cm);

% The "regression" as a planar fit is now accomplished by SVD.
% This presumes errors in all three variables. In fact, it makes
% presumptions that the noise variance is the same for all the
% variables. Be very careful, as this fact is built into the model.
% If your goal is merely to fit z(x,y), where x and y were known
% and only z had errors in it, then this is the wrong way
% to do the fit.
[U,S,V] = svd(XYZ0,0);

% The singular values are ordered in decreasing order for svd.
% The vector corresponding to the zero singular value tells us
% the direction of the normal vector to the plane. Note that if
% your data actually fell on a straight line, this will be a problem
% as then there are two vectors normal to your data, so no plane fit.
% LOOK at the values on the diagonal of S. If the last one is NOT
% essentially small compared to the others, then you have a problem
% here. (If it is numerically zero, then the points fell exactly in
% a plane, with no noise.)
diag(S);

% Assuming that S(3,3) was small compared to S(1,1), AND that S(2,2)
% is significantly larger than S(3,3), then we are ok to proceed.
% See that if the second and third singular values are roughly equal,
% this would indicate a points on a line, not a plane.
% You do need to check these numbers, as they will be indicative of a
% potential problem.
% Finally, the magnitude of S(3,3) would be a measure of the noise
% in your regression. It is a measure of the deviations from your
% fitted plane.

% The normal vector is given by the third singular vector, so the
% third (well, last in general) column of V. I'll call the normal
% vector P to be consistent with the question notation.
P = V(:,3);

% The equation of the plane for ANY point on the plane [x,y,z]
% is given by
%
%    dot(P,[x,y,z] - cm) == 0
%
% Essentially this means we subtract off the column mean from our
% point, and then take the dot product with the normal vector. That
% must yield zero for a point on the plane. We can also think of it
% in a different way, that if a point were to lie OFF the plane,
% then this dot product would see some projection along the normal
% vector.
%
% So if your goal is now to predict Z, as a function of X and Y,
% we simply expand that dot product.
% 
%   dot(P,[x,y,z]) - dot(P,cm) = 0
%
%   P(1)*X + P(2)*Y + P(3)*Z - dot(P,cm) = 0
%
% or simply (assuming P(3), the coefficient of Z) is not zero...


ZHat = (dot(P,cm) - P(1)*X - P(2)*Y)/P(3)
4

0 回答 0