首页 文章

3D曲线拟合

提问于
浏览
15

我有 a,b 点的离散规则网格及其相应的 c 值,我进一步插值得到一条平滑的曲线 . 现在从插值数据,我还想创建一个曲线拟合的多项式方程 . 如何在多项式中拟合三维图?

我试着在MATLAB中这样做 . 我在MATLAB(r2010a)中使用了Surface拟合工具箱来曲线拟合三维数据 . 但是,如何在MATLAB / MAPLE或任何其他软件中找到适合一组数据的公式以获得最佳效果 . 有什么建议?最有用的还有一些真实的代码示例,可以在Web上查看PDF文件等 .

这只是我数据的一小部分 .

a = [ 0.001 .. 0.011];

b = [1, .. 10];

c = [ -.304860225, .. .379710865];

提前致谢 .

3 回答

  • 17

    要将曲线拟合到一组点上,我们可以使用ordinary least-squares回归 . MathWorks有一个solution page描述了这个过程 .

    举个例子,让我们从一些随机数据开始:

    % some 3d points
    data = mvnrnd([0 0 0], [1 -0.5 0.8; -0.5 1.1 0; 0.8 0 1], 50);
    

    正如@BasSwinckels所示,通过构造所需的design matrix,您可以使用mldividepinv来表示为 Ax=b Ax=b

    % best-fit plane
    C = [data(:,1) data(:,2) ones(size(data,1),1)] \ data(:,3);    % coefficients
    
    % evaluate it on a regular grid covering the domain of the data
    [xx,yy] = meshgrid(-3:.5:3, -3:.5:3);
    zz = C(1)*xx + C(2)*yy + C(3);
    
    % or expressed using matrix/vector product
    %zz = reshape([xx(:) yy(:) ones(numel(xx),1)] * C, size(xx));
    

    接下来我们可视化结果:

    % plot points and surface
    figure('Renderer','opengl')
    line(data(:,1), data(:,2), data(:,3), 'LineStyle','none', ...
        'Marker','.', 'MarkerSize',25, 'Color','r')
    surface(xx, yy, zz, ...
        'FaceColor','interp', 'EdgeColor','b', 'FaceAlpha',0.2)
    grid on; axis tight equal;
    view(9,9);
    xlabel x; ylabel y; zlabel z;
    colormap(cool(64))
    

    1st_order_polynomial


    如上所述,我们可以通过向自变量矩阵( Ax=b 中的 A )添加更多项来获得更高阶的多项式拟合 .

    假设我们想要使用常数,线性,相互作用和平方项(1,x,y,xy,x ^ 2,y ^ 2)拟合二次模型 . 我们可以手动执行此操作:

    % best-fit quadratic curve
    C = [ones(50,1) data(:,1:2) prod(data(:,1:2),2) data(:,1:2).^2] \ data(:,3);
    zz = [ones(numel(xx),1) xx(:) yy(:) xx(:).*yy(:) xx(:).^2 yy(:).^2] * C;
    zz = reshape(zz, size(xx));
    

    统计工具箱中还有一个帮助函数x2fx,它有助于为几个模型订单构建设计矩阵:

    C = x2fx(data(:,1:2), 'quadratic') \ data(:,3);
    zz = x2fx([xx(:) yy(:)], 'quadratic') * C;
    zz = reshape(zz, size(xx));
    

    最后,John D'Errico在文件交换中有一个很好的函数polyfitn,允许您指定所涉及的各种多项式阶数和术语:

    model = polyfitn(data(:,1:2), data(:,3), 2);
    zz = polyvaln(model, [xx(:) yy(:)]);
    zz = reshape(zz, size(xx));
    

    2nd_order_polynomial

  • 3

    文件交换可能有一些更好的功能,但手动执行的一种方法是:

    x = a(:); %make column vectors
    y = b(:);
    z = c(:);
    
    %first order fit
    M = [ones(size(x)), x, y];
    k1 = M\z; 
    %least square solution of z = M * k1, so z = k1(1) + k1(2) * x + k1(3) * y
    

    同样,您可以进行二阶拟合:

    %second order fit
    M = [ones(size(x)), x, y, x.^2, x.*y, y.^2];
    k2 = M\z;
    

    这似乎对您提供的有限数据集有数值问题 . 输入 help mldivide 了解更多详情 .

    要在某个常规网格上进行插值,您可以这样做:

    ngrid = 20;
    [A,B] = meshgrid(linspace(min(a), max(a), ngrid), ...
                     linspace(min(b), max(b), ngrid));
    M = [ones(numel(A),1), A(:), B(:), A(:).^2, A(:).*B(:), B(:).^2];
    C2_fit = reshape(M * k2, size(A)); % = k2(1) + k2(2)*A + k2(3)*B + k2(4)*A.^2 + ...
    
    %plot to compare fit with original data
    surfl(A,B,C2_fit);shading flat;colormap gray
    hold on
    plot3(a,b,c, '.r')
    

    可以使用下面的TryHard给出的公式来完成三阶拟合,但是当订单增加时,公式很快变得乏味 . 如果必须多次执行此操作,最好编写一个可以构造 x 给定 xyorder 的函数 .

  • 2

    这听起来像是一个哲学问题,而不是特定的实现,特别是位 - "how does one find a formula that fits a set of data to the best advantage?"在我的经验中,你必须根据你想要达到的目标做出选择 .

    什么为您定义“最佳”?对于数据拟合问题,您可以继续添加越来越多的多项式系数并获得更好的R ^ 2值......但最终会“过度拟合”数据 . 高阶多项式的缺点是样本数据范围之外的行为,你已经习惯了它适合你的响应面 - 它可以在一些狂野的方向快速消失,这可能不适合你想要模拟的任何东西 .

    您是否了解了您所适合的系统/数据的物理行为?这可以用作用于创建数学模型的方程组的基础 . 我的建议是使用你可以逃脱的最经济(简单)的模型 .

相关问题