[1]Sohum Thakkar

[3]André J. Ferreira–Martins

1]QC Ware Corp, Palo Alto, USA and Paris, France 2]IRIF, Université Paris Cité
and CNRS, France 3]Itaú Unibanco, São Paulo, Brazil *]Corresponding authors.
Emails<EMAIL_ADDRESS><EMAIL_ADDRESS>
# Improved Financial Forecasting via Quantum Machine Learning

Skander Kazdaghli  Natansh Mathur  Iordanis Kerenidis  Samurai Brito  [  [  [
[

###### Abstract

Quantum algorithms have the potential to enhance machine learning across a
variety of domains and applications. In this work, we show how quantum machine
learning can be used to improve financial forecasting. First, we use classical
and quantum Determinantal Point Processes to enhance Random Forest models for
churn prediction, improving precision by almost 6%. Second, we design quantum
neural network architectures with orthogonal and compound layers for credit
risk assessment, which match classical performance with significantly fewer
parameters. Our results demonstrate that leveraging quantum ideas can
effectively enhance the performance of machine learning, both today as
quantum-inspired classical ML solutions, and even more in the future, with the
advent of better quantum hardware.

###### keywords:

Computational Finance, Machine Learning, Quantum Computing, Credit Risk, Churn
Prediction

##  1 Introduction

Quantum computing is a rapidly evolving field that promises to revolutionize
various domains, and finance is no exception. There is a variety of
computationally hard financial problems for which quantum algorithms can
potentially offer advantages [1, 2, 3, 4], for example in combinatorial
optimization [5, 6], convex optimization [7, 8], monte carlo simulations [9,
10, 11], and machine learning [12, 13, 14, 15].

In this work, we explore the potential of quantum machine learning methods in
improving the performance of forecasting in finance, specifically focusing on
two use cases within the business of Itaú Unibanco, the largest bank in Latin
America.

In the first use case, we aim to improve the performance of Random Forest
methods for churn prediction. We introduce quantum algorithms for
Determinantal Point Processes (DPP) sampling [16], and develop a method of DPP
sampling to enhance Random Forest models. We evaluate our model on the churn
dataset using classical DPP sampling algorithms and perform experiments on a
scaled-down version of the dataset using quantum algorithms. Our results
demonstrate that, in the classical setting, the proposed algorithms outperform
the baseline Random Forest in precision, efficiency, and bottom line, and also
offer a precise understanding of how quantum computing can impact this kind of
problem in the future. The quantum algorithm run on an IBM quantum processor
gives similar results as the classical DPP on small batch dimensions but
falters as the dimensions grow bigger due to hardware noise.

In the second use case, we aim to explore the performance of neural network
models for credit risk assessment by incorporating ideas from quantum compound
neural networks [17]. We start by using quantum orthogonal neural networks
[17], which add the property of orthogonality for the trained model weights to
avoid redundancy in the learned features [18]. These orthogonal layers, which
can be trained efficiently on a classical computer, are the simplest case of
what we call compound neural networks, which explore an exponential space in a
structured way. For our use case, we design compound neural network
architectures that are appropriate for financial data. We evaluate their
performance on a real-world dataset and show that the quantum compound neural
network models both have far fewer parameters and achieve better accuracy and
generalization than classical fully-connected neural networks.

This paper is organized as follows: In sections 2-4, we focus on the churn
prediction use case and present the DPP-based quantum machine learning
methods. In sections 5-7, we present quantum neural network models for risk
assessment. Finally, in section 8, we conclude the paper and discuss potential
future research directions.

##  2 DPP-Enhanced Random Forest Models for Churn Prediction

###  2.1 DPP-Random Forest Model

The Random Forest algorithm was introduced in 2001 by Leo Brieman [19] and has
easily become one of the most popular supervised machine learning algorithms
in use. It consists of an ensemble of decision trees, each trained on a
uniform subsample of rows and columns from the dataset.

In this section, we propose an extension of the Random Forest, called the DPP-
Random Forest (DPP-RF), which utilizes Determinantal Point Processes (DPPs)
instead of uniform sampling to subsample rows and columns for individual
decision trees. In the original RF algorithm, subsampling makes the model more
robust to variance in the training data; however, the use of uniform sampling
runs the risk of improperly representing the dataset and missing under-sampled
areas. DPP sampling better preserves the diversity of the dataset, and
corrects for sampling bias [20]. We first introduce the theory and techniques
of DPP sampling, then present the algorithm.

###  2.2 Determinantal Point Processes

We will now introduce the Determinantal Point Process (DPP), which lies at the
core of the methodology behind our solution to the churn problem. DPPs are a
class of probabilistic models that can be used to sample diverse subsets of
items from a larger set. They were first formalized by Macchi in 1975 as a way
to model fermions in quantum mechanics [21]. More recently, these models are
showing increasing promise in the context of machine learning [20], where they
can be used for a variety of tasks, such as building unbiased estimators for
linear regression [22], performing monte carlo estimation [23], and promoting
diversity in model outputs [24].

####  2.2.1 Definitions

A point process $P$ on a set $Y$ is a probability measure over the subsets of
$Y$. Sampling from a point process on $Y$ will produce some subset $S\subseteq
Y$ with probability $P(S)$. A repulsive point process is a point process in
which points that are more similar to each other are less likely to be
selected together.

A determinantal point process (DPP) is a particular case of a repulsive point
process, in which the selection probability of a subset of items $T\subseteq
Y$ is given by a determinant. Given a real, symmetric $n\times n$ matrix $K$
indexed by the elements of $Y$:

$P\\{T\subseteq S\\}=\det(K_{T,T})\ ,$

where $K_{T,T}$ denotes the $|T|\times|T|$ submatrix indexed by the set $T$
and $n$ is the cardinality of $Y$. In other words, the marginal distribution
$P\\{T\subseteq S\\}$ is defined by the subdeterminants of K.

The above is the most general definition, but in machine learning, we
typically focus on a slightly more restrictive class of DPPs called
$L$-ensembles. In $L$-ensembles, the whole distribution, not just the
marginals, is given by the subdeterminant of a real, symmetric $n\times n$
matrix $L$.

$P\\{S\\}\propto\det(L_{S,S})\ .$

Just like $K$, $L$ is indexed by the elements of $Y$. Because of some
convenient properties of the determinant [20], we can explicitly write down
the distribution of an $L$-ensemble:

$P\\{S\\}=\frac{\det(L_{S,S})}{\det(L+I)}\ .$

In machine learning literature, DPPs are typically defined over a set of
points $\mathbf{X}$, with each item $\mathbf{x}_{i}$ a row in the data matrix
$\mathbf{X}$. If we preprocess $\mathbf{X}$ such that its columns are
orthonormal and choose $\mathbf{L}$ to be the inner-product similarity matrix,
i.e. $L=\mathbf{X}\mathbf{X}^{T}$, then the distribution becomes even simpler
to write down. Instead of explicitly computing the $\mathbf{L}$ matrix, we can
write the distribution in terms of the data matrix $\mathbf{X}$ itself,
courtesy of the Cauchy-Binet formula,

$P\\{S\\}=\frac{\det(\mathbf{X}_{S})^{2}}{\det(\mathbf{X}\mathbf{X}^{T}+\mathbf{I})}\
.$ (1)

Moreover, the distribution will almost surely produce samples with size $d$,
the rank of the orthogonalized data matrix $X$. This kind of DPP is denoted
$d$-DPP. We will focus here on an application of sampling from a $d$-DPP from
a data matrix X.

####  2.2.2 Unbiased Least-Squares Regression

One unique feature of the DPP compared to i.i.d sampling techniques is that it
can lead to provably unbiased estimators for least-squares linear regression
[25] [26]. Given an $n\times d$ data matrix $\mathbf{X}$ and a target vector
$\mathbf{y}\in\mathbb{R}^{n}$, where $n\gg d$, we wish to approximate the
least squares solution
$\mathbf{w}^{\ast}=\mathrm{argmin}_{\mathbf{w}}||\mathbf{X}\mathbf{w}-\mathbf{y}||$.
$\mathbf{w}^{\ast}$ represents the best-fit parameters to a linear model to
predict $\mathbf{y}$.

Surprisingly, if we sample $d$ points $S$ from DPP($\mathbf{X}\mathbf{X}^{T}$)
and solve the reduced system of equations
$\mathbf{y}_{S}=\mathbf{X}_{S}\mathbf{w}$, we get an unbiased estimate of
$\mathbf{w}^{\ast}$. Formally, if $S\sim
d\text{-DPP}_{\mathbf{L}}(\mathbf{X}\mathbf{X}^{T})$,

$\mathop{\mathbb{E}}[\mathbf{X}_{S}^{-1}\mathbf{y}_{S}]=\mathrm{argmin}_{\mathbf{w}}||\mathbf{X}\mathbf{w}-\mathbf{y}||=\mathbf{w}^{\ast}\
.$ (2)

This allows us to create an ”ensemble” of unbiased linear regressors, each
trained on a DPP sample. In some regard, this was the inspiration for trying
an ensemble of _decision trees_ trained on DPP samples, as detailed in Sec.
2.1.

####  2.2.3 Algorithms for Sampling

There are several efficient algorithms for sampling from DPPs and computing
their properties. The naive sampling method — calculating all subdeterminants
and performing $l2$ sampling — takes exponential time. The first major leap in
making DPP sampling feasible on today’s computers was the ”spectral method”
[27, 28]. This algorithm performs an eigendecomposition of the kernel matrix
before applying a projection-based iterative sampling approach. Thus, the
first sample takes $O(nd^{2})$ time, and subsequent samples take $O(d^{3})$.

Monte Carlo methods have been proposed to approximate the DPP distribution
[29, 30], though they are not exact, and are often still prohibitively slow
with a runtime of $O(n\text{poly}(d))$ per sample.

In a counter-intuitive result, [31] and [32] proposed methods that avoid
performing the full DPP sampling procedure on large parts of the basis set.
This approach resulted in a significant reduction in runtime, making DPPs more
practical for mid-to-large-scale datasets. These techniques allow exact
sampling of subsequent $d$-DPP samples in $O(\text{poly}(d))$, independent of
the size of the full basis set $n$. Many of these algorithms are implemented
in the open-source DPPy library[33], which we used in the experiments in this
paper.

Recent work has shown that quantum computers are in principle able to sample
from DPPs in even lower complexity in some cases. We describe this quantum
algorithm in Sec. 4. This and several other algorithms which arise from the
techniques introduced in [16] are a budding area of research in the quantum
computing space, and will hopefully inspire more applications like the one we
describe in this paper. For example, in [34], DPPs and deterministic DPPs were
used to improve the methods for the imputation of clinical data.

###  2.3 DPP-RF Algorithm Outline

In principle, the DPP-RF uses DPP sampling on the whole dataset to select
diverse subsets of data on which to train decision trees. However, sampling
from a DPP on large datasets (like the entire churn dataset of 174,000 points)
can take copious time, especially when using current open-source
implementations of DPP sampling. To be able to test these techniques quickly,
a novel sampling procedure was developed which preserves many of the benefits
of DPPs, but does not require sampling from the full dataset. The procedure
can be summarized as follows:

  1. 1.

Divide the training set uniformly into smaller batches;

  2. 2.

Sample $S_{1}\sim d$-DPP$(\mathbf{X}_{batch}\mathbf{X}_{batch}^{T})$ data
points from every batch;

  3. 3.

Sample $S_{2}\sim d$-DPP$(\mathbf{X}_{S_{1}}^{T}\mathbf{X}_{S_{1}})$ features;

  4. 4.

Train a first group $G_{1}$ of $N_{1}$ decision trees on these small patches
of data;

  5. 5.

Aggregate the patches of data resulting from step 2 to create a larger dataset
$\mathbf{X}_{agg}$;

  6. 6.

Repeat for $N_{2}$ times: sample $S_{3}\sim
d$-DPP$(\mathbf{X}_{agg}^{T}\mathbf{X}_{agg})$ features to create a long
matrix;

  7. 7.

Train a second group $G_{2}$ of $N_{2}$ decision trees on these new datasets;

  8. 8.

Combine $G_{1}$ and $G_{2}$ by aggregating them to make predictions (similar
to the classical Random Forest algorithm).

Figure 1: Steps 1 to 7 of the DPP-Random Forest algorithm.

##  3 Classical DPP-RF Results

The DPP-RF algorithm was designed for the purpose of predicting customer churn
in the bank. In this section, we define this use case and present the results.
In addition, we benchmark our proposed DPP-RF method by constructing models on
public tabular classification tasks.

###  3.1 Use Case Introduction: Churn Prediction

Churn, defined as a customer withdrawing more than a certain amount of money
in a single month, is a significant concern for retail banks. Our objective is
to predict which customers are most likely to churn in the next three months
using customer data from the previous six months.

The primary dataset used in this study consists of 304,000 datapoints, with
153 features for each datapoint. Each datapoint represents a banking customer
at a particular month in time, with the features representing various aspects
of their activity over the previous six months. The target variable is a
binary flag indicating whether or not the customer churned in at least one of
the following three months. The data was anonymized and standardized before
being split into training and test sets based on time period, with 130,000
datapoints being set aside as the test set and 174,000 datapoints used for
training. The data was split in a way that did not produce any significant
covariate shift between the train and test sets.

With the end goal of preventing churn, the model works by flagging customers
with the highest risk of potential churn. For these flagged customers, the
bank can deploy a representative to intervene and better understand their
needs. However, resource limitations make it necessary to flag a relatively
small number of customers with high confidence. The focus of this exploration
was to reduce false positives in the flagged customers to increase the
efficiency of bank interventions. In terms of the precision-recall trade-off,
our model should be tuned to provide the highest possible precision for low
recall values. Despite this simplification to a classification problem, the
primary business KPI is the amount of withdrawal money correctly captured by
the model, as discussed more in Sec. 3.4.

This use case already had a solution in production: a Random Forest classifier
[19], whose performance was used as a benchmark. The model in production
already captured a significant amount of churn, but there was clear room for
improvement in the amount of withdrawals captured (see Fig. 3). Moreover,
given the large number of customers in the dataset and the relative
homogeneity of the population of interest, there existed an opportunity to
employ techniques that explicitly try to explore diversity in the data.

We focused on three key performance indicators (KPIs): the precision-recall
curve, the training time and the bottom line.

###  3.2 Precision-recall

To evaluate the performance of our proposed method, we optimized
hyperparameters and measured the precision for a low fixed recall (6% in this
case). As seen in Figure 2, our method showed an improvement in precision from
$71.6\%$ for the benchmark model to $77.5\%$ with the new model. Our method
also provided similar improvements in precision for the relevant range of
small recall.

Figure 2: Precision-recall curve for the test set. Using DPP with the Random
Forest algorithm shows an improvement of $5.9\%$.

###  3.3 Training Time

The DPP-RF model has a longer training time compared to the traditional random
forest on a classical computer: it took 54 minutes to train the model with the
best hyperparameters using, compared to 311 seconds for the benchmark model.
The models were trained on a computer with an Intel© Core™ i5-8350U CPU
running at 1.70 GHz, 24 GB of RAM and Windows 10 version 21H2, compilation
19044.2604.

The computational bottleneck in this algorithm is the DPP sampling. Instead of
simulating quantum DPP circuits (which is infeasible for large datasets), we
used a classical SVD-based sampling algorithm [28] implemented in the dppy
library [33]. We believe that improved classical sampling techniques [32] and
future quantum techniques (Sec. 4) can reduce the runtime dramatically.

Hyperparameters were selected using a grid search with 5-fold cross-validation
over the typical RF parameters n_estimators, max_depth, min_samples_leaf,
min_samples_split, max_features, max_samples, and the DPP-RF-specific
parameter batch_size. The training time of a DPP-RF model depends heavily on
this batch_size parameter, which is the size of the batches from which we take
DPP samples. Choosing a batch size higher than 1000 can increase runtime
dramatically. Thus, in our hyperparameter search, we limited the batch size to
less than 1000.

Within the bank, the churn model is retrained just once every few months, so
the training time was not prohibitive. However, faster sampling algorithms
still serve to increase the range of feasible hyperparameters (especially the
batch size).

###  3.4 Bottom Line — Withdrawals Captured

From a business perspective, the most direct indicator of the success of the
model is the amount of assets under management (AUM) that can be salvaged via
interventions. Thus, we evaluated the amount of money withdrawn every month by
the 500 customers flagged by the model, i.e., the 500 customers which had the
highest predicted probability of churning in one of the following 3 months. As
seen in figures 3, 4 and 5, our model showed substantial overall improvements.
The true financial impact of these predictions is dependent on the success of
the interventions as well as the bank’s profit-per-dollar-AUM.

Figure 3: Classical benchmark (BM) vs DPP-RF solution: money withdrawn per
month by the flagged 500 customers, comparing the benchmark model (blue line)
to the DPP-RF one (orange line). On the y axis, we have monetary values (not
shown). The green line represents the total amount of money withdrawn by all
customers in each month. The purple line is the sum of the 500 largest
withdrawals, which is the maximum value that the model could capture. The red
line represents the withdrawals captured by randomly flagging 500
observations. The y-axis units are omitted for confidentiality. Figure 4:
Classical benchmark vs DPP-RF solution - percentage of total withdrawals
captured per month, that is, relative to the green line in Fig. 3. On average
over the 11 test months, the BM model captures 61.42% of the total, whilst the
DPP-RF model captures 62.77% — an improvement of 1.35%. Figure 5: Classical
benchmark vs DPP-RF solution - the percentage of maximum money possible to be
captured (given $\text{n\\_flags}=500$ customers flagged every month), that
is, relative to the purple line in Fig. 3. On average over the 11 test months,
the BM model captures 69.18% of the total, whilst the DPP-RF model captures
70.72% — an improvement of 1.54%.

###  3.5 Summary of Results

The proposed DPP Random Forest model provides significant improvements in
precision and bottom line, while taking significantly longer to train. The
results are summarized in the following table.

Table 1: Summarized comparison between models. Metric |  Benchmark model |  Proposed model
---|---|---
Precision |  71.6% |  77.5%
% total withdrawals captured |  61.42% |  62.77%
% maximum possible withdrawals captured |  69.18% |  70.72%
Train time |  311s |  54 min
|  |

###  3.6 Further Benchmarks

We further benchmarked our model on various classification datasets from
OpenML. All except one (madelon) of these datasets were used in [35] and
preprocessed accordingly. They were chosen to be representative of a wide
variety of classification tasks. Each dataset was split into train,
validation, and test sets. For each model, 400 sets of hyperparameters were
randomly chosen and evaluated on the validation set. Both models used the same
hyperparameter space, except for the addition of the batch_size parameter for
the DPP-RF. The hyperparameters which gave the best results on the validation
set were evaluated on the test set, and the results are reported below. Models
were evaluated with the ROC-AUC metric111The area under the receiver operating
characteristic curve (ROC-AUC) is a common metric for two-class classification
tasks, and evaluates the ability of the model to produce a proper ranking of
datapoints by likelihood of being class 1..

Table 2: Comparison of DPP Random Forest and Random Forest models for different datasets. The results are reported via the ROC-AUC metric. Dataset |  Random Forest |  DPP Random Forest
---|---|---
madelon |  0.916 |  0.941
credit-default |  0.856 |  0.856
house-pricing |  0.948 |  0.939
jannis |  0.866 |  0.870
eye movements |  0.704 |  0.710
bank-marketing |  0.881 |  0.881
wine |  0.901 |  0.906
california |  0.962 |  0.963
|  |

##  4 Quantum DPP-RF

###  4.1 Quantum Circuits for Determinantal Point Processes

Classical DPP sampling algorithms have improved significantly since their
inception, but it still remains infeasible to sample from large datasets like
ours. Recent work by Kerenidis and Prakash [16] has shown that a quantum
computer can more natively perform DPP sampling, achieving a gate complexity
of $O(nd)$ and a circuit depth of $O(d\log(n))$ for an orthogonal matrix of
size $n\times d$. The classical time complexity for sampling is $O(d^{3})$
[28]. Note that when $n$ is very large, then one can reduce the number of rows
to $O(d^{2})$ before performing the sampling [32].

For a thorough review of the quantum methods and circuits, we refer the reader
to [16]. The circuit is described in brief below.

Given an orthogonal matrix
$\mathbf{X}=(\mathbf{x}^{1},\mathbf{x}^{2},\dots,\mathbf{x}^{n})\in\mathbb{R}^{n\times
d}$, the quantum DPP circuit applied on $\mathbf{X}$ performs the following
operation:

$\mathcal{D}(\mathbf{X})|0^{n}\rangle=\sum_{\begin{subarray}{c}|S|=d\\\
S\in\\{0,1\\}^{n}\end{subarray}}\det(\mathbf{X}_{S})|1_{S}\rangle\ ,$ (3)

where $\mathbf{X}_{S}$ is the $\mathbb{R}^{d\times d}$ submatrix obtained
after sampling the rows of $\mathbf{X}$ indexed by $S$; $1_{S}$ is the
characteristic vector of $S$ (with 1’s in the positions indexed by the
elements of S) and $\mathcal{D}(\mathbf{X})$ represents the quantum $d$-DPP
circuit, as detailed below.

Figure 6: Clifford Loader circuit $\mathcal{C}(x)$ for $x\in\mathbb{R}^{8}$.

Figure 7: DPP circuit as a series of Clifford Loaders.

Thus, the probability of sampling $S$, i.e., of measuring $|1_{S}\rangle$, is:
$Pr(S)=\det(\mathbf{X}_{S})^{2}=\det(L_{S,S})$, where
$L=\mathbf{X}\mathbf{X}^{T}$. This draws the link between the quantum
determinantal sampling circuit and the classical $d$-DPP model as seen in eq.
1.

To construct the quantum $d$-DPP circuit, we need to first introduce a circuit
known as a Clifford loader, which performs the following operation:

$\mathcal{C}(\mathbf{x})=\sum_{i=1}^{n}x_{i}Z^{i-1}XI^{n-i},\quad\text{for}\quad\mathbf{x}\in\mathbb{R}^{n}\
.$ (4)

The Clifford loader was shown to have a log-depth circuit in [16], and is
shown for $n=8$ in fig 7, in which the gates represented by vertical lines are
RBS gates — parameterized, hamming weight preserving two-qubit gates.

The full quantum $d$-DPP circuit is a series of $d$ Clifford loaders, one for
each orthogonal column of $\mathbf{X}$:

$\mathcal{D}(\mathbf{X})=\mathcal{C}(\mathbf{x}^{1})\mathcal{C}(\mathbf{x}^{2})\dots\mathcal{C}(\mathbf{x}^{d})\
.$ (5)

An example of a $d$-DPP circuit as a series of Cliffords for $n=4$ is shown in
fig. 7.

###  4.2 Hardware Experiment Results

As a hardware experiment, we aimed to implement a simplified version of our
algorithm on a quantum processor. We chose to use the ”ibmq_guadelupe”
16-qubit chip, which is only capable of running small quantum DPP circuits for
matrices of certain dimensions, such as $(4,2),(5,2),(5,3),(6,2),(8,2)$. As a
result, we had to reduce the size of our problem.

To accomplish this, we defined reduced train/test sets: a train set of
$\sim$1000 points from 03/2019 and a test set of $\sim$10000 points from
04/2019. The quantum hardware-ready simplified algorithm is outlined in Figure
8. It includes the following steps:

  1. 1.

Applying PCA to reduce the number of columns from 153 to $d=2,3$;

  2. 2.

Dividing the dataset into batches of $n=4,5,6,8$ points;

  3. 3.

Sampling $S\sim d$-DPP$(\mathbf{X}_{batch}\mathbf{X}_{batch}^{T})$ rows from
each batch, resulting in small $d\times d$ patches of data;

  4. 4.

Aggregating these patches to form a larger dataset, then training one decision
tree on this dataset.

Figure 8: Quantum hardware-ready procedure for DPP sampling.

We repeated this process for a number of trees and estimated the F1222We chose
the F1 score as the evaluation metric for two reasons. Firstly, a single
decision tree, unlike the random forest, does not provide an estimate of the
likelihood, which is required for the computation of ROC-AUC as we had used
before. Secondly, we had an imbalanced dataset and thus needed a metric that
balanced precision and recall. The F1 score is particularly effective when
used in these scenarios. score for every tree. We then compared the results
for different sampling methods: uniform sampling, quantum DPP sampling using a
simulator, and quantum DPP sampling using a quantum processor.

Figure 9: Diagonal Clifford loader.

Figure 10: Semi-diagonal Clifford loader.

The IBM quantum computer only allows using RBS gates on adjacent qubits, so we
cannot use the circuit described in section 4.1. Instead, we use two different
Clifford loader architectures which only use adjacent-qubit connectivity,
visualized in fig 10. The diagonal Clifford loader is explained in [16], and
the semi-diagonal loader is a modification that halves the circuit depth. As
an error mitigation measure, we disregarded all results that did not have the
expected hamming weight ($d$). The results are shown in the violin plots in
Figures 11 and 12.

Figure 11: Decision trees performance using Quantum DPP sampling with diagonal
Clifford loaders. Figure 12: Decision trees performance using Quantum DPP
sampling with semi-diagonal Clifford loaders.

The results indicate that for small matrix dimensions — up to (6,2) — the IBM
quantum processor gives results similar to the ones achieved with the
simulator. However, as the dimensions grow bigger, the samples from the
quantum DPP circuits lead to worse classifier performance. This highlights the
limitations of the available quantum hardware, which is prone to errors.

##  5 Quantum Neural Networks for credit risk assessment

###  5.1 Quantum Neural Networks with Orthogonal and Compound Layers

In recent years, variational/parameterized quantum circuits [36] have become
very prominent as NISQ-friendly QML techniques. When applied to classification
problems, they are commonly known as Variational Quantum Classifiers (VQC)
[37]. The quantum circuits associated with VQCs may be schematically thought
of as composed of three layers: the feature map $\mathcal{U}_{\Phi(\vec{x})}$,
which encodes classical data $\vec{x}$ into quantum states; the variational
layer $W(\bm{\theta})$, which is the part of the circuit parameterized by a
set of parameters $\bm{\theta}$ which are learned in the training process; and
finally, the measurement layer, which measures the quantum registers and
produces classical information used in training and inference.

The feature map and variational layers can take different forms, called
_ansätze_ , consisting of many possible different quantum gates in different
configurations. Such immense freedom raises an important question: how should
one choose an architecture for a given problem, and can it be expected to
yield a quantum advantage? This question is of major practical importance, and
although benchmark results have been shown for very particular datasets [37,
38], there is little consensus on which ansätze are good choices for machine
learning.

In our work, we use quantum neural networks with orthogonal and compound
layers. Although these neural networks roughly match the general VQC
construction, they produce well-defined linear algebraic operations, which not
only makes them much more interpretable but gives us the ability to analyze
their complexity and scalability. Because we understand the actions of these
layers precisely, we are able to identify instances for which we can design
efficient classical simulators, allowing us to classically train and test the
models on real-scale datasets.

A standard feed-forward neural network layer modifies an input vector by first
multiplying it by a weight matrix and then applying a non-linearity to the
result. Feed-forward neural networks usually use many such layers and learn to
predict a target variable by optimizing the weight matrices to minimize a loss
function. Enforcing the orthogonality of these weight matrices, as proposed in
[39], brings theoretical and practical benefits: it reduces the redundancy in
the trained weights and can avoid the age-old problem of vanishing gradients.
However, the overhead of typical projection-based methods to enforce
orthogonality prevents mainstream adoption.

In [17], an improved method of constructing orthogonal neural networks using
quantum ideas was developed. We describe it below in brief.

###  5.2 Data Loaders

In order to perform a machine learning task with a quantum computer, we need
to first load classical data into the quantum circuit.

#### Unary data loading circuits

Figure 13: Three possible unary data loaders for $d$-dimensional vectors
($d=8$). From left to right: the parallel, diagonal, and semi-diagonal
circuits have respectively a circuit depth of $log(d)$, $d$, and $d/2$. The X
gate represents the Pauli X gate, and the vertical lines represent $RBS$ gates
with tunable parameters.

The first way we will load classical data is an example of _amplitude
encoding_ , which means that we load the (normalised) vector elements as the
amplitudes of a quantum state. In [40], three different circuits to load a
vector $\bm{x}\in\mathbb{R}^{d}$ using $d-1$ gates are proposed. The circuits
range in depth from $O(log(d))$ to $O(d)$, with varying qubit connectivity
(see Fig.13). They use the _unary_ amplitude encoding, where a vector
$\bm{x}=(x_{1},\cdots,x_{d})$ is loaded in the quantum state
$\ket{\bm{x}}=\frac{1}{\|x\|}\sum_{i=1}^{d}x_{i}\ket{e_{i}}$, where
$\ket{e_{i}}$ is the quantum state with all qubits in $\ket{0}$ except the
$i^{th}$ qubit in state $\ket{1}$ (e.g., $\ket{e_{3}}=\ket{00100000}$). The
circuit uses $RBS$ gates: a parameterized two-qubit hamming weight-preserving
gate implementing the unitary given by Eq.6:

$RBS(\theta)=\left(\begin{array}[]{cccc}1&0&0&0\\\
0&\cos\theta&\sin\theta&0\\\ 0&-\sin\theta&\cos\theta&0\\\
0&0&0&1\end{array}\right)\ .$ (6)

The parameters $\theta_{i}:i\in\\{1,...,d-1\\}$ of the $d-1$ RBS gates are
classically pre-computed to ensure they encode the correct vector
$\ket{\mathbf{x}}$.

####  $RY$-loading circuits

We will also use data loading procedures beyond the unary basis. In
particular, for a normalized input vector $\bm{x}\in\mathbb{R}^{d}$, we use
$d$ qubits, where on each of the qubits, we apply an $RY(\theta)$ rotation
gate where the angle parameter on the $i^{th}$ qubit is $\theta_{i}=2\pi
x_{i}$, according to Eq. 7. Multiplication with $2\pi$ allows us to cover the
entire range of the $\sin$ and $\cos$ functions. This technique loads the data
in the entire $2^{d}$-dimensional Hilbert space encompassing all the hamming
weights from $0$ to $d$. This loading technique has constant depth independent
of $d$, and we refer to it as the _RY loading_ , whose circuit for $d=8$ is
illustrated in Fig. 14.

$RY(\theta)\ket{0}=\cos{\frac{\theta}{2}}\ket{0}+\sin{\frac{\theta}{2}}\ket{1}$
(7)

####  $H$-loading circuits

Lastly, we define a different technique for loading the data in the entire
$2^{d}$-dimensional Hilbert space, which loads the vector in the unary basis
and then applies a Hadamard gate on each qubit. This operation applies a
Fourier transform on $\mathbb{Z}_{2}$ and gives us a state encompassing all
the hamming weights from $0$ to $d$ at no additional cost to the circuit
depth. We call this the $H$-_loading_ , whose circuit for $d=8$ is illustrated
in Fig. 14.

$H\ket{0}=\frac{\ket{0}+\ket{1}}{\sqrt{2}}\hskip
28.45274ptH\ket{1}=\frac{\ket{0}-\ket{1}}{\sqrt{2}}$ (8)

####

The $RY$ and $H$ loading circuits spread the data over all $2^{n}$ bases to
allow future RBS-based neural network layers to utilize an exponential space.
In contrast, unary loaders spread the data across the $n$ bases with hamming
weight 1, and since RBS gates are hamming-weight-preserving (known as match
gates [41]), they cannot change this. Their action when using unary loaders is
thus restricted to a much smaller space.

(a) RY Loader

(b) Hadamard Loader

Figure 14: Non-Unary Loaders.

###  5.3 Quantum Orthogonal and Compound Layers

(a) Butterfly

(b) Pyramid

(c) X

Figure 15: Parameterized Quantum Circuits for Orthogonal and Compound Layers.
Vertical lines represent two-qubit RBS gates, parameterized with independent
angles $\theta$, which are shown as the same in the figures for simplicity.

Quantum orthogonal layers consist of a unary data loader plus a parametrized
quantum circuit made of $RBS$ gates, while quantum compound layers consist of
a general data loader plus a parametrized quantum circuit made of $RBS$ gates.

$RBS$ gates and circuits preserve the hamming weight of the input state, and
thus if we use a unary data loader, then the output of the layer will be
another vector in unary amplitude encoding. Similarly, if the loaded quantum
state is a superposition of only basis states of hamming weight $k$, so is the
output state. More generally, we can think of such hamming-weight preserving
circuits with $n$ qubits as block-diagonal unitaries that act separately on
$n+1$ subspaces, where the $k^{th}$ subspace is defined by all computational
basis states with hamming weight equal to $k$. The dimension of these
subspaces is equal to $n\choose k$. The first block of this unitary is an
${n\times n}$ orthogonal matrix, such that when a vector is loaded in the
unary basis, this circuit simply performs orthogonal matrix multiplication. In
general, the $k$-th block of this unitary applies a _compound matrix of order
k_ of the $n\times n$ unary matrix. The dimension of this $k$-th order
compound matrix is ${n\choose k}\times{n\choose k}$. We refer to the layers
that use bases beyond the unary as _compound layers._

There exist many possibilities for building a parametrized quantum circuit
made of $RBS$ gates which can be used in a quantum orthogonal or compound
layer, each with different properties.

The Pyramid circuit (Fig.15), proposed in [17], is a parameterised quantum
circuit composed of exactly $n(n-1)/2$ $RBS$ gates. This circuit requires only
adjacent qubit connectivity, which makes it suitable for most superconducting
qubit hardware. In addition, when restricted to the unary basis, the pyramid
circuit expresses exactly the Special Orthogonal Group, i.e. orthogonal
matrices with the determinant equal to $+1$. To allow this circuit to express
the entire orthogonal group, we can add a final $Z$ gate on the last qubit.
This allows us to express orthogonal matrices with a $-1$ determinant as well.
The pyramid circuit is, therefore, very general and covers all the possible
orthogonal matrices of size $n\times n$.

The $\mathbf{X}$ circuit (Fig.15), introduced in [42], uses just $O(n)$ gates
and has nearest-neighbor connectivity. Due to reduced depth and gate
complexity, it accumulates less hardware noise.

The Butterfly circuit (Fig.15) is inspired by the classical fast Fourier
transform algorithm, and uses $O(nlog(n))$ gates. It was also introduced in
[42], and despite having reduced expressivity compared to the Pyramid circuit,
it often performs just as well.

In [17], a method is proposed to train orthogonal layers for the unary basis
by computing the gradient of each parameter $\theta_{i}$ using
backpropagation. This backpropagation method for the pyramid circuit (which is
the same for any circuit with $RBS$ gates) takes time $O(n^{2})$,
corresponding to the number of gates, and provides a polynomial improvement in
runtime compared to the previously known orthogonal neural network training
algorithms which relied on an $O(n^{3})$ SVD operation [39]. Since the runtime
corresponds to the number of gates, it is lower for the butterfly and $X$
circuits. See Table 3 for full details on the comparison between the three
types of circuits. For the compound layers, we need to consider the entire
$2^{n}\times 2^{n}$ space and thus train an exponential size weight matrix,
which takes exponential time on a classical computer. In principle, a compound
layer can also be trained using the parameter shift rule for quantum circuits,
which can be more efficient since the number of parameters is polynomial in
the input size, though noise in current quantum hardware makes this
impractical for the time being.

Table 3: Comparison of different parameterized quantum circuits for orthogonal and compound layers with $n$ qubits. Circuit |  Hardware Connectivity |  Depth |  # Gates
---|---|---|---
Pyramid |  Nearest-Neighbor |  $2n-3$ |  $\frac{n(n-1)}{2}$
X |  Nearest-Neighbor |  $n-1$ |  $2n-3$
Butterfly |  All-to-all |  $log(n)$ |  $\frac{n}{2}log(n)$
|  |  |

###  5.4 Expectation-per-Subspace Compound Layer

Figure 16: Expectation-per-subspace Compound Layer. In the final step, we
combine the $n\choose 0$ and the $n\choose n$ subspaces and calculate their
overall expectation.

We describe here a compound layer that we call the Expectation-per-subspace
compound layer. This layer involves loading the input vector using a non-unary
basis which could be done either via the $RY$-loading or the $H$-loading
circuit as previously defined. Then, we apply a parameterized quantum circuit
with RBS gates, e.g. a pyramid circuit, which performs the compound matrix
operation on all the fixed hamming weight subspaces. More precisely, we can
think of the operation as performing the matrix-vector multiplication of an
${n\choose k}\times{n\choose k}$ matrix with an ${n\choose k}$-dimensional
vector for each hamming weight $k$ from $0$ to $n$. Note that for $0$ and $n$,
the dimension is $1$, and hence the unitary acts as identity.

If we look at the output quantum state, it defines a distribution over a
domain of size $2^{n}$. Given the exponential size of the distribution, it is
not advisable to try and train the entire distribution, since that would take
exponential time. However, one can still try to use a loss function that
contains some information about the distribution. For example, one can use the
expectation of the distribution, which is what normally happens in variational
quantum algorithms where one approximates this expectation by using a number
of measurement outcomes. Given the fact that our unitary is block-diagonal,
one can try to define a more complex loss function that contains more
information about the distribution. In particular, one can split the domain of
the distribution into $n+1$ subdomains, one for each subspace, and then train
on all these expectations.

This is what we do in the Expectation-per-subspace compound layer, where for
each $k$ from $0$ to $n$, we take the outputs corresponding to the hamming
weight $k$ strings and sort them. Now, for each $k$, we assign values which
are equally spaced between two bounds $a$ and $b$ (which are $0$ and $10$, in
our models) to the $n\choose k$ strings. We normalize the outputs using the
$L1$-norm to correspond to a probability distribution over the $n\choose k$
values between $a$ and $b$, and then we calculate the expectation value for
that hamming weight. This gives us a set of $n+1$ values corresponding to each
hamming weight. Since for hamming weight $0$ and $n$ the dimension of the
subspace is 1 (the all-zero and all-one strings), we combine them and
calculate the expectation for these two together and make the layer have $n$
outputs. The entire operation is illustrated in Fig. 16.

While these compound layers do increase classical simulation complexity, they
do not increase quantum complexity. The advent of better quantum hardware will
allow us to test larger compound layers that explore much larger portions of
the Hilbert space.

##  6 QNN Results with Classical Simulation

In this second study, we focus on the problem of credit-default prediction
associated with credit applications from Small and Medium Enterprises (SMEs).
We report the results of our neural network architectures on this use case and
on public datasets.

###  6.1 Use Case Introduction: Credit Risk Prediction

The credit operation is one of the largest and most important operations in a
retail banking institution. Throughout the credit journey (life-cycle) of a
customer within the bank, several different models are used at different
points of the journey and for different purposes, such as the determination of
interest rates, offering of different products, etc.

The credit granting model is a particularly important one since it determines
whether or not a credit relationship will be established. It is also
particularly challenging in the case of SMEs, where the relationship with the
bank often starts only when the SME submits an application for credit, so very
little data is available.

Given these challenges, we propose the use of quantum techniques aiming at
improving the predictive performance of the credit granting model.

The credit granting decision may be seen as a binary classification problem,
in which the objective is to predict if the SME will default on credit. More
specifically, we are interested in calculating the so-called probability of
default (PD), which is given by $P(\hat{y}=1|\bm{x})$. The PD information is
used internally for other pipelines primarily concerned with the determination
of credit ratings for the SMEs (though in this study, we focus solely on the
PD model), so the PD distribution is the main output of interest from the
model. For this reason, we do not threshold the probability outputs from the
model — thus, we use threshold-independent classification metrics to evaluate
its predictive performance. Namely, the main Key Performance Indicator (KPI)
that we use is the Gini score, constructed from the Area Under the Curve (AUC)
of the ROC (Receiver Operating Characteristic) curve as
$\text{Gini}=2\times\text{AUC}-1$. The Gini score is easily interpreted by the
business team and allows for a holistic estimation of the model’s impact.

In this study, we chose to focus on the development of an ”internal model” of
credit default, which only uses features collected by Itaú, without
considering any external information (from credit bureaus, for instance). The
dataset used consists of $\approx 141,500$ observations, each one represented
by 32 features: 31 numerical and 1 categorical. Each observation represents a
given SME customer in a specified reference month, whose observed target
indicates its default behavior, and whose features consist of internal
information about the company. The data was anonymized, standardized, and
split into training and test sets based on the time period: the training set
consists of $\approx 74,700$ observations covering 12 months of data, while in
the test set, we have $\approx 66,800$ observations covering the subsequent 8
months.

###  6.2 Neural Network Architectures for Credit Risk

To compare the performance of the orthogonal and compound layers to the
classical baseline, we designed three neural network architectures. Each
architecture had three layers: an encoding layer, an experimental layer, and a
classification head. The encoding layer was a standard linear layer of size
$32\times 8$ followed by a $\tanh$ activation. Its purpose is to bring the
dimension of the features down to 8, which is a reasonable simulation size for
both proposed quantum layers. The second layer was the experimental layer of
size $8\times 8$ (described below). Finally, the third layer, the
classification head, was a linear layer of size $8\times 2$ followed by
softmax to predict the probabilities.

The first quantum neural network architecture, named OrthoResNN, uses an
$8\times 8$ orthogonal experimental layer implemented with a semi-diagonal
loader and $X$ circuit. Note that the final output of the layer is provided by
measurements. We add a skip connection by adding the input of the orthogonal
layer to the output 333The model cannot learn to ”ignore” the orthogonal layer
via the skip connection, because the orthogonal layer cannot change the
magnitude of the input.. The layer is followed by a $\tanh$ activation
function. This architecture is illustrated in Fig. 17.

Figure 17: Architecture of the OrthoResNN model.

Our second architecture, ExpResNN, replaces the experimental layer with an
$8\times 8$ Expectation-per-subspace compound layer. We use the the $H$-loader
to encode our data. The layer is again followed by a $\tanh$ activation
function. Fig. 18 illustrates the ExpResNN architecture.

Figure 18: Architecture of the ExpResNN model using the $H$-loader.

And finally, the classical architecture, ResNN, used an $8\times 8$ linear
residual layer followed by $\tanh$.

###  6.3 Methods and Training

The training of the networks was performed using the JAX package by Google. We
train our models for $500$ epochs. To identify suitable hyperparameters, we
performed a search over learning rate, learning-rate-halving points, and batch
size. The hyperparameter search was performed with the ray-tune library.

The dataset contains a large number of missing values, which motivated the
experimentation of different imputation techniques such as zero-filling,
round-robin imputation (implemented in the Python package scikit-learn), and
MICE [43]. The best results were achieved with round-robin imputation using
scikit-learn’s IterativeImputer with Bayesian ridge regression. This was the
pre-processing employed in all the results on the SME dataset.

###  6.4 Results

In our experimental setup, we consider the fully connected residual layer
(ResNN) as the classical benchmark. We performed the same experiment with an
orthogonal layer using the semi-diagonal loader and the X circuit
(OrthoResNN). Finally, we tried the expectation-per-subspace compound layer
with the Hadamard loader and X circuit (ExpResNN). While the performance of
the OrthoResNN and ExpResNN remained nearly the same as the FNN layer, these
new layers learn the angles of $2n$ RBS gates instead of $n^{2}$ elements of
the weight matrix, dramatically reducing the number of parameters needed. The
results are shown in Table 4.

Table 4: Comparison between different architectures. Model |  Trainable parameters |  Test Gini
---|---|---
ResNN |  64 |  54.20 %
OrthoResNN |  13 |  54.29%
ExpResNN |  13 |  53.95%
|  |

The results show that quantum orthogonal and compound layers can preserve the
performance of fully-connected layers on this dataset while using a fraction
of the trainable parameters. We note that the ExpResNN did not show advantages
over the OrthoResFNN.

###  6.5 Further Benchmarks

We compared the orthogonal and linear layers on public classification datasets
from OpenML as we did for the Random Forests (Sec. 3.6). Again, we used a
three-layer architecture: a linear encoding layer to $16$ dimensions with GeLU
activation, a $16\times 16$ experiment layer with tanh activation, and a
binary classification head. These architectures did not have residual
connections.

We compared two models: OrthoFNN and FNN. OrthoFNN used a $16\times 16$
pyramid circuit as the experiment layer, for a total of 120 trainable
parameters. FNN used a feed-forward linear layer with 256 trainable
parameters.

For all the datasets and models, we use a batch size of $128$ and a learning
rate of $10^{-4}$. Each is trained for $500$ epochs and evaluated with the
ROC-AUC metric. The results are summarised in Table 5.

Table 5: Comparison of OrthoFNN and FNN for different datasets. The results are reported via the ROC-AUC metric. Dataset |  FNN |  OrthoFNN
---|---|---
madelon |  0.624 |  0.636
credit-default |  0.647 |  0.634
house-pricing |  0.780 |  0.781
jannis |  0.791 |  0.797
eye movements |  0.580 |  0.578
bank-marketing |  0.847 |  0.850
wine |  0.553 |  0.640
california |  0.634 |  0.668
|  |

##  7 QNN Results on Quantum Hardware

###  7.1 Implementation of Quantum Circuits

Using a classical computer, inference using an orthogonal layer takes time
$O(n^{2})$, while for a general compound layer, this time is exponential in
$n$. Using a quantum computer, inference with an orthogonal or compound layer
uses a quantum circuit that has depth $O(n)$ (Pyramid or X) or $O(\log(n))$
(Butterfly), and $O(n^{2})$ gates. Therefore, one may find a further advantage
if the inference is performed on a quantum computer. This motivated us to test
the inference step for classically trained OrthoResNN and ExpFNN models
(ExpResNN from the classical experiments without a residual connection) on
currently available quantum hardware.

The data loader and orthogonal/compound layer circuits employed in our model
architectures are NISQ-friendly and particularly suitable for superconducting
qubits, with low depth and nearest-neighbors qubit connectivity. Thus, we
chose to use IBM’s 27-qubit machine ibm_hanoi (see Fig.19).

Figure 19: Topology graph of the 27-qubit _ibm_hanoi_ machine used to perform
our hardware experiments. The colors in the qubits indicate readout assignment
error; and in the connections the CNOT error — dark blue is low, purple is
high.

To perform inference on _ibm_hanoi_ , we used the semi-diagonal data loader
and X circuit to implement the OrthoResNN model; and the Hadamard loader and X
circuit for the ExpFNN model — the same architectures described in Sec. 6.
Both neural networks were trained classically, and the trained parameters were
used to construct our quantum circuits for inference.

Given the large size of the test dataset ($66,750$ data points), we decided to
perform inference using the trained models on a small test subsample of $300$
test points, corresponding to the maximum number of circuits we could send in
one job to the IBM machine. After testing different subsamples with the
OrthoResNN model444We made sure to pick a representative subsample of
observations, for which we would have neither a heavily under nor
overestimated Gini score compared to the one for the full test set. By
selecting 50 different subsamples of randomly chosen 300 observations and
performing the classical inference, we found Gini scores between 45.9% and
64.3%, with an average of 54.44%. This further supports the fact that the
subsample that we chose correctly yields the (approximate) expected value for
this Gini score distribution, which then yields an unbiased subsample Gini
score. The same procedure was performed for the ExpFNN experiment., we
selected one for which we achieved a subsample test Gini score of 54.19% using
a noiseless simulator (blue ROC curve in Figure 20). The same was done for the
ExpFNN experiment, yielding a subsample test Gini of 53.90% with the noiseless
simulator. These values were taken as the best possible Gini scores if the
inference was performed on noiseless quantum hardware, which could then be
compared with the values actually achieved with ibm_hanoi.

The circuits were then run on the quantum processor. Due to its limited
Hilbert space of size $n$, the OrthoResNN has a natural error-rejection
procedure: any measurements outside of the unary basis can be disregarded as
errors. As a result, the inference yielded a Gini score of $50.19\%$, as shown
in the orange ROC curve in Fig. 20. The achieved Gini was not too far from the
noise-free simulation result (54.19%), but there was clearly room for
improvement in order to close the 4 pp difference. We also attempted inference
with the more complex ExpFNN, which yielded a Gini score of 40.20%, much
farther from the noiseless simulation Gini of 53.90%. Since the ExpFNN uses
the entire $2^{n}$-dimensional Hilbert space, it is more prone to errors due
to noise, as the error-rejection procedure used for the OrthoResNN cannot be
employed.

###  7.2 Improving the Hardware Results with Error Mitigation Techniques

Error mitigation and error suppression techniques undoubtedly play a very
important role in NISQ-era quantum computing. While these techniques alone may
not be sufficient to fully overcome the imperfections of current quantum
systems, they can push the practical limits of what can be achieved. As a next
step, for the OrthoResNN model, we experimented with various error mitigation
and suppression approaches, going beyond the simple hamming-weight
postselection procedure, in an attempt to close the gap of 4 pp between the
Gini score from the noisy simulation and the one from hardware execution.

The first approach that we tried was a correlated readout mitigator. This is a
purely classical post-processing technique which demands the construction of a
calibration circuit for each one of the possible $2^{N}$ states of the full
$N$ qubits Hilbert space. The calibration circuits’ execution (simulated using
ibm_hanoi’s backend information, in our case) yields a $2^{N}\times 2^{N}$
assignment matrix, which is used to understand how errors might occur during
readout. One can see that this method rapidly becomes intractable as the
number of qubits $N$ increases. In our case, for $N=8$, the Gini score
improved to 50.24%, a small improvement of only 0.05 pp.

Thus, in order to investigate the effect of more robust error suppression and
mitigation techniques in our results, we moved on to a new round of hardware
experiments, performing the inference by executing the exact same OrthoResNN
circuits via the Qiskit Runtime [44] service using the Sampler primitive,
which allows one to use circuit optimization as well as error suppression and
mitigation techniques, as detailed below.

Firstly, we used circuit optimization at the point of circuit transpilation
and compilation by setting the optimization_level parameter to the highest
possible value, 3. This performs the following circuit optimization routines:
Layout selection and routing (VF2 layout pass and SABRE layout search
heuristics [45]; 1 qubit gate optimization (chains of single-qubit u1, u2, u3
gates are combined into a single gate); commutative cancellation (cancelling
of redundant self-adjoint gates); 2 qubit KAK optimization (decomposition of
2-qubit unitaries into the minimal number of uses of 2-qubit basis gates).

Secondly, we used the Dynamical Decoupling error suppression technique [46,
47]. This technique works as a pulse schedule by inserting a DD pulse sequence
into periods of time in which qubits are idle. The DD pulses effectively
behave as an identity gate, thus not altering the logical action of the
circuit, but having the effect of mitigating decoherence in the idle periods,
reducing the impact of errors.

Thirdly, we used the M3 (Matrix-free Measurement Mitigation) error mitigation
technique [48] by setting the Sampler resilience_level parameter to 1 (the
only option available for the Sampler primitive). This provides mitigated
quasi-probability distributions after the measurement. M3 works in a reduced
subspace defined by the noisy input bitstrings supposed to be corrected, which
is often much smaller than the full $N$ qubits Hilbert space. For this reason,
this method is much more efficient than the matrix-based readout mitigator
technique mentioned above. M3 provides a matrix-free preconditioned iterative
solution method, which removes the need to construct the full reduced
assignment matrix but rather computes individual matrix elements, which uses
orders of magnitude less memory than direct factorization.

By employing these three techniques, we were able to achieve a Gini score of
53.68% for the OrthoResNN (Fig. 20). This is a 3.49 pp improvement from the
initial 50.19% Gini of the unmitigated run, falling only 0.53 pp behind the
ideal noiseless execution (54.21% Gini)! This remarkable result underscores
the NISQ-friendliness of the orthogonal layer and highlights the importance of
error suppression and mitigation techniques in the NISQ era.

It is important to note that circuit optimization, error suppression, and
mitigation techniques typically result in some classical/quantum pre/post-
processing overhead to the overall circuit runtime. Some of these techniques
are based on heuristics and/or do not have efficient scaling at larger circuit
sizes. It is important to balance the desired levels of optimization and
resilience with the required time for the full execution, especially as the
circuit sizes increase.

Figure 20: ROC Curves with Gini score for the ideal simulation, hardware
execution, and the error-mitigated hardware execution.

##  8 Conclusion

In this work, we have explored the potential of quantum machine learning
methods in improving forecasting in finance, with a focus on two specific use
cases within the Itaú business: churn prediction and credit risk assessment.
Our results demonstrate that the proposed algorithms, which leverage quantum
ideas, can effectively enhance the performance of Random Forest and neural
network models, achieving better accuracy and training with fewer parameters.

In the present day, quantum hardware is not powerful enough to provide real
improvements or conclusive large-scale benchmarks. Performance enhancements
can be achieved today by turning these quantum ideas into classical ML
solutions run on GPUs. However, with the advent of better quantum hardware, we
expect these methods to run faster and produce even better results when run on
quantum computers.

The general nature of the proposed methods makes them applicable to other use
cases in finance and beyond, although they must be tuned to specific datasets
and tasks. We hope this work inspires confidence that QML research holds
promise both for today as well as for the coming era of scaled, fault-tolerant
quantum hardware.

Disclaimer This paper is a research collaboration between Itaú Unibanco and QC
Ware. Any opinions, findings, conclusions or recommendations expressed in this
material are those of the authors and do not necessarily reflect the views of
Itaú Unibanco. This paper is not and does not constitute or intend to
constitute investment advice or any investment service. It is not and should
not be deemed to be an offer to purchase or sell, or a solicitation of an
offer to purchase or sell, or a recommendation to purchase or sell any
securities or other financial instruments. Moreover, all data used in this
study is compliant with the Brazilian General Data Protection Law.

## Declarations

### Ethical Approval

Not applicable as this study did not involve any human or animal subjects.

### Competing interests

S.T., S.K., N.M., and I.K. are employed by QC Ware, while A.F. and S.B. are
employed by Itaú, both of which may stand to benefit from the results of this
research. The authors declare no other competing interests.

### Authors’ contributions

S.T., S.K. N.M., and A.F. wrote the manuscript body and performed the
classical experiments. N.M., S.K., and A.F. performed the quantum hardware
experiments. All authors reviewed the manuscript.

### Funding

This research did not receive any specific grant from funding agencies in the
public, commercial, or not-for-profit sectors.

### Availability of data and materials

The primary datasets and model used in this study are proprietary to Itau and
QC Ware and are not publicly available. Additional benchmark datasets are
available via openML.

## References

  *   * Herman et al. [2022] Herman, D., Googin, C., Liu, X., Galda, A., Safro, I., Sun, Y., Pistoia, M., Alexeev, Y.: A Survey of Quantum Computing for Finance. Papers 2201.02773, arXiv.org (January 2022). https://ideas.repec.org/p/arx/papers/2201.02773.html
  * Egger et al. [2020] Egger, D.J., Gambella, C., Marecek, J., McFaddin, S., Mevissen, M., Raymond, R., Simonetto, A., Woerner, S., Yndurain, E.: Quantum computing for finance: State-of-the-art and future prospects. IEEE Transactions on Quantum Engineering 1, 1–24 (2020) https://doi.org/10.1109/TQE.2020.3030314
  * McKinsey & Company [2021] McKinsey & Company: Quantum computing: An emerging ecosystem and industry use cases. Accessed: February 16, 2023 (2021)
  * Bouland et al. [2020] Bouland, A., Dam, W., Joorati, H., Kerenidis, I., Prakash, A.: Prospects and challenges of quantum finance. arXiv (2020). https://doi.org/10.48550/ARXIV.2011.06492 . https://arxiv.org/abs/2011.06492
  * Leclerc et al. [2022] Leclerc, L., Ortiz-Guitierrez, L., Grijalva, S., Albrecht, B., Cline, J.R.K., Elfving, V.E., Signoles, A., Henriet, L., Del Bimbo, G., Sheikh, U.A., Shah, M., Andrea, L., Ishtiaq, F., Duarte, A., Mugel, S., Caceres, I., Kurek, M., Orus, R., Seddik, A., Hammammi, O., Isselnane, H., M’tamon, D.: Financial Risk Management on a Neutral Atom Quantum Processor. arXiv (2022). https://doi.org/10.48550/ARXIV.2212.03223 . https://arxiv.org/abs/2212.03223
  * Rebentrost and Lloyd [2018] Rebentrost, P., Lloyd, S.: Quantum computational finance: quantum algorithm for portfolio optimization. arXiv preprint arXiv:1811.03975 (2018)
  * Kerenidis et al. [2019] Kerenidis, I., Prakash, A., Szilágyi, D.: Quantum algorithms for portfolio optimization. In: Proceedings of the 1st ACM Conference on Advances in Financial Technologies. AFT ’19, pp. 147–155. Association for Computing Machinery, New York, NY, USA (2019). https://doi.org/10.1145/3318041.3355465 . https://doi.org/10.1145/3318041.3355465
  * Rebentrost et al. [2022] Rebentrost, P., Luongo, A., Bosch, S., Lloyd, S.: Quantum computational finance: martingale asset pricing for incomplete markets. arXiv (2022). https://doi.org/10.48550/ARXIV.2209.08867 . https://arxiv.org/abs/2209.08867
  * Doriguello et al. [2022] Doriguello, J.a.F., Luongo, A., Bao, J., Rebentrost, P., Santha, M.: Quantum Algorithm for Stochastic Optimal Stopping Problems with Applications in Finance. In: Le Gall, F., Morimae, T. (eds.) 17th Conference on the Theory of Quantum Computation, Communication and Cryptography (TQC 2022). Leibniz International Proceedings in Informatics (LIPIcs), vol. 232, pp. 2–1224. Schloss Dagstuhl – Leibniz-Zentrum für Informatik, Dagstuhl, Germany (2022). https://doi.org/10.4230/LIPIcs.TQC.2022.2 . https://drops.dagstuhl.de/opus/volltexte/2022/16509
  * Suzuki et al. [2020] Suzuki, Y., Uno, S., Raymond, R., Tanaka, T., Onodera, T., Yamamoto, N.: Amplitude estimation without phase estimation. Quantum Information Processing 19(2), 75 (2020) https://doi.org/10.1007/s11128-019-2565-2
  * Giurgica-Tiron et al. [2022] Giurgica-Tiron, T., Kerenidis, I., Labib, F., Prakash, A., Zeng, W.: Low depth algorithms for quantum amplitude estimation. Quantum 6, 745 (2022) https://doi.org/10.22331/q-2022-06-27-745
  * Pistoia et al. [2021] Pistoia, M., Ahmad, S.F., Ajagekar, A., Buts, A., Chakrabarti, S., Herman, D., Hu, S., Jena, A., Minssen, P., Niroula, P., Rattew, A., Sun, Y., Yalovetzky, R.: Quantum Machine Learning for Finance. arXiv (2021). https://doi.org/10.48550/ARXIV.2109.04298 . https://arxiv.org/abs/2109.04298
  * Emmanoulopoulos and Dimoska [2022] Emmanoulopoulos, D., Dimoska, S.: Quantum machine learning in finance: Time series forecasting. arXiv e-prints, 2202 (2022)
  * Alcazar et al. [2020] Alcazar, J., Leyton-Ortega, V., Perdomo-Ortiz, A.: Classical versus quantum models in machine learning: insights from a finance application. Machine Learning: Science and Technology 1(3), 035003 (2020) https://doi.org/10.1088/2632-2153/ab9009
  * Nguyen and Chen [2022] Nguyen, N., Chen, K.-C.: Bayesian quantum neural networks. IEEE Access 10, 54110–54122 (2022) https://doi.org/10.1109/ACCESS.2022.3168675
  * Kerenidis and Prakash [2022] Kerenidis, I., Prakash, A.: Quantum machine learning with subspace states. arXiv preprint arXiv:2202.00054 (2022)
  * Landman et al. [2022] Landman, J., Mathur, N., Li, Y.Y., Strahm, M., Kazdaghli, S., Prakash, A., Kerenidis, I.: Quantum Methods for Neural Networks and Application to Medical Image Classification. Quantum 6, 881 (2022) https://doi.org/10.22331/q-2022-12-22-881
  * Arjovsky et al. [2016] Arjovsky, M., Shah, A., Bengio, Y.: Unitary evolution recurrent neural networks. In: International Conference on Machine Learning, pp. 1120–1128 (2016). PMLR
  * Breiman [2001] Breiman, L.: Random forests. Machine Learning 45(1), 5–32 (2001) https://doi.org/10.1023/A:1010933404324
  * Kulesza and Taskar [2012] Kulesza, A., Taskar, B.: Determinantal point processes for machine learning. Foundations and Trends® in Machine Learning 5(2–3), 123–286 (2012) https://doi.org/10.1561/2200000044
  * Macchi [1975] Macchi, O.: The coincidence approach to stochastic point processes. Advances in Applied Probability 7(1), 83–122 (1975) https://doi.org/10.2307/1425855
  * Derezinski and Mahoney [2021] Derezinski, M., Mahoney, M.W.: Determinantal point processes in randomized numerical linear algebra. Notices of the American Mathematical Society 68(1), 34–45 (2021)
  * Bardenet and Hardy [2020] Bardenet, R., Hardy, A.: Monte carlo with determinantal point processes. The Annals of Applied Probability 30(1), 368–417 (2020)
  * Elfeki et al. [2019] Elfeki, M., Couprie, C., Riviere, M., Elhoseiny, M.: Gdpp: Learning diverse generations using determinantal point processes. In: International Conference on Machine Learning, pp. 1774–1783 (2019). PMLR
  * Derezinski [2018] Derezinski, M.: Volume sampling for linear regression. UC Santa Cruz Electronic Theses and Dissertations 36 (2018)
  * Dereziński et al. [2018] Dereziński, M., Warmuth, M.K., Hsu, D.: Leveraged volume sampling for linear regression. In: Proceedings of the 32nd International Conference on Neural Information Processing Systems. NIPS’18, pp. 2510–2519. Curran Associates Inc., Red Hook, NY, USA (2018)
  * Kulesza and Taskar [2011] Kulesza, A., Taskar, B.: K-dpps: Fixed-size determinantal point processes. In: Proceedings of the 28th International Conference on International Conference on Machine Learning. ICML’11, pp. 1193–1200. Omnipress, Madison, WI, USA (2011)
  * Hough et al. [2006] Hough, J.B., Krishnapur, M., Peres, Y., Virág, B.: Determinantal Processes and Independence. Probability Surveys 3(none), 206–229 (2006) https://doi.org/10.1214/154957806000000078
  * Anari et al. [2016] Anari, N., Oveis Gharan, S., Rezaei, A.: Monte carlo markov chain algorithms for sampling strongly rayleigh distributions and determinantal point processes. In: Feldman, V., Rakhlin, A., Shamir, O. (eds.) 29th Annual Conference on Learning Theory. Proceedings of Machine Learning Research, vol. 49, pp. 103–115. PMLR, Columbia University, New York, New York, USA (2016). https://proceedings.mlr.press/v49/anari16.html
  * Li et al. [2016] Li, C., Sra, S., Jegelka, S.: Fast mixing markov chains for strongly rayleigh measures, dpps, and constrained sampling. In: Lee, D.D., Sugiyama, M., Luxburg, U., Guyon, I., Garnett, R. (eds.) Advances in Neural Information Processing Systems 29: Annual Conference on Neural Information Processing Systems 2016, December 5-10, 2016, Barcelona, Spain, pp. 4188–4196 (2016). https://proceedings.neurips.cc/paper/2016/hash/850af92f8d9903e7a4e0559a98ecc857-Abstract.html
  * Derezinski et al. [2019] Derezinski, M., Calandriello, D., Valko, M.: Exact sampling of determinantal point processes with sublinear time preprocessing. In: Wallach, H., Larochelle, H., Beygelzimer, A., Alché-Buc, F., Fox, E., Garnett, R. (eds.) Advances in Neural Information Processing Systems, vol. 32. Curran Associates, Inc., Vancouver, Canada (2019). https://proceedings.neurips.cc/paper/2019/file/fa3060edb66e6ff4507886f9912e1ab9-Paper.pdf
  * Calandriello et al. [2020] Calandriello, D., Derezinski, M., Valko, M.: Sampling from a k-dpp without looking at all items. In: Larochelle, H., Ranzato, M., Hadsell, R., Balcan, M.F., Lin, H. (eds.) Advances in Neural Information Processing Systems, vol. 33, pp. 6889–6899. Curran Associates, Inc., Online (2020). https://proceedings.neurips.cc/paper/2020/file/4d410063822cd9be28f86701c0bc3a31-Paper.pdf
  * Gautier et al. [2019] Gautier, G., Polito, G., Bardenet, R., Valko, M.: Dppy: Dpp sampling with python. Journal of Machine Learning Research 20(180), 1–7 (2019)
  * Kazdaghli et al. [2023] Kazdaghli, S., Kerenidis, I., Kieckbusch, J., Teare, P.: Improved clinical data imputation via classical and quantum determinantal point processes. arXiv preprint arXiv:2303.17893 (2023)
  * Grinsztajn et al. [2022] Grinsztajn, L., Oyallon, E., Varoquaux, G.: Why do tree-based models still outperform deep learning on typical tabular data? In: Thirty-sixth Conference on Neural Information Processing Systems Datasets and Benchmarks Track (2022). https://openreview.net/forum?id=Fp7__phQszn
  * Benedetti et al. [2019] Benedetti, M., Lloyd, E., Sack, S., Fiorentini, M.: Parameterized quantum circuits as machine learning models. Quantum Science and Technology 4(4), 043001 (2019) https://doi.org/10.1088/2058-9565/ab4eb5 arXiv:1906.07682 [quant-ph]
  * Havlíček et al. [2018] Havlíček, V., Córcoles, A.D., Temme, K., Harrow, A.W., Kandala, A., Chow, J.M., Gambetta, J.M.: Supervised learning with quantum-enhanced feature spaces. Nature 567, 209–212 (2018)
  * Liu et al. [2021] Liu, Y., Arunachalam, S., Temme, K.: A rigorous and robust quantum speed-up in supervised machine learning. Nature Physics 17(9), 1013–1017 (2021) https://doi.org/10.1038/s41567-021-01287-z arXiv:2010.02174 [quant-ph]
  * Jia et al. [2019] Jia, K., Li, S., Wen, Y., Liu, T., Tao, D.: Orthogonal deep neural networks. IEEE transactions on pattern analysis and machine intelligence (2019)
  * Johri et al. [2021] Johri, S., Debnath, S., Mocherla, A., Singh, A., Prakash, A., Kim, J., Kerenidis, I.: Nearest centroid classification on a trapped ion quantum computer. npj Quantum Information (to appear), arXiv:2012.04145 (2021)
  * Jozsa and Miyake [2008] Jozsa, R., Miyake, A.: Matchgates and classical simulation of quantum circuits. Proceedings of the Royal Society A: Mathematical, Physical and Engineering Sciences 464(2100), 3089–3106 (2008) https://doi.org/10.1098/rspa.2008.0189
  * Cherrat et al. [2022] Cherrat, E.A., Kerenidis, I., Mathur, N., Landman, J., Strahm, M., Li, Y.Y.: Quantum vision transformers. arXiv preprint arXiv:2209.08167 (2022)
  * van Buuren and Groothuis-Oudshoorn [2011] Buuren, S., Groothuis-Oudshoorn, K.: mice: Multivariate imputation by chained equations in r. Journal of Statistical Software 45(3), 1–67 (2011) https://doi.org/10.18637/jss.v045.i03
  * contributors [2021] contributors, Q.: Qiskit: An open-source framework for quantum computing (2021) https://doi.org/10.5281/zenodo.2573505
  * Li et al. [2018] Li, G., Ding, Y., Xie, Y.: Tackling the Qubit Mapping Problem for NISQ-Era Quantum Devices. arXiv e-prints, 1809–02573 (2018) https://doi.org/10.48550/arXiv.1809.02573 arXiv:1809.02573 [cs.ET]
  * Viola and Lloyd [1998] Viola, L., Lloyd, S.: Dynamical suppression of decoherence in two state quantum systems. Phys. Rev. A 58, 2733 (1998) https://doi.org/10.1103/PhysRevA.58.2733 arXiv:quant-ph/9803057
  * Ezzell et al. [2022] Ezzell, N., Pokharel, B., Tewala, L., Quiroz, G., Lidar, D.A.: Dynamical decoupling for superconducting qubits: a performance survey (2022) arXiv:2207.03670 [quant-ph]
  * Nation et al. [2021] Nation, P.D., Kang, H., Sundaresan, N., Gambetta, J.M.: Scalable mitigation of measurement errors on quantum computers. PRX Quantum 2, 040326 (2021) https://doi.org/%****␣Main.bbl␣Line␣850␣****10.1103/PRXQuantum.2.040326 arXiv:2108.12518 [quant-ph]


* Chael et al. (2018b) —. 2018b, MNRAS, 478, 5209
  * Chael et al. (2019) Chael, A., Bouman, K., Johnson, M., et al. 2019, eht-imaging: v1.1.0: Imaging interferometric data with regularized maximum likelihood, Zenodo
  * Chan et al. (2009) Chan, C.-k., Liu, S., Fryer, C. L., et al. 2009, ApJ, 701, 521
  * Chan et al. (2015a) Chan, C.-k., Psaltis, D., Özel, F., et al. 2015a, ApJ, 812, 103
  * Chan et al. (2015b) Chan, C.-K., Psaltis, D., Özel, F., Narayan, R., & Skadowski, A. 2015b, ApJ, 799, 1
  * Chandra et al. (2017) Chandra, M., Foucart, F., & Gammie, C. F. 2017, ApJ, 837, 92
  * Chandra et al. (2015) Chandra, M., Gammie, C. F., Foucart, F., & Quataert, E. 2015, ApJ, 810, 162
  * Chashkina et al. (2021) Chashkina, A., Bromberg, O., & Levinson, A. 2021, MNRAS, 508, 1241
  * Chatterjee et al. (2020) Chatterjee, K., Younsi, Z., Liska, M., et al. 2020, MNRAS, 499, 362
  * Chatterjee et al. (2021) Chatterjee, K., Markoff, S., Neilsen, J., et al. 2021, MNRAS, 507, 5281
  * Chen & Yuan (2020) Chen, A. Y., & Yuan, Y. 2020, ApJ, 895, 121
  * Chen et al. (2018) Chen, A. Y., Yuan, Y., & Yang, H. 2018, ApJ, 863, L31
  * Chernoglazov et al. (2021) Chernoglazov, A., Ripperda, B., & Philippov, A. 2021, ApJ, 923, L13
  * Comisso & Sironi (2018) Comisso, L., & Sironi, L. 2018, Phys. Rev. Lett., 121, 255101
  * Crinquand et al. (2021) Crinquand, B., Cerutti, B., Dubus, G., Parfrey, K., & Philippov, A. 2021, A&A, 650, A163
  * Crinquand et al. (2020) Crinquand, B., Cerutti, B., Philippov, A., Parfrey, K., & Dubus, G. 2020, Phys. Rev. Lett., 124, 145101
  * Crumley et al. (2019) Crumley, P., Caprioli, D., Markoff, S., & Spitkovsky, A. 2019, MNRAS, 485, 5105
  * Cruz-Osorio et al. (2021) Cruz-Osorio, A., Fromm, C. M., Mizuno, Y., et al. 2021, Nature Astronomy, arXiv:2111.02517 [astro-ph.HE]
  * Davelaar et al. (2018) Davelaar, J., Mościbrodzka, M., Bronzwaer, T., & Falcke, H. 2018, A&A, 612, A34
  * De Villiers & Hawley (2003) De Villiers, J.-P., & Hawley, J. F. 2003, ApJ, 589, 458
  * De Villiers et al. (2003) De Villiers, J.-P., Hawley, J. F., & Krolik, J. H. 2003, ApJ, 599, 1238
  * Del Zanna et al. (2007) Del Zanna, L., Zanotti, O., Bucciantini, N., & Londrillo, P. 2007, A&A, 473, 11
  * Dexter et al. (2009) Dexter, J., Agol, E., & Fragile, P. C. 2009, ApJ, 703, L142
  * Dexter et al. (2010) Dexter, J., Agol, E., Fragile, P. C., & McKinney, J. C. 2010, ApJ, 717, 1092
  * Dexter & Fragile (2013) Dexter, J., & Fragile, P. C. 2013, MNRAS, 432, 2252
  * Dexter et al. (2020) Dexter, J., Jiménez-Rosales, A., Ressler, S. M., et al. 2020, MNRAS, 494, 4168
  * Do et al. (2019) Do, T., Hees, A., Ghez, A., et al. 2019, Science, 365, 664
  * Doeleman et al. (2008) Doeleman, S. S., Weintroub, J., Rogers, A. E. E., et al. 2008, Nature, 455, 78
  * Dolence et al. (2009) Dolence, J. C., Gammie, C. F., Mościbrodzka, M., & Leung, P. K. 2009, ApJS, 184, 387
  * El Mellah et al. (2021) El Mellah, I., Cerutti, B., Crinquand, B., & Parfrey, K. 2021, arXiv e-prints, arXiv:2112.03933
  * Emami et al. (2021) Emami, R., Anantua, R., Chael, A. A., & Loeb, A. 2021, ApJ, 923, 272
  * Event Horizon Telescope Collaboration et al. (2019a) Event Horizon Telescope Collaboration, Akiyama, K., Alberdi, A., et al. 2019a, ApJ, 875, L1, (M87Paper I)
  * Event Horizon Telescope Collaboration et al. (2019b) —. 2019b, ApJ, 875, L2, (M87Paper II)
  * Event Horizon Telescope Collaboration et al. (2019c) —. 2019c, ApJ, 875, L3, (M87Paper III)
  * Event Horizon Telescope Collaboration et al. (2019d) —. 2019d, ApJ, 875, L4, (M87Paper IV)
  * Event Horizon Telescope Collaboration et al. (2019e) —. 2019e, ApJ, 875, L5, (M87Paper V)
  * Event Horizon Telescope Collaboration et al. (2019f) —. 2019f, ApJ, 875, L6, (M87Paper VI)
  * Event Horizon Telescope Collaboration et al. (2021a) Event Horizon Telescope Collaboration, Akiyama, K., Algaba, J. C., et al. 2021a, ApJ, 910, L12
  * Event Horizon Telescope Collaboration et al. (2021b) —. 2021b, ApJ, 910, L13
  * Event Horizon Telescope Collaboration et al. (2022a) Event Horizon Telescope Collaboration, Akiyama, K., Alberdi, A., et al. 2022a, ApJ, 2
  * Event Horizon Telescope Collaboration et al. (2022b) —. 2022b, ApJ, 3
  * Event Horizon Telescope Collaboration et al. (2022c) —. 2022c, ApJ, 4
  * Event Horizon Telescope Collaboration et al. (2022d) —. 2022d, ApJ, 6
  * Event Horizon Telescope Collaboration et al. (2022e) —. 2022e, ApJ, 5
  * Event Horizon Telescope Collaboration et al. (2022f) —. 2022f, ApJ, 1
  * Falcke (1996) Falcke, H. 1996, in Unsolved Problems of the Milky Way, ed. L. Blitz & P. J. Teuben, Vol. 169, 169
  * Falcke et al. (1993) Falcke, H., Mannheim, K., & Biermann, P. L. 1993, A&A, 278, L1
  * Falcke & Markoff (2000) Falcke, H., & Markoff, S. 2000, A&A, 362, 113
  * Falcke et al. (2000) Falcke, H., Melia, F., & Agol, E. 2000, ApJ, 528, L13
  * Farah et al. (2022) Farah, J. R., et al. 2022, in prep., 00, 0
  * Fishbone & Moncrief (1976) Fishbone, L. G., & Moncrief, V. 1976, ApJ, 207, 962
  * Foucart et al. (2016) Foucart, F., Chandra, M., Gammie, C. F., & Quataert, E. 2016, MNRAS, 456, 1332
  * Foucart et al. (2017) Foucart, F., Chandra, M., Gammie, C. F., Quataert, E., & Tchekhovskoy, A. 2017, MNRAS, 470, 2240
  * Fragile et al. (2007) Fragile, P. C., Blaes, O. M., Anninos, P., & Salmonson, J. D. 2007, ApJ, 668, 417
  * Fragile et al. (2012) Fragile, P. C., Gillespie, A., Monahan, T., Rodriguez, M., & Anninos, P. 2012, ApJS, 201, 9
  * Fragile et al. (2014) Fragile, P. C., Olejar, A., & Anninos, P. 2014, ApJ, 796, 22
  * Fromm et al. (2022) Fromm, C. M., Cruz-Osorio, A., Mizuno, Y., et al. 2022, A&A, 660, A107
  * Gammie et al. (2003) Gammie, C. F., McKinney, J. C., & Tóth, G. 2003, ApJ, 589, 444
  * Generozov et al. (2014) Generozov, A., Blaes, O., Fragile, P. C., & Henisey, K. B. 2014, ApJ, 780, 81
  * Georgiev et al. (2022) Georgiev, B., et al. 2022, in prep., 00, 0
  * Giacomazzo & Rezzolla (2007) Giacomazzo, B., & Rezzolla, L. 2007, Classical and Quantum Gravity, 24, S235
  * Goddi et al. (2021) Goddi, C., Martí-Vidal, I., Messias, H., et al. 2021, ApJ, 910, L14
  * Gold et al. (2017) Gold, R., McKinney, J. C., Johnson, M. D., & Doeleman, S. S. 2017, ApJ, 837, 180
  * Gold et al. (2020) Gold, R., Broderick, A. E., Younsi, Z., et al. 2020, ApJ, 897, 148
  * Gravity Collaboration et al. (2018) Gravity Collaboration, Abuter, R., Amorim, A., et al. 2018, A&A, 618, L10
  * Gravity Collaboration et al. (2019) —. 2019, A&A, 625, L10
  * Gravity Collaboration et al. (2020a) Gravity Collaboration, Jiménez-Rosales, A., Dexter, J., et al. 2020a, A&A, 643, A56
  * Gravity Collaboration et al. (2020b) Gravity Collaboration, Abuter, R., Amorim, A., et al. 2020b, A&A, 638, A2
  * Hawley (2000) Hawley, J. F. 2000, ApJ, 528, 462
  * Hornstein et al. (2007) Hornstein, S. D., Matthews, K., Ghez, A. M., et al. 2007, ApJ, 667, 900
  * Hoshino (2013) Hoshino, M. 2013, ApJ, 773, 118
  * Hoshino (2015) —. 2015, Phys. Rev. Lett., 114, 061101
  * Howes (2010) Howes, G. G. 2010, MNRAS, 409, L104
  * Huang et al. (2009) Huang, L., Takahashi, R., & Shen, Z.-Q. 2009, ApJ, 706, 960
  * Hunter (2007) Hunter, J. D. 2007, Computing in Science and Engineering, 9, 90
  * Ichimaru (1977) Ichimaru, S. 1977, ApJ, 214, 840
  * Igumenshchev et al. (2003) Igumenshchev, I. V., Narayan, R., & Abramowicz, M. A. 2003, ApJ, 592, 1042
  * Inchingolo et al. (2018) Inchingolo, G., Grismayer, T., Loureiro, N. F., Fonseca, R. A., & Silva, L. O. 2018, ApJ, 859, 149
  * Issaoun et al. (2019) Issaoun, S., Johnson, M. D., Blackburn, L., et al. 2019, ApJ, 871, 30
  * Issaoun et al. (2021) —. 2021, ApJ, 915, 99
  * Johnson et al. (2015) Johnson, M. D., Fish, V. L., Doeleman, S. S., et al. 2015, Science, 350, 1242
  * Johnson et al. (2018a) Johnson, M. D., Narayan, R., Psaltis, D., et al. 2018a, ApJ, 865, 104
  * Johnson et al. (2018b) —. 2018b, ApJ, 865, 104, (J18)
  * Kerr (1963) Kerr, R. P. 1963, Phys. Rev. Lett., 11, 237
  * Kisaka et al. (2020) Kisaka, S., Levinson, A., & Toma, K. 2020, ApJ, 902, 80
  * Koide et al. (1999) Koide, S., Shibata, K., & Kudoh, T. 1999, ApJ, 522, 727
  * Komissarov (2001) Komissarov, S. S. 2001, MNRAS, 326, L41
  * Kunz et al. (2015) Kunz, M. W., Schekochihin, A. A., Chen, C. H. K., Abel, I. G., & Cowley, S. C. 2015, Journal of Plasma Physics, 81, 325810501
  * Kunz et al. (2016) Kunz, M. W., Stone, J. M., & Quataert, E. 2016, Phys. Rev. Lett., 117, 235101
  * Lee & Gammie (2021) Lee, D., & Gammie, C. F. 2021, ApJ, 906, 39
  * Lee et al. (2022) Lee, D., et al. 2022, in prep.
  * Leung et al. (2011) Leung, P. K., Gammie, C. F., & Noble, S. C. 2011, ApJ, 737, 21
  * Levinson & Cerutti (2018) Levinson, A., & Cerutti, B. 2018, A&A, 616, A184
  * Li et al. (2013) Li, Z., Morris, M. R., & Baganoff, F. K. 2013, ApJ, 779, 154
  * Liska et al. (2018) Liska, M., Hesp, C., Tchekhovskoy, A., et al. 2018, MNRAS, 474, L81
  * Liska et al. (2019) Liska, M., Chatterjee, K., Tchekhovskoy, A., et al. 2019, arXiv e-prints, arXiv:1912.10192
  * Macquart et al. (2006) Macquart, J.-P., Bower, G. C., Wright, M. C. H., Backer, D. C., & Falcke, H. 2006, ApJ, 646, L111
  * Marrone et al. (2006a) Marrone, D. P., Moran, J. M., Zhao, J.-H., & Rao, R. 2006a, ApJ, 640, 308
  * Marrone et al. (2006b) Marrone, D. P., Moran, J. M., Zhao, J.-H., & Rao, R. 2006b, in Journal of Physics Conference Series, Vol. 54, Journal of Physics Conference Series, 354
  * Marshall et al. (2018) Marshall, M. D., Avara, M. J., & McKinney, J. C. 2018, MNRAS, 478, 1837
  * Marszewski et al. (2021) Marszewski, A., Prather, B. S., Joshi, A. V., Pandya, A., & Gammie, C. F. 2021, ApJ, 921, 17
  * Medeiros et al. (2018) Medeiros, L., Chan, C.-k., Özel, F., et al. 2018, ApJ, 856, 163
  * Mizuno et al. (2021) Mizuno, Y., Fromm, C. M., Younsi, Z., et al. 2021, MNRAS, 506, 741
  * Mościbrodzka & Falcke (2013) Mościbrodzka, M., & Falcke, H. 2013, A&A, 559, L3
  * Mościbrodzka et al. (2016) Mościbrodzka, M., Falcke, H., & Shiokawa, H. 2016, A&A, 586, A38
  * Mościbrodzka et al. (2014) Mościbrodzka, M., Falcke, H., Shiokawa, H., & Gammie, C. F. 2014, A&A, 570, A7
  * Mościbrodzka & Gammie (2018) Mościbrodzka, M., & Gammie, C. F. 2018, MNRAS, 475, 43
  * Mościbrodzka et al. (2011) Mościbrodzka, M., Gammie, C. F., Dolence, J. C., & Shiokawa, H. 2011, ApJ, 735, 9
  * Mościbrodzka et al. (2009) Mościbrodzka, M., Gammie, C. F., Dolence, J. C., Shiokawa, H., & Leung, P. K. 2009, ApJ, 706, 497
  * Mościbrodzka et al. (2021) Mościbrodzka, M., Janiuk, A., & De Laurentis, M. 2021, MNRAS, 508, 4282
  * Most & Noronha (2021) Most, E. R., & Noronha, J. 2021, Phys. Rev. D, 104, 103028
  * Most et al. (2021) Most, E. R., Noronha, J., & Philippov, A. A. 2021, arXiv e-prints, arXiv:2111.05752
  * Narayan et al. (2022) Narayan, R., Chael, A., Chatterjee, K., Ricarte, A., & Curd, B. 2022, MNRAS, 511, 3795
  * Narayan et al. (2003) Narayan, R., Igumenshchev, I. V., & Abramowicz, M. A. 2003, PASJ, 55, L69
  * Narayan et al. (1998) Narayan, R., Mahadevan, R., Grindlay, J. E., Popham, R. G., & Gammie, C. 1998, ApJ, 492, 554
  * Narayan et al. (2012) Narayan, R., SÄ dowski, A., Penna, R. F., & Kulkarni, A. K. 2012, MNRAS, 426, 3241
  * Narayan & Yi (1994) Narayan, R., & Yi, I. 1994, ApJ, 428, L13
  * Narayan & Yi (1995a) —. 1995a, ApJ, 444, 231
  * Narayan & Yi (1995b) —. 1995b, ApJ, 452, 710
  * Narayan et al. (1995) Narayan, R., Yi, I., & Mahadevan, R. 1995, Nature, 374, 623
  * Narayan et al. (1996) —. 1996, A&AS, 120, 287
  * Nathanail et al. (2020) Nathanail, A., Fromm, C. M., Porth, O., et al. 2020, MNRAS, 495, 1549
  * Nathanail et al. (2021) Nathanail, A., Mpisketzis, V., Porth, O., Fromm, C. M., & Rezzolla, L. 2021, arXiv e-prints, arXiv:2111.03689
  * Nättilä & Beloborodov (2021) Nättilä, J., & Beloborodov, A. M. 2021, ApJ, 921, 87
  * Neilsen et al. (2013) Neilsen, J., Nowak, M. A., Gammie, C., et al. 2013, ApJ, 774, 42
  * Oliphant (2007) Oliphant, T. E. 2007, Computing in Science & Engineering, 9, 10
  * Olivares et al. (2019) Olivares, H., Porth, O., Davelaar, J., et al. 2019, A&A, 629, A61
  * Olivares et al. (2022) Olivares, H., et al. 2022, in prep., 00, 0
  * Olivares Sánchez et al. (2018) Olivares Sánchez, H., Porth, O., & Mizuno, Y. 2018, in Journal of Physics Conference Series, Vol. 1031, Journal of Physics Conference Series, 012008
  * Özel et al. (2000) Özel, F., Psaltis, D., & Narayan, R. 2000, ApJ, 541, 234
  * Pandya et al. (2016) Pandya, A., Zhang, Z., Chandra, M., & Gammie, C. F. 2016, ApJ, 822, 34
  * Parfrey et al. (2019) Parfrey, K., Philippov, A., & Cerutti, B. 2019, Phys. Rev. Lett., 122, 035101
  * Porth et al. (2017) Porth, O., Olivares, H., Mizuno, Y., et al. 2017, Computational Astrophysics and Cosmology, 4, 1
  * Porth et al. (2019) Porth, O., Chatterjee, K., Narayan, R., et al. 2019, ApJS, 243, 26
  * Prather et al. (2021) Prather, B., Wong, G., Dhruv, V., et al. 2021, The Journal of Open Source Software, 6, 3336
  * Prather et al. (2022) Prather, B., et al. 2022, in prep.
  * Psaltis et al. (2018) Psaltis, D., Johnson, M., Narayan, R., et al. 2018, arXiv e-prints, arXiv:1805.01242
  * Pu & Broderick (2018) Pu, H.-Y., & Broderick, A. E. 2018, ApJ, 863, 148
  * Qian et al. (2018) Qian, Q., Fendt, C., & Vourellis, C. 2018, ApJ, 859, 28
  * Quataert (1998) Quataert, E. 1998, The Astrophysical Journal, 500, 978
  * Quataert & Gruzinov (2000) Quataert, E., & Gruzinov, A. 2000, ApJ, 545, 842
  * Rees et al. (1982) Rees, M. J., Begelman, M. C., Blandford, R. D., & Phinney, E. S. 1982, Nature, 295, 17
  * Reid et al. (2019) Reid, M. J., Menten, K. M., Brunthaler, A., et al. 2019, ApJ, 885, 131
  * Ressler et al. (2018) Ressler, S. M., Quataert, E., & Stone, J. M. 2018, MNRAS, 478, 3544
  * Ressler et al. (2019) —. 2019, MNRAS, 482, L123
  * Ressler et al. (2020a) —. 2020a, MNRAS, 492, 3272
  * Ressler et al. (2015) Ressler, S. M., Tchekhovskoy, A., Quataert, E., Chandra, M., & Gammie, C. F. 2015, MNRAS, 454, 1848
  * Ressler et al. (2017) Ressler, S. M., Tchekhovskoy, A., Quataert, E., & Gammie, C. F. 2017, MNRAS, 467, 3604
  * Ressler et al. (2020b) Ressler, S. M., White, C. J., Quataert, E., & Stone, J. M. 2020b, ApJ, 896, L6
  * Ripperda et al. (2020) Ripperda, B., Bacchini, F., & Philippov, A. A. 2020, ApJ, 900, 100
  * Ripperda et al. (2022) Ripperda, B., Liska, M., Chatterjee, K., et al. 2022, ApJ, 924, L32
  * Ripperda et al. (2019) Ripperda, B., Bacchini, F., Porth, O., et al. 2019, ApJS, 244, 10
  * Riquelme et al. (2012) Riquelme, M. A., Quataert, E., Sharma, P., & Spitkovsky, A. 2012, ApJ, 755, 50
  * Rowan et al. (2017) Rowan, M. E., Sironi, L., & Narayan, R. 2017, The Astrophysical Journal, 850, 29
  * Ryan et al. (2018) Ryan, B. R., Ressler, S. M., Dolence, J. C., Gammie, C., & Quataert, E. 2018, The Astrophysical Journal, 864, 126
  * Scepi et al. (2021) Scepi, N., Dexter, J., & Begelman, M. C. 2021, arXiv e-prints, arXiv:2107.08056
  * Shakura & Sunyaev (1973) Shakura, N. I., & Sunyaev, R. A. 1973, A&A, 24, 337
  * Shapiro et al. (1976) Shapiro, S. L., Lightman, A. P., & Eardley, D. M. 1976, ApJ, 204, 187
  * Shepherd (1997) Shepherd, M. C. 1997, in Astronomical Society of the Pacific Conference Series, Vol. 125, Astronomical Data Analysis Software and Systems VI, ed. G. Hunt & H. Payne, 77
  * Sironi et al. (2015) Sironi, L., Petropoulou, M., & Giannios, D. 2015, MNRAS, 450, 183
  * Sironi et al. (2021) Sironi, L., Rowan, M. E., & Narayan, R. 2021, ApJ, 907, L44
  * Sądowski et al. (2013) Sądowski, A., Narayan, R., Tchekhovskoy, A., & Zhu, Y. 2013, MNRAS, 429, 3533
  * Stone et al. (2016) Stone, J. M., Marrone, D. P., Dowell, C. D., et al. 2016, ApJ, 825, 32
  * Sądowski et al. (2016) Sądowski, A., Wielgus, M., Narayan, R., et al. 2016, Monthly Notices of the Royal Astronomical Society, 466, 705
  * Tchekhovskoy et al. (2011) Tchekhovskoy, A., Narayan, R., & McKinney, J. C. 2011, MNRAS, 418, L79
  * Tiede (2022) Tiede, P. 2022, in prep.
  * van der Walt et al. (2011) van der Walt, S., Colbert, S. C., & Varoquaux, G. 2011, Computing in Science and Engineering, 13, 22
  * von Fellenberg et al. (2018) von Fellenberg, S. D., Gillessen, S., Graciá-Carpio, J., et al. 2018, ApJ, 862, 129
  * Vourellis et al. (2019) Vourellis, C., Fendt, C., Qian, Q., & Noble, S. C. 2019, ApJ, 882, 2
  * Wang et al. (2013) Wang, Q. D., Nowak, M. A., Markoff, S. B., et al. 2013, Science, 341, 981
  * Werner et al. (2017) Werner, G. R., Uzdensky, D. A., Begelman, M. C., Cerutti, B., & Nalewajko, K. 2017, Monthly Notices of the Royal Astronomical Society, 473, 4840
  * White & Quataert (2022) White, C. J., & Quataert, E. 2022, ApJ, 926, 136
  * White et al. (2019) White, C. J., Quataert, E., & Blaes, O. 2019, ApJ, 878, 51
  * White et al. (2016) White, C. J., Stone, J. M., & Gammie, C. F. 2016, ApJS, 225, 22
  * Wielgus et al. (2022) Wielgus, M., Marchili, N., Martí-Vidal, I., et al. 2022, in prep., 00, 0
  * Witzel et al. (2018) Witzel, G., Martinez, G., Hora, J., et al. 2018, ApJ, 863, 15
  * Wong & Gammie (2022) Wong, G. N., & Gammie, C. F. 2022, in prep.
  * Wong et al. (2021) Wong, G. N., Ryan, B. R., & Gammie, C. F. 2021, ApJ, 907, 73
  * Wong et al. (2022) Wong, G. N., Prather, B. S., Dhruv, V., et al. 2022, ApJS, 259, 64
  * Yarza et al. (2020) Yarza, R., Wong, G. N., Ryan, B. R., & Gammie, C. F. 2020, ApJ, 898, 50
  * Yoon et al. (2020) Yoon, D., Chatterjee, K., Markoff, S. B., et al. 2020, MNRAS, 499, 3178
  * Younsi et al. (2020) Younsi, Z., Porth, O., Mizuno, Y., Fromm, C. M., & Olivares, H. 2020, in Perseus in Sicily: From Black Hole to Cluster Outskirts, ed. K. Asada, E. de Gouveia Dal Pino, M. Giroletti, H. Nagai, & R. Nemmen, Vol. 342, 9
  * Younsi et al. (2012a) Younsi, Z., Wu, K., & Fuerst, S. V. 2012a, A&A, 545, A13
  * Younsi et al. (2012b) —. 2012b, A&A, 545, A13
  * Yuan et al. (2002) Yuan, F., Markoff, S., & Falcke, H. 2002, A&A, 383, 854
  * Yuan & Narayan (2014) Yuan, F., & Narayan, R. 2014, ARA&A, 52, 529
  * Yuan et al. (2003) Yuan, F., Quataert, E., & Narayan, R. 2003, The Astrophysical Journal, 598, 301
  * Yuan et al. (2018) Yuan, Q., Wang, Q. D., Liu, S., & Wu, K. 2018, MNRAS, 473, 306
  * Zhdankin et al. (2019) Zhdankin, V., Uzdensky, D. A., Werner, G. R., & Begelman, M. C. 2019, Phys. Rev. Lett., 122, 055101
  * Zhdankin et al. (2017) Zhdankin, V., Werner, G. R., Uzdensky, D. A., & Begelman, M. C. 2017, Phys. Rev. Lett., 118, 055103