π

# MPOGames: Efficient Multimodal Partially Observable Dynamic Games

Oswin So$1,2,*$, Paul Drews$2$, Thomas Balch$2$, Velin Dimitrov$2$, Guy
Rosman$2$, Evangelos A. Theodorou$1$ $1$ Georgia Institute of Technology.$2$
Toyota Research Institute.$*$ Corresponding Author<EMAIL_ADDRESS>
This work was supported by the Toyota Research Institute (TRI). This article
solely reflects the opinions and conclusions of its authors and not TRI or any
other Toyota entity. Their support is gratefully acknowledged.

###### Abstract

Game theoretic methods have become popular for planning and prediction in
situations involving rich multi-agent interactions. However, these methods
often assume the existence of a single local Nash equilibria and are hence
unable to handle uncertainty in the intentions of different agents. While
maximum entropy (MaxEnt) dynamic games try to address this issue, practical
approaches solve for MaxEnt Nash equilibria using linear-quadratic
approximations which are restricted to unimodal responses and unsuitable for
scenarios with multiple local Nash equilibria. By reformulating the problem as
a POMDP, we propose MPOGames, a method for efficiently solving MaxEnt dynamic
games that captures the interactions between local Nash equilibria. We show
the importance of uncertainty-aware game theoretic methods via a two-agent
merge case study. Finally, we prove the real-time capabilities of our approach
with hardware experiments on a 1/10th scale car platform.

[b]2023 IEEE International Conference on Robotics and Automation (ICRA)

Figure 1:  The ego agent (left) seeks to merge into a lane safely without
collisions. With multiple local minima present, MPOGames predicts likely plans
for other agents, infers their probabilities, then hedges against all
possibilities. Methods that fail to do so may result in catastrophic
collisions.

##  I Introduction and Related Work

In recent years, robots have been seen increasing use in applications with
multiple interacting agents. In these applications, actions taken by the robot
cannot be considered in isolation and must take the responses of other agents
into account, especially in areas such as autonomous driving [1], robotic arms
[2] and surgical robotics [3] where failures can be catastrophic. For example,
in the autonomous driving setting, two agents at a merge must negotiate and
agree on which agent merges first (Fig. 1). Failure to consider the responses
of the other agent could result in serious collisions.

Game-theoretic (GT) approaches have become increasingly popular for planning
in multi-agent scenarios [4, 5, 6, 7, 8, 9, 10, 11, 12]. By assuming that
agents act rationally, these approaches decouple the problem of planning and
prediction into the problems of finding objective functions for each agent and
solving for Nash equilibrium of the resulting non-cooperative game. This
provides a more principled and interpretable alternative to the prediction of
other agents when compared to black-box approaches that use deep neural
networks [13, 14].

In practice, assuming rational agents is too strict. Humans face cognitive
limitations and make irrational decisions that are satisfactory but
suboptimal, an idea known as “bounded rationality” or “noisy rationality” [15,
16]. Even autonomous agents are rarely exactly optimal and often make
suboptimal decisions due to finite computational resources. One method of
addressing this disconnect is the Maximum Entropy (MaxEnt) framework, which
has been applied to diverse areas such as inverse reinforcement learning [17,
7], forecasting [18] and biology [19]. In particular, this combination of
MaxEnt and GT has been proposed for learning objectives for agents from a
traffic dataset [7]. However, since the exact solution is computationally
intractable, the authors use a linear quadratic (LQ) approximation of the game
for which a closed-form solution can be efficiently obtained. The LQ
approximation, popular among recent approaches (e.g., [5, 7]), only considers
a single local minimum and consequently produces unimodal predictions for each
agent. This fails to capture complex uncertainty-dependent hedging behaviors
which mediate between other agents’ possible plans in the true solution of the
MaxEnt game.

These hedging behaviors are closely related to partially observable Markov
decision process (POMDP) in the single-player setting, where solutions must
consider the distribution of states when solving for the optimal controls.
Similar multi-agent scenarios have been considered via POMDPs in [20, 21, 22].
Also related are partially observable stochastic games which are usually
intractable [1] and mostly limited to the discrete setting [23, 24, 25, 26].

One big challenge for methods that seek to address partial observability and
multimodality is the problem of computational efficiency and scalability
(e.g., in the case of POMDPs). Furthermore, only a few of the existing GT
methods perform experiments on hardware [10, 11, 27], where the robustness of
the method to noise and delays in the system is crucial. To tackle both of
these issues, we present MPOGames, a framework for efficiently solving
multimodal MaxEnt dynamic games by reformulating the problem as a game of
incomplete information and efficiently finding approximate solutions to the
equivalent POMDP problem. We showcase the real-time potential of our method
via extensive hardware comparisons.

Our contributions in this work are as follows.

  1. 1.

A computationally efficient framework for GT planning in multi-agent scenarios
that produces more accurate solutions in the presence of multiple local minima
compared to previous works.

  2. 2.

A case-study in simulation comparing how GT methods handle multimodality that
show the importance of inference and handling uncertainty.

  3. 3.

Analysis of the performance of GT methods in multimodal situations with a
1/10th scale car platform. These results demonstrate the robustness and
tractability of MPOGames on hardware.

##  II Preliminaries

###  II-A Discrete Dynamic Games and Nash Equilibria

We consider a discrete dynamic game with $N$ players. Let
${\mathbf{u}}_{t}=[u_{t}^{1},\dots,u_{t}^{N}]=[u_{t}^{i},{\mathbf{u}}_{t}^{\neg
i}]\in\mathbb{R}^{n_{{\mathbf{u}}}}$ denote the joint control vector for all
players, where we use $(\cdot)^{i}$ and $(\cdot)^{\neg i}$ to denote the
partition of $(\cdot)$ into parts belonging to agent $i$ and other agents
respectively. We denote by ${\mathbf{x}}_{t}\in\mathbb{R}^{n_{{\mathbf{x}}}}$
the full state of the system at timestep $t$ with nonlinear dynamics

${\mathbf{x}}_{t+1}=f({\mathbf{x}}_{t},u_{t}^{1},\dots,u_{t}^{N})=f({\mathbf{x}}_{t},{\mathbf{u}}_{t}).$
(1)

Time indices are dropped when they are clear from the context. Each agent’s
objective is to minimize a corresponding finite-horizon cost function $J^{i}$
with running cost $l^{i}$ and terminal cost $\Phi^{i}$ with respect to the
control trajectory ${\mathbf{U}}=[{\mathbf{u}}_{1},\dots,{\mathbf{u}}_{T-1}]$

$\min_{U^{i}\in\mathcal{U}^{i}}J^{i}(U^{i},{\mathbf{U}}^{{\neg
i}})\;=\;\min_{U^{i}\in\mathcal{U}^{i}}\Phi^{i}({\mathbf{x}}_{T})+\sum_{t=1}^{T-1}l^{i}({\mathbf{x}}_{t},{\mathbf{u}}_{t}),$
(2)

where $\mathcal{U}^{i}$ denotes the set of feasible control trajectories for
agent $i$. Note that the objective function $J^{i}$ is a function of both
$U^{i}$ and ${\mathbf{U}}^{\neg i}$, but the optimization is only over
$U^{i}$. Hence, the quantity of interest here is the Nash Equilibrium (NE),
i.e., find controls ${\mathbf{U}}^{*}$ such that for each agent $i$,

$J^{i}(U^{i*},{\mathbf{U}}^{{\neg i}*})\leq J^{i}(U^{i},{\mathbf{U}}^{{\neg
i}*}),\quad\forall U^{i}\in\mathcal{U}^{i}.$ (3)

Figure 2: Optimal policies for the ego agent P${\mathfrak{1}}$ (left) and non-
ego agent P${\mathfrak{2}}$ (right) on the toy problem with $\epsilon=0.1$.
For ${\color[rgb]{0.0234375,0.5546875,0.83203125}\pi^{\mathfrak{1}}}$, only
MaxEnt and POMDP hedges against both local NE of $\pi^{\mathfrak{2}}$. The
deterministic policy of NE is denoted by the Dirac delta.

###  II-B Maximum Entropy Dynamic Games

Agents rarely act exactly optimally due to cognitive or computational
limitations and settle for satisfactory decisions [16]. We can model
suboptimality via the noisy rationality model [15], where controls are
stochastic according to the MaxEnt framework [28, 17]. Consider a stochastic
control policy $\pi^{i}(u^{i},{\mathbf{x}})$ and introduce an entropy
regularization term to the original objective (2). We then consider the
expected cost under all agents’ policies $\vpi$

$J^{i}(\vpi)=\operatorname*{\mathbb{E}}_{{\mathbf{u}}_{t}\sim\vpi}\left[\Phi^{i}({\mathbf{x}}_{T})+\sum_{t=1}^{T-1}l^{i}({\mathbf{x}}_{t},{\mathbf{u}}_{t})-\frac{1}{\beta}H[\pi^{i}(\cdot|{\mathbf{x}}_{t})]\right],$
(4)

where $\beta>0$ denotes inverse temperature or the rationality coefficient,
and $H[\pi]$ is the Shannon entropy of $\pi$, defined as

$H[\pi]\coloneqq-\operatorname*{\mathbb{E}}_{u\sim\pi}[\log\pi(u)]=-\int\pi(u)\log\pi(u)\differential{u}.$
(5)

The resulting MaxEnt NE problem (referred to as “Entropic Cost Equilibrium” in
[7]) is then to find stochastic policies $\vpi^{*}$ such that the following
holds for each agent

$J^{i}(\pi^{i*},\pi^{{\neg i}*})\leq J^{i}(\pi^{i},\pi^{{\neg
i}*}),\quad\forall\pi^{i}\in\Pi^{i},$ (6)

for feasible control policies $\Pi^{i}$ from any initial state. Let the value
function $V^{i}$ for agent $i$ given the policies of other agents $\vpi^{{\neg
i}}$ be

$V^{i}_{t}({\mathbf{x}}_{t})\coloneqq\inf_{\pi^{i}}\\{\,J^{i}({\mathbf{x}},\pi^{i},\vpi^{{\neg
i}})\,\\}.$ (7)

The optimal policy $\pi^{i*}$ takes the form [7, 9]

$\pi^{i*}(u^{i}|{\mathbf{x}}_{t})={Z^{i}({\mathbf{x}}_{t})}^{-1}\exp\left(-\beta
Q({\mathbf{x}}_{t},u^{i})\right),$ (8)

where

$Q^{i}({\mathbf{x}}_{t},u^{i}_{t})\coloneqq\operatorname*{\mathbb{E}}_{\vpi^{{\neg
i}}}\left[{V^{i}}_{t+1}({\mathbf{x}}_{t+1})+l^{i}({\mathbf{x}}_{t},{\mathbf{u}})\right],$
(9)

and $Z^{i}$ denotes the partition function

$Z^{i}({\mathbf{x}}_{t})\coloneqq\int\exp\left(-\beta\operatorname*{\mathbb{E}}_{\vpi^{{\neg
i}}}\left[{V^{i}}_{t+1}({\mathbf{x}}_{t+1})+l^{i}({\mathbf{x}}_{t},{\mathbf{u}})\right]\right)\differential{u^{i}}.$
(10)

Then, the value function $V^{i}$ takes the form [29, Appendix A]

$V^{i}_{t}({\mathbf{x}}_{t})=-\ln Z^{i}({\mathbf{x}}_{t}).$ (11)

While this gives us an expression for $\pi^{i*}$ and $V^{i}$, it is generally
intractable to compute the integral (10) and solve for $\vpi^{*}$ and $Z^{i}$
under general nonlinear costs and dynamics.

One way of resolving this problem is to approximate the problem as a linear-
quadratic (LQ) game [7, 9] in a manner similar to iLQR [30], an approach also
taken in the iLQGames [5, 31] family of NE solvers. By doing so, the
approximate game can be solved exactly [9, Lemma 1 and Lemma 2].

###  II-C LQ is problematic for Multiodal MaxEnt NE

While taking LQ approximations yields a computationally efficient closed-form
simulation for both the original NE and the MaxEnt NE problem, this
approximation can be especially problematic in the MaxEnt NE setting. While
the value function in the deterministic case does not depend on the cost
function at states away from the NE, this is not the case for MaxEnt NE due to
the integral in (8) and (11).

Illustrative Toy Example.  We illustrate this using a toy example where the
MaxEnt NE can be solved for in closed-form. Consider the following two-agent
single timestep dynamic game with
${\mathbf{x}}=[{\color[rgb]{0.0234375,0.5546875,0.83203125}x^{\mathfrak{1}}},{\color[rgb]{0.98828125,0.30859375,0.1875}x^{\mathfrak{2}}}]\in\mathbb{R}^{2}$
and
${\mathbf{u}}=[{\color[rgb]{0.0234375,0.5546875,0.83203125}u^{\mathfrak{1}}},{\color[rgb]{0.98828125,0.30859375,0.1875}u^{\mathfrak{2}}}]\in\mathbb{R}^{2}$
with single-integrator dynamics and costs defined by

$\displaystyle{\color[rgb]{0.0234375,0.5546875,0.83203125}J^{\mathfrak{1}}}$
$\displaystyle=\frac{1}{2}({\color[rgb]{0.0234375,0.5546875,0.83203125}u_{0}^{\mathfrak{1}}})^{2}+\frac{3}{2}({\color[rgb]{0.0234375,0.5546875,0.83203125}x_{1}^{\mathfrak{1}}}-{\color[rgb]{0.98828125,0.30859375,0.1875}x_{1}^{\mathfrak{2}}})^{2},$
(12)
$\displaystyle{\color[rgb]{0.98828125,0.30859375,0.1875}J^{\mathfrak{2}}}$
$\displaystyle=\frac{1}{2}({\color[rgb]{0.98828125,0.30859375,0.1875}u_{0}^{\mathfrak{2}}})^{2}+\Phi^{\mathfrak{2}}({\color[rgb]{0.98828125,0.30859375,0.1875}x_{1}^{\mathfrak{2}}}),$
(13)
$\displaystyle\Phi^{\mathfrak{2}}({\color[rgb]{0.98828125,0.30859375,0.1875}x_{1}^{\mathfrak{2}}})$
$\displaystyle=\operatorname*{\sigma_{\textrm{SM}}}\left(\frac{3}{2}({\color[rgb]{0.98828125,0.30859375,0.1875}x_{1}^{\mathfrak{2}}}-1)^{2},\;\frac{3}{2}({\color[rgb]{0.98828125,0.30859375,0.1875}x_{1}^{\mathfrak{2}}}+1)^{2}+\epsilon\right),$
(14) $\displaystyle\operatorname*{\sigma_{\textrm{SM}}}(a,b)$
$\displaystyle\coloneqq-\ln(e^{-a}+e^{-b}),$ (15)

where $\epsilon>0$ controls the “height” of the left local minima (Fig. 2 top
right) and
${\color[rgb]{0.0234375,0.5546875,0.83203125}x_{0}^{\mathfrak{1}}}={\color[rgb]{0.98828125,0.30859375,0.1875}x_{0}^{\mathfrak{2}}}=0$.
The $\operatorname*{\sigma_{\textrm{SM}}}$ in (14) acts as a smooth minimum of
the two quadratic state costs. In other words, P${\mathfrak{2}}$ wants to
reach either $-1$ or $1$ (with a small preference for $1$), while
P${\mathfrak{1}}$ wants to track P${\mathfrak{2}}$.

NE:  In the regular NE case, the optimization problem to solve for
${\color[rgb]{0.0234375,0.5546875,0.83203125}u_{0}^{\mathfrak{1}}}$ is exactly
quadratic and gives
${\color[rgb]{0.0234375,0.5546875,0.83203125}u_{0}^{{\mathfrak{1}}*}}=\frac{3}{4}{\color[rgb]{0.98828125,0.30859375,0.1875}u_{0}^{\mathfrak{2}}}$.
Note that solution is independent of $\epsilon$. We can also solve for
${\color[rgb]{0.98828125,0.30859375,0.1875}u_{0}^{{\mathfrak{2}}*}}$
numerically to obtain
${\color[rgb]{0.98828125,0.30859375,0.1875}u_{0}^{{\mathfrak{2}}*}}={\color[rgb]{0.98828125,0.30859375,0.1875}x_{1}^{{\mathfrak{2}}*}}\approx
0.73$.

Exact MaxEnt NE:  When solving the MaxEnt NE case,
${\color[rgb]{0.98828125,0.30859375,0.1875}\pi^{{\mathfrak{2}}*}}$ is a
Gaussian mixture. Note the “shape” of the cost function
${\color[rgb]{0.98828125,0.30859375,0.1875}J^{\mathfrak{2}}}$ determines the
form of ${\color[rgb]{0.98828125,0.30859375,0.1875}\pi^{\mathfrak{2}}}$.
Consequently, the optimal
${\color[rgb]{0.0234375,0.5546875,0.83203125}\pi^{{\mathfrak{1}}*}}$ is
Gaussian with mean $0.25\tanh(\epsilon/2)$. In comparison to the NE case, both
${\color[rgb]{0.0234375,0.5546875,0.83203125}\pi^{{\mathfrak{1}}*}}$ and
${\color[rgb]{0.98828125,0.30859375,0.1875}\pi^{{\mathfrak{2}}*}}$ depend on
$\epsilon$. Moreover, for $\epsilon\to 0$, the mean of
${\color[rgb]{0.0234375,0.5546875,0.83203125}\pi^{{\mathfrak{1}}*}}$
approaches zero, which is qualitatively different behavior to the previous
case.

LQ Approx MaxEnt NE:  Using the LQ approximation, we obtain the following two
sets of local LQ NE:

$\displaystyle{\color[rgb]{0.0234375,0.5546875,0.83203125}\pi^{{\mathfrak{1}}*,(1)}}$
$\displaystyle\approx\mathcal{N}(+0.55,0.5),$
$\displaystyle{\color[rgb]{0.98828125,0.30859375,0.1875}\pi^{{\mathfrak{2}}*,(1)}}$
$\displaystyle\approx\mathcal{N}\left(+0.73,0.53\right),$
$\displaystyle{\color[rgb]{0.0234375,0.5546875,0.83203125}\pi^{{\mathfrak{1}}*,(2)}}$
$\displaystyle\approx\mathcal{N}(-0.55,0.5),$
$\displaystyle{\color[rgb]{0.98828125,0.30859375,0.1875}\pi^{{\mathfrak{2}}*,(2)}}$
$\displaystyle\approx\mathcal{N}\left(-0.73,0.53\right).$

For each method, we plot $\pi$ for both agents in Fig. 2. While the true
MaxEnt ${\color[rgb]{0.0234375,0.5546875,0.83203125}\pi^{\mathfrak{1}}}$
hedges against both modes of
${\color[rgb]{0.98828125,0.30859375,0.1875}\pi^{\mathfrak{2}}}$ with near zero
mean, neither LQ approximations capture this.

##  III Solving MaxEnt Dynamic Games as POMDPs

Figure 3: Diagram of the proposed algorithm. In this work, we use QMDP for the
POMDP solver. Belief updates and QMDP solving are performed at a much faster
rate than the computation of local NE.

###  III-A POMDPs: Better Approximations for MaxEnt NE

The main problem with LQ approximations for MaxEnt NE is that they work poorly
when the non-ego agents’ policies are not well approximated by a unimodal
Gaussian. In the previous example, both LQ approximations alone are terrible
approximations of true MaxEnt policy for $\pi^{\mathfrak{2}}$ (Fig. 2). On the
other hand, if we use a more complex distribution such as a mixture model to
model the optimal policy for each timestep, the number of mixture components
will grow exponentially in the time horizon and is not computationally
tractable.

To resolve this issue, we consider the following reformulation of the problem
proposed in [9]. Taking inspirations from bounded rationality, we assume the
non-ego agents are only capable of solving MaxEnt NE with LQ structure, but
the ego agent does not know which LQ approximation is used. To model this,
consider a set of $M$ different LQ NE
$\\{({\mathbf{X}}^{(z)},{\mathbf{U}}^{(z)})_{z=1}^{Z}\\}$, and let
$z\in\\{1,\dots,Z\\}$ be a discrete latent random variable that determines
which approximation is used by the non-ego agents. We now consider an
extension of the MaxEnt Dynamic Game (6) where all non-ego agents have full
knowledge of the value of $z$ and are playing the corresponding optimal policy
$\pi^{(z)}$. However, the non-ego agents incorrectly assume the ego agent
knows what the true mode is. This is now a dynamic game with incomplete
information — the ego agent only has a belief $b$ of what the true game is but
does not know what the true dynamics are [32, 9]. As shown in [9], the problem
of computing the Bayesian equilibrium [32, Chapter 10] is equivalent to
solving the following POMDP problem with state
$\tilde{{\mathbf{x}}}\coloneqq[{\mathbf{x}},b]$ (where we assume the ego agent
is agent ${\mathfrak{1}}$)

$V_{t}^{\mathfrak{1}}(\tilde{{\mathbf{x}}}_{t})=\inf_{\pi^{\mathfrak{1}}}\operatorname*{\mathbb{E}}_{z\sim
b_{t}}\operatorname*{\mathbb{E}}_{{\mathbf{u}}_{t}}\left[l^{\mathfrak{1}}({\mathbf{x}}_{t},{\mathbf{u}}_{t})+V_{t+1}^{\mathfrak{1}}(\tilde{{\mathbf{x}}}_{t+1})\right],$
(16)

where ${\mathbf{u}}^{{\neg\mathfrak{1}}}\sim\vpi^{{\neg\mathfrak{1}},z}$ and
the dynamics of the belief $b$ follow the Bayesian filtering update law

$b_{t+1}(z)\propto\pi^{{\neg\mathfrak{1}},z}_{t}({\mathbf{u}}^{{\neg\mathfrak{1}}}_{t}|z)\,b_{t}(z).$
(17)

Since computing $b_{t}$ requires $\pi^{{\neg\mathfrak{1}},z}_{t-1}$, the LQ NE
is solved from ${\mathbf{x}}_{t-1}$. We follow [9] and choose the prior belief
$b_{0}$ as the MaxEnt distribution over the sum of agent value functions

$b_{0}=\inf_{b}\operatorname*{\mathbb{E}}_{z\sim
b}\left[\sum_{i=1}^{N}V^{i,z}({\mathbf{x}}_{0})\right]-\frac{1}{\beta}H[b].$
(18)

######  Remark.

Unlike the usual POMDP setup in literature with discrete state and action
spaces [33, 34, 35], the above setup can be considered a POMDP with hybrid
states $[{\mathbf{x}},z]$ (${\mathbf{x}}$ is continuous but $z$ is discrete),
infinite-dimensional controls $\pi^{\mathfrak{1}}$ and continuous observations
${\mathbf{x}}_{t}$. Since ${\mathbf{x}}$ is fully observable, the belief
$\tilde{b}$ is degenerate on the ${\mathbf{x}}$ coordinates.

###  III-B MPOGames

Leveraging this POMDP reformulation, we now present the MPOGames algorithm.
First, we solve for a set of different LQ MaxEnt NE, which can be done in
parallel. Since we assume each non-ego agent has picked one of these local NE
and executes $\vpi^{{\neg\mathfrak{1}},z}$, we perform Bayesian inference
using (17) to form a belief $b$ over $z$. This belief is used to solve a POMDP
to give an optimal policy $\pi^{{\mathfrak{1}}}$ for the ego agent. In
practice, only the mean ego policy is used. A diagram of MPOGames is
summarized in Fig. 3. Intuitively, this reformulation improves the fidelity of
$\pi^{{\neg\mathfrak{1}}}$ from a single to a mixture of Gaussians. This
enables a more accurate expectation over non-ego agents’ policies
$\pi^{{\neg\mathfrak{1}}}$ in (9), which is key to the behavior of the MaxEnt
solution.

Returning to the toy example, MPOGames recovers similar hedging behavior to
the exact MaxEnt case (ego P${\mathfrak{1}}$).

MPOGames:  The value functions for P${\mathfrak{1}}$ take on the values
$V^{{\mathfrak{1}},(1)}=0.17$ and $V^{{\mathfrak{1}},(2)}=0.09$. Consequently,
for $\epsilon=0.1$,

$b_{0}=[0.52,0.48],\quad{\color[rgb]{0.0234375,0.5546875,0.83203125}\pi^{{\mathfrak{1}}}}=\mathcal{N}\left(0.0,0.5\right),$
(19)

which displays similar hedging behavior as the exact MaxEnt NE solution (see
Fig. 2 left).

###  III-C Approximate Solutions for POMDPs

Note that MPOGames does not prescribe a POMDP solver. However, exact solutions
for POMDPs are generally intractable [33]. Two popular approximations for
solving POMDPs are $Q_{\text{MDP}}$ [35] and fast informed bound (FIB) [34].
$Q_{\text{MDP}}$ has been used successfully in practice [36] and performs well
when the particular action has little impact on the reduction in state
uncertainty [33]. It can be shown that when the state ${\mathbf{x}}_{t}$ of
all agents is directly observable, the $Q_{\text{MDP}}$ and FIB approximations
are identical and equal to the MaxEnt NE value function. Moreover, if the
value function $V$ for each mode is known, then the $Q_{\text{MDP}}$ can be
computed easily.

An alternative approach for solving POMDPs is to use point-based methods which
yield more accurate solutions [33]. Instead of changing the information
structure (i.e., assuming the state is known), these methods work with the
full value function on beliefs. For example, [20, 37, 38] consider a scenario
tree based structure for solving the POMDP. However, these methods still scale
exponentially in the time horizon and have not been demonstrated to run in
real-time in existing works [20, 21]. In comparison, QMDP with known $V$ has
constant complexity in the time horizon.

In this work, we implement MPOGames using the QMDP approximation for its
effectiveness at low computational costs. While the value function
$V^{{\mathfrak{1}},z}$ for each mode is not known exactly, we take the
quadratic approximation of $V$ obtained from solving the (feedback) NE for
each mode. The POMDP is solved with QMDP using the new belief and value
functions, i.e.,

$\min_{\pi^{{\mathfrak{1}}}_{t}}\operatorname*{\mathbb{E}}_{u^{{\mathfrak{1}}}_{t}}\operatorname*{\mathbb{E}}_{z\sim
b_{t}}\left[l^{{\mathfrak{1}},z}+V^{{\mathfrak{1}},z}_{t+1}({\mathbf{x}}^{z}_{t+1})\right],$
(20)

and the resulting solution is sent to a lower level trajectory tracker. Due to
the efficiency of QMDP, we are able to solve it much faster than we can solve
for local NE. Hence, we decouple the two parts in our algorithm to enable more
responsive belief updates. The low level trajectory tracking is also decoupled
to enable fast disturbance rejection from model errors in practice. The three
decoupled nodes and their respective runtimes are presented in Fig. 3.

##  IV Simulation experiments

###  IV-A Autonomous Vehicle Merge

Figure 4:  Top Left, Top Right: Simulation snapshots of a merge scenario. The
ego agent is on the left and in gold. The circles denote the radius for each
robot. Red circles indicate collision. Being unable to consider multiple modes
can lead to catastrophic failures such as collisions or freezing if the
inferred mode is incorrect. Bottom Left: Visualization of the two local Nash
equilibria for this problem. Bottom Right: Inferred probabilities for each
mode via MaxEnt inference.

We consider a merge scenario with two agents merging into the same lane while
respecting track boundaries and avoiding collisions and denote the ego agent
by P0.

Agents are modeled as kinematic bicycles. To easily handle both types of state
constraints, we consider a hybrid Cartesian-Curvilinear state representation
similar to [39]. Let
$x_{i}=[p_{x},p_{y},v,\theta,\zeta,s,n,\xi]\in\mathbb{R}^{8}$ and
$u^{i}=[\differential{\zeta},a]$ denote the state and control for agent $i$.
The tracks are represented as a constant width lane around a centerline. Both
the signed curvature $\kappa$ and the Cartesian coordinates of the centerline
are represented as a cubic spline in the arclength $s$.

The cost for each agent consists of a quadratic cost for staying close to the
centerline, maintaining a target velocity, and a soft constraint for staying
inside the road boundaries. We handle collisions by treating them as soft-
costs using a quadratic cost function on the Cartesian coordinates. Since this
is concave in the position, there are two local NE — either the ego agent goes
first, or the non-ego agent goes first (Fig. 4). Consequently, we do not know
which NE the non-ego agent is considering and must perform inference.

We consider the following game-theoretic methods in this scenario:

  1. 1.

(NoInf) No inference is done for the non-ego agent’s considered mode, and we
assume that the non-ego agent considers the same mode as the one we are
considering. This is the approach taken by methods such [5, 4] which only
consider a single local NE.

  2. 2.

(ML) Inference is done. However, instead of solving the POMDP, the controls
corresponding to the maximum-likelihood (ML) mode are used as in [12].

  3. 3.

(MPOGames) Inference is done, and the resulting POMDP is solved using the QMDP
approximation.

Although both ML and MPOGames perform inference, ML naively takes the most
probable mode and ignores the rest. If the inferred probabilities are
$b=[0.5+\epsilon,0.5-\epsilon]$ for small $\epsilon$, ML applies the policy
for the first mode, while MPOGames will consider both modes almost equally.

Each method is run in a receding-horizon fashion with the results shown in
Fig. 4 and Fig. 5. While NoInf performs well when the considered mode is the
same as the non-ego agent, it fails when the mode is different — either both
agents believe the other agent will yield and crash into each other, or both
agents yield resulting in both agents stopping. Despite modeling the response
of other agents, the “freezing” behavior here happens because the predicted
response is incorrect. Methods that perform inference (ML, MPOGames) infer the
correct mode in the noiseless case and successfully perform the merge. Note
that the likelihood of the wrong hypothesis when using ML actually increases
at the end of the trajectory in Fig. 5. This occurs because the trajectory
prefix of both modes becomes increasingly similar at the merge point despite
the cost of the wrong mode being much larger.

Figure 5:  In contrast to the MaxEnt probabilities in Fig. 4, naive inference
results in both modes having similar probabilities near the end (eg. dashed
line) due to both modes having similar prefixes at the merge point, despite
the wrong mode having much higher cost. Top: Inferred probabilities for each
mode via MaxEnt inference. Left: Prediction error for each mode. Right: Cost
of each predicted mode.

We now consider the case where the non-ego agent’s plans are not perfectly
rational and perturb the non-ego agent’s controls slightly such that there is
large uncertainty regarding the non-ego agent’s mode (Fig. 6). The high
uncertainty of the true mode means that the maximum-likelihood mode oscillates
between the two modes. Consequently, the discontinuity of ML results in a
discontinuous, jerky ego control. In comparison, MPOGames hedges between the
two possibilities and outputs a smooth, continuous control.

Choice of inference technique.  Additionally, we investigate the difference
between a naive method of inference via a user-specified Gaussian measurement
model [12] and the MaxEnt inference presented in this paper. The results for
naive inference are shown in Fig. 5 (cf. Fig. 4). Despite inferring the
correct mode, this probability actually decreases towards the end despite
larger costs.

Figure 6: Comparison between the ML and MPOGames methods when the non-ego
agent uses suboptimal controls. The discontinuous nature of ML solution means
that small perturbations in the non-ego agent’s controls (top) results in a
discontinuous, jerky ego control (center) when the uncertainty oscillates
around $0.5$ (bottom). On the other hand, MPOGames is a continuous function of
the belief with smoother controls.

Implementation Details. The symmetric structure of the coupling costs in the
problems we consider makes this a potential game. Consequently, we solve for
the (feedback) NE by reformulating the problem as an optimal control problem
in a manner similar to [27, 40]. However, we note that the NE can be solved
with other “backends” such as iLQGames [5], ALGames [4] or SQP [39], with the
value function and policy recovered by performing an additional backwards pass
of iLQGames. Next, although the belief update (17) requires access to the
control inputs of other players, this information may not be directly
available. Instead, we perform a least-squares estimate of
${\mathbf{u}}_{t}^{{\neg\mathfrak{1}}}$ using ${\mathbf{x}}_{t}$ and
${\mathbf{x}}_{t+1}$. We also assume prior knowledge on how to find a
comprehensive set of NE. In this work, this is done using a heuristic
initialization process and fixed to two NE. Determining how many NE exist and
how to find this set of NE is beyond the scope of this work. We use CasADi
[41] and acados [42] to solve the resulting optimal control problem with
$\Delta t=$0.05\text{\,}\mathrm{s}$$ and horizon $T=60$.

##  V Hardware experiments

Finally, we validate our approach on hardware in a merge scenario. We use the
TRIKart [43], a modified version of the F1Tenth platform [44]. Optitrack [45]
is used for localization.

###  V-A No Perturbations

We first investigate how each method performs in the cases where the non-ego
agent is cooperative (Yield), selfish (NoYield), or using the ML or MPOGames
strategies. The first two cases (Yield, NoYield) are realized by performing
NoInf with an initialization strategy that is biased towards yielding or not
yielding respectively. The success rate of merging over ten trials is shown in
Table I, where success is defined loosely as both agents reaching the center
lane without colliding or stopping.

Yield, NoYield — Mismatch of Local NE:  As predicted from Section IV-A,
failures occur when there is a mismatch in local NE between the two ages
(i.e., Yield-Yield and NoYield-NoYield). In the Yield-Yield case, all failures
arose from both agents coming to a complete stop before the merge point. Due
to noise, one of the agents eventually merges, but this only happens after a
long delay. The two successes in the NoYield-NoYield example arose when one
car was much closer to the merge point than the other car due to hardware
differences, resolving the ambiguity. However, NoYield-NoYield was also the
only combination that resulted in collisions (see Fig. 7).

ML — Vulnerablility to system latency:  The simulation experiments in Section
IV-A make the standard assumption that there is no delay between receiving
observations and outputting the controls. However, this is not true in the
real world. In the ML-ML setup, although all trials result in successful
merges, both agents exhibit large oscillatory controls (Fig. 10). This
behavior is inherent to the discontinuous behavior of ML. Despite having the
same system latency, MPOGames-MPOGames handles this situation more gracefully
due to accounting for uncertainty (Fig. 10).

TABLE I:  Merge success rate on the scale car platform over ten trials from
the same initial conditions.

|  | Non-ego Agent
---|---|---
|  | Yield | NoYield | ML | MPOGames
Ego Agent | Yield | 0.0 |  |  |
NoYield | 1.0 | 0.2 |  |
ML | 1.0 | 1.0 | 1.0 |
MPOGames | 1.0 | 1.0 | 1.0 | 1.0

Figure 7: NoYield-NoYield setup. Collisions occur when agents do not agree on
the same local NE. Left is ego agent.

Figure 8: Results on the scale car hardware when both agents use ML. The
processing latency leads to an unstable, oscillatory system. Top: Output
accelerations. Bottom: inferred probabilities of $P0$ (Ego) going first.

Figure 9: Same setup as Fig. 10, but both agents use MPOGames.

Figure 10: Comparison of accelerations between ML and MPOGames on hardware
when the non-ego agent applies high uncertainty controls. Note the oscillation
for ML.

###  V-B Controls with High Uncertainty

Finally, we compare how each method performs when the non-ego agent has
explicitly suboptimal controls. Specifically, we add sinusoidal noise to the
non-ego agent’s controls such that it is difficult to determine the non-ego
agent’s intentions. Here, we only compare ML and MPOGames and show the results
in Fig. 10. This is similar to the scenario in Fig. 6, and indeed we see the
same oscillations that ML displays on hardware. However, in this case, the
non-ego agent is intentionally suboptimal. Again, MPOGames handles this
uncertainty gracefully and correctly hedges against both NE.

##  VI Discussion

Summary:  We have proposed a MaxEnt game-theoretic framework for planning and
inference in multi-agent situations with multiple local Nash equilibrium.
Through simulation and hardware case studies, we have shown that MPOGames
provides a real-time method for hedging against the uncertainty of other
agents and improves performance compared to previous game-theoretic
approaches.

Limitations and future work: Using a GT planner for joint planning and
prediction assumes the knowledge of the cost function and dynamics for all
agents. However, this is rarely true and consequently must be estimated using
techniques such as inverse reinforcement learning [17, 46]. While the proposed
reformulation into a POMDP relies on the key assumption that the non-ego
agents know the true mode, we demonstrate that this restriction may be relaxed
(e.g., MPOGames-MPOGames in Section V). Also, our framework assumes that non-
ego agents are only capable of solving LQ games. Relaxing this assumption to
model more intelligent non-ego agents is left as future work. Moreover, we
have assumed prior knowledge on finding all NE. Doing so efficiently for
general games is left as future work. Finally, the current work assumes that
there are a discrete set of modes, but there exist scenarios when there are a
continuum of NE. We leave this direction as future work.

## References

  * [1] D. Sadigh, N. Landolfi, S. S. Sastry, S. A. Seshia, and A. D. Dragan, “Planning for cars that coordinate with people: leveraging effects on human actions for planning and active information gathering over human internal state,” _Autonomous Robots_ , vol. 42, no. 7, pp. 1405–1426, 2018.
  * [2] H. B. Amor, G. Neumann, S. Kamthe, O. Kroemer, and J. Peters, “Interaction primitives for human-robot cooperation tasks,” in _2014 IEEE international conference on robotics and automation (ICRA)_. IEEE, 2014, pp. 2831–2837.
  * [3] T. Zhou and J. P. Wachs, “Early prediction for physical human robot collaboration in the operating room,” _Autonomous Robots_ , vol. 42, no. 5, pp. 977–995, 2018.
  * [4] S. L. Cleac’h, M. Schwager, and Z. Manchester, “Algames: A fast solver for constrained dynamic games,” _arXiv preprint arXiv:1910.09713_ , 2019. [Online]. Available: https://arxiv.org/abs/1910.09713
  * [5] D. Fridovich-Keil, E. Ratner, L. Peters, A. D. Dragan, and C. J. Tomlin, “Efficient iterative linear-quadratic approximations for nonlinear multi-player general-sum differential games,” in _2020 IEEE International Conference on Robotics and Automation (ICRA)_. IEEE, 2020, pp. 1475–1481.
  * [6] F. Laine, D. Fridovich-Keil, C.-Y. Chiu, and C. Tomlin, “Multi-hypothesis interactions in game-theoretic motion planning,” in _2021 IEEE International Conference on Robotics and Automation (ICRA)_. IEEE, 2021, pp. 8016–8023.
  * [7] N. Mehr, M. Wang, and M. Schwager, “Maximum-entropy multi-agent dynamic games: Forward and inverse solutions,” _arXiv preprint arXiv:2110.01027_ , 2021\. [Online]. Available: https://arxiv.org/abs/2110.01027
  * [8] W. Schwarting, A. Pierson, S. Karaman, and D. Rus, “Stochastic dynamic games in belief space,” _IEEE Transactions on Robotics_ , 2021.
  * [9] O. So, K. Stachowicz, and E. A. Theodorou, “Multimodal maximum entropy dynamic games,” _arXiv preprint arXiv:2201.12925_ , 2022. [Online]. Available: https://arxiv.org/abs/2201.12925
  * [10] Z. Wang, R. Spica, and M. Schwager, “Game theoretic motion planning for multi-robot racing,” in _Distributed Autonomous Robotic Systems_. Springer, 2019, pp. 225–238.
  * [11] Z. Wang, T. Taubner, and M. Schwager, “Multi-agent sensitivity enhanced iterative best response: A real-time game theoretic planner for drone racing in 3d environments,” _Robotics and Autonomous Systems_ , vol. 125, p. 103410, 2020.
  * [12] L. Peters, D. Fridovich-Keil, C. J. Tomlin, and Z. N. Sunberg, “Inference-based strategy alignment for general-sum differential games,” _arXiv preprint arXiv:2002.04354_ , 2020. [Online]. Available: https://arxiv.org/abs/2002.04354
  * [13] T. Salzmann, B. Ivanovic, P. Chakravarty, and M. Pavone, “Trajectron++: Dynamically-feasible trajectory forecasting with heterogeneous data,” in _European Conference on Computer Vision_. Springer, 2020, pp. 683–700.
  * [14] Y. Ban, X. Li, G. Rosman, I. Gilitschenski, O. Meireles, S. Karaman, and D. Rus, “A deep concept graph network for interaction-aware trajectory prediction,” in _2022 International Conference on Robotics and Automation (ICRA)_. IEEE, 2022, pp. 8992–8998.
  * [15] R. D. Luce, “Individual choice behavior, a theoretical analysis,” _Bull. Amer. Math. Soc_ , vol. 66, no. 1960, pp. 259–260, 1960.
  * [16] H. A. Simon, “Theories of bounded rationality,” _Decision and organization_ , vol. 1, no. 1, pp. 161–176, 1972.
  * [17] B. D. Ziebart, A. L. Maas, J. A. Bagnell, A. K. Dey, _et al._ , “Maximum entropy inverse reinforcement learning.” in _Aaai_ , vol. 8. Chicago, IL, USA, 2008, pp. 1433–1438.
  * [18] B. P. Evans and M. Prokopenko, “A maximum entropy model of bounded rational decision-making with prior beliefs and market feedback,” _Entropy_ , vol. 23, no. 6, p. 669, 2021.
  * [19] A. De Martino and D. De Martino, “An introduction to the maximum entropy approach and its application to inference problems in biology,” _Heliyon_ , vol. 4, no. 4, p. e00596, 2018.
  * [20] D. Qiu, Y. Zhao, and C. L. Baker, “Latent belief space motion planning under cost, dynamics, and intent uncertainty,” in _Proceedings of Robotics: Science and Systems (RSS)_ , 2020.
  * [21] H. Hu and J. F. Fisac, “Active uncertainty learning for human-robot interaction: An implicit dual control approach,” _arXiv preprint arXiv:2202.07720_ , 2022. [Online]. Available: https://arxiv.org/abs/2202.0772
  * [22] H. Hu, K. Nakamura, and J. F. Fisac, “Sharp: Shielding-aware robust planning for safe and efficient human-robot interaction,” _IEEE Robotics and Automation Letters_ , vol. 7, no. 2, pp. 5591–5598, 2022.
  * [23] R. Emery-Montemerlo, G. Gordon, J. Schneider, and S. Thrun, “Approximate Solutions For Partially Observable Stochastic Games with Common Payoffs,” p. 8.
  * [24] E. A. Hansen, D. S. Bernstein, and S. Zilberstein, “Dynamic Programming for Partially Observable Stochastic Games,” p. 6.
  * [25] K. Horák and B. Bošanský, “Solving Partially Observable Stochastic Games with Public Observations,” vol. 33, pp. 2029–2036. [Online]. Available: https://www.aaai.org/ojs/index.php/AAAI/article/view/4032
  * [26] A. Kumar and S. Zilberstein, “Dynamic Programming Approximations for Partially Observable Stochastic Games,” p. 6.
  * [27] T. Kavuncu, A. Yaraneri, and N. Mehr, “Potential ilqr: A potential-minimizing controller for planning multi-agent interactive trajectories,” _arXiv preprint arXiv:2107.04926_ , 2021. [Online]. Available: https://arxiv.org/abs/2107.04926
  * [28] E. T. Jaynes, “On the rationale of maximum-entropy methods,” _Proceedings of the IEEE_ , vol. 70, no. 9, pp. 939–952, 1982.
  * [29] O. So, Z. Wang, and E. A. Theodorou, “Maximum entropy differential dynamic programming,” _arXiv preprint arXiv:2110.06451_ , 2021. [Online]. Available: https://arxiv.org/abs/2110.06451
  * [30] W. Li and E. Todorov, “Iterative linear quadratic regulator design for nonlinear biological movement systems.” in _ICINCO (1)_. Citeseer, 2004, pp. 222–229.
  * [31] D. Fridovich-Keil, V. Rubies-Royo, and C. J. Tomlin, “An iterative quadratic method for general-sum differential games with feedback linearizable dynamics,” in _2020 IEEE International Conference on Robotics and Automation (ICRA)_. IEEE, 2020, pp. 2216–2222.
  * [32] M. Maschler, S. Zamir, and E. Solan, _Game theory_. Cambridge University Press, 2013.
  * [33] M. J. Kochenderfer, _Decision making under uncertainty: theory and application_. MIT press, 2015.
  * [34] M. Hauskrecht, “Value-function approximations for partially observable markov decision processes,” _Journal of artificial intelligence research_ , vol. 13, pp. 33–94, 2000.
  * [35] M. L. Littman, A. R. Cassandra, and L. P. Kaelbling, “Learning policies for partially observable environments: Scaling up,” in _Machine Learning Proceedings 1995_. Elsevier, 1995, pp. 362–370.
  * [36] K. Hauser, “Online planning in continuous pomdps with open-loop information-gathering plans,” in _Proc. of the Int. Conf. on Machine Learning (ICML)_ , 2011.
  * [37] D. Bernardini and A. Bemporad, “Stabilizing Model Predictive Control of Stochastic Constrained Linear Systems,” vol. 57, no. 6, pp. 1468–1480. [Online]. Available: http://ieeexplore.ieee.org/document/6082380/
  * [38] A. Mesbah, “Stochastic model predictive control with active uncertainty learning: A Survey on dual control,” vol. 45, pp. 107–117. [Online]. Available: https://linkinghub.elsevier.com/retrieve/pii/S1367578817301232
  * [39] E. L. Zhu and F. Borrelli, “A sequential quadratic programming approach to the solution of open-loop generalized nash equilibria,” _arXiv preprint arXiv:2203.16478_ , 2022. [Online]. Available: https://arxiv.org/abs/2203.16478
  * [40] M. Bhatt, A. Yaraneri, and N. Mehr, “Efficient constrained multi-agent interactive planning using constrained dynamic potential games,” _arXiv preprint arXiv:2206.08963_ , 2022. [Online]. Available: https://arxiv.org/abs/2206.08963
  * [41] J. A. Andersson, J. Gillis, G. Horn, J. B. Rawlings, and M. Diehl, “Casadi: a software framework for nonlinear optimization and optimal control,” _Mathematical Programming Computation_ , vol. 11, no. 1, pp. 1–36, 2019.
  * [42] R. Verschueren, G. Frison, D. Kouzoupis, J. Frey, N. v. Duijkeren, A. Zanelli, B. Novoselnik, T. Albin, R. Quirynen, and M. Diehl, “acados—a modular open-source framework for fast embedded optimal control,” _Mathematical Programming Computation_ , vol. 14, no. 1, pp. 147–183, 2022.
  * [43] V. Dimitrov, P. Drews, T. Balch, X. Cui, A. A. Allaban, G. Rosman, and S. McGill. (2022, April) Rapid iterative design and testing with trikart. [Online]. Available: https://medium.com/toyotaresearch/rapid-iterative-design-and-testing-with-trikart-f5a3f59f2dd9
  * [44] M. O’Kelly, H. Zheng, D. Karthik, and R. Mangharam, “F1tenth: An open-source evaluation environment for continuous control and reinforcement learning,” _Proceedings of Machine Learning Research_ , vol. 123, 2020.
  * [45] Optitrack systems. https://optitrack.com. (accessed on 4 September 2022).
  * [46] S. L. Cleac’h, M. Schwager, and Z. Manchester, “Lucidgames: Online unscented inverse dynamic games for adaptive trajectory prediction and planning,” _arXiv preprint arXiv:2011.08152_ , 2020.


# Consensus Maximisation Using Influences of Monotone Boolean Functions

Ruwan Tennakoon1 David Suter2 Erchuan Zhang2 Tat-Jun Chin3 and Alireza Bab-
Hadiashar1
1RMIT University  Melbourne Australia. 2Edith Cowen University  Perth
Australia.
3University of Adelaide  Adelaide Australia.
{ruwan.tennakoon<EMAIL_ADDRESS>{d.suter<EMAIL_ADDRESS>
<EMAIL_ADDRESS>
###### Abstract

Consensus maximisation (MaxCon), which is widely used for robust fitting in
computer vision, aims to find the largest subset of data that fits the model
within some tolerance level. In this paper, we outline the connection between
MaxCon problem and the abstract problem of finding the maximum upper zero of a
Monotone Boolean Function (MBF) defined over the Boolean Cube. Then, we link
the concept of influences (in a MBF) to the concept of outlier (in MaxCon) and
show that influences of points belonging to the largest structure in data
would generally be smaller under certain conditions. Based on this
observation, we present an iterative algorithm to perform consensus
maximisation. Results for both synthetic and real visual data experiments show
that the MBF based algorithm is capable of generating a near optimal solution
relatively quickly. This is particularly important where there are large
number of outliers (gross or pseudo) in the observed data.

##  1 Introduction

The popular Maximum Consensus (MaxCon) criterion for robust fitting (as
typified by that of RANSAC [10]), seeks the maximum sized feasible set. Here
feasible means that all data points belonging to a “structure” (the inlier
set) fits its model within a tolerance level.

Given a set of $n$ data points
$\mathcal{D}=\left\\{\mathbf{p}_{i}\right\\}_{i=1}^{n}$ and a tolerance level
$\epsilon$, the MaxCon criterion for robust fitting can be written as:

$\displaystyle\underset{\theta,~{}\mathcal{I}\subseteq\mathcal{D}}{\max}\left|\mathcal{I}\right|$
$\displaystyle\textrm{subject
to}~{}~{}r_{\mathbf{p}_{i}}\left(\theta\right)\leq\epsilon\quad\forall\mathbf{p}_{i}\in\mathcal{I}$

where $r_{\mathbf{p}_{i}}\left(\theta\right)$ is the distance of
$\mathbf{p}_{i}$ from the model $\theta$.

A subset $\mathcal{I}$ can be represented by length-$n$ bit-vector,
$\mathbf{x}$, where the $i$’th position of the bit vector denote the inclusion
($x_{i}=1$) or exclusion ($x_{i}=0$) of the data point $i$. The above shows
that, each subset can be represented by a vertex of the $n$-dimensional
Boolean Cube. Therefore, any statement about which of the subsets are
feasible, is a statement on the evaluation of a Boolean function over the
$n$-dimensional Boolean Cube. This Boolean function
$f:\left\\{0,1\right\\}^{n}\rightarrow\left\\{0,1\right\\}$, outputs 1 (for
infeasible), or 0 (for feasible) for any vertex of the Boolean Cube. Since
MaxCon is a search for the maximum sized feasible subset of the data, it is
inherently a search on the Boolean Cube.

Such a view immediately opens the huge theory and very many associated
mathematical tools, developed during the study of Boolean functions; for the
purposes of the analysis of the MaxCon problem and for devising algorithms to
solve this problem.

Also, there is an additional significant observation. The Boolean function
associated with a MaxCon problem belongs to a special class of Boolean
functions called the Monotone Boolean Functions (MBF) [14].

######  Definition 1.1.

A MBF of $n$ variables is a mapping
$f:\left\\{0,1\right\\}^{n}\rightarrow\left\\{0,1\right\\}$ such that
$\alpha\prec\beta$ implies that $f\left(\alpha\right)\leq
f\left(\beta\right)$. Here we apply the natural ordering relation on
$n$-dimensional Boolean Cube: $\beta$ follow $\alpha$ ($\alpha\prec\beta$) if
for any $i$, the equality $\alpha_{i}\leqslant\beta_{i}$ is satisfied.

This is easy to see: If a subset is feasible, then adding more data to that
subset can only move the function towards infeasibility and once infeasible
adding more points will not change the subset back to feasible. Likewise,
deleting points from a subset can move only towards feasible. MaxCon, when
viewed in the above sketched framework, is nothing more than the search for
the maximum upper zero (see the definition below) of the above-mentioned
monotone infeasibility function.

######  Definition 1.2 (Upper zero of a MBF).

Upper zero of the MBF $f\left(\cdot\right)$ is a vertex $\alpha$ for which
$f\left(\alpha\right)=0$ and, for all $\alpha\prec\beta$ the relation
$f\left(\beta\right)=1$ is satisfied.

######  Definition 1.3 (Maximum Upper zero of a MBF [15]).

Maximum Upper zero of the MBF $f\left(\cdot\right)$ is a vertex $\alpha$ for
which $f\left(\alpha\right)=0$ and, for all
$\left\|\beta\right\|_{1}>\left\|\alpha\right\|_{1}$ the relation
$f\left(\beta\right)=1$ is satisfied. Here,
$\left\|\alpha\right\|_{1}=\sum_{i}\alpha_{i}$.

A simple example in Figure 1 illustrates the concepts involved. The figure
shows that, along paths going up the Hasse diagram (a directed version of the
Boolean cube), the function can only stay constant or increase (never
decrease). The vertex “11110” (representing the subset
$\left\\{\mathbf{p}_{1},\mathbf{p}_{2},\mathbf{p}_{3},\mathbf{p}_{4}\right\\}$)
is the MaxCon solution for this problem as it is the highest feasible subset
in the Hasse diagram (maximum upper zero). The vertex “11001” is an example of
an upper zero as there are no feasible nodes that follows it
($\succ\textrm{`11001'}$) on the Hasse diagram (can be seen as a local
optimum).

(a) 2D line fitting problem

(b) Hasse Diagram

Figure 1: (a) An example of a 2D line fitting problem with 5 data points
together with (b) the associated Boolean Cube and MBF - represented in Hasse
diagram format. The 5D Boolean Cube is drawn flattened onto 2D and oriented so
that the higher up a vertex appears, the larger is the number of 1’s in the
coordinates. A Boolean Function maps the Boolean cube to 0 or 1. We illustrate
by colouring: “white” nodes map to 1, coloured nodes map to 0. This example is
a Monotone Boolean Function moving up the picture, the value of the function
only ever increases, it never translate in the opposite direction. Red node is
the maximum upper zero in this example.

Whilst Monotone Boolean Functions (MBFs) have been extensively studied for a
variety of application domains, including learning theory [1] (and computer
vision is, these days, highly dominated by learning style approaches), there
appears to be relatively little attention to MBFs in computer vision: and more
specifically in model-based computer vision. Reference [25] appears to be a
recent exception - but even this is tackling very different considerations
from our main areas of interest (namely, [25] is concerned with efficient
Conditional Random Field (CRF) calculations, and CRF modelling is very
different to the geometric modelling we refer to).

In this paper, we present a novel view point on the MaxCon problem using the
Monotone Boolean Functions theory and use that to develop some efficient
solutions. For this purpose, we concentrate on a property of all Boolean
Functions, called Influence [21] (definition in Section 2.1); a property that
has a particular relationship with the Fourier Transform of a Boolean Function
when that Boolean Function is Monotone. In summary our main contributions
include:

  1. 1.

We provide a precise (and abstract) definition of the MaxCon problem in terms
of finding the maximum upper zero of a Monotone Boolean Function defined over
the Boolean Cube where, vertices correspond to subsets of data and the output
of the MBF correspond to whether the subset can be feasibly covered by a model
with set tolerance.

  2. 2.

We link the concept of influences in a MBF to the concept of outlier in MaxCon
and provide theoretical analysis to show that influences of points belonging
to the largest structure in data would be smallest under “ideal” conditions
(defined in Section 3).

  3. 3.

Based on the above analysis, we derive a greedy algorithm that search for the
maximum upper zero of a MBF to efficiently solve the MaxCon problem.

##  2 Background

###  2.1 Boolean Functions and Influences

For any Boolean function
$f:\left\\{0,1\right\\}^{n}\rightarrow\left\\{0,1\right\\}$, the probability
that the $i$’th coordinate in a random length-$n$ binary vector, $\mathbf{x}$,
“affects” the output of the function is defined as the influence [21]:

$\mathrm{Inf}_{i}\left[f\right]=\underset{\mathbf{x}\sim\left\\{0,1\right\\}^{n}}{\mathrm{\mathbf{Pr}}}\left[f\left(\mathbf{x}\right)\neq
f\left(\mathbf{x}^{\oplus i}\right)\right]$ (1)

where $\mathbf{x}^{\oplus i}$ is equal to $\mathbf{x}$ with the $i$’th
coordinate flipped. For a MBF, the influence of the $i$’th coordinate is equal
to the degree-1 Fourier coefficient [21].

The most obvious strategy for influence approximation would be to uniformly
sample the cube.

$\mathrm{Inf}_{i}\left[f\right]=\underset{\mathbf{x}\sim\left\\{0,1\right\\}^{n}}{\mathbb{E}}\delta\left[f\left(\mathbf{x}\right)\neq
f\left(\mathbf{x}^{\oplus i}\right)\right]$ (2)

where $\delta\left[\cdot\right]=1$ if the condition inside is true (0
otherwise). However, in Section 4.1, we empirically show that this may not be
the most efficient when it comes to solving MaxCon. A better approach would be
to use the notion of weighted influence [8]. In such an approach one typically
defines a Bernoulli measure $\mu_{q}(\mathbf{x})$ over the vertices of the
Boolean Cube. Operationally, this can correspond to sampling by selecting to
“turn bits $i$” independently and with probability $q$. Uniform sampling
corresponds to $q=0.5$. Sampling with low $q$ concentrates the measure towards
the bottom of the cube (small Hadamard norm) and high $q$ towards the top. The
influences under $\mu_{q}$ distribution can now be defined as:

$\mathrm{Inf}_{i}^{(q)}\left[f\right]=\underset{\mathbf{x}\sim\mu_{q}(\mathbf{x})}{\mathbb{E}}\delta\left[f\left(\mathbf{x}\right)\neq
f\left(\mathbf{x}^{\oplus i}\right)\right].$ (3)

The influences under different $q$ values are analysed further in Section 4.1.

###  2.2 Determining feasibility and the concept of a basis

Searching for the upper zeros of a MBF will require the evaluation of the
function on a given subset. In the case of MaxCon, this involves finding if a
given subset is feasible or infeasible. The feasibility/infeasibility of a
subset can be obtained efficiently via solving the following minmax problem
and checking if the resulting objective value is within $\epsilon$ [13, 5]:

$\underset{\theta}{\min}~{}\underset{i}{\max}~{}r_{\mathbf{p}_{i}}\left(\theta\right)$
(4)

The above can be solved exactly for cases where
$r_{\mathbf{p}_{i}}\left(\theta\right)$ is quasiconvex using algorithms based
on bisection [9]. The solution of equation (4) for a subset
$\mathcal{I}^{(j)}$, will return three quantities: the optimal objective value
$g^{(j)}$, basis $\mathcal{B}^{(j)}$ and the fitted model $\theta^{(j)}$. A
basis for a set, is a small subset such that the value of the final objective
on that subset is equal to the final objective on the whole set [7].
Intuitively, a basis in the MaxCon setting is a subset of points “that
prevents the enclosing structure from shrinking further”. The cardinality of
the basis in MaxCon is connected to the concept of “combinatorial dimension”
and we denote it by $p+1$ ($p$ is the minimum number of data points required
for the unique determination of the model).

###  2.3 Related work in computer vision

Consensus maximisation is widely used for robust fitting in computer vision.
Early work mostly focused on solving the consensus maximisation using
randomised hypothesise-and-test techniques [10, 17, 24]. More recently, some
focus has shifted towards finding globally optimal solutions [5, 31, 2, 18, 4,
19]. As MaxCon is known to be NP-hard, the run time of the global methods
scale exponentially in the general case [6]. This has lead to the search of
deterministic and/or near optimal algorithms [3, 23, 16, 30].

To the best of our knowledge, ours is the first work in computer vision that
explicitly poses the MaxCon problem in terms of finding the maximum upper zero
of a MBF. However, there are several works that are related, in some sense, to
the proposed framework. Though not expressed that way in the original works,
the $\textrm{A}^{*}$ “tree searches” of [2] and [5] are in fact searches on
this cube (where nodes reached by different paths become repeated nodes in the
tree constructed by starting from the top of the Hasse Diagram). In this
context, it is interesting to note that RANSAC (and it’s many derivatives)
search amongst subsets towards the bottom of the Hasse Diagram (minimal sized
subsets or $p$-sized); and use these to “index” up to higher subsets by
greedily including all data points the fit within tolerance of the model found
on the $p$-subset. Furthermore, the methods that use the
$\textrm{L}_{\infty}$-Norm (equation 4) for outlier removal in the context of
geometric fitting, can also be seen as traversing the Boolean cube [26, 22,
19]. However, the above methods do not explicitly utilize the theory and tools
associated with MBF.

##  3 Proposed Method

As described in the introduction, MaxCon can be seen as finding the maximum
upper zero of a MBF. This paper is devoted to the exploitation of information
embedded in the Influence function to find the maximum upper zero,
efficiently. To explain the method, we first note that the influence of an
inlier data point is likely to be smaller than the influence of an outlier
data point. An example of this relationship for a 2D line-fitting problem is
shown in Figure 2. The essential and intuitive reason is that inclusion of an
outlier, into a feasible subset, most likely turns that subset infeasible
(thereby “influencing” the function). In contrast adding an inlier leaves the
set feasible. It is the “breaking” or “creation” of an infeasible $p+1$ sized
basis that is responsible for flipping the outcome of the MBF and therefore
the addition/deletion of an outlier triggers this more than that of an inlier.
Below, we prove that the above property of influences hold under “ideal”
conditions. The proofs of theorems in the following text are available in the
supplementary materials.

Figure 2: An example 2D line-fitting problem and the corresponding influences.
The data generation procedure is described in detail in Section 4.1.
Influences are estimated using equation (3) with $m=1000$ and $q=0.5$. Data
points are sorted such that inliers comes before outliers.

###  3.1 Ideal single structure case

We first formally define an ideal structure in data. Let
$L_{k}:=\left\\{\mathbf{x}\in\left\\{0,1\right\\}^{n}:\left\|\mathbf{x}\right\|_{1}=k\right\\}$
be the level $k$ in the $n$-dimensional Boolean cube and, $L_{\leq k}$ the
levels below $k+1$.

######  Definition 3.1.

Given a monotone Boolean function $f$, for $\mathbf{x}^{k}\in L_{k}$ ($p<k\leq
n$), $f$ is called ideal with respect to $\mathbf{x}^{k}$ if

$f\left(\mathbf{x}\right)=\begin{cases}0&\forall\mathbf{x}\in
B_{\mathbf{x}^{k}}\cup L_{\leq p}\\\ 1&\text{others}\end{cases}$ (5)

where
$B_{\mathbf{x}^{k}}=\left\\{\mathbf{x}\in\left\\{0,1\right\\}^{n}:d\left(\mathbf{x},\mathbf{x}^{k}\right)=l,\mathbf{x}\in
L_{k-l}\right\\}$ for all $0\leq l\leq k-p-1$, is the Boolean sub-cube
determined by $\mathbf{x}^{k}$. Here, $d\left(\cdot,\cdot\right)$ is the
hamming distance.

Intuitively, an ideal structure has only one structure: points are either
inlier or outliers to that single structure, and no other structure exists
other than the trivial with $p$ or less points.

######  Theorem 3.1.

If a monotone Boolean function $f$ is ideal with respect to $\mathbf{x}^{k}\in
L_{k}$, then

$\mathrm{Inf}_{i}\left[f\right]=\begin{cases}\frac{C^{n-1}_{p}-C^{k-1}_{p}}{2^{n}}&\text{if}~{}i~{}\text{is
inlier}\\\
\frac{C^{n-1}_{p}+\sum_{l=p+1}^{k}C^{k}_{l}}{2^{n}}&\text{if}~{}i~{}\text{is
outlier}\end{cases}$ (6)

######  Corollary 3.1.1.

$\mathrm{Inf}_{j}\left[f\right]-\mathrm{Inf}_{i}\left[f\right]=\frac{1}{2^{n}}\sum_{l=p+1}^{k}C_{l}^{k}+C_{p}^{k-1}>0$
(7)

where $i$ is any inlier and $j$ is any outlier.

The above shows that for an ideal structure case the influence of an inlier is
strictly smaller than the influence of an outlier.

###  3.2 Ideal $K$-structure case

In what follows, we generalize the above result to ideal $K$-structure ($K\geq
1$), namely there are $K$ upper zeros only and the Boolean sub-cubes
determined by these upper zeros are disjoint above level $p$.

######  Definition 3.2.

Let $f$ be a monotone Boolean function and
$\left\\{\mathbf{x}^{k_{r}}\right\\}_{r=1}^{K}$ are upper zeros, where
$p<k_{1}\leq k_{2}\leq\cdots\leq k_{K}\leq n$, then $f$ is called $K$-ideal
with respect to $\left\\{\mathbf{x}^{k_{r}}\right\\}_{r=1}^{K}$ if

$\begin{matrix}1)&d\left(B_{\mathbf{x}^{k_{i}}}\setminus L_{\leq
p},B_{\mathbf{x}^{k_{j}}}\setminus L_{\leq p}\right)>0\quad\forall k_{i}\neq
k_{j}&\\\
2)&f\left(\mathbf{x}\right)=\begin{cases}0&\forall\mathbf{x}\in\bigcup_{r=1}^{K}B_{\mathbf{x}^{k_{r}}}\cup
L_{\leq p}\\\ 1&\text{others}\end{cases}\end{matrix}$ (8)

The first condition states that the sub-cubes determined by
$\mathbf{x}^{k_{r}}$ have no elements in common. Let
$\mathcal{S}_{\mathbf{x}^{k_{r}}}^{c}$ be the set of inlier (if $c=1$) or the
set of outliers (if $c=0$) to $\mathbf{x}^{k_{r}}$. Then we can define

$\mathcal{S}_{c_{1}c_{2}\cdots
c_{K}}=\bigcap_{r=1}^{K}\mathcal{S}_{\mathbf{x}^{k_{r}}}^{c_{r}}\quad
c_{r}\in\left\\{0,1\right\\}$ (9)

which represent the index set of inlier to structures where bit string
$c_{1}c_{2}\cdots c_{K}$ is one. For example, $\mathcal{S}_{11\cdots 1}$ is
the index set of points that are inliers with respect to all
$\mathbf{x}^{k_{r}}$, and $\mathcal{S}_{00\cdots 0}$ is the index set of
points that are outliers with respect to all $\mathbf{x}^{k_{r}}$.

######  Theorem 3.2.

$2^{n}\mathrm{Inf}_{\mathcal{S}_{c_{1}\cdots
c_{K}}}\left[f\right]=C_{p}^{n-1}+\sum_{\underset{1\leq r\leq
K}{c_{r}=0}}\sum_{l=p+1}^{k_{r}}C_{l}^{k_{r}}-\sum_{\underset{1\leq r\leq
K}{c_{r}=1}}C_{p}^{k_{r}-1}$ (10)

where $\mathrm{Inf}_{\mathcal{S}_{c_{1}\cdots c_{K}}}\left[f\right]$ denote
the influence $\mathrm{Inf}_{i}\left[f\right]$ of
$i\in\mathcal{S}_{c_{1}\cdots c_{K}}$.

######  Corollary 3.2.1.

The influences have the following ordered relationship

$\forall\alpha,\beta\in\left\\{0,1\right\\}^{K},\alpha>\beta\Rightarrow\mathrm{Inf}_{\mathcal{S}_{\alpha}}\left[f\right]<\mathrm{Inf}_{\mathcal{S}_{\beta}}\left[f\right]$
(11)

which means $\mathrm{Inf}_{x}\left[f\right]$ is a real-valued monotone
decreasing Boolean function. If any $\mathcal{S}_{\bullet}=\emptyset$ then
$\mathrm{Inf}_{\mathcal{S}_{\bullet}}\left[f\right]$ is not defined.

Following Corollary 11 and the Definition 8 we can see that the influence of a
point that belongs to a structure with a higher upper zero is smaller that
that of a point belongs to a structure with a lower upper zero.

###  3.3 Finding the maximum upper zero

The above theoretical analysis shows that influences of points belonging to
the largest structure in data would be smallest under “ideal” conditions.
Thus, in the ideal case, a very simple process of estimation of the
influences, followed by thresholding, would immediately lead to the sought
after MaxCon solution. However, we need to recognise that:

  * •

The non-ideal case (which applies to all realistic data sets) will be a
perturbation away from this, where the influences of inlier and outlier come
closer together.

  * •

We can only work with estimated influences, not the actual influences. Thus we
need to assume that the estimated influences largely follow the ordering given
mathematically in our derivations and according with the above intuition.

Nonetheless, continuity of behaviour arguments suggest the departure from
ideal will be only partial for many data sets; and that relatively standard
ways for attacking such scenarios, should still be viable. We empirically
verify that this is often true in our experiments.

Using the above intuition (“influences of points belonging to the largest
structure in data would be smallest”) we formulate an algorithm - MBF-MaxCon
(Algorithm 1) that starts with the set of all data points and then gradually
removes one data point at a time (data point with the largest influence) until
the remaining set of points is within the $\epsilon$ band. In the algorithm,
we also use the notion that: If there are outliers in a subset of data, then
at least one of them should belong to the basis returned by solving equation
(4) [26]. This additional constraint enables us to compute only $p+1$
influences per elimination of a data point (rather than for all points in
subset) which leads to a more efficient algorithm.

The estimation of Influence introduces some noise to the proposed algorithms,
and the solution returned by them may not include all the inlier points of a
given structure. To partially overcome this, we introduce a local expansion
step (Algorithm 2). In this step, starting from the initial solution, at each
iteration, the distance-1 upper neighborhood (Hasse Diagram) of the current
solution is explored and the current solution is updated if there is any
feasible set. This process is repeated until there are no feasible subsets in
the distance-1 upper neighborhood. This local update guarantees that the
algorithm will find an upper zero (local optimal).

Algorithm 1 MBF-MaxCon - algorithm for finding the maximum consensus set using
influences of BMFs.

0: $\left\\{\mathbf{p}_{i}\right\\}_{i=1}^{n}$, $\epsilon$, $m$, $q$.

1: $\mathbf{x}\leftarrow\left[1,\dots,1\right]_{\left[1\times n\right]}$

2: repeat

3: $\mathcal{I}^{(t)}\leftarrow\left\\{i:x_{i}=1\right\\}$

4: Solve equation (4) for subset $\mathcal{I}^{(t)}$ and obtain
$\mathcal{B}^{(t)}$

5: Estimate $\mathrm{Inf}_{i}^{(q)}\left[f\right]\quad\forall
i\in\mathcal{B}^{(t)}$ using equation (3)111When estimating influences, one
can use the monotonic nature of $f$ to save some computations (In a MBF, the
function value before the bit flip has some information regarding the value
after). The algorithm used for estimating the influences is available in
supplementary materials..

6:
$e\leftarrow\underset{i\in\mathcal{B}^{(t)}}{\textrm{argmax}}~{}\mathrm{Inf}_{i}^{(q)}\left[f\right]$

7: $x_{e}\leftarrow 0$

8: until $f(\mathbf{x})=0$

9: $\mathbf{x}\leftarrow$ Run local expansion step in (Algorithm 2)

10: return $\mathcal{I}\leftarrow\left\\{i:x_{i}=1\right\\}$

Algorithm 2 Local expansion step.

0: $\left\\{\mathbf{p}_{i}\right\\}_{i=1}^{n}$, $m$, initial feasible set
$\mathbf{x}$.

1: repeat

2: updated $\leftarrow$ false

3: for all $i:x_{i}=0$ do

4: $\bar{\mathbf{x}}\leftarrow\mathbf{x}$; $\bar{{x}}_{i}\leftarrow 1$

5: if$~{}f(\bar{\mathbf{x}})=0$ then$~{}\mathbf{x}\leftarrow\bar{\mathbf{x}}$;
updated $\leftarrow$ true; break;

6: end for

7: until updated=false

8: return $\mathbf{x}$

It is important to note that one needs to re-estimate the influences at each
iteration, $t$ of algorithm 1. If $\mathrm{Inf}_{i}^{[t]}\left[f\right]$ is
the influence of point $i$ at iteration $t$ of the algorithm, then
$\mathrm{Inf}_{i}^{[t]}\left[f\right]\neq\mathrm{Inf}_{i}^{[t-1]}\left[f\right]$.
This is because function at level $t$ is a restricted version of the function
at level $t-1$ and the relationship between influences at different levels is
derived in supplementary materials. Another allied intuition is that though,
as mentioned before, noise (from both the estimation process and from the
departure of the data from that of being ideal) will raise the level of the
influence of some inliers (and decrease the values of some outliers) to the
point where the estimated influences of some inliers will be above those of
some outliers: at each stage we only remove the largest influence data point
which will be away from the “polluted” data division; and that re-estimation
afterwards allows the possibility for the re-estimated influences to be
“cleaner”.

##  4 Results

We evaluated the performance of the proposed algorithms, on both synthetic and
real data experiments, and compared those with the state-of-the-art
techniques. All experiments were executed in MATLAB on a computer with Intel
Xeon E5-1650 CPU, 16GB RAM and Ubuntu 16.04 OS. The publicly available codes
were used to obtain the results for improved $\textrm{A}^{*}$ tree
search222Code: MaxConTreeSearch ($\textrm{A}^{\ast}$-NAPA-DIBP) [2] and Lo-
RANSAC333Code: Consensus-maximization-with-biconvex-programming [17]. Our
implementation is publicly available at Link removed to preserve anonymity.

###  4.1 Influence estimation

Estimating the influences via equation (3) requires two hyper-parameters: The
number of randomly sampled bit-vectors $m$ and, the sampling probability $q$
of the Bernoulli measure $\mu_{q}(x)$. In this section we explore the effects
of the two hyper-parameters on influence estimation using a 2D-line fitting
problem. Here, the number of data points ($n$) is set to $15$, out of which
25% would be outliers ($n_{o}$). The number of points is chosen to be
relatively small as the exact influence calculation time increases
exponentially with $n$. A subset of ($n-n_{o}$) randomly selected points
(inliers) were then perturbed with uniformly distributed noise in the range
$\left[-0.1,0.1\right]$. The remaining $n_{o}$ data points (outliers) were
then perturbed with uniformly distributed noise from
$\left[-5,-0.1\right)\cup\left(0.1,5\right]$. The inlier threshold $\epsilon$
was set to $0.1$ for all the experiments in this section.

First, we compare the estimated influences with the exact influences for
different combinations of $m$ and $q$. The exact influences,
$\overline{\mathrm{Inf}}_{i}\left[f\right]$, are computed by taking the
expectation in equation (2) over all $2^{n}$ vertices of the Boolean cube. The
estimation error of all influence values can then be calculated as:
$\frac{1}{n}\sum_{i\in\left[n\right]}\left({\mathrm{Inf}}_{i}^{(q)}\left[f\right]-\overline{\mathrm{Inf}}_{i}\left[f\right]\right)^{2}$.

The results in Figure 3 shows that, for all $q$ values, the influence
estimation error is high when only a few samples are used. However, the error
decreases exponentially with increasing the number of samples. We can also see
that the estimated influences are closest to their exact values when $q=0.5$.
When $q$ is varied in either direction, the error increases. This is because
at $q=0.5$, the computed values are an unbiased estimate of the influences.
Changing the sampling distribution changes the definition of the orthonormal
basis and hence introduces a bias [8].

However, in solving MaxCon we do not seek for an unbiased estimate of the
influences. What we are after is a definition of influences where the
separation between influences of inliers and ourliers are maximum. We used the
same synthetic experiment to analyse the separation. Here we define the
separation as:

$\textrm{Separation}=\underset{i\in\mathcal{D}_{out}}{\min}~{}\mathrm{Inf}_{i}^{(q)}\left[f\right]-\underset{j\in\mathcal{D}_{in}}{\min}~{}\mathrm{Inf}_{j}^{(q)}\left[f\right]$
(12)

where $\mathcal{D}_{in}$ is the set of inlier data points and
$\mathcal{D}_{out}=\mathcal{D}\setminus\mathcal{D}_{in}$. The results in
Figure 4 show that the separation increase with $m$. However, in relation to
$q$, the separation is maximum when the $q$ value is small. Using $q<(p+1)/n$
will also decrease the separation as most samples at this probability will be
trivially feasible. In summary, changing $q$ leads to different measures of
influence where the order (influences of inliers are smaller than outliers) is
much the same, though the separation can vary. In our experiments, we use a
$q$ value between $(p+1)/n$ and $0.3$.

Figure 3: The variation of the influence estimation error with the number of
samples used ($m$) and the sampling probability ($q$). The figure show the
statistics across 100 random runs for each combination.

Figure 4: The variation in the separation (between influences of inliers and
outliers) with the number of samples used ($m$) and the sampling probability
($q$). The figure show the statistics across 100 random runs for each
combination.

###  4.2 Controlled experiments with synthetic data

To study the behaviour of the proposed algorithm under a controlled setting,
similar to [2], we conducted experiments on an 8-dimensional robust linear
regression problem with synthetically generated data. First, a set of $n=200$
data points on a randomly instantiated model $\theta\in\mathbb{R}^{8}$ was
generated. The data set was then perturbed using the same process described in
Section 4.1 to get a data set that is corrupted by noise and outliers. In our
experiments, the number of outliers, $n_{o}$, were varied in the range of
$\left[5,40\right]$ (upper bound determined by computation time of
$\textrm{A}^{\ast}$-NAPA-DIBP used for obtaining ground truth
inliers/outliers). The error of a method is computed by comparing the
cardinality of the MaxCon solution found by that method
($\left|\mathcal{I}_{\bullet}\right|$) to the ground truth cardinality
obtained using $\textrm{A}^{\ast}$-NAPA-DIBP
($\left|\mathcal{I}_{\mathrm{A}^{\ast}}\right|$).

Ablation study: To identify the importance of each component in our overall
algorithm, we conduct an ablation study using the above data. Here, we analyse
four variants of our algorithm: 1) MBF-MaxCon-nR: Simple algorithm where all
the influences are computed at start and the data point with largest influence
is removed iteratively until the remaining subset is feasible (no re-
estimation of influences at each iteration). 2) MBF-MaxCon-nB: Same as MBF-
MaxCon-nR but the influences of all remaining points are recomputed (not just
points in basis) at the end of each iteration. 3) MBF-MaxCon-nL: Same as
algorithm 1 without the local expansion in line 9. 4) MBF-MaxCon: Proposed
algorithm 1.

The error of each variant and the computation times are shown in Fig. 5. The
results show that re-estimating influences at the end of each iteration has a
significant effect on the final outcome (see section 3.3 for explanation).
Furthermore, the results show that, on average, the use of both local
expansion and basis has helped to move the solutions closer to the global
optimal. The results also show that the most significant contribution in terms
of the use of basis is in computational efficiently.

(a)

(b)

Figure 5: Results of the ablation study for 8-dimensional robust linear
regression with synthetic data (a) Number of inliers found compared with the
global optimal (obtained using $A^{*}$) and (b) Variation of computational
time with number of outliers. The experiments were repeated 100 times and the
error-bars indicate the $0.05$-th and $0.95$-th percentile.

Comparative analysis: Next, we compare the performance of the proposed method
with the most relevant methods in literature: $\textrm{A}^{\ast}$-NAPA-DIBP
[2], RANSAC [10] and Lo-RANSAC [17]. Both RANSAC variants were run with the
number of RANSAC iterations set to match the computation time of MBF-MaxCon.

The computation time for the above methods are shown in Figure 6(b). The
figure shows that when the number of outliers are low ($<30$)
$\textrm{A}^{\ast}$-NAPA-DIBP converges to a solution relatively quickly.
However, the computational time of $\textrm{A}^{\ast}$-NAPA-DIBP increases
exponentially with the number of outliers. On the other hand, the
computational time of the proposed algorithms increase linearly with the
number of outliers. This is clearly predictable as our algorithms take one
step across each level and the deeper down is the MaxCon solution,
proportionally longer is the “search”. $\textrm{A}^{\ast}$-NAPA-DIBP has a
much more sophisticated search that allows backtracking of routes explored and
this causes the exponential behaviour when that is heavily exercised. Figure
6(a) shows the difference between the number of inliers returned by
$\textrm{A}^{\ast}$-NAPA-DIBP and other methods. On average the proposed
method MBF-MaxCon returns a solution with usually close to the same number of
inliers as the A* method. The figure also shows the 0.05th and 0.95th
percentile distances from $A^{*}$ solution over 100 random runs. This shows
that in few cases the solution returned by MBF-MaxCon can be up to 2-4 inliers
away from the optimal solution (around 1% error). The main summary is that the
proposed methods never “go exponential”, unlike $A^{*}$, in runtime, but
compares often favourably in terms of accuracy.

(a)

(b)

Figure 6: Results for 8 dimensional robust linear regression with synthetic
data (a) Number of inliers found compared with the global optimal (obtained
using $A^{*}$) and (b) Variation of computational time with number of
outliers. The experiments were repeated 100 times and the error-bars indicate
the $0.05$-th and $0.95$-th percentile.

###  4.3 Linearized fundamental matrix estimation

In this section we examine the performances of the proposed methods for
linearized fundamental matrix estimation. Provided that the point matches
between two views are given as $[\textbf{p}_{1},\textbf{p}_{2}]$ where
$\textbf{p}_{j}=(x_{j},y_{j},1)^{\top}$ is a coordinate of a point in view
$j$, each rigid motion in the scene can be modeled using the fundamental
matrix $F\in\mathcal{R}^{3\times 3}$ as [27]:
$\textbf{p}_{1}^{\top}F\textbf{p}_{2}=0$. In our experiments we use the
linearized version presented in [5] together with the algebraic error [12]
(chapter 11). For each image pair, the input is a set of SIFT [20] feature
matches generated using VLFeat [28].

Single dominant motion: Following [2], we used the first five crossroads image
pairs from the sequence “00” of the KITTI Odometry dataset [11] in our
experiments. The inlier threshold $\epsilon$ is set to $0.03$ for all image
pairs. The Number of inliers returned by each method ($n_{i}$) and the
computation times are shown in the first part of Table 1. The results reported
for the probabilistic methods are the mean (min, max) over 100 random runs.
The results show that the proposed methods on average have produced solutions
that are close to the optimum solution returned by $\textrm{A}^{\ast}$-NAPA-
DIBP [2].

The distribution of errors by each algorithm over 100 repeated runs for all
the frames in sequence “00” of the KITTI Odometry data set is shown in Figure
7. The main message is that we operate in a time cost regime a little better
than $\textrm{A}^{\ast}$-NAPA-DIBP and around the same as we allowed for Lo-
RANSAC but we generally get much closer to $\textrm{A}^{\ast}$-NAPA-DIBP
performance - including often finding the optimal, which RANSAC or Lo-RANSAC
rarely do in this experiment.

Table 1: Linearized fundamental matrix estimation result. “Same Comp.” refers to running RANSAC with the same time budget as MBF-MaxCon and, “sp=0.99” refers to running RANSAC with the success probability $0.99$.  |  |  | A*-NAPA
---
-DIBP
MBF-MaxCon |  | RANSAC
---
Same Comp.
| Lo-RANSAC
---
Same Comp.
| RANSAC
---
sp=0.99
| Lo-RANSAC
---
sp=0.99
104-108 | $n_{i}$ | 289 |  | 288.52
---
(289, 285)
| 282.03
---
(286, 277)
| 284.54
---
(287, 283)
| 271.18
---
(282, 254)
| 281.41
---
(284, 276)
Time (s) | 10.96 | 1.78 | 1.78 | 1.78 | 0.004 | 0.04
198-201 | $n_{i}$ | 296 |  | 293.10
---
(296, 291)
| 291.88
---
(294, 290)
| 292.87
---
(294, 291)
| 287.00
---
(293, 272)
| 290.35
---
(293, 287)
Time (s) | 4.04 | 2.05 | 2.05 | 2.05 | 0.004 | 0.04
417-420 | $n_{i}$ | 366 |  | 364.44
---
(366, 359)
| 363.54
---
(365, 361)
| 364.02
---
(365, 363)
| 357.38
---
(364, 343)
| 362.33
---
(364, 359)
Time (s) | 7.93 | 2.59 | 2.59 | 2.59 | 0.004 | 0.06
579-582 | $n_{i}$ | 523 |  | 520.68
---
(523, 514)
| 518.04
---
(521, 511)
| 520.89
---
(522, 520)
| 502.31
---
(519, 463)
| 517.30
---
(512, 497)
Time (s) | 4.28 | 3.22 | 3.22 | 3.22 | 0.004 | 0.11
738-742 | $n_{i}$ | 462 |  | 460.97
---
(462, 457)
| 455.02
---
(459, 451)
| 457.24
---
(460, 455)
| 438.55
---
(455, 409)
| 451.35
---
(459, 435)
Time (s) | 2.97 | 2.72 | 2.72 | 2.72 | 0.005 | 0.1
breadcube | $n_{i}$ | $\sim$ |  | 64.36
---
(68, 58)
| 61.77
---
(65, 59)
| 63.66
---
(66, 61)
| 59.04
---
(65, 55)
| 62.17
---
(66, 59)
Time (s) | >3600 | 19.07 | 19.07 | 19.07 | 0.993 | 1.07
breadtoy | $n_{i}$ | $\sim$ |  | 107.48
---
(115, 102)
| 105.24
---
(111, 102)
| 107.23
---
(111, 104)
| 100.31
---
(107, 93)
| 105.21
---
(109, 101)
Time (s) | >3600 | 38.37 | 38.37 | 38.37 | 0.526 | 0.68
cubetoy | $n_{i}$ | $\sim$ |  | 58.51
---
(61, 52)
| 54.54
---
(58, 52)
| 56.14
---
(58, 55)
| 52.20
---
(56, 48)
| 55.24
---
(57, 52)
Time (s) | >3600 | 15.38 | 15.38 | 15.38 | 0.45 | 0.65
Figure 7: The distribution of errors by each algorithm over 100 repeated runs
for all the frames in sequence “00” of the KITTI Odometry data set.

Multiple motions: The above data set contains a single dominant motion and the
number of outliers are limited (around 13-22). To study the behaviour of the
proposed algorithm in the presence of multiple structures, we conducted
experiments on three sequences from the AdelaideRMF data-set [29]. The inlier
threshold $\epsilon$ was set to $0.015$ for all image pairs. In these
sequences, there are two rigidly moving objects. The results in the second
part of Table 1 show that the method $\textrm{A}^{\ast}$-NAPA-DIBP did not
find a solution after $3600$s (1 hour), where as the proposed methods found
solutions in around 15-40 seconds. Once again, when compared to Lo-RANSAC, we
generally obtain higher consensus sets (and thus implicitly closer to what
$A^{*}$ is capable of, but on these data sets, would require astronomically
more computation.

##  5 Conclusion

We have applied a new perspective to the long standing problem of MaxCon. This
perspective recognises that the underlying mathematical object is a Monotone
Boolean Infeasibility function, defined over the Boolean Cube. Such a
perspective immediately identifies a rich mathematical theory that can be
applied. Very probably, we have only scratched that surface here. But we have
been able to take at least one element of the theory (i.e. Influence) and link
that concept to the concept of outlier (in MaxCon) and shown that already,
without borrowing further from the rich theory, that we can derive algorithms
that are already at least competitive, in some aspects. Specifically:

  1. 1.

The approach sometimes achieves the true MaxCon, whereas RANSAC rarely
achieves the MaxCon solution. Indeed, it is well known that this is a feature
of RANSAC based methods in general.

  2. 2.

The approach can (mostly) achieve close to $A^{*}$ (provably optimal) in a
similar time budget or faster.

  3. 3.

Unlike $A^{*}$ the approach will never go exponential in runtime or memory
requirements, and - as above - unlike RANSAC variants, it can more reliably
obtain the optimal or close to optimal result - given a similar time budget.

The proposed algorithm can be used on model fitting problems where the
residuals, $r_{\mathbf{p}_{1}}\left(\theta\right)$, are quasi-convex and the
run-time of the algorithm would depend on the existence of the an efficient
oracle to evaluate the feasibility/infeasibility of a subset of points.

## References

  * [1] Nader H Bshouty and Christino Tamon.  On the fourier spectrum of monotone functions.  Journal of the ACM (JACM), 43(4):747–770, 1996.
  * [2] Z. Cai, T. Chin, and V. Koltun.  Consensus maximization tree search revisited.  In 2019 IEEE/CVF International Conference on Computer Vision (ICCV), pages 1637–1645, 2019.
  * [3] Zhipeng Cai, Tat-Jun Chin, Huu Le, and David Suter.  Deterministic consensus maximization with biconvex programming.  In Proceedings of the European Conference on Computer Vision (ECCV), pages 685–700, 2018.
  * [4] Dylan Campbell, Lars Petersson, Laurent Kneip, and Hongdong Li.  Globally-optimal inlier set maximisation for simultaneous camera pose and feature correspondence.  In Proceedings of the IEEE International Conference on Computer Vision, pages 1–10, 2017.
  * [5] T. Chin, P. Purkait, A. Eriksson, and D. Suter.  Efficient globally optimal consensus maximisation with tree search.  In 2015 IEEE Conference on Computer Vision and Pattern Recognition (CVPR), pages 2413–2421, 2015.
  * [6] Tat-Jun Chin, Zhipeng Cai, and Frank Neumann.  Robust fitting in computer vision: Easy or hard?  In Proceedings of the European Conference on Computer Vision (ECCV), pages 701–716, 2018.
  * [7] Tat-Jun Chin and David Suter.  The maximum consensus problem: recent algorithmic advances.  Synthesis Lectures on Computer Vision, 7(2):1–194, 2017.
  * [8] Irit Dinur and Samuel Safra.  On the hardness of approximating minimum vertex cover.  Annals of mathematics, pages 439–485, 2005.
  * [9] David Eppstein.  Quasiconvex programming.  Combinatorial and Computational Geometry, 52(287-331):3, 2005.
  * [10] Martin A. Fischler and Robert C. Bolles.  Random sample consensus: A paradigm for model fitting with applications to image analysis and automated cartography.  Commun. ACM, 24(6):381–395, June 1981.
  * [11] Andreas Geiger, Philip Lenz, and Raquel Urtasun.  Are we ready for autonomous driving? the kitti vision benchmark suite.  In 2012 IEEE Conference on Computer Vision and Pattern Recognition, pages 3354–3361. IEEE, 2012.
  * [12] Richard Hartley and Andrew Zisserman.  Multiple view geometry in computer vision.  Cambridge university press, 2003.
  * [13] Q. Ke and T. Kanade.  Quasiconvex optimization for robust geometric reconstruction.  In Tenth IEEE International Conference on Computer Vision (ICCV’05) Volume 1, volume 2, pages 986–993 Vol. 2, 2005.
  * [14] A D Korshunov.  Monotone boolean functions.  Russian Mathematical Surveys, 58(5):929–1001, oct 2003.
  * [15] E.G. Kul’yanov.  An algorithm for finding the maximal upper zero of an arbitrary monotonic function of the algebra of logic.  USSR Computational Mathematics and Mathematical Physics, 15(4):267 – 269, 1975.
  * [16] Huu Le, Tat-Jun Chin, and David Suter.  An exact penalty method for locally convergent maximum consensus.  In Proceedings of the IEEE Conference on Computer Vision and Pattern Recognition, pages 1888–1896, 2017.
  * [17] K Lebeda, J Matas, and O Chum.  Fixing the locally optimized ransac.  In British machine vision conference, 2012, 2012.
  * [18] Hongdong Li.  Consensus set maximization with guaranteed global optimality for robust geometry estimation.  In 2009 IEEE 12th International Conference on Computer Vision, pages 1074–1080. IEEE.
  * [19] Hongdong Li.  A practical algorithm for l triangulation with outliers.  In 2007 IEEE Conference on Computer Vision and Pattern Recognition, pages 1–8. IEEE, 2007.
  * [20] David G Lowe.  Object recognition from local scale-invariant features.  In Proceedings of the seventh IEEE international conference on computer vision, volume 2, pages 1150–1157. Ieee, 1999.
  * [21] Ryan O’Donnell.  Analysis of Boolean Functions.  Cambridge University Press, 2014.
  * [22] Carl Olsson, Olof Enqvist, and Fredrik Kahl.  A polynomial-time bound for matching and registration with outliers.  In 2008 IEEE Conference on Computer Vision and Pattern Recognition, pages 1–8. IEEE, 2008.
  * [23] Thomas Probst, Danda Pani Paudel, Ajad Chhatkuli, and Luc Van Gool.  Convex relaxations for consensus and non-minimal problems in 3d vision.  In Proceedings of the IEEE International Conference on Computer Vision, pages 10233–10242, 2019.
  * [24] Rahul Raguram, Ondrej Chum, Marc Pollefeys, Jiri Matas, and Jan-Michael Frahm.  Usac: a universal framework for random sample consensus.  IEEE transactions on pattern analysis and machine intelligence, 35(8):2022–2038, 2012.
  * [25] Srikumar Ramalingam, Chris Russell, Lubor Ladicky, and Philip H.S. Torr.  Efficient minimization of higher order submodular functions using monotonic boolean functions.  Discrete Applied Mathematics, 220:1 – 19, 2017.
  * [26] Kristy Sim and Richard Hartley.  Removing outliers using the l$\backslash$infty norm.  In 2006 IEEE Computer Society Conference on Computer Vision and Pattern Recognition (CVPR’06), volume 1, pages 485–494. IEEE, 2006.
  * [27] P. H. S. Torr and D. W. Murray.  The development and comparison of robust methods for estimating the fundamental matrix.  International Journal of Computer Vision, 24(3):271–300, 1997.
  * [28] A. Vedaldi and B. Fulkerson.  VLFeat: An open and portable library of computer vision algorithms.  http://www.vlfeat.org/, 2008.
  * [29] Hoi Sim Wong, Tat-Jun Chin, Jin Yu, and David Suter.  Dynamic and hierarchical multi-structure geometric model fitting.  In 2011 International Conference on Computer Vision, pages 1044–1051. IEEE, 2011.
  * [30] Heng Yang and Luca Carlone.  A quaternion-based certifiably optimal solution to the wahba problem with outliers.  In Proceedings of the IEEE International Conference on Computer Vision, pages 1665–1674, 2019.
  * [31] Y. Zheng, S. Sugimoto, and M. Okutomi.  Deterministically maximizing feasible subsystem for robust model fitting with unit norm constraint.  In CVPR 2011, pages 1825–1832, 2011.