WGMwgm QEqe EPep PMSpms BECbec DEde

[1]

[cor1]Corresponding<EMAIL_ADDRESS>
# ContrastNet: Unsupervised Feature Learning by Autoencoder and Prototypical
Contrastive Learning for Hyperspectral Imagery Classification

Zeyu Cao  Xiaorun Li  Zhejiang University, College of Electrical Engineering,
38 Zheda Rd, Hangzhou, China, 310027  Yueming Feng  State Grid Jiaxing Power
Supply Company, Jiaxing Zhejiang, China, 314100  Shuhan Chen  Chaoqun Xia
Liaoying Zhao  HangZhou Dianzi University, China Institute of Computer
Application Technology, Xiasha Higher Education Zone, Hangzhou, China, 310018

###### Abstract

Hyperspectral classification is a fundamental problem for applying
hyperspectral technology, and unsupervised learning is a promising direction
to address the issue. Unsupervised feature learning algorithm extracts
representative features, which help efficiently classify the hyperspectral
pixels. This paper combines the popular contrastive learning method
(prototypical contrastive learning) and the classic representation learning
method (autoencoder) to design an unsupervised feature learning network for
hyperspectral classification. First, two different encoders are used to
extract different features as augmentation function in a contrastive network
(ContrastNet). Second, the prototypical contrastive learning method is adopted
to train a contrastive network. Third, features extracted by the contrastive
network are used for classification. Experiments have proved that our two
proposed autoencoder networks show good feature learning capabilities by
themselves and the designed contrastive learning network can further learn
better representative features from the two modules. The three modules compose
an efficient framework for unsupervised feature learning. Our method also
performs reasonably fast in the testing phase, which implies that it is
applicable in practice blackunder specific conditions.

###### keywords:

Unsupervised learning, autoencoder, contrastive learning, hyperspectral
imagery classification

##  1 Introduction

With the rapid development of deep learning, much progress has been made in
the hyperspectral classification field. Feature extraction is an essential
part in the hyperspectral classification. And recent years, many great
supervised deep learning methods for hyperspectral classification were good at
feature extraction, such as one-dimensional convolutional neural networks
(1DCNN) [1D-CNN], three-dimensional convolutional neural networks (3DCNN)
[3DCNN], deep recurrent neural networks (RNN) [RNN], deep feature fusion
network (DFFN) [DFFN], spectral-spatial residual network (SSRN) [SSRN],
blackadaptive spectral-spatial multiscale contextual feature extraction
[adaptive], fully contextual network (FullyContNet) [contextual], 3D attention
denoising network [Attention]. redBesides, some machine learning methods are
also good at feature extraction, such as multitask learning-based reliability
analysis [multitask], joint sparse representation and multitask learning
[jointsparse]. Nowadays, it is widely admitted that the fusion of spatial
information and spectral information is the best strategy to extract features
for hyperspectral classification. Furthermore, deep learning models show great
power in finding the two kinds of information. With enough labeled samples,
many supervised learning methods can effectively learn the spatial and
spectral information simultaneously (e.g., SSRN) and update the models’
parameters, thus getting good hyperspectral classification results.

However, there are some disadvantages to supervised learning methods. First,
severe data dependence is the notorious disadvantage of supervised deep
learning. Without enough labeled samples, many supervised deep learning
methods work badly. Second, many supervised deep learning methods only rely on
the labels to supervise the extraction of features. As a result, with the
increase of labels, the models need to be retrained entirely, which seems
expensive in the practical circumstances.

Researchers have paid their attention to unsupervised learning methods to
avoid the disadvantages of supervised deep learning in the past years.
Unsupervised learning methods are mostly feature extraction algorithms, with a
suitable feature extractor, we can use a simple classifier to get good
classification results. To some extent, unsupervised learning methods are more
applicable in practical scenarios. For unsupervised learning, there are some
classic methods widely applied to all kinds of fields, such as principal
component analysis (PCA) [PCA] and independent component analysis (ICA) [ICA].
These traditional algorithms remain efficient in the pre-process or feature
extraction field.

Generally speaking, unsupervised learning aims to dig out information from
data other than labels; there are two major ways to achieve the objective:
representative learning and discriminative learning. There are two famous
representative learning models: autoencoder (AE) [AE] and generative
adversarial network (GAN) [GAN]. The two methods both use the self-similarity
of data and then construct a model that can map the original data to some
certain distributions. In this way, we can get generated distributions which
are similar to the distributions of real samples. As a result, we can get
forgery samples that are similar to real samples or get good feature
extractors. Inspired by autoencoder, stacked sparse autoencoder (SSAE) [SSAE]
used autoencoder for sparse spectral feature learning and multiscale spatial
feature learning, respectively and then fused these two kinds of feature for
hyperspectral classification. Furthermore, GAN [GANforC] has also been proved
to be useful in extracting spectral-spatial features for hyperspectral
classification. Except for these methods, three-dimensional convolutional
autoencoder (3DCAE) [3DCAE], semi-supervised GAN [semi-GAN], GAN-Assisted
CapsNet (TripleGAN) [triple-gan] and many other papers also originated from
the representative learning. Generally, many researchers have acknowledged the
utility of the representative learning methods in the remote sensing field.

On the other hand, discriminative learning is also an original idea in
exploring the information of data. blackSupport vector machine (SVM) [svm] may
be the most classic discriminative learning algorithm, and decision tree
[dtree], multilayer perceptron (MLP) [mlp], conditional random field (CRF)
[crf] are all methods based on discriminative learning. Among various
discriminative learning models, an impressive method is the contrastive
learning method. Unlike autoencoder and GAN, the contrastive learning method
does not care about the details of the data. It cares about the difference
between the different samples. In other words, contrastive learning aims to
discriminate different samples and does not aim to find samples’ exact
distributions. By comparing different samples, a discriminator is trained and
can be used for classification or other tasks. By finding the difference
between two randomly selected samples, the discriminator learns to extract
features that are helpful for downstream tasks.black Contrastive learning can
be applied in supervised [sucontra] as well as unsupervised learning
[uncontra].

For unsupervised contrastive learning, there are no label information can be
utilized at all. Therefore a paradigm is widely used for unsupervised
learning: we use two kinds of data augmentation methods to transform a single
sample into two transformed samples because contrastive learning requires the
contrast between two samples. In this way, we introduce a kind of information:
the two samples originate from the same sample. Then with a proper model
structure and a well-designed objective function, we can train an encoder that
can map the samples to latent encoding space, and the features from the
encoding space can well represent the original samples. When the encoder acts
as a feature extractor, many downstream tasks can be done. Much work has been
done in the computer vision field, such as momentum contrast [Moco] and a
simple framework for contrastive learning [simclr]. They proposed some
original structures for unsupervised contrastive learning and got promising
results in many visual datasets. Furthermore, prototypical contrastive
learning (PCL) [PCL] combined cluster algorithm with contrastive learning, and
got better performances than other methods, even outperformed some supervised
methods. blackBesides, contrastive learning is also adopted for self-
supervised learning from unlabeled video [timecontra] and for speaker
recognition [speakcontra]. Generally, unsupervised contrastive learning has
shown its great power in the computer vision field.

To some extent, representative learning method cares about the accurate
details of data, while discriminative learning method pays attention to the
information in a higher level. Naturally, we hope to use the advantages of
both methods in feature extraction. However, as far as we know, few
contrastive learning methods have been reported in the hyperspectral
classification field. There are major problems to be addressed to transfer the
unsupervised contrastive learning method to the process of hyperspectral
image. First, the computer vision field’s typical data augmentation methods
seem not reasonable in the hyperspectral classification field. For example,
color distortion to the normal images will not change objects’ spatial
features, so it is an acceptable transformation. When distorting hyperspectral
images’ spectral value, the spectral information is destroyed, and the
spectral information is essential for hyperspectral classification, so this is
not an acceptable transformation for hyperspectral images. Second,
unsupervised contrastive learning for normal images demands large batch size
and lots of GPUs, because the larger the batch size is, the more information
can be utilized from negative samples. However, if it is directly applied to
hyperspectral images, the requirement for more computing resources will be too
expensive.

PCL showed great power of unsupervised learning in processing normal images.
However, as mentioned above, it is hard to apply PCL to hyperspectral images
directly. So we must make some adjustments to PCL in order to apply it to
hyperspectral classification. With the blossom of AE-based methods in the
hyperspectral images field, we find that the encoder in an autoencoder
structure can be seen as a transformation, just like other data augmentation
methods. They both map the original images from the original distribution to
another distribution. The difference is, data augmentations usually do not
change the shape of images, while encoder usually encodes the images to
vectors with less dimensions. Another essential truth is that data
augmentation methods usually keep most information of the original images to
make them able to be correctly classified. Likewise, the encoder in a trained
autoencoder also can keep most information in the original images.

As a result, we find it reasonable to treat the encoder as a transformation
function that plays the role of data augmentation. Then the question is to
find two encoders with different distributions, just like two kinds of data
augmentations. So we choose designed adversarial autoencoder (AAE) and
variational autoencoder (VAE) modules to act as two kinds of data augmentation
because AAE and VAE can generate relatively fixed and slightly different
distributions. Then the PCL only needs to learn from two kinds of vectors
other than two entire images.

In this way, we can not only apply PCL to hyperspectral images but also lower
the demand for deep networks, because the network for PCL only needs to
process vectors with low dimensions. By encoding the images to shape fixed
vectors (e.g., 1024-d), we can assign the feature extraction and contrastive
learning into two parts, and each part will not be too complicated. Through
this method, we make contrastive learning methods affordable for those who
have limited computing resources.

To apply the unsupervised contrastive learning to hyperspectral
classification, we addressed the two problems and proposed ContrastNet, a deep
learning model for unsupervised feature learning of hyperspectral images. Our
innovations can be concluded as follows:

  1. 1)

We introduced discriminative learning methods into remote sensing field.
Specifically, we applied the prototypical contrastive learning to unsupervised
feature learning for hyperspectral classification.

  2. 2)

Enlightened by the autoencoder structure, we redintroduced two autoencoder-
based modules: adversarial autoencoder module(AAE module), variational
autoencoder module(VAE module). The two modules can work alone for feature
extraction. When using the extracted features for hyperspectral
classification, they outperformed many baselines, even supervised methods.

  3. 3)

We proposed an excellent pipeline to process the hyperspectral images, which
combines two kinds of unsupervised learning methods (representative learning
and contrastive learning). And this pipeline shows good performance in
unsupervised feature learning for hyperspectral classification. Furthermore,
the pipeline only demands a single GPU, which means it is more affordable than
the original prototypical contrastive learning method.

##  2 Related Work

###  2.1 Variational Autoencoder

AE structure [AE] is a widely used unsupervised feature extraction model. The
standard AE is composed of an encoder and a decoder. The encoder encodes the
input images to features with fixed shape. Moreover, the features are called
latent code. The decoder decodes the latent code to generate forgery images
similar to the input images. The AE’s objective function is usually
reconstruction loss(e.g., Euclidean distances between output and input).
Furthermore, the whole model is optimized by the backpropagation algorithm.

However, original AE does not restrict the distribution of latent code. Thus
variational autoencoder (VAE) [VAE] is proposed to optimize the distribution
of latent code. VAE used KL-divergence and the reparameterization trick to
restrict latent code distribution so that randomly sampled value from a normal
distribution can be decoded to an image with similar distribution to the
training samples. The latent code generated by VAE performs better in the
classification task, so VAE is widely used for unsupervised feature
extraction.

###  2.2 Adversarial Autoencoder

Adversarial autoencoder (AAE) [AAE] is another AE-based structure which adopts
adversarial restriction on the latent code. Unlike VAE, AAE does not use KL-
divergence to measure the distance between latent code distribution and normal
distribution. AAE chooses to use the idea of adversarial to restrict the
distribution of latent code. By training a generative adversarial network
(GAN) [GAN] with the latent code and randomly selected samples from a normal
distribution, the distribution of latent code will become an expected normal
distribution. As a result, the latent code generated by AAE will be helpful
for classification, and its distribution will be similar but different from
the VAE’s latent code distribution.

###  2.3 Prototypical Contrastive Learning

Usual unsupervised contrastive learning methods aim to find an embedding space
where every sample can be discriminated among other instances in the whole
dataset, such as momentum contrast (MoCo) [Moco]. But this objective may be
too hard to achieve without enough computing resources or sufficient
iteration. Although this kind of distribution may be beneficial in some
downstream tasks, it is not very fit for classification. Excellent
distribution for classification demands for clear clusters composed of the
data in the same classes. And typical unsupervised contrastive learning can’t
satisfy this demand.

PCL [PCL] is an unsupervised representation learning method based on
contrastive learning. Unlike most contrastive learning methods, PCL introduces
prototypes as latent variables to help find the maximum-likelihood estimation
of the network parameters in an expectation-maximization framework. PCL uses
ProtoNCE loss, a generalized version of the InfoNCE [InfoNCE] loss for
contrastive learning by encouraging representations to be closer to their
assigned prototypes. The experiments show that PCL can get better performances
than other baselines in the classification tasks, so we choose PCL as our
method’s base structure.

##  3 Proposed Method

###  3.1 VAE Module

To process the hyperspectral images, we design a VAE network to extract
features from the original hyperspectral images. The structure of the VAE
module is shown in Fig.1. As the figure shows, we first use PCA to preprocess
the whole hyperspectral image to reduce the number of the channels. Then
sliding windows are adopted to split the image into small patches with the
shape black$s\times s\times k$. $s$ is the size of sliding windows, and $k$ is
the channel of the preprocessed hyperspectral image. The structure of the
encoder network is enlightened by the idea of hybrid convolution in HybridSN
[HybridSN]. We used several 3D convolutional layers and 2D convolutional
layers to extract the input patches’ spectral and spatial information. We also
used several 3D transposed convolutional layers and 2D transposed
convolutional layers to build the decoder network.

Furthermore, to avoid the risk of overfitting, we adopted batch normalization
[bn] after every convolutional layer. Every linear layer is combined with a
rectified linear unit (ReLU) activation layer except the output layer. We use
a global pooling layer to fix the shape of the feature map, thus fix the shape
of extracted features (e.g., 1024-d) after reshaping operation.

According to the original VAE, we let the encoder output two variables $\mu$
and $\sigma$. Then the latent code can be computed by Eq. 1. $z$ is the latent
code, and $\varepsilon$ is a random value drawn from the normal distribution.
The loss function for VAE network can be divided into two parts:
$\mathcal{L}_{\mu,\sigma^{2}}$ and $\mathcal{L}_{recon}$.
$\mathcal{L}_{\mu,\sigma^{2}}$ controls the distribution of latent code, and
$\mathcal{L}_{recon}$ ensures the preservation of the input image’s
information. As Eq. 2 shows, $N$ is the number of samples,
$\mathcal{L}_{\mu,\sigma^{2}}$ computes the KL divergence between latent code
and normal distribution. Moreover, Eq. 3 shows the detail of reconstruction
loss $\mathcal{L}_{recon}$, $N$ is the number of samples, the input image
patch is $I^{x,y,z}$ with the shape of $s\times s\times k$, and
$\hat{I}^{x,y,z}$ is the output image patch of the decoder. The whole loss
function can be represented by Eq. 4.

Figure 1: Structure of VAE module.

${z=\mu+\varepsilon\times\sigma}$ (1)

${\mathcal{L}_{\mu,\sigma^{2}}=\frac{1}{2}\sum_{i=1}^{N}\left(\mu_{(i)}^{2}+\sigma_{(i)}^{2}-\log\sigma_{(i)}^{2}-1\right)}$
(2)

$\mathcal{L}_{recon}=\sum_{i=1}^{N}\left(\frac{1}{s\times s\times
k}\sum_{x=0}^{s-1}\sum_{y=0}^{s-1}\sum_{z=0}^{k-1}\left(I^{x,y,z}-\hat{I}^{x,y,z}\right)^{2}\right)$
(3)

$\mathcal{L}=\mathcal{L}_{\mu,\sigma^{2}}+\mathcal{L}_{recon}$ (4)

When we finished the VAE module training, we can map all the patches to VAE
features. Thus we can use these VAE features to build the training dataset for
ContrastNet.

###  3.2 AAE Module

Similar to the VAE module, the whole hyperspectral image is preprocessed by
PCA in the AAE module. The structure of the AAE network is shown in Fig.2. The
encoder and the decoder in the AAE module are the same with that in the VAE
module except for the linear layers before the latent code. The encoder in the
AAE module directly outputs the latent code and acts as a GAN generator.

According to the original AAE, the latent code is the output of the generator
and the discriminant input. Another input of the discriminant is a randomly
drawn variable from the normal distribution. The training of the AAE can be
divided into two phases: reconstruction phase and regularization phase. In the
reconstruction phase, we use reconstruction loss to optimize the parameters in
the network. Reconstruction loss is shown in Eq. 5, it is the same as the
reconstruction loss in the VAE module. $\mathcal{L}_{recon}^{AAE}$ is
reconstruction loss, its definition is the same as $\mathcal{L}_{recon}$.

In the regularization phase, the discriminator and the encoder are optimized.
At first, we optimize the discriminator and then the generator (encoder). In
order to stablize the training process of GAN, we adopt Wasserstein GAN (WGAN)
[WGAN] loss as our optimization function for the discriminator and the
generator. Discriminant loss and generator loss are shown in Eq. 6 and 7,
where $\mathcal{L}_{D}$ represents the loss of the discriminant, and
$\mathcal{L}_{G}$ represents the loss of the generator. $P_{g}$ is the
distribution of the generated samples(latent code), and $P_{r}$ is the
distribution of the real samples(samples from a normal distribution). x
represents random samples from $P_{g}$ and $P_{r}$.

Figure 2: Structure of AAE module.

$\mathcal{L}_{recon}^{AAE}=\sum_{i=1}^{N}\left(\frac{1}{s\times s\times
k}\sum_{x=0}^{s-1}\sum_{y=0}^{s-1}\sum_{z=0}^{k-1}\left(I^{x,y,z}-\hat{I}^{x,y,z}\right)^{2}\right)$
(5) ${\mathcal{L}_{D}=\mathbb{E}_{x\sim
P_{g}}\left[\mathbf{D}(x)\right]-\mathbb{E}_{x\sim
P_{r}}\left[\mathbf{D}(x)\right]}$ (6)

${\mathcal{L}_{G}=-\mathbb{E}_{x\sim P_{g}}\left[\mathbf{D}(x)\right]}$ (7)

Similarly, when we finished the training of the AAE module, we can map all the
patches to AAE features. Thus we can use these AAE features to build the
training dataset for ContrastNet.

###  3.3 ContrastNet

Before we introduce the ContrastNet, the specific description of typical
unsupervised contrast learning should be clarified. Let the blacktraining set
$X=\\{{x_{1},x_{2},...,x_{n}}\\}$ contain $n$ images. And our goal is to find
an embedding function $f_{\theta}$(realized via a deep neural network) that
maps $X$ to another space $V=\\{{v_{1},v_{2},...,v_{n}}\\}$ where
$v_{i}=f_{\theta}(x_{i})$, and $v_{i}$ can describe $x_{i}$ best. Usually the
objective is achieved by the InfoNCE loss function [cpc] as Eq. 8 shows.
$v^{\prime}_{i}$ is a positive embedding for instance $i$, and
$v^{\prime}_{j}$ includes one positive embedding (assume
$v^{\prime}_{0}=v^{\prime}_{i}$ ) and r negative embeddings for other
instances, and $\tau$ is a temperature hyper-parameter. In MoCo [Moco], these
embeddings are obtained by feeding $x_{i}$ to a momentum encoder parametrized
by $\theta^{\prime}$, $v^{\prime}_{i}=f_{\theta^{\prime}}(x_{i})$, where
$\theta^{\prime}$ is a moving average of $\theta$.

${\mathcal{L}_{\mathrm{InfoNCE}}=\sum_{i=1}^{n}-\log\frac{\exp\left(v_{i}\cdot
v_{i}^{\prime}/\tau\right)}{\sum_{j=0}^{r}\exp\left(v_{i}\cdot
v_{j}^{\prime}/\tau\right)}}$ (8)

When minimizing the Eq. 8, the distance between $v_{i}$ and $v^{\prime}_{i}$
are becoming close, and the distance between $v_{i}$ and $r$ negative
embeddings are becoming far. However, original InfoNCE loss uses a fixed
$\tau$, which means the same concentration on a certain scale. Moreover, it
computes the distance between $v_{i}$ and all kinds of negative embeddings,
some of which are not representative.

$\mathcal{L}_{Proto}=\sum_{i=1}^{n}-\left(\frac{1}{M}\sum_{m=1}^{M}\log\frac{\exp\left(v_{i}\cdot
c_{s}^{m}/\phi_{s}^{m}\right)}{\sum_{j=0}^{r}\exp\left(v_{i}\cdot
c_{j}^{m}/\phi_{j}^{m}\right)}\right)$ (9)

As a result, PCL proposed a new loss function named ProtoNCE. As Eq. 9 and Eq.
10 show, based on InfoNCE, ProtoNCE contains a second part
$\mathcal{L}_{Proto}$. It use prototypes $c$ instead of $v^{\prime}$, and
replace the fixed temperature $\tau$ with a per-prototype concentration
estimation $\phi$.

$\mathcal{L}_{ProtoNCE}=\mathcal{L}_{Proto}+\mathcal{L}_{\mathrm{InfoNCE}}$
(10)

The detailed algorithm of PCL can be found in Algorithm 1. There are two
encoders in the PCL algorithm, query encoder maps $x_{i}$ to $v_{i}$, and the
momentum encoder maps $x_{i}$ to $v^{\prime}_{i}$. Practically, we do not use
the same $x_{i}$ to get $v_{i}$ and $v^{\prime}_{i}$, two images with proper
data augmentation will be the input of two encoders for better performance.
blackFor example, we can generate $x^{\prime}_{i}$ by applying color
distortion to $x_{i}$, and get $x^{\prime\prime}_{i}$ by applying random crop
to $x_{i}$. Then we use $x^{\prime}_{i}$ and $x^{\prime\prime}_{i}$ to get
$v_{i}$ and $v^{\prime}_{i}$. The query encoder’s parameters can be
represented by $\theta$, and the momentum encoder’s parameters can be
represented by $\theta^{\prime}$ which meets Eq. 11. In algorithm 1,
$k-means(\cdot,\cdot)$ [faiss] is the GPU implementation of a widely used
clustering algorithm. And $C^{m}$ is the sets of $k_{m}$ prototypes(center of
clusters). $Concentration(\cdot,\cdot)$ is the function of concentration
estimation $\phi$, as Eq. 12 shows, $\phi$ is calculated by the momentum
features $\\{v^{\prime}_{z}\\}^{Z}_{z=1}$ that are within the same clusters as
a prototype $c$, $\alpha$ is a smooth parameter to ensure that small clusters
do not have overly-large $\phi$, and the typical value is 10.
$SGD(\cdot,\cdot)$ [sgd] is a kind of optimization function to update the
parameters in the query encoder.

Algorithm 1 Prototypical Contrastive Learning.

1:encoder $f_{\theta}$, training dataset $X$, number of clusters
$K=\\{k_{m}\\}^{M}_{m=1}$

2:$\theta^{\prime}=\theta$

3:while not MaxEpoch do

4: $V^{\prime}=f_{\theta}^{\prime}(X)$

5: for $m=1$ to $M$ do

6: $C^{m}=k-means(V^{\prime},k_{m})$

7: $\phi^{m}=Concentration(C^{m},V^{\prime})$

8: end for

9: for $x$ in $Dataloader(X)$ do

10: $v=f_{\theta}(x),v^{\prime}=f^{\prime}_{\theta}(x)$

11:
$\mathcal{L}_{ProtoNCE}(v,v^{\prime},\\{C^{m}\\}^{M}_{m=1},\\{\phi^{m}\\}^{M}_{m=1})$

12: $\theta=SGD(\mathcal{L}_{ProtoNCE},\theta)$

13: $\theta^{\prime}=0.999*\theta^{\prime}+0.001*\theta$

14: end for

15:end while

$\theta^{\prime}=0.999*\theta^{\prime}+0.001*\theta$ (11)

$\phi=\frac{\sum_{z=1}^{Z}\left\|v_{z}^{\prime}-c\right\|_{2}}{Z\log(Z+\alpha)}$
(12)

The proposed ContrastNet is similar to the original PCL algorithm, but there
are several differences. Firstly, ContrastNet uses AAE and VAE features as
inputs, not two perturbed images. Secondly, the structure of ContrastNet is
more straightforward and specialized in the ”contrastive” part, not the
”encoding” part. Thirdly, we adopt the projection head trick in sim-CLR
[simclr] to improve the performance.

Figure 3:  Structure of ContrastNet. Clustering means k-means algorithm, and
ProtoNCE is $\mathcal{L}_{ProtoNCE}$ in Eq. 10. M-step and E-step are the
phases of Expectation-Maximization algorithm.

Fig. 3 is an overview of the ContrastNet. Expectation-Maximization (EM)
algorithm is adopted to train the ContrastNet. In the E-step, VAE features are
fed into the momentum encoder, then $C^{m}$ and $\phi^{m}$ are calculated. In
the M-step, $\mathcal{L}_{ProtoNCE}$ is calculated based on the updated
features and variables in the E-step, then the query encoder is updated via
backpropagation, the momentum encoder is updated via moving average (Eq. 11).
Similarly, all convolutions in the ContrastNet are composed of a convolutional
layer and a batch normalization layer, and every linear layer is combined with
a ReLU layer. The 1024-d AAE and VAE features will be reshaped to $4\times
4\times 64$ feature maps first, in this way, the spatial information extracted
can be preserved. And then, the feature maps are fed into the ContrastNet. The
structure of ContrastNet is light and straightforward because the previous
process has achieved the feature extraction. The ContrastNet only needs to pay
attention to the contrastive learning. We adopted the projection head
structure in the ContrastNet, that is, we use the output of the projection
head for training, and use the features before the head for testing.
Projection head acts as a nonlinear trasformation $g()$, and embedded variable
$v$ is the transformation of contrast feature $z$, just as Eq. 13 shows.

$v=g(z)$ (13)

We use $v$ to supervise the update of ContrastNet in the training phase, and
use $z$ for downstream work in the testing phase. In this way, the extracted
contrast features will contain more original information. In the testing
phase, we only need the query encoder to extract contrast features in Fig. 3,
thus can save much testing time for the ContrastNet.

##  4 Experiments and Analysis

To illustrate the great feature learning ability of ContrastNet, we use
classification as application to validate the performance. As a simple and
commonly used classifier, SVM was used for highlighting the effectiveness of
the extracted features, which were learned by ContrastNet. Classification
experiments were conducted using three publicly available HSI datasets: Indian
Pines(IP), University of Pavia(PU), and Salinas Scene(SA). Eleven other
methods were adopted as baselines for comparison.

They are linear discriminant analysis (LDA) [LDA], locality-preserving
dimensionality reduction (LFDA) [LFDA], sparse graph-based discriminant
analysis (SGDA) [SGDA], sparse and low-rank graph for discriminant analysis
(SLGDA) [SLGDA], deep convolutional neural networks (1D-CNN) [1D-CNN],
supervised deep feature extraction (S-CNN) [S-CNN], tensor principal component
analysis (TPCA) [TPCA], stacked sparse autoencoder (SSAE) [SSAE], unsupervised
deep feature extraction (EPLS) [EPLS], and 3D convolutional autoencoder
(3DCAE) [3DCAE]. blackAll experiments are repeated three times, and we
recorded the mean value of all results.

###  4.1 Datasets

Experiments are conducted using three publicly available HSI datasets: Indian
Pines (IP), University of Pavia (PU), and Salinas Scene (SA). The IP dataset
was gathered using the Airborne Visible / Infrared Imaging Spectrometer
(AVIRIS) hyperspectral sensor. It consisted of 145$\times$145 pixels and 224
spectral reflectance bands in the wavelength range of $0.4\sim 2.5\times
10^{-6}$ meters. After removing bands covering the region of water absorption,
200 bands were left, and the ground truth was composed of 16 classes. In this
dataset, there is a serious problem of uneven sample size. The class with the
largest number of samples has more than two thousand samples, while the class
with the smaller number of samples has only 20 samples. The details of the IP
dataset are shown in Table 1.

Table 1: Details of the IP dataset.  IP Dataset
---
Class Name | Index | Number
Alfalfa | 0 | 46
Corn-notill | 1 | 1428
Corn-mintill | 2 | 830
Corn | 3 | 237
Grass-pasture | 4 | 483
Grass-trees | 5 | 730
Grass-pasture-mowed | 6 | 28
Hay-windrowed | 7 | 478
Oats | 8 | 20
Soybean-notill | 9 | 972
Soybean-mintill | 10 | 2455
Soybean-clean | 11 | 593
Wheat | 12 | 205
Woods | 13 | 1265
Buildings-Grass-Trees-Drives | 14 | 386
Stone-Steel-Towers | 15 | 93
Total |  | 10249

The PU dataset was a hyperspectral image with 610$\times$340 pixels and 103
spectral bands in the wavelength range of $0.43\sim 0.86\times 10^{-6}$
meters. The ground truth contained nine classes. The characteristics of this
dataset are that the background pixels are complex, the area of the connected
region is small, and there are many discontinuous pixels. The details of the
PU dataset are shown in Table 2.

Table 2: Details of the PU dataset.  UP Dataset
---
Class Name | Index | Number
Asphalt | 0 | 6631
Meadows | 1 | 18649
Gravel | 2 | 2099
Trees | 3 | 3064
Painted metal sheets | 4 | 1345
Bare Soil | 5 | 5029
Bitumen | 6 | 1330
Self-Blocking Bricks | 7 | 3682
Shadows | 8 | 947
Total |  | 42776

The SA dataset was an image with 512$\times$217 pixels and 224 spectral bands.
This dataset was in the wavelength range of $0.36\sim 2.5\times 10^{-6}$
meters. Twenty water-absorbing bands were discarded from the SA dataset, so it
was a dataset with 204 spectral bands and 16 landcover classes in the ground
truth. The background pixels of this dataset are not complex, and there is no
problem of uneven number of samples. Compared with the other two datasets, it
is easier to classify. The details of the SA dataset are shown in Table 3.

Table 3: Details of the SA dataset.  SA Dataset
---
Class Name | Index | Number
Brocoli_green_weeds_1 | 0 | 2009
Brocoli_green_weeds_2 | 1 | 3726
Fallow | 2 | 1976
Fallow_rough_plow | 3 | 1394
Fallow_smooth | 4 | 2678
Stubble | 5 | 3959
Celery | 6 | 3579
Grapes_untrained | 7 | 11271
Soil_vinyard_develop | 8 | 6203
Corn_senesced_green_weeds | 9 | 3278
Lettuce_romaine_4wk | 10 | 1068
Lettuce_romaine_5wk | 11 | 1927
Lettuce_romaine_6wk | 12 | 916
Lettuce_romaine_7wk | 13 | 1070
Vinyard_untrained | 14 | 7268
Vinyard_vertical_trellis | 15 | 1807
Total |  | 54129

###  4.2 Parameters Setting

####  4.2.1 Training sample and testing sample selection

10% samples of each class are randomly selected for training SVM in the Indian
Pines and University of Pavia datasets, and 5% samples of each class are used
for training in the Salinas dataset, the rest in each dataset is used for
testing.

####  4.2.2 Parameters in model

All experiments were implemented on a PC with 16GB RAM and a GTX 1080 GPU. The
coding environment was pytorch [pytorch]. The shape of patches is $27\times
27\times 15$ in PU and SA dataset, and $27\times 27\times 30$ in the IP
dataset. When we use $27\times 27\times 15$ patches, our networks’ details are
shown as follows.

Table 4 and Table 5 show the structure of AAE. The structure of the
discriminant is simple, so it is omitted. The latent code dimension is 128 and
extracted AAE features are flattened outputs of the pooling layer in the
encoder(1024-d). In the training phase, we use two Adam [adam] optimizers to
optimize the encoder and the decoder, and two SGD [sgd] optimizers to optimize
the generator and the discriminant. Two Adam optimizers are set with the
learning rate of 0.001, and weight decay is set to 0.0005. The optimizer for
the generator is set with a 0.0001 learning rate and no weight decay. The
optimizer for the discriminant is set with a 0.00005 learning rate and no
weight decay. The batch size is set to 128. We train the AAE module for 20
epochs, and only save the parameters of the encoder.

Table 4: Structure of AAE Encoder. The shape of data is defined in pytorch
style. -1 means batchsize in the shape array. AAE Encoder
---
Layer(type) | Output Shape
input | [-1,1,15,27,27]
Conv3d | [-1,8,9,25,25]
BatchNorm3d | [-1,8,9,25,25]
ReLU | [-1,8,9,25,25]
Conv3d | [-1,16,5,23,23]
BatchNorm3d | [-1,16,5,23,23]
ReLU | [-1,16,5,23,23]
Conv3d | [-1,32,3,21,21]
BatchNorm3d | [-1,32,3,21,21]
ReLU | [-1,32,3,21,21]
Conv2d | [-1,64,19,19]
BatchNorm2d | [-1,64,19,19]
ReLU | [-1,64,19,19]
AdaptiveAvgPool2d | [-1,64,4,4]
Linear | [-1,512]
ReLU | [-1,512]
Linear | [-1,128]
Table 5: Structure of AAE Decoder. The shape of data is defined in pytorch
style. -1 means batchsize in the shape array. AAE Decoder
---
Layer(type) | Output Shape
input | [-1,128]
Linear | [-1,256]
ReLU | [-1,256]
Linear | [-1,23104]
ReLU | [-1,23104]
ConvTransposed2d | [-1,96,21,21]
BatchNorm2d | [-1,96,21,21]
ReLU | [-1,96,21,21]
ConvTransposed3d | [-1,16,5,23,23]
BatchNorm3d | [-1,16,5,23,23]
ReLU | [-1,16,5,23,23]
ConvTransposed3d | [-1,8,9,25,25]
BatchNorm3d | [-1,8,9,25,25]
ReLU | [-1,8,9,25,25]
ConvTransposed3d | [-1,1,15,27,27]
BatchNorm3d | [-1,1,15,27,27]

Table 6 and Table 7 show the structure of VAE. The dimension of latent code is
also 128 and extracted VAE features are flattened outputs of the pooling layer
in the encoder(1024-d) too. In the training phase, we use two Adam [adam]
optimizers to optimize the encoder and the decoder. Two Adam optimizers are
set with the learning rate of 0.001, and weight decay is set to 0.0005. The
two outputs of VAE encoder are used to compute the latent code. Then the
latent code is fed into the decoder. The batch size is set to 128. We train
the VAE module for 30 epochs, and only save the parameters of the encoder.

Table 6: Structure of VAE Encoder. The shape of data is defined in pytorch
style. -1 means batchsize in the shape array. VAE Encoder
---
input | [-1,1,15,27,27]
Conv3d | [-1,8,9,25,25]
BatchNorm3d | [-1,8,9,25,25]
ReLU | [-1,8,9,25,25]
Conv3d | [-1,16,5,23,23]
BatchNorm3d | [-1,16,5,23,23]
ReLU | [-1,16,5,23,23]
Conv3d | [-1,32,3,21,21]
BatchNorm3d | [-1,32,3,21,21]
ReLU | [-1,32,3,21,21]
Conv2d | [-1,64,19,19]
BatchNorm2d | [-1,64,19,19]
ReLU | [-1,64,19,19]
AdaptiveAvgPool2d | [-1,64,4,4]
Linear | [-1,512]
ReLU | [-1,512]
Linear/Linear | [-1,128]/[-1,128]
Table 7: Structure of VAE Decoder. The shape of data is defined in pytorch
style. -1 means batchsize in the shape array. VAE Decoder
---
Layer(type) | Output Shape
input | [-1,128]
Linear | [-1,256]
ReLU | [-1,256]
Linear | [-1,23104]
ReLU | [-1,23104]
ConvTransposed2d | [-1,96,21,21]
BatchNorm2d | [-1,96,21,21]
ReLU | [-1,96,21,21]
ConvTransposed3d | [-1,16,5,23,23]
BatchNorm3d | [-1,16,5,23,23]
ReLU | [-1,16,5,23,23]
ConvTransposed3d | [-1,8,9,25,25]
BatchNorm3d | [-1,8,9,25,25]
ReLU | [-1,8,9,25,25]
ConvTransposed3d | [-1,1,15,27,27]
BatchNorm3d | [-1,1,15,27,27]

The two encoders in ContrastNet have the same structure, which is shown in
Table 8. Two inputs are 1024-d vectors, and the dimension of contrast features
and projected features are both 128-d. We train the ContrastNet with an SGD
optimizer, and the learning rate is 0.003, weight decay is set to 0.001.
Number of negative samples $r=640$, and batch size is 128, $\tau$ in Eq. 8 is
0.01, number of clusters $K=[1000,1500,2500]$. According to the original PCL
method, we only train the network with Eq. 8 for 30 epochs, then train the
network with Eq. 10 for 170 epochs. The learning rate will be multiplied with
0.1 when the trained epoch is more than 120 and 160.

Table 8: Structure of ContrastNet Encoder. The shape of data is defined in
pytorch style. -1 means batchsize in the shape array. ContrastNet Encoder
---
Layer(type) | Output Shape
input | [-1,1024]
ConvTransposed2d | [-1,64,6,6]
BatchNorm2d | [-1,64,6,6]
ReLU | [-1,64,6,6]
ConvTransposed2d | [-1,64,8,8]
BatchNorm2d | [-1,64,8,8]
ReLU | [-1,64,8,8]
Conv2d | [-1,128,6,6]
BatchNorm2d | [-1,128,6,6]
ReLU | [-1,128,6,6]
Conv2d | [-1,64,4,4]
BatchNorm2d | [-1,64,4,4]
ReLU | [-1,64,4,4]
Conv2d | [-1,32,2,2]
ReLU | [-1,32,2,2]
Linear | [-1,128]
ReLU | [-1,128]
Linear | [-1,128]

###  4.3 Criteria

The quantitative results are evaluated by overall accuracy (OA) and average
accuracy (AA).

As Eq. 14 shows, overall accuracy is represented as $\mathbf{S}_{OA}$,
$N_{correct}$ means the sum of samples that are classified correctly, and
$N_{total}$ means all labeled samples. When applied to binary classification,
$\mathbf{S}_{OA}$ means accuracy. $\mathbf{S}_{OA}$ ranges from 0 to 1, and
the bigger the $\mathbf{S}_{OA}$ is, the better the classifier is.
$\mathbf{S}_{OA}$ treats different classes by weighting the amount of data in
each class, but it tends to be dominated by classes that contain more samples
than other classes.

$\mathbf{S}_{OA}=\frac{N_{correct}}{N_{total}}$ (14)

Eq. 15 is the definition of average accuracy, $\mathbf{S}_{AA}$ means average
accuracy, $n$ is the number of classes, $\text{accuracy}_{i}$ is the
classification accuracy of $\mathbf{i}$th class($i\leq n$). Similar to
$\mathbf{S}_{OA}$, $\mathbf{S}_{AA}$ also ranges from 0 to 1, and it is also
expected to be close to 1. However, when there are some classes with few
samples, the statistical variable deviation of accuracies of these classes
will be high, thus leading to distortion of the average accuracy.

$\mathbf{S}_{AA}=\frac{\sum\limits_{i=1}^{n}\text{accuracy}_{i}}{n}$ (15)

###  4.4 Results and Analysis

####  4.4.1 Classification Experiments

The results in the Indian Pines dataset are shown in Table 9. It is clear that
the proposed structures performed best in many classes. AAE got the best
results in five classes while ContrastNet got the best results in six classes,
and VAE performed best in 3 classes.

To some extent, AAE and VAE are two better feature extractors than 3DCAE in
term of AA, and AAE is better than VAE in both AA and OA. However, VAE shows
great ability in extracting features from the classes with few numbers, such
as class 0, 6, 8 and 15. The number of samples in these classes are less than
100, while other classes may have over a thousand samples. In these classes,
VAE performs better than AAE, though in other classes it is not as well as AAE
is. As for ContrastNet, it shows a significant increase in OA than AAE and
VAE, thus outperform all methods. This can be an evidence indicating that
ContrastNet can extract some core information from two distributions by
contrastive learning. However, ContrastNet did not perform well enough in AA
in the Indian Pines dataset. According to the table, the accuracy of
ContrastNet are lower than AAE’s and VAE’s in class 0, 3, 5, 6, 7, 8, 12, 14,
15. And this problem is severe in the classes containing few samples, such as
class 0, 6, 8 and 15. We guess the reason for it may be the imbalance of the
number of samples in the Indian Pines dataset. Classes with large numbers
dominated the location of clusters. Thus the ContrastNet tended to perform
well in these classes, which decreased accuracy in other classes. However,
ContrastNet still shows the best OA, and its AA is very close to 3DCAE, which
outperforms other baselines. In other words, ContrastNet makes a trade-off
between AA and OA. It tends to increase OA at the expense of lowering AA. As a
result, ContrastNet gets the best OA and relatively good AA in the IP dataset.

Generally speaking, the AAE and VAE modules we designed show excellent ability
in extracting features unsupervisely, and VAE is better at extracting
information from classes with a small sample size. The ContrastNet we proposed
utilized the features generated by AAE and VAE, further improved the ability
of extracting features. However, limited by the clustering idea it introduces,
it improves the OA while reducing the focus on the few sample classes, thus
leading to a decline in AA.

Table 9: Comparison in the Indian Pines dataset. 10% labeled samples are used
for training, and the rest are used for testing. AA menas average accuracy and
OA means overall accuracy. The best values in each row are in bold. And the
methods in bold are structures proposed in this paper. Some of the data in the
table is quoted from [3DCAE].

Class | Supervised Feature Extraction | Unsupervised Feature Extraction
---|---|---
LDA | LFDA | SGDA | SLGDA | 1D-CNN | S-CNN | PCA | TPCA | SSAE | EPLS | 3DCAE | AAE | VAE | ContrastNet
0 | 58.54 | 29.63 | 42.59 | 39.62 | 43.33 | 83.33 | 39.02 | 60.97 | 56.25 | 58.72 | 90.48 | 100.00 | 100.00 | 85.37
1 | 69.88 | 75.59 | 80.89 | 85.56 | 73.13 | 81.41 | 72.30 | 87.00 | 69.58 | 59.91 | 92.49 | 81.63 | 78.78 | 97.15
2 | 65.86 | 75.42 | 65.71 | 74.82 | 65.52 | 74.02 | 72.02 | 94.51 | 75.36 | 71.34 | 90.37 | 95.27 | 92.37 | 97.95
3 | 73.71 | 58.12 | 64.10 | 49.32 | 51.31 | 71.49 | 55.87 | 79.34 | 64.58 | 74.31 | 86.90 | 99.22 | 97.34 | 95.62
4 | 90.32 | 95.17 | 94.57 | 95.35 | 87.70 | 90.11 | 93.09 | 93.08 | 88.81 | 97.95 | 94.25 | 95.17 | 93.87 | 96.09
5 | 92.09 | 96.12 | 98.39 | 95.59 | 95.10 | 94.06 | 94.67 | 96.34 | 87.00 | 96.44 | 97.07 | 98.73 | 98.27 | 96.80
6 | 96.00 | 11.53 | 50.00 | 36.92 | 56.92 | 84.61 | 80.00 | 76.00 | 90.00 | 54.02 | 91.26 | 96.00 | 98.67 | 70.67
7 | 98.14 | 93.87 | 99.80 | 99.75 | 96.64 | 98.37 | 98.37 | 99.76 | 89.72 | 88.99 | 97.79 | 99.84 | 99.77 | 98.68
8 | 11.11 | 0.00 | 0.00 | 5.00 | 28.89 | 33.33 | 88.89 | 100.00 | 100.00 | 58.89 | 75.91 | 96.30 | 98.15 | 70.37
9 | 73.80 | 80.89 | 59.30 | 69.11 | 75.12 | 86.05 | 74.49 | 79.51 | 77.19 | 73.10 | 87.34 | 87.01 | 78.86 | 97.45
10 | 55.41 | 83.02 | 84.44 | 89.91 | 83.49 | 82.98 | 69.58 | 85.42 | 77.58 | 70.78 | 90.24 | 89.08 | 81.75 | 98.40
11 | 76.92 | 86.32 | 77.04 | 86.78 | 67.55 | 73.40 | 65.29 | 84.24 | 72.00 | 57.51 | 95.76 | 93.51 | 90.64 | 93.57
12 | 91.30 | 79.25 | 99.09 | 99.51 | 96.86 | 87.02 | 98.37 | 98.91 | 87.80 | 99.25 | 97.49 | 98.56 | 98.56 | 95.32
13 | 93.32 | 88.87 | 92.50 | 96.45 | 96.51 | 94.38 | 91.39 | 98.06 | 93.48 | 95.07 | 96.03 | 95.73 | 93.24 | 98.51
14 | 67.72 | 60.00 | 68.68 | 61.79 | 39.08 | 75.57 | 48.99 | 87.31 | 72.36 | 91.26 | 90.48 | 97.31 | 97.02 | 96.73
15 | 90.36 | 53.68 | 85.26 | 84.16 | 89.40 | 79.76 | 87.95 | 96.38 | 97.22 | 91.27 | 98.82 | 98.02 | 98.81 | 79.76
AA(%) | 76.89 | 66.72 | 72.65 | 73.10 | 71.66 | 84.44 | 76.89 | 89.31 | 81.18 | 77.43 | 92.04 | 95.09 | 93.51 | 91.78
OA(%) | 76.88 | 81.79 | 80.05 | 85.19 | 79.66 | 80.72 | 76.88 | 88.55 | 79.78 | 77.18 | 92.35 | 91.80 | 88.03 | 97.08

The results in the SA dataset are shown in Table 10. Similar to the IP dataset
results, AAE performed best in five classes, VAE and ContrastNet both
performed best in three classes. Since the background of the dataset is simple
and there is no problem of sample size imbalance, many methods perform well on
this dataset. AAE and VAE still get higher AA than other baselines, which
indicates their great feature extraction ability. However, there are some
methods outperform AAE and VAE in OA. S-CNN, TPCA, 3DCAE get higher OA than
VAE, and S-CNN also gets higher OA than AAE. For VAE, it basically shows
similar performance to AAE except in class 7 and 14. The two classes have
11271 and 7268 samples respectively, while other classes basically have less
than 5000 samples. It seems the performance of VAE will be affected by the
ratio of the number of samples between each class. According to the AA and OA,
the performance of VAE is close to 3DCAE. For AAE, it has less OA than S-CNN
but close, it is still an efficient feature extraction method. In this simpler
dataset, AAE and VAE also keep relatively better performances than many
methods, and this show their good generalization ability.

Though all methods got good performances in this dataset, ContrastNet still
got the highest OA and got the highest AA. Unlike the IP dataset, the number
of samples in the SA dataset is relatively balanced, so ContrastNet got a
higher AA than AAE and VAE, and there is no obvious decrease of accuracy in
each class. On the contrary, the accuracy of several classes has been
significantly improved after combining the characteristics of AAE and VAE,
such as class 7 and 14. There are two reasons for this phenomenon. On the one
hand, ContrastNet has utilized the information shared by AAE and VAE. On the
other hand, the sample size of these two classes is relatively large, thus may
contain more cluster centers. Thanks to the introduction of clustering ideas,
ContrastNet is easier to extract features from such classes. As a result,
ContrastNet outperforms all methods in AA and OA.

Table 10: Comparison in the Salinas dataset. 5% labeled samples are used for
training, and the rest are used for testing. AA menas average accuracy and OA
means overall accuracy. The best values in each row are in bold. And the
methods in bold are structures proposed in this paper. Some of the data in the
table is quoted from [3DCAE].

Class | Supervised Feature Extraction | Unsupervised Feature Extraction
---|---|---
| LDA | LFDA | SGDA | SLGDA | 1D-CNN | S-CNN | PCA | TPCA | SSAE | EPLS | 3DCAE | AAE | VAE | ContrastNet
0 | 99.16 | 99.20 | 99.65 | 98.14 | 97.98 | 99.55 | 97.48 | 99.88 | 100.00 | 99.99 | 100.00 | 99.98 | 99.74 | 99.93
1 | 99.94 | 99.19 | 99.33 | 99.44 | 99.25 | 99.43 | 99.52 | 99.49 | 99.52 | 99.92 | 99.29 | 100.00 | 99.79 | 99.80
2 | 99.79 | 99.75 | 99.30 | 99.29 | 94.43 | 98.81 | 99.41 | 99.04 | 94.24 | 98.75 | 97.13 | 100.00 | 100.00 | 99.95
3 | 99.77 | 99.71 | 99.21 | 99.57 | 99.42 | 97.45 | 99.77 | 99.84 | 99.17 | 98.52 | 97.91 | 99.09 | 99.57 | 98.01
4 | 98.98 | 98.47 | 99.07 | 98.06 | 96.60 | 97.96 | 98.70 | 98.96 | 98.82 | 98.33 | 98.26 | 99.42 | 99.71 | 99.48
5 | 99.89 | 99.09 | 99.57 | 99.32 | 99.51 | 99.83 | 99.65 | 99.80 | 100.00 | 99.92 | 99.98 | 99.97 | 99.86 | 99.94
6 | 99.97 | 99.66 | 99.27 | 99.33 | 99.27 | 99.59 | 99.94 | 99.84 | 99.94 | 97.69 | 99.64 | 99.93 | 99.96 | 99.79
7 | 81.84 | 86.89 | 89.78 | 89.48 | 86.79 | 94.40 | 83.90 | 84.11 | 80.73 | 78.86 | 91.58 | 91.21 | 86.91 | 99.53
8 | 99.90 | 97.34 | 100.00 | 99.65 | 99.08 | 98.85 | 99.97 | 99.60 | 99.47 | 99.54 | 99.28 | 99.69 | 99.99 | 99.71
9 | 96.31 | 95.85 | 97.99 | 97.94 | 93.71 | 97.35 | 96.89 | 95.76 | 92.12 | 95.98 | 96.65 | 98.46 | 96.54 | 99.80
10 | 99.61 | 97.94 | 99.53 | 99.06 | 94.55 | 97.71 | 96.84 | 96.14 | 96.62 | 98.60 | 97.74 | 99.57 | 100.00 | 99.80
11 | 99.67 | 99.64 | 100.00 | 100.00 | 99.59 | 98.73 | 99.95 | 99.07 | 97.75 | 99.44 | 98.84 | 100.00 | 99.44 | 99.98
12 | 99.20 | 97.38 | 98.47 | 97.82 | 97.50 | 96.72 | 99.54 | 100.00 | 95.81 | 98.85 | 99.26 | 99.89 | 97.24 | 98.20
13 | 96.56 | 93.18 | 96.07 | 90.47 | 94.08 | 95.22 | 97.24 | 95.74 | 96.65 | 98.56 | 97.49 | 99.08 | 96.49 | 98.62
14 | 73.60 | 63.04 | 65.40 | 69.51 | 66.52 | 95.61 | 76.68 | 79.54 | 79.73 | 83.13 | 87.85 | 93.69 | 87.90 | 99.53
15 | 98.48 | 99.00 | 99.34 | 99.00 | 97.48 | 99.44 | 97.90 | 98.40 | 99.12 | 99.50 | 98.34 | 99.73 | 99.61 | 99.57
AA(%) | 96.42 | 95.33 | 96.37 | 96.01 | 94.73 | 97.39 | 96.46 | 93.24 | 95.61 | 96.55 | 97.45 | 98.73 | 97.67 | 99.48
OA(%) | 92.18 | 91.22 | 91.82 | 93.31 | 91.30 | 97.92 | 92.87 | 96.57 | 92.11 | 92.35 | 95.81 | 97.10 | 95.23 | 99.60

Experiment results in the PU dataset are shown in Table 11. Similarly, VAE
shows a relatively good performance and AAE outperforms all baselines in AA
and OA. VAE shows a surprisingly low accuracy in class 0, this may caused by
the comlex background of the PU dataset. Class 0 in the PU dataset is asphalt,
and it is often distributed in elongated connected areas, which makes the
pixels around the center pixels more complicated, and it is difficult for VAE
to encode such patches efficiently. However, AAE still performs well on such a
complex dataset, especially in class 6.

Although the background is complex, ContrastNet still outperforms all methods
in this dataset. In the PU dataset, because of the use of the characteristics
of AAE and VAE for contrastive learning, ContrastNet’s accuracy on most
classes surpasses AAE and VAE. Furthermore, it surpasses all methods in class
0, 1, 2 and 7. This result shows great power of ContrastNet in feature
learning.

Table 11: Comparison in the University of Pavia dataset. 10% labeled samples
are used for training, and the rest are used for testing. AA menas average
accuracy and OA means overall accuracy. The best values in each row are in
bold. And the methods in bold are structures proposed in this paper. Some of
the data in the table is quoted from [3DCAE].

Class | Supervised Feature Extraction | Unsupervised Feature Extraction
---|---|---
| LDA | LFDA | SGDA | SLGDA | 1D-CNN | S-CNN | PCA | TPCA | SSAE | EPLS | 3DCAE | AAE | VAE | ContrastNet
0 | 78.41 | 93.45 | 91.37 | 94.66 | 90.93 | 95.40 | 90.89 | 96.17 | 95.72 | 95.95 | 95.21 | 95.39 | 81.31 | 99.49
1 | 83.69 | 97.36 | 97.23 | 97.83 | 96.94 | 97.31 | 93.27 | 97.95 | 94.13 | 95.91 | 96.06 | 98.96 | 97.65 | 99.98
2 | 73.02 | 71.41 | 66.08 | 77.27 | 69.43 | 81.21 | 82.60 | 86.50 | 87.47 | 94.33 | 91.32 | 97.49 | 91.00 | 99.06
3 | 93.68 | 91.00 | 91.19 | 93.18 | 90.32 | 95.83 | 92.41 | 94.84 | 96.91 | 99.28 | 98.28 | 96.94 | 94.79 | 97.75
4 | 100.00 | 97.99 | 99.33 | 98.51 | 99.44 | 99.91 | 98.98 | 100.00 | 99.76 | 99.92 | 95.55 | 99.94 | 99.94 | 99.81
5 | 88.51 | 87.55 | 80.31 | 90.08 | 73.69 | 95.29 | 92.00 | 94.76 | 95.76 | 93.57 | 95.30 | 99.29 | 99.91 | 99.90
6 | 85.75 | 80.23 | 75.26 | 85.34 | 83.42 | 87.05 | 85.83 | 91.89 | 91.18 | 98.17 | 95.14 | 100.00 | 99.42 | 99.83
7 | 74.49 | 86.98 | 84.11 | 90.49 | 83.65 | 87.35 | 82.96 | 89.04 | 82.47 | 91.23 | 91.38 | 97.90 | 94.77 | 98.79
8 | 99.11 | 99.26 | 95.14 | 99.37 | 98.23 | 95.66 | 100.00 | 98.94 | 100.00 | 99.78 | 99.96 | 96.32 | 88.62 | 94.84
AA(%) | 86.29 | 89.47 | 86.67 | 91.86 | 87.34 | 94.75 | 90.99 | 95.64 | 93.71 | 96.33 | 95.36 | 98.03 | 94.16 | 98.83
OA(%) | 83.75 | 92.77 | 90.58 | 94.15 | 89.99 | 92.78 | 91.37 | 94.45 | 93.51 | 95.13 | 95.39 | 98.14 | 94.53 | 99.46

Generally speaking, the proposed AAE, VAE, and ContrastNet in this paper are
good models that outperform many supervised and unsupervised methods.

blackAt first glance, it seems counterintuitive that unsupervised methods can
surpass supervised methods. However, as far as this experiment is concerned,
all methods have the same number of labels, and the accuracy of classification
majorly depends on feature extraction. The supervised method directly uses the
label information and the backpropagation algorithm to extract features (10%
or fewer samples are used), while the unsupervised method extracts the
features based on the characteristics of the data itself (all samples are
used). In terms of the sample size used when extracting features, unsupervised
methods are much larger than supervised methods. Therefore, using unsupervised
methods for feature learning is easier to learn representative and robust
generalization features.

In contrast, the features learned by supervised methods are more susceptible
to the difference between the distribution of labeled samples and the
distribution of the real world. Moreover, using different network structures
for feature learning will have different learning efficiency. Therefore,
unsupervised feature learning may surpass supervised methods in feature
extraction.

Among the three methods, autoencoder-based AAE and VAE are models that make
use of reconstruction information and concentrate on the single sample itself.
ContrastNet is a model that learns the similarity between two different
distributions of a single sample, and also learn the dissimilarity between
positive samples and negative samples. According to the experiment results,
AAE is a better feature extraction method than VAE under most circumstances;
VAE shows potential in extracting features from the class with small sample
size; ContrastNet performs better than AAE and VAE in OA and AA when the
number of samples is balanced. When the number of samples is imbalanced,
ContrastNet tends to get higher OA and lower AA than AAE and VAE.

####  4.4.2 Visualization Experiments

To intuitively show the clustering effect of ContrastNet, we use t-SNE [TSNE]
method to visualize the features extracted by AAE, VAE, and ContrastNet. The
visualizations are shown in Fig. 4, 5, 6. In Fig. 4, it is a little hard to
tell which feature distribution is better for classification in the IP
dataset, and our experiments also proved that. blackFeatures in the PU dataset
are shown in Fig. 5, clustering effect is obvious in the ContrastNet feature
distribution, such as the features in class 0 and 8. Features in the same
class tend to be close to each other in the ContrastNet, which is due to the
clustering algorithm in ContrastNet. Similarly, features in the SA dataset are
shown in Fig. 6. blackAAE features in class 7 and VAE features in class 1 are
mixed with other classes. However, they are both separated from other classes
when they become contrast features. This is a good illustration of our
method’s advantage. As a result, we can roughly make a conclusion that
ContrastNet has the ability to decrease the inner-class distances and
increasing the inter-class distances.

(a) AAE

(b) VAE

(c) ContrastNet

Figure 4: Visualization of extracted features in the Indian Pines dataset.
Each class is corresponding to a kind of color.

(a) AAE

(b) VAE

(c) ContrastNet

Figure 5: Visualization of extracted features in the University of Pavia
dataset. Each class is corresponding to a kind of color.

(a) AAE

(b) VAE

(c) ContrastNet

Figure 6: Visualization of extracted features in the Salinas dataset. Each
class is corresponding to a kind of color. Figure 7: Parameter analysis based
on the IP dataset. When altering the temperature value and window size, the
performance of different modules are different. In this figure, $\tau$ means
the temperature parameter in Eq. 8, and horizontal axis is the window size for
input patches.

####  4.4.3 Parameter Analysis

blackThe window size and temperature value $\tau$ may affect the performance
of our models, so we conducted some parameter analysis experiments. We vary
the value of $\tau$ in $[0.01,0.05,0.1]$, and vary the window size in
$[19,21,23,25,27]$. With different parameters, we recorded the OA of
classification results in the IP dataset. Similar to the classification
experiments before, we use 10% labels and features for training the SVM and
the rest for testing. As Fig. 7 shows, classification results of contrast
features are significantly better than AAE features and VAE features with all
kinds of parameters. The value of $\tau$ only affects the contrast features,
so VAE features and AAE features do not show relevance with $\tau$. The
performance of contrast features is relevant to $\tau$ when training the
ContrastNet. It is not hard to tell from the figure that when $\tau=0.01$,
contrast features can get the best classification results. Moreover, when the
window size increases, performances of VAE features and AAE features become
better. It may owe to the increase of spatial information in patches. However,
the change of window size does not make a tremendous difference in contrast
features’ performance. Because ContrastNet can get valuable information from
contrastive learning, and the contribution of increased spatial information is
relatively limited. These experiments show that ContrastNet is a robust and
parameter-insensitive model.

####  4.4.4 Computational Aalysis

The computational analysis is shown in Table 12. Because the training of
ContrastNet demands extracted AAE and VAE features, the training time of
ContrastNet is considerably longer than other algorithms. However, when
extracting features, the time consumption of ContrastNet is fairly fast than
others. Because for AAE and VAE, only encoders are used in the feature
extraction phase. Moreover, the structure of AAE and VAE is relatively simple
and efficient, so the ContrastNet performs fast in the feature extraction
phase.

Table 12: Computational analysis of different feature extration algorithms.
The best value in each row is in bold. And the methods in bold are structures
proposed in this paper. Some of the data in the table is quoted from [3DCAE].
Indian Pines dataset
---
| 1D-CNN | S-CNN | TPCA | SSAE | EPLS | 3DCAE | ContrastNet
Training(s) | 20.6 | 119 | 44 | 240.1 | 111.9 | 1156 | 2073
Inference (s) | 17.4 | 3.2 | 150 | 2.76 | 47.3 | 5.22 | 12.77
Salinas dataset
| 1D-CNN | S-CNN | TPCA | SSAE | EPLS | 3DCAE | ContrastNet
Training(s) | 43.1 | 134 | 241 | 513 | 103.5 | 1159 | 4784
Inference (s) | 83.8 | 62 | 816 | 16.3 | 192.2 | 26.4 | 36.65
University of Pavia dataset
| 1D-CNN | S-CNN | TPCA | SSAE | EPLS | 3DCAE | ContrastNet
Training(s) | 32.1 | 332 | 44 | 491 | 127.13 | 1168 | 3778
Inference (s) | 150.07 | 106 | 150 | 31.9 | 141.5 | 32.04 | 27.82

##  5 Conclusion

In this paper, we proposed an unsupervised feature learning method based on
autoencoder and contrastive learning. This method combines unsupervised
representative methods and unsupervised discriminative methods, learning to
extract better features for hyperspectral classification than other baseline
methods. In the proposed method, we designed two efficient autoencoder
structure: VAE and AAE, and designed a model named ContrastNet for contrastive
learning, which reduces the computing resources demand of contrastive
learning. In fact, the VAE and AAE modules can be replaced by any other
feature extraction algorithms. Our experiments show that the proposed method
can extract more representative features and keep a high feature extraction
speed in the testing phase. Our work shows that unsupervised learning still
has great potential in the remote sensing field, and we hope others can get
more exciting ideas through our exploration.

## Acknowledgment

The work was supported in part by the Joint Fund of the Ministry of Education
of China under Grant 6141A02022362, and in part by the Project of Science and
Technology Plan of State Grid Jiaxing Power Supply Company under Grant
KM20180179.


# On Bounds for Greedy Schemes in String Optimization based on Greedy
Curvatures

Bowen Li, Brandon Van Over, Edwin K. P. Chong, and Ali Pezeshki This work is
supported in part by the AFOSR under award FA8750-20- 2-0504 and by the NSF
under award CNS-2229469.B. Li, B. Van Over, E. K. P. Chong and A. Pezeshki are
with the Department of Electrical and Computer Engineering, Colorado State
University, Fort Collins, CO 80523, USA. {bowen.li, b.van_over, edwin.chong,
<EMAIL_ADDRESS>
###### Abstract

We consider the celebrated bound introduced by Conforti and Cornuéjols (1984)
for greedy schemes in submodular optimization. The bound assumes a submodular
function defined on a collection of sets forming a matroid and is based on
greedy curvature. We show that the bound holds for a very general class of
string problems that includes maximizing submodular functions over set
matroids as a special case. We also derive a bound that is computable in the
sense that they depend only on quantities along the greedy trajectory. We
prove that our bound is superior to the greedy curvature bound of Conforti and
Cornuéjols. In addition, our bound holds under a condition that is weaker than
submodularity.

##  I INTRODUCTION

In many sequential decision-making or machine learning problems, we encounter
the problem of optimally choosing a _string_ (ordered set) of actions over a
finite horizon to maximize a given objective function. String optimization
problems have the added complexity relative to _set_ optimization problems in
that the objective function depends on both the set of actions taken and the
order in which the actions are taken. In such problems, determining the
optimal solution (optimal string of actions) can become computationally
intractable with increasing size of state/action space and optimization
horizon. Therefore, we often have to resort to approximate solutions. One of
the most common approximation schemes is the _greedy scheme_ , in which we
sequentially select the action that maximizes the increment in the objective
function at each step. A natural question that arises is how good is the
greedy solution to a _string optimization_ problem relative to the optimal
solution?

In this paper, we derive a _ratio bound_ for the performance of greedy
solutions to string optimization problems relative to that of the optional
solution. By a ratio bound, we mean a bound of the form

$\frac{f(G_{K})}{f(O_{K})}\geq\beta,$

where $f(G_{K})$ and $f(O_{K})$ are the objective function values of the
greedy solution and optimal solution, respectively, and $K$ is the horizon.
The bound guarantees that the greedy solution achieves at least a factor of
$\beta$ of the optimal solution. Our focus is on establishing a bound in which
the factor $\beta$ is easily computable.

Prior work: Bounding greedy solutions to set optimization problems has a rich
history, especially in the context of _submodular_ set functions. The most
celebrated results are from Fisher et al. [1] and Nemhauser et al. [2]. They
showed that $\beta\geq 1/2$ over a finite general set matroid [1] and
$\beta\geq 1-((K-1)/K)^{K}>(1-e^{-1})$ over a finite uniform set matroid with
$K$ being the horizon [2]. Improved bounds (with $\beta$ larger than
$1-e^{-1}$) have been developed by Conforti and Cornuéjols [3], Vondrak [4,
5], and Wang et al. [6], by introducing various notions of _curvature_ for set
submodular functions under uniform and/or general matroid settings.

However, for most of these bounds, computing the value of $\beta$ is
intractable and hence the bound is effectively not computable in problems with
large action spaces and decision horizons. In addition, the value of $\beta$
depends on the values of the objective function beyond the optimization
horizon. One notable exception is the bound derived by Conforti and Cornuéjols
[3] that is based on a quantity called greedy curvature. This bound is easily
computable. Another notable exception is the bound derived by Welikala et al.
[7], which is also computable and in some problems has shown to be larger than
other curvature bounds. However, the value of $\beta$ for this bound depends
on values of the objective function beyond the optimization horizon.

Recently, Zhang et al. [8] extended the concept of set submodularity to string
submodularity and established a $\beta=(1-e^{-1})$ bound for greedy solutions
of string optimization problems. They also derived improved bounds involving
various notions of curvature for string submodular functions. However, these
curvatures have the same computational intractability issue as mentioned
above. Alaei et al. [9] proved $(1-e^{-1})$ and $(1-1/e^{(1-1/e)})$ bounds for
online advertising and query rewriting problems, both of which can be
formulated under the framework of string submodular optimization.

Main contributions: The contributions of this paper are as follows:

  1. 1.

We first derive a computable bound for greedy solutions to string optimization
problems by extending the notion of greedy curvature and the proof technique
used by Conforti and Cornuéjols [3] from set functions to string functions.
Our bound coincides with the greedy curvature bound in [3] if reducing our
string objective function to a set objective function. However, in deriving
our bound, we don’t use submodularity. Rather we rely on weaker conditions on
the string objective function. Therefore, our bounding result is stronger.
This result is reported in Theorem 1.

  2. 2.

We then establish an even stronger result. We derive another computable bound
that relies on weaker assumptions than those used in deriving our first bound.
This is done in Theorem 2.

  3. 3.

We then show that the value $\beta_{2}$ of the bound in Theorem 2 is larger
than the value $\beta_{1}$ of the bound in Theorem 1 (i.e., the generalization
of the Conforti and Cornuéjols bound). This is shown in Theorem 3.

Organization: The paper is organized as follows. Section II introduces all the
mathematical preliminaries and formulates the string optimization problem.
Section III presents the mains results regarding our computable bounds.
Applications of our theoretical results to task scheduling and sensor coverage
problems are demonstrated in Section IV. Finally, conclusions are given in
Section V.

##  II DEFINITIONS FOR STRING OPTIMIZATION

In this section, we introduce notation, terminology, and definitions that we
use in our paper and present a formulation for a general string optimization
problem.

######  Definition 1.

Let $\mathbb{S}$ be a set. In our context, $\mathbb{S}$ is called the ground
set and its elements are called symbols.

  1. 1.

Let $s_{1},s_{2},\ldots,s_{k}\in\mathbb{S}$. Then $S:=s_{1}s_{2}\cdots s_{k}$
is a string with length $|S|=k$.

  2. 2.

Let $K$ be a positive integer, called the horizon.

  3. 3.

Let $\mathbb{S}_{K}$ be the set of all positive strings of length up to $K$,
including the empty string $\varnothing$. It is also called the uniform string
matroid of rank $K$.

######  Definition 2.

Consider two strings $S=s_{1}s_{2}\cdots s_{m}$ and $T=t_{1}t_{2}\cdots t_{n}$
in $\mathbb{S}_{K}$.

  1. 1.

Define the concatenation of $S$ and $T$ as $ST:=s_{1}s_{2}\cdots
s_{m}t_{1}t_{2}\cdots t_{n}$

  2. 2.

We say that $P$ is a prefix of $S$ if $S=PU$ for some $U\in\mathbb{S}_{K}$, in
which we write $P\preceq S$.

  3. 3.

We say that $\mathbb{T}\subseteq\mathbb{S}_{K}$ is prefix-closed if for all
$S\in\mathbb{T}$ and $P\preceq S$, $P\in\mathbb{T}$.

  4. 4.

Let $f:\mathbb{S}_{K}\rightarrow\mathbb{R}$ be the objective function.

######  Definition 3.

Let $\mathbb{X}\subseteq\mathbb{S}_{K}$. Then $(\mathbb{X},\mathbb{S})$ is a
finite rank $K$ string matroid if

  1. 1.

$|A|\leq K$ for all $A\in\mathbb{X}$.

  2. 2.

If $B\in\mathbb{X}$ and $A\preceq B$, then $A\in\mathbb{X}$.

  3. 3.

For every $A,B\in\mathbb{X}$ where $|A|+1=|B|$, there exists $b$ that is a
symbol in $B$ and $Ab\in\mathbb{X}$.

######  Definition 4.

$f:\mathbb{S}_{K}\rightarrow\mathbb{R}$ is string submodular on
$(\mathbb{X},\mathbb{S})$ if

  1. 1.

$\forall A\preceq B\in\mathbb{X}$, $f(A)\leq f(B)$.

  2. 2.

$\forall A\preceq B\in\mathbb{X}\text{ and }\forall j\in\mathbb{S}$ such that
$Aj,Bj\in\mathbb{X}$, $f(Aj)-f(A)\geq f(Bj)-f(B)$.

######  Remark 1.

  1. 1.

Unlike the permutation invariance in the set case, the order of the symbols in
a string matters. Different orders of the same set of symbols represent
different strings.

  2. 2.

Note that $ST\in\mathbb{S}_{K}$ if and only if $m+n\leq K$.

  3. 3.

In the matroid literature, a prefix-closed collection is also said to satisfy
the hereditary property and is called as independence system. An independence
system does not require condition $(3)$ in Definition 3.

######  Remark 2.

The definitions of finite rank string matroid and string submodular function
are introduced in [10]. Their formulations are inspired by the definitions of
set matroid and set submodular function in [3]. Most of the previous work is
established in the set case. Here, we extend the theoretical results to the
more general string case.

Consider an objective function $f$ and a prefix-closed $\mathbb{T}$ with
$\max_{S\in\mathbb{T}}|S|=K$. The general string-optimization is given by

$\displaystyle\text{maximize }f(S)$ (1) $\displaystyle\text{subject to
}S\in\mathbb{T}$

######  Remark 3.

  1. 1.

The stipulation that $f(\varnothing)=0$ is without loss of generality because
if $f(\varnothing)\neq 0$, we can subtract $f(\varnothing)$ from all values of
$f$ without changing the maximizer of the optimization problem (1).

  2. 2.

The constraint set $\mathbb{T}$ here is not necessarily the uniform matroid
$\mathbb{S}_{K}$ or a finite rank $K$ string matroid
$(\mathbb{X},\mathbb{S})$. In particular, our analysis applies to general
prefix-closed sets as constraint sets subject to certain assumptions.

  3. 3.

The stipulation that $\max_{S\in\mathbb{T}}|S|=K$ is without loss of
generality because we can always define $K:=\max_{S\in\mathbb{T}}|S|$.

######  Definition 5.

Any solution to the optimization problem (1) is said to be optimal, denoted by
$O_{L}=o_{1}o_{2}\cdots o_{L}$, where $L\in\\{1,\ldots,K\\}$.

######  Remark 4.

  1. 1.

The length $L$ of an optimal solution $O_{L}$ could be anything from $1$ to
$K$.

  2. 2.

There may be multiple optimal solutions.

  3. 3.

If $\mathbb{S}$ is finite, then $\mathbb{S}_{K}$ and $\mathbb{T}$ are finite,
and an optimal solution always exists.

######  Definition 6.

We define $G_{K}=g_{1}g_{2}\cdots g_{K}$ to be a greedy solution if for all
$k\in\\{1,2,\ldots,K\\}$,

$g_{k}=\operatorname*{arg\,max}_{s\in\mathbb{S}:g_{1}\cdots
g_{k-1}s\in\mathbb{T}}f(g_{1}\cdots g_{k-1}s).$

Note that we are going all the way to horizon $K$. This implicitly assumes
that $f$ is nondecreasing along $G_{K}$ and $f(g_{1})>0$. Theses implicit
conditions are reflected in assumptions $\mathbf{A_{1}}\text{ and
}\mathbf{A_{3}}$ in Section III.

######  Remark 5.

  1. 1.

The length of a greedy solution $G_{K}$ is always $K$.

  2. 2.

A greedy scheme is the one that adds a symbol to the existing string at each
time $k\text{ for }k\in\\{1,\ldots,K-1\\}$, so that the resulting string
produces the highest value of $f$ without regard to the future times. In other
words, a greedy solution is generated by a greedy scheme.

  3. 3.

The $\operatorname*{arg\,max}$ above could be nonunique, in which case there
would be multiple greedy solutions.

######  Definition 7.

We define $\beta:=f(G_{K})/f(O_{K})$ as the performance bound of the greedy
solution.

######  Remark 6.

$f(G_{K})$ can be computed exactly, but $f(O_{K})$ is computational
intractable. Finding a valid lower bound of $\beta$ is equivalent to finding
an upper bound of $f(O_{K})$.

##  III Main Results

In this section, we present our bounding results. We first introduce the
assumptions that our bounding results rely on. Then we generalize the greedy
curvature bound of Conforti and Cornuéjols [3] from set optimization problems
to string optimization problems. Then, we establish a stronger bounding result
that relies on fewer assumptions and has a larger bound value.

Given any string $S=s_{1}s_{2}\cdots s_{|S|}\in\mathbb{S}_{K}$ and
$i,j\in\\{1,2,\ldots,|S|\\}$, denote $S_{i:j}:=s_{i}\cdots s_{j}$ if $i\leq j$
and $S_{i:j}=\varnothing$ if $i>j$. For simplicity, we use the abbreviation
$S_{k}:=S_{1:k}\;(k\in\\{0,\ldots,|S|\\})$. Note that by definition,
$S_{0}=\varnothing$.

Let $S\in\mathbb{T}$. Define the $\Delta$ notation as follows: for each
$k\in\\{1,\ldots,|S|\\}$, $\Delta(S_{k}):=f(S_{k})-f(S_{k-1})$, called the
$k^{\text{th}}$ increment of $f(S)$. In [3], $\Delta(S_{k})$ is called a
discrete derivative.

For any $S\in\mathbb{T}$, let
$\mathbb{S}(S):=\\{s\in\mathbb{S}:Ss\in\mathbb{T}\\}$. The set $\mathbb{S}(S)$
contains the actions that are feasible with respect to $\mathbb{T}$ following
string $S$. The set
$\mathbb{S}(G_{k-1})=\\{s\in\mathbb{S}:G_{k-1}s\in\mathbb{T}\\}$ is frequently
referred to in our derivations. When $k=1$,
$\mathbb{S}(G_{k-1})=\mathbb{S}(\varnothing)$. For $k>1$,
$\mathbb{S}(G_{k-1})=\\{s\in\mathbb{S}:G_{k-1}s\in\mathbb{T}\\}$ is nonemepty
because $g_{k}\in\mathbb{S}(G_{k-1})$ by definition.

Assumptions: We introduce the following three key assumptions regarding
feasibility along the greedy sequence and the diminishing return property
along the optimal and greedy paths:

$\displaystyle\mathbf{A_{1}:}$ $\displaystyle\text{ For each
}k\in\\{1,\ldots,K\\},o_{k}\in\mathbb{S}(G_{k-1})\text{ and }$
$\displaystyle\text{ }\Delta(G_{k-1}s)>0\text{ for }s\in\mathbb{S}(G_{k-1});$
$\displaystyle\mathbf{A_{2}:}$ $\displaystyle\text{ For each
}k\in\\{1,\ldots,K\\},\Delta(O_{k})\leq f(o_{k});$
$\displaystyle\mathbf{A_{3}:}$ $\displaystyle\text{ For each
}k\in\\{1,\ldots,K\\},\Delta(G_{k})\leq f(g_{k}).$

Our bounding results (Theorem 1-3) rely on a subset or all of the the above
three assumptions. We note that $\mathbf{A_{2}}$ and $\mathbf{A_{3}}$ are
weaker conditions than string submodularity, they only involve diminishing
return along the optimal and greedy sequences instead of any sequence in
$\mathbb{T}$. Any string submodular function satisfies $\mathbf{A_{2}}$ and
$\mathbf{A_{3}}$, as stated in the following lemma.

######  Lemma 1.

If $f$ is string submodular, then $f$ satisfies $\mathbf{A_{2}}$ and
$\mathbf{A_{3}}$.

###### Proof.

By condition (1) of string submodularity in Definition 4, it is easy to obtain
$\Delta(G_{k})\geq 0\text{ for }k\in\\{1,\ldots,K\\}.$

By condition (2) of string submodularity in Definition 4, we can set
$A=\varnothing\text{ and }B=O_{k-1}\text{ for }k\in\\{1,\ldots,K\\}$. Then
$Ao_{k}=o_{k}$ and $Bo_{k}=O_{k}$ are both feasible in the domain of $f$ and
$f(O_{k})-f(O_{k-1})=\Delta(O_{k})\leq f(o_{k})\text{ for
}k\in\\{1,\ldots,K\\}$. Similarly, we can obtain $\Delta(G_{k})\leq
f(g_{k})\text{ for }k\in\\{1,\ldots,K\\}$. ∎

We now establish our first bounding result. We define the greedy curvature for
a string function as

$\displaystyle\alpha_{G}$
$\displaystyle:=\max_{k\in\\{2,\ldots,K\\}}\max_{s\in\\{s|\Delta(G_{k-1}s)>0,s\in\mathbb{S}(G_{k-1})\\}}\frac{f(s)}{\Delta(G_{k-1}s)}.$
(2)

We have the following theorem.

######  Theorem 1.

Assuming $\mathbf{A_{1}},\mathbf{A_{2}}\text{ and }\mathbf{A_{3}}$,
$f(G_{K})/f(O_{K})\geq\beta_{1}$, where

$\beta_{1}=\frac{1}{K}+\frac{1}{\alpha_{G}}*\frac{K-1}{K}.$

###### Proof.

Rewriting $\alpha_{G}$ according to its definition gives us:

$\displaystyle\alpha_{G}\geq\frac{f(s)}{\Delta(G_{k-1}s)}$

for every $k\in\\{2,\ldots,K\\}$ and
$s\in\\{s|\Delta(G_{k-1}s)>0,s\in\mathbb{S}(G_{k-1})\\}.$

By assumption $\mathbf{A_{1}}$,

$\displaystyle f(o_{k})\leq\alpha_{G}\Delta(G_{k-1}o_{k})\text{ for
}k\in\\{2,\ldots,K\\}.$ (3)

Moreover, by assumption $\mathbf{A_{2}}$,

$f(O_{K})=\sum_{k=1}^{K}\Delta(O_{k})\leq\sum_{k=1}^{K}f(o_{k})\stackrel{{\scriptstyle(a)}}{{\leq}}Kf(g_{1}),$
(4)

in which $(a)$ holds since $f(g_{1})$, by the definition of greedy solution,
dominates each term in $\sum_{k=1}^{K}f(o_{k})$.

Then,

$\displaystyle f(O_{K})$
$\displaystyle=\sum_{k=1}^{K}(f(O_{k})-f(O_{k-1}))=\sum_{k=1}^{K}\Delta(O_{k})$
(5)
$\displaystyle\stackrel{{\scriptstyle(b)}}{{\leq}}\sum_{k=1}^{K}f(o_{k})\stackrel{{\scriptstyle(c)}}{{\leq}}f(g_{1})+\sum_{k=2}^{K}\alpha_{G}\Delta(G_{k-1}o_{k})$
$\displaystyle\stackrel{{\scriptstyle(d)}}{{\leq}}f(g_{1})+\sum_{k=2}^{K}\alpha_{G}\Delta(G_{k})$
$\displaystyle=\alpha_{G}(f(g_{1})+\sum_{k=2}^{K}\Delta(G_{K}))+(1-\alpha_{G})f(g_{1})$
$\displaystyle=\alpha_{G}f(G_{K})+(1-\alpha_{G})f(g_{1}),$

where $(b)$ is due to assumption $\mathbf{A_{2}}$, $(c)$ is by the relation in
inequality $(\ref{beta_greedy_1})$, and $(d)$ is by the definition of $g_{k}$.

Dividing by $f(O_{K})$ on both sides of $(\ref{beta_greedy_3})$ yields:

$\displaystyle\alpha_{G}\frac{f(G_{K})}{f(O_{K})}$ $\displaystyle\geq
1-(1-\alpha_{G})\frac{f(g_{1})}{f(O_{K})}\stackrel{{\scriptstyle(d)}}{{\geq}}1-\frac{1-\alpha_{G}}{K}$
(6) $\displaystyle\frac{f(G_{K})}{f(O_{K})}$
$\displaystyle\geq\frac{1}{K}+\frac{1}{\alpha_{G}}*\frac{K-1}{K}.$

where $(d)$ is due to inequality $(\ref{beta_greedy_2})$ and the fact that, by
assumption $\mathbf{A_{3}}$, $\alpha_{G}\geq
f(g_{k})/\Delta(G_{k-1}g_{k})=f(g_{k})/\Delta(G_{k})\geq 1$ for any
$k\in\\{2,\ldots,K\\}$. ∎

We note that the above proof technique is similar to the proof technique used
in Theorem 3.1 in [3] by Conforti and Cornuéjols, in deriving their greedy
curvature bound for set optimization problems. However, Theorem 1 is more
general in that it applies to bounding greedy solutions in string optimization
problems, which subsume set optimization problems, and it does not rely on
having submodularity. The weaker assumptions A2 and A3 suffice. For a set
submodular function, the value of the bound in Theorem 1 coincides with the
value of the bound of Conforti and Cornuéjols. The result in Theorem 1 also
generalizes our prior result for bounding string optimizations [10].

We now establish a stronger result.

######  Theorem 2.

Assuming $\mathbf{A_{1}}$ and $\mathbf{A_{2}}$,
$f(G_{K})/f(O_{K})\geq\beta_{2}$, where

$\beta_{2}=\frac{f(G_{K})}{\sum_{k=1}^{K}\max_{s\in\mathbb{S}(G_{k-1})}f(s)}.$
(7)

###### Proof.

$\displaystyle
f(O_{K})=\sum_{k=1}^{K}(f(O_{k})-f(O_{k-1}))=\sum_{k=1}^{K}\Delta(O_{k})$

By $\mathbf{A_{2}}$, $\Delta(O_{k})\leq f(o_{k})$ for $k=1,\ldots,K$.
Therefore,

$\sum_{k=1}^{K}\Delta(O_{k})\stackrel{{\scriptstyle(a)}}{{\leq}}\sum_{k=1}^{K}f(o_{k})\stackrel{{\scriptstyle(b)}}{{\leq}}\sum_{k=1}^{K}\max_{s\in\mathbb{S}(G_{k-1})}f(s),$

where $(a)$ is by $\mathbf{A_{2}}$ and $(b)$ is by $\mathbf{A_{1}}$ and the
definition of the $\max$ operator. ∎

Theorem 2 essentially shows that the sum of the objective values under the
first $K$ greedy and feasible actions serves a valid upper bound for
$f(O_{K})$. This upper bound is easy to compute along the greedy trajectory
and is obtained under weaker assumptions than those in Theorem 1.

We now show that the bound in Theorem 2 is better than the bound in Theorem 1.

######  Theorem 3.

Assuming $\mathbf{A_{1}},\mathbf{A_{2}}\text{ and }\mathbf{A_{3}}$,
$\beta_{2}\geq\beta_{1}$.

###### Proof.

Recall the improved greedy curvature proposed in [10]:

$\alpha_{k}:=\max_{s\in\\{s|\Delta(G_{k-1}s)>0,s\in\mathbb{S}(G_{k-1})\\}}f(s)/\Delta(G_{k-1}s)$

for $k\in\\{2,\ldots,K\\}$. By definition of $\alpha_{k}$ and $\alpha_{G}$, it
is not hard to see that:

$1\leq\alpha_{k}\leq\alpha_{G}\text{ for }k\in\\{2,\ldots,K\\},$ (8)

because the definition of $\alpha_{G}$ includes an additional $\max$ operator
with respect to all $k\in\\{2,\ldots,K\\}$ while $\alpha_{k}$ only includes
the $\max$ operator with respect to $s$.

Rewriting the definition of $\alpha_{k}$ gives:

$f(s)\leq\alpha_{k}\Delta(G_{k-1}s)\leq\alpha_{k}\Delta(G_{k})$ (9)

for every $k\in\\{2,\ldots,K\\}$ and
$s\in\\{s|\Delta(G_{k-1}s)>0,s\in\mathbb{S}(G_{k-1})\\}.$

Hence,

$\displaystyle\sum_{k=1}^{K}\max_{s\in\mathbb{S}(G_{k-1})}f(s)$
$\displaystyle\stackrel{{\scriptstyle(a)}}{{\leq}}f(g_{1})+\sum_{k=2}^{K}\alpha_{k}\Delta(G_{k})$
(10)
$\displaystyle\stackrel{{\scriptstyle(b)}}{{\leq}}\alpha_{G}f(G_{K})+(1-\alpha_{G})f(g_{1})$
$\displaystyle\stackrel{{\scriptstyle(c)}}{{\leq}}\alpha_{G}f(G_{K})+(1-\alpha_{G})\frac{f(O_{K})}{K}$

where $(a)$ is by inequality $(\ref{comparison_2})$, $(b)$ is by inequality
$(\ref{comparison_1})$, and $(c)$ is due to inequality $(\ref{beta_greedy_2})$
and the fact that $\alpha_{G}\geq 1$.

Dividing $\sum_{k=1}^{K}\max_{s\in\mathbb{S}(G_{k-1})}f(s)$ on both sides of
$(\ref{comparison_3})$ yields:

$\displaystyle\frac{\alpha_{G}f(G_{K})}{\sum_{k=1}^{K}\max_{s\in\mathbb{S}(G_{k-1})}f(s)}=\alpha_{G}\beta_{\text{
Univ }}$ $\displaystyle\geq
1-(1-\alpha_{G})\frac{f(O_{K})}{K\sum_{k=1}^{K}\max_{s\in\mathbb{S}(G_{k-1})}}$
$\displaystyle\geq 1-\frac{(1-\alpha_{G})}{K},$

which further implies

$\beta_{2}\geq\frac{1}{K}+\frac{1}{\alpha_{G}}*\frac{K-1}{K}=\beta_{1}.$

∎

##  IV APPLICATIONS

###  IV-A Task Scheduling

As a canonical optimization problem in operations research, task scheduling
was studied in [11] and further investigated in [8] and [12]. In task
scheduling, we aim to assign agents at multiple stages throughout the task to
maximize the probability of successful completion. The detailed mathematical
formulation is as follows.

Assume a task is comprised of $K$ stages, and an agent needs to be assigned at
stage $k\;(1\leq k\leq K)$ to accomplish the task. We have a pool of $N$
agents, denoted by $M_{1},\ldots,M_{N}$, and no agent is allowed to be
repeatedly selected. For each agent $M_{i}\;(1\leq i\leq N)$, the probability
of accomplishing the task at each stage $k\;(1\leq k\leq K)$ is denoted by
$P(m_{i}^{k})$. Therefore, the probability of accomplishment after selecting
$M_{i_{1}},M_{i_{2}},\ldots,M_{i_{K}}$ in the end is:

$f(M_{i_{1}}M_{i_{2}}\cdots
M_{i_{K}}):=1-\prod_{k=1}^{K}(1-p(m_{i_{k}}^{k})),$ (11)

where the subscript $i_{*}$ in $M$ and $m$ represents the selected agent
number. Note that the probability of accomplishing the task for each agent is
stage-variant. $f(M_{i_{1}}M_{i_{2}}\cdots M_{i_{K}})$ in
$(\ref{task_scheduling})$ is henceforth a string function and order dependent.
Streeter and Golovin proved in [11] that the greedy schemes yields a
$\beta_{0}=(1-e^{-1})$ performance bound when $f$ bears submodularity and
monotonicity.

Consider the following example with $K=3$, $N=5$, and the probabilities of
successful completion in Table. I. The optimal sequence $(O_{3})$ and greedy
sequence $(G_{3})$ are equal in this example with
$G_{3}=O_{3}=M_{1}M_{2}M_{3}$. We can verify that assumption $\mathbf{A_{1}}$
is satisfied, and $f$ is string submodular. The three performance bounds are
as follows:

$\displaystyle\beta_{2}=\frac{f(M_{1}M_{2}M_{3})}{f(M_{1})+f(M_{2})+f(M_{3})}=0.7816;$
$\displaystyle\beta_{1}=\frac{1}{K}+\frac{1}{\alpha_{G}}*\frac{K-1}{K}=0.5893;$
$\displaystyle\beta_{0}=1-e^{-1}=0.6321.$

We observe that the $\beta_{2}$ bound outperforms the other two bounds. This
example shows the curvature bounds are not necessarily better than the
original $(1-e^{-1})$ bound. On some problems they are better, in others they
are not. Therefore, we need to choose the largest of the three bounds. But
this is easy, because all three bounds are easily compuatble. In the next
example, we show a plot that compares the three bounds under different
parameters of another optimization problem.

|  Stage 1 |  Stage 2 |  Stage 3
---|---|---|---
$M_{1}$ |  0.2 |  0.16 |  0.14
$M_{2}$ |  0.18 |  0.16 |  0.14
$M_{3}$ |  0.16 |  0.14 |  0.14
$M_{4}$ |  0.14 |  0.12 |  0.10
$M_{5}$ |  0.12 |  0.1 |  0.08
TABLE I: Probability of successful completion for each agent at each stage.

###  IV-B Multi-agent Sensor Coverage

The multi-agent sensor coverage problem was originally studied in [13] and
further analyzed in [14] and [7]. In a given mission space, we need to find a
placement of a set of homogeneous sensors to maximize the probability of
detecting randomly occurring events. This problem can be formulated as a set
optimization problem, which is a special case in the string setting. In this
subsection, we demonstrate our theoretical results by applying them to a
discrete version of the coverage problem. Our simplified version can be easily
generalized to more complicated settings within the same theoretical
framework.

The mission space $\Omega\subset\mathbb{R}^{2}$ is modeled as a non-self-
intersecting polygon where $K$ homogeneous sensors will be placed to detect a
randomly occurring event in $\Omega$. For simplicity of calculation, we assume
both the sensors and the random event can only be placed and occur at lattice
points. We denote the feasible space for sensor placement and event occurrence
as $\Omega^{F}$. Our goal is to maximize the overall likelihood of successful
detection in the mission space, as illustrated in Fig. 1.

Figure 1: Sensor Coverage in a Mission Space

The likelihood of event occurrence over $\Omega^{F}$ is given by an event mass
function $R:\Omega^{F}\xrightarrow{}\mathbb{R}_{\geq 0}$, and we assume that
$\sum_{\mathbf{x}\in\Omega^{F}}R(\mathbf{x})<\infty$. $R(\mathbf{x})$ may
reflect a particular distribution if some prior information is available.
Otherwise, $R(\mathbf{x})=1$ when no prior information is obtained. The
locations of all the sensors are represented as
$\mathbf{s}=(\mathbf{s}_{1},\mathbf{s}_{2},\ldots,\mathbf{s}_{K})\in(\Omega^{F})^{K}$,
where $\mathbf{s}_{i}\;(1\leq i\leq K)$ are coordinates of the placed sensors.
Each sensor placed at $\mathbf{s_{i}}$ can detect any occurring event at
location $\mathbf{x}$ with probability
$p(\mathbf{x},\mathbf{s}_{i})=e^{-\lambda\|\mathbf{x}-\mathbf{s}_{i}\|}$,
where $\lambda$ is the decay rate characterizing how quick the sensing
capability decays along the distance.

Assuming all the sensors are working independently, the probability of
detecting an occurring event at location $\mathbf{x}\in\Omega^{F}$ after
placing $K$ homogeneous sensors at locations $\mathbf{s}$ is
$p(\mathbf{x},\mathbf{s})=1-\prod_{i=1}^{K}\left(1-p(\mathbf{x},\mathbf{s}_{i})\right)$.
In order to consider the whole feasible space for event occurrence, we need to
incorporate the event mass function. Our objective function becomes
$H(\mathbf{s})=\sum_{\mathbf{x}\in\Omega^{F}}R(\mathbf{x})p(\mathbf{x},\mathbf{s})$.
The goal is to maximize $H(\mathbf{s})$:

$\displaystyle\text{maximize }H(\mathbf{s})$ (12) $\displaystyle\text{subject
to }\mathbf{s}\in(\Omega^{F})^{K}.$

If $n$ lattice points in $\Omega^{F}$ are feasible for sensor placement, we
need to select $K$ out of $n$ locations with its complexity being
$n!/\left(K!(n-K)!\right)$. This becomes a set optimization problem, and
exhaustive search is computationally intractable when $n$ is large. Therefore,
greedy algorithm is an approach for an approximate solution in polynomial
time. It was proved that the continuous version of $H(\mathbf{s})$ is
submodular in [14], and it is not difficult to verify that its discrete
version is also submodular.

In our experiment, we consider a rectangular mission space $\Omega$ of size
$40\times 30$, and $K=4$ homogeneous sensors are waiting to be deployed on
those integer coordinates within $\Omega$, denoted by $\Omega^{F}$. For a
point $\mathbf{p}=(s_{x},s_{y})$, the event mass function is given by
$R(\mathbf{p})=\left(s_{x}+s_{y}\right)/\left(s_{x\text{ max}}+y_{y\text{
max}}\right)$, where $s_{x\text{ max}}=40$ and $s_{y\text{ max}}=30$ are the
largest values of the $x$ and $y$ coordinates in the mission space. This event
mass function implies that the random events are more likely to occur around
the top right corner of the rectangular mission space.

A comparison of the performance bounds under different decay rate is shown in
Fig. 2. Small decay rates imply good sensing capability and strong
submodularity, under which the greedy scheme does not produce nearly optimal
objective function value. In Fig. 2, the $\beta_{2}$ bound (red line) always
exceeds the $\beta_{1}$ bound (blue line), illustrating Theorem 3. In
addition, we can observe instances in which the $\beta_{2}$ bound is larger
than the $\beta_{0}=(1-e^{-1})$, while the $\beta_{1}$ bound is below this
value.

Figure 2: Performance bound comparison under different decay rates

##  V CONCLUSION AND FUTURE WORK

We derived a computable bound $\beta_{1}$ for greedy solutions to string
optimization problems by extending the notion of greedy curvature and the
proof technique used by Conforti and Cornuéjols [3] from set functions to
string functions. However, in deriving our bound we did not use submodularity.
Rather we relied on weaker conditions on the string objective function, and
hence produced a stronger bounding result. We then derived another computable
bound $\beta_{2}$ that relies on even weaker assumptions than those used in
deriving our first bound $\beta_{1}$, further strengthening the bounding
result. We also showed that this second bound $\beta_{2}$ has a larger value
than the first bound $\beta_{2}$. We demonstrated the superiority of our bound
on two applications, task scheduling and multi-agent sensor coverage.

## References

  * [1] M. L. Fisher, G. L. Nemhauser, and L. A. Wolsey, _An analysis of approximations for maximizing submodular set functions—II_. Springer Berlin Heidelberg, 1978.
  * [2] G. L. Nemhauser, L. A. Wolsey, and M. L. Fisher, “An analysis of approximations for maximizing submodular set functions—I,” _Math. Program._ , vol. 14, pp. 265–294, 1978.
  * [3] M. Conforti and G. Cornuéjols, “Submodular set functions, matroids and the greedy algorithm: Tight worst-case bounds and some generalizations of the Rado-Edmonds theorem,” _Discrete. Appl. Math._ , vol. 7, no. 3, pp. 251–274, 1984.
  * [4] J. Vondrák, “Submodularity and curvature: The optimal algorithm (combinatorial optimization and discrete algorithms),” _RIMS Kokyuroku Bessatsu B_ , vol. 23, pp. 253–266, 2010.
  * [5] G. Calinescu, C. Chekuri, M. Pal, and J. Vondrák, “Maximizing a monotone submodular function subject to a matroid constraint,” _SIAM J. Comput._ , vol. 40, no. 6, pp. 1740–1766, 2011.
  * [6] Z. Wang, B. Moran, X. Wang, and Q. Pan, “Approximation for maximizing monotone non-decreasing set functions with a greedy method,” _J. Comb. Optim._ , vol. 31, pp. 29–43, 2016.
  * [7] S. Welikala, C. G. Cassandras, H. Lin, and P. J. Antsaklis, “A new performance bound for submodular maximization problems and its application to multi-agent optimal coverage problems,” _Automatica_ , vol. 144, p. 110493, 2022.
  * [8] Z. Zhang, E. K. P. Chong, A. Pezeshki, and W. Moran, “String submodular functions with curvature constraints,” _IEEE Trans. Autom. Control_ , vol. 61, no. 3, pp. 601–616, 2015.
  * [9] S. Alaei, A. Makhdoumi, and A. Malekian, “Maximizing sequence-submodular functions and its application to online advertising,” _Manage. Sci._ , vol. 67, no. 10, pp. 6030–6054, 2021.
  * [10] B. Van Over, B. Li, E. K. Chong, and A. Pezeshki, “An improved greedy curvature bound in finite-horizon string optimization with an application to a sensor coverage problem,” in _Proc. 62nd IEEE Conf. Decision Control (CDC)_ , Singapore, Dec. 2023, pp. 1257–1262.
  * [11] M. Streeter and D. Golovin, “An online algorithm for maximizing submodular functions,” in _Proc. Adv. in Neural Inf. Process. Syst. 21 (NIPS 2008)_ , vol. 21, Vancouver, BC, Canada, Dec. 2008, pp. 1577–1584.
  * [12] Y. Liu, Z. Zhang, E. K. P. Chong, and A. Pezeshki, “Performance bounds with curvature for batched greedy optimization,” _J. Optim. Theory Appl._ , vol. 177, pp. 535–562, 2018.
  * [13] M. Zhong and C. G. Cassandras, “Distributed coverage control and data collection with mobile sensor networks,” _IEEE Trans. Autom. Control_ , vol. 56, no. 10, pp. 2445–2455, 2011.
  * [14] X. Sun, C. G. Cassandras, and X. Meng, “Exploiting submodularity to quantify near-optimality in multi-agent coverage problems,” _Automatica_ , vol. 100, pp. 349–359, 2019.