Date of publication xxxx 00, 0000, date of current version xxxx 00, 0000.
10.1109/TQE.2020.DOI

Corresponding author: Soronzonbold Otgonbaatar (email:
soronzonbold.otgonbaatar@dlr.de).

# Exploiting the Quantum Advantage for Satellite Image Processing: Review and
Assessment

SORONZONBOLD OTGONBAATAR1,2  DIETER KRANZLMÜLLER2  Remote Sensing Technology
Institute, German Aerospace Center DLR, Weßling, 82234, Germany (email:
<EMAIL_ADDRESS> Institut für Informatik, Ludwig-Maximilians-
Universität München, Oettingenstr. 67, München, 80538, Germany (email:
<EMAIL_ADDRESS>
###### Abstract

This article examines the current status of quantum computing in Earth
observation (EO) and satellite imagery. We analyze the potential limitations
and applications of quantum learning models when dealing with satellite data,
considering the persistent challenges of profiting from quantum advantage and
finding the optimal sharing between high-performance computing (HPC) and
quantum computing (QC). We then assess some parameterized quantum circuit
models transpiled into a Clifford+T universal gate set. The T-gates shed light
on the quantum resources required to deploy quantum models, either on an HPC
system or several QC systems. In particular, if the T-gates cannot be
simulated efficiently on an HPC system, we can apply a quantum computer and
its computational power over conventional techniques. Our quantum resource
estimation showed that quantum machine learning (QML) models, with a
sufficient number of T-gates, provide the quantum advantage if and only if
they generalize on unseen data points better than their classical counterparts
deployed on the HPC system and they break the symmetry in their weights at
each learning iteration like in conventional deep neural networks. We also
estimated the quantum resources required for some QML models as an initial
innovation. Lastly, we defined the optimal sharing between an HPC+QC system
for executing QML models for hyperspectral satellite images. These are a
unique dataset compared to other satellite images since they have a limited
number of input qubits and a small number of labeled benchmark images, making
them less challenging to deploy on quantum computers.

###### Index Terms:

Earth observation, hyperspectral images, image classification, quantum machine
learning, quantum computers, quantum resource estimation, remote sensing.

=-15pt

##  I Introduction

###  I-A Why quantum computing for Earth observation?

Earth observation (EO) methodologies tackle optimization and artificial
intelligence (AI) problems involving big datasets obtained from instruments
mounted on space-borne and airborne platforms. Some optimization and AI
problems combined with big EO datasets are intractable computational problems
for conventional high-performance computing (HPC) systems. In addition, EO
datasets themselves are complex heterogeneous image datasets, compared with
conventional red-green-blue (RGB) images, characterized by so-called 4V
features comprising _volume_ , _variety_ , _velocity_ , and _veracity_ [1];
here, _volume_ refers to big EO datasets (e.g., Terabytes of data per day
collected, for instance, by the European Space Agency), _variety_ refers to
distinct spectral data such as multispectral, and hyperspectral pixel data,
_velocity_ refers to the speed of change on the Earth’s surface, and
_veracity_ refers to imperfect datasets such as noisy images or remotely-
sensed images partly covered by clouds [2]. In general, EO problems also
include calibration and integer optimization problems in synthetic aperture
radar (SAR) applications [3, 4], a Bayesian paradigm (e.g., Gaussian
processes) for retrieving physical parameters from remotely-sensed datasets
[5, 6], uncertainty estimates for EO predictions [7], solving partial
differential equations (PDEs) for climate modeling and digital twin Earth
paradigms [8], and identifying objects on the Earth’s surface [9].
Furthermore, some computational problems like AI training architectures are
computationally expensive and inherently intractable problems or
$\mathbf{NP}$-$\mathbf{hard}$ problems (see Fig. I-B) [10]; Non-deterministic
($\mathbf{NP}$) polynomial problems are computational problems where there are
no known efficient commonly-used algorithms for finding their solutions in a
reasonable polynomial time (_i.e._ a polynomial number of steps) but can be
verified in a polynomial time given their solutions, and
$\mathbf{NP}$-$\mathbf{hard}$ problems are computational problems harder than
$\mathbf{NP}$ problems. On the other hand, quantum machines harnessing quantum
physics phenomena like entanglement can solve some challenging problems faster
and more efficiently than their counterpart conventional machines ranging from
integer optimization problems [11, 12, 13] to AI techniques [14, 15, 16, 17,
18] and PDEs, [19, 20], and even quantum-inspired algorithms for solving PDEs
[21]. Thus, quantum algorithms’ computational advantages (or quantum
advantage) over conventional algorithms inspire enough to examine and identify
computationally intractable problems with EO methodologies and hard EO
datasets for near- and far-term quantum machines.

###  I-B Do we really need quantum machines?

Quantum machines can be generally divided into three families comprising
quantum annealers [22], quantum simulators [23, 24], and universal quantum
computers [25]. These quantum machines promise computational advantage for
computing notoriously difficult problems over conventional computers according
to computational complexity theorems/conjectures [26, 27]; computational
complexity theorems draw boundaries between computational problems according
to their hardness for finding their solutions (see Fig. I-B) [10]. At the
moment, quantum machines are designed to tackle specific forms and kinds of
intractable computational problems, e.g., quantum annealers for quadratic
unconstrained binary optimization (QUBO) problems or simulating the Ising
Hamiltonian [11], and quantum simulators for mimicking some physical
Hamiltonian [28, 29]. Research communities ranging from high energy physics
[24], condensed-matter physics [29], AI [15] to EO [30] are in the exploration
phase of identifying and investigating their hard problems for quantum
platforms. Furthermore, classical computational methods for intractable
computational problems reach their limitations and potential accuracy due to
the classical computational resource required and the complexity of both EO
challenges and datasets. As stated earlier, some computational techniques are
intractable problems on conventional machines and computationally expensive,
even on the HPC system. Thus, to go beyond current computational methods
integrated with large-scale datasets to find a better solution and utilize low
computational cost, it is inevitable to examine and identify computationally
demanding problems in EO applications for novel near- and long-term quantum
machines. More importantly, gaining insight into programming these novel
computing machines and their potential advantages and imperfections for
computational problems is vital.

[t!](topskip=0pt, botskip=0pt, midskip=0pt)[width=0.8]figures/soron1.pdf The
computational complexity conjecture draws boundaries between computational
problems according to their hardness based on the required classical and
quantum computational resources. In particular, the computational problem
denoted by the green star is easy to solve for both quantum machines and
classical computers, the computational problem denoted by the orange star is
easy for quantum machines but hard for classical machines, and the
computational problem denoted by the black star is hard for classical
computers. Still, no known efficient quantum algorithmic approaches exist for
quantum machines.

###  I-C State of the art of quantum computing for Earth observation

Quantum computing is a novel computing paradigm that promises to find
solutions to some intractable computational problems more efficiently and
faster by exploiting quantum superposition and entanglement than conventional
computing techniques if and only if one considers ideal quantum complexity
measures without overhead considerations like a distillation of Toffoli gates
in the real quantum machines, e.g., the classical versions of the Toffoli
gates are transistors in a conventional computer [31]. Quantum machines are a
kind of computer constructed using the primitives of a quantum computing
method, such as quantum bits (qubits) and quantum gates, in contrast to
traditional classical bits and transistors. Digital quantum machines can be
decomposed into three layers [32]:

  1. 1.

a _quantum state preparation_ or a _quantum data encoding_ layer,

  2. 2.

a _quantum unitary evolution_ or a _parametrized quantum gate_ layer,

  3. 3.

a _quantum measurement_ layer.

For gaining insight into computing EO problems involving big datasets on
quantum machines, some studies already exist for processing a _variety_ of EO
datasets to tackle EO challenges using hybrid classical-quantum approaches
(see Fig. I-C); hybrid classical-quantum approaches are exchangeable with
quantum artificial intelligence (QAI) and quantum machine learning (QML). A
_variety_ of datasets includes hyperspectral, multispectral, and polarimetric
EO images.

[t!](topskip=0pt, botskip=0pt, midskip=0pt)[width=]figures/soron2.pdf A hybrid
classical-quantum approach for computational and machine learning tasks. A
quantum layer includes implicitly quantum data encoding, parametrized quantum
gates, and quantum measurement layers.

####  I-C1 Earth observation images

We can generalize that EO images are third-order tensors regardless of a
_variety_. Furthermore, a hyperspectral image is a remotely-sensed image
denoted by $\mathbb{R}^{I\times J\times K}$ where $I$ and $J$ are its spatial
dimensionality, and $K$ means hundreds of its narrow-spaced spectral bands (or
features), e.g., the Pavia University hyperspectral image described by
$\mathbb{R}^{610\times 340\times 103}$ tensor. Multispectral images are a
third-order tensor $\mathbb{R}^{I\times J\times K}$ with at most $K=12$
spectral bands. The main difference between them is the spectral bands’ number
and spacing. In contrast, polarimetric images are characterized by the
scattering property $S$ of ground targets; each pixel is described by a
$2\times 2$ scattering matrix but not by spectral bands as in hyperspectral
and multispectral images. Hence, we could assume that polarimetric images have
$K=3$ informative features if the scattering matrix is symmetric and $K=4$
informative features otherwise (see Figure 1) [33].

####  I-C2 Quantum machine learning for Earth observation images

Climate AI tasks involve analyzing satellite images that consist of thousands
of pixels and hundreds of spectral bands. For example, Eurosat multispectral
images have a size of $64\times 64$ pixels and $12$ spectral bands, which can
be represented as $\mathbb{R}^{64\times 64\times 12}$ [34]. In contrast, the
digital quantum machines currently available on the market have less than a
hundred noisy qubits and around depth-five of faulty quantum gates [35].
Therefore, the main challenge is to embed satellite images in a _quantum data
encoding_ layer, regardless of the size of quantum machines and their quantum
errors. To address this challenge, the authors of the articles [36, 37, 38,
39, 40] proposed and utilized a two-level embedding scheme. This scheme
comprises a classical layer for dimensionality reduction and a _quantum data
encoding_ layer for dimensionally-reduced images. In other words, they used a
hybrid classical-quantum approach, embedding classical datasets in a _quantum
data encoding_ layer and optimizing a _parametrized quantum gate_ layer of
digital quantum computers with the help of a conventional classical computer.
However, the Eurosat dataset they used is a large dataset consisting of low-
dimensional and easy-to-classify images and thus has low _veracity_. Most EO
datasets, on the other hand, are small datasets containing high-dimensional
and hard-to-classify images or high _veracity_ images. For example, the
multispectral UC Merced Land Use dataset has a size of $245\times 245$ pixels
and $3$ spectral bands, which can be represented as $\mathbb{R}^{245\times
245\times 3}$ [41]. To investigate the performance of quantum machines with
varying depths of a _parametrized quantum gate_ layer, the authors of the
article [42] utilized this dataset and polarimetric EO images for natural
embedding in input qubits without a dimensionality reduction technique [43].
It is important to note that the quality of the datasets used plays a crucial
role in data-driven tasks for hybrid classical-quantum approaches [44].
Therefore, the article’s authors [45] analyzed the power of EO image datasets
for training digital quantum machines.

Furthermore, a quantum annealer is a type of quantum simulator that is
designed to simulate an Ising Hamiltonian equivalent to QUBO problems [22]. In
recent articles’ authors [46, 47], they analyzed classification problems posed
as QUBO problems, belonging to $\mathbf{NP}$-$\mathbf{hard}$ problems, on a
D-Wave quantum annealer. They employed binary hyperspectral EO images since a
D-Wave quantum annealer promises to converge to a better solution to
$\mathbf{NP}$-$\mathbf{hard}$ problems. Some studies also transformed a
support vector machine (SVM) into a QUBO problem [48] and optimized it on a
D-Wave quantum annealer when analyzing EO image datasets [49, 50, 33]. To
embed large EO datasets in a D-Wave quantum annealer, the authors of [51] used
a K-fold technique and the concept of a coreset since a D-Wave quantum
annealer has around $5,000$ qubits arranged according to an expressly limited
topology. A D-Wave quantum annealer was also proposed for a notoriously hard
feature selection task and multi-label SVM for remotely-sensed hyperspectral
images [52].

Lastly, quantum-inspired algorithms are becoming increasingly popular in both
academic and industrial circles due to their energy and computational
efficiency. These algorithms are inspired by the potential advantages of
quantum algorithms, such as the quantum-inspired quantum Fourier
transformation [53], quantum-inspired AI/ML [54], and the use of tensor
networks to compress deep neural networks (DNNs) [55]. Tensor networks are
designed to compute quantum many-body systems efficiently [56], and they are
currently being used to simulate quantum circuits on modern GPU tensor cores
[57]. Thanks to these advancements, quantum tensor networks have been
successfully utilized to decrease the weights of physics-informed neural
networks (PINNs) and increase the resolution of hyperspectral images [58].

Figure 1: [Top] example hyperspectral, multispectral, and polarimetric images,
[Bottom Left] their third-order tensor representation, and [Bottom Right] each
pixel/target in polarimetric images is characterized by the complex numbered
scattering matrix in contrast to hyperspectral and multispectral images. Here,
$s_{ij}$ denotes a scattering element given sent/reflected horizontal $H$ or
vertical $V$ polarized beam.

####  I-C3 Selecting Earth observation data for quantum machines

When working with quantum machines in EO challenges, it is vital to choose
remotely sensed datasets based on the principle that “the more features in the
dataset, the less quantum resources required.” Studies have shown that
processing multispectral images requires more quantum gates and qubits than
hyperspectral and polarimetric images [36, 43]. This is because multispectral
images need global feature capturing, with each pixel dependent on its
neighbors, making processing more resource-intensive. On the other hand,
hyperspectral and polarimetric images contain informative spectral information
for each pixel. They can be embedded in qubits without the constraint of their
neighbors, making processing less resource-intensive [43]. For instance, one
QML model known as a quantum convolutional neural network (QCNN) requires
approximately $4,000$ quantum gates only to embed the element
$\mathbb{R}^{64\times 64\times 12}$ in the Eurosat dataset and roughly
$60,000$ quantum gates for embedding the multispectral image
$\mathbb{R}^{300\times 290\times 3}$ illustrated in Figure 1 in the input
qubits [59]. Hence, multispectral images are not viable for deploying QCNNs on
today’s quantum machines, even on future quantum machines. However, a hybrid
classical-quantum model requires fewer quantum resources than QCNNs. The
authors of the article [36] used only $16$ quantum gates for encoding the
Eurosat and the multispectral image $\mathbb{R}^{300\times 290\times 3}$
depending on the compressing quality. In contrast, we can embed the pixels of
a hyperspectral image, e.g., the Pavia University hyperspectral image, in the
input qubits using only at least three and, at most, about $103$ quantum gates
thanks to their abundant spectral bands [33]. As for polarimetric images, we
need at most five quantum gates due to their doppelgänger feature to qubits or
the one-to-one mapping between polarimetric images and qubits [43].

Based on the above analysis, hyperspectral satellite images are much more
appropriate for designing and assessing QML models and tackling climate
challenges than multispectral and polarimetric images since they have abundant
spectral information and fewer quantum resources required than other remotely
sensed datasets. More importantly, QML models generalize better on small-scale
datasets than their classical alternatives [60], whereas a hyperspectral
dataset has limited labeled images (or small-scale datasets) compared to
multispectral datasets and has more features than both multispectral and
polarimetric datasets.

###  I-D How and when do quantum machines outperform conventional computers?

It is becoming increasingly clear that quantum processing units (QPUs) will
soon be working alongside conventional classical computers, like how central
processing units (CPUs) and general processing units (GPUs) are used in
heterogeneous computing [30]. We are currently in the era of high-performance
computing (HPC), and the emergence of quantum computing (QC) is a new and
exciting concept in heterogeneous computing. It involves integrating a CPU+GPU
with QPUs designed to handle specific computational problems (see Fig. I-D).
For instance, a quantum annealer is designed to tackle only QUBO problems, and
neutral atom platforms can simulate certain chemical Hamiltonians. Depending
on the difficulty level of the computational problems, we may need to program
a challenging heterogeneous computing environment (i.e., CPU+GPU+QPUs) or a
conventional one (i.e., CPU+GPU).

[t!](topskip=0pt, botskip=0pt, midskip=0pt)[width=0.9]figures/soron4.pdf Novel
heterogeneous computing: a high performance and quantum computing paradigm.
Here, conventional heterogeneous computing refers to the programming of CPU
and GPU, whereas we call novel heterogeneous computing when integrating QPUs
with CPUs and GPUs. QPUs can be several parallel quantum machines based on
different quantum technologies such as quantum annealing, neutral atoms,
superconducting, and photonic.

QPUs, except for quantum annealers, currently consist of around 100 error-
prone qubits and low-depth, faulty quantum gates. The authors of the article
[61] coined these devices as ”noisy intermediate-scale quantum (NISQ)
devices.” However, for practical computational problems, there is no
demonstration of the computational advantage of NISQ devices over a
conventional classical computer. Therefore, estimating the quantum resources
required to tackle hard computational and ML problems is vital to achieving a
quantum advantage in EO. It is worth noting that some quantum algorithms can
be simulated efficiently using a conventional classical computer. For this
reason, any reasonable quantum resource estimation of a quantum algorithm
should consider non-Clifford T-gates, error rates of qubits and quantum gates,
and the execution time of single- and two-qubit quantum gates [62].

Figure 2: A visual representation of traditional NNs.

Non-Clifford T-gates are the most resource-expensive part of implementing a
quantum algorithm, compared to Clifford quantum gates or CNOT, Hadamard,
Phase, and measurement gates. Even the Gottesman-Knill theorem states
(informally) that non-Clifford T-gates cannot be efficiently simulated on a
conventional classical computer. In contrast, Clifford quantum gates can be
simulated in polynomial time using a conventional classical computer without
any restriction on entanglement [62, 63]. Specifically, quantum algorithms
consisting only of Clifford quantum gates can be simulated in
$\mathcal{O}(n^{2}m)$ polynomial steps with $n$ qubits and $m$ Clifford
operations. However, quantum algorithms consisting of Clifford+T gates take
exponential steps $\mathcal{O}(\kappa t^{3}\epsilon^{-2})$, with the number of
T-gates known as T count ($t$), stabilizer state ($k$) growing exponentially
$\mathcal{O}(2^{t})$, and an error rate ($\epsilon$) [62]. We note that some
quantum algorithms can be efficiently simulated using a sophisticated
classical technique like a tensor network on GPU tensor cores [64].

The Clifford+T gate set $\\{S,H,\text{CNOT},T\\}$ is considered a universal
gate set for digital QPUs. This is due to the feasibility of quantum error-
correcting, known as a surface code [65]. More importantly, the surface code
enables the creation of fault-tolerant digital quantum computers that surpass
the NISQ-era computers [35]. In contrast to NISQ computers, fault-tolerant
quantum computers are made up of error-free qubits and quantum gates that are
transpiled into the Clifford+T gate set. Therefore, this shows that for
quantum advantage in EO applications to be reached if and only if our quantum
learning models have a sufficiently high number $\mathcal{O}(10^{12})$ of
T-gates and generalize on unseen data points [66]. Otherwise, we can simulate
them efficiently using conventional classical computing resources.

Further, a hybrid classical-quantum approach for computational EO problems is
embedding high-dimensional classical data in a limited number of qubits and
optimizing the weights of a parameterized quantum model [67, 36]. There is yet
another challenging question: how notoriously difficult computational problems
can take advantage of both HPC and QC systems or when we should execute them
on an HPC instead of a QC system and vice versa. We decompose the
parameterized quantum model into the Clifford+T gate set at each learning
iteration to tackle these issues. If the parameterized quantum model only
includes Clifford gates and a small number of T-gates [68], then we execute it
on the HPC system instead of the QC machines since we already know that
Clifford gates and hundreds of T-gates can be simulated efficiently using a
conventional classical computer. We re-emphasize that quantum learning models
can be simulated efficiently using a classical computer without the need for
quantum computers if they do not have a high number of T-gates.

####  I-D1 Quantum machine learning: symmetry-breaking

Symmetry-breaking refers to asymmetric tunable weights of traditional ML
models such that the weights capture and rank the dataset’s features during
training. Consider a neural network with a single hidden layer.
Mathematically, it is defined as illustrated in Figure 2:

$\displaystyle h^{1}_{z}$
$\displaystyle=f\left(w_{0,z}+\sum_{k=0}^{2}w_{0,k,z}x_{k}\right),\quad
z=1,\dots,5,$ (1) $\displaystyle\tilde{y}_{t}$
$\displaystyle=f\left(h_{0,t}+\sum_{k=1}^{5}w_{1,k,t}h^{1}_{k}\right),\quad
t=1,2,$ (2)

where $f\left(\cdot\right)$ is a non-linear activation function, $w$’s denotes
a tuneable weight, and $x_{k}$ is the dataset’s feature. We note that $w$’s
must have different values identical to a linear regression model
$\tilde{y}\sim w_{0}+w_{1}\cdot x_{0}+w_{2}\cdot x_{1}$. If the model weights
are symmetric $w_{1}=w_{2}$, it has not learned the dataset’s feature. To
capture the dataset’s feature, the learning model must have asymmetric weights
$w_{1}\neq w_{2}$, or the learning model must break the symmetry in its
weights. Identical to the symmetry-breaking in conventional ML, the authors of
the article [18] implicitly demonstrated that QML models also must break
symmetry in their weights, resulting in better generalizability or more
expressive power and higher effective dimension than their classical
counterparts. In particular, they identified and disregarded some redundant
weights in their quantum models that are symmetric (e.g., the same digital
values) and do not simultaneously increase the QML model’s expressive power.
They, however, did not estimate the hardness of their QML models characterized
by non-Clifford T-gates that can be implemented efficiently on quantum
machines and otherwise difficult on conventional HPC systems.

Furthermore, to outperform classical learning models deployed on an HPC
system, we should invent and design QML models having thousands of T-gates,
and their expressive power (signaling the symmetry-breaking in QML models) is
higher than their classical counterparts [16]. There is (still) no such QML
model with thousands of T-gates and higher expressive power on unseen data
points than its classical counterpart.

##  II Quantum resource estimation for hyperspectral images

A hyperspectral imaging satellite, such as the EnMAP satellite, is a type of
imaging instrument mounted on a satellite and used to sense spectral
reflectances [69]. The mission of this satellite is to collect hyperspectral
imaging data that provides crucial information for scientific inquiries,
societal grand challenges, and key stakeholders and decision-makers. This
information pertains to various topics, such as climate change impact and
interventions, hazard and risk assessment, biodiversity and ecosystem
processes, land cover changes, and surface processes.

We already have seen that hyperspectral images require less quantum resources
than other remotely sensed datasets. They also have limited label information,
and there is limited availability of benchmark hyperspectral images compared
to conventional benchmark remote-sensing datasets, such as multispectral
images [70, 71]. When training QML models on limited benchmark-oriented
labeled hyperspectral image datasets, a classical layer can reduce the
dimensionality of the hyperspectral image dataset’s spectral bands due to the
limited number of input qubits. However, the degree of dimensionality
reduction required for the given hyperspectral image dataset depends on the
utilized quantum machines. Regardless of their error, this means whether we
can access a quantum machine with qubits $\leq 100$ or $>100$. The role of
classical machines in pre-processing the hyperspectral image dataset is
reduced as we can feed many informative features to a quantum machine with
less dimensionality reduction, especially as the number of qubits of quantum
machines increases. We assume we used EnMAP hyperspectral images with $103$
spectral bands and $610\times 340$ spatial dimensions. The EnMAP hyperspectral
images also have $207,400$ data points and $103$ features, which are small-
scale image datasets compared to conventional multispectral images. To execute
the QML model on the quantum machine having $\leq 100$ input qubits, we can
either reduce the spectral bands of the EnMAP hyperspectral images from $103$
to at most $100$ or select the most informative $100$ bands to be compatible
with the input qubits by utilizing a classical machine. Instead, for quantum
machines with more than $100$ input qubits, we can use a classical machine to
persevere more spectral bands of the EnMAP hyperspectral images when
performing the dimensionality reduction or the feature selection technique in
the spectral bands.

Toward quantum resource estimation, we assessed four different PQC models
expressed by the Clifford+T gate set (see Figs. 3-6). The Clifford+T gate set
is defined by $U_{1}$, $U_{2}$, $U_{3}$ and CNOT gates:

$\begin{split}U_{1}(\lambda)&=\begin{pmatrix}1&0\\\
0&e^{i\lambda}\end{pmatrix},\quad
U_{2}(\lambda,\phi)=\frac{1}{\sqrt{2}}\begin{pmatrix}1&-e^{i\phi}\\\
e^{i\lambda}&e^{i(\lambda+\phi)}\end{pmatrix},\\\
U_{3}(\lambda,\phi,\gamma)&=\begin{pmatrix}\cos(\lambda/2)&-e^{i\gamma}\sin(\lambda/2)\\\
-e^{i\phi}\sin(\lambda/2)&e^{i(\phi+\gamma)}\cos(\lambda/2)\end{pmatrix},\end{split}$
(3)

where, for example, $U_{1}(\pi/4)=T$, $U_{1}(\pi/2)=S$, $U_{2}(0,\pi)=H$.
Hence, the Clifford+T gate set is {$U_{1}(\pi/2)$, $U_{2}(0,\pi)$, CNOT,
$U_{1}(\pi/4)$}.

We have chosen the PQC models in Figs. 3-6 as benchmark QML models identical
to conventional benchmark deep learning (DL) models, such as Resnet [72]. The
quantum resource required for executing them on the quantum machine is
$\mathcal{O}(1)$ (constant time) if there is either no sign of T-gates or a
low number of T-gates. In particular, we will deploy them on either the HPC
system or the quantum machines depending on the existence and the number of
T-gates in their configuration during the training phase. Furthermore, the
number of T-gates defines the quantum resource required for deploying QML
models on quantum computers.

Figure 3: We transpiled a real-amplitude quantum circuit having depth-one into
the Clifford+T gate set. It is used to demonstrate the power of a PQC model by
the authors of the article [16].

Figure 4: We transpiled an energy-based quantum circuit having depth-one into
the Clifford+T gate set. This PQC model is proposed for the NISQ device by the
authors of the article [73].

We used the symmetry-breaking concept inherited from conventional neural
networks to determine the number of T-gates in our four PQCs [74]. Again, we
strongly emphasize that QML models break the symmetry in their weights to
decrease their redundant parameterized quantum gates, resulting in better
generalization on unseen data points than conventional neural networks [18].
Namely, each weight within a parameterized quantum layer must have different
digital values for capturing unique features. Therefore, we assumed that each
layer of the QML models must have, at most, a single T-gate at each learning
iteration, and our QML models having depth-one can only have one T-gate.

As for the quantum resource required for executing them on the quantum
hardware, we assumed also:

  1. 1.

If our PQCs have $10^{8}$ T-gates and $5$ logical qubits then we need
$158,431$ physical qubits (_i.e._ $9,375$ state distillation qubits, and
$149,056$ physical qubits) with a surface code distance of $d=25$, and our QML
models then take around $5$ hours per shot.

  2. 2.

If our PQCs have three T-gates and $5$ logical qubits then we need $50,700$
physical qubits (_i.e._ $14,400$ state distillation qubits, and $36,300$
physical qubits) with a surface code distance of $d=11$, and our QML models
then take around $8.12^{-8}$ hours per shot.

  3. 3.

If our PQCs have one T-gate and five logical qubits, then we need $15,135$
physical qubits (_i.e._ $14,400$ state distillation qubits, and $735$ physical
qubits) with a surface code distance of $d=7$, and our QML models then take
around $2.07^{-8}$ hours per shot.

Based on the study of the article’s authors [75, 76], we estimated the quantum
resources required for deploying QML models on error-correcting quantum
machines known as surface code quantum computers. Our estimation considers
that the quantum gate error is about $p=10^{-3}$, and the single round of the
surface code takes around $10^{-6}$ seconds. Here, the hours refer to T-gates
preparation; the article’s authors [75] provided a detailed spreadsheet for
the quantum resource estimation. The quantum resource estimation demonstrates
whether the QML models have to be deployed on quantum computers or not [68,
77], and it also generates the number of physical qubits required for
deploying quantum algorithms on the surface code quantum computers.

Figure 5: We transpiled a strongly-entangling quantum circuit having depth-one
transpiled into the Clifford+T gate set. This PQC model is proposed to build a
powerful quantum learning model in the article [78].  Figure 6: We transpiled
a hardware-efficient quantum circuit having depth-one into the Clifford+T gate
set. This PQC is used for quantum variational inference in the article [79].

##  III Conclusion

We assessed the quantum resource required to execute QML models on a digital
quantum computer to obtain a quantum advantage. We demonstrated that some
quantum advantage can only be obtained if and only if QML models have a
sufficient number of T-gates and generalize better on unseen data points than
their classical counterparts. To count the T-gates of a particular QML model,
we used the strong assumption that the QML models must break the symmetry in
their weights $-$ identical to the symmetry-breaking in conventional deep
learning models $-$ so that they become a more powerful model than their
counterpart classical learning models. Based on the number of T-gates, we
proposed a new HPC+QC paradigm (novel heterogeneous computing). In particular,
we can simulate QML models on an HPC system (_i.e._ CPU+GPU) if they comprise
a few hundred T-gates.

Toward quantum advantage in Earth observation, we focused on QML models for
hyperspectral images acquired by the EnMAP satellite since QML models can be
trained on a limited labeled dataset, and our hyperspectral images have
limited label information compared with multispectral images. For QML models,
we utilized four parameterized quantum circuits and estimated the quantum
resources required for deploying them on digital quantum machines. We found
that we can deploy our QML models on an HPC system instead of a QC system
since they only have a single T-gate due to the symmetry-breaking assumption.
To design QML models with around $\mathcal{O}(10^{8})$ that cannot be executed
on an HPC system, they must have almost a depth of $\mathcal{O}(10^{8})$,
which is impractical for current and future quantum computers. Toward quantum
advantage, it seems, therefore, reasonable to build, first, a special-purpose
digital quantum computer for some practically significant computational task
instead of a universal digital quantum computer similar to a D-Wave quantum
annealer.

As future and ongoing work, we will invent and design a QML model with a
reasonable depth that cannot be simulated on HPC systems but can be executed
efficiently on QC systems and simultaneously has more expressive power over
classical learning models.

## Acknowledgment

It is a pleasure to thank Gottfried Schwarz for reading, commenting, and
improving the quality of this paper.

## References

  * [1] M. Reichstein, G. Camps-Valls, B. Stevens, M. Jung, J. Denzler, N. Carvalhais, and Prabhat, “Deep learning and process understanding for data-driven earth system science,” _Nature_ , vol. 566, no. 7743, pp. 195–204, Feb 2019\. [Online]. Available: https://doi.org/10.1038/s41586-019-0912-1
  * [2] P. Ebel, Y. Xu, M. Schmitt, and X. X. Zhu, “Sen12ms-cr-ts: A remote-sensing data set for multimodal multitemporal cloud removal,” _IEEE Transactions on Geoscience and Remote Sensing_ , vol. 60, pp. 1–14, 2022.
  * [3] Y. Shi, X. x. Zhu, and R. Bamler, “Nonlocal compressive sensing-based sar tomography,” _IEEE Transactions on Geoscience and Remote Sensing_ , vol. 57, no. 5, pp. 3015–3024, 2019.
  * [4] X. X. Zhu, S. Montazeri, M. Ali, Y. Hua, Y. Wang, L. Mou, Y. Shi, F. Xu, and R. Bamler, “Deep learning meets SAR: Concepts, models, pitfalls, and perspectives,” _IEEE Geoscience and Remote Sensing Magazine_ , vol. 9, no. 4, pp. 143–172, 2021.
  * [5] G. Camps-Valls, J. Verrelst, J. Munoz-Mari, V. Laparra, F. Mateo-Jimenez, and J. Gomez-Dans, “A survey on gaussian processes for earth-observation data analysis: A comprehensive investigation,” _IEEE Geoscience and Remote Sensing Magazine_ , vol. 4, no. 2, pp. 58–78, 2016.
  * [6] D. Narmandakh, C. Butscher, F. Doulati Ardejani, H. Yang, T. Nagel, and R. Taherdangkoo, “The use of feed-forward and cascade-forward neural networks to determine swelling potential of clayey soils,” _Computers and Geotechnics_ , vol. 157, p. 105319, 2023. [Online]. Available: https://www.sciencedirect.com/science/article/pii/S0266352X23000769
  * [7] A. Stritih, P. Bebi, and A. Grêt-Regamey, “Quantifying uncertainties in earth observation-based ecosystem service assessments,” _Environmental Modelling and Software_ , vol. 111, pp. 300–310, 2019. [Online]. Available: https://www.sciencedirect.com/science/article/pii/S1364815218300884
  * [8] J. Pathak, S. Subramanian, P. Harrington, S. Raja, A. Chattopadhyay, M. Mardani, T. Kurth, D. Hall, Z. Li, K. Azizzadenesheli, P. Hassanzadeh, K. Kashinath, and A. Anandkumar, “Fourcastnet: A global data-driven high-resolution weather model using adaptive fourier neural operators,” 2022\.
  * [9] G. Cheng, X. Xie, J. Han, L. Guo, and G. S. Xia, “Remote sensing image scene classification meets deep learning: Challenges, methods, benchmarks, and opportunities,” _IEEE Journal of Selected Topics in Applied Earth Observations and Remote Sensing_ , vol. 13, pp. 3735–3756, 2020\.
  * [10] S. Arora and B. Barak, _Computational Complexity: A Modern Approach_. Cambridge University Press, 2009\.
  * [11] E. Farhi, J. Goldstone, S. Gutmann, and M. Sipser, “Quantum computation by adiabatic evolution,” _arXiv preprint quant-ph/0001106_ , 2000.
  * [12] A. Lucas, “Ising formulations of many np problems,” _Frontiers in Physics_ , vol. 2, p. 5, 2014. [Online]. Available: https://www.frontiersin.org/article/10.3389/fphy.2014.00005
  * [13] P. Rebentrost, M. Mohseni, and S. Lloyd, “Quantum support vector machine for big data classification,” _Phys. Rev. Lett._ , vol. 113, p. 130503, Sep 2014\. [Online]. Available: https://link.aps.org/doi/10.1103/PhysRevLett.113.130503
  * [14] J. Allcock, C.-Y. Hsieh, I. Kerenidis, and S. Zhang, “Quantum algorithms for feedforward neural networks,” _ACM Transactions on Quantum Computing_ , vol. 1, no. 1, oct 2020. [Online]. Available: https://doi.org/10.1145/3411466
  * [15] J. Biamonte, P. Wittek, N. Pancotti, P. Rebentrost, N. Wiebe, and S. Lloyd, “Quantum machine learning,” _Nature_ , vol. 549, no. 7671, pp. 195–202, Sep 2017. [Online]. Available: https://doi.org/10.1038/nature23474
  * [16] A. Abbas, D. Sutter, C. Zoufal, A. Lucchi, A. Figalli, and S. Woerner, “The power of quantum neural networks,” _Nature Computational Science_ , vol. 1, no. 6, pp. 403–409, Jun 2021. [Online]. Available: https://doi.org/10.1038/s43588-021-00084-1
  * [17] C. Gyurik and V. Dunjko, “On establishing learning separations between classical and quantum machine learning with classical data,” Aug. 2022, arXiv:2208.06339 [quant-ph]. [Online]. Available: http://arxiv.org/abs/2208.06339
  * [18] T. Haug, K. Bharti, and M. Kim, “Capacity and quantum geometry of parametrized quantum circuits,” _PRX Quantum_ , vol. 2, p. 040309, Oct 2021. [Online]. Available: https://link.aps.org/doi/10.1103/PRXQuantum.2.040309
  * [19] A. M. Childs, J.-P. Liu, and A. Ostrander, “High-precision quantum algorithms for partial differential equations,” _Quantum_ , vol. 5, p. 574, nov 2021\. [Online]. Available: https://doi.org/10.22331%2Fq-2021-11-10-574
  * [20] A. J. Pool, A. D. Somoza, M. Lubasch, and B. Horstmann, “Solving partial differential equations using a quantum computer,” in _2022 IEEE International Conference on Quantum Computing and Engineering (QCE)_ , 2022, pp. 864–866.
  * [21] N. Gourianov, M. Lubasch, S. Dolgov, Q. Y. van den Berg, H. Babaee, P. Givi, M. Kiffner, and D. Jaksch, “A quantum-inspired approach to exploit turbulence structures,” _Nature Computational Science_ , vol. 2, no. 1, pp. 30–37, 2022.
  * [22] D-Wave Systems, “D-Wave Quantum Computing,” 2023.
  * [23] C. Gross and I. Bloch, “Quantum simulations with ultracold atoms in optical lattices,” _Science_ , vol. 357, no. 6355, pp. 995–1001, 2017. [Online]. Available: https://www.science.org/doi/abs/10.1126/science.aal3837
  * [24] L. Funcke, T. Hartung, K. Jansen, and S. Kühn, “Review on Quantum Computing for Lattice Field Theory,” _PoS_ , vol. LATTICE2022, p. 228, 2023.
  * [25] IBM, “Ibm quantum computing,” 2023.
  * [26] A. Acín, I. Bloch, H. Buhrman, T. Calarco, C. Eichler, J. Eisert, D. Esteve, N. Gisin, S. J. Glaser, F. Jelezko, S. Kuhr, M. Lewenstein, M. F. Riedel, P. O. Schmidt, R. Thew, A. Wallraff, I. Walmsley, and F. K. Wilhelm, “The quantum technologies roadmap: a european community view,” _New Journal of Physics_ , vol. 20, no. 8, p. 080201, aug 2018. [Online]. Available: https://dx.doi.org/10.1088/1367-2630/aad1ea
  * [27] S. Aaronson, “How much structure is needed for huge quantum speedups?” 2022.
  * [28] M. Dalmonte, B. Vermersch, and P. Zoller, “Quantum simulation and spectroscopy of entanglement hamiltonians,” _Nature Physics_ , vol. 14, no. 8, pp. 827–831, Aug 2018. [Online]. Available: https://doi.org/10.1038/s41567-018-0151-7
  * [29] S. Lu, M. C. Bañuls, and J. I. Cirac, “Algorithms for quantum simulation at finite energies,” _PRX Quantum_ , vol. 2, p. 020321, May 2021. [Online]. Available: https://link.aps.org/doi/10.1103/PRXQuantum.2.020321
  * [30] S. Otgonbaatar, O. Nurmi, M. Johansson, J. Mäkelä, T. Kocman, P. Gawron, Z. Puchala, J. Mielczarek, A. Miroszewski, and C. O. Dumitru, “Quantum computing for climate change detection, climate modeling, and climate digital twin,” Tech. Rep., November 2023. [Online]. Available: https://elib.dlr.de/198760/
  * [31] R. Babbush, J. R. McClean, M. Newman, C. Gidney, S. Boixo, and H. Neven, “Focus beyond quadratic speedups for error-corrected quantum advantage,” _PRX Quantum_ , vol. 2, no. 1, mar 2021. [Online]. Available: https://doi.org/10.1103%2Fprxquantum.2.010103
  * [32] M. A. Nielsen, I. Chuang, and L. K. Grover, “Quantum Computation and Quantum Information,” _American Journal of Physics_ , vol. 70, no. 5, pp. 558–559, 05 2002. [Online]. Available: https://doi.org/10.1119/1.1463744
  * [33] S. Otgonbaatar and M. Datcu, “Assembly of a coreset of earth observation images on a small quantum computer,” _Electronics_ , vol. 10, no. 20, 2021\. [Online]. Available: https://www.mdpi.com/2079-9292/10/20/2482
  * [34] P. Helber, B. Bischke, A. Dengel, and D. Borth, “Eurosat: A novel dataset and deep learning benchmark for land use and land cover classification,” _IEEE Journal of Selected Topics in Applied Earth Observations and Remote Sensing_ , vol. 12, no. 7, pp. 2217–2226, 2019.
  * [35] R. Acharya and et al, “Suppressing quantum errors by scaling a surface code logical qubit,” 2022.
  * [36] S. Otgonbaatar and M. Datcu, “Classification of remote sensing images with parameterized quantum gates,” _IEEE Geoscience and Remote Sensing Letters_ , pp. 1–5, 2021.
  * [37] P. Gawron and S. Lewiński, “Multi-spectral image classification with quantum neural network,” in _IGARSS 2020-2020 IEEE International Geoscience and Remote Sensing Symposium_. IEEE, 2020, pp. 3513–3516.
  * [38] D. A. Zaidenberg, A. Sebastianelli, D. Spiller, B. Le Saux, and S. L. Ullo, “Advantages and bottlenecks of quantum machine learning for remote sensing,” in _2021 IEEE International Geoscience and Remote Sensing Symposium IGARSS_ , 2021, pp. 5680–5683.
  * [39] A. Sebastianelli, D. A. Zaidenberg, D. Spiller, B. Le Saux, and S. L. Ullo, “On circuit-based hybrid quantum neural networks for remote sensing imagery classification,” _IEEE Journal of Selected Topics in Applied Earth Observations and Remote Sensing_ , pp. 1–1, 2021.
  * [40] M. K. Gupta, M. Romaszewski, and P. Gawron, “Potential of quantum machine learning for processing multispectral earth observation data,” _TechRxiv_ , 2023.
  * [41] Y. Yang and S. Newsam, “Bag-of-visual-words and spatial extensions for land-use classification,” in _Proceedings of the 18th SIGSPATIAL International Conference on Advances in Geographic Information Systems_ , ser. GIS ’10. New York, NY, USA: Association for Computing Machinery, 2010, p. 270–279. [Online]. Available: https://doi.org/10.1145/1869790.1869829
  * [42] S. Otgonbaatar, G. Schwarz, M. Datcu, and D. Kranzlmüller, “Quantum transfer learning for real-world, small, and high-dimensional remotely sensed datasets,” _IEEE Journal of Selected Topics in Applied Earth Observations and Remote Sensing_ , vol. 16, pp. 9223–9230, 2023.
  * [43] S. Otgonbaatar and M. Datcu, “Natural embedding of the stokes parameters of polarimetric synthetic aperture radar images in a gate-based quantum computer,” _IEEE Transactions on Geoscience and Remote Sensing_ , pp. 1–8, 2021.
  * [44] H.-Y. Huang, M. Broughton, M. Mohseni, R. Babbush, S. Boixo, H. Neven, and J. R. McClean, “Power of data in quantum machine learning,” 2021.
  * [45] M. K. Gupta, M. Beseda, and P. Gawron, “How quantum computing-friendly multispectral data can be?” in _IGARSS 2022-2022 IEEE International Geoscience and Remote Sensing Symposium_. IEEE, 2022, pp. 4153–4156.
  * [46] E. Boyda, S. Basu, S. Ganguly, A. Michaelis, S. Mukhopadhyay, and R. R. Nemani, “Deploying a quantum annealing processor to detect tree cover in aerial imagery of california,” _PLOS ONE_ , vol. 12, no. 2, pp. 1–22, 02 2017. [Online]. Available: https://doi.org/10.1371/journal.pone.0172505
  * [47] S. Otgonbaatar and M. Datcu, “Quantum annealing approach: Feature extraction and segmentation of synthetic aperture radar image,” in _IGARSS 2020 - 2020 IEEE International Geoscience and Remote Sensing Symposium_ , 2020, pp. 3692–3695.
  * [48] D. Willsch, M. Willsch, H. De Raedt, and K. Michielsen, “Support vector machines on the d-wave quantum annealer,” _Computer Physics Communications_ , vol. 248, p. 107006, 2020. [Online]. Available: https://www.sciencedirect.com/science/article/pii/S001046551930342X
  * [49] G. Cavallaro, M. Riedel, M. Richerzhagen, J. A. Benediktsson, and A. Plaza, “On understanding big data impacts in remotely sensed image classification using support vector machine methods,” _IEEE Journal of Selected Topics in Applied Earth Observations and Remote Sensing_ , vol. 8, no. 10, pp. 4634–4646, 2015.
  * [50] A. Delilbasic, B. L. Saux, M. Riedel, K. Michielsen, and G. Cavallaro, “A single-step multiclass svm based on quantum annealing for remote sensing data classification,” 2023.
  * [51] S. Otgonbaatar, M. Datcu, and B. Demir, “Coreset of hyperspectral images on a small quantum computer,” in _IGARSS 2022 - 2022 IEEE International Geoscience and Remote Sensing Symposium_ , 2022, pp. 4923–4926.
  * [52] S. Otgonbaatar and M. Datcu, “A quantum annealer for subset feature selection and the classification of hyperspectral images,” _IEEE Journal of Selected Topics in Applied Earth Observations and Remote Sensing_ , vol. 14, pp. 7057–7065, 2021.
  * [53] J. Chen, E. M. Stoudenmire, and S. R. White, “The quantum fourier transform has small entanglement,” _arXiv: 2210.08468_ , 2022.
  * [54] E. M. Stoudenmire and D. J. Schwab, “Supervised learning with quantum-inspired tensor networks,” _arXiv_ , 2016. [Online]. Available: https://arxiv.org/abs/1605.05775
  * [55] Z.-F. Gao, S. Cheng, R.-Q. He, Z. Y. Xie, H.-H. Zhao, Z.-Y. Lu, and T. Xiang, “Compressing deep neural networks by matrix product operators,” _Phys. Rev. Res._ , vol. 2, p. 023300, Jun 2020. [Online]. Available: https://link.aps.org/doi/10.1103/PhysRevResearch.2.023300
  * [56] F. Verstraete, T. Nishino, U. Schollwöck, M. C. Bañuls, G. K. Chan, and M. E. Stoudenmire, “Density matrix renormalization group, 30 years on,” _Nature Reviews Physics_ , Apr 2023. [Online]. Available: https://doi.org/10.1038/s42254-023-00572-5
  * [57] H. Huang, X.-Y. Liu, W. Tong, T. Zhang, A. Walid, and X. Wang, “High performance hierarchical tucker tensor learning using gpu tensor cores,” _IEEE Transactions on Computers_ , pp. 1–1, 2022.
  * [58] S. Otgonbaatar and D. Kranzlmüller, “Quantum-inspired tensor network for earth science,” in _IGARSS 2023 - 2023 IEEE International Geoscience and Remote Sensing Symposium_ , 2023, pp. 788–791.
  * [59] F. Fan, Y. Shi, T. Guggemos, and X. X. Zhu, “Hybrid quantum-classical convolutional neural network model for image classification,” _IEEE Transactions on Neural Networks and Learning Systems_ , pp. 1–15, 2023.
  * [60] M. C. Caro, H.-Y. Huang, M. Cerezo, K. Sharma, A. Sornborger, L. Cincio, and P. J. Coles, “Generalization in quantum machine learning from few training data,” _Nature Communications_ , vol. 13, no. 1, p. 4919, Aug 2022. [Online]. Available: https://doi.org/10.1038/s41467-022-32550-3
  * [61] J. Preskill, “Quantum computing in the nisq era and beyond,” _Quantum_ , vol. 2, p. 79, 2018.
  * [62] S. Bravyi and D. Gosset, “Improved classical simulation of quantum circuits dominated by clifford gates,” _Phys. Rev. Lett._ , vol. 116, p. 250501, Jun 2016. [Online]. Available: https://link.aps.org/doi/10.1103/PhysRevLett.116.250501
  * [63] S. Aaronson and D. Gottesman, “Improved simulation of stabilizer circuits,” _Physical Review A_ , vol. 70, no. 5, nov 2004. [Online]. Available: https://doi.org/10.1103%2Fphysreva.70.052328
  * [64] J. Tindall, M. Fishman, M. Stoudenmire, and D. Sels, “Efficient tensor network simulation of ibm’s kicked ising experiment,” 2023.
  * [65] D. Litinski, “A game of surface codes: Large-scale quantum computing with lattice surgery,” _Quantum_ , vol. 3, p. 128, mar 2019. [Online]. Available: https://doi.org/10.22331%2Fq-2019-03-05-128
  * [66] M. Hinsche, M. Ioannou, A. Nietner, J. Haferkamp, Y. Quek, D. Hangleiter, J.-P. Seifert, J. Eisert, and R. Sweke, “A single T-gate makes distribution learning hard,” _arXiv preprint arXiv:2207.03140_ , 2022.
  * [67] A. Mari, T. R. Bromley, J. Izaac, M. Schuld, and N. Killoran, “Transfer learning in hybrid classical-quantum neural networks,” _Quantum_ , vol. 4, p. 340, Oct. 2020. [Online]. Available: https://doi.org/10.22331/q-2020-10-09-340
  * [68] M. E. Beverland, P. Murali, M. Troyer, K. M. Svore, T. Hoefler, V. Kliuchnikov, G. H. Low, M. Soeken, A. Sundaram, and A. Vaschillo, “Assessing requirements to scale to practical quantum advantage,” 2022.
  * [69] (2022) German hyperspectral satellite: EnMAP. [Online]. Available: https://www.enmap.org/mission/
  * [70] G. Cheng, J. Han, and X. Lu, “Remote sensing image scene classification: Benchmark and state of the art,” _Proceedings of the IEEE_ , vol. 105, no. 10, pp. 1865–1883, 2017.
  * [71] M. Paoletti, J. Haut, J. Plaza, and A. Plaza, “Deep learning classifiers for hyperspectral imaging: A review,” _ISPRS Journal of Photogrammetry and Remote Sensing_ , vol. 158, pp. 279–317, 2019. [Online]. Available: https://www.sciencedirect.com/science/article/pii/S0924271619302187
  * [72] K. He, X. Zhang, S. Ren, and J. Sun, “Deep residual learning for image recognition,” 2015.
  * [73] E. Farhi and H. Neven, “Classification with quantum neural networks on near term processors,” 2018.
  * [74] R. Fok, A. An, and X. Wang, “Spontaneous symmetry breaking in neural networks,” 2017.
  * [75] A. G. Fowler and C. Gidney, “Low overhead quantum computation using lattice surgery,” 2019.
  * [76] C. Gidney, “Stim: a fast stabilizer circuit simulator,” _Quantum_ , vol. 5, p. 497, jul 2021. [Online]. Available: https://doi.org/10.22331%2Fq-2021-07-06-497
  * [77] M. Reiher, N. Wiebe, K. M. Svore, D. Wecker, and M. Troyer, “Elucidating reaction mechanisms on quantum computers,” _Proceedings of the National Academy of Sciences_ , vol. 114, no. 29, pp. 7555–7560, 2017. [Online]. Available: https://www.pnas.org/doi/abs/10.1073/pnas.1619152114
  * [78] M. Schuld, A. Bocharov, K. M. Svore, and N. Wiebe, “Circuit-centric quantum classifiers,” _Physical Review A_ , vol. 101, no. 3, Mar 2020. [Online]. Available: http://dx.doi.org/10.1103/PhysRevA.101.032308
  * [79] M. Benedetti, B. Coyle, M. Fiorentini, M. Lubasch, and M. Rosenkranz, “Variational inference with a quantum computer,” _Phys. Rev. Appl._ , vol. 16, p. 044057, Oct 2021. [Online]. Available: https://link.aps.org/doi/10.1103/PhysRevApplied.16.044057


# A Unitary Transform Based Generalized Approximate Message Passing

###### Abstract

We consider the problem of recovering an unknown signal
${\mathbf{x}}\in{\mathbb{R}}^{n}$ from general nonlinear measurements obtained
through a generalized linear model (GLM), i.e.,
${\mathbf{y}}=f\left({\mathbf{A}}{\mathbf{x}}+{\mathbf{w}}\right)$, where
$f(\cdot)$ is a componentwise nonlinear function. Based on the unitary
transform approximate message passing (UAMP) and expectation propagation, a
unitary transform based generalized approximate message passing (GUAMP)
algorithm is proposed for general measurement matrices $\bf{A}$, in particular
highly correlated matrices. Experimental results on quantized compressed
sensing demonstrate that the proposed GUAMP significantly outperforms state-
of-the-art GAMP and GVAMP under correlated matrices $\bf{A}$.

Index Terms—  GLM, AMP, GAMP, message passing, quantized compressed sensing

##  1 Introduction

We consider the general problem of inference on generalized linear models
(GLM) [1], i.e., recovering an unknown signal
${\mathbf{x}}\in{\mathbb{R}}^{n}$ from a noisy linear transform followed by
componentwise nonlinear measurements

$\displaystyle{\mathbf{y}}=f\left({\mathbf{A}}{\mathbf{x}}+{\mathbf{w}}\right),$
(1)

where ${\bf{A}}\in\mathbb{R}^{m\times n}$ is a known linear mixing matrix,
${\bf{w}}\sim\mathcal{N}({\bf{w}};0,\sigma^{2}{\bf{I}}_{m})$ is an i.i.d.
Gaussian noise with known variance, and ${f}(\cdot):\mathbb{R}^{m\times
1}\to\mathcal{Q}^{m\times 1}$ is an componentwise nonlinear link function.
Denote ${\mathbf{z}}\triangleq{\mathbf{A}}{\mathbf{x}}\in{\mathbb{R}}^{m}$ as
the hidden linear transform outputs, the componentwise nonlinear function
$f(\cdot)$ can be equivalently described in a probabilistic way using a fully
factorized likelihood function as follows

$\displaystyle
p({\mathbf{y}}|{\mathbf{z}})=\prod\limits_{i=1}^{m}p(y_{i}|z_{i}),$ (2)

where $p(y_{i}|z_{i})$ is determined by the specific nonlinear function
$f(\cdot)$. The prior distribution $p({\mathbf{x}})$ of $\bf{x}$ is also
assumed to be fully factorized
$p({\mathbf{x}})=\prod\limits_{j=1}^{n}p(x_{j})$ for simplicity. GLM inference
has wide applications in science and engineering such as wireless
communications, signal processing, and machine learning. In the special case
of identity function $f(\cdot)$, the nonlinear GLM in (1) will reduce to the
popular standard linear model (SLM) as follows

${\mathbf{y}}={\mathbf{A}}{\mathbf{x}}+{\mathbf{w}}.$ (3)

Thus, GLM is actually an extention of SLM from linear measurements to
nonlinear measurements, which are prevalent in some real-world applications
such as quantized compressed sensing, pattern classification, phase retrieval,
etc.

A variety of algorithms have been proposed for inference over GLMs and SLMs.
Among them, the past decade has witnessed an advent of one distinguished
family of probabilistic algorithms called message passing algorithm. Among
them, the most famous ones are the approximate message passing (AMP) algorithm
[2, 3] for SLM and generalized approximate message passing (GAMP) [4] for GLM,
which have been proved to be optimal under i.i.d. Gaussian matrices $\bf{A}$.
However, both AMP and GAMP diverge for general $\bf{A}$. In [5, 6, 7], the
GAMP is incorporated into the sparse Bayesian learning (SBL) to improve the
robustness of GAMP and reduce the computation complexity of SBL. Vector AMP
(VAMP) [8] (or OAMP [9], one similar algorithm to VAMP) and generalized VAMP
(GVAMP) [10] have been proposed to improve the performance with general
$\bf{A}$ for SLM and GLM, respectively. The VAMP is derived from expectation
propagation (EP) [11], another powerful message passing algorithm. Remarkably,
it has been demonstrated in [12, 13] that all the AMP, GAMP, and GVAMP, can be
derived concisely as special instances of EP under different assumptions and
thus be unified within a single EP framework. Despite significant improvement
in robustness, in particular for right-orthogonally invariant matrices, VAMP
and GVAMP still suffer poor convergence in some more challenging scenarios,
e.g., the measurement matrix is highly correlated [14, 15].

In this work, we focus on the AMP variant: unitary approximate message passing
(UAMP) [16], which was formerly called UTAMP [16, 17, 18]. The motivation
behind UAMP is this: since AMP has proven to work well when the elements of
$\mathbf{A}$ in (3) are uncorrelated, one can artificially constructs such a
“good” and equivalent measurement model simply by performing a singular value
decomposition (SVD) on
${\mathbf{A}}={\mathbf{U}}{\bm{\Sigma}}{\mathbf{V}}^{\rm T}$, where
${\mathbf{U}}\in{\mathbb{R}}^{m\times r}$,
${\bm{\Sigma}}\in{\mathbb{R}}^{r\times r}$,
${\mathbf{V}}\in{\mathbb{R}}^{n\times r}$ and $r={\rm rank}({\mathbf{A}})$,
and by left multiplying ${\mathbf{U}}^{\rm T}$ on the original SLM (3), thus
leading to an equivalent SLM as [16]

$\displaystyle\tilde{\mathbf{b}}={\mathbf{Q}}{\mathbf{x}}+\tilde{\mathbf{w}},$
(4)

where $\tilde{\mathbf{b}}={\mathbf{U}}^{\rm
T}{\mathbf{y}}\in{\mathbb{R}}^{r}$,
${\mathbf{Q}}={\bm{\Sigma}}{\mathbf{V}}^{\rm T}\in{\mathbb{R}}^{r\times n}$
and $\tilde{\mathbf{w}}={\mathbf{U}}^{\rm T}{\mathbf{w}}\in{\mathbb{R}}^{r}$
correspond to unitary-transformed pseudo linear measurements, linear mixing
matrix, and additive noise, respectively. Subsequently, one can readily run
standard AMP on the equivalent model (4), which leads to one version of UAMP.
Further, two averaging operations can be performed to obtain the other version
of UAMP, where the number of matrix-vector products per iteration is reduced
from 4 to 2 [16, 17, 18]. Despite its simplicity, UAMP has been shown to be
more robust than AMP and VAMP under some types of “tough” measurement
matrices, e.g., highly correlated matrices [17, 18]. Nevertheless, the
extension of UAMP to the case of GLM has remained lacking.

In this paper, we extend UAMP from SLM to GLM and propose the generalized UAMP
(GUAMP). The key idea is to utilize the unified Bayesian inference framework
in [13, 19] and EP [11] to iteratively decompose the original nonlinear
measurement model (1) into a series of pseudo SLMs, whereby UAMP could be
conducted. Conceptually, GUAMP consists of two modules, namely AMP module and
GAMP module, as shown in Fig. 1. Extrinsic information [20] is exchanged per
iteration between the two modules before GUAMP finally converge. Experimental
results demonstrate that GUAMP significantly outperforms both GAMP and GVAMP
under highly correlated matrices $\bf{A}$.

##  2 The GUAMP Algorithm

The key idea is that, inspired by UAMP[16], we introduce an additional hidden
variable ${\mathbf{b}}\in{\mathbb{R}}^{r}$ and thus an equivalent
representation of ${\mathbf{z}}\triangleq{\mathbf{A}}{\mathbf{x}}$ as follows

$\displaystyle{\mathbf{b}}={\bm{\Sigma}}{\mathbf{V}}^{\rm
T}{\mathbf{x}}\triangleq{\mathbf{Q}}{\mathbf{x}},$ (5a)
$\displaystyle{\mathbf{z}}={\mathbf{U}}{\mathbf{b}},$ (5b)

where, as in the original UAMP [16], ${\mathbf{U}}\in{\mathbb{R}}^{m\times
r}$, ${\bm{\Sigma}}\in{\mathbb{R}}^{r\times r}$,
${\mathbf{V}}\in{\mathbb{R}}^{n\times r}$ are obtained from SVD
${\mathbf{A}}={\mathbf{U}}{\bm{\Sigma}}{\mathbf{V}}^{\rm T}$ when ${\rm
rank}({\mathbf{A}})=r$. Note that the underdetermined ($m\leq n$) and
overdetermined ($m\leq n$) cases are treated in a unified way.

The corresponding factor graph of the original GLM (1) can be equivalently
shown in Fig.1 (a). Subsequently, using the unified modular framework in [13],
the inference on the factor graph in Fig.1 (a) can be decomposed into two
modules, namely one AMP module and one GAMP module, as shown in Fig.1 (b). At
a high level, in each iteration, module A performs standard AMP w.r.t.
$\bf{x}$ on an equivalent SLM with a pseudo measurement matrix $\bf{Q}$ while
module B performs GAMP w.r.t. $\bf{b}$ on an equivalent GLM with a pseudo
measurement matrix $\bf{U}$, rather than the original measurement matrix
$\bf{A}$. Module A and Module B exchange extrinsic information with each other
in the same way as [13] and this process proceeds before convergence. Please
refer to [13] for further details of the modular representation perspective
via EP [11].

Fig. 1: (a). Equivalent factor graph representation of GLM. (b) Modular
representation of the GUAMP algorithm.

In the following, we describe the implementation details of GUAMP. First of
all, we initialize ${\bm{\tau}}^{b}(0)\in{\mathbb{R}}^{r}$,
$\hat{\mathbf{s}}_{\rm B}(-1)\in{\mathbb{R}}^{m}$, ${\mathbf{b}}_{{\rm
A}}^{\rm ext}(1)\in{\mathbb{R}}^{r}$ and variances ${\mathbf{v}}_{{\rm
A}}^{{\rm ext},b}(1)\in{\mathbb{R}}^{r}$ in module B, and initialize
$\hat{\mathbf{p}}_{\rm A}(0)\in{\mathbb{R}}^{r}$, $\hat{\bm{\tau}}_{\rm
A}^{p}(0)\in{\mathbb{R}}^{r}$, $\hat{\mathbf{x}}(0)\in{\mathbb{R}}^{n}$. The
number of outer iterations between module A and module B is set as $\rm T_{\rm
max}$, while the number of inner iterations in module A and B during an outer
iteration are set as $T_{\rm A}$ and $T_{\rm B}$, respectively. Next, we
present the details of the operation in the two modules.

###  2.1 The GAMP Module

The extrinsic means ${\mathbf{b}}_{{\rm A}}^{\rm ext}(t)\in{\mathbb{R}}^{r}$
and variances ${\mathbf{v}}_{{\rm A}}^{{\rm ext},b}(t)\in{\mathbb{R}}^{r}$
transmitting from module A to module B can be viewed as the prior means and
variances of $\mathbf{b}$. i.e.,

$\displaystyle p({\mathbf{b}})={\mathcal{N}}({\mathbf{b}};{\mathbf{b}}_{{\rm
A}}^{\rm ext}(t),{\rm diag}({\mathbf{v}}_{{\rm A}}^{{\rm ext},b}(t))).$ (6)

In addition, ${\mathbf{z}}={\mathbf{U}}{\mathbf{b}}$ and
${\mathbf{y}}|{\mathbf{z}}\sim p({\mathbf{y}}|{\mathbf{z}})$. Thus we could
run the GAMP algorithm treating $\bf{b}$ as the unknown signal and $\bf{U}$ as
the measurement matrix as follows.

  * •

[Step 1 (B)] Perform output linear step to obtain $\hat{\mathbf{p}}_{{\rm
B}}(t)\in{\mathbb{R}}^{m}$ and ${\bm{\tau}}^{p}_{{\rm
B}}(t)\in{\mathbb{R}}^{m}$ as

$\displaystyle{\bm{\tau}}^{p}_{{\rm B}}(t)$
$\displaystyle=|{\mathbf{U}}|^{2}{\bm{\tau}}^{b}(t),$ (7a)
$\displaystyle\hat{\bm{p}}_{{\rm B}}(t)$
$\displaystyle={\mathbf{U}}\hat{\mathbf{b}}(t)-{\bm{\tau}}^{p}_{{\rm
B}}(t)\hat{\mathbf{s}}_{{\rm B}}(t-1).$ (7b)

  * •

[Step 2 (B)] Perform output linear step to obtain $\hat{s}_{{\rm B},i}(t)$ and
${\tau}_{{\rm B},i}(t)$ as

$\displaystyle\hat{s}_{{\rm B},i}(t)=g_{\rm out}(\hat{p}_{{\rm
B},i}(t),y_{i},\tau^{p}_{{\rm B},i}(t)),$ (8a) $\displaystyle{\tau}_{{\rm
B},i}(t)=-\frac{\partial g_{\rm out}(\hat{p}_{{\rm
B},i}(t),y_{i},\tau^{p}_{{\rm B},i}(t))}{\partial\hat{p}},$ (8b)

for $i=1,2,\cdots,m$, where

$\displaystyle g_{\rm out}(\hat{p},y,\tau^{p})$
$\displaystyle=\frac{\hat{z}^{0}-\hat{p}}{\tau_{p}},$ (9)
$\displaystyle-\frac{\partial g_{\rm
out}(\hat{p},y,\tau^{p})}{\partial\hat{p}}$ $\displaystyle=\frac{\tau_{p}-{\rm
Var}(z|\hat{p},y)}{(\tau_{p})^{2}},$ (10)

and the posterior means $\hat{z}^{0}$ and variances ${\rm Var}(z|\hat{p},y)$
are computed w.r.t. the posterior
$\propto{\mathcal{CN}}(z;\hat{p},\tau_{p})p(y|z)$. See [4] for further
details.

  * •

[Step 3 (B)] Perform input linear step as

$\displaystyle{\bm{\tau}}^{r}_{{\rm B}}(t)$
$\displaystyle=\left(|{\mathbf{U}}^{\rm T}|^{2}{\bm{\tau}}^{s}_{{\rm
B}}(t)\right)^{-1},$ (11a) $\displaystyle\hat{\mathbf{r}}_{{\rm B}}(t)$
$\displaystyle=\hat{\mathbf{b}}(t)+{\bm{\tau}}^{r}_{{\rm
B}}(t)\odot\left({\mathbf{U}}^{\rm T}\hat{\mathbf{s}}_{{\rm B}}(t-1)\right).$
(11b)

  * •

[Step 4 (B)] Perform input nonlinear step in Module B to obtain the posterior
means and variances of variable $\mathbf{b}$ as: For $j=1,2,\cdots,r$,

$\displaystyle\hat{b}_{j}(t)$ $\displaystyle=\frac{\hat{r}_{{\rm
B},j}(t){v}_{{\rm A},j}^{{\rm ext},b}(t)+{b}_{{\rm A},j}^{\rm
ext}(t)\tau^{r}_{{\rm B},j}(t)}{\tau^{r}_{{\rm B},j}(t)+{v}_{{\rm A},j}^{{\rm
ext},b}(t)},$ (12a) $\displaystyle\tau_{j}^{b}(t)$
$\displaystyle=\frac{{v}_{{\rm A},j}^{{\rm ext},b}(t)\tau^{r}_{{\rm
B},j}(t)}{\tau^{r}_{{\rm B},j}(t)+{v}_{{\rm A},j}^{{\rm ext},b}(t)}.$ (12b)

After running $T_{\rm B}$ iterations, the extrinsic means ${\mathbf{b}}_{{\rm
B}}^{\rm ext}(t)$ and variances ${\mathbf{v}}_{{\rm B}}^{{\rm ext},b}(t)$ of
$\mathbf{b}$ from module B to module A can be calculated as follows [13]

$\displaystyle{\mathbf{b}}_{{\rm B}}^{\rm ext}(t)=\hat{\mathbf{r}}_{{\rm
B}}(t),$ (13) $\displaystyle{\mathbf{v}}_{{\rm B}}^{{\rm
ext},b}(t)={\bm{\tau}}^{r}_{{\rm B}}(t).$ (14)

Remark: It is worth noting that if the GLM (1) degenerates to the SLM (3), it
can be verified that the extrinsic information $\bf{b}_{{\rm B}}^{\rm ext}$
and ${\mathbf{v}}_{{\rm B}}^{{\rm ext},b}(t)$ from module B to module A are
always ${\mathbf{b}}_{{\rm B}}^{\rm ext}={\mathbf{U}}^{\rm T}{\mathbf{y}}$ and
${\mathbf{v}}_{{\rm B}}^{{\rm ext},b}=\sigma^{2}{\mathbf{1}}_{r}$.
Consequently, the GUAMP reduces to the UAMP in the special case of SLM.

###  2.2 The AMP Module

Algorithm 1 GUAMP Algorithm

1: Initialize ${\bm{\tau}}^{b}(0)\in{\mathbb{R}}^{r}$, $\hat{\mathbf{s}}_{\rm
B}\in{\mathbb{R}}^{m}$, ${\mathbf{b}}_{{\rm A}}^{\rm
ext}(1)\in{\mathbb{R}}^{r}$ and variances ${\mathbf{v}}_{{\rm A}}^{{\rm
ext},b}(1)\in{\mathbb{R}}^{r}$ in module B, initialize $\hat{\mathbf{p}}_{\rm
A}(0)\in{\mathbb{R}}^{r}$, $\hat{\bm{\tau}}_{\rm
A}^{p}(0)\in{\mathbb{R}}^{r}$, $\hat{\mathbf{x}}(0)\in{\mathbb{R}}^{n}$. Set
$T_{\rm{max}}$, $T_{\rm A}$ and $T_{\rm B}$.

2: for $t=1,\cdots,T_{\rm{max}}$ do

3: for $t_{\rm B}=1,\cdots,T_{\rm B}$ do

4: // [Step 1 (B)] Perform output linear step

5: Calculate $\hat{\mathbf{p}}_{{\rm B}}(t)$ (7b) and ${\bm{\tau}}^{p}_{{\rm
B}}(t)$ (7a).

6: // [Step 2 (B)] Perform output linear step

7: Calculate $\hat{s}_{{\rm B},i}(t)$ (8a) and ${\tau}_{{\rm B},i}(t)$ (8b).

8: // [Step 3 (B)] Perform input linear step

9: Calculate the posterior means $\hat{\mathbf{r}}_{{\rm B}}(t)$ (11b) and
variances ${\bm{\tau}}^{r}_{{\rm B}}(t)$ (11a) of $\mathbf{b}$.

10: // [Step 4 (B)] Perform input nonlinear step

11: Calculate the posterior means (12a) and variances (12b) of $\mathbf{b}$.

12: end for

13: Set ${\mathbf{b}}_{{\rm B}}^{\rm ext}(t)=\hat{\mathbf{r}}_{{\rm B}}(t)$
and ${\mathbf{v}}_{{\rm B},b}^{{\rm ext}}(t)={\bm{\tau}}^{r}_{{\rm B}}(t)$.

14: for $t_{\rm A}=1,\cdots,T_{\rm A}$ do

15: // [Step 1 (A)] Perform output nonlinear step.

16: Calculate $\hat{\mathbf{s}}_{{\rm A}}(t)\in{\mathbb{R}}^{m}$ (16a) and
${\bm{\tau}}^{s}_{{\rm A}}(t)\in{\mathbb{R}}^{m}$ (16b).

17:
${\mathbf{z}}_{\text{B}}^{\text{post}}(t)={\text{E}}[{\mathbf{z}}|{\mathbf{y}},{\mathbf{z}}_{\text{A}}^{\text{ext}}(t),{\mathbf{v}}_{\text{A}}^{\text{ext}}(t)]$

18: // [Step 2 (A)] Perform input linear step

19: Calculate $\hat{\mathbf{r}}_{{\rm A}}(t)\in{\mathbb{R}}^{m}$ (17b) and
${\bm{\tau}}_{{\rm A}}^{r}(t)\in{\mathbb{R}}^{m}$ (17a).

20: // [Step 3 (A)] Perform input nonlinear step

21: Calculate the posterior means (19a) and variances (19b) of $\mathbf{x}$,

22: // [Step 4 (A)] Perform output linear step

23: Calculate ${\bm{\tau}}^{p}_{{\rm A}}(t)$ (19a) and variances
$\hat{\mathbf{p}}_{{\rm A}}(t)$ (20b)

24: end for

25: Set ${\mathbf{b}}_{{\rm A}}^{\rm ext}(t)=\hat{\mathbf{p}}_{{\rm A}}(t)$
and ${\mathbf{v}}_{{\rm A}}^{{\rm ext},b}(t)={\bm{\tau}}^{p}_{{\rm A}}(t)$.

26: end for

27: Return $\hat{\mathbf{x}}$.

As shown in [13], the extrinsic means ${\mathbf{b}}_{{\rm B}}^{\rm ext}(t)$
and variances ${\mathbf{v}}_{{\rm B}}^{{\rm ext},b}(t)$ can be regarded as the
pseudo observations and variances of $\mathbf{b}$ in module A, i.e.,

$\displaystyle\tilde{\mathbf{b}}(t)={\mathbf{Q}}{\mathbf{x}}+\tilde{\bm{\xi}}(t),$
(15)

where $\tilde{\mathbf{b}}(t)\triangleq{\mathbf{b}}_{{\rm B}}^{\rm ext}(t)$,
$\tilde{\bm{\xi}}(t)\sim\mathcal{N}({\mathbf{0}},{\rm diag}({\mathbf{v}}_{{\rm
B}}^{{\rm ext},b}(t)))$. Consequently, we could run standard AMP with $T_{\rm
A}$ iterations on this pseudo-linear model (15), as shown below.

  * •

[Step 1 (A)] For $i=1,2,\cdots,r$, $\hat{s}_{{\rm A},i}(t)$ and
$\tau^{s}_{{\rm A},i}(t)$ can be calculated as

$\displaystyle\hat{s}_{{\rm A},i}(t)=\frac{b_{{\rm B},i}^{\rm
ext}(t)-\hat{p}_{{\rm A},i}(t-1)}{v^{{\rm ext},b}_{{\rm B},i}(t)+\tau_{{\rm
A},i}^{p}(t-1)},$ (16a) $\displaystyle\tau^{s}_{{\rm A},i}(t)=\frac{1}{v_{{\rm
B},i}^{{\rm ext},b}(t)+\tau^{p}_{{\rm A},i}(t-1)}.$ (16b)

  * •

[Step 2 (A)] Perform input linear step to obtain $\hat{\mathbf{r}}_{{\rm
A}}(t)\in{\mathbb{R}}^{n}$ and ${\bm{\tau}}_{{\rm
A}}^{r}(t)\in{\mathbb{R}}^{n}$ as

$\displaystyle\bm{\tau}^{r}_{\rm A}(t)$ $\displaystyle=\left(|\mathbf{Q}^{\rm
T}|^{2}\bm{\tau}^{s}_{\rm A}(t)\right)^{-1},$ (17a)
$\displaystyle\hat{\mathbf{r}}_{\rm A}(t)$
$\displaystyle=\hat{\mathbf{x}}(t-1)+\bm{\tau}^{r}_{\rm
A}(t)\odot\left(\mathbf{Q}^{\rm T}\hat{\mathbf{s}}_{\rm A}(t-1)\right),$ (17b)

which can be characterized by a pseudo model

$\displaystyle\hat{\bm{r}}_{\rm A}(t)=\mathbf{x}+\tilde{\mathbf{w}}(t),$ (18)

where $\tilde{\mathbf{w}}(t)\sim\mathcal{N}(\mathbf{0},{\rm
diag}(\bm{\tau}^{r}_{\rm A}(t)))$.

  * •

[Step 3 (A)] Perform input nonlinear step in Module A to obtain the posterior
means and of variances of $\mathbf{x}$.

$\displaystyle\hat{x}_{j}(t)$ $\displaystyle={\rm E}[x_{j}|\hat{r}_{{\rm
A},j}(t),\tau^{r}_{{\rm A},j}(t)],$ (19a) $\displaystyle\tau^{x}_{j}(t)$
$\displaystyle={\rm Var}[x_{j}|\hat{r}_{{\rm A},j}(t),\tau^{r}_{{\rm
A},j}(t)]$ (19b)

for $j=1,2,\cdots,n$, where ${\rm E}[\cdot]$ and ${\rm Var}[\cdot]$ denotes
the posterior means and variances with respect to the likelihood calculated in
(18) and the prior $p({\mathbf{x}})$. Even if the nuisance parameters in
$p({\mathbf{x}})$ are unknown, EM algorithm [21] can be easily incorporated
into the GUAMP to learn them similarly as [22].

  * •

[Step 4 (A)] Perform output linear step to obtain ${\bm{\tau}}^{p}_{{\rm
A}}(t)\in{\mathbb{R}}^{r}$ and $\hat{\mathbf{p}}_{{\rm
A}}(t)\in{\mathbb{R}}^{r}$, i.e.,

$\displaystyle{\bm{\tau}}^{p}_{{\rm A}}(t)$
$\displaystyle=|{\mathbf{Q}}|^{2}{\bm{\tau}}^{x}(t),$ (20a)
$\displaystyle\hat{\mathbf{p}}_{{\rm A}}(t)$
$\displaystyle={\mathbf{Q}}\hat{\mathbf{x}}(t)-{\bm{\tau}}^{p}_{\rm
A}(t)\odot\hat{\mathbf{s}}_{{\rm A}}(t-1).$ (20b)

After running $T_{\rm A}$ iterations of AMP on module A, according to [13],
the extrinsic mean ${\mathbf{b}}_{{\rm A}}^{\rm ext}(t)$ and variance
${\mathbf{v}}_{{\rm A}}^{{\rm ext},b}(t)$ are

$\displaystyle{\mathbf{b}}_{{\rm A}}^{\rm ext}(t)=\hat{\mathbf{p}}_{{\rm
A}}(t),{\mathbf{v}}_{{\rm A}}^{{\rm ext},b}(t)={\bm{\tau}}^{p}_{{\rm A}}(t).$
(21)

Now we could again run the GAMP on module B. The pseudo code of GUAMP is
summarized as Algorithm 1.

##  3 Numerical Simulation

Fig. 2: The ${\rm dNMSE}(\hat{\mathbf{z}})$ ($1$ bit) ${\rm
NMSE}(\hat{\mathbf{z}})$ ($2$ bit) and versus the correlation coefficient
$\rho$. (a): ${\rm SNR}=20$ dB, (b): ${\rm SNR}=40$ dB. Fig. 3: The ${\rm
dNMSE}(\hat{\mathbf{z}})$ ($1$ bit) and ${\rm NMSE}(\hat{\mathbf{z}})$ ($2$
bit) versus the the number of iterations at ${\rm SNR}=20$ dB.

We verify the efficacy of GUAMP on quantized compressed sensing problem with
correlated measurement matrices $\bf{A}$. Specifically, consider the
correlated case where ${\mathbf{A}}$ is constructed as
${\mathbf{A}}={\mathbf{R}}_{L}{\mathbf{H}}{\mathbf{R}}_{R}$ [23], where
${\mathbf{R}}_{L}={\mathbf{R}}_{1}^{\frac{1}{2}}\in{\mathbb{R}}^{m\times m}$
and ${\mathbf{R}}_{U}={\mathbf{R}}_{2}^{\frac{1}{2}}\in{\mathbb{R}}^{n\times
n}$, the $(i,j)$th element of both ${\mathbf{R}}_{1}$ and ${\mathbf{R}}_{2}$
is $\rho^{|i-j|}$ and $\rho$ is termed as the correlation coefficient here,
${\mathbf{H}}\in{\mathbb{R}}^{m\times n}$ is a random matrix whose elements
are drawn i.i.d. from ${\mathcal{N}}(0,1)$. The elements $x_{i}$ of
${\mathbf{x}}$ are drawn from a Bernoulli Gaussian distribution
$p(x_{i})=(1-\lambda)\delta(x_{i})+\lambda{\mathcal{N}}(x_{i};0,1/\lambda)$
and $\lambda=0.1$. The $1$ bit (or probit) model ${\mathbf{y}}={\rm
sign}({\mathbf{z}}+{\mathbf{w}})$ and $2$ bit quantization models are
considered where ${\mathbf{z}}={\mathbf{A}}{\mathbf{x}}$ and
${\mathbf{w}}\sim{\mathcal{N}}({\mathbf{0}},\sigma^{2}{\mathbf{I}}_{m})$. The
number of measurements and the number of unknowns are $m=2048$ and $n=512$.
For $2$ bit quantization, the dynamic range of the quantizer is restricted to
$[-3\sigma_{z},3\sigma_{z}]$, where $\sigma_{z}^{2}$ is the variance of $z$
and is calculated to be $\sigma_{z}^{2}=\|{\mathbf{A}}\|_{\rm F}/m=n/m=0.25$
as we normalize the measurement matrix $\mathbf{A}$ such that
$\|{\mathbf{A}}\|_{\rm F}=n$, and the thresholds for $2$ bit quantization are
$\\{-1.5,0,1.5\\}$. The SNR is defined as ${\rm
SNR}=\|{\mathbf{A}}{\mathbf{x}}\|_{2}^{2}/(m\sigma^{2})$. The debiased
normalized mean squared error (dNMSE) ${{\rm
dNMSE}}(\hat{\mathbf{z}})\triangleq\underset{c}{\operatorname{min}}~{}{\|{\mathbf{z}}-c\hat{\mathbf{z}}\|_{2}^{2}}/{\|{\mathbf{z}}\|_{2}^{2}}$
and NMSE ${{\rm
NMSE}}(\hat{\mathbf{z}})\triangleq{\|{\mathbf{z}}-\hat{\mathbf{z}}\|_{2}^{2}}/{\|{\mathbf{z}}\|_{2}^{2}}$
are used for $1$ and $2$ bit quantization, respectively. The GUAMP is compared
with the GAMP and GVAMP. We set $T_{\rm B}=1$ and $T_{\rm A}=4$ (it is
numerically found $T_{\rm A}>1$ improves the robustness of GUAMP). All the
results are averaged over $500$ Monte Carlo (MC) trials. Our code is available
at https://github.com/RiverZhu/GUAMP.

Fig. 2 shows the average dNMSE ($1$ bit) and NMSE ($2$ bit) versus the
correlation coefficient $\rho$. It can be seen that when $\rho$ is small and
near zero, GUAMP performs the same as GVAMP, and GAMP. As $\rho$ increases,
both GVAMP and GAMP tend to diverge for moderate $\rho>0.05$. By contrast, the
proposed GUAMP achieves significantly better performances for large values of
$\rho$.

Fig. 3 plots the average dNMSE ($1$ bit) and NMSE ($2$ bit) versus iteration
for different correlation coefficients $\rho\in\\{0,0.15,0.35\\}$ at ${\rm
SNR}=20$ dB. It can be seen that at $\rho=0$, all the algorithms accurately
recover $\mathbf{x}$ and converge in about tens of iterations. For
$\rho=0.15$, GAMP fails and GVAMP starts to diverge, while at a larger value
of $\rho=0.35$, both GAMP and GVAMP fail completely while GUAMP converges
remarkably well, showing that GUAMP is more robust than both GVAMP and GAMP
under correlated matrices.

##  4 Conclusion

In this paper, inspired from the UAMP [16] and a unified inference framework
[13], we propose a novel algorithm called generalized UAMP (GUAMP) to address
the GLM inference problem under general measurement matrices. Numerical
experiments demonstrate that GUAMP significantly outperforms GAMP and GVAMP
under correlated measurement matrices. In current implementation of GUAMP, no
damping is used although it is believed to further improve its robustness. The
adaptive step scheme [24] is also expected to improve GUAMP’s convergence and
thus it is interesting to formulate the cost function and unveil the relations
between GUAMP and the alternating direction method of multipliers (ADMM),
which is left as an important future work.

## References

  * [1] P. McCullagh and J. A. Nelder, Generalized linear models.  Routledge, 2019.
  * [2] Y. Kabashima, “A cdma multiuser detection algorithm on the basis of belief propagation,” Journal of Physics A: Mathematical and General, vol. 36, no. 43, p. 11111, 2003.
  * [3] D. L. Donoho, A. Maleki, and A. Montanari, “Message-passing algorithms for compressed sensing,” Proceedings of the National Academy of Sciences, vol. 106, no. 45, pp. 18914–18919, 2009.
  * [4] S. Rangan, “Generalized approximate message passing for estimation with random linear mixing,” in 2011 IEEE International Symposium on Information Theory Proceedings, pp. 2168–2172, IEEE, 2011.
  * [5] M. Al-Shoukairi, P. Schniter, and B. D. Rao, “A gamp-based low complexity sparse bayesian learning algorithm,” IEEE Transactions on Signal Processing, vol. 66, no. 2, pp. 294–308, 2017.
  * [6] C. K. Thomas and D. Slock, “Posterior variance predictions in sparse bayesian learning under approximate inference techniques,” in 2020 54th Asilomar Conference on Signals, Systems, and Computers, pp. 1375–1379, IEEE, 2020.
  * [7] C. K. Thomas and D. Slock, “Save-space alternating variational estimation for sparse bayesian learning,” in 2018 IEEE Data Science Workshop (DSW), pp. 11–15, IEEE, 2018.
  * [8] S. Rangan, P. Schniter, and A. K. Fletcher, “Vector approximate message passing,” IEEE Transactions on Information Theory, vol. 65, no. 10, pp. 6664–6684, 2019.
  * [9] J. Ma and L. Ping, “Orthogonal amp,” IEEE Access, vol. 5, pp. 2020–2033, 2017.
  * [10] P. Schniter, S. Rangan, and A. K. Fletcher, “Vector approximate message passing for the generalized linear model,” in 2016 50th Asilomar Conference on Signals, Systems and Computers, pp. 1525–1529, IEEE, 2016.
  * [11] T. P. Minka, A family of algorithms for approximate Bayesian inference.  PhD thesis, Massachusetts Institute of Technology, 2001.
  * [12] X. Meng, S. Wu, L. Kuang, and J. Lu, “An expectation propagation perspective on approximate message passing,” IEEE Signal Processing Letters, vol. 22, no. 8, pp. 1194–1197, 2015.
  * [13] X. Meng, S. Wu, and J. Zhu, “A unified bayesian inference framework for generalized linear models,” IEEE Signal Processing Letters, vol. 25, no. 3, pp. 398–402, 2018.
  * [14] M. T. Ivrlac, W. Utschick, and J. A. Nossek, “Fading correlations in wireless mimo communication systems,” IEEE Journal on selected areas in communications, vol. 21, no. 5, pp. 819–828, 2003.
  * [15] D. Chizhik, J. Ling, P. W. Wolniansky, R. A. Valenzuela, N. Costa, and K. Huber, “Multiple-input-multiple-output measurements and modeling in manhattan,” IEEE Journal on Selected Areas in Communications, vol. 21, no. 3, pp. 321–331, 2003.
  * [16] Q. Guo and J. Xi, “Approximate message passing with unitary transformation,” arXiv preprint arXiv:1504.04799, 2015.
  * [17] M. Luo, Q. Guo, M. Jin, Y. C. Eldar, D. Huang, and X. Meng, “Unitary approximate message passing for sparse bayesian learning,” IEEE transactions on signal processing, vol. 69, pp. 6023–6039, 2021.
  * [18] Z. Yuan, Q. Guo, and M. Luo, “Approximate message passing with unitary transformation for robust bilinear recovery,” IEEE Transactions on Signal Processing, vol. 69, pp. 617–630, 2020.
  * [19] J. Zhu, “A comment on the” a unified bayesian inference framework for generalized linear models”,” arXiv preprint arXiv:1904.04485, 2019.
  * [20] Q. Guo and D. D. Huang, “A concise representation for the soft-in soft-out lmmse detector,” IEEE communications letters, vol. 15, no. 5, pp. 566–568, 2011.
  * [21] A. P. Dempster, N. M. Laird, and D. B. Rubin, “Maximum likelihood from incomplete data via the em algorithm,” Journal of the Royal Statistical Society: Series B (Methodological), vol. 39, no. 1, pp. 1–22, 1977.
  * [22] J. P. Vila and P. Schniter, “Expectation-maximization gaussian-mixture approximate message passing,” IEEE Transactions on Signal Processing, vol. 61, no. 19, pp. 4658–4672, 2013.
  * [23] D.-S. Shiu, G. J. Foschini, M. J. Gans, and J. M. Kahn, “Fading correlation and its effect on the capacity of multielement antenna systems,” IEEE Transactions on communications, vol. 48, no. 3, pp. 502–513, 2000.
  * [24] S. Rangan, P. Schniter, E. Riegler, A. K. Fletcher, and V. Cevher, “Fixed points of generalized approximate message passing with arbitrary matrices,” IEEE Transactions on Information Theory, vol. 62, no. 12, pp. 7464–7474, 2016.