# Selective Cascade of Residual ExtraTrees

Qimin Liu
Department of Psychology and Human Development
Vanderbilt University
Nashville, TN
<EMAIL_ADDRESS>
Fang Liu, Ph.D.
Department of Applied and Computational Mathematics & Statistics
University of Notre Dame
Notre Dame, IN
<EMAIL_ADDRESS>
###### Abstract

We propose a novel tree-based ensemble method named Selective Cascade of
Residual ExtraTrees (SCORE). SCORE draws inspiration from representation
learning, incorporates regularized regression with variable selection
features, and utilizes boosting to improve prediction and reduce
generalization errors. We also develop a variable importance measure to
increase the explainability of SCORE. Our computer experiments show that SCORE
provides comparable or superior performance in prediction against ExtraTrees,
random forest, gradient boosting machine, and neural networks; and the
proposed variable importance measure for SCORE is comparable to studied
benchmark methods. Finally, the predictive performance of SCORE remains stable
across hyper-parameter values, suggesting potential robustness to hyper-
parameter specification.

Conflicts of Interest: The authors declare that they have no conflicts of
interest.

keywords: Extremely randomized trees, boosting, ensemble learning, regularized
regression, variable importance measure, explainability

##  1 Introduction

Ensemble learning methods combine multiple learning algorithms to achieve
better learning performance than that from individual learning algorithms. The
success of ensemble methods, in part, can be attributed to the diversity among
the constituent members, which helps mitigate over-fitting and reduce the
generalization error[1].

We focus on tree-based ensemble methods for regression. Tree-based ensemble
methods construct more than one decision tree and, by appreciating the
diversity among the trees, increase the generalizability of the ensemble. Such
diversity often comes from perturbation in the optimization process of
individual trees. For example, tree bagging generates bootstrap samples, from
which decision trees are obtained and averaged [2]. Random subspace selects a
pseudo-random subset of features when building trees [3]. Random forest (RF)
selects a random subset of the features for each candidate split, and trains
the trees with bootstrap samples [4].

The computational costs increase drastically as the number of trees increase
in tree-based ensemble approaches given both the optimization and the
perturbation processes. One way to reduce the computational burden is to
replace optimization with randomization processes. C4.5, for example, randomly
chooses from the best 20 splits in constructing individual trees [1];
extremely randomized trees (ExtraTrees) randomly selects cut-point for each
candidate split variable [5]; and [6] proposes “randomly” choosing a non-
tested feature at each level of the tree without using any training data.

However, the extreme randomness and lack of optimization can lead to
prediction bias for individual tree. [7] suggests that a continuum of tree
randomization exists between two extremes: deterministic algorithms with
injected randomization and the complete randomization paradigm; and proposes
an algorithm to generate a range of models between the two extremes. Another
approach to reduce bias is through gradient boosting machine (GBM) [8, 9].
Rather than combining trees independently in a parallel fashion, GBM builds
one tree at a time in sequential orders and each new tree aims to predict the
residuals from the previous tree (i.e., functional gradient descent). GBM
performs better than RF if parameters are tuned carefully. However,
disadvantages exist for GBM: First, hyperparameters in GBM is arguably harder
to tune than RF; Second, GBM can be prone to overfitting; Third, the model
training can be computationally burdensome given the sequential order in tree
construction.

Recent works have emerged in embedding decision trees in complex structures,
such as neural networks (NNs) [10], to further reduce bias. [11] stacks the RF
onto deep NNs to improve the prediction accuracy in classification for the
$n<p$ case. [12] embeds RF in a NN architecture, referred to as the “deep
forest” approach. [13] adopts a multi-layer architecture where gradient
boosting decision trees serve as base learners. With output from prior layers
functioning as input for next layers, these methods optimize base learners in
each layer by improving the output of the prior layer so that the quality of
input increases for the subsequent layers. The effectiveness of these
approaches highlights that optimized base learners associate with predictive
accuracy for the current layer and serve as useful inputs in the multi-layer
context.

As one of multi-layer learning methods, NNs have gained vast popularity in
both machine learning and applied sciences. However, for multi-layer methods
with deep structures, difficulty arises in tracing the contributions of the
input features. That is, the deep structure can generate black-box predictions
that lack interpretability or explainability. The deep structure in NNs also
requires a relatively large amount of training data and careful tuning of
hyper-parameters. NN models often 1) involve a large of amount of parameters,
2) can be subject to overfitting were inappropriate or insufficient
regularizations employed, 3) rely on high-performance computational resources
to yield fast training and learning. For these reasons, NNs do not necessarily
offer the best performance for all types of learning tasks on all types of
data. Regardless, representation learning, as one of the key properties of
NNs, has received wide attention for its role in enhancing predictive
performance [14]. We propose to regard decision trees as representation
features. Specifically, we develop an ensemble learning approach that
incorporates representation learning, sparsity regularization, and boosting,
to achieve accurate predictions.

We term our proposed method Selective Cascade of Residual ExtraTrees (SCORE).
Our proposed tree-based ensemble method can offer high prediction accuracy,
explainability regardless of deep structures, and low computational costs.
Methodologically, SCORE incorporates several effective machine learning
strategies, including extremely randomly trees [5], boosting [15, 16],
sparsity regularization [17], and representation learning [18]. The reasons
for choosing each of the specific learning techniques to build SCORE are
provided below.

  * •

To maximize the diversity of the base learners and reduce computational costs,
SCORE employs the extremely randomized trees (ExtraTrees) as the base
learners. In our experiments, the computational cost of SCORE is much lower
than the benchmarks methods.

  * •

To reduce prediction bias in each individual tree in the ensemble and set up a
framework to measure “importance” of the input features for predictions, SCORE
filters out irrelevant trees for prediction via sparsity-regularized
regression, referred to as the TreeSelection step. The selected trees are
regarded as representation features. The TreeSelection is a critical step for
SCORE and differentiates it from other tree ensemble approaches. Existing tree
ensemble approaches often include irrelevant features in building trees or
constructing layers, which can associate with low prediction accuracy. In our
experiments, the predictions from SCORE are significantly better in general
than its competitors without a tree selection step.

  * •

To further enhance prediction accuracy, we employ the boosting technique by
forming multiple layers of ExtraTrees + TreeSelection, where the outcome in
each layer is the residuals given the predictions up to that layer. Boosting
itself is an ensemble method for improving prediction, by training weak
learners sequentially and combining them to form a strong learner. In the
context of SCORE, each layer of ExtraTrees + TreeSelection can be regarded as
a weak learner.

  * •

To improve the explainability of the predictions by SCORE, we develop a
conceptually and computationally simple yet effective variable importance
measure (VIM). The VIM correctly identifies the top predictors in our
experiments.

In summary, SCORE leverages the advantages of the existing tree ensemble
approaches but improves on computational cost by using a random set of base
learners, on outcome prediction via boosting and regularization, and on
explainability with a VIM. The ensemble is reflected not only horizontally
through building random tress but also vertically via boosting. In addition,
our experiments suggests that SCORE does not require precise tuning on the
hyperparameters compared to other state-of-the-art ensemble methods to yield
comparable or superior results.

##  2 Selective Cascade of Residual ExtraTrees (SCORE)

Let $y$ represent the observed outcome and $\bar{f}$ be the prediction of $y$
from an ensemble method. [19] decomposes the generalization error
$E(\bar{f}-y)^{2}$ of the ensemble method as follows

$\begin{split}\textstyle
M^{-2}\left[(\sum_{i}(E(f_{i})-y))^{2}+\sum_{i}E(f_{i}-E(f_{i}))^{2}+\sum_{i}\sum_{j\neq
i}E\left((f_{i}-E(f_{i}))(f_{j}-E(f_{j}))\right)\right]\end{split}$ (1)

where $M$ is the number of individual base learners in the ensemble and
$f_{i}$ refers to prediction from the $i$-th learner in the ensemble for
$i=1,\ldots,M$. The terms in Eqn (1) are respectively the squared average
bias, average variance, and average covariance across the ensemble members.
The covariance serves as an indicator of the ensemble diversity. The smaller
the covariance, the greater diversity there is among the ensemble members.

Our proposed SCORE procedure aims at 1) maximizing the diversity and thus
reducing the covariance term among the ensemble trees in Eqn (1) by employing
the extremely randomized tree procedure; 2) reducing the mean squared error
term (the sum of the squared bias and the variance terms in Eqn (1) by
implementing a selection step via sparsity-promoting regularized regression on
individual trees as the representation features; 3) further reducing the bias
and the variance of the ensemble by using weighted averages of predictions
from the multi-layered boosting. Algorithm 1 lists the steps for the SCORE
procedure. Figure 1 depicts an example of a 2-layer SCORE with 4 original
attributes ($X_{1}$ to $X_{4}$), along with the algorithmic flowchart.
$T_{i}^{(l)}$ refers to the $i$-th individual randomized tree constructed in
layer $l$ for $l=1,2$. $Y^{(l)}$ is the outcome in layer $l$, which is the
residual from the regression in layer $l-1$ for $l>1$ (that is,
$Y^{(l)}=Y^{(l-1)}-\hat{Y}^{(l-1)}$), and is the original $Y$ when $l=1$.

Algorithm 1 The SCORE Algorithm

1:input: training set $\mathcal{D}_{t}$; validation set $\mathcal{D}_{v}$;
stopping threshold $\tau_{0}$

2:output: trained SCORE ensemble leaner and VIM.

3:Employ ExtraTrees to construct an initial set of randomized trees (base
learners) on $\mathcal{D}_{t}$. Calculate the mean squared prediction error
MSE(1) on $\mathcal{D}_{v}$.

4:Let $\Delta^{(1)}$= MSE(1); $l\leftarrow 2$

5:while $\Delta^{(l-1)}>\tau_{0}$ do

6: TreeSelection: Apply sparsity-regularized regression to outcome $Y^{(l-1)}$
to select relevant representation features and calculate residuals $Y^{(l)}$.

7: ExtraTrees: Employ ExtraTrees to predict $Y^{(l)}$ with either the
global-$X$ or local-$X$ input feature option.

8: Calculate MSE(l) on $\mathcal{D}_{v}$, and set $\Delta^{(l)}\\!=\\!$
MSE${}^{(l-1)}\\!-$MSE(l).

9: $l\leftarrow 1+1$

10:end while

11:Calculate VIM.

Figure 1: An example ensemble structure in SCORE and algorithm flowchart
(green circles with solid edges and red circles with dashed edges represent
selected and non-selected trees via TreeSelection, respectively.)

In what follows, we describe the ExtraTrees and TreeSelection steps of the
SCORE procedure in more details, and propose an approach to measure the
importance of the original input feature which aids the interpretability of
the SCORE procedure and its output. We use $\mathbf{X}^{(l)}$ to denote the
predictor attributes in layer $l$, and provide two options. The first option,
referred to as the _local $X$ option_, limits $X^{(l)}$ only to the
representation attributes learned from layer $l-1$ and selected by the
TreeSelection step (see Sec 2.2). The second option, referred to as the
_global $X$ option_, allows $\mathbf{X}^{(l)}$ to include not only the
selected representation attributes from layer $l-1$, but also the original
input features as well as the learned representation features from all
previous layers.

###  2.1 Construction of Extremely Randomized Trees

The first step in each layer is the construction of extremely randomized trees
as base learners. We employ the ExtraTrees procedure to achieve that goal.
ExtraTrees offers prediction accuracy comparable to RF and higher than tree
bagging due to the “extreme” randomness. The randomization process in
ExtraTrees drastically increases diversity among trees at a lower
computational cost compared to RF. The reduced computational costs come from
decreased optimization in ExtraTrees as optimization process can be a limiting
factor for computation speed.

ExtraTrees starts by randomly selecting $K$ attributes out of all attributes.
Then, ExtraTrees chooses a random cut-point from the observed range of each of
the $K$ attributes. After that, ExtraTrees selects the attribute resulting in
the smallest MSE as the final split variable. A pre-set minimal node size can
be used as a stopping rule for the splitting process. Denote the number of
randomized trees in layer $l$ by $M$ (for notation simplicity, we set $M$
fixed across all layers in SCORE, but $M$ can vary by layers), and let
$T_{i}^{(l)}$ be the $i$-th randomized tree for $i=1,\ldots,M$. We express the
tree-based regression for each tree $T_{i}^{(l)}$ in layer $l$ as

$Y^{(l)}=T_{i}^{(l)}(\mathbf{X}^{(l)})+\epsilon_{i}^{(l)},$ (2)

where $\epsilon_{i}^{(l)}$ is the error term, the distribution of which does
not need to be specified.

###  2.2 Tree Selection via Sparsity Regularization

The predicted outcome $\hat{T}_{i}^{(l)}$ from each randomized tree can be
regarded as a representation feature. Not every representation feature is a
good predictor of $Y^{(l)}$, especially given the way how trees are built in
ExtraTrees. Including irrelevant features in a regression model can lead to
increased bias in prediction and generalization errors. To mitigate this
problem, SCORE applies a sparsity-promoting regularization to select more
relevant representation features from the $M$ randomized trees. The
regularized loss function in layer $l$ is

$L\left(Y^{(l)},\boldsymbol{\beta}^{(l)}\mathbf{\hat{T}}^{(l)}\right)+\lambda^{(l)}R\left(\boldsymbol{\beta}^{(l)}\right),$
(3)

where $\boldsymbol{\beta}^{(l)}$ quantifies the “weights” for the trees in
$\mathbf{\hat{T}}^{(l)}=(T_{1}^{(l)},\ldots,T_{M}^{(l)})$, and
$R(\boldsymbol{\beta}^{(l)})$ is the regularizer on $\boldsymbol{\beta}^{(l)}$
with tuning parameter $\lambda^{(l)}$. Different types of loss function can be
used on $L$, such as $l_{1},l_{2}$ or negative log-likelihood. If the least-
squared regression ($l_{2}$ loss function) with the $l_{1}$ regularizer [20]
is used, then Eqn (3) becomes

$\big{|}\big{|}Y^{(l)}-\boldsymbol{\beta}^{(l)}\mathbf{\hat{T}}^{(l)}\big{|}\big{|}_{2}^{2}+\lambda^{(l)}||\boldsymbol{\beta}^{(l)}||_{1}.$
(4)

Minimizing Eqn (4) leads to estimate
$\hat{\boldsymbol{\beta}}^{(l)}=(\hat{\beta}_{1}^{(l)},\ldots,\hat{\beta}_{1}^{(M)})$,
some of which are exactly 0 with the $l_{1}$ regularization. The trees
associated with non-zero $\hat{\beta}$ make the set (local $X$ option) or a
subset (global $X$ option) of the input features for the next layer. The
residual
$Y^{(l)}-\hat{Y}^{(l)}=Y^{(l)}-\hat{\boldsymbol{\beta}}^{(l)}\mathbf{\hat{T}}^{(l)}$
serves as the outcome $Y^{(l+1)}$ for the next layer.

###  2.3 Cascade Layers and Final Outcome Prediction

Multiple layers with the dual ExtraTrees+TreeSelection step are built, each
with a newly calculated outcome and a new set of input features (newly
identified representation tree features from the last layer, plus or not the
input features from all previous layers). The layers are stacked on until
there is no meaningful improvement on the prediction MSE for the validation
set. Suppose there are totally $L$ layers after running the SCORE algorithm,
the overall prediction model can be written as

$\displaystyle
Y=\textstyle\sum^{L}_{l=1}\left(\beta_{0}^{(l)}+\boldsymbol{\beta}^{(l)}\mathbf{T}^{(l)}(\mathbf{X}^{(l)})\right)+\varepsilon^{(L)}.$
(5)

Denote the final output prediction by $\hat{Y}$. If each layer is weighted
equally, then

$\displaystyle\hat{Y}=\textstyle\sum^{L}_{l=1}\left(\hat{\beta}_{0}^{(l)}+\hat{\boldsymbol{\beta}}^{(l)}\hat{\mathbf{T}}^{(l)}(\mathbf{X}^{(l)})\right).$
(6)

A generalized version of the SCORE prediction is to apply different weights to
different layers, say by incorporating a learning rate $\gamma\in(0,1]$ in a
similar manner as in GBM,

$\displaystyle\hat{Y}=\textstyle\sum^{L}_{l=1}\gamma^{l}\left(\hat{\beta}_{0}^{(l)}+\hat{\boldsymbol{\beta}}^{(l)}\mathbf{T}^{(l)}(\mathbf{X}^{(l)})\right).$
(7)

When $\gamma=1$, it reduces to equation (6). For $\gamma<1$, the deeper a
layer is, the less contribution it will have toward the final prediction of
$Y$. The weighting scheme can help to mitigate overfitting with an oversized
cascade structure.

###  2.4 Variable Importance Measure (VIM)

The VIM we propose for SCORE is based on the weighted frequency of each
variable being selected for node splitting in all the randomized trees, where
the weights are the estimated coefficients $\hat{\boldsymbol{\beta}}$ from the
the regularized regression in Eqn (3). If a regression coefficient
$\beta^{(l)}$ is associated with the original input features ($l=1$ with local
X option; and $l\geq 1$ with global X option), then it can be used directly to
quantify the contributions of that feature in the outcome prediction. For
$\beta^{(l)}$ associated with a representation feature, its value will be
traced back, layer by layer, to the original input features used for building
these trees. To achieve this, we define two frequency vectors
$\mathbf{c}_{X,i}^{(l)}$ and $\mathbf{c}_{T,i}^{(k,l)}$.

Denote the dimension of the original $\mathbf{X}$ by $p$. Let
$\mathbf{c}_{X,i}^{(l)}$ be a column vector of length $p$ for the $i$-th
randomized tree in the $l$-th layer. $\mathbf{c}_{X,i}^{(l)}$ contains the
frequency on each input attribute used as a split variable in that tree. For
example, suppose $\mathbf{X}=(X_{1},X_{2},X_{3},X_{4},X_{5})$ ($p=5$) and the
1st tree in the 1st layer uses $X_{3},X_{4}$ and $X_{5}$ as split variables,
each once at some nodes in the tree, then
$\mathbf{c}^{(1)}_{X,1}=(0,0,1,1,1)$. If $\mathbf{X}^{(l)}$ employs the global
X option, then the original $\mathbf{X}$ participates in constructing trees in
each layer. For example, $X_{3}$ is used twice as the split variable in two
nodes, and $X_{1}$ and $X_{5}$ once each as the split variable in the 10-th
tree in the 2nd layer, then $\mathbf{c}^{(2)}_{X,10}=(1,0,2,0,1)$. If the
local X option is used, then the original $\mathbf{X}$ no longer directly
contributes to the construction of trees for layer $l>1$ (that is,
$\mathbf{c}^{(l)}_{X,i}\equiv\mathbf{0}$ for $l>2$ in the local X option
case), but its contribution will be acknowledged through the selected trees in
the previous layer.

We also define $\mathbf{c}_{T,i}^{(k,l)}$ for the layer pair $(l,k)$: $l\geq
2$ and $1\leq k\leq l-1$, which is a column vector of length $m^{(k)}$,
containing the frequency of the learned representation attributes from a
previous layer $k$ that are selected as split variables for the $i$-th
randomized tree in the $l$-th layer. If the local X option is employed, then
$\mathbf{c}_{T,i}^{(k,l)}\equiv 0$ for $k<l-1$. Imagine that TreeSelection
chooses $m^{(1)}=50$ representation features out of $M=500$ randomized trees
from layer 1. The 5-th randomized tree in the 2nd layer of SCORE uses 15
representation features, once each, out of the $50$ features. Then
$\mathbf{c}_{T,5}^{(1,2)}$ is a vector of length 50, with 1 in the positions
corresponding to the used 15 representation features and 0 elsewhere.

We are now ready to define the VIM for the original input features
$\mathbf{X}$. If the local X option is employed, then the VIM for $\mathbf{X}$
is defined as follows:

$\displaystyle\\!\mbox{VIM}\\!=\\!\textstyle\sum_{l=1}^{L}\\!\mbox{VIM}^{(l)}\\!=\\!\sum_{l=1}^{L}\\!\bigg{\\{}\\!\sum_{i_{1}=1}^{M^{(1)}}\\!\sum_{i_{2}=1}^{M^{(2)}}\\!\ldots\\!\sum_{i_{l}=1}^{M^{(l)}}\\!|\beta_{i_{l}}^{(l)}|c^{(l-1,l)}_{T,i_{l},i_{l-1}}\\!c^{(l-2,l-1)}_{T,i_{l-1},i_{l-2}}\ldots
c^{(1,2)}_{T,i_{2},i_{1}}\\!\mathbf{c}^{(1)}_{X,i_{1}}\\!\\!\bigg{\\}}.$ (8)

$\mbox{VIM}^{(l)}$ quantifies the importance of the original $\mathbf{X}$ in
learning the $l$-th layer of the SCORE structure, $c^{(l-1,l)}_{T,i,j}$ denote
the $j$-th element in $\mathbf{c}^{(l-1,l)}_{T,i}$, and $M^{(l)}$ refers to
the number of representation features selected via the regularized regression
in layer $l$. If the global X option is employed, the calculation of VIM is
more complicated as the original and learned representation features not only
make the input attributes for the next immediate layer, but possibly also for
each layer after that. Hence, the defined VIM has to account for this
multiplicity effect. Specifically, we calculate

  * •

the direct contribution of $\mathbf{X}$ to each layer $l$ for $1\leq l\leq L$

$\displaystyle\\!\\!\\!\\!\mbox{VIM:direct}=\textstyle\sum_{l=1}^{L}\mbox{VIM}^{(l)}\\!=\\!\sum_{l=1}^{L}\sum_{i=1}^{M^{(l)}}|\beta_{i}^{(l)}|\mathbf{c}^{(l)}_{X,i};$
(9)

  * •

the indirect contribution of $\mathbf{X}$ through the selected representation
features in layer $l$ for $2\leq l\leq L$

$\displaystyle\textstyle\mbox{VIM:indirect}=\sum_{l=2}^{L}\mbox{VIM}^{(l)}$
$\displaystyle=$
$\displaystyle\sum_{l=2}^{L}\sum_{i_{1}=1}^{M^{(1)}}\\!\sum_{i_{2}=1}^{M^{(2)}}\ldots\sum_{i_{l-1}=1}^{M^{(l-1)}}\\!\sum_{i_{l}=1}^{M^{(l)}}\left\\{|\beta_{i_{l}}^{(l)}|c^{(l,l-1)}_{T,i_{l},i_{l-1}}c^{(l-1,l-2)}_{T,i_{l-1},i_{l-2}}\ldots
c^{(2,1)}_{T,i_{2},i_{1}}\mathbf{c}^{(1)}_{X,i_{1}}\right\\}.$ (10)

The final VIM for the original $\mathbf{X}$ is defined as the sum of Eqns (9)
and (• ‣ 2.4),

VIM = VIM:direct + VIM:indirect (11)

##  3 Experiments

We compare SCORE with three ensemble methods – ExtraTrees, RF, GBM – in
prediction accuracy, variable importance, and computational efficiency. In
addition, we include NN as another benchmark method due to our use of
representation learning. The reasons for choosing the ExtraTrees, RF and GBM
for comparison are as follows. First, SCORE builds parallel trees similar to
RF and ExtraTrees . Second, SCORE incorporates boosting using functional
gradients similar to GBM. Third, ExtraTrees, as one of the base steps for
SCORE, epitomizes the use of extreme randomization to maximize diversity among
base learners. All compared methods enjoy vast popularity among practitioners
for prediction tasks. In addition, all included methods have accompanying
VIMs.

###  3.1 Experiment Setup

We run five experiments. The first three experiments use the Friedman 1, 2,
and 3 data sets, the fourth runs on the real-life Boston Housing data, and the
fifth on the real-life World Happiness Report data [21].

The outcome $Y$ in Friedman 1 data is simulated as $Y\\!=\\!10\sin(\pi
X_{1}X_{2})\\!+\\!20(X_{3}-0.5)^{2}+10X_{4}+5X_{5}+\varepsilon$, where
$X_{j}\sim U(0,1)$ for $j=1,\ldots,5$ and $\varepsilon\sim N(0,1)$. To examine
the efficiency of each method in measuring variable important, we add
irrelevant features to the training and testing data. Specifically, we examine
two scenarios: 1) 5 additional irrelevant input attributes sampled from
$U(0,1)$, totaling to $p=10$; 2) 45 irrelevant input attributes, of which 5
are sampled from $U(0,1)$, 20 from $N(0,1)$, and 20 from Lognormal$(0,1)$,
totaling to $p=50$. We examine different training sizes at
$10p,50p,100p,200p,400p$, respectively; and the size of the testing data is
fixed at 5,000. The outcomes in the Friedman 2 and Friedman 3 experiments are
simulated from $Y\\!=\\!\displaystyle
X_{1}^{2}+(X_{2}X_{3}-(X_{2}X_{4})^{-2})^{0.5}+\varepsilon$ with
$\varepsilon\sim N(0,125)$, and
$Y=\displaystyle\mbox{atan}((X_{2}X_{3}-(X_{2}X_{4})^{-1}))/X_{1})+\varepsilon$
with $\varepsilon\sim N(0,0.1)$, respectively, where $X_{1}\sim U[0,100],\
X_{2}\sim U[40\pi,560\pi],\ X_{3}\sim U[0,1]$, $X_{4}\sim U[1,11]$ in both
cases. Five irrelevant $X$’s are simulated from $N(0,1)$ and added to the
training and testing data ($p=9$) in Friedman 2 and Friedman 3 data. Both the
training and the test data are of size $n=5,000$. The Boston housing data
contain 506 census tracts of Boston from the 1970 census. The outcome variable
is housing price. There are 17 input attributes including housing locations,
location characteristics, etc (see Appendix for the list of variables). The
world happiness report data (https://worldhappiness.report/ed/2018/) contain
the survey data on the state of happiness. The outcome is the happiness index,
with 10 being the best and 0 being the worst possible life. After excluding
data with missingness, we include 382 cases. There are 17 input attributes
which measure the economic, quality-of-life, and political aspects of
different countries (see Appendix for the list of variables).

R packages randomForest, ranger, gbm were employed to run RF, ExtraTrees, and
GBM, respectively; and the R package caret was used to tune the hyper-
parameters with the 10-fold cross validation (for GBM, the maximum number of
trees was fixed at 500 in each experiment).

For NN, the number of hidden layers was set at 2 to match the depth of SCORE
for all 5 experiment (for the Friedman 1 experiment, NNs with one hidden layer
were also examined but the performance was not good). R package neuralnet was
used to compare the number of hidden nodes per hidden layer, early stopping
rules, and learning rates. The hyper-parameter settings that lead to the
highest prediction accuracy on the test data were used: 2 hidden layers with
10 hidden nodes per layer in the Friedman 1 and 2 experiments and the Boston
Housing data, 2 hidden layers with 5 hidden nodes per layer in the Friedman 3
experiment and the World Happiness Report data.

For SCORE, we examine 9 hyper-parameter settings at each $n$ in the Friedman 1
experiment (Table 1) by varying the minimal node size to training size ratio
($1/5,1/4,1/3$), the number of layers and the number of trees per layer ($L=1$
with 500 trees, $L=2$ with 250 trees per layer), and the global $X$ vs. local
$X$ options for $L=2$, in the ExtraTrees step. For the TreeSelection step in
SCORE, we employed the R package glmnet to run the lasso regression with the
built-in procedure to tune the hyper-parameter.

Table 1: Nine SCORE hyper-parameter scenarios in the Friedman 1 experiment minimal node size to | $L=1$ | $L=2$
---|---|---
training size ratio |  | local X option | global X option
1/5 | 1 | 4 | 7
1/4 | 2 | 5 | 8
1/3 | 3 | 6 | 9

In the other 4 experiments, rather than tuning the hyper-parameters in each
experiment, we examined findings from the Friedman 1 experiment to evaluate
how the hyper-parameters impact the performance in SCORE and pre-specified the
tuning parameters in other experiments based on it. Specifically, we used a
two-layer SCORE structure with the global $X$ option and 5 trees as the
minimal node size for Friedman 2 and 3 and the Boston housing data, 25 trees
as the minimal node size for the world happiness data.

For the VIM measures in ExtraTrees, RF, and GBM, the R function varImp from
package caret was applied. Though NNs, especially those with deep structures,
are known for lacking explainability, there exists some VIMs for fully-
connected feedforward NNs [22], which are the types of NNs we applied to the
five experiments. We applied the Olden method via R package NeuralNetTools to
calculate the VIM for the NNs, which measures the relative importance of the
input features as the sum of the product of the raw input-to-hidden, hidden-
to-output weights.

Across all 5 methods, the prediction, VIM and computational time results were
averaged across 100 repeats. For the Friedman 1, 2, and 3 experiments, we
directly simulated 100 sets of training data from the models. For the Boston
housing data, we randomly selected 400 cases as the training samples without
replacement, and used the rest (106) as testing samples. The process was
repeated 100 times. For the World happiness data, we randomly selected 250
cases as the training sample and used the rest (132) as the test samples; and
the process was repeated 100 times to generate 100 repetitions.

###  3.2 Results on Prediction Accuracy

The prediction results for the Friedman 1 data are displayed in Figure 2.

Figure 2: Average ($\pm$ SD) prediction MSE over 100 repeats in Friedman 1
data

First, the prediction improves in all methods as $n$ increases. Compared to
the 3 ensemble methods, SCORE shows smaller prediction MSE and better
performance over ExtraTrees and RF, and was generally comparable to GBM across
all the experimental settings. Compared to NN, SCORE is universally better
when $p=50$, and for $n=100;500;1,000$ when $p=10$; NN is slightly better than
SCORE at $n=4,000$ for $p=10$ and $p=50$, at $2,000$ when $p=10$. In sum, good
prediction of NN requires a large $n$. In addition, the prediction based on NN
is less stable compared to other methods (much larger SD of MSE for all
scenarios). There are some differences among the 9 SCORE cases with varying
hyper-parameter settings, but they are generally comparable. Specifically,
having more layers (SCORE cases 4 to 9) improves the prediction accuracy
compared to the single-layered SCORE cases (1 to 3). The minimal node size to
the training size ratio seems to affect the predictive accuracy: increasing
ratio seems to lead to higher predictive accuracy. This is not unexpected as
larger trees tend to over-fit whereas the higher ratio helps control over-
fitting. The global $X$ option (cases 7 to 9) offers slight advantages over
the local $X$ option (cases 4 to 6). Overall, the small differences among
different hyper-parameter settings imply that SCORE is relatively robust to
the hyper-parameter specification and can deliver comparable prediction
without much tuning.

Table 2 presents the prediction accuracy for the Friedman 2, Friedman 3, the
Boston housing, and the world happiness experiments. SCORE delivers the
smallest MSE in general compared to the three 3 approaches.

Table 2: Mean Prediction MSE (SD) in Experiments 2 to 5 | Friedman 2 | Friedman 3 | Boston | happiness
---|---|---|---|---
|  |  | housing | report
SCORE |  $\mathbf{16.6\times 10^{3}}$ (418.3) |  $\mathbf{0.0119}\;(3.443\times 10^{-4}$) |  9.0 (2.1) |  0.07 (0.015)
ExtraTrees |  $17.1\times 10^{3}$ (367.8) |  $0.0120\;(3.235\times 10^{-4}$) | 13.0 (4.0) | 0.08 (0.014)
RF |  $17.1\times 10^{3}$ (360.7) |  $0.0120\;(3.197\times 10^{-4}$) | 10.3 (3.5) |  0.07 (0.016)
GBM |  $35.4\times 10^{3}$ (728.2) |  $0.0236\;(6.031\times 10^{-4}$) | 11.8 (3.2) | 0.11 (0.018)
NN |  $21.5\times 10^{3}$ (4277.0) |  $0.0127\ (2.280\times 10^{-3}$) | 16.2 (5.) | 0.10 (0.034)

###  3.3 Results on Variable Importance

The list of the identified important variables are presented in Table 3, based
on the average normalized VIM across the 100 repetitions. We applied $1/p$ as
the cutoff to identify the relevant attributes in all methods.

Table 3: Identified Relevant Input Attributes

Data | Friedman | Boston | happiness
---|---|---|---
$(1/p)^{\ddagger}$ | 1 ($p\\!=\\!50;n\\!=\\!5,000;2\%$) | 2 (11.1%) | 3 (11.1%) | housing (5.88%) | report (5.88%)
SCORE | $X_{2},X_{1},X_{3},X_{5},X_{4}$ | $X_{3},X_{2}$ | $X_{3},X_{1},X_{2}$ | rm, lstat, age, crim, lon, | sdml, sdl, lGDPpc,
|  |  |  | b, dis, tract, lat, nox | life, posaf, del, sup
ET | $X_{4},X_{2},X_{1},X_{5},X_{3}$ | $X_{3},X_{2}$ | $X_{3},X_{2},X_{1}$ | rm, lstat | sdml, del, life,
|  |  |  |  | lGDPpc, posaf
RF | $X_{4},X_{2},X_{1},X_{5},X_{3}$ | $X_{3},X_{2}$ | $X_{3},X_{1},X_{2}$ | lstat, rm | sdml, life, sdl
GBM | $X_{4},X_{2},X_{1},X_{3},X_{5}$ | $X_{3},X_{2}$ | $X_{3},X_{2},X_{1}$ | lstat, rm | sdml, life, posaf,
|  |  |  |  | lGDPpc
NN | $X_{2},X_{1},X_{3},X_{4},X_{5},X_{8}$ | $X_{1}$ | $X_{3},X_{1},X_{2}$ | lstat, dis, rm, nox, | sdml, sdl
|  |  |  | crim, tax, tract |

The attributes are listed in the order of descending VIM.. .
---
‡ $>1/p$ is used as the cutoff for normalized VIM to identify relevant
attributes.

In the Friedman 1 experiment, all approaches succeeded in identifying $X_{1}$
to $X_{5}$ as the five most important features and the other 45 variables
(which we added artificially) as irrelevant. For SCORE and ExtraTrees, $X_{1}$
to $X_{5}$ had considerably larger weights than the 45 irrelevant features,
but ExtraTrees put much less weight on $X_{3}$ (3.6%) compared to SCORE
(11.5%).

In the Friedman 2 experiment, all the methods except for NN selected $X_{2}$
and $X_{3}$ as the top two attributes while NN picked $X_{2}$. The sub-optimal
performance is likely due to the difficulty associated with the underlying
model. The true underlying model that is simulated to simulate the data shows
that $X_{2}$ and $X_{4}$ appear together as in the product term
$(X_{2}X_{4})^{-2}$, which can be very small given the large magnitude of
$X_{2}$, resulting in a smaller contribution toward $Y$ relative to others
terms. In other words, the large VIM for $X_{2}$ likely comes from the
$X_{2}X_{3}$ term. The failure to select $X_{1}$ is likely due to the the
large spread of $X^{2}_{1}$ ($0$ to $10^{4}$) compared to the standard
deviation of the error term, making the signals from $X_{1}$ indistinguishable
from the background noise.

In the Friedman 3 experiment, all methods successfully identified
$X_{1},X_{2},X_{3}$ as relevant variables but missed $X_{4}$. Though $X_{4}$
is ranked in the 4-th place in SCORE and ExtraTrees, the corresponding VIM is
only slightly higher than for the rest. The failure for identifying $X_{4}$ in
this case is likely due to the same reason listed for Friedman 2 data.

For the Boston housing data, all methods agree that rm (the number of rooms)
and lstat (proportion of people with lower socioeconomic status) are two
important variables for predicting the housing price. For the next 3 important
variables, SCORE identified age (proportion of owner-occupied units built
prior to 1940), crim (per capita crime rate by town), and lon (longitude of
census tract); ExtraTrees identified – not listed in Table 3 – ptratio (pupil-
teacher ratio by town), tax(full-value property-tax rate per USD 10,000), and
indus(proportion of non-retail business acres per town); RF identified nox
(nitric oxides concentration), crim (per capita crime rate by town), and dis
(weighted distances to 5 Boston employment centres); GBM identified tract
(census tract), town (name of the town), dis; and NN identified dist, nox, and
crim. The 10 important predictors identified by SCORE included almost all the
top-five predictors identified by RF and NN. Note that tract (census tract),
town (name of the town), and (lon, lat) (longitude & latitude of a census
tract) containing about the same information, so identifying one or the other
does not matter much from an interpretation perspective.

For the world happiness report data, all methods agree that sdml (standard
deviation/mean of the ladder index) is the most important predictor for
happiness. For the next 4 important variables, SCORE identified sdl (standard
deviation of the ladder index), lGDPpc (log GDP per capita), life (healthy
life expectancy at birth), and posaf (positive affect); ExtraTrees identified
del (delivery quality), life (healthy life expectancy at birth), lGDPpc (log
GDP per capita), posaf (positive affect); RF identified sdml (standard
deviation/mean of the ladder index), life (healthy life expectancy at birth),
sdl (standard deviation of the ladder index), posaf (positive affect), and
lGDPpc (log GDP per capita); GBM identified, sdml (standard deviation/mean of
the ladder index), life (healthy life expectancy at birth), posaf (positive
affect), lGDPpc (log GDP per capita), and sup (social support); NN identifies
sdl (standard deviation of the ladder index), gini2 (GINI index Word Bank
estimate), gini1(GINI index Word Bank estimate, average 2000-2015), del
(delivery quality). SCORE and RF agree on all top five important variables.

In summary, the VIM procedure associated with SCORE shows good performance in
identifying important predictors and in general agrees with RF and GBM on the
important variable list in these 5 experiments.

###  3.4 Computational Time

We present the computational time from the Friedman 1 experiment (Table 4)
(averaged over 100 repeats). The observations are similar in other
experiments. The time spent on tuning hyper-parameters via the 10-fold CV is
included in all methods except for SCORE and NNs, where the hyper-parameters
were pre-specified by comparing a limited set of scenarios. All methods were
run in R version 3.5.3 on 12 core Intel(R) Haswell processors with 256 GB RAM.
In summary, SCORE was the fastest in general. The speed of NN is greatly
affected by the stopping criterion used (early stopping can save a significant
amount of time). The computational advantage for SCORE is particularly evident
for the larger $n$ case. Moreover, the computational time for two-layered
SCORE (cases 4 to 9) with 250 trees per layer is less than the singled-layered
SCORE (cases 1 to 3) with 500 trees. The computational time also decreased as
the minimal node size to training size ratio increases as the tree stops
growing earlier with a larger ratio. The global $X$ options, as using more
input attributes, was more time-consuming than the local options.

Table 4: Average Computational Time (seconds) in Friedman 1 Data

. . | $n$ | 100 | 500 | 1,000 | 2,000 | 4,000 | . .
---|---|---|---|---|---|---|---
| SCORE case 1 | 0.75 | 2.09 | 2.28 | 4.25 | 9.61 |
| SCORE case 2 | 0.74 | 2.20 | 2.32 | 4.30 | 10.39 |
| SCORE case 3 | 0.77 | 2.46 | 2.44 | 4.28 | 10.11 |
| SCORE case 4 | 0.55 | 1.92 | 6.55 | 19.21 | 43.78 |
$p=10$ | SCORE case 5 | 0.51 | 1.92 | 7.00 | 21.77 | 49.54 |
| SCORE case 6 | 0.50 | 1.97 | 7.96 | 24.03 | 59.87 |
| SCORE case 7 | 0.77 | 2.12 | 2.34 | 4.34 | 9.48 |
| SCORE case 8 | 0.75 | 2.21 | 2.36 | 4.34 | 10.34 |
| SCORE case 9 | 0.77 | 2.44 | 2.45 | 4.34 | 10.34 |
| ExtraTrees | 12.89 | 14.52 | 26.29 | 79.40 | 168.72 |
| RF | 10.07 | 79.26 | 204.00 | 754.71 | 1830.51 |
| GBM | 51.71 | 64.88 | 68.54 | 126.27 | 145.81 |
| NN | 0.38 | 8.43 |  21.23 |  75.20 |  91.44 |
| $n$ | 500 | 2,500 | 5,000 | 10,000 | 20,000 |
| SCORE case 1 | 2.96 | 5.63 | 14.94 | 42.7 | 117.10 |
| SCORE case 2 | 3.23 | 5.50 | 14.94 | 48.58 | 149.09 |
| SCORE case 3 | 3.71 | 5.55 | 15.26 | 63.14 | 180.27 |
| SCORE case 4 | 2.78 | 20.13 | 60.37 | 143.4 | 387.52 |
$p=50$ | SCORE case 5 | 2.85 | 21.85 | 63.99 | 165.95 | 521.57 |
| SCORE case 6 | 2.98 | 23.46 | 71.08 | 215.01 | 721.38 |
| SCORE case 7 | 2.92 | 5.66 | 15.11 | 42.82 | 116.46 |
| SCORE case 8 | 3.10 | 5.62 | 15.56 | 48.74 | 149.47 |
| SCORE case 9 | 3.42 | 5.62 | 15.63 | 58.77 | 183.54 |
| ExtraTrees | 79.18 | 228.78 | 470.48 | 1862.43 | 4462.34 |
| RF | 461.02 | 3967.89 | 13,084.36 | 32,237.32 | 10,6692.39 |
| GBM | 1,023.23 | 1,458.58 | 291.24 | 528.53 | 807.90 |
| NN∗ | 1.96 |  73.17 | 75.67 |  50.45 | 127.47 |

$L=1$ in SCORE 1 to 3 with 500 trees; $L=2$ in SCORE 4 to 9 with 250 trees per
layer; the minimal node size to training size ratio is 1/5 for SCORE 1, 4, 7,
1/4 for SCORE 2, 5, 8; and 1/3 for SCORE 3, 6, 8; SCORE 4 to 6 use the local X
option; and SCORE 7 to 9 employ the global X option. Refer to Table 1.
---
∗ The computational time does not increase with $n$ for NN because the
stopping criteria, which is the threshold for the partial derivatives of the
error function, are different across different $n$. Larger errors were allowed
for larger $n$ due to non-convergence issues. Specifically, the error was set
at 0.01 (default in the neuralnet function), 0.025, 0.04, 0.1, 0.3, and 0.4
for the 5 $n$, respectively.

##  4 Discussion

We have proposed a novel tree ensemble method, Selective Cascade of Residual
ExtraTrees (SCORE), for regression prediction. SCORE incorporates a
TreeSelection step that weighs each tree based on its relevance to the outcome
prediction. This step filters out the irrelevant trees, which helps reduce not
only the bias but also the variance for prediction without hurting the
diversity among the ensemble members. Moreover, the weights obtained through
the TreeSelection step can be used to evaluate the importance of predictors
involved in building each tree, on which our new VIM is based. The boosting
step in SCORE helps to further reduce the bias from the potential over-
randomization in ExtraTrees. The experiments show that SCORE offers comparable
or superior performance in prediction and predictor selection compared to
ExtraTrees, RF, and GBM at significantly less computational time, owing to its
robustness to hyper-parameter specifications. Similar observations are
obtained comparing SCORE to NNs, suggesting NNs are not universally better
than other leaning approaches, depending on the specific data as well as
training data sizes. Moreover, compared to other tree-based ensemble methods
[12, 13, 11], SCORE’s ability to measure variable importance as relevant to
the prediction process offers explainability on its prediction, and provides
transparency of and insight into the modeling processes.

Future research may extend the framework of SCORE to classification tasks. The
investigation on the robustness of SCORE to hyper-parameter specifications in
terms of prediction performance also awaits a comprehensive investigation.

## Acknowledgement

We thank Dr. Gitta Lubke and two anonymous referees for their useful and
constructive comments on the project and the manuscript.

## References

  * [1] T. G. Dietterich, “An Experimental comparison of three methods for constructing ensembles of decision trees: bagging, boosting, and randomization,” _Machine Learning_ , vol. 40, pp. 139–157, 2000.
  * [2] L. Breiman, “Bagging predictors,” _Machine Learning_ , vol. 24, no. 2, pp. 123–140, 8 1996.
  * [3] T. K. Ho, “The random subspace method for constructing decision forests,” _IEEE Transactions on Pattern Analysis and Machine Intelligence_ , vol. 20, no. 8, pp. 832–844, 1998.
  * [4] ——, “Random decision forest,” in _Proceedings of the 3rd International Conference on Document Analysis and Recognition_ , 1995, pp. 278–282.
  * [5] P. Geurts, D. Ernst, and L. Wehenkel, “Extremely randomized trees,” _Machine Learning_ , vol. 63, no. 1, pp. 3–42, 2006.
  * [6] W. Fan, H. Wang, P. Yu, and S. Ma, “Is random model better? On its accuracy and efficiency,” in _Third IEEE International Conference on Data Mining_ , no. October 2015. IEEE Comput. Soc, 2003, pp. 51–58.
  * [7] F. T. Liu, K. M. Ting, Y. Yu, and Z.-H. Zhou, “Spectrum of Variable-Random Trees,” _Journal of Artificial Intelligence Research_ , vol. 32, pp. 355–384, 2008.
  * [8] J. H. Friedman, “Greedy function approximation: A gradient boosting machine,” _Annals of Statistics_ , vol. 29, no. 5, pp. 1189–1232, 2001\.
  * [9] L. Mason, J. Baxter, P. Bartlett, and M. Frean, “Boosting algorithms as gradient descent,” in _Advances in Neural Information Processing Systems_ , 2000.
  * [10] G. E. Hinton, S. Osindero, and Y.-W. Teh, “A fast learning algorithm for deep belief nets,” _Neural computation_ , vol. 18, no. 7, pp. 1527–1554, 2006\.
  * [11] Y. Kong and T. Yu, “A Deep Neural Network Model using Random Forest to Extract Feature Representation for Gene Expression Data Classification,” _Scientific Reports_ , vol. 8, 2018.
  * [12] Z. H. Zhou and J. Feng, “Deep forest: Towards an alternative to deep neural networks,” _IJCAI International Joint Conference on Artificial Intelligence_ , pp. 3553–3559, 2017.
  * [13] J. Feng, Y. Yu, and Z. H. Zhou, “Multi-layered gradient boosting decision trees,” _Advances in Neural Information Processing Systems_ , pp. 3551–3561, 2018.
  * [14] Y. Bengio, A. Courville, and P. Vincent, “Representation learning: A review and new perspectives,” _IEEE transactions on pattern analysis and machine intelligence_ , vol. 35, no. 8, pp. 1798–1828, 2013.
  * [15] Y. Freund and R. E. Schapire, “A decision-theoretic generalization of on-line learning and an application to boosting,” _Journal of Computer and System Sciences_ , vol. 55, no. 1, pp. 119–139, 8 1997.
  * [16] R. E. Schapire, “The Strength of Weak Learnability,” _Machine Learning_ , vol. 5, no. 2, p. 197–227, 1990.
  * [17] N. Li and Z. H. Zhou, “Selective ensemble under regularization framework,” in _Lecture Notes in Computer Science (including subseries Lecture Notes in Artificial Intelligence and Lecture Notes in Bioinformatics)_ , 2009.
  * [18] Y. Bengio, A. Courville, and P. Vincent, “Representation learning: A review and new perspectives,” _IEEE Transactions on Pattern Analysis and Machine Intelligence_ , 2013.
  * [19] G. Brown, J. L. Wyatt, and P. Tino, “Managing diversity in regression ensembles,” _Journal of Machine Learning Research_ , vol. 6, pp. 1621–2650, 2005.
  * [20] R. Tibshirani, “Regression shrinkage and selection via the Lasso,” _Journal of the Royal Statistical Society. Series B (Methodological)_ , vol. 58, no. 1, pp. 267–288, 1996.
  * [21] J. Helliwell, R. Layard, and J. Sachs, “World Happiness Report 2019,” Sustainable Development Solutions Network, New York, Tech. Rep., 2019. [Online]. Available: https://worldhappiness.report/ed/2019/
  * [22] J. D. Olden, M. K. Joy, and R. G. Death, “An accurate comparison of methods for quantifying variable importance in artificial neural networks using simulated data,” _Ecological modelling_ , vol. 178, no. 3-4, pp. 389–397, 2004.

## Appendix

Input Attributes in Boston Housing Data  label attribute crim per capita crime
rate by town zn proportion of residential land zoned for lots over 25,000
sq.ft indus proportion of non-retail business acres per town chas tract bounds
Charles River? nox nitric oxides concentration (parts per 10 million) rm
average number of rooms per dwelling age proportion of owner-occupied units
built prior to 1940 dis weighted distances to five Boston employment centres
rad index of accessibility to radial highway tax full-value property-tax rate
per USD $10,000$ ptratio pupil-teacher ratio by town b $1000(B-0.63)^{2}$,
where $B$ is the proportion of blacks by town lstat percentage of lower status
of the population town name of town tract census tract lon longitude of census
tract lat latitude of census tract

Input Attributes in World Happiness Report Data  label attribute yr year
lGPDpc logged GDP per Capita sup social support life healthy life expectancy
at birth free freedom to make life choices genr generosity corr perceptions of
corruption posaf positive affect negaf negative affect govcon confidence in
national government dem democratic quality del delivery quality sdl standard
deviation of ladder by country-year sdml standard deviation/mean of ladder by
country-year gini GINI index (World Bank estimate) gini2 GINI index (World
Bank estimate), average 2000-2015 gini3 GINI of household income reported in
Gallup world poll, by country-year


# Learning Social Cost Functions for
Human-Aware Path Planning

Andrea Eirale
Department of Electronics and Telecommunications
Politecnico di Torino
Torino, TO, 10129, Italy
<EMAIL_ADDRESS>
Matteo Leonetti
Department of Informatics
King’s College London
London, UK
<EMAIL_ADDRESS>
Marcello Chiaberge
Department of Electronics and Telecommunications
Politecnico di Torino
Torino, TO, 10129, Italy
<EMAIL_ADDRESS>

###### Abstract

Achieving social acceptance is one of the main goals of Social Robotic
Navigation. Despite this topic has received increasing interest in recent
years, most of the research has focused on driving the robotic agent along
obstacle-free trajectories, planning around estimates of future human motion
to respect personal distances and optimize navigation. However, social
interactions in everyday life are also dictated by norms that do not strictly
depend on movement, such as when standing at the end of a queue rather than
cutting it. In this paper, we propose a novel method to recognize common
social scenarios and modify a traditional planner’s cost function to adapt to
them. This solution enables the robot to carry out different social navigation
behaviors that would not arise otherwise, maintaining the robustness of
traditional navigation. Our approach allows the robot to learn different
social norms with a single learned model, rather than having different modules
for each task. As a proof of concept, we consider the tasks of queuing and
respect interaction spaces of groups of people talking to one another, but the
method can be extended to other human activities that do not involve motion.

_K_ eywords Human-Aware Navigation $\cdot$ Path Planning $\cdot$ Social
Robotics

##  1 Introduction

Robot navigation, especially indoor, has made tremendous progress, enabling
the first commercial robots to enter a large number of homes as well as public
spaces, such as shopping centers, airports, and restaurants. Most such robots
are based on _geometric_ navigation, treating people as dynamic obstacles. On
the other hand, _social_ navigation aims to achieve socially acceptable robot
behavior—a central goal of sharing spaces with service robots.

Most of the work in social navigation is focused on the direction and velocity
of the people surrounding the robot, to respect personal space and act
naturally in scenarios that require coordination, such as narrow corridors.
Our work is orthogonal to such efforts and considers how robot path planning
should be affected not by people’s velocity but by their activity. Indeed,
even when people do not move at all, such as in a queue, the robot should not
just consider them as static obstacles and cut in front of them.

We propose a learning-based approach that leverages existing widespread
planners. We focus the learning effort on the social aspect, rather than on
areas of navigation in which current planners are very effective, such as
obstacle avoidance. For this reason, contrary to end-to-end learning, we
propose to learn an additional term of the cost function to use in planning,
while retaining the default behavior everywhere else. Planners, to be
efficient, work on a cost map, while social navigation depends on many other
aspects, such as whether two obstacles are, in fact, people, and whether they
are talking to each other or not. Our model takes in input a representation of
the position of the people in the vicinity of the robot and the position of
the goal, and outputs a social cost map in addition to the usual obstacle one.
We present an architecture and a methodology to train such a neural network,
which effectively projects the social navigation problem into classical 2D
navigation.

As a proof of concept, we use two common situations that appear frequently in
public spaces: people queuing and people talking in small groups. Our
methodology is not tied to these particular scenarios and can be extended to
any other setting in which practitioners can provide a dataset. Both scenarios
have been extensively studied separately. We show how the same method can
effectively plan for both scenarios at the same time, with the potential to
incorporate many more. We evaluate the generalization abilities of the
proposed architecture; we analyze the generated plans in a simulated
environment in Gazebo [1], and we demonstrate the resulting behaviors in the
real world on a Pal Robotics TIAGo.

##  2 Related Work

The first examples of designing a behavioral policy able to consider social
factors date back more than twenty years, with the deployment of RHINO [2] and
MINERVA [3] as tour guides in museums. In those cases, the robotic system
perceived people as dynamic, non-responsive obstacles. In the following years,
researchers focused on systems able to distinguish humans from inanimate
objects.

###  2.1 Motion-based Human-Aware Navigation

Many works focused on proxemics, inflating and reshaping the cost function of
people in the scenery, depending on the velocity and the facing direction of
each person [4, 5, 6, 7]. One of the most well-known solutions in this
category is the Robot Operating System (ROS) social navigation layers
111http://wiki.ros.org/social_navigation_layers, which alter the cost map with
a Gaussian distribution around people, with an increased cost in the direction
of motion. Mateus et al. [8] exploit deep learning with hybrid asymmetric
Gaussian functions to represent proxemics for modeling comfort distance with
pedestrians. Truong and Ngo [9] took states (position, orientation, motion,
and hand poses) and social interaction information relative to the robot into
account to model extended personal space and social interaction space using
two-dimensional Gaussian functions.

Social Force Model (SFM) [10] is another popular technique to model human and
robot motions. Ferrer et al. [11] extend the SFM to present a proactive
navigation approach, where the robot can produce actions with a minimum impact
on surrounding pedestrians. Yang et al. [12] employ the SFM on an
omnidirectional robot to achieve natural and human-like movements. More recent
works consider autonomous agents able to plan an optimal trajectory, avoiding
collisions with people and obstacles, by planning around estimates of future
human motion. These predictions are often achieved with deep neural networks,
like generative adversarial networks [13], convolutional neural networks [14,
15], and attention transformers [16].

Intention-aware navigation is the direct evolution of the previous methods. In
this case, by exploiting behavior prediction models, it is possible not only
to predict the future movements of humans but also to estimate their final
goals. Bai et al. [17] model the uncertainty of human intent in the Partially
Observable Markov Decision Process (POMDP) framework. Navigation in dense
human crowds is another interesting task, where solutions focused on avoiding
blockage due to human activities [18] and cooperating with people through
interacting Gaussian processes [19]. Mavrogiannis et al. employ geometric [20]
and topological invariant [21] representations to model the coupling among
trajectories of multiple navigating agents.

Deep Reinforcement Learning has been used in many works for prediction in
crowd navigation domains. Chen et al. [22] apply CADRL, a deep reinforcement
learning framework for socially aware multi-agent collision avoidance. Everett
et al. [23] exploit an actor-critic variant to relax prior assumptions and
learn policies and agent motion models simultaneously. Furthermore, Tai et al.
[24] train a generative adversarial imitation learning model on a dataset
generated using the social force model [10]. Finally, Chen et al. [25] use
attention-based reinforcement learning to produce interaction-aware collision
avoidance behaviors.

The works cited above, as well as the majority of scientific literature on
social navigation, focus on people’s movements: their trajectory, intention,
and destination. However, several social norms are independent of movement but
dependent on people’s activity. For instance, walking between two people
talking to each other is considered rude, even if they are just standing
still. Such scenarios are currently under-represented in the scientific
literature and form the basis of this work.

###  2.2 Task-specific Social Navigation

Other works have been more focused on solving specific social tasks during
navigation. Xiao et al. [26] train a Performer architecture [27] in an
imitation-learning fashion to learn a cost function and improve an MPC
controller. This system is then used to solve different social tasks, such as
moving around corners and respecting comfort distance in a human-aware manner.
However, unlike our approach, Performer-MPC requires the cost function to be
learned individually for each navigation scenario.

Within the queue following, one of the social navigation scenarios we consider
in this paper, Nakauchi et al. [28] designed a dedicated pipeline. Their work
focuses on how a line of people can be defined and then perceived by the
autonomous agent. The navigation system generates a series of goals (depending
on the number of people in the queue) until the back of the line is reached,
and does not provide full path planning. In [29, 30], Banisetty et al. exploit
simple geometric reasoning to compute a social goal at the beginning of the
line. The true goal of the navigation, placed at the end of the line, is then
replaced with the social goal. However, this system only allows the robot to
join the queue and is unable to follow the line until the robot reaches the
final goal.

Many works focused on detecting groups of people, based on the definition of
F-formations provided by Kendon [31]. Cristiani et al. [32] and Setti et al.
[33] proposed Hough Voting with the people’s positions and head orientations
to detect F-formations, while Hedayati et al. [34, 35] treated the problem as
a binary classification. A few works focused on unsupervised group detection
in dynamic, egocentric views [36, 37], while Schmuck et al. [38, 39] proposed
GROup detection With Link prediction (GROWL), demonstrating the effectiveness
of Graph Neural Network in the detection of people interaction. While our work
does not strictly focus on recognizing the complex connections within groups
of people, our experimentation shows that, for social navigation purposes,
implicit knowledge of the interaction between people arising from simple
spatial information is often sufficient to recognize and avoid sparse groups
of people. However, group detection could be integrated within our methodology
to further improve the navigation, especially through dense crowds.

Differently from all the work described above, our approach is targeted at all
scenarios in which people are static, and yet their position or activity
dictates socially acceptable robot paths.

##  3 Background

Grid-based path planners, such as Dijkstra’s or A*, work with a discrete
representation of the environment, which can be formalized as a two-
dimensional, undirected graph $M=(V,E)$, usually referred to as a _grid map_.
In this graph, $V$ is a vertex set of size $n\times m$, where each vertex
corresponds to a cell in the map, and $E$ is an edge set of size $k$. Edges
connect adjacent map cells. A path is a sequence of vertices
$P=(v_{1},v_{2},...,v_{n})$ such that there exists an edge $e_{i,i+1}\in E$
connecting $v_{i}$ and $v_{i+1}$, with an associated distance $d_{i,i+1}$ for
$1\leq i<n-1$.

Given a starting vertex $v_{s}\in V$ and a goal $v_{g}\in V$, a shortest
feasible path from $v_{s}$ to $v_{g}$ is defined as
$P^{*}=(v_{1},v_{2},...,v_{n})$, where $v_{1}=v_{s}$ and $v_{n}=v_{g}$. This
path is a minimizer of:

$P^{*}=arg\,min_{v1,...,vn\in V}\sum_{i=1}^{n-1}{f_{c}(v_{i},v_{i+1})},$ (1)

where $f_{c}:V\times V\rightarrow\mathbb{R^{+}}$ is a cost function.

Cost functions can be defined with several terms, encoding both the distance
between two nodes and the desirability of each node. A _cost map_ $C$
corresponding to a grid map $M$ is an $n\times m$ matrix, with a cost
$c_{i}\in\mathcal{R^{+}}$ for each element $v_{i}\in V$, representing the cost
of entering the node. Usually, a priori information about the environment,
like walls and other static obstacles, is gathered in a global cost map
$C_{G}$. Sensor data, containing the detection of obstacles around the agent,
are instead included in a local cost map $C_{L}$ centered on the robot’s
reference frame. The total cost of a transition is, therefore:

$f_{c}(v_{i},v_{i+1})=d_{i,i+1}+c_{G,i+1}+c_{L,i+1}\;,$ (2)

where $c_{G,i+1}$ and $c_{L,i+1}$ are the costs of $v_{i+1}$ respectively in
$C_{G}$ and $C_{L}$.

Figure 1: Representation of the pipeline used to integrate the social costs
into the traditional navigation system cost function. The social grid map
$M_{S}$ is obtained from the position of the goal and people. The planner will
compute the shortest path towards the goal considering the social cost map
$C_{S}$ in addition to all the other obstacle cost maps.

##  4 Methodology

We want to retain all the desirable properties of grid-based navigation
systems, such as path optimality and obstacle avoidance, while introducing
socially acceptable behaviors. All social aspects that require the robot to
avoid a particular node in the map can be encoded in an additional cost map.
We learn the social cost function $f_{S}:M_{S}\rightarrow C_{S}$ through a
deep neural network, able to generate the social cost map $C_{S}$ encoding
_social obstacles_.

Figure 1 shows the pipeline of the methodology. Positions of the goal and
people are gathered into a social local grid map $M_{S}$, centered on the
robot’s reference frame. A learned encoder-decoder produces the social cost
map, which is added to the other navigation costs to produce the total cost
map. This process repeats continually, as the social cost map is regularly
updated to reflect the latest position of the people in the scene.

The total cost expressed in Equation 2 can then be rewritten as:

$f_{c}(v_{i},v_{i+1})=d_{i,i+1}+c_{G,i+1}+c_{L,i+1}+c_{S,i+1}\;,$ (3)

where $c_{S,i+1}$ is the costs added by the social cost map $C_{S}$.

In the rest of this section, we will describe the network’s architecture of
the social cost function $f_{S}$, how the training dataset is created, and how
the network is deployed in a real system.

###  4.1 Network Architecture

Figure 2: Representation of the Encoder-Decoder neural network used in this
work.

The network is an encoder-decoder consisting of convolutional, maxpooling, and
upscaling layers with $342,049$ (learnable) parameters. The network structure
is shown in Figure 2. The input discretizes a 24x24m area with a resolution of
0.2m, represented in a grid map with 120x120 cells. These measures are
compatible with a typical indoor Lidar sensor, usually in the 8–18m range.

###  4.2 Dataset creation

We generate the dataset with a series of RGB images for each human activity we
want the system to react to. The dataset consists of pairs $\\{(M_{S}^{(i)}$,
$C_{S}^{(i)})\\}$ of grid map and corresponding cost map.

The maps for a particular human activity only have elements pertaining to that
activity, excluding other elements, such as obstacles. The label image, that
is, the intended output of the network, is a grey-scale cost map with the same
dimensions as the input image, representing only the costs corresponding to
the input map.

As a proof of concept, in this work, we consider the social scenarios of a
queue of people and that of small groups of people talking. In principle, any
scenario that only requires the robot to avoid certain areas can be added.

####  4.2.1 Queuing

In this scenario, we expect the robot to follow the queue to the goal without
cutting the line. We create grid maps with a line of people behind a goal,
represented with two different colors, so that the network can distinguish the
people from the goal. Each sample we generate is randomized through a series
of varying parameters: the number of people in the queue (which also
influences the total length of the queue), the distance between two people,
the distance between the first person and the goal, and the deviation of each
person from the center axis of the queue. The label image represents a
U-shaped obstacle surrounding the queue. The thickness of the U-shaped
obstacle walls should be tuned according to the robot’s footprint, in order
for the planner to plan around the queue and leave the only available path
toward the goal behind the queue (cf. Figure 1). This “social corridor"
towards the goal should be narrow enough to force the robot to wait behind the
last person.

####  4.2.2 Groups of People

In this scenario, we expect the robot to go around groups of people rather
than through them. The input grid map contains the small groups of people. In
this case, the randomized parameters are: the number of people for each group,
their position in the map, and the number of groups. The label image consists
of a variable number of virtual obstacles (depending on the position of the
people), filling the gap between them.

####  4.2.3 Further refinements

Each training image contains at least one instance for each considered social
scenario, positioned randomly within the grid map space. To improve
generalization, we also diversify the dataset according to the following
criteria. In the queuing scenario, with a small probability, the training
image does not include the queue. Instead, it has the goal in a random
position, with the label image not containing any virtual obstacle. This
variation prevents the network from associating the social obstacle with the
goal alone. Similarly, to prevent the network from associating the obstacle
with each person instead of a group of people, a certain number of isolated
people have a probability of spawning in the image. Each generated sample pair
is oriented differently to obtain rotation invariance in the prediction.

Every sample pair of grid map and label is created programmatically with an
automatic generator, which may raise the question: Can such a generator be
used online—do we need the neural network at all? The generator for each
scenario can produce pairs $(M_{S}^{(i)},C_{S}^{(i)})$, but it is not
equivalent to the function $f_{S}$ that translates a general social grid map
into a cost map. The generator first creates the social obstacles in
randomized locations and then places people and the goal. This is easier than
the inverse operation, and works on a single scenario. For this reason,
generators are used only to produce samples for the function $f_{S}$, which is
then learned to classify a given grid map into its corresponding cost map and
to do so with possibly different scenarios simultaneously.

###  4.3 Deployment

We deployed the method in simulation and on a real robot by integrating it
into a ROS2 [40] node. Periodically, this node retrieves the position of the
goal and people and converts this information into the input grid map image
representation $M_{S}$ used by the network. The generated image is then
provided to the neural network, and the output prediction of the social cost
map $C_{S}$ is published on an OccupancyGrid ROS2 topic. This map is used as a
local cost map by the Nav2 [41] planner, in addition to all the other
traditional cost maps provided to the navigation system. From Nav2 parameters,
we select a wavefront Dijkstra expanded holonomic planner in combination with
a Model Predictive Path Integral Controller (MPPI). In order to wait in line
behind each person, we use a simple Nav2 behavior tree with Wait nodes in
combination with a short obstacle maximum detection range, which enables the
robot to reach the queue and then wait for the following person to proceed,
freeing the computed path, before advancing.

##  5 Experimental Validation

In this section, we first evaluate the generalization of the learned social
cost function, and then demonstrate the resulting behaviors in a number of
scenarios both in simulation and on a real robot.

###  5.1 Learning Evaluation

The network was trained with a generated dataset containing $5\times 10^{5}$
samples for 5 epochs on a machine with an Intel Core i7-9700K CPU and an
Nvidia GeForce RTX 2080 Ti, taking approximately 30 minutes. Generated cost
maps may differ from the expected test label in a few cells but still produce
a socially correct path. For this reason, instead of the accuracy of the cost
map on the test set, we measure the success rate of the generated cost maps.
We consider a cost map successful if the path induced by the map is correct;
that is, it follows the queue appropriately and does not navigate across
groups.

We produced three different test sets, containing respectively images only
with queues, only with groups, and with queues and groups mixed. Table 1 shows
the parameters for training, testing, and success rates. In every dataset we
generated for training and testing, where the queue appears, the distance
between two people in the line and from the first person and the goal is in
the range 0.5-1.5 meters, while the deviation of each person from the center
axis of the queue is in the range 0-0.5 meters from each side 222The source
code we used to generate each dataset and to train the network is available at
https://github.com/PIC4SeR/SocialCostFunction. In all cases the success rate
is above $95\%$, which we consider sufficiently reliable.

Table 1: Parameters used for training and test datasets. The "Success" column represents the percentage of samples in which the network could correctly recognize and delimit each social scenario. Dataset |  Samples |  People in queue |  People in group |  Number of groups |  Success
---|---|---|---|---|---
Train |  5e5 |  2 to 5 |  2 |  1 |  n/a
Test queue only |  400 |  2 to 15 |  n/a |  n/a |  95.75%
Test groups only |  300 |  n/a |  2 to 7 |  1 to 5 |  98.67%
Test queue and groups |  300 |  2 to 9 |  2 to 4 |  1 to 3 |  97.33%

(a)

(b)

(c)

Figure 3: Simulation of navigation in a demo environment with three different
groups of people interacting with each other and one queue, visualized in
$(a)$ with Gazebo. Magenta patches in the cost map represent areas with high
costs corresponding to obstacles, while the cyan-to-blue gradient is a
decreasing cost that inflates the obstacles. The green line is the computed
path towards the goal, positioned at the end of the queue. In $(b)$, the
social cost map is included, allowing the navigation algorithm to plan a path
that drives the robot at the beginning of the queue without crossing any group
of people. On the contrary, in $(c)$ the social cost map is not considered;
the robot cuts the line directly to the goal and passes between a group of
people talking.

(a)

(b)

(c)

Figure 4: Simulation of the robot approaching the queue from different
directions: from the front $(a)$, the side $(b)$, and the back $(c)$ of the
line of people. As can be seen, with the social cost map, the planner is
always able to set the robot in the queue correctly.

(a)

(b)

Figure 5: Simulation of the robot moving between multiple groups of people
$(a)$ and approaching a long queue of people $(b)$. The network can recognize
and delimit the whole queue and each group, assigning high costs to the
interaction areas between the people. Figure 6: Tests carried out on the robot
in a real environment. With the social cost map enabled, the robot can
approach the queue from the beginning, stay behind the last person as the line
proceeds, and wait for its turn until the goal at the end of the queue is
reached.

###  5.2 Simulated Scenarios

We now demonstrate the system in several practical scenarios in simulation,
comparing the generated plan with and without the learned social costs. We
built Gazebo [1] simulation environments containing a queue of people, a
variable number of groups of people interacting, and other people passing by
(Figure 3a). The robot mounts a lidar sensor to detect physical obstacles
around it, while the information regarding the position of people is obtained
directly from the simulator.

In Figure 3, the robot is given a goal in front of a queue. The network can
recognize and delimit the queue as well as the area between each group of
people talking. Given the social obstacles, the planner computes a path around
the groups and the line of people, queuing correctly behind the last person in
line. In comparison, without the social layer, the robot plans between the
people talking and directly to the goal, cutting the line (Figure 3c).

We further demonstrate the queue task by approaching the line of people from
different angles. The robot is positioned in front, on the side, and behind
different queues containing a variable number of people. As seen from Figure
4, the network can correctly recognize the queue regardless of where the robot
comes from. Furthermore, the provided social cost map allows the planner to
compute an appropriate path for the robot to queue correctly, position behind
the last person, and follow the line of people until the goal is reached.
Although the maximum number of people appearing in training set queues is 5,
the network can correctly recognize and delimit a queue with any number of
people, as long as the whole line is sufficiently aligned and contained within
the input image (Figure 5$(b)$).

We further demonstrate the talking groups scenario in Figure 5$(a)$. The
network can correctly identify different groups formed by a variable number of
people. The network’s virtual obstacles often have irregular shapes, but they
are sufficient to prevent the robot from crossing the interaction zone.

A series of limitations of the network are dictated by the dataset used for
its training. For example, the network cannot correctly recognize the queue if
it is not settled along a straight line, since our dataset only contained
straight queues. Our tests demonstrated that as long as each person’s
deviation from the center of the image does not exceed 1m (from both sides),
the network can classify the queue correctly, but it often fails otherwise.
Varying the maximum distance between each person in the queue could also lead
to instability of the recognition. We noted that recognition starts failing
when the mutual distance between people in the line is around 3m. Furthermore,
while groups of people not in a queue but near one do not usually interfere,
they may be interpreted as a second queue if placed near the goal. Similarly,
if a goal is given too close to a particularly aligned group, it could be
misinterpreted as a line of people.

Other incorrect results could arise in the presence of a crowd of people. The
network can correctly recognize multiple groups of people and set the
appropriate costs for each interaction space, even when such groups are close
to each other. However, if the distances between people of different groups
are smaller than the distances of people within a group, the network can
produce incorrect virtual obstacles. Further refinement of the dataset can
overcome some of these limitations, but if the environment is extremely
crowded, other strategies that require the robot to ask people to move become
necessary.

###  5.3 Real-robot deployment

The same system used in the simulation environment was deployed on a PAL
Robotics TIAGo robot. Similarly to the simulation, the Lidar sensor detects
obstacles around the robot. The information regarding the number and position
of people in the scene are extracted from YoloV8 [42], using an RGB-D camera
mounted on TIAGo’s head. We recreated the tests made in simulation in the real
environment (one example is shown in Figure 6).

For the queue task, once the people in line were correctly recognized and
localized, the robot was given different starting points: in front, on the
side, and behind the queue. Also in the real-life setting, the network
recognized and delimited the queue correctly. The robot entered the line of
people and followed it until the goal was reached, independently from the
direction it came from.

We recreated a scenario with multiple groups with different numbers of people.
A goal was provided at the end of the map. The network delimited the
interaction areas, and the robot reached the goal without crossing any group
of people.

The reliability of the perception and tracking system may limit the
correctness of the computed social cost. In the queuing task, when the robot
enters the line of people, the last person in the queue can obstruct the
vision of the people behind. For this reason, the social grid map cannot be
just generated from perceptions each time. To mitigate this problem, we keep
track of every person seen in the past time instants, similarly to how the
navigation stack remembers obstacles that are not presently in view. We place
a marker in the last position where a person was detected and remove it when
the area is in view and the laser sees no obstacle.

##  6 Conclusions

This paper proposes a novel method to enhance Social Robotic Navigation by
modifying the planner cost function. Differently from past literature, we do
not focus on proxemics or estimates of people’s movement based on kinematic
information. Instead, we consider areas the robot should plan to avoid based
on the position of people in the environment, in accordance to specific social
norms. We develop a proof of concept based on two common social tasks: queuing
in a line of people and respecting interaction spaces of people talking with
each other. We demonstrate the effectiveness of our solution both in
simulation and in a real-world scenario, comparing it with traditional goal-
based navigation. We focused on learning and navigation, and used off-the-
shelf libraries for computer vision. The social maps could be further refined
if the vision module was able to detect more than people’s positions, for
instance, whether or not they were talking, and with whom.

Adapting the robot’s behavior to what people are doing in the environment has
not received as much attention as how people are moving. Some social norms
affect navigation despite people not moving at all. We expect this line of
research to receive more attention within social navigation in the near
future.

## References

  * [1] Design and use paradigms for gazebo, an open-source multi-robot simulator.  In 2004 IEEE/RSJ international conference on intelligent robots and systems (IROS)(IEEE Cat. No. 04CH37566), volume 3, pages 2149–2154. IEEE, 2004.
  * [2] Wolfram Burgard, Armin B Cremers, Dieter Fox, Dirk Hähnel, Gerhard Lakemeyer, Dirk Schulz, Walter Steiner, and Sebastian Thrun.  The museum tour-guide robot rhino.  In Autonome Mobile Systeme 1998: 14. Fachgespräch Karlsruhe, 30. November–1. Dezember 1998, pages 245–254. Springer, 1999.
  * [3] Sebastian Thrun, Michael Beetz, Maren Bennewitz, Wolfram Burgard, Armin B Cremers, Frank Dellaert, Dieter Fox, Dirk Haehnel, Chuck Rosenberg, Nicholas Roy, et al.  Probabilistic algorithms and the interactive museum tour-guide robot minerva.  The international journal of robotics research, 19(11):972–999, 2000.
  * [4] Rachel Kirby.  Social robot navigation.  Carnegie Mellon University, 2010.
  * [5] Marina Kollmitz, Kaijen Hsiao, Johannes Gaa, and Wolfram Burgard.  Time dependent planning on a layered social cost map for human-aware robot navigation.  In 2015 European Conference on Mobile Robots (ECMR), pages 1–6. IEEE, 2015.
  * [6] Leonardo Scandolo and Thierry Fraichard.  An anthropomorphic navigation scheme for dynamic scenarios.  In 2011 IEEE International Conference on Robotics and Automation, pages 809–814. IEEE, 2011.
  * [7] Fang Fang, Manxiang Shi, Kun Qian, Bo Zhou, and Yahui Gan.  A human-aware navigation method for social robot based on multi-layer cost map.  International Journal of Intelligent Robotics and Applications, 4:308–318, 2020.
  * [8] Andre Mateus, David Ribeiro, Pedro Miraldo, and Jacinto C Nascimento.  Efficient and robust pedestrian detection using deep learning for human-aware navigation.  Robotics and Autonomous Systems, 113:23–37, 2019.
  * [9] Xuan-Tung Truong and Trung-Dung Ngo.  Dynamic social zone based mobile robot navigation for human comfortable safety in social environments.  International Journal of Social Robotics, 8:663–684, 2016.
  * [10] Dirk Helbing and Peter Molnar.  Social force model for pedestrian dynamics.  Physical review E, 51(5):4282, 1995.
  * [11] Gonzalo Ferrer and Alberto Sanfeliu.  Proactive kinodynamic planning using the extended social force model and human motion prediction in urban environments.  In 2014 IEEE/RSJ International Conference on Intelligent Robots and Systems, pages 1730–1735. IEEE, 2014.
  * [12] Chun-Tang Yang, Tianshi Zhang, Li-Pu Chen, and Li-Chen Fu.  Socially-aware navigation of omnidirectional mobile robot with extended social force model in multi-human environment.  In 2019 IEEE International Conference on Systems, Man and Cybernetics (SMC), pages 1963–1968. IEEE, 2019.
  * [13] Agrim Gupta, Justin Johnson, Li Fei-Fei, Silvio Savarese, and Alexandre Alahi.  Social gan: Socially acceptable trajectories with generative adversarial networks.  In Proceedings of the IEEE conference on computer vision and pattern recognition, pages 2255–2264, 2018.
  * [14] Abduallah Mohamed, Kun Qian, Mohamed Elhoseiny, and Christian Claudel.  Social-stgcnn: A social spatio-temporal graph convolutional neural network for human trajectory prediction.  In Proceedings of the IEEE/CVF conference on computer vision and pattern recognition, pages 14424–14432, 2020.
  * [15] Dapeng Zhao and Jean Oh.  Noticing motion patterns: A temporal cnn with a novel convolution operator for human trajectory prediction.  IEEE Robotics and Automation Letters, 6(2):628–634, 2020.
  * [16] Anirudh Vemula, Katharina Muelling, and Jean Oh.  Social attention: Modeling attention in human crowds.  In 2018 IEEE international Conference on Robotics and Automation (ICRA), pages 4601–4607. IEEE, 2018.
  * [17] Haoyu Bai, Shaojun Cai, Nan Ye, David Hsu, and Wee Sun Lee.  Intention-aware online pomdp planning for autonomous driving in a crowd.  In 2015 ieee international conference on robotics and automation (icra), pages 454–460. IEEE, 2015.
  * [18] Chonhyon Park, Jan Ondřej, Max Gilbert, Kyle Freeman, and Carol O’Sullivan.  Hi robot: Human intention-aware robot planning for safe and efficient navigation in crowds.  In 2016 IEEE/RSJ International Conference on Intelligent Robots and Systems (IROS), pages 3320–3326. IEEE, 2016.
  * [19] Pete Trautman, Jeremy Ma, Richard M Murray, and Andreas Krause.  Robot navigation in dense human crowds: Statistical models and experimental studies of human–robot cooperation.  The International Journal of Robotics Research, 34(3):335–356, 2015.
  * [20] Christoforos I Mavrogiannis and Ross A Knepper.  Multi-agent path topology in support of socially competent navigation planning.  The International Journal of Robotics Research, 38(2-3):338–356, 2019.
  * [21] Christoforos I Mavrogiannis and Ross A Knepper.  Multi-agent trajectory prediction and generation with topological invariants enforced by hamiltonian dynamics.  In Algorithmic Foundations of Robotics XIII: Proceedings of the 13th Workshop on the Algorithmic Foundations of Robotics 13, pages 744–761. Springer, 2020.
  * [22] Yu Fan Chen, Miao Liu, Michael Everett, and Jonathan P How.  Decentralized non-communicating multiagent collision avoidance with deep reinforcement learning.  In 2017 IEEE international conference on robotics and automation (ICRA), pages 285–292. IEEE, 2017.
  * [23] Michael Everett, Yu Fan Chen, and Jonathan P How.  Motion planning among dynamic, decision-making agents with deep reinforcement learning.  In 2018 IEEE/RSJ International Conference on Intelligent Robots and Systems (IROS), pages 3052–3059. IEEE, 2018.
  * [24] Lei Tai, Jingwei Zhang, Ming Liu, and Wolfram Burgard.  Socially compliant navigation through raw depth inputs with generative adversarial imitation learning.  In 2018 IEEE international conference on robotics and automation (ICRA), pages 1111–1117. IEEE, 2018.
  * [25] Changan Chen, Yuejiang Liu, Sven Kreiss, and Alexandre Alahi.  Crowd-robot interaction: Crowd-aware robot navigation with attention-based deep reinforcement learning.  In 2019 international conference on robotics and automation (ICRA), pages 6015–6022. IEEE, 2019.
  * [26] Xuesu Xiao, Tingnan Zhang, Krzysztof Marcin Choromanski, Tsang-Wei Edward Lee, Anthony Francis, Jake Varley, Stephen Tu, Sumeet Singh, Peng Xu, Fei Xia, Sven Mikael Persson, Dmitry Kalashnikov, Leila Takayama, Roy Frostig, Jie Tan, Carolina Parada, and Vikas Sindhwani.  Learning model predictive controllers with real-time attention for real-world navigation.  In Karen Liu, Dana Kulic, and Jeff Ichnowski, editors, Proceedings of The 6th Conference on Robot Learning, volume 205 of Proceedings of Machine Learning Research, pages 1708–1721. PMLR, 14–18 Dec 2023.
  * [27] Krzysztof Choromanski, Valerii Likhosherstov, David Dohan, Xingyou Song, Andreea Gane, Tamas Sarlos, Peter Hawkins, Jared Davis, Afroz Mohiuddin, Lukasz Kaiser, et al.  Rethinking attention with performers.  arXiv preprint arXiv:2009.14794, 2020.
  * [28] Yasushi Nakauchi and Reid Simmons.  A social robot that stands in line.  Autonomous Robots, 12:313–324, 2002.
  * [29] Santosh Balajee Banisetty, Scott Forer, Logan Yliniemi, Monica Nicolescu, and David Feil-Seifer.  Socially aware navigation: A non-linear multi-objective optimization approach.  11(2), jul 2021.
  * [30] Santosh Balajee Banisetty, Vineeth Rajamohan, Fausto Vega, and David Feil-Seifer.  A deep learning approach to multi-context socially-aware navigation.  In 2021 30th IEEE International Conference on Robot & Human Interactive Communication (RO-MAN), pages 23–30, 2021.
  * [31] Adam Kendon.  Conducting interaction: Patterns of behavior in focused encounters, volume 7.  CUP Archive, 1990.
  * [32] Marco Cristani, Loris Bazzani, Giulia Paggetti, Andrea Fossati, Diego Tosato, Alessio Del Bue, Gloria Menegaz, and Vittorio Murino.  Social interaction discovery by statistical analysis of f-formations.  In BMVC, volume 2, pages 10–5244, 2011.
  * [33] Francesco Setti, Chris Russell, Chiara Bassetti, and Marco Cristani.  F-formation detection: Individuating free-standing conversational groups in images.  PloS one, 10(5):e0123783, 2015.
  * [34] Hooman Hedayati, Annika Muehlbradt, Daniel J Szafir, and Sean Andrist.  Reform: Recognizing f-formations for social robots.  In 2020 IEEE/RSJ International Conference on Intelligent Robots and Systems (IROS), pages 11181–11188. IEEE, 2020.
  * [35] Hooman Hedayati, Daniel Szafir, and Sean Andrist.  Recognizing f-formations in the open world.  In 2019 14th ACM/IEEE International Conference on Human-Robot Interaction (HRI), pages 558–559. IEEE, 2019.
  * [36] Angelique Taylor, Darren M Chan, and Laurel D Riek.  Robot-centric perception of human groups.  ACM Transactions on Human-Robot Interaction (THRI), 9(3):1–21, 2020.
  * [37] Viktor Schmuck and Oya Celiktutan.  Rica: Robocentric indoor crowd analysis dataset.  IMU, 127(74,234):31–172, 2020.
  * [38] Viktor Schmuck and Oya Celiktutan.  Growl: Group detection with link prediction.  In 2021 16th IEEE International Conference on Automatic Face and Gesture Recognition (FG 2021), pages 1–8. IEEE, 2021.
  * [39] Viktor Schmuck and Oya Celiktutan.  igrowl: Improved group detection with link prediction.  IEEE Transactions on Biometrics, Behavior, and Identity Science, 2022.
  * [40] Steven Macenski, Tully Foote, Brian Gerkey, Chris Lalancette, and William Woodall.  Robot operating system 2: Design, architecture, and uses in the wild.  Science Robotics, 7(66):eabm6074, 2022.
  * [41] Steven Macenski, Francisco Martin, Ruffin White, and Jonatan Ginés Clavero.  The marathon 2: A navigation system.  In 2020 IEEE/RSJ International Conference on Intelligent Robots and Systems (IROS), 2020.
  * [42] Glenn Jocher, Ayush Chaurasia, and Jing Qiu.  Ultralytics YOLO, January 2023.