Robust estimation 시리즈 (1) - RANSAC이란?

논문 소개

RANSAC은 1981년에 Fischler와 Bolles가 제안한 Model fitting 방법론이다. RANSAC의 이름은 RANdom SAmple Consensus를 줄여서 만들어졌다. 논문이 처음 제안했던 방법도 컴퓨터 비전 분야에서의 사용 방법을 제안하였지만, 크게 주목을 받지 못하다가 2004년 Nister가 제안하는 이미지 프레임간 relative motion을 robust하게 찾아내는 5-point algorithm + RANSAC 프레임워크를 통해 주목을 받게 되었다.

그 후, 많은 양의 데이터로부터 outlier를 거르고 model fitting을 하는 방식으로 RANSAC이 많이 사용되고 있다. Fundamental matrix, Essential matrix, Homography, PnP, Rigid point cloud estimation 등등에서 사용된다.

오리지널 논문 링크: Fischler 1981 - Random Sample Consensus: A Paradigm for Model Fitting with Applications to Image Analysis and Automated Cartography

 


RANSAC 알고리즘 이해하기

RANSAC의 알고리즘을 사실 굉장히 간단하다.

 

간단한 Linear regression 예시

우선 간단하게 Linear regression으로 예시를 들어보자. Linear regression이 무엇인지 잘 모르겠다면, 어떤 점들이 흩뿌려져 있고, 그 점들의 중앙을 가장 잘 통과하는 라인이 어떤 라인인지 찾는 것이라고 보면 된다.

RANSAC에서는 우선 제일 먼저 Minimal set 크기의 무작위 데이터 샘플링과 함께 시작한다. Minimal set 데이터, 우리가 만들고 싶어하는 모델을 만드는데 필요한 최소한의 데이터 갯수를 의미한다. Linear regression에서는 ‘선’ 이라는 모델을 만들려고 하니까, ‘점’ 두개가 필요하겠다. 이번 예시에서는 점 하나는 무조건 그래프의 (0,0)으로 고정시켜서 문제를 간단하게 만들기로 한다.

무작위로 데이터를 뽑은 결과, 오른쪽 아래 있는 데이터가 뽑혔다. 이 데이터와 (0,0)에 고정된 위치를 이어서 ‘선’을 그렸다. 이때 우리는 이 데이터를 inlier라고 가정한다. 기존의 RANSAC 방법을 그대로 따르면 이 선 위에 올라가는 다른 데이터의 수를 고르면 된다. 하지만 로보틱스 기술 특성상, 모든 센서에는 노이즈가 포함되기 때문에, 이 선 위에 곧바로 올라가는 데이터는 거의 없을 것이다. 그러므로 어느정도 데이터에 여지를 줘서, 만큼의 거리 안에 들어오는 데이터의 수를 센다. 3개의 데이터가 이 공간에서 발견되었다.

그러면 노트에다가 ‘현재까지 가장 잘 나온 모델은… 1. 샘플링 된 데이터 (i.e. 모델 파라미터로 쓰인 데이터), 2. inlier의 수’를 기록한다.

Attempt 1 - 3개의 inlier가 발견됨.

이 후, 다른 데이터를 다시 무작위 샘플링하고, 아까와 같이 모델 (i.e. ‘선’)을 그린 다음에, 만큼의 공간 안에서 데이터의 수를 센다. 이때, 새로 센 inlier의 수가 이전에 찾은 inlier의 수 보다 많으면, ‘현재까지 제일 잘 나온 모델’의 정보에서 이전의 기록을 삭제하고 새 데이터로 덮어씌운다.

Attempt 2 - 6개의 inlier가 발견됨. Attempt 1의 3개 inlier보다 많으므로, Attempt 2가 ‘현재까지 제일 잘 나온 모델’이 된다.

이 작업을 계속해서 반복하다보면, ‘가장 많은 inlier를 가진 모델‘을 계산할 수 있다.

논문에서도 보여주듯이, 수학적으로도 몇번의 샘플을 돌려야 가장 좋은 모델을 찾을 수 있는지에 대한 수식이 있다.

  • 는 샘플링의 횟수
  • 는 우리가 고른 데이터가 inlier일 확률
  • 는 전체 데이터의 inlier : outlier 비율
  • 는 minimal set을 고르기 위한 데이터의 갯수

이 식은 두가지 방법으로 사용할 수 있다.

  1. Inlier를 얻어내고 싶은 확률을 설정해서, 샘플링의 횟수를 알아내기
  2. 샘플링의 횟수를 설정해서, Inlier를 얻어낼 확률을 설정하기.

정확도를 위해서는 전자를, Real-time 환경에서의 속도 안정성을 위해서는 후자를 선택하면 된다.

를 대략 95% 또는 98%가 나올 수 있게 샘플링의 수를 맞춰두면, 왠만하면 실패하지 않고 잘 추정해내는 것을 볼 수 있다.

 

컴퓨터 비전 예시

이제 컴퓨터 비전에 RANSAC을 적용해본 예시를 들어본다.

Homography matrix를 구하기 위해서는, 두개의 이미지에서 공유하는 4쌍의 feature pair가 필요하다. 이 때 각각의 feature pair는 삼각측량을 했을 때 3D point를 구상해야하므로 총 4개의 3D point가 생기는데, 이 4개의 3D point 들은 하나의 plane 위에 있어야한다. 여기서 Homography matrix를 사용하면 하나의 이미지 위 4개의 feature들의 위치를 다른 이미지 위 4개의 feature들의 위치를 알 수 있다. 또, 여기서는 사용하지 않지만 이 homography matrix를 SVD 등을 통해 decompose하면 Rotation가 Translation을 얻을 수 있다.

Homography matrix를 구하기 위해서는, 계산에 들어가는 feature pair에 outlier가 없다는 것이 전제이다.

Homography 의 원리. 4쌍의 correspondence가 만드는 plane 위의 4개의 3D point.

두개의 이미지 사이의 homography를 구해보기로 하자. 두개의 이미지를 준비했고, 양쪽 이미지에서 feature detection 및 descriptor 추출을 했다. FLANN을 사용해서 descriptor matching을 했고, 총 100개의 feature matching이 만들어졌다. 하지만, 이들 중 분명 몇개는 잘못된 feature matching이 있을 것이며, 이런 feature pair가 homography 계산에 들어가면 완전히 잘못된 결과가 나올 것이다.

5-point algorithm을 위해 두 이미지에서 feature를 미리 뽑아놓음

사실 RANSAC을 쓰지 않고도 homography 계산을 하는 방법은 있다. 100개의 피쳐매칭 중 어떤 피쳐매칭을 써야할지 모르겠을 때는, 존재하는 모든 경우의 수를 전부 계산한 후 가장 에러가 적은 방식을 고르는 것도 한가지 방법이다. 하지만 이 방법은 너무 오래 걸리기 때문에 실시간 컴퓨터 비전에서 돌리기에는 적절한 방법이 아니다.

존재하는 모든 경우의 수를 전부 계산하는 것도 하나의 방법이지만 너무 오래걸린다.

그러므로, 우리는 RANSAC을 사용해 볼 것이다. 우선, 랜덤하게 minimal set 데이터를 (i.e. feature pair 4개를) 고른다. 그리고, 이 데이터로 모델 추정을 한다 (i.e. Homography 계산을 한다).

Random sampling

Homography matrix가 나왔으니, Image 1에 있는 feature들의 위치에 해당 homography matrix를 곱해주면, image 2에서 해당 feature들이 어디에서 보일지 예측할 수 있다. Image 1에 보이는 100개의 feature들의 위치에 homography matrix를 전부 곱해줘서, image 2에서 어디에 나올지 예상 위치를 찾는다.

Image 1의 feature 위치에 Homography를 곱하면 Image 2의 feature 위치를 알 수 있다.

그리고 그 예상되는 위치와 실제 feature matching이 된 위치의 차이를 (i.e. residual) 구한다. 이 때, 각각의 feature 위치마다의 에러 값을 residual 로 고려하며, SSE를 이용할 경우 을 계산한다.

또, 설계자의 취향에 따라 특정 threshold 미만의 residual 값은 무시하는 방법도 있다. 이는 센서 noise 값을 고려해서 noise 제거와 outlier 제거를 명확히 구분하는 방법 중 하나이다. 이 residual 값들을 전부 더해서 SSE 값을 구하고, ‘현재까지 가장 좋은 모델 값’으로 저장한다.

‘현재까지 가장 좋은 모델 값’에 대한 업데이트가 끝나면, 다시 랜덤하게 또 샘플링을 한 후, 위의 과정을 반복한다.

Random sampling

반복하는 과정에서 좋은 모델이 나오면 SSE의 값이 적을 것이다. 현재의 모델이 이전의 모델보다 SSE가 작을 경우, 이전의 모델은 잊어버리고 현재의 모델로 ‘현재까지 가장 좋은 모델 값’을 덮어씌운다.

좋은 모델이 추정된 경우 (i.e. outlier가 없이 추정된 경우), residual이 작게 나오고 SSE가 작게 나온다.

RANSAC이 끝나는 조건은 여러가지 방법이 있지만, 많이 쓰이는 방법은 아래와 같다.

  1. 정해놓은 iteration 수가 전부 돌았을 때 (e.g. 100번의 iteration을 돌라고 설계했고, 100번을 다 돌았을 때)
  2. 정해놓은 residual threshold보다 더 낮은 에러가 나왔을 때 (e.g. pixel RMSE가 2.0 미만이면 iteration을 중단)

하지만 1번과 2번 방법 둘 다 장단점이 있다.

1번의 경우, 운이 좋아서 RANSAC 초반에 좋은 모델을 찾게 되어도 무조건 정해놓은 iteration이 다 돌 때 까지 계산을 돌려야한다. 즉, early-termination이 불가능하고, 쓸데없는 계산을 하게 된다.

2번의 경우, 운이 좋아서 RANSAC 초반에 좋은 모델을 찾게 되면 RANSAC을 빨리 종료할 수 있다. 하지만 outlier들끼리의 패턴 (i.e. local minima)에 빠지게 되면 헤어나올 방법이 없다.

어쨌든 끝나면서

  • Model과 Model을 이루는 inlier 데이터 (i.e. Outlier가 제거된 모델 + 데이터)
  • Sum of Squared Error (SSE) 값

을 가지고 끝난다.

 


RANSAC의 코드

RANSAC의 알고리즘은 사실 이렇게 간단하지만, 의외로 RANSAC을 직접 짜서 쓰시는 분들을 주변에서 보기 힘들다. 이는 사실 코드를 구현하시는 분들 중에서 스크래치로 시작하시는 분들보다, 이미 잘 만들어진 라이브러리에서 가져와서 변형하시는 분들이 많은데, 유명한 OpenCV나 PCL 등에는 딱 RANSAC() 이라는 함수가 없다. OpenCV에서 구현된 방식은 Fundamental Matrix 를 구하는 함수 내부에 RANSAC이 들어가있거나, solvePnPRANSAC() 같이 특정 목적을 가진 함수 내부에 구현이 되어있다.

그도 그럴만한게 RANSAC의 경우, model fitting 프레임워크로 봐야지, 실질적인 내부 계산은 우리가 어떤 계산을 하고 싶냐에 따라 달라지기 때문에, template class로 놔두는게 맞다고 본다. 하지만 OpenCV와 같은 오픈소스 라이브러리에서 template class로 구현하면 C나 Python, Java, javascript 등 해당 라이브러리가 지원하는 다른 언어에서의 함수 이름의 통일성이라던지가 어려워진다. 그리고 무엇보다 접근성이 많이 떨어지기 때문에, 지금까지 RANSAC을 자주 사용하는 함수 내부에 심어놓은 것으로 보인다.

하지만 이번 GSoC 2020을 통해서 RANSAC이 따로 떨어져 나올 수 있는 가능성이 보였다.

아래는 내가 RANSAC template class를 구현한다면 어떻게 만들지 간단하게 만들어본 디자인이다.

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
// A rough design of RANSAC template class

template <typename Precision, typename MODEL>
class RANSAC
{
public:
// constructor
RANSAC(MODEL& model)
{
mMinimalSet = MODEL.minimalSet();
mMaxIteration = 100; // user-defined
mModel = model;
mResidualThreshold = static_cast<Precision>(3);
mModelThreshold = static_cast<Precision>(2);
}

// run
MODEL run()
{
std::random_device randDev;
std::mt19937 mtGenerator(randDev());

int32_t iterNum = 0;

while true
{
// Terminate when max iteration is reached
if (iterNum == mMaxIteration)
return mModel;

if (iterate(mtGenerator));
{
if (mModel.residual < mModelTreshold)
return mModel;
}

iterNum++;
}
}

private:
bool iterate(std::mt19937 randNum)
{
// Calculate residual based on random selection of data points.
// Also somehow saves model params in model.
std::vector<Precision> residuals = mModel.run(randNum);

Precision sumResidual = 0;

// Residual thresholding - ignoring low residuals
for (const auto& res : residuals )
{
res = res < mDataResidualThreshold? 0 : res;
sumResidual += res;
}

// Update best model if residual is low
if(sumResidual < mBestModelRes)
{
mModel.saveResidualSum(sumResidual);
return true;
}

return false;
}

MODEL mModel;
int32_t mMinimalSet = 0;
int32_t mMaxIteration = 0;
Precision mBestModelRes = 0;
Precision mDataResidualThreshold = 0;
Precision mModelThreshold = 0;
}