# Provable superior accuracy in machine learned quantum models

Chengran Yang<EMAIL_ADDRESS>Nanyang Quantum Hub, School of Physical and
Mathematical Sciences,
Nanyang Technological University, Singapore 637371  Complexity Institute,
Nanyang Technological University, Singapore 637335  Centre for Quantum
Technologies, National University of Singapore, 3 Science Drive 2, Singapore
117543  Andrew J. P. Garner  Nanyang Quantum Hub, School of Physical and
Mathematical Sciences,
Nanyang Technological University, Singapore 637371  Complexity Institute,
Nanyang Technological University, Singapore 637335  Institute for Quantum
Optics and Quantum Information,
Austrian Academy of Sciences, Boltzmanngasse 3, Vienna 1090, Austria  Feiyang
Liu  Shenzhen Institute for Quantum Science and Engineering and Department of
Physics,
Southern University of Science and Technology, Shenzhen 518055, China  Nora
Tischler  Centre for Quantum Computation and Communication Technology,
Australian Research Council,
Centre for Quantum Dynamics, Griffith University, Brisbane, Queensland 4111,
Australia  Dahlem Center for Complex Quantum Systems, Freie Universität
Berlin, 14195 Berlin, Germany  Jayne Thompson  Centre for Quantum
Technologies, National University of Singapore, 3 Science Drive 2, Singapore
117543  Man-Hong Yung  Shenzhen Institute for Quantum Science and Engineering
and Department of Physics,
Southern University of Science and Technology, Shenzhen 518055, China  Mile Gu
<EMAIL_ADDRESS>Nanyang Quantum Hub, School of Physical and
Mathematical Sciences,
Nanyang Technological University, Singapore 637371  Complexity Institute,
Nanyang Technological University, Singapore 637335  Centre for Quantum
Technologies, National University of Singapore, 3 Science Drive 2, Singapore
117543  Oscar Dahlsten<EMAIL_ADDRESS>Shenzhen Institute for Quantum
Science and Engineering and Department of Physics,
Southern University of Science and Technology, Shenzhen 518055, China  Wolfson
College, University of Oxford, Oxford, OX2 6UD, United Kingdom

###### Abstract

In modelling complex processes, the potential past data that influence future
expectations are immense. Models that track all this data are not only
computationally wasteful but also shed little light on what past data most
influence the future. There is thus enormous interest in dimensional reduction
– finding automated means to reduce the memory dimension of our models while
minimizing its impact on its predictive accuracy. Here we construct
dimensionally reduced quantum models by machine learning methods that can
achieve greater accuracy than provably optimal classical counterparts. We
demonstrate this advantage on present-day quantum computing hardware. Our
algorithm works directly off classical time-series data and can thus be
deployed in real-world settings. These techniques illustrate the immediate
relevance of quantum technologies to time-series analysis and offer a rare
instance where the resulting quantum advantage can be provably established.

## Introduction

Predictive modelling permeates quantitative science, from the rise and fall of
stock markets to the dynamics of neural spike trains. The making of
quantitative predictions often relies on learning the potential causal
dynamics that underlie vast amounts of observational data. In this vein, there
is a significant interest to explain more with less – enabling us to pinpoint
key observations from the past that have maximal influence on future
behaviour. This motivates efficient models, which generate faithful future
predictions while retaining only a small fraction of information from the
past. The benefits of such models are two-fold. Firstly, they help us minimize
the memory we need to make predictions about the future. Secondly, they enable
us to sample a wide array of possible futures by adjusting a few key initial
parameters. Indeed, much of machine learning, from dimensional reduction and
feature detection to auto-encoders, adheres to such philosophy [1, 2, 3].

Figure 1:  Scheme for our quantum causal-inference algorithm. The algorithm
takes as input time-series data drawn from some stochastic process, together
with a specification of the number of dimensions of available quantum memory
(the number of mutually distinguishable quantum states). The algorithm
identifies a quantum machine that best predicts future outcomes of the
process, subject to such memory constraints. The resulting machines can
achieve significantly greater accuracy than classical counterparts.

Even efficient models can be highly memory intensive. Complexity scientists in
computational mechanics, for example, have extensively studied the inference
of predictive models for stochastic processes [4, 5, 6, 7]. These processes
govern the behaviour of a system observed at discrete points in time. The goal
is then to infer a systematic method of encoding past observations into some
machine, whose dynamics allow it to sequentially generate statistically
faithful conditional future behaviour at all future time-steps. Sophisticated
techniques have been developed for realistic settings where predictive models
with memory constraints are trained on finite time-series data – techniques
that have successfully inferred structure in diverse contexts, from
distinguishing states of brain activity to causal inference in geomagnetic
time-series [8, 9, 10]. Meanwhile, theory indicated that as processes become
progressively non-Markovian, the number of distinct configurations a machine
needs to store relevant past information can grow without bound [11, 12].
Thus, memory constraints on such models invariably lead to distortions (i.e.
statistical errors) in predicted futures.

Here we ask: _Can memory-limited quantum models – quantum machines armed with
quantum memories – reduce such distortions?_ To answer this, we develop a
quantum model discovery algorithm that infers such models directly from raw
data (see Figure 1). Applying this algorithm, we show that our inferred
quantum models can achieve accuracies in prediction impossible in provably
optimal classical limits and illustrate that these advantages can grow with
the non-Markovianity of the processes under study. These results complement
recent studies showing how quantum mechanics can reduce memory requirements
for the idealistic limiting case of exact modelling [13, 14, 15, 16, 17, 18,
19, 11], which left an open question on whether quantum advantages prevail in
more realistic memory-limited scenarios where there is a trade-off between
memory and accuracy. Our results prove that the answer is in the affirmative,
and additionally give an algorithmic means to infer such dimensionality-
reduced quantum models directly from experimental data. We implement a model
found by our algorithm on the IBM quantum computer (“ibmq athens”). Despite
the noise, a statistically significant accuracy advantage remains.

## Results

Figure 2: Comparison of classical and quantum models that produce statistical
samples of a stochastic process’s future behaviour. A classical model stores
classical states in its memory. In quantum models, the quantum memory and the
output register – initialized in state $|0\rangle$ – interact through the
unitary transformation $U$, and a measurement of the output register in the
computational basis yields $x_{i}$. This can be repeated to sample $x_{0:L}$
for as big $L$ as desired. $U(L)$ denotes the unitary coupling the memory with
$L$ consecutive output registers.

Classical models and distortion – Consider a physical system observed at
discrete times $t$ with corresponding outcomes $x_{t}\in\mathcal{X}$ over some
finite alphabet $\mathcal{X}$. We assume that the system is stationary and
stochastic, such that each $x_{t}$ is drawn from random variable $X_{t}$, all
of which are governed by a time-translation invariant distribution
$P(\overleftarrow{X},\overrightarrow{X})$. Here, $\overleftarrow{X}:=\cdots
X_{-2}X_{-1}$ and $\overrightarrow{X}:=X_{0}X_{1}X_{2}\cdots$ denote random
variables governing stochastic outcomes in the past and future.

An (exact) predictive model is a systematic algorithm that takes each possible
past $\overleftarrow{x}$, encodes it via the deterministic map $\mathcal{E}$
to some suitable state $\mathcal{E}(\overleftarrow{x})$ within a physical
memory $M$, such that the same systematic action $\mathcal{P}$ on $M$ at each
subsequent time-step sequentially outputs $x_{0},x_{1},x_{2}\cdots$ with
probability $P(\overrightarrow{X}=\overrightarrow{x}|\overleftarrow{x})$, as
shown in Figure 2. Knowing the state of the memory $M$ is thus as useful as
$\overleftarrow{x}$ for purposes of inferring future statistics. Together,
these properties moreover imply that the model is _unifilar_ – if a predictive
model in state $r_{\overleftarrow{x}}=\mathcal{E}(\overleftarrow{x})$ emits
some output $x_{0}$, we know with certainty that the model will transition to
state $r_{\overleftarrow{x},x_{0}}=\mathcal{E}(\overleftarrow{x},x_{0})$ [5].
This guarantees that regardless of how many steps $L$ of the process are
produced, it is always possible to work out what state we should update $M$
to, such that the model can continue the simulation process.

The classical predictive models of provably minimal memory dimension are
$\varepsilon$-machines [5, 20]. $\varepsilon$-machines have been extensively
studied in the context of complex systems [8, 21, 22, 23, 24]. They capitalize
on the fact that if two pasts, $\overleftarrow{x}$ and
$\overleftarrow{x}^{\prime}$ have coinciding future statistics, then it is
unnecessary to distinguish them for purposes of future prediction. This then
divides the set of all pasts into equivalence classes, according to the
relation $\overleftarrow{x}\sim\overleftarrow{x}^{\prime}$ if and only if
$P(\overrightarrow{X}|\overleftarrow{x})=P(\overrightarrow{X}|\overleftarrow{x}^{{}^{\prime}})$.
The $\varepsilon$-machine then allocates in memory one configuration state for
each equivalence class. These states are known as _causal states_ , denoted by
$\mathcal{E}(\overleftarrow{x})$.

The number of causal states, $d_{c}$, is known as a process’s _topological
complexity_ – capturing the minimal number of configurations a predictive
model requires, whilst exhibiting correct conditional future behaviour. For
Markovian processes, for example, the map
$\mathcal{E}({\overleftarrow{x}})=x_{-1}$ is predictive, since knowing the
very last output is sufficient for purposes of prediction without distortion.
Thus for Markovian processes, $d_{c}$ needs never exceed the size of the
alphabet $\mathcal{X}$. When a process becomes highly non-Markovian, its
topological complexity can scale without bound [25]. The simplest example is a
discrete-time clock, which ‘ticks’ by emitting ‘$1$’ every $T$ time-steps and
‘$0$’ at all other times.

In such scenarios of diverging $d_{c}$, finite memory means a necessary
distortion of predictive accuracy. Specifically, given the constraint that our
memory $M$ has a fixed number of dimensions $\hat{d}_{c}\leq d_{c}$, we may
define a candidate model as one that encodes $\mathcal{E}(\overleftarrow{x})$
within $M$ – such that systematic action $\mathcal{P}$ on $M$ generates a
sample from
$\hat{P}\left(\overrightarrow{X}|\mathcal{E}(\overleftarrow{x})\right)$.

Our goal is then to find a model such that the deviation of $\hat{P}$ from
$P$– the distortion–is minimal. There are a number of possible ways to
quantify this distortion. Here, we make use of the normalized conditional
Kullback–Leibler (KL) divergence

$D_{L}(P,\hat{P}|\overleftarrow{x}):=\frac{1}{L}\mathcal{D}_{\rm
KL}(P(X_{0:L}|\overleftarrow{x})~{}||~{}\hat{P}(X_{0:L}|\mathcal{E}(\overleftarrow{x}))),$
(1)

where $X_{0:L}=X_{0}X_{1}\ldots X_{L}$ and $\mathcal{D}_{\rm
KL}(P||Q)=\sum_{x}P(x)\log(P(x)/Q(x))$. The KL divergence holds operational
significance in diverse contexts including machine learning, hypothesis
testing and thermodynamics [26, 27]. The distortion of a model $\hat{P}$
relative to the process $P$ is then quantified as:

$D_{e}(P,\hat{P})=\lim_{L\to\infty}\sum_{\overleftarrow{x}}P(\overleftarrow{x})\,D_{L}(P,\hat{P}|\overleftarrow{x}),$
(2)

i.e. the mean KL divergence over all pasts.

Quantum Models – In conventional models, the memory $M$ is assumed to be
classical. Quantum mechanics, however, enables memories configured in
superposition states. In this scenario, a $d$-dimensional quantum model is one
where each past $\overleftarrow{x}$ is encoded into a suitable quantum state
$|\sigma_{\overleftarrow{x}}\rangle$ within a $d$\- dimensional quantum memory
$M$. Meanwhile, there must exist some quantum process whose repeated
application allows sampling of the appropriate conditional futures. Figure 2
illustrates the quantum circuit picture. At each time-step, $M$ is coupled via
a time-step independent unitary $U$ to an output register initially in
$|0\rangle$. The output register is then measured, resulting in emitted
outcomes $x_{1},x_{2},\ldots$ governed by
$\hat{P}(\overrightarrow{X}|\sigma_{\overleftarrow{x}})$.

There are processes where quantum models can achieve zero distortion (i.e.
$\hat{P}(\overrightarrow{X}|\sigma_{\overleftarrow{x}})=P(\overrightarrow{X}|\overleftarrow{x})$)
with a lower memory dimension than classically possible [19]. Experiments
demonstrate that stochastic processes with three causal states can be
accurately modelled by encoding the past within a single qubit [28]. Thus the
quantum topological complexity $d_{q}$ can be strictly lower than its
classical counterpart. However, systematic means to construct such models
require optimisation over an exponentially growing parameter space [29, 18].

Zero-distortion is an idealized limit. In practice, we are unlikely to know
the exact process $P$ we wish to model. Instead, we are given a finite
sequence of observations drawn from this process. In such scenarios,
statistical fluctuations imply that aiming for zero distortion is likely
meaningless and often represents harmful overfitting. Typically, we start with
general finite and fixed memory resources. Thus, minimizing distortion when
given finite memory $\hat{d}_{q}<d_{q}$ and finite data is of both of
conceptual and practical importance. In such settings, the capabilities of
quantum models to minimize this distortion are not known. Indeed, the task
remains under active study even in classical settings [30, 31, 32]. As such,
our goal here is to find an algorithm for inferring quantum models of some
given dimension $\hat{d}_{q}$ that attempts to minimize distortion, and use it
to establish that such models exhibit less distortion than any classical
counterpart.

Quantum model discovery algorithm – (see Box 1). The algorithm takes two
inputs: (i) a data sequence $x_{0:L}$, assumed to be drawn from some
stationary stochastic process $P$ and (ii) our desired memory dimension
$\hat{d}_{q}$. From this data, the algorithm learns a quantum model of the
process, including (i) an encoding map $\mathcal{E}$ which specifies how to
encode each past into a quantum memory of dimension $\hat{d}_{q}$, and (ii)
the relevant physical process – described by a unitary $U$ as in Figure 2 –
whose repeated application produces entangled output states. Measurements on
such output states emit outputs that approximate the process’s conditional
future behaviour.

Box 1: Quantum model discovery algorithm.  Inputs:
$\hat{d}_{q}$ – the desired model memory dimension. $x_{0:L}$ – a length $L$
sequence from the process. Outputs:
$U$ – the unitary dynamics of the quantum model. $\mathcal{E}$ – the encoding
map from observed histories to $\hat{d}_{q}$–dimensional quantum memory
states. Algorithm: 1. To learn the unitary operator $U$: (a) Randomly
initialize parameter set $A$ corresponding up to completeness to the set of
Kraus operator matrices $\\{A_{x}:=\langle x|U|0\rangle\\}$. (b) Evaluate the
cost function $C:=-\log P(x_{0:L};A)$ and its gradient $\nabla C$. (c) Update
the Kraus operators $\\{A_{x}\\}$ using gradient-descent based method (such as
Adam). (d) Repeat (b)-(c) until the cost function decrease is sufficiently
small. (e) Save final $A$ as candidate parameters. (f) Repeat (a)-(e) as
desired, to minimize the impact of initial choice of $A$. Choose the final $A$
with the lowest C. (g) For that $A$, recover the completeness relation for
$\\{A_{x}\\}$. (h) Construct a unitary operator $U$ from these $\\{A_{x}\\}$.
2. To compute the encoding $\mathcal{E}$: (a) Compute the leading eigenvector
$|\sigma_{0}\rangle$ of one of the Kraus operators ($A_{0}$). (b)
$\mathcal{E}$ is then defined as
$\mathcal{E}\\!\left(\overleftarrow{x}\right):=A_{\overleftarrow{x}}|\sigma_{0}\rangle/\lVert
A_{\overleftarrow{x}}|\sigma_{0}\rangle\rVert_{2}$.

To operate the learning algorithm, we cast the problem as minimizing some cost
function over a parameter set $B$. This involves finding (a) an effective
parametrisation of all such models to optimize over and (b) a computable cost
function that proxies the mean KL-divergence $D_{e}$.

To tackle (a), first note that our model’s output behaviour is entirely
defined by a family of Kraus operators $A=\\{A_{x}\\}$, where $A_{x}=\langle
x|U|0\rangle$ captures the action of the memory when it interacts with an
output register that is subsequently observed to have output $x$. Observe also
that $A_{x}$ informs the encoding map from pasts onto quantum memory states.
To see this, note that a memory initially in $\rho$ is mapped to $A_{x}\rho
A_{x}^{\dagger}$ after observation of $x$. Thus, observation of a past
$\overleftarrow{x}$ implies applying a sequence of Kraus operators
$A_{\overleftarrow{x}}:=A_{x_{0}}A_{x_{-1}}\ldots$. We can therefore encode
the $\overleftarrow{x}$ as

$\mathcal{E}\\!\left(\overleftarrow{x}\right):=|\sigma_{\overleftarrow{x}}\rangle:=cA_{\overleftarrow{x}}|\sigma_{0}\rangle,$
(3)

where $c$ is the constant of normalization. Here, the initial quantum state
$|\sigma_{0}\rangle$ can be the leading eigenvector of any Kraus operator
$A_{x}$ (see Appendix C). This naïvely suggests optimizing over $A$, but this
is generally cumbersome as $A$ is constrained by the completeness relation
$\sum_{x}(A_{x})^{\dagger}A_{x}=\mathbb{I}$. In the Methods, we instead devise
a way to optimize over a set of unconstrained $\hat{d}_{q}\times\hat{d}_{q}$
complex matrices $B=\\{B^{x}\\}$, whose value enables the deduction of $A$ via
tensor network techniques [33].

To tackle (b), the obvious candidate cost function is $D_{e}$ itself. However,
computation of this quantity is costly – especially as we are given only a
sample $x_{0:L}$ of the desired process $P$. We thus employ a more efficiently
computable proxy,

$C(B)=-\log_{2}P_{B}(x_{0:L}),$ (4)

representing the negative log-likelihood of producing $x_{0:L}$ with a model
parameterized by $B$. In the limit of $L\to\infty$, a model minimising the
negative log-likelihood of Eq. (4) exactly replicates the desired stochastic
behavior [34, 35]. With both the cost function and parameter spaces
established, the use of optimisation algorithms such as Adam is enabled (see
the Methods). In the Methods, we also show that $C$ can be computed directly
from $B$ without first deducing $A$, boosting optimisation efficiency. Once
optimisation concludes, the corresponding $U$ for generating predictions is
retrieved from $B$ (see the Methods).

(a)

(b)

Figure 3: Example: The asymmetric process. (a) $\varepsilon$-machine for the
asymmetric process, parametrised by $p,q\in(0,1)$. The nodes are the causal
states $s_{0}$, $s_{1}$ and $s_{2}$. The arrows indicate transitions between
states, where each label $x|P(x)$ indicates output $x$ is produced with
probability $P(x)$. (b) Classical vs quantum distortion. Distortion in the
$2$-dimensional models for the asymmetric process across a range of process
parameters $p$, $q$. The left-hand-side shows the optimal approximate
classical model’s distortion; the right-hand-side shows the distortion of the
most accurate quantum model found by our algorithm.

Provable Quantum Advantage – We are now in a position to answer the open
question: _can memory-constrained quantum models exhibit improved accuracy_?
Answering this question in full generality is non-trivial as almost all
existing analytical results for predictive modelling apply only to zero
distortion. While sophisticated numerical methods for inferring classical
models exist [36, 7], they do not tell us whether the resulting models are
provably optimal in minimizing KL-divergence. In the Methods we establish
provable upper bounds on the quality of classical methods of fixed dimension.
We will show that our inferred quantum models can exceed these optimal bounds.

The first example is the family of simple asymmetric processes. For any
process in this family, exact modelling requires three causal states
classically (see its optimal predictor in Figure 3a), but can be achieved
quantumly on a qubit [17, 19]. Thus, if we set the desired dimensions to $2$,
our quantum model could theoretically achieve a KL-divergence of $0$, while
some distortion for classical models is inevitable. Indeed, in the methods, we
evaluate the minimal error achieved by $2$-dimensional classical models for
process parameters; and compare them with errors in our trained quantum models
(see Figure 3b). We observed that errors in our trained quantum models are
negligible in comparison to the best possible classical counterparts.

Our second example involves a family of generalised $3$-state quasi-cycles
(see Figure 4a). Quasi-cycles are for example of interest in
thermodynamics—certain quasi-cycles are valid thermal operations in that they
preserve the Gibb’s thermal state, but interestingly without respecting
detailed balance [37]. They constitute a family of processes where neither
quantum nor classical models can generally achieve perfect accuracy when
limited to $2$ dimensional memories [18]. The resulting error of our inferred
quantum model (see Appendix E) is compared with that of the classical limit in
Figure 4b. Clearly, a broad area of process parameters $(p,\delta)$ can be
accurately approximated by a compressed quantum model, whose error is up to an
order of magnitude better than the classical models. The algorithm also
achieves near perfect accuracy over a one-dimensional subspace of the quasi-
cycles – retrieving a class of quantum models that was only recently
discovered using exhaustive numerical search [18]).

(a)

(b)

Figure 4: Example: Quasi-cycle process. (a) $\varepsilon$-machine for the
$3$-state quasicycle, parametrised by $p\in(0,1),\delta\in(0,p]$. The nodes
are the causal states $s_{0}$, $s_{1}$ and $s_{2}$. The arrows indicate
transitions between states, where each label $x|P(x)$ indicates output $x$ is
produced with probability $P(x)$. (b) Classical vs quantum error comparison.
Error in the $2$-dimensional models for the quasi-cycle, across a range of
process parameters $p$, $\delta$. The left-hand-side shows the optimal
approximate classical model’s error; the right-hand-side shows the error of
the most accurate quantum model found by our algorithm.

We now apply our inference algorithm to an iconic non-Markovian process - the
family of discrete-time renewal processes. Renewal processes represent an
interpolation between Poissonian processes and ticking clocks, characterising
many natural processes including models of lifetimes [38] and neural spike
trains [39]. At each time-step, the process may emit one of two outputs: $0$,
representing no tick, and $1$ representing a tick, and their probability of
ticking at each time-step depends only on the number of time-steps since their
previous tick [11, 40, 12]. In the event that this probability is constant, we
have a Poissonian process whose topological complexity is zero; as its action
at each time is independent of past events (and thus considered memory-less).

Here we consider a class of period-$N$ uniform renewal processes
$\mathcal{U}_{N}$. These represent processes where the number of $0$s between
each two neighboring ticks is uniformly distributed between $0$ and $N-1$. As
$N$ scales, the process becomes progressively more non-Markovian. Indeed, the
classical topological complexity for $\mathcal{U}_{N}$ is $N$, the intuition
being that a maximally faithful model would need to store the number of time-
steps since the last tick to decide with what probability it should tick at
each time-step (see corresponding $\varepsilon$-machine in Figure 6a).
Reducing the number of distinct memory configurations to $\hat{d}<N$ would
force the model to lose track of this information – resulting in unavoidable
distortion.

To study the magnitude of this distortion, we applied our algorithm to infer
memory-constrained models of $\mathcal{U}_{N}$ for cases where $N$ ranged
between 2 and 5 inclusive, for a variety of memory limits $\hat{d}_{q}<N$. We
then compare the error with a lower bound of optimal classical models with
similar memory constraints (see the Methods). The resulting distortions are
illustrated in the upper graph of Figure 6b. Even though the learned quantum
models have the additional handicap of being trained on finite-length data, we
see that quantum models outperform optimal classical counterparts away from
the zero-error regime. Models with larger $\hat{d}_{q}$ or $\hat{d}_{c}$
maintain smaller error but at the cost of a larger memory, and requiring a
larger dataset to train. In the case of $N=4$, there is a significant gap
between the errors of $\hat{d}_{q}=3,4$ and the error of $\hat{d}_{q}=2$.
Thus, the algorithm could here be used to identify a sweet spot for the choice
of $\hat{d}_{q}$: setting $\hat{d}_{q}=3$ gives a good accuracy for a small
dimension. To examine the resilience of the quantum advantage on present-day
noisy devices, we implement the quantum model on one of IBM’s quantum systems
(“ibmq athens”) for the case of $\hat{d}_{q}=2,N=3$. Such quantum models still
exhibit statistically significant accuracy advantage over the classical bound,
as is shown in Figure 5.

Figure 5: Accuracy advantage on present-day quantum hardware. Despite the
experimental noise, the quantum model (orange) retains an accuracy advantage
of $2.12$ standard deviations over the classical bound (yellow) for the given
dimension. For comparison, the blue bar shows the expected distortion of the
learned quantum model without experimental noise. See Methods for details.

(a)

(b)

Figure 6: Example: The uniform renewal process. (a) $\varepsilon$-machine for
an $N$-step discretization of the uniform renewal process (in which events are
spaced by a uniform random interval). The nodes are the causal states $s_{0}$,
…$s_{N-1}$. The arrows indicate transitions between states, where each label
$x|P(x)$ indicates that output $x$ is produced with probability $P(x)$. (b)
Comparison of classical and quantum errors. The left-hand-side shows the lower
bound on the optimal approximate classical model’s error; the right-hand-side
shows the error of the most accurate quantum model found by our algorithm. The
accuracy of the resulting quantum model is restricted by finite data whereas
the classical models are assumed with full knowledge of the stochastic
processes.

## Discussion

Inferring models to predict future behaviour based on past observations is a
key task permeating quantitative science. Ideally, when two candidate models
require storing the same amount of data from the past, the one that offers
predictions with greater statistical fidelity is preferred. Here we
demonstrate that quantum models have a provable advantage, and present an
algorithm to infer such advantaged models directly from data. When given
training data in terms of a data sequence and a memory of limited dimension,
our algorithm is able to encode relevant past information into the memory and
discover a quantum process that can sequentially output predictions. The
accuracy of these predictions – as measured by their statistical fidelity to
the true process generating the training data – can greatly exceed that of
provably optimal classical counterparts. The accuracy advantage is realized on
an IBM quantum computer, illustrating potential near-term applicability. Our
work substantially expands upon prior work in memory-enhanced quantum models,
which focus on the memory advantages in the exact modelling of a given
stochastic process [17, 12, 19]. The results here illustrate that when the
memory dimension is fixed, quantum models can exhibit a significant accuracy
advantage. This establishes the advantage of quantum models in real-world
settings, where there is a trade-off between model complexity and accuracy. In
this context, our results provide a quantum counterpart to classical tools for
structural inference [6, 41, 7] and dimensional reduction in stochastic
modelling [42].

The resulting circuits for executing such models simply involve sequential
application of the same unitary process $U$ at each time-step. The
dimensionality of this unitary is then bounded by the memory dimension of the
quantum model. As such, a fixed memory dimension guarantees that the
simulation of a process for $T$ time-steps has gate complexity that scales
linearly with $T$. Moreover, unitarity means such quantum models can generate
a quantum superposition of all possible futures if the outputs remain
unmeasured. Thus the quantum model produces a key resource for various quantum
data analytics algorithms, such as amplitude estimation, Grover’s search,
value at risk, and importance sampling [43, 44, 45, 46]. The capacity for our
models to generate data by sequential application of the same unitary makes
them particularly amenable to certain contemporary hardware architectures
(e.g., loop-based photonics where only a single optical circuit implementing
$U$ needs to be engineered [47]). Our techniques thus pave the way for
quantum-enhanced tools to analyse time-series data from diverse fields such as
financial forecasting [22, 21], neuroscience [8, 10], seismology [48, 49], and
natural language processing [50].

There are a number of promising avenues for future research. One caveat of our
present learning algorithms is that they remain classical. They can thus
suffer the curse of dimensionality, and perform best when studying processes
that are close to those with low quantum topological complexity. While some
processes within this class have an arbitrary number of causal states and are
already sufficient to demonstrate scaling advantage [51], there is a clear
motivation to develop alternative algorithms to learn quantum models of higher
memory dimension. One possibility is to adopt variational circuit techniques,
where candidate quantum models are parameterized by repeated application of
some low-depth circuit. One can then leverage potential interferometric means
to efficiently ascertain distortion in candidate models [52]. There are also a
number of connections to neighbouring fields. One is the dimensional
advantages of quantum clocks [53, 54], as a clock can be considered a
continuum limit of a particular discrete renewal process. Meanwhile prediction
also plays a key role in information engines that extract free energy from
structured patterns [55, 24, 56, 57, 58], where the KL divergence between
expectation and reality induces unavoidable energetic waste [59, 60]. Thus our
work may herald more effective quantum engines in memory-limited regimes, even
when extracting free energy from purely classical correlations.

## Methods

Model Parametrization – The candidate quantum models of dimension
$\hat{d}_{q}$ for a process with output alphabet $\mathcal{X}$ can be
parametrised by a complete set of $|\mathcal{X}|$ Kraus operators
$A=\\{A_{x}\\}_{x\in\mathcal{X}}$, where $A_{x}$ denotes the Kraus operator
describing how the model updates upon emission of output $x$. However, this
parametrisation is non-ideal for optimisation due to the completeness
constraint.

Here, we demonstrate an alternative parametrisation using
$B=\\{B_{x}\\}_{x\in\mathcal{X}}$, where each $B_{x}$ is a general
$\hat{d}_{q}\times\hat{d}_{q}$ complex matrix. Notably, we show that given any
such $B$, it is possible to recover a corresponding set $A$ via the following
process: First, consider the linear map
$\mathcal{E}^{B^{*}}\\!\left(\cdot\right):=\sum_{x}B_{x}^{\dagger}\cdot
B_{x}$. By construction, this map is completely positive, so its leading
eigenvalue $\lambda$ is real and positive, and the associated eigenmatrix $V$
is positive semi-definite. Thus, $V$ admits a decomposition $V=W^{\dagger}W$
(where $W$ is invertible). We can then set $A_{x}$ such that

$A_{x}:=WB_{x}{W}^{-1}/\sqrt{\lambda}.$ (5)

Every set $\\{A_{x}\\}$ formed in this way satisfies the completeness relation
(see Appendix A). Thus, each $B_{x}$ can be used directly to construct a valid
set of corresponding Kraus operators $A_{x}$. This further results in valid
unitary operator $U$ (see Appendix B). Likewise, via Eq. (5) one can also
infer the encoding map $\mathcal{E}$ from $B$ (see Appendix C).

Computing the Cost Function – To see how to compute the cost function
$C(B)=-\log_{2}P_{B}(x_{0:L})$ directly from $B$, first consider a model with
Kraus operators $A=\\{A_{x}\\}$ initialized in state $\rho_{0}$ at $t=0$. The
probability it outputs $x$ at $t=1$ is then given by
$P(x|\rho_{0})|_{A}=\mathrm{Tr}(A_{x}\rho_{0}(A_{x})^{\dagger})$, whereby the
state transitions to $\rho_{1}=A_{x}\rho_{0}(A_{x})^{\dagger}$. As such,
repeated iterations of the model will output $x_{0:L}$ with probability

$P(x_{0:L}|\rho_{0})|_{A}=\mathrm{Tr}(A_{x_{0:L}}\rho_{0}A_{x_{0:L}}^{\dagger}),$
(6)

where $A_{x_{0:L}}=A_{x_{L-1}}\cdots A_{x_{0}}$. In addition, if $\rho_{0}$ is
the stationary memory state averaged over all histories of the process, we
obtain the probability of output sequence $x_{0:L}$ when averaged over all
pasts. We can then write this likelihood directly in terms of $B^{x}$ by
applying Eq. (5),

$P_{B}(x_{0:L})=\mathrm{Tr}(B_{x_{0:L}}\tilde{\rho}{B_{x_{0:L}}}^{\dagger}V)/\lambda^{L},$
(7)

where $\tilde{\rho}=W^{-1}\bar{\rho}{W^{-1}}^{\dagger}$ is the leading
eigenmatrix of $\mathcal{E}^{B}$.

Training process – Here we specify the details of the training process to
minimise the cost function $C(B)=-\log_{2}P_{B}(x_{0:L})$. At this stage, any
number of different optimisation techniques could be employed. Specifically,
we used Adam optimisation [61].

The Adam method is a sophisticated form of gradient descent. Recall that
standard gradient descent involves computing $\nabla C$, the partial
derivatives with respect to each degree of freedom in $B$ (referred to as a
free parameter), and updating $B$ according to $B^{\prime}=B-\eta\nabla C$,
where $\eta>0$ is an adjustable learning rate. The Adam method fine-tunes this
by using individual adaptive learning rates for each free parameter, computed
using estimates of the first and second moments of the gradients $\nabla C$
[61]. Due to the capacity to relate the cost function directly to $B$ (see Eq.
(7)), this can be done without ever needing to distill the quantum circuit
itself.

In each run of the training process, we begin with a random $B$. The update
process is then repeated until the decrease of updated cost function is below
a chosen threshold. Like nearly all gradient descent variants, Adam can
converge to local minima. However, in our experiments we found roughly three
quarters of all runs give the same minimal value for $C$, and no further gains
were made after taking the optimal of three runs (see Appendix D and E).

Box 2: Finding a lower bound on distortion
of unifilar models Inputs:
$P$ – the process $P$ with Markov order $\kappa$ $\hat{d}_{c}$ – the dimension
of the classical model. Outputs:
A lower bound on distortion $D_{e}$ of all predictive models. Algorithm: (a)
Pick one of the encoding $\mathcal{E}$ merging causal states of process $P$.
(b) Optimize the future statistics
$\hat{P}(X_{0:\kappa}|\mathcal{E}(\overleftarrow{x}))$ to minimize
$D_{\kappa}(P,\hat{P})=\sum_{\overleftarrow{x}}P(\overleftarrow{x})D_{\kappa}(P,\hat{P}|\overleftarrow{x}).$
(c) Pick a new encoding $\mathcal{E}$ that merges causal states and repeat
(a)-(c) until having gone through all encoding that merge causal states. (d)
Return the minimal $D_{\kappa}(P,\hat{P})$ among all encodings.

The result of this algorithm for the examples discussed in this article is
given in Appendix G.

Lower bound on distortion of classical models – Here we outline the details of
evaluating a lower bound on distortion of classical unifilar models with
finite dimension $\hat{d}_{c}$. We prove that a bound on the minimum
distortion of approximate unifilar models can be obtained by merging causal
states, and then optimizing over the statistics of the next $\kappa$ outputs
for a process with Markov order $\kappa$ (see Appendix F), as summarized in
Box 2. These analytical results reduce the search for bounds on the optimal
classical model from an infinite–dimensional problem to a finite–dimensional
search, which in general can be exhaustively completed numerically.
Furthermore, the bound becomes tight when $\kappa=1$ (Theorem 2 of Appendix
F).

Implementation of quantum models for the discrete renewal process on IBM’s
devices – To investigate the real-world performance of the models found by the
algorithm, we implement one of the discovered quantum models for the discrete
renewal process on IBM’s superconducting cloud quantum computer. We consider
the case of $\hat{d}_{q}=2,N=3$ and decompose the quantum model into a quantum
circuit consisting of single qubit gates and CNOT gates, as shown in Figure 7.
The IBM device used was “ibmq athens”. We ran the circuit $40,000$ times for
each input quantum state and obtained the output probabilities by measuring
the output register (see Appendix H).

Figure 7: Circuit diagram for a single step output. The upper line represents
the memory qubit while the lower line represents the output register. The
circuit involves $8$ single qubit gates ($U(\theta_{i})$) and $3$ C-NOT gates.
Each single qubit gate is specified by three real parameters. $U(\theta_{1})$
encodes the input memory state.

Figure 5 compares the distortion realized on the IBM device with the optimal
classical bound for $\hat{d}_{c}=2$ and the ideal realization of the learned
quantum model. The 40,000 runs were divided into 50 batches, and the
distortion was calculated for each batch. The error bar depicts the standard
deviation in distortion of those 50 batches. Our quantum model maintains a
statistically significant accuracy advantage even in a noisy environment:
there is a $2.12$ standard deviations gap between the experimental quantum
distortion and our lower bound on the best classical distortion.

## Acknowledgements

We are grateful for discussions with Thomas Elliott, Jim Crutchfield and
Kelvin Onggadinata. This research is supported by the National Research
Foundation (NRF), Singapore, under its NRFF Fellow programme (Award No. NRF-
NRFF2016-02), the Singapore Ministry of Education Tier 1 Grant No.
2018-T1-002-043 (RG162/19) and Grant No. FQXI R-710-000-146-720 (Are quantum
agents more energetically efficient at making predictions?) from the
Foundational Questions Institute and Fetzer Franklin Fund (a donor–advised
fund of Silicon Valley Community Foundation), the Quantum Engineering Program
QEP-SP3. N.T. acknowledges support by the Griffith University Postdoctoral
Fellowship Scheme.

## References

  * [1] Laurens Van Der Maaten, Eric Postma, and Jaap Van den Herik.  Dimensionality reduction: a comparative.  J Mach Learn Res, 10(66-71):13, 2009.
  * [2] Ethem Alpaydin.  Introduction to machine learning.  MIT press, 2020.
  * [3] Andrew Ng et al.  Sparse autoencoder.  CS294A Lecture notes, 72(2011):1–19, 2011.
  * [4] Judea Pearl et al.  Models, reasoning and inference.  Cambridge, UK: CambridgeUniversityPress, 2000.
  * [5] Cosma Rohilla Shalizi and James P. Crutchfield.  Computational mechanics: Pattern and prediction, structure and simplicity.  Journal of Statistical Physics, 104(3-4):817–879, 2001.
  * [6] Cosma Rohilla Shalizi, Kristina Lisa Shalizi, and James P Crutchfield.  An algorithm for pattern discovery in time series.  arXiv preprint cs/0210025, 2002.
  * [7] Cosma Shalizi and Kristina Lisa Klinkner.  Blind construction of optimal nonlinear recursive predictors for discrete sequences.  arXiv preprint arXiv:1408.2025, 2014.
  * [8] Robert Haslinger, Kristina Lisa Klinkner, and Cosma Rohilla Shalizi.  The Computational Structure of Spike Trains.  Neural Computation, 22(1):121–157, jan 2010.
  * [9] R. W. Clarke, M. P. Freeman, and N. W. Watkins.  Application of computational mechanics to the analysis of natural data: An example in geomagnetism.  Physical Review E, 67(1):016203, 2003.
  * [10] Roberto N Muñoz, Angus Leung, Aidan Zecevik, Felix A Pollock, Dror Cohen, Bruno van Swinderen, Naotsugu Tsuchiya, and Kavan Modi.  General anesthesia reduces complexity and temporal asymmetry of the informational structures derived from neural recordings in drosophila.  Physical Review Research, 2(2):023219, 2020.
  * [11] Thomas J. Elliott and Mile Gu.  Superior memory efficiency of quantum devices for the simulation of continuous-time stochastic processes.  npj Quantum Information, 4(1):18, 2018.
  * [12] T. J. Elliott, C. Yang, F. C. Binder, A. J. P. Garner, J. Thompson, and M. Gu.  Extreme Dimensionality Reduction with Quantum Modeling.  Physical Review Letters, 125(26):260501, dec 2020.
  * [13] Mile Gu, Karoline Wiesner, Elisabeth Rieper, and Vlatko Vedral.  Quantum mechanics can reduce the complexity of classical models.  Nature Communications, 3:762, 2012.
  * [14] C. Aghamohammadi, S. P. Loomis, J. R. Mahoney, and J. P. Crutchfield.  Extreme Quantum Advantage for Rare-Event Sampling.  Physical Review X, 8(1):11025, 2017.
  * [15] Felix C Binder, Jayne Thompson, and Mile Gu.  Practical Unitary Simulator for Non-Markovian Complex Processes.  Physical Review Letters, 120(24):240502, 2018.
  * [16] Chengran Yang, Felix C. Binder, Varun Narasimhachar, and Mile Gu.  Matrix product states for quantum stochastic modeling.  Phys. Rev. Lett., 121:260602, Dec 2018.
  * [17] Farzad Ghafari, Nora Tischler, Jayne Thompson, Mile Gu, Lynden K. Shalm, Varun B. Verma, Sae Woo Nam, Raj B. Patel, Howard M. Wiseman, and Geoff J. Pryde.  Dimensional quantum memory advantage in the simulation of stochastic processes.  Phys. Rev. X, 9:041013, Oct 2019.
  * [18] Qing Liu, Thomas J Elliott, Felix C Binder, Carlo Di Franco, and Mile Gu.  Optimal stochastic modeling with unitary quantum dynamics.  Physical Review A, 99(6):062110, 2019.
  * [19] Jayne Thompson, Andrew J. P. Garner, John R. Mahoney, James P. Crutchfield, Vlatko Vedral, and Mile Gu.  Causal Asymmetry in a Quantum World.  Physical Review X, 8(3):31013, 2018.
  * [20] James P. Crutchfield and Karl Young.  Inferring statistical complexity.  Physical Review Letters, 63(2):105–108, 1989.
  * [21] Jae Suk Yang, Wooseop Kwak, Taisei Kaizoji, and In Mook Kim.  Increasing market efficiency in the stock markets.  European Physical Journal B, 61(2):241–246, 2008.
  * [22] Joongwoo Brian Park, Jeong Won Lee, Jae Suk Yang, Hang Hyun Jo, and Hie Tae Moon.  Complexity analysis of the stock market.  Physica A: Statistical Mechanics and its Applications, 379(1):179–187, 2007.
  * [23] Alexander B. Boyd, Dibyendu Mandal, and James P. Crutchfield.  Identifying Functional Thermodynamics in Autonomous Maxwellian Ratchets.  New Journal of Physics, 18(2):023049, jul 2015.
  * [24] Andrew JP Garner, Jayne Thompson, Vlatko Vedral, and Mile Gu.  Thermodynamics of complexity and pattern manipulation.  Physical Review E, 95(4):042140, 2017.
  * [25] Sarah E Marzen and James P Crutchfield.  Informational and causal architecture of discrete-time renewal processes.  Entropy, 17(7):4891–4917, 2015.
  * [26] Art Hobson.  Concepts in statistical mechanics.  CRC Press, 1987.
  * [27] Thomas M Cover and Joy A Thomas.  Elements of information theory.  John Wiley & Sons, 2012.
  * [28] Farzad Ghafari, Nora Tischler, Jayne Thompson, Mile Gu, Lynden K. Shalm, Varun B. Verma, Sae Woo Nam, Raj B. Patel, Howard M. Wiseman, and Geoff J. Pryde.  Dimensional quantum memory advantage in the simulation of stochastic processes.  Phys. Rev. X, 9:041013, Oct 2019.
  * [29] Samuel P Loomis and James P Crutchfield.  Strong and weak optimizations in classical and quantum models of stochastic processes.  Journal of Statistical Physics, 176(6):1317–1342, 2019.
  * [30] Max Tegmark and Tailin Wu.  Pareto-optimal data compression for binary classification tasks.  Entropy, 22(1):7, 2020.
  * [31] Sarah E Marzen and James P Crutchfield.  Predictive rate-distortion for infinite-order markov processes.  Journal of Statistical Physics, 163(6):1312–1338, 2016.
  * [32] Sarah E Marzen and James P Crutchfield.  Nearly maximally predictive features and their dimensions.  Physical Review E, 95(5):051301, 2017.
  * [33] G. Vidal.  Entanglement renormalization.  Physical Review Letters, 99(22):220405, 2007.
  * [34] Ted Petrie.  Probabilistic functions of finite state markov chains.  The Annals of Mathematical Statistics, 40(1):97–115, 1969.
  * [35] B. Juang and L. R. Rabiner.  A probabilistic distance measure for hidden markov models.  AT T Technical Journal, 64(2):391–408, Feb 1985.
  * [36] David Kelly, Mark Dillingham, Andrew Hudson, and Karoline Wiesner.  A new method for inferring hidden markov models from noisy time sequences.  PloS one, 7(1):e29703, 2012.
  * [37] M. Horodecki and J. Oppenheim.  Fundamental limitations for quantum and nanoscale thermodynamics.  Nat Commun, (4):2059, 2013.
  * [38] Joseph L Doob.  Renewal theory from the point of view of the theory of probability.  Transactions of the American Mathematical Society, 63(3):422–438, 1948.
  * [39] James P Crutchfield, Michael Robert DeWeese, and Sarah E Marzen.  Time resolution dependence of information measures for spiking neurons: Scaling and universality.  Frontiers in computational neuroscience, 9:105, 2015.
  * [40] Sarah Marzen and James P Crutchfield.  Informational and Causal Architecture of Continuous-time Renewal Processes.  Journal of Statistical Physics, 168(1):109–127, 2017.
  * [41] Geoffrey E Hinton, Terrence Joseph Sejnowski, Tomaso A Poggio, et al.  Unsupervised learning: foundations of neural computation.  MIT press, 1999.
  * [42] L. van der Maaten, E. Postma, and J. van den Herik.  Dimensionality Reduction: A Comparative Review.  Technical Report, Tilburg University, 2007.
  * [43] Lov K Grover.  Quantum mechanics helps in searching for a needle in a haystack.  Physical review letters, 79(2):325, 1997.
  * [44] Gilles Brassard, Peter Hoyer, Michele Mosca, and Alain Tapp.  Quantum amplitude amplification and estimation.  Contemporary Mathematics, 305:53–74, 2002.
  * [45] Stefan Woerner and Daniel J Egger.  Quantum risk analysis.  npj Quantum Information, 5(1):1–8, 2019.
  * [46] Carsten Blank, Daniel K Park, and Francesco Petruccione.  Quantum-enhanced analysis of discrete stochastic processes.  arXiv preprint arXiv:2008.06443, 2020.
  * [47] A. Schreiber, K. N. Cassemiro, V. Potoček, A. Gábris, P. J. Mosley, E. Andersson, I. Jex, and Ch. Silberhorn.  Photons walking the line: A quantum walk with adjustable coin operations.  Phys. Rev. Lett., 104:050502, Feb 2010.
  * [48] G. Bertello, P. J. Arduin, F. Boschetti, and D. Weatherley.  Application of Computational Mechanics to the analysis of seismic time-series via numerical optimisation.  New Generation Computing, 27(1):1–23, nov 2008.
  * [49] T. Chelidze and T. Matcharashvili.  Complexity of seismic process; measuring and applications - A review.  Tectonophysics, 431(1-4):49–60, feb 2007.
  * [50] Lawrence R. Rabiner.  A Tutorial on Hidden Markov Models and Selected Applications in Speech Recognition.  Proceedings of the IEEE, 77(2):257–286, 1989.
  * [51] Thomas J Elliott, Chengran Yang, Felix C Binder, Andrew JP Garner, Jayne Thompson, and Mile Gu.  Extreme dimensionality reduction with quantum modeling.  Physical Review Letters, 125(26):260501, 2020.
  * [52] Farzad Ghafari, Nora Tischler, Carlo Di Franco, Jayne Thompson, Mile Gu, and Geoff J Pryde.  Interfering trajectories in experimental quantum-enhanced stochastic simulation.  Nature communications, 10(1):1–8, 2019.
  * [53] Yuxiang Yang, Lennart Baumgärtner, Ralph Silva, and Renato Renner.  Accuracy enhancing protocols for quantum clocks.  arXiv preprint arXiv:1905.09707, 2019.
  * [54] Mischa P Woods, Ralph Silva, Gilles Pütz, Sandra Stupar, and Renato Renner.  Quantum clocks are more accurate than classical ones.  arXiv preprint arXiv:1806.00491, 2018.
  * [55] S. Still, D. A. Sivak, A. J. Bell, and G. E. Crooks.  Thermodynamics of Prediction.  Physical Review Letters, 109(12):120604, sep 2012.
  * [56] A. J. P. Garner.  The fundamental thermodynamic bounds on finite models.  arXiv preprint arXiv:1912.03217, dec 2019.
  * [57] Alexander B. Boyd, Dibyendu Mandal, and James P. Crutchfield.  Thermodynamics of modularity: Structural costs beyond the landauer bound.  Phys. Rev. X, 8:031036, Aug 2018.
  * [58] Alexander B Boyd, James P Crutchfield, and Mile Gu.  Thermodynamic machine learning through maximum work production.  arXiv preprint arXiv:2006.15416, 2020.
  * [59] Artemy Kolchinsky and David H Wolpert.  Dependence of dissipation on the initial distribution over states.  Journal of Statistical Mechanics: Theory and Experiment, 2017(8):083202, 2017.
  * [60] Paul M Riechers and Mile Gu.  Initial-state dependence of thermodynamic dissipation for any quantum process.  Physical Review E, 103(4):042145, 2021.
  * [61] Diederik P Kingma and Jimmy Ba.  Adam: A method for stochastic optimization.  arXiv preprint arXiv:1412.6980, 2014.
  * [62] Roger A Horn and Charles R Johnson.  Matrix analysis.  Cambridge university press, 2012.
  * [63] Sarah E. Marzen and James P. Crutchfield.  Nearly maximally predictive features and their dimensions.  Phys. Rev. E, 95:051301, May 2017.

##  Appendix A Proof of Completeness

Here we prove that for each model parametrised by $B=\\{B_{x}\\}_{x}$, the set
$A=\\{A_{x}\\}_{x}$ is complete. This follows from

$\displaystyle\sum_{x}A_{x}^{\dagger}A_{x}$
$\displaystyle=\sum_{x}\frac{1}{\lambda}({W}^{-1})^{\dagger}{B_{x}}^{\dagger}W^{\dagger}WB_{x}W^{-1}$
$\displaystyle=\frac{1}{\lambda}({W}^{\dagger})^{-1}\left(\sum_{x}{B_{x}}^{\dagger}W^{\dagger}WB_{x}\right)W^{-1}$
$\displaystyle=({W}^{\dagger})^{-1}W^{\dagger}WW^{-1}=\mathbb{I},$ (8)

where we use $({W}^{\dagger})^{-1}=({W}^{-1})^{\dagger}$.

##  Appendix B Building Unitary Circuits

Once we have the Kraus operators $A_{x}$ that describe a model, standard
techniques enable construction of a full unitary circuit. Specifically, we
need to find a $|\mathcal{X}|\hat{d}_{q}\times|\mathcal{X}|\hat{d}_{q}$
unitary operator $U$ such that

$\langle x|U|0\rangle=A_{x}.$ (9)

As a result, some of the elements of the unitary operators are predefined,
i.e.,

$\langle j|\langle x|U|k\rangle|0\rangle=A_{x,jk}.$ (10)

These predefined elements forms $\hat{d}_{q}$ columns of the unitary
operators, i.e., $|\phi_{k0}\rangle:=U|k\rangle|0\rangle$. Such columns
$|\phi_{k0}\rangle$ are mutually orthogonal quantum states, due to the
completeness relation of $A_{x}$. Since an operator $U$ is unitary if and only
if its columns form an orthogonal basis, the remaining task is to find the
complementary quantum states $|\phi_{kx}\rangle$ so that they form orthogonal
states. This can be done by the Gram-Schmidt process [62].

##  Appendix C Quantum encoding maps

The encoding map $\mathcal{E}$ should map the past $\overleftarrow{x}$ onto a
pure quantum state $|\sigma_{\overleftarrow{x}}\rangle$, which is used as the
initial configuration of the quantum model. Assume first that the process has
finite Markov order $\kappa$ such that only the latest $\kappa$ steps of a
given past $\overleftarrow{x}$ affect future statistics, i.e.,
$P(\overrightarrow{X}|\overleftarrow{x})\equiv
P(\overrightarrow{X}|x_{-\kappa:0})$. Writing $0$ as the element of
$\mathcal{X}$ that also coincides with the initial state of the output
register, for the purpose of the model, any given past $\overleftarrow{x}$ is
equivalent to the concatenation $0^{\infty},x_{-\kappa:0}$, where
$0^{\infty}:=\cdots 0,0,0$. That is,

$P(\overrightarrow{x}|\overleftarrow{x})=P(\overrightarrow{x}|0^{\infty},x_{-\kappa:0}).$
(11)

Recall that an output of $0^{\infty}$ corresponds to applying Kraus operator
$A_{0}$ infinitely many times on the initial quantum state. As a result, any
initial quantum state converges to its leading eigenvector
$|\sigma_{0}\rangle$ of $A_{0}$ after applying $A_{0}$ an infinite number of
times. The quantum state transitions to
$|\sigma_{x_{-\kappa:0}}\rangle:=A_{x_{-\kappa:0}}|\sigma_{0}\rangle$ after
outputting $x_{-\kappa:0}$. The choice of output symbol here was arbitrary, as
such, the argument can be repeated for any output $x$ in place of $0$.

When the process has infinite Markov order, perfect initialization of the
machine would require an observation of the entire process history. However,
one can approximate by taking a long but finite string, which will have
exponentially diminishing error with the length of the string [63]. As such,
the above choice of $|\sigma_{0}\rangle$ remains valid.

##  Appendix D Hyperparameters

The length of the sample in the first two examples was $1000$, while that of
the third example was $5000$. The third example – the renewal process – has
more causal states, and training quantum models with higher $\hat{d}_{q}$ also
requires more data. For all of the three examples, the learning rate $\eta$ is
$0.1$ and the optimizer is Adam [61]. We repeat the training three times and
choose the quantum model of minimal cost.

The evaluation of the error $D_{e}(P,\hat{P})$ relies on computing the KL-
divergence
$\mathcal{D}_{\rm
KL}(P(X_{0:L}|\overleftarrow{x})||\hat{P}(X_{0:L}|\epsilon(\overleftarrow{x}))$.
Here, the length of the past $\overleftarrow{x}$ is $5$ and the length of the
future $L$ is $1$ (see Appendix G). This amounts to assuming that the Markov
order of the process is less than or equal to 5, which is true by construction
for our examples. The update process is then repeated until $\Delta C$ is less
than $0.1$.

##  Appendix E Learned Quantum Models

Here we present details of the quantum models produced by our model inference
algorithm. We show the Kraus operators together with the corresponding
encoding quantum states. As $\hat{d}_{q}$ is set to $2$, each encoding quantum
state is represented by a point in a Bloch sphere.

Case 1: the asymmetric process – For a $\hat{d}_{q}=2$ model of the asymmetric
process (recall Figure 3a) with $p=0.3,q=0.8$, we found the Kraus operators:

$\begin{split}&A_{0}=\begin{bmatrix}0.676&0.317\\\
0.316&0.150\end{bmatrix},\\\ &A_{1}=\begin{bmatrix}-0.264&0.534\\\
-0.336&0.729\end{bmatrix},\\\ &A_{2}=\begin{bmatrix}-0.241&-0.095\\\
0.449&0.227\end{bmatrix}.\\\ \end{split}$ (12)

For the asymmetric processes, real amplitudes are sufficient to train a
quantum model with $\hat{d}_{q}=2$ to good accuracy [19].

The quantum states $|\sigma_{\overleftarrow{x}}\rangle$ associated with
encodings of pasts of length $5$ are plotted in Figure 8. The more likely
quantum states lie in three clusters. This coincides with our expectation from
the classical $\varepsilon$-machine, which has $3$ causal states– but here,
this has been discovered directly from the sample during training. Moreover,
the location of these states aligns with our expectation from the theoretical
optimal quantum states ascertained in [19]: two clusters (pink and blue)
approximately correspond to orthogonal quantum states, while the third (green)
lies biased between the two.

Figure 8:  Example: Asymmetric process – encoded states. Quantum states
produced by the learned $\hat{d}_{q}=2$ encoding map for the asymmetric
process with $p=0.3$, $q=0.8$, represented on a Bloch sphere (viewed sideways
and top-down). Each point represents a pure state mapped to by one history (of
length $5$), and the colour represents the history’s associated causal state.
For larger points (representing the most probable $99\%$ of histories), the
area is proportional to the probability of that history. For smaller points
(representing the next $0.999\%$ of histories), the opacity is proportional to
the probability.

Case 2: the quasi-cycle – For the quasi-cycle, we show the case where
$p=0.5,\delta=0.1$. The Kraus operators found for $\hat{d}_{q}=2$ are

$\begin{split}&A_{0}=\begin{bmatrix}0.021-0.077i&-0.104+0.084i\\\
-0.277+0.236i&0.656+0.098i\\\ \end{bmatrix},\\\
&A_{1}=\begin{bmatrix}0.508+0.135i&-0.233+0.465i\\\
0.279-0.171i&0.073+0.277i\\\ \end{bmatrix},\\\
&A_{2}=\begin{bmatrix}0.371+0.271i&0.304+0.010i\\\
-0.209+0.473i&0.055+0.308i\\\ \end{bmatrix}.\end{split}$ (13)

The quantum states $|\sigma_{\overleftarrow{x}}\rangle$ are plotted in Figure
9. As shown in Figure 9, most of the states (weighted by probability of
occurrence) lie in three distinct clusters – but here, no pair of clusters is
orthogonal.

Figure 9: Example: Quasicycle – encoded states. Quantum states produced by the
learned $\hat{d}_{q}=2$ encoding map for a quasicycle with $p=0.5,\delta=0.1$,
represented on a Bloch sphere. The diagram is interpreted as in Figure 8.

Case 3: the uniform renewal process – For the discrete renewal process (Figure
6a) we show the results for $N=3,\hat{d}_{q}=2$. The Kraus operators found are

$\begin{split}&A_{0}=\begin{bmatrix}0.064+0.170i&0.043+0.246i\\\
-0.196+0.825i&0.499-0.079i\end{bmatrix},\\\
&A_{1}=\begin{bmatrix}0.490-0.053i&0.304+0.753i\\\
0.005-0.068i&0.048-0.142i\end{bmatrix}\end{split}$ (14)

The associated encoded states are plotted in Figure 10.

Figure 10: Example: Renewal process – encoded states. Quantum states produced
by the learned $\hat{d}_{q}=2$ encoding map for an $N=3$ discretization of the
uniform renewal process, represented on a Bloch sphere. The diagram is
interpreted as in Figure 8.

##  Appendix F Bounding the minimum-distortion classical model

###  F.1 Past partitions, future morphs, and pre-models

Let $\overleftarrow{X}$ be the set of all possible observable histories of
some process $P$, and $\mathcal{R}$ be a countable set of possible model
memory states. A partition of histories or encoding $\mathcal{E}$ is then a
function $\mathcal{E}:\overleftarrow{X}\to\mathcal{R}$. We refer to the number
of distinct values in $\mathcal{R}$ taken by $\mathcal{E}$ (excluding those
that occur with measure zero with respect to $P(\overleftarrow{X})$) as the
model dimension of $\mathcal{E}$, which we write as $D_{c}$. That is, each
$\mathcal{E}$ effectively divides the set of histories into $D_{c}$ mutually
exclusive partitions. To use this encoding as a model, we must also supply a
second map $\mathcal{F}:\mathcal{R}\to\overrightarrow{X}$ from memory states
to a set of statistics over the future of the process known as the future
morphs of the memory states. Together, the composition of $\mathcal{E}$ and
$\mathcal{F}$ defines a conditional probability distribution
$\hat{P}(\overrightarrow{X}|\mathcal{E}(\overleftarrow{x}))$, which can be
directly compared with a process’s conditional probability distribution
$P(\overrightarrow{X}|\overleftarrow{x})$ for each history
$\overleftarrow{x}$.

As a clarifying example, let us cast the $\varepsilon$–machine into this
language. Here $\mathcal{E}=\varepsilon$, the partition of histories according
to the equivalence relation $\overleftarrow{x}\sim\overleftarrow{x}^{\prime}$
if and only if
$P(\overrightarrow{X}|\overleftarrow{x})=P(\overrightarrow{X}|\overleftarrow{x}^{{}^{\prime}})$.
We write $\mathcal{R}\equiv\mathcal{S}$, where the elements of $\mathcal{S}$
are known as the process’s causal states. Meanwhile, the future morph
associated with any casual state is the same as any constituent pasts within
the class, i.e. for all $s_{i}\in\mathcal{S}$,
$\hat{P}(\overrightarrow{X}|s_{i})=P(\overrightarrow{X}|\overleftarrow{x})$
where $\overleftarrow{x}$ is an arbitrary choice of the past such that
$\varepsilon(\overleftarrow{x})=s_{i}$.

In general, a couple $(\mathcal{E},\mathcal{F})$ does not always correspond to
a stationary hidden Markov model, since there might not be a systematic action
$\mathcal{P}$ that can be repeatedly applied to $\mathcal{R}$ to sequentially
produce a sample from the morph. As such, we call the couple
$(\mathcal{E},\mathcal{F})$ a pre-model:

######  Definition 1 (Pre-model).

A pre-model is the couple $(\mathcal{E},\mathcal{F})$ where $\mathcal{E}$ is a
function from histories ($\overleftarrow{X}$) to memory states
($\mathcal{R}$), and $\mathcal{F}$ is a function from memory states
($\mathcal{R}$) to distributions over future statistics
($P(\overrightarrow{X}|r),\forall r\in\mathcal{R}$).

We write the set of pre-models with maximum dimension $d$ as
$\mathcal{M}^{d}_{\infty}$. This set will be used to derive bounds on the
properties of unifilar hidden Markov models, as we shall subsequently show.
First, we define a particular class of pre-model:

######  Definition 2 ($K$-unifilar pre-model).

For $K\in\mathbb{Z}^{+}$, a $K$-unifilar pre-model is the triple
$(\mathcal{E},\mathcal{F},\mathcal{P}_{K})$ consisting of a pre-model
$(\mathcal{E},\mathcal{F})$ with a deterministic map
$\mathcal{P}_{K}:\mathcal{R}\times\mathcal{X}^{\otimes K}\to\mathcal{R}$ such
that
$\mathcal{P}_{K}\\!\left(\mathcal{E}(\overleftarrow{x}),x_{0:K}\right)=\mathcal{E}(\overleftarrow{x}x_{0:K})$
for all $x_{0:K}$ and all $\overleftarrow{x}$.

We write the set of all $K$-unifilar pre-models of maximum dimension $d$ as
$\mathcal{M}^{d}_{K}$. Crucially, the valid morphs of $K$-unifilar pre-models
are highly constrained:

######  Lemma 1.

A $K$-unifilar pre-model $(\mathcal{E},\mathcal{F},\mathcal{P}_{K})$ can
always be specified by the triple
$(\mathcal{E},\mathcal{F}_{K},\mathcal{P}_{K})$ where $\mathcal{F}_{K}$ is a
map from $\mathcal{R}$ to probability distributions over words of length $K$.

###### Proof.

First, let $r:=\mathcal{E}(\overleftarrow{x})$ and note that the definitions
of $\mathcal{P}_{K}$ and $\mathcal{F}$ imply that
$\mathcal{F}\left(\mathcal{E}\\!\left(\overleftarrow{x}x_{0:K}\right)\right)=\mathcal{F}\left(\mathcal{P}_{K}\\!\left(\mathcal{E}\\!\left(\overleftarrow{x}\right),x_{0:K}\right)\right)$
allowing for any substitutions of the form (for any $L>K$):

$P(X_{K:L}\\!=\\!x_{K:L}|rx_{0:K})=P(X_{0:L\\!-\\!K}\\!=\\!x_{K:L}|\mathcal{P}_{K}(r,x_{0:K})),$
(15)

since if this was not the case then $\mathcal{F}$ would not consistently
assign the correct future morph for some memory states.

For notational brevity, we recursively define the set of functions
$\\{r_{n}:\mathcal{R}\times\mathcal{X}^{nK}\to\mathcal{R}\\}_{n}$ as:

$r_{n}\\!\left(r,x_{0:nK}\right)=\mathcal{P}_{K}\\!\left(r_{n\\!-\\!1}\\!\left(r,x_{0:(n\\!-\\!1)K}\right),x_{(n\\!-\\!1)K:nK}\right)$
(16)

for $n\geq 1$, and $r_{0}=r$. For any $M\in\mathbb{Z}^{+}$, we can thus expand
the probability distribution

$\displaystyle P(X_{0:MK}\\!=\\!x_{0:MK}\,|\,r)$
$\displaystyle=P(X_{0:K}\\!=\\!x_{0:K}\,|\,r)\,P(X_{K:KM}\\!=\\!x_{K:MK}\,|\,rx_{0:K})$
$\displaystyle=P(X_{0:K}\\!=\\!x_{0:K}\,|\,r_{0})\,P(X_{0:(M\\!-\\!1)K}\\!=\\!x_{K:MK}\,|\,r_{1}\\!\left(r,x_{0:K}\right))$
$\displaystyle=P(X_{0:K}\\!=\\!x_{0:K}\,|\,r_{0})\,P(X_{0:K}\\!=\\!x_{K:2K}\,|\,r_{1}\\!\left(r,x_{0:K}\right))$
$\displaystyle\quad\cdot
P(X_{K:(M\\!-\\!1)K}\\!=\\!x_{2K:MK}\,|\,r_{1}\\!\left(r,x_{0:K}\right)x_{K:2K})$
$\displaystyle=P(X_{0:K}\\!=\\!x_{0:K}\,|\,r_{0})\,P(X_{0:K}\\!=\\!x_{K:2K}\,|\,r_{1}\\!\left(r,x_{0:K}\right))$
$\displaystyle\quad\cdot
P(X_{0:(M\\!-\\!2)K}\\!=\\!x_{2K:MK}\,|\,r_{2}\\!\left(r,x_{0:2K}\right))$
$\displaystyle=\ldots$
$\displaystyle=\prod_{i=0}^{M-1}P(X_{0:K}\\!=\\!x_{iK:(i+1)K}\,|\,r_{i}(r,x_{0:iK})).$
(17)

Here, Bayesian expansion allows us to make the first and third equalities, and
$K$-unifilarity (via Eqs. 15 and 16) allows us to make the substitutions for
the second and fourth equalities. Thus, we can use this to generate a future
morph for words of any length $KM$, as a function of the probability
distribution over the next $K$ symbols. (Probabilities over words of length
that are not multiples of $K$ can always be found by taking marginals of a
longer word that is a multiple of $K$.)

The contrapositive implies that if $\mathcal{F}$ does not assign probabilities
in this way, then it cannot satisfy Eq. 15, and hence is not a $K$-unifilar
pre-model. ∎

An immediate corollary of this is that a $1$-unifilar pre-model is exactly a
unifilar hidden Markov model. It also immediately follows that any
$K$-unifilar pre-model for $K\in\mathbb{Z}^{+}$ is also an $NK$-unifilar pre-
model for $N\in\mathbb{Z}^{+}$ – and particularly, unifilar models are
$K$-unifilar pre-models for all $K\in\mathbb{Z}^{+}$. For any maximum
dimension $d$, this sets up the hierarchy
$\mathcal{M}^{d}_{1}\subseteq\mathcal{M}^{d}_{K}\subseteq\mathcal{M}^{d}_{\infty}$.

###  F.2 Calculating the distortion $D_{e}$

For the purposes of this article, we are interested in minimizing the
distortion (Eq. 2) of models. Particularly, consider two processes $P$ and
$\hat{P}$, and for calculational convenience, define the arguments within the
limit of Eq. 1 and Eq. 2 as

$\displaystyle D_{L}(P,\hat{P})$
$\displaystyle=\sum_{\overleftarrow{x}}P(\overleftarrow{x})\,D_{L}(P,\hat{P}|\overleftarrow{x}),$
(18)

such that (Eq. 2) may be written

$D_{e}(P,\hat{P})=\lim_{L\to\infty}D_{L}(P,\hat{P}).$ (19)

######  Lemma 2.

For stationary processes $P$ and $\hat{P}$:

$D_{e}(P,\hat{P})=D_{K}(P,\hat{P})=D_{1}(P,\hat{P})\qquad\forall
K\in\mathbb{Z}^{+}.$ (20)

###### Proof.

According to Bayesian rules, for $L>1$:

$P(x_{0:L}|\overleftarrow{x})=P(x_{0}|\overleftarrow{x})P(x_{1:L}|\overleftarrow{x}x_{0}).$
(21)

Then (expanding the definition of KL-divergence),

$\displaystyle LD_{L}(P,\hat{P})$ $\displaystyle\hskip
30.00005pt=\sum_{\overleftarrow{x}}P(\overleftarrow{x})\sum_{x_{0:L}}P(x_{0:L}|\overleftarrow{x})\log\frac{P(x_{0:L}|\overleftarrow{x})}{\hat{P}(x_{0:L}|\overleftarrow{x})}$
$\displaystyle=\sum_{\overleftarrow{x}}P(\overleftarrow{x})\sum_{x_{0}}P(x_{0}|\overleftarrow{x})\log\frac{P(x_{0}|\overleftarrow{x})}{\hat{P}(x_{0}|\overleftarrow{x})}$
$\displaystyle\qquad+\sum_{\overleftarrow{x}}P(\overleftarrow{x})\sum_{x_{0:L}}P(x_{0:L}|\overleftarrow{x})\log\frac{P(x_{1:L}|\overleftarrow{x}x_{0})}{\hat{P}(x_{1:L}|\overleftarrow{x}x_{0})}.$
(22)

Since $P$ and $\hat{P}$ are stationary, we can substitute
$P(x_{1:L}|\overleftarrow{x}x_{0})=P(x_{0:L-1}|\overleftarrow{x})$, and hence
for $L>1$:

$LD_{L}(P,\hat{P})=D_{1}(P,\hat{P})+(L-1)D_{L-1}(P,\hat{P}).$ (23)

Then, inductively

$D_{L}(P,\hat{P})=D_{1}(P,\hat{P})\qquad\forall L\geq 1,$ (24)

and hence also the limit

$D_{e}(P,\hat{P})=\lim_{L\to\infty}D_{L}(P,\hat{P})=D_{1}(P,\hat{P}).$ (25)

∎

Similarly, this enables a shortcut to calculate the distortion between a
process and a $K$-unifilar pre-model:

######  Lemma 3.

Let $P$ be some stationary process, and
$(\mathcal{E},\mathcal{F}_{K},\mathcal{P}_{K})$ be a $K$-unifilar pre-model
whose future morph is generated by
$\hat{P}(X_{0:K}|\mathcal{E}(\overleftarrow{x}))$. Then,

$D_{e}(P,\hat{P})=D_{K}(P,\hat{P}).$ (26)

###### Proof.

We can group each word of $K$ contiguous symbols from process $P$ over
alphabet $\mathcal{X}$ into a single symbol of process $Q$ with alphabet
$\mathcal{Y}:=\mathcal{X}^{K}$, such that the distributions over the two
processes are related by

$Q(\cdots,y_{0},y_{1},\cdots)=P(\cdots,x_{0:R},x_{R:2R},\cdots).$ (27)

We similarly group $\hat{P}$ into $\hat{Q}$. Since
$(\mathcal{E},\mathcal{F},\mathcal{P}_{K})$ is a $K$-unifilar pre-model for
$\hat{P}$, it immediately defines a unifilar model for $\hat{Q}$ with encoding
$\mathcal{E}_{Q}:\overleftarrow{Y}\to\mathcal{R}$ satisfying
$\mathcal{E}_{Q}(\overleftarrow{y})=\mathcal{E}(\overleftarrow{x})$ where
$\overleftarrow{y}$ is the grouping of $\overleftarrow{x}$. Such combination
has no effect on the KL–divergence and hence the distortion satisfies

$D_{e}(Q,\hat{Q})=KD_{e}(P,\hat{P}).$ (28)

Meanwhile, using Lemma 2):

$\displaystyle D_{e}(Q,\hat{Q})$ $\displaystyle=D_{1}(Q,\hat{Q})$
$\displaystyle=\sum_{\overleftarrow{y}}Q(\overleftarrow{y})D_{1}\\!\left(Q,\hat{Q}|\overleftarrow{y}\right)$
$\displaystyle=\sum_{\overleftarrow{x}}P(\overleftarrow{x})KD_{K}\\!\left(P,\hat{P}|\overleftarrow{x}\right)$
$\displaystyle=KD_{K}(P,\hat{P}).$ (29)

And hence $D_{e}(P,\hat{P})=D_{K}(P,\hat{P})$. ∎

###  F.3 Bounding the minimum distortion

A consequence of the above lemmas is to make it easier to calculate bounds
that minimize $K$-unifilar pre-models and unifilar models:

######  Lemma 4.

Let $P$ be a stationary process, and let $K,N\in\mathbb{Z}^{+}$. Then for
every dimension $d$

$\min_{M\in\mathcal{M}^{d}_{K}}D_{e}(P,M)\geq\min_{M\in\mathcal{M}^{d}_{KN}}D_{K}(P,M).$
(30)

###### Proof.

First, Lemma 3 tells us that for any $K$-unifilar pre-model
$D_{e}(P,M)=D_{K}(P,M)$ and hence
$\min_{M\in\mathcal{M}_{K}}D_{e}(P,M)=\min_{M\in\mathcal{M}_{K}}D_{K}(P,M)$.
Next, the hierarchy of models
$\mathcal{M}^{d}_{K}\subseteq\mathcal{M}^{d}_{KN}$, and hence a minimization
over $\mathcal{M}^{d}_{KN}$ lower bounds a minimization over
$\mathcal{M}^{d}_{K}$. This proves Eq. 30. ∎

Why might we choose some finite $K$ rather than minimizing $D_{e}$ over
generic ($K=\infty$) pre-models? First, $\mathcal{M}_{\infty}$ is too
permissive, and $D_{e}$ is too forgiving. Consider any (stationary, aperiodic)
process, and let us assign a $1$–dimensional model, but whose single future
morph defines probabilities arbitrarily far into the future, with statistics
that increasingly match those of $P$ (e.g. by assigning a weighted average of
$P$’s causal state’s future morphs). For such a (highly non-unfilar) model,
$D_{e}\to 0$, saturating the arithmetically obvious bound. Second, the
parameter space of
$P\left(\overrightarrow{X}|\mathcal{E}(\overleftarrow{x})\right)$ is a priori
infinite, and the above Lemmas allow us to greatly reduce the relevant
parameter space of future morphs in our calculations to probability
distributions over words of finite length $K$.

Although the above Lemmas allow us to reduce the dimensionality of
$\mathcal{F}$ assigned to each memory state, there still remains the problem
that there are infinitely many possible partitions of the past $\mathcal{E}$.
We will argue, with the next set of Lemmas, that for exploring minimum bounds,
we can restrict ourselves to the finite set of coarse–grainings over the
causal states.

First, we note that the problem of finding a minimum distortion model is
trivial if we allow models with high–enough dimension:

######  Lemma 5.

For a process $P$ with topological complexity $d_{c}$, then for all model
dimensions $\hat{d}_{c}\geq d_{c}$, there is a zero–distortion unifilar model
$M$ such that

$D_{K}(P,M)=0\qquad\forall K\in\mathbb{Z}^{+}.$ (31)

###### Proof.

Such a minimum distortion model can be realized by the $\varepsilon$–machine
of $P$, which has zero distortion. By definition of causal states,
$P(X_{0:K}|\overleftarrow{x})=P(X_{0:K}|\varepsilon(\overleftarrow{x}))$ for
all $K$ and $\overleftarrow{x}$, and hence $\mathcal{D}_{\rm
KL}(P(X_{0:K}|\overleftarrow{x})\,||\,P(X_{0:K}|\varepsilon(\overleftarrow{x})))=0$
for all $L$ and $\overleftarrow{x}$. This implies that all
$D_{K}(P,\hat{P})=0$. ∎

######  Lemma 6.

For any process $P$, for every dimension $d$, and for every
$K\in\mathbb{Z}^{+}$, if two pasts $\overleftarrow{y}$ and $\overleftarrow{z}$
have identical future morphs in $P$, then there exists a pre-model
$\left(\mathcal{E}^{d}_{\min},\mathcal{F}\right)$ that minimizes $D_{K}$ and
satisfies

$\mathcal{E}^{d}_{\min}(\overleftarrow{y})=\mathcal{E}^{d}_{\min}(\overleftarrow{z}).$
(32)

###### Proof.

If $d\geq d_{c}$, the $\varepsilon$–machine has the minimum distortion (see
Lemma 5) and the causal state partition $\varepsilon$ satisfies Eq. 32 by
definition. Likewise, if $d_{c}=1$, there is only one encoding onto a single
state and so Eq. 32 is trivially satisfied. Thus, it remains only to prove the
cases where $1<d<d_{c}$.

Suppose $\overleftarrow{y}$ and $\overleftarrow{z}$ have identical future
morphs in $P$, and we have a pre-model $(\mathcal{E},\mathcal{F})$ where
encoding $\mathcal{E}(\overleftarrow{y})\neq\mathcal{E}(\overleftarrow{z})$.
Then without loss of generality (by switching the labels of
$\overleftarrow{y}$ and $\overleftarrow{z}$ if necessary) the implied
statistics $\hat{P}$ satisfy

$D_{K}(P,\hat{P}|\overleftarrow{y})\leq D_{K}(P,\hat{P}|\overleftarrow{z}),$
(33)

where $D_{K}(P,\hat{P}|\overleftarrow{x})$ is defined in Eq. 18.

We can then construct a new pre-model with identical $\mathcal{F}$, but with a
new mapping $\mathcal{E}^{{}^{\prime}}$ identical to $\mathcal{E}$ in every
way, except it now maps $\overleftarrow{z}$ to
$\mathcal{E}\\!\left(\overleftarrow{y}\right)$ instead of
$\mathcal{E}\\!\left(\overleftarrow{z}\right)$. This has implied statistics
$\hat{P}^{\prime}$, and

$D_{K}(P,\hat{P}^{\prime}|\overleftarrow{z})=D_{K}(P,\hat{P}|\overleftarrow{y})\leq
D_{K}(P,\hat{P}|\overleftarrow{z}).$ (34)

Hence $D_{K}$ of $(\mathcal{E}^{\prime},\mathcal{F})$ satisfies

$\displaystyle D_{K}(P,\hat{P}^{\prime})$
$\displaystyle=P(\overleftarrow{z})D_{K}(P,\hat{P}^{\prime}|\overleftarrow{z})+\sum_{\overleftarrow{x}\neq\overleftarrow{z}}P(\overleftarrow{x})D_{K}(P,\hat{P}^{\prime}|\overleftarrow{x})$
$\displaystyle=P(\overleftarrow{z})D_{K}(P,\hat{P}^{\prime}|\overleftarrow{z})+\sum_{\overleftarrow{x}\neq\overleftarrow{z}}P(\overleftarrow{x})D_{K}(P,\hat{P}|\overleftarrow{x})$
$\displaystyle\leq
P(\overleftarrow{z})D_{K}(P,\hat{P}|\overleftarrow{z})+\sum_{\overleftarrow{x}\neq\overleftarrow{z}}P(\overleftarrow{x})D_{K}(P,\hat{P}|\overleftarrow{x})$
$\displaystyle=D_{K}(P,\hat{P}).$ (35)

It then follows that for any pre-model where $\overleftarrow{y}$ and
$\overleftarrow{z}$ with identical future morphs map to different memory
states, there is another pre-model of the same (or lower) dimension such that
$\overleftarrow{y}$ and $\overleftarrow{z}$ map to the same memory state, and
this pre-model has the same or lower distortion. Hence, among the minimum
distortion pre-models of a given process, there will always be a pre-model
where $\overleftarrow{y}$ and $\overleftarrow{z}$ map to the same memory
state. ∎

The above lemma implies the following:

######  Lemma 7.

For a process $P$ with topological complexity $d_{c}$ For every dimension
$d<d_{c}$, and every $K\in\mathbb{Z}^{+}$ there is a pre-model $M$ that
minimizes $D_{K}(P,M)$ whose encoding $\mathcal{E}^{d}_{\min}$ is a map onto a
coarse graining of the causal states $\mathcal{S}$.

###### Proof.

This follows from Lemma 6 by noting that pasts in the same causal state have
the same future morph by definition. ∎

Such coarse–grained models admit a computational shortcut for calculating
their distortion: we can define $\pi_{i}:=\sum_{\overleftarrow{x}\in
s_{i}}P(\overleftarrow{x})$ for each causal state $s_{i}\in\mathcal{S}$, such
that the distortion is

$D_{K}(P,\hat{P})=\sum_{i=1}^{d_{c}}\pi_{i}\,D_{K}(P,\hat{P}|s_{i}),$ (36)

where for each $s_{i}$,
$D_{K}(P,\hat{P}|s_{i}):=D_{K}(P,\hat{P}|\overleftarrow{x})$ for one arbitrary
choice of $\overleftarrow{x}\in s_{i}$ (since by definition the value is the
same for all such choices).

Crucially, the above Lemma enables an exhaustive search for $D_{K}$–optimal
pre-models for any process with finite topological complexity – instead of
having to consider an infinite number of possible partitions, we can iterate
through the various combinations of causal states.

######  Lemma 8.

For a process $P$, for every dimension $d$, and every $K\in\mathbb{Z}^{+}$
there is pre-model $M=(\mathcal{E}^{d}_{\min},\mathcal{F}^{d}_{\min})\in
M^{d}_{\infty}$ that minimizes $D_{K}(P,M)$ such that every memory state
$r_{i}$ whose entire pre-image (except for a measure zero subset with respect
to $P(\overleftarrow{x})$) is mapped to the causal state $s_{i^{\prime}}$, has
the future morph
$\mathcal{F}^{d}_{\min}(r_{i})=P(\overrightarrow{X}|s_{i^{\prime}})$.

###### Proof.

Again, we consider the case where model dimension $d<d_{c}$ since the case
$d\geq d_{c}$ is trivially satisfied by the $\varepsilon$–machine (Lemma 5).
From Lemma 7 we may restrict ourselves to coarse-grainings of
$\varepsilon$–machines.

Let $s_{i}$ be the causal states of the process, indexed such that $i=1\ldots
k$ correspond to states that are not merged, and $i=k+1\ldots d_{c}$ are
merged into states $r_{1}\ldots r_{d-k}$ (i.e. the states of the model are
$\mathcal{R}=\\{s_{1},\ldots s_{k},r_{1},\ldots r_{d-k}\\}$). As the Lemma’s
claim is vacuously true if $k=0$, we consider the cases when $k\geq 1$. Then,
using the coarse–graining structure, and choosing an arbitrary past
$\overleftarrow{x}_{i}\in s_{i}$ for each causal state,

$D_{K}(P,\hat{P})=\sum_{i=1}^{k}\pi_{i}D_{K}(P,\hat{P}|\overleftarrow{x}_{i})+\sum_{i=k+1}^{d}\pi_{i}D_{K}(P,\hat{P}|\overleftarrow{x}_{i}).$
(37)

If for any $i\in[1,k]$ we have $D_{K}(P,\hat{P}|\overleftarrow{x}_{i})>0$, we
can instead form a new pre-model that assigns the future morph
$\mathcal{F}(s_{i})=P(\overrightarrow{X}|s_{i})$. Then, in this new pre-model
$D_{k}(P,\hat{P}^{\prime}|\overleftarrow{x}_{i})=0$ for all such $s_{i}$.
Thus,

$D_{K}(P,\hat{P}^{\prime})=\sum_{i=k+1}^{d}\pi_{i}D_{K}(P,\hat{P}|s_{i})\leq
D_{K}(P,\hat{P}).$ (38)

∎

In simpler words: when we look for a $D_{K}$–optimal pre-model by merging
states of an $\varepsilon$–machine, we should not alter the future morphs of
any “unmerged” states from how they are in the $\varepsilon$–machine.

######  Theorem 1.

Let $P$ be any process. Then, for every model dimension $d$, the algorithm in
Box.2 produces a bound on the minimum distortion $D_{e}$ of a unifilar model
by considering all possible mergers of causal states, and searching over the
next $K$ output statistics.

###### Proof.

Lemma 7 tells us that a $D_{K}$–optimal pre-model exists that is formed by
merging causal states, but says nothing as to whether such a pre-models is
unifilar. However, Lemma 4 implies that such a $D_{K}$-optimized pre-model
bounds the $D_{e}$–optimized unifilar models. Moreover, $D_{K}$, as a function
of $\hat{P}(X_{0:K}|\mathcal{E}(\overleftarrow{x}))$, depends only on the next
$K$ output statistics. As such, to find this minimum bound, we can iterate
through the various coarse–grainings of casual states of a given dimension,
and search through finite parameter space of $\hat{P}$ to minimize $D_{K}$.
The minimum found here will lower bound the lowest value of $D_{e}$ achievable
by a unifilar model. The proof is illustrated in Figure 11.

Figure 11: $P$ represents the given stochastic process. For a given model
dimension $d$, $\mathcal{M}_{\infty}^{d}$ represent all pre-models,
$\mathcal{M}_{K}^{d}$ represents all K-unifilar pre-models, and
$\mathcal{M}_{1}^{d}$ represents all unifilar models.

∎

###  F.4 Minimum–distortion approximations of Markov processes

For Markov processes, this bound is tight. First:

######  Lemma 9.

Let $P$ be a process with Markov order $\kappa$. Then, any encoding map
$\mathcal{E}$ formed by coarse–graining the causal states of $P$ admits an
$\kappa$-unifilar pre-model. An $\kappa$-unifilar pre-model that minimizes
$D_{\kappa}(P,\hat{P})$ is in this set.

###### Proof.

Let $\varepsilon$ be the encoding map onto causal states, and
$\mathcal{C}:\mathcal{S}\to\mathcal{R}$ be the coarse–graining map from causal
states $\mathcal{S}$ to memory states $\mathcal{R}$, such that
$\mathcal{E}=\mathcal{C}\cdot\varepsilon$. We must demonstrate that there
exists a $\mathcal{P}_{\kappa}$ such that
$\mathcal{P}_{\kappa}\left(\mathcal{E}\\!\left(\overleftarrow{x}\right),x_{0:\kappa}\right)=\mathcal{E}\\!\left(\overleftarrow{x}x_{0:\kappa}\right)$
for all $\overleftarrow{x}$ and all $x_{0:\kappa}$. Since $P$ has Markov order
$\kappa$, there exists an
$\tilde{\varepsilon}:\mathcal{X}^{\otimes\kappa}\to\mathcal{S}$ that acts on
only the final $\kappa$ symbols of the history to identify the causal states,
and for every $\overleftarrow{x}$ that ends in the same $x_{-\kappa:0}$,
$\varepsilon(\overleftarrow{x})=\tilde{\varepsilon}(x_{-\kappa:0})$. We thus
have an $\tilde{\mathcal{E}}=\mathcal{C}\cdot\tilde{\varepsilon}$ that maps to
the same memory state as $\mathcal{E}$ for every $x_{-\kappa:0}$ and every
$\overleftarrow{x}$ that ends in $x_{0:\kappa}$.

Thus, now consider
$\mathcal{P}_{\kappa}\left(\mathcal{E}\\!\left(\overleftarrow{x}\right),x_{0:\kappa}\right)=\mathcal{E}\\!\left(\overleftarrow{x}x_{0:\kappa}\right)=\tilde{\mathcal{E}}(x_{0:\kappa})$.
A candidate $\mathcal{P}_{\kappa}$ can be defined such that it is completely
independent of its first argument (the current machine state
$r\in\mathcal{R}$) and instead takes exactly the same value as
$\tilde{\mathcal{E}}$ applied to its second input (the recently output word
$x_{0:\kappa}\in\mathcal{X}^{\kappa}$). Thus, there exists an
$\kappa$-unifilar pre-model for any such $\mathcal{E}$.

Now, from Lemma 7, we know that the $D_{\kappa}$ optimal ($\infty$-unifilar)
pre-model $(\mathcal{E},\mathcal{F})$ is formed by merging causal states, but
a priori we have not demonstrated that such a model is $\kappa$-unifilar. To
calculate $D_{\kappa}$, we need only evaluate the probabilities associated
with the first $\kappa$ outputs. Thus, we can define $\mathcal{F}_{\kappa}$ as
per Lemma 1 such that $\mathcal{F}_{\kappa}$ and $\mathcal{F}$ perfectly agree
on the statistics of the first word of length $\kappa$ in the morph. Then,
taking the encoding map of $\mathcal{E}$ from this pre-model, we form a
$\kappa$-unifilar model
$(\mathcal{E},\mathcal{F}_{\kappa},\mathcal{P}_{\kappa})$ where
$\mathcal{P}_{\kappa}$ is defined as above. Since $D_{\kappa}$ only depends on
the first $\kappa$ steps of the future, minimizing the value of
$D_{\kappa}(P,\hat{P})$ provides the optimal $\kappa$-unifilar pre-model. As
such, we have formed a minimum–distortion $\kappa$-unifilar pre-model. ∎

When the process is Markovian, these lemmas then give us a systematic method
for finding the minimum distortion unifilar hidden Markov model.

######  Theorem 2.

Let $P$ be a Markov process. Then, for every model dimension $d$, we can find
a minimum distortion unifilar model by merging causal states, and searching
over the next–output statistics.

###### Proof.

If $d\geq d_{c}$, this encoding is the $\varepsilon$-machine (Lemma 5). In the
case $d<d_{c}$, we use Lemma 7, to find a $\mathcal{E}^{d}_{\min}$ that is a
coarse–graining of causal states. Then Lemma 9 specialized to $\kappa=1$
implies that this model can be made unifilar. ∎

##  Appendix G Classical minimum distortion models: Examples

As illustrative examples, we present the minimum–distortion approximate
classical model for asymmetric process with $p=0.3,q=0.8$ in Figure 12 and for
quasi cycle with $p=0.5,\delta=0.1$ in Figure 13.

Figure 12: Approximate model for asymmetric process. The approximate classical
model for the asymmetric example with $p=0.3$ and $q=0.8$. $s_{0}$ and $s_{1}$
of the original asymmetric process are merged into $r_{0}$, with transition
probabilities identified by a minimization. $s_{2}$ is mapped to $r_{1}$ with
an unchanged future morph. Figure 13: Approximate model for quasicycle. The
approximate classical model for the quasicycle example with $p=0.5$ and
$\delta=0.1$. $s_{0}$ and $s_{1}$ of the original quasicycle process are
merged into $r_{0}$. The transition probability is minimized by going through
the probability vector space. $s_{2}$ is mapped to $r_{1}$ with an unchanged
future morph.

For the third example, we evaluate $D_{N}$ as a lower bound on the error of
the optimal classical model for any discrete renewal process with $N$ states,
as shown in Figure 6b.

##  Appendix H Details of experimental quantum probability distribution.

Figure 14 compares the conditional probability distribution realized on “ibmq
athens” with the corresponding noiseless distribution of the learned quantum
model.

Figure 14:  Comparison of experimental (orange) and theoretical (blue)
realizations of the learned quantum model to the true process.  Each datum
represents the difference between the model’s conditional probability for a
given past, and the true process’s conditional probability.