2015-05-04 87 views
0

我有一组三维数据(300点),创建一个看起来像两个锥体或椭圆体相互连接的表面。我想要一种方法来找到这个数据集的最佳拟合椭球或圆锥的方程。回归方法并不重要,越容易越好。我基本上需要一种方法,一个代码或一个matlab函数来计算这些数据的椭圆方程的常量。如何将椭圆锥体拟合成一组数据?

回答

0

你也可以尝试fminsearch,但为了避免落在局部最小值上,考虑到系数的数量(试图消除其中一些),你将需要一个好的起点。

这里是一个2D椭圆的例子:

% implicit equation 
fxyc = @(x, y, c_) c_(1)*x.^2 + c_(2).*y.^2 + c_(3)*x.*y + c_(4)*x + c_(5).*y - 1; % free term locked to -1 

% solution (ellipse) 
c_ = [1, 2, 1, 0, 0]; % x^2, y^2, x*y, x, y (free term is locked to -1) 

[X,Y] = meshgrid(-2:0.01:2); 
figure(1); 
fxy = @(x, y) fxyc(x, y, c_); 
c = contour(X, Y, fxy(X, Y), [0, 0], 'b'); 
axis equal; 
grid on; 
xlabel('x'); 
ylabel('y');   
title('solution');   

% we sample the solution to have some data to fit 
N = 100; % samples 
sample = unique(2 + floor((length(c) - 2)*rand(1, N))); 
x = c(1, sample).'; 
y = c(2, sample).'; 

x = x + 5e-2*rand(size(x)); % add some noise 
y = y + 5e-2*rand(size(y)); 

fc = @(c_) fxyc(x, y, c_); % function in terms of the coefficients 
e = @(c) fc(c).' * fc(c); % squared error function 

% we start with a circle 
c0 = [1, 1, 0, 0, 0]; 

copt = fminsearch(e, c0) 

figure(2); 
plot(x, y, 'rx'); 
hold on 
fxy = @(x, y) fxyc(x, y, copt); 
contour(X, Y, fxy(X, Y), [0, 0], 'b'); 
hold off; 
axis equal; 
grid on; 
legend('data', 'fit'); 
xlabel('x');     %# Add an x label 
ylabel('y');   
title('fitted solution'); 

enter image description here

0

该matlab函数fit可以采取任意适合的表达式。它需要弄清楚参数,但可以完成。

您将首先创建一个fittype对象,该对象具有代表您的预期形式的字符串。你需要锻炼自己的表达是最适合你期待什么,我要带a cone expression from the Mathworld site为例,并重新安排它针对z

ft = fittype('sqrt((x^2 + y^2)/c^2) + z_0', ... 
    'independent', {'x', 'y'}, 'coeff', {'c', 'z_0'}); 

如果它是一个简单的表格MATLAB可以计算出哪些是变量,哪些是系数,但是有些更复杂的因素,您可能需要帮助。

的“fitoptions”对象持有的方法配置:根据您的数据集,你可能需要花费一些时间,指定的上限和下限,起始值等

fo = fitoptions('Upper', [one, for, each, of, your, coeffs, in, the, order, they, appear, in, the, string], ... 
    'Lower', [...], `StartPoint', [...]); 

然后得到输出

[fitted, gof] = fit([xvals, yvals], zvals, ft, fo); 

警告:我已经用2D数据集完成了这么多工作,并且docs状态可以工作三次,但是我自己没有这样做,所以上面的代码可能无法工作,请检查文档以确保您已经得到你的语法正确。

这可能值得从一个简单的拟合表达式开始,这是线性的,以便您可以使您的代码正常工作。然后将表达式换成圆锥体并玩耍,直到你得到看起来像你期望的东西。

当你合适的时候,一个好的诀窍是你可以在你的fit中使用的字符串表达式上使用eval函数来评估字符串的内容,就好像它是一个matlab表达式一样。这意味着您需要使用与您的字符串表达式中的变量和系数名称相同的工作空间变量。