# Metric Learning as a Service with Covariance Embedding

Imam Mustafa Kamal, Hyerim Bae, and Ling Liu IM. Kamal was with the Institute
of Intelligent Logistics Big Data, Pusan National University, 30-Jan-jeon
Dong, Geum-Jeong Gu, 609-753, Busan, South Korea, Email:
<EMAIL_ADDRESS>H. Bae was with a major in Industrial Data Science &
Engineering, Department of Industrial Engineering, Pusan National University,
30-Jan-jeon Dong, Geum-Jeong Gu, 609-753, Busan, South Korea, Email:
<EMAIL_ADDRESS>(Corresponding author). L. Liu was with the College of
Computing, Georgia Institute of Technology, 266 Ferst Drive, Atlanta,
30332-0765, GA, USA, Email<EMAIL_ADDRESS>
###### Abstract

With the emergence of deep learning, metric learning has gained significant
popularity in numerous machine learning tasks dealing with complex and large-
scale datasets, such as information retrieval, object recognition and
recommendation systems. Metric learning aims to maximize and minimize inter-
and intra-class similarities. However, existing models mainly rely on distance
measures to obtain a separable embedding space and implicitly maximize the
intra-class similarity while neglecting the inter-class relationship. We argue
that to enable metric learning as a service for high-performance deep learning
applications, we should also wisely deal with inter-class relationships to
obtain a more advanced and meaningful embedding space representation. In this
paper, a novel metric learning is presented as a service methodology that
incorporates covariance to signify the direction of the linear relationship
between data points in an embedding space. Unlike conventional metric
learning, our covariance-embedding-enhanced approach enables metric learning
as a service to be more expressive for computing similar or dissimilar
measures and can capture positive, negative, or neutral relationships.
Extensive experiments conducted using various benchmark datasets, including
natural, biomedical, and facial images, demonstrate that the proposed model as
a service with covariance-embedding optimizations can obtain higher-quality,
more separable, and more expressive embedding representations than existing
models.

###### Index Terms:

Metric learning, semantic similarity, siamese network, covariance metric, deep
learning.

This work has been submitted to the IEEE for possible publication. Copyright
may be transferred without notice, after which this version may no longer be
accessible.

##  1 Introduction

Deep learning has penetrated many businesses, science, and engineering domains
with a remarkable generalization performance, including self-driving cars,
smart manufacturing, smart cities, healthcare diagnostics, spam detection, and
cloud resource management. The notion of a distance or similarity measure
plays a crucial role in deep-learning algorithms. Conventionally, standard
distance metrics, such as Euclidean, cosine, and Mahalanobis, are applied by
incorporating a priori knowledge of the domain. However, it is often
challenging to devise metrics that are pertinent to particular data and tasks
of interest.

Metric learning aims to automatically construct task-specific distance or
similarity metrics using weakly supervised data in a machine-learning fashion.
The outcome can be used to perform various tasks, such as $k$-nearest neighbor
($k$-NN), clustering, and information retrieval. Moreover, metric learning can
offer a natural solution to various machine learning problems because it has
numerous benefits, including robustness to noisy data, a high generalization
to unseen categories, the capability of a dimensionality-reduction model,
reliability as a feature extraction model (for fine-tuning classification),
and the ability to operate with small sample datasets [1]. Accordingly,
research into metric learning has received increasing attention over the past
few years and has been applied to many different fields, such as information
retrieval [2], recommender systems [3], social media mining [4], face
recognition [5], and speech recognition [6], to name a few. Siamese, Triplet,
and N-pair networks are prominent metric learning frameworks. A Siamese
network is trained using two tandem networks with two different input data to
compute their similarity. A triplet network is trained using three networks to
define the similarity between three images, two of which will be similar
(anchor and positive samples), and the third will be unrelated (a negative
example). The N-pair network is trained using the cosine similarity to
calculate the pairwise distance of N samples, where N $>$ 3\. Each of these
networks is able to project high-dimensional data into a low-dimensional
embedding space while preserving data separability. The samples of different
categories must be dissimilar (far from each other); otherwise, the samples
from the same categories must have a similar latent representation. Defining
the closeness and remoteness among data or samples in complex data, such as
natural images, is a non-trivial task. Accordingly, the similarity metric is
critical for obtaining a meaningful embedding space representation.

When similarity measures are not given a priori, although a generic function,
such as Euclidean distance, can be adopted, doing so frequently produces
unsatisfactory results [7]. If two data vectors have no attribute values in
common, they may have a smaller distance than the other pair of data vectors
containing the same attribute values. In addition, as one of the main
drawbacks of the cosine similarity, the magnitude of the vectors is not
considered, and only their direction is considered [8]. In practice, this
means that the differences in values are not fully considered. To the best of
our knowledge, existing metric learning models only explicitly maximize intra-
class similarity and implicitly neglect the interclass relationship. Thus, the
latent space representation is possibly less meaningful because it cannot
capture inter-class connections. Moreover, existing metric learning models
such as Siamese, Triplet, and N-pair networks often suffer from a low
performance under a particular condition. The Siamese network results in fewer
semantic outcomes because it neglects the inter-class relationship. The
triplet network can diverge when the positive and negative samples have
similar features. However, the N-pair network requires a large number of
batches to obtain an appropriate result that is computationally expensive.

In this paper, we propose a novel metric learning method called CovNet with
three original contributions. First, CovNet employs covariance to determine
the relationship between two pairs of inputs. The covariance signifies the
direction of the linear relationship between the two vectors. By direction, we
are referring to whether the variables are directly or inversely proportional
to each other. Increasing the value of one variable might have a positive or
negative impact on the value of the other variable. In addition, unlike the
cosine similarity, the covariance subtracts the means before taking the dot
product, making it invariant to shifts. Second, we provide three types of
mapping functions to embody the covariance embedding: inner-class mapping
(IM), inner-class with semi-inter-class mapping (ISIM), and inner- and inter-
class mapping (IIM). Finally, we conducted extensive experiments on several
benchmark datasets and demonstrated that CovNet outperforms existing deep
metric learning models, such as Siamese, Triplet, and N-pair networks. The
covariance metric is more expressive than the Euclidean and cosine similarity
metrics because it can capture three possible relationships between two
variables: a positive, neutral, or negative correlation. Compared with
existing metric learning models, as it can obtain accurate inner and inter-
class relationships, we show that CovNet is more semantic and expressive.

The remainder of this paper is organized as follows. Section 2 reviews the
literature on metric learning models. Section 3 presents the proposed method.
An extensive experiment along with a discussion is presented in Section 4, and
finally, we outline some conclusions and future directions for this research
in Section 5.

##  2 Related works

Metric learning is an approach based directly on a similarity measure that
aims to establish a relationship between data points. The relation provides a
semantic similarity such that close data points will be considered similar,
whereas remote data points will be considered dissimilar. It can also be
categorized as a dimensionality reduction because it maps high-dimensional
data points to a low-dimensional space while preserving its separable features
[9]. Hence, it can also be useful for large-scale data or multimedia
applications, which are ubiquitous in the modern era. The popularity of metric
learning emerged in 2002 with the pioneering work of Xing et al. [10], who
formulated it as a convex optimization problem. Therefore, several metrics
learning applications in various domains have been introduced by some
scholars, such as medicine, security, social media mining, speech recognition,
information retrieval, recommender systems, and computer vision. In terms of
business, the usage of metric learning as a service platform, such as
recommendation system and image search similarity, makes consumption and
buying decisions more effective and the user experience more comfortable by
recommending only relevant items. It allows service providers to predict the
customer’s usage behavior by reusing professional data mining services.

According to the type of supervision applied during training, the metric
learning model mainly falls into two main categories: supervised and weakly
supervised learning. In supervised learning, the model is trained using data,
where each sample has label information as a standard discriminative model or
classifier. In summary, the model learns to map data points based on the
label. Thus, the data points in the same label are considered similar and are
placed in a close space representation; otherwise, they are mapped in a remote
space representation. Popular supervised metric learning models, such as a
linear discriminant analysis (LDA), margin maximizing discriminant analysis
(MMDA), learning with side information (LSI), relevant component analysis
(RCA), and neighborhood component analysis (NCA). Because supervised metric
learning directly employs a label in defining the similarity among the data
points, the embedding space becomes rigid to the label, is prone to an
overfitting similarly to a standard (deep learning) classifier, has difficulty
capturing the semantic relationship across data points lying within different
categories, and has problems dealing with new (unlisted) categories. In weakly
supervised learning, the model has access to a set of data points with
supervision at the tuple level (particularly pairs, triplets, quadruplets, or
N pairs of data points). Thus, its mechanism is called weakly supervised
because the model does not directly map the similarity of data points based on
the corresponding label, i.e., the prominent models of weakly supervised
metric learning such as Siamese, Triplet, and N-pair networks.

The Siamese network algorithm was first introduced by Bromley et al. to two
handwritten signatures in 1994. The model was defined as a binary classifier
to distinguish whether the pair of data points belonged to the same class or
originated in a different class. Owing to its simplicity and applicability,
the Siamese network has become the most widely used network by scholars. Zhang
et al. [11] developed a content-based image retrieval framework using a
Siamese network. In addition, Qiao et al. [12] employed a Siamese network to
define user identity linkages through web browsing. Moreover, Yang et al. [13]
presented terahertz image verification using a symmetric Siamese network. As
represented by a binary classification problem, the Siamese network solely
overlooks the inner-class relationship. Thus, the connection between inter-
class data points can be less semantic. To address this issue, the triplet
network was introduced by Schroff et al. [14] and was originally devised for
face recognition. Herein, the data points are mapped in a triplet consisting
of an anchor and positive and negative data points; accordingly, the model can
learn the similarity of all possible pairs of data points. Some previous
triplet network applications are as follows. Zhang et al. [15] adopted a
triplet loss for remote-sensing image retrieval. Boutros et al. [16]
implemented a triplet loss for mask-image recognition. Lue et al. [17] also
devised a visible-thermal-person re-identification framework using triplet
loss. Learning the similarity in a triplet fashion is a promising technique
for enhancing the semantic similarity of inter-class structures. Nonetheless,
this mechanism is challenging to implement. Many scholars show that a triplet
loss often suffers from slow convergence, partially because they employ only
one negative example while not interacting with the other negative classes in
each update [18]. Moreover, if some positive and negative data points are
similar, the performance can be degraded. Consequently, Sohn and Kihyuk
introduced N-pair networks to solve this problem. In an N-pair network, the
model learns the pairwise distance or similarity between all possible pairs in
N data points. Commonly, N is represented as the number of batches. Therefore,
the network is updated by calculating the pairwise similarity within a batch.
Among the latest N-pair network applications, Chen and Deng [19] utilized
N-pair loss for image retrieval and clustering. In addition, Pal et al. [20]
employed N-pair loss in biomedical image classification. Moreover, Espejo et
al. [21] extended the N-pair loss with a $(C_{N,2}+1)$-pair loss function for
keyword spotting. Indeed, an N-pair network is a promising metric learning
model for obtaining reliable semantic representations. Nevertheless, the
N-pair loss tends to have inferior performance when the number of batches (N)
is small. Moreover, employing a large number of batches requires a
significantly high computational cost.

In metric learning, the distance or similarity measure also plays a crucial
role in defining the embedding space. The base network of a metric learning
model commonly produces a vector. Its value must represent the original input
both accurately and semantically. The standard distance measure used to define
the distance or similarity between vectors is Euclidean and cosine.
Nevertheless, these often fail to capture the idiosyncrasies of the data of
interest [22]. Thus, Norouzi et al. [7] introduced a framework for learning a
broad class of binary hash functions based on a triplet ranking loss designed
to preserve the relative similarity. Zhang et al. [8] proposed a spherical
embedding constraint (SEC) to regularize the distribution of the norms.
Moreover, some scholars have recently used distance measures to accomplish a
better semantic embedding space representation. Xu et al. [23] used a bi-level
distance metric to enhance the similarity accuracy, and Ye et al. [24]
incorporated multi-metric learning to capture multi-perspective data
relationships. Based on the above approaches, we extended this research in
another direction. Unlike the previous studies, we devised an expressive but
straightforward method for determining the relationship between embedding
vectors using covariance metrics.

Figure 1: Covariance between two embedding vectors capturing negative,
neutral, and positive correlations.

$cov(z,z^{\prime})=\frac{\sum_{i=1}^{s}(z_{i}-\bar{z})(z_{i}^{\prime}-\bar{z}^{\prime})}{s-1}$
(1)
$\overrightarrow{cov}(z,z^{\prime})=(z_{i}-\bar{z})(z_{i}^{\prime}-\bar{z}^{\prime})$
(2)
$\begin{split}\sum\overrightarrow{cov}(z,z^{\prime})\begin{cases}>0,&\textnormal{positive
correlation}\\\ \approx 0,&\textnormal{uncorrelated}\\\
<0,&\textnormal{negative correlation}\\\ \end{cases}\end{split}$ (3)

##  3 Methodology

Rather than Euclidean distance and cosine similarity, the covariance metric is
more expressive because it can simultaneously capture three possible relations
between two variables, as depicted in Fig. 1. Given two random variables, $z$
and $z^{\prime}$, the relationship between them becomes negative when the
$cov(z,z^{\prime})<0$, neutral or uncorrelated when $cov(z,z^{\prime})\approx
0$, and positive when $cov(z,z^{\prime})>0$. This value can be extracted from
the original covariance equation denoted in Eq. 1, where $s$ represents the
total number of variables, $\bar{z}$ is the mean value of $z$ derived from
$\bar{z}=\frac{\sum_{i=1}^{s}z_{i}}{s}$ and $\bar{z}^{\prime}$ is equal to
$\bar{z}^{\prime}=\frac{\sum_{i=1}^{s}z_{i}^{\prime}}{s}$. Therefore, unlike
Euclidean and cosine similarity metrics, the covariance metric is invariant to
shifts. If $z$ and $z^{\prime}$ is a vector with $q$ dimension. The covariance
between them can be formulated as
$cov(z,z^{\prime})=\sum_{i=1}^{q}(z_{i}-\bar{z})(z_{i}^{\prime}-\bar{z}^{\prime})$,
where $cov(z,z^{\prime})$ is a scalar value. We represent the covariance
vector $\overrightarrow{cov}(z,z^{\prime})$, as shown in Eq. 2, instead of a
scalar because we cannot directly minimize $cov(z,z^{\prime})$ using an
optimization algorithm. Note that, the summation of
$\overrightarrow{cov}(z,z^{\prime})$ will be positive, neutral (uncorrelated),
or negative, as denoted in Eq. 3.

###  3.1 Problem definition

Let us denote $X=\\{x_{0},x_{1},...,x_{N-1}\\}$ as real value $N$ with a
$p$-dimensional space, $X\in\mathbb{R}^{p}$, and
$Y=\\{y_{0},y_{1},...,y_{N-1}\\}$ as a label of $X$. The embedding network
($F$) will project $X$ into $q$-dimensional space ($z=F(x)$), where
$z\in\mathbb{R}^{q}$ and $q\ll p$. We define a pair input
$\\{x,x^{\prime}\\}$, where $x=x_{i}$, $x^{\prime}=x_{j}$, and $i\neq j$.
$z^{\prime}$ can also be obtained after $x^{\prime}$ is processed using $F$.
Subsequently, we can obtain the covariance value, represented as a vector,
between $z$ and $z^{\prime}$ using Eq. 2. By using a mapping function (IM,
IIM, or ISIM), we can obtain a label of $\mathcal{X}=\\{x,x^{\prime}\\}$,
which is represented as $\mathcal{Y}$. Note that $Y$ is a raw label whereas
$\mathcal{Y}$ represents the final label generated by the mapping function.
The covariance vector $\overrightarrow{cov}(z,z^{\prime})$ is directly
projected to become $\mathcal{Y}$ through the function $G$
($G(\overrightarrow{cov}(z,z^{\prime}))=\mathcal{Y}$), where $G$ represents a
dense layer with the number of units $|\mathcal{Y}|$. When the label is binary
$|\mathcal{Y}|=1$, and when the label is categorical $|\mathcal{Y}|$ is equal
to the total number of categories. Unlike the Euclidean and cosine similarity
metrics, the covariance operation may result in a range value of
$[-\infty,\infty]$; hence we cannot directly minimize it by using an
optimization algorithm. Thus, we employ a classification approach to
implicitly obtain the covariance value between $z$ and $z^{\prime}$.

(a) IM

(b) IIM

(c) ISIM

Figure 2: Mapping function comparison: IM, ISIM, and IIM. Herein there are
four categories $\\{y_{0},y_{1},y_{2},y_{3}\\}$ and all possible generated
labels by mapping function are denoted as $\mathcal{Y}$.

###  3.2 Data processing and network model

Before the network model learns the data, we apply a data processing technique
which is defined as a mapping function. In this study, we present three types
of mapping functions to determine a pair of inputs
($\mathcal{X}=\\{x,x^{\prime}\\}$) and its label ($\mathcal{Y}$), namely, IM,
IIM, and ISIM, as illustrated in Fig. 2. Note that $\mathcal{X}$ is an image
and $\mathcal{Y}$ is a text. In the IM, the data are paired if they lie in the
same class or category. For instance, in the CIFAR-10 dataset, if $x$ is a
bulldog, $x^{\prime}$ can be a puddle
($\\{x=\text{bulldog},x=\text{puddle}\\},\mathcal{Y}=\text{dog}$) and if $x$
is Persian, $x^{\prime}$ can be a ragdoll
($\\{x=\text{Persian},x=\text{ragdoll}\\},\mathcal{Y}=\text{cat}$), and so on.
Thus, suppose $C$ is the total number of categories or classes in a dataset,
and the total number of classes in IM ($n_{class}$) is equal to the total
number of categories in a dataset ($n_{class}$ = 10 in the CIFAR-10 dataset).
The detailed IM mapping algorithm is presented in Algorithm 1. It indicates
that the IM explicitly maximizes the inner-class similarity. Nevertheless, it
can also implicitly minimize the inter-class similarity after the convergence
(when all data points successfully fit with their corresponding categories).
In the IIM, the data are paired based on all pair combinations among the
categories. For example, ($\\{x$=dog, $x^{\prime}$=another dog$\\}$,
$\mathcal{Y}=\text{dog-dog}$), ($\\{x$=dog, $x^{\prime}$=cat$\\}$,
$\mathcal{Y}=\text{dog-cat}$), ($\\{x$=dog, $x^{\prime}$=airplane$\\}$,
$\mathcal{Y}=\text{dog-airplane}$), ($\\{x$=dog, $x^{\prime}$=automobile$\\}$,
$\mathcal{Y}=\text{dog-automobile}$), ($\\{x$=dog, $x^{\prime}$=truck$\\}$,
$\mathcal{Y}=\text{dog-truck}$), ($\\{x$=dog, $x^{\prime}$=ship$\\}$,
$\mathcal{Y}=\text{dog-ship}$), ($\\{x$=dog, $x^{\prime}$=deer$\\}$,
$\mathcal{Y}=\text{dog-deer}$), ($\\{x$=dog, $x^{\prime}$=bird$\\}$,
$\mathcal{Y}=\text{dog-bird}$), $\\{x$=dog, $x^{\prime}$=frog$\\}$,
$\mathcal{Y}=\text{dog-frog}$), ($\\{x$=dog, $x^{\prime}$=horse$\\}$,
$\mathcal{Y}=\text{dog-horse}$), and so forth for other categories. The total
number of classes in IIM can be denoted as $n_{class}=C+\binom{C}{2}$. Because
all possible combinations in the class categories are explicitly observed, IIM
can learn to maximize and minimize inner- and inter-class similarities
simultaneously. The detailed IIM mapping algorithm is presented in Algorithm
2. In addition, compared with IM and ISIM, IIM mapping can be robust in
extracting the semantic similarity in both inner-class and inter-class
relationships because it learns similarity of all combinations of categories.
However, because the number of pair combinations can be significantly
increased when the number of categories in a dataset is high, it will become a
classification with a large number of categories. In which, this is still an
active area of research in machine learning [25]. In ISIM, the data are paired
in the same as in the original Siamese network mapping. Therefore, there are
only two possible labels (binary) for ISIM, namely, a pair containing the same
class category (e.g. $\\{x$=dog, $x^{\prime}$=another dog$\\},\mathcal{Y}=1$),
and a pair consisting of a different class category (e.g. $\\{x$=dog,
$x^{\prime}$=cat$\\},\mathcal{Y}=0$). Because it is assumed that the
similarity between a dog and a cat is equal to that between a dog and an
airplane, the capturing of the inter-class relationship may be explicitly
neglected. The detailed ISIM mapping is presented in Algorithm 3.

Figure 3: Covariance network (CovNet). Herein the input $\\{x,x^{\prime}\\}$
and output $\mathcal{Y}$ are obtained from the mapping function.

Input: $X=${$x_{0},x_{1},...,x_{N-1}$}, $Y=${$y_{0},y_{1},...,x_{N-1}$}

Output: $\mathcal{X}$, $\mathcal{Y}$

$CG$ ;

// In $CG$, $x$ is grouped and indexed by $y$

$\mathcal{X}$ = [] ;

// initialize $\mathcal{X}$ as an empty list

$\mathcal{Y}$ = [] ;

// initialize $\mathcal{Y}$ as an empty list

1 for  _$i\leftarrow 0$ to $N-1$_ do

2 do

$x_{j}\leftarrow$ random($CG_{y_{i}}$) ;

// get random $x_{j}$ in $CG$, where $y_{j}=y_{i}$

3

4 while  _$x_{i}=x_{j}$_ ;

$\mathcal{X}.$append($\\{x_{i},x_{j}\\}$) ;

// pair $x$ and $x^{\prime}$

$\mathcal{Y}.$append($y_{i}$) ;

// label, note that $y_{i}=y_{j}$

5

6 end for

Algorithm 1 IM

Input: $X=${$x_{0},x_{1},...,x_{N-1}$}, $Y=${$y_{0},y_{1},...,x_{N-1}$}

Output: $\mathcal{X}$, $\mathcal{Y}$

$CG$ ;

// In $CG$, $x$ is grouped and indexed by $y$

$\mathcal{X}$ = [] ;

// initialize $\mathcal{X}$ as an empty list

$\mathcal{Y}$ = [] ;

// initialize $\mathcal{Y}$ as an empty list

1 for  _$i\leftarrow 0$ to $N-1$_ do

2 for  _$j\leftarrow 0$ to $C-1$_ do

3 if  _$y_{i}=y_{j}$_ then

4 do

$x_{j}\leftarrow$ random($CG_{y_{i}}$) ;

// get random $x_{j}$ in $CG$,where $y_{j}=y_{i}$

5

6 while  _$x_{i}=x_{j}$_ ;

7 else

$x_{j}\leftarrow$ random($CG_{y_{j}}$) ;

// get random $x_{j}$ in $CG$, where $y_{j}\neq y_{i}$

8

9 end if

$\mathcal{X}.$append($\\{x_{i},x_{j}\\}$) ;

// pair $x$ and $x^{\prime}$

10 if  _$i\leq j$_ then

$\mathcal{Y}.$append(str($y_{i}$)+str($y_{j}$)) ;

// label as a string

11

12 else

$\mathcal{Y}.$append(str($y_{j}$)+str($y_{i}$)) ;

// label as a string in ascending order only

13

14 end if

15

16 end for

17

18 end for

Algorithm 2 IIM

Input: $X=${$x_{0},x_{1},...,x_{N-1}$}, $Y=${$y_{0},y_{1},...,x_{N-1}$}

Output: $\mathcal{X}$, $\mathcal{Y}$

$CG$ ;

// In $CG$, $x$ is grouped and indexed by $y$

$\mathcal{X}$ = [] ;

// initialize $\mathcal{X}$ as an empty list

$\mathcal{Y}$ = [] ;

// initialize $\mathcal{Y}$ as an empty list

1 for  _$i\leftarrow 0$ to $N-1$_ do

/* add a matching pair */

2 do

$x_{j}\leftarrow$ random($CG_{y_{i}}$) ;

// get random $x_{j}$ in $CG$, where $y_{j}=y_{i}$

3

4 while  _$x_{i}=x_{j}$_ ;

$\mathcal{X}.$append($\\{x_{i},x_{j}\\}$) ;

// pair $x$ and $x^{\prime}$

$\mathcal{Y}.$append(1) ;

// a matching label, $\mathcal{Y}=1$

/* add a non-matching pair */

5 do

$y_{j}\leftarrow$ random_int($0,C-1$) ;

// to make sure $y_{i}\neq y_{j}$

6

7 while  _$y_{j}=y_{i}$_ ;

$x_{j}\leftarrow$ random($CG_{y_{j}}$) ;

// get random $x_{j}$ in $CG$ where $y_{j}\neq y_{i}$

$\mathcal{X}.$append($\\{x_{i},x_{j}\\}$) ;

// pair $x$ and $x^{\prime}$

$\mathcal{Y}.$append(0) ;

// a non-matching label, $\mathcal{Y}=0$

8

9 end for

Algorithm 3 ISIM

The data format obtained from the processing technique is a pair input,
$\\{x,x^{\prime}\\}$ with a single label, $\mathcal{Y}$. They were trained in
a supervised manner using a network model called CovNet. Thus,
$\\{x,x^{\prime}\\}$ is processed to estimate the label $\mathcal{\hat{Y}}$.
The overall architecture of CovNet is illustrated in Fig. 3. We define $F$ as
an embedding network represented as a convolutional neural network (CNN). We
visually present two embedding networks. However, because they share their
weights, they are physically a single network ($F$). As a standard CNN, the
model consists of convolutional, max-pooling, batch normalization, and dropout
layers. In the last layer, we employ a global average-pooling layer to capture
the extracted feature generated from the convolutional layer. Subsequently, a
dense layer with L2 normalization (L2 norm.) is applied to obtain the
embedding vector ($z$). The L2 norm. has also been commonly used in previous
metric learning models because it can result in stability during training.
Thus, because $z$ and $z^{\prime}$ are already normalized through the L2 norm,
if the $z$ and $z^{\prime}$ are centered (have zero means), it (Eq. 2) will
have the same result as the cosine similarity. In the tail layer, we
incorporate it with a dense layer ($G$) using a softmax ($S$) or sigmoid
function to determine the class probability.

###  3.3 Learning mechanism

CovNet was trained as a standard neural network with a back-propagation
algorithm. The learning mechanism of CovNet is described in Algorithm 4. The
inputs are $\mathcal{X}$ and $\mathcal{Y}$, the values of which are obtained
from the mapping function used (IM, IIM, or ISIM, as defined in Algorithms 1,
2, and 3, respectively). Note that each of the samples in $\mathcal{X}$
consists of tuples $\\{x,x^{\prime}\\}$, and $\mathcal{Y}$ is represented as a
one-hot vector (in IM and IIM mapping) and a binary vector (in ISIM mapping).
The output of the algorithm is the optimum embedding network model
($F^{*}_{{\phi}}$). Initially, $F_{\phi}$ produces $z$ and $z^{\prime}$ given
input $x$ and $x^{\prime}$, respectively. Subsequently, we can obtain
$\overrightarrow{cov}(z,z^{\prime})$, which is formulated in Eq. 2. Then,
$G_{\psi}$ will classify whether $\overrightarrow{cov}(z,z^{\prime})$ agrees
on the same label. A predefined loss_function calculates the discrepancy loss
between $\mathcal{Y}$ and $\mathcal{\hat{Y}}$. If we employ IM and IIM as our
mapping function, the loss_function is the categorical cross-entropy
($\mathcal{L}_{CE}$) formulated in Eqs. 4 and 5. However, if we utilize ISIM
as our mapping function, the loss_function is a binary cross-entropy
($\mathcal{L}_{BE}$), denoted in Eq. 6. For simplicity, the model parameters
$F$($\phi$) and $G$ ($\psi$) are wrapped as a single parameter $\Omega$, and
will be updated using a gradient-descent based optimization algorithm. All
aforementioned processes were repeated until reaching the predefined number of
epochs ($NoEpoch$). Note that, in practice, the model is trained and updated
using a predefined number of batches.

$S(\mathcal{\hat{Y}}_{i})=\frac{{\exp}^{\mathcal{\hat{Y}}_{i}}}{\sum_{i=1}^{n_{class}}\exp^{\mathcal{\hat{Y}}_{i}}}$
(4)
$\mathcal{L}_{CE}=-\sum_{i=1}^{n_{class}}{\mathcal{Y}}_{i}\log(S(\mathcal{\hat{Y}}_{i}))$
(5)
$\mathcal{L}_{BE}=-[\mathcal{Y}\log(\mathcal{\hat{Y}})+(1-\mathcal{Y})\log(1-\mathcal{\hat{Y}})]$
(6)

Input: $\mathcal{X}$ =
{$\\{x_{0},x^{\prime}_{0}\\},\\{x_{1},x^{\prime}_{1}\\},...,\\{x_{N-1},x^{\prime}_{N-1}\\}$},
$\mathcal{Y}$ = {$\mathcal{Y}_{0},\mathcal{Y}_{1},...,\mathcal{Y}_{N-1}$}

Output: $F^{*}_{{\phi}}$

1 for  _$i\leftarrow 1$ to $NoEpoch$_ do

$z=F_{{\phi}_{i}}(x)$ ;

// embedding network extracts $z$

$z^{\prime}=F_{{\phi}_{i}}(x^{\prime})$ ;

// embedding network extracts $z^{\prime}$

2

$\overrightarrow{cov}(z,z^{\prime})=$ Eq. 2 ;

// extract covariance vector

$\mathcal{\hat{Y}}=G_{{\psi}_{i}}(\overrightarrow{cov}(z,z^{\prime}))$;

// estimate $\mathcal{Y}$

$\mathcal{L}_{{\Omega}_{i}}({\phi}_{i},{\psi}_{i})=$loss_function($\mathcal{Y},\mathcal{\hat{Y}}$)
;

// calculate supervised loss betwen $\mathcal{Y}$ and $\mathcal{\hat{Y}}$

${\Omega}_{i}\leftarrow{\Omega}_{i}-\eta\nabla_{{\Omega}_{i}}{\mathcal{L}}_{{\Omega}_{i}}$
;

// update model parameter, $\Omega=\\{\phi,\psi\\}$

3

4 end for

Algorithm 4 training procedure of CovNet Figure 4: Dataset overview

##  4 Experiments

This section presents the performance evaluation of our proposed model
(CovNets), followed by a comparison and discussion between CovNets and
existing metric learning models.

TABLE I: Embedding network ($F$) of CovNet, Siamese, Triplet, and N-pair
networks in all datasets.

CIFAR-10:
---
$Conv(64,3,Relu)-BN-Conv(64,3,Relu)-MP(3)-BN-DO(0.5)-Conv(128,3,Relu)-BN-$
$Conv(128,3,Relu)-MP(3)-BN-DO(0.5)-Conv(256,3,Relu)-BN-Conv(256,3,Relu)-$
$Conv(256,3,Relu)-MP(3)-BN-DO(0.5)-AP-Dense(256,Relu)-BN-
DO(0.5)-Dense(100,Tanh)-L2Norm.$
Food-11 & Colorectal:
$Conv(64,3,Relu)-BN-Conv(64,3,Relu)-MP(3)-BN-DO(0.5)-Conv(128,3,Relu)-BN-$
$Conv(128,3,Relu)-MP(3)-BN-DO(0.5)-Conv(256,3,Relu)-BN-Conv(256,3,Relu)-$
$Conv(512,3,Relu)-MP(3)-BN-DO(0.5)-AP-Dense(256,Relu)-BN-
DO(0.5)-Dense(100,Tanh)-L2Norm.$
Adience:
$ZP(5)-FaceNet-Dense(512,Relu)-BN-DO(0.25)-Dense(100,Tanh)-L2Norm.$
Yale:
$FaceNet-Dense(512,Relu)-BN-DO(0.25)-Dense(100,Tanh)-L2Norm.$

###  4.1 Experimental setting

We evaluated our model on various well-known datasets including those with
natural images (CIFAR-10 [26] and Food-11 [27]), biomedical images (Colorectal
[28]), and facial images (Adience [29] and Yale [30]). Sample images from each
dataset are illustrated in Fig. 4. The CIFAR-10 dataset consists of 50,000
$32\times 32$ pixels color training images and 10,000 test images, labeled
over ten categories: airplanes, automobiles, birds, cats, deer, dogs, frogs,
horses, ships, and trucks. We obtained 10% of the training set for the
validation set. Food-11 is a dataset containing 16,643 food images grouped
into 11 major food categories: bread, dairy products, desserts, eggs, fried
foods, meat, noodles/pasta, rice, seafood, soup, and vegetable/fruit. The
original images had various dimensions and were split into 9,866 training,
3,430 validation, and 3,347 testing sets. For simplicity, we uniformly resized
the image dimensions to 128$\times$128 pixels. We employ zero-contrast
normalization and ZCA whitening as data preprocessing on Food-11 and CIFAR-10.
Colorectal cancer is a biomedical dataset containing histological tiles from
patients with colorectal cancer. It is made up of 150×150 pixels color images
from eight classes, i.e., debris, mucosa, tumor, adipose, stroma, lympho,
complex, and empty. The original dataset consisted of 5,000 samples. We
randomly split the datasets into training, validation, and testing sets at a
ratio of 70%, 10%, and 20%, respectively. We reshaped the images into
128$\times$128 pixels, standardized, and normalized to a value range of 0–1.
Adience is composed of face images scraped from Flickr.com albums that were
labeled for age and gender. The benchmark uses eight classes for age groups
(0–2, 4–6, 8–13, 15–20, 25–32, 38–43, 48–53, 60+). A total of 38,740 images
were split into five groups of 4,484, 3,730, 3,894, 3,446, 3,816, and 19,370.
The Yale dataset contains 165 grayscale images of 15 individuals in GIF
format. There were 11 images per subject, 1 for each different facial
expression or configuration: center-light, w/glasses, happy, left-light,
without glasses, normal, right-light, sad, sleepy, surprised, and wink. We
preprocessed the Yale and Adience datasets by cropping and centering on the
faces and resized them into 150$\times$150 pixels and 160$\times$160 pixels
color images for the Adience and Yale datasets, respectively. In addition, we
conducted five cross-validations to assess the model performance.

TABLE II: Comparison of network architectures.

Model | Input mapping | Embeding net. | Merging layer | Tail layer | Loss function
---|---|---|---|---|---
CovNet v1 | IM | $F_{dataset}$ | Covariance | Dense($n_{class}$, | Categorical
Softmax) | crossentropy
CovNet v2 | IIM | $F_{dataset}$ | Covariance | Dense($n_{class}$, | Categorical
Softmax) | crossentropy
CovNet v3 | ISIM | $F_{dataset}$ | Covariance | BN + Dense(1, | Binary
Sigmoid) | crossentropy
Siamese | ISIM | $F_{dataset}$ | Euclidean dist. | BN + Dense(1, | Contrastive
Sigmoid) | loss
Triplet | TM* | $F_{dataset}$ | Triplet dist. | – | Triplet loss
N-pair | IM | $F_{dataset}$ | – | – | N-pair cosine
similarity

  * *

Triplet mapping.

We compared our experimental results with the state-of-art metric learning
models, such as Siamese, Triplet, and N-pair networks, to assess the
effectiveness of our proposed model. We employed the same embedding network
($F$) in all models for a fair comparison, as shown in Table I. Here,
$Conv(i,j,k)$ denotes a convolutional layer with $i$ number of filters, a
kernel size of $j\times j$, and $k$ activation functions. In addition, $MP(i)$
is a 2D max-pooling layer with a size of $i\times i$, $BN$ is a batch
normalization layer, and $DO(i)$ corresponds to the dropout layer with
probability $i$. In addition, $AP$ is a 2D average-pooling layer, $Dense(i,j)$
represents a dense layer with $i$ neurons and $j$ activation functions.
Moreover, $ZP(i)$ indicates zero-padding with size $i$. For the Adience and
Yale datasets, we employed a pre-trained network (FaceNet) [14], which is an
inception model trained on the MS-Celeb-1M dataset, as a backbone network.

The overall network architecture of all models is summarized in Table II. The
value of $F$ can vary based on the dataset, which is referred to in Table I.
Here, CovNet v1 and CovNet v2 have different properties in both the input
mapping and tail layer. In this case, CovNet v1 utilizes IM, whereas CovNet v2
employs IIM, which is notably more complex than IM. Because $n_{class}$ is
defined by all pair combinations of the class label in CovNet v2, the
$n_{class}$ of the CovNet v2 value is significantly higher than that of CovNet
v1. The $n_{class}$ of CovNet v1 is equal to the total number of class
categories. In addition, CovNet v3 has the same architecture as the Siamese
network except in the merging layer, which is represented as a binary
classification. Both triplet and N-pair networks are more straightforward than
the others because they do not have a tail layer. Nevertheless, their learning
mechanism can be trickier because triplet and N-pair losses are directly
optimized during the training phase. We employed batch sizes of 128, 128, 32,
64, and 16 for the CIFAR-10, Food-11, Colorectal, Adience, and Yale datasets,
respectively. Each network was trained for 200 epochs using Adam optimization.
An early stopping mechanism was executed when the best model was obtained
based on the validation sets.

$accuracy=\frac{\sum_{i=1}^{N}\sum_{j=1}^{k}c_{i}^{j}}{N\times k}$ (7)
$\begin{split}c_{i}^{j}\begin{cases}=1,&\mathcal{\hat{Y}}_{i}^{j}=\mathcal{Y}_{i}\\\
0,&\text{otherwise}\\\ \end{cases}\end{split}$ (8)

CovNet-v1CovNet-v2CovNet-v3SiameseTripletN-pairN-
pair(1000)$0.6$$0.62$$0.64$$0.66$$0.68$$0.7$$0.72$$0.74$$0.76$$0.78$$0.8$$0.82$AccuracyCIFAR-10$k$=10$k$=40$k$=70$k$=100

CovNet-v1CovNet-v2CovNet-v3SiameseTripletN-pairN-
pair(1000)$0.6$$0.62$$0.64$$0.66$$0.68$$0.7$$0.72$$0.74$$0.76$$0.78$$0.8$$0.82$CovNet-v1CovNet-v2CovNet-v3SiameseTripletN-
pairN-
pair(1000)$0.2$$0.3$$0.4$$0.5$$0.6$$0.7$AccuracyFood-11$k$=10$k$=40$k$=70$k$=100$k$=10$k$=40$k$=70$k$=100

CovNet-v1CovNet-v2CovNet-v3SiameseTripletN-pairN-
pair(1000)$0.6$$0.62$$0.64$$0.66$$0.68$$0.7$$0.72$$0.74$$0.76$$0.78$$0.8$$0.82$CovNet-v1CovNet-v2CovNet-v3SiameseTripletN-
pairN-
pair(1000)$0.2$$0.3$$0.4$$0.5$$0.6$$0.7$CovNet-v1CovNet-v2CovNet-v3SiameseTripletN-
pairN-
pair(1000)$0.55$$0.6$$0.65$$0.7$$0.75$$0.8$$0.85$$0.9$$0.95$$1$AccuracyColorectal$k$=10$k$=40$k$=70$k$=100$k$=10$k$=40$k$=70$k$=100$k$=10$k$=40$k$=70$k$=100

CovNet-v1CovNet-v2CovNet-v3SiameseTripletN-pairN-
pair(1000)$0.6$$0.62$$0.64$$0.66$$0.68$$0.7$$0.72$$0.74$$0.76$$0.78$$0.8$$0.82$CovNet-v1CovNet-v2CovNet-v3SiameseTripletN-
pairN-
pair(1000)$0.2$$0.3$$0.4$$0.5$$0.6$$0.7$CovNet-v1CovNet-v2CovNet-v3SiameseTripletN-
pairN-
pair(1000)$0.55$$0.6$$0.65$$0.7$$0.75$$0.8$$0.85$$0.9$$0.95$$1$CovNet-v1CovNet-v2CovNet-v3SiameseTripletN-
pairN-
pair(1000)$0.64$$0.66$$0.68$$0.7$$0.72$$0.74$$0.76$$0.78$$0.8$$0.82$$0.84$AccuracyAdience$k$=10$k$=40$k$=70$k$=100$k$=10$k$=40$k$=70$k$=100$k$=10$k$=40$k$=70$k$=100$k$=2$k$=6$k$=10$k$=14

CovNet-v1CovNet-v2CovNet-v3SiameseTripletN-pairN-
pair(1000)$0.6$$0.62$$0.64$$0.66$$0.68$$0.7$$0.72$$0.74$$0.76$$0.78$$0.8$$0.82$CovNet-v1CovNet-v2CovNet-v3SiameseTripletN-
pairN-
pair(1000)$0.2$$0.3$$0.4$$0.5$$0.6$$0.7$CovNet-v1CovNet-v2CovNet-v3SiameseTripletN-
pairN-
pair(1000)$0.55$$0.6$$0.65$$0.7$$0.75$$0.8$$0.85$$0.9$$0.95$$1$CovNet-v1CovNet-v2CovNet-v3SiameseTripletN-
pairN-
pair(1000)$0.64$$0.66$$0.68$$0.7$$0.72$$0.74$$0.76$$0.78$$0.8$$0.82$$0.84$CovNet-v1CovNet-v2CovNet-v3SiameseTripletN-
pairN-
pair(1000)$0.7$$0.75$$0.8$$0.85$$0.9$$0.95$$1$AccuracyYale$k$=10$k$=40$k$=70$k$=100$k$=10$k$=40$k$=70$k$=100$k$=10$k$=40$k$=70$k$=100$k$=2$k$=6$k$=10$k$=14$k$=1$k$=2$k$=3$k$=4

Figure 5: Effect of varying $k$ (nearest neighbor) on image search accuracyy.

###  4.2 Near neighbor analysis

The outcome of an embedding network is commonly represented as a vector having
a significantly lower dimension than the original input. The main objective of
the embedding method is to preserve the data separability in the embedding
space to be the same as possible as in the original space. Accordingly, any
machine learning task with complex and large dimensional datasets, such as
images, can be simply solved using the near neighbor algorithm. Equation 7
measures the accuracy of the near neighbors of each point in the embedding
space based on the corresponding label, where $k$ represents the number of
neighbors, and $N$ is the total number of data. The $c_{i}^{j}$ value is 1
when the data point and the predicted neighbor have the same label; otherwise,
it is zero, as shown in Eq. 8. Note that the denominator of Eq. 7 is $N\times
k$, which represents the total number of neighbors of each sample in the
dataset. In this study, we employed a simple $k$-nearest neighbor algorithm
with cosine similarity to measure the closeness among the points (represented
as vectors) in the embedding space.

The $k$ nearest neighbor performance of the embedding network for all models
is shown in Fig. 5. In general, the accuracy decreases slightly with
increasing values of $k$ (except for Adience, in which the accuracy is
significantly decreased). This indicates that the higher the value of $k$ is,
the more diverse the neighbors are. Nevertheless, in real applications, such
as search engines and recommendation systems, users are commonly interested
only in the top-10 neighbors. CovNet v2 outperforms other methods in all
datasets because it explicitly minimizes the inner-class separability and
maximizes the inter-class diversity using IIM. A Siamese network is the most
competitive existing method for the proposed models. However, it has a
slightly lower performance compared to the other approaches for the Adience
and Yale datasets, which used a pre-trained network. In comparison, both
Siamese and CovNet v3 utilized binary cross-entropy, which directly neglected
the inter-class relationship. However, CovNet v3 is generally better than a
Siamese network. Thus, we can reveal that employing covariance in the merging
layer is more robust than the Euclidean distance. The triplet network for both
Food-11 and Colorectal has worse performance than the other models. This is
because some samples in different classes can be nearly similar, as shown in
Fig. 4; for example, images labeled 3 and 6 in Food-11 are on the same plate,
and those labeled 6 and 7 in colorectal have nearly similar characteristics.
Accordingly, the distance between anchor-positive and anchor-negative can be
inconstant. In addition, the N-pair with a small batch achieves a low
performance. By contrast, if we increase the number of batches to 1000
(N-pair(1000)), its performance significantly increases. However, the
computational cost is considerably increased. Finally, we can observe that if
we utilize pre-trained networks (Adience and Yale) as a backbone, the
performances of all models are nearly competitive with each other.

As a service, the goal of computing technology is to perform business services
more efficiently and effectively. Nowadays, many real-world applications must
deal with large-scale datasets, such as image similarity searches and content-
based recommendation systems. Thus, metric learning can be efficient because
it can project high-dimensional images into a small vector. A vector
representation must also reflect the original data effectively. One of the
interesting applications of near-neighbor (or similarity) problems is image
search similarity. In this case, given a query image, the model retrieves $k$
similar images as the query image. The visualization of the image similarity
is illustrated in Figs. 6 and 7. The first column (yellow box) represents the
query image, whereas the remaining columns correspond to the top-10 similar
images. The red box represents irrelevant images based on the label of the
query image. In Fig. 6, the irrelevant images mainly occur on dogs-cats and
deer-horses, which share similar features. In Fig. 7, because the label (face
age) naturally has an ordinal relationship, the error commonly occurs in the
nearest age; for instance, 15–20 year-olds are often predicted as 25–32 year-
olds. In general, we can infer that CovNet v2 outperforms the other models.
Here, CovNet v1 and v3 slightly outperformed the Siamese network. Compared
with the triplet and N-pair networks, the Siamese network is more competitive
with our proposed models.

(a) CovNet v1

(b) CovNet v2

(c) CovNet v3

(d) Siamese

(e) Triplet

(f) N-Pair

Figure 6: Image search application on CIFAR-10: The first column (yellow box)
represents the query image, the remaining columns are the 10 nearest images of
the query image, and the irrelevant image result is marked by the red box.

(a) CovNet v1

(b) CovNet v2

(c) CovNet v3

(d) Siamese

(e) Triplet

(f) N-Pair

Figure 7: Image search application on Adience: The first column (yellow box)
represents the query image, the remaining columns are the 10 nearest images of
the query image, and the irrelevant image result is marked by the red box.

###  4.3 Semantic analysis

In this section, the semantic relationships between both inter-class
separability and inter-class relationships are observed. We present Fig. 8 as
a visualization of the embedding space to assess inter-class separability. The
output of the embedding network was a vector with 100 dimensions. Accordingly,
we utilized t-stochastic embedding (t-SNE) to project them into a two-
dimensional space. As illustrated in Fig. 8, our proposed methods have a
better embedding space separability compared with Siamese, Triplet, and N-pair
networks. Therefore, we can reveal that covariance embedding generates
separable outcomes more effectively than the Euclidean and cosine similarity
metrics. Note that, unlike the Euclidean and cosine similarity metrics,
covariance embedding is more expressive in capturing the relationship between
two inputs because it can reflect the positive, negative, or even neutral
correlation simultaneously. In addition, for the merging layer (see Table II),
CovNet obtains a vector that is directly assessed using categorical or binary
cross-entropy. By contrast, other models result in a scalar in the merging
layer. Thus, this mechanism can be more effective in preserving the data
separability in the embedding space because assessing a vector (multivariate)
as a feature is more effective for guaranteeing separability than representing
it as a scalar (univariate).

(a) CovNet v1

(b) CovNet v2

(c) CovNet v3

(d) Siamese

(e) Triplet

(f) N-pair

Figure 8: Visualization of latent separability on CIFAR-10 by projecting an
embedding space to a two-dimensional space using t-SNE.

(a) CovNet v1

(b) CovNet v2

(c) CovNet v3

(d) Siamese

(e) Triplet

(f) N-pair

Figure 9: Pearson correlation matrix to capture the semantic relationship
among categories on CIFAR-10.

In addition, we argue that a reliable embedding-space representation can be
assessed not only through inter-class separability but also by inter-class
relationships. For instance, semantically, we agree that dog-cat and truck-
automobile share similar shapes and properties. Therefore, a good embedding-
space representation must also preserve this special proximity property. This
property is useful for content-based recommendation systems and search engine
applications. For better services, in addition to recommending similar
products, it is better if the system also offers different products that have
similar characteristics as the user query. Consequently, we employed the
Pearson correlation method to assess the relationships among the inter-class
images. The results are summarized in a correlation matrix, as shown in Fig.
9. Compared with other models, our CovNet is more expressive in describing an
inter-class relationship. Unlike the existing models, the inter-class
correlation value is still high, whereas the value is still lower than the
inner-class correlation value. Inherently, a well-known similar object in
CIFAR-10 is in a different class, such as cat-dog, airplane-bird, deer-horse,
and automobile-truck (red box in Fig. 9). Nevertheless, unlike our CovNet
(particularly v1 and v2), they have low correlation values in Siamese,
Triplet, and N-pair networks. The correlation matrix developed using Siamese
is slightly similar to that of CovNet v3 because they have a similar
architecture. However, semantically, CovNet v3 has a slightly higher
correlation in terms of inter-class relationships than Siamese. In the triplet
network, the inter-class correlation value is too aggressive, and it has a
higher magnitude than its inner-class correlation. Moreover, the N-pair
network is too fierce in defining inter-class relationships, and thus some
inner-class correlation values also become low.

TABLE III: Inspection of embedding space separability using a standard
classifier (support vector machine).

Model | Accuracy % (mean $\pm$ standard deviation)
---|---
CIFAR-10 | Food-11 | Colorectal | Adience | Yale
CovNet v1 | $85.4\pm 0.48$ | $81.2\pm 0.72$ | $84.7\pm 0.53$ | $65.3\pm 0.49$ | $100.0\pm 0$
CovNet v2 | $86.1\pm 0.39$ | $82.3\pm 0.65$ | $86.5\pm 0.46$ | $66.0\pm 0.37$ | $100.0\pm 0$
CovNet v3 | $85.2\pm 0.41$ | $79.9\pm 0.70$ | $83.6\pm 0.67$ | $64.4\pm 0.43$ | $100.0\pm 0$
Siamese | $85.2\pm 0.42$ | $79.0\pm 0.69$ | $82.1\pm 0.50$ | $64.4\pm 0.47$ | $100.0\pm 0$
Triplet | $74.1\pm 0.61$ | $36.6\pm 0.94$ | $77.0\pm 0.73$ | $65.0\pm 0.68$ | $100.0\pm 0$
N-pair | $71.7\pm 0.58$ | $48.9\pm 0.74$ | $48.5\pm 0.65$ | $62.1\pm 0.65$ | $100.0\pm 0$
N-pair(1000) | $83.5\pm 0.35$ | $75.1\pm 0.67$ | $67.2\pm 0.59$ | $64.6\pm 0.37$ | $100.0\pm 0$

###  4.4 Classification and parameter

Because of their effectiveness, several metric learning models, such as
Siamese, Triplet, and N-pair networks, are commonly used for fine-tuning
classification. In this study, a metric learning model acts as a feature
extraction network (base model). Here, we train only a few layers on top of
it. Meanwhile, the weights of the pre-trained network were not updated during
the training. Commonly, we use a deep neural network (DNN) as the top model to
classify the data. However, because of its stochastic nature, and for the sake
of a fair comparison, we employ a standard classifier, i.e., an SVM with a
radial basis function (RBF) kernel, to classify the embedding outcome. The
classification performances are presented in Table III. In line with the
previous analysis results in the previous sections, CovNet v2 has the best
performance compared with other metric learning models. Siamese and CovNet v3
have a competitive performance because their architectures are quite similar.
The triplet network is slightly better than N-Pair, except in Food-11. In
addition, increasing the number of batches in the N-pair network can increase
the classification performance, except in the Yale dataset, which is a
relatively small dataset (165 instances). Thus, it was unaffected by the
number of batches. Moreover, all models become competitive with each other if
we utilize a pre-trained network, as shown in the Adience and Yale datasets.

TABLE IV: Comparison of the number of trained parameters in each model.

Model | Number of trainable parameters
---|---
CIFAR-10 | Food-11 | Colorectal | Adience* | Yale*
CovNet v1 | 1,830,294 | 2,486,523 | 2,486,523 | 119,180 | 119,887
CovNet v2 | 1,834,839 | 2,492,078 | 2,489,078 | 122,008 | 130,492
CovNet v3 | 1,829,585 | 2,485,713 | 2,485,713 | 118,673 | 118,673
Siamese | 1,829,585 | 2,485,713 | 2,485,713 | 118,673 | 118,673
Triplet | 1,829,284 | 2,485,412 | 2,485,412 | 970,340 | 128,472
N-pair | 1,829,284 | 2,485,412 | 2,485,412 | 970,340 | 128,472

  * *

Employing FaceNet [14] (pre-trained network) as a base model.

Finally, we present a comparison of the number of trainable parameters in each
model, as described in Table II. In general, the number of parameters in each
model is not significantly different. CovNet v2 has a higher trained parameter
because the tail layer has more neuron units compared with the others, as
described in Section 4.1. CovNet v3 and Siamese have the same numbers of
parameters because they have the same architecture except for the merging
layer. The merging layer in Table II has no parameters. The triplet and N-pair
networks have the least number of parameters because the triplet loss and
cosine similarity loss can be directly optimized such that a tail layer is not
required. Nevertheless, their performance is significantly lower than that of
the other models. The numbers of parameters in the Adience and Yale datasets
are less than the others because they employ a pre-trained network, FaceNet,
as a base model.

##  5 Conclusion

A desirable metric learning model should not only separate the dataset based
on its categories it should also maintain the inter-class semantic
relationship. We present CovNets, a novel metric learning method, as a service
approach with covariance embedding optimization. Through the development of
CovNet, we show two important properties: (1) A desirable metric learning
model should not only separate the dataset based on its categories but also
maintain the inter-class semantic relationship. (2) Covariance-enhanced
embedding can make complex machine learning tasks more expressive, more
explainable, and more efficient. For example, by producing covariance
embedding, face verification simply involves thresholding the distance between
two embeddings. Similarly, object recognition, recommender systems, and an
image search similarity when producing covariance embedding with CovNets can
be reduced to a $k$-NN classification problem. High-dimensional data
clustering with CovNets feature-engineering through covariance embedding can
be reduced to a linear-space problem solvable using simple bottom-up
techniques, including agglomerative clustering. Extensive empirical
experiments demonstrate the effectiveness of CovNets over representative
state-of-art metric learning models, including Siamese, Triplet, and N-pair
networks.

## Acknowledgment

This research was supported by the National Research Foundation of Korea(NRF)
grant funded by the Korea government(MSIT) (No. 2020R1A2C1102294).

## References

  * [1] P. Moutafis, M. Leng, and I. A. Kakadiaris, “An overview and empirical comparison of distance metric learning methods,” _IEEE Transactions on Cybernetics_ , vol. 47, no. 3, pp. 612–625, 2017.
  * [2] Q. Liu, W. Li, Z. Chen, and B. Hua, “Deep metric learning for image retrieval in smart city development,” _Sustainable Cities and Society_ , vol. 73, p. 103067, 2021.
  * [3] H. Wu, Q. Zhou, R. Nie, and J. Cao, “Effective metric learning with co-occurrence embedding for collaborative recommendations,” _Neural Networks_ , vol. 124, pp. 308–318, 2020.
  * [4] Y. Liu, D. Pi, and L. Cui, “Metric learning combining with boosting for user distance measure in multiple social networks,” _IEEE Access_ , vol. 5, pp. 19 342–19 351, 2017.
  * [5] J. Lu, G. Wang, W. Deng, and K. Jia, “Reconstruction-based metric learning for unconstrained face verification,” _IEEE Transactions on Information Forensics and Security_ , vol. 10, no. 1, pp. 79–89, 2015.
  * [6] P.-J. Last, H. A. Engelbrecht, and H. Kamper, “Unsupervised feature learning for speech using correspondence and siamese networks,” _IEEE Signal Processing Letters_ , vol. 27, pp. 421–425, 2020.
  * [7] M. Norouzi, D. J. Fleet, and R. R. Salakhutdinov, “Hamming distance metric learning,” in _Advances in Neural Information Processing Systems_ , F. Pereira, C. J. C. Burges, L. Bottou, and K. Q. Weinberger, Eds., vol. 25. Curran Associates, Inc., 2012\.
  * [8] D. Zhang, Y. Li, and Z. Zhang, “Deep metric learning with spherical embedding,” in _Advances in Neural Information Processing Systems_ , H. Larochelle, M. Ranzato, R. Hadsell, M. F. Balcan, and H. Lin, Eds., vol. 33. Curran Associates, Inc., 2020, pp. 18 772–18 783.
  * [9] I. M. Kamal and H. Bae, “Super-encoder with cooperative autoencoder networks,” _Pattern Recognition_ , vol. 126, p. 108562, 2022.
  * [10] E. Xing, M. Jordan, S. J. Russell, and A. Ng, “Distance metric learning with application to clustering with side-information,” in _Advances in Neural Information Processing Systems_ , S. Becker, S. Thrun, and K. Obermayer, Eds., vol. 15. MIT Press, 2003.
  * [11] K. Zhang, S. Qi, J. Cai, D. Zhao, T. Yu, Y. Yue, Y. Yao, and W. Qian, “Content-based image retrieval with a convolutional siamese neural network: Distinguishing lung cancer and tuberculosis in ct images,” _Computers in Biology and Medicine_ , vol. 140, p. 105096, 2022.
  * [12] Y. Qiao, Y. Wu, F. Duo, W. Lin, and J. Yang, “Siamese neural networks for user identity linkage through web browsing,” _IEEE Transactions on Neural Networks and Learning Systems_ , vol. 31, no. 8, pp. 2741–2751, 2020.
  * [13] X. Yang, H. Guo, N. Wang, B. Song, and X. Gao, “A novel symmetry driven siamese network for thz concealed object verification,” _IEEE Transactions on Image Processing_ , vol. 29, pp. 5447–5456, 2020.
  * [14] F. Schroff, D. Kalenichenko, and J. Philbin, “Facenet: A unified embedding for face recognition and clustering,” in _2015 IEEE Conference on Computer Vision and Pattern Recognition (CVPR)_ , 2015, pp. 815–823.
  * [15] M. Zhang, Q. Cheng, F. Luo, and L. Ye, “A triplet nonlocal neural network with dual-anchor triplet loss for high-resolution remote sensing image retrieval,” _IEEE Journal of Selected Topics in Applied Earth Observations and Remote Sensing_ , vol. 14, pp. 2711–2723, 2021.
  * [16] F. Boutros, N. Damer, F. Kirchbuchner, and A. Kuijper, “Self-restrained triplet loss for accurate masked face recognition,” _Pattern Recognition_ , vol. 124, p. 108473, 2022.
  * [17] H. Liu, Y. Chai, X. Tan, D. Li, and X. Zhou, “Strong but simple baseline with dual-granularity triplet loss for visible-thermal person re-identification,” _IEEE Signal Processing Letters_ , vol. 28, pp. 653–657, 2021.
  * [18] K. Sohn, “Improved deep metric learning with multi-class n-pair loss objective,” in _Advances in Neural Information Processing Systems_ , D. Lee, M. Sugiyama, U. Luxburg, I. Guyon, and R. Garnett, Eds., vol. 29. Curran Associates, Inc., 2016\.
  * [19] B. Chen and W. Deng, “Deep embedding learning with adaptive large margin n-pair loss for image retrieval and clustering,” _Pattern Recognition_ , vol. 93, pp. 353–364, 2019.
  * [20] A. Pal, Z. Xue, B. Befano, A. C. Rodriguez, L. R. Long, M. Schiffman, and S. Antani, “Deep metric learning for cervical image classification,” _IEEE Access_ , vol. 9, pp. 53 266–53 275, 2021.
  * [21] I. López-Espejo, Z.-H. Tan, and J. Jensen, “A novel loss function and training strategy for noise-robust keyword spotting,” _IEEE/ACM Transactions on Audio, Speech, and Language Processing_ , vol. 29, pp. 2254–2266, 2021.
  * [22] A. Bellet, A. Habrard, and M. Sebban, “A survey on metric learning for feature vectors and structured data,” _CoRR_ , vol. abs/1306.6709, 2013.
  * [23] J. Xu, L. Luo, C. Deng, and H. Huang, “Bilevel distance metric learning for robust image recognition,” in _Advances in Neural Information Processing Systems_ , S. Bengio, H. Wallach, H. Larochelle, K. Grauman, N. Cesa-Bianchi, and R. Garnett, Eds., vol. 31. Curran Associates, Inc., 2018.
  * [24] H.-J. Ye, D.-C. Zhan, X.-M. Si, Y. Jiang, and Z.-H. Zhou, “What makes objects similar: A unified multi-metric learning approach,” in _Advances in Neural Information Processing Systems_ , D. Lee, M. Sugiyama, U. Luxburg, I. Guyon, and R. Garnett, Eds., vol. 29. Curran Associates, Inc., 2016.
  * [25] F. Abramovich and M. Pensky, “Classification with many classes: Challenges and pluses,” _Journal of Multivariate Analysis_ , vol. 174, p. 104536, 2019.
  * [26] A. Krizhevsky, “Learning multiple layers of features from tiny images,” Tech. Rep., 2009.
  * [27] A. Singla, L. Yuan, and T. Ebrahimi, “Food/non-food image classification and food categorization using pre-trained googlenet model,” in _Proceedings of the 2nd International Workshop on Multimedia Assisted Dietary Management_ , ser. MADiMa ’16. New York, NY, USA: Association for Computing Machinery, 2016, p. 3–11.
  * [28] J. N. Kather, C.-A. Weis, F. Bianconi, S. M. Melchers, L. R. Schad, T. Gaiser, A. Marx, and F. G. Z”ollner, “Multi-class texture analysis in colorectal cancer histology,” _Scientific reports_ , vol. 6, p. 27988, 2016.
  * [29] E. Eidinger, R. Enbar, and T. Hassner, “Age and gender estimation of unfiltered faces,” _IEEE Transactions on Information Forensics and Security_ , vol. 9, no. 12, pp. 2170–2179, 2014.
  * [30] Y. Ling, X. Yin, and S. Bhandarkar, “Sirface vs. fisherface: recognition using class specific linear projection,” in _Proceedings 2003 International Conference on Image Processing (Cat. No.03CH37429)_ , vol. 3, 2003, pp. III–885.


interval, such that a cell $i$ is guaranteed to change its position. The total
set of speeds for $n_{cells}$ is $S:=\cup_{i=1}^{n_{cells}}S_{i}$. The average
speed $s_{av}$ is then defined as arithmetic average:

$\displaystyle s_{av}:=\frac{1}{\lvert S\rvert}\sum_{s\in S}s.$

The set of normalized velocities $V_{i}$ (or, alternatively, displacements) of
cell $i$ is given by:

$\displaystyle
V_{i}:=\left\\{\mathbf{v}\in\mathbb{R}^{2}:\mathbf{v}=\frac{\mathbf{x}^{i}(t_{(k+1)l_{i}})-\mathbf{x}^{i}(t_{kl_{i}})}{\lVert\mathbf{x}^{i}(t_{(k+1)l_{i}})-\mathbf{x}^{i}(t_{kl_{i}})\rVert},k\in\mathbb{N},(k+1)l_{i}\leq
n_{time}\right\\},$

and the total set of normalized velocities $V$ is given by
$V:=\cup_{i=1}^{n_{cells}}V_{i}$.

The directionality ratio $r_{i}(t_{j})$ of cell $i$ over a time interval of
length $t_{j}$ is given by:

$\displaystyle
r_{i}(t_{j}):=\frac{\sum_{k=1}^{j}\lVert\mathbf{x}^{i}(t_{k})-\mathbf{x}^{i}(t_{k-1})\rVert}{\lVert\mathbf{x}^{i}(t_{j})-\mathbf{x}^{i}(t_{0})\rVert}.$

The population and the time averages of the directionality ratio are given by,
respectively:

$\displaystyle r(t_{j})$
$\displaystyle=\frac{1}{n_{cell}}\sum_{i=1}^{n_{cell}}r_{i}(t_{j})$
$\displaystyle\bar{r}$
$\displaystyle=\frac{1}{n_{time}}\sum_{j=1}^{n_{time}}r(t_{j}).$

Velocity autocorrelation $v^{i}_{ac}(t_{j})$ of cell $i$ over a time interval
of length $t_{j}$ is given by:

$\displaystyle
v^{i}_{ac}(t_{j}):=\frac{1}{n_{time}-j}\sum_{k=1}^{n_{time}-j}\frac{\mathbf{v}^{i}(t_{j+k})\cdot\mathbf{v}^{i}(t_{k})}{\lVert\mathbf{v}^{i}(t_{j+k})\rVert\phantom{.}\lVert\mathbf{v}^{i}(t_{k})\rVert},$

where $j=1,\ldots,n_{time}-1$ and
$\mathbf{v}^{i}(t_{k})=\left(\mathbf{x}^{i}(t_{k})-\mathbf{x}^{i}(t_{k-1})\right)/\Delta
t$. Velocity autocorrelation of the population $v_{ac}(t)$ is simply the
arithmetic average of each cell’s velocity autocorrelation. To compute
$v_{ac}$ we used the time step of $12min$, whereas for all other quantities
involving time dependence (e.g. $msd$) we used the time step of $0.12min$.

We define the guidance parameter $G\in[0,1]$ similarly as in [67]:

$\displaystyle
G:=\frac{1}{\lvert\Theta\rvert}\sum_{\theta\in\Theta}g(\theta),$

where $\Theta:=\cup_{i=1}^{n_{cells}}\Theta_{i}$. The set $\Theta_{i}$ of
angles between a displacement vector of cell $i$ and the ECM fibers is defined
by

$\displaystyle\Theta_{i}:=\left\\{\theta\in[-\frac{\pi}{2},\frac{\pi}{2}]:\theta=arcsin\left(\frac{x^{i}_{2}(t_{(k+1)l_{i}})-x^{i}_{2}(t_{kl_{i}})}{\lVert\mathbf{x}^{i}(t_{(k+1)l_{i}})-\mathbf{x}^{i}(t_{kl_{i}})\rVert}\right),k\in\mathbb{N},(k+1)l_{i}\leq
n_{time}\right\\},$

where $x^{i}_{2}$ is the $y$-component of $\mathbf{x}^{i}$. The function
$g:[-\frac{\pi}{2},\frac{\pi}{2}]\rightarrow\left\\{0,1\right\\}$ is given by

$\displaystyle g(\theta)=\begin{cases*}1,\text{ if
}\lvert\theta\rvert\leq\pi/4\\\ 0,\text{ else}\end{cases*}.$

Thus, $G$ increases when displacements are aligned with the horizontal axis.

##  Appendix D Simulation of the PDMP

In order to compute the trajectories of the cell motility process, the
following algorithm is used:

  1. 1.

Set $(\nu_{0},\mathbf{X}_{0})\in A\times\Gamma$ and $t=T_{0}=0$

  2. 2.

For $k=0,1,\ldots$
Generate interarrival time $\Delta_{k}=T_{k+1}-T_{k}$, whose distribution
function is given by:

$\displaystyle\mathbb{P}\left(\Delta_{k}\leq\tau\right)=1-\exp\left(-\int_{t}^{t+\tau}a_{0}(\nu_{t},\phi_{\nu}(s,\mathbf{X}_{t}))ds\right)$
(D.1)

Compute $\mathbf{X}_{t+s^{-}}:=\phi_{\nu}(s,\mathbf{X}_{t})$
Set $t=T_{k+1}=T_{k}+\Delta_{k}$
Generate $(\nu_{t},\mathbf{X}_{t})\sim
Q(\cdot;(\nu_{t^{-}},\mathbf{X}_{t^{-}}))$

Algorithm 1 Simulation of the PDMP

To generate the interarrival time $\Delta_{k}$, we need to solve for $\tau$ in
the following equation:

$\displaystyle
f(\tau):=\int_{t}^{t+\tau}a_{0}(\nu_{t},\phi_{\nu_{t}}(s,\mathbf{X}_{t}))ds+\ln(1-u)=0,$
(D.2)

where $u$ is uniformly distributed on the interval $(0,1)$. Notice that the
evaluation of the integral by a quadrature rule requires computing the
solution $\mathbf{X}_{t+s}=\phi_{\nu_{t}}(s,\mathbf{X}_{t})$ up to time $s$,
where $s$ is a quadrature point. Moreover, using an iterative method to solve
(D.2) requires computing the integral at each iteration. Therefore, it is
imperative to devise an efficient method to sample from distribution (D.1). In
the following, we propose a general method to generate the next event time.

###  D.1 Simulation of the next event time

Let $T_{k}\leq t<T_{k+1}$ and let
$\mathbf{G}(\cdot;h):\Gamma\rightarrow\Gamma$ represent a numerical method to
solve the ODE system

$\displaystyle\frac{d}{dt}\mathbf{X}_{t}=\mathbf{H}_{\nu_{t}}(\mathbf{X}_{t})$

for one time step $h$. That is,
$\mathbf{X}_{t+h}=\mathbf{G}(\mathbf{X}_{t};h)$ is the numerical solution of
the above ODE system at time $t+h$.

Let $[T_{k},T_{k+1})^{2}\ni(s^{\prime},t^{\prime})\mapsto
A_{0}(t^{\prime},s^{\prime}):=\int_{s^{\prime}}^{t^{\prime}}a_{0}(\nu_{t},\phi_{\nu_{t}}(u,\mathbf{X}_{t}))du$
denote the integrated rate function.

The method to find the root $\tau$ of equation (D.2) is given in Algorithm 2.
First, in steps (1-22), we find the upper bound $\tau_{max}$ by solving the
ODE system for $n$ steps with step size $h$ and store the solution, the
computed rate $a_{0}$, and the integrated rate $A_{0}$ at these time steps.
Then, for any $\tau\leq\tau_{max}$ we can compute $A_{0}(t+\tau,t)$ upon using
the stored value of $A_{0}$ at time $t+\tau_{i}$, where
$i=\left\lfloor{\frac{\tau}{h}}\right\rfloor$ and $\tau_{i}=ih$:

$\displaystyle A_{0}(t+\tau,t)$
$\displaystyle=\int_{t}^{t+\tau}a_{0}(\nu_{t},\phi_{\nu_{t}}(s,\mathbf{X}_{t}))ds$
$\displaystyle=\int_{t}^{t+\tau_{i}}a_{0}(\nu_{t},\mathbf{X}_{t+s})ds+\int_{t+\tau_{i}}^{t+\tau_{i}+\tau-\tau_{i}}a_{0}(\nu_{t},\phi_{\nu_{t}}(s,\mathbf{X}_{t}))ds$
$\displaystyle=A_{0}(t+\tau_{i},t)+\int_{t+\tau_{i}}^{t+\tau_{i}+h_{i}}a_{0}(\nu_{t},\phi_{\nu_{t}}(s,\mathbf{X}_{t}))ds.$

Algorithm 2 Event time computation

1:procedure Initialization

2: Input: Time $t=T_{k}$; $(\nu_{t},\mathbf{X}_{t})$; time step $h$ and ODE
method $\mathbf{G}$; $n,m\in\mathbb{N}$.

3: Set $s_{0}:=0,n_{0}:=0$, create $List_{a_{0}}$, $List_{A_{0}}$, and
$List_{\mathbf{X}}$.

4: Append $List_{a_{0}}\leftarrow a_{0}(\nu_{t},\mathbf{X}_{t})$,
$List_{A_{0}}\leftarrow 0$, $List_{\mathbf{X}}\leftarrow\mathbf{X}_{t}$.

5: Set $\tau_{max}:=nh$, $s_{0}:=n_{0}h$.

6: Set initial condition $\mathbf{X}_{s}\leftarrow List_{\mathbf{X}}[last]$.

7: Set $A_{0}^{0}:=List_{A_{0}}[last]$.

8: for $i=1$ to $n$ do

9: $s_{i}:=t+s_{0}+ih$.

10: Compute $\mathbf{X}_{s_{i}}:=\mathbf{G}(\mathbf{X}_{s};h)$ and
$a_{0}(\nu_{t},\mathbf{X}_{s_{i}})$.

11: Compute $A_{0}(s_{i};t+s_{0})$ with quadrature points $s_{j}$,
$j=0,\ldots,i$.

12: $A_{0}(s_{i};t):=A_{0}^{0}+A_{0}(s_{i};t+s_{0})$.

13: Append $List_{a_{0}}\leftarrow a_{0}(\nu_{t},\mathbf{X}_{s_{i}})$,
$List_{A_{0}}\leftarrow A_{0}(s_{i};t)$,
$List_{\mathbf{X}}\leftarrow\mathbf{X}_{s_{i}}$.

14: end for

15: Generate $u\sim U(0,1)$.

16: if $List_{A_{0}}[last]<-ln(1-u)$ then

17: $n:=n+m$.

18: $n_{0}=n_{0}+m$.

19: go to 5.

20: end if

21: Output: $\tau_{max}$, $List_{a_{0}}$, $List_{A_{0}}$, and
$List_{\mathbf{X}}$.

22:end procedure

23:procedure Evalution of $f$

24: Input: $\tau$; time step $h$; $List_{a_{0}}$, $List_{A_{0}}$; Integrated
interpolation method $I$.

25: Set $i:=\left\lfloor{\frac{\tau}{h}}\right\rfloor$ and $h_{i}=\tau-ih$.

26: Output:
$f(\tau)=List_{A_{0}}[i]+I(List_{a_{0}}[i],List_{a_{0}}[i+1];h_{i},h)+\ln(1-u)$

27:end procedure

28:procedure Event time

29: Find the root $\tau$ of $f(\tau)=0$ using 23 and a root finding method.

30: $i:=\left\lfloor{\frac{\tau}{h}}\right\rfloor$.

31: Compute $\mathbf{X}_{t+\tau}:=\mathbf{G}(List_{\mathbf{X}}[i];\tau-ih)$.

32: Output: $\tau$, $\mathbf{X}_{t+\tau}$

33:end procedure

To compute the last integral in the expression above, we interpolate the
integrand using the stored values $a_{0}(\nu_{t},\mathbf{X}_{t+\tau_{i}})$ and
$a_{0}(\nu_{t},\mathbf{X}_{t+\tau_{i+1}})$. Let

$\displaystyle
I\left(a_{0}(\nu_{t},\mathbf{X}_{t+\tau_{i}}),a_{0}(\nu_{t},\mathbf{X}_{t+\tau_{i+1}});h_{i},h\right):=\int_{t+\tau_{i}}^{t+\tau_{i}+h_{i}}a_{0}(\nu_{t},\phi_{\nu_{t}}(s,\mathbf{X}_{t}))$

denote the approximation of the integral using interpolated integrand. We can
use the following interpolation methods for $t+\tau_{i}<s<t+\tau_{i+1}$:

  1. 1.

Piecewise constant
Forward:
$a_{0}(\nu_{t},\phi_{\nu_{t}}(s,\mathbf{X}_{t}))=a_{0}(\nu_{t},\phi_{\nu_{t}}(\tau_{i},\mathbf{X}_{t}))=a_{0}(\nu_{t},\mathbf{X}_{t+\tau_{i}})$.

$\displaystyle
I\left(a_{0}(\nu_{t},\mathbf{X}_{t+\tau_{i}}),a_{0}(\nu_{t},\mathbf{X}_{t+\tau_{i+1}});h_{i},h\right)=h_{i}a_{0}(\nu_{t},\mathbf{X}_{t+\tau_{i}}).$
(D.3)

Backward:
$a_{0}(\nu_{t},\phi_{\nu_{t}}(s,\mathbf{X}_{t}))=a_{0}(\nu_{t},\phi_{\nu_{t}}(\tau_{i+1},\mathbf{X}_{t}))=a_{0}(\nu_{t},\mathbf{X}_{t+\tau_{i+1}})$

$\displaystyle
I\left(a_{0}(\nu_{t},\mathbf{X}_{t+\tau_{i}}),a_{0}(\nu_{t},\mathbf{X}_{t+\tau_{i+1}});h_{i},h\right)=h_{i}a_{0}(\nu_{t},\mathbf{X}_{t+\tau_{i+1}}).$
(D.4)

  2. 2.

Average:
$a_{0}(\nu_{t},\phi_{\nu_{t}}(s,\mathbf{X}_{t}))=\frac{1}{2}\left(a_{0}(\nu_{t},\mathbf{X}_{t+\tau_{i}})+a_{0}(\nu_{t},\mathbf{X}_{t+\tau_{i+1}})\right)$.

$\displaystyle
I\left(a_{0}(\nu_{t},\mathbf{X}_{t+\tau_{i}}),a_{0}(\nu_{t},\mathbf{X}_{t+\tau_{i+1}});h_{i},h\right)=\frac{1}{2}h_{i}\left(a_{0}(\nu_{t},\mathbf{X}_{t+\tau_{i}})+a_{0}(\nu_{t},\mathbf{X}_{t+\tau_{i+1}})\right).$
(D.5)

  3. 3.

Piecewise linear:

$\displaystyle
a_{0}(\nu_{t},\phi_{\nu_{t}}(s,\mathbf{X}_{t}))=(s-t-\tau_{i})\frac{a_{0}(\nu_{t},\mathbf{X}_{t+\tau_{i+1}})-a_{0}(\nu_{t},\mathbf{X}_{t+\tau_{i}})}{h}+a_{0}(\nu_{t},\mathbf{X}_{t+\tau_{i}}).$
$\displaystyle I\big{(}a_{0}(\nu_{t},\mathbf{X}_{t+\tau_{i}}),a_{0}(\nu_{t},$
$\displaystyle\mathbf{X}_{t+\tau_{i+1}});h_{i},h\big{)}.$
$\displaystyle=h_{i}a_{0}(\nu_{t},\mathbf{X}_{t+\tau_{i}})\left(1-\frac{h_{i}}{h}\right)+\frac{h^{2}_{i}}{2h}a_{0}(\nu_{t},\mathbf{X}_{t+\tau_{i+1}}).$
(D.6)

Thus, $f(\tau)$ can be evaluated using equations (D.3)-(3):

$\displaystyle
f(\tau)=A_{0}(t+\tau_{i},t)+I\left(a_{0}(\nu_{t},\mathbf{X}_{t+\tau_{i}}),a_{0}(\nu_{t},\mathbf{X}_{t+\tau_{i+1}});h_{i},h\right)+\ln(1-u).$

Using the interpolations above, we can now employ, for example, Newton’s
method to find the root of equation (D.2):

$\displaystyle\tau_{l+1}=\tau_{l}-\frac{f(\tau_{l})}{a_{0}(\nu_{t},\phi_{\nu_{t}}(\tau_{l},\mathbf{X}_{t}))},$

or any other root-finding method (note that $a_{0}>0$).

Once the root is found, we simply advance the ODE system for one time step as
described in Steps (28-33) of the Algorithm 2.

Note that we solve the ODE system for $n+1=\tau_{max}/h+1$ steps and we also
solve for the interarrival time $\tau$ primarily by using a look-up table.
Moreover, obtaining a relatively sharp upper bound $\tau_{max}$ does not yield
a large computational overhead, since one simply can start the Algorithm 2
with a small $n,m$. Consequently, choosing an initial guess close to the sharp
upper bound for Newton’s method results in a faster convergence. In case of
thinning methods (see [45] or [8] for adaptive method), after each rejection
one needs to solve the ODE system for time period that is, on average,
approximately the same as $\tau_{max}$ (in the best case scenarios for both
methods, i.e. when the bound $\tau_{max}$ and the bound for the rate function
in thinning methods are sharp). Of course, these arguments hold when the
computational cost of solving the ODE system is relatively large.

###  D.2 Sampling from the transition measure

Given the time $t$ of the next event and $\mathbf{X}_{t}$ we need to sample
from the transition distribution $Q(\cdot,(\nu_{t^{-}},\mathbf{X}_{t^{-}}))$.
Recalling Section 3.2 and the proof of Proposition 4.3, in order to sample
from the transition measure it is sufficient to simulate the index
$j\in\left\\{1,\ldots,2M\right\\}$ of the next reaction, since the continuous
component of the process does not jump. The discrete distribution of the next
reaction index is given by equation (4.7):

$\displaystyle\mathbb{P}(j|\bm{\alpha}_{2}(\nu_{t^{-}}),\mathbf{X}_{t^{-}})=\frac{a_{j}(\bm{\alpha}_{2}(\nu_{t^{-}}),\mathbf{X}_{t^{-}})}{a_{0}(\bm{\alpha}_{2}(\nu_{t^{-}}),\mathbf{X}_{t^{-}})}.$

To simulate from the discrete distribution one can use the fairly efficient
Vose Alias Method [86].

## References

  * [1] M. Abercrombie.  The croonian lecture, 1978 - the crawling movement of metazoan cells.  Proceedings of the Royal Society of London B: Biological Sciences, 207(1167):129–147, 1980.
  * [2] A. Y. Alexandrova, K. Arnold, S. Schaub, J. M. Vasiliev, J.-J. Meister, A. D. Bershadsky, and A. B. Verkhovsky.  Comparative dynamics of retrograde actin flow and focal adhesions: formation of nascent adhesions triggers transition from fast to slow flow.  PLoS ONE, 3(9):e3234, 2008.
  * [3] Y. Aratyn-Schaus, P. W. Oakes, and M. L. Gardel.  Dynamic and structural signatures of lamellar actomyosin force generation.  Molecular Biology of the Cell, 22(8):1330–1339, 2011.
  * [4] D. Arcizet, B. Meier, E. Sackmann, J. O. Rädler, and D. Heinrich.  Temporal analysis of active and passive transport in living cells.  Physical Review Letters, 101:248103, Dec 2008.
  * [5] C. Arrieumerlou and T. Meyer.  A local coupling model and compass parameter for eukaryotic chemotaxis.  Developmental Cell, 8(2):215–227, 2005.
  * [6] N. Q. Balaban, U. S. Schwarz, D. Riveline, P. Goichberg, G. Tzur, I. Sabanay, D. Mahalu, S. Safran, A. Bershadsky, L. Addadi, et al.  Force and focal adhesion assembly: a close relationship studied using elastic micropatterned substrates.  Nature Cell Biology, 3(5):466–472, 2001.
  * [7] G. I. Bell.  Models for the specific adhesion of cells to cells.  Science, 200(4342):618–627, 1978.
  * [8] A. Bouchard-Côté, S. J. Vollmer, and A. Doucet.  The bouncy particle sampler: A non-reversible rejection-free Markov chain Monte Carlo method.  ArXiv e-prints, Oct. 2015.
  * [9] K. Burridge and M. Chrzanowska-Wodnicka.  Focal adhesions, contractility, and signaling.  Annual Review of Cell and Developmental Biology, 12(1):463–519, 1996.
  * [10] B. A. Camley, Y. Zhao, B. Li, H. Levine, and W.-J. Rappel.  Crawling and turning in a minimal reaction-diffusion cell motility model: coupling cell shape and biochemistry.  Physical Review E, 95:012401, Jan 2017.
  * [11] C. E. Chan and D. J. Odde.  Traction dynamics of filopodia on compliant substrates.  Science, 322(5908):1687–1691, 2008.
  * [12] W. T. Chen.  Mechanism of retraction of the trailing edge during fibroblast movement.  Journal of Cell Biology, 90(1):187–200, 1981.
  * [13] C. K. Choi, M. Vicente-Manzanares, J. Zareno, L. A. Whitmore, A. Mogilner, and A. R. Horwitz.  Actin and $\alpha$-actinin orchestrate the assembly and maturation of nascent adhesions in a myosin II motor-independent manner.  Nature Cell Biology, 10(9):1039–1050, 2008.
  * [14] C. A. Copos, S. Walcott, J. C. del Álamo, E. Bastounis, A. Mogilner, and R. D. Guy.  Mechanosensitive adhesion explains stepping motility in amoeboid cells.  Biophysical Journal, 112(12):2672 – 2682, 2017.
  * [15] L. P. Cramer.  Forming the cell rear first: breaking cell symmetry to trigger directed cell migration.  Nature Cell Biology, 12:628 – 632, 2010.
  * [16] M. H. Davis.  Piecewise-deterministic Markov processes: A general class of non-diffusion stochastic models.  Journal of the Royal Statistical Society. Series B (Methodological), 46(3):353–388, 1984.
  * [17] M. H. A. Davis.  Markov models and optimization.  Chapman and Hall, 1993.
  * [18] S. Deguchi, T. Ohashi, and M. Sato.  Tensile properties of single stress fiber isolated from cultured vascular smooth muscle cells.  Journal of Biomechanics, 39(14):2603–2610, 2006.
  * [19] R. B. Dickinson and R. T. Tranquillo.  A stochastic model for adhesion-mediated cell random motility and haptotaxis.  Journal of Mathematical Biology, 31(6):563–600, 1993.
  * [20] P. Dieterich, R. Klages, R. Preuss, and A. Schwab.  Anomalous dynamics of cell migration.  Proceedings of the National Academy of Sciences, 105(2):459–463, 2008.
  * [21] C. Engwer, T. Hillen, M. Knappitsch, and C. Surulescu.  Glioma follow white matter tracts: a multiscale DTI-based model.  Journal of Mathematical Biology, 71(3):551–582, 2015.
  * [22] T. Erdmann and U. S. Schwarz.  Stability of adhesion clusters under constant force.  Physical Review Letters, 92(4):108102, 2004.
  * [23] J. T. Finer, R. M. Simmons, and J. A. Spudich.  Single myosin molecule mechanics: piconewton forces and nanometre steps.  Nature, 368(6467):113–119, 1994.
  * [24] A. Fletcher, M. Osterfield, R. Baker, and S. Shvartsman.  Vertex models of epithelial morphogenesis.  Biophysical Journal, 106(11):2291 – 2304, 2014.
  * [25] P. Friedl, P. Noble, E. Shields, and K. Zänker.  Locomotor phenotypes of unstimulated CD45RAhigh and CD45ROhigh CD4+ and CD8+ lymphocytes in three-dimensional collagen lattices.  Immunology, 82(4):617, 1994.
  * [26] M. L. Gardel, I. C. Schneider, Y. Aratyn-Schaus, and C. M. Waterman.  Mechanical integration of actin and adhesion dynamics in cell migration.  Annual review of cell and developmental biology, 26:315–333, 2010\.
  * [27] M. L. Gardel, M. T. Valentine, J. C. Crocker, A. R. Bausch, and D. A. Weitz.  Microrheology of entangled F-actin solutions.  Physical Review Letters, 91:158302, Oct 2003.
  * [28] B. Geiger, J. P. Spatz, and A. D. Bershadsky.  Environmental sensing through focal adhesions.  Nature Reviews. Molecular Cell Biology, 10(1):21, 2009.
  * [29] D. T. Gillespie.  A rigorous derivation of the chemical master equation.  Physica A: Statistical Mechanics and its Applications, 188(1-3):404–425, 1992.
  * [30] C. Giverso and L. Preziosi.  Mechanical perspective on chemotaxis.  Physical Review E, 98(6):062402, 2018.
  * [31] R. Gorelik and A. Gautreau.  Quantitative and unbiased analysis of directional persistence in cell migration.  Nature protocols, 9(8):1931, 2014.
  * [32] A. Groh and A. K. Louis.  Stochastic modelling of biased cell migration and collagen matrix modification.  Journal of Mathematical Biology, 61(5):617–647, 2010.
  * [33] P. Guthardt Torres, I. B. Bischofs, and U. S. Schwarz.  Contractile network models for adherent cells.  Physical Review Letters E, 85:011913, Jan 2012.
  * [34] W. R. Holmes and L. Edelstein-Keshet.  A comparison of computational models for eukaryotic cell shape and motility.  PLoS Computational Biology, 8(12):e1002793, 2012.
  * [35] W. R. Holmes, B. Lin, A. Levchenko, and L. Edelstein-Keshet.  Modelling cell polarization driven by synthetic spatially graded Rac activation.  PLoS Computational Biology, 8(6):e1002366, 2012.
  * [36] P. Hotulainen and P. Lappalainen.  Stress fibers are generated by two distinct actin assembly mechanisms in motile cells.  The Journal of Cell Biology, 173(3):383–394, 2006.
  * [37] T. Kalwarczyk, N. Ziȩbacz, A. Bielejewska, E. Zaboklicka, K. Koynov, J. Szymański, A. Wilk, A. Patkowski, J. Gapiński, H.-J. Butt, and R. Hołyst.  Comparative analysis of viscosity of complex liquids and cytoplasm of mammalian cells at the nanoscale.  Nano Letters, 11(5):2157–2163, 2011.
  * [38] E. Kassianidou, A. C. Brand, U. S. Schwarz, and S. Kumar.  Geometry and network connectivity govern the mechanics of stress fibers.  Proceedings of the National Academy of Sciences, 114(10):2622–2627, 2017.
  * [39] E. Kassianidou and S. Kumar.  A biomechanical perspective on stress fiber structure and function.  Biochimica et Biophysica Acta (BBA)-Molecular Cell Research, 1853(11):3065–3074, 2015.
  * [40] J. Kelkel and C. Surulescu.  A multiscale approach to cell migration in tissue networks.  Mathematical Models and Methods in Applied Sciences, 22(03):1150017, 2012.
  * [41] K. E. Kubow, V. D. Shuklis, D. J. Sales, and A. R. Horwitz.  Contact guidance persists under myosin inhibition due to the local alignment of adhesions and individual protrusions.  Scientific Reports, 7(1):14380, 2017.
  * [42] S. Kumar, I. Z. Maxwell, A. Heisterkamp, T. R. Polte, T. P. Lele, M. Salanga, E. Mazur, and D. E. Ingber.  Viscoelastic retraction of single living stress fibers and its impact on cell shape, cytoskeletal organization, and extracellular matrix mechanics.  Biophysics Journal, 90(10):3762–3773, 2006.
  * [43] M. Kuragano, T. Q. Uyeda, K. Kamijo, Y. Murakami, and M. Takahashi.  Different contributions of nonmuscle myosin iia and iib to the organization of stress fiber subtypes in fibroblasts.  Molecular biology of the cell, 29(8):911–922, 2018.
  * [44] D. A. Lauffenburger and A. F. Horwitz.  Cell migration: a physically integrated molecular process.  Cell, 84(3):359–369, 1996.
  * [45] P. A. Lewis and G. S. Shedler.  Simulation of nonhomogeneous Poisson processes by thinning.  Naval Research Logistics (NRL), 26(3):403–413, 1979.
  * [46] F. Li, S. D. Redick, H. P. Erickson, and V. T. Moy.  Force measurements of the $\alpha$5$\beta$1 integrin–fibronectin interaction.  Biophysical Journal, 84(2):1252–1262, 2003.
  * [47] L. Li, S. F. Nørrelykke, and E. C. Cox.  Persistent cell motion in the absence of external signals: A search strategy for eukaryotic cells.  PLoS ONE, 3(5):1–11, 05 2008.
  * [48] Y.-J. Liu, M. Le Berre, F. Lautenschlaeger, P. Maiuri, A. Callan-Jones, M. Heuzé, T. Takaki, R. Voituriez, and M. Piel.  Confinement and low adhesion induce fast amoeboid migration of slow mesenchymal cells.  Cell, 160(4):659–672, 2015.
  * [49] A. Manhart, D. Oelz, C. Schmeiser, and N. Sfakianakis.  An extended filament based lamellipodium model produces various moving cell shapes in the presence of chemotactic signals.  Journal of Theoretical Biology, 382:244–258, 2015.
  * [50] A. F. Marée, A. Jilkine, A. Dawes, V. A. Grieneisen, and L. Edelstein-Keshet.  Polarization and movement of keratocytes: a multiscale modelling approach.  Bulletin of Mathematical Biology, 68(5):1169–1211, 2006.
  * [51] A. F. M. Marée, V. A. Grieneisen, and L. Edelstein-Keshet.  How cells integrate complex stimuli: The effect of feedback from phosphoinositides and cell shape on cell polarization and motility.  PLoS Computational Biology, 8(3):1–20, 03 2012.
  * [52] A. M. Mastro, M. A. Babich, W. D. Taylor, and A. D. Keith.  Diffusion of a small molecule in the cytoplasm of mammalian cells.  Proceedings of the National Academy of Sciences, 81(11):3414–3418, 1984.
  * [53] A. Mogilner.  Mathematics of cell motility: have we got its number?  Journal of Mathematical Biology, 58(1-2):105, 2009.
  * [54] C. Möhl, N. Kirchgessner, C. Schäfer, B. Hoffmann, and R. Merkel.  Quantitative mapping of averaged focal adhesion dynamics in migrating cells by shape normalization.  Journal of Cell Science, 125(1):155–165, 2012.
  * [55] J. E. Molloy, J. E. Burns, J. Kendrick-Jones, R. T. Tregear, and D. C. White.  Movement and force produced by a single myosin head.  Nature, 378(6553):209–212, 1995.
  * [56] Y. Mori, A. Jilkine, and L. Edelstein-Keshet.  Wave-pinning and cell polarity from a bistable reaction-diffusion system.  Biophysical Journal, 94(9):3684–3697, 2008.
  * [57] M. Murrell, P. W. Oakes, M. Lenz, and M. L. Gardel.  Forcing cells into shape: the mechanics of actomyosin contractility.  Nature Reviews Molecular Cell Biology, 16(8):486–498, 2015.
  * [58] A. Narang.  Spontaneous polarization in eukaryotic gradient sensing: a mathematical model based on mutual inhibition of frontness and backness pathways.  Journal of Theoretical Biology, 240(4):538–553, 2006.
  * [59] P. W. Oakes, E. Wagner, C. A. Brand, D. Probst, M. Linke, U. S. Schwarz, M. Glotzer, and M. L. Gardel.  Optogenetic control of RhoA reveals zyxin-mediated elasticity of stress fibers.  Nature Communications, 8:15817, 2017.
  * [60] D. Oelz and C. Schmeiser.  How do cells move? Mathematical modeling of cytoskeleton dynamics and cell migration. Cell mechanics: from single scale-based models to multiscale modeling. Chapman and Hall, 2010.
  * [61] H. G. Othmer, S. R. Dunbar, and W. Alt.  Models of dispersal in biological systems.  Journal of Mathematical Biology, 26(3):263–298, 1988.
  * [62] K. Paňková, D. Rösel, M. Novotnỳ, and J. Brábek.  The molecular mechanisms of transition between mesenchymal and amoeboid invasiveness in tumor cells.  Cellular and molecular life sciences, 67(1):63–71, 2010.
  * [63] R. Paul, P. Heil, J. P. Spatz, and U. S. Schwarz.  Propagation of mechanical stress through the actin cytoskeleton toward focal adhesions: model and experiment.  Biophysical Journal, 94(4):1470–1482, 2008.
  * [64] R. J. Petrie, A. D. Doyle, and K. M. Yamada.  Random versus directionally persistent cell migration.  Nature Reviews. Molecular Cell Biology, 10(8):538–549, 2009.
  * [65] L. Preziosi and G. Vitale.  A multiphase model of tumor and tissue growth including cell adhesion and plastic reorganization.  Mathematical Models and Methods in Applied Sciences, 21(09):1901–1932, 2011.
  * [66] M. Raftopoulou and A. Hall.  Cell migration: Rho GTPases lead the way.  Developmental Biology, 265(1):23–32, 2004.
  * [67] G. R. Ramirez-San Juan, P. W. Oakes, and M. L. Gardel.  Contact guidance requires spatial control of leading-edge protrusion.  Molecular Biology of the Cell, 28(8):1043–1053, 2017.
  * [68] A. J. Ridley.  Rho GTPases and cell migration.  Journal of Cell Science, 114(15):2713–2722, 2001.
  * [69] A. J. Ridley, M. A. Schwartz, K. Burridge, R. A. Firtel, M. H. Ginsberg, G. Borisy, J. T. Parsons, and A. R. Horwitz.  Cell migration: integrating signals from front to back.  Science, 302(5651):1704–1709, 2003.
  * [70] A. Roycroft and R. Mayor.  Molecular basis of contact inhibition of locomotion.  Cellular and Molecular Life Sciences, 73(6):1119–1130, 2016\.
  * [71] B. Rubinstein, K. Jacobson, and A. Mogilner.  Multiscale two-dimensional modeling of a motile simple-shaped cell.  Multiscale Modeling & Simulation, 3(2):413–439, 2005.
  * [72] R. J. Russell, S.-L. Xia, R. B. Dickinson, and T. P. Lele.  Sarcomere mechanics in capillary endothelial cells.  Biophysical Journal, 97(6):1578–1585, 2009.
  * [73] J. Satulovsky, R. Lui, and Y.-l. Wang.  Exploring the control circuit of cell migration by mathematical modeling.  Biophysical Journal, 94(9):3671–3683, 2008.
  * [74] D. Selmeczi, S. Mosler, P. H. Hagedorn, N. B. Larsen, and H. Flyvbjerg.  Cell motility as persistent random motion: Theories from experiments.  Biophysical Journal, 89(2):912 – 931, 2005.
  * [75] D. Shao, W.-J. Rappel, and H. Levine.  Computational model for cell morphodynamics.  Physical Review Letters, 105(10):108104, 2010.
  * [76] A. Shellard, A. Szabó, X. Trepat, and R. Mayor.  Supracellular contraction at the rear of neural crest cell groups drives collective chemotaxis.  Science, 362(6412):339–343, 2018.
  * [77] J. R. D. Soiné, C. A. Brand, J. Stricker, P. W. Oakes, M. L. Gardel, and U. S. Schwarz.  Model-based traction force microscopy reveals differential tension in cellular actin bundles.  PLoS Computational Biology, 11(3):e1004076, 2015.
  * [78] J. Stricker, Y. Aratyn-Schaus, P. W. Oakes, and M. L. Gardel.  Spatiotemporal constraints on the force-dependent growth of focal adhesions.  Biophysical Journal, 100(12):2883–2893, 2011.
  * [79] T. M. Svitkina, I. G. Surguchova, A. B. Verkhovsky, V. I. Gelfand, M. Moeremans, and J. De Mey.  Direct visualization of bipolar myosin filaments in stress fibers of cultured fibroblasts.  Cell Motility and the Cytoskeleton, 12(3):150–156, 1989.
  * [80] R. T. Tranquillo, D. A. Lauffenburger, and S. Zigmond.  A stochastic model for leukocyte random motility and chemotaxis based on receptor binding fluctuations.  The Journal of Cell Biology, 106(2):303–309, 1988.
  * [81] A. Uatay.  A mathematical model of contact inhibition of locomotion: coupling contractility and focal adhesions.  arXiv:1901.02386, 2019.
  * [82] U.S.Schwarz and M.L.Gardel.  United we stand integrating the actin cytoskeleton and cell matrix adhesions in cellular mechanotransduction.  Journal of Cell Science, 125(13):3051–3060, 2012.
  * [83] A. B. Verkhovsky and G. G. Borisy.  Non-sarcomeric mode of myosin ii organization in the fibroblast lamellum.  The Journal of Cell Biology, 123(3):637–652, 1993.
  * [84] A. B. Verkhovsky, T. M. Svitkina, and G. G. Borisy.  Self-polarization and directional motility of cytoplasm.  Current Biology, 9(1):11–S1, 1999.
  * [85] M. Vicente-Manzanares, J. Zareno, L. Whitmore, C. K. Choi, and A. F. Horwitz.  Regulation of protrusion, adhesion dynamics, and polarity by myosins IIA and IIB in migrating cells.  The Journal Of Cell Biology, 176(5):573–580, 2007.
  * [86] M. D. Vose.  A linear algorithm for generating random numbers with a given distribution.  IEEE Transactions on Software Engineering, 17(9):972–975, 1991\.
  * [87] P. T. Yam, C. A. Wilson, L. Ji, B. Hebert, E. L. Barnhart, N. A. Dye, P. W. Wiseman, G. Danuser, and J. A. Theriot.  Actin–myosin network reorganization breaks symmetry at the cell rear to spontaneously initiate polarized cell motility.  Journal of Cell Biology, 178(7):1207–1221, 2007.
  * [88] R. Zaidel-Bar, M. Cohen, L. Addadi, and B. Geiger.  Hierarchical assembly of cell–matrix adhesion complexes.  Biochemical Society Transactions, 32(3):416–420, 2004.
  * [89] J. Zhu and A. Mogilner.  Comparison of cell migration mechanical strategies in three-dimensional matrices: a computational study.  Interface focus, 6(5):20160040, 2016.
  * [90] F. Ziebert and I. S. Aranson.  Computational approaches to substrate-based cell motility.  npj Computational Materials, 2:16019, 2016.