根據空間點座標擬合平面和直線(線性迴歸和svd分解)
阿新 • • 發佈:2019-02-12
第一種:如果在測量得到的資料中,x,y值都是確認沒有誤差的,而誤差只是出現在z值上,則可以使用線性迴歸的方法,此方法最小二乘的目標是在z方向上de殘差
Matlab 程式碼
- % 隨機生成一組(x,y,z),這些點的座標離一個空間平面比較近
- x0=1,L1=2;
- y0=1,L2=2;
- x=x0+rand(20,1)*L1;
- y=y0+rand(20,1)*L2;
- z=1+2*x+3*y;
- scatter3(x,y,z,'filled')
- hold on;
- X = [ones(length(x),1) x y];
- % 擬合,其實是線性迴歸,但可以用來擬合平面
- % 輸出為 b = [b(1) b(2) b(3)] 表示 z = b(1) + b(2)*x + b(3)*y 是擬合出來的平面的方程
- [b,bint,r,rint,stats] = regress(z,X,95);
- % 圖形繪製
- xfit = min(x):0.1:max(x);
- yfit = min(y):0.1:max(y);
- [XFIT,YFIT]= meshgrid (xfit,yfit);
- ZFIT = b(1) + b(2) * XFIT + b(3) * YFIT;
- mesh(XFIT,YFIT,ZFIT);
第二中: 如果在測量得到的資料中,x,y,z都存在誤差,則最小化的目標應該是測量點到平面距離的殘差
Matlab程式碼
- % 隨機生成一組(x,y,z),這些點的座標離一個空間平面比較近
- x0=1,L1=2;
- y0=1,L2=2;
- x=x0+rand(20,1)*L1;
- y=y0+rand(20,1)*L2;
- z=1+2*x+3*y;
- scatter3(x,y,z,'filled')
- hold on;
- planeData=[x,y,z];
- % 協方差矩陣的SVD變換中,最小奇異值對應的奇異向量就是平面的方向
- xyz0=mean(planeData,1);
- centeredPlane=bsxfun(@minus,planeData,xyz0);
- [U,S,V]=svd(centeredPlane);
- a=V(1,3);
- b=V(2,3);
- c=V(3,3);
- d=-dot([a b c],xyz0);
- % 圖形繪製
- xfit = min(x):0.1:max(x);
- yfit = min(y):0.1:max(y);
- [XFIT,YFIT]= meshgrid (xfit,yfit);
- ZFIT = -(d + a * XFIT + b * YFIT)/c;
- mesh(XFIT,YFIT,ZFIT);
而根據空間點擬合一條空間直線的思路比較直接,就是最小化這些點到直線的距離
- % 隨機生成一組點,這寫點距離直線l比較近,l的過點[1,1,1],方向向量為[1,2,3]
- lineData=bsxfun(@plus, [1,1,1], (-1:.1:1).'*[1,2,3]);
- Noise=rand(size(lineData))*.1;
- lineData=lineData+Noise;
- scatter3(lineData(:,1), lineData(:,2), lineData(:,3),'filled')
- hold on;
- % 擬合的直線必過所有座標的算數平均值
- xyz0=mean(lineData,1),
- % 協方差矩陣奇異變換,與擬合平面不同的是
- % 所得直線的方向實際上與最大奇異值對應的奇異向量相同
- centeredLine=bsxfun(@minus,lineData,xyz0);
- [U,S,V]=svd(centeredLine);
- direction=V(:,1);
- % 畫圖
- t=-8:0.1:8;
- xx=xyz0(1)+direction(1)*t;
- yy=xyz0(2)+direction(2)*t;
- zz=xyz0(3)+direction(3)*t;
- plot3(xx,yy,zz)