Deterministic Hypothesis Generation for

Robust Fitting of Multiple Structures

Kwang Hee Lee, Chanki Yu, and Sang Wook Lee, Member, IEEE

Abstract— We present a novel algorithm for generating robust and consistent hypotheses for multiple-structure model fitting.
Most of the existing methods utilize random sampling which produce varying results especially when outlier ratio is high. For a
structure where a model is fitted, the inliers of other structures are regarded as outliers when multiple structures are present.
Global optimization has recently been investigated to provide stable and unique solutions, but the computational cost of the
algorithms is prohibitively high for most image data with reasonable sizes. The algorithm presented in this paper uses a maximum
feasible subsystem (MaxFS) algorithm to generate consistent initial hypotheses only from partial datasets in spatially overlapping
local image regions. Our assumption is that each genuine structure will exist as a dominant structure in at least one of the local
regions. To refine initial hypotheses estimated from partial datasets and to remove residual tolerance dependency of the MaxFS
algorithm, iterative re-weighted L1 (IRL1) minimization is performed for all the image data. Initial weights of IRL1 framework are
determined from the initial hypotheses generated in local regions. Our approach is significantly more efficient than those that use
only global optimization for all the image data. Experimental results demonstrate that the presented method can generate more
reliable and consistent hypotheses than random-sampling methods for estimating single and multiple structures from data with a
large amount of outliers. We clearly expose the influence of algorithm parameter settings on the results in our experiments.

Index Terms—Fitting of multiple structures, hypothesis generation, maximum feasible subsystem (MaxFS), iterative re-
weighted L1(IRL1) minimization

—————————— u ——————————

1 INTRODUCTION

I N many computer vision problems observations or

measurements are frequently contaminated with outliers
noise [18]. Furthermore, the existence of multiple struc-
tures makes the problem more difficult since the inliers be-
and noise, and thus robust estimation is needed for model longing to other structures are regarded as outliers
fitting. The “hypothesize-and-verify” framework is the (pseudo-outliers). To the best of our knowledge, there is no
core of many robust geometric fitting methods. The Ran- stopping criterion to guarantee consistent and reliable re-
dom Sample Consensus (RANSAC) algorithm [24] is a sults for multi-structure data.
widely used robust estimation technique, and most of the When the number of iterations is insufficient, the ran-
state-of-the-art methods are based on random sampling. dom sampling-based techniques provide varying results
They are comprised of two steps: (1) random hypotheses for the same data and parameter settings. Despite their ro-
generation and (2) verification. These steps are performed bustness, the random sampling-based methods provide no
iteratively. Many hypotheses of geometric model are ran- guarantee of consistency in their solutions due to the ran-
domly generated from minimal subsets of the data. The domized nature [2]. Since many extensions of random
goal of random sampling is to generate many putative hy- sampling-based frameworks also follow the same heuristic
potheses for given geometric model. In the verification of random sampling, none of them can guarantee deter-
step, the hypotheses are evaluated according to robustness ministic solutions [2].
criterion to find the best model. The various approaches to improve the efficiency of the
Random sampling-based methods have some draw- random hypothesis generation for the estimation of single
backs. In the majority of random sampling-based methods, structure have been investigated [18, 29, 21, 22, 19, 20, 27].
the number of iterations needed to guarantee a desired They have been developed to increase the frequency of hit-
confidence can be determined by a priori knowledge such ting all-inlier samples. Unfortunately, these methods are
as inlier ratio and inlier scale. For the single-structure data limited to single structure problem. In recent years, to deal
with unknown inlier ratio, it is crucial to determine the ad- with multiple structure data, some guided sampling meth-
equate number of iteration. The more heavily contami- ods [15, 16, 23] have been proposed. To perform guided
nated the data are, the lower the probability of hitting an sampling, a series of tentative hypotheses are generated
all-inlier subset is. However, since the true inlier ratio is a from minimal subsets of the data in advance. Then, guided
priori unknown in many practical situations, it is necessary sampling based on preference analysis is performed. The
to be determined by users. The standard stopping criterion quality of the initial hypotheses depends on the proportion
in RANSAC is based on an assumption that a model gen- of pseudo-outliers and gross outliers. Furthermore, since
erated from an uncontaminated minimal subset is con- [15, 16, 23] have no clear stopping criterion, it is difficult to
sistent with all-inliers. However, in practice, this assump- determine optimal number of iterations needed. Although
tion may be discrepant either increasing runtime or esti- the methods mentioned above improve the probability of
mating incorrect solutions since inliers are perturbed by hitting an all-inlier subsets in some ways, consistent per-
formance cannot be guaranteed for the data with unknown
2

inlier ratio or the number of structures. Many multiple- fitting method for multi-structure data has been pro-
structure model fitting methods also start with random hy- posed. To reduce high computation cost, the MaxFS al-
pothesis generation [9, 10, 11, 13, 12, 14, 17, 8]. Some of gorithm is performed for small subset like the method
these algorithms classify dataset based on randomly gen- presented in this paper. However, while initial subset is
erated hypotheses and find model parameters [9, 10, 11]. obtained by sorting keypoint matching scores in [36],
Due to the same nature of random sampling, however, var- the presented method does not require application-spe-
ying results may be produced from the same dataset. cific knowledge. In [36], there exist dependencies be-
Global optimization has recently been actively investi- tween hypotheses generated, since a sequential “fitting-
gated for model fitting problems in computer vision [6, 2, and-removing” procedure is used. Therefore, it is im-
3, 25, 26]. Li has developed a global optimization method possible to generate hypotheses in parallel. On the other
for the algebraic DLT (Direct Linear Transformation) prob- hand, the presented method can be immediately paral-
lem that has fixed bounded variables [2]. He suggested an lelized since there is no dependency between all of the
exact bilinear-MIP (Mixed Integer Programming) formula- hypotheses generated.
tion and obtained globally optimal solutions using an LP The rest of the paper is organized as follows: Section
(Linear Programming)-based BnB (Branch-and-Bound) al- 2 introduces two deterministic methods for geometric
gorithm. In [3], Yu et al. directly solved Big-M based MILP fitting. Section 3 describes our algorithm based on
problem. While these methods guarantee globally optimal MaxFS and IRL1 frameworks. Section 4 shows the ex-
solution, high computational cost is required in general. perimental results on synthetic and real data, and we
For a large dataset, the global optimization methods re- conclude in Section 5.
quire a great deal more running time than RANSAC. Fur-
thermore, the presence of image features from multiple
2 DETERMINISTIC METHODS FOR ROBUST
structures makes their computation cost even higher.
In this paper, we present a novel approach to reliable GEOMETRIC FITTING
and consistent hypothesis generation for multiple- In this section, we briefly describe two main optimiza-
structure model fitting. Unlike previous random sam- tion techniques that we employ in our method.
pling methods, our method generates hypotheses using
deterministic optimization techniques, and thus pro- 2.1 Maximum Feasible Subsystem (MaxFS)
duces consistent results given a set of images. A whole Problem
image dataset is split into spatially overlapping circular The aim of a MaxFS algorithm is to find the largest car-
regions of subsets, and a maximum feasible subsystem dinality set with constraints that are feasible [4, 2, 3]. In
(MaxFS) problem is solved to generate consistent initial other words, it finds a feasible subsystem containing the
hypotheses in each region. Because of the reduction of largest number of inequalities for an infeasible linear sys-
data size using local regions, the MaxFS algorithm can tem Ax ≤ b with real matrix A∈Â ´ , real vector b∈Â and n k k

generate the initial hypothesis from each local region variable x∈Â . The objective of the MaxFS and RANSAC
n

with reasonable efficiency. Since the MaxFS algorithm are the same. However, unlike the RANSAC, the MaxFS
yields a globally optimal solution only for the image guarantees a global solution. The MaxFS problem admits
subset in a spatially local region and the result is influ- the following mixed integer linear programming (MILP)
enced by a residual tolerance value, an iterative re- formulation by introducing a binary variable y for each of i

weighted L1 (IRL1) minimization is carried out using all the inequalities:

the data in the image to compensate for fitting errors
from the subset and to get rid of the residual tolerance k
dependency. Our method is developed to find multiple min
x,y
åy i
structures under the assumption that a good hypothesis i =1

for each genuine structure will be found in at least one n

of the spatially local regions. subject to åa

j =1
ij x j £ bi + M i yi , "i (1)
It may be noted that the use of spatial restriction for
hypothesis generation is not unprecedented. There x Î Â n , yi Î {0,1}, i = 1,...,k
have been approaches to using spatial coherence in lo-
cal regions for estimating single and multiple struc- where Mi is a large positive number that converts an infea-
tures/motions [27, 28, 12]. They are based on random sible inequality into a feasible one when yi =1. The case
sampling and exploit spatial coherence mainly to in- where yi =0 indicates that the inequality is feasible. Note
crease the chance of finding all-inlier samples. Those that if yi =1, then the ith constraint is automatically deac-
methods therefore have the limitations of random sam-
tivated. This MILP formulation is known as the Big-M
pling when the inlier ratio is unavailable. Our algo-
rithm uses spatially overlapping local regions to gener- method [1].
ate the initial hypothesis efficiently and to exhaustively Generally, MILP problems are solved using the LP-
search for all genuine structures. It provides con-sistent based BnB (Linear Programming-based Branch and Bound)
hypotheses regardless of inlier ratio, noise effect and the or LP-based BnC (LP-based Branch and Cut) methods. LP-
number of structures. based BnB and BnC guarantee global optimality of their
Recently, in our previous work [36], a deterministic solutions [3, 30]. The MILP problem is expensive in terms
 3

(a) (b) (c)

Fig. 1. Overview of our approach. (a) Overlapping circular regions for hypothesis generation. (b) Detection of dominant structures in local
regions. (c) IRL1 minimization for refining initial hypothesis.

of computational cost, which is known to be NP-hard.

Thus, only relatively small problems can be solved practi- Algorithm 1. IRL1 algorithm [31]
cally. While the exact Big-M MILP formulation is useful for Input: y = Ax
small models, it is not effective on larger models for rea-
sons of computational inefficiency [1]. To solve the geo- Output: xt
metric fitting problem including heavily contaminated
multi-structure data, MaxFS demands a much higher com- 1: Initialize: Set the weights wi0=1 for i = 1…d
putational cost than RANSAC style methods.

2.2 Iterative Re-weighted L1 (IRL1) Minimization 2: Repeat

3: t = t+1
The most common deterministic method for geometric
fitting is the least-squares estimation. The best fit mini-
4: xt = arg min | Wt x |1 , s.t.y = Ax
x
mizes the sum of the squared residuals (L2-norm). This
5: 1
method is optimal for Gaussian noise, but the algorithm wit +1 =
mostly fails in the presence of outliers. Using the sum of | xit | +a
the absolute residuals (L1-norm) relatively brings about
better results than using L2-norm in the presence of outliers, 6: Until convergence or a fixed number of times
since L1 minimization puts less weight on large residuals
than L2 minimization. Nevertheless, L1 minimization still
cannot guarantee robustness in the case of severely con- 3 HYPOTHESIS GENERATION USING
taminated data or multiple structured data with outliers DETERMINISTIC OPTIMIZATION TECHNIQUES
because the influence function has no cut off point. Alt- We present a deterministic algorithm to generate a reli-
hough these methods always have a global solution, relia- able and consistent hypothesis for robust fitting of multiple
ble output cannot be guaranteed in the presence of severe structures. The whole data space is split into densely over-
outliers. lapping circular regions (see Fig. 1(a)). In each region, an
Iterative re-weighted L1 (IRL1) minimization has been initial hypothesis is generated from the maximum inliers.
presented by Candѐs, Wakin and Boyd [31]. The IRL1 algo- However, the initial hypothesis may be slightly biased
rithm solves a sequence of weighted L1-norm problems even if it is generated from pure inliers, since the estimated
where the weights are determined according to the esti- inliers are from the local region. To estimate the best fit for
mated coefficient magnitude. The IRL1 minimization algo- each genuine structure from the initial hypothesis, IRL1
rithm is summarized in Algorithm 1. W is a diagonal minimization is performed for all of the data in the image
weighting matrix from the tth iteration with ith diagonal el- (see Fig. 1(c)). The algorithm of our method is summarized
ement wit and a is a stability parameter which affects the in Algorithm 2.
stability of the algorithm.
In [31], experimental results show that it often outper-
forms standard L1 minimization. Although each iteration
of the algorithm solves a convex optimization problem, the
overall algorithm does not. Therefore, one cannot expect
this algorithm to always find a global minimum. Conse-
quentially, it is important to determine a good starting
point for the algorithm. In [31], initializing with the solu-
tion to standard L1 minimization has been introduced.
However, it is unsuitable for problems with heavily con-
taminated multi-structure data.
4

Algorithm 2. MaxFS-IRL1 algorithm

Input: input data D, the number of data (in each subset) Rj s

xjc
k, the fractional constant α, Mi (Big-M), the resid-
ual tolerance ε, the variance of weight function σ
s

Output: hypothesis sets {Θ*j}j=1,…l

(a) (b)
1: Split D into l subset Dj

2: For j = 1 to l

3: Estimate initial parameter Q̂ MaxFS using MaxFS rmaxdensity

from Eq. (6) (from Dj)
4: Initialize weights from Q̂ MaxFS from Eq. (7)
5: Refine parameter Θ*j using IRL1 minimization
(from D) : (c) (d)

6: Repeat
xjc s
7: Solve the weighted L1 minimization prob- s
lem from Eq. (8)
8: Update the weights from Eq. (9)

9: Until convergence or a fixed number of times

(e) (f)
10: End
Fig. 2. (a) Center xjc and radius Rj of circular region j and step size s.
(b) Input image and data points. (c) Data density map estimated from
3.1 Determination of Spatially Overlapping Circular 2D KDE method. (d) Estimation of rmaxdensity. (e) Estimation of s and xjc.
Regions (f) Spatially overlapping circular regions.
Our algorithm splits a whole image region into many
overlapping small circular regions. The input data (D) is
over the regions, the circular windows include approxi-
the union of l subsets
mately the same number of data points. Thus, the compu-
l
tational cost for initial hypothesis generation is about the
D = {xi }iN=1 = ! D j , (2) same.
j =1 Our algorithm places the centers of circular windows
(xjc s) at a regular interval s in the horizontal and the verti-
where Dj is the set of data in the circular region j. Neigh- cal directions of image as illustrated in Fig. 2 (a), but the
boring subsets, Dm and Dn, share common data in the over- radii of the circles (Rj s) vary to keep the number of data
lapping region (Dm∩Dn≠Ø). points approximately the same. Once k is set, the smallest
It is assumed that every structure in the image appears circle that contains k data points should appear where the
as a dominant structure in at least one of the local regions, data density is highest and we take its radius r as a ref-
maxdensity

and thus it suffices for our algorithm to find only one dom- erence for determining the step size s and Rj s for individual
inant structure in a circular region. The remaining struc- circles. We compute the density of data using the 2D Kernel
tures are missed due to the dominant structure’s initial de- Density Estimate (KDE) described in [7], and Figure 2 (c)
tection, yet they can be found in the other local regions. shows the data density map for the image data points
Figure 1 (b) shows an example where an undetected struc- shown in Fig. 2(b). In Fig. 2 (d), the red circle represents the
ture in a window (Fig. 1 (b), middle) becomes a dominant highest-density region with k-nearest neighbors and the
structure in a neighboring window (Fig. 1 (b), bottom). dotted red line does r . maxdensity

The circles’ sizes and positions should be determined The step size s is determined as follows:
depending on the number of data points to be included in
the circular regions. The number of data points k that a cir- s = rmaxdensity α, (3)
cle covers should be larger than the minimum number re-
quired for fitting a desired model for hypothesis genera- where α is a fractional constant. The subset Dj around the
tion. If k becomes larger, the result is more reliable but the circle center xjc is defined as:
computational cost is higher. The smaller the interval be-
tween the circles is, the slimmer the chance of missing
small structures becomes. To maintain even performance
 5

D j = {xi Î D | |xi - xcj | £ R j }, If yi=1, the i data is an outlier and the corresponding con-
th

(4) straint is deactivated automatically. We use a linear con-

R j = min(rj , 2rmaxdensity ),
straint cTΘ=1, rather than the commonly used ||Θ||=1, where
c is a problem dependent vector determined by the user [5].
where rj is the minimum raidus for the circle that contains
Note that our MaxFS algorithm solves Eq. (6) for every
k-nearest neighbors. Figure 2 (e) shows the center points of
subset Dj. Then, a series of hypotheses { Θ̂1MaxFS ,…,
the circular regions (blue crosses), and an example of spa-
Θ̂l MaxFS } are generated from the maximum inlier-set
tially overlapping circular regions (green circles) is shown
where l is the number of hypotheses.
in Fig. 2 (f).
In the hypothesis generation stage, although the MaxFS
Since r is the shortest possible radius, rj is always
maxdensity
algorithm obtains a global solution for the local subset, the
longer than r and thus substantial overlap between the
maxdensity
parameter vectors estimated in spatially restricted regions
adjacent circles is guaranteed. The constant α controls the
can be biased againt the true structures. To refine initial hy-
extra degree of overlap. The maximum Rj is set to 2r maxdensity
potheses estimated from local datasets, IRL1 minimization
since rj becomes unmeaningfully long where data density
is performed for all data in the image.
is very low. The number of datapoints k and the constant α
are the most important parameters in our algorithm and
we experimentally investigate the influence of their set- 3.3 Hypothesis Refinement Using the IRL1
tings on the results in Secion 4. Minimization
In each local region, given an initial hypothesis gener-
ated from the MaxFS algorithm, initial weights of all the
3.2 Initial Hypothesis Generation using the MaxFS
data can be determined as
Algorithm
We determine the maximum feasible inliers in each local ˆ MaxF S
- | riQ |
region. An initial hypothesis is generated from the maxi- wi(0) = exp( ), i = 1,..., N ,
mum feasible inliers in each local region. A simple way of s
deterministically solving this problem would be to per- (7)
form an exhaustive search. This is intractable due to the MaxFS
combinatorial explosion of subsets. On the other hand, ran- where riQ̂ is the residual of data xi to the initial hypoth-
dom-sampling methods cannot guarantee consistent re- esis and σ controls the width of the weight function. Unlike
sults due to their randomized nature and stopping crite- the original IRL1 [31] minimization algorithm, our IRL1
rion depends on the prior knowledge of information such minimization algorithm uses the MaxFS algorithm to gen-
as inlier ratio. erate the initial weights. Since the initial hypothesis gener-
Unlike random-sampling methods, the MaxFS algo- ated from the MaxFS algorithm is robust to outliers, it may
rithm can guarantee maximum feasible inliers, even provide much better initial weights than the original IRL1
though the inlier ratio is unknown. Moreover, by splitting [31] minimization that uses the standard L1 minimization
the problem into many small parts, the MaxFS algorithm as the method to generate initial weights.
can be performed quickly and efficiently. We use the alge- After initial weights are set, the IRL1 minimization iter-
braic Direct Linear Transformation (DLT) to estimate hy- atively performs the two-step algorithm shown below. The
pothesis parameters [5]. We can then formulate the DLT- first step solves the following weighted version of L1 min-
based geometric fitting problem as a MaxFS problem. Each imization:
subset Dj is partitioned into the inlier-set DjI and the outlier-
N
set DjO with DjI⊆Dj, DjO⊆Dj, DjI∪DjO=Dj and DjI∩DjO=Ø.
A maximum residual tolerance, ε>0, provides a bound
ˆ , yˆ } = argmin
{Θ t
Θ, y
∑w y t
i i
i =1
for the algebraic residual di=|aiTΘ| at point i, where aiT is T (8)
each row vector of A in the homogeneous equation AΘ=0: subject to | ai Θ | ≤ yi , ∀ i
cTΘ = 1,
di = |aiTΘ| ≤ ε, ε>0. (5)
Θ Î Â n , yi = [0, ¥], i = 1,...,N ,
The MaxFS formulation of Eq. (5) is as follows:
where t is the index of current iteration.
k
ˆ MaxFS , yˆ } = argmin In the next step, the weights of all the data are updated
{Θ
Θ, y ∑y i as follows:
i =1
(6)
subject to | ai T Θ | ≤ e + M i yi , ∀ i - | r Θt |
ˆ

c T Θ = 1, w ( t +1)
i = exp( i ), i = 1,..., N ,
σ
Θ Î Ân , yi ∈{0,1}, i = 1,...,k . (9)

where Mi is a large positive number (called Big-M value). where riΘ̂ t is the residual of data xi to the current hypothe-
The case where yi =0 indicates that the i data is an inlier. th
sis and σ controls the width of the weight function. In each
iteration, Eq. (8) and Eq. (9) are alternately performed (see
6

Fig. 3. Examples for line fitting. (Top row) Input data with different outlier ratio (50%, 66%, 75%, 20%, 33.3%, 42.8%). (Middle row) The gener-
ated hypotheses using the proposed method. (Bottom row) Fitting results.

Algorithm 2, 6~9). In other words, the current parameter

vector Θ̂ t is estimated by solving Eq. (8) and then weights4.1 Experiments with Synthetic Datasets: Single
of all the data are updated based on the current parameter and Multiple Line Fitting
vector by Eq. (9). The algorithm is repeated until conver- The first set of results we show is produced from single
gence or for a fixed number of times. and multiple 2D line fitting. We performed the DLT-based
Note that an advantage of IRL1 minimization is that the MaxFS and IRL1 algorithms for each data subset and the
results are not influenced by residual tolerance, unlike results are shown in Figs. 3 and 4. The residual tolerance ε
with RANSAC or MaxFS. Therefore, our IRL1 minimiza- was set to 1 and the Big-M value of Eq. (6) was fixed to
tion not only refines a hypothesis biased by local fitting, 10000. The variance of weight function σ was fixed to 10.
but also eliminates the residual tolerance dependency. The number of data points in the circular region k was fixed
to 40. In the single line test, each line includes 100 inliers
3.4 Complementary Role of MaxFS and IRL1 with Gaussian noise and various gross outliers. Noise level
The hypotheses estimated from only MaxFS may deviate was fixed at 0.03, and the number of gross outliers varies
from the true structures since they are estimated in local from 100 to 300. In multiple line fitting tests, the same ratio
regions. We use the IRL1 procedure to globally refine these of gross outliers and Gaussian noise was added. Final fit-
hypotheses. The initial weights of the conventional IRL1 ting from generated candidate hypotheses was performed
procedure are usually generated from the standard L1 min- using Li’s method [8].
imization which often fails when there are a large amount Figure 3 shows input synthetic data (Top row), corre-
of outliers. We employ the MaxFS algorithm to generate sponding hypotheses results (Middle row) and the final fit-
much more reliable initial weights than the standard L1 ting results (Bottom row). The results show our algorithm
minimization. It is our intention to carefully combine the generates good hypotheses from true structures. Even
two algorithms to complement each other’s limitations. when the ratio of gross outliers increases, the proportion of
The contribution of each component can be seen in Fig. 11 good hypotheses is fairly consistent.
and Fig. 12. Figure 4 shows (First column) initial hypotheses from the
MaxFS algorithm, (Second column) hypotheses after the
first reweighted iteration, (Third column) hypotheses after
4 EXPERIMENTAL RESULTS the third reweighted iteration and (Last column) fitting re-
We have implemented our algorithm in MATLAB using sults. We used the IRL1 minimization from all of the data
the LP/MILP solver GUROBI [30] which provides func- to reduce the errors from subset data. In the first row, noise
tions for the LP/MILP and a desktop with Intel i5-2500 level was fixed at 0.01 and the number of gross outliers was
3.30GHz (4 cores) and 3GB RAM was used for experiments. fixed at 200. In the second row, noise level is fixed at 0.03
We used 4 cores only for solving mixed integer linear pro- and the number of gross outliers was the same. With the
gramming in each MaxFS problem. We measured the ac- residuals of the whole dataset being used for estimation,
tual elapsed computation time and tested our proposed hypotheses are refined for genuine structures.
method on several synthetic and real datasets.
 7

Fig. 4. Hypothesis refinement through reweighted L1 iterations. (First column) Initial hypotheses from Max-FS framework. (Second column)
Hypotheses after the first reweighted iteration. (Third column) Hypotheses after the third reweighted iteration. (Last column) Fitting results.

for the homography estimation and 20 and above for the

4.2 Experiments with Real Datasets affine fundamental matrix estimation, and the computa-
We tested four methods including our MaxFS-IRL1 al- tion time gradually increases with k. For the MaxFS-IRL1
gorithm on several real datasets. For performance evalua- algorithm, we set k to 30 for the homography estimation
tion and comparison, we measured the actual elapsed com- and 20 for the affine fundamental matrix estimation to at-
putation time. Images and keypoint correspondences were tain both accuracy and computational efficiency. Secondly,
acquired from the Oxford VGG dataset [32] and the Ade- we investigated the effects of α in the range of [0.25, 1.5]
laide RMF dataset [34, 35]. The image pairs used in the ex- with predetermined k. Figs. 7(c) and 7(d) show the re-pro-
periment are shown in Figs. 5 and 6. Yellow crosses indi- jection errors and computation time for the homography
cate the gross outliers randomly generated and other col- estimation with the MaxFS-IRL1 method, respectively.
ored squares indicate the inliers of each structure. Each da- Figs. 8(c) and 8(d) show the re-projection errors and the
taset included the various rates of outliers. We used man- computation time for affine fundamental matrix estima-
ually labeled keypoint correspondences which were ob- tion with the MaxFS-IRL1 method, respectively. Based on
tained by SIFT matching. Given true inlier-sets and the these results, we set α to 1 for the homography estimation
number of structures, for each structure, a single hypothe- and 0.5 for the affine fundamental matrix estimation by
sis having minimum re-projection errors between ground considering both accuracy and computational efficiency.
truth inliers and estimated hypothesis parameters was se- We empirically examined the effects of parameters ε (re-
lected. The quality of hypotheses estimated was evaluated sidual tolerance) and σ (variance of weight function).
by averaging re-projection errors for all structures. Firstly, we investigated the effect of ε with fixed σ=1 in the
range of [0.1, 0.3, … 2.9] for the image pairs for homogra-
4.2.1 Analysis: MaxFS-IRL1 Framework phy estimation (Fig. 5) and in the range of [0.0005, 0.001,
We performed the MaxFS-IRL1 method to estimate pla- …, 0.005] for the image pairs for affine fundamental matrix
nar homography and affine fundamental matrix. For our estimation (Fig. 6). Figs. 9(a) and 9(b) show the re-projec-
MaxFS-IRL1 algorithm, the Big-M value in Eq. (6) was set tion errors from the MaxFS algorithm and MaxFS-IRL1
to 10000 for both applications. For estimation we used Di- method for one of several structures in the CollegeIII data
rect Linear Transformation (DLT) and the residual was with 60% of the outliers included and the Carchipscube
taken as the Sampson distance [5]. data with 40% of the outliers included. The results clearly
We experimentally examined the effects of parameter k show that our algorithm is stable within a wide range of
(subset size) and α (fractional constant) on different da- the inlier tolerance ε. We set ε=0.5 for homography estima-
tasets with 70% of outliers. Firstly, we investigated the ef- tion and ε=0.001 for affine fundamental matrix estimation
fects of k in the range of [10, 50] with fixed s=100 for the based on our test results. Secondly, while σ is varied in the
homography estimation and s=50 for the affine fundamen- range of [1, 2, … 5] with fixed ε, we recorded the re-projec-
tal matrix estimation. Figs. 7(a) and 7(b) show the re-pro- tion error in each iteration. Figs. 9(c) and 9(d) shows the
jection errors and computation time for the homography results for the Neem data with 50% of the outliers included
estimation with the MaxFS-IRL1 method, respectively. and the Cubebreadtoychips data with 30% of the outlier
Only the results from three datasets are shown in the plots. included. In most experiments, when σ<5, the results were
Figs. 8(a) and 8(b) show the re-projection errors and the stable and accurate. Considering the convergence speed
computation time for the affine fundamental matrix esti- and accuracy, we set σ=3 by default for both applications.
mation, respectively. It can be seen in Figs. 7 and 8 that
high accuracy is achieved for the k from about 30 and above
8

(a) Johnsona (b) Neem (c) CollegeIII

(d) CollegeII (e) Unionhouse (f) Elderhalla

(g) Sene

Fig. 5. Test image-sets used for homography estimation experiment. Yellow crosses indicate the gross outliers randomly generated and other
colored squares indicate the inliers of each structure.

(a) Biscuitbook (b) Carchipscube (c) Breadcartoychips

(d) Gamebiscuit (e) Cubetoy (f) Dinobooks

(g) Cube (h) Biscuitb

okbox
Fig. 6. Test image-sets used for affine fundamental matrix estimation experiment. Yellow crosses indicate the gross outliers randomly gener-
ated and other colored squares indicate the inliers of each structure.

elapsed computation time for our method was not limited

4.2.2 Comparison with Other Methods but measured. To ensure a fair comparison, three compet-
Our algorithm was compared with uniform random ing methods were performed for similar periods of time
sampling (Random Sampling) [24, 33], NAPSAC [27] and (elapsed computation time) along with our algorithm. In
the state-of-the-art algorithm Multi-GS [23, 34]. We imple- this experiment, we selected all of the parameters based on
mented NAPSAC in MATLAB. For performance evalua- section 4.2.1 for our MaxFS-IRL1 framework.
tion, we measured elapsed computation time and the num- Homography Estimation: We tested the performance of
ber of generated hypotheses (L) and computed the re-pro- our proposed method for estimating planar homographies
jection errors (maximum error and standard deviation). on real image data. For the homography estimation, the
The results for the four algorithms are summarized in Ta- Big-M value of Eq. (6) was fixed to 10000, the residual tol-
bles 1 and 2 with the best results shown in bold. erance ε was fixed to 0.5, the variance of weight function σ
For Random Sampling, NAPSAC and Multi-GS, 50 ran- was fixed to 3, the number of iterations in IRL1 was fixed
dom runs were performed. Unlike these competitors, since
our method runs untill completion of the algorithm, the
 9

to 5, fractional constant α was fixed to 1 and the number of hypothesis is reliable but may not yield the minimum er-
data points in the circular region k was fixed to 30. rors for true inliers. Note that random sampling-based
Table 1 summarizes the performance of four methods for methods produced large variations in their results.
estimating planar homography for seven datasets (see Fig.
5. (a-g)) with various outlier ratios. The results demon- 4.2.3 Performance under Increasing Outlier Rates
strate that our method yields more reliable and consistent We compared the performance of four algorithms un-
results with reasonable computational efficiency. der different outlier rates. We increased the outlier ratio
Affine Fundamental Matrix Estimation: We also tested from 10% to 80% for the CollegeIII data and the Biscuit-
the performance of our proposed method for estimating an bookbox data. Figs. 10(a) and 10(c) show the re-projection
affine fundamental matrix on real image data. For the af- errors produced by the four methods on the two test da-
fine fundamental matrix estimation, the Big-M value was tasets as the outlier ratio increases and Figs. 10(b) and 10(d)
fixed to 10000, the residual tolerance ε was fixed to 0.001, show the corresponding standard deviations for the re-
the variance of weight function σ was fixed to 3, the num- projection errors. Our algorithm outperforms the other al-
ber of iterations in IRL1 also was fixed to 5, fractional con- gorithms when outlier ratio is high over similar periods of
stant α was fixed to 0.5 and the number of data points in time (elapsed computation time). Since the probability of
the circular region k was fixed to 20. producing an all-inlier subset decreases with random sam-

(a) (b) (c) (d)

Fig. 8. The effect of parameter k (subset size) and α (fractional constant) on different datasets with 70% of outliers for the affine fundamental
matrix estimation with MaxFS-IRL1 method. (a) The re-projection errors obtained from different k. (b) The computation time measured from
different k. (c) The re-projection errors obtained from different α. (d) The computation time measured from different α .

(a) 1) (Ga (c) (d)

meb
Fig. 7. The effect of parameter k (subset size) and α (fractional constant) on different datasets with 70% of outliers for the homography
estimation with MaxFS-IRL1 method. (a) The re-projection errors obtained from different k. (b) The computation time measured from different
k. (c) The re-projection errors obtained from different α. (d) The computation time measured from different α.

Table 2 shows the performance of the algorithms for esti- pling-based approaches as the outlier ratio increases, the
mating the affine fundamental matrix for seven datasets maximum errors and standard deviations increase
(see Fig. 6. (a-g)) with varied outlier ratios. The results substantially. On the other hand, more robust results were
show that our algorithm generates high-quality hypothe- yielded with our MaxFS-IRL1 method since an initial hy-
ses with reasonable efficiency for the datasets with high pothesis generated from the MaxFS was rarely influenced
outlier ratios and thus finds all the true structures stably by the outlier ratio.
from all of the datasets. In most cases, our method results
in the smallest errors on all the test datasets except for the 4.2.4 Combination of MaxFS and IRL1
Biscuit book dataset with a relatively low outlier ratio. In Figs. 11 and 12, we compared the results using the
Since our method performs IRL1 minimization for all of conventional IRL1 method and using our MaxFS-IRL1
the data in order to generate a hypothesis while other method. We used the CollegeIII dataset. For the conven-
methods generate a hypothesis from minimal subsets, its tional IRL1 method, the initial weights were generated
10

from standard L1 minimization from local datasets to com- model fitting. For reliable hypothesis generation with rea-
pete on par with the MaxFS-IRL1 method. We can see the sonable computational efficiency, we employ a MaxFS
contribution of each stage from Figs. 11 and 12. (maxium feasible subsystem) algorithm, a global optimiza-
Figures 11(a) and 11(d) show initial fitting results from tion technique, only in spatially localized image regions
standard L1 minimization and the MaxFS algorithm for the and refine the hypotheses using an IRL1 (re-weighted L1)
CollegeIII data with no outliers included. Figures 11(b) and minimization method. To search out all of the structures
11(e) show the fitting results after the first reweighted iter- thoroughly, the local circular regions spatially overlap with
ation. Figs. 11(c) and 11(f) show the final fitting results the neighboring ones, and the number of data in each local
from the conventional IRL1 method and the MaxFS-IRL1 region is distributed evenly for computational efficiency.
method. The model parameters of major structure are estimated in
Figures 12(a-c) show fitting results from the conven- each local image region, and those of the remainder can be
tional IRL1 method in each iteration step and Figures 12(d- found in one of the neighboring regions. The IRL1 minimi-
f) show those from the MaxFS-IRL1 method for the Col- zation is performed over all the image data to refine initial
legeIII data with 30% of the outliers included. hypotheses estimated from subsets and to get rid of resid-
As appears by these results, when the outlier ratio is low, ual tolerance dependency of the MaxFS algorithm.
the standard L1 minimization can provide good results. Our experiments show that without prior knowledge of
When there are severe outliers, however, standard L1 min- the inlier ratio, inlier scale and the number of structures,

(a) (b) (c) (d)

Fig. 10. The performance of four algorithms under different outlier rates. Graphs (a) and (c) show the re-projection errors produced by the four
methods on the two test datasts as outlier ratio increases. (b) and (d) show the corresponding standard deviations for the re-projection errors.

(a) (b) (c) (d)

Fig. 9. (a) and (b) show the re-projection errors from the MaxFS algorithm and MaxFS-IRL1 method for one of several structures in the
CollegeIII data with 60% of the outliers included and the Carchipscube data with 40% of the outliers included. (c) and (d) show the re-projection
errors obtained from different σ for each iteration for the Neem data with 50% of the outliers included and the Cubebreadtoychips data with
30% of the outliers included.

imization frequently fails. To obtain good fitting results us- our method generates consistent hypotheses which are
ing the IRL1 method, it is important to determine good in- more reliable than the random sampling-based methods as
itial weights for the algorithm. Therefore, the MaxFS and outlier ratios increase.
the IRL1 algorithms are complementary and we show that
their combination is highly effective for the task of robust
multiple-structure fitting.

5 CONCLUSION
We present a new deterministic approach to reliable and
consistent hypothesis generation for multiple-structure
 11

(a) (b) (c)

(d) (e) (f)

Fig. 11. (a) and (d) show initial fitting results from the standard L1 minimization and the MaxFS algorithm for the CollegeIII data with no outlier
included. (b) and (e) show the fitting results after the first reweighted iteration. (c) and (f) show the final fitting results from the conventional
IRL1 method and the MaxFS-IRL1 method.

(a) (b) (c)

(d) (e) (f)

Fig. 12. (a) and (d) show initial fitting results from the standard L1 minimization and the MaxFS algorithm for the CollegeIII data with 30% of
the outlier included. (b) and (e) show the fitting results after the first reweighted iteration. (c) and (f) show the final fitting results from the
conventional IRL1 method and the MaxFS-IRL1 method.

[3] C. Yu, Y. Seo, and S. W. Lee. Photometric stereo from maximum

feasible lambertian reflec-tions. In ECCV, 2010.
[4] Y. Zheng, S. Sugimoto and M. Okutomi. Deterministically max-
REFERENCES
imizing feasible subsystem for robust model fitting with unit
[1] J. W. Chinneck. Feasibility and Infeasibility in Optimization: Al- norm constraint. In CVPR 2011: 1825-1832
gorithms and Computational Methods, 1st edition, Springer, [5] R. Hartley and A. Zisserman. Multiple view geometry in com-
Heidelberg 2007. puter vision 2nd ed. Cambridge University Press, 2004.
[2] H. Li. Consensus Set Maximization with guaranteed global opti- [6] O. Enqvist, F. Kahl. Two view geometry estimation with outliers.
mality for robust geometry estimation. In ICCV 2009. In BMVC 2009.
12

TABLE 1
PERFORMANCE OF VARIOUS METHODS ON HOMOGRAPHY ESTIMATION FOR THE SEVERAL REAL DATASETS.

Data Outlier Ratio[%] Method Random NAPSAC Multi-GS MaxFS-IRL1

Elapsed time [sec] 15 15 15 13.71
Johnsona 10 Max Error 1.0871 1.5245 1.7234 0.8512
Std 0.0563 0.1351 0.1973 0
L 12028 10926 1855 24
Elapsed time [sec] 10 10 10 9.18
Neem 20 Max Error 1.178 1.6387 1.0993 0.9294
Std 0.0508 0.1324 0.0352 0
L 8386 7529 1496 12
Elapsed time [sec] 35 35 35 33.15
CollegeIII 30 Max Error 1.9025 2.2848 1.9946 1.4536
Std 0.0643 0.1105 0.0983 0
L 25043 23068 2628 51
Elapsed time [sec] 50 50 50 49.19
CollegeII 40 Max Error 1.5704 1.8473 1.574 1.4851
Std 0.0136 0.0668 0.0135 0
L 35775 33819 2589 35
Elapsed time [sec] 10 10 10 7.61
Unionhouse 50 Max Error 0.6507 0.8214 0.65 0.5817
Std 0.0148 0.0515 0.0165 0
L 8559 7868 2239 9
Elapsed time [sec] 20 20 20 18.62
Elderhalla 60 Max Error 2.9963 2.1821 1.9426 1.6893
Std 0.2685 0.1278 0.0711 0
L 16922 16080 2721 20
Elapsed time [sec] 40 40 40 39.09
Sene 70 Max Error 1.0424 0.9275 0.6983 0.491
Std 0.0953 0.1023 0.0331 0
L 29208 28875 3061 35

TABLE 2
PERFORMANCE OF VARIOUS METHODS ON AFFINE FUNDAMENTAL MATRIX ESTIMATION FOR THE SEVERAL REAL DATASETS.

Data Outlier Ratio[%] Method Random NAPSAC Multi-GS MaxFS-IRL1

Elapsed time [sec] 20 20 20 21.7
Biscuitbook 20 Max Error 1.0147 1.0141 1.0529 1.0591
Std 0.0047 0.0056 0.0138 0
L 20444 21041 1839 108
Elapsed time [sec] 20 20 20 17.08
Carchipscube 30 Max Error 0.5542 0.5023 0.5356 0.4927
Std 0.0171 0.0036 0.0112 0
L 29832 27005 3295 74
Elapsed time [sec] 30 30 30 28.8
Breadcartoychips 40 Max Error 0.789 0.7454 0.7344 0.7114
Std 0.0205 0.0084 0.0101 0
L 31450 28929 3023 159
Elapsed time [sec] 35 35 35 35.53
Gamebiscuit 50 Max Error 0.7196 0.6752 0.6863 0.6744
Std 0.0085 0.0038 0.0066 0
L 31028 28927 2914 178
Elapsed time [sec] 30 30 30 28.58
Cubetoy 50 Max Error 0.6469 0.6317 0.6257 0.6076
Std 0.0085 0.0058 0.0048 0
L 24468 26301 2328 194
Elapsed time [sec] 40 40 40 40.39
Dinobooks 60 Max Error 2.8726 2.3352 2.2649 2.2078
Std 0.1693 0.0567 0.04 0
L 25307 23109 2599 156
Elapsed time [sec] 30 30 30 32.4
Cube 70 Max Error 0.5301 0.5486 0.5197 0.5053
Std 0.0079 0.0086 0.004 0
L 26585 27652 2558 159
 13

[7] Z. I. Botev, J. F. Grotowski and D. P. Kroese Kernel density esti- [35] http://cs.adelaid.edu.au/~hwong/doku.php?id=data
mation via diffusion. Annals of Statistics 2010. Volume 38, Num- [36] Kwang Hee Lee, Sang Wook Lee.: Deterministic Fitting of Mul-
ber 5, Pages 2916—2957 tiple Structures using Iterative MaxFS with Inlier Scale Estima-
[8] H. Li. Two-view motion segmentation from linear programming tion. In ICCV 2013.
relaxation. In CVPR, 2007.
[9] R. Toldo, A. Fusiello. Robust Multiple Structures Estimation
with JLinkage. In ECCV 2008: 537-547.
[10] T.-J. Chin, H. Wang, D. Suter. Robust Fitting of Multiple Struc-
tures: The Statistical Learning Approach. In ICCV 2009 : 413-420.
[11] W. Zhang, J. Kosecka. Nonparametric Estimation of Multiple
Structures with Outliers. LNCS, DVW 2006: 4358. 60-74.
[12] Y. Kanazawa, H. Kawakami. Detection of planar regions with
uncalibrated stereo using distributions of feature points. In:
BMVC 2004.
[13] M. Zuliani, C. Kenney and B. Manjunath. The multiransac algo-
rithm and its application to detect planar homographies. In ICIP
2005 : Volume 3.
[14] H. S. Wong, T.–J. Chin, J.Yu and D. Suter. Efficient Multi-struc-
ture Robust Fitting with Incremental Top-k Lists Comparison. In
ACCV 2010: 553-564.
[15] T.-J. Chin, J. Yu and D. Suter. Accelerated Hypothesis Genera-
tion for Multistructure Data via Preference Analysis. IEEE Trans.
Pattern Anal. Mach. Intell. 34(4): 2012 625-638
[16] H. S. Wong, T.-J. Chin, J. Yu and D. Suter. Dynamic and hierar-
chical multi-structure geometric model fitting. In ICCV 2011:
1044-1051
[17] A. Delong, A. Osokin, H. Isack, and Y. Boykov. Fast approximate
energy minimization with label costs. In CVPR, 2010.
[18] O. Chum, J. Matas and J. Kittler. Locally optimized RANSAC. In:
DAGM 2003.
[19] B. J. Tordo and D. W. Murray. Guided-MLESAC: Faster image
transform estimation by using matching priors. IEEE Trans. Pat-
tern Anal. Mach. Intell. (27) 2005:1523-1535
[20] O. Chum and J. Matas. Matching with PROSAC- progressive
sample consensus. In CVPR.2005.
[21] K. Ni, H. Jin and F. Dellaert. GroupSAC: Efficient consensus in
the presence of groupings. In ICCV. 2009.
[22] T. Sattler, B. Leibe and L. Kobbelt. SCRAMSAC: Improving
RANSAC's efficiency with a spatial consistency filter. In: ICCV
2009.
[23] T.-J. Chin, J. Yu, and D. Suter. Accelerated hypothesis generation
for multi-structure robust fitting. In ECCV 2010.
[24] M. A. Fischler and R. C. Bolles. RANSAC: A paradigm for model
fitting with applications to image analysis and automated car-
tography. Comm. of the ACM 24 1981: 381-395
[25] R. Hartley and F. Kahl. Global Optimization through Rotation
Space Search. IJCV, 82(1):64-79, 2009.
[26] F. Kahl. Multiple view geometry and the L-infinity norm. In
ICCV 2005: 1002–1009.
[27] D. R. Myatt, Philip H. S. Torr, Slawomir J. Nasuto, J. Mark Bishop,
R. Craddock. NAPSAC: High Noise, High Dimensional Robust
Estimation - it's in the Bag. In BMVC 2002.
[28] K. Schindler and D. Suter. Two-View Multibody Structure-and-
Motion with Outliers through Model Selection. IEEE Trans. Pat-
tern Anal. Mach. Intell. 28(6): 983-995, 2006.
[29] R. Raguram, J. M. Frahm and M. Pollefeys. Exploiting uncer-
tainty in random sample consensus. In ICCV 2009: 2074-2081.
[30] Gurobi optimization. http://www.gurobi.com/, 2010
[31] E. J. Cand`es, M. B. Wakin and S. P. Boyd. Enhancing sparsity by
reweighted L1 minimization, The Journal of Fourier Analysis
and Applications, vol. 14, no. 5, pp. 877–905, 2008.
[32] http://www.robots.ox.ac.uk/~vgg/data/
[33] http://www.csse.uwa.edu.au/~pk/research/matlabfns/
[34] http://cs.adelaide.edu.au/~tjchin/doku.php