# Quantum Reservoir Computing and Risk Bounds

Naomi Mona Chmielewski, Nina Amini, Joseph Mikael

► To cite this version:

Naomi Mona Chmielewski, Nina Amini, Joseph Mikael. Quantum Reservoir Computing and Risk Bounds. 2025. hal-04884902

HAL Id: hal-04884902

https://hal.science/hal-04884902v1

Preprint submitted on 14 Jan 2025

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# Quantum Reservoir Computing and Risk Bounds

**Naomi Mona Chmielewski**

L2S, CNRS-CentraleSupélec-
Université Paris-Sud
Université Paris-Saclay
3 Rue Joliot Curie
91190 Gif-sur-Yvette, France;
EDF Lab
Palaiseau, France

NAOMI.CHMIELEWSKI@EDF.FR

**Nina Amini**

L2S, CNRS-CentraleSupélec-
Université Paris-Sud
Université Paris-Saclay
3 Rue Joliot Curie, 91190 Gif-sur-Yvette, France

NINA.AMINI@CENTRALESUPELEC.FR

**Joseph Mikael**

EDF Lab
Palaiseau, France

JOSEPH.MIKAEL@EDF.FR

## Abstract

We propose a way to bound the generalisation errors of several classes of quantum reservoirs using the Rademacher complexity. We give specific, parameter-dependent bounds for two particular quantum reservoir classes. We analyse how the generalisation bounds scale with growing numbers of qubits. Applying our results to classes with polynomial readout functions, we find that the risk bounds converge in the number of training samples. The explicit dependence on the quantum reservoir and readout parameters in our bounds can be used to control the generalisation error to a certain extent. It should be noted that the bounds scale exponentially with the number of qubits $n$. The upper bounds on the Rademacher complexity can be applied to other reservoir classes that fulfill a few hypotheses on the quantum dynamics and the readout function.

**Keywords:** Quantum Machine Learning, Quantum Reservoir Computing, Rademacher Complexity, Risk Bounds, Generalisation Error

## 1 Introduction

As an alternative to classical machine learning models, which tend to require a high number of parameters and consequently large data sets, quantum machine learning (QML) models for near-term quantum devices present a promising lead for understanding the potential of quantum technologies. Many current near-term QML methods are based on some form of variational quantum circuits and can be used for different learning tasks, such as classification (Schuld et al. (2020)), autoencoding (Romero et al. (2017)) and generative models (Romero et al. (2019)). Phenomena like barren plateaus (McClean et al. (2018)) and general poor trainability (Anschuetz and Kiani (2022)) of variational quantum methods have moti-
vated the search for circuit designs and quantum machine learning models that circumvent
these issues.

*Reservoir computing (RC)* has recently garnered some attention as an energy and data efficient alternative to more mainstream variants of recurrent neural networks (RNN) for time series forecasting problems (Tanaka et al. (2019); Gilpin (2023)). Instead of a large neural network with trainable weights, RC uses a fixed-weight reservoir described by its input-dependent dynamics and a final readout layer, the latter of which is usually a simple model such as a linear or polynomial regression. The RC models that we are concerned with in this paper possess a memory of past inputs and are able to adapt to the dynamics of the input process. While initially conceived as a type of neural network (Jaeger (2001)), physical implementations of RC such as the use of silicone octopus arms have become a novel way to introduce different reservoir dynamics (Nakajima (2020)).

Two important properties that we typically require of a class of RC models are the so called *echo state property (ESP)* and the *fading memory property (FMP)*, both of which are characterisations of how the reservoir's "memory" of past inputs should behave over time. The ESP ensures that the reservoir forgets its initial conditions in the infinite time limit. On the other hand, a reservoir with the FMP will yield similar outputs for two input processes that behaved differently in the distant past but similarly in the recent past.

Recent propositions for *quantum reservoir computing (QRC)* propose to leverage the complex dynamics of quantum systems to reproduce the behaviour of complex input data dynamics. Since the only trainable part of a quantum reservoir is in the classical linear regression, the trainability issues of other QML models mentioned earlier are no longer of concern. As pointed out by Chen and Nurdin (2019), for the quantum system to have ESP, it must be dissipative (see e.g. Fujii and Nakajima (2016)). As an example, Suzuki et al. (2022) make use of the so called depolarising error.

One of the properties that is often studied in machine learning is a model's ability to generalise on unseen data. We call the *generalisation error* or *risk* of a hypothesis class the expected value of the loss over all possible samples according to some data distribution. Since this distribution is generally unknown, it is bounded through the maximal difference within the class of hypothesis functions between the generalisation error and the empirical error, i.e. the average of the model's errors evaluated on a finite number of training samples. Essentially, an upper bound on this quantity can give an idea of how different the loss on a random data sample will be, compared to the empirical risk. This upper bound is often called a risk bound. Ideally, the risk bound should decrease as the number of samples in the empirical risk increases. The risk bound can also give some information on how to scale the parameters of the learning model in order to improve its ability to generalise. A common way to bound the risk of a hypothesis class is through the *Rademacher complexity*, which is essentially a measure of the “richness” of the class by quantifying how well the class can adapt to random noise. The Rademacher complexity is defined on a finite number of samples, and for a hypothesis class to exhibit good generalisation, we expect the Rademacher complexity to go down as the number of samples goes up.

If the hypothesis class is universal, the difference between the in-class minimum and the minimum over all measurable functions can be made arbitrarily small. Formally, a universal class of RC is defined to be a class such that for any process characterised by a dynamical map with fading memory, there exists a member of the reservoir class that
can approximate the fading memory map arbitrarily well. Some popular reservoir classes have been proven to be universal by Grigoryeva and Ortega (2018a,b). Universal reservoir classes based on quantum dynamics have been analysed in Chen and Nurdin (2019); Chen et al. (2020); Sannia et al. (2024). Monzania and Pratia (2024) give universality conditions for reservoir classes based on quantum dynamics. For our analysis, we restrict ourselves to quantum reservoir classes with a fixed number of qubits, similarly to Bu et al. (2022), which is a more realistic setting. This is in contrast to the standard procedure when analysing universality properties of quantum reservoir classes where the number of qubits is generally unbounded.

Gonon et al. (2020) propose a certain number of hypotheses which allow them to establish risk bounds for reservoir classes equipped with a general readout class. In this paper, we propose to bound the Rademacher complexity of a general quantum reservoir class with finite parameter space equipped with different readout classes. We show that, under certain constraints, the dynamics of the quantum reservoir, which are governed by completely positive trace preserving (CPTP) maps, verify the aforementioned hypotheses.

We establish risk bounds for different quantum reservoirs and provide an analysis of the scaling of the risk bound in the number of training samples. We show how the choice of the readout class can impact the scaling of the risk bound in the number of qubits. As stated above, the specific classes of quantum reservoirs that we study employ polynomial readouts to prove universality. We show that this choice leads to an unfavourable scaling in the number of qubits. This problem can be somewhat alleviated by choosing either a linear readout function, or a polynomial readout of fixed (small) degree, in a fixed number of variables. Another method to introduce the structure of a polynomial algebra that has been suggested by Sannia et al. (2024) and Monzania and Pratia (2024) is spatial multiplexing. We compare the scaling of a polynomial readout and spatial multiplexing.

The generalisation abilities of certain parameterised quantum circuit configurations have been studied for example by Bu et al. (2022) and Caro et al. (2022), using the Rademacher complexity and covering numbers respectively. The learning problem in these cases is based on independent and identically distributed samples and is thus not immediately applicable to reservoir computing. Gonon and Jacquier (2023) establish approximation error bounds for trainable variational quantum circuits as well as quantum extreme learning machines (a form of reservoir computing that does not take into account past input states) using the Fourier transform. To our knowledge, risk bounds for the specific case of reservoir computing with quantum dynamics for time series forecasting have not yet been studied.

This paper is organised as follows: In Section 2 we introduce the mathematical background and recall some definitions. In Section 3 we establish some general conditions on the quantum reservoir classes as well as the input and target processes. In Section 4 we establish bounds on the Rademacher complexity for any quantum reservoir that verifies the conditions established in Section 3 when equipped with a polynomial readout class, as well as for linear and spatial multiplexing readout classes. In Section 5 we define slightly modified versions of the QRC classes introduced by Chen and Nurdin (2019) and Chen et al. (2020) and show that they verify the hypotheses from Section 3. Finally, in Section 6 we establish risk bounds for the reservoir classes from Section 5 that explicitly depend on the reservoir and readout parameters.
## 2 Framework

### 2.1 General definitions

We first define some of the notation that will be used in the remainder of the paper, namely the reservoir system and its filter and functional, the input-output processes, the echo state and fading memory property, which are two important properties of a reservoir system, and some norms that we use in the remainder of the document.

**Hilbert Space.** For the Hilbert space $\mathcal{H} = \mathbb{C}^{2n}$, $n \in \mathbb{N}$, define the set $\mathcal{B}(\mathcal{H})$ (which we abbreviate as $\mathcal{B}$ in the following) of bounded self-adjoint trace-class linear operators that act on $\mathcal{H}$, as well as the dual space $\mathcal{B}^*(\mathcal{H})$ of all possible linear transformations $T(\mathcal{B}) : \mathcal{B} \to \mathcal{B}$ which we abbreviate as $\mathcal{B}^*$ and which forms a finite Banach space with induced norm $\|T\|_1 = \sup_{A \in \mathcal{B}, A \neq 0} \frac{\|T(A)\|_1}{\|A\|_1}$.

**Density Matrices.** Define the set of quantum density matrices

$\bar{\mathcal{B}}^+(\mathcal{H}) := \{\rho \in \mathcal{B} : \rho \ge 0, \text{tr}[\rho] = 1\}$, as well as the hyperplane of traceless operators $\mathcal{B}_0(\mathcal{H}) := \{A \in \mathcal{B} : \text{tr}[A] = 0\}$ which we abbreviate as $\bar{\mathcal{B}}^+$ and $\mathcal{B}_0$ respectively. In the following we may also call elements in $\bar{\mathcal{B}}^+$ quantum states.

**Input and Output Sequences.** We now define the classical input and output sequences as well as the quantum reservoir states. Let $\mathbb{Z}_- = \{\dots, -2, -1, 0\}$ be the set of all negative integers. We consider the semi-infinite input sequences

$\mathbf{v} = (\dots, v_{-1}, v_0) \in (D_v)^{\mathbb{Z}_-} \subset \mathbb{R}^{\mathbb{Z}_-}$, where $D_v$ is a compact subset of $\mathbb{R}$, and output sequences

$\mathbf{y} = (\dots, y_{-1}, y_0) \in \mathbb{R}^{\mathbb{Z}_-}$ as well as reservoir states $\rho = (\dots, \rho_{-1}, \rho_0) \in (\bar{\mathcal{B}}^+)^{\mathbb{Z}_-}$.

**Input and Output Processes.** To model the stochastic processes of the input and output sequences, we consider random variable inputs $\mathbf{V} = (V_t)_{t \in \mathbb{Z}_-}$ and outputs $\mathbf{Y} = (Y_t)_{t \in \mathbb{Z}_-}$ with a causal Bernoulli shift structure, i.e. for $q_v, q_y \in \mathbb{N}^*$ and measurable functions $G^v : (\mathbb{R})^{q_v} \to D_v$ and $G^y : (\mathbb{R})^{q_y} \to \mathbb{R}$ of independent and identically distributed random variables $\xi = (\xi_t^v, \xi_t^y)_{t \in \mathbb{Z}_-}$ with values in $\mathbb{R}^{q_v} \times \mathbb{R}^{q_y}$ such that for any $t \in \mathbb{Z}_-$ we have

$$ 
\begin{aligned}
V_t &= G^v(\dots, \xi_{t-1}^v, \xi_t^v) \\
Y_t &= G^y(\dots, \xi_{t-1}^y, \xi_t^y).
\end{aligned} 
\quad (1) 
$$

See (Dedecker et al., 2007, Section 3.1) for more details on Bernoulli shifts.

**Quantum Reservoir System.** The reservoir system is defined by the system of equations

$$ 
\begin{cases} 
\rho_t = T(v_t, \rho_{t-1}) \\ 
y_t = h(\rho_t) 
\end{cases} 
\quad (2) 
$$

where $h: \mathcal{B} \to \mathbb{R}$ is called the readout function and $T: D_v \times \mathcal{B} \to \mathcal{B}$ is the completely positive trace preserving (CPTP) map that describes the dynamics of the quantum reservoir. A CPTP map is a completely positive linear mapping that preserves the trace of any operator.
Conceptually, the reservoir at time $t$ is described by its reservoir state $\rho_t$ and the subsequent input $v_{t+1}$ is injected into the dynamics of the reservoir, which is left to evolve for some unit of time. At time $t+1$, the new reservoir state $\rho_{t+1}$ encodes all of the previously seen inputs, transformed through the reservoir dynamics (see Fig. 1).

<figure><img src="image_2.png" /><figcaption>Figure 1: Schematic of the temporal evolution of a quantum reservoir. The purple dots designate qubits, the outline determines the reservoir. At time $t$, the reservoir is in state $\rho_t$. A new input $v_{t+1}$ is injected and the full reservoir is made to evolve according to the dynamics induced by the CPTP map $T$. After all the data has been injected in this way, the final reservoir state $\rho_0$ is processed in a readout function $h$ which produces the prediction.</figcaption></figure>

**Echo State Property.** A reservoir system satisfies the *echo state property (ESP)* if for any $v \in (D_v)^{\mathbb{Z}_-}$, (2) has a unique solution $\rho \in (\tilde{\mathbb{B}}^+)^{\mathbb{Z}_-}$ when restricted to the space $\tilde{\mathbb{B}}^+$ of density matrices. This translates the idea that the reservoir dynamics are independent of its initial conditions.

**Reservoir Filter and Functional.** In the case where the system satisfies the ESP, it induces a causal and time invariant reservoir filter $M_h^T$ as shown by Grigoryeva and Ortega (2018a) which allows us to describe the input transformation of all past inputs, as opposed to the CPTP map which only describes the transformation between time steps $t-1$ and $t$:

$$M_h^T: (D_v)^{\mathbb{Z}_-} \to \mathbb{R}^{\mathbb{Z}_-} \\ \mathbf{v} \mapsto h \circ M^T(\mathbf{v})$$

where $M^T(\mathbf{v}) := (\dots, T(v_{-1}, \rho_{-2}), T(v_0, \rho_{-1}))$, with the corresponding functional

$$H_h^T: (D_v)^{\mathbb{Z}_-} \to \mathbb{R} \\ \mathbf{v} \mapsto h \circ H^T(\mathbf{v})$$

where we define $H^T(\mathbf{v}) := M^T(\mathbf{v})_0$, and the subscript 0 indicates the last element in the semi-infinite sequence, i.e. $M^T(\mathbf{v})_0 = T(v_0, \rho_{-1}) = \rho_0$. The subscript $h$ and superscript $T$ are meant to highlight the dependence of the filter and functional on the CPTP map $T$ and the readout function $h$. In the following we will omit this subscript and superscript of the reservoir functional whenever the dependence on $T$ and $h$ is clear from the definitions.

**Fading Memory.** Let $\mathbf{w} = \{w_t\}_{t \ge 0}$ be a decreasing positive-valued sequence such that $\lim_{t \to \infty} w_t = 0$. Then we call the norm $\|\cdot\|_{1,\mathbf{w}}$ defined as $\|\mathbf{v}\|_{1,\mathbf{w}} := \sum_{t \in \mathbb{Z}_-} |v_t|w_{-t}$ for $\mathbf{v} \in (D_v)^{\mathbb{Z}_-}$ the $(1, \mathbf{w})$-norm and define the space $\ell_-^{1,w} := \{\mathbf{v} \in \mathbb{R}^{\mathbb{Z}_-} | \|\mathbf{v}\|_{1,\mathbf{w}} < \infty\}$. We say that
a reservoir map has **w-fading memory** if for any $v, s \in (D_v)^{\mathbb{Z}_-}$, for all $\varepsilon > 0$ there exists $\delta > 0$ such that $|M_h^T(v)_0 - M_h^T(s)_0| < \varepsilon$ whenever $\sup_{t \in \mathbb{Z}_-} |v_t - s_t| w_{-t} < \delta$ (alternatively, we say that a reservoir map has the fading memory property (FMP) if there exists a decreasing positive-valued weighting sequence **w** such that the reservoir map has **w**-fading memory). This formalises the idea that a reservoir map should focus more on recent inputs and forget the distant past over time.

**Norms.** For an operator $\Phi$ we write $\|\Phi\|_2 = \sqrt{\text{tr}[\Phi^\dagger\Phi]}$ for the Schatten-2 norm. Note that for a matrix $A$, the Schatten-2 norm $\|A\|_2$ is the same as the Frobenius norm $\|A\|_F = \sqrt{\sum_{i,j} A_{ij}^2}$ as well as the $L^2$-norm of the vectorised matrix $\|\text{vec}(A)\|_2 = \sqrt{\sum_k (\text{vec}(A))_k^2}$. For the remainder of the text, we will use the notation $\|\cdot\|_2$ to denote all three of these norms. This is justified by the fact that any CPTP map $T$ is a linear map, meaning that there exists a matrix $\hat{T}$ such that we can write the linear application as a matrix multiplication $T(A) = \hat{T}A$. We also denote by $\|\cdot\|_2$ the operator norm induced by the Schatten-2 norm, i.e. $\|T\|_2 := \sup_{A \neq 0} \frac{\|T(A)\|_2}{\|A\|_2} = \sup_{A \neq 0} \frac{\|T(A)\|_2}{\|A\|_2} = \sup_{A \in \mathcal{B}_0} \|T(A)\|_2$.

## 2.2 Empirical Risk and Rademacher Complexity

In this subsection we formally define the generalisation error, which measures the ability of the reservoir to generalise on unseen data, and then introduce the empirical risk. We also define the Rademacher complexity, the measure that we use to bound the generalisation error.

We consider the loss functions $\ell(\hat{y}, y) = f(\hat{y} - y)$ where $f: \mathbb{R} \to \mathbb{R}$ is $L_\ell$-Lipschitz-continuous and such that $f(0) = 0$. Once again, this condition is necessary for the proofs of the theorems in Section 6. For a hypothesis functional $H$, the generalisation error is given by

$$R(H) := \mathbb{E}_{(V,Y) \sim P}[\ell(H(V), Y_0)]$$

where $P$ is the joint distribution of $V$ and $Y$. As this quantity is generally intractable since $P$ is unknown, the standard procedure to bound the generalisation error is to uniformly bound the maximal difference between the true generalisation error and the empirical risk error for some hypothesis class. The empirical risk error $\hat{R}_m(H)$ of the hypothesis functional $H$ is typically defined on $m$ independent and identically distributed samples. In the context of time series, we do not usually have such i.i.d. samples; instead, a common way to define the empirical risk is to consider the $m$ subsequences $(V_{-m+1}, \dots, V_{-i-1}, V_i)$ for $i = 0$ to $m-1$ of the input process to construct $m$ samples. Note that the reservoir filter is defined on an infinite sequence, so that we must append the subsequence with zeros. We thus define $V_{-m+1:-i} := (\dots, 0, 0, V_{-m+1}, \dots, V_{-i-1}, V_i)$ and the empirical risk of hypothesis $H$

$$\hat{R}_m(H) := \frac{1}{m} \sum_{i=0}^{m-1} \ell(H(V_{-m+1:-i}), Y_{-i}).$$

The goal of this paper is to bound the quantity

$$\sup_{H \in \mathcal{H}} \left| R(H) - \hat{R}_m(H) \right| \qquad (3)$$
for specific classes $\mathcal{H}$ of quantum reservoirs with dependencies in the parameters of the CPTP and readout maps.

One common way to bound (3) is to find a bound that depends on the Rademacher complexity, and to subsequently bound the Rademacher complexity. In the case of time series, the standard definition of the Rademacher complexity needs to be adapted. Usually, the Rademacher complexity is defined on $k$ independent identically distributed samples. In our case, we only have one sample (the time series). We thus make use of the version defined in (Gonon et al., 2020, Section 4.1) which introduces ghost samples, which are “imagined copies” of the input process: For $k \in \mathbb{N}$, denote by $\tilde{\mathbf{V}}^{(0)}, \dots, \tilde{\mathbf{V}}^{(k-1)}$ independent copies of $\mathbf{V}$, then the Rademacher complexity of the hypothesis class $\mathcal{H}$ is defined as

$$ \mathcal{R}_k(\mathcal{H}) = \frac{1}{k} \mathbb{E} \left[ \sup_{H \in \mathcal{H}} \left\| \sum_{j=0}^{k-1} \varepsilon_j H \left( \tilde{\mathbf{V}}^{(j)} \right) \right\|_2 \right] $$

where $\varepsilon_0, \dots, \varepsilon_{k-1}$ are independent and identically distributed Rademacher random variables, independent of $\tilde{\mathbf{V}}^{(0)}, \dots, \tilde{\mathbf{V}}^{(k-1)}$.

## 3 General Assumptions on the Data Distribution and Quantum Reservoir Classes

In this section, we establish general conditions for the quantum reservoir classes that we wish to analyse. The conditions that we require of the reservoir map are notably contractivity in the space of quantum density matrices and Lipschitz-continuity in the input space, whereas the conditions on the readout map are mostly related to boundedness and Lipschitz-continuity. We also require the parameter space to be finite.

To be able to bound the generalisation error of the reservoir classes that we describe in Section 5, we first recall the assumptions on the data distribution introduced by Gonon et al. (2020):

(A$_1^{IO}$) For $I = y, v$ the functional $G^I$ as defined in (1) is $L_I$-Lipschitz continuous when restricted to $(\ell_{-}^{1, w^I}(\mathbb{R}^{q_I}), \|\cdot\|_{1, w^I})$ for some strictly decreasing weighting sequence $w^I$: $\mathbb{N} \to (0, 1]$ with finite mean, that is, $\sum_{j \in \mathbb{N}} j w_j^I < \infty$. More specifically, there exists $L_I > 0$ such that for all $\mathbf{u}^I = (u_t^I)_{t \in \mathbb{Z}_-} \in \ell_{-}^{1, w^I}(\mathbb{R}^{q_I})$ and $\mathbf{s}^I = (s_t^I)_{t \in \mathbb{Z}_-} \in \ell_{-}^{1, w^I}(\mathbb{R}^{q_I})$ it holds that

$$ \|G^I(\mathbf{u}^I) - G^I(\mathbf{s}^I)\|_2 \le L_I \|\mathbf{u}^I - \mathbf{s}^I\|_{1, w^I}. $$

(A$_2^{IO}$) Additionally, let the innovations in (1) satisfy $\mathbb{E}[\|\xi_0^I\|_2] < \infty$ for $I = y, v$.

(A$_3^{IO}$) We also assume that $D_{w^v} := \sup_{i \in \mathbb{N}} \frac{w_i^v}{w_i^1} < 1$ and $D_{w^v} := \sup_{i \in \mathbb{N}} \frac{w_i^v}{w_i^1} < 1$ and that there exists $M_\xi > 0$ such that for all $t \in \mathbb{Z}_-$, $\|\xi_t^I\|_2 < M_\xi$, where $I = v, y$.

We designate the set of conditions on the input-output processes as $A^{IO} = \{(A_1^{IO}), (A_2^{IO}), (A_3^{IO})\}$.
## 3.1 A General Class of Quantum Reservoirs

We begin by defining a general class of quantum reservoirs and readout maps.

First, we suppose that the parameter space $\Theta$ which parameterises the quantum reservoir maps is of finite cardinality:

(A$^{\Theta}$) The size of the parameter space $\Theta$ is finite, i.e. $|\Theta|:= \text{Card}(\Theta) < \infty$.

This condition can always be verified by restricting the parameter space appropriately. We then define, for a fixed number $n$ of qubits, the general class $\mathcal{F}_n^{\text{QRC}}(\Theta)$ of parameterised quantum reservoir maps on an $n$-qubit system with parameters in $\Theta$ as the set of CPTP maps $T$ that map a real-valued input $v$ as well as a self-adjoint trace class linear operator $A$ to another self-adjoint trace class linear operator, and that meets the two following conditions:

(A$_1^T$) We require all CPTP maps to be $L_R$-Lipschitz-continuous on the real-valued closed input interval $D_v = [a, b] \subset \mathbb{R}$, i.e. there exists a constant $0 < L_R$ such that for all $T^\theta \in \mathcal{F}_n^{\text{QRC}}(\Theta)$, for any fixed density operator $\rho \in \bar{\mathbb{B}}^+$ and any fixed parameter $\theta \in \Theta$, for all inputs $v_1, v_2 \in D_v$ we have $\|T^\theta(v_1, \rho) - T^\theta(v_2, \rho)\|_2 \le L_R|v_1 - v_2|$.

(A$_2^T$) We further require all CPTP maps to be $r$-contractions in the space of density operators, i.e. there exists a constant $0 < r < 1$ such that for all $T^\theta \in \mathcal{F}_n^{\text{QRC}}(\Theta)$, for any fixed $v \in D_v$ and for any fixed $\theta \in \Theta$, for all $\rho_1, \rho_2 \in \bar{\mathbb{B}}^+$, we have $\|T^\theta(v, \rho_1) - T^\theta(v, \rho_2)\|_2 \le r\|\rho_1 - \rho_2\|_2$. Note that this condition implies that the reservoir system has ESP and FMP according to (Martinez-Pena and Ortega, 2023, Proposition 1).

We designate the set of conditions on the CPTP maps as $A^T = \{(A_1^T), (A_2^T)\}$.

**Remark 1** In general, we also require the class $\mathcal{F}_n^{\text{QRC}}(\Theta)$ to be separable in the space of bounded continuous functions when equipped with the supremum norm so that we can conclude that the supremum over the reservoir class is measurable. In the case of quantum reservoirs, the classes that we consider are always subsets of the set of bounded operators that act on the finite-dimensional Hilbert space $\mathcal{H}$ with fixed number of qubits $n$, which itself is a vector space and becomes a finite-dimensional metric space when equipped with the supremum norm.

We denote by $\mathcal{F}_n^O$ a general class of readout functions $h: \bar{\mathbb{B}}^+ \cup 0 \to \mathbb{R}$ on an $n$-qubit system (where $0$ denotes the $2^n \times 2^n$ matrix with all zeros) that verifies the following conditions:

(A$_1^h$) $\mathcal{F}_n^O$ contains the zero function.

(A$_2^h$) $\mathcal{F}_n^O$ is separable in the space of bounded continuous functions when equipped with the supremum norm.

(A$_3^h$) The maps $h$ in $\mathcal{F}_n^O$ are $L_h$-Lipschitz continuous where there exists $\bar{L}_h > 0$ such that $L_h \le \bar{L}_h$ for all $h \in \mathcal{F}_n^O$.
$(A_4^h)$ There exists $L_{h,0} > 0$ such that $|h(0)| \le L_{h,0}$.

We designate the set of conditions on the readout maps as $A^h = \{(A_1^h), (A_2^h), (A_3^h), (A_4^h)\}$

For a CPTP map $T^\theta$, we define

$\prod_{i=0}^\infty T^\theta(v_{-i})\rho_{-\infty} := \lim_{N \to \infty} T^\theta(v_0, T^\theta(v_{-1}, \dots, T^\theta(v_{-N}, \rho_{-N})))$. Then we can write the infinite composition of reservoir maps $T^\theta$ as

$$H^{T^\theta}(\mathbf{v}) := \prod_{i=0}^{\infty} T^\theta(v_{-i})\rho_{-\infty}. \qquad (4)$$

We then define the general class of reservoir functionals as

$$\mathcal{H}_n^{\text{QRC}} := \left\{ H : (D_v)^{\mathbb{Z}_-} \to \mathbb{R} \mid H(\mathbf{v}) = h(H^{T^\theta}(\mathbf{v})), h \in \mathcal{F}_n^O, T^\theta \in \mathcal{F}_n^{\text{QRC}}(\Theta) \right\}. \quad (5)$$

Note that in general these expressions depend on the initial condition $\rho_{-\infty}$. In the following lemma we show that when condition $(A_2^T)$ is met, the expressions are no longer dependent on the initial condition:

**Lemma 2** Let $\mathcal{H}_n^{\text{QRC}}$ be the general class of reservoir functionals defined in (5). Then, all filters $H \in \mathcal{H}_n^{\text{QRC}}$ are independent of the initial state.

**Proof** The contractivity condition $(A_2^T)$ necessarily implies the condition on $T$ in (Chen and Nurdin, 2019, Theorem 3), namely that for all inputs $v \in D_v$, there exists $\epsilon \in (0, 1)$ such that the input-dependent CPTP map $T(v, \cdot)$ restricted to the hyperplane $\mathcal{B}_0$ of $2^n \times 2^n$ traceless Hermitian operators satisfies $\|T\|_{\mathcal{B}_0} \|2 \le 1 - \epsilon$. Lemma 2 is a straightforward consequence of (Chen and Nurdin, 2019, Theorem 3) and its following remark, which states that any quantum system that satisfies the aforementioned condition on $T$, applied an infinite number of times, maps any initial state $\rho_{-\infty}$ to the state $\prod_{i=0}^\infty T^\theta(v_{-i}) (\frac{I}{2^n})$. ■

We can thus set the initial reservoir state to be the completely mixed state, i.e. $\rho_{-\infty} = \frac{I}{2^n}$, which makes the reservoir maps in (4) independent of the initial condition. This means that our analyses are independent of the initial state, and we can establish general results.

# 4 Bounds on the Rademacher Complexity for Specific Readout Functions

We now move on to analyse the Rademacher complexity of quantum reservoirs with specific classes of readout functions. We first begin by defining the polynomial readout class since it has been used multiple times to prove universality of quantum reservoirs (see e.g. Chen and Nurdin (2019); Chen et al. (2020)). We find a bound on the Lipschitz constants of the members of the class, which allows us to establish a general bound on the Rademacher complexity of all reservoir classes which employ this readout class, and where the reservoir maps verify the hypotheses highlighted in Section 3. Additionally, we analyse the Rademacher complexity of two alternative readout classes, namely a linear readout and spatial multiplexing and compare them to the Rademacher complexity of the polynomial readout.
## 4.1 Polynomial Readout Map

The polynomial nature of the readout class used in the reservoir classes by Chen and Nurdin (2019) and Chen et al. (2020) is used to establish universality of the quantum reservoir classes. In fact, a common way to prove universality in the case of reservoirs is to evoke the Stone-Weierstrass theorem, which, amongst other things, requires the class to be an algebra. The polynomial readout class generates a polynomial algebra.

We define a class of multivariate polynomial readout functions with a slight modification to the one mentioned above, specifically with a maximal degree and a bound on the absolute value of the constant bias:

$$ \mathcal{F}_{n,R_{\max},C_{\max}}^{\text{poly}} = \left\{ h : \bar{\mathcal{B}}^{+} \cup 0 \to \mathbb{R} \middle| \begin{array}{l} h(\rho) = C + \sum_{d=1}^{R} \sum_{i_1=1}^{n} \sum_{i_2=i_1+1}^{n} \cdots \sum_{i_n=i_{n-1}+1}^{n} \sum_{r_{i_1}+\cdots+r_{i_n}=d} w_{i_1, \dots, i_n}^{r_{i_1}, \dots, r_{i_n}} \operatorname{tr}[\sigma_Z^{(i_1)} \rho]^{r_{i_1}} \cdots \operatorname{tr}[\sigma_Z^{(i_n)} \rho]^{r_{i_n}}, \\ 0 \le |C| \le C_{\max}, w_{i_1, \dots, i_n}^{r_{i_1}, \dots, r_{i_n}} \in [0, 1], 1 \le R \le R_{\max} \end{array} \right\} \quad (6) $$

where $\operatorname{tr}[\sigma_Z^{(i)} \rho]$ denotes the expectation of the projective measurement along the Z-axis on the i-th qubit; $C_{\max}$ is a bound on the maximal magnitude of the constant bias $C$; $R$, $R_{\max}$ and the $r_i$ are positive valued integers; $R_{\max}$ is a bound on the degree $R$ of the polynomials; and the $w_{i_1, \dots, i_n}^{r_{i_1}, \dots, r_{i_n}}$ are real valued weight parameters. The reason for this modification of the readout function is analogous to the reason for the restriction in the number of qubits.

Note that this is a class of polynomials with a constant bound on the maximal degree, meaning that it is a finite dimensional vector space so that it is separable and condition ($A_2^h$) in Section 3.1 is immediately verified. Additionally, it is clear that this class contains the zero function by setting all parameters to zero, so that condition ($A_1^h$) is also verified. We also have $|h(0)| = |C| \le C_{\max}$ for all $h \in \mathcal{F}_{n,R_{\max},C_{\max}}^{\text{poly}}$ meaning that ($A_4^h$) is verified. It remains to show condition ($A_3^h$), and to find an explicit expression for the Lipschitz constant, which we do in the following proposition:

**Proposition 3** The polynomial readout class $\mathcal{F}_{n,R_{\max},C_{\max}}^{\text{poly}}$ introduced in (6), with a fixed number $n \in \mathbb{N}$ of qubits, a maximal degree $R_{\max}$ of the polynomial readout, and a bound $C_{\max}$ on the maximal absolute value of the constant term in the readout verifies assumption ($A_3^h$), i.e. all maps $h \in \mathcal{F}_{n,R_{\max},C_{\max}}^{\text{poly}}$ are Lipschitz-continuous with bounded Lipschitz constants $L_h \le \bar{L}_h^{\text{poly}}$, where we write

$$ \bar{L}_h^{\text{poly}} := n\sqrt{2^n} R_{\max} \cdot \left( \binom{n+R_{\max}}{R_{\max}} - 1 \right) . \quad (7) $$

The proof is provided in Appendix A. This proposition is used in Section 6 to establish risk bounds for classes of quantum reservoir functionals that employ this readout function. In the following theorem we bound the Rademacher complexity of the general class of reservoir maps that fulfills the general hypotheses laid out in Section 3.1, equipped with the polynomial readout class (6):
**Theorem 4** Consider the class of quantum reservoirs

$\mathcal{H}_{n,R_{\max},C_{\max}}^{\text{poly}} := \{ H : (D_v)^{\mathbb{Z}_-} \to \mathbb{R} | H(\mathbf{v}) = h(H^{T^{\theta}}(\mathbf{v})), h \in \mathcal{F}_{n,R_{\max},C_{\max}}^{\text{poly}}, T^{\theta} \in \mathcal{F}_{n}^{\text{QRC}}(\Theta) \}$, i.e. the general class of reservoir maps equipped with the polynomial readout class introduced in (6), with a fixed number $n \in \mathbb{N}$ of qubits, a maximal degree $R_{\max}$ of the polynomial readout, and a bound $C_{\max}$ on the maximal absolute value of the constant term in the readout. Then the Rademacher complexity is bounded as

$$ \mathcal{R}_k(\mathcal{H}_{n,R_{\max},C_{\max}}^{\text{poly}}) \le \frac{1}{\sqrt{k}} \left( |\Theta| \cdot R_{\max} \cdot \left( \binom{n+R_{\max}}{R_{\max}} - 1 \right) + C_{\max} \right) \quad (8) $$

for any $k \ge 2$.

The proof is provided in Appendix B.

The size of the parameter space contributes linearly to the Rademacher complexity and can be chosen to fit the desired complexity.

## 4.2 Alternative Readout Functions

### 4.2.1 LINEAR AND FIXED-POLYNOMIAL READOUT

Note that the bound in (8) is unfavourable in both the number of qubits $n$ and the maximal degree of the polynomial $R_{\max}$. One way to alleviate this problem is to restrict the readout class to a low-degree polynomial, or even a linear readout, similar to Yasuda et al. (2023). This would reduce the upper bound on the Rademacher complexity by a factor $(\binom{n+R_{\max}}{R_{\max}} - 1)$:

**Corollary 5** Let $\mathcal{F}_{n,C_{\max}}^{\text{lin}} = \{h \in \mathcal{F}_{n,R_{\max},C_{\max}}^{\text{poly}} | R = R_{\max} = 1\}$ be the class of linear readout functions obtained by restricting the maximal polynomial degree in $\mathcal{F}_{n,C_{\max}}^{\text{lin}}$ to one. Consider the class of quantum reservoirs

$\mathcal{H}_{n,C_{\max}}^{\text{lin}} := \{ H : (D_v)^{\mathbb{Z}_-} \to \mathbb{R} | H(\mathbf{v}) = h(H^{T^{\theta}}(\mathbf{v})), h \in \mathcal{F}_{n,C_{\max}}^{\text{lin}}, T^{\theta} \in \mathcal{F}_{n}^{\text{QRC}}(\Theta) \}$, i.e. the general class of reservoir maps equipped with the linear readout class. Then the Rademacher complexity is bounded as

$$ \mathcal{R}_k(\mathcal{H}_{n,C_{\max}}^{\text{lin}}) \le \frac{1}{\sqrt{k}} (|\Theta| \cdot R_{\max} + C_{\max}) \quad (9) $$

for any $k \ge 2$.

**Proof** This bound is obtained by simply replacing the number $(\binom{n+R}{R} - 1)$ of sums in the proof of Theorem 4 by one. ■

While this significantly reduces the scaling with $n$, this class no longer induces any non-linearity in the measurement outcomes. An intermediary solution between the two previous readout classes that remains non-linear in the measurements but with a more favourable Rademacher complexity than that of (8) can be obtained through spatial multiplexing, which we explore in the next section.
### 4.2.2 SPATIAL MULTIPLEXING

Using spatial multiplexing to introduce a polynomial algebra as suggested by Monzania and Pratia (2024) can yield a similar Rademacher upper bound to the bound of the linear readout in (9), under the correct hypotheses: For a total number of qubits $n$, suppose we have a finite number $\ell$ of independent reservoir systems, each with $n/\ell$ qubits (where $\ell$ is such that $n/\ell$ is an integer). After injecting the input into each reservoir system, we consider the resulting reservoir state $\rho$ to be the joint state $\rho = \rho_1 \otimes \cdots \otimes \rho_\ell$ of each individual reservoir state $\rho_1, \dots, \rho_\ell$. Then the readout function of the complete state is obtained by multiplying linear readouts of the individual states, and adding a constant bias $C$, i.e.

$$h^{\text{SM}}(\rho) = \left( \sum_{i_1=1}^{n/\ell} w_{i_1} \text{tr} [\sigma_Z^{(i_1)} \rho] \right) \cdots \left( \sum_{i_\ell=(\ell-1)\frac{n}{\ell}+1}^{n} w_{i_\ell} \text{tr} [\sigma_Z^{(i_\ell)} \rho] \right) + C.$$

This induces the readout class

$$\begin{aligned}
\mathcal{F}_{n, \ell_{\max}, C_{\max}}^{\text{SM}} = & \{ h : \bar{\mathcal{B}}^{+} \cup 0 \to \mathbb{R} \mid \\
& h^{\text{SM}}(\rho) = C + \sum_{i_1=1}^{n/\ell} \sum_{i_2=n/\ell+1}^{2n/\ell} \cdots \sum_{i_\ell=(\ell-1)\frac{n}{\ell}+1}^{n} w_{i_1} w_{i_2} \cdots w_{i_\ell} \text{tr} [\sigma_Z^{(i_1)} \rho] \text{tr} [\sigma_Z^{(i_2)} \rho] \cdots \text{tr} [\sigma_Z^{(i_\ell)} \rho], \\
& 0 \le |C| \le C_{\max}, w_{i_j} \in [0, 1], 1 \le \ell \le \ell_{\max}, n/\ell \in \mathbb{N} \}
\end{aligned}
\tag{10}$$

where $\ell_{\max}$ denotes the maximal number of independent reservoir systems and $C_{\max}$ is an upper bound on the absolute value of the constant bias. Conditions $(A_1^h)$, $(A_2^h)$ and $(A_4^h)$ are verified in the same way as for the readout class $\mathcal{F}_{n, R_{\max}, C_{\max}}^{\text{poly}}$ discussed in Section 4.1. Condition $(A_3^h)$ is verified in the following proposition:

**Proposition 6** The spatial multiplexing readout class $\mathcal{F}_{n, \ell_{\max}, C_{\max}}^{\text{SM}}$ introduced in (10), with a fixed number $n \in \mathbb{N}$ of qubits, a maximal partition $\ell_{\max}$ of the $n$ qubits, and a bound $C_{\max}$ on the maximal absolute value of the constant term in the readout verifies assumption $(A_3^h)$, i.e. all maps $h \in \mathcal{F}_{n, \ell_{\max}, C_{\max}}^{\text{SM}}$ are Lipschitz-continuous with bounded Lipchitz constants $L_h \le \bar{L}_h^{\text{SM}}$, where we write

$$\bar{L}_h^{\text{SM}} := \ell_{\max}^{5/2} \sqrt{2^n}.$$

The proof is provided in Appendix C.

We can bound the Rademacher complexity in a similar way to Theorem 4:

**Corollary 7** Consider the class of quantum reservoirs

$$\mathcal{H}_{n, \ell_{\max}, C_{\max}}^{\text{SM}} := \left\{ H : (D_v)^{\mathbb{Z}_-} \to \mathbb{R} \mid H(v) = h(H^{T^\theta}(v)), h \in \mathcal{F}_{n, \ell_{\max}, C_{\max}}^{\text{SM}}, T^\theta \in \mathcal{F}_n^{\text{QRC}}(\Theta) \right\}, \\ \text{i.e. the general class of reservoir maps equipped with the spatial multiplexing readout class as defined in (10). Then the Rademacher complexity is bounded as}$$

$$\mathcal{R}_k(\mathcal{H}_{n, \ell_{\max}, C_{\max}}^{\text{SM}}) \le \frac{1}{\sqrt{k}} (|\Theta| \cdot n + C_{\max})$$
for any $k \ge 2$.

**Proof** The proof proceeds analogously to that of Theorem 4, only replacing the sums and indices, which leads to $\ell$ sums of $n/\ell$ indices. ■

We can interpret $\ell_{max}$ to be analogous to $R_{max}$ in the class $\mathcal{F}_{n,R_{max},C_{max}}^{\text{poly}}$, as $\ell_{max}$ limits the degree of the polynomial in the measurements of the reservoir state. Note that, compared to Theorem 4, there is no dependence in the maximal degree of the polynomial in the bound on the Rademacher complexity.

The scaling of this bound being significantly more favourable than that of the polynomial readout class, it seems that the spatial multiplexing readout is a more interesting readout class in terms of the Rademacher complexity. In general, the reservoir systems of the individual subsystems are governed by different dynamics, so that the remaining results of the paper are not directly applicable to the general spatial multiplexing case. In particular, the tensor product of Lipschitz-continuous maps is not necessarily Lipschitz-continuous, which means that an ensemble of reservoir classes that individually verify the hypotheses $A^T$ does not necessarily verify the hypotheses $A^T$ when taken as a tensor product. A more thorough analysis is required to establish general risk bounds when using spatial multiplexing. For the remainder of the document we thus consider the polynomial readout class $\mathcal{F}_{n,R_{max},C_{max}}^{\text{poly}}$.

**Remark 8** In practice, the condition ($A^\Theta$) is likely verified for quantum systems, as the accessible parameters are discretised through the limitations on the physical operations that are performed in order to apply the CPTP map. Depending on the chosen quantum system and the achievable precision of the physical operations, the size of the accessible parameter space might be very large, and the upper bound in Theorem 4 may become meaningless. A simple solution consists in restricting the parameter space artificially; the generalisation bounds in Section 6 can be used to choose an appropriately restricted parameter space.

# 5 Two Specialised Subclasses of Quantum Reservoirs

We now introduce two more specific classes of quantum reservoirs, namely the *Partial Trace Reservoir (PTR)* and the *Random Reinitialisation Reservoir (RRR)*, which are adaptations of the reservoir classes by Chen and Nurdin (2019) and Chen et al. (2020), respectively. We use the properties of these reservoirs to establish risk bounds in Section 6 that depend on their parameters as well as the parameters of the polynomial readout class $\mathcal{F}_{n,R_{max},C_{max}}^{\text{poly}}$ and the number of training samples.

## 5.1 Partial Trace Reservoir

We now define the PTR subclass of the reservoir map from Chen and Nurdin (2019), which slightly modifies the input injection and readout function from the first proposal of a QRC introduced by Fujii and Nakajima (2016). Modified versions of this class are frequently considered in quantum reservoir computing (see e.g. Yasuda et al. (2023), Sannia et al. (2024)). Conceptually, we inject the input $v$ into the first qubit in the mixed state $(v|0\rangle\langle0| + (1-v)|1\rangle\langle1|)$, then apply a global Hamiltonian on the full system during a time
$\tau$, and subsequently reset the first qubit in order to inject the next input; see Fig. 2 for a visual representation.

<figure><img src="image_3.png" /><figcaption>Figure 2: Schematic of the CPTP map of a partial trace reservoir. The purple dots designate the reservoir qubits, the solid outline determines the reservoir. The single pink dot designates the ancillary “input qubit”, and the dashed line designates the system on which we apply the unitary evolution determined by the XY-Hamiltonian. At time $t$, the reservoir is in state $\rho_t$. A new input $v_{t+1}$ is injected in the ancillary qubit by setting it to the mixed state $v_{t+1}|0\rangle\langle0| + (1 - v_{t+1})|1\rangle\langle1|$ and the reservoir qubits along with the ancillary qubit are made to evolve according to the unitary map induced by the Hamiltonian. After time $\tau$, the new reservoir state $\rho_{t+1}$ is obtained by tracing out the ancilla qubit. Finally, the next input is injected into the ancilla qubit.</figcaption></figure>

More concretely, define the input space

$$D_v^{\text{PTR}} = \left[ \frac{1}{2} \left( 1 - \sqrt{2} - 1 \right) + \epsilon_{\text{PTR}} \right] \left( \frac{1}{2} \left( 1 + \sqrt{2} - 1 \right) - \epsilon_{\text{PTR}} \right) \text{ for some fixed } 0 < \epsilon_{\text{PTR}} < \sqrt{2} - 1. \text{ We consider an } n\text{-qubit system and define the CPTP map}$$

$$T^{J,\gamma,\tau} : D_v^{\text{PTR}} \times \mathcal{B} \to \mathcal{B} \\ (v, A) \mapsto \text{tr}_0 \left[ e^{-iH(J,\gamma)\tau} A \otimes (v|0\rangle\langle0| + (1-v)|1\rangle\langle1|) e^{iH(J,\gamma)\tau} \right]$$

where $\text{tr}_0$ denotes the partial trace with respect to the first qubit, $\tau \in \mathbb{R}_+$ is a time parameter, and $H(J, \gamma)$ is an XY-Hamiltonian of the form

$$H(J, \gamma) = \sum_{i=0}^{n} \sum_{j=i+1}^{n} J^{i,j} \left( \sigma_X^{(i)} \sigma_X^{(j)} + \sigma_Y^{(i)} \sigma_Y^{(j)} \right) + \gamma \sum_{i=0}^{n} \sigma_Z^{(i)} \quad (11)$$

where $\sigma_X^{(i)}, \sigma_Y^{(i)}$ and $\sigma_Z^{(i)}$ are Pauli-X, -Y and -Z matrices applied to the $i$-th qubit (i.e. $\sigma_X^{(i)} = I \otimes \cdots \otimes I \otimes \sigma_X \otimes I \otimes \cdots \otimes I$ and so on), and the parameters $J = \{J^{ij}\}$ and $\gamma$ are real-valued constants.

Define $\Theta_{\text{PTR}}$ as the space of such reservoir parameters $J, \gamma$ as well as $\tau \in \mathbb{R}_+$. For a fixed number $n$ of qubits, we can then define the class of PTR quantum reservoir maps
$$ \mathcal{F}_n^{\text{PTR}} := \{ T^{J,\gamma,\tau} | (J,\gamma,\tau) \in \Theta_{\text{PTR}} \} . $$

In the following Lemma, we show that the class of reservoir maps $\mathcal{F}_n^{\text{PTR}}$ verifies Assumptions $A^T$. We use this Lemma in the next section to prove a generalisation bound for the PTR class of quantum reservoirs.

**Lemma 9** For any fixed parameters $(J, \gamma, \tau) \in \Theta_{\text{PTR}}$, the CPTP map

$T^{J,\gamma,\tau} : D_v^{\text{PTR}} \times \mathcal{B} \to \mathcal{B}$, $T^{J,\gamma,\tau}(v, A) \mapsto \text{tr}_0 [e^{-iH(J,\gamma)\tau} A \otimes (v|0\rangle\langle 0| + (1-v)|1\rangle\langle 1|) e^{iH(J,\gamma)\tau}]$ with Hamiltonian $H(J, \gamma)$ is strictly contractive in the space of density operators whenever $0 < \epsilon_{\text{PTR}} < \sqrt{2}-1$, with contractivity constant

$r_{\text{PTR}}(\epsilon_{\text{PTR}}) = 2\sqrt{2}(\epsilon_{\text{PTR}}^2 - \epsilon_{\text{PTR}}\sqrt{2}-1)+1$. Additionally, the map is Lipschitz-continuous in the space of inputs $D_v^{\text{PTR}}$ with Lipschitz-constant $L_R^{\text{PTR}} = 2$.

The proof is provided in Appendix D. Note that we can somewhat control the strength of the contractivity constant by choosing a convenient value for $\epsilon_{\text{PTR}}$.

In the following, we consider the parameter space $\Theta_{\text{PTR}}$ to be discrete and finite. We then define the class of reservoir functionals associated with the class of PTR reservoir maps as

$$ \mathcal{H}_{n,R_{\text{max}},C_{\text{max}}}^{\text{PTR}} := \{ H : (D_v)^{\mathbb{Z}_-} \to \mathbb{R} \mid H(v) = h(H^{T^{J,\gamma,\tau}}(v)), h \in \mathcal{F}_{n,R_{\text{max}},C_{\text{max}}}^{\text{poly}}, T^{J,\gamma,\tau} \in \mathcal{F}_n^{\text{PTR}} \}. \quad (12) $$

## 5.2 Random Reinitialisation Reservoir

We now define the RRR subclass of the reservoir map from Chen et al. (2020), which provides an interesting, more general alternative to the PTR class for quantum platforms on which the Hamiltonian (11) or the input injection in a mixed state are not available. The RRR map can be interpreted as a probabilistic injection of the input $v$, where for some parameter $\alpha \in (0, 1)$, with probability $(1-\alpha)v$ we apply the CPTP map $T_0$ to the reservoir, and with probability $(1-\alpha)(1-v)$ we apply the CPTP map $T_1$ (see Fig. 3 for a visual representation).

Mathematically, define the input space $D_v^{\text{RRR}} = [\epsilon_{\text{RRR}}, \frac{1}{2} - \epsilon_{\text{RRR}}]$ for some fixed $0 < \epsilon_{\text{RRR}} < \frac{1}{4}$. We consider an $n$-qubit system and define the CPTP map

$$ T^{\alpha, \sigma} : D_v^{\text{RRR}} \times \mathcal{B} \to \mathcal{B} \\ (v, A) \mapsto (1-\alpha)(vT_0 + (1-v)T_1)(A) + \alpha\sigma $$

where $\sigma$ is a fixed arbitrary quantum state, $\alpha$ is the parameter that determines the probability of reinitialising the reservoir state to $\sigma$ and $T_0 \neq T_1$ are fixed CPTP maps, contractive in the space of density matrices, with contraction constants $r_0$ and $r_1$ such that $r_0 + r_1 < 1$. If $r_0 = r_1$, we write $r_0 = r_1 = r_{0,1}$. If $r_0 \neq r_1$, we impose $r_0 > r_1$.

For a fixed number $n$ of qubits, we can then define the class of RRR quantum reservoir maps
<figure><img src="image_4.png" /><figcaption>Figure 3: Schematic of the CPTP map of a Random Reintialisation Reservoir. The purple dots designate the reservoir qubits, the solid outline determines the reservoir. At time t, the reservoir is in state $\rho_t$. With probability $v_{t+1}(1-\alpha)$ we apply the map $T_0$ to the reservoir to obtain the new reservoir state $\rho_{t+1} = T_0(\rho_t)$, with probability $(1-v_{t+1})(1-\alpha)$ we apply the map $T_1$ to the reservoir to obtain the new reservoir state $\rho_{t+1} = T_1(\rho_t)$, and with probability $\alpha$ we reintialise the reservoir to some fixed quantum state $\sigma$.</figcaption></figure>

$$ \mathcal{F}_n^{RRR} := \{T^{\alpha, \sigma} | (\alpha, \sigma) \in \Theta_{RRR}\} \qquad (13) $$

where $\Theta_{RRR} := [\alpha_{min}, \alpha_{max}] \times \Sigma$ is some compact subset of $(0, 1) \times \mathbb{B}^+$.

In the following Lemma we show that the class of reservoir maps $\mathcal{F}_n^{RRR}$ verifies Assumptions $A^T$. We use this Lemma in the next section to prove a generalisation bound for the RRR class of quantum reservoirs.

**Lemma 10** The class of RRR quantum reservoir maps $\mathcal{F}_n^{RRR}$ defined in (13) verifies assumptions $A^T$. More specifically, any CPTP map $T^{\alpha, \sigma} \in \mathcal{F}_n^{RRR}$ with arbitrary $n$ is strictly contractive in the space of density operators. If $r_0 = r_1 = r_{0,1}$, the contractivity constant is given by $r_{RRR} = (1 - \alpha_{min})r_{0,1}$. If $r_0 > r_1$, the contractivity constant is given by $r_{RRR} = 1 - \epsilon_{RRR}$. Furthermore, any CPTP map $T^{\alpha, \sigma} \in \mathcal{F}_n^{RRR}$ with arbitrary $n$ is Lipschitz-continuous in the space of inputs $D_v^{RRR}$ with Lipschitz-constant $L_R^{RRR} = 2(1 - \alpha_{min})$.

The proof is provided in Appendix E.

In the following, we consider the parameter space $\Theta_{RRR}$ to be discrete and finite. Finally, we define the class of reservoir functionals associated with the class of RRR reservoir maps as
$$ \mathcal{H}_{n,R_{\max},C_{\max}}^{\text{RRR}} := \left\{ H : (D_v^{\text{RRR}})^{\mathbb{Z}_-} \to \mathbb{R} \mid H(\mathbf{v}) = h(H^{T^{\alpha,\sigma}}(\mathbf{v})), h \in \mathcal{F}_{n,R_{\max},C_{\max}}^{\text{poly}}, T^{\alpha,\sigma} \in \mathcal{F}_n^{\text{RRR}} \right\}. \quad (14) $$

Using the results of this section as well as the previous sections, we can derive an upper bound on the generalisation error of the RRR and PTR quantum reservoir classes.

## 6 Generalisation Bounds of the Quantum Reservoir Classes

We are now ready to present the main results. We can apply (Gonon et al., 2020, Theorem 14) to the RRR and PTR reservoirs, with the appropriate constants. To highlight the dependence in the reservoir parameters, and to simplify notation, we omit the dependencies on parameters that are not determined by the reservoir map or the readout map. The explicit expressions can be found in Appendix H.

**Theorem 11** Let $n \in \mathbb{N}$ be a fixed number of qubits. Consider the PTR class of reservoir functionals $\mathcal{H}_{n,R_{\max},C_{\max}}^{\text{PTR}}$ defined in (12) as well as the input and target processes defined in Section 2, fulfilling the hypotheses $A^{\text{IO}}$. Assume that Assumption (A$^{\Theta}$) is verified, with $|\Theta_{\text{PTR}}|$ the size of the parameter space. Define $\zeta_{\max} := \max(r_{\text{PTR}}(\epsilon_{\text{PTR}}), D_{w^y}, D_{w^v})$. Then for all $m \in \mathbb{N}$ such that $\log(m) < m \log(\zeta_{\max}^{-1})$ and for all $\delta \in (0, 1)$, with probability at least $1 - \delta$ we have

$$ \sup_{H \in \mathcal{H}_{n,R_{\max},C_{\max}}^{\text{PTR}}} \left| R(H) - \hat{R}_m(H) \right| \le \mathcal{O} \left( \bar{L}_h^{\text{poly}} \max \left\{ \frac{P_1^p}{m}, P_2^p \frac{\log m}{m}, P_3^p \sqrt{\frac{\log m}{m}}, P_4^p \sqrt{\frac{\log 4/\delta}{2m}} \right\} \right) $$

where the reservoir parameter dependent expressions are given by

$$ P_1^p = \max \left\{ \frac{r_{\text{PTR}}(\epsilon_{\text{PTR}})}{1 - r_{\text{PTR}}(\epsilon_{\text{PTR}})}, 2L_{\ell} \zeta_{\max}^{-1}, L_{\ell} C_y \zeta_{\max}^{-1} \right\} $$

$$ P_2^p = \max \{1, C_{\max}, 2\zeta_{\max}^{-1}C_v\} $$

$$ P_3^p = \frac{1}{\bar{L}_h^{\text{poly}} \sqrt{\log \zeta_{\max}^{-1}}} \max \left\{ |\Theta_{\text{PTR}}| \cdot R_{\max} \cdot \left( \binom{n+R_{\max}}{R_{\max}} - 1 \right), C_{\max} \right\} $$

$$ P_4^p = \frac{1}{1 - r_{\text{PTR}}(\epsilon_{\text{PTR}})} \max \{r_{\text{PTR}}(\epsilon_{\text{PTR}}), 2M_{\xi}L_v \|w^v\|_1\} $$

and $\bar{L}_h^{\text{poly}}$ is as in (7), i.e. $\bar{L}_h^{\text{poly}} = R_{\max} \cdot n\sqrt{2^n} \left( \binom{n+R_{\max}}{R_{\max}} - 1 \right)$.

The proof can be found in Appendix F.

We can proceed analogously for the RRR reservoir class, distinguishing the case where $r_0 = r_1$ from the case where $r_0 > r_1$ to obtain the following result:
**Theorem 12** Let $n \in \mathbb{N}$ be a fixed number of qubits. Consider the RRR class of reservoir functionals $\mathcal{H}_{n,R_{\max},C_{\max}}^{RRR}$ defined in (14) as well as the input and target processes defined in Section 2, fulfilling the hypotheses $A^{IO}$. Assume that Assumption (A$^\Theta$) is verified, with $|\Theta_{RRR}|$ the size of the parameter space. Define $\zeta_{\max} := \max(r_{RRR}, D_{w^y}, D_{w^v})$. Then for all $m \in \mathbb{N}$ such that $\log(m) < m \log(\zeta_{\max}^{-1})$ and for all $\delta \in (0, 1)$, with probability at least $1 - \delta$ we have

$$ \sup_{H \in \mathcal{H}_{n,R_{\max},C_{\max}}^{RRR}} |R(H) - \hat{R}_m(H)| \\ \le \mathcal{O} \left( \bar{L}_h^{\text{poly}} \max \left\{ \frac{P_1^r}{m}, P_2^r \frac{\log m}{m}, P_3^r \sqrt{\frac{\log m}{m}}, P_4^r \sqrt{\frac{\log 4/\delta}{2m}} \right\} \right) $$

where the reservoir parameter dependent expressions are given by

$$ P_1^r = \begin{cases} \max \left\{ \frac{(1-\alpha_{\min})r_{0,1}}{1-(1-\alpha_{\min})r_{0,1}}, 1 \right\} & \text{if } r_0 = r_1 = r_{0,1} \\ \max \left\{ \frac{1-\epsilon_{RRR}}{\epsilon_{RRR}}, 1 \right\} & \text{if } r_0 > r_1 \end{cases} $$

$$ P_2^r = \max \{1, 2(1 - \alpha_{\min}) \zeta_{\max}^{-1} L_{\ell} C_v, C_{\max} L_{\ell}\} $$

$$ P_3^r = \frac{1}{\bar{L}_h^{\text{poly}} \sqrt{\log \zeta_{\max}^{-1}}} \max \left\{ |\Theta_{RRR}| \cdot R_{\max} \cdot \left( \binom{n + R_{\max}}{R_{\max}} - 1 \right), C_{\max} \right\} $$

$$ P_4^r = \begin{cases} \frac{(1-\alpha_{\min})}{1-(\alpha_{\min})r_{0,1}} \max \{r_{0,1}, 2M_{\xi}L_v \|w^v\|_1\} & \text{if } r_0 = r_1 = r_{0,1} \\ \frac{1}{\epsilon_{RRR}} \max \{1 - \epsilon_{RRR}, 2(1 - \alpha_{\min}) M_{\xi}L_v \|w^v\|_1\} & \text{if } r_0 > r_1 \end{cases} $$

and $\bar{L}_h^{\text{poly}}$ is as in (7), i.e. $\bar{L}_h^{\text{poly}} = R_{\max} \cdot n\sqrt{2^n} \left( \binom{n+R_{\max}}{R_{\max}} - 1 \right)$.

The proof can be found in Appendix G.

From the above theorems we see that the main scaling issues of the generalisation error with respect to the number of qubits comes from the Lipschitz constant of the polynomial readout.

**Remark 13** Note that we have analysed the generalisation error in an idealised setting, without taking into account the estimation of the observables. In reality, multiple shots are necessary to average over the measurement results. The risk bound presented in this paper can help prevent overfitting in the idealised setting of stationary time series without noise.

## 7 Discussion

In this paper, we have established bounds on the Rademacher complexity of a general class of quantum reservoir maps with different readout classes, as well as risk bounds for the more specific quantum reservoir classes which we have introduced in Section 5, based on the classes introduced by Chen and Nurdin (2019) and Chen et al. (2020).
If one wants to use the aforementioned classes for a forecasting task, our risk bounds can give an idea on how to choose a convenient set of parameters that are more likely to generalise well.

On the one hand, we find risk bounds that scale as $O(\sqrt{\log m/m})$ in the number of training samples $m$, which suggests a convergence rate for the generalisation error. On the other hand, our analysis shows that the choice of the class of readout maps can significantly influence the way that the risk bound scales with the number of qubits. As shown in Section 6, the choice of a polynomial readout class leads to the risk bounds scaling as $O(n\sqrt{2^n}(n+R_{max}/R_{max}))$ for a fixed number of qubits $n$, where $R_{max}$ is an upper bound on the degree of the polynomial. As discussed in Section 4.2, this problem can be somewhat alleviated by choosing either a linear readout function, or a polynomial readout of fixed (small) degree, in a fixed number of variables.

The risk bounds we provide are not directly comparable to the risk bounds established for universal classical reservoir classes, since the modified quantum reservoir classes studied here are a priori no longer universal.

The bad scaling of the risk bound in the number of qubits motivates the search for universal reservoir classes that do not require a polynomial readout class. Grigoryeva and Ortega (2018a) have shown that the reservoir class called State Affine System (SAS) is universal even with a linear readout. Martinez-Pena and Ortega (2023) showed that all input-dependent quantum reservoir maps can be written in a similar form, and thus universality results from Grigoryeva and Ortega (2018a) can be applied to reservoir systems with linear readouts when the necessary conditions are verified. However, verifying these conditions is not trivial; apart from the difficulty of finding explicit matrix representations of the CPTP maps, which is necessary to apply the results from Martinez-Pena and Ortega (2023), establishing risk bounds on those representations is not an easy undertaking, as these matrices necessarily depend on the input variables, making it difficult to separate inputs and reservoir states in the decomposition of the reservoir functional.

Another potential avenue is to utilise the rather general results in Sannia et al. (2024) which establishes conditions on the Lindbladian for the reservoir class to be universal when using spatial multiplexing. In forthcoming work, we develop device-specific classes that can harness the natural dynamics of specific quantum systems, including experimental verification of the theoretical results.

## Acknowledgments and Disclosure of Funding

This work was supported by the European Union's Horizon 2020 research and innovation programs under grant agreement No. 817482 (PASQuanS) and No. 101079862 (PASQuanS2) as well as the ANR projects Q-COAST (ANR-19-CE48-0003) and IGNITION (ANR-21-CE47-0015). Part of this research was conducted during a visit to the Institute for Mathematical and Statistical Innovation (IMSI), which is supported by the National Science Foundation (Grant No. DMS-1929348). The authors would like to thank Paulin Jacquot for useful discussions and technical advice.
Appendix A. Proof of Proposition 3

In the following we show that all readout functions in the class $\mathcal{F}_{n,R_{\max},C_{\max}}^{\text{poly}}$ are Lipschitz-continuous in the space of square $2^n \times 2^n$ matrices $\rho$ with $\text{tr}[\sigma_Z^{(i)}\rho] \in [-1, 1]$ for all $i$, with Lipschitz-constants bounded by $\bar{L}_h^{\text{poly}}$ as defined in (7).

Let $\rho_1, \rho_2 \in \bar{\mathbb{B}}^+ \cup 0$ be two density matrices or the square $2^n \times 2^n$ matrix with only zero entries. Denote by $g_i(\rho) := \text{tr}[\sigma_Z^{(i)}\rho] \in [-1, 1]$ the expectation of operator $\sigma_Z^{(i)}$ (i.e. of the measurement operator of the $i$-th qubit on the Z-axis). The function $g_i$ is a Lipschitz-continuous function with Lipschitz-constant $\sqrt{2^n}$ for all $i$, since

$$ \begin{aligned} \left| \text{tr}(\sigma_Z^{(i)} \rho_1) - \text{tr}(\sigma_Z^{(i)} \rho_2) \right| &= \left| \text{tr}(\sigma_Z^{(i)} (\rho_1 - \rho_2)) \right| \\ &\le \sqrt{\text{tr}(\sigma_Z^{(i)\dagger} \sigma_Z^{(i)})} \sqrt{\text{tr}((\rho_1 - \rho_2)^\dagger (\rho_1 - \rho_2))} \\ &= \sqrt{2^n} \|\rho_1 - \rho_2\|_2 \end{aligned} $$

where we have used linearity of the trace operator and matrix multiplication in the first step, the Cauchy-Schwarz inequality in the second step, and the fact that $\sigma_Z^{(i)}$ is hermitian and unitary in the last step.

Next, define

$$ f(x_{i_1}, \dots, x_{i_n}) = w_{i_1, \dots, i_n}^{r_{i_1}, \dots, r_{i_n}} (x_{i_1})^{r_{i_1}} \dots (x_{i_n})^{r_{i_n}} $$

where $x_{i_k} \in [-1, 1]$ for all $i_k$. Then for each $x_{i_k}$, the partial derivative of $f(x_{i_1}, \dots, x_{i_n})$ with respect to $x_{i_k}$ is

$$ D_{x_{i_k}} f(x_{i_1}, \dots, x_{i_n}) = w_{i_1, \dots, i_n}^{r_{i_1}, \dots, r_{i_n}} r_{i_k} (x_{i_1})^{r_{i_1}} \dots (x_{i_k})^{r_{i_k-1}} \dots (x_{i_n})^{r_{i_n}} . $$

It follows that for each $x_{i_k}$, $\|D_{x_{i_k}} f(x_{i_1}, \dots, x_{i_n})\|_2 \le \max_{r_{i_1}, \dots, r_{i_n}} |w_{i_1, \dots, i_n}^{r_{i_1}, \dots, r_{i_n}}| \cdot R \le R$. Recall that an everywhere differentiable multivariate function is Lipschitz continuous if it has bounded partial derivatives, in which case the Lipschitz constant is given by the square root of the sum of the square of the upper bounds of the partial derivatives. That means that the Lipschitz constant $L_f$ of $f$ is given by

$$ L_f = \sqrt{S_{i_1}^2 + \dots + S_{i_n}^2} $$

where $S_{i_k}$ is the upper bound of the $i$th partial derivative, i.e. $R$ in our case. This can be shown by applying the mean-valued theorem in multiple variables as well as the Cauchy-Schwarz inequality, see for example (Eriksson et al., 2003, Theorem 54.2 and following discussion). Thus, the function $f$ is Lipschitz-continuous with Lipschitz constant $L_f = \sqrt{n} \cdot R$. Then we can write

$$ \begin{aligned} & \|f(g_1(\rho_1), \dots, g_n(\rho_1)) - f(g_1(\rho_2), \dots, g_n(\rho_2))\|_2 && (15) \\ & \le L_f \|(g_1(\rho_1) - g_1(\rho_2), \dots, g_n(\rho_1) - g_n(\rho_2))\|_2 \\ & = L_f \sqrt{(g_1(\rho_1) - g_1(\rho_2))^2 + \dots + (g_n(\rho_1) - g_n(\rho_2))^2} \\ & \le L_f \sqrt{n \cdot 2^n \|\rho_1 - \rho_2\|_2^2} = n \cdot R \cdot \sqrt{2^n} \|\rho_1 - \rho_2\|_2. \end{aligned} $$
Note that the sum of Lipschitz-continuous functions is Lipschitz-continuous, where the Lipschitz constant is the sum of the Lipschitz constants of the individual summands. In our case, all summands have the same Lipschitz-constant $L_{f,g} = R \cdot n\sqrt{2^n}$ so that the Lipschitz-constant $L_h$ of $h$ is simply $L_{f,g}$ multiplied by the number of summands. To find this number, note that the way the $r_{i_1}, \dots, r_{i_n}$ are partitioned is in fact a known combinatorics problem, called stars and bars or multichOOSE. The problem is formulated as follows: For a fixed number $d$ of stars, and a fixed number $n$ of bins (delimited by $n-1$ bars), how many ways are there to distribute the $d$ stars into the $n$ bins (allowing empty bins). Then the number of ways to partition the integer $d$ into an $n$-tuple is given by $\binom{n+d-1}{d}$, see for example (Feller, 1968, Chapter II.5). Summing over $d = 1, \dots, R$, we can then bound the complete Lipschitz-constant of $h$ by

$$L_h \le R \cdot n\sqrt{2^n} \left( \binom{n+R}{R} - 1 \right).$$

In particular, this is bounded by

$$\bar{L}_h^{\text{poly}} = R_{\max} \cdot n\sqrt{2^n} \left( \binom{n+R_{\max}}{R_{\max}} - 1 \right).$$

## Appendix B. Proof of Theorem 4

For the class $\mathcal{H}_{n,R_{\max},C_{\max}}^{\text{poly}}$ of quantum reservoir functionals, begin by writing

$$ \begin{aligned} & \mathbb{E} \left[ \sup_{H \in \mathcal{H}_{n, R_{\max}, C_{\max}}^{\text{poly}}} \left| \sum_{j=0}^{k-1} \varepsilon_j H(\tilde{\mathbf{V}}^{(j)}) \right| \right] \\ &= \mathbb{E} \left[ \sup_{\substack{h \in \mathcal{H}_{n, R_{\max}, C_{\max}}^{\text{poly}} \\ T^{\theta} \in \mathcal{H}_{n}^{\text{RC}}(\Theta)}} \left| \sum_{j=0}^{k-1} \varepsilon_j h \left( \prod_{i=0}^{\infty} T^{\theta} \left( \tilde{V}_{-i}^{(j)} \right) \rho_{-\infty} \right) \right| \right] \\ &= \mathbb{E} \left[ \sup_{\substack{w_{i_1, \dots, i_n}^{r_{i_1}, \dots, r_{i_n}} \in [0,1] \\ 1 \le R \le R_{\max}, \ 0 \le C \le C_{\max}}} \sup_{\theta \in \Theta} \left| \sum_{j=0}^{k-1} \varepsilon_j \left( C \sum_{d=1}^{R} \sum_{i_1=1}^{n} \sum_{i_2=i_1+1}^{n} \dots \sum_{i_n=i_{n-1}+1}^{n} \sum_{r_{i_1}+\dots+r_{i_n}=d} w_{i_1, \dots, i_n}^{r_{i_1}, \dots, r_{i_n}} \operatorname{tr} \left[ \sigma_Z^{(i_1)} \rho_0^{\theta, j} \right]^{r_{i_1}} \dots \operatorname{tr} \left[ \sigma_Z^{(i_n)} \rho_0^{\theta, j} \right]^{r_{i_n}} \right) \right| \right] \end{aligned} $$

where we have written $\rho_0^{\theta,j} := \prod_{i=0}^{\infty} T^{\theta}(\tilde{V}_{-i}^{(j)}) \rho_{-\infty}$. Using the triangle inequality as well as the subadditivity of the supremum, we can bound the last term as follows:
$$
\begin{align}
& \mathbb{E} \left[ \sup_{\substack{w^{r_{i_1}, \dots, r_{i_n}}_{i_1, \dots, i_n} \in [0,1] \\ 1 \le R \le R_{\max}, \ 0 \le C \le C_{\max}}} \sup_{\theta \in \Theta} \left| \sum_{j=0}^{k-1} \varepsilon_j \left( C + \sum_{d=1}^{R} \sum_{i_1=1}^{n} \sum_{i_2=i_1+1}^{n} \dots \sum_{i_n=i_{n-1}+1}^{n} \sum_{r_{i_1}+\dots+r_{i_n}=d} w^{r_{i_1}, \dots, r_{i_n}}_{i_1, \dots, i_n} \operatorname{tr} [\sigma_Z^{(i_1)} \rho_0^{\theta,j}]^{r_{i_1}} \dots \operatorname{tr} [\sigma_Z^{(i_n)} \rho_0^{\theta,j}]^{r_{i_n}} \right) \right| \right. \nonumber \\
& \le \mathbb{E} \left[ \sup_{0 \le C \le C_{\max}} \left| \sum_{j=0}^{k-1} \varepsilon_j C \right| \right] \nonumber \\
& \quad + \mathbb{E} \left[ \sup_{1 \le R \le R_{\max}} \sum_{d=1}^{R} \sum_{i_1=1}^{n} \dots \sum_{i_n=i_{n-1}+1}^{n} \sum_{r_{i_1}+\dots+r_{i_n}=d} \sup_{w^{r_{i_1}, \dots, r_{i_n}}_{i_1, \dots, i_n} \in [0,1]} \left| \sum_{j=0}^{k-1} \varepsilon_j w^{r_{i_1}, \dots, r_{i_n}}_{i_1, \dots, i_n} \operatorname{tr} [\sigma_Z^{(i_1)} \rho_0^{\theta,j}]^{r_{i_1}} \dots \operatorname{tr} [\sigma_Z^{(i_n)} \rho_0^{\theta,j}]^{r_{i_n}} \right| \right]. \tag{16}
\end{align}
$$

For the first term, we can write

$$
\left( \mathbb{E} \left[ \sup_{0 \le |C| \le C_{\max}} \left| \sum_{j=0}^{k-1} \varepsilon_j C \right| \right] \right)^2 \le \mathbb{E} \left[ \sup_{0 \le |C| \le C_{\max}} |C|^2 \left| \sum_{j=0}^{k-1} \varepsilon_j \right|^2 \right] = \sum_{j=0}^{k-1} \mathbb{E} [\varepsilon_j^2] C_{\max}^2 = k C_{\max}^2 \quad (17)
$$

where we used Jensen's inequality twice in the first inequality, as well as the fact that $C$ is
independent of $j$ and the submultiplicativity of the norm; and the fact that $\mathbb{E}[\varepsilon_j \varepsilon_m] = 0$ if
$m \neq j$ in the first equality.

For the second term (16), we write

$$
\begin{align}
& \mathbb{E} \left[ \sup_{1 \le R \le R_{\max}} \sum_{d=1}^{R} \sum_{i_1=1}^{n} \right. \nonumber \\
& \qquad \left. \dots \sum_{i_n=i_{n-1}+1}^{n} \sum_{r_{i_1}+\dots+r_{i_n}=d} \sup_{w^{r_{i_1}, \dots, r_{i_n}}_{i_1, \dots, i_n} \in [0,1]} \left| \sum_{j=0}^{k-1} \varepsilon_j w^{r_{i_1}, \dots, r_{i_n}}_{i_1, \dots, i_n} \operatorname{tr} [\sigma_Z^{(i_1)} \rho_0^{\theta,j}]^{r_{i_1}} \dots \operatorname{tr} [\sigma_Z^{(i_n)} \rho_0^{\theta,j}]^{r_{i_n}} \right| \right] \nonumber \\
& \le \sum_{d=1}^{R_{\max}} \sum_{i_1=1}^{n} \dots \sum_{i_n=i_{n-1}+1}^{n} \sum_{r_{i_1}+\dots+r_{i_n}=d} \mathbb{E} \left[ \sup_{w^{r_{i_1}, \dots, r_{i_n}}_{i_1, \dots, i_n} \in [0,1]} \left| \sum_{j=0}^{k-1} \varepsilon_j w^{r_{i_1}, \dots, r_{i_n}}_{i_1, \dots, i_n} \operatorname{tr} [\sigma_Z^{(i_1)} \rho_0^{\theta,j}]^{r_{i_1}} \dots \operatorname{tr} [\sigma_Z^{(i_n)} \rho_0^{\theta,j}]^{r_{i_n}} \right| \right]. \tag{18}
\end{align}
$$

We can once again take arguments that are independent of *j* out of the sum and use the
submultiplicativity of the norm to bound the term behind the sums:
$$
\begin{align*}
& \mathbb{E} \left[ \sup_{\substack{w^{r_{i_1}, \dots, r_{i_n}}_{i_1, \dots, i_n} \in [0,1] \\ \theta \in \Theta}} \left| \sum_{j=0}^{k-1} \varepsilon_j w^{r_{i_1}, \dots, r_{i_n}}_{i_1, \dots, i_n} \operatorname{tr} [\sigma_Z^{(i_1)} \rho_0^{\theta, j}]^{r_{i_1}} \cdots \operatorname{tr} [\sigma_Z^{(i_n)} \rho_0^{\theta, j}]^{r_{i_n}} \right| \right] \\
&= \mathbb{E} \left[ \sup_{\substack{w^{r_{i_1}, \dots, r_{i_n}}_{i_1, \dots, i_n} \in [0,1] \\ \theta \in \Theta}} \left| w^{r_{i_1}, \dots, r_{i_n}}_{i_1, \dots, i_n} \sum_{j=0}^{k-1} \varepsilon_j \operatorname{tr} [\sigma_Z^{(i_1)} \rho_0^{\theta, j}]^{r_{i_1}} \cdots \operatorname{tr} [\sigma_Z^{(i_n)} \rho_0^{\theta, j}]^{r_{i_n}} \right| \right] \\
&\le \mathbb{E} \left[ \sup_{\substack{w^{r_{i_1}, \dots, r_{i_n}}_{i_1, \dots, i_n} \in [0,1] \\ \theta \in \Theta}} \left| w^{r_{i_1}, \dots, r_{i_n}}_{i_1, \dots, i_n} \sum_{j=0}^{k-1} \varepsilon_j \operatorname{tr} [\sigma_Z^{(i_1)} \rho_0^{\theta, j}]^{r_{i_1}} \cdots \operatorname{tr} [\sigma_Z^{(i_n)} \rho_0^{\theta, j}]^{r_{i_n}} \right| \right] \\
&\le \mathbb{E} \left[ \sup_{\theta \in \Theta} \left| \sum_{j=0}^{k-1} \varepsilon_j \operatorname{tr} [\sigma_Z^{(i_1)} \rho_0^{\theta, j}]^{r_{i_1}} \cdots \operatorname{tr} [\sigma_Z^{(i_n)} \rho_0^{\theta, j}]^{r_{i_n}} \right| \right].
\end{align*}
$$

Using the fact that the parameter space is finite, we write

$$
\begin{align}
& \mathbb{E} \left[ \sup_{\theta \in \Theta} \left| \sum_{j=0}^{k-1} \varepsilon_j \operatorname{tr} [\sigma_Z^{(i_1)} \rho_0^{\theta,j}]^{r_{i_1}} \cdots \operatorname{tr} [\sigma_Z^{(i_n)} \rho_0^{\theta,j}]^{r_{i_n}} \right| \right] \tag{19} \\
& \le \sum_{\theta \in \Theta} \mathbb{E} \left[ \left| \sum_{j=0}^{k-1} \varepsilon_j \operatorname{tr} [\sigma_Z^{(i_1)} \rho_0^{\theta,j}]^{r_{i_1}} \cdots \operatorname{tr} [\sigma_Z^{(i_n)} \rho_0^{\theta,j}]^{r_{i_n}} \right| \right]. \nonumber
\end{align}
$$

We can bound the term behind the sum over the parameter space as follows:

$$
\begin{align}
& \left( \mathbb{E} \left[ \left| \sum_{j=0}^{k-1} \varepsilon_j \operatorname{tr} [\sigma_Z^{(i_1)} \rho_0^{\theta,j}]^{r_{i_1}} \cdots \operatorname{tr} [\sigma_Z^{(i_n)} \rho_0^{\theta,j}]^{r_{i_n}} \right| \right] \right)^2 \tag{20} \\
&\le \mathbb{E} \left[ \left| \sum_{j=0}^{k-1} \varepsilon_j \operatorname{tr} [\sigma_Z^{(i_1)} \rho_0^{\theta,j}]^{r_{i_1}} \cdots \operatorname{tr} [\sigma_Z^{(i_n)} \rho_0^{\theta,j}]^{r_{i_n}} \right|^2 \right] \nonumber \\
&= \sum_{j=0}^{k-1} \mathbb{E} [\varepsilon_j^2] \mathbb{E} \left[ \left| \operatorname{tr} [\sigma_Z^{(i_1)} \rho_0^{\theta,j}]^{r_{i_1}} \cdots \operatorname{tr} [\sigma_Z^{(i_n)} \rho_0^{\theta,j}]^{r_{i_n}} \right|^2 \right] \nonumber \\
&\le k \nonumber
\end{align}
$$

where we used Jensen's inequality in the first inequality, independence between the Rademacher variables and the input variables and the fact that $\mathbb{E}[\varepsilon_j \varepsilon_m] = 0$ if $m \neq j$ in the first equality, as well as the fact that $\operatorname{tr}[\sigma_Z^{(i)} \rho_0^{\theta,j}]^{r_i} \in [-1, 1]$.
Finally, combining (19) and (20), and using the same combinatorial argument as in Appendix A to bound the number of terms of sums in (18), we have

$$ \sum_{d=1}^{R_{\max}} \sum_{i_1=1}^{n} \cdots \sum_{i_n=i_{n-1}+1}^{n} \sum_{r_{i_1}+\cdots+r_{i_n}=d} \mathbb{E} \left[ \sup_{w_{i_1, \dots, i_n} \in [0,1] \atop \theta \in \Theta} \left| \sum_{j=0}^{k-1} \varepsilon_j w_{i_1, \dots, i_n}^{r_{i_1}, \dots, r_{i_n}} \operatorname{tr} \left[ \sigma_Z^{(i_1)} \rho_0^{\theta, j} \right]^{r_{i_1}} \cdots \operatorname{tr} \left[ \sigma_Z^{(i_n)} \rho_0^{\theta, j} \right]^{r_{i_n}} \right| \right] \\ \le R_{\max} \left( \binom{n+R_{\max}}{R_{\max}} - 1 \right) |\Theta| \sqrt{k}. $$

Combining this with (17), we get the final result.

## Appendix C. Proof of Proposition 6

The proof of Proposition 6 proceeds similarly to the proof of Proposition 3, with a few modifications in the constants.

First, define the multivariate function

$$ f(x_{i_1}, \dots, x_{i_\ell}) = w_{i_1} \cdots w_{i_\ell} x_{i_1} \cdots x_{i_\ell}. $$

Then, the partial derivatives of $f$ with respect to each variable $x_{i_k}$ is given by

$$ D_{x_{i_k}} f(x_{i_1}, \dots, x_{i_\ell}) = w_{i_1} \cdots w_{i_\ell} x_{i_1} \cdots x_{i_{k-1}} x_{i_{k+1}} \cdots x_{i_\ell} $$

and for each variable $x_{i_k}$ the partial derivative is bounded by

$$ |D_{x_{i_k}} f(x_{i_1}, \dots, x_{i_\ell})| \le |w_{i_1}| \cdots |w_{i_\ell}| \le 1 $$

where the last inequality follows from the fact that $|w_{i_j}| \le 1$ by design. Using the same argument as in Appendix A, we conclude that the function $f$ is Lipschitz-continuous with Lipschitz constant $L_f = \ell$.

Then, with the functions $g_i(\rho)$ as defined in Appendix A and following a similar argument as in (15), we can write

$$ \begin{aligned} & \|f(g_1(\rho_1), \dots, g_\ell(\rho_1)) - f(g_1(\rho_2), \dots, g_\ell(\rho_2))\|_2 \\ & \le L_f \|(g_1(\rho_1) - g_1(\rho_2), \dots, g_\ell(\rho_1) - g_\ell(\rho_2))\|_2 \\ & = L_f \sqrt{(g_1(\rho_1) - g_1(\rho_2))^2 + \dots + (g_\ell(\rho_1) - g_\ell(\rho_2))^2} \\ & \le L_f \sqrt{\ell \cdot 2^n \|\rho_1 - \rho_2\|_2^2} = \ell \cdot \sqrt{\ell \cdot 2^n} \|\rho_1 - \rho_2\|_2. \end{aligned} $$

Finally, by noting that there are $\ell$ sums in the readout function $h^{\mathrm{SM}}$, we conclude that $h^{\mathrm{SM}}$ is Lipschitz-continuous w.r.t. $\rho$, with Lipschitz constant $L_h^{\mathrm{SM}} = \ell^{5/2} \sqrt{2^n}$, and an upper bound on the Lischitz constants of the class $\mathcal{F}_{n, R_{\max}, C_{\max}}^{\mathrm{SM}}$ is given by

$$ \bar{L}_h^{\mathrm{SM}} = \ell^{5/2} \max \sqrt{2^n}. $$
## Appendix D. Proof of Lemma 9

We begin by proving contractivity in the space of density operators. The proof makes use of (Rastegin, 2012, Proposition 1), which states that for an operator $\tilde{Q}$ defined on a composite system $\mathcal{H}_A \otimes \mathcal{H}_B$ with $\dim(\mathcal{H}_A) = m$ and $\dim(\mathcal{H}_B) = n$ we have

$$ \|\mathrm{tr}_B[\tilde{Q}]\|_2 \le \sqrt{n} \|\tilde{Q}\|_2 . \qquad (21) $$

To show that there exists an $r_{PTR} \in (0, 1)$ such that for fixed parameters $(J, \alpha, \gamma)$, for any two quantum states $\rho_1, \rho_2$ and for all $v \in D_v$

$$ \|T^{J,\gamma,\tau}(v, \rho_1) - T^{J,\gamma,\tau}(v, \rho_2)\|_2 \le r_{PTR} \|\rho_1 - \rho_2\|_2, $$

note that $T^{J,\gamma,\tau}$ is a linear map, and $\rho_1, \rho_2$ being quantum states, we have $\mathrm{tr}[\rho_1 - \rho_2] = 0$ so that we can restrict $T^{J,\gamma,\tau}$ to the hyperplane $\mathcal{B}_0$ of $2^n \times 2^n$ traceless Hermitian operators. In the remainder of this proof we write $T = T^{J,\gamma,\tau}$ and $H = H(J, \gamma)\tau$ for ease of notation, as the parameters do not intervene in the calculations. Define $\mathcal{V} = \begin{pmatrix} v & 0 \\ 0 & 1-v \end{pmatrix}$. Using (21) we have for all $A \in \mathcal{B}_0$,

$$ \|T(v, A)\|_2 = \|\mathrm{tr}_0[e^{-iH}A \otimes \mathcal{V}e^{iH}]\|_2 \le \sqrt{2} \|e^{-iH}A \otimes \mathcal{V}e^{iH}\|_2. $$

Notice also that

$$ \begin{aligned} \|e^{-iH}A \otimes \mathcal{V}e^{iH}\|_2 &= \|A \otimes \mathcal{V}\|_2 = \sqrt{\mathrm{tr}[[A \otimes \mathcal{V}]^* [A \otimes \mathcal{V}]]} \\ &= \sqrt{\mathrm{tr}[A^*A \otimes \mathcal{V}^*\mathcal{V}]} = \sqrt{\mathrm{tr}[A^*A] \mathrm{tr}[\mathcal{V}^*\mathcal{V}]} = \|A\|_2 (v^2 + (1-v)^2). \end{aligned} $$

We thus have, for $v \in [\frac{1}{2}(1 - \sqrt{2} - 1) + \epsilon_{PTR}, \frac{1}{2}(1 + \sqrt{2} - 1) - \epsilon_{PTR}]$,

$$ \sup_{\substack{A \in \mathcal{B}_0 \\ \|A\|_2 = 1}} \|T(v, A)\|_2 \le r_{PTR}(\epsilon_{PTR}) $$

with

$$ \begin{aligned} r_{PTR}(\epsilon_{PTR}) &= \sqrt{2} \left( \left( \frac{1}{2}(1 + \sqrt{2} - 1) - \epsilon_{PTR} \right)^2 + \left( 1 - \frac{1}{2}(1 + \sqrt{2} - 1) + \epsilon_{PTR} \right)^2 \right) \\ &= 2\sqrt{2} \left( \epsilon_{PTR}^2 - \epsilon_{PTR}\sqrt{2} - 1 \right) + 1 \end{aligned} $$

where we indicate the dependence of the contractivity constant on the choice of $\epsilon_{PTR}$. We proceed similarly to prove Lipschitz-continuity: For any inputs $v_1, v_2 \in D_v$ and any density matrix $\rho \in \tilde{\mathcal{B}}^+$ we have
$$
\begin{align*}
& \|T(v_1, \rho) - T(v_2, \rho)\|_2 \\
&= \left\| \operatorname{tr}_0 \left[ e^{-iH} \rho \otimes \begin{pmatrix} v_1 & 0 \\ 0 & 1-v_1 \end{pmatrix} e^{iH} \right] - \operatorname{tr}_0 \left[ e^{-iH} \rho \otimes \begin{pmatrix} v_2 & 0 \\ 0 & 1-v_2 \end{pmatrix} e^{iH} \right] \right\|_2 \\
&= \left\| \operatorname{tr}_0 \left[ e^{-iH} \rho \otimes \begin{pmatrix} v_1 & 0 \\ 0 & 1-v_1 \end{pmatrix} e^{iH} - e^{-iH} \rho \otimes \begin{pmatrix} v_2 & 0 \\ 0 & 1-v_2 \end{pmatrix} e^{iH} \right] \right\|_2 \\
&= \left\| \operatorname{tr}_0 \left[ e^{-iH} \left( \rho \otimes \begin{pmatrix} v_1 & 0 \\ 0 & 1-v_1 \end{pmatrix} - \rho \otimes \begin{pmatrix} v_2 & 0 \\ 0 & 1-v_2 \end{pmatrix} \right) e^{iH} \right] \right\|_2 \\
&= \left\| \operatorname{tr}_0 \left[ e^{-iH} \rho \otimes \begin{pmatrix} v_1 - v_2 & 0 \\ 0 & v_2 - v_1 \end{pmatrix} e^{iH} \right] \right\|_2 \\
&\le \sqrt{2} \left\| \rho \otimes \begin{pmatrix} v_1 - v_2 & 0 \\ 0 & v_2 - v_1 \end{pmatrix} \right\|_2 \\
&= \sqrt{2 \operatorname{tr} \left[ \left( \rho^\dagger \otimes \begin{pmatrix} v_1 - v_2 & 0 \\ 0 & v_2 - v_1 \end{pmatrix} \right) \left( \rho \otimes \begin{pmatrix} v_1 - v_2 & 0 \\ 0 & v_2 - v_1 \end{pmatrix} \right) \right]} \\
&= \sqrt{2 \operatorname{tr} (\rho^\dagger \rho) \operatorname{tr} \left( \begin{pmatrix} v_1 - v_2 & 0 \\ 0 & v_2 - v_1 \end{pmatrix}^2 \right)} \\
&\le \sqrt{2 \cdot 2(v_1 - v_2)^2} = 2\|v_1 - v_2\|_2
\end{align*}
$$

where the second equality follows from the linearity of the partial trace operator, the first inequality follows from (21), and the second inequality follows from the fact that the purity of any quantum state is bounded by one.

Appendix E. Proof of Lemma 10

We begin by proving contractivity in the space of density operators: For any $\rho_1, \rho_2 \in \bar{\mathcal{B}}^+$ and for any fixed $v \in D_v^{RRR}$ we have

$$
\begin{align*}
& \|T^{\alpha,\sigma}(v, \rho_1) - T^{\alpha,\sigma}(v, \rho_2)\|_2 \\
&= (1-\alpha) \|(vT_0 + (1-v)T_1)\rho_1 - (vT_0 + (1-v)T_1)\rho_2\|_2 \\
&\le (1-\alpha)(v \|T_0(\rho_1) - T_0(\rho_2)\|_2 + (1-v) \|T_1(\rho_1) - T_1(\rho_2)\|_2) \\
&\le (1-\alpha)(vr_0 + (1-v)r_1) \|\rho_1 - \rho_2\|_2 \\
&= (1-\alpha)(v(r_0 - r_1) + r_1) \|\rho_1 - \rho_2\|_2
\end{align*}
$$

where we used the fact that $T_0$ and $T_1$ are $r_0$- and $r_1$-contractions, respectively. If $r_0 = r_1 = r_{0,1}$, the result is immediate noting that $\alpha \in [\alpha_{\min}, \alpha_{\max}] \subset (0,1)$. If $r_0 > r_1$, we have $1 > r_0 + r_1$ which implies $0 < r_0 - r_1 < 1 - 2r_1$, thus

$$v(r_0 - r_1) + r_1 < v(1 - 2r_1) + r_1 = r_1(1 - 2v) + v < 1 - v < 1 - \epsilon_{RRR}$$
where the second to last inequality follows from the fact that $r_1 < 1$ and $v < 1/2$.

The proof of the Lipschitz-continuity is a straightforward application of the triangle inequality. For any $v_1, v_2 \in D_v^{RRR}$ and for any fixed $\rho \in \tilde{B}^+$ we have

$$ 
\begin{aligned}
& \|T^{\alpha, \sigma}(v_1, \rho) - T^{\alpha, \sigma}(v_2, \rho)\|_2 \\
&= \|(1-\alpha)(v_1 T_0 + (1-v_1)T_1)\rho + \alpha\sigma - (1-\alpha)(v_2 T_0 + (1-v_2)T_1)\rho - \alpha\sigma\|_2 \\
&\le (1-\alpha) \|(v_1 - v_2)T_0(\rho) + (v_2 - v_1)T_1(\rho)\|_2 \\
&\le (1-\alpha) |v_1 - v_2| \|T_0(\rho) + T_1(\rho)\|_2 \\
&\le (1-\alpha) |v_1 - v_2| (\|T_0(\rho)\|_2 + \|T_1(\rho_1)\|_2) \\
&\le 2(1-\alpha_{\min}) |v_1 - v_2|
\end{aligned}
$$

where the last inequality comes from the fact that $T_0$ and $T_1$ are CPTP maps and thus produce a density matrix with purity bounded by 1.

## Appendix F. Proof of Theorem 11

The proof of this theorem hinges on the application of (Gonon et al., 2020, Theorem 14) which we restate here while adapting to our particular case:

**Theorem 14** Let $\mathcal{H}_n^{QRC}$ be the hypothesis class of quantum reservoir functionals specified in (5) associated to a class $\mathcal{F}_n^{QRC}(\Theta)$ of reservoir maps verifying the conditions $A^T$ and a class $\mathcal{F}_n^O$ of readout maps verifying conditions $A^h$. Suppose that both the input $V$ and the target $Y$ processes verify the hypotheses $A^{IO}$. Assume additionally that there exists a constant $C_{QRC}$ such that the Rademacher complexity satisfies $\mathcal{R}_k(\mathcal{H}_n^{QRC}) \le C_{QRC}/\sqrt{k}$. Furthermore, define $\zeta_{\max} := \max(r, D_{wy}, D_{w^v})$. Then there exist constants $C_0, C_1, C_2, C_3, C_{bd} > 0$ such that for all $m \in \mathbb{N}^+$ satisfying $\log m < m \log \zeta_{\max}^{-1}$ and for all $\delta \in (0, 1)$ with probability at least $1 - \delta$ we have

$$ \sup_{H \in \mathcal{H}_n^{QRC}} \left| R(H) - \hat{R}_m(H) \right| \le \frac{(1-r^m) C_0^{QRC} + C_1}{m} + \frac{C_2 \log m}{m} + \frac{C_3 \sqrt{\log m}}{\sqrt{m}} + \frac{C_{bd} \sqrt{\log \frac{4}{\delta}}}{\sqrt{2m}} \quad (22) $$

with constants

$$ 
\begin{aligned}
C_1 &= \frac{L_{\ell}(2\bar{L}_h + C_y)}{\zeta_{\max}} \\
C_2 &= \frac{2S}{\log \zeta_{\max}^{-1}} + \frac{L_{\ell}L_R\bar{L}_h C_v}{\zeta_{\max} \log \zeta_{\max}^{-1}} \\
C_3 &= \frac{2L_{\ell}C_{QRC}}{\sqrt{\log \zeta_{\max}^{-1}}} \\
C_{bd} &= 2L_{\ell} \left( \frac{\bar{L}_h}{1-r} (r + L_R M_{\xi} L_v \|w^v\|_1) + M_{\xi} L_y \|w^y\|_1 \right)
\end{aligned}
$$
where

$$S = L_{\ell} \bar{L}_{h} + \mathbb{E}[|\ell(0, Y_0)|] + L_{h,0} L_{\ell} \\ C_I = \frac{2L_I \mathbb{E}[\|\xi_0^I\|_2]}{1 - D_{w^I}}, \quad I = v, y$$

and

$$C_0^{QRC} = \frac{2rL_{\ell}\bar{L}_h}{1-r}.$$

Note that in the adaptation of the theorem we have used the fact that any convergent reservoir functional is bounded by one in the space of density matrices, and have replaced the constant $M_f$ by one. The dependence in the upper bound of the Rademacher complexity intervenes in the third summand in (22).

Then we can calculate the constant factor in the Rademacher complexity $C_0^{PTR}$ by plugging the constants from Lemma 9 into Theorem 14:

$$C_0^{PTR} = \frac{2r_{PTR}(\epsilon_{PTR})L_{\ell}\bar{L}_h^{\text{poly}}}{1 - r_{PTR}(\epsilon_{PTR})} \\ = \frac{(2\sqrt{2}(\epsilon_{PTR}^2 - \epsilon_{PTR}\sqrt{\sqrt{2}-1}) + 1)L_{\ell}R_{\text{max}} \cdot n\sqrt{2^n}\binom{n+R_{\text{max}}}{R_{\text{max}}} - 1}{\sqrt{2}(\epsilon_{PTR}^2 - \epsilon_{PTR}\sqrt{\sqrt{2}-1})}.$$

The constant factor in the bound on the Rademacher complexity is given by $C_{PTR} = |\Theta_{PTR}| \cdot R_{\text{max}} \cdot \binom{n+R_{\text{max}}}{R_{\text{max}}} - 1 + C_{\text{max}}$ as shown in Theorem 4.

We can then calculate the explicit bound which we state in Appendix H by injecting the above constants and all relevant constants derived in this document into the bound in (22).

## Appendix G. Proof of Theorem 12

The proof of this theorem is analogous to that of Theorem 11, but distinguishing the case where $r_0 = r_1$ from the case where $r_0 > r_1$.

As in Appendix F, we calculate $C_0^{RRR}$ by plugging the constants from Lemma 10 into Theorem 14 while distinguishing the case when $r_0 = r_1 = r_{0,1}$ from the case when $r_0 > r_1$:

$$C_0^{RRR} = \begin{cases} \frac{2(1-\alpha_{\min})r_{0,1}L_{\ell}R_{\text{max}} \cdot n\sqrt{2^n}\binom{n+R_{\text{max}}}{R_{\text{max}}} - 1}{1-(1-\alpha_{\min})r_{0,1}} & \text{if } r_0 = r_1 = r_{0,1} \\ \frac{2(1-\epsilon_{RRR})L_{\ell}R_{\text{max}} \cdot n\sqrt{2^n}\binom{n+R_{\text{max}}}{R_{\text{max}}} - 1}{\epsilon_{RRR}} & \text{if } r_0 > r_1 \end{cases}$$

By Theorem 4, we have $C_{RRR} = |\Theta_{RRR}| \cdot R_{\text{max}} \cdot \binom{n+R_{\text{max}}}{R_{\text{max}}} - 1 + C_{\text{max}}$.

We can then calculate the explicit bound which we state in Appendix H by injecting the above constants and all relevant constants derived in this document into the bound in (22).
Appendix H. Explicit risk bounds for $\mathcal{H}_{n,R_{\max},C_{\max}}^{\text{PTR}}$ and $\mathcal{H}_{n,R_{\max},C_{\max}}^{\text{RRR}}$

Replacing the general constants $r$, $L_R$ and $C_0^{\text{QRC}}$ in Theorem 14 with their more explicit forms for the PTR reservoir class $\mathcal{H}_{n,R_{\max},C_{\max}}^{\text{PTR}}$ as well as the RRR reservoir class $\mathcal{H}_{n,R_{\max},C_{\max}}^{\text{RRR}}$ that we establish in Appendix F and Appendix G respectively, we obtain the generalisation bounds, which we state in the following with the explicit constants.

Under the same hypotheses as in Theorem 11, for all $\delta \in (0, 1)$, with probability at least $1 - \delta$ we have

$$ \sup_{H \in \mathcal{H}_{n,R_{\max},C_{\max}}^{\text{PTR}}} \left|R(H) - \hat{R}_m(H)\right| \\ \le (C_0^P(n, R_{\max}, \epsilon_{\text{PTR}}) + C_1(n, R_{\max}, \zeta_{\max})) \frac{1}{m} \\ + C_2^P(n, R_{\max}, \zeta_{\max}, C_{\max}) \frac{\log m}{m} \\ + C_3^P(n, R_{\max}, \zeta_{\max}, |\Theta_{\text{PTR}}|, C_{\max}) \sqrt{\frac{\log m}{m}} \\ + C_4^P(n, R_{\max}, \zeta_{\max}, \epsilon_{\text{PTR}}) \sqrt{\frac{\log 4/\delta}{2m}} $$

where

$$ C_0^P(n, R_{\max}, \epsilon_{\text{PTR}}) = \frac{2(1 - r_{\text{PTR}}(\epsilon_{\text{PTR}})^m) r_{\text{PTR}}(\epsilon_{\text{PTR}}) L_{\ell} \bar{L}_h^{\text{poly}}}{1 - r_{\text{PTR}}(\epsilon_{\text{PTR}})} \\ C_1(n, R_{\max}, \zeta_{\max}) = L_{\ell}(2\bar{L}_h^{\text{poly}} + C_y) \zeta_{\max}^{-1} \quad (23) $$

$$ C_2^P(n, R_{\max}, \zeta_{\max}, C_{\max}) = \frac{2(L_{\ell}\bar{L}_h^{\text{poly}} + \mathbb{E}[|\ell(0, Y_0)|] + C_{\max}L_{\ell})}{\log \zeta_{\max}^{-1}} + \frac{2L_{\ell}\bar{L}_h^{\text{poly}}C_v}{\zeta_{\max}\log \zeta_{\max}^{-1}} \\ C_3^P(n, R_{\max}, \zeta_{\max}, |\Theta_{\text{PTR}}|, C_{\max}) = \frac{2L_{\ell}(|\Theta_{\text{PTR}}| \cdot R_{\max} \cdot \binom{n+R_{\max}}{R_{\max}} - 1) + C_{\max}}{\sqrt{\log \zeta_{\max}^{-1}}} \\ C_4^P(n, R_{\max}, \zeta_{\max}, \epsilon_{\text{PTR}}) = 2L_{\ell}\left(\frac{\bar{L}_h^{\text{poly}}}{1 - r_{\text{PTR}}(\epsilon_{\text{PTR}})}(r_{\text{PTR}}(\epsilon_{\text{PTR}}) + 2M_{\xi}L_v\|w^v\|_1) + M_{\xi}L_y\|w^y\|_1\right) $$

and $\bar{L}_h^{\text{poly}}$ is as in (7), i.e. $\bar{L}_h^{\text{poly}} = R_{\max} \cdot n\sqrt{2^n} \binom{n+R_{\max}}{R_{\max}} - 1$, and $r_{\text{PTR}}(\epsilon_{\text{PTR}})$ as in Lemma 9, i.e. $r_{\text{PTR}}(\epsilon_{\text{PTR}}) = 2\sqrt{2}(\epsilon_{\text{PTR}}^2 - \epsilon_{\text{PTR}}\sqrt{2} - 1) + 1$.

Here we have highlighted the dependence in the parameters of the reservoir class in the constants $C_0^P, C_1, C_2^P, C_3^P$ and $C_4^P$.

Similarly, under the same hypotheses as in Theorem 12, for all $\delta \in (0, 1)$, with probability at least $1 - \delta$ we have
$$
\begin{align*}
\sup_{H \in \mathcal{H}_{n, R_{\max}, C_{\max}}^{\text{RRR}}}
& \left| R(H) - \hat{R}_m(H) \right| \\
& \le \left( C_0^{a,R}(n, R_{\max}, \alpha_{\min}, r_{0,1}) + C_1(n, R_{\max}, \zeta_{\max}) \right) \frac{1}{m} \\
& \quad + C_2^R(n, R_{\max}, \zeta_{\max}, C_{\max}, \alpha_{\min}) \frac{\log m}{m} \\
& \quad + C_3^R(n, R_{\max}, \zeta_{\max}, |\Theta_{\text{RRR}}|, C_{\max}) \sqrt{\frac{\log m}{m}} \\
& \quad + C_4^{a,R}(n, R_{\max}, \zeta_{\max}, \alpha_{\min}, r_{0,1}) \sqrt{\frac{\log 4/\delta}{2m}}
\end{align*}
$$

if $r_0 = r_1 = r_{0,1}$, and

$$
\begin{align*}
\sup_{H \in \mathcal{H}_{n, R_{\max}, C_{\max}}^{\text{RRR}}}
& \left| R(H) - \hat{R}_m(H) \right| \\
& \le \left( C_0^{b,R}(n, R_{\max}, \epsilon_{\text{RRR}}) + C_1(n, R_{\max}, \zeta_{\max}) \right) \frac{1}{m} \\
& \quad + C_2^R(n, R_{\max}, \zeta_{\max}, C_{\max}, \alpha_{\min}) \frac{\log m}{m} \\
& \quad + C_3^R(n, R_{\max}, \zeta_{\max}, |\Theta_{\text{RRR}}|, C_{\max}) \sqrt{\frac{\log m}{m}} \\
& \quad + C_4^{b,R}(n, R_{\max}, \zeta_{\max}, \alpha_{\min}, \epsilon_{\text{RRR}}) \sqrt{\frac{\log 4/\delta}{2m}}
\end{align*}
$$

if $r_0 > r_1$, where
$$C_0^{a,R}(n, R_{\max}, \alpha_{\min}, r_{0,1}) = 2L_{\ell} \bar{L}_h^{\text{poly}} \frac{(1 - (1 - \alpha_{\min})^m r_{0,1}^m) (1 - \alpha_{\min}) r_{0,1}}{1 - (1 - \alpha_{\min}) r_{0,1}}$$

$$C_0^{b,R}(n, R_{\max}, \epsilon_{RRR}) = 2L_{\ell} \bar{L}_h^{\text{poly}} \frac{(1 - (1 - \epsilon_{RRR})^m) (1 - \epsilon_{RRR})}{\epsilon_{RRR}}$$

$$C_2^R(n, R_{\max}, \zeta_{\max}, C_{\max}, \alpha_{\min})$$

$$= \frac{2 \left( L_{\ell} \bar{L}_h^{\text{poly}} + \mathbb{E}[|\ell(0, Y_0)|] + C_{\max} L_{\ell} \right)}{\log \zeta_{\max}^{-1}} + \frac{2(1 - \alpha_{\min}) L_{\ell} \bar{L}_h^{\text{poly}} C_v}{\zeta_{\max} \log \zeta_{\max}^{-1}}$$

$$C_3^R(n, R_{\max}, \zeta_{\max}, |\Theta_{RRR}|, C_{\max}) = \frac{2L_{\ell} \left( |\Theta_{RRR}| \cdot R_{\max} \cdot \left( \binom{n+R_{\max}}{R_{\max}} - 1 \right) + C_{\max} \right)}{\sqrt{\log \zeta_{\max}^{-1}}}$$

$$C_4^{a,R}(n, R_{\max}, \zeta_{\max}, \alpha_{\min}, r_{0,1})$$

$$= 2L_{\ell} \frac{\bar{L}_h^{\text{poly}}}{1 - (1 - \alpha_{\min}) r_{0,1}} \left( (1 - \alpha_{\min}) r_{0,1} + 2(1 - \alpha_{\min}) M_{\xi} L_v \|w^v\|_1 \right) + M_{\xi} L_y \|w^y\|_1$$

$$C_4^{b,R}(n, R_{\max}, \zeta_{\max}, \alpha_{\min}, \epsilon_{RRR})$$

$$= 2L_{\ell} \frac{\bar{L}_h^{\text{poly}}}{\epsilon_{RRR}} (1 - \epsilon_{RRR} + 2(1 - \alpha_{\min}) M_{\xi} L_v \|w^v\|_1) + M_{\xi} L_y \|w^y\|_1$$

and $\bar{L}_h^{\text{poly}}$ is as in (7), i.e. $\bar{L}_h^{\text{poly}} = R_{\max} \cdot n\sqrt{2^n} \left( \binom{n+R_{\max}}{R_{\max}} - 1 \right)$, and the constant $C_1$ is as in (23).

The Big-O bounds follows straightforwardly by considering the parameters that do not explicitly depend on the reservoir choice as constants (i.e. constants related to the Input-Output distribution as well as the choice of loss function).
References

Eric R. Anschuetz and Bobak T. Kiani. Quantum variational algorithms are swamped with traps. *Nature Communications*, 13(1), 2022.

Kaifeng Bu, Dax Enshan Koh, Lu Li, Qingxian Luo, and Yaobo Zhang. On the statistical complexity of quantum circuits. *Physical Review A*, 105(6), 2022.

Matthias C. Caro, Hsin-Yuan Huang, M. Cerezo, Kunal Sharma, Andrew Sornborger, Lukasz Cincio, and Patrick J. Coles. Generalization in quantum machine learning from few training data. *Nature Communications*, 13(1), 2022.

Jiayin Chen and Hendra I. Nurdin. Learning nonlinear input-output maps with dissipative quantum systems. *Quantum Information Processing*, 18(7), 2019.

Jiayin Chen, Hendra I. Nurdin, and Naoki Yamamoto. Temporal information processing on noisy quantum computers. *Physical Review Applied*, 14(2), 2020.

Jérôme Dedecker, Paul Doukhan, Gabreil Land, José Rafael Leon R., Sana Louhichi, and Clémentine Prieur. *Weak Dependence With Examples and Applications*. Springer New York, NY, 2007.

K. Eriksson, D. Estep, and C. Johnston. *Applied Mathematics: Body and Soul*. Springer, 2003.

William Feller. *An Introduction to Probability Theory and Its Applications*. John Wiley and Sons, Inc., 1968.

Keisuke Fujii and Kohei Nakajima. Harnessing disordered ensemble quantum dynamics for machine learning. *Physical Review Applied*, 8(2), 2016.

William Gilpin. Model scale versus domain knowledge in statistical forecasting of chaotic systems. *Physical Review Research*, 5(4), 2023.

Lukas Gonon and Antoine Jacquier. Universal approximation theorem and error bounds for quantum neural networks and quantum reservoirs. *arXiv preprint*, 2023.

Lukas Gonon, Lyudmila Grigoryeva, and Juan-Pablo Ortega. Risk bounds for reservoir computing. *Journal of Machine Learning Research*, 21(240), 2020.

Lyudmila Grigoryeva and Juan-Pablo Ortega. Echo state networks are universal. *Neural Networks*, 108, 2018a.

Lyudmila Grigoryeva and Juan-Pablo Ortega. Universal discrete-time reservoir computers with stochastic inputs and linear readouts using non-homogeneous state-affine systems. *Journal of Machine Learning Researchs*, 19, 2018b.

Herbert Jaeger. The “echo state” approach to analysing and training recurrent neural networks. *GMD Report, German National Research Institute for Computer Science*, 148, 2001.
Rodrigo Martinez-Pena and Juan-Pablo Ortega. Quantum reservoir computing in finite dimensions. *Physical Review E*, 107(3), 2023.

Jarrod R. McClean, Sergio Boixo, Vadim N. Smelyanskiy, Ryan Babbush, and Hartmut Neven. Barren plateaus in quantum neural network training landscape. *Nature Communications*, 9(1), 2018.

Francesco Monzania and Enrico Pratia. Universality conditions of unified classical and quantum reservoir computing. *arXiv preprint*, 2024.

Kohei Nakajima. Physical reservoir computing—an introductory perspective. *Japanese Journal of Applied Physics*, 59(6), 2020.

Alexey E. Rastegin. Relations for certain symmetric norms and anti-norms before and after partial trace. *Journal of Statistical Physics*, 148(6), 2012.

Jonathan Romero, Jonathan P. Olson, and Alan Aspuru-Guzik. Quantum autoencoders for efficient compression of quantum data. *Quantum Science and Technology*, 5(4), 2017.

Jonathan Romero, Jonathan P. Olson, and Alan Aspuru-Guzik. A generative modeling approach for benchmarking and training shallow quantum circuits. *npj Quantum Information*, 5(1), 2019.

Antonio Sannia, Rodrigo Martinez-Pena, Miguel C. Soriano, Gian Luca Giorgi, and Roberta Zambrini. Dissipation as a resource for quantum reservoir computing. *Quantum*, 8, 2024.

Maria Schuld, Alex Bocharov, Krysta Svore, and Nathan Wiebe. Circuit-centric quantum classifiers. *Physical Review A*, 101(3), 2020.

Yudai Suzuki, Qi Gao, Ken C. Pradel, Kenji Yasuoka, and Naoki Yamamoto. Natural quantum reservoir computing for temporal information processing. *Scientific Reports*, 12 (1), 2022.

Gouhei Tanaka, Toshiyuki Yamane, Jean Benoit Héroux, Ryosho Nakane, Naoki Kanazawa, Seiji Takeda, Hidetoshi Numata, Daiju Nakano, and Akira Hirose. Recent advances in physical reservoir computing: A review. *Neural Networks*, 115, 2019.

Toshiki Yasuda, Yudai Suzuki, Tomoyuki Kubota, Kohei Nakajima, Qi Gao, Wenlong Zhang, Satoshi Shimono, Hendra I. Nurdin, and Naoki Yamamoto. Quantum reservoir computing with repeated measurements on superconducting devices. *arXiv preprint*, 2023.