# In shadow tomography, how is the state reconstructed from its shadows?

I'm reading Huang et al. (2020) (nature physics), where the authors present a version of Aaronson's shadow tomography scheme as follows (see page 11 in the arXiv version):

We want to estimate a state $$\rho$$. We apply a number of random unitary evolutions, $$\rho\mapsto U\rho U^\dagger$$, picking $$U$$ from an ensemble $$\mathcal U$$. For each choice of unitary $$U$$, we perform a measurement, observing a state $$|b\rangle$$. We then apply the inverse evolution to this state, obtaining $$U^\dagger|b\rangle$$. On average, this procedure leaves us with the state $$\mathcal M(\rho) =\mathbb E_{U\sim\mathcal U}\Bigg[ \sum_b \underbrace{\langle b|U\rho U^\dagger|b\rangle}_{\text{prob of observing |b\rangle}} \!\!\!\!\!\!\! \overbrace{(U^\dagger |b\rangle\!\langle b|U)}^{\text{post-measurement state}} \Bigg].$$ The claim is then that $$\mathcal M$$ can be inverted to obtain, on average, the original state: $$\hat\rho=\mathcal M^{-1}(U^\dagger |\hat b\rangle\!\langle \hat b| U) \,\,\text{ is such that } \,\, \mathbb E[\hat \rho]=\rho.$$

Is there an easy way to see how one would go in performing such inversion? The authors mention that we are thinking here of $$\mathcal M$$ as a linear map, so I suppose we represent $$\rho$$ as a vector in some operatorial basis of Hermitians, which is fine, but to then perform the inversion of this linear map we would need to characterise it, which in this case I think would mean to know the values of $$\mathcal M(\sigma_i)$$ for some complete basis of Hermitians $$\{\sigma_i\}_i$$.

• In practice, I think what you described is the best way. Just compute $\mathcal M(|i\rangle\!\langle j|)$ and find the inverse by standard linear algebra methods. Dec 2, 2020 at 10:19

Under the assumption that the ensemble $$\mathcal{U}$$ faithfully produces the Haar expectations at least to the second moment, the inversion can be performed as suggested in the last paragraph of the question:

Define:

$$\theta_b = U^{\dagger}|b\rangle\langle b| U$$

and replace the averaging over the ensemble by Haar averaging. (This step is is done only to compute the resulting quantity. In practice we average over the ensemble):

$$\mathcal{M}(\rho) = \int d\mu_{\small{HARR}} \sum_b \mathrm{Tr}(\rho \theta_b) \theta_b$$

($$d\mu_{\small{HARR}}$$ is normalized to a unit volume). Suppose that we need to measure the expectation of an operator $$O$$ in the state $$\rho$$. We measure this operator in the final state average the expectations over the ensemble and the computational basis projectors, we obtain the following result:

$$\mathrm{Tr}(\mathcal{M}(\rho) O) = \int d\mu_{\small{HARR}} \sum_b \mathrm{Tr}(\rho \theta_b) \mathrm{Tr}(O \theta_b)$$

Using the following identity:

$$\int d\mu_{\small{HARR}} \mathrm{Tr}(\rho \theta_b) (O \theta_b) = \frac{ \mathrm{Tr}(\rho) \mathrm{Tr}(O)+ \mathrm{Tr}(\rho O) }{N(N+1)}$$

($$N$$ is the dimension of the Hilbert space), thus we obtain:

$$\mathrm{Tr}(\mathcal{M}(\rho) O) = \frac{ \mathrm{Tr}(O)+ \mathrm{Tr}(\rho O) }{(N+1)}$$

From which the required expectation $$\mathrm{Tr}(\rho O)$$ can be computed.

Of course, we can repeat the procedure for a complete basis of the Lie algebra $$\mathfrak{u}(N)$$ for a full tomography of $$\rho$$.

Explanation of the integration formula

$$\theta_b = U^{\dagger}|b\rangle\langle b| U$$ is a one-dimensional projector can be written as:

$$\theta_b = |\Psi \rangle\langle\Psi|$$

where $$|\Psi\rangle$$ is a unit vector

When $$U$$ uniformly scans the unitary group manifold, $$|\Psi\rangle$$ will uniformly scan a $$2N-1$$ dimensional round sphere. This is because the metric of a round sphere is Haar. Thus, we need to evaluate an integral of the form:

$$\int_{S^{2N-1} }\langle\Psi| A |\Psi \rangle \langle\Psi| B |\Psi\rangle d\mu_{\small{ROUND}}(S^{2N-1})$$

The integrand is a homogeneous polynomial of the components of the unit vector $$|\Psi\rangle = \sum_{a=1}^N \psi_a |a\rangle$$:

$$\langle\Psi| A |\Psi \rangle \langle\Psi| B |\Psi\rangle = \sum_{a,b,c,d} A_{ab} B_{cd} \bar{\psi_a}\bar{\psi_c}\psi_b \psi_d$$

Integrals of homogeneous polynomials on the sphere can be evaluated as Gaussian integrals up to an overall normalization, which can be computed from a trivial particular example. In this case we can use the Wick's theorem: The constant can be evaluated by taking $$A=B=I_{N\times N}$$, we get:

$$\mathrm{Const.} = \frac{1}{N(N+1)}$$

• great! thanks. Can you add a source, or sketch the derivation, for that identity?
– glS
Dec 1, 2020 at 11:40
• Another way to get this is by symmetry. What we actually compute is the Haar average $P_H(\psi)=\int \mathrm{d}\psi |\psi\rangle\langle \psi|^{\otimes 2}$. This commutes with permutations, thus it has to be proportional to the projector onto the symmetric subspace by Schur's lemma. Thus $\mathrm{tr}(A\otimes B P_H(\psi)) = C \mathrm{tr}(A\otimes B P_{\mathrm{Sym}_2})$, now apply to swap trick to get the result. $C$ can be computed by setting $A$ and $B$ to the identity. (BTW the requirement on the ensemble is that it should be a unitary 2-design) Dec 2, 2020 at 9:32

This is not a complete answer but a couple of comments that should help clarify these details:

(i) $$\mathcal{M}$$ is not invertible as a quantum channel but as a linear map; namely, even though the inverse exists, $$\mathcal{M}^{-1}$$ is not CP, even though it is a linear map. Assuming the input-output dimensions are the same, the only CP maps that are invertible are unitaries. This is the same problem that arises when studying dynamical semigroups for Markovian evolutions, they are described as a family $$\{ e^{t \mathcal{L}} \}_{t \geq 0}$$ where $$\mathcal{L}$$ is the Lindbladian; clearly $$e^{-t \mathcal{L}}$$ is the inverse, but it is not CP.

(ii) The authors assume tomographic completeness,'' that is, for each $$\sigma \neq \rho, \exists U \in \mathcal{U} \text{ s.t. } \langle b| U \sigma U^{\dagger} | b \rangle \neq \langle b| U \rho U^{\dagger} | b \rangle$$ which ensures that $$\mathcal{M}^{-1}$$ exists and is unique.

(iii) The classical shadow, $$\hat{\rho}$$ is not necessarily positive semidefinite since it originates from the application of a not CP map on some state.

(iv) $$\mathcal{M}^{-1}$$ is applied to the classically stored measurement outcomes, $$U^{\dagger} | \hat{b} \rangle \langle \hat{b} | U$$ and not the quantum states themselves.

• thanks, this is useful info. For me the main question remains how one does actually calculate the inverse $\mathcal M^{-1}$, be it a physical channel or not. Which I suppose amounts to asking how to characterise $\mathcal M(\sigma_i)$ on a complete basis of Hermitians $\{\sigma_i\}_i$. I guess this must be done beforehand somehow, and that it only needs to be done once to be then able to use the result for any state in need for characterisation. Regarding the non-positivity of $\hat\rho$, I understood that such non-positivity was to go away when we average over many unitaries/outcomes
– glS
Dec 1, 2020 at 9:11
• @glS The easiest way I know to compute the inverse is to use the Choi-Jamiolkowski isomorphism, which sends $\mathcal{M} \mapsto \rho_{\mathcal{M}}$, then invert the matrix $\rho_{\mathcal{M}}$ and then do the inverse CJ transformation to get $\mathcal{M}^{-1}$ back. Dec 1, 2020 at 10:32
• don't you still need to know the matrix representing $\mathcal M$ (thought of as a linear operator) to do that? Computing the Choi essentially amounts to computing $\mathcal M(|i\rangle\!\langle j|)$
– glS
Dec 1, 2020 at 10:45