# Graph Convolutional Neural Networks for analysis of EEG signals, BCI
application

Mirfarid Musavian Ghazani
Skolkovo Institute of Science and Technology
(SKOLTECH), CDISE, Moscow, Russia
<EMAIL_ADDRESS>
&Anh Huy Phan
Skolkovo Institute of Science and Technology
(SKOLTECH), CDISE, Moscow, Russia,

###### Abstract

Decoding brain signals has gained many attention and has found much
applications in recent years such as Brain Computer Interfaces, communicating
with controlling external devices using the user’s intentions, occupies an
emerging field with the potential of changing the world, with diverse
applications from rehabilitation to human augmentation. This being said brain
signal analysis, EEG brain signal analysis in particular, is a challenging
task. With the advances and achievements in the field of deep learning in
problem solving with using only raw data, few attempts has been carried in
recent years, to apply deep learning to tackle EEG among other types of brain
signals. In this study, we propose a novel loss function, called DeepCSP to
extend the classical Common Spatial Patterns to a non linear, differentiable
module to serve as the loss function to enforce linearly separable latent
representations of EEG signals belonging to different classes in an end to end
manner on raw signals without the need to perform extensive feature
engineering. With recent generalizations of deep learning methods to work on
arbitrarily structured graphs and the introduced loss we have proposed two
light weight models to decode EEG signals and carried experiments to show
their performance.

_K_ eywords DeepCSP $\cdot$ Graph convolution $\cdot$ Deep learning $\cdot$
Motor Imaginary

##  1 Introduction

A Brain Computer Interfaces (BCI) is "a combination of hardware and software
that allows brain activities to control external devices or even computers".
Decoding and classification of mental intents using the data gathered from the
brain activities recorded from subjects in a very carefully designed
experiment setting is the main goal in these applications.

Among the various data devices and technologies used in recording brain
activities, EEG recording modalities are of special interest. In the case of
this kind of data, two main parameters are spatial resolution, referring to
the spatial proximity of the recorded data point and temporal resolution,
referring to the sampling frequency. The high temporal resolution of EEG
recordings combined with relatively much lower costs compared to other methods
and the fact that it is non invasive, has been the most popular tool in BCI
applications and therefore is chosen as the focus of this study. There are
variety of applications for BCIs, from control of robotic arms and various
other prosthetic devices [1] to speech synthesizers [2] and other
communication devices [3]. In general.

It can be argued that however attempting to decode the brain signals obtained
among the various experiment settings faces lots of challenges, decoding Motor
Imagery (MI) signals is relatively more troublesome.

BCI’s can be defined as "a combination of hardware and software that allows
brain activities to control external devices or even computers"[4]. It can be
argued that one of the most challenging and crucial parts of these interfaces,
is designing the experiments to record data from brain. Motor Imagery (MI) is
defined as "mental rehearsal of movements without the movements being
executed" [5]. An example of MI can be the intent of hand movement versus the
physical act of hand movement. When movement is imagined, the activity in the
brain is in the motor cortex, like the real movement.

Despite vast and successful applications and developments in deep learning in
recent years and their generalizations to Graph Convolutional Neural Networks,
the amount of studies that have focused on the BCI Motor Imagery applications
are limited.

In this study, we develop a new loss function, DeepCSP, inspired by the
ubiquitous Common Spatial Patterns method used in BCI application. The
applicability of non linear extension of this method and the effect of
building a differentiable module of one of the most commonly used signal
processing, on the adjustment and optimizing the parameters of the deep
convolutional architectures to get separable latent features of EEG signals
can be power full addition to end to end BCI applications.

By representing EEG signals as graphs we investigate the positive effect of
the features learned from GCNs, powerful paradigm linking deep learning to
analysis and classification of various types of graphs, on the performance of
the Motor Imagery EEG signals classification.

##  2 Literature Review

In BCI systems, machine learning has a major role in the classification of EEG
signals. Supervised algorithms such as SVM (with different kernels) [6, 7],
LDA [8, 9] and Bayesian methods [7] have been among the most popular
algorithms, mainly due to their simplicity and fast training [2]. One main
drawback of these methods is that they require data preprocessing and manual
feature extraction in order to be able to learn the data.

Most of the published techniques using the mentioned class of machine learning
algorithms use some form of feature extraction followed by classification with
a relatively simple model. This fact causes these techniques to be fast and
relatively less likely to suffer from overfitting, in the meantime their
accuracy and efficiency heavily depend on the quality of the feature
extraction methods used. They function by transforming the training data to a
space, mainly such that the different classes are more easily distinguishable
and then find a hyperplane in that space that tries to separate the classes of
data. The Linear Discriminant Analysis (LDA) is one of the commonly used
classifiers of this group of approaches which assumes the normality and
homoscedasticity (meaning that the classes have the same variance) to create a
linear decision surface [10]. One of the other commonly used classifiers is
Naive Bayes, a simple probabilistic classifier which uses the concept of
maximum likelihood to make decisions [11], and K Nearest Neighbors, which
queries the training data using a distance metric to find the K samples
closest to the one being classified to determine the label [2]. The most
common feature extraction methods include: Fourier transforms and its
variants, wavelet transforms, principal component analysis (PCA), independent
component analysis (ICA), autoregressive methods, or combinations of those
techniques. The Fourier transform is likely the most common feature extraction
method. It provides information on the frequency spectrum of the signals,
extracting the spectral information but failing to capture the temporal
information of the signals. One of the methods used to overcome this
shortcoming heavily deployed in the literature is to use the short time
Fourier transforms. The Fourier transform is used to calculate features, such
as power spectra, which match nicely with EEG literature on brain activities
that occur within different frequency bands [12, 2].

As mentioned before new approaches to tackle the challenge decoding EEG
signals, and minimizing the need to perform feature engineering steps on the
signals, some of the recently published studies have attempted to use deep
learning architectures, while maintaining the competitive performance levels.

One of the early introduced deep learning architectures presenting successful
results at classifying Steady state visually evoked potential, SSVEP signals,
used series of Spatial convolutional layers to mix the input EEG signals,
followed by 1D temporal convolutional layers [13]. before feeding the
extracted features from these two layers, the features have been switched into
the frequency domain using fast Fourier transform.

In another study both the raw form of the EEG signals and extracted spectral
features have been used as input for doing classification on the dataset of
recorded EEG signals of subjects, while imagining fragments from a span of
music, selected from different [14]. The results showed the capacity of CNNs
to classify imagery based EEG. In a yet another study, the same authors
employed a convolutional autoencoder (CAE) to train a CNN on the same dataset,
using cross trial similarity constraint encoding in order to learn
individually adapted components that are shared in between the subjects, in
order to capture invariance between trials within and across subjects [15].
The use of autoencoders in this study is one of the earliest successful
attempts to use deep learning architectures to do data augmentation to tackle
the issue of the low number of samples of EEG data.

A novel representation of EEG signals, where the channels of the EEG signals
were transformed into two dimensional space and then by using Fourier
transform across several time bins, the data was reduced to fewer frames [16].
By doing so, a sequence of images is generated for the whole EEG trial and
then fed into a combined CNN and Long Short Term Memory (LSTM), and then a
fully connected layer was used for classification. This method was able to
outperform the classical methods by considerable margins.

In [17] use of fully convolutional neural networks was explored by Lawhern et
al, on a variety of BCI paradigms, including P300 visual evoked potentials,
error related negativity responses (ERN), movement related cortical potentials
(MRCP), and sensory motor rhythms (SMR). The network proposed in this study
includes temporal convolution layer to learn frequency filters, then a depth
wise convolution layer, connected to each feature map individually, enabling
the model to learn frequency specific spatial filters, followed by a separable
convolution layer, which attempts to learns a temporal summary for each
feature map individually, followed by a point wise convolution, which learns
how to optimally mix the feature maps together. The achieved state of the art
results of this study in many BCI paradigms can be linked to the low number of
parameters of the model and well designed architecture employed.

In [18] proposed a neural network based on the CSP algorithm, consisting of
three steps, first a temporal convolutional for each channel, then a spatio
temporal projection of the filtered signal, which is implemented as a single
fully connected feedforward layer or two divided layers, one for spatial
projection and another for temporal projection. Finally, the extracted
features are fed into the classification layer. The weights of the spatial
layer are initialized by CSP. In a similar paper, [19] also focused on
optimizing a neural network for spatial filter optimization by using only two
fully connected layers.

In the paper [20], an augmented CSP (ACSP) algorithm is suggested by using
overlapping frequency bands. to extract features and use them as an input for
the CNN. This study suffers from the need to do a fair amount of manual
feature engineering and even does not take into account the temporal
characteristics of the EEG signals which can be accounted for relatively poor
results achieved in the motor imagery paradigm

[21] used a stacked autoencoder (SAE) and a CNN to classify the two class,
left/right hand MI task, three electrodes C3, Cz and C4, which are located on
the motor area of the brain. The signals from these selected electrodes are
transformed using a short time Fourier transform (STFT). The output of the
convolutional layer is passed to SAE and then to a classifier. This study has
outperformed SVM.

[22] proposes and experiment with three deep learning architectures, Deep
ConvNet, Shallow ConvNet, and Hybrid Convnet for classifying motor imagery
signals and outperforms the accuracy of the widely used filter bank common
spatial patterns (FBCSP) with minimal preprocessing applied to EEG signals.
Deep ConvNet consists of four convolution max pooling blocks, with a special
first block designed to handle EEG input, followed by three standard
convolution max pooling blocks and softmax classification layer. Shallow
ConvNet architecture consists only of two layers, one temporal convolutional
layer with bigger kernel size compared to Depp ConvNet and then a spatial
convolutional layer, followed by standard non linearity and mean pooling
layer. Hybrid ConvNet just stacks both networks and fuses the features learned
from them after the final layer. The results of this study show that the
shallow networks if designed in an efficient manner, can achieve better
results than more deep models and even standard classifiers used in the BCI
applications.

Few attempts have been carried in the literature to apply Graph Convolutional
architectures to decode brain signals, even fewer attempts have been aimed at
BCI motor imagery applications. In [23], authors used phase locking value
(PLV) connectivity estimation to construct the graph and used manually
extracted time frequency domain feature univariate features, including the
power spectral density (PSD) feature, the differential asymmetry (DASM)
feature, the rational asymmetry (RASM) feature and the differential caudality
(DCAU) features in order to decode the various emotional states of the
subjects. The architecture employed consists of a graph convolutional layer,
two ReLu layers, and a fully connected layer followed by the softmax layer to
perform the classification. The published results in this study, outperforms
SVM, Deep Belief Networks and GCNN.

In another study [24], the EEG signals have been represented as graph based on
the euclidean distance of the electrodes in sensor space. This matrix has been
used to embed the signals in a CSP fashion and then the embedded graph has
been cut into several temporal slices. Two models has been introduced to
extract features from these slices, one CNN model including a convolutional
layer followed by a MaxPooling layer and the other architecture using Long
Short Term Memory (LSTM) units to explore the temporal attentive temporal
dynamics of the slices. These introduced models achieved superior results in
decoding Motor Imagery EEG signals.

In s very recent study [25], authors have proposed a novel architecture
combining the BiLSTM with the GCN to decode EEG tasks which used Pearson
Correlation as an estimate for the graph adjacency matrix, reporting to have
achieved the state of the art performance on Motor Imagery EEG classification
on Physionet dataset.

##  3 Proposed Methods

###  3.1 DeepCSP, novel EEG tailored loss

Inspired by [26], where the authors introduced the formulation of DeepLDA
version of Linear Discriminant Analysis (LDA) enabling the deep learning
networks to impose a separation of classes in the latent space produced by the
network by maximizing the ratio of the between and within class scatter by
means of end to end backpropagation, we introduce the same concept to Common
Spatial patterns (CSP).

CSP is the most common method used in motor imagery signals to achieve
accurate spatial features and easily distinguishable classes. Motor Imagery
experiment deals with predicting the intent of the subjects, commanded to them
by a cue, such as on the screen. In this study, we only investigate
classifying two class cases.

Let $C_{1}$ and $C_{2}$ be different intents, imagined by the subjects, for
example imagining left hand vs right hand movement. Each trial in which
subjects are requested to perform a task, is called an epoch. The EEG signal
recorded for each epoch can be represented as
$\mathbf{X}\in\mathbb{R}^{D\times T}$, $D$ denoting electrodes (channels), and
$T$ the number of time points sampled in each trial. $y\in\\{0,1\\}$ being the
class labels, the objective of the classification is to learn the relations in
$\\{\mathbf{X}_{i},y_{i}\\}_{i=1}^{N}$, where $N$ is the total number of
training trials, such that it can perform well predicting the unseen trials in
test sets.

CSP method reduces the dimensionality of each trial by extracting spatial
components that are common to both classes, in a way that the signals after
being projected onto these components have high variance in between the
classes. This method uses simultaneous diagonalization of the covariance
matrices of the classes. Covariance matrices for trials belonging to each
class. $C_{1}$ and $C_{2}$ are given by:

$C_{1}={{(X_{1}X_{1}^{T})}\over{trace(X_{1}X_{1}^{T})}}$ (1)

$C_{2}={{(X_{2}X_{2}^{T})}\over{trace(X_{2}X_{2}^{T})}}$ (2)

Then these covariances are averaged to get averaged normalized covariance
matrices $\overline{C_{1}}$ and $\overline{C_{2}}$, and add to obtain the
composite covariance matrix C.

$`C=\overline{C_{1}}+\overline{C_{2}}$ (3)

Since $C$ is symmetric, its eigen decomposition has the form

$\mathbf{C}=\mathbf{U}{\mathbf{\Lambda}}\mathbf{U}^{\top}$ (4)

from which we can obtain the whitening transformation matrix $P$ as follows:

$\mathbf{P}={\mathbf{\Lambda}}^{-{{1}/{2}}}\mathbf{U}^{\top}$ (5)

The transformed covariance matrices
${\mathbf{\Sigma}}_{1}=\mathbf{P}\mathbf{C}_{1}\mathbf{P}^{\top}$ and
${\mathbf{\Sigma}}_{2}=\mathbf{P}\mathbf{C}_{2}\mathbf{P}^{\top}$ share common
eigen vectors and the sum of them has the following property:

${\mathbf{\Sigma}}_{1}+{\mathbf{\Sigma}}_{2}=\mathbf{I}_{D\times D}.$ (6)

The spatial projection matrix $W$ is

$\mathbf{W}=\mathbf{P}\mathbf{V}^{\top}$ (7)

Due to the equation eq. 6 top $n$ eigenvectors of ${\mathbf{\Sigma}}_{1}$,
when sorted in decreasing order, both account for the maximum variance in the
first class and also the minimum variance in the second class and vice versa.
Studies shows that $n=2$ or $n=3$ helps to achieve good empirical performance
[27]. The CSP features of the EEG data are computed as

$\mathbf{Z}=\log\left(\mathbf{W}_{2n}^{\top}\mathbf{X}\mathbf{X}^{\top}\mathbf{W}_{2n}\right)$
(8)

We define the new loss objective as

$J(W)={{trace(\mathbf{W}_{2n}^{\top}C_{1}\mathbf{W}_{2n})\over{trace(\mathbf{W}_{2n}^{\top}\left(C_{1}+C_{2}\right)\mathbf{W}_{2n})}}}$
(9)

imposing separation onto the latent representation over the output of the
hidden layers of the deep learning architecture of choice, by not only
learning the spatial filters but also optimizing the parameters of the network
to ensure more separation using backpropagation.

###  3.2 Proposed shallow networks using DeepCSP loss

As discussed, the results of [22] shows that shallow networks can outperform
classic methods like FBCSP and deep network architectures in BCI motor imagery
applications.

We propose a light weight model composed of two blocks. First block is
composed of 3 layers of 1D temporal convolutional filters with multi scale
temporal filters, with kernel sizes of respectively half, one third and one
forth of the sampling frequency of the devices used in the recording of EEG
signals, in an attempt to extract temporal patterns of the input data. The
output of these layers are concatenated and fed into the DeepCSP to obtain the
spatial filters and reduce the dimension of the output of temporal block
further as discussed in 3.1. Second block consists of two layers of fully
connected layers to do classification. A sketch of this model can be found in
fig. 1

Figure 1: Shallow CNN with DeepCSP loss

###  3.3 Proposed shallow GCN networks using DeepCSP loss

As mentioned, although combining CNN with graph theory, the Graph
Convolutional Neural Network (GCN) has shown competitive results in
classifying EEG signals, but very few research has been done to utilize these
architectures in MI BCI applications. One main issue in using these class of
models is not having a clear cut consensus on how to represents the EEG
signals as graphs, in other words the connections and relationship in between
the channels.

In an attempt to address this issue, we investigate two approaches. First we
explore the effect of constructing the graphs from EEG signals based on
statistical methods developed in graph theoretical analysis of EEG signals.
Then we explore constructing the graph dynamically.

###  3.4 Statistical methods to construct the graph from EEG

Here we list the methods widely used in the literature in estimation of
various connectivity in between the channels, that we have utilized to
construct the graph from EEG signals.

  * •

Coherence (Coh) is one of the most commonly utilized connectivity estimators;
it is a measurement of the linear relationship of two signals at a specific
frequency [28]. Given two time series $x$ and $y$, coherence is given by:

$coh^{2}_{xy}(f)=\frac{|G_{xy}(f)^{2}|}{G_{xx}(f)G_{yy}(f)}$ (10)

where the $G_{xy}(f)$ is the estimated cross spectral density between $x$ and
$y$, while $G_{xx}(f)$ and $G_{yy}(f)$ are the autospectrum of $x$ and $y$
respectively. The result is a symmetric matrix of size $[n\\_channels\times
n\\_channels]$ bearing no information about the directionality of the
interaction, with values within the range $[0,1]$.

  * •

Directed Phase Lag Index (dPLI) was introduced in [29] to capture the phase
and lag relationship as a measure of directed functional connectivity.

  * •

Phase Locking Value (PLV) is one of the pioneer methods introduced in [30]; it
utilizes the Hilbert representation of an EEG time series and quantifies their
interaction based on their instantaneous phase in a specific band frequency.

  * •

Imaginary part of Phase Locking Value (IPLV) was proposed to resolve PLV’s
sensitivity to volume conduction and common reference effects. IPLV is
computed similarly as PLV, but taking the imaginary part of the summation
[31].

  * •

Phase Lag Index (PLI) [32] proposed as an alternative (to PLV) phase
synchronization estimator that is less prone to the effects of common sources
(namely, volume conduction and active reference electrodes). These effects can
artificially generate functional connectivity as the same signal signal is
measured at different electrodes [33]. PLI estimates the asymmetry in the
distribution of two time series’ instantaneous phase differences.

  * •

Weighted Phase Lag Index and Debiased Weighted Phase Lag Index (wPLI) and
(dwPLI) were introduced in [34] to achieve better estimations

After defining the connectivity between channels, we use it as prior in form
of the weighted adjacency matrix and use a layer of GraphSage convolution [35]
to extract graph features and concatenate them with the temporal features as
in fig. 1. Our experiments show that defining the temporal features as node
feature inputs to GraphSage yields better results.

##  4 Results and Discussion

In order to test the validity of proposed methods, we have used two Motor
Imagery datasets and compared the results with some of the state of the art
deep learning architectures. First, we give a brief overview of the datasets
that we have run the models, then compare the results subject dependent
classification results then analyze the proposed methods.

###  4.1 Datasets

The first Motor Imagery dataset is the lla dataset from BCI competition 4
[36]. This dataset consists of EEG recordings of 9 subjects. Four cue based
BCI paradigms in these experiments is the imagination of left hand, right
hand, both feet and tongue. 12 trials were conducted for each class, in 6 runs
separated by short breaks on two different days, a total of 576 trails. We
have selected only the trials belonging to the intent of the left hand vs
right hand movement of this dataset, since at this point the DeepSCP only
supports two classes.

Each trial was started by a fixation mark on a computer screen for the
subjects along with a short warning sound. after 2 seconds subjects were cued
by an arrow to perform the associated intent depending on the direction of the
arrow, which stayed on the screen for 1.25 seconds, without any feedback. The
subjects were asked to repeat the intent until the fixation mark on the screen
disappeared. The signals were sampled with. $250Hz$ and bandpass filtered
between $0.5Hz$ and $100Hz$. An additional $50Hz$ notch filter was enabled to
suppress line noise.

The second Motor Imagery dataset [37], consisted of eight runs, five of them
for training and three with feedback for validation each run having 20 trials,
50 trials per class for training and 30 trials per class for validation.
Subjects were given the task of performing a cued motor imagery (MI) of the
right hand and of the for 5 seconds. The Experiment has been carried as the
one mentioned in the first dataset, with a sampling rate of $512Hz$. The
electrodes placement was designed for obtaining three Laplacian derivations.
Center electrodes at positions C3, Cz, and C4 and four additional electrodes
around each center electrode with a distance of 2.5 cm, 15 electrodes total.
The participants were aged between 20 and 30 years, 8 naive to the task, and
had no known medical or neurological diseases.

##  5 Experiments and Results

We have used the ShallowConv and DeepConv CNN’s proposed in [22] and EEGNet
fourth version [17]. We have used the official implementations provided by the
authors and used the recommended training settings by the authors for each
one. The optimizer being pytorch implementation of AdamW with $0.00625$
learning rate and zero wight decay and for the DeepConv Net the learning rate
of $0.01$ and weight decay of $0.0005$ with a batch size of 64 and Adam
optimizer with a learning rate of $0.001$. We have used the EEGNet
implementation provided in the BrainCode, a deep learning toolbox to decode
raw time domain EEG 111https://braindecode.org/.

We have trained the Shallow CNN DeepCSP and Shallow GCN with DeepCSP with
separate SGD optimizers, one with a learning rate of $0.01$ for updating the
parameters of the CNN and GCN part with respect to proposed DeepSCP loss, and
one with $0.1$ for updating the parameters for the linear layers of the
classification block in consecutive order. Our experiments showed that setting
the number of components of the spatial filters learned by the loss to 4,
yields better results. So here we include only those experiments. One thing
worth mentioning here is that in training, we feed the whole training set to
get more accurate spatial filters. Small scale model and data set makes this
doable even when not training the model on GPU.

Table 1: Comparison of results on the BCI competition dataset [36] Table 2:
Comparison of results on the Motor Imagery dataset [37]

The results presented in table 1 and table 2, the proposed DeepCSP adds
considerable amount to the accuracy of the classification. Using graph based
extracted features, increases the accuracy a small amount too. This comparison
of the different approaches to define connectivity matrices of EEG signals and
represent them as graphs, effects the model accuracy dramatically. We confirm
the findings in [23], that PLV method for connectivity estimation, is a good
place to start.

One explanation for these increases in the accuracy can be found if we take a
look at the output features of the temporal filters in Shallow DeepCSP. In
order to see this, in fig. 2, we plot the transformed outputs of the temporal
block, the number of components is set to two. The separation of the classes
can be easily noticed, which makes it each for the network to optimize the
parameters of the classification block to perform classifying. One has to keep
in mind that the spatial filters are learned in the training and then is used
to transform the unseen trials.

Figure 2: The output of the temporal layer in Shallow DeepCSP after being
projected by the spatial filters learned by DeepCSP loss for subject number 3
in BCI competition dataset from [36] . $C_{0}$ and $C_{1}$ denotes the EEG
signals of intent of right hand and left hand movement, respectively. For
visualization reasons we have chosen number of components of spatial filters
as two

The spatial filters learned from DeepCSP loss, can be used to plot EEG
topographic maps (topomap) to find the brain activation. In the case of
performing intent of right or left hand movement, we should see activation in
the motor cortex region of the brain which corresponds to electrodes $C3$ and
$C4$. "An event with the high feature for the first component (channel C3) and
low feature for the second component (channel C4) relates to the right hand
imagery movement. On the contrary, this event relates to the left hand imagery
movement" [38]. The plot for the component corresponding to 0 in the far left
of section C of fig. 3 fails to illustrate this. Mainly because one of the
shortcomings of CSP method is it requires carefully chosen filter bands. On
the other hand, we can clearly see two patterns with different signs over
regions of $C3$ and $C4$ in the Topomaps obtained from DeepCSP in an end to
end manner.

Figure 3: A) Position of electrodes in [36] B) Topomaps generated from spatial
filters obtained from DeepCSP C) Topomaps generated from spatial filters
obtained from CSP

The results for our limited attempts to employ dynamic graph convolutional
networks did not yield good results and have not been presented in here. The
reason for the poor results can mainly be contributed to the large increase in
the number of parameters in the architecture to update the graph dynamically
at each iteration, which causes the model to over fit very quickly in early
epochs, which is mainly due to the large the noise to signal ratio in the EEG
signals.

##  6 Conclusion and future research plans

Our preliminary experiments provided in the previous section, the proposed non
linear, differentiable generalization of common spatial filters as a new loss,
DeepCSP, shows promising results in Motor Imagery BCI applications by forcing
separation on the output feature space of the CNN and GCN layers. Despite the
good performance of the model proposed, it is designed for two class
experiments. This can be solved by extending the DeepCSP for multi class case,
benefiting from the existing literature of CSP method for multi class
settings. We have seen that representing EEG signals as graphs and extracting
graph features using recent graph convolutions, adds to the performance of the
model. However, it is noteworthy to mention that among various kinds of graph
convolutions available, we have only experimented with GraphSage, only because
of the availability of the source code and ease of use of it. The main goal of
future research is to extend this model to subject independent applications.

## References

  * [1] Dennis J McFarland and Jonathan R Wolpaw.  Brain-computer interface operation of robotic and prosthetic devices.  Computer, 41(10):52–56, 2008.
  * [2] Fabien Lotte, Marco Congedo, Anatole Lécuyer, Fabrice Lamarche, and Bruno Arnaldi.  A review of classification algorithms for eeg-based brain–computer interfaces.  Journal of neural engineering, 4(2):R1, 2007.
  * [3] Cuntai Guan, Manoj Thulasidas, and Jiankang Wu.  High performance p300 speller for brain-computer interface.  In IEEE International Workshop on Biomedical Circuits and Systems, 2004., pages S3–5. IEEE, 2004.
  * [4] Rabie A Ramadan and Athanasios V Vasilakos.  Brain computer interface: control signals review.  Neurocomputing, 223:26–44, 2017.
  * [5] Charles Spence and Ophelia Deroy.  Crossmodal mental imagery.  In Multisensory imagery, pages 157–183. Springer, 2013.
  * [6] Hyohyeong Kang, Yunjun Nam, and Seungjin Choi.  Composite common spatial pattern for subject-to-subject transfer.  IEEE Signal Processing Letters, 16(8):683–686, 2009.
  * [7] Mahnaz Arvaneh, Cuntai Guan, Kai Keng Ang, and Chai Quek.  Optimizing the channel selection and classification accuracy in eeg-based bci.  IEEE Transactions on Biomedical Engineering, 58(6):1865–1873, 2011\.
  * [8] Wojciech Samek, Carmen Vidaurre, Klaus-Robert Müller, and Motoaki Kawanabe.  Stationary common spatial patterns for brain–computer interfacing.  Journal of neural engineering, 9(2):026013, 2012.
  * [9] Mahnaz Arvaneh, Cuntai Guan, Kai Keng Ang, and Chai Quek.  Eeg data space adaptation to reduce intersession nonstationarity in brain-computer interface.  Neural computation, 25(8):2146–2171, 2013.
  * [10] Peng Xu, Ping Yang, Xu Lei, and Dezhong Yao.  An enhanced probabilistic lda for multi-class brain computer interface.  PloS one, 6(1), 2011.
  * [11] Juliano Machado, Alexandre Balbinot, and Adalberto Schuck.  A study of the naive bayes classifier for analyzing imaginary movement eeg signals using the periodogram as spectral estimator.  In 2013 ISSNIP Biosignals and Biorobotics Conference: Biosignals and Robotics for Better and Safer Living (BRC), pages 1–4. IEEE, 2013.
  * [12] Mu Li and Bao-Liang Lu.  Emotion classification based on gamma-band eeg.  In 2009 Annual International Conference of the IEEE Engineering in medicine and biology society, pages 1223–1226. IEEE, 2009.
  * [13] Hubert Cecotti.  A time–frequency convolutional neural network for the offline classification of steady-state visual evoked potential responses.  Pattern Recognition Letters, 32(8):1145–1153, 2011.
  * [14] Sebastian Stober, Daniel J Cameron, and Jessica A Grahn.  Using convolutional neural networks to recognize rhythm stimuli from electroencephalography recordings.  In Advances in neural information processing systems, pages 1449–1457, 2014.
  * [15] Sebastian Stober, Avital Sternin, Adrian M Owen, and Jessica A Grahn.  Deep feature learning for eeg recordings.  arXiv preprint arXiv:1511.04306, 2015.
  * [16] Pouya Bashivan, Irina Rish, Mohammed Yeasin, and Noel Codella.  Learning representations from eeg with deep recurrent-convolutional neural networks.  arXiv preprint arXiv:1511.06448, 2015.
  * [17] Vernon J Lawhern, Amelia J Solon, Nicholas R Waytowich, Stephen M Gordon, Chou P Hung, and Brent J Lance.  Eegnet: a compact convolutional neural network for eeg-based brain–computer interfaces.  Journal of neural engineering, 15(5):056013, 2018.
  * [18] Eder Santana, Austin J Brockmeier, and Jose C Principe.  Joint optimization of algorithmic suites for eeg analysis.  In 2014 36th Annual International Conference of the IEEE Engineering in Medicine and Biology Society, pages 2997–3000. IEEE, 2014.
  * [19] Ayhan Yuksel and Tamer Olmez.  A neural network-based optimal spatial filter design method for motor imagery classification.  PloS one, 10(5), 2015.
  * [20] Huijuan Yang, Siavash Sakhavi, Kai Keng Ang, and Cuntai Guan.  On the use of convolutional neural networks and augmented csp features for multi-class motor imagery of eeg signals classification.  In 2015 37th Annual International Conference of the IEEE Engineering in Medicine and Biology Society (EMBC), pages 2620–2623. IEEE, 2015\.
  * [21] Yousef Rezaei Tabar and Ugur Halici.  A novel deep learning approach for classification of eeg motor imagery signals.  Journal of neural engineering, 14(1):016003, 2016.
  * [22] Robin Tibor Schirrmeister, Jost Tobias Springenberg, Lukas Dominique Josef Fiederer, Martin Glasstetter, Katharina Eggensperger, Michael Tangermann, Frank Hutter, Wolfram Burgard, and Tonio Ball.  Deep learning with convolutional neural networks for eeg decoding and visualization.  Human brain mapping, 38(11):5391–5420, 2017.
  * [23] Z. Wang, Y. Tong, and X. Heng.  Phase-locking value based graph convolutional neural networks for emotion recognition.  IEEE Access, 7:93711–93722, 2019.
  * [24] Dalin Zhang, Kaixuan Chen, Debao Jian, and Lina Yao.  Motor imagery classification via temporal attention cues of graph embedded eeg signals.  IEEE Journal of Biomedical and Health Informatics, 2020.
  * [25] Yimin Hou, Shuyue Jia, Shu Zhang, Xiangmin Lun, Yan Shi, Yang Li, Hanrui Yang, Rui Zeng, and Jinglei Lv.  Deep feature mining via attention-based bilstm-gcn for human motor imagery recognition.  arXiv preprint arXiv:2005.00777, 2020.
  * [26] Matthias Dorfer, Rainer Kelz, and Gerhard Widmer.  Deep linear discriminant analysis.  arXiv preprint arXiv:1511.04707, 2015.
  * [27] Johannes Müller-Gerking, Gert Pfurtscheller, and Henrik Flyvbjerg.  Designing optimal spatial filters for single-trial eeg classification in a movement task.  Clinical neurophysiology, 110(5):787–798, 1999.
  * [28] Guido Nolte, Ou Bai, Lewis Wheaton, Zoltan Mari, Sherry Vorbach, and Mark Hallett.  Identifying true brain interaction from eeg data using the imaginary part of coherency.  Clinical neurophysiology, 115(10):2292–2307, 2004.
  * [29] Cornelis J Stam and Elisabeth CW van Straaten.  Go with the flow: use of a directed phase lag index (dpli) to characterize patterns of phase relations in a large-scale model of brain dynamics.  Neuroimage, 62(3):1415–1428, 2012.
  * [30] Jean-Philippe Lachaux, Eugenio Rodriguez, Jacques Martinerie, and Francisco J Varela.  Measuring phase synchrony in brain signals.  Human brain mapping, 8(4):194–208, 1999.
  * [31] Sepideh Sadaghiani, René Scheeringa, Katia Lehongre, Benjamin Morillon, Anne-Lise Giraud, Mark D’Esposito, and Andreas Kleinschmidt.  Alpha-band phase synchrony is related to activity in the fronto-parietal adaptive control network.  Journal of Neuroscience, 32(41):14305–14310, 2012.
  * [32] Cornelis J Stam, Guido Nolte, and Andreas Daffertshofer.  Phase lag index: assessment of functional connectivity from multi channel eeg and meg with diminished bias from common sources.  Human brain mapping, 28(11):1178–1193, 2007.
  * [33] Martin Hardmeier, Florian Hatz, Habib Bousleiman, Christian Schindler, Cornelis Jan Stam, and Peter Fuhr.  Reproducibility of functional connectivity and graph measures based on the phase lag index (pli) and weighted phase lag index (wpli) derived from high resolution eeg.  PloS one, 9(10), 2014.
  * [34] Martin Vinck, Robert Oostenveld, Marijn Van Wingerden, Franscesco Battaglia, and Cyriel MA Pennartz.  An improved index of phase-synchronization for electrophysiological data in the presence of volume-conduction, noise and sample-size bias.  Neuroimage, 55(4):1548–1565, 2011.
  * [35] Will Hamilton, Zhitao Ying, and Jure Leskovec.  Inductive representation learning on large graphs.  In Advances in neural information processing systems, pages 1024–1034, 2017.
  * [36] Michael Tangermann, Klaus-Robert Müller, Ad Aertsen, Niels Birbaumer, Christoph Braun, Clemens Brunner, Robert Leeb, Carsten Mehring, Kai J Miller, Gernot Mueller-Putz, et al.  Review of the bci competition iv.  Frontiers in neuroscience, 6:55, 2012.
  * [37] David Steyrl, Reinhold Scherer, Josef Faller, and Gernot R Müller-Putz.  Random forests in non-invasive sensorimotor rhythm brain-computer interfaces: a practical and convenient non-linear classifier.  Biomedical Engineering/Biomedizinische Technik, 61(1):77–86, 2016\.
  * [38] Anh Huy Phan and Andrzej Cichocki.  Tensor decompositions for feature extraction and classification of high dimensional datasets.  Nonlinear theory and its applications, IEICE, 1(1):37–68, 2010\.


# Interaction-Aware Sampling-Based MPC
with Learned Local Goal Predictions

Walter Jansma, Elia Trevisan, Álvaro Serra-Gómez and Javier Alonso-Mora This
research is supported by the project “Sustainable Transportation and Logistics
over Water: Electrification, Automation and Optimization (TRiLOGy)” of the
Netherlands Organization for Scientific Research (NWO), domain Science (ENW),
and the Amsterdam Institute for Advanced Metropolitan Solutions (AMS) in the
Netherlands.The authors are with the Cognitive Robotics Department, TU Delft,
<EMAIL_ADDRESS>{e.trevisan, a.serragomez<EMAIL_ADDRESS>
###### Abstract

Motion planning for autonomous robots in tight, interaction-rich, and mixed
human-robot environments is challenging. State-of-the-art methods typically
separate prediction and planning, predicting other agents’ trajectories first
and then planning the ego agent’s motion in the remaining free space. However,
agents’ lack of awareness of their influence on others can lead to the
freezing robot problem. We build upon Interaction-Aware Model Predictive Path
Integral (IA-MPPI) control and combine it with learning-based trajectory
predictions, thereby relaxing its reliance on communicated short-term goals
for other agents. We apply this framework to Autonomous Surface Vessels (ASVs)
navigating urban canals. By generating an artificial dataset in real sections
of Amsterdam’s canals, adapting and training a prediction model for our
domain, and proposing heuristics to extract local goals, we enable effective
cooperation in planning. Our approach improves autonomous robot navigation in
complex, crowded environments, with potential implications for multi-agent
systems and human-robot interaction.

Dataset, Prediction Model, Video and Code available at:
https://sites.google.com/view/ia-mppi-lbm/

##  I Introduction

Cities characterized by dense networks of urban canals, such as Amsterdam,
could greatly benefit from deploying Autonomous Surface Vessels (ASVs) for
various tasks including deliveries, transportation of people, and garbage
collection [1]. However, navigating autonomously in urban canals amidst mixed
human-robot crowds presents a significant challenge. Urban canals are
typically narrow, frequently congested, and lack the structured nature of
roads. While not as strictly enforced as on roads, navigation principles like
right-of-way and right-hand conventions should still be considered. Thus, akin
to autonomous ground robots among pedestrian crowds, successful navigation in
urban canals relies on cooperation and awareness of interactions [2].

Recently, a sampling-based Model Predictive Control (MPC) called Interaction-
Aware Model Predictive Path Integral (IA-MPPI) control has been developed for
generating cooperative motion plans in urban canals among multiple non-
communicating vessels while maintaining awareness of navigation rules [3].
This algorithm assumes rational and homogeneous agents, exact sensing of
states, and knowledge of local goals. In real-time, the algorithm samples
thousands of input sequences to approximate the optimal input sequence that
enables all agents to progress toward their goals cooperatively. In scenarios
where the local goals of other vessels are unavailable, such as in mixed
human-robot environments or due to lack of communication, this previous
approach has approximated these goals using a constant velocity model over a
given horizon. However, in narrow and crowded environments, vessels often need
to execute complex maneuvers to navigate tight intersections and avoid
collisions while adhering to navigation rules. In such situations, relying
solely on a constant velocity approximation can lead to inaccurate
predictions, which can adversely affect the performance of the motion planner
in terms of deadlocks, collisions, navigation rule violations, traveled
distance, and travel time.

In this paper, we present a framework (see Fig. 1) that utilizes a learning-
based trajectory prediction method to improve the estimation of agents’
intended destinations. We introduce heuristics to extract local goals from the
predicted trajectories and provide the motion planner with the flexibility to
influence the behavior of other agents while expecting cooperation in
collision avoidance.

Figure 1: Overview of the proposed framework. Firstly, the prediction model
utilizes information from all elements in the scene to forecast trajectories
for obstacle agents. Meanwhile, the global planner, equipped with the map,
start, and goal positions, generates a path for the ego agent. Subsequently,
the local goal extractor leverages this information to determine appropriate
local goals for the motion planner. With inputs derived from the scene and the
local goals, the Interaction-Aware Model Predictive Path Integral (IA-MPPI)
algorithm simultaneously plans and predicts input sequences for all agents in
the scene. The first input of the sequence is then assigned to the ego agent,
and the algorithm iterates.

###  I-A Related Work

Robot motion planning in dynamic environments is a challenging problem for
which a series of classical and heuristic-based approaches have been developed
[4], such as the Dynamic Window Approach [5] or Reciprocal Velocity Obstacles
[6, 7]. Despite their successful applications, e.g. to non-holonomic robots
[8] or vessels in open waters [9], the motions planned by this class of
methods are often reactive. This, especially in crowded environments, can lead
to unsafe and unpredictable behaviors.

Model Predictive Control (MPC) has become a popular approach to trajectory
planning for autonomous vehicles [10] because of its ability to optimize
accounting for the system’s dynamics and constraints. Moreover, by planning
over a sufficiently large horizon, MPC can anticipate dynamic obstacles
resulting in trajectories that are less reactive. To anticipate other agents,
however, the free space over the entire planning horizon needs to be computed
[11], which requires knowledge about other agents’ positions in the future. If
all the agents in the environment are autonomous, communication and
distributed optimization can be used to plan trajectories in multi-agent
environments [12]. In mixed human-robot environments, however, such
communication is not possible and predictions of the future motion of the
other agents have to be employed. For instance, recent work on MPC for rule-
aware navigation in urban canals uses constant velocity to model the future
behavior of other vessels [13].

In interaction-rich scenarios, however, constant velocity can be an inaccurate
approximation which may lead to unsafe motion plans [14]. Therefore, several
works rely on learning-based models to predict the future motion of other
agents [15] and can include prediction confidence [16] and multimodality [17].
These methods, however, decouple prediction and planning which, in high-
interaction environments, may lead the ego agent to wrongly assume that no
collision-free path exists [18]. To avoid the so-called freezing robot problem
the robot has to expect cooperation in collision avoidance from the other
agents [19]. Coupled prediction and planning can be done with MPC by modeling
the interacting agents as a system, but it quickly becomes expensive to solve
via constrained optimization leading to long computation times and short
planning horizons [20, 21].

Building upon a novel sampling-based Model Predictive Control (MPC)
framework[22], Interaction-Aware Model Predictive Path Integral (IA-MPPI)
control [3] has successfully demonstrated decentralized coupled predictions
and planning in real-time, accommodating long prediction horizons, nonlinear
dynamics, and discontinuous cost functions in multi-agent environments. While
IA-MPPI has exhibited superior performance compared to optimization-based MPC
approaches that rely on fixed predictions of other agents’ motion, it
necessitates knowledge of their near-term local goals, which can either be
communicated or estimated.

###  I-B Contribution

This paper presents a novel framework for interaction-aware decentralized
motion planning in urban canals without relying on communication. Our
framework encompasses the following contributions:

  * •

Realistic Dataset: We generate and publish a realistic dataset of simulated
rule-abiding vessel trajectories in real sections of Amsterdam’s urban canals.

  * •

Learning-Based Trajectory Prediction: We adapt a pedestrian prediction model
[23] to vessels and train it specifically for urban canals. This approach
enables us to generate trajectory predictions for other agents.

  * •

Local Goal Extraction: We propose heuristics to extract local goals from the
predicted trajectories, thereby providing the motion planner with information
about where agents intend to go.

  * •

Communication-Free Coupled Prediction and Planning: By combining the local
goal extraction with the IA-MPPI control [3], we achieve coupled prediction
and planning without the need for communication. This approach ensures that
the ego agent can influence the behavior of other agents while anticipating
cooperation in collision avoidance.

We validate our planning framework through extensive simulated experiments,
comparing it against baseline approaches and providing insights into the
benefits of coupled prediction and planning over decoupled methods. The
framework can be adapted to other robot types beyond vessels.

##  II Interaction-Aware MPPI

In this section, we introduce the main ideas of IA-MPPI [3], upon which our
proposed framework is built. For details on the method, models used and cost
function please refer to the original paper. For insights on the underlying
sampling-based MPC, one can refer to the work on Information-Theoretic MPC
[22]. In short, IA-MPPI assumes that all the agents are homogenous and
rational, i.e. have the same model and cost function. Under this assumption,
we can create a large multi-agent system and plan input sequences resulting in
cooperative trajectories for the ego agent as well as all the obstacle agents.
This being a decentralized planning framework, we then apply the first input
of the sequence to our ego agent, observe the environment and plan again. In
more detail, IA-MPPI models the ego-agent $i$ as a discrete-time dynamical
system,

$\mathbf{q}_{i,t+1}=\mathcal{F}(\mathbf{q}_{i,t},\mathbf{u}_{i,t})$ (1)

where $\mathbf{q}_{i,t}$ and $\mathbf{u}_{i,t}$ are, respectively, the state
and the input of the ego-agent at timestep $t$. The state
$\mathbf{q}_{i,t}=[\mathbf{p}_{i,t},\mathbf{v}_{i,t}]$ contains the position
and velocity of the agent. IA-MPPI assumes that all agents in the environment
are homogenous. The state and the input of the multi-agent system consisting
of the ego-agent and the obstacle agents can therefore be stacked, resulting
in,

$\begin{aligned}
\mathbf{q}&=\begin{bmatrix}\mathbf{q}_{i}^{\top}&\mathbf{q}^{\top}_{j}\end{bmatrix}^{\top},\\\
\mathbf{u}&=\begin{bmatrix}\mathbf{u}_{i}^{\top}&\mathbf{u}^{\top}_{j}\end{bmatrix}^{\top},\end{aligned}\quad\forall
j\in\mathcal{M}\setminus i,$ (2)

where $\left(.\right)_{j}$ is a variable that the ego-agent $i$ estimates of
agent $j$ and $\mathcal{M}=\\{0,1,...,m\\}$ is the set of all agents in the
scene. By also stacking the state transition functions $\mathcal{F}$ over all
agents, we obtain a model for the multi-agent system
$\mathbf{q}_{t+1}=\mathcal{G}(\mathbf{q}_{t},\mathbf{u}_{t})$. Given a
planning horizon $T$ and a prior input sequence
$\mathbf{U}=[\mathbf{u}_{0},\mathbf{u}_{1},\dots,\mathbf{u}_{T-1}]$, IA-MPPI
samples $K$ input sequences for the entire multi-agent system,

$\tilde{\mathbf{U}}_{k}=[\tilde{\mathbf{u}}_{0,k},\tilde{\mathbf{u}}_{1,k},\dots,\tilde{\mathbf{u}}_{T-1,k}],\quad\tilde{\mathbf{u}}_{t,k}=\mathcal{N}(\mathbf{u}_{t},\nu\mathbf{\Sigma})$
(3)

with $k=1,\dots,K$, variance $\Sigma$ and scaling parameter $\nu$. At the
first iteration, the prior input sequence $\mathbf{U}$ is initialized at zero.
By the end of this section, it will become clear how this prior input sequence
is updated in subsequent iterations. Having a model for the multi-agent
system, we can forward simulate the $K$ input sequences into $K$ state
trajectories $\mathbf{Q}_{k}$ for the multi-agent system,

$\mathbf{Q}_{k}=\big{[}\mathbf{q}_{0},\,\mathcal{G}(\mathbf{q}_{0},\tilde{\mathbf{u}}_{k,0}),\,\dots,\,\mathcal{G}(\mathbf{q}_{k,T-1},\,\tilde{\mathbf{u}}_{k,T-1})\big{]}.$
(4)

Each of the resulting state trajectories is evaluated with respect to both an
agent-centric cost as well as a system-wide cost, resulting in a total sample
cost $S_{k}$. The reader can refer to the original publication for details on
the cost function [3]. For the scope of our paper, it is important to know
that the agent-centric cost includes a tracking cost to encourage progress
towards a local goal $p_{g}$ computed as,

$C_{\text{tracking}}=k_{\text{tracking}}\frac{||\mathbf{p}_{g}-\mathbf{p}_{t}||_{2}}{||\mathbf{p}_{g}-\mathbf{p}_{t_{0}}||_{2}},$
(5)

where $p_{t}$ is the position of the agent at timestep $t$, $p_{t_{0}}$ is the
position of the agent at the beginning of the planning horizon and
$k_{tracking}$ is a tuning parameter. Notice that we need to know the position
of the local goal of each agent. For the ego agent, the local goal is
extracted from a global plan. For all the other agents, the local goal has to
be either communicated or estimated. We propose in the following section how
this goal can be estimated. Once $S_{k}$, $\forall k\in[1,\dots,K]$ has been
computed, importance sampling weights $w_{k}$ can be calculated as,

$w_{k}=\frac{1}{\eta}\exp\biggl{(}\frac{-1}{\lambda}(S_{k}-S_{min})\biggr{)},\quad\sum_{k=0}^{K-1}w_{k}=1,$
(6)

where $S_{min}$ is the minimum sampled cost, $\eta$ a normalization factor and
$\lambda$ a tuning parameter. We then compute an approximation of the optimal
control sequence through a weighted average of the sampled control sequences,

$\mathbf{U}^{*}=\sum_{k=0}^{K-1}w_{k}\tilde{\mathbf{U}}_{k}$ (7)

and apply the first input $\mathbf{u}^{*}_{i,0}$ to the ego-agent. We can now
use a time-shifted version of $\mathbf{U}^{*}$ as the prior input sequence
$\mathbf{U}$ to warm-start the sampling strategy at the next iteration.

##  III Predicting goal positions

In Fig. 1 we provide an overview of the proposed framework. In Section III-A,
we outline the prediction model. In Section III-B, we describe the dataset we
have collected to train a prediction model that is interaction and rule-aware.
In Section III-C, we present the steps taken to port the prediction model to
urban vessel environments. In Section III-D, we propose a heuristic to extract
a local goal suitable for IA-MPPI using the predicted trajectories.

###  III-A Interaction-aware trajectory prediction method

Our approach leverages interaction-aware trajectory prediction for goal
estimation. We employ an adapted version of Social-VRNN [23], which was
originally designed for pedestrians, to obtain trajectory predictions.
However, we remark that our framework is agnostic to the choice of trajectory
predictor as long as it accounts for obstacles and interactions between agents
in the environment.

Social-VRNN [23] is an interaction-aware trajectory prediction method that
leverages a generative model based on Variational Recurrent Neural Networks
(VRNNs) [24]. The model combines three types of contextual cues to define a
joint representation of an agent’s current state: information on the past
trajectory of the agent of interest, environment context, and agent-agent
interactions. The input to predict the trajectory of agent $i$ is denoted as:

$\textbf{x}=\\{\textbf{v}^{i}_{-T_{0}:0},\textbf{O}^{i}_{env},\textbf{O}^{-i}_{int}\\},$
(8)

where $\textbf{v}^{i}_{-T_{0}:0}$ corresponds to the sequence of velocity
states over the previous observed horizon $T_{O}$ of the agent of interest
$i$. The environment information $\textbf{O}^{i}_{env}$ is represented in the
form of a grid map extracted around the agent of interest. Then,
$\textbf{O}^{-i}_{int}$ represents the information on agent-agent
interactions. It is a vector with the relative positions and velocities of all
other agents from agent $i$’s perspective, listed in ascending order based on
the absolute distance to it. The output of the model is a sequence of velocity
probability distributions represented by $T_{H}$ diagonal gaussian
distributions
$\mathcal{N}(\mu_{\textbf{v},k},\text{diag}(\sigma^{2}_{\textbf{v},k}))$. For
details on the method and its architecture, please refer to the original paper
[23].

###  III-B Artificial Dataset

In the absence of a publicly available dataset for short-term vessel
trajectory prediction, an artificial dataset of vessel interactions is
collected in a simulation environment. In order to obtain trajectories that
resemble those of real vessels in urban canals, four real canal section maps
in Amsterdam: the Herengracht (HG), the Prinsengracht (PG) and the Bloemgracht
(BG) are used to collect data. The Open Crossing (OC) environment is created
to collect vessel interactions in open water. Data on an additional
environment, the Amstel (AM), is included only for testing our framework’s
generalization to environments not seen during training. Figure 2 depicts two
of these canal sections. The yellow rectangles correspond to the areas in
which start and goal locations are randomly initialized. These areas are
placed around the entire map and in each canal section to improve the
diversity of the trajectories and interactions.

To collect the data, more than four thousand experiments are conducted by
initializing up to four vessels simultaneously in the mentioned environments.
Each vessel is assigned a randomized start and goal location in one of the
predefined areas. All sampled locations are ensured to be collision-free. The
vessels run a centralized IA-MPPI to sail toward their respective goals while
accounting for navigation rules. This ensures that the recorded trajectories
are safe, interaction-aware, and mostly rule-abiding.

For each experiment and vessel in the environment we record the current
timestamp, the vessel ID, its position and velocity in the global frame. Each
timestamp is unique across timesteps and experiments, which enables to
identify vessels belonging to the same scene. The specifications of the
artificial vessel dataset can be found in Table I. In order to evaluate the
prediction model, 10% of the dataset is used as the test set. The remaining
data is used for training and is split into a training set (72%) and a
validation set (18%). The distribution of data from each scenario is ensured
to be equal in all splits.

TABLE I: Specifications of the artificial dataset. Exp. refers to the number of experiments done in each scenario. Frames and Vessels refer to the total number of frames and vessels present in the data set, respectively. All data is recorded at a rate of 10. Scenario | Exp. | Frames | Vessels
---|---|---|---
Herengracht | 1000 | 406229 | 2499
Prinsengracht | 1247 | 420285 | 4122
Bloemgracht | 1188 | 372173 | 3564
Open Crossing | 1182 | 417515 | 3544
Amstel | 79 | 23468 | 316
Total | 4696 | 1639670 | 14045
Figure 2: Canal sections of the Bloemgracht and Prinsengracht. The black areas
are the canals. The yellow rectangles correspond to the initialization areas
in which goals and starting locations were randomly initialized for each agent
during the simulations.

###  III-C Model Training and Adaptation

We adapt the variational inference architecture presented in [23] to generate
unimodal trajectory probability predictions of vessels. In contrast to humans,
vessels are slower and have lower-order dynamics, which results in less
reactive behaviors and smoother trajectories. To take this into account and
avoid overfitting to the dataset, we reduce the dimensionality of the method’s
latent space. We also add an L2-regularization term to the loss function and
weight it with a hyperparameter we define as $\gamma$.

####  III-C1 Hyperparameters

The model is trained using backpropagation through time and the RMSProp [25]
optimizer. With a time step of $\Delta T=0.4$ seconds, the prediction horizon
is set to $T_{H}=24$ steps (9.6 seconds) and the previous horizon to
$T_{O}=14$ steps (5.6 seconds). Furthermore, we employ learning rate starting
at $\alpha$ = $1\mathrm{e}{-4}$ that decays by a factor of 0.9 after every
gradient step. The regularization weight is kept at $\gamma=0.0001$. Finally,
the model is trained for $4\mathrm{e}{4}$ training steps, using early
stopping.

###  III-D Local Goal Extraction

In eq. (5) we show that the IA-MPPI needs to know the local goal $p_{g}$ of
each agent. There are two requirements for a goal to be suitable: it has to
lie within a radius $r_{p_{g}}$ from the agent it corresponds to and cannot be
in space occupied by static obstacles. Therefore, we first search the
predicted trajectory backward until we obtain a position $p_{\leq r_{p_{g}}}$
within the desired radius. If $p_{\leq r_{p_{g}}}$ is in collision with a
static obstacle, we construct a circle centered on the agent’s position
$p_{a}$ with radius $p_{a}-p_{\leq r_{p_{g}}}$ and find the point on the
circle closest to $p_{\leq r_{p_{g}}}$ which is not in collision with static
obstacles. This goal extraction method is illustrated in Fig. 3. Once the
goals for all agents are predicted, IA-MPPI can plan interaction-aware
trajectories in a decentralized fashion.

Figure 3: A visual illustration of how the local goal is extracted from a
colliding trajectory prediction. Figure 4: The displacement error of the
predictions from CVM and the LBM over the prediction horizon for each canal
section. The solid line represents the mean error and the shaded area
represents 30% of the standard deviation.

##  IV Experiments

The experiments are conducted in real maps of Amsterdam’s canals, namely the
Herengracht (HG), Bloemgracht (BG), Prinsengracht (PG), and the Amstel (AM).
In addition, experiments are conducted in an Open Crossing (OC) map without
static obstacles. In Section IV-A we evaluate the prediction model, in Section
IV-B we show the performances of the proposed framework for motion planning,
and in Section IV-C we highlight the benefits of coupled prediction and
planning with respect to a decoupled approach.

###  IV-A Prediction Accuracy

In Fig. 4 we compare the proposed Learning-Based Model (LBM) to a Constant
Velocity Model (CVM) on test data. We evaluate the methods against the
displacement error at each prediction step, which is defined as the Euclidean
distance between a prediction and the ground truth. In all maps the LBM
outperforms the CVM, showing a lower average displacement error and a smaller
standard deviation. Note that the Amstel map was previously unseen during
training, demonstrating generalization capabilities.

Figure 5: Examples of experiments in the low-interaction scenario (left) and
high-interaction scenario (right).

###  IV-B Interaction-Aware Motion Planning with Predictions

In this study, we evaluate the performance of the proposed decentralized
framework that uses a Learning-Based prediction Model to extract local goals
(IA-MPPI-LBM), by comparing it against a decentralized approach that extracts
local goals from a Constant Velocity Model (IA-MPPI-CVM) and decentralized
with communication (IA-MPPI-w/comm.), which assumes perfect knowledge of other
agents’ local goals. It is important to stress that, in similar experiments,
the IA-MPPI-CVM which serves as the communication-free baseline in our
comparisons has already been demonstrated to outperform an optimization-based
Model Predictive Control (MPC) approach that relies on fixed predictions [3].

In the simulated experiments taking place in real sections of the canals of
Amsterdam, we randomize the initial positions and goals of four interacting
agents, all running the same algorithm. To challenge each method, we design
regions within which each agent’s start and goal position are randomly
initialized in a way that forces all four agents to interact in a narrow
section of the map. These high-interaction scenarios are discussed in Section
IV-B1.

For completeness, we also design experiments where agents’ starting and goal
positions are randomized across much larger spaces. In these experiments,
however, vessels don’t often interact and usually have larger free spaces to
avoid each other. These low-interaction scenarios are discussed in Section
IV-B2.

An example of experiments in low- and high-interaction scenarios is shown in
Fig. 5. The IA-MPPI plans with a time horizon $T$ of 100 time steps with step
size $\delta T=0.1s$ and $K=4500$ samples. Each method is evaluated on the
same set of randomly initialized experiments. For fairness, metrics such as
rule violations, goal displacement error, total traveled distance, and time
are only displayed for experiments that ended successfully with all methods.

TABLE II: Successes (Succ.), Deadlocks (Deadl.), Collisions (Coll.), Rule
Violations (Rule Viol.) and Goal Displacement Error (Goal DE) for all methods
in high-interaction scenarios per canal sections.

| Method | Succ. / Deadl. | Rule Viol. | Goal DE
---|---|---|---|---
| / Coll.
HG | IA-MPPI-CVM | 18 / 0 / 2 | 16 | 5.82 m
IA-MPPI-LBM (ours) | 19 / 1 / 0 | 16 | 5.26 m
IA-MPPI-w/comm. | 20 / 0 / 0 | 16 |
PG | IA-MPPI-CVM | 19 / 0 / 1 | 11 | 7.30 m
IA-MPPI-LBM (ours) | 19 / 1 / 0 | 5 | 4.11 m
IA-MPPI-w/comm. | 20 / 0 / 0 | 5 |
BG | IA-MPPI-CVM | 17 / 0 / 3 | 7 | 4.98 m
IA-MPPI-LBM (ours) | 20 / 0 / 0 | 4 | 4.13 m
IA-MPPI-w/comm. | 20 / 0 / 0 | 3 |

####  IV-B1 High-Interaction Scenario

The experiments in high-interaction scenarios are conducted in narrow
intersections in the Bloemgracht, Herengracht, and Prinsengracht. Since the
Amstel canal is very wide and the Open Crossing has no static map constraints,
it is difficult to generate experiments with high-interactions, and thus these
two maps are excluded from this experiment section. The results of the
experiments are summarized in Table II and Figure 6. It can be seen that in
these high-interaction scenarios, the LBM consistently outperforms the CVM in
terms of the goal displacement error (Goal DE). As a consequence, the motion
planning framework that estimates other agents’ local goals using predictions
from the LBM outperforms the framework that uses the CVM on all the metrics.
Moreover, we demonstrate our framework with the LBM has negligible performance
losses compared to the method with perfect communication.

Figure 6: This figure displays the distribution of the total traveled distance
and total traveled time of the vessels during the experiments in the high-
interaction scenarios. The results are displayed per map and for each method.

####  IV-B2 Low-Interaction Scenarios

Table III and Fig. 7 summarize the results in low-interaction scenarios. Note
that we here also test on the Open Crossing maps and the Amstel, which the LBM
has not previously seen in training. The results show that also when the start
and goal positions of all agents are randomly initialized over large areas,
our proposed communication-free framework with the LBM performs just as well
as the baseline with full communication, even in a map unseen in training.
However, perhaps unsurprisingly, the framework that approximates the local
goals with a CVM can also achieve the same performance as the framework with
full communication. Intuitively, in low-interaction scenarios where agents
mostly navigate straight to their goal, CVM is a reasonably good approximator.

Figure 7: This figure displays the distribution of the total traveled distance
and total traveled time by vessels during the experiments in the low-
interaction scenarios. The results are displayed per map and for each method.

###  IV-C Decoupled Prediction and Planning

The framework we proposed utilizes a Learning-Based Model (LBM) to predict
trajectories for obstacle agents and extract local goals while employing
Interaction-Aware Model Predictive Path Integral (IA-MPPI) for coupled
predictions and planning. To assess the advantages of this framework, we
compare it to a planner without interaction awareness (MPPI-LBM), which
decouples prediction and planning. Like other state-of-the-art methods, MPPI-
LBM treats the predicted future trajectories of obstacle agents as occupied
space and plans the ego agent’s motion without considering interaction
awareness. This approach reduces the system size and computational burden by
minimizing the space to be sampled. However, apart from this difference, MPPI-
LBM shares the same sampling strategy and cost function as the proposed IA-
MPPI-LBM. We conducted 100 low-interaction experiments across the Amstel,
Bloemgracht, Herengracht, Open Crossing, and Prinsengracht, comparing
different methods. Table IV presents the outcomes, including total successes,
deadlocks, collisions, and rule violations. Again, the proposed IA-MPPI-LBM
shows similar performances to the method with communication (IA-MPPI-w/comm).

However, MPPI-LBM exhibited a significantly lower success rate and a higher
number of rule violations. The LBM, while trained to be somewhat rule- and
interaction-aware in its predictions, occasionally struggles to capture
complex reciprocal collision avoidance maneuvers when agents are in close
proximity. This, combined with the motion planner’s unawareness of the ego
agent’s influence on other agents’ motion and their cooperation in collision
avoidance, often led the MPPI-LBM to wrongly assume that no feasible solution
existed. Consequently, this resulted in agents drifting into collisions due to
their large inertia.

TABLE III: Successes (Succ.), Deadlocks (Deadl.), Collisions (Coll.), Rule
Violations (Rule Viol.), and Goal Displacement Error (Goal DE) for all methods
in the various canal sections.

| Method | Succ. / Deadl. | Rule Viol. | Goal DE
---|---|---|---|---
| / Coll.
HG | IA-MPPI-CVM | 38 / 0 / 2 | 22 | 5.31 m
IA-MPPI-LBM (ours) | 37 / 0 / 3 | 26 | 5.42 m
IA-MPPI-w/comm. | 37 / 2 / 1 | 26 |
PG | IA-MPPI-CVM | 38 / 0 / 2 | 14 | 3.12 m
IA-MPPI-LBM (ours) | 38 / 0 / 2 | 13 | 2.94 m
IA-MPPI-w/comm. | 38 / 0 / 2 | 14 |
BG | IA-MPPI-CVM | 38 / 0 / 2 | 17 | 5.36 m
IA-MPPI-LBM (ours) | 39 / 0 / 1 | 20 | 4.93 m
IA-MPPI-w/comm. | 40 / 0 / 0 | 20 |
OC | IA-MPPI-CVM | 40 / 0 / 0 | 27 | 3.98 m
IA-MPPI-LBM (ours) | 40 / 0 / 0 | 29 | 4.04 m
IA-MPPI-w/comm. | 40 / 0 / 0 | 27 |
AM | IA-MPPI-CVM | 40 / 0 / 0 | 22 | 3.12 m
IA-MPPI-LBM (ours) | 39 / 1 / 0 | 18 | 3.49 m
IA-MPPI-w/comm. | 40 / 0 / 0 | 21 |

TABLE IV: Successes (Succ.), Deadlocks (Deadl.), Collisions (Coll.), and Rule
Violations (Rule Viol.) for the non-interactive MPPI and the IA-MPPI baseline
in low-interaction scenarios.

Method | Succ. / Deadl. / Coll. | Rule Viol.
---|---|---
MPPI-LBM. | 72 / 1 / 27 | 40
IA-MPPI-LBM (ours) | 97 / 1 / 2 | 30
IA-MPPI-w/comm. | 98 / 0 / 2 | 28

##  V Conclusions

In this paper, we introduced a framework that combines a learning-based
trajectory prediction model with Interaction Aware MPPI, enabling
decentralized and communication-free coupled prediction and planning. Our
experimental results demonstrated the superiority of our Learning-Based Model
(LBM) over the Constant Velocity Model (CVM) in accurately predicting the
trajectories of interacting vessels, even in unseen maps. Through simulated
experiments in Amsterdam’s canals, we showed that our motion planning
framework achieved comparable performance to a method with ground truth
knowledge of local goals, which was shown to outperform classical
optimization-based MPC approaches with decoupled prediction and planning in
previous work [3]. Additionally, we highlighted the limitations of the CVM in
tight environments with multiple interacting agents. Finally, by comparing our
approach with a non-interactive planner, we emphasized the advantages of
coupled planning and predictions.

## References

  * [1] W. Wang, T. Shan, P. Leoni, D. Fernández-Gutiérrez, D. Meyers, C. Ratti, and D. Rus, “Roboat II: A Novel Autonomous Surface Vessel for Urban Environments,” in _2020 IEEE/RSJ International Conference on Intelligent Robots and Systems (IROS)_ , Oct. 2020, pp. 1740–1747.
  * [2] P. Trautman and A. Krause, “Unfreezing the robot: Navigation in dense, interacting crowds,” in _2010 IEEE/RSJ International Conference on Intelligent Robots and Systems_ , Oct. 2010, pp. 797–803.
  * [3] L. Streichenberg, E. Trevisan, J. J. Chung, R. Siegwart, and J. Alonso-Mora, “Multi-Agent Path Integral Control for Interaction-Aware Motion Planning in Urban Canals,” in _2023 IEEE International Conference on Robotics and Automation (ICRA)_ , May 2023, pp. 1379–1385.
  * [4] M. G. Mohanan and A. Salgoankar, “A survey of robotic motion planning in dynamic environments,” _Robotics and Autonomous Systems_ , vol. 100, pp. 171–185, Feb. 2018.
  * [5] D. Fox, W. Burgard, and S. Thrun, “The dynamic window approach to collision avoidance,” _IEEE Robotics & Automation Magazine_, vol. 4, no. 1, pp. 23–33, Mar. 1997.
  * [6] J. van den Berg, M. Lin, and D. Manocha, “Reciprocal velocity obstacles for real-time multi-agent navigation,” in _2008 IEEE International Conference on Robotics and Automation (ICRA)_ , 2008, pp. 1928–1935.
  * [7] J. van den Berg, S. J. Guy, M. Lin, and D. Manocha, “Reciprocal n-Body Collision Avoidance,” in _Robotics Research_ , ser. Springer Tracts in Advanced Robotics, C. Pradalier, R. Siegwart, and G. Hirzinger, Eds. Springer, 2011, pp. 3–19.
  * [8] J. Alonso-Mora, A. Breitenmoser, M. Rufli, P. Beardsley, and R. Siegwart, “Optimal Reciprocal Collision Avoidance for Multiple Non-Holonomic Robots,” in _Distributed Autonomous Robotic Systems_ , ser. Springer Tracts in Advanced Robotics, A. Martinoli, F. Mondada, N. Correll, G. Mermoud, M. Egerstedt, M. A. Hsieh, L. E. Parker, and K. Støy, Eds. Berlin, Heidelberg: Springer, 2013, vol. 83, pp. 203–216.
  * [9] Y. Kuwata, M. T. Wolf, D. Zarzhitsky, and T. L. Huntsberger, “Safe Maritime Autonomous Navigation With COLREGS, Using Velocity Obstacles,” _IEEE Journal of Oceanic Engineering_ , vol. 39, no. 1, pp. 110–119, Jan. 2014, conference Name: IEEE Journal of Oceanic Engineering.
  * [10] M. Zanon, J. V. Frasch, M. Vukov, S. Sager, and M. Diehl, “Model Predictive Control of Autonomous Vehicles,” in _Optimization and Optimal Control in Automotive Systems_ , ser. Lecture Notes in Control and Information Sciences, H. Waschl, I. Kolmanovsky, M. Steinbuch, and L. del Re, Eds. Cham: Springer International Publishing, 2014, pp. 41–57.
  * [11] B. Paden, M. Čáp, S. Z. Yong, D. Yershov, and E. Frazzoli, “A survey of motion planning and control techniques for self-driving urban vehicles,” _IEEE Transactions on Intelligent Vehicles_ , vol. 1, no. 1, pp. 33–55, Apr. 2016, publisher: Institute of Electrical and Electronics Engineers Inc.
  * [12] C. E. Luis, M. Vukosavljev, and A. P. Schoellig, “Online Trajectory Generation With Distributed Model Predictive Control for Multi-Robot Motion Planning,” _IEEE Robotics and Automation Letters_ , vol. 5, no. 2, pp. 604–611, Apr. 2020.
  * [13] J. de Vries, E. Trevisan, J. van der Toorn, T. Das, B. Brito, and J. Alonso-Mora, “Regulations Aware Motion Planning for Autonomous Surface Vessels in Urban Canals,” in _2022 IEEE International Conference on Robotics and Automation_ , May 2022, pp. 3291–3297.
  * [14] H. Zhu and J. Alonso-Mora, “Chance-Constrained Collision Avoidance for MAVs in Dynamic Environments,” _IEEE Robotics and Automation Letters_ , vol. 4, no. 2, pp. 776–783, 2019.
  * [15] H. Zhu, F. M. Claramunt, B. Brito, and J. Alonso-Mora, “Learning Interaction-Aware Trajectory Predictions for Decentralized Multi-Robot Motion Planning in Dynamic Environments,” _IEEE Robotics and Automation Letters_ , vol. 6, no. 2, pp. 2256–2263, Apr. 2021.
  * [16] D. Fridovich-Keil, A. Bajcsy, J. F. Fisac, S. L. Herbert, S. Wang, A. D. Dragan, and C. J. Tomlin, “Confidence-aware motion prediction for real-time collision avoidance1,” _The International Journal of Robotics Research_ , vol. 39, no. 2-3, pp. 250–265, Mar. 2020, publisher: SAGE Publications Ltd STM.
  * [17] E. Schmerling, K. Leung, W. Vollprecht, and M. Pavone, “Multimodal Probabilistic Model-Based Planning for Human-Robot Interaction,” in _2018 IEEE International Conference on Robotics and Automation (ICRA)_ , May 2018, pp. 3399–3406, iSSN: 2577-087X.
  * [18] C. Mavrogiannis, F. Baldini, A. Wang, D. Zhao, P. Trautman, A. Steinfeld, and J. Oh, “Core Challenges of Social Robot Navigation: A Survey,” _ACM Transactions on Human-Robot Interaction_ , vol. 12, no. 3, pp. 36:1–36:39, Apr. 2023.
  * [19] P. Trautman, J. Ma, R. M. Murray, and A. Krause, “Robot navigation in dense human crowds: Statistical models and experimental studies of human–robot cooperation,” _The International Journal of Robotics Research_ , vol. 34, no. 3, pp. 335–356, Mar. 2015, publisher: SAGE Publications Ltd STM.
  * [20] K. Liu, N. Li, H. E. Tseng, I. Kolmanovsky, and A. Girard, “Interaction-Aware Trajectory Prediction and Planning for Autonomous Vehicles in Forced Merge Scenarios,” _IEEE Transactions on Intelligent Transportation Systems_ , vol. 24, no. 1, pp. 474–488, Jan. 2023, conference Name: IEEE Transactions on Intelligent Transportation Systems.
  * [21] D. Sadigh, S. Sastry, S. A. Seshia, and A. D. Dragan, “Planning for Autonomous Cars that Leverage Effects on Human Actions,” in _Proceedings of Robotics: Science and Systems_ , Jun. 2016.
  * [22] G. Williams, P. Drews, B. Goldfain, J. M. Rehg, and I. A. Theodorou, “Information-Theoretic Model Predictive Control: Theory and Applications to Autonomous Driving,” _IEEE Transactions on Robotics_ , vol. 34, no. 6, pp. 1603–1622, Dec. 2018.
  * [23] B. F. d. Brito, H. Zhu, W. Pan, and J. Alonso-Mora, “Social-VRNN: One-Shot Multi-modal Trajectory Prediction for Interacting Pedestrians,” in _Proceedings of the 2020 Conference on Robot Learning_. PMLR, Oct. 2021, pp. 862–872.
  * [24] J. Chung, K. Kastner, L. Dinh, K. Goel, A. C. Courville, and Y. Bengio, “A Recurrent Latent Variable Model for Sequential Data,” in _Advances in Neural Information Processing Systems_ , vol. 28. Curran Associates, Inc., 2015\.
  * [25] T. Tieleman and G. Hinton, “Lecture 6.5-rmsprop: Divide the gradient by a running average of its recent magnitude,” _COURSERA: Neural networks for machine learning_ , vol. 4, no. 2, pp. 26–31, 2012.