# WasteNet: Waste Classification at the Edge for Smart Bins

Gary White, Christian Cabrera, Andrei Palade, Fan Li, Siobhán Clarke  School
of Computer Science and Statistics,
Trinity College Dublin, Dublin, Ireland
{whiteg5, cabrerac, paladea, fali<EMAIL_ADDRESS>
###### Abstract

Smart Bins have become popular in smart cities and campuses around the world.
These bins have a compaction mechanism that increases the bins’ capacity as
well as automated real-time collection notifications. In this paper, we
propose WasteNet, a waste classification model based on convolutional neural
networks that can be deployed on a low power device at the edge of the
network, such as a Jetson Nano. The problem of segregating waste is a big
challenge for many countries around the world. Automated waste classification
at the edge allows for fast intelligent decisions in smart bins without
needing access to the cloud. Waste is classified into six categories: paper,
cardboard, glass, metal, plastic and other. Our model achieves a 97%
prediction accuracy on the test dataset. This level of classification accuracy
will help to alleviate some common smart bin problems, such as recycling
contamination, where different types of waste become mixed with recycling
waste causing the bin to be contaminated. It also makes the bins more user
friendly as citizens do not have to worry about disposing their rubbish in the
correct bin as the smart bin will be able to make the decision for them.

##  I Introduction

Increases in urbanisation and economic development around the world has led to
increased waste generation [1, 2]. In the past century waste production has
risen tenfold and by 2025 it is expected to double again [3]. This waste is
often not recycled properly and is usually either dumped in a landfill or
incinerated [4]. Landfills, such as Laogang in Shanghai, Sudokwon in Seoul,
Jardim Gramacho in Rio de Janeiro and Bordo Poniente in Mexico City are among
the world’s largest and show that this is a global problem [5]. Waste
incinerators are being used to deal with the waste problem, but there are a
number of concerns over fine particle emission, which can lead to a number of
associated health risks as well as contributing to global warming [6, 7]. This
can lead to the accumulation of polluting materials that threaten the
environment [8, 9, 10] as well as human [11, 12] and animal health [13, 14].
Waste products can be classified based on various factors, such as
consumption, production, chemical and physical properties to allow for
effective reuse and recycling. The ability to reduce and recycle waste
material more effectively not only reduces the impact on the environment, but
is also more cost effective [15, 16].

Waste classification can happen at any stage in the waste management pipeline,
but is it better to happen when the waste is initially being disposed of to
avoid potential recyclable materials becoming contaminated. Solar powered
compaction bins have become increasingly popular in smart cities and campuses
as they can allow for eight times the capacity of a traditional bin due to the
compaction and can also send an alert when the bin is full [17]. However,
these bins do not have any waste classification mechanism, which can lead to a
recycling bin becoming contaminated [18]. General waste can be classified into
the categories shown in Figure 1 to allow for more efficient recycling. The
most common contaminants in an urban environment are napkins, plastic food
wrapping, plastic bags, coffee cups, coffee sleeves, rubber gloves and plastic
medical waste [19, 20]. A number of schemes have been tried to inform and
better educate the public on the categories and benefits of waste
classification in different demographics [21]. An automated classification
mechanism would solve this problem and encourage user engagement as they would
not be worried about placing an item in a recycling bin that could possible
contaminate it.

Figure 1: Waste Categories

Deep neural networks have shown increased accuracy in a number of recent
benchmark competitions in machine learning and pattern recognition [22]. We
can frame the waste classification problem as an image classification problem
where we have a camera that takes an image of the waste and uses a deep neural
network to identify what class the waste belongs to. For image-based tasks
standard deep neural networks are hard to train and encounter problems when
trying to scale [23]. Traditional dense connections do not have the
translation-invariance property, so any slight change in the size or placement
of an image would be difficult to detect. Convolutional neural networks (CNNs)
solve this problem by using multiple hidden convolutional layers that extract
high level features. Convolutional layers consist of learnable parameters,
called filters that go through the whole image i.e., it steps across the
width, length and colour channels of the image and calculates the inner
product of the input and the filter. This leads to a feature map of the
filter, which is used to classify the images.

Another building block of convolution networks, is a pooling layer, which
operates on the feature maps. The pooling layer reduces the spatial size of
the representation to reduce the chance of overfitting and cuts down the
number of parameters. Max pooling is an approach that divides the space into
individual regions and picks the maximum value for each region. CNNs also make
use of the ReLu activation function of the form $f(x)=max(0,x)$ as this leads
to faster training without effecting the accuracy of the network [24]. The
combination of these recent advancements has allowed CNNs to be used in a
range of computer vision tasks, such as: face detection [25], license plate
recognition [26], multi-label classification [27], noise recognition [28] and
object detection [29].

In this paper we present WasteNet, a deep neural network model that improves
waste classification accuracy. This model is tested on a Jetson Nano edge
device that can easily be deployed at the edge to allow waste classification
in smart bins. Our model achieves an improved state of the art accuracy score
of 97% on the TrashNet dataset [30]. This is a large improvement on the
original SVM approach that achieved an accuracy of 63% and CNN approach that
achieved an accuracy of 22% [30]. A recent survey of automated image-based
waste classification papers, has shown a state of the art accuracy on the
TrashNet dataset is 88.42% [31].

The paper is organised as follows, Section II presents the related work that
has been conducted on waste classification. Section III presents the design of
the deep learning model and the categories of waste that it uses when making
the classification decision. Section IV presents the experimental setup that
was used to evaluate our proposed deep learning model against other state of
the art approaches and Section V presents the results of those experiments.
Section VI concludes the paper and presents some future work.

##  II Related Work

Traditional approaches to waste classification have focused on the physical
properties, such as the weight, form (solid, liquid, aqueous or gaseous) and
the type of process that has generated the waste [32]. These properties can
then be used by a machine learning [33] or fuzzy set theory [34] algorithm to
classify the waste according to the set categories. For incinerators, the
classification of waste can be for the thermal conversion in waste-to-energy
research [35]. In vitro testing, such as Microtex, which uses bioluminescent
bacteria to detect hazardous waste has also been conducted [36]. Alternative
techniques have also been used to estimate waste generation using GIS-based
modelling [37] and multi-agent systems to simulate household behaviour [38].

In recent years, visual methods using either images in the visible [39] or
infrared [33] spectrum have become popular. Image classification has become a
major research area thanks to the release of large publicly available
datasets, such as ImageNet [40], which has led to the development of large
deep neural network models. Recent computer vision approaches have been used
to localise and classify waste on the streets [41, 42]. This allows street
cleaning equipment to find the areas that have the most waste and to focus on
cleaning those areas. Computer vision approaches have also been applied to
internal cleaning robots to identify and detect the type of waste that the
robot is about to clean on the floor [43, 44, 45]. This can also be applied at
a larger urban scale with Urban Swarms for autonomous waste management [46].
At the moment these urban swarms are proposed to allow for the efficient
collection of waste throughout the city from bins, but could easily be
extended to contain waste localisation and classification technology to allow
the agents in the swarm to automatically pick up and classify waste as they
move through the streets.

Previous image-based waste classification systems have focused on specific
subsections of waste classification, such as classifying plastic into:
polyethylene terephthalate, high-density polyethylene, polypropylene and
polystyrene [47] or focusing on household waste [48]. Our approach is designed
to be deployed at the edge and so would serve as a first pass classification
mechanism that first divides the waste into the classes defined in the
TrashNet dataset: plastic, glass, paper, metal, cardboard and other. This
first pass then allows for further classification of specific classes, such as
plastic at dedicated plastic recycling centres. Other waste classification
approaches make use of other devices, such as inductors and bridge sensors to
influence the classification of waste [49]. In this paper we focus on the use
of image classification for the classification of the waste at the edge of the
network.

Figure 2: Urban Smart Bin

##  III WasteNet Design

In urban environments self-compacting bins have become popular because due to
a compression mechanism the can hold eight time more waste than a typical bin
[17]. However, these bins compact all the waste together without
classification, which can lead to recyclable material becoming contaminated.
In this paper we propose using a deep neural network-based model to classify
the waste as it is added to the bin. Figure 2 shows how these bins accept
waste in a separate container at the top, which could easily be fitted with a
camera to allow for images of the waste to be taken. These images would then
be fed to our model, which would select the correct container for the waste to
be disposed in to ensure that it would be recycled properly. These bins are
also fitted with solar panels allowing for a low power edge device, such as a
Nvidia Jetson Nano to run the model locally. The transparent section at the
bottom of the bin shows how plastic, cardboard and paper are currently mixed
and squashed together at the bottom of the bin, without waste classification.
The remainder of this section describes the design of the WasteNet neural
network model.

###  III-A Transfer Learning

The design of the WasteNet model uses transfer learning to leverage knowledge
from a source task [44, 50]. We leverage knowledge from models trained for a
general image classification task on the ImageNet dataset [40]. Transfer
learning can also provide a number of benefits, such as improving baseline
performance, speeding up overall model development and training time and also
getting overall improved model performance compared to building the model from
scratch. This is especially important in deep learning where models can have a
very long and intensive training time.

Figure 3 shows the different categories of transfer learning based on the
relationship between the source and target distributions. The most simple case
is regular learning where the source and target have the same distributions
and are required to perform the same tasks. When the source and target have
the same distribution or are in the same domain but the tasks that they are
required to perform are different this is called inductive transfer learning.
This category can further be broken down depending upon whether the source
domains contain labelled data or not: if a lot of labelled data in the source
domain are available then it is multi-task learning and if there are no
labelled data from the source domain then it is self-taught learning [51].
When the source and target distribution are not the same but the tasks are
similar it is called transductive transfer learning. In this situation, no
labelled data in the target domain are available, while a lot of data in the
source domains are available. The final category is unsupervised transfer
learning, where there is a difference in both the source and target
distribution and tasks. This category focuses on solving unsupervised tasks in
the target domain such as clustering [52] and dimensionality reduction [53],
with no labelled data available in the source and target domains in training.
A linear cost function can be used to minimise the difference between the
source and target domain distribution in unsupervised transfer learning [54,
55]. In this paper we focus on inductive transfer learning, where the original
model is trained on the same source and target distribution, but performs a
different tasks of waste classification.

Figure 3: Transfer Learning Categories

####  III-A1 Feature Extraction

Deep learning models are layered architectures that learn different features
at different layers. These layers are finally connected to a final layer,
which is usually fully connected in the case of classification to get the
final output class. This layered architecture allows us to utilize pre-trained
networks without a final layer as a fixed feature-extractor for other tasks.
Figure 4 shows the different transfer learning approaches that can be applied
to this task. Figure 4a shows the VGG-16 model that can be downloaded directly
to an edge device [56]. Figure 4b shows the VGG model as a feature extractor
where we freeze (fix weights and don’t train) all the blocks of convolutions
layers and the flattening layer in the dashed blue box. We only update the
fully connected classifier block at the end of the model. This allows the new
model in this case to transform the image from a new domain task into a large
dimension vector based on hidden states, thus enabling us to extract features
from a new domain task using the source domain. This is one of the most widely
used methods of transfer learning for deep neural networks.

(a)

(b)

(c)

Figure 4: Deep Transfer Learning Approaches

####  III-A2 Fine Tuning

In fine tuning, the weights of the last few layers of the network are updated
as shown in Figure 4c with a pink box and trained as well as the fully
connected layers at the end of the model, for the classification task. This
means that this method is slightly more resource intensive as we have to train
some of the previous layers. As deep neural networks are layered with the
initial layers capturing the most basic features, such as edges and the later
layers capturing more specific details about the task, we can freeze some of
the first blocks and update the later ones. In Figure 4c we freeze the first
three layers in the dashed blue box, while fine tuning the last two blocks to
suit the task. This allows us to use the knowledge in terms of the overall
architecture of the network and use its states as the starting point for our
retraining step allowing us to achieve better performance in less time.

One of the problems with updating a model using fine-tuning is that some of
the parameters in the non-frozen layers have to be updated to solve the new
problem. This can overwrite parameters that the network has learned before
leading to “catastrophic forgetting” of the knowledge that was previously
acquired [57, 58]. Previous approaches to preventing catastrophic forgetting
have used an ensemble of neural networks. When there is a new task, the
algorithm creates a new network and shares the representation between the
tasks [59, 60]. However, these approaches are not suitable for the edge due to
the space and complexity restrictions as the number of networks increases
linearly with the number of new tasks to be learned.

####  III-A3 Hybrid Tuning

In this paper we propose a hybrid transfer learning approach that combines the
benefits of feature extraction and fine-tuning. Our hybrid approach has two
stages, in the first pre-training stage the base network is used as a feature
extractor by freezing the lower layers of the network and only updating the
weights of the top layer. Once the loss function begins to stabilise and the
network has reached a high level of accuracy with the new top layer of the
network, the remaining layers of the network are gradually unfrozen. The
remaining layers of the network are trained using discriminative layer
training, where a different learning rate is applied to each layer. This
allows us to apply a lower learning rate to the low-level layer
representations and adapt the weights of the higher-level layers faster as
they contain more domain specific information.

This gradual unfreezing rather than fine-tuning all the layers at once reduces
the risk of catastrophic forgetting. We first unfreeze the top layer as this
contains the least general knowledge and update the weights [61]. We then
unfreeze the next lower frozen layer and repeat using the updated learning
rate, until we fine-tune all layers with convergence at the last iteration. To
generate the WasteNet model we apply our hybrid tuning approach to a DenseNet
model that has been trained on ImageNet [24].

###  III-B Data Augmentation

Generating artificial data based on existing observations is a technique in
machine learning to control overfitting, improve model accuracy and
generalisation [62]. The idea behind this technique also known as data
augmentation is that we follow a set process, taking existing data, such as
images from our training set and applying some image transformation operations
on them, such as translation, zooming, shearing and rotation to produce a new
set of alternative images. The randomness of the process means that we do not
get the same images each time. This helps to stop the deep learning model from
overfitting on the local training data. In our experiments we use the
ImageDataGenerator
class111https://keras.io/preprocessing/image/#imagedatagenerator-class from
Keras to provide a number of transformations for generating new images, such
as: zooming, rotation, translation, randomly flipping images horizontally, and
filling new pixels with their nearest surround pixels.

##  IV Experimental Setup

###  IV-A Dataset

To evaluate our model we use the TrashNet dataset [30]. The dataset contains
multiple classes: plastic, paper, glass, metal, cardboard and other. It also
contains a range of images of waste in different orientations and positions as
seen in Figure 5. The exposure and lighting for each of the photos is also
varied in the dataset. Each image is resized to 512x384 with 3 colour
channels. The original dataset is 3.5GB in size, containing 2527 images in
total. The number of images in each class varies with 594 paper, 501 glass,
410 metal, 482 plastic, 493 cardboard and 137 for other materials. For the
evaluation, the images are split into a train, validation and test set with a
50, 25, 25 split.

(a)

(b)

(c)

(d)

(e)

(f)

(g)

(h)

(i)

Figure 5: Images of Waste Categories

###  IV-B Metrics

We use a number of classification metrics to comprehensively show how the
classification models perform. The metrics are based on the number of True
Positives (TP), False Positives (FP), True Negatives (TN) and False Negatives
(FN) and are defined as:

$Accuracy=\frac{True\>Positive+True\>Negative}{TP+TN+FP+FN}$ (1)
$Precision=\frac{True\>Positive}{True\>Positive+False\>Positive}$ (2)
$TPR/Recall=\frac{True\>Positive}{True\>Positive+False\>Negative}$ (3)
$F1=\frac{2\cdot Precision\cdot Recall}{Precision+Recall}$ (4)

###  IV-C Performance Comparison

We compare the results of our classification approach against a range of other
state of the art deep learning models. For each model we also use a specific
instance of the models with a different numbers of layers. This can be seen in
the ResNet model, where we use ResNet34, Resnet50 and ResNet101 to evaluate
the impact that model size has on the final accuracy metrics [63]. We also
evaluate against other state of the art deep neural network models, such as
DenseNet [64], SqeezeNet [65], AlexNet [66] and VGG [56].

##  V Results

###  V-A Classification Accuracy

Table I shows the results for each of the trash classification methods. We
evaluate our WasteNet approach against 10 other state of the art deep neural
networks using 4 different metrics for a comprehensive analysis. We can see
that the WasteNet approach shows an increase in accuracy compared to existing
state of the art approaches, such as DenseNet169 and ResNet101. There is a
0.017 difference in accuracy from WasteNet compared to the next most accurate
model. There is a large difference in accuracy between the WasteNet and
AlexNet model, with the AlexNet model achieving an accuracy of 0.784, which is
0.186 less than WasteNet. The precision and recall metrics follow a similar
pattern with WasteNet achieving the best results. The F1 Score, which is a
combination of the precision and recall values are also highest for the
WasteNet model.

Model | Accuracy | Precision | Recall | F1 Score
---|---|---|---|---
WasteNet | 0.970 | 0.970 | 0.970 | 0.970
DenseNet169 | 0.953 | 0.954 | 0.953 | 0.953
DenseNet121 | 0.942 | 0.942 | 0.942 | 0.941
ResNet101 | 0.939 | 0.939 | 0.939 | 0.938
ResNet50 | 0.937 | 0.937 | 0.937 | 0.937
VGG16 | 0.928 | 0.927 | 0.928 | 0.927
ResNet34 | 0.918 | 0.918 | 0.918 | 0.918
DenseNet161 | 0.913 | 0.917 | 0.913 | 0.912
VGG19 | 0.907 | 0.908 | 0.907 | 0.907
SqueezeNet | 0.811 | 0.814 | 0.811 | 0.812
AlexNet | 0.784 | 0.786 | 0.784 | 0.782

Table I: Waste Classification Accuracy

Figure 6 shows the training and validation loss of the network as it is being
trained. The orange and blue lines show the pre-training and pre-validation
loss. This is when only the top layer of the network is being trained and the
remaining layers are frozen for the first 100 epochs. We can see that during
this period the loss function is high as the new top layer of the network is
learning to classify the waste material.

Figure 6: Loss Over Epochs

At around 100 epochs the loss value begins to stabilise and falls to a lower
amount. This is when the second stage of training begins, where we gradually
unfreeze the remaining layers of the network and begin training the layers at
different learning rates. The training loss at this stage is shown by the
green line and the validation loss is shown be the red line. We can see how
gradually unfreezing the layers reintroduced some variation especially in the
validation loss function. This can be helpful to avoid the network becoming
stuck in local maximums or minimums. After 1000 epochs both the validation and
training loss remain quite stable.

Figure 7 shows the pre-training and training accuracy. The blue line shows the
network during pre-training, when the top layer of the network is trained and
all the others remain frozen. The orange line shows the network when the
remaining layers are being gradually unfrozen and trained at different
learning rates. Similar to Figure 6 there is a lot of variation in the pre-
training stage, when the network is beginning to classify the waste. We can
see that accuracy in this stage varies a lot starting at around 0.9 and
dropping to 0.75. Once the network begins to stabilise the accuracy at around
100 epochs the remaining layers are gradually unfrozen to allow for further
accuracy improvement. The accuracy still varies a lot with a noticeable
increase around 700 epochs, when there is a lot of variation in the validation
loss function in Figure 6. After this increase the network accuracy remains
quite stable.

Figure 7: Accuracy Over Epochs

###  V-B Classification Confusion Matrix

Figure 8 shows the confusion matrix for our WasteNet approach on the test
dataset. The y-axis shows the true labels and the x-axis shows the predicted
labels. We can see that most of the labels have been correctly predicted as
the middle top left to bottom right diagonal has most of the results. Overall,
the model has achieved a good prediction accuracy of 0.970. However, we can
also see the classes that have been mislabelled in the confusion matrix. The
most confused classification is predicting metal when the true label is glass,
which happens 5 times. The second most confused classification is also when
the material is glass, but the prediction is plastic.

Figure 8: Waste Classification Confusion Matrix

To get further insight into how the deep neural network is making the
predictions we can use gradient-based localisation to identify the areas in
the image that caused the model to make a prediction for a particular class
[67]. Figure 9 shows a heat map of the regions of the image that cause the
network to fire when making the classification decision. These images are
taken from the training dataset and are the images that were originally
misclassified by the network. Figure 9a shows an image where the network has
made the prediction metal, when the image was actually glass. The network has
assigned a probability of 0.00 so was not confident about the prediction. This
caused a big change in the network when the true label was revealed, leading
to a loss of 9.28. Figure 9e shows an example where the network has predicted
that the bottle is plastic, when it turned out to be glass. We can see that
the network has focused on the specific areas to identify the bottle such as
the neck and base of the cylinder. As the network is trained on more images in
the training dataset it becomes better at classifying the difference between
plastic and glass bottles.

(a) test

(b)

(c)

(d)

(e)

(f)

Figure 9: Prediction/Actual Loss/Probability

We can see in all these cases that the model had a low probability value close
to 0.00. This could be used to avoid these misclassifications in the real-
world, where the model would only classify waste when it had a probability of
greater than 0.9. Waste with a probability less than this could be separated
into a different bin to allow for human sorting on that subset.

##  VI Conclusion and Future Work

The results of the experimentation have provided a number of interesting
insights for future waste classification. The very high accuracy of 97% that
we are able to achieve shows the potential of these methods. As these models
can now be deployed on devices at the edge of the network, such as a Jetson
Nano, it allows decision making and artificial intelligence at the edge of the
network. The models can be readily deployed in the next versions of smart bins
to ensure that the compressed waste is all of the same type. This would
increase the amount of waste being recycled, while also reducing the human
workload. The method can also be applied to large sorting distribution centres
to automatically sort waste into different categories to allow for easier
recycling.

For future work we plan on focusing on waste classification at larger
distributions centres, which typically use video feeds. We are planning to
updated our model and collect additional data to allow for the automatic
classification of waste at a large industrial scale where objects may be
overlapping and moving on a conveyor belt.

## Acknowledgment

This research has been funded by Science Foundation Ireland (SFI) under grant
13/IA/1885. The Titan V used for this research was donated by the NVIDIA
Corporation.

## References

  * [1] N. J. Bandara, J. P. A. Hettiaratchi, S. Wirasinghe, and S. Pilapiiya, “Relation of waste generation and composition to socio-economic factors: a case study,” _Environmental monitoring and assessment_ , vol. 135, no. 1-3, pp. 31–39, 2007.
  * [2] A. D. Vethaak and H. A. Leslie, “Plastic debris is a human health issue,” 2016\.
  * [3] D. Hoornweg, P. Bhada-Tata _et al._ , _What a waste: a global review of solid waste management_ , 2012, vol. 15.
  * [4] A. M. King, S. C. Burgess, W. Ijomah, and C. A. McMahon, “Reducing waste: repair, recondition, remanufacture or recycle?” _Sustainable development_ , vol. 14, no. 4, pp. 257–267, 2006.
  * [5] H. Robinson _et al._ , “The composition of leachates from very large landfills: an international review,” _Commun Waste Res Manag_ , vol. 8, no. 1, pp. 19–32, 2007.
  * [6] J. Zhou, S. Wu, Y. Pan, L. Zhang, Z. Cao, X. Zhang, S. Yonemochi, S. Hosono, Y. Wang, K. Oh _et al._ , “Enrichment of heavy metals in fine particles of municipal solid waste incinerator (mswi) fly ash and associated health risk,” _Waste Management_ , vol. 43, pp. 239–246, 2015.
  * [7] I. Naroznova, J. Møller, and C. Scheutz, “Global warming potential of material fractions occurring in source-separated organic household waste treated by anaerobic digestion or incineration under different framework conditions,” _Waste management_ , vol. 58, pp. 397–407, 2016.
  * [8] J. P. da Costa, P. S. Santos, A. C. Duarte, and T. Rocha-Santos, “(nano) plastics in the environment–sources, fates and effects,” _Science of the Total Environment_ , vol. 566, pp. 15–26, 2016.
  * [9] J. Hopewell, R. Dvorak, and E. Kosior, “Plastics recycling: challenges and opportunities,” _Philosophical Transactions of the Royal Society B: Biological Sciences_ , vol. 364, no. 1526, pp. 2115–2126, 2009.
  * [10] A. Cózar, M. Sanz-Martín, E. Martí, J. I. González-Gordillo, B. Ubeda, J. Á. Gálvez, X. Irigoien, and C. M. Duarte, “Plastic accumulation in the mediterranean sea,” _PLoS One_ , vol. 10, no. 4, p. e0121762, 2015.
  * [11] K.-H. Kim, E. Kabir, and S. Kabir, “A review on the human health impact of airborne particulate matter,” _Environment international_ , vol. 74, pp. 136–143, 2015.
  * [12] S. L. Wright and F. J. Kelly, “Plastic and human health: a micro issue?” _Environmental science & technology_, vol. 51, no. 12, pp. 6634–6647, 2017\.
  * [13] K. Mattsson, S. Jocic, I. Doverbratt, and L.-A. Hansson, “Nanoplastics in the aquatic environment,” in _Microplastic Contamination in Aquatic Environments_. Elsevier, 2018, pp. 379–399.
  * [14] J. G. Derraik, “The pollution of the marine environment by plastic debris: a review,” _Marine pollution bulletin_ , vol. 44, no. 9, pp. 842–852, 2002\.
  * [15] A. U. Zaman, “A comprehensive study of the environmental and economic benefits of resource recovery from global waste management systems,” _Journal of cleaner production_ , vol. 124, pp. 41–50, 2016.
  * [16] Y. Xu, J. Li, Q. Tan, A. L. Peters, and C. Yang, “Global status of recycling waste solar panels: A review,” _Waste Management_ , vol. 75, pp. 450–458, 2018.
  * [17] J. Poss, J. Satwicz, B. Richmond, and M. Taylor, “Solar powered compaction apparatus,” Oct. 24 2006, uS Patent 7,124,680.
  * [18] A. Andrews, M. Gregoire, H. Rasmussen, and G. Witowich, “Comparison of recycling outcomes in three types of recycling collection units,” _Waste management_ , vol. 33, no. 3, pp. 530–535, 2013.
  * [19] A. Heathcote, T. Wilson, K. Milnes, J. Thompson, J. Woods, and D. Zipursky, “Conducting a waste audit in the killam library at dalhousie university to establish waste habits and locate problem areas.”
  * [20] N. Sajjad, S. Hassan, J. Qadir, R. Ali, and D. Shah, “Analyses of the recycling potential of medical plastic wastes,” in _Handbook of Research on Environmental and Human Health Impacts of Plastic Pollution_. IGI Global, 2020, pp. 178–199.
  * [21] M. Mourad, “Recycling, recovering and preventing “food waste”: competing solutions for food systems sustainability in the united states and france,” _Journal of Cleaner Production_ , vol. 126, pp. 461–477, 2016.
  * [22] G. White and S. Clarke, “Urban intelligence with deep edges,” _IEEE Access_ , vol. 8, pp. 7518–7530, 2020.
  * [23] D. C. Ciresan, U. Meier, J. Masci, L. M. Gambardella, and J. Schmidhuber, “Flexible, high performance convolutional neural networks for image classification,” in _Twenty-Second International Joint Conference on Artificial Intelligence_ , 2011.
  * [24] A. Krizhevsky, I. Sutskever, and G. E. Hinton, “Imagenet classification with deep convolutional neural networks,” in _Advances in neural information processing systems_ , 2012, pp. 1097–1105.
  * [25] H. Jiang and E. Learned-Miller, “Face detection with the faster r-cnn,” in _2017 12th IEEE International Conference on Automatic Face & Gesture Recognition (FG 2017)_. IEEE, 2017, pp. 650–657.
  * [26] R. Polishetty, M. Roopaei, and P. Rad, “A next-generation secure cloud-based deep learning license plate recognition for smart cities,” in _2016 15th IEEE International Conference on Machine Learning and Applications (ICMLA)_. IEEE, 2016, pp. 286–293.
  * [27] J. Wang, Y. Yang, J. Mao, Z. Huang, C. Huang, and W. Xu, “Cnn-rnn: A unified framework for multi-label image classification,” in _Proceedings of the IEEE conference on computer vision and pattern recognition_ , 2016, pp. 2285–2294.
  * [28] J. Cao, M. Cao, J. Wang, C. Yin, D. Wang, and P.-P. Vidal, “Urban noise recognition with convolutional neural network,” _Multimedia Tools and Applications_ , pp. 1–21, 2018.
  * [29] X. Chen, H. Ma, J. Wan, B. Li, and T. Xia, “Multi-view 3d object detection network for autonomous driving,” in _Proceedings of the IEEE Conference on Computer Vision and Pattern Recognition_ , 2017, pp. 1907–1915.
  * [30] M. Yang and G. Thung, “Classification of trash for recyclability status.” [Online]. Available: http://cs229.stanford.edu/proj2016/report/ThungYang-ClassificationOfTrashForRecyclabilityStatus-report.pdf
  * [31] V. Ruiz, Á. Sánchez, J. F. Vélez, and B. Raducanu, “Automatic image-based waste classification,” in _From Bioinspired Systems and Biomedical Applications to Machine Learning_ , J. M. Ferrández Vicente, J. R. Álvarez-Sánchez, F. de la Paz López, J. Toledo Moreo, and H. Adeli, Eds. Cham: Springer International Publishing, 2019, pp. 422–431.
  * [32] D. Zekkos, E. Kavazanjian Jr, J. D. Bray, N. Matasovic, and M. F. Riemer, “Physical characterization of municipal solid waste for geotechnical purposes,” _Journal of Geotechnical and Geoenvironmental Engineering_ , vol. 136, no. 9, pp. 1231–1241, 2010.
  * [33] P. Kuritcyn, K. Anding, E. Linß, and S. Latyev, “Increasing the safety in recycling of construction and demolition waste by using supervised machine learning,” in _Journal of Physics: Conference Series_ , vol. 588, no. 1. IOP Publishing, 2015, p. 012035\.
  * [34] N. Musee, L. Lorenzen, and C. Aldrich, “New methodology for hazardous waste classification using fuzzy set theory: Part i. knowledge acquisition,” _Journal of hazardous materials_ , vol. 154, no. 1-3, pp. 1040–1051, 2008\.
  * [35] H. Zhou, Y. Long, A. Meng, Q. Li, and Y. Zhang, “Classification of municipal solid waste components for thermal conversion in waste-to-energy research,” _Fuel_ , vol. 145, pp. 151–157, 2015.
  * [36] G. Bitton, K. Jung, and B. Koopman, “Evaluation of a microplate assay specific for heavy metal toxicity,” _Archives of Environmental Contamination and Toxicology_ , vol. 27, no. 1, pp. 25–28, 1994.
  * [37] N. V. Karadimas and V. G. Loumos, “Gis-based modelling for the estimation of municipal solid waste generation and collection,” _Waste Management & Research_, vol. 26, no. 4, pp. 337–346, 2008.
  * [38] X. Meng, Z. Wen, and Y. Qian, “Multi-agent based simulation for household solid waste recycling behavior,” _Resources, Conservation and Recycling_ , vol. 128, pp. 535–545, 2018.
  * [39] G. E. Sakr, M. Mokbel, A. Darwich, M. N. Khneisser, and A. Hadi, “Comparing deep learning and support vector machines for autonomous waste sorting,” in _2016 IEEE International Multidisciplinary Conference on Engineering Technology (IMCET)_. IEEE, 2016, pp. 207–212.
  * [40] J. Deng, W. Dong, R. Socher, L.-J. Li, K. Li, and L. Fei-Fei, “ImageNet: A Large-Scale Hierarchical Image Database,” in _CVPR09_ , 2009.
  * [41] M. S. Rad, A. von Kaenel, A. Droux, F. Tieche, N. Ouerhani, H. K. Ekenel, and J.-P. Thiran, “A computer vision system to localize and classify wastes on the streets,” in _International Conference on Computer Vision Systems_. Springer, 2017, pp. 195–204.
  * [42] P. Zhang, Q. Zhao, J. Gao, W. Li, and J. Lu, “Urban street cleanliness assessment using mobile edge computing and deep learning,” _IEEE Access_ , vol. 7, pp. 63 550–63 563, 2019.
  * [43] B. Ramalingam, A. Lakshmanan, M. Ilyas, A. Le, and M. Elara, “Cascaded machine-learning technique for debris classification in floor-cleaning robot application,” _Applied Sciences_ , vol. 8, no. 12, p. 2649, 2018.
  * [44] L. Torrey and J. Shavlik, “Transfer learning,” in _Handbook of research on machine learning applications and trends: algorithms, methods, and techniques_. IGI Global, 2010, pp. 242–264.
  * [45] B. Ramalingam, A. Lakshmanan, M. Ilyas, A. Le, and M. Elara, “Cascaded machine-learning technique for debris classification in floor-cleaning robot application,” _Applied Sciences_ , vol. 8, no. 12, p. 2649, 2018.
  * [46] A. L. Alfeo, E. C. Ferrer, Y. L. Carrillo, A. Grignard, L. A. Pastor, D. T. Sleeper, M. G. Cimino, B. Lepri, G. Vaglini, K. Larson _et al._ , “Urban swarms: A new approach for autonomous waste management,” in _2019 International Conference on Robotics and Automation (ICRA)_. IEEE, 2019, pp. 4233–4240.
  * [47] J. Bobulski and M. Kubanek, “Waste classification system using image processing and convolutional neural networks,” in _Advances in Computational Intelligence_ , I. Rojas, G. Joya, and A. Catala, Eds. Cham: Springer International Publishing, 2019, pp. 350–361.
  * [48] C. Srinilta and S. Kanharattanachai, “Municipal solid waste segregation with cnn,” in _2019 5th International Conference on Engineering, Applied Sciences and Technology (ICEAST)_ , July 2019, pp. 1–4.
  * [49] Y. Chu, C. Huang, X. Xie, B. Tan, S. Kamal, and X. Xiong, “Multilayer hybrid deep-learning method for waste classification and recycling,” _Computational Intelligence and Neuroscience_ , vol. 2018, 2018.
  * [50] G. White and S. Clarke, “Short-term qos forecasting at the edge for reliable service applications,” _IEEE Transactions on Services Computing_ , pp. 1–1, 2020.
  * [51] J. Han, J. Pei, and M. Kamber, _Data mining: concepts and techniques_. Elsevier, 2011.
  * [52] W. Dai, Q. Yang, G.-R. Xue, and Y. Yu, “Self-taught clustering,” in _Proceedings of the 25th international conference on Machine learning_. ACM, 2008, pp. 200–207.
  * [53] Z. Wang, Y. Song, and C. Zhang, “Transferred dimensionality reduction,” in _Joint European Conference on Machine Learning and Knowledge Discovery in Databases_. Springer, 2008, pp. 550–565.
  * [54] S. Nejatian, V. Rezaie, H. Parvin, M. Pirbonyeh, K. Bagherifard, and S. K. Syed Yusof, “An innovative linear unsupervised space adjustment by keeping low-level spatial data structure,” _Knowledge and Information Systems_ , vol. 59, no. 2, pp. 437–464, May 2019. [Online]. Available: https://doi.org/10.1007/s10115-018-1216-8
  * [55] A. Pirbonyeh, V. Rezaie, H. Parvin, S. Nejatian, and M. Mehrabi, “A linear unsupervised transfer learning by preservation of cluster-and-neighborhood data organization,” _Pattern Analysis and Applications_ , vol. 22, no. 3, pp. 1149–1160, Aug 2019. [Online]. Available: https://doi.org/10.1007/s10044-018-0753-9
  * [56] K. Simonyan and A. Zisserman, “Very Deep Convolutional Networks for Large-Scale Image Recognition,” _arXiv e-prints_ , p. arXiv:1409.1556, Sep 2014.
  * [57] R. M. French, “Catastrophic forgetting in connectionist networks,” _Trends in cognitive sciences_ , vol. 3, no. 4, pp. 128–135, 1999.
  * [58] R. Kemker, M. McClure, A. Abitino, T. Hayes, and C. Kanan, “Measuring Catastrophic Forgetting in Neural Networks,” _arXiv e-prints_ , p. arXiv:1708.02072, Aug 2017.
  * [59] S.-W. Lee, C.-Y. Lee, D.-H. Kwak, J. Kim, J. Kim, and B.-T. Zhang, “Dual-memory deep learning architectures for lifelong learning of everyday human behaviors,” in _Proceedings of the Twenty-Fifth International Joint Conference on Artificial Intelligence_ , ser. IJCAI’16. AAAI Press, 2016, pp. 1669–1675. [Online]. Available: http://dl.acm.org/citation.cfm?id=3060832.3060854
  * [60] A. A. Rusu, N. C. Rabinowitz, G. Desjardins, H. Soyer, J. Kirkpatrick, K. Kavukcuoglu, R. Pascanu, and R. Hadsell, “Progressive neural networks,” _arXiv preprint arXiv:1606.04671_ , 2016.
  * [61] J. Yosinski, J. Clune, Y. Bengio, and H. Lipson, “How transferable are features in deep neural networks?” in _Advances in neural information processing systems_ , 2014, pp. 3320–3328.
  * [62] M. D. Bloice, C. Stocker, and A. Holzinger, “Augmentor: an image augmentation library for machine learning,” _arXiv preprint arXiv:1708.04680_ , 2017.
  * [63] K. He, X. Zhang, S. Ren, and J. Sun, “Deep residual learning for image recognition,” _CoRR_ , vol. abs/1512.03385, 2015. [Online]. Available: http://arxiv.org/abs/1512.03385
  * [64] G. Huang, Z. Liu, L. van der Maaten, and K. Q. Weinberger, “Densely Connected Convolutional Networks,” _arXiv e-prints_ , p. arXiv:1608.06993, Aug 2016.
  * [65] F. N. Iandola, M. W. Moskewicz, K. Ashraf, S. Han, W. J. Dally, and K. Keutzer, “Squeezenet: Alexnet-level accuracy with 50x fewer parameters and <1mb model size,” _CoRR_ , vol. abs/1602.07360, 2016. [Online]. Available: http://arxiv.org/abs/1602.07360
  * [66] A. Krizhevsky, “One weird trick for parallelizing convolutional neural networks,” _CoRR_ , vol. abs/1404.5997, 2014. [Online]. Available: http://arxiv.org/abs/1404.5997
  * [67] R. R. Selvaraju, A. Das, R. Vedantam, M. Cogswell, D. Parikh, and D. Batra, “Grad-cam: Why did you say that? visual explanations from deep networks via gradient-based localization,” _CoRR_ , vol. abs/1610.02391, 2016. [Online]. Available: http://arxiv.org/abs/1610.02391


# Data Analysis using Riemannian Geometry
and Applications to Chemical Engineering

Alexander Smith¶, Benjamin Laubach†, Ivan Castillo†, and Victor M. Zavala¶
¶Department of Chemical and Biological Engineering
University of Wisconsin-Madison, 1415 Engineering Dr, Madison, WI 53706, USA
†Chemometrics, AI and Statistics
The Dow Chemical Company, 332 SH 332 E, Lake Jackson, TX 77566, Corresponding
Author<EMAIL_ADDRESS>
###### Abstract

We explore the use of tools from Riemannian geometry for the analysis of
symmetric positive definite matrices (SPD). An SPD matrix is a versatile data
representation that is commonly used in chemical engineering (e.g.,
covariance/correlation/Hessian matrices and images) and powerful techniques
are available for its analysis (e.g., principal component analysis). A key
observation that motivates this work is that SPD matrices live on a Riemannian
manifold and that implementing techniques that exploit this basic property can
yield significant benefits in data-centric tasks such classification and
dimensionality reduction. We demonstrate this via a couple of case studies
that conduct anomaly detection in the context of process monitoring and image
analysis.

##  1 Introduction

The assumption that data lies in a Euclidean space is pervasive throughout
science and engineering and is the basis of diverse data analysis techniques
used in these domains. Making this blanket assumption, however, is not always
appropriate and can affect the accuracy/interpretability of such techniques or
even break fundamental physical laws.

Recognizing that data can live in spaces that are governed by _non-Euclidean_
geometry is critical to appropriately representing, manipulating, and
analyzing certain data objects. A simple example of this arises when computing
travel paths between a couple of points that are located on the surface of the
Earth; when computing the distance between such points, the _elliptic_
geometry of the Earth surface must be taken into consideration. If Euclidean
geometry is assumed, travel paths between antipodal locations (e.g., United
States and China) can require traversing the Earth (not physically-realizable
paths).

Elliptic geometry is a non-Euclidean geometry in which one of the postulates
of Eucliean geometry (sum of interior angles of a triangle is equal to 180o)
no longer holds, because of the presence of positive curvature in a surface
[7]. Elliptic geometry was first proposed by Bernhard Riemann in the 19th
century, and he further developed these ideas into the field that is now known
as Riemannian geometry [30].

An important example in which assuming Euclidean geometry can lead to spurious
results is in the analysis of _symmetric positive definite_ (SPD) matrices
(e.g., covariance/correlation matrices). SPD matrices lie on a high-
dimensional space which is governed by Riemannian geometry (known as a
Riemannian manifold) [19]. Standard techniques for the analysis of SPD
matrices (e.g., PCA or basic matrix norms) do not take this property into
consideration and can lead to misleading results. For instance, computing the
distance between SPD matrices in Euclidean geometry (e.g., via the Frobenius
norm) ignores the fact that such matrices live on a Riemannian manifold, and
this can yield misleading results [12, 2, 34]. Specifically, the so-called
_swelling effect_ can occur when applying operations in Euclidean geometry to
SPD matrices [2]. This effect introduces spurious results by inflating the
determinants of SPD matrices and can also distort the results of commonly used
methods [34]. Computing interpolations and averages of SPD matrices, which is
key in understanding physical systems (e.g. Brownian motion), can also break
physical conservation laws if performed under Euclidean geometry [24].

In this work, we focus our attention on the use of techniques from Riemannian
geometry for the analysis of data objects that can be represented as SPD
matrices. An SPD matrix is a simple but versatile data representation that is
widely used in multivariate analysis techniques such as PCA [29, 39, 42]. SPD
representations are also used in process control, monitoring, and anomaly
detection [11, 16, 20, 40, 31, 33], in the study of functional brain networks
[28, 35, 38, 14, 33], in object detection [36, 25, 41], in biomedical image
analysis [24, 27], in the analysis of Laplacian matrices in graph theory, and
in the analysis of Hessian matrices in optimization [22]. In applications such
as image analysis, an SPD representation can be obtained by applying
transformations to a raw data object (e.g., smoothing via kernel functions
and/or combination of image features). The defining feature of an SPD matrix
is that all its eigenvalues are real and positive. In the context of
optimization, it is well-known that an SPD Hessian matrix defines positive
curvature of a multi-dimensional quadratic function and is key in defining the
geometry of objective and constraint functions (e.g., convex or non-convex).
In the context of statistics, it is well-known that an SPD covariance matrix
defines a multi-dimensional ellipsoid (a surface with positive curvature) and
that the level sets of a multivariate Gaussian probability density function
are ellipsoids.

This paper provides a practical introduction to the Riemannian geometry of SPD
matrices and demonstrates applications of interest to the process systems
engineering community. Specifically, we illustrate the benefits of exploiting
the Riemannian geometry of SPD matrices and discuss how these tools can be
incorporated into common dimensionality reduction and classification
workflows. We also provide example applications of interest in science and
engineering. The first application focuses on analysis of covariance matrices
derived from multivariate time series, which is a common task in process
monitoring. Our application focuses on the so-called Tennessee Eastman Process
(TEP) [9]. The TEP is a process where anomalies/faults are systematically
introduced which shifts the relationships between the measured variables.
Covariance matrices encode these changing relationships and are then used to
predict what type of anomaly the process is experiencing. The second
application is in defect/anomaly detection of textiles taken from the MVTEC AD
dataset [4]. Here, grayscale images of textiles are represented as covariance
matrices by characterizing the relationships between the original image and
multiple transformations of the image (e.g., smoothing with different
kernels). The image transformations emphasize different features of the
original image (e.g., edges and fibers). Subsequently, when an anomaly/defect
is introduced (e.g., a cut or discoloration of the textile) these
relationships will change, impacting the covariance matrix and allowing us to
detect anomalies. This type of analysis can also be applied to other relevant
image/field datasets such as those arising in microscopy and flow cytometry
[33, 32]. All data and scripts needed to conduct such analyses is provided as
open-source code in https://github.com/zavalab/ML/tree/master/RiemannianSPD.
With this, we aim to provide a concise and easy introduction to non-experts to
the field of Riemannian geometry.

##  2 Riemannian Manifolds

The key observation driving this work is that SPD matrices lie on a Riemannian
manifold; we thus begin our discussion by characterizing such manifolds. In
this section, we aim to provide an intuitive understanding of manifolds and
equip the reader with knowledge of their key properties.

An $n$-dimensional manifold is a topological space (a space where closeness
and connectedness are defined but not directly measurable) that locally
resembles an $n$-dimensional Euclidean space. Specifically, an $n$-dimensional
manifold $M$ is a set where every point $p\in M$ has an open neighborhood
$U\subset M$ (known as a _chart_) that can be mapped to an open set of
$n$-dimensional Euclidean space $V\subset\mathbb{R}^{n}$ via a one-to-one,
onto, and continuous mapping $f:U\rightarrow V$ (i.e., a homeomorphism). The
set of charts whose union covers the manifold is known as an _atlas_
$\bigcup_{i=1}^{n}U_{i}=M$. The chart/atlas nomenclature is derived from
navigation along the surface of the Earth (a $2$D manifold); here, charts are
flat (Euclidean) maps of the Earth that are collected in an atlas.

In Figure 1 we illustrate multiple topological spaces embedded in a 3D
Euclidean space. The spiked (b) and smooth (c) hollow spheres are examples of
$2$D manifolds. The space in (a) represents a couple of cones that intersect
at a single point and is not manifold; this is because a neighborhood drawn
around the intersecting point will look like a smaller version of the
intersecting cones, which cannot be mapped to $2$D Euclidean space through a
homeomorphism.

(a) Non-Manifold

(b) Non-Differentiable Manifold

(c) Differentiable Manifold

Figure 1: (a) Space composed of a couple of cones intersecting at a single
point. This is a non-manifold space because any neighborhood formed around the
intersecting point is not homeomorphic to 2D Euclidean space (the neighborhood
is a smaller version of the two intersecting cones). (b) Represents a 2D
manifold (all points and associated neighborhoods can be mapped to 2D
Euclidean space) but is not a differentiable manifold because of the cusps
occurring at the edges of the manifold (differential is not defined
everywhere). (c) A smooth sphere is a 2D manifold that is also differentiable
(curves on the surface can be differentiated everywhere).

Manifolds can also be endowed with geometric structure; for our analysis, we
are particularly interested in whether or not a given manifold is
_differentiable_. In simple terms, a differentiable manifold is a manifold for
which calculus (e.g., computing derivatives and integrals) can be performed on
the charts that make up the manifold atlas [18]. This also means that _curves_
on the surface of the manifold can be analyzed from a geometric perspective
using calculus. A curve is defined as a continuous function
$\gamma:[a,b]\rightarrow M$ mapping the interval $[a,b]\in\mathbb{R}$ to the
manifold $M$. We will not cover the specific mathematical requirements that
make a manifold differentiable in general, but refer interested readers to the
following reference for details [18]. Examples of a differentiable and non-
differentiable manifold are shown in Figure 1. We can see that the spiked
sphere in (b) has multiple cusps where a derivative cannot be defined for a
curve, making it non-differentiable. For the smooth sphere in (c), a
derivative can be taken anywhere, allowing for more complex
operations/transformations to be performed on the manifold [23].

We now restrict our attention to a special class of manifolds known as
_Riemannian manifolds_. A Riemannian manifold is a differentiable manifold
equipped with a defined _tangent space_ at each point in the manifold $p\in
M$, denoted as $T_{p}M$ [19]. The tangent space is the set of tangent vectors
of all curves passing through point $p\in M$. The tangent space $T_{p}M$ is a
vector (linear) space that is of the same dimension as the manifold itself. An
illustration of a tangent space is shown in Figure 2 for a smooth sphere (a
$2$D Riemannian manifold). For a Riemannian manifold, the tangent space
$T_{p}M$ is equipped with an inner product $g_{p}:T_{p}M\times
T_{p}M\rightarrow\mathbb{R}$, along with a norm metric
$|\cdot|_{p}:T_{p}M\rightarrow\mathbb{R}$ defined by
$|v|_{p}=\sqrt{g_{p}(v,v)}$ for any vector $v\in T_{p}M$. These properties
allow us to define the _length_ of a curve on the manifold surface. A
differentiable curve $\gamma:[a,b]\rightarrow M$ assigns to each $t\in(a,b)$ a
tangent vector $\gamma^{\prime}(t)\in T_{\gamma(t)}M$; thus, to obtain the
length of the curve $L(\gamma)$, we integrate the norm of the tangent vectors
along the curve (i.e., arc length):

$L(\gamma):=\int_{a}^{b}|\gamma^{\prime}(t)|_{\gamma(t)}\ dt$ (2.1)

In our analysis, we are primarily interested in measuring the _shortest_ curve
between a couple points on the manifold (known as a _geodesic_). Given a
couple of points on a Riemannian manifold $p,q\in M$ and the set of all curves
$\gamma:[a,b]\rightarrow M$ such that $\gamma(a)=p$ and $\gamma(b)=q$, the
geodesic $\bar{\gamma}$ is the curve with the shortest total length
$L(\bar{\gamma})$:

$L(\bar{\gamma}):=\text{inf}\\{L(\gamma)\ |\ \gamma:[a,b]\rightarrow
M,\text{with}\ \gamma(a)=p,\gamma(b)=q\\}$ (2.2)

An illustration of the geodesic between a couple of points on the smooth
sphere is shown in Figure 2. Geodesics are a powerful tool for quantifying the
relationship between points on a manifold surface and can be used to compute
summarizing statistics for points on the surface (such as means and variances)
[23].

(a)

(b)

Figure 2: (a) Illustration of a Riemannian manifold ($M$) and the associated
tangent space at a point $p\in M$. (b) Illustration of a geodesic $\gamma(t)$
constructed between two points $p,q\in M$, along with the associated tangent
space vector $v\in T_{p}M$. The exponential map
($\text{exp}_{p}(v):v\rightarrow q$) and the logarithmic map
($\text{log}_{p}(q):q\rightarrow v$) are also shown.

There are direct relationships between the tangent space of a Riemannian
manifold and geodesics, such as the _exponential map_ and the _logarithmic
map_. For a tangent vector $v\in T_{p}M$ constructed at point $p\in M$, there
exists a unique geodesic $\gamma:[0,1]\rightarrow M$ such that $\gamma(0)=p$
and $\gamma^{\prime}(0)=v$. The vector $v\in T_{p}M$ is mapped to the endpoint
of the geodesic $\gamma(1)\in M$ through the exponential map (see Figure 2):

$\displaystyle\text{exp}_{p}(v)=\gamma(1).$ (2.3)

The inverse of the exponential map is the logarithmic map, which maps the
point $\gamma(1)$ in the neighborhood of $p\in M$ to a vector in the tangent
space $v\in T_{p}M$ (see Figure 2):

$\displaystyle\text{log}_{p}(\gamma(1))=v$ (2.4)

These functions provide a mapping from the surface of the manifold to the
tangent space of a given point. The tangent space, which is a (linear) vector
space, can be analyzed using standard techniques designed for Euclidean spaces
(e.g., classification, regression, and dimensionality reduction) while
properly capturing the relationships between points defined by the geometry of
the manifold (e.g., geodesics). For instance, we can use these constructs to
analyze and operate on the space of SPD matrices, as we discuss next.

##  3 Riemannian Geometry of SPD Matrices

This section will provide the mathematical background needed to understand the
geometry of the space defined by SPD matrices. We first introduce the needed
notation and define the properties of different matrix spaces that are used in
our analysis. We then demonstrate how the properties of these matrix spaces
can be used to quantify the geometry and structure of the space defined by SPD
matrices. The main message is that SPD matrices lie on a Riemannian manifold
and that important computations (e.g., matrix operations, summarizing
statistics, classification, regression, and dimensionality reduction) can be
performed by respecting the geometry of this manifold by conducting these on
the tangent space; respecting such properties can lead to important
improvements in efficiency and interpretability.

###  3.1 Matrix Spaces, Properties, and Notation

We define spaces and properties of matrices that reflect the structure of SPD
matrices. We denote
$\mathcal{S}(n):=\\{\mathbf{S}\in\mathcal{M}(n),\mathbf{S}=\mathbf{S}^{T}\\}$
as the set of symmetric $n\times n$ matrices in the space of square, real
matrices $\mathcal{M}(n)$ and the set
$\mathcal{P}(n):=\\{\mathbf{P}\in\mathcal{S}(n),u^{T}\mathbf{P}u>0,\forall\
u\in\mathbb{R}^{n}\\}$ as the set of all $n\times n$ SPD matrices.

We also define the exponential and logarithmic mappings used in our analysis.
The matrix exponential $\text{exp}(\mathbf{P})$, where
$\mathbf{P}\in\mathcal{S}(n)$, is defined as:

$\text{exp}(\mathbf{P}):=\mathbf{U}\
\text{diag}(\text{exp}(\lambda_{1}),...,\text{exp}(\lambda_{n}))\
\mathbf{U}^{T}$ (3.5)

where $\mathbf{U}$ represents the matrix of eigenvectors of $\mathbf{P}$ and
$\lambda_{1}>...>\lambda_{n}$ represent the eigenvalues of $\mathbf{P}$ (also
denoted as $\lambda_{i}(\mathbf{P})$). We also define the inverse operation;
the matrix logarithm $\text{log}(\mathbf{P})$ as:

$\text{log}(\mathbf{P}):=\mathbf{U}\
\text{diag}(\text{log}(\lambda_{1}),...,\text{log}(\lambda_{n}))\
\mathbf{U}^{T}.$ (3.6)

The following properties should also be considered in the analysis [21]:

  * •

$\forall\ \mathbf{P}\in\mathcal{P}(n)$ we have that $\text{det}(\mathbf{P})>0$

  * •

$\forall\ \mathbf{P}\in\mathcal{P}(n)$ we have that
$\mathbf{P}^{-1}\in\mathcal{P}(n)$

  * •

$\forall\ \mathbf{P}\in\mathcal{P}(n)$ we have that
$\text{log}(\mathbf{P})\in\mathcal{S}(n)$

  * •

$\forall\ \mathbf{S}\in\mathcal{S}(n)$ we have that
$\text{exp}(\mathbf{S})\in\mathcal{P}(n)$

We also define the Frobenius inner product for matrices $\mathbf{A}$ and
$\mathbf{B}$ as:

$\langle\mathbf{A},\mathbf{B}\rangle_{F}:=\text{Tr}(\mathbf{A}^{T}\mathbf{B}),$
(3.7)

where Tr($\cdot$) represents the matrix trace operator. The Frobenius norm for
a matrix $\mathbf{A}$ is:

$\displaystyle||\mathbf{A}||_{F}=\sqrt{\text{Tr}(\mathbf{A}^{T}\mathbf{A})}=\sqrt{\sum_{i=1}^{n}\lambda_{i}(\mathbf{A})}.$
(3.8)

###  3.2 Manifold of SPD Matrices

To understand the Riemannian geometry of the space of SPD Matrices
($\mathcal{P}(n)$), we first need to construct a Riemannian metric, which will
allow us to computed distances and other relationships between points in a
manifold. The metric we consider for Riemannian manifolds is known as the
Affine Invariant Riemannian Metric (AIRM); a detailed derivation of the metric
is found in the work of Bhatia [5]. We aim to provide an intuitive
understanding of this metric, and its meaning.

Given SPD matrices $\mathbf{A}\in\mathcal{P}(n)$ and
$\mathbf{B}\in\mathcal{P}(n)$, we construct a geodesic
$\gamma(t):[0,1]\rightarrow\mathcal{P}(n)$ parameterized as:

$\displaystyle\gamma(t)=\text{exp}\big{(}(1-t)\cdot\text{log}(\mathbf{A})+t\cdot\text{log}(\mathbf{B})\big{)}$
(3.9)

where we have log maps $\mathcal{P}(n)\rightarrow\mathcal{S}(n)$ and exp maps
$\mathcal{S}(n)\rightarrow\mathcal{P}(n)$. Here, we are using the $n\times n$
identity matrix $\mathbf{I}\in\mathcal{P}(n)$ as our tangent space basis.
Informally, we are mapping our matrices from our SPD manifold $\mathcal{P}(n)$
to the tangent vector space $\mathcal{S}(n)$ via the logarithmic map,
constructing a line between these points in the tangent space, and then
projecting the constructed line back to the manifold via the exponential map.
This is guaranteed to be the shortest length path between the points in
$\mathcal{P}(n)$. The corresponding geodesic distance between matrices
$\mathbf{A}$ and $\mathbf{B}$ is given by:

$\displaystyle
d_{g}(\mathbf{A},\mathbf{B}):=||\text{log}(\mathbf{A})-\text{log}(\mathbf{B})||_{F}.$
(3.10)

Here, note that we are simply projecting the matrices from the Riemannian
manifold to the tangent vector space prior to measuring their distance using
the Frobenius norm [5]. However, an important consideration must be made when
using the $n\times n$ identity matrix $\mathbf{I}$ as the tangent space basis.
In many applications, the data may lie within a particular neighborhood that
is far from $\mathbf{I}$ on the manifold. Thus, projections to the tangent
space $T_{I}M$ can result in distortions of the data [23]. Intuitively, one
can think of this as similar to an analysis of a projection of the Earth
surface onto a plane tangent to the 0’ latitude and 0’ longitude point (as
many maps are represented). In this projection, landmasses near the edge of
the projection are highly distorted, whereas points near (0’,0’) have almost
no distortion. Thus, our aim is to identify a distance with respect to a
tangent space defined by the matrices of concern. This distance can be
constructed using a critical property of this Riemannian manifold and metric:
_congruence invariance_.

Congruence invariance states that, for any $n\times n$ invertible matrix
$\mathbf{X}$ and matrices $\mathbf{A},\mathbf{B}\in\mathcal{P}(n)$:

$\displaystyle
d_{g}(\mathbf{X}^{T}\mathbf{A}\mathbf{X},\mathbf{X}^{T}\mathbf{B}\mathbf{X})=d_{g}(\mathbf{A},\mathbf{B})$
(3.11)

We thus have that linear transformations of the given matrices do not impact
the geodesic distance on the manifold. This property allows us to redefine the
geodesic distance as:

$\displaystyle d_{g}(\mathbf{A},\mathbf{B})$
$\displaystyle=d_{g}(\mathbf{I},\mathbf{A}^{-1/2}\mathbf{B}\mathbf{A}^{-1/2})$
(3.12a)
$\displaystyle=||\text{log}(\mathbf{I})-\text{log}(\mathbf{A}^{-1/2}\mathbf{B}\mathbf{A}^{-1/2})||_{F}$
(3.12b)
$\displaystyle=||\text{log}(\mathbf{A}^{-1/2}\mathbf{B}\mathbf{A}^{-1/2})||_{F}$
(3.12c)
$\displaystyle=\bigg{(}\sum_{i=1}^{n}\text{log}^{2}\lambda_{i}(\mathbf{A}^{-1}\mathbf{B})\bigg{)}^{1/2}$
(3.12d)

where $\mathbf{I}$ is the $n\times n$ identity matrix,
$\mathbf{A}=\mathbf{A}^{1/2}\mathbf{A}^{1/2}$, and
$\lambda_{i}(\mathbf{A}^{-1}\mathbf{B})$ are the eigenvalues of
$\mathbf{A}^{-1}\mathbf{B}$.

We can simplify the geodesic distance further; we have that the geodesic
between matrices $\mathbf{I}$ and
$\mathbf{A}^{-1/2}\mathbf{B}\mathbf{A}^{-1/2}$ is:

$\displaystyle\gamma_{0}(t)=\text{exp}\big{(}\text{log}(\mathbf{A}^{-1/2}\mathbf{B}\mathbf{A}^{-1/2})t\big{)}=(\mathbf{A}^{-1/2}\mathbf{B}\mathbf{A}^{-1/2})^{t}$
(3.13)

We can then leverage congruence invariance to shift this geodesic to our
matrices of interest as:

$\displaystyle\gamma(t)=\mathbf{A}^{1/2}\big{(}\gamma_{0}(t)\big{)}\mathbf{A}^{1/2}=\mathbf{A}^{1/2}(\mathbf{A}^{-1/2}\mathbf{B}\mathbf{A}^{-1/2})^{t}\mathbf{A}^{1/2}$
(3.14)

where $\gamma(0)=\mathbf{A}$ and $\gamma(1)=\mathbf{B}$, which provides a
geodesic that is independent of the tangent space basis. Essentially, we are
leveraging congruence invariance to translate our points to a neighborhood of
$\mathbf{I}$, which allows us to compute distances with minimal distortion,
and then translate the points back through these linear transformations. We
can apply this same logic to the exponential map (and the logarithmic map);
for matrices $\mathbf{A},\mathbf{B}\in\mathcal{P}(n)$ and $\mathbf{T}_{B}\in
T_{A}M$ where $T_{A}M\subset\mathcal{S}(n)$:

$\displaystyle\mathbf{B}=\text{exp}_{A}(\mathbf{T_{B}})=\mathbf{A}^{1/2}\text{exp}(\mathbf{A}^{-1/2}\mathbf{T_{B}}\mathbf{A}^{-1/2})\mathbf{A}^{1/2}$
(3.15)

$\displaystyle\mathbf{T_{B}}=\text{log}_{A}(\mathbf{B})=\mathbf{A}^{1/2}\text{log}(\mathbf{A}^{-1/2}\mathbf{B}\mathbf{A}^{-1/2})\mathbf{A}^{1/2}$
(3.16)

Here, we are mapping the matrix $\mathbf{B}\in\mathcal{P}(n)$ to the tangent
space vector $\mathbf{T}_{B}\in T_{A}M$ which is centered at
$\mathbf{A}\in\mathcal{P}(n)$ via $\text{log}_{A}(\mathbf{B})$ and inversely
through $\text{exp}_{A}(\mathbf{B})$. Thus, with these newly defined geodesics
and mappings, we are able to compute relationships between SPD matrices with
minimal distortions in the tangent space. However, when given a large dataset
with multiple SPD matrices, the choice of a tangent space basis may not be
immediately clear. In this case, the geometric mean of the matrices on the
manifold is typically identified and used as reference point.

###  3.3 SPD Matrix Means and Tangent Spaces

In the analysis of a set of SPD matrices, we often need to identify a center
point on the SPD manifold that will minimize the distortion of all geometric
relationships between the matrices of the dataset when mapped to the tangent
space. This matrix is the (Riemannian) geometric mean of the matrices [44].
For a set of SPD matrices $\mathbf{A}_{i}$, the geometric mean (see Figure 3)
is the matrix $\bar{\mathbf{A}}$ that minimizes the sum of squared geodesic
distances to all other matrices in the set:

$\bar{\mathbf{A}}:=\mathop{\text{argmin}}_{\mathbf{M}}\
\sum_{i=1}^{n}||\text{log}(\mathbf{M}^{-1/2}\mathbf{A}_{i}\mathbf{M}^{-1/2})||_{F}$
(3.17)

We can see that the geometric mean is obtained by solving a matrix
optimization problem. For the SPD manifold this problem is geodesically convex
(similar to Euclidean convexity) [1]. This optimization problem can be solved
by using classical optimization algorithms that have been adapted to geometric
setting (e.g., gradient descent) [1]. A detailed review on these approaches
can be found in the work by Absil, Mahony, and Sepulchre [1].

(a) Geometric mean.

(b) Tangent space projection.

Figure 3: (a) Illustration of the geometric mean $\bar{\mathbf{A}}$ of a set
of matrices
$\\{\mathbf{A}_{1},\mathbf{A}_{2},\mathbf{A}_{3},\mathbf{A}_{4}\\}\in\mathcal{P}(n)$.
The geometric mean represents a point on the manifold that minimizes the
geodesic to all other matrices in the set. (b) Representation of the tangent
space $T_{\bar{A}}M$ at the geometric mean. The set of matrices
$\mathbf{A}_{i}$ are projected (through the logarithmic map) onto this tangent
space with minimal geometric distortion $\mathbf{T_{A_{i}}}\in T_{\bar{A}}M$.

Given a set of SPD matrices $\mathbf{A}_{i}$ and a geometric mean
$\bar{\mathbf{A}}$, we can construct a tangent space at the geometric mean
$T_{\bar{A}}M$, and project the SPD matrices onto the tangent space
$\mathbf{T_{A_{i}}}=\text{log}_{\bar{A}}(\mathbf{A}_{i})$, as shown in Figure
3. The matrices are now represented in a vector (linear) space that reflects
the geometry of the SPD manifold. Projecting the data into a vector space
allows us to apply common matrix analysis methods such as PCA.

##  4 Case Study : Process Monitoring

Figure 4: Simplified illustration of the Tennessee Eastman Process (TEP) and
resulting multivariate time series sensor data. Process sensors measure values
such as temperature, pressure, flow, and level. From the multivariate time
series data it is difficult to distinguish whether there is a fault occurring,
or what type of fault may be occurring. Thus, a simplified and informative
representation of the data is needed to be able to distinguish when the
process is behaving normally or is experiencing a particular fault.

We focus on data obtained from a simulated industrial process known as the
Tennessee Eastman Process (TEP) [9]. This dataset is a widely used benchmark
dataset for testing and comparing various anomaly (i.e. fault) detection
methods [17, 43, 8]. Figure 4 provides a high-level illustration of the
process along with the multivariate time series data that is produced by the
sensors monitoring the process. The process has a total of 52 measurements, 41
are process variables, 11 are manipulated variables. There are 20 different
potential faults, which are defined in Table 1 (for further details see
Appendix A in [9]). Our aim is to use geometric methods to detect and classify
the presence and type of fault using only the multivariate sensor data.

[!h]  Fault ID Fault Name Type Fault 1 A/C feed ratio, B composition constant
(stream 4) Step Fault 2 B composition, A/C ratio constant (stream 4) Step
Fault 3 D feed temperature (stream 2) Step Fault 4 Reactor cooling water inlet
temperature Step Fault 5 Condenser cooling water inlet temperature Step Fault
6 A feed loss (stream 1) Step Fault 7 C header pressure loss - reduced
availability (stream 4) Step Fault 8 A, B, C feed composition (stream 4)
Random variation Fault 9 D feed temperature (stream 2) Random variation Fault
10 C feed temperature (stream 4) Random variation Fault 11 Reactor cooling
water inlet temperature Random variation Fault 12 Condenser cooling water
inlet temperature Random variation Fault 13 Reaction kinetics Slow drift Fault
14 Reactor cooling water valve Sticking Fault 15 Condenser cooling water valve
Sticking Fault 16-20 Unknown Unknown

Table 1: Types of Faults for Tennessee Eastman Process [9].

The difficulty in distinguishing potential faults in the TEP is illustrated in
Figure 4, where multivariate time series sensor signals for the 52 monitored
variables are shown for the TEP without faults, and in the presence of fault 7
(step change in the component ”C” header pressure) and fault 12 (random
variation in the condenser cooling water inlet). The complexity of the
multivariate process dynamics and of the number of sensor measurements make it
difficult to reliably identify if a fault is occurring and to distinguish
between fault types. Our method focuses on simplifying the data by quantifying
the relationships between the 52 measured variables through covariance
matrices and leveraging the geometry of if covariance matrices to detect and
distinguish faults in the TEP.

###  4.1 Data Pre-Processing

(a) Covariance matrices.

(b) Geometric mean.

(c) Tangent space projection.

Figure 5: Representation of the data pre-processing workflow for the geometric
analysis of the TEP data. (a) Covariance matrices are constructed from the
sensor data multivariate time series forming a set of SPD matrices
($\mathbf{P}_{i}\in\mathcal{P}(n)$). (b) The geometric mean of the matrices
($\bar{\mathbf{P}}\in\mathcal{P}(n)$) is identified as the point that
minimizes the squared geodesic distance to all other points. (c) All derived
covariance matrices are mapped to the tangent space $T_{\bar{P}}M$ through the
logarithmic mapping. This maps the matrices into a vector space that reflects
the manifold geometry. The mapped data can then be analyzed with commonly
defined dimensionality reduction and classification/regression methods.

To begin our analysis of the TEP, we must first pre-process the TEP data into
multiple covariance matrices. To accomplish this, we represent each of the 52
measured variables as a univariate random variable $x_{i}$ where
$i=1,2,...,52$ and we denote the collection of signals as a multivariate
random vector $\mathbf{X}=(x_{1},...,x_{n})$ where $n=52$. We denote the
observations of each signal at time $t=1,2,..,m$ as
$x_{i}(t)\in\mathbb{R}^{m}$. We use this representation to construct the
sample covariance matrix for the process data $\mathbf{P}\in\mathcal{P}(n)$
as:

$\displaystyle\mathbf{P}:=\frac{1}{m-1}\mathbf{X}\mathbf{X}^{T}$ (4.18)

The TEP dataset consists of multiple separate simulations of the process, both
with and without faults. Thus, for each simulation we construct a sample
covariance matrix $\mathbf{P}_{i}$. Our goal is to pair each simulation sample
covariance matrix with the fault occurring in the simulation. Figure 5
provides examples of the sample covariance matrices that are constructed from
simulations containing no faults, fault 7, and fault 12. We note that there
are no obvious differences between the covariance matrices that would identify
a given fault. These covariance matrices lie on the SPD manifold, and can be
integrated into our geometric framework. An illustration of this computational
workflow is found in Figure 5. We compute the geometric mean
$\bar{\mathbf{P}}$ for our set of covariance matrices
$\mathbf{P}_{i}\in\mathcal{P}(n)$, the matrices are then projected to the
tangent space $T_{\bar{P}}M$ centered at the geometric mean via the
logarithmic mapping $\mathbf{T_{P_{i}}}=\text{log}_{\bar{P}}(\mathbf{P}_{i})$.
The data is now projected into a vector space that retains the geometric
characteristics of the SPD matrices with minimal distortion, and can be
integrated in dimensionality reduction and classification algorithms to
perform analysis.

###  4.2 Principal Geodesic Analysis

Mapping the process data covariance matrices to the tangent space provides an
avenue for the application of common dimensionality reduction techniques.
Here, we apply PCA to the matrices mapped to the tangent (vector) space. PCA
applied on the tangent space of the SPD manifold is commonly known as
Principal Geodesic Analysis (PGA), as it identifies the geodesics that capture
the most variance in the data [15]. An example comparison of PGA versus PCA
(directly on the covariance matrices) is presented in Figure 6. The
simulations with no faults are colored in red and the faulty simulations are
represented by different grayscale values. We can see that using only the
first two components in PGA, we are able to perfectly separate the faulty and
non-faulty simulations; on the other hand, when applying PCA directly on the
covariance matrices we can see that there is significant overlap between the
faulty and non-faulty simulations. The comparison of these projections
demonstrates that capturing the geometry of the Riemannian manifold in the
analysis of covariance matrices can improve the performance with minimal added
complexity.

This improvement in separation of the data through the geometric approach is
due in part to the congruence invariance of our defined metric on the
Riemannian manifold. As previously stated, congruence invariance means that
any $n\times n$ invertable matrix $\mathbf{X}$ applied to a set of covariance
matrices $\mathbf{P}_{i}\in\mathcal{P}(n)$ does not impact the geodesic
distance between the two matrices:

$\displaystyle
d_{g}(\mathbf{X}^{T}\mathbf{P}_{i}\mathbf{X},\mathbf{X}^{T}\mathbf{P}_{j}\mathbf{X})=d_{g}(\mathbf{P}_{i},\mathbf{P}_{j})$
(4.19)

Operations such as re-scaling and normalization, which can be represented
algebraically as invertable square matrices, have no impact on the geodesic
distance between the matrices [3]. Therefore, there is no need to select
specific scaling or normalization strategies when applying our geometry based
analysis of the data (e.g., PGA). However, this is not true when ignoring the
data geometry, making methods such as PCA susceptible to the chosen framework
(or lack of) for normalization/scaling. We perform no scaling prior to PCA in
this case to ensure a direct comparison between PCA and PGA.

(a) Principal component analysis.

(b) Principal geodesic analysis.

Figure 6: Comparison of PCA applied to the raw covariance matrices and PCA
applied to data mapped onto the tangent space (PGA). The red points represent
the simulations where no fault is occurring and the grayscale points represent
simulations with different faults. (a) PCA on the raw covariance matrices
shows minimal separation in the data; faultless simulations are overlapped
with faulty simulations. (b) PCA performed in the tangent space provides
perfect separation between the faulty and faultless simulations, and also
shows clustering of the faulty systems into separate groups. This demonstrates
that simple considerations for the geometry of the data can yield improved
results.

###  4.3 Classification and Clustering Results

PGA analysis of the covariance matrices also reveals that there is definite
clustering of the data with respect to the different fault types within the
TEP dataset. This suggests that classification of the different fault types
can be done directly using the tangent space of the SPD manifold. To
investigate this we use a simple linear (ridge) classifier. We compare the
prediction accuracy of the linear classifier using coefficients of the tangent
space projected matrices versus the coefficients of the non-transformed
covariance matrices as input. In the analysis, we perform a simple train-test
split of the data, where $30\%$ of the data is used for testing and $70\%$ of
the data is use for training. Figure 7 illustrates the dramatic increase in
accuracy when the model incorporates geometric information, which is reflected
in the normalized confusion matrices. Here, a value $x\in[0,1]$ on the
diagonal indicates an accuracy of $x*100\%$ when classifying a particular
fault. All values in the off diagonal (e.g., row $i$, column $j$, where $i\neq
j$) represent the percentage of covariance matrices associated with fault $i$
that have been incorrectly labeled as experiencing fault $j$. When the SPD
manifold is accounted for via the tangent space projection, there is perfect
classification of the data (with the exception of faults $3,4,9$ and $15$).
When the manifold geometry is ignored, there are few instances where high
classification accuracy is achieved. The faults $3,4,9$ and $15$ have been
shown in prior work to be difficult to classify [43]. We also note that these
faults are only misclassified within their group (are never classified as
having no fault), which suggests that there is limited quantifiable difference
in the covariance matrices for these faults. The inclusion of more information
around these particular faults may correct this issue and further increase
accuracy.

(a) Classification raw covariance matrices.

(b) Classification in tangent space ($T_{\bar{P}}M$).

Figure 7: Comparison of linear classification on the raw covariance matrices
versus matrices mapped to the tangent space $T_{\bar{P}}M$. (a) Classification
of the covariance matrices without regard for the data geometry results in
poor classification accuracy. (b) Simple mapping of the data to the
appropriate tangent space provides a dramatic improvement in classification
accuracy.

##  5 Case Study - Image Analysis

Another important application of covariance matrices is in the analysis of
images [26, 10, 6]. Here, we focus on an analysis of real images from textile
manufacturing that contain non-defective and defective woven textiles taken.
The images were obtained from the public MVTEC AD dataset [4]. Example images
of non-defective and defective textiles are found in Figure 8. Our goal is to
classify textiles using a linear classifier and Riemannian geometry.

(a) Normal textiles.

(b) Defective textiles.

Figure 8: Example images from the woven textile dataset [4]. (a)
Representative sample of woven textiles with no defects. (b) Representative
sample of woven textiles that are considered defective.

###  5.1 Data Pre-Processing

Covariance matrices are useful data representations for images because they
are invariant to translations and rotations [36]. These matrices can also be
used to combine multiple image features that can be quantified through filters
and kernel methods [37]. Here, we select nine image features, as shown in
Figure 9. We implement both the Frangi and Hessian filters, which are designed
to detect edges and fiber structures of the woven textures [13]. We apply
these filters to the image; multiple transformed images that are smoothed
through a Gaussian filter of varying strength are shown in Figure 9. We use
Gaussian filters to emphasize features of different scale within the image
[45]. This yields nine total feature images (including the original image),
which we use to construct a sample covariance matrix. We do this by taking
each feature image (which we treat as a $64\times 64$ matrix) and vectorizing
each image $x_{i}\in\mathbb{R}^{m}$, where $i=1,2,...,9$ and $m=4096$. Each
feature image can therefore be represented as a realization of a multivariate
random vector $\mathbf{X}=(x_{1},x_{2},...,x_{9})$. We use this observation to
construct a sample covariance matrix $\mathbf{P}\in\mathcal{P}(n)$ (see Figure
9) where $n=9$:

$\mathbf{P}:=\frac{1}{m-1}\mathbf{X}\mathbf{X}^{T}.$ (5.20)

The covariance matrices are SPD, and can be directly integrated into our
geometric analysis framework.

Figure 9: Workflow for the construction of an image covariance matrix. (left)
Image filters and transformations emphasize specific characteristics of an
image. Gaussian filters emphasize features of different scale within an image,
and the Frangi and Hessian filters capture important fiber and edge features
of an image. (right) The covariance between each image representation can be
computed and used to form a covariance matrix. Importantly, the covariance
matrix representation is invariant to transformations such as rotation and
translation which are present in the textile images.

###  5.2 Classification Results

We first apply PGA to the image covariance matrices (see Figure 10). PGA
reveals a distinct clustering and separation of the data with only two
components. The PGA suggests that a simple classifier can be used to separate
the defective and non-defective textile samples (after projecting to the
tangent space). Thus, we construct a simple linear support vector machine
classifier, which we train with $70\%$ of the image data and test our trained
model on the remaining $30\%$ of the data. The covariance representation of
the images provides a simple characterization of the data and its various
features, while also imbuing the data with the inherent geometry of SPD
matrices. This results in the trained model being able to separate the
defective and non-defective samples in the testing data with $92\%$ accuracy.
We also compare our PGA analysis of the covariance matrix representation to
PCA applied to the raw image data in Figure 10. The PCA analysis reveals
almost no separation of the data, and results in poor classification accuracy
when used as input to a linear SVM (50 % accuracy). This is likely due to
sensitivity of the analysis on the raw images to rotations and translations of
the textiles.

(a) PGA on covariance matrices.

(b) PCA on raw images.

Figure 10: (a) PGA of the covariance matrices derived from the set of woven
textile images. The analysis reveals clustering of the data into groups
representing the defective and non-defective textiles. A simple linear SVM
classifier is able to separate the defective and non-defective textile images
with $92\%$ accuracy. (b) PCA analysis of the raw images; without the
covariance matrix representation and considerations for the data geometry, we
see there is almost no separation in the data.

##  6 Conclusions and Future work

We have presented an introduction to concepts of Riemannian geometry for
symmetric positive definite (SPD) matrices and show how these concepts can be
used in applications of interest to chemical engineers. Specifically, we
discussed approaches to capture the geometry of the SPD manifold (a Riemannian
manifold) in dimensionality reduction and classification tasks. Through a
couple of case studies, we demonstrated that capturing such geometry can lead
to significant improvements in accuracy and interpretability. This work
focused primarily on the manifold of SPD matrices, but we note that similar
concepts can be applied to over matrix representations such as fixed-rank
matrices, orthogonal matrices, and shape spaces as well as to other datasets
where a Riemannian manifold is known to exist (e.g., points on a sphere).
Furthermore, geometry can be leveraged in optimization and control to relax
certain constraints that are already enforced by the underlying system
geometry (e.g., optimization on the SPD manifold ensures that all solutions
are SPD).

## Acknowledgments

We acknowledge funding from the U.S. National Science Foundation (NSF) under
BIGDATA grant IIS-1837812.

## References

  * [1] P-A Absil, Robert Mahony, and Rodolphe Sepulchre.  Optimization algorithms on matrix manifolds.  Princeton University Press, 2009.
  * [2] Vincent Arsigny, Pierre Fillard, Xavier Pennec, and Nicholas Ayache.  Geometric means in a novel vector space structure on symmetric positive-definite matrices.  SIAM journal on matrix analysis and applications, 29(1):328–347, 2007.
  * [3] Alexandre Barachant and Marco Congedo.  A plug&play p300 bci using information geometry.  arXiv preprint arXiv:1409.0107, 2014.
  * [4] Paul Bergmann, Michael Fauser, David Sattlegger, and Carsten Steger.  Mvtec ad–a comprehensive real-world dataset for unsupervised anomaly detection.  In Proceedings of the IEEE/CVF Conference on Computer Vision and Pattern Recognition, pages 9592–9600, 2019.
  * [5] Rajendra Bhatia.  Positive definite matrices.  Princeton university press, 2009.
  * [6] Charlene E Caefer, Jerry Silverman, Oded Orthal, Dani Antonelli, Yaron Sharoni, and Stanley R Rotman.  Improved covariance matrices for point target detection in hyperspectral data.  Optical Engineering, 47(7):076402, 2008.
  * [7] Isaac Chavel.  Riemannian geometry: a modern introduction, volume 98.  Cambridge university press, 2006.
  * [8] Leo H Chiang, Evan L Russell, and Richard D Braatz.  Fault diagnosis in chemical processes using fisher discriminant analysis, discriminant partial least squares, and principal component analysis.  Chemometrics and intelligent laboratory systems, 50(2):243–252, 2000.
  * [9] James J Downs and Ernest F Vogel.  A plant-wide industrial process control problem.  Computers & chemical engineering, 17(3):245–255, 1993.
  * [10] Ian L Dryden, Alexey Koloydenko, and Diwei Zhou.  Non-Euclidean statistics for covariance matrices, with applications to diffusion tensor imaging.  The Annals of Applied Statistics, 3(3):1102–1123, 2009.
  * [11] Bo Feng, Mengyin Fu, Hongbin Ma, Yuanqing Xia, and Bo Wang.  Kalman filter with recursive covariance estimation—sequentially estimating process noise covariance.  IEEE Transactions on Industrial Electronics, 61(11):6253–6263, 2014\.
  * [12] P Thomas Fletcher, Conglin Lu, Stephen M Pizer, and Sarang Joshi.  Principal geodesic analysis for the study of nonlinear statistics of shape.  IEEE transactions on medical imaging, 23(8):995–1005, 2004.
  * [13] Alejandro F Frangi, Wiro J Niessen, Koen L Vincken, and Max A Viergever.  Multiscale vessel enhancement filtering.  In International conference on medical image computing and computer-assisted intervention, pages 130–137. Springer, 1998.
  * [14] Joaquín Goñi, Martijn P Van Den Heuvel, Andrea Avena-Koenigsberger, Nieves Velez De Mendizabal, Richard F Betzel, Alessandra Griffa, Patric Hagmann, Bernat Corominas-Murtra, Jean-Philippe Thiran, and Olaf Sporns.  Resting-brain functional connectivity predicted by analytic measures of network communication.  Proceedings of the National Academy of Sciences, 111(2):833–838, 2014.
  * [15] Mehrtash T Harandi, Mathieu Salzmann, and Richard Hartley.  From manifold to manifold: Geometry-aware dimensionality reduction for spd matrices.  In European conference on computer vision, pages 17–32. Springer, 2014.
  * [16] M Willjuice Iruthayarajan and S Baskar.  Covariance matrix adaptation evolution strategy based design of centralized pid controller.  Expert systems with Applications, 37(8):5775–5781, 2010.
  * [17] Wenfu Ku, Robert H Storer, and Christos Georgakis.  Disturbance detection and isolation by dynamic principal component analysis.  Chemometrics and intelligent laboratory systems, 30(1):179–196, 1995.
  * [18] Serge Lang.  Fundamentals of differential geometry, volume 191.  Springer Science & Business Media, 2012.
  * [19] John M Lee.  Riemannian manifolds: an introduction to curvature, volume 176.  Springer Science & Business Media, 2006.
  * [20] M Mansouri, MZ Sheriff, R Baklouti, M Nounou, H Nounou, A Ben Hamida, and N Karim.  Statistical fault detection of chemical process-comparative studies.  Journal of Chemical Engineering & Process Technology, 7(1):282–291, 2016.
  * [21] Maher Moakher and Philipp G Batchelor.  Symmetric positive-definite matrices: From geometry to applications and visualization.  In Visualization and Processing of Tensor Fields, pages 285–298. Springer, 2006.
  * [22] Jorge Nocedal and Stephen Wright.  Numerical optimization.  Springer Science & Business Media, 2006.
  * [23] Xavier Pennec.  Intrinsic statistics on Riemannian manifolds: Basic tools for geometric measurements.  Journal of Mathematical Imaging and Vision, 25(1):127–154, 2006\.
  * [24] Xavier Pennec.  Manifold-valued image processing with spd matrices.  In Riemannian geometric statistics in medical image analysis, pages 75–134. Elsevier, 2020.
  * [25] Fatih Porikli.  Achieving real-time object detection and tracking under extreme conditions.  Journal of Real-Time Image Processing, 1(1):33–40, 2006.
  * [26] Fatih Porikli and Oncel Tuzel.  Fast construction of covariance matrices for arbitrary size image windows.  In 2006 International Conference on Image Processing, pages 1581–1584. IEEE, 2006.
  * [27] Nara M Portela, George DC Cavalcanti, and Tsang Ing Ren.  Semi-supervised clustering for mr brain image segmentation.  Expert Systems with Applications, 41(4):1492–1497, 2014.
  * [28] Anqi Qiu, Annie Lee, Mingzhen Tan, and Moo K Chung.  Manifold learning on brain functional networks in aging.  Medical image analysis, 20(1):52–60, 2015.
  * [29] Gregory C Reinsel.  Elements of multivariate time series analysis.  Springer Science & Business Media, 2003.
  * [30] Bernhard Riemann.  On the Hypotheses which lie at the Bases of Geometry.  Birkhäuser, 2016.
  * [31] Evan L Russell, Leo H Chiang, and Richard D Braatz.  Fault detection in industrial processes using canonical variate analysis and dynamic principal component analysis.  Chemometrics and intelligent laboratory systems, 51(1):81–93, 2000\.
  * [32] Alexander D Smith, Paweł Dłotko, and Victor M Zavala.  Topological data analysis: concepts, computation, and applications in chemical engineering.  Computers & Chemical Engineering, 146:107202, 2021.
  * [33] Alexander D Smith and Victor M Zavala.  The euler characteristic: A general topological descriptor for complex data.  Computers & Chemical Engineering, 154:107463, 2021.
  * [34] Stefan Sommer, François Lauze, Søren Hauberg, and Mads Nielsen.  Manifold valued statistics, exact principal geodesic analysis and the effect of linear approximations.  In European conference on computer vision, pages 43–56. Springer, 2010.
  * [35] Olaf Sporns.  Network analysis, complexity, and brain function.  Complexity, 8(1):56–60, 2002.
  * [36] Oncel Tuzel, Fatih Porikli, and Peter Meer.  Region covariance: A fast descriptor for detection and classification.  In European conference on computer vision, pages 589–600. Springer, 2006.
  * [37] Oncel Tuzel, Fatih Porikli, and Peter Meer.  Pedestrian detection via classification on Riemannian manifolds.  IEEE transactions on pattern analysis and machine intelligence, 30(10):1713–1727, 2008.
  * [38] Gaël Varoquaux, Flore Baronnet, Andreas Kleinschmidt, Pierre Fillard, and Bertrand Thirion.  Detection of brain functional-connectivity difference in post-stroke patients using group-level covariance modeling.  In International Conference on Medical Image Computing and Computer-Assisted Intervention, pages 200–208. Springer, 2010.
  * [39] Larry Wasserman.  All of statistics: a concise course in statistical inference, volume 26.  Springer, 2004.
  * [40] Barry M Wise and Neal B Gallagher.  The process chemometrics approach to process monitoring and fault detection.  Journal of Process Control, 6(6):329–348, 1996.
  * [41] Yang Xu, Zebin Wu, Jun Li, Antonio Plaza, and Zhihui Wei.  Anomaly detection in hyperspectral images based on low-rank and sparse representation.  IEEE Transactions on Geoscience and Remote Sensing, 54(4):1990–2000, 2015.
  * [42] Jianfeng Yao, Shurong Zheng, and ZD Bai.  Sample covariance matrices and high-dimensional data analysis.  Cambridge University Press Cambridge, 2015.
  * [43] Shen Yin, Steven X Ding, Adel Haghani, Haiyang Hao, and Ping Zhang.  A comparison study of basic data-driven fault diagnosis and process monitoring methods on the benchmark tennessee eastman process.  Journal of process control, 22(9):1567–1581, 2012.
  * [44] Kisung You and Hae-Jeong Park.  Re-visiting Riemannian geometry of symmetric positive definite matrices for the analysis of functional connectivity.  Neuroimage, 225:117464, 2021.
  * [45] Ian T Young and Lucas J Van Vliet.  Recursive implementation of the Gaussian filter.  Signal processing, 44(2):139–151, 1995.