RANSAC算法

NSAC（Random Sample Consensus）是一种鲁棒的参数估计方法，特别适用于数据中包含大量离群点（outliers）的场景。其核心思想是通过反复随机采样一小部分数据来生成候选模型，然后评估该模型支持的内点数量，最终选择内点最多的模型。下面给出RANSAC的数学推导，涵盖迭代次数计算、内点判定准则以及算法收敛性分析。

​ RANSAC的核心思想是：随机选取一小部分点（最小子集）来生成候选模型，然后用整个数据集验证该模型，统计符合该模型的点（内点）数量。重复多次，选择内点最多的模型作为最终估计。

1. 问题定义与数学符号

​ 给定样本散点数据集D={xi}Ni=1�={��}�=1�，xi�� 可以表示空间点坐标或对应关系等。我们希望估计一个参数为θ∈Rd�∈��的模型（例如：直线、平面、基本矩阵等）。模型可由最小的m�个数据点唯一确定（例如直线模型2点、平面3个点）。数据集D�包含两类点：

• 内点（inlier points）：符合真实模型（在允许误差内）的数据点，假设其噪声服从某种分布（如高斯分布）。
• 离群点（outlier points）：不符合模型的点，可能来自错误测量、其他结构等。
• 内点率ϵ� ：数据中内点所占比例，未知。

RANSAC������ 的目标是找到一组参数θ∗�∗，使得内点数量最大（或某种代价函数最小）。

​ RANSAC������的核心思想是：随机选取一小部分点（最小子集）来生成候选模型，然后用整个数据集验证该模型，统计符合该模型的点（内点）数量。重复多次，选择内点最多的模型作为最终估计。

2. RANSAC算法步骤

输 入：样本数据集D� ，最小采样点数m�，距离阈值τ�，置信度p� (如0.99)，最大迭代次数Kmax���� 。

输 出：最优化模型参数θfinal������ 及内点集I∗�∗。

1. 初始化：令最优点数nbest=0�����=0, 迭代次数 k=0�=0。
2. 随机采样：从D� 中均匀随机抽取m� 个不同点构成子集Sk��。
3. 模型假设：由Sk�� 计算模型参数θk�� (通常通过最小二乘，PCA等其他方法或直接求解优化)。
4. 内点计数： 对每个xi∈D��∈�，计算距离di=d(xi,θk)��=�(��,��) 。若di≤τ��≤�，则将xi��判为内点。记内点集Ik={xi:di<τ}��={��:��<�}，内点数集nk=|Ik|��=|��|.
5. 更新最优：若nk>nbest��>����� ，则令θ∗=θk�∗=��,I∗=Ik�∗=��，nbest=nk�����=��；并根据当前nbest����� 更新所需迭代次数K�。
6. 迭代优化：k←k+1�←�+1，若k<K�<�且未达到其他终止条件，返回步骤2;
7. 最终优化：使用I∗�∗ 中的所有内点重新估计模型，得到最终参数θfinal������。

3.迭代次数K� 的概率推导

3.1 基本公式

​ 设数据中内点的真实比例为ε�。一次采用中，m� 个点全部为内点的概率为：

p=εm(1)(1)�=��

若每次采用独立，则经过K�次采用后，至少有一次全内点采样的概率（即成功概率）Psuccess��������为：

Psuccess=1−(1−p)K(2)(2)��������=1−(1−�)�

我们希望Psuccess�������� 不低于某个置信度pconf����� （如0.99）。解得所需最少迭代次数：

K=[log(1−pconf)log(1−εm)](3)(3)�=[log⁡(1−�����)log⁡(1−��)]

此公式是RANSAC������ 理论的核心。注意ε� 通常未知，需估计。

3.2 自适用迭代次数

​ 实际中，我们不知道ε� ，但可以在迭代过程中根据当前找到的最佳内点数nbest����� 来估计内点率：

^ε=nbestN(4)(4)�^=������

那么一次采样迭代次数全内点的概率估计为^p=^εm�^=�^�。为保证置信度pconf����� ，所需迭代次数为：

Knew=[log(1−pconf)log(1−^p)](5)(5)����=[log⁡(1−�����)log⁡(1−�^)]

若当前已迭代次数k<Knew�<����，则继续；否则可提前终止，这种方法避免了固定K�的盲目性。

3.3 概率下界分析

当ε� 很小，K� 可能性极大。例如，ε=0.1,m=3�=0.1,�=3，则p=0.001�=0.001,K≈log(0.01)/log(0.999)≈4603�≈log⁡(0.01)/log⁡(0.999)≈4603 。因此RANSAC������ 适用于内点率不太低的问题。

4. 模型质量度量

​ RANSAC������ 用内点数量作为模型质量的指标，这等价于一个0-1损失函数：

L(xi,θ)={0,d(xi,θ)<τ1,d(xi,θ)≥τ(6)(6)�(��,�)={0,�(��,�)<�1,�(��,�)≥�

总损失为离群点个数，最大化内点数量即最小化离群点个数。这种硬阈值损失虽然简单，但忽略了距离的具体数值。

6.最终模型优化

​ 找到最佳内点集I∗�∗ 后，为了获取更精确的参数估计，通常使用使用内点重新进行模型拟合（如最小二乘）。若内点噪声服从独立同分布的高斯分布，则最大似然估计等价于最小化平方距离和：

θfinal=argminθ∑xi∈I∗d2(xi,θ)(7)(7)������=arg⁡min�∑��∈�∗�2(��,�)

这一步显著提高了参数精度，因为最初只有m� 个点，而现在利用了所有内点。

7. RANSAC 平面拟合案例

​ RANSAC通过反复随机采样一小部分点（最小子集）来生成候选平面，然后用整个点云评估该平面的内点数量（即到平面距离小于阈值的点数）。重复多次后，选择内点最多的平面作为最佳模型，最后用所有内点重新优化平面参数（如使用PCA）。对于平面拟合，最小子集大小为 m=3（因为三点确定一个平面）。内点判定基于点到平面的垂直距离。

7.1 算法步骤

输 入：散点点云 P={pi}Ni=1�={��}�=1� ，距离阈值τ�，置信度pconf�����（通常0.99）。

输 出：最优平面参数（法向量n� 和常数d�）及内点集的索引。

1. 初始化：令最优内点数nbest=0�����=0，迭代次数k=0�=0，最大迭代次数K=∞�=∞（或是一个大数）。

2. while (k<K) do

   a.随机采样：从P�中随机抽取3个不共线的点p1,p2,p3�1,�2,�3 。

   b.计算平面模型：由三点确定平面法平面nk=(p2−p1)×(p3−p1)��=(�2−�1)×(�3−�1) ，归一化; 平面方程 nTkx+dk=0����+��=0，其中dk=−nkp1��=−���1。

   c.内点计数：对每个点pi�� ,计算距离di=|nTkpi+dk|��=|�����+��| （因为nk�� 是单位向量，此即垂直距离）。若di≤τ��≤� ，则判为内点。记内点集Ik��，内点数nk=|Ik|��=|��|。

   d.更新最优：若nk>nbest��>����� ,则更新最优模型(n∗,d∗)=(nk,dk)(�∗,�∗)=(��,��)，I∗=Ik,nbest=nk�∗=��,�����=��；并根据当前内点率^ε=nbestN�^=������ 重新计算所需的迭代次数：

   K=log(1−pconf)log(1−^ε3)(8)(8)�=log⁡(1−�����)log⁡(1−�^3)

   e.k=k+1�=�+1

3. 最终优化：使用最优内点集I∗�∗ 重新拟合平面(例如通过PCA���或最小二乘 )、得到更精确的参数(nfinal,dfinal)(������,������) 。

4. 返回（nfinal������,dfinal������）及I∗�∗。

注意：在实际应用中：需要考虑设置最大迭代次数Kmax���� 。

7.2 RANSAC拟合平面的 matlab 仿真代码

​ 以下代码生成一个带有大量离群点的平面点云，用RANSAC拟合平面，并与标准PCA拟合结果对比。

​文件一：RANSAC的平面拟合算法函数

function [a,b,c,d]=ransac_fitplane(data,N,t,T) figure;plot3(data(1,:),data(2,:),data(3,:),'o');hold on; % 显示数据 iter = N; %抽样次数N number = size(data,2); % 总数 maxNum=0; %符合拟合模型的数据的个数 for i=1:iter %循环次数 sampleidx = randperm(number); sampleidx =sampleidx(1,1:3); sample = data(:,sampleidx); %取样本 [a1,a2,a3,a4]=get_nice_plane(sample);%拟合直线方程 z=ax+by+c plane = [-a1/a3,-a2/a3,-1,-a4/a3];%建模型 mask=abs(plane*[data; ones(1,size(data,2))]); %求每个数据到拟合平面的距离 total=sum(mask<t); %计算数据距离平面小于一定阈值的数据的个数 index= mask<t; if total>T nsample=data(:,index); [a1,a2,a3,a4]=get_nice_plane(nsample); plane = [-a1/a3,-a2/a3,-1,-a4/a3]; %z=ax+by+c mask=abs(plane*[data; ones(1,size(data,2))]); total=sum(mask<t); %计算数据距离平面小于一定阈值的数据的个数 end; if total>maxNum %记录最大total maxNum=total; bestplane=plane;%最好的拟合平面 bestindex=index; bestplane2=[a1,a2,a3,a4]; end end %显示符合最佳拟合的数据 maxNum %最大一致集 avgerror=abs(bestplane*[data; ones(1,size(data,2))]); avgerror=sum(avgerror(find(avgerror<t)))/maxNum %计算一致集内平均误差 a=bestplane2(1);b=bestplane2(2);c=bestplane2(3);d=bestplane2(4); % 图形绘制 temp1=data(1,bestindex); temp2=data(2,bestindex); xfit = min(temp1):0.2:max(temp1); yfit = min(temp2):0.2:max(temp2); [XFIT,YFIT]= meshgrid (xfit,yfit); ZFIT = bestplane(1)*XFIT+bestplane(2)*YFIT+bestplane(4); mesh(XFIT,YFIT,ZFIT);grid on; xlabel('X'); ylabel('Y'); end

​文件二：利用PCA拟合平面的函数

function [a,b,c,d]=get_nice_plane(data) planeData=data'; % 协方差矩阵的SVD变换中，最小奇异值对应的奇异向量就是平面的方向 xyz0=mean(planeData,1); centeredPlane=bsxfun(@minus,planeData,xyz0); [~,~,V]=svd(centeredPlane); a=V(1,3); b=V(2,3); c=V(3,3); d=-dot([a b c],xyz0); end

​文件三：测试主程序

mu=[0 0 0]; S=[2 0 4;0 4 0;4 0 8]; data1=mvnrnd(mu,S,400); %外点 mu=[2 2 2]; S=[8 1 4;1 8 2;4 2 8]; data2=mvnrnd(mu,S,500); data=[data1',data2'];%% 相当于原始数据 [a,b,c,d]=ransac_fitplane(data,3000,0.5,300)

代码说明：

• 数据生成：内点来自真实平面加小高斯噪声，离群点在大范围内均匀随机分布，比例较高。
• RANSAC主循环：
  • 每次随机选取抽取3个点，计算平面。
  • 用所有点计算到平面的距离，统计内点。
  • 若内点创新点创新高，则更新最优模型，并根据新内点率动态调整所需迭代次数K�。
  • 循环直到达到K� 或最大迭代次数。
• 最终优化：使用找到的内点集，通过PCA��� 重新拟合平面（因为内点中有小噪声，PCA可给出更精准的法向量）。