Efficient quantum state tracking in noisy environments

The quantum states are encoded in the transverse spatial mode—the shape—of single photons. We describe these modes in the Laguerre–Gaussian basis $\{\vert l_i,p_i\rangle\}$ [24], where each randomly chosen state is given by the superposition $\vert \psi\rangle { = } \sum_i c_i \vert l_i,p_i\rangle$ (with $\sum_i |c_i|^2 = 1$). For the investigated qutrits in the experiment we use the three logical basis states $\{\vert 0\rangle\equiv\vert -1,0\rangle, \vert 1\rangle\equiv\vert 0,2\rangle, \vert 2\rangle\equiv\vert 1,0\rangle\}$, keeping the Gouy phase of the same order and therefore preventing the shape from rotating as the photons propagate [25]. The experiment is conducted with highly attenuated CW laser light (Thorlabs, mean photon number $|\alpha|^2 = 0.01$) together with two spatial light modulators (SLM, Meadowlark Optics) and a single photon detector (Perkin-Elmer, ${\sim} 100$ Hz dark count rate). The first SLM is used to prepare the unknown quantum state $\vert \psi_t\rangle$ via a hologram that changes phase and amplitude of the initial state $\vert l = 0,p = 0\rangle$ from the laser source. A second hologram displayed on the other SLM is then measuring the state $\vert \psi_m\rangle$ from one randomly chosen basis set (see section 2.3 for details). Finally, a single-mode fibre and the detector are acting as a local filter with the number of detected photons proportional to the overlap $\left| \langle \psi_m|\psi_t\rangle \right|^2$ between the prepared and measured mode. We note that, although all experiments were performed on a qutrit (d = 3) as proof of principle, the setup is capable of preparing and measuring much higher dimensions as shown in previous work [26]. The schematic, details, and characterisation of this simple but effective prepare-and-measure experiment are described in [26].

MEG tomography [13, 14] is an online algorithm adapted from machine learning techniques that can estimate and track quantum states. The idea is to construct an iterative procedure, in which an estimate $\hat{\rho}_t$ of the (mixed or pure) quantum state at the iteration t is updated to a more accurate estimate $\hat{\rho}_{t+1}$, given a new measurement performed on the unknown quantum state. The MEG update rule [14] for the estimate $\hat{\rho}_{t+1}$ is given by:

$$\begin{align} \hat{\rho}_{t+1} = \frac{\exp(\log(\hat{\rho_t}) - \eta_t\nabla L_t)}{\mathrm{tr} \ {\exp(\log(\hat{\rho_t}) - \eta_t\nabla L_t)}}, \end{align} \tag{ 1 }$$

where η_t is the learning rate, and L_t is the loss function,

$$\begin{align} L_t = (\mathrm{tr}(\hat{\rho}_{t}X_t) - y_t)^2, \end{align} \tag{ 2 }$$

and its gradient is given by:

$$\begin{align} \nabla L_t = 2(\mathrm{tr}(\hat{\rho}_{t}X_t) - y_t)X_t.\end{align} \tag{ 3 }$$

The pair $(X_t, y_t)$ denotes the measurement record where X_t is the measurement operator and y_t is the experimentally obtained average value of measuring the state with operator X_t . The measurements have to be informationally complete as a requirement for performing any quantum tomography procedure. Therefore, at each iteration, a measurement operator X_t is chosen at random from a complete set of bases, and together with the measurement result y_t , the estimate of the state $\hat{\rho}_{t}$ is updated to $\hat{\rho}_{t+1}$. The learning rate η_t determines the weight given to the new obtained information at each iteration. For tomography on stationary states, η_t can be chosen to decrease continuously with the number of iterations while for online state tracking it needs to be constant or adaptive, depending on former measurement outcomes.

We generalised the originally proposed MEG in [14] to make it applicable to qudits (d > 2), and also introduced modifications to make the algorithm more suitable for our experiment.

2.3.1. Measurement schemes

We extended MEG to work with general higher-dimensional systems d > 2 and local Pauli measurements. We use two measurement schemes to perform our state estimation and online tracking. The first scheme is utilising Mutually Unbiased Bases (MUBs) [27], measurements that are optimal for extracting information [28]. For a d-dimensional quantum system, we have a maximum of d + 1 basis sets, each set consistent of d states. However, complete sets of MUBs are not know for every dimension and so it is important to also introduce a second scheme using generalised Pauli operators [29]. These operators are given by the set $\{u_{jk}, v_{jk}, w_{l}\}$:

$$\begin{align} \ u_{jk} & = \vert j\rangle \langle k\vert + \vert k\rangle \langle \,j\vert,\qquad \qquad && 0\leqslant \,j\leqslant k\leqslant d \qquad \end{align} \tag{ 4a }$$

$$\begin{align} \ \ \ \ v_{jk} & = -i\,\vert j\rangle \langle k\vert +i\,\vert k\rangle \langle\, j\vert,\qquad \quad && \ 0\leqslant j\leqslant k\leqslant d \qquad \ \ \ \end{align} \tag{ 4b }$$

$$\begin{align}w_l & = \sqrt{\frac{2}{l+1}}\left(-l\,\vert l\rangle \langle l\vert + \sum_{j = 0}^{l}{\vert j\rangle \langle\, j\vert} \right) && 0\leqslant l \leqslant d. \qquad \end{align} \tag{ 4c }$$

The operators reduce to the Pauli and Gell-Mann matrices for d = 2 and d = 3, respectively. The total number of operators is $d^2 -1$, and for each operator, there are d eigenstates that can be measured. The extension to generalised Pauli operators gives us eigenbases which are orthonormal and informationally complete, thus making the algorithm applicable to any dimension. Projective measurements of MUBs and eigenbases of Pauli operators are conveniently implemented in the experiment, thus making our setup ideal for testing the MEG algorithm.

2.3.2. Modified update rule

The learning rate is usually chosen to be adaptive to guarantee convergence in the case of noisy measurements as shown in [14]. However, upon exploring the effect of choosing different learning rates, we found that a constant learning rate is sufficient to obtain an estimate with the maximum possible accuracy obtainable given our experimental capabilities. Thus, we fix the learning rate to be $\eta_t = \eta = 5$.

The measurement schemes we use provide us with more information at each iteration. In the case of MUB measurements, we randomly select a basis at each iteration. In order to calculate the probabilities, we need to projectively measure all states of that basis for normalisation of the measurement counts. The same holds when we use generalised Pauli measurements, we still need to measure all the eigenstates of a given Pauli operator. Therefore, we modify the loss function to include multiple measurement results obtained at a given iteration as follows:

$$\begin{align} L_t = \sum_{i = 1}^{d}{\left(\mathrm{tr}(\hat{\rho}_{t}X_t^{(i)}) - y_t^{(i)}\right)^2}, \end{align} \tag{ 5 }$$

with gradient:

$$\begin{align} \nabla L_t = 2\sum_{i = 1}^{{d}}\left(\mathrm{tr}(\hat{\rho}_{t}X_t^{(i)}) - y_t^{(i)}\right)X_t^{(i)}, \end{align} \tag{ 6 }$$

where $X_t^{(i)} = \vert \psi_i\rangle \langle \psi_i\vert$ is the projector of the ith state of the MUB (or eigenstate of the Pauli operator) basis that was randomly chosen at iteration t out of d + 1 (or $d^2-1$) possible bases, and $y_t^{(i)}$ is the measured probability. Although this modification is motivated by how we perform our measurements, it is beneficial for many experimental platforms, where often several projective measurements are necessary to calculate probabilities.

Since we are restricting the experiments to pure states in our investigations—true for systems that have very high control over their preparation, we can increase the estimation accuracy using this prior information. The method we use is to project the state into the subspace of its largest eigenvalue, or equivalently, finding the pure state closest to our estimate in fidelity. We illustrate this method with an example: we can write the estimated density matrix $\hat{\rho}_t$ of a pure quantum state $\vert \phi_t\rangle$ in the form:

$$\begin{align} \hat{\rho}_t = \lambda\vert \phi_t\rangle \langle \phi_t\vert + (1-\lambda) \hat{\rho}_{\textrm{noise}}, \end{align} \tag{ 7 }$$

where the parameter $\lambda \in [0,1]$ indicates the noise strength, and $\hat{\rho}_{\textrm{noise}} = I_d/d$ is the completely mixed state representing the noise in the estimate, with I_d the d dimensional identity matrix. The special cases λ = 0 and λ = 1 correspond to the noiseless and completely noisy states respectively. Now, let $\vert \gamma_t\rangle$ be an eigenstate of $\hat{\rho}_t$ with eigenvalue γ_t , i.e. $\hat{\rho}_t \vert \gamma_t\rangle = \gamma_t\vert \gamma_t\rangle$. Then, we can write:

$$\begin{align} \gamma_t & = \langle \gamma_t|\hat{\rho}_t|\gamma_t\rangle \qquad \qquad \qquad \qquad \qquad \quad \end{align} \tag{ 8a }$$

$$\begin{align} & = \lambda\langle \gamma_t|\phi_t\rangle\langle \phi_t|\gamma_t\rangle + (1-\lambda) \frac{\langle \gamma_t|I_d|\gamma_t\rangle}{d} \end{align} \tag{ 8b }$$

$$\begin{align} & = \lambda\left|\langle \gamma_t|\phi_t\rangle\right|^2 + \frac{1-\lambda}{d}.\qquad \qquad \qquad\end{align} \tag{ 8c }$$

From there, we can see that the maximum eigenvalue is:

$$\begin{align} \gamma_{t, \textrm{max}} = \lambda\ + \frac{1-\lambda}{d}, \end{align} \tag{ 9 }$$

which is obtainable when there is maximum overlap between γ_t and φ_t , in the best case $\left|\langle \gamma_t|\phi_t\rangle\right|^2 = 1$, i.e. $\vert \gamma_t\rangle = \vert \phi_t\rangle$. This shows that the true pure state $\vert \phi_t\rangle$ is an eigenvector of the density matrix $\hat{\rho_t}$ corresponding to the maximum eigenvalue, in other words, closest to our estimator in fidelity. So, after each iteration we calculate the eigenvalues of the state estimate and select the corresponding eigenvector to be the estimate of the pure state $\vert \phi_t\rangle$ with density matrix $\hat{\rho}_{\phi}$. Finally, the infidelity $1 - f(\hat{\rho}_{\phi},\hat{\Omega}_{\psi}) = 1 - \left( \mathrm{tr} \sqrt{\sqrt{\hat{\rho}_{\phi}} \, \hat{\Omega}_{\psi} \, \sqrt{\hat{\rho}_{\phi}}} \right)^2$ is used for benchmarking of the algorithm. Here, $f(\hat{\rho}_{\phi},\hat{\Omega}_{\psi})$ is the fidelity and $\hat{\Omega}_{\psi}$ is the density matrix of the theoretical prepared state $\vert \psi_t\rangle$. Our results show convergence to the unknown state within experimental limitations, most prominently mode-dependent losses, in all investigated situations and even for excessive experimental noise. We note that MEG estimates a full density matrix and hence related parameters of interest (as sometimes used in other tomography schemes, e.g. shadow tomography [30, 31]) can be calculated with high accuracy.

We investigated the qutrit (d = 3) state estimation and online tracking performance of MEG using mutually unbiased bases (MUBs) and generalised Pauli operator measurements, which makes the algorithm generic and applicable to higher dimensions. We also investigated the MEG performance in the presence of small and large statistical noise, and also environmental noise. In all cases we continuously evolve the prepared state in time following:

$$\begin{equation} \vert \psi_t\rangle = \exp{\left( -i \sigma \omega t \right)} \vert \psi_0\rangle, \end{equation} \tag{ 10 }$$

where σ is a Hermitian matrix and $\omega = \tfrac{1.3}{t_{\textrm{tot}}}$ the rate of change dependent on the total amount of iterations t_tot. The rate of change was chosen to allow the state to evolve to a minimum fidelity compared to the initial state and then back to its initial state. We studied three possible cases of σ: $0_{3,3}$ which is the 3-dimensional zero matrix, σ_z which is the generalised 3-dimensional Pauli Z matrix given by [29, 32]:

$$\begin{align} \sigma_z = \frac{1}{\sqrt{3}} \begin{pmatrix}1 & 0& 0 \\ 0 & 1 & 0\\ 0 & 0 & -2 \end{pmatrix}, \end{align} \tag{ 11 }$$

and σ_r which is a general random Hermitian matrix. The zero matrix results in the identity evolution leaving the initial state stationary, while the other two matrices simulate a rotation in the Hilbert space of a qutrit. To ensure a fair comparison, the learning rate—a crucial parameter for performance and convergence of the algorithm—is kept constant at η = 5, which shows good performance in all cases albeit not optimal for $\sigma = 0_{3,3}$. The learning rate can be optimised individually, dependent on investigated scenario and physical system.

An initial baseline measurement using MUBs and high signal rates ($N \approx 10^6$ photons per measurement) is shown in figure 1(a). The high count rate ensures that statistical Poissonian distributed counting noise ($\Delta N = \sqrt{N} \approx 10^3$), background ($N_{\mathrm{back}} \approx 50$) and dark counts ($N_{\mathrm{dark}} \approx 100$), are negligible and do not deteriorate the performance. For this case, the algorithm finds density matrices with a median purity above 99.6% in all investigated scenarios, as to be expected for our system where we have a high degree of control over the state preparation. Since the fidelities in the experiments are very high, we plot the infidelity, $1 - f(\hat{\rho}_{\phi},\hat{\Omega}_{\psi})$, for clarity. The solid lines in figure 1(a) are the median infidelities of 50 randomly chosen states of qutrits (according to Haar measure), with the shaded regions bounded by the upper and lower quartile $(50 \pm 25)\%$.

Zoom In Zoom Out Reset image size

In theory, MEG can estimate any quantum state with arbitrary precision and accuracy given an infinite amount of iterations. However, experimental imperfections such as misalignment and mode-dependent loss significantly limit the minimum infidelities achievable by QST algorithms (not just MEG). The infidelities for MEG listed in table 1 are similar to those obtained using other QST algorithms applied to our system, e.g. root-approach and maximum-likelihood estimators [33] where we set the rank to 1.

Table 1. MEG performance indicators. Top: iterations needed to reach $10\%$ infidelity. Bottom: mean infidelity. All values achieved with a constant learning rate η = 5. Uncertainties are the boundaries of the interquartile ranges, i.e. a mean infidelity of $\left(5.3~\substack{+3.2 \\ -2.2}\right)$% means that the middle 50%—between 25%–75%, sorted by their infidelity—of our investigated states showed infidelities between 3.1% and 8.5%.

	MUB			Pauli operators
Signal rate (Hz)	$\sigma = 0_{3,3}$	$\sigma = \sigma_z$	$\sigma = \sigma_{r}$	$\sigma = 0_{3,3}$	$\sigma = \sigma_z$	$\sigma = \sigma_{r}$
	Iterations for infidelity $\lt10\%$
102	$14~\substack{+17 \\ -10}$	$13~\substack{+26 \\ -6}$	$16~\substack{+26 \\ -12}$	$30~\substack{+32 \\ -21}$	$23~\substack{+68 \\ -14}$	$24~\substack{+38 \\ -17}$
106	$4~\substack{+2 \\ -1}$	$4~\substack{+3 \\ -1}$	$4~\substack{+4 \\ -1}$
	Mean infidelity (%)
102	$5.3~\substack{+3.2 \\ -2.2}$	$5.4~\substack{+3.4 \\ -2.4}$	$5.6~\substack{+3.4 \\ -2.5}$	$4.9~\substack{+2.8 \\ -2.0}$	$5.4~\substack{+3.5 \\ -2.4}$	$5.0~\substack{+3.3 \\ -2.1}$
106	$3.4~\substack{+1.8 \\ -1.2}$	$5.1~\substack{+2.2 \\ -1.9}$	$4.7~\substack{+2.3 \\ -1.9}$

We next show the robustness of MEG to noise using MUBs, reducing the count rate per iteration to $N \approx 10^2$. This means that the values for Poissonian, dark, and background noise become comparable and play a significant role in state estimation and tracking. Robustness to statistical noise in the limit of a small number of repeated preparations of the unknown state—here the number of photons—is especially interesting for systems like ions and superconducting qubits, as the preparation is considerably more time consuming compared to our system. Again, we benchmark the performance of MEG by calculating the infidelity in each iteration, see figure 1(b). State estimation and online tracking are in excellent agreement with an overall mean infidelity ${\sim}5.5\%$, slightly higher than in the high count rate case, but well within the uncertainty bars. This already shows the robustness of the modified algorithm to statistical noise. We would like to highlight here that the unitary state evolution (following equation (10)) is done in preparation. The noise is added to the measurement, not the preparation. Infidelities below $10\%$ are achieved within around 15 iterations, four times more than in the baseline measurement with high count rates (see table 1 for details).

We furthermore studied the performance of MEG under high statistical noise, but using generalised Pauli operator measurements rather than MUBs, shown in figure 1(c). As generalised Pauli operators can be mathematically defined, using them instead of MUBs is a powerful tool to describe states far beyond the qubit. We observe that state estimation and online tracking are again in excellent agreement and the overall mean infidelity ${\sim}5.1\%$ is slightly better than in the MUB case, but overlapping within uncertainty bars. This demonstrates that generalised Pauli measurements are a viable path to extend MEG towards arbitrary high dimensions in realistic experimental settings. The initial convergence towards the unknown state is slower with around 25 iterations necessary to reach an infidelity of $10\%$. This is to be expected as MUBs span the qutrit state space more efficiently compared to the Pauli operator measurements. Table 1 shows details of the performance.

Since a changing state leads to changing probabilities of finding a system in a particular state, we also show that we are able to follow these probabilities. Figure 2 illustrates one example of these probabilities p$_i = \left| \langle i|\psi_t\rangle \right|^2$ for the logical basis states ${\vert i\rangle,i = 0,1,2}$, out of the 50 randomly chosen qutrits analysed in figure 1 using MUBs. The theoretically expected probabilities in the experiment (dashed lines) are followed nicely by online state tracking of the deterministic evolution in the experiment (solid lines). This confirms that the low infidelities directly correspond to the capability to predict measurement outcomes on the states well. We also highlight that MEG quickly converges to a useful estimate (after 3-4 iterations), even though the initial guess is a completely mixed qutrit (p$_i = 1/3 \forall i$) to avoid bias in the results. The small offset between theory and experiment could be further decreased by using an adaptive learning rate.

Zoom In Zoom Out Reset image size

We tested the efficacy of MEG under significant environmental noise experimentally—by adding a light bulb close to the detector while keeping the signal rate at $N \approx 10^2$ photons per measurement—and in simulations—adding Poissonian distributed noise to the raw data. The strength of the noise in the experiment is controlled via the distance and angle of the light bulb to the single photon detector. As a proof of principle demonstration, we only looked at stationary states using MUBs in our experiment, yet the simulations show similar behaviour for all estimation and online tracking scenarios investigated in the previous section. The median infidelities of 50 randomly chosen states for different methods (solid—experiment, dashed—simulation) and environmental noise strengths (colours) are shown in figure 3. Shaded uncertainty ranges similar to figure 1 are omitted for clarity. As expected, the achievable infidelity and number of iterations to get there both become higher as the environmental noise increases, however, MEG is showing exceptional robustness.

Zoom In Zoom Out Reset image size

We observe that up to $N_{back} \approx 2$ kHz of experimental and simulated noise (not shown explicitly in the figure), the infidelities are comparable to the no noise case within uncertainty bars. This corresponds to a SNR ≈ 0.05, so only one in 20 photons carry the actual information on the unknown quantum state prepared in the experiment. Here and in figure 3, N_back relates purely to environmental noise from the light bulb, excluding other noise sources like e.g. dark counts. It is however noteworthy that even at 0 kHz of artificially introduced noise, the naturally occurring counting noise together with background and dark count rates are almost double the signal count rate, i.e. ${\sim}170$ compared to ${\sim}100$ photons per measurement. Experimental and simulated data for $N_{back} = \left[1,2.5\right]$ kHz show good overlap within uncertainty bars, indicating the validity of the noise model in simulation. Therefore we expect to achieve infidelities of around $10\%$ up to $\textrm{SNR} \approx 0.01$ (orange dashed line). The seemingly smoother performance of the simulations stems from 20 runs of randomly added shot-to-shot noise for each of the 50 investigated states. In other words, the results from simulations reflect the median infidelity of $50\times 20 = 1000$ traces, compared to 50 in the experiment.