RANSAC
1 概述
引用RANSAC[1]提出者的介绍“RANSAC过程与通常的平滑技术相反：不是用尽可能多的点去获得一个初始解并在以后消除无效点，RANSAC使用满足可行条件的尽量少的初始数据集并在可能时用一致性数据集扩大它”，这段话说明了RANSAC的根本思想。RANSAC提出的目的就是应对数据中的大比例外点，这在最小二乘法中是不可能的，因为，外点的存在会严重影响结果的精度。下面我就算法基本流程及参数确定、算法实现及实验结果来介绍自己对RANSAC的理解。
2 算法基本流程及参数确定
算法的步骤可以说是非常“简单”，大概有以下四步：
1 ）：确定求解问题需要的最小数据个数
2 ）：确定求解问题的最多迭代次数
3 ）：根据求解得到的模型参数，判断在此参数下的内点个数
4 ）：根据内点个数判断是否满足停止条件
由上面这看起来非常简单的步骤我们可以判断出，这里面大概需要提供4个参数，分别是求解问题的最少数据个数min_set、最大迭代次数times、内点判断阈值t 和停止判断条件 stop_num。
首先看min_set，对于不同的模型求解，这个参数是不同的。比如求解直线方程，需要知道两个点，所以min_set = 2;如果求解平面方程则需要不在同一直线上的三个点; 如果求解基本矩阵的话需要不退化的8个点 ...
其次是times，这里只要迭代次数足够大，以确保在迭代过程中由最小点集组成的样本中至少有一次没有外点的概率p足够大就行（一般取这个概率为0.99）。但是问题这就需要知道数据中内点的概率w了，这样我们才能根据前面条件的逆，也就是由最小点集组成的样本中都存在外点的概率为1-p来计算最大迭代次数。在一次迭代中含有外点的概率为1-w^min_set，求解times可由下面等式计算：
（1 - w ^ min_set）^ times = 1 - p
times = log(1 - p) / log(1 - w ^ min_set）
对于内点阈值t，我们需要知道更多的信息，比如测量误差的分布。假设我们知道了这个分布，就有办法根据假设检验计算出这个阈值t。这个地方也没搞明白，想要研究的可以参考[2]或者[1]或者自行百度吧。
最后是停止条件stop_num，这个也需要知道数据中内点的概率w，根据这个概率可以计算出数据中内点个数，当模型参数计算出的内点个数与数据集中内点个数差不多的时候就可以停止了。当然因为前面的最大迭代次数已经在大概率下保证了迭代过程中一定有一个正确的解，所以也可以计算出所有模型参数下的内点个数，取最大的就行。
由上面的分析可以看出来，RANSAC其实并不像看起来那么简单，它需要我们知道数据的分布情况，而这又是非常困难的。不过现在也有了许多改进的RANSAC算法，参数更少，但是还需要知道部分数据分布情况。对于这个缺点，有想法的朋友欢迎讨论！
3 算法实现及实验结果
在这里，我先给定一条直线y = 3 * x + 4，生成一些准确的数据点，然后再生成一些外点，最后使用RANSAC来计算模型参数。这里的算法参数都是随便给的，也试过按照[2]，噪声点是在原数据的基础上加上均值为0，标准方差为1的高斯噪声，然后计算阈值。但是这样计算出的效果并不明显，基本上每次都是在所有的模型参数下数据集中的点都被判定为内点。后来是随便生成了一些外点，具体可以看程序。
#include <iostream>
#include <vector>
#include <math.h>
#include <stdlib.h> //生成随机点
#include <time.h> //根据系统时间生成随机点种子
#include <random> //高斯噪声生成
#include <fstream>
#define random(a) (rand()%a) //0-a
struct DATA{
int x;
int y;
int flag;
};
class RANSAC{
private:
std::vector<DATA> _data; //输入数据
int _min_size; //最小子集
float _ratio; //数据外点比例，不知道的情况下默认0.5
float _p; //至少一次采样，没有外点的概率,默认0.99
std::vector<float> infor_entro; //信息熵
std::vector<DATA>minset; //每次采样的最小子集
public:
int times; //迭代次数
std::vector<std::vector<float> > param; //计算得到的模型参数
std::vector<int> inliners; //每个模型参数的内点个数
RANSAC();
RANSAC(std::vector<DATA>& data, int min_size, float p = 0.99, float ratio = 0.5){
_data = data;
_min_size = min_size;
_ratio = ratio;
_p = p;
}
void iterator();
void compute_infor_entro();
};
/**
* 计算迭代次数，最小集采样，计算模型参数，计算内点个数
* 这里默认当进行距离阈值判断时，点为内点的概率为0.95,并且默认噪声服从均值为0,标准差为1的高斯分布
*/
void RANSAC::iterator() {
times = (int)(log(1-_p)/log(1-pow((1-_ratio),_min_size)));//迭代次数
minset.resize(_min_size);
std::vector<DATA> temdata;
param.resize(times);
for(int j = 0; j < times; j++)
{
temdata = _data;
srand((unsigned)time(0) + j * 10);
for(int i = 0; i < _min_size; i ++)//选取最小点集
{
int local = random(temdata.size());
minset[i] = temdata[local];
temdata[local] = temdata.back();
temdata.pop_back();
}
//计算模型参数
param[j].push_back((float)(minset[1].y - minset[0].y) / (minset[1].x - minset[0].x));
param[j].push_back(minset[0].y - param[j][0] * minset[0].x) ;
//计算每个参数对应的内点个数
int inlinepoint = 0;
for(int i = 0; i < _data.size(); i ++)
{
if(fabs(param[j][0]) > 100000 || fabs(param[j][1]) > 10000)//防止出现无穷大数
break;
if(abs(_data[i].x * param[j][0] - _data[i].y + param[j][0])/sqrt(pow(param[j][0],2) + 1) <= 3.84)//距离阈值判断，按照高斯分布且判断内点概率0.95的卡方检验设置
inlinepoint ++;
}
inliners.push_back(inlinepoint);
}
}
/**
* 根据给定模型生成数据
* @param data 存放生成数据
* @param good 准确的数据个数
* @param bad 加噪声的数据个数
*/
void generatordata(std::vector<DATA>& data, int good = 50, int bad = 50)
{
data.resize(good + bad);
srand((unsigned)time(NULL));
for(int i = 0; i < good; i ++) //生成准确点
{
data[i].x = random(good);
data[i].y = 3 * data[i].x + 4;
data[i].flag = 1;
}
double mean = 0;
double stddev = 1;
std::default_random_engine generator;
std::normal_distribution<double> dist(mean, stddev);
//srand((unsigned)time(0) + good*2);
for(int i = good; i < good + bad; i ++)//噪点
{
data[i].x = random(good);
data[i].y = random(good * 3);
//data[i].y = 3 * data[i].x + 4 + dist(generator);
///如果按照测量误差为均值0 标准差1来生成数据的话，大多数情况下所有点都被判别为内点，这样就和处理大量外点情况不太相符
data[i].flag = 0;
}
}
int main ()
{
std::vector<DATA> data;
generatordata(data);
RANSAC it1(data,2);
it1.iterator();
//存放数据
std::ofstream of("Ransac.txt");
of << "data: " << std::endl;
for(int i = 0; i < data.size(); i ++)
{
of << data[i].x << /*" " << data[i].y <<*/ std::endl;
}
for(int i = 0; i < data.size(); i ++)
{
of /*<< data[i].x << " " */<< data[i].y << std::endl;
}
of << "拟合参数： " <<std::endl;
for(int i = 0; i < it1.times; i ++)
{
of << it1.param[i][0] /*<< " " << it1.param[i][1] */<< std::endl;
}
for(int i = 0; i < it1.times; i ++)
{
of /*<< it1.param[i][0] << " "*/ << it1.param[i][1] << std::endl;
}
//输出参数及内点个数
for(int i = 0; i < it1.times; i ++)
{
std::cout << "y = " << it1.param[i][0] << "x + " << it1.param[i][1] << std::endl;
std::cout << "内点个数： " << it1.inliners[i] << std::endl << std::endl;
}
return 0;
}测试发现，大部分情况下计算还是比较准确的，但是有个别情况有少许偏差，下面是其中一个存在偏差的测试图。
其中蓝色点是由直线y = 3 * x + 4 生成的，姜黄色点是随即生成的外点。数据中内点、外点各50个，灰色点是RANSAC计算出来的内点数最多的模型 y = 3.09091 * x 生成的。它的内点数比y = 3 * x + 4 多几个。
4 参考资料
[1] Fischler M A, Bolles R C. Random sample consensus: a paradigm for model fitting with applications to image analysis and automated cartography[J]. Communications of the ACM, 1981, 24(6): 381-395.
[2] 计算机视觉中的多视图几何 P73（3.7 鲁棒估计）
5293
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
