3D Coordinates Line of Fit

Question

David Robie 2018-6-20

0
链接

此问题的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/406619-3d-coordinates-line-of-fit

评论： Wes Anderson 2019-12-3

Hello, I have an Nx3 matrix which represents sets of coordinates in 3D space. Is there a way to calculate a line of best fit (or any type of regression) to generate an equation for approximating expected data points? Thank you!

0 个评论
显示 -2更早的评论隐藏 -2更早的评论

请先登录，再进行评论。

请先登录，再回答此问题。

Answer 1

Anton Semechko 2018-6-20

3
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/406619-3d-coordinates-line-of-fit#answer_325562

在 MATLAB Online 中打开

Hey, David. Here is a demo you may find useful:

function best_fit_3D_line_demo
% Generate sample dataset
% -------------------------------------------------------------------------
% Standard deviation of noise
s=1;
% Random direction in space
r=randn(3,1);
r=r/(norm(r)+eps);
% Random data points along r
N=1E3; % number of point samples
t=(10*s)*(2*rand(N,1)-1);
Xo=bsxfun(@times,t,r');
% Add (isotropic) Gaussian noise to Xo
X=Xo+s*randn(N,3);
% Offset X (relative to the origin) by a random amount
X=bsxfun(@plus,X,5*s*randn(1,3));
% Find line of best fit (in least-squares sense) through X
% -------------------------------------------------------------------------
X_ave=mean(X,1);            % mean; line of best fit will pass through this point  
dX=bsxfun(@minus,X,X_ave);  % residuals
C=(dX'*dX)/(N-1);           % variance-covariance matrix of X
[R,D]=svd(C,0);             % singular value decomposition of C; C=R*D*R'
% NOTES:
% 1) Direction of best fit line corresponds to R(:,1)
% 2) R(:,1) is the direction of maximum variances of dX 
% 3) D(1,1) is the variance of dX after projection on R(:,1)
% 4) Parametric equation of best fit line: L(t)=X_ave+t*R(:,1)', where t is a real number
% 5) Total variance of X = trace(D)
% Coefficient of determineation; R^2 = (explained variance)/(total variance)
D=diag(D);
R2=D(1)/sum(D);
% Visualize X and line of best fit
% -------------------------------------------------------------------------
% End-points of a best-fit line (segment); used for visualization only 
x=dX*R(:,1);    % project residuals on R(:,1) 
x_min=min(x);
x_max=max(x);
dx=x_max-x_min;
Xa=(x_min-0.05*dx)*R(:,1)' + X_ave;
Xb=(x_max+0.05*dx)*R(:,1)' + X_ave;
X_end=[Xa;Xb];
figure('color','w')
axis equal 
hold on
plot3(X_end(:,1),X_end(:,2),X_end(:,3),'-r','LineWidth',3) % best fit line 
plot3(X(:,1),X(:,2),X(:,3),'.k','MarkerSize',13)           % simulated noisy data
set(get(gca,'Title'),'String',sprintf('R^2 = %.3f',R2),'FontSize',25,'FontWeight','normal')
xlabel('X','FontSize',20,'Color','k')
ylabel('Y','FontSize',20,'Color','k')
zlabel('Z','FontSize',20,'Color','k')
view([20 20])
drawnow

9 个评论
显示 7更早的评论隐藏 7更早的评论

Anton Semechko 2018-6-21

编辑：Anton Semechko 2018-6-21

在 MATLAB Online 中打开

I modified your script to speed-up generation of guesses.

david_script;

Output:

    Number of plausible guesses = 2450
    Best fit line : L(t) = Xo + t*r, where
    Xo  = [ 60.6388 584.6318 133.9739 ]
    r   = [ -0.1015 -0.9698 -0.2216 ]
    R^2 = 0.9999

function guesses=david_script
% Receiver positions & times
p1=[-1  0  0];
p2=[ 0 -1  0];
p3=[ 1  0  1];
p4=[ 0  1  0];
p5=[ 0  0  1];
t1=0.0031;
t2=0.0056;
t3=0.0019;
t5=0.0022;
% Generate a list of plausible guesses
percent = 0.01;
[x,y,z]=deal(cell(201,1));
[x_i,y_i]=meshgrid(single(-800:800),single(-800:800));
x_i=x_i(:);
y_i=y_i(:);
z_i=ones(size(x_i),'single');
fprintf('Generating plausible guesses\n')
for i=0:200    
    [x{i+1},y{i+1},z{i+1}]=get_plausible_guesses(x_i(:),y_i(:),single(i)*z_i(:));
    fprintf('Iteration %3u/201\n',i+1)
end
x=cell2mat(x);
y=cell2mat(y);
z=cell2mat(z);
guesses=[x y z];
guessNumber=size(guesses,1); 
fprintf('\n\nNumber of plausible guesses = %u\n\n',guessNumber)
% Get parameters of best-fit-thought the guesses and visulize
% -------------------------------------------------------------------------
[Hg,Hl]=best_fit_3D_line(guesses);
% Also visulize receiver positions
% -------------------------------------------------------------------------
P={p1 p2 p3 p4 p5};
Hp=zeros(1,5);
for i=1:5
    Hp(i)=plot3(P{i}(1),P{i}(2),P{i}(3),'.','MarkerSize',30);
end
legend([Hg Hl Hp],{'guesses' 'best-fit line' 'receiver 1' 'receiver 2' 'receiver 3' 'receiver 4' 'receiver 5'})
if nargout<1, clear guesses; end
      function [x,y,z]=get_plausible_guesses(x,y,z)
      % Remove grid points based on specified contraints
      % ---------------------------------------------------------------------
      % Constraint # 1
      idx=x>=y;
      x(idx)=[];
      y(idx)=[];
      z(idx)=[];
      % Constraint # 2
      idx=z<=2*x;
      x(idx)=[];
      y(idx)=[];
      z(idx)=[];
      % Distance of the two receivers with the lowest tdoa
      D4=(x-p4(1)).^2 + (y-p4(2)).^2 + (z-p4(3)).^2;
      D3=(x-p3(1)).^2 + (y-p3(2)).^2 + (z-p3(3)).^2;
      idx=D4>=D3;
      x(idx)=[];
      y(idx)=[];
      z(idx)=[];
      D3(idx)=[];
      D4(idx)=[];
      % Distance of the receiver with the next lowest tdoa
      D5=(x-p5(1)).^2 + (y-p5(2)).^2 + (z-p5(3)).^2;
      idx=D3>=D5;
      x(idx)=[];
      y(idx)=[];
      z(idx)=[];
      D3(idx)=[];
      D4(idx)=[];
      D5(idx)=[];
      % Distance of the receiver with the next lowest tdoa
      D1=(x-p1(1)).^2 + (y-p1(2)).^2 + (z-p1(3)).^2;
      idx=D5>=D1;
      x(idx)=[];
      y(idx)=[];
      z(idx)=[];
      D1(idx)=[];
      D3(idx)=[];
      D4(idx)=[];
      D5(idx)=[];
      % Distance of the receiver with the next lowest tdoa
      D2=(x-p2(1)).^2 + (y-p2(2)).^2 + (z-p2(3)).^2;
      idx=D1>=D2;
      x(idx)=[];
      y(idx)=[];
      z(idx)=[];
      D1(idx)=[];
      D2(idx)=[];
      D3(idx)=[];
      D4(idx)=[];
      D5(idx)=[];
      % Enforce additional constraints
      % ---------------------------------------------------------------------
      [D1,D2,D3,D4,D5]=deal(sqrt(D1),sqrt(D2),sqrt(D3),sqrt(D4),sqrt(D5));
      tdoa1 = (D1-D4)/343;
      tdoa2 = (D2-D4)/343;
      tdoa3 = (D3-D4)/343;
      tdoa5 = (D5-D4)/343;
      clear D1 D2 D3 D4 D5
      % If the calculated value of any tdoa deviates more than 5% of the ideal,
      % then this is not a plausible solution
      idx= (tdoa1>=(t1+(percent*t1))) | (tdoa1<=(t1-(percent*t1))) | ...
          (tdoa2>=(t2+(percent*t2))) | (tdoa2<=(t2-(percent*t2))) | ...
          (tdoa3>=(t3+(percent*t3))) | (tdoa3<=(t3-(percent*t3))) | ...
          (tdoa5>=(t5+(percent*t5))) | (tdoa5<=(t5-(percent*t5)));
      % Plausible guesses
      x(idx)=[];
      y(idx)=[];
      z(idx)=[];
      end
end
function [Hg,Hl]=best_fit_3D_line(X)
N=size(X,1);
% Find line of best fit (in least-squares sense) through X
% -------------------------------------------------------------------------
X_ave=mean(X,1);            % mean; line of best fit will pass through this point
dX=bsxfun(@minus,X,X_ave);  % residuals
C=(dX'*dX)/(N-1);           % variance-covariance matrix of X
[R,D]=svd(C,0);             % singular value decomposition of C; C=R*D*R'
% Coefficient of determination; R^2 = (explained variance)/(total variance)
D=diag(D);
R2=D(1)/sum(D);
% Visualize X and line of best fit
% -------------------------------------------------------------------------
% End-points of a best-fit line (segment); used for visualization only
x=dX*R(:,1);    % project residuals on R(:,1)
x_min=min(x);
x_max=max(x);
dx=x_max-x_min;
Xa=(x_min-0.05*dx)*R(:,1)' + X_ave;
Xb=(x_max+0.05*dx)*R(:,1)' + X_ave;
X_end=[Xa;Xb];
figure('color','w')
axis equal
hold on
Hl=plot3(X_end(:,1),X_end(:,2),X_end(:,3),'-r','LineWidth',3); % best fit line
Hg=plot3(X(:,1),X(:,2),X(:,3),'.k','MarkerSize',13);         
set(get(gca,'Title'),'String',sprintf('R^2 = %.3f',R2),'FontSize',25,'FontWeight','normal')
xlabel('X','FontSize',20,'Color','k')
ylabel('Y','FontSize',20,'Color','k')
zlabel('Z','FontSize',20,'Color','k')
view([20 20])
drawnow
% Display line parameters
% -------------------------------------------------------------------------
fprintf('Best fit line : L(t) = Xo + t*r, where\n')
fprintf('Xo  = [ '); fprintf('%.4f ',X_ave);   fprintf(']\n') 
fprintf('r   = [ '); fprintf('%.4f ',R(:,1)'); fprintf(']\n') 
fprintf('R^2 = %.4f\n',R2)
end

Wes Anderson 2019-12-3

Hi,

I have the same issue, but I'd like my regression to go from the axis to those data points, and to get an equation out of it...any idea how to do that?

Thank you

Wes Anderson 2019-12-3

*from the origin, not the axis. Sorry.

请先登录，再进行评论。

3D Coordinates Line of Fit

0 个评论
显示 -2更早的评论隐藏 -2更早的评论

采纳的回答

9 个评论
显示 7更早的评论隐藏 7更早的评论

更多回答（0 个）

另请参阅

类别

标签

产品

版本

Community Treasure Hunt

3D Coordinates Line of Fit

0 个评论 显示 -2更早的评论隐藏 -2更早的评论

采纳的回答

9 个评论 显示 7更早的评论隐藏 7更早的评论

更多回答（0 个）

另请参阅

类别

标签

产品

版本

Community Treasure Hunt

0 个评论
显示 -2更早的评论隐藏 -2更早的评论

9 个评论
显示 7更早的评论隐藏 7更早的评论