{'fraction_non_alphanumeric': 0.04082468946087847, 'fraction_numerical': 0.014700986041746703, 'mean_word_length': 4.357574094401756, 'pattern_counts': {'":': 0, '<': 0, '<?xml version=': 0, '>': 1, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 3, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
The Hintons in your Neural Network: a Quantum Field Theory View of Deep Learning

In this work we develop a quantum field theory formalism for deep learning, where input signals are encoded in Gaussian states, a generalization of Gaussian processes which encode the agent's uncertainty about the input signal. We show how to represent linear and non-linear layers as unitary quantum gates, and interpret the fundamental excitations of the quantum model as particles, dubbed ``Hintons''. On top of opening a new perspective and techniques for studying neural networks, the quantum formulation is well suited for optical quantum computing, and provides quantum deformations of neural networks that can be run efficiently on those devices. Finally, we discuss a semi-classical limit of the quantum deformed models which is amenable to classical simulation.


Introduction and Contributions
Since its inception more than 10 years ago, deep learning has achieved some stunning successes, transforming entire fields such as speech recognition, automated translation, computer vision and protein folding. It's impact on industry and the economy is so large that dedicated chips are being developed just to process neural network workloads. A natural question arises, what will drive the next wave of innovation in this field?
One candidate for the next disruption in AI is quantum computing. Quantum computing is in its nascent phase with only a modest number of qubits or qumodes (i.e. states of light) available to perform computations. Moreover, it is not at all clear how to map neural networks onto a quantum computer in way that provide clear benefits. Yet, the confluence of deep learning and quantum computing holds promise and deserves exploration.
In this paper we propose a direct mapping of a deep neural network onto a optical quantum computer through the language of quantum field theory. We consider data, e.g. an image, as a finite sampling of an underlying continuous signal. We model this using quantum optical Gaussian states and show that they represent a generalization of a Gaussian process. We show how to compute the posterior of this Quantum Gaussian Process given input data. Then we show how both linear layers and nonlinear layers can be represented as unitary transformations on these Gaussian states. In the process we define a new quantum nonlinearity analogous to the classical softplus.
Because nonlinearities map us outside of the scope of Gaussian states we discuss a number of special cases which are tractable to simulate classically. The full architecture can only be simulated on an optical quantum computer.
We also provide details on how to map our model onto an optical quantum computer, including quantum Gaussian state preparation, and the implementation of linear and nonlinear layers in terms of elementary universal gates.
While we have implemented and tested our tractable approximations we emphasize that our contribution is theoretical; we do not expect that these tractable approximations perform better than their classical cousins. By describing deep learning in the language of quantum field theory we pave the way for the development of novel quantum neural network architectures in the future. We find the match between neural networks and quantum field theory very natural indeed. And perhaps amusingly, in the process of formulating deep learning in the language of quantum field theory we discovered a new particle: The "Hinton" is the elementary excitation of the quantum field from which optical quantum neural networks are made.

Related Work
Several works have recently investigated quantum neural networks as a variational quantum circuit using qubit architectures, e.g. (Farhi & Neven, 2018;Verdon et al., 2018;Beer et al., 2019;Cong et al., 2019;Huggins et al., 2019;Bondesan & Welling, 2020). Quantum optical neural networks similar to those described in this work have been arXiv:2103.04913v1 [quant-ph] 8 Mar 2021 discussed in (Killoran et al., 2019;Steinbrecher et al., 2019;Lau et al., 2017;Das et al., 2018b). Optical devices have also been considered for implementations of classical neural networks on photonic hardware due to the efficiency of matrix multiplication with optical instruments in (Shen et al., 2017). While the implementation of linear layers considered in all these works is the same as ours, our work differs from previous literature in the following two main aspects: 1) we employ Gaussian states to reason about uncertainty due to discretization errors; 2) we show how to implement quantum analogs of popular nonlinearities such as softplus using unitary gates. Technical details of these differences will be discussed after we introduce our framework below. Using Gaussian states for interpolating data that is further processed by a neural network is inspired by the recent works (Li & Marlin, 2016;Finzi et al., 2020), which however are concerned only with classical networks.
Quantum algorithms for Gaussian processes have been studied in (Zhao et al., 2019;Das et al., 2018a). However, the focus there is speeding up GP regression and not devising a quantum neural network. Further, none of these works discusses the connection between quantum optical Gaussian states and Gaussian processes.

Probabilistic Numeric Neural Networks
We start by reviewing the framework of probabilistic numeric neural networks for classifying an input signal with missing data (Finzi et al., 2020). This model uses a GP on a continuous space X = R D to interpolate the input signal and defines a neural network on this GP representation. To simplify the discussion we consider here a finite dimensional input space. This also allows us to present a simple numerical procedure to study the models, and does not change the conceptual findings of the present paper. Note that our approach is also related to (Li & Marlin, 2016) but we shall keep the term probability numeric to emphasize the connection with the field of probabilistic numerics (Cockayne et al., 2019).
corresponding to the observation y i of a signal (field) ϕ : X → R at location x i ∈ X . We assume an underlying finite input grid X = {1, . . . , l} ×d and denote by ϕ the vector with components ϕ x . We consider a prior GP with zero mean and kernel k, and compute the posterior GP(µ , k ) to interpolate the signal (Rasmussen et al., 2006): Here n is measurement noise. We then apply a sequence of linear layers and nonlinearities on the input random vector ϕ (1) ∼ GP(µ , k ): (2) In the intermediate layers, ϕ has shape N C × |X |, with N C being the number of channels, and components ϕ a,x . σ a point-wise non-linearity, B adds a bias, and A is linear.
After L transformations of type (2), the output feature undergoes a global average pooling, This implements an invertible version of the usual operation. Finally, the first element of ϕ (L+1) is passed through a linear layer that produces an output whose mean is interpreted as logits for classification. Denoting L = B • A, the chain of operations of the PNCNN is

Quantum Mechanics
In this section we review some basics of quantum mechanics, in particular quantum fields. We refer the reader to e.g. (Sakurai & Napolitano, 2017) for a standard introduction. In the quantum formalism we associate to every classical configuration of a random field ϕ = {ϕ x } x∈X ∈ R |X | (for simplicity we ignore the channels in this section) a vector |ϕ . The span of all these vectors is a vector space H, with elements given by superpositions: Note that ψ(ϕ) are the coefficients ∈ C that are used to combine vectors |ϕ . The bra-ket notation might feel a little strange to readers unfamiliar with quantum mechanics, and is used to represent abstract elements in H.
H is a Hilbert space equipped with a scalar product ϕ|ϕ = δ(ϕ − ϕ ), so that ϕ|Ψ = ψ(ϕ). We call a linear operator on H a quantum field (indicated with a hat). The algebra of quantum fields is generated by the pairs { ϕ x , π x } x∈X : They are self-adjoint and satisfy the canonical commutation relations Quantum dynamics needs to preserve the norm of quantum states and acts by unitary operators U t = e −it H , where H is a self-adjoint Hamiltonian. Instead of evolving states we can equivalently evolve the observables A we are going to measure: which can be computed using the Baker-Campbell-Hausdorff formula: A(t) satisfies the Heisenberg equation of motion: Measurements reduce a quantum superposition to a classical configuration 1 . This projection |Ψ → |φ occurs with probability |ψ(φ)| 2 . We shall make use below of the following more general fact about partial measurements (see e.g. (Watrous, 2018)).
We note that Bayes rule follows from proposition 3.1: a subsequent measurement of an observable on H 2 with outcome 1 One can consider measurements of non-diagonal self-adjoint fields O, which project onto an eigenstate of O, but in this work we shall restrict to measurements of diagonal operators. ϕ π m Figure 1. The covariance ellipse 1 2 (z − m) T C(z − m) ≤ 1, also known as "quantum blob" (de Gosson, 2006), in a 2d phase space z = (ϕ, π). The area is proportional to Det(C) > 1 (uncertainty principle) and is preserved by (20), (21). y 2 on the state |ζ y1 will give outcome probability: which coincides with p(y 1 , y 2 )/p(y 1 ).

Gaussian States
Define the 2|X | dimensional vector of operators: After introducing the symplectic form J, (10) reads as: Gaussian states are specified uniquely by their mean and covariance defined as: Here and below we shall denote by 1, 2 the first and second |X | components related to ϕ, π sectors: Proposition 3.2. The covariance matrix C satisfies: See Appendix A for a proof. The condition C + iJ ≥ 0 encodes the uncertainty principle and distinguishes quantum Gaussian states from classical Gaussian distributions on phase space (Bartlett et al., 2012). Figure 1 gives a visualization of this fact.
We denote a Gaussian state by |m, C 2 . The wave-function ϕ|m, C of a Gaussian state is a Gaussian function of ϕ, albeit with complex quadratic form (Hudson, 1974).
We now show that the following unitary transformations, whose Hamiltonians are at most quadratic in R i , implement the most general transformations among Gaussian states: The Lie group Sp 2|X | (R) is the symplectic group of matrices satisfying SJS T = J, which include rotations and scaling. An element X in its Lie algebra, S = e X , satisfies JX + ensuring unitarity of ω(S). D and ω implement symmetry transformations (resp. translations and linear symplectic transforms 3 ) on H. The following proposition shows that R intertwines their action on H and the fundamental action on R 2|X | : Proposition 3.3. The unitaries of (20) and (21) represent symplectic affine transformations of the canonical operators This result follows from formula (11) and the commutation relations (15). See Appendix A for details. This related proposition gives the effect on the mean and covariance of Gaussian states and is also proved in the Appendix A. Proposition 3.4. Under the unitaries of (20) and (21), the Gaussian states transform as

Quantum Extensions of Probabilistic Numeric NNs
We introduce in this section a series of quantum operations that generalize the classical layers of a probabilistic numeric NN.

State Preparation
We start by showing how to perform Bayesian inference with Gaussian states. From 3.1 we have the following result. Lemma 4.1 (Quantum GP Inference). Let |0, C be a Gaussian prior state such that C 11 we consider the posterior: 3 More precisely, ω is a representation of the metaplectic group, a double cover of the symplectic group (de Gosson, 2006).
where |y x is an eigenstate of ϕ x . We have: with µ , k as in (1).
Proof. It follows directly from Prop. 3.1 with |Ψ = |0, C , H 1 the space corresponding to locations x i and the formulas for the Gaussian conditionals (Rasmussen et al., 2006).
The quantum GP inference step allows one to encode a classical signal in a quantum state in such a way that quantum entanglement represents an agent's uncertainty about discretization errors.

Quantum Linear Layers
Next, we show how to perform the quantum equivalent of a linear layer that acts on the quantum fields R in the same way as a classical linear layer acts on a classical field R. From Prop. 3.3 we have the following natural definition. Definition 4.1. A quantum linear layer is the unitary: where D(ξ) and ω(S) generalize the bias and multiplication by the weight matrix respectively.

Quantum Non-Linearity
The definition of a unitary operator whose time evolution correspond to the action of a non-linearity is more involved, and in fact has remained elusive in the previous quantum neural network literature. Similarly to classical non-linearities, a quantum non-linearity acts pointwise on the quantum fields. This restricts the associated Hamiltonian to be x,a H x,a , where H x,a acts non-trivially only on the quantum fields at x, a. As a design principle, we consider the following class of time evolutions which map ϕ x,a to a function σ( ϕ x,a ). (Recall that a function of an operator is defined by its Taylor series.) This will provide a simple way to embed classical neural networks in our framework. Proposition 4.2. Under the time evolution generated by the fields evolve according to the equations of motion: where h.c. means the Hermitian conjugate of the expression preceding it.
Proof. The equations of motion follow from (12) upon us- The following proposition relates σ to f (the proof checks the time derivatives and is in Appendix B). Proposition 4.3. The ODEs (32), (33) have solutions where Prop. 4.3 gives a general, albeit implicit, solution to the problem of constructing a quantum non-linearity. Next, we give an explicit solution for the case of softplus, a popular smooth version of ReLU. Lemma 4.4 (Quantum Softplus). The softplus non-linearity with temperature β, corresponds to time evolution from time 0 to time 1 under Proof. This Hamiltonian is of the general form (30) with The antiderivative of its inverse, βe +βx , is F (x) = e βx . So F −1 (y) = β −1 log(y) and the result then follows immediately from Prop. 4.3.
We now move on to discussing architectures and symmetries of quantum neural networks. We shall resume the study of non-linearities in section 7 where we will look at a semiclassical limit of the quantum models.

Architecture
We define a quantum neural network by putting these pieces together: where U σ and U lin are as in Prop. 4.2 and Def. 4.1, while is a global average pooling operator, where P is as in (5).
To make a prediction we proceed similarly to (Finzi et al., 2020). We discard the spatial locations that have not been aggregated over by averaging, act with a final linear classifier and finally measure the means ϕ c for the c = 1, . . . , C channels, C being the number of classes: These are interpreted as the logits for the classification task at hand. For simplicity we introduced the state |ζ D in lemma 4.1 without referring to the channel dimension. To make sense of the action of U NN on it, we add extra registers for the channel dimension which are initialized to the vacuum state, which is a Gaussian state with zero mean and unit covariance (Adesso et al., 2014). While we considered here only a global average pooling at the end following (Finzi et al., 2020), it is possible to pool features in intermediate layers as well by simply discarding registers associated to the modes to be discarded. Computing the logits of (40) is in general intractable classically and requires a quantum computer.
Our definition of a quantum neural network is similar to that of (Killoran et al., 2019). However, w.r.t. that work we introduce the following two main novelties: 1) we use Gaussian states for data interpolation (lemma 4.1); 2) we use a unitary gate implementing the non-linearity (see Sec. 4.3). (Killoran et al., 2019) discusses two sets of non-linearities: the first uses a quantum channel to implement a non-linearity of the type ϕ x,a → σ( ϕ x,a ), but, w.r.t. our setting, that implementation requires to double the number of quantum registers and then discard half of them. 4 The second strategy is to employ Hamiltonians such as H cubic = ϕ 3 a,x or H Kerr = ( ϕ 2 a,x + π 2 a,x ) 2 as non-linearities, since either of them, together with the unitaries of (20) and (21), form a simple set of universal gates for quantum computation with continuous variables (Lloyd & Braunstein, 1999), i. e. any other quantum gate can be expressed in terms of those. 5 While requiring minimal resources, these choices are likely to perform worse than our choice. Indeed in the semiclassical limit discussed in section 7, these correspond to low degree polynomial non-linearities, which are not efficient for classical neural network approximation (Pinkus, 1999).

Symmetries
Generally, symmetries in classical neural networks are realized as linear maps g ∈ G that acts on the activations ϕ as ρ(g)ϕ, where ρ is a representation matrix. On top of translations, prominent examples of G in ML are rotations (Cohen et al., 2018) and permutations (Maron et al., 2018). Having replaced the linear action on activations with ω, we define unitary representations of G on quantum states by ω(S g := ρ(g)⊕ρ * (g)), where ρ * (g) = ρ(g −1 ) T is the dual representation, ensuring symplecticity of S g . We have For example, in case of translations along the µ ∈ {1, . . . , d} axis, (τ µ ϕ) a,x = ϕ a,x+eµ , we have S τµ = τ µ ⊕ τ µ , which translates both ϕ, π variables equally.
Equivariance of a quantum linear layer ω(S) now amounts to the commutation relations: where the second formula follows from the group homomorphism property: ω(S) ω(S ) = ω(SS ). The characterization of symmetries presented here completely solves the problem of designing equivariant quantum linear layers by reducing the problem to designing equivariant classical linear layers with symplectic weight matrices S which commute with ρ(g) ⊕ ρ * (g). As in the classical case, since the non-linearities act pointwise, they will be invariant under operations that permute the coordinates ( ϕ x,a , π x,a ) → ( ϕ x ,a , π x ,a ), such as spatial symmetries, ensuring equivariance of the whole architecture.
As an illustration, the condition to be a convolution. Similarly, for rotations G = SO (3), irreps are self-dual, so S g = ρ(g) ⊕ ρ(g) and each block of S is a group convolution (Cohen et al., 2018).

Hintons
We discuss now a particle interpretation of the formalism introduced for neural networks. We introduce the operators: We can then use a particle basis (a.k.a. Fock space) for H, where we identify a zero mean and unit covariance Gaussian state with no particles (vacuum) |Ω = |0, 1 such that b |Ω = 0 and create an orthogonal basis by acting with different monomials ( b † a1,x1 ) ni · · · ( b † am,xm ) nm on |Ω , n i being the number of particles at channel a i and location x i . In the quantum optical setting, the particles are called photons. In the neural network context we dub these fundamental excitations "Hintons" after G. Hinton, a founding father of the field of deep learning.
In the next sections, we will analyze tractable limits of the quantum neural network defined so far to get more insights into the architectural design.

The Case of Classical Probabilistic NNs
Here we show how a classical probabilistic NN can be embedded in the quantum model of Sec. 4.4. First we prove the following quantum representation of the push forward of a GP under a generic classical (invertible) map (proof in the Appendix C).
Lemma 6.1. Let |m, C be a Gaussian state and U a unitary such that: Then: where f #p denotes the push forward of p under f .
We already know that U σ has this property, so we only need to constrain the quantum linear layers so that they do not mix ϕ with π.
Theorem 6.2. Consider the quantum network (38) with: We have the quantum-classical duality: where the rhs is the push forward of the GP posterior (1) under the map Φ of (6). That is, the logits of Eq. (40) computed by the quantum neural network coincide exactly with those computed by the probabilistic numeric NN of section 3.1 with weights and biases A ( ) , ξ ( ) .
Proof. Recalling the linear transformation law of Prop. 3.3 and the non-linear action of U σ , the proof then follows from lemma 6.1 with F = Φ.

The Semiclassical Limit
Let us denote by R = (ϕ, π) the classical fields corresponding to the quantum operators introduced above. Recall that notationally we distinguish between classical and quantum fields by the absence or presence of a hat.
Note that despite producing entangled states, the quantum linear layers acting on Gaussian states can be simulated efficiently on a classical computer as the action amounts to the matrix multiplications of Prop. 3.3. See also (Bartlett et al., 2002). In fact, at the linear level the only difference between a quantum evolution and a probabilistic classical evolution of a Gaussian Liouville measure in phase space is the covariance condition C + iJ ≥ 0 (Prop. 3.2) coming from the non-commutativity of position and momenta in quantum mechanics. We remark also that this condition is preserved by classical evolution thanks to the symplectic nature of classical mechanics (De Gosson, 2009). In Sec. 6 we showed how restricting the linear layers to block diagonal matrices led to the classical model, which corresponds to the push forward of an initial Gaussian Liouville distribution under a neural network, but this time only involving the ϕ field. In this section, as an intermediate step towards studying the full quantum model, we change the non-linearity in such a way that the modified model corresponds to the push forward of a initial Gaussian measure under a neural network, involving both the ϕ and the π fields. We can interpret the resulting model as a semi-classical limit of the quantum model since it uses elements of quantum mechanics (uncertainty relation for the covariance) as well as classical mechanics (for the non-linearity). Conceptually, we have the inclusion of models as special cases depicted in figure 2.

Quantum
Semiclassical Classical classical nonlinearity block diag. linear We define the semiclassical model by simply replacing U σ with a classical Hamiltonian evolution under which the phase space measure evolves into a new classical phase space measure, which corresponds to the classical limit of the equation of motion (Sakurai & Napolitano, 2017).
Proposition 7.1. Under the classical time evolution generated by the Hamiltonian the fields transform as: where F (x) = 1/f (x).
Proof. The proof follows by simply checking the classical equations of motion: We note that classical and quantum equations of motions and solutions ( In particular the following is the classical counterpart of lemma 4.4. Proposition 7.2 (Symplectic Softplus). Replacing operators with classical variables in lemma 4.4 we have We then define a neural network that pushes forward the input Gaussian Liouville distribution on phase space to an output distribution p out by alternating linear layers with non-linear classical layers U σ . Analogously to the original construction of (40), its mean is then used as the logits for classification: We implemented and tested the semiclassical neural network with symplectic softplus non-linearity. The experiments are performed on a simple classification task for irregularly sampled time series. As already remarked, we do not claim that the semi-classical network performs better than a classical probabilistic numeric neural network with similar capacity. We therefore emphasize that the purpose of our experiments is simply to check that after adding momenta and using a different non-linearity than standard, the model performs similarly to a classical baseline. In the implementation we saw that taking β = 0.1 for the non-linearity avoids numerical instabilities for negative ϕ due to the presence of e −βϕ in (52). We refer to Appendix D for details of architecture, task and learning algorithm tested. We leave the study of the fully quantum case as an outstanding open problem which will most likely be tackled only when error corrected quantum optical computers will become available.

Quantum Optical Implementation
We here describe practical considerations for implementing our model on an optical quantum computer. Optical devices are attractive since one can implement matrix multiplication  Fig. 3. While the implementation of linear layers is common to these works, the quantum nonlinearities and the data embedding in Gaussian states are original to our work.

State preparation
The state preparation step encodes the input signal D = {(x i , y i )} N i=1 into quantum registers. We compute the GP posterior as in (1) for a set of points x ∈ X and then we create an input Gaussian state by acting on the vacuum state |Ω defined in Sec. 5 with the linear layer D(ξ = (µ , 0)) ω(S = A ⊕ A −1 ), where A is a square root of k . The vacuum state can be created using lasers (Nielsen & Chuang, 2000) and we defer to the next section the implementation of the linear layer. This procedure incurs a complexity similar to the classical GP inference, that is O(N 3 ). In Appendix E.1 we give further comments on the state preparation step. We see developing a more efficient quantum implementation of quantum GP inference as in interesting future direction.

Linear layer
The implementation of a unitary ω(S) on a quantum optical computer is a well studied problem and we limit ourselves here to simply say that it can be done by a sequence of beamsplitters and phase shifters. Appendix E.2 provides more details of the procedure and its complexity.

Non linearity
Quantum computers can perform arbitrary computations if given a set of universal gates. For quantum optical com-puters, one can take the quadratic Hamiltonians and the cubic Hamiltonian ϕ 3 (Lloyd & Braunstein, 1999). We can implement a non-linearity with Hamiltonian (30) by approximating it with the truncation of the Taylor series of the function f to order k: We can then use the standard procedures for quantum simulation (Nielsen & Chuang, 2000) to implement e i H (k) . This is detailed in Appendix E.3, where we derive an explicit decomposition in terms of gates which use only the universal Hamiltonians π x,a , π 2 x,a and ϕ 3 x,a . This fully specifies a protocol to implement our model on a quantum optical device.

Conclusion and Outlook
In this paper we introduced a new formalism for deep learning in terms of quantum fields. We made the following contributions: 1) we showed how to use Gaussian states for Bayesian inference, 2) we devised unitary operators that implement standard non-linearities, 3) we presented tractable limits of the quantum network, 4) we discussed how to implement our models on a quantum computer.
We implemented the semi-classical architecture to check that it performed as expected (which it did). But because we do not expect that these tractable limits perform any better than classical models we choose not to benchmark them against each other.
Exciting directions for the near future are: studying approximate solutions that get closer to the full quantum model; finding efficient ways to do quantum GP inference on quantum hardware; developing further quantum non-linearities and the quantum formalism for classical models.

A. Proofs of Section 3
Proposition A.1. The covariance matrix C satisfies: Proof. We give the proof for |X | = 1 and refer the reader to (Simon et al., 1987) for the general case. Suppressing the x index and denoting · = ψ| · |ψ , we have: where σ 2 ϕ = ϕ 2 − ϕ 2 and similarly for π. For |X | = 1, This is the uncertainty relation in the stronger form due to Robertson-Schroedinger, proving the statement for |X | = 1.
Proposition A.2. The unitaries of (20) and (21) represent symplectic affine transformations of the canonical operators Proof. The result then follows from the Baker-Campbell-Hausdorff formula and the commutators: Under the unitaries of (20) and (21), the Gaussian states transform as Proof. From proposition 3.3 and the definitions (16) and (17), we have for any |Ψ : The result follows by specifying |Ψ = |m, C .
Proof. We can check directly the formulas by differentiating w.r.t. t to show that the time evolved fields satisfy the equation of motions. Rewriting the first equation as F ( ϕ x,a (t)) = F ( ϕ x,a (0)) + t and differentiating the l.h.s.: which equals ∂ t (F ( ϕ x,a (0)) + t) = 1, showing that ϕ x,a (t) satisfies (32). For the second equation we differentiate the first term in the parenthesis: which shows that π x,a (t) satisfies (33).

C. Proofs of Section 6
Lemma C.1. Let |m, C be a Gaussian state and U a unitary such that: Then: where f #p denotes the push forward of p under f .
Now we note the representation of a projector as the average: a formula which can be proved by comparing matrix elements of the two operators in arbitrary states. Then we have: In going from the third to the fourth line we used the following representation, valid for any f : and we have used the definition of push forward of a distribution, i.e. that if x ∼ p X , then p Z = f #p X is the distribution of the random variable z = f (x) which can be obtained explicitly via the change of variable formula: where the second equality follows upon changing variables to x = f (x ).
We also present an alternative proof the lemma, which relies on: The fact that U |x = c(x) |F (x) for some c(x) follows from and the delta normalization fixes the proportionality factor: using that for a g(x) having the unique solution g(x) = 0 at x 0 , one has Given equation (81), the lemma is immediate: which coincides with the definition of push forward.

D. Experiments
In this section we discuss numerical experiments to show that the semi-classical neural networks introduced in Sec. 7 performs on par with its classical counterpart for a simple irregular time series classification task.

D.1. Datasets
Similarly to (Li & Marlin, 2015), to perform controllable experiments with irregular data, we select the UCR time series dataset (Chen et al., 2015). We assume that the input time series is uniformly sampled in [0, 1] and then sample randomly a fraction f of the input time series, where f varies from 1 (no subsampling) to .1 (only 10% of the data is retained). This gives us coordinates x i and values y i that can be used to prepare an input Gaussian state as discussed below in the training algorithm section. We select the same random subsampling mask for the the training data and a different one for the test data, which we expect is a good model of a realistic irregular measurement and a more challenging machine learning task.
Since our algorithm has cubic complexity in number of points of the input signal and the purpose of the experiments is only to validate that the semiclassical architecture can perform on par with a classical probabilistic numeric baseline, we restrict the focus here only to the following dataset which has small number of data points: SyntheticControl (train: 300, test: 300, classes: 6, length time series 60).

D.2. Training Algorithm
The training algorithm is composed of two steps. First, we train the GP using Alg. 1. Then we use the resulting posterior mean and covariance -more precisely the posterior mean µ , the cholesky decomposition of the posterior covariance L and its inverse M as inputs to the PNN and SPNN. The PNN discards M , which is used for the momenta sector, which is decoupled for PNN and does not intervene in the prediction. Then we use µ , L, M to sample from a normal distribution that is defined as in Alg. 2, i.e. we set all the other means to zero and the covariance is k −1 = M M T for the momenta sector and channel 0 while identity elsewhere. This returns a batch of samples z i that are then passed through the neural network as in Alg. 3 to finally be averaged to produce the logits used for classification. The next section describes the architectures in more details.
Note that GP inference is the bottleneck here and we could use the methods of (Li & Marlin, 2016) to improve on that. Our interest is however not to have a state of the art classifier, . . , B (B being data size), grid X , GP kernel parameters Θ. Train GP: Compute Θ by maximizing marginal likelihood of data in mini-batches: Input:µ , L, M (output of Alg. 1), number of channels N C , number of sample N samples . Create Normal: N (m, C) with : N c × |X | × 2 and C: Sample: z i from N (m, C) using the reparametrization trick for i = 1, . . . , N samples Output: {z i } so we test the model in the simplest setting of O(N 3 ) exact GP inference and simply validate that the performance of the semi-classical neural network architecture is on par with a classical counterpart. Note that separating the GP training from the neural network training allows us to amortize the GP inference as a preprocessing step and have an O(N 2 ) algorithm for network propagation as in BNN case.

D.3. Architecture
For the GP to interpolate the data we take Matern kernel with ν = 0.5 as implemented in gpytorch (Gardner et al., 2018). We compare three models: 1) a baseline MLP (BNN), 2) a probabilistic numeric MLP (PNN) and 3) a semiclassical version of the PNN (SPNN). All models share the same basic architecture, which is a MLP with three layers. In all cases the hidden sizes are equal to the number of points in Output: Logits l c the original grid of the data, |X |, times the number of channels N C . For BNN we parameterize the weight matrices directly i.e. the learnable parameters are the weight matrices. Instead for PNN and SPNN we parametrize the logarithm of the weight matrix as the learnable parameters. We did check that for PNN this did not impact perform comparing to the case where we use the same parameterization of BNN. We do this to ensure simply that the weight matrices of the SPNN are symplectic. Indeed recalling the definitions from the main text, see Sec. 3.2, a linear layer with input/output dimensions 2N has learnable parameters A, B, C ∈ R N ×N and the weight matrix S is then constructed as: A, B, C are denoted h in Alg. 3. For PNN we simply take B = C = 0 and discard the lower right block. We have also tried to restrict S to be a free matrix and use the parametrization from (de Gosson, 2006) but found that the matrix exponential (which allow use to parametrize more general symplectic matrices) worked best. We used pytorch (Paszke et al., 2019)

D.4. Results
We present results in table 1. We pretrained the GPs using Adam optimizer with learning rate 0.1 from (Paszke et al., 2019) for 20 epochs using batch size of 50. We trained the neural network using default Adam optimizer for 1000 epochs with batch size 50 as well. We used N samples = 100 and N C = 2.
The BNN baseline as well as PNN and SPNN perform on the original unsampled dataset similarly to the MLP used in (Wang et al., 2016) (95 for SyntheticControl), which also has three layers but additionally dropout and uses ReLU instead of softplus. This validates the choice of architecture we made.
The results presented are on the average over three random seeds for the random sampling. The error bars are considerable in all cases due to the fact that different sampling masks lead to different data points that are retained and this choice affects the classification of the signal. When we start to subsample the data (using the same procedure for all models), both the PNN and SPNN n average perform better than the BNN, confirming the findings of (Li & Marlin, 2015;2016;Finzi et al., 2020) showing the the usefulness of the GP model for irregularly sampled data. (Li & Marlin, 2015) report accuracies averaged over all the 43 UCR datasets and not directly comparable with our restricted setting, while (Li & Marlin, 2016;Finzi et al., 2020) use different datasets for which scalable methods for GP inference are required. Finally, we note that the performance of PNN and SPNN have a considerable overlap within the confidence region of the results, validating the claim that they perform similarly.

E. Details of quantum optical implementation E.1. State preparation strategies
The state preparation step encodes the input signal D = {(x i , y i )} N i=1 into quantum registers. There a several options for that. First we note that while the prescription of lemma 4.1, i.e. to create a prior state and to project onto y i the registers associated to locations x i , is appealing conceptually and allows us to formulate inference entirely in quantum language, it is not straightforward to implement. This is because projections are enacted by partial measurements and those give outcome y i only with a small probability of success. To overcome this, we consider the following alternatives. The simplest is the one explained in the main text, that we repeat here. We use classical hardware to compute the GP posterior as in (1) for a set of points x ∈ X . Then we create a an input Gaussian state by acting on the vacuum state |Ω defined in Sec. 5 with the linear layer D(ξ = (µ , 0)) ω(S = A ⊕ (A) −1 ), where A is a square root of k . The vacuum state can be created using lasers (Nielsen & Chuang, 2000) and we defer to the next section the implementation of the linear layer. While this procedure is straightforward it incurs a complexity similar to the classical GP inference, that is O(N 3 ). Another alternative is to use the fact that on quantum hardware we can perform matrix inversion exponentially faster (Harrow et al., 2009). Under some conditions (in particular access to QRAM) one can compute posterior mean and covariance in polylog time (Zhao et al., 2019;Das et al., 2018a). Under similar conditions, quantum singular value decomposition can compute A in polylog time (Kerenidis & Prakash, 2016). The catch in our setting however is that to read out the values of µ and A we would still need O(N 2 ) operations. We leave developing a more efficient quantum implementation of quantum GP inference as in interesting future challenge.

E.2. Linear layer
The implementation of a unitary ω(S) on a quantum optical computer is a well studied problem and we limit here to give a high level description, see (Weedbrook et al., 2012) for more details. We first decompose the unitary in terms of elementary linear optical gates, we use the group homomorphism property ω(S) ω(S ) = ω(SS ) together with the Bloch-Messiah decomposition S = KΣL with K, L symplectic and orthogonal and Σ = diag(e r1 , . . . , e r M , e −r1 , . . . , e −r M ). Σ can be implemented directly using optical parametric amplifiers. The orthogonal matrices K, L can be further decomposed using Givens rotations as product of rotations of two components and implemented in terms of beamsplitters and phase shifters. Note that the number of gates required to decompose an arbitrary matrix can grow quadratically with the dimension (Reck et al., 1994). The bias D can be also easily