# Training embedding quantum kernels with data re-uploading quantum neural
networks

Pablo Rodriguez-Grasa  Department of Physical Chemistry, University of the
Basque Country UPV/EHU, Apartado 644, 48080 Bilbao, Spain  EHU Quantum Center,
University of the Basque Country UPV/EHU, Apartado 644, 48080 Bilbao, Spain
TECNALIA, Basque Research and Technology Alliance (BRTA), 48160 Derio, Spain
<EMAIL_ADDRESS>Yue Ban  Departamento de Física, Universidad
Carlos III de Madrid, Avda. de la Universidad 30, 28911 Leganés, Spain
TECNALIA, Basque Research and Technology Alliance (BRTA), 48160 Derio, Spain
Mikel Sanz  Department of Physical Chemistry, University of the Basque Country
UPV/EHU, Apartado 644, 48080 Bilbao, Spain  EHU Quantum Center, University of
the Basque Country UPV/EHU, Apartado 644, 48080 Bilbao, Spain  IKERBASQUE,
Basque Foundation for Science, Plaza Euskadi 5, 48009, Bilbao, Spain  Basque
Center for Applied Mathematics (BCAM), Alameda de Mazarredo, 14, 48009 Bilbao,
Spain<EMAIL_ADDRESS>
###### Abstract

Kernel methods play a crucial role in machine learning and the Embedding
Quantum Kernels (EQKs), an extension to quantum systems, have shown very
promising performance. However, choosing the right embedding for EQKs is
challenging. We address this by proposing a $p$-qubit Quantum Neural Network
(QNN) based on data re-uploading to identify the optimal $q$-qubit EQK for a
task ($p$-to-$q$). This method requires constructing the kernel matrix only
once, offering improved efficiency. In particular, we focus on two cases:
$n$-to-$n$, where we propose a scalable approach to train an $n$-qubit QNN,
and $1$-to-$n$, demonstrating that the training of a single-qubit QNN can be
leveraged to construct powerful EQKs.

##  I Introduction

Quantum computing is a promising computational paradigm for tackling some
complex computational problems which are classically intractable. In
particular, its potential in enhancing machine learning tasks has garnered
significant attention Biamonte  _et al._ (2017); Carleo  _et al._ (2019);
Dunjko and Briegel (2018); Dalzell  _et al._ (2023) with parametrized quantum
circuits as the most common approach Benedetti  _et al._ (2019); Skolik  _et
al._ (2022); Schuld  _et al._ (2020); Farhi and Neven (2018). Although there
are evidences of quantum advantage in some tailored problems Sweke  _et al._
(2021); Jerbi  _et al._ (2021); Pirnay  _et al._ (2023); Gyurik and Dunjko
(2023), advantage over classical counterparts for practical applications
remains as an area of active research.

Amidst the exploration of quantum machine learning models, previous studies,
particularly in Refs. Havlíček  _et al._ (2019); Schuld and Killoran (2019),
have delineated a categorization into explicit and implicit models. In
explicit models, data undergoes encoding into a quantum state, and then a
parametrized measurement is performed. In contrast, implicit or kernel models
are based on a weighted summation of inner products between encoded data
points. A specialized category within parametrized quantum circuits, which can
be considered separately, consists of data re-uploading models Pérez-Salinas
_et al._ (2020). This architecture features an alternation between encoding
and processing unitaries, yielding expressive models that have found extensive
use Schuld  _et al._ (2021a); Caro  _et al._ (2021); Ono  _et al._ (2023).
However, Ref. Jerbi  _et al._ (2023) shows that these models can be mapped
onto an explicit model, thereby unifying these three approaches within the
framework of linear models within Hilbert space.

Beyond this overarching theoretical framework, significant effort has been
devoted to studying the performance of quantum kernel methods Huang  _et al._
(2021); Peters  _et al._ (2021); Bartkiewicz  _et al._ (2020); Kusumoto  _et
al._ (2021). Furthermore, theorizing about them as a way to explain quantum
machine learning models has also been pursued Schuld and Killoran (2019);
Schuld (2021); Schuld  _et al._ (2021b). This interest is substantiated by
several compelling factors. Firstly, the act of embedding data into quantum
states provides direct access to the exponentially vast quantum Hilbert space,
where inner products can be efficiently computed. Secondly, the creation of
embedding quantum kernels that pose classical intractability yet offers the
potential for quantum advantages Liu  _et al._ (2021). Finally, in both
classical machine learning and kernel methods, the representer theorem
Schölkopf and Smola (2001) ensures that they consistently achieve a training
loss lower than or equal to that of explicit models, provided the same
encoding and training set are considered. However, it is worth noting that
this enhanced expressivity may be accompanied by reduced generalization
capability, a point explored and discussed in Ref. Jerbi  _et al._ (2023).

While the arguments above might suggest a preference for kernel or implicit
models over explicit ones, practical challenges persist. Firstly, the
computational complexity associated with constructing the kernel matrix scales
quadratically with the number of training samples. Additionally, selecting the
appropriate embedding that defines the kernel function is contingent on the
specific problem at hand, demanding careful consideration. Regarding the last
one, two approaches to train kernels customized for the dataset in question
have been considered: multiple kernel learning and kernel target alignment. In
the former, a combination of kernel functions is trained to minimize an
empirical risk function Vedaie  _et al._ (2020); Ghukasyan  _et al._ (2023).
In the latter, it is defined a kernel using an embedding ansatz which is
trained to align with an ideal kernel Hubregtsen  _et al._ (2022). However,
these techniques require the computation of the kernel matrix at each training
iteration, incurring substantial computational costs.

In this Article, we employ a quantum neural network (QNN) to create an
embedding quantum kernel (EQK). This strategy utilizes QNN training to
construct the corresponding kernel matrix, significantly reducing the
computational cost when compared to constructing the kernel matrix at each
training step. In particular, we utilize the data re-uploading architecture
for the QNN and propose a method to straightforwardly scale the training up to
an $p-$qubit QNN. Following training on a given dataset, this QNN is used to
generate a $q-$qubit EQK. This construction is denoted as the $p$-to-$q$
approach. Here, we will focus on two specific cases: the $n$-to-$n$ and the
$1$-to-$n$. The former is fundamentally relevant, since we can prove that the
classification accuracy of the QNN never decreases with $n$ and that
entanglement plays a role in this accuracy. Additionally, we have empirically
observed that the classification accuracy of the EQK not only is higher than
the one of the QNN, as expected, but also that it does not decrease with $n$.
We also propose the use of the training of a single-qubit QNN to create an
$n$-qubit EQK with entanglement, referred to as the $1$-to-$n$ approach. We
substantiate our proposals with numerical evidence, demonstrating that the
training of a small-sized QNN, potentially even simulatable classically,
adeptly selects parameters to form a potent EQK tailored for a specific task.
Our results consider up to 10 qubits, showing that our approach circumvents
typical training challenges associated with scaling parametrized quantum
circuits.

The article is organized into the following sections. We begin by introducing
quantum kernel methods and EQKs. Afterwards, we present the data re-uploading
architecture and our approach to scaling it up to an $n$-qubit QNN.
Subsequently, we detail the construction of EQKs derived from the training of
these QNNs. Finally, we present numerical results supporting the
practicability of our approach.

##  II Quantum kernel methods

Quantum machine learning models can be categorized into explicit and implicit
models, as discussed in Refs. Schuld and Killoran (2019); Havlíček  _et al._
(2019). Explicit models are based on parametrized quantum circuits Benedetti
_et al._ (2019); Abbas  _et al._ (2021); Cerezo  _et al._ (2021); Bharti  _et
al._ (2022), which are tailored to process classical information $\bm{x}$
encoded into a quantum state as $\rho(\bm{x})=S(\bm{x})|0\rangle\langle
0|^{\otimes n}S^{\dagger}(\bm{x})$, where $S$ refers to some quantum
embedding. These models have proven resilient to noise Sharma  _et al._
(2020), making them attractive for seeking quantum advantages in near-term
quantum processors. However, they face scalability challenges due to the
barren plateau problem during training McClean  _et al._ (2018); Cerezo  _et
al._ (2021); Holmes  _et al._ (2022); Anschuetz and Kiani (2022); Ragone  _et
al._ (2023). Implicit or kernel models, on the other hand, are defined by the
linear combination

$f_{\bm{\alpha},X}=\sum_{i=1}^{M}\alpha_{i}\;k(\bm{x},\bm{x}_{i}),$ (1)

where $k(\bm{x},\bm{x}_{i})=\mathrm{tr}(\rho(\bm{x})\rho(\bm{x}_{i}))$ is the
kernel function and $\\{\bm{x}_{i}\\}_{i=1}^{M}$ denotes training points from
a training set $X$. An important insight from classical machine learning
theory is captured by the representer theorem. This theorem asserts that, when
provided with a feature encoding and a training set, implicit models in the
form of Eq. 1 consistently attain training losses that are either equal to or
lower than those of explicit quantum models.

Determining the optimal $\bm{\alpha}$ parameters for implicit models involves
solving a convex optimization problem, necessitating the construction of the
kernel matrix $K$ with entries $k_{ij}=k(\bm{x}_{i},\bm{x}_{j})$. To
accomplish this, inner products between all training points must be evaluated
on a quantum computer Buhrman  _et al._ (2001); Fanizza  _et al._ (2020);
Cincio  _et al._ (2018), involving $\mathcal{O}(M^{2})$ evaluations. Although
there are instances where these inner products can be computed without
explicitly constructing the feature map, our primary focus here is on
embedding quantum kernels (EQKs). EQKs, prevalent in quantum kernel literature
and proven to be universal in Ref. Gil-Fuster  _et al._ (2023), are
constructed as follows:

$k(\bm{x}_{i},\bm{x}_{j})=|\langle\phi(\bm{x}_{i})|\phi(\bm{x}_{j})\rangle|^{2}=|\langle\bm{0}|S^{\dagger}(\bm{x}_{i})S(\bm{x}_{j})|\bm{0}\rangle|^{2},$
(2)

where $|\bm{0}\rangle:=|0\rangle^{\otimes n}$. These overlap evaluations serve
as input for a support vector machine (SVM) algorithm. Although quantum SVM
approaches have been explored Rebentrost  _et al._ (2014); Park  _et al._
(2023), we assume classical computation for this step in this work, incurring
a time complexity of $\mathcal{O}(M^{3})$.

Given the representer theorem and the inherent capability of quantum devices
to access exponentially large feature spaces, quantum kernel methods emerge as
promising candidates for optimal quantum machine learning models Schuld
(2021). However, they come with their challenges, such as the exponentially
vanishing of kernel values Thanasilp  _et al._ (2022); Kübler  _et al._
(2021), and, as addressed in this work, the critical issue of selecting the
appropriate kernel for a specific task.

To identify the optimal kernel, the common strategy involves employing a
parametrized quantum embedding $S_{\bm{\theta}}(\bm{x})$ Lloyd  _et al._
(2020), where $\bm{\theta}$ represents the trainable parameters. This
embedding defines a trainable quantum kernel $k_{\bm{\theta}}(\bm{x})$, which
is then optimized based on a specific figure of merit. In the study by
Hubregtsen et al. Hubregtsen  _et al._ (2022), the optimization problem is
formulated using kernel target alignment. Another avenue explored in works
such as Vedaie  _et al._ (2020); Ghukasyan  _et al._ (2023) involves multiple
kernel learning, where the goal is to determine the optimal combination of
different kernels for a specific task. However, these approaches require
constructing the kernel matrix at each training step, which implies a high
computational cost. In our work, we propose a novel method for training
embedding quantum kernels that only necessitates constructing the kernel
matrix once.

Figure 1: Iterative training of a two-qubit data re-uploading QNN. In Step 1,
a single-qubit QNN is trained to obtain optimal model parameters denoted as
$\bm{\theta}^{*}$. Moving to Step 2, training for the two-qubit QNN is
initiated, initializing new extra parameters to 0, while the parameters of the
first qubit are set to $\bm{\theta}^{*}$. This iterative approach can scale
the QNN up to $n\leq L+1$ qubits, ensuring that the $n$-qubit QNN performs at
least as effectively as the $n-1$ version.

##  III Scaling data re-uploading for $n-$qubit QNN

Figure 2: An schematic illustration of how a single-qubit EQK, constructed
from a single-qubit QNN, can enhance classification outcomes. The QNN aims to
group points of the same class while fixing the decision plane (left Bloch
sphere). The SVM, employing the single-qubit EQK, fine-tunes the decision
boundary parameters to find the optimal hyperplane (right Bloch sphere). As a
result, data points that were previously misclassified can now be correctly
assigned to their respective labels. We also provide an example demonstrating
a scenario in which the QNN accuracy plateaus out for a specific dataset. By
using the corresponding EQK, we obtain however an increase in accuracy,
denoted as $\Delta\mathrm{Acc}$.

As reported by Pérez-Salinas et al. in Ref. Pérez-Salinas  _et al._ (2020),
the data re-uploading model incorporates layers composed of data-encoding and
training unitaries. This approach effectively introduces non-linearities to
the model allowing to capture complex patterns on data Schuld  _et al._
(2021a); Casas and Cervera-Lierta (2023); Caro  _et al._ (2021); Barthe and
Pérez-Salinas (2023). In fact, it has been demonstrated that a single qubit
quantum classifier possesses universal capabilities Pérez-Salinas  _et al._
(2021).

While various encoding strategies could be considered for this architecture,
we specifically adopt the easiest one defining

$\mathrm{QNN}_{\bm{\theta}}(\bm{x})\equiv\prod_{l=1}^{L}U(\bm{\theta}_{l})\;U(\bm{x})=U(\bm{\theta}_{L})\;U(\bm{x})\dots
U(\bm{\theta}_{1})\;U(\bm{x}).$ (3)

Here, $L$ denotes the number of layers, $U$ represents a generic
$\mathrm{SU}(2)$ unitary, and the vector
$\bm{\theta}=\\{\bm{\theta}_{1},...,\bm{\theta}_{L}\\}$ encompasses the
trainable parameters. To leverage this model to construct a binary classifier,
one must select two label states that are maximally separated in Hilbert
space. The training objective involves instructing the model to collectively
rotate points belonging to the same class, bringing them closer to their
corresponding label state.

Starting with the data re-uploading single-qubit QNN architecture, we can
naturally extend it to create multi-qubit QNNs. The introduction of more
qubits enhances the model’s expressivity by increasing the number of trainable
parameters per layer and offering the potential for entanglement between
qubits which enhances the expressivity of the model Panadero  _et al._ (2023).

In this work, we propose an iterative training approach for multi-qubit QNNs.
In our construction, the $n$-qubit QNN is defined as

$\mathrm{QNN}_{\bm{\theta},\bm{\varphi}}(\bm{x})=\prod_{l=1}^{L}\left(\prod_{s=1}^{n-1}\mathrm{CU}^{s}_{s+1}(\bm{\varphi}_{l}^{(s)})\;\left(\bigotimes_{r=1}^{n}U(\bm{\theta}^{(r)}_{l})\right)\;U(\bm{x})^{\otimes
n}\right),$ (4)

where $\mathrm{CU}^{s}_{s+1}$ denotes the controlled version of the general
$\mathrm{SU}(2)$ unitary with control in the $(s+1)$-th qubit and target in
the $s$-th qubit, and $\bm{\theta}$ and $\bm{\varphi}$ refer to the trainable
parameters of single-qubit and two-qubit gates, respectively. The total number
of trainable parameters in this architecture is $3(2n-1)L$.

For the training of the $n$-qubit QNN, we propose an iterative construction
starting from a single-qubit QNN. The initialization process is depicted in
Figure 1. Initially, we train a single-qubit QNN and utilize its parameters to
initialize the two-qubit QNN. During the initialization of the two-qubit QNN,
we set $\bm{\varphi}^{(1)}_{l}=\bm{0}$ for all $l\in[1,L]$, initializing the
parameters of the first qubit with those obtained from the training of the
single-qubit step. Consequently, the entangling layers do not have any action,
ensuring that, with a local measurement on the first qubit, we commence in the
output state of the single-qubit QNN training.

This process can be employed to scale up the QNN architecture, allowing the
construction of QNNs with up to $n=L+1$ qubits. When adding an extra qubit,
the supplementary entangling gates are initialized as identities, and the
training begins with the optimal parameters obtained from the previous step.
Essentially, this formalized approach signifies a systematic and scalable
improvement in the QNN’s performance with the incorporation of each additional
qubit.

Figure 3: Embedding quantum kernels generated from quantum neural networks
training. The kernel matrix element $k_{ij}$ is defined as the probability of
measuring all qubits in the state $|0\rangle$, denoted as $P_{\bm{0}}$. On the
left, we have the $n$-to-$n$ proposal, constructed by directly utilizing the
trained data re-uploading $n$-qubit QNN as a quantum feature map. On the
right, we show the construction of an Embedded EQK from the training of a
single-qubit QNN, named as $1$-to-$n$.

##  IV Constructing EQKs from QNN training

In this section, we propose a method for constructing trained embedding
quantum kernels (EQK) using data re-uploading quantum neural networks (QNN).
The idea is to train the QNN for a classification task and leverage its
architecture to generate EQKs tailored to the specific task, thereby enhancing
performance on the given dataset. The motivation for combining these two
binary classification approaches, QNN and kernel methods, arises from two
perspectives.

Firstly, we investigate whether the QNN can effectively select a suitable
embedding kernel for a specific task. This approach may lead to more efficient
kernel training compared to previous methods, as we only need to construct the
kernel matrix once. Secondly, the QNN’s performance is contingent on its
training. Utilizing the corresponding EQK construction might produce superior
results compared to relying solely on the QNN, even in cases where the
training process has not been optimal. Numerical results illustrating this can
be found in Appendix E.

The underlying concept of this method is formalized in the Appendix and
clarified in Figure 2, specifically for a single-qubit scenario. Initially,
the QNN rotates data points of the same class near their label state while
maintaining the decision hyperplane fixed. This hyperplane is taken as the
equator of the Bloch sphere, i.e., $\langle\hat{\sigma}_{z}\rangle=0$. After
training the QNN, the resulting feature map is obtained by preserving the
parameters acquired during training. This feature map is then utilized to
construct an EQK. The subsequent application of the SVM algorithm, using this
kernel, aims to identify the optimal separation hyperplane in the feature
space. In this context, the optimization is equivalent to adjusting the
optimal measurement $\mathcal{M}(\bm{\gamma})$ while keeping the data points
fixed.

Building upon this insight, we propose constructing EQKs using QNN training in
two ways: the $n$-to-$n$ approach and the $1$-to-$n$. In the first one, an
$n$-qubit QNN is trained using the iterative method proposed earlier to
directly construct the corresponding EQK of $n$ qubits. As depicted in Figure
3, a multi-qubit QNN of the form of Eq. 4 is trained while fixing the
trainable parameters of the embedding $\bm{\theta}^{*}$ and
$\bm{\varphi}^{*}$. Then, these parameters are subsequently used to construct
an EQK which is defined by

$k_{ij}=|\langle\bm{0}|\mathrm{QNN}_{\bm{\theta}^{*},\bm{\varphi}^{*}}(\bm{x}_{i})^{\dagger}\;\mathrm{QNN}_{\bm{\theta}^{*},\bm{\varphi}^{*}}(\bm{x}_{j})|\bm{0}\rangle|^{2}.$
(5)

This method allows us to scale the QNN as much as possible during training,
and when it reaches a performance plateau, we can utilize the trained feature
map to construct an EQK.

In the $1$-to-$n$ construction, we train a single-qubit
$\mathrm{QNN}_{\bm{\theta}}$ from Eq. 3, fixing the $\bm{\theta}^{*}$ and we
leverage this training to construct the kernel

$\begin{split}k_{ij}=|\langle\bm{0}|&\prod_{l=1}^{L-1}\left(U(\bm{x}_{i})^{\dagger\otimes
n}\;U(\bm{\theta}^{*}_{l})^{\dagger\otimes
n}\;E\right)U(\bm{x}_{i})^{\dagger\otimes n}\\\ &U(\bm{x}_{j})^{\otimes
n}\left(\prod_{l=1}^{L-1}E\;U(\bm{\theta}^{*}_{l})\;U(\bm{x}_{j})\right)|\bm{0}\rangle|^{2},\end{split}$
(6)

where $E$ denotes an entangling operation, such as a cascade of
$\mathrm{CNOT}$ or $\mathrm{CZ}$ gates, among other possibilities. It is worth
noting that the training is conducted for a single-qubit QNN and does not
explicitly consider entanglement. Nevertheless, we will present numerical
results demonstrating that this training alone is sufficient to select
parameters for constructing a customized multi-qubit kernel tailored to a
specific task.

Certainly, one can combine both the $n$-to-$n$ and the $1$-to-$n$
architectures and generalize it to $n$-to-$m\cdot n$, where $m$ represents
some integer. In this scenario, an $n$-qubit QNN is trained and utilized to
implement the same design as in the $1$-to-$n$ construction. However, in this
case, each qubit of the QNN is embedded into $m$ qubits, introducing
entanglement between layers.

##  V Numerical results

Let us now explore the practical implementation of our approaches and evaluate
their performance across various datasets detailed in Appendix F. Both the
training and test sets consist of 500 data points each.

For the QNN training in this section, we employ the Adam optimizer with a
batch size of 24. In the initial step for the QNN training ($n=1$), we use a
learning rate of 0.05 and consider 30 epochs. For $n>1$, we use a learning
rate of 0.005 and 10 epochs.

The results presented in Figure 4 demonstrate the test accuracies for two
distinct datasets—the corners dataset and the circles dataset—for the
$1$-to-$n$ and the $n$-to-$n$ constructions as functions of the number of
qubits $n$. These experiments consider $L=7$ layers, scaling up the QNN to
$n=L+1=8$. In the $1$-to-$n$ construction, we introduce entanglement between
layers using a cascade of CNOT gates, defined as

$E=\prod_{s=1}^{n-1}\mathrm{CNOT}^{s+1}_{s},$ (7)

where the subscript refers to the control qubit, and the superscript refers to
the target qubit.

The numerical experiments illustrate a consistent improvement in the
performance of the $n$-qubit QNN as more qubits are introduced for both
datasets. Additionally, the corresponding kernel, represented by the
$n$-to-$n$ construction, exhibits enhanced performance compared to the QNN
alone. Interestingly, the $1$-to-$n$ construction performs remarkably well,
showing a rapid increase in performance with a small number of qubits,
although it eventually plateaus. While the improvement in the $n$-to-$n$
scenario is slower, its performance consistently surpasses that of the QNN,
and the scaling of this model is consistent and rigorous. For additional
numerical experiments, refer to Appendix G.

These results demonstrate, on one hand, the scalability of a data re-uploading
QNN and the improved performance achieved by constructing the corresponding
EQK. On the other hand, they highlight the capability to create complex
kernels through the training of simple QNN architectures. Despite the small
QNNs being trainable using a classical computer, the obtained parameters could
be utilized to construct a classically hard embedding quantum kernel.

Figure 4: Test accuracies for the two EQK architectures presented in this
work, namely $1$-to-$n$ and $n$-to-$n$, plotted against the number of qubits
$n$. Additionally, we include the test accuracy of the $n$-qubit QNN using the
iterative training proposed. Results are showcased for the corners dataset in
the top plot and the circles dataset in the bottom plot.

##  VI Conclusions

We have addressed the challenge of training a $q$-qubit embedding quantum
kernel by introducing a practical approach that constructs them from a data
re-uploading $p$-qubit quantum neural network. This method significantly
reduces computational costs compared to previous techniques in which the
kernel matrix was constructed at each training step. We focus on two specific
cases for constructing embedding quantum kernels, namely, the $n$-to-$n$ and
the $1$-to-$n$ cases. The first involves using the output of the training of
an $n$-qubit QNN directly as a feature map (referred to as the $n$-to-$n$
approach), for which we have proposed a scalable training method. The second
approach enables the creation of multi-qubit EQKs with entanglement from the
training of a single-qubit QNN (referred to as the $1$-to-$n$ approach). Given
the current limitations in scaling parametrized quantum circuits, our combined
approach offers an alternative for creating trainable scalable kernels as
quantum machine learning models, mitigating the challenges associated with
their scalability and reducing the dependence on the training for model
success. Moreover, the $1$-to-$n$ approach can be considered a quantum
inspired proposal, paving the way for exploring the construction of
classically hard-to-estimate kernels from the training of small architectures,
even those that can be simulated with a classical computer.

##  VII Acknowledgements

The authors would like to thank Maria Schuld for her valuable comments which
led to a more solid presentation of the results and Sofiene Yerbi and Adrián
Pérez-Salinas for their worthy insights during QTML 2023. PR and YB
acknowledge financial support from the CDTI within the Misiones 2021 program
and the Ministry of Science and Innovation under the Recovery, Transformation
and Resilience Plan—Next Generation EU under the project “CUCO: Quantum
Computing and its Application to Strategic Industries”. PR and MS acknowledge
support from EU FET Open project EPIQUS (899368) and HORIZON-CL4-
2022-QUANTUM01-SGA project 101113946 OpenSuperQPlus100 of the EU Flagship on
Quantum Technologies, the Spanish Ramón y Cajal Grant RYC-2020-030503-I,
project Grant No. PID2021-125823NA-I00 funded by
MCIN/AEI/10.13039/501100011033 and by “ERDF A way of making Europe” and “ERDF
Invest in your Future”, and from the IKUR Strategy under the collaboration
agreement between Ikerbasque Foundation and BCAM on behalf of the Department
of Education of the Basque Government. This project has also received support
from the Spanish Ministry of Economic Affairs and Digital Transformation
through the QUANTUM ENIA project call - Quantum Spain, and by the EU through
the Recovery, Transformation and Resilience Plan - NextGenerationEU. YB
acknowledges support from the Spanish Government via the project
PID2021-126694NA- C22 (MCIU/AEI/FEDER, EU).

## References

  * Biamonte  _et al._ (2017) J. Biamonte, P. Wittek, N. Pancotti, P. Rebentrost, N. Wiebe, and S. Lloyd, Nature 549, 195 (2017).
  * Carleo  _et al._ (2019) G. Carleo, I. Cirac, K. Cranmer, L. Daudet, M. Schuld, N. Tishby, L. Vogt-Maranto, and L. Zdeborová , Reviews of Modern Physics 91 (2019), 10.1103/revmodphys.91.045002.
  * Dunjko and Briegel (2018) V. Dunjko and H. J. Briegel, Reports on Progress in Physics 81, 074001 (2018).
  * Dalzell  _et al._ (2023) A. M. Dalzell, S. McArdle, M. Berta, P. Bienias, C.-F. Chen, A. Gilyén, C. T. Hann, M. J. Kastoryano, E. T. Khabiboulline, A. Kubica, G. Salton, S. Wang, and F. G. S. L. Brandão, “Quantum algorithms: A survey of applications and end-to-end complexities,” (2023), arXiv:2310.03011 [quant-ph] .
  * Benedetti  _et al._ (2019) M. Benedetti, E. Lloyd, S. Sack, and M. Fiorentini, Quantum Science and Technology 4, 043001 (2019).
  * Skolik  _et al._ (2022) A. Skolik, S. Jerbi, and V. Dunjko, Quantum 6, 720 (2022).
  * Schuld  _et al._ (2020) M. Schuld, A. Bocharov, K. M. Svore, and N. Wiebe, Physical Review A 101 (2020), 10.1103/physreva.101.032308.
  * Farhi and Neven (2018) E. Farhi and H. Neven, “Classification with quantum neural networks on near term processors,” (2018), arXiv:1802.06002 [quant-ph] .
  * Sweke  _et al._ (2021) R. Sweke, J.-P. Seifert, D. Hangleiter, and J. Eisert, Quantum 5, 417 (2021).
  * Jerbi  _et al._ (2021) S. Jerbi, L. M. Trenkwalder, H. Poulsen Nautrup, H. J. Briegel, and V. Dunjko, PRX Quantum 2 (2021), 10.1103/prxquantum.2.010328.
  * Pirnay  _et al._ (2023) N. Pirnay, R. Sweke, J. Eisert, and J.-P. Seifert, Physical Review A 107 (2023), 10.1103/physreva.107.042416.
  * Gyurik and Dunjko (2023) C. Gyurik and V. Dunjko, “On establishing learning separations between classical and quantum machine learning with classical data,” (2023), arXiv:2208.06339 [quant-ph] .
  * Havlíček  _et al._ (2019) V. Havlíček, A. D. Córcoles, K. Temme, A. W. Harrow, A. Kandala, J. M. Chow, and J. M. Gambetta, Nature 567, 209 (2019).
  * Schuld and Killoran (2019) M. Schuld and N. Killoran, Physical review letters 122, 040504 (2019).
  * Pérez-Salinas  _et al._ (2020) A. Pérez-Salinas, A. Cervera-Lierta, E. Gil-Fuster, and J. I. Latorre, Quantum 4, 226 (2020).
  * Schuld  _et al._ (2021a) M. Schuld, R. Sweke, and J. J. Meyer, Physical Review A 103 (2021a), 10.1103/physreva.103.032430.
  * Caro  _et al._ (2021) M. C. Caro, E. Gil-Fuster, J. J. Meyer, J. Eisert, and R. Sweke, Quantum 5, 582 (2021).
  * Ono  _et al._ (2023) T. Ono, W. Roga, K. Wakui, M. Fujiwara, S. Miki, H. Terai, and M. Takeoka, Physical Review Letters 131 (2023), 10.1103/physrevlett.131.013601.
  * Jerbi  _et al._ (2023) S. Jerbi, L. J. Fiderer, H. Poulsen Nautrup, J. M. Kübler, H. J. Briegel, and V. Dunjko, Nature Communications 14, 517 (2023).
  * Huang  _et al._ (2021) H.-Y. Huang, M. Broughton, M. Mohseni, R. Babbush, S. Boixo, H. Neven, and J. R. McClean, Nature Communications 12 (2021), 10.1038/s41467-021-22539-9.
  * Peters  _et al._ (2021) E. Peters, J. Caldeira, A. Ho, S. Leichenauer, M. Mohseni, H. Neven, P. Spentzouris, D. Strain, and G. N. Perdue, “Machine learning of high dimensional data on a noisy quantum processor,” (2021), arXiv:2101.09581 [quant-ph] .
  * Bartkiewicz  _et al._ (2020) K. Bartkiewicz, C. Gneiting, A. Černoch, K. Jiráková, K. Lemr, and F. Nori, Scientific Reports 10 (2020), 10.1038/s41598-020-68911-5.
  * Kusumoto  _et al._ (2021) T. Kusumoto, K. Mitarai, K. Fujii, M. Kitagawa, and M. Negoro, npj Quantum Information 7 (2021), 10.1038/s41534-021-00423-0.
  * Schuld (2021) M. Schuld, arXiv preprint arXiv:2101.11020 (2021).
  * Schuld  _et al._ (2021b) M. Schuld, F. Petruccione, M. Schuld, and F. Petruccione, Machine Learning with Quantum Computers , 217 (2021b).
  * Liu  _et al._ (2021) Y. Liu, S. Arunachalam, and K. Temme, Nature Physics 17, 1013 (2021).
  * Schölkopf and Smola (2001) B. Schölkopf and A. Smola,  _Smola, A.: Learning with Kernels - Support Vector Machines, Regularization, Optimization and Beyond. MIT Press, Cambridge, MA_ , Vol. 98 (2001).
  * Vedaie  _et al._ (2020) S. S. Vedaie, M. Noori, J. S. Oberoi, B. C. Sanders, and E. Zahedinejad, “Quantum multiple kernel learning,” (2020), arXiv:2011.09694 [quant-ph] .
  * Ghukasyan  _et al._ (2023) A. Ghukasyan, J. S. Baker, O. Goktas, J. Carrasquilla, and S. K. Radha, “Quantum-classical multiple kernel learning,” (2023), arXiv:2305.17707 [quant-ph] .
  * Hubregtsen  _et al._ (2022) T. Hubregtsen, D. Wierichs, E. Gil-Fuster, P.-J. H. S. Derks, P. K. Faehrmann, and J. J. Meyer, Physical Review A 106 (2022), 10.1103/physreva.106.042431.
  * Abbas  _et al._ (2021) A. Abbas, D. Sutter, C. Zoufal, A. Lucchi, A. Figalli, and S. Woerner, Nature Computational Science 1, 403–409 (2021).
  * Cerezo  _et al._ (2021) M. Cerezo, A. Sone, T. Volkoff, L. Cincio, and P. J. Coles, Nature Communications 12 (2021), 10.1038/s41467-021-21728-w.
  * Bharti  _et al._ (2022) K. Bharti, A. Cervera-Lierta, T. H. Kyaw, T. Haug, S. Alperin-Lea, A. Anand, M. Degroote, H. Heimonen, J. S. Kottmann, T. Menke, _et al._ , Reviews of Modern Physics 94, 015004 (2022).
  * Sharma  _et al._ (2020) K. Sharma, S. Khatri, M. Cerezo, and P. J. Coles, New Journal of Physics 22, 043006 (2020).
  * McClean  _et al._ (2018) J. R. McClean, S. Boixo, V. N. Smelyanskiy, R. Babbush, and H. Neven, Nature Communications 9 (2018), 10.1038/s41467-018-07090-4.
  * Holmes  _et al._ (2022) Z. Holmes, K. Sharma, M. Cerezo, and P. J. Coles, PRX Quantum 3 (2022), 10.1103/prxquantum.3.010313.
  * Anschuetz and Kiani (2022) E. R. Anschuetz and B. T. Kiani, Nature Communications 13 (2022), 10.1038/s41467-022-35364-5.
  * Ragone  _et al._ (2023) M. Ragone, B. N. Bakalov, F. Sauvage, A. F. Kemper, C. O. Marrero, M. Larocca, and M. Cerezo, “A unified theory of barren plateaus for deep parametrized quantum circuits,” (2023), arXiv:2309.09342 [quant-ph] .
  * Buhrman  _et al._ (2001) H. Buhrman, R. Cleve, J. Watrous, and R. de Wolf, Physical Review Letters 87 (2001), 10.1103/physrevlett.87.167902.
  * Fanizza  _et al._ (2020) M. Fanizza, M. Rosati, M. Skotiniotis, J. Calsamiglia, and V. Giovannetti, Physical Review Letters 124 (2020), 10.1103/physrevlett.124.060503.
  * Cincio  _et al._ (2018) L. Cincio, Y. Subaşı, A. T. Sornborger, and P. J. Coles, New Journal of Physics 20, 113022 (2018).
  * Gil-Fuster  _et al._ (2023) E. Gil-Fuster, J. Eisert, and V. Dunjko, arXiv preprint arXiv:2309.14419 (2023).
  * Rebentrost  _et al._ (2014) P. Rebentrost, M. Mohseni, and S. Lloyd, Physical review letters 113, 130503 (2014).
  * Park  _et al._ (2023) S. Park, D. K. Park, and J.-K. K. Rhee, Scientific Reports 13, 3288 (2023).
  * Thanasilp  _et al._ (2022) S. Thanasilp, S. Wang, M. Cerezo, and Z. Holmes, “Exponential concentration and untrainability in quantum kernel methods,” (2022), arXiv:2208.11060 [quant-ph] .
  * Kübler  _et al._ (2021) J. M. Kübler, S. Buchholz, and B. Schölkopf, “The inductive bias of quantum kernels,” (2021), arXiv:2106.03747 [quant-ph] .
  * Lloyd  _et al._ (2020) S. Lloyd, M. Schuld, A. Ijaz, J. Izaac, and N. Killoran, “Quantum embeddings for machine learning,” (2020), arXiv:2001.03622 [quant-ph] .
  * Casas and Cervera-Lierta (2023) B. Casas and A. Cervera-Lierta, Physical Review A 107 (2023), 10.1103/physreva.107.062612.
  * Barthe and Pérez-Salinas (2023) A. Barthe and A. Pérez-Salinas, “Gradients and frequency profiles of quantum re-uploading models,” (2023), arXiv:2311.10822 [quant-ph] .
  * Pérez-Salinas  _et al._ (2021) A. Pérez-Salinas, D. López-Núñez, A. García-Sáez, P. Forn-Díaz, and J. I. Latorre, Physical Review A 104 (2021), 10.1103/physreva.104.012405.
  * Panadero  _et al._ (2023) I. Panadero, Y. Ban, H. Espinós, R. Puebla, J. Casanova, and E. Torrontegui, “Regressions on quantum neural networks at maximal expressivity,” (2023), arXiv:2311.06090 [quant-ph] .
  * Holmes and Jain (2006) D. E. Holmes and L. C. Jain,  _Innovations in machine learning_ (Springer, 2006).

##  Appendix A Related work on quantum kernel training

One of the main drawbacks of quantum kernel methods is determining the best
kernel which depends on the specific problem. In cases where no prior
knowledge about the input data is available, approaches have been developed to
construct parametrized embedding quantum kernels that are trained for a
specific task. These kernel training methods are based on multiple kernel
learning Vedaie  _et al._ (2020); Ghukasyan  _et al._ (2023) and kernel target
alignment Hubregtsen  _et al._ (2022).

###  A.1 Multiple kernel learning

Multiple kernel learning methods involve the creation of a combined kernel

$k_{\bm{\varphi},\bm{\omega}}(\bm{x}_{i},\bm{x}_{j})=f_{\bm{\omega}}\left(\\{k_{\bm{\varphi}}^{(r)}(\bm{x}_{i},\bm{x}_{j})\\}_{r=1}^{R}\right),$
(8)

where $k_{\bm{\varphi}}^{(r)}$ represents a set of parametrized embedding
quantum kernels, and $f_{\bm{\omega}}$ is a combination function. Two common
choices for this function result in either a linear kernel combination

$k_{\bm{\varphi},\bm{\omega}}(\bm{x}_{i},\bm{x}_{j})=\sum_{r=1}^{R}\omega_{r}k_{\bm{\varphi}}^{(r)}(\bm{x}_{i},\bm{x}_{j}),$
(9)

or a multiplicative kernel combination

$k_{\bm{\varphi}}(\bm{x}_{i},\bm{x}_{j})=\prod_{r=1}^{R}k_{\bm{\varphi}}^{(r)}(\bm{x}_{i},\bm{x}_{j}).$
(10)

When it comes to choosing model parameters, in Ref. Vedaie  _et al._ (2020)
they use an empirical risk function to minimize, while in Ref. Ghukasyan  _et
al._ (2023) the authors opt for a convex minimization problem. However, in
both constructions a complete kernel matrix is computed at each optimization
step.

###  A.2 Kernel target alignment

Training kernels using kernel target alignment is based on the similarity
measure between two kernel matrices, denoted as $K_{A}$ and $K_{B}$. This
similarity measure, known as kernel alignment Holmes and Jain (2006), is
defined as

$KA(K_{A},K_{B})=\frac{\mathrm{tr}(K_{A}K_{B})}{\sqrt{\mathrm{tr}(K_{A}^{2})\cdot\mathrm{tr}(K_{B}^{2})}}.$
(11)

and corresponds to the cosine of the angle between kernel matrices if we see
them as vectors in the space of matrices with the Hilbert-Schmidt inner
product. The use of kernel alignment to train kernels is by defining the ideal
kernel matrix $K^{*}$ whose entries are
$k^{*}_{ij}=k^{*}(\bm{x}_{i},\bm{x}_{j})=y_{i}y_{j}$. Given a parametrized
quantum embedding kernel
$k_{\bm{\theta}}(\bm{x}_{i},\bm{x}_{j})=|\langle\phi_{\bm{\theta}}(\bm{x}_{i})|\phi_{\bm{\theta}}(\bm{x}_{j})\rangle|^{2}$
whose kernel matrix is named as $K_{\bm{\theta}}$, the kernel-target alignment
is defined as the kernel alignment between $K$ and the ideal kernel, i.e.

$TA(K_{\bm{\theta}})=KA(K_{\bm{\theta}},K^{*})=\frac{\sum_{ij}y_{i}y_{j}\;k_{\bm{\theta}}(\bm{x}_{i},\bm{x}_{j})}{M\sqrt{\sum_{ij}k_{\bm{\theta}}(\bm{x}_{i},\bm{x}_{j})^{2}}},$
(12)

where we used that $y_{i}^{2}=1$ independently from the label and $M$ is the
number of training points. Thus, authors in Ref. Hubregtsen  _et al._ (2022)
use this quantity as cost function which is minimized using gradient descent
over the $\bm{\theta}$ parameters. This strategy is closely related to the
approach outlined in Ref. Lloyd  _et al._ (2020), where the authors explore
the construction of quantum feature maps with the aim of maximizing the
separation between different classes in Hilbert space. However, in the case of
the kernel target alignment strategy, similar to multiple kernel learning, the
construction of the kernel matrix (or at least a subset of it) at each
training step is required.

##  Appendix B Data re-uploading

To construct a data re-uploading QNN for binary classification, each label is
associated with a unique quantum state, aiming to maximize the separation in
the Bloch sphere. For the single qubit quantum classifier, the labels +1 and
-1 are represented by the computational basis states $|0\rangle$ and
$|1\rangle$, respectively. The objective is to appropriately tune the
parameters $\\{\bm{\theta}_{l}\\}_{l=1}^{L}$ that define the state

$|\phi_{\bm{\theta}}(\bm{x}_{i})\rangle=U(\bm{\theta}_{L})\;U(\bm{x})\dots
U(\bm{\theta}_{1})\;U(\bm{x})|0\rangle,$ (13)

to rotate the data points of the same class close to their corresponding label
state. To achieve this, we use the fidelity cost function

$f_{\mathrm{cost}}=\frac{1}{M}\sum_{i=1}^{M}\Bigl{(}1-|\langle\phi^{i}_{l}|\phi_{\bm{\theta}}(\bm{x}_{i})\rangle|^{2}\Bigr{)},$
(14)

where $|\phi_{l}^{i}\rangle$ represents the correct label state for the data
point $\bm{x}_{i}$. This optimization process occurs on a classical processor.
Once the model is trained, the quantum circuit is applied to a test data point
$\bm{x}_{t}$, and the probability of obtaining one of the label states is
measured. If this probability surpasses a certain threshold (set as 1/2 in
this case), the data point is classified into the corresponding label state
class. Formally, the decision rule can be expressed as

$\hat{y}[\bm{x}_{t}]=\begin{cases}+1&\text{if }|\langle
0|\phi(\bm{x}_{t})\rangle|^{2}\geq 1/2,\\\ -1&\text{if }|\langle
0|\phi(\bm{x}_{t})\rangle|^{2}<1/2.\end{cases}$ (15)

When consider a multi-qubit architecture, the idea is the same but we need to
properly choose the corresponding label states. In this work, for a $n$-qubit
QNN we consider as label states the ones defined by the projectors
$|0\rangle\langle 0|\otimes\mathds{1}^{(n-1)}$ and $|1\rangle\langle
1|\otimes\mathds{1}^{(n-1)}$, which corresponds to a local measurement in the
first qubit.

##  Appendix C Connection with the representer theorem

In this section, we aim to leverage the representer theorem to demonstrate
that the kernels derived from the data re-uploading Quantum Neural Network
(QNN) will be at least as effective as the corresponding QNN. In our
construction, the $n$-qubit QNN is defined in Eq. 4 in the main text. There,
we can extract the last layer to separate it into a variational part, which
will be absorbed into the measurement when we define the corresponding quantum
model, and the encoding part which will define the quantum feature map for the
construction of the kernel. This corresponds to

$\mathrm{QNN}_{\bm{\theta},\bm{\varphi}}(\bm{x})=\underbrace{\prod_{s=1}^{n-1}\mathrm{CU}^{s}_{s+1}(\bm{\lambda}^{(s)})\;\left(\bigotimes_{r=1}^{n}U(\bm{\omega}^{(r)})\right)}_{=V(\bm{\lambda},\bm{\omega})}\underbrace{\prod_{l=1}^{L-1}\left(\prod_{s=1}^{n-1}\mathrm{CU}^{s}_{s+1}(\bm{\varphi}_{l}^{(s)})\;\left(\bigotimes_{r=1}^{n}U(\bm{\theta}^{(r)}_{l})\right)\;U(\bm{x})^{\otimes
n}\right)}_{=S_{\bm{\theta},\bm{\varphi}}(\bm{x})},$ (16)

where we defined $\bm{\lambda}\equiv\bm{\varphi}_{L}$ and
$\bm{\omega}\equiv\bm{\theta}_{L}$. This formulation aligns with the
definition of a quantum model from Ref. Schuld (2021),

$f(\bm{x})=\mathrm{tr}(\rho(\bm{x})\mathcal{M}).$ (17)

In our case,
$\rho(\bm{x}):=\rho_{\bm{\theta},\bm{\varphi}}(\bm{x})=S_{\bm{\theta},\bm{\varphi}}(\bm{x})|0\rangle\langle
0|^{\otimes n}S_{\bm{\theta},\bm{\varphi}}(\bm{x})^{\dagger}$, and the
variational measurement
$\mathcal{M}:=\mathcal{M}_{\bm{\lambda},\bm{\omega}}=V(\bm{\lambda},\bm{\omega})\;(\hat{\sigma}_{z}\otimes\mathds{1}^{(n-1)})\;V(\bm{\lambda},\bm{\omega})^{\dagger}$.

Once this model is trained over a dataset, we determine the parameters that
minimize the $f_{\mathrm{cost}}$ defined in the previous section, fixing
$\bm{\beta}^{*},\bm{\theta}^{*},\bm{\lambda}^{*}$, and $\bm{\omega}^{*}$. If
we now use the corresponding feature map
$S_{\bm{\theta}^{*},\bm{\varphi}^{*}}$ to construct an embedding quantum
kernel, we are effectively replacing the measurement from the optimization
$\mathcal{M}_{\bm{\lambda}^{*},\bm{\omega}^{*}}$ by the optimal measurement

$\mathcal{M}_{\mathrm{opt}}=\sum_{m=1}^{M}\alpha_{m}\;\rho_{\bm{\theta}^{*},\bm{\varphi}^{*}}(\bm{x}_{m})$
(18)

which, by the representer theorem, defines the optimal quantum model

$f_{\mathrm{opt}}(\bm{x})=\sum_{m=1}^{M}\alpha_{m}\mathrm{tr}\left(\rho_{\bm{\theta}^{*},\bm{\varphi}^{*}}(\bm{x})\;\rho_{\bm{\theta}^{*},\bm{\varphi}^{*}}(\bm{x}_{m})\right)$
(19)

that minimizes the regularized empirical risk function. Thus, we proved that
constructing the embedding quantum kernel from QNN training will perform
equally or better in terms of training loss than the QNN alone.

##  Appendix D Hyperplane defined in the Bloch sphere

The objective of the data re-uploading quantum neural network (QNN) is to map
all points with label +1 to the $|0\rangle$ state on the Bloch sphere, while
mapping points with label -1 to the $|1\rangle$ state. In an ideal scenario,
all +1 points would be rotated to the $|0\rangle$ state and all -1 points to
the $|1\rangle$ state. Using this feature map to construct the embedding
quantum kernel, all points would be support vectors associated to the same
Lagrange multiplier $\alpha$. The SVM generates a separating hyperplane
defined by the equation

$\sum_{i\in\mathrm{SV}}\alpha_{i}\;y_{i}\;k(\bm{x}_{i},\bm{x})+b=0.$ (20)

In this scenario, assuming a balanced dataset translates to an equal number of
support vectors, denoted as $N_{\mathrm{SV}}$, for each class, and results in
$b=0$ due to symmetry. Working with this equation by considering the sum of
support vectors separately for the +1 class ($\mathrm{SV}_{+1}$) and -1 class
($\mathrm{SV}_{-1}$), we obtain

$\begin{split}\sum_{i\in\mathrm{SV}_{+1}}\alpha_{i}\;y_{i}\;|\langle\phi(\bm{x}_{i})|\phi(\bm{x})\rangle|^{2}+\sum_{i\in\mathrm{SV}_{+1}}\alpha_{i}\;y_{i}\;|\langle\phi(\bm{x}_{i})|\phi(\bm{x})\rangle|^{2}=N_{\mathrm{SV}}\alpha\;|\langle
0|\phi(\bm{x}\rangle|^{2}-N_{\mathrm{SV}}\alpha\;|\langle
1|\phi(\bm{x}\rangle|^{2}=0,\end{split}$ (21)

which is equivalent to the hyperplane defined by

$|\langle 0|\phi(\bm{x}\rangle|^{2}=|\langle 1|\phi(\bm{x}\rangle|^{2},$ (22)

mirroring the decision boundary of the QNN part. Therefore, in the case of a
perfect training, the SVM becomes redundant. Even if the points are not
perfectly mapped to their label states and the separating hyperplane of the
SVM is not exactly the one from Eq. 22, as long as all points are correctly
classified by the QNN, the SVM will yield the same results. Thus, the SVM part
is only meaningful for the single-qubit QNN to construct a single-qubit EQK
case when the QNN training is suboptimal and requires to adjust the
measurement of the decision boundary.

##  Appendix E Noisy simulations

We have introduced a combined protocol for binary classification. It is worth
noting that the kernel estimation part involves the utilization of a quantum
circuit with twice the depth of the QNN part. In practical implementations on
current noisy intermediate-scale quantum devices, the consideration of a
larger circuit warrants careful attention. This is due to the susceptibility
of larger circuits to elevated noise levels, potentially affecting the overall
performance.

As previously mentioned, increasing the number of layers in the model enhances
its expressivity. However, when constructing the kernel using a high number of
layers, noise can have a detrimental impact on the performance, resulting in
poorer results for the combined protocol compared to using only the QNN.
Therefore, in this section, our objective is to perform simulations to
visualize the trade-off between the number of layers in the QNN and using the
combined protocol.

To characterize the noise incorporated into these simulations, we employ the
operator sum representation of the noise channel $\varepsilon$ acting on a
quantum state $\rho$,

$\varepsilon(\rho)=\sum_{k}E_{k}\rho E_{k}^{\dagger},$ (23)

where $E_{k}$ represents the respective Kraus operators. In our simulations,
we account for two single-qubit noise channels that are applied after each
quantum gate. Firstly, we consider amplitude damping error which is described
by the Kraus operators

$\displaystyle K_{0}$ $\displaystyle=\begin{pmatrix}1&0\\\
0&\sqrt{1-\gamma}\end{pmatrix},$ (24) $\displaystyle K_{1}$
$\displaystyle=\begin{pmatrix}0&\sqrt{\gamma}\\\ 0&0\end{pmatrix},$ (25)

Here, $\gamma\in[0,1]$ is the amplitude damping probability, which can be
defined as $p=1-e^{-\Delta t/T_{1}}$, with $T_{1}$ representing the thermal
relaxation time and $\Delta t$ denoting the duration of application of a
quantum gate. The second noise source under consideration is the phase flip
error, which can be described by the following Kraus operators

$\displaystyle K_{0}$ $\displaystyle=\sqrt{1-\alpha}\begin{pmatrix}1&0\\\
0&1\end{pmatrix},$ (26) $\displaystyle K_{1}$
$\displaystyle=\sqrt{\alpha}\begin{pmatrix}1&0\\\ 0&-1\end{pmatrix}.$ (27)

In this case, $\alpha\in[0,1]$ represents the probability of a phase flip
error, which can be defined as $2\alpha=1-e^{-\Delta t/(2T_{2})}$, with
$T_{2}$ denoting the spin-spin relaxation time or dephasing time. The value of
$\alpha$ also falls within the interval $[0,1]$.

The current state-of-the-art experimental values for the noise parameters of a
superconducting quantum processor are as follows: $T_{1}$ falls within the
range of $50-150\;\mu s$, $T_{2}$ is in the range of $25-75\;\mu s$, and
$\Delta t$ varies from $10-50\;ns$. It is important to note that noise becomes
more pronounced when $T_{1}$ and $T_{2}$ decrease and/or when $\Delta t$
increases. To explore the scenario with the worst noise, we consider the
extreme values within these ranges, leading to $\gamma^{*}=0.001$ and
$\alpha^{*}=0.0005$.

In our simulations, we simplify the noise characterization by defining a noise
strength parameter $\tau\equiv\alpha=\gamma$ for the sake of simplicity. We
choose the single-qubit data re-uploading quantum neural network (QNN) as the
kernel selection part. To observe a transition where it is no longer
advantageous to include the support vector machine (SVM) part, we explore
values of $\tau$ ranging from $0.005$ to $0.030$ in steps of $0.005$, while
also including the noiseless case where $\tau=0$. These values correspond to
examining the range $5\gamma^{*}\leq\gamma\leq 30\gamma^{*}$ and
$10\alpha^{*}\leq\alpha\leq 60\alpha^{*}$. Notably, even in these highly
adverse conditions, which are significantly worse than those of current real
hardware quantum devices, we find that the combined protocol remains suitable.
This conclusion aligns with expectations since the protocol utilizes only a
small number of qubits and quantum gates.

Figure E.1: Relative improvement, as a function of the number of layers and
noise strength, comparing the single qubit QNN accuracy performance with the
combined protocol utilizing a two-qubit embedding kernel. Results are shown
for the sinus, corners, and spiral datasets.

In Figure E.1, we present the relative improvement, denoted as

$\text{Relative
improvement}=\frac{\mathrm{Acc}_{\mathrm{QNN+SVM}}-\mathrm{Acc}_{\mathrm{QNN}}}{\mathrm{Acc}_{\mathrm{QNN}}},$
(28)

plotted as a function of the noise strength parameter $\tau$ and the number of
layers $L$. The architecture employed corresponds to the $1$-to-$2$ case,
utilizing $\mathrm{CNOT}$ gates for entanglement between layers. For QNN
training in this instance, we limit it to two epochs with a learning rate of
0.05. The objective is to demonstrate that the resulting EQK is not highly
dependent on the training specifics of the corresponding QNN and even a sub-
optimal QNN training can lead to a powerful EQK for the specific task.

The relative improvements vary depending on the dataset under consideration
and range from $-40\%$ to $40\%$. As expected, the worst performance of the
combined protocol is observed at higher noise levels and for a greater number
of layers, corresponding to a longer-depth quantum circuit. It is crucial to
emphasize once more that these high noise strength values significantly
surpass the noise parameters of current quantum devices. Therefore, even
though the combined protocol necessitates running a double-depth circuit with
more qubits, considering the noise model described earlier, it remains
suitable for actual quantum devices.

##  Appendix F Datasets

Figure F.1:  From left to right: sinus, corners, spiral and circles datasets
generated and used in this work. These figures are generated with 2000 data
points.

In this work we utilized four artificial datasets for evaluating the
performance of the proposed embedded quantum kernels. Here we explain how they
were created:

  * •

Sinus Dataset: The sinus dataset is defined by a sinusoidal function,
specifically $f(x_{1})=-0.8\;\sin(\pi x_{1})$. Points located above this
sinusoidal curve are categorized as class -1, while points below it are
assigned to class +1.

  * •

Corners Dataset: The corners dataset comprises four quarters of a
circumference with a radius of 0.75, positioned at the corners of a square.
Points located inside these circular regions are labeled as class -1, while
points outside are classified as class +1.

  * •

Spiral Dataset: The spiral dataset features two spirals formed by points
arranged along a trajectory defined by polar coordinates. The first spiral,
denoted as class +1, originates at the origin (0, 0) and spirals outward in a
counter-clockwise direction, forming a curve. The second spiral, labeled class
-1, mirrors the first spiral but spirals inward in a clockwise direction.
These spirals are generated by varying the polar angle, selected randomly to
create the data points. The radial distance from the origin for each point
depends on the angle, creating the characteristic spiral shape. Noise is added
to the data points by introducing random perturbations to ensure they do not
align perfectly along the spirals.

  * •

Circles Dataset: The circles dataset is created using two concentric circles
that define an annular region. The inner circle has a radius of
$\sqrt{2/\pi}$, while the outer circle has a radius of $0.5\sqrt{2/\pi}$. Data
points located within the annular region are labeled as -1, while those
outside the region are labeled as +1.

##  Appendix G Numerical experiments

In the main we provide numerical results for two datasets for the $n$-to-$n$
and $1$-to-$n$ constructions up to 8 qubits. Here we provide the tables of
more numerical experiments with the same hyperparameters: learning rate of
0.05 and 30 epochs for the $n=1$ QNN and learning rate of 0.005 and 10 epochs
for $n>1$. Here we provide more results for different datasets and considering
also CZ as source of entanglement. All accuracies refer to test accuracies.

Dataset | QNN accuracy | EQK type | EQK accuracy
---|---|---|---
Sinus | 0.890 |  $3q$ $\mathrm{CNOT}$ | 0.960
Sinus |  |  $3q$ $\mathrm{CZ}$ | 0.948
Corners | 0.886 |  $3q$ $\mathrm{CNOT}$ | 0.954
Corners |  |  $3q$ $\mathrm{CZ}$ | 0.940
Spiral | 0.800 |  $3q$ $\mathrm{CNOT}$ | 0.994
Spiral |  |  $3q$ $\mathrm{CZ}$ | 0.866
Circles | 0.698 |  $3q$ $\mathrm{CNOT}$ | 0.866
Circles |  |  $3q$ $\mathrm{CZ}$ | 0.864
Table 1: Numerical results for the $1$-to-$3$ architecture, considering two
types of entangling gates ($\mathrm{CNOT}$ and $\mathrm{CZ}$), on four
distinct datasets.

In Table 1, the outcomes demonstrate striking similarity when introducing
entanglement through $\mathrm{CZ}$ gates or $\mathrm{CNOT}$ gates, with the
exception of the spiral dataset. Given the ambiguity surrounding the method of
entanglement introduction, one might contemplate utilizing controlled
rotations with random parameters as a potential source of entanglement.
Notably, it is observed that even starting from a single-qubit QNN achieving
accuracies of less than $90\%$, the construction of EQKs yields accuracies
surpassing $95\%$.

Figure G.1: Explicit quantum circuit for constructing the $1$-to-$3$ embedding quantum kernel with cascade of $CNOT$ gates as source of entanglement. The parameters $\bm{\theta}^{*}$ are the ones obtained during the training of the single-qubit QNN. Note that the parameters of the last layer are irrelevant. Dataset | n=2 | n=3 | n=4 | n=5 | n=6 | n=7 | n=8 | n=9 | n=10
---|---|---|---|---|---|---|---|---|---
Corners | 0.890 | 0.954 | 0.974 | 0.978 | 0.982 | 0.980 | 0.978 | 0.978 | 0.978
Circles | 0.832 | 0.866 | 0.950 | 0.970 | 0.974 | 0.970 | 0.978 | 0.966 | 0.962
Table 2: Numerical results for the $1$-to-$n$ architecture for up to $n=10$
qubits. These results up to 8 qubits are depicted in the main text.

In the $1$-to-$n$ construction presented in Table 2 for the corners and
circles datasets, we observe that the accuracy rises rapidly as qubits are
added, reaching a peak at $n=4$. Subsequently, the accuracy plateaus,
attaining a maximum at $n=6$ for the corners dataset and $n=8$ for the circles
dataset. After reaching these points, the accuracy gradually decreases with
additional qubits.

Dataset | n=1 | n=2 | n=3 | n=4 | n=5 | n=6 | n=7 | n=8
---|---|---|---|---|---|---|---|---
Corners (QNN) | 0.886 | 0.934 | 0.948 | 0.952 | 0.954 | 0.956 | 0.960 | 0.962
Corners (EQK) | 0.898 | 0.948 | 0.960 | 0.966 | 0.970 | 0.972 | 0.974 | 0.974
Circles (QNN) | 0.698 | 0.792 | 0.820 | 0.858 | 0.934 | 0.944 | 0.944 | 0.946
Circles (EQK) | 0.796 | 0.820 | 0.906 | 0.968 | 0.974 | 0.974 | 0.976 | 0.980
Table 3: Numerical results for the $n$-to-$n$ architecture for up to $n=L+1=8$
qubits. These results are depicted in the main text.

In the $n$-to-$n$ approach, as presented in Table 3, the accuracy consistently
increases with the addition of qubits for both the QNN and the EQK.
Additionally, as discussed in the main text, the EQK consistently outperforms
the corresponding QNN architecture for the same value of $n$.

Dataset | Layers | QNN accuracy | EQK accuracy
---|---|---|---
Sinus | $L=5$ | 0.948 | 0.970
Sinus | $L=6$ | 0.956 | 0.964
Sinus | $L=7$ | 0.966 | 0.972
Sinus | $L=8$ | 0.958 | 0.964
Corners | $L=5$ | 0.948 | 0.970
Corners | $L=6$ | 0.936 | 0.950
Corners | $L=7$ | 0.934 | 0.948
Corners | $L=8$ | 0.916 | 0.920
Spiral | $L=5$ | 0.952 | 0.996
Spiral | $L=6$ | 0.974 | 0.998
Spiral | $L=7$ | 0.978 | 1.000
Spiral | $L=8$ | 0.980 | 0.998
Circles | $L=5$ | 0.786 | 0.812
Circles | $L=6$ | 0.808 | 0.814
Circles | $L=7$ | 0.792 | 0.820
Circles | $L=8$ | 0.844 | 0.902
Table 4: Numerical results for the $2$-to-$2$ architecture for different
number of layers. For the experiments in the main text we choose $L=7$.

Finally, in Table 4, we present results for the four datasets considering
different numbers of layers for the $n$-to-$n$ architecture with $n=2$. Adding
layers increases the expressivity of the QNN but does not guarantee better
accuracies, as we can observe. Again, we see that the EQK consistently
outperforms the QNN.


# Three–dimensional Gaussian fluctuations of spectra of overlapping stochastic
Wishart matrices

Jeffrey Kuan and Zhengye Zhou

###### Abstract

In [DP18], the authors consider eigenvalues of overlapping Wishart matrices
and prove that its fluctuations asymptotically convergence to the Gaussian
free field. In this brief note, their result is extended to show that when the
matrix entries undergo stochastic evolution, the fluctuations asymptotically
converge to a three-dimensional Gaussian field, which has an explicit contour
integral formula. This is analogous to the result of [Bor14] for stochastic
Wigner matrices.

##  1 Introduction and Main Result

###  1.1 Motivation

The two–dimensional Gaussian free field appears as the asymptotic fluctuations
of the height functions of several two–dimensional models, such as eigenvalues
of sub–matrices of Wigner random matrices and interacting particle systems on
a two–dimensional lattice. For example, in [BF14], the authors prove that an
interacting particle system on $\mathbb{Z}\times\mathbb{Z}_{+}$ has asymptotic
fluctuations described by the Gaussian free field on the upper half–plane
$\mathbb{H}$. This interacting particle system can be constructed through the
representation theory of the infinite–dimensional unitary group $U(\infty)$
[BK08]. Similarly, in [Bor14] it is proven that for Wigner matrices whose
entries undergo stochastic evolution, the fluctuations of the eigenvalues of
sub–matrices are described by a three–dimensional Gaussian field which
generalizes the Gaussian free field on $\mathbb{H}$. This same
three–dimensional field occurs in a non–commutative random walk [Kua16].

In [BK10], the representation theory of the infinite–dimensional orthogonal
group $O(\infty)$ is used to construct an interacting particle system, this
time on the quadrant $\mathbb{Z}_{\geq 0}\times\mathbb{Z}_{+}$. The vertical
boundary $\\{0\\}\times\mathbb{Z}_{+}$ acts as a reflecting barrier for the
particles. In [Kua14], it is proven that the fixed–time asymptotic
fluctuations of the interacting particle system are the Gaussian free field on
the upper half–plane with the unit disk removed, $\mathbb{H}-\mathbb{D}$.
Upcoming work by the second author generalizes this result to multiple time
parameters, wherein a three–dimensional Gaussian field appears in the
asymptotic limit. Thus, it is natural to ask if the same three–dimensional
Gaussian field occurs in the context of random matrices. An appropriate analog
is Wishart matrices, because the eigenvalues are always non–negative,
corresponding to the interacting particles being restricted to be to the right
of the reflecting barrier. In [DP18], it is shown that the eigenvalues of
overlapping Wishart matrices have fluctuations that converge to the Gaussian
free field. In the present work, their latter result is generalized to the
case when the entries of the Wishart matrices undergo stochastic evolution.

###  1.2 Main Result

Consider an infinite array $\\{Z_{ij}(t):i,j\geq 1\\}$ of stochastic processes
(not necessarily Markov). The $Z_{ij}(t)$ are allowed to be real, complex, or
quaternionic, corresponding to the usual $\beta=1,2,4$ values. In all cases,
assume that all real and imaginary components of $Z_{ij}(t)$ are mutually
independent. Assume the following moment conditions for all fixed times $t\geq
0$:

$\mathbb{E}Z_{ij}(t)=0,\quad\mathbb{E}|Z_{ij}(t)|^{2}=1,\quad\mathbb{E}|Z_{ij}(t)|^{4}=1+\frac{2}{\beta},\quad\sup_{i,j\geq
1}\mathbb{E}|Z_{ij}(t)|^{4+\epsilon}<\infty.$

The first, second, and fourth moments agree with the standard real, complex
and quaternion Gaussian random variables for all $t\geq 0$. Further suppose
that for all $s,t\geq 0$

$\mathbb{E}\left|Z_{ij}(s)\right|^{2}\left|Z_{ij}(t)\right|^{2}=c_{2}(s,t),\quad\quad\mathbb{E}|Z_{ij}(s)Z_{ij}(t)]|=c_{1}(s,t).$

From this infinite array, consider rectangular submatrices $B_{i}(s)$ for
$1\leq i\leq k$ with sizes $m_{i}L\times n_{i}L$ that overlap at $m_{ij}L$
rows and $n_{ij}L$ columns. Define the Wishart matrices
$W_{i}(s)=B_{i}(s)^{*}B_{i}(s)/L$. The statistics of interest are

$X_{i}^{(p)}(s):=\mathrm{Tr}(W_{i}^{p}(s)))-\mathbb{E}\mathrm{Tr}(W_{i}^{p}(s))$

for positive integers $p$. The result of this paper is that as
$L\rightarrow\infty$, the random variables $X_{i}^{(p)}$ converge to a
Gaussian field with an explicit covariance.

######  Theorem 1.1.

Suppose that $m_{i},n_{j},m_{ij},n_{ij}$ all depend on $L$ in such a way that

$\frac{m_{i}}{L}\rightarrow\mu_{i},\quad\quad\frac{n_{j}}{L}\rightarrow\nu_{j},\quad\quad\frac{m_{ij}n_{ij}L^{2}}{m_{i}n_{i}m_{j}n_{j}}\rightarrow\theta_{ij}$

as $L\rightarrow\infty$. Then
$(X_{1}^{(p_{1})}(s_{1}),\ldots,X_{k}^{(p_{k}}(s_{k}))$ converges in
distribution to a centered Gaussian vector $(\xi_{1},\ldots,\xi_{k})$ with
covariance

$\mathbb{E}[\xi_{i}\xi_{j}]=\frac{4c_{1}(s_{i},s_{j})p_{i}p_{j}}{\beta(i\pi)^{2}}\oint\oint(\mu_{i}+\nu_{i}+2\Re\zeta_{i})^{p_{i}-1}(\mu_{j}+\nu_{j}+2\Re\zeta_{j})^{p_{i}-1}\log\left|\frac{\theta_{ij}^{-1}-{\zeta_{i}}\zeta_{j}}{\theta_{ij}^{-1}-{\zeta_{i}}\bar{\zeta_{j}}}\right|\frac{\Im\zeta_{i}}{\zeta_{i}}\frac{\Im\zeta_{j}}{\zeta_{j}}d\zeta_{i}d\zeta_{j}\\\
+\frac{4(c_{2}(s_{i},s_{j})-1-2\beta^{-1}c_{1}(s_{i},s_{j}))\theta_{ij}p_{i}p_{j}}{(\pi
i)^{2}}\oint\oint(\mu_{i}+\nu_{i}+2\Re\zeta_{i})^{p_{i}-1}(\mu_{j}+\nu_{j}+2\Re\zeta_{j})^{p_{i}-1}\frac{\left(\Im\zeta_{i}\right)^{2}}{\zeta_{i}}\frac{\left(\Im\zeta_{j}\right)^{2}}{\zeta_{j}}d\zeta_{i}d\zeta_{j}$

where the $\zeta_{i},\zeta_{j}$ contours are semicircles over the upper half
plane $\mathbb{H}$ with radii $\sqrt{\mu_{i}\nu_{i}},\sqrt{\mu_{j}\nu_{j}}.,$
respectively.

Remark. When $s_{i}=s_{j}$, this reduces to Lemma 2.2 of [DP18]. Also, by the
Cauchy–Schwarz inequality,
$c_{2}(s,t)\leq\sqrt{\mathbb{E}|Z_{ij}(s)|^{4}|Z_{ij}(t)|^{4}}=1+2/\beta$, so
the constant in the second line is always non–positive.

Acknowledgments. The authors would like to thank Alexei Borodin, Ioana
Dumitriu and Elliot Paquette for helpful discussions. The simulations done in
section 3 were aided by Texas A&M University High Performance Research
Computing.

##  2 Proof

The proof is in essence identical to [DP18]. As explained in [DP18], the
centered trace $X^{(p)}_{i}(t)$ can be written as a sum over closed walks of
length $2p$ on the complete bipartite graph $K(\mathbb{N},\mathbb{N}).$ More
explicitly, represent a walk as a function
$w:\\{1,\ldots,2p\\}\rightarrow\mathbb{N}$ and define

$Z_{w}(t)=\prod_{i=1}^{p}Z_{w(2p-1),w(2p)}(t)\overline{Z_{w(2p-1),w(2p)}(t)}-\mathbb{E}\left[\prod_{i=1}^{p}Z_{w(2p-1),w(2p)}(t)\overline{Z_{w(2p-1),w(2p)}(t)}\right].$

With this definition,

$X_{i}^{(p)}(t)=\sum_{w\in\mathcal{S}_{p}}Z_{w}(t),$

where $\mathcal{S}_{p}$ is the set of all closed walks of length $2p$ on the
complete bipartite graph $K(\mathbb{N},\mathbb{N})$.

Thus, when computing the covariance of $X_{i}^{(k)}$ and $X_{j}^{(l)}$, the
sum is over pairs of walks. As explained on page 8 of [DP18], the only pairs
of walks that contribute in the limit are those that cover exactly $k+l$
vertices and traverse a common edge. There are two types of pairs of walks
that can arise this way. The first type is that both walks are depth–first
search walks of two trees that cover a common edge. For example, the first
walk can be $1\rightarrow 2\rightarrow 1\rightarrow 3\rightarrow 1$ while the
second walk is $1\rightarrow 2\rightarrow 1\rightarrow 4\rightarrow 1$. The
contributions of these walks is then

$\frac{4(\mathbb{E}[Z_{11}^{2}(s)Z_{11}^{2}(t)]-1)k\cdot
l\cdot\theta_{ij}}{(\pi i)^{2}}\oint\oint\
(\mu_{1}+\nu_{1}+2\Re\zeta_{1})^{k-1}(\mu_{2}+\nu_{2}+2\Re\zeta_{2})^{l-1}\frac{(\Im\zeta_{1})^{2}}{\zeta_{1}}\frac{(\Im\zeta_{2})^{2}}{\zeta_{2}}d\zeta_{1}d\zeta_{2}.$

(compare to equation (8)) of [DP18]. The second type of pairs of walks are
when both walks cover a unicyclic graph, traversing the cycle once and making
excursions along trees attached to this cycle. The cycle must be common to
both walks so that it is traversed twice by the union of the walks. This type
of walk contributes

$-\frac{8\mathbb{E}|Z_{ij}(s)Z_{ij}(t)|k\cdot l\cdot\theta_{ij}}{\beta(\pi
i)^{2}}\oint\oint\
(\mu_{1}+\nu_{1}+2\Re\zeta_{1})^{k-1}(\mu_{2}+\nu_{2}+2\Re\zeta_{2})^{l-1}\frac{(\Im\zeta_{1})^{2}}{\zeta_{1}}\frac{(\Im\zeta_{2})^{2}}{\zeta_{2}}d\zeta_{1}d\zeta_{2}\\\
\frac{4\mathbb{E}[Z_{ij}(s)Z_{ij}(t)]kl}{\beta(\pi
i)^{2}}\oint\oint(\mu_{1}+\nu_{1}+2\Re\zeta_{1})^{k-1}(\mu_{2}+\nu_{2}+2\Re\zeta_{2})^{l-1}\log\left|\frac{\theta_{ij}^{-1}-{\zeta_{1}}\zeta_{2}}{\theta_{ij}^{-1}-{\zeta_{1}}\bar{\zeta_{2}}}\right|\frac{\Im\zeta_{1}}{\zeta_{1}}\frac{\Im\zeta_{2}}{\zeta_{2}}d\zeta_{1}d\zeta_{2}$

Combining the two sums yields the expression in the statement of the Theorem.

Finally, the asymptotic fluctuations are asymptotically Gaussian because of
the same arguments as in section 3 of [DP18].

##  3 Examples

Integrating by parts twice yields, for the main term,

$\frac{c_{1}(s_{i},s_{j})\theta_{ij}^{-1}}{2\beta\pi^{2}}\oint\oint\frac{(\mu_{i}+\nu_{i}+2\Re(\zeta_{i}))^{p_{i}}(\mu_{j}+\nu_{j}+2\Re(\zeta_{j}))^{p_{j}}}{(\theta_{ij}^{-1}-\zeta_{i}\zeta_{j})^{2}}d\zeta_{i}d\zeta_{j}$

where the $\zeta_{i},\zeta_{j}$ contours are now full circles in $\mathbb{C}$
with radii $\sqrt{\mu_{i}\nu_{i}},\sqrt{\mu_{j}\nu_{j}},$ respectively. The
error term can be evaluated directly with residues without integrating by
parts. For example, if $p_{i}=p_{j}=1$ then, using the Taylor series
$(\theta_{ij}^{-1}-\zeta_{i}\zeta_{j})^{-2}=\theta_{ij}^{2}(1+2\theta_{ij}\zeta_{i}\zeta_{j}+\cdots),$
the covariance equals

$2c_{1}(s_{i},s_{j})\beta^{-1}\theta_{ij}\mu_{i}\nu_{i}\mu_{j}\nu_{j}+\theta_{ij}\mu_{i}\mu_{j}\nu_{i}\nu_{j}(c_{2}(s_{i},s_{j})-1-2\beta^{-1}c_{1}(s_{i},s_{j})).$

This can be checked with Python: for example, if $Z_{11}(t)$ is a real
Ornstein–Uhlenbeck process with $\mathbb{E}[Z_{11}(s_{i})Z_{11}(s_{j})]=1/2$,
then $\mathbb{E}[Z_{11}(s_{i})^{2}Z_{11}(s_{j})^{2}]=1+2\cdot(1/2)^{2}=3/2$,
and if $\mu_{i}=4,\nu_{i}=2,\mu_{j}=\nu_{j}=1=m_{ij}=n_{ij}$, then
$\theta_{ij}=1/8$ and the prediction for the covariance is $2\cdot 1/2\cdot
1/8\cdot 8+1/8\cdot 8\cdot(3/2-1-2\cdot 1/2)=1/2$, and a simulation of 16,800
samples of Wishart matrices at $L=400$ yields $\approx 0.51310$ for the
covariance. Higher moments are difficult to simulate, because the Gaussian
free field has logarithmic fluctuations, requiring essentially exponential
computation time. However, with the aid of Texas A&M High Performance Research
Computing, 48,000 samples at size $L=400$ predicts a covariance of $\approx
181$ when $p_{i}=p_{j}=2$ and
$(\mu_{i},\nu_{i},\mu_{j},\nu_{j},m_{ij},n_{ij})=(5,2,3,1,1,1),$ and
$Z_{11}(t)$ as above. The formula yields

$2c_{1}(s_{i},s_{j})\theta_{ij}\mu_{i}\nu_{i}\mu_{j}\nu_{j}\left(4(\mu_{i}+\nu_{i})(\mu_{j}+\nu_{j})+2\theta_{ij}\mu_{i}\nu_{i}\mu_{j}\nu_{j}\right)\\\
+(c_{2}(s_{i},s_{j})-1-2\beta^{-1}c_{1}(s_{i},s_{j}))\theta_{ij}\mu_{i}\nu_{i}\mu_{j}\nu_{j}\cdot
4(\mu_{i}+\nu_{i})(\mu_{j}+\nu_{j}),$

which numerically evaluates to $354-168=186$. Both simulations are within
$3\%$ of the correct answer.

The Python code is available from the first author upon request.

## References

  * [BF14] Alexei Borodin and Patrik L. Ferrari.  Anisotropic growth of random surfaces in $2+1$ dimensions.  Comm. Math. Phys., 325(2):603–684, 2014.
  * [BK08] Alexei Borodin and Jeffrey Kuan.  Asymptotics of Plancherel measures for the infinite-dimensional unitary group.  Adv. Math., 219(3):894–931, 2008.
  * [BK10] Alexei Borodin and Jeffrey Kuan.  Random surface growth with a wall and Plancherel measures for $\rm O(\infty)$.  Comm. Pure Appl. Math., 63(7):831–894, 2010.
  * [Bor14] Alexei Borodin.  CLT for spectra of submatrices of Wigner random matrices.  Moscow Mathematical Journal, 14(1), March 2014.
  * [DP18] Ioana Dumitriu and Elliot Paquette.  Spectra of overlapping wishart matrices and the gaussian free field.  Random Matrices: Theory and Applications, 07(02):1850003, 2018.
  * [Kua14] Jeffrey Kuan.  The Gaussian free field in interlacing particle systems.  Electron. J. Probab., 19:no. 72, 31, 2014.
  * [Kua16] Jeffrey Kuan.  Three–dimensional Gaussian fluctuations of non–commutative random surfaces along time–like paths.  Advances in Mathematics, 303(5):716–744, 2016.