文章目录
- RANSAC算法简介
- RANSAC算法基本思想和流程
- 迭代次数推导
- RANSAC与最小二乘区别
- RANSAC直线拟合代码（C++及Python版本）
- C++版本代码
- Python版本代码如下：
- RANSAC优缺点
- 参考
RANSAC算法简介
RANSAC(RANdom SAmple Consensus,随机采样一致)算法是从一组含有“外点”(outliers)的数据中正确估计数学模型参数的迭代算法。“外点”一般指的的数据中的噪声，比如说匹配中的误匹配和估计曲线中的离群点。所以，RANSAC也是一种“外点”检测算法。RANSAC是一个非确定性算法，在某种意义上说，它会产生一个在一定概率下合理的结果，其允许使用更多次的迭代来使其概率增加。RANSAC算最早是由Fischler和Bolles在SRI上提出用来解决LDP(Location Determination Problem，位置确定问题)问题的。
对于RANSAC算法来说一个基本的假设就是数据是由“内点”和“外点”组成的。“内点”就是组成模型参数的数据，“外点”就是不适合模型的数据。同时RANSAC假设：在给定一组含有少部分“内点”的数据，存在一个程序可以估计出符合“内点”的模型。
RANSAC算法基本思想和流程
RANSAC的思想比较简单，主要有以下几步：
- 给定一个数据集S,从中选择建立模型所需的最小样本数（空间直线最少可以由两个点确定，所以最小样本数是2，空间平面可以根据不共线三点确定，所以最小样本数为3，拟一个圆时，最小样本数是3），记选择数据集为S1
- 使用选择的数据集S1计算得到一个数学模型M1
- 用计算的模型M1去测试数据集中剩余的点，如果测试的数据点在误差允许的范围内，则将该数据点判为内点(inlier)，否则判为外点(outlier)，记所有内点组成的数据集为S1*，S1* 称作 S1的一致性集合
- 比较当前模型和之前推出的最好的模型的“内点”的数量，记录最大“内点”数量时模型参数和“内点”数量
- 重复1-4步，直到迭代结束或者当前模型已经足够好了(“内点数目大于设定的阈值”)；每次产生的模型要么因为内点太少而被舍弃，要么因为比现有的模型更好而被选用
以RANSAC直线拟合为例，如下图所示，首先在点集中随机选择两个点，求解它们构成的直线参数，再计算点集中剩余点到该直线的距离，距离小于设置的距离阈值的点为内点，统计内点个数；接下来再随机选取两个点，同样统计内点个数，以此类推；其中内点个数最多的点集即为最大一致集，最后将该最大一致集里面的点利用最小二乘拟合出一条直线。
迭代次数推导
根据上面RANSAC基本原理的介绍，在这算法流程中存在两个重要的参数需要设置，迭代次数（采样次数）和距离阈值。
迭代的次数我们应该选择多大呢？这个值是否可以事先知道应该设为多少呢？还是只能凭经验决定呢？ 这个值其实是可以估算出来的。下面我们就来推算一下。
假设内点在数据中的占比为t:
t = n i n l i e r s n i n l i e r s + n o u t l i e r s t=\frac{n_{inliers}} {n_{inliers}+n_{outliers}} t=ninliers+noutliersninliers
那么我们每次计算模型使用n个点的情况下，选取的点至少有一个外点的概率为：
1 − t n 1-t^n 1−tn
那么，在迭代k次计算模型都至少采样的一个外点的概率为：
$ (1-tn)k$
那么，能采样到正确的n个内点去计算出正模型的概率就是：
P = 1 − ( 1 − t n ) k P=1-(1-t^n)^k P=1−(1−tn)k
上式等效为：
1 − P = ( 1 − t n ) k 1-P=(1-t^n)^k 1−P=(1−tn)k
然后，对上式的两边取对数，得到：
k = l o g ( 1 − P ) l o g ( 1 − t n ) k = \frac{log(1-P)}{log(1-t^n)} k=log(1−tn)log(1−P)
内点的概率t通常是一个先验值。然后P 是我们希望RANSAC得到正确模型的概率。如果事先不知道t 的值，可以使用自适应迭代次数的方法。也就是一开始设定一个无穷大的迭代次数，然后每次更新模型参数估计的时候，用当前的内点比值当成t 来估算出迭代次数。
RANSAC与最小二乘区别
最小二乘法尽量去适应包括外点在内的所有点。因此，最小二乘法只适合与误差较小的情况。试想一下这种情况，假使需要从一个噪音较大的数据集中提取模型（比方说只有20%的数据时符合模型的）时，最小二乘法就显得力不从心了。例如下图，肉眼可以很轻易地看出一条直线（模式），但算法却找错了。
相反，RANSAC能得出一个仅仅用内点计算出模型，并且概率还足够高。但是，RANSAC并不能保证结果一定正确，为了保证算法有足够高的合理概率，必须小心的选择算法的参数（迭代次数，点与模型之间误差阈值）。经实验验证，对于包含80%误差的数据集，RANSAC的效果远优于直接的最小二乘法。
RANSAC直线拟合代码（C++及Python版本）
C++版本代码
结合OpenCV库实现了RANSAC直线拟合，在主函数中给出了一个例子，仅供大家参考。
//====================================================================//
//Program:RANSAC直线拟合，并与最小二乘法结果进行对比
//Data:2020.2.29
//Author:liheng
//Version:V1.0
//====================================================================//
#include <iostream>
#include <opencv2/opencv.hpp>//RANSAC 拟合2D 直线
//输入参数：points--输入点集
// iterations--迭代次数
// sigma--数据和模型之间可接受的差值,车道线像素宽带一般为10左右
// （Parameter use to compute the fitting score）
// k_min/k_max--拟合的直线斜率的取值范围.
// 考虑到左右车道线在图像中的斜率位于一定范围内，
// 添加此参数，同时可以避免检测垂线和水平线
//输出参数:line--拟合的直线参数,It is a vector of 4 floats
// (vx, vy, x0, y0) where (vx, vy) is a normalized
// vector collinear to the line and (x0, y0) is some
// point on the line.
//返回值：无
void fitLineRansac(const std::vector<cv::Point2f>& points,cv::Vec4f &line,int iterations = 1000,double sigma = 1.,double k_min = -7.,double k_max = 7.)
{unsigned int n = points.size();if(n<2){return;}cv::RNG rng;double bestScore = -1.;for(int k=0; k<iterations; k++){int i1=0, i2=0;while(i1==i2){i1 = rng(n);i2 = rng(n);}const cv::Point2f& p1 = points[i1];const cv::Point2f& p2 = points[i2];cv::Point2f dp = p2-p1;//直线的方向向量dp *= 1./norm(dp);double score = 0;if(dp.y/dp.x<=k_max && dp.y/dp.x>=k_min ){for(int i=0; i<n; i++){cv::Point2f v = points[i]-p1;double d = v.y*dp.x - v.x*dp.y;//向量a与b叉乘/向量b的摸.||b||=1./norm(dp)//score += exp(-0.5*d*d/(sigma*sigma));//误差定义方式的一种if( fabs(d)<sigma )score += 1;}}if(score > bestScore){line = cv::Vec4f(dp.x, dp.y, p1.x, p1.y);bestScore = score;}}
}int main()
{cv::Mat image(720,1280,CV_8UC3,cv::Scalar(125,125,125));//以车道线参数为(0.7657,-0.6432,534,548)生成一系列点double k = -0.6432/0.7657;double b = 548 - k*534;std::vector<cv::Point2f> points;for (int i = 360; i < 720; i+=10){cv::Point2f point(int((i-b)/k),i);points.emplace_back(point);}//加入直线的随机噪声cv::RNG rng((unsigned)time(NULL));for (int i = 360; i < 720; i+=10){int x = int((i-b)/k);x = rng.uniform(x-10,x+10);int y = i;y = rng.uniform(y-30,y+30);cv::Point2f point(x,y);points.emplace_back(point);}//加入噪声for (int i = 0; i < 720; i+=20){int x = rng.uniform(1,640);int y = rng.uniform(1,360);cv::Point2f point(x,y);points.emplace_back(point);}int n = points.size();for (int j = 0; j < n; ++j){cv::circle(image,points[j],5,cv::Scalar(0,0,0),-1);}//RANSAC 拟合if(1){cv::Vec4f lineParam;fitLineRansac(points,lineParam,1000,10);double k = lineParam[1] / lineParam[0];double b = lineParam[3] - k*lineParam[2];cv::Point p1,p2;p1.y = 720;p1.x = ( p1.y - b) / k;p2.y = 360;p2.x = (p2.y-b) / k;cv::line(image,p1,p2,cv::Scalar(0,255,0),2);}//最小二乘法拟合if(1){cv::Vec4f lineParam;cv::fitLine(points,lineParam,cv::DIST_L2,0,0.01,0.01);double k = lineParam[1] / lineParam[0];double b = lineParam[3] - k*lineParam[2];cv::Point p1,p2;p1.y = 720;p1.x = ( p1.y - b) / k;p2.y = 360;p2.x = (p2.y-b) / k;cv::line(image,p1,p2,cv::Scalar(0,0,255),2);}cv::imshow("image",image);cv::waitKey(0);return 0;
}
代码运行结果如下：
其中绿色直线为RANSAC直线拟合结果，红色直线为最小二乘法拟合直线结果。
Python版本代码如下：
#!/usr/bin/env python3
#coding=utf-8#============================#
#Program:RANSAC_Line.py
#
#Date:20-2-29
#Author:liheng
#Version:V1.0
#============================#import numpy as np
import random
import mathimport cv2def fitLineRansac(points,iterations=1000,sigma=1.0,k_min=-7,k_max=7):"""RANSAC 拟合2D 直线:param points:输入点集,numpy [points_num,1,2],np.float32:param iterations:迭代次数:param sigma:数据和模型之间可接受的差值,车道线像素宽带一般为10左右（Parameter use to compute the fitting score）:param k_min::param k_max:k_min/k_max--拟合的直线斜率的取值范围.考虑到左右车道线在图像中的斜率位于一定范围内，添加此参数，同时可以避免检测垂线和水平线:return:拟合的直线参数,It is a vector of 4 floats(vx, vy, x0, y0) where (vx, vy) is a normalizedvector collinear to the line and (x0, y0) is somepoint on the line."""line = [0,0,0,0]points_num = points.shape[0]if points_num<2:return linebestScore = -1for k in range(iterations):i1,i2 = random.sample(range(points_num), 2)p1 = points[i1][0]p2 = points[i2][0]dp = p1 - p2 #直线的方向向量dp *= 1./np.linalg.norm(dp) # 除以模长，进行归一化score = 0a = dp[1]/dp[0]if a <= k_max and a>=k_min:for i in range(points_num):v = points[i][0] - p1dis = v[1]*dp[0] - v[0]*dp[1]#向量a与b叉乘/向量b的摸.||b||=1./norm(dp)# score += math.exp(-0.5*dis*dis/(sigma*sigma))误差定义方式的一种if math.fabs(dis)<sigma:score += 1if score > bestScore:line = [dp[0],dp[1],p1[0],p1[1]]bestScore = scorereturn lineif __name__ == '__main__':image = np.ones([720,1280,3],dtype=np.ubyte)*125# 以车道线参数为(0.7657, -0.6432, 534, 548)生成一系列点k = -0.6432 / 0.7657b = 548 - k * 534points = []for i in range(360,720,10):point = (int((i-b)/k),i)points.append(point)# 加入直线的随机噪声for i in range(360,720,10):x = int((i-b)/k)x = random.sample(range(x-10,x+10),1)y = iy = random.sample(range(y - 30, y + 30),1)point = (x[0],y[0])points.append(point)# 加入噪声for i in range(0,720,20):x = random.sample(range(1, 640), 1)y = random.sample(range(1, 360), 1)point = (x[0], y[0])points.append(point)for point in points:cv2.circle(image,point,5,(0,0,0),-1)points = np.array(points).astype(np.float32)points = points[:,np.newaxis,:]# RANSAC 拟合if 1:[vx, vy, x, y] = fitLineRansac(points,1000,10)k = float(vy) / float(vx) # 直线斜率b = -k * x + yp1_y = 720p1_x = (p1_y-b) / kp2_y = 360p2_x = (p2_y-b) / kp1 = (int(p1_x),int(p1_y))p2 = (int(p2_x), int(p2_y))cv2.line(image,p1,p2,(0,255,0),2)# 最小二乘法拟合if 1:[vx, vy, x, y] = cv2.fitLine(points, cv2.DIST_L2, 0, 0.1, 0.01)k = float(vy) / float(vx) # 直线斜率b = -k * x + yp1_y = 720p1_x = (p1_y - b) / kp2_y = 360p2_x = (p2_y - b) / kp1 = (int(p1_x), int(p1_y))p2 = (int(p2_x), int(p2_y))cv2.line(image, p1, p2, (0, 0, 255), 2)cv2.imshow('image',image)cv2.waitKey(0)上述代码可以改进的地方：迭代求取直线参数过程中，将每次迭代的直线的内点保存下来，当迭代结束时，可以对这些内点采用最小二乘法拟合出最佳的直线参数
RANSAC优缺点
RANSAC的优点是它能鲁棒的估计模型参数。例如，它能从包含大量局外点的数据集中估计出高精度的参数。
RANSAC的缺点是它计算参数的迭代次数没有上限；如果设置迭代次数的上限，得到的结果可能不是最优的结果，甚至可能得到错误的结果。RANSAC只有一定的概率得到可信的模型，概率与迭代次数成正比。RANSAC的另一个缺点是它要求设置跟问题相关的阀值。
RANSAC只能从特定的数据集中估计出一个模型，如果存在两个（或多个）模型，RANSAC不能找到别的模型。
RANSAC可以对局外点进行剔除，这一点是比较好的，但它也并不是完美的，当它对于拟合两条平近似平行直线分布的点时，拟合出来的结果并不是最好的，甚至可以说拟合的结果是不对的。因为最终的正确结果有可能并不是经过所给的数据点的，如下图所示：
绿色线：RANSA拟合结果;红色线:最小二乘法拟合结果;蓝色线:期望的理想结果
当我们要从如图所示的数据点中拟合出一条直线的时候（即把两条离得较近的近似的平行的线当做一条线来拟合），理想的情况是蓝色线所示，但实际情况会如绿色线所示，因为蓝色线所示的点并没有经过数据点中的两个点，RANSAC是从给的数据点中进行拟合选择内点最多的情况，而实际情况可能理想直线并不经过数据点，这是这个算法的局限所在。可以考虑将RANSAC和最小二乘结合进行，这样可以结合二者的优势，得到较为理想的结果。
参考
RANSAC算法详解(附Python拟合直线模型代码)
RANSAC估计——以直线拟合为例
RANSAC拟合直线
