Dimensionality Reduction with Variational Encoders Based on Subsystem Purification

Abstract

Efficient methods for encoding and compression are likely to pave the way toward the problem of efficient trainability on higher-dimensional Hilbert spaces, overcoming issues of barren plateaus. Here, we propose an alternative approach to variational autoencoders to reduce the dimensionality of states represented in higher dimensional Hilbert spaces. To this end, we build a variational algorithm-based autoencoder circuit that takes as input a dataset and optimizes the parameters of a Parameterized Quantum Circuit (PQC) ansatz to produce an output state that can be represented as a tensor product of two subsystems by minimizing $Tr({\rho }^2)$. The output of this circuit is passed through a series of controlled swap gates and measurements to output a state with half the number of qubits while retaining the features of the starting state in the same spirit as any dimension-reduction technique used in classical algorithms. The output obtained is used for supervised learning to guarantee the working of the encoding procedure thus developed. We make use of the Bars and Stripes (BAS) dataset for an 8 × 8 grid to create efficient encoding states and report a classification accuracy of 95% on the same. Thus, the demonstrated example provides proof for the working of the method in reducing states represented in large Hilbert spaces while maintaining the features required for any further machine learning algorithm that follows.

1. Introduction

Variational quantum algorithms in the NISQ [1] era provide a promising route toward developing useful algorithms that allow for optimizing states in higher dimensional spaces by tuning the polynomial number of parameters. The most prominent techniques within variational methods include the Variational Quantum Eigensolver (VQE) [2], Quantum Approximate Optimization Algorithm (QAOA) [3], and other classical machine learning-inspired ones. We ask the readers to refer to [4] for an exhaustive study on quantum machine learning with applications in chemistry [5], physics [6], supervised learning and optimization [7]. Within the context of optimization and machine learning in general, some of the major problems that need to be addressed include encoding classical data, finding an expressible enough ansatz (Expressibility) [8], and efficiently computing gradients (Trainability) [9], generalizability [10]. These problems are interlinked and thus not treated independently in general.

As we move away from the NISQ era toward deep parameterized quantum circuits (PQC), one of the major problems with regard to trainability that needs addressing is the problem of vanishing gradients referred to as barren plateaus [11]. This might be an effect of working with a large number of qubits [11], expressive circuit ansatz [12], noise inducement [13] or the use of global cost functions in the learning [14]. Having efficient procedures to reduce the dimensionality the representation of the input quantum state helps create efficient encoding schemes that could later be used as inputs to other machine learning algorithms, where the cost functions on higher dimensional spaces with expressive ansatz are less likely to be trainable. To this end, we develop machine learning techniques that allow for compact representations of a given input quantum state.

Within the classical machine learning community, autoencoders have been effectively used to develop low-dimensional representation of samples generated from a given probability distribution [15]. Inspired by these techniques, work on Quantum Autoencoders [16,17] have allowed for people to develop compact representations against a fixed finite state. It is not clear that such tensor product states with a fixed finite state are always possible while retaining the maximal possible information. Here, we show that if one were to relax the condition toward maintaining a fixed finite state, a better compact representation can be generated that can be post-processed toward classification. We develop techniques to create subsystem purifications for a given set of inputs and follow it by creating superpositions of these purifications indexed using the subsystem number. This representation is further used for doing classification achieved by applying variational methods over parameterized quantum circuits restricted to this compact representation and showing the learning of the method. We apply an ansatz to create subsystem purification on the Bars and Stripes (BAS) dataset and show that one can reduce the number of qubits required to represent the data by half and achieve a 95% classification accuracy on the Bars and Stripes (BAS) dataset. The demonstrated example shows proof for the working of the method in reducing states represented in large Hilbert spaces while maintaining the features required for any further machine learning algorithm that follows. The scheme thus proposed can be extended to problems with states in large Hilbert spaces where dimensionality reduction plays a key role with regard to the trainability of the parameterized quantum circuit.

2. Method

Given an ensemble of input states, $E=\left\{|{\psi }_i\rangle \right\}$, the objective is to construct a low-dimensional representation of states sampled from this distribution E. Let $|{\psi }_i\rangle$ be a state over $n_A+n_B$ qubits. We design a protocol that allows for us to create an equivalent compact representation of $|\psi \rangle$ with $max(n_A,n_B)+1$ qubits. To simplify the discussion, let us assume that $n_A=n_B$, and thus, we create a representation using half the qubits. We do this in 2 stages.

Stage 1:

In the first stage, we apply a unitary $U(\overrightarrow{\theta })$ that decomposes ${|{\psi }_i\rangle }_{A,B}$ into ${|\alpha {(\theta )}_i\rangle }_A\otimes {|\beta {(\theta )}_i\rangle }_B$. To produce such a tensor product structure, we could minimize the entropy on either subsystem A or B till a zero entropy is achieved. Thus, we could optimize the cost function,

$$C_B(\overrightarrow{\theta })={〈S\left(t{r}_A\left[U(\overrightarrow{\theta }){|\psi \rangle }_{AB}{\langle \psi |}_{AB}{U}^{†}(\overrightarrow{\theta })\right]\right)〉}_{\left\{|\psi \rangle \right\}}$$

(1)

where $t{r}_A$ represents the tracing operation over the qubits of subsystem A, $<.{>}_{\left\{|\psi \rangle \right\}}$ represents the averaging over the $\left\{|\psi \rangle \right\}$, and $S(\rho )=tr(\rho log(\rho ))$ is the entropy of a given density matrix $\rho$. The cost function $C_B(\overrightarrow{\theta })$ attains a maximum value equal to $log(n_B)$ when ${\rho }_B$ is maximally mixed and equal to 0 when ${\rho }_B$ is in a pure state. Figure 1 shows a schematic representation of the ansatz used for $U(\theta )$. The hardware efficient ansatz used here is restricted to $R_y$ gates on each qubit followed by a ladder of CNOT operators that couples adjacent qubits. This is sufficient for our purpose as the input data vectors used here, associated with bars and stripes, are vectors with real amplitudes.

Variational quantum algorithms have been studied in the past to create thermal systems by minimizing the free energy of the output state [18,19]. The main problem tackled in these papers involves developing techniques that allow one to compute the gradients of Entropy required to be optimized over the training. The issue arises from not having exact representations that can compute the logarithm of a given density matrix efficiently. Furthermore, to avoid numerical instabilities in the entropy function arising from the density matrix of pure states being singular, here, we alternatively maximize over the cost function,

$$C_{AB}(\overrightarrow{\theta })=\frac{1}{2}{〈T{r}_A({{\rho }_A}^2)+T{r}_B({{\rho }_B}^2)〉}_{\left\{|\psi \rangle \right\}}$$

(2)

where ${\rho }_A=T{r}_B(U(\overrightarrow{\theta }){|\psi \rangle }_{AB}{\langle \psi |}_{AB}{U}^{†}(\overrightarrow{\theta }))$ and ${\rho }_B=T{r}_A(U(\overrightarrow{\theta }){|\psi \rangle }_{AB}{\langle \psi |}_{AB}{U}^{†}(\overrightarrow{\theta }))$. $C_{AB}$ attains a maximum value of 1 when ${\rho }_A$ or ${\rho }_B$ are pure states resulting in $Tr({\rho }_{A/B}^2)=Tr({\rho }_{A/B})=1$ and attains a least value $2/2^n$. $C_{AB}$ is a convex function as $Tr({\rho }_{A/B}^2)$ is a convex function of ${\rho }_{A/B}$ [18]. Figure 2 shows a schematic representation to compute $Tr({\rho }^2)$ using a destructive swap test. The optimization landscape thus has one local minimum dictated by the expressivity of the ansatz used to capture it.

The parameters $\overrightarrow{\theta }$ are variationally optimized to obtain ${\overrightarrow{\theta }}^{*}={\mathrm{argmax}}_{\overrightarrow{\theta }}{C}_{AB}(\overrightarrow{\theta })$. If $C_{AB}(\overrightarrow{\theta })$ reaches an optimal value of zero, we can express ${|\psi \rangle }_{AB}={|\varphi \rangle }_A⨂{|\eta \rangle }_B$, thus expressing a state with $2^{n+m}$ degrees of freedom, effectively using $2^n+2^m$ degrees of freedom. Having expressed the input state as a tensor product of subsystems, we now move to stage 2 of the algorithm.

Stage 2:

Note that the above representation still makes use of $2n$ qubits to capture the features of $|\psi \rangle$. If the subsystems are not equal in size, additional ancillary qubits are included to match the system size. We now show how this representation can be compressed using $n+1$ qubits. To carry this out, we apply an ancillary CSWAP (controlled swap/Fredkin) gate acting on the qubits of systems A and B. The additional ancillary qubit works as an index register for these states. Thus, we obtain $|0\rangle |\varphi \rangle |\eta \rangle +|1\rangle |\eta \rangle |\varphi \rangle$. If $|\eta \rangle$ and $|\varphi \rangle$ are not orthogonal states, then there exists at least one basis element $|g\rangle$ in the computational basis with a nonzero coefficient in both these states. Without a loss of generality, let us assume that the measurement collapses onto $|g\rangle$, giving rise to $\frac{1}{\sqrt{1+c^2}}(|0\rangle |\varphi \rangle +c{e}^{i\alpha }|1\rangle |\eta \rangle )\otimes |g\rangle$, where c and $\alpha$ are real numbers. The factor $c{e}^{i\alpha }$ is generated from the relative difference in the coefficients of the state corresponding to $|g\rangle$. Figure 3 shows a schematic representation of the main steps involved in creating a superposition, with the ancilla register being used as an index to the subsystem outputs of Stage 1.

Output:

Thus, we have successfully managed to map the input state $|\psi \rangle$ to $\frac{1}{\sqrt{1+c^2}}(|0\rangle |\varphi \rangle +c{e}^{i\alpha }|1\rangle |\eta \rangle )$, apart from an arbitrary relative phase $\varphi$ and amplitude c in the representation. Note that this procedure is reversible, thus preserving all information content encoded into input state $|\psi \rangle$. To show that it is reversible, one just needs to take two copies of the output state $\frac{1}{\sqrt{1+c^2}}(|0\rangle |\varphi \rangle +c{e}^{i\alpha }|1\rangle |\eta \rangle )$, measure the corresponding ancilla to project out $|\varphi \rangle |\eta \rangle$, and then apply the inverse of $U(\overrightarrow{\theta })$, giving back $|\psi \rangle$. Thus, the encoding scheme allows for us to create a representation of the input state $|\psi \rangle$ with $2n$ qubits into only $n+1$ qubits. We now show that the presence of the arbitrary phase and relative coefficient can be ignored if we work with an ansatz with a specific structure for an L2 cost function, which can be generalized to other cost functions as well.

Let $V(\overrightarrow{\alpha })$ be the ansatz used for classification after creating the compact state representation. We calculate the label corresponding to each state by averaging the expectation value across an ensemble of representative states for each datapoint obtained via projection of the second state. The expected label for a datapoint indexed by t is thus given by

$$\mathrm{Expected}\mathrm{Label}=\sum _{i}P(i)\langle {\tilde{\psi }}_{i,t}|V^{†}(\overrightarrow{\alpha })[Z\otimes I^{\otimes n}]V(\overrightarrow{\alpha })|{\tilde{\psi }}_{i,t}\rangle$$

(3)

where i indexes the projection of the second subsystem and $P(i)$ is the probability of that projection. Using the L2 norm, the classification cost function with respect to this ansatz is given by

$$\mathrm{Classification}\mathrm{Cost}=\sum _t{(l_t-\sum _{i}P(i)\langle {\tilde{\psi }}_{i,t}|V^{†}(\overrightarrow{\alpha })[Z\otimes I^{\otimes n}]V(\overrightarrow{\alpha })|{\tilde{\psi }}_{i,t}\rangle )}^2$$

(4)

where $|\psi \rangle =\frac{1}{\sqrt{1+c^2}}(|0\rangle |\varphi \rangle +c{e}^{i\alpha }|1\rangle |\eta \rangle )$. Note here $V(\overrightarrow{\alpha })$ is an ansatz that acts on $n+1$ qubits. The cost expression chosen abive can be re-expressed on the non-projected state as follows:

$$\mathrm{Classification}\mathrm{Cost}=\sum _t{(l_t-\langle \psi |V^{†}(\overrightarrow{\alpha })\otimes I^{\otimes n}[Z\otimes I^{\otimes n}\otimes I^{\otimes n}]V(\overrightarrow{\alpha })\otimes I^{\otimes n}|\psi \rangle )}^2$$

(5)

where $|\psi \rangle =|0\rangle |\varphi \rangle |\eta \rangle +|1\rangle |\eta \rangle |\varphi \rangle$. The above extension is supported by the presence of $P(i)$ in the original expression. We now choose an ansatz that allows us to get rid of the effects of arbitrary relative phase and amplitude via the projection of the second subsystem in our original expression. Let $V(\overrightarrow{{\alpha }_1},\overrightarrow{{\alpha }_2})=|0\rangle \langle 0|\otimes V_0({\overrightarrow{\alpha }}_0)+|1\rangle \langle 1|\otimes V_1({\overrightarrow{\alpha }}_1)$, where $\overrightarrow{\alpha }=(\overrightarrow{{\alpha }_0},\overrightarrow{{\alpha }_1})$. Thus we obtain

$$\mathrm{Classification}\mathrm{Cost}=\sum _t{(l_t-\langle \varphi |V_0^{†}(\overrightarrow{{\alpha }_0})[Z\otimes I^{\otimes n-1}]V_0(\overrightarrow{{\alpha }_0})|\varphi \rangle -\langle \eta |V_1^{†}(\overrightarrow{{\alpha }_1})[Z\otimes I^{\otimes n-1}]V_1(\overrightarrow{{\alpha }_1})|\eta \rangle )}^2$$

(6)

The above expression thus obtained is akin to the averaging over the ensemble for each datapoint producing a state that is oblivious to the relative amplitude and phase factor produced via the projection. For both $V_0(\overrightarrow{{\alpha }_0})$ and $V_1(\overrightarrow{{\alpha }_1})$, we use the hardware efficient ansatz, as shown in Figure 1, which we shall employ for supervised learning in the next section.

3. Results

To demonstrate the working of the method described above, we pick a toy dataset with images of Bars and Stripes (BAS) and build a compact representation of it. The BAS dataset we consider is a square grid with either some columns being only vertically filled (Bars) or some rows being horizontally filled (Stripes) [21]. One can easily generate such a supervised dataset and realize that the distribution from which these images are sampled has a low entropy characterization. We randomly sample 1000 data points from a grid size of $16\times 16$ from the BAS dataset consisting of 131,068 datapoints represented using amplitude encoding on eight qubits. Note that, despite other classically trained and more efficient encodings that could be straightforwardly employed, here, we use amplitude encoding, as the representation is efficient in the number of qubits while remaining data-agnostic.

Applying the protocol described above, we reduce the representation of the state into a tensor product of two subsystems of equal size. Figure 4 shows the learning of optimal parameters $\overrightarrow{\theta }$ as the cost function falls. Note that the dataset has been factorized from a non-trivially entangled state. We use the standard gradient descent [22] approach in performing the training. Note that the cost function drops to zero, implying that the representation thus created is exact with a lossless transformation created by $U(\overrightarrow{\theta })$. For the 16 × 16 grid case, the ansatz $U(\overrightarrow{\theta })$ is made of D = 5 layers, while that for the 8 × 8 grid is made of D = 3 layers. At this point, we apply a layer of swap gates to reduce the eight-qubit representation of 16 × 16 grid samples into five qubits and the six-qubit representation of 8 × 8 grid samples into four qubits.

We now use this as input for performing supervised classification. We use approximately 80% of the samples from the output of the encoded samples for training and keep the remaining 20% of the samples for testing. An ansatz $V(\overrightarrow{{\alpha }_1},\overrightarrow{{\alpha }_2})=|0\rangle \langle 0|\otimes V_0({\overrightarrow{\alpha }}_0)+|1\rangle \langle 1|\otimes V_1({\overrightarrow{\alpha }}_1)$ with the same number of qubits as that of the input samples is trained, with the expectation value of pauli-Z operator being used as a label for differentiating between bars and stripes. The input image is classified as an image of bars if the expectation value is positive and as stripes image if it is negative. We use the sum of two norm errors over the dataset labels (1 for bars and $-1$ for stripes) as the cost function (6) to be minimized over, i.e,

$$\mathrm{Classification}\mathrm{Cost}=\sum _t{(l_t-\langle \varphi |V_0^{†}(\overrightarrow{{\alpha }_0})[Z\otimes I^{\otimes n-1}]V_0(\overrightarrow{{\alpha }_0})|\varphi \rangle -\langle \eta |V_1^{†}(\overrightarrow{{\alpha }_1})[Z\otimes I^{\otimes n-1}]V_1(\overrightarrow{{\alpha }_1})|\eta \rangle )}^2$$

(7)

where the summation index i labels the dataset, $l_i$ refers to the labels corresponding to the sample input and $|{\psi }_i\rangle$ is used to denote the compact representation of the state that the above encoding scheme provides. For the 8 × 8 grid, a total of 508 bars and stripe images are produced, with half of them belonging to each category. We use 400 of these samples for training and 108 samples for testing. By default, we have used 8192 shots for each measurement during the training for creating subsystem purifications and carrying out the classification. We use 1000 representative samples corresponding to each datapoint, which are characterized by an arbitrary relative phase and amplitude. Figure 5 shows the cost of optimizing the parameters of $V(\overrightarrow{\theta })$ as a function of the number of iterations. We achieved a 95% accuracy on the testing data, showing that the method used to generate the compact representation did not destroy the features of the input state.

4. Runtime Analysis of Encoding Scheme

Here we shall analytically compute the required runtime for the protocol described above. Let us assume that the input ensemble of N quantum states over n qubits supports a compact representation, allowing us to use the above protocol to encode with half the number of qubits. Let our ansatz to be optimized be made of d layers. Thus, stage 1 involves optimization over $2ndN$ parameters for N samples. Using a destructive swap test to compute fidelity with an error $ϵ$, we would require $O(1/{ϵ}^2)$ samples. Thus, the runtime complexity scales evaluate $O(ndN/{ϵ}^2)$ quantum circuits per iteration for Stage 1. Stage 2 involves sampling output states after measuring the second system in a fixed basis. If K representative samples characterized by arbitrary phase and amplitude are used for each data point, then the overall runtime is bounded by $O(KNTnd/{ϵ}^2)$.

5. Discussion and Conclusions

We discuss a scheme that allows for a compact representation of states in higher-dimensional Hilbert spaces using half the number of qubits. The output thus created serves as a good starting state for any further machine learning algorithm that might follow. The protocol is based on designing a quantum circuit that allows creating a tensor product subsystems and demonstrates results on bars and stripes datasets for 8 × 8 grid and 16 × 16 grid. We further use this output to create compact representations with half the number of qubits as compared to the starting state. To show that this representation is a lossless encoding, we use it to perform supervised learning using variational circuits on the entire dataset of an 8 × 8 grid and reproduce a 95% accuracy on the training dataset (consisting of 108 samples). Unlike quantum autoencoders where the compact representations rely on being able to optimize against a fixed garbage state, here, the relaxed restriction on the tensor product helps provide compact representations in cases where a fixed garbage state would not be feasible. The protocol described works in the absence of noise, as it requires creating pure tensor product subsystems. Further future investigation with regard to how one might overcome this hurdle is required across all quantum autoencoder protocols that aspire to create efficient codes. One might also be interested in carrying out machine learning by using weighted quantum circuits that run on the subsystems independently and compare their performance against the compact representations created thereby. One can also imagine using low entropic entangled states that stage 1 protocol outputs as input states for entanglement forging [23] and look for useful applications with them. We would like to conclude by saying that efficient methods for encoding and compression are likely to pave the way toward the problem of efficient trainability on higher-dimensional Hilbert spaces, and this work serves as a step in that direction.