首页 文章

是否可以通过样本数据的单次传递来进行代数曲线拟合?

提问于
浏览
0

我想对2D数据点进行曲线拟合,但由于各种原因 - 实际上不可能同时在内存中包含大量样本数据,并且遍历所有这些数据是一个昂贵的过程 .

(这样做的原因是实际上我需要根据我正在读取磁盘的千兆字节数据同时适应数千条曲线,因此它是懒散的) .

请注意,多项式系数的数量将受到限制(可能是5-10),因此精确拟合将极不可能,但这是可以的,因为我试图在随机噪声的 lot 数据中找到基础模式 . 我理解如何使用遗传算法将曲线拟合到数据集,但这需要多次通过样本数据,因此对我的应用程序来说不实用 .

有没有办法通过单次传递数据来拟合曲线,其中从样本到样本必须维持的状态是最小的?

我应该补充一点,数据的性质是点可能位于X轴上0.0到1.0之间的任何位置,但Y值总是1.0或0.0 .

所以,在Java中,我正在寻找一个具有以下接口的类:

public interface CurveFit {
   public void addData(double x, double y);
   public List<Double> getBestFit(); // Returns the polynomial coefficients
}

实现它的类不必在其实例字段中保留大量数据,即使对于数百万个数据点也不超过一千字节 . 这意味着您不能只是存储数据,因为您可以在以后通过它进行多次传递 .

编辑:有些人建议在单次通过中找到最佳曲线可能是不可能的,但是不需要最佳拟合,就像我们可以在单次通过中获得最佳曲线一样 .

如果我们有一种方法可以从一条曲线开始,然后一种方法来修改它以使它稍微靠近新数据点,这可能是一种方法的基础 - 实际上是一种梯度下降的形式 . 希望有足够的数据(并且数据充足),我们得到一个非常好的曲线 . 也许这会激励某人找到解决方案 .

8 回答

  • 2

    是的,这是一个预测 . 对于

    y = X beta + error
    

    其中小写的术语是向量,X是矩阵,你有解决方案向量

    \hat{beta} = inverse(X'X) X' y
    

    根据OLS页面 . 你 almost never 想直接计算这个,而是使用LR,QR或SVD分解 . 统计学文献中有很多参考文献 .

    如果你的问题只有一个参数(而x也是一个向量),那么这只会减少到y和x之间的交叉积的总和 .

  • 0

    如果你不介意你会得到一条直线“曲线”,那么你只需要六个变量用于任何数据量 . 这是我即将出版的书中的源代码;我确信您可以弄清楚DataPoint类的工作原理:

    Interpolation.h:

    #ifndef __INTERPOLATION_H
    #define __INTERPOLATION_H
    
    #include "DataPoint.h"
    
    class Interpolation
    {
    private:
      int m_count;
      double m_sumX;
      double m_sumXX;  /* sum of X*X */
      double m_sumXY;  /* sum of X*Y */
      double m_sumY;
      double m_sumYY;  /* sum of Y*Y */
    
    public:
      Interpolation();
    
      void addData(const DataPoint& dp);
    
      double slope() const;
      double intercept() const;
    
      double interpolate(double x) const;
      double correlate() const;
    };
    
    #endif // __INTERPOLATION_H
    

    Interpolation.cpp:

    #include <cmath>
    
    #include "Interpolation.h"
    
    Interpolation::Interpolation()
    {
      m_count = 0;
      m_sumX = 0.0;
      m_sumXX = 0.0;
      m_sumXY = 0.0;
      m_sumY = 0.0;
      m_sumYY = 0.0;
    }
    
    void Interpolation::addData(const DataPoint& dp)
    {
      m_count++;
      m_sumX += dp.getX();
      m_sumXX += dp.getX() * dp.getX();
      m_sumXY += dp.getX() * dp.getY();
      m_sumY += dp.getY();
      m_sumYY += dp.getY() * dp.getY();
    }
    
    double Interpolation::slope() const
    {
      return (m_sumXY - (m_sumX * m_sumY / m_count)) /
        (m_sumXX - (m_sumX * m_sumX / m_count));
    }
    
    double Interpolation::intercept() const
    {
      return (m_sumY / m_count) - slope() * (m_sumX / m_count);
    }
    
    
    double Interpolation::interpolate(double X) const
    {
      return intercept() + slope() * X;
    }
    
    
    double Interpolation::correlate() const
    {
      return m_sumXY / sqrt(m_sumXX * m_sumYY);
    }
    
  • 2

    为什么不使用一些固定大小的环形缓冲区(比如最后1000个点)并对缓冲数据进行基于标准QR分解的最小二乘法?一旦缓冲区填满,每次获得一个新点时,您将替换最旧的点并重新调整 . 这样你就拥有了一个有限的工作集,它仍然具有一些数据局部性,没有实时流(无记忆)处理的所有挑战 .

  • 0

    您是否限制多项式系数的数量(即拟合多项式中x的最大幂)?

    如果没有,那么你不需要“最佳拟合”算法 - 你总是可以将N个数据点完全拟合到N个系数的多项式 .

    只需使用矩阵求解N个未知数的N个联立方程(多项式的N个系数) .

    如果您限制系数的最大数量,那么您的最大值是多少?

    Following your comments and edit:

    你想要的是一个低通滤波器来滤除噪声,不适合噪声的多项式 .

  • 0

    鉴于您的数据的性质:

    这些点可能位于X轴上0.0到1.0之间的任何位置,但Y值始终为1.0或0.0 .

    然后你甚至不需要一次通过,因为这两行将完全通过每个点:

    • X = [0.0 ... 1.0],Y = 0.0

    • X = [0.0 ... 1.0],Y = 1.0

    两个短线段,单位长度和每个点落在一条线上或另一条线上 .

    不可否认,在一次通过中找到适合任意点的良好曲线的算法很有意思,但是(基于你的问题),这不是你需要的 .

  • 0

    假设您不知道哪个点应该属于哪条曲线,那么像Hough Transform这样的东西可能会提供您所需要的 .

    霍夫变换是一种允许您识别数据集内的结构的技术 . 一种用途是用于计算机视觉,其允许容易地识别视野内的线和边界 .

    这种情况的好处:

    • 每个点只需要考虑一次

    • 您不需要为每个候选行保留数据结构,只需要一个(复杂的,多维的)结构

    • 每条线的处理很简单

    • 您可以在任何时候停止并输出一组商品火柴

    • 您永远不会丢弃任何数据,因此它不依赖于任何偶然的引用位置

    • 您可以在准确度和内存要求之间进行权衡

    • 不仅限于完全匹配,也会突出显示部分匹配 .

    一种方法

    要找到三次拟合,您需要构建一个4维Hough空间,您可以在其中投影每个数据点 . 霍夫空间内的热点将通过这些点为您提供立方体的参数 .

  • 0

    您需要一个超定线性系统的解决方案 . 常用的方法是正规方程(通常不推荐),QR分解和奇异值分解(SVD) . 维基百科有不错的解释,Trefethen and Bau非常好 . 你的选择:

    • 通过正规方程实现核外实现 . 这需要产品 A'A ,其中 A 的行数多于列数(因此结果非常小) . 矩阵 A 完全由样本位置定义,因此您不必存储它,因此计算 A'A 相当便宜(如果您不需要为节点位置命中内存,则非常便宜) . 一旦计算出 A'A ,您就可以通过输入数据获得解决方案,但该方法可能不稳定 .

    • 实现核外QR分解 . 古典格兰施密特将是最快的,但你必须小心稳定性 .

    • 使用分布式内存进行内核处理(如果您有可用的硬件) . 像PLAPACK和SCALAPACK这样的库可以做到这一点,性能应该比1好得多 . 并行可伸缩性并不是很好,但是如果它是一个你甚至会考虑串行的问题大小就会很好 .

    • 使用迭代方法计算SVD . 根据系统的光谱特性(可能在预处理之后),这可以非常快速地收敛,并且不需要存储矩阵(在您的情况下,矩阵有5-10列,每列都是输入数据的大小 . 一个好的库因为这是SLEPc,你只需要找到带有矢量的Vandermonde矩阵的乘积(所以你只需要存储样本位置) . 这是非常可扩展的并行 .

  • 0

    我相信我根据this代码的修改版本找到了我自己的问题的答案 . 对于那些感兴趣的人,我的Java代码是here .

相关问题