我使用MATLAB工作,而对于问题我提出我可以使用P = polyfit(X,Y,1)来估计最佳拟合线用于在板的散射数据。我想知道我可以依靠其资源,以实现符合C ++拟合算法。我知道有很多的算法这个问题,对我来说我想到的算法应该是速度快,同时能获得polyfit功能的MATLAB中的可比的精度。
我会建议从头开始对其进行编码。这是C ++一个非常简单的实现。您可以编写了两个从您的数据直接从这里公式最小二乘(同样的方法,polyfit
)截距和梯度
http://en.wikipedia.org/wiki/Simple_linear_regression#Fitting_the_regression_line
这些封闭形式的公式,你可以使用循环很容易地评价自己。如果您正在使用更高程度的千篇一律然后你描述上述这我会建议一个矩阵库或更复杂的算法,但对于简单线性回归是你所需要的。矩阵和线性代数程序会矫枉过正这样的问题(在我看来)。
此页面描述了算法的维基百科相比更简单,无需额外的步骤来计算方式等:http://faculty.cs.niu.edu/~hutchins/csci230/best-fit.htm。从那里几乎引述,在C ++中它是:
#include <vector>
#include <cmath>
struct Point {
double _x, _y;
};
struct Line {
double _slope, _yInt;
double getYforX(double x) {
return _slope*x + _yInt;
}
// Construct line from points
bool fitPoints(const std::vector<Point> &pts) {
int nPoints = pts.size();
if( nPoints < 2 ) {
// Fail: infinitely many lines passing through this single point
return false;
}
double sumX=0, sumY=0, sumXY=0, sumX2=0;
for(int i=0; i<nPoints; i++) {
sumX += pts[i]._x;
sumY += pts[i]._y;
sumXY += pts[i]._x * pts[i]._y;
sumX2 += pts[i]._x * pts[i]._x;
}
double xMean = sumX / nPoints;
double yMean = sumY / nPoints;
double denominator = sumX2 - sumX * xMean;
// You can tune the eps (1e-7) below for your specific task
if( std::fabs(denominator) < 1e-7 ) {
// Fail: it seems a vertical line
return false;
}
_slope = (sumXY - sumX * yMean) / denominator;
_yInt = yMean - _slope * xMean;
return true;
}
};
请注意,这两种算法和维基百科(http://en.wikipedia.org/wiki/Simple_linear_regression#Fitting_the_regression_line)的算法失效的情况下,点的“最佳”的描述是一条垂直线。他们失败,因为他们使用
y = k*x + b
线方程,其固有不能够描述垂直线。如果您还想覆盖情况下,当数据点是“最好”的垂直线所描述的,你需要一个线拟合算法,它使用
A*x + B*y + C = 0
线方程。您仍然可以修改当前的算法来产生等式:
y = k*x + b <=>
y - k*x - b = 0 <=>
B=1, A=-k, C=-b
在上面的代码的术语:
B=1, A=-_slope, C=-_yInt
和“然后”所述if
的块为分母等于0检查代替// Fail: it seems a vertical line
,产生下面的行式:
x = xMean <=>
x - xMean = 0 <=>
A=1, B=0, C=-xMean
我刚刚注意到,原来的文章中,我指的是已被删除。而这个网页提出了线拟合一点点不同的公式:http://hotmath.com/hotmath_help/topics/line-of-best-fit.html
double denominator = sumX2 - 2 * sumX * xMean + nPoints * xMean * xMean;
...
_slope = (sumXY - sumY*xMean - sumX * yMean + nPoints * xMean * yMean) / denominator;
该公式是因为nPoints*xMean == sumX
和nPoints*xMean*yMean == sumX * yMean == sumY * xMean
相同。
线的方程为AX + + C = 0。
因此,它可以很容易地(当B是不那么接近零)转换为Y =(-A / B)* X +(-C / B)
typedef double scalar_type;
typedef std::array< scalar_type, 2 > point_type;
typedef std::vector< point_type > cloud_type;
bool fit( scalar_type & A, scalar_type & B, scalar_type & C, cloud_type const& cloud )
{
if( cloud.size() < 2 ){ return false; }
scalar_type X=0, Y=0, XY=0, X2=0, Y2=0;
for( auto const& point: cloud )
{ // Do all calculation symmetric regarding X and Y
X += point[0];
Y += point[1];
XY += point[0] * point[1];
X2 += point[0] * point[0];
Y2 += point[1] * point[1];
}
X /= cloud.size();
Y /= cloud.size();
XY /= cloud.size();
X2 /= cloud.size();
Y2 /= cloud.size();
A = - ( XY - X * Y ); //!< Common for both solution
scalar_type Bx = X2 - X * X;
scalar_type By = Y2 - Y * Y;
if( fabs( Bx ) < fabs( By ) ) //!< Test verticality/horizontality
{ // Line is more Vertical.
B = By;
std::swap(A,B);
}
else
{ // Line is more Horizontal.
// Classical solution, when we expect more horizontal-like line
B = Bx;
}
C = - ( A * X + B * Y );
//Optional normalization:
// scalar_type D = sqrt( A*A + B*B );
// A /= D;
// B /= D;
// C /= D;
return true;
}
您还可以使用或过去this implementation还有documentation here。
直线拟合可以以不同的方式acomplished。最小二乘方法最大限度地减少平方距离的总和。但你可以采取另一种成本函数例如(未平方)的距离。但normaly使用squred距离(最小二乘)。还有一种可能性,以不同方式定义的距离。 Normaly你只需要使用“Y”轴的距离。但你也可以使用总/垂直距离。那里的距离在x方向和y方向计算的。这可能是更好的选择,如果你有同样的错误在X方向(让它成为计量数据的时间),你没有启动您在数据保存的确切时间measurment。对于最小二乘和总体最小二乘法拟合线存在封闭形式的算法。所以,如果你安装了其中的一个,你会得到与该数据点的平方距离的总和最小的线。你不能在你的defenition的培训就业处安装一个更好的线。你可以只修改定义为例服用另一种成本函数或以另一种方式定义的距离。
有很多关于进入的数据,你能想到的拟合模型的东西,但normaly它们都使用了“最小二乘法拟合线”,你应该罚款最次。但是,如果你有一个特殊的情况下,可能有必要考虑你做什么。以最小二乘法在也许几分钟内完成。什么方法最适合你的问题的思考envolves理解数学,它可以:-)采取indefinit时间。
注意:这个答案是不是回答这个问题,而是已被标记为“复制”这一项(不正确地在我看来),没办法了新的答案添加到它的这一“Line closest to a set of points”。
这个问题问的:
查找其距离所有的点是最低的线?按距离我的意思是点和线之间的最短距离。
“点和线之间的”距离的最通常的解释是欧几里德距离和“从所有点”是距离的总和中最常见的解释(以绝对值或平方值)。
当目标是最小化的平方欧氏距离的总和,线性回归(LST)是不使用的算法。另外,线性回归不会导致垂直线。要使用的算法是“总最小二乘”。举例wikipedia的问题说明,并在this answer有关详细配方数学堆栈交换见。
以适应线y=param[0]x+param[1]
简单地做到这一点:
// loop over data:
{
sum_x += x[i];
sum_y += y[i];
sum_xy += x[i] * y[i];
sum_x2 += x[i] * x[i];
}
// means
double mean_x = sum_x / ninliers;
double mean_y = sum_y / ninliers;
float varx = sum_x2 - sum_x * mean_x;
float cov = sum_xy - sum_x * mean_y;
//检查零VARx前提
param[0] = cov / varx;
param[1] = mean_y - param[0] * mean_x;
更多关于该主题http://easycalculation.com/statistics/learn-regression.php(式是相同的,只是它们相乘并通过N,样品SZ划分)。如果你想3D数据使用类似的方法,以适应平面 - http://www.mymathforum.com/viewtopic.php?f=13&t=8793
免责声明:所有的二次拟合是线性和最优的意义上,他们减少参数的噪音。但是,您可能感兴趣的不是数据的减少噪音。您可能还需要忽略异常,因为他们可以BIA是你的解决方案,极大地。这两个问题都可以用RANSAC来解决。请参阅我的文章在: