Complexity of quantum tomography from genuine non-Gaussian entanglement

Domain: physics
Published: 2025-12-05
Editorial outcome: minor_revision

Abstract
Quantum state tomography typically requires exponentially many copies of a quantum state, due to the complex correlations present in large systems. We show that, for bosonic systems, the scaling is completely determined by the nature of these correlations. Motivated by the Hong-Ou-Mandel effect and boson sampling, we define Gaussian-entanglable (GE) states, produced by generalized interference between separable bosonic modes. GE states greatly extend the Gaussian family, encompassing separable states, multi-mode Gottesman-Kitaev-Preskill codes, entangled cat states, and boson-sampling outputs—resources for error correction and quantum advantage. We prove that any pure GE state of m modes can be learned efficiently, requiring only poly(m) copies, via a protocol based on Gaussian unitaries, local tomography, and classical post-processing; for boson-sampling states, no Gaussian unitaries are needed. For states outside GE, we define an operational monotone—the minimal number of ancillary modes needed to make them GE—which exactly characterizes the exponential tomography overhead. We also show that deterministic generation of NOON states with N≥3 via two-mode interference is impossible.

Introduction
Quantum state tomography1–3 aims to provide a classical description of quantum states from experimental data. As a fundamental task in physics and quantum information processing, it not only facilitates the tests of quantum physics, but also enables quantum device certification and benchmarking4. However, quantum state tomography is a challenging task in general—each measurement requires a new copy of the state as measurement collapses the original quantum state. To learn an approximate description of a quantum system with n qubits, the required number of samples is exponential in n5,6. In the case of bosonic systems, the infinite-dimensional Hilbert space casts additional challenges to the problem of tomography7.

The challenge in quantum state tomography arises from the complex nature of quantum correlations between subsystems: highly entangled states involving superpositions of a large number of states are generally hard to learn, while separable states with classical correlations are naturally easier to learn. At the same time, one may identify subsets of non-separable quantum states that are easier to learn. For instance, refs. 7–9 developed efficient learning protocols for Gaussian bosonic quantum states10,11—a family of experimentally friendly states—and its generalization of t-doped Gaussian states. Thus, the interplay between properties of quantum correlation and learning sample complexity in bosonic quantum systems is unclear. In this work, we approach the problem by considering how quantum correlations are generated from separable states.

We begin by considering the example of the Hong–Ou–Mandel (HOM) effect12, where two incoming identical and separable single photons interfere on a beam-splitter to deterministically generate an entangled output, manifesting a bunching effect between indistinguishable photons. The output state |20⟩+|02⟩ is the N=2 case of the family of NOON states |N0⟩+|0N⟩, which are maximally entangled N-photon states useful to enhance quantum metrology13,14. Similar to HOM, in boson sampling experiments15, multiple modes of photons interfere on a network of beam-splitters creating highly multipartite entangled states, whose measurement statistics are hard to simulate classically, creating an opportunity to demonstrate quantum advantage16,17.

Inspired by the HOM effect12 and boson sampling experiments15, we define Gaussian-entanglable (GE) states as states produced by generalized interference described by Gaussian processes10 between separable bosonic modes (see Fig. 1). Despite the fact that GE states enable quantum computational advantage, we prove that an m-mode pure GE state is learnable with only poly(m) samples, by providing an explicit protocol—Gaussian disentangling (GDE) tomography—involving only Gaussian unitary, local tomography (e.g., homodyne or heterodyne) and classical processing. The protocol is robust to experimental imperfections, including weak loss and noise, and arbitrary homogeneous loss. We provide numerical simulations for learning the states generated from passing number states, cat states and thermal states over a lossy beam-splitter, demonstrating GDE’s advantage over conventional approaches.

For non-GE states, we derive an operation-driven monotone—the number of additional modes necessary to reduce them to GE states. Moreover, we prove that this monotone exactly captures the exponential overhead in sample complexity of learning the states. Besides their clear implications for the sample complexity of learning, our results also settle a long-standing open problem, namely that NOON states with N≥3 photons cannot be generated deterministically by interfering an arbitrary product state, as they are not in the GE class. Earlier studies addressed only restricted cases involving number state inputs in linear optical settings18–20.

Results
To characterize quantum correlations beyond HOM and boson sampling in bosonic systems, we introduce non-Gaussian-entanglable (NGE) states: these are the quantum states which cannot be generated via a generalized interference experiment (see Fig. 1), where separable quantum states are input to Gaussian protocols10,11. Conversely, quantum states generated via generalized interference experiments are GE states. GE states are operationally meaningful, especially considering the recent progress of multi-mode quantum state engineering21,22. By definition, GE states include all Gaussian states and their convex mixtures23–28. In addition, they include states that are neither separable nor convex mixtures of Gaussian states, such as the output states of boson sampling experiments15, which involve identical photons scattering through a linear optical interferometer and whose measurement statistics is computationally hard to classically simulate. Moreover, the recently proposed multi-mode Gottesman–Kitaev–Preskill (GKP) states, with superior error correction performance over single-mode versions29–32, are GE—they can be engineered via applying Gaussian unitaries on single-mode GKP states31,32, reducing the challenge of multi-mode state engineering to preparing single-mode GKP states. On the contrary, NGE states will generally require a multi-mode non-Gaussian operation—or equivalently realized by Gaussian operations assisted with a sufficient amount of non-Gaussian ancilla33, creating a large overhead in quantum information processing.

To sharpen the distinction between GE and NGE states, we develop a convex resource theory34 of genuine non-Gaussian entanglement. With GE states as the free states, the free operations are Gaussian protocols. Besides proving that NOON states with N≥3 are NGE, we provide other examples of NGE states, such as superposition of two-mode squeezed vacuum states (TMSV). Furthermore, we introduce and analyze two monotones to quantify genuine non-Gaussian entanglement, one based on entanglement entropy (the NG entanglement entropy) and the other one based on the minimal number of ancillary modes necessary to convert an NGE state into a GE state (the GE cost). Moreover, we propose the GDE tomography protocol based on Gaussian unitary, local tomography (e.g., homodyne or heterodyne) and classical processing, and show its efficient learning in pure GE states and special classes of mixed GE states. We also establish a continuity proof for the robustness of our protocol for general states close to pure GE states. We show that the sample complexity of tomography for pure GE states grows polynomially with the number of modes, while NGE states generally require an exponential overhead, quantified by their GE cost. This implies that output states in boson sampling experiments are easy to learn, despite their measurement statistics being hard to sample from with a classical computer15. We also provide numerical verification of the tomography process. The main results of this work are summarized in Table 1.

Continuous-variable systems and state learning
Continuous-variable quantum systems model continuous degrees of freedom, such as optical fields35, mechanical motions of ions36, and microwave cavities37. Owing to their information-rich infinite-dimensional Hilbert space, these systems are crucial for quantum communication38, quantum sensing39,40, quantum error correction22,33, and quantum computing41 (see Supplemental Note 1A).

While general states of continuous-variable quantum systems are unbounded in energy, we are interested in practically relevant states with energy constraints. In particular, we consider m-mode quantum states with finite energy moment constraint, [Tr(E^k_m ρ)]^(1/k) ≤ mE, where the energy operator E^m = 1/4 Σ^m_ℓ=1 (q^2_ℓ + p^2_ℓ) is a quadratic function of the position and momentum quadratures q^ℓ and p^ℓ in natural units ℏ=2. We denote the set of m-mode pure states satisfying this k-th energy moment constraint as H^k_{m,E}. We denote the set of mixed states satisfying the same energy constraint as S^k_{m,E}(H^2).

Learning a state within a trace distance error of ϵ is equivalent to pinpointing it in state space within an ϵ-ball with respect to the trace distance metric. The number of such ϵ-balls—the volume—of the state space that the unknown state belong to can manifest the challenge of state learning. It has been shown in ref. 7 that, for pure states subject to k-th moment constraint, the volume of the state space grows double-exponentially with the system size, scaling as exp(m^k/ϵ^2).

Gaussian state and protocols
Among bosonic states and operations, Gaussian states (e.g., coherent states and squeezed states), and operations (e.g., beam-splitters, phase rotations, and squeezing) are not only analytically tractable but also experimentally friendly10,11.

Gaussian states are specific types of quantum states with their Wigner characteristic function in a Gaussian form42,43; therefore, a Gaussian state ρ^ is uniquely specified by its mean ξ := Tr[r^ρ^] and covariance matrix V := 1/2 Tr[{r^−ξ, r^−ξ}^T ρ^]10,11, where r^=[q^1, p^1, ..., q^m, p^m]^T is the vector of quadrature operators. For instance, the TMSV, defined as |ζ_{r,ϕ}⟩_{AB} = (cosh r)^−1 Σ^∞_{k=0} e^{ikϕ}(tanh r)^k |k⟩_A |k⟩_B, have zero-mean and are uniquely determined by a covariance matrix that depends on the squeezing strength r and phase angle ϕ. Here |k⟩ is the Fock state with photon number k. Gaussian quantum channels are completely positive trace-preserving maps that transform all Gaussian states into Gaussian states (see Supplemental Note 1B). In the following text, we shall use the notation G(H) and GU(H) to represent the set of all Gaussian quantum channels and unitaries on the Hilbert space H. More generally, a Gaussian protocol27 includes the following set of quantum operations: (1) Gaussian unitaries on the system or part of the system, (2) partial trace, (3) composition with ancillary vacuum states, (4) probabilistically applying the above operations. Note that we exclude Gaussian measurements followed by the conditional operations (1)–(4) from the set of Gaussian protocols we consider, because these are more challenging from an experimental standpoint and the resulting feed-forward scheme might not reduce to a convex mixture of Gaussian channels.

Gaussian protocols (1)–(4) alone, without non-Gaussian ancilla, are limited in power—they are not universal for various tasks44–51. As it turns out, even when allowing for non-Gaussian input states, the generalized Gaussian interference process in Fig. 1 can only generate a limited class of states, as we define below.

Gaussian-entanglable state
Consider a bi-partition of the Hilbert space H=H_A ⊗ H_B. A quantum state ρ^_{AB} is (bi-partite) GE if it can be obtained by probabilistically applying poly-many Gaussian channels G_x on a separable state σ^_{AB,y} of H_A ⊗ H_B, ρ^_{AB} = Σ_{x,y} p_{x,y} G_x(σ^_{AB,y}) with p_{x,y} being a joint probability distribution. Equivalently, a GE state can be expressed as:
ρ^_{AB} = Σ_{x,y} p_{x,y} G_x(σ^_{AB,y}), (1)
where p refers to a joint probability distribution and G_x ∈ G(H_A ⊗ H_B) denotes a Gaussian channel.

In practice, the GE state in Eq. (1) can be obtained by sampling from the joint probability distribution p_{x,y}, preparing the separable input state σ^_{AB,y}, and implementing the channel G_x. When G_x is unitary, it corresponds to randomly applying Gaussian unitaries on separable states. By definition, the set of GE states is convex—a mixture of GE states is still GE.

As schematically illustrated in Fig. 1, GE states can be non-Gaussian, and even not a convex mixture of Gaussian states. For instance, all entangled states generated by beam-splitter interference of non-Gaussian states, such as number states, are GE, e.g., the state (|0⟩_A |2⟩_B + |2⟩_A |0⟩_B)/√2 from the HOM effect12. Similarly, the entangled cat state 1/√2(|α⟩_A |α⟩_B + |−α⟩_A |−α⟩_B) is GE, as it can be generated using a balanced beam-splitter on a product of cat state and vacuum, 1/√2|cat⟩_{α,A} ⊗ |0⟩_B. Here |α⟩ is a coherent state with amplitude α, and |cat⟩_α ∝ |α⟩+|−α⟩ denotes a cat state. (Here, we have introduced the underline to distinguish between coherent states and number states, omitted in the special case of the vacuum state |0⟩.) Moreover, one can obtain a state 1/√2(|cat⟩_{α,A} |0⟩_B + |0⟩_A |cat⟩_{α,B}) via interference between two identical cat states |cat⟩_α. Other examples include GE states generated by the multi-mode GKP encoding process: |ψ⟩ = U^g ( |ψ⟩_m ⊗ |GKP⟩_m0), where U^g is a multi-mode Gaussian unitary, |ψ⟩_m is an arbitrary m-mode state, |GKP⟩ is the canonical GKP state. Any multi-mode GKP state is also GE29–32.

As we show, one can obtain any pure GE states from Gaussian unitary on a product state, as shown below.

Theorem 1. (Decomposition of pure GE states): Any pure GE state can be expressed in the following form:
|ψ⟩_{AB} = U^g |ϕ⟩_A |ϕ'⟩_B, (2)
where U^g is a Gaussian unitary, |ϕ⟩_A and |ϕ'⟩_B are local states of systems A and B.

The proof, as elaborated in Supplemental Note 2, is based on the concavity of entropy, Gaussian unitary decompositions, and more importantly the fact that the only states that remain pure under multi-mode pure loss are coherent states, which generalizes a single-mode result that is crucial for solving the minimum entropy output conjecture and the Gaussian state majorization conjecture52,53 (see Lemma 7 in “Methods”).

GE states possess several distinctive properties (see details in Supplemental Note 3). For a pure Gaussian state with a non-degenerate symplectic spectrum, the condition that its covariance matrix is block-diagonal (linear uncorrelated) is equivalent to the state being a product state (independent). This equivalence fails when the symplectic eigenvalues are degenerate. In addition, Gaussian protocols map GE states to GE states. At last, one can extend the definition of bi-partite GE states in Eq. (1) and Theorem 1 to a multipartite m-mode case.

Furthermore, one can denote the covariance matrix of |ψ⟩_{AB} as V_{AB}. For any symplectic matrix S such that S V_{AB} S^T = Λ is the diagonal form, the corresponding Gaussian unitary U^g is equivalent to U^g = D^ξ U^S U^O, where D^ξ is a displacement, U^S is a squeezing, and U^O is a passive unitary. In Eq. (2), U^g = D^ξ U^S U^O. If U^O is a tensor product of phase shifters, then by first measuring the mean vector and covariance matrix, then applying the partial inverse U^y_S D^y_ξ of the Gaussian unitary, any GE state can be transformed into a passive-separable state (GE states with only passive beam-splitters as the entangling operation). The next section outlines a tomography protocol for GE states that exploits this counter-rotation step.

Theorem 2. (Sample complexity of pure GE states): Consider the set of m-mode pure states {|ψ⟩} obtained by applying a Gaussian unitary to a tensor product of single-mode states ⊗^m_{j=1} |ψ_j⟩, where the states satisfy energy moment constraint, |ψ⟩ ∈ H^2_{m,E}. Then the tomography of such states requires M copies such that:
Θ(poly(m, E, 1/ϵ, log(1/δ))) ≤ M ≤ O(poly(m, E, 1/ϵ, log(1/δ))), (3)
where δ denotes the failure probability and ϵ refers to the trace distance between the input and the reconstructed state.

The upper bound in Theorem 2 is demonstrated in Supplemental Note 7D, while the lower bound is determined by sample complexity of product states with moment constraints7. Therefore, the tomography of pure GE states is sample-efficient. For the case with a moderate number of modes in which constant prefactors in the sample complexity become important, we provide an explicit upper bound in Supplemental Note 7F to offer more practical guidance for GE-state learning.

Below, we present the GDE learning protocol (see Fig. 2), involving a Gaussian unitary and then local state tomography. The local state tomography can be performed via different approaches, for example using heterodyne detection and classical processing. We begin with the simple example of passive separable states, where the physical Gaussian unitary is not needed. Then, we proceed to general pure GE states. In the next section, we generalize towards mixed GE states.

As a concrete example, consider the output state of an m-mode boson sampling process, known as passive-separable states { |ψ_{ps}⟩ }54, generated by applying a passive Gaussian unitary U^O to a product state { ⊗^n_{j=1} |ψ_j⟩ }, where |ψ_j⟩ is a single-mode state (for multi-mode subsystem, we refer to Supplemental Note 7A). We present a tomography algorithm that takes as input M copies of the unknown state |ψ_{ps}⟩ = U^O ⊗^n_{j=1} |ψ_j⟩, output the classical description of a boson-sampling output state |ψ_{ps}⟩ with smaller than ϵ^2 fidelity error with the input state |ψ_{ps}⟩ (leading to a trace distance error smaller than ϵ by the Fuchs–van de Graaf inequality55) and a failure probability lower than δ. As shown in Fig. 3, without relying on any adaptive strategy, the algorithm first performs heterodyne measurement on all M copies of the quantum states and then performs classical post-processing to obtain an estimate of the state.

To describe the post-processing algorithm, we denote the heterodyne measurement result q^{(ℓ)} + i p^{(ℓ)} for the k-th mode in the ℓ-th copy. We collect all results as vectors γ^{(ℓ)} = (q^{(ℓ)}_1, p^{(ℓ)}_1, ..., q^{(ℓ)}_m, p^{(ℓ)}_m)^T, 1≤ℓ≤M. Given the measurement results γ(ℓ), the algorithm first performs a lower bound estimate W. The exhaustive search ends when the fidelity lower bound estimate satisfies W≥1–ϵ^2, outputting the description of the corresponding estimate |ψ_{ps}⟩ = U^O0 ⊗^n_{j=1} |ψ_j^0⟩. Otherwise, a failure is declared.

As sketched in Fig. 3, to obtain the fidelity lower bound estimate, we first calculate α^{(ℓ)} = O^0 γ^{(ℓ)} and reconstruct all local states {ρ^_{j,O0}}^m_{j=1}, whose matrix element in the number basis is
(ρ^_{j,O0})_{kk'} = M^−1 Σ^M_{ℓ=1} ∫_{C^2} d^2u/π e^{|u|^2/2 + i h α^{(ℓ)}, u i} Ω_u χ_{kk'}(u), (4)
where Ω denotes the symplectic form (see “Methods”) and χ_{kk'}(u) refer to the characteristic function of the operator |k⟩⟨k'|. The above reconstruction is a standard heterodyne tomography56, but other reconstruction algorithms can also be adopted57,58. Finally, following ref. 59, we can compute the fidelity witness estimate
W = 1 − Σ^n_{j=1} (1 − F_j), (5)
where each F_j = ⟨ψ_j^0| ρ^_{j,O0} |ψ_j^0⟩ is the local state fidelity for the current estimate.

A challenge in the classical processing is that it need to go over all possible guesses of U^O0 and local states, which are parameterized by continuous parameters with infinite possibilities. The way out of this dilemma is that we only require a finite error in the state learning: the finite-error relaxation allows one to establish a discrete set of ϵ-balls of the possible U^O0 and local states. Essentially, learning the parameters up to a discrete set of values suffices for meeting the ϵ error requirement. Accordingly, for each of the ϵ-balls, we take the mean values as the choice and perform a fidelity test, which also has a failure probability of not satisfying the needed error tolerance. The failure probability then builds up with the number of ϵ-balls—generally a huge number. On the other hand, the failure probability decreases exponentially with the number of samples M; therefore, the required number of samples is logarithmic of the number of ϵ-balls. As a result, the key part of the proof of Theorem 2 is to upper bound the number of ϵ-balls (the volume) required for boson-sampling states, as we analyze in the “Methods” section in detail.

With the sample-efficient learning of passive-separable states introduced, the generalization to arbitrary GE states can be obtained by applying a Gaussian unitary to reduce a GE state to a passive-separable state, hence, we name the protocol Gaussian-disentangling (GDE) learning. To do so, we first estimate the mean and covariance matrix of the input state to estimate a symplectic matrix S which diagonalizes the covariance matrix. Here, one can use either homodyne or heterodyne for covariance matrix estimation. Based on S, we apply a Gaussian unitary to reduce the input GE state to a zero-mean passive-separable state, which then can be sample-efficiently learned via the algorithm sketched in Fig. 3. The details can be found in Supplemental Note 7.

Toward efficient learning of mixed GE states
In practice, preparing ideal pure GE states may be infeasible. Noise, loss, and decoherence can lead to mixed states, in which case the counter-rotation step may fail to map it to a passive separable form. Nevertheless, we find that some specific mixed GE states remain learnable via GDE tomography, with sample complexity polynomial in the number of modes.

First, we demonstrate that our GDE learning algorithm is robust against small perturbations caused by general (potentially non-Gaussian) noise. Even if the unknown state is a slightly perturbed GE state—technically, a state constructed by sending an input state, which is ϵ_{channel}-close to a Gaussian unitary state, into a channel—GDE learning algorithm remains effective and tomography remains efficient. In “Methods”, we show that the corresponding estimation errors for moments have fixed floors of order O(√ϵ_{channel} + ϵ_{channel}) and O(ϵ_{channel}) + O(ϵ_{channel}^{1/8}) for the approximated passive-separable state after counter-rotation. If we assume that the input error ϵ_{channel} vanishes, we have the following theorem:

Theorem 3. Consider an approximated pure GE state ρ^ acting on H^2_{m,E} whose entangling channel is ϵ_{channel}-close to an arbitrary Gaussian unitary. The sample complexity of learning this class of state is
M = O(poly(m, E, 1/ϵ, log(1/δ))), (6)
for an error ϵ > ϵ_{channel} in trace distance and failure probability δ.

The inequality ϵ > ϵ_{channel} arises because the measured state is known to be ϵ_{channel}-close to the target, which induces an irreducible error floor for the estimators. A detailed discussion appears in Supplemental Note 8.

Second, we extend our discussion to GE states generated by applying a Gaussian unitary to a product of arbitrarily mixed states. In particular, this class includes the output of boson sampling under uniform photon loss and situations with detectors of equal efficiency. In the case where the symplectic eigenvalues of the input state has a finite eigen-gap Δ, only moment estimation, Gaussian unitary, and local tomography will be used. Then, the following corollary holds:

Corollary 4. (Sample complexity of states entangled by a Gaussian unitary): Consider the set of m-mode states {ρ^} obtained by applying a Gaussian unitary to a tensor product of single-mode states ⊗^m_{j=1} ρ^_j, subject to an energy-moment constraint ρ^ ∈ S^2_{m,E}. Meanwhile, the states’ symplectic eigenvalues are non-degenerate, up to a finite eigen-gap Δ. The GDE protocol needs M = O(poly(m, E, Δ^−1, ϵ^−1, log(δ^−1))) copies for this state family, for an error ϵ in trace distance and failure probability δ.

In Supplemental Note 8, we show that GDE protocol is able to achieve the sample complexity in Corollary 4, with local state tomography switched from fidelity witness in the pure-state case to heterodyne local tomography. The case of degenerate symplectic eigenvalues is not covered by the above Corollary. In Method, we show that there exist protocols relying on general measurements that in-principle guarantees sample-efficient learning of such degenerate states. Additionally, we provide a generalized GDE protocol with the hope of efficient learning such states with practical Gaussian operations, however, subject to an unproven conjecture.

Below, we apply Corollary 4 and Theorem 3 to practical experimental scenarios. In this regard, detector inefficiency is a ubiquitous and practically relevant imperfection that leads to mixed NGE states. Here, it can be modeled as ideal detection preceded by pure-loss channels. As uniform loss commutes with beam-splitters, it maps boson-sampling outputs to a passive Gaussian unitary acting on a product of mixed states. Consequently, Corollary 4 applies in the situation with non-degenerate symplectic eigenvalues, where both the sample and time complexity stay polynomial in the number of modes.

For non-uniform loss, the degraded channel remains close to a Gaussian unitary in the low-loss regime. Thus, Lemma 11 in “Methods” guarantees a vanishing reconstruction error floor, and Theorem 3 shows that the sample complexity increases only polynomially with the degradation level.

To evaluate the performance of learning GE states, we numerically simulate the GDE tomography process for two-mode GE states with non-degenerate symplectic eigenvalues. As shown in Fig. 4, the GDE method, optimized over the copies allocated to estimating the mean vector and the covariance matrix (red and pink), achieves lower errors and markedly better scaling than direct tomography (black and gray)2. Details of the simulation are presented in Supplemental Note 7B.

Non-Gaussian-entanglable states
While we have shown that GE states are efficient learnable, we know that general quantum states are hard to learn7. This immediately indicates the existence of states that do not belong to the GE class defined in Eq. (1), which we refer to as NGE states. By definition, these states feature genuine non-Gaussian entanglement, which cannot be produced by Gaussian protocols. In practice, we show that NGE states can be produced by a multi-mode non-Gaussian unitary, or by a local non-Gaussian unitary acting on a non-trivial GE state. Furthermore, entanglement-breaking channels60, followed by Gaussian protocols, can destroy genuine non-Gaussian entanglement. Moreover, it can be verified that Gaussian unitary maps NGE states to NGE states.

From Theorem 1, we have an efficient protocol to identify pure NGE states (see “Methods”). As a result, we provide a few examples of NGE states (see Fig. 1), whose membership proofs can be found in Supplemental Note 4. Firstly, we prove that the NOON states |NOON⟩ ∝ |N⟩_A |0⟩_B + |0⟩_A |N⟩_B are NGE for N>2. This means that, aside from the trivial N=1 case, the N=2 case from the HOM effect is the only nontrivial instance of a beam-splitter generating NOON state—those with N>2 cannot be generated deterministically by any Gaussian protocol acting on two-mode separable states, including linear optical operations. Note that, for this specific problem, previous approaches have only established such a no-go result for number state inputs in a linear optical setting18–20, where the photon number is preserved and therefore it suffices to disprove the case with two-mode Fock state. Our approach completes the general case by reducing general Gaussian unitaries to passive unitaries.

In addition, we show that the superposition of a TMSV state and a π-phase rotated TMSV state, |Ψ_{sTMSV}⟩ ∝ |ζ_{r,0}⟩ + |ζ_{r,π}⟩ with r>0, is also NGE. Such states are along the line of recent single-mode experiments generating non-classical superposition of Gaussian states in quantum systems with universal control61.

Resource measures of genuine non-Gaussian entanglement
To quantify genuine non-Gaussian entanglement, we further employ the quantum resource theory framework34, which has been successful in quantifying quantum entanglement, quantum coherence and quantum non-Gaussianity. Here, we model all GE states as free states. To be general and also to provide insight into multi-mode state learning, we consider the case of multi-partite GE states, with n partitions A={A_1, ⋯, A_n}, and define the Hilbert space H_A accordingly. On this account, Gaussian protocols are free operations as they map free states to free states (see Supplemental Note 3C). However, due to the multipartite nature of entanglement, the general quantification of genuine non-Gaussian entanglement faces similar challenges to the case of multipartite entanglement62.

To tackle this issue, we define the non-Gaussian entanglement entropy, abbreviated simply as the NG entropy (see Supplemental Note 5). Given an arbitrary n-partite pure state |ψ⟩_A ∈ H_A, denote ρ^_{A}(U^g) = U^g |ψ⟩⟨ψ|_A U^{g†} and its reduced state as ρ^_{A_j}(U^g) = Tr_{k≠j} [ρ^_{A}(U^g)]; the NG entropy is defined as:
E_{NG}(|ψ⟩⟨ψ|_A) := min_{U^g ∈ GU(H_A)} 1/n Σ^n_{j=1} H(ρ^_{A_j}(U^g)), (7)
with H(ρ^) = − Tr[ρ^ log ρ^] being the von Neumann entropy of the state ρ^. Here, we minimize over all Gaussian unitaries. Given an arbitrary pure GE state |ψ⟩, the resource measure E_{NG} satisfies the following criteria: (1) E_{NG}(|ψ⟩⟨ψ|) = 0 if and only if ψ is GE; (2) E_{NG}(|ψ⟩⟨ψ|) is invariant under Gaussian unitaries. A numerical evaluation of the NG entropy for pure NGE states is shown in Fig. 5. For operational purposes, we may also define a measure as the maximum overlap between a given pure state and an arbitrary pure GE state (see “Methods”). This measure distinctly separates GE and NGE states: it is zero for every pure GE state and non-zero for any NGE state, and invariant under Gaussian unitaries. As a “golden unit” of multipartite entanglement exists due to the complex structure of multi-partite entanglement, we further introduce a vector resource measure to characterize multipartite genuine non-Gaussian entanglement. This resource measure is based on the number of ancillary modes necessary to make an NGE state GE: the GE cost defined below. To get a clearer picture, let us begin with a simple bi-partite system of two modes A and B. As shown in Fig. 6a, any NGE state ρ^_{AB} can extended into a tri-partite state ρ^_{AB} ⊗ σ^_{C} by introducing an ancillary mode C, such that the state is GE between A and BC. By applying Gaussian unitary operations swapping systems B and C. Likewise, for m single-mode subsystems, any NGE state can be produced by introducing an (m–1)—mode non-Gaussian ancillary state entangled with the first mode.

For simplicity of analysis, let us consider an n-mode state with n-partition, while the more general case can be found in Supplemental Note 6. The GE states are those generated by applying n-mode Gaussian protocols on a n-mode fully separable state. Then, the corresponding GE cost is given as follows:
Definition 5. (GE cost vector and GE cost function): Given an arbitrary n-mode density matrix ρ over H_A, a GE cost vector is defined as a list of ancillary-mode numbers in a decreasing order:
R(ρ) = {|A'_1|, ..., |A'_n|}, (|A'_1| ≥ |A'_2| ≥ ... ≥ |A'_n|) (8)
of an extension ρ_{ext} satisfying ρ := Tr_{A'_1, ..., A'_n} [ρ_{ext}], such that ρ_{ext} is a GE state that can be obtained by applying a Gaussian protocol on a separable state between AA'_1, A_2A'_2, ..., and A_n A'_n. We will denote the k-th entry of R as R_k and the number of nonzero entries of R as ℓ.
The GE cost function is then defined as:
R_f(ρ) = min_{ρ_{ext}} R_1 (9)
where the minimization is taken over all possible GE extensions ρ_{ext} of the original state ρ.

Figure 6b shows an example of six modes, where ancilla modes A'_1 and A'_3 are introduced to extend the state to a GE state. The definition of a GE cost vector is closely related to the concept of Young diagram, as it is connected to irreducible representations of the permutation group. Note that although we have introduced the extension procedure with swap operations, general Gaussian protocols can be applied.

Note that the GE cost vector is not unique, but the GE cost function is. Moreover, the GE cost function R_f (as well as the one-norm of the GE cost vector) is bounded as kRk_1/ℓ ≤ R_f ≤ kRk_1. In the next section, we will show that the GE cost function provides useful prior information for learning, as it guides how resources must be allocated in a tomography protocol.

Learning pure NGE states
The GE cost is operationally meaningful, as it directly connects to the generation procedure of an NGE state. Moreover, as the clear structure of entanglement indicates, the cost of learning NGE states can be quantified by the GE cost, as we show in the following. In particular, to leverage the Gaussian entangling procedure in defining the GE cost, we consider a learning algorithm based on the one for GE states, with the local states corresponding to the minimal subsystems of modes {R_k + 1} that cannot be further disentangled by any Gaussian unitary (see details on “generalized GE states” in Supplemental Note 7).

When the number of modes m is large, the sample complexity is shown in Proposition 14 in “Methods”, which is an exponential function of max_k R_k. It is worth noting that, for any family of states, the GE cost plays the role of prior information that guides the choice between two strategies: (i) applying a counter-rotation followed by local tomography, or (ii) performing tomography directly on the GE state space. The GE cost determines which strategy is feasible and, in turn, sets the overall sample complexity.

In general, one can have a tight scaling with respect to the GE cost for the non-degenerate scenario:
Theorem 6. (Sample complexity for non-degenerate case): Consider an energy-constrained m−mode pure state in H^2_{m,E}, whose GE cost is specified by an l-element vector R and the function R_f. Meanwhile, its covariance matrix has non-degenerate symplectic eigenvalues with an eigen-gap Δ. The required sample complexity for tomography is:
M = poly(ℓ, E, 1/Δ, 1/ϵ, log(1/δ)) * Θ(E^{R_f}) (10)
where poly(ℓ, E, 1/ϵ, log(1/δ)) refers to a polynomial function of ℓ, E, the error ϵ in trace distance, and failure probability. Θ denotes the theta notation.

A detailed proof of Theorem 6 can be found in Supplemental Note 7E. In particular, we prove that the exponent of E/ϵ in the upper bound achieved by our protocol matches that of the lower bound demonstrated in Proposition 14. However, a residual gap remains, scaling as ℓ^2 R_f. Notably, the GE cost can be connected with the recently proposed t-doped bosonic Gaussian states7—states generated from t-doped Gaussian unitaries. Our multi-mode Gaussian unitaries are the NGE-class non-Gaussian states. In addition, the GE cost and the tomography overhead is directly related to the GE cost. In particular, any t-doped state can be extended into a GE state with (κ−1) ancillary modes. Furthermore, the GE cost is R=(κ^t−1, 0, ...).

Methods
Preliminary
An m-mode bosonic quantum system is completely specified by the statistical moments of its canonical quadrature operator vector r^=[q^1, p^1, ..., q^m, p^m]^T. Its first and second moments are captured by the mean vector ξ := Tr[r^ρ^] and the covariance matrix V := 1/2 Tr[{r^−ξ, r^−ξ}^T ρ^]10,11. The first-moment changes are modeled by the displacement operation D^ξ=exp(−iξ^T Ω r^), where ξ=(ξ_1, ..., ξ_{2m}) denotes the vector of displacements and Ω=⊕^m_{j=1} (0 1; −1 0) refers to the symplectic form. Then, given an arbitrary operator X^, one can define the Wigner characteristic function χ_X(ξ)=Tr[X^ D^ξ]. For notational simplicity, we will let χ_{kk'}(ξ) denote the characteristic function of the operator |k⟩⟨k'|, where {|k⟩, k=0, ..., 1} is the Fock basis.

A useful lemma for the proof of Theorem 1
The proof of Theorem 1 is based on the following lemma:
Lemma 7. (Generalizing Lemma 2 in ref. 52): The only states that remain pure under multi-mode pure loss are coherent states.
The detailed proof of this lemma is given in Supplemental Note 2.

Verification of pure NGE states
From Theorem 1, we have the following corollary for the identification of NGE states:
Corollary 8. (Verification of pure NGE states): Let |ψ⟩_{AB} be a bi-partite pure state and let S be a symplectic matrix that diagonalizes its covariance matrix, and U^g be the corresponding Gaussian unitary. (1) If the symplectic eigenvalues S are non-degenerate, and U^g |ψ⟩_{AB} is not a product state, then |ψ⟩_{AB} is NGE. (2) If there is no beam-splitter network U^O such that U^O U^g S^{−1} |ψ⟩_{AB} is a product state, then |ψ⟩_{AB} is NGE.

Maximal fidelity in state preparation
Consider the maximal fidelity of state preparation of the state |ψ⟩_A ∈ H_A, the state |ψ^0⟩_A ∈ H_A can be obtained as follows:
F_{NGE}(|ψ⟩⟨ψ|_A) = max_{|ϕ⟩_j ∈ H_{A_j}} |⟨ψ|_A U^g ⊗^n_{j=1} |ϕ_j⟩|^2 (11)
= max_{U^g ∈ GU(H_A)} |⟨ψ|_A U^g |ψ⟩⟨ψ|_A U^{g†}| × (12)
where |X|× = max_{|ϕ_j⟩ ∈ H_{A_j}} ⟨ϕ_j| X |ϕ_j⟩ refers to the cross norm. Note that the cross norm can be bounded by the operator norm as |X|× ≤ |X^T A'_0|_(A ⊗ A), which can be computed efficiently.

Volume of the state space
A crucial step in the proofs of Theorem 2 and Proposition 14 is to show that Algorithm 2 in Supplemental Note 7A requires fidelity-witnesses on exp(poly(m)/ϵ) passive-separable states, so that applying the union bound56 yields a polynomial sample complexity (see Supplemental Note 7C,D).

Intuitively, given that the state is known to lie within a specific class, the problem of determining how many states can be distinguished by tomography reduces to quantifying the amount of classical information required to specify any state in that class within an error ϵ. The associated quantity can also be referred to as the state description complexity or volume of the state space, whose logarithm determines the number of classical bits required to transmit the information. More rigorously, the following definition is given:
Definition 9. (ϵ-covering net): Let S denote a set of m-mode pure states. Then, there exists a ϵ-covering net S_{net} ⊆ S, such that for every |ψ⟩ ∈ S, there exists |ψ^0⟩ ∈ S_{net} satisfying |⟨ψ|ψ^0⟩|^2 ≥ 1−ϵ. Thus, the number of target states that require certification in m modes is given by the size of the ϵ-covering net. To visualize the upper bound, let us look at a more general case with GE states:
Corollary 10. (Volume of GE state space): Consider the set of m-mode pure GE states {|ψ⟩ ∈ H^1_{m,E}} produced by applying a Gaussian unitary on the tensor product of single-mode states. The size of the corresponding ϵ-covering net is bounded by:
|S_{net}| ≤ exp(O(poly(m^3, E, 1/ϵ^2))) (13)

The proof of Corollary 10 is shown in Supplemental Note 7. In contrast, the size of ϵ-covering net for all pure states in H^1_{m,E} is lower bounded as |S_{net}| ≥ exp(Θ(m^2 E / ϵ^2))7.

In our algorithm, we propose local fidelity witnesses to approximate the global fidelity witness for passive-separable states (see Fig. 3). This requires slightly fewer observables than the “volume” of passive-separable states, which is already smaller than that of GE states in Eq. (13)—at the cost of a controllable global fidelity witness error (see Supplemental Note 7A for details).

In our discussion of mixed-state learning, we use an ϵ_2-covering net of symplectic orthogonal matrices. Formally, let
S_{net,O0} ⊆ U(m) (14)
be a set of 2m×2m symplectic orthogonal matrices such that, for every symplectic orthogonal matrix O, there exists O^0 ∈ S_{net,O0} with |O_{jk} − O^0_{jk}| ≤ ϵ_2 ∀j,k ∈ {1, ..., 2m}, where {O_{jk}} denotes the matrix elements of O.

A useful lemma for mixed GE-state learning
In the main text, we mention a continuity proof for the learning protocol for pure GE states. Although the details are presented in Supplemental Note 8B, we present a crucial lemma here that may be of interest by itself.
First, the state of interest can be defined as follows:
ρ^_{mix} = C_A (⊗^m_{j=1} ρ^_j) (15)
where C_A is an arbitrary m-mode quantum channel acting on the subsystems {A_1, ..., A_m}, U^g(·)=U^g · U^{g†} refers to an arbitrary Gaussian unitary channel, |C_A|_E,2 = sup_{|Ψ⟩_{AB}} |(C_A ⊗ I_B) (|Ψ⟩⟨Ψ|_{AB})|_1 is the diamond norm with |Ψ⟩_{AB} ∈ S^2_{m,E} second moment of input energy, |ψ⟩ = U^g (⊗^m_{j=1} |ψ_j⟩) is the ideal state, and we have ϵ_{total} = ϵ_{mix} + ϵ_{channel} + ϵ_{input} as the total error. Then, we have a lemma concerning all errors.

Lemma 11. Let ρ^_{mix} be an approximate GE state described by Eq. (15). The estimation errors of the GDE protocol, designed for pure states, obey
|ρ^_{mix} − |ψ⟩⟨ψ||_1/2 ≤ ϵ_{total}, ϵ ≤ f( |ξ−ξ~|_2, |V_{pure}−V_e|_1/4 ) + ϵ_{ps}^0, (16)
|ξ−ξ~|_2 ≤ 8 √E ϵ_{mix} + 20 √ϵ_{mix} E, (17)
|V_{pure}−V_e|_1 ≤ 40 √E √ϵ_{mix} + ϵ_{ps}^0, (18)
ϵ_{ps}^0 ≤ ϵ_{channel} + f( |ξ−ξ~|_2, |V_{pure}−V_e|_1/4 ) + ϵ_{ps,net}, (19)
where ϵ denotes the reconstruction error in tomography, ξ and V_{pure} denote the displacement and covariance matrix of the state |ψ⟩, f(|ξ−ξ~|_2, |V_{pure}−V_e|_1/4) = O(|ξ−ξ~|_2) + O(|V_{pure}−V_e|_1/4) is pure error of estimating moments, the first terms in Eqs. (17)–(19) constitute a fixed error floor, the second terms of Eqs. (17) and (18) give the errors for an ideal pure input, ϵ_{ps}^0 is the error of estimating the approximated passive separable state after counter-rotation based on estimating moments, ϵ_{ps,net} denotes the corresponding error where we have an ideal pure GE state.

In this regard, the GDE learning protocol remains robust provided that the approximation errors ϵ_{channel} and ϵ_{input} vanish. When the input error ϵ_{input} vanishes, we have Theorem 3 in the main text. A detailed derivation of Lemma 11 is shown in Supplemental Note 8.

Learning degenerate mixed GE states
To go beyond Corollary 4 which only addresses the non-degenerate case, we begin with a complexity upper bound that guarantee efficiency, achieved by general joint measurements.
Theorem 12. (Application of Theorem 1.5 in ref. 63): Consider the set of m-mode states {ρ^} obtained by applying a Gaussian unitary to a tensor product of single-mode states ⊗^m_{j=1} ρ^_j, subject to an energy-moment constraint ρ^ ∈ S^2_{m,E}. There exists an tomography algorithm for this class of states that achieves the sample complexity M = O(poly(m, E, 1/ϵ, log(1/δ))).

A concrete proof of Theorem 12 is shown in Supplemental Note 8.

In general, Theorem 12 implies the existence of an efficient learning algorithm for mixed Gaussian-unitary-entanglable states, possibly using non-Gaussian joint measurements. Whether homodyne or heterodyne measurements with classical post processing suffice remains open and is left to future work.

In this regard, we propose in Supplemental Note 8 a generalized-GDE tomography protocol for arbitrary mixed states entangled by a Gaussian unitary. The protocol reduces back to the GDE protocol for pure states, where the input channel is a Gaussian unitary. The first few steps of generalized-GDE are identical to the GDE protocol, where we estimate the covariance matrix and then perform the disentangling Gaussian unitary. When the symplectic eigenvalues are degenerate, we continue with the next step of generalized-GDE. We classically build an ϵ-covering net of all 2m×2m symplectic orthogonal matrices, so that for any O there exists O^0 in the net with entrywise deviation at most ϵ_2.

ρ^_{ps,O0} = U^O0 (⊗^m_{j=1} ρ^_{j,O0}) U^{O0†}, where ρ^_{j,O0} is the reconstructed state from heterodyne measurement in Eq. (4). To precisely obtain all of these reconstructed states, we can build up L/m(1/ϵ_2)^{4m^2} local observables and rely on shadow tomography56 to estimate them from heterodyne data on O(log L) * O(poly(m)) copies.

For each of the ρ^_{ps,O0}, we perform a reduced-state test, which is the replacement of fidelity witness that only works for pure states. In this test, we apply virtual passive rotations U^y_{Of} classically, compute reduced states of U^y_{Of} ρ^_{ps,O0} U^y_{Of}, and compare them with the locally learned ρ^ in Eq. (4). In principle, scanning O ∈ U(m) determines whether O^0 is correct. In this regard, we have the following Lemma 13.

Lemma 13. Given two arbitrary m-mode states ρ^_1 and ρ^_2 that satisfy the following condition:
Tr[|0⟩⟨0| ⊗ ... ⊗ |0⟩⟨0| U^y_{Of} ρ^_1 U^y_{Of}] = Tr[|0⟩⟨0| ⊗ ... ⊗ |0⟩⟨0| U^y_{Of} ρ^_2 U^y_{Of}], (20)
(∀j=1, ..., m; ∀O ∈ U(m))
where U(m) denotes the set of all symplectic orthogonal matrices corresponding to arbitrary m-mode passive Gaussian operations. Then, we have:
ρ^_1 = ρ^_2. (21)

The concrete proof of Lemma 13 is shown in Supplemental Note 8.

Note that Lemma 13 requires checking all passive transforms with O ∈ U(m); however, we only have access to the reconstructed reduced states in the r.h.s. of Eq. (20) following a search over all possible O ∈ S_{net,O0} from the early steps of the generalized GDE protocol. Nevertheless, as O^0 has only m^2 independent parameters, we conjecture that going over the O_f’s among S_{net,O0} in Eq. (14) is sufficient for guaranteeing small errors in the local reduced-state test. This conjecture is subject to future study.

Learning states with the same GE cost
Regarding the states with the same GE cost, we have the following proposition:
Proposition 14. (Sample complexity of NGE states): Consider the set of m-mode pure states H^2_{m,E} with a GE cost vector R=(R_1, ..., R_ℓ)^T satisfying R_f > 7. Then, the tomography of such states will take M copies that satisfy:
M ≥ Θ((1−δ) E^{R_f+1} / m log(E/ϵ)) (22)
M ≤ O(poly(ℓ^{5R_f+9}, E^{2R_f+3}, 1/ϵ^{2R_f+6}, log(1/δ)^{32R_f+2})) (23)
where δ denotes the failure probability and ϵ refers to the trace distance between the input and the reconstructed state. Note that we ℓ=O(m) by definition.

A detailed proof of Proposition 14 is shown in Supplemental Note 7C. In addition, it is shown that the sample-complexity scaling is dominated by covariance-matrix estimation when R_f is small (e.g., R_f ≤ 7).

Peer review rounds
[{"round": 1, "reviewer_comments": {"reviewer_1": "The manuscript “Complexity of quantum tomography from genuine non-Gaussian entanglement” by X. Zhao et al. deals with\nthe relevant topic of quantum state tomography focusing on relation between the quantum correlation of the state under\ninvestigation and the scaling of the number of copies required by the tomographic process. The authors prove that the\ntomography of a multimode pure “Gaussian-entanglable” state requires a number of copies that scales polynomially with\nrespect to the number of modes. To address the non Gaussian-entanglable states they introduce an operation-driven\nmonotone corresponding to the number of ancillary modes that one should add to reduce the states to Gaussian-entanglable\nones.\nOverall, the manuscript is indeed technical and it is of interest for the specific field of research, namely, quantum state\ntomography, but also for the quantum information community, since it points out the relation between the quantum\ncorrelations of the considered systems and the complexity of quantum tomography.\nThe manuscript is clear, well-written and well-organized. The reference list is rather complete and also the method section is\nwell-supporting the main text.\nAs a minor comment, at pag. 3, in the definition of the two-mode squeezed vacuum state, there should be “-1/2” instead of\n“-1”.\nI recommend the publication of the manuscript in its present form.", "reviewer_2": "The authors theoretically investigate quantum tomography of pure non-Gaussian quantum states. The authors show that\nefficient reconstruction can be performed if the state is so-called Gaussian-entanglable (GE), i.e. if it can be generated from a\nproduct state by Gaussian operation. The authors derive bounds on sample complexity and they also obtain and present\nseveral other interesting results concerning GE states and their transformations, and the non-GE states. Intuitively, efficient\ntomography of the GE states is possible because for m-mode (pure) GE states the number of parameters that characterize\nthe state is given by the number of parameters that characterize the disentangling Gaussian operation and the number of\nparameters that characterize m-mode product state, which results in polynomial scaling. The concept of GE states appears\nto be very useful and represents an important addition to theory of quantum resources for bosonic quantum systems. Given\nthe current wide interest in continuous-variable quantum information processing at various platforms (light, microwaves\ncombined with superconducting systems, optomechanics, etc.) which typically requires engineering of various highly\nnonclassical and entangled states as a crucial resource, the present work is timely, and potentially of broad impact and\ninterest.\nThe authors should consider and address the following points:\n1) The authors mostly focus on reconstruction of pure states. However, in practice one cannot assume that the\nexperimentally generated state is exactly pure. It would be useful to provide more details on the possibility to extend the\npresented results to mixed states, considering especially the cases when the state is close to a pure state, but it cannot be\nconverted to a product of (generally mixed) single-mode states by a single unitary Gaussian operation.\n2) In the text in the Supplementary information below Eq. (18) the symbol omega is used, while it appears that a symbol w is\nused instead before, c.f. e.g. Eq. (17). Please check the notation.\n3) Proof that the NOON states with N>2 are NGE: The authors first prove that if this state were GE it would have to be\npossible to disentangle it by a passive Gaussian unitary operation. Since passive Gaussian unitary operations preserve the\ntotal number of photons, this implies that it suffices to prove that a NOON with N>2 state cannot be generated from any two-\nmode Fock state |k>|N-k> by a passive unitary Gaussian operation, i.e. by interference at a beam splitter. This is what the\nauthors essentially do in their calculations on page 11 of the Supplementary Information but it might be useful to emphasize\nthe above described simple physical picture there.\n4) Algorithm 1 specifies a procedure for estimation of the covariance matrix and displacement of the m-mode state. The\nauthors propose to utilize homodyne and heterodyne measurements and also photon counting. It would be helpful to clarify\nwhether the photon counting is really necessary, since the first and second moments of quadrature operators can be\nestimated from homodyne and/or heterodyne measurements. In fact, by measuring simultaneously each mode with a\nseparate heterodyne detector, one can reconstruct full covariance matrix and vector of displacements from this\nmeasurement. Perhaps there is some need to ensure independence of the estimators of each matrix element but additional\nclarification would be very helpful. Note that from the experimental point of view determination of the second quadrature\nmoments by photon counting may be much more challenging than the homodyne or heterodyne measurements.\n5) Another issue related to the previous point: the detectors utilized in experiments need not be perfect and in particular their\ndetection efficiency can be smaller than 1. Can the proposed tomography be performed also with such realistic imperfect\ndetectors?\n6) In Figure 5 the authors present results of numerical simulations which demonstrate efficiency of their tomographic\nprotocol. The results are presented for a thermal state correlated at a beam splitter. It would be interesting to see such results\nalso for other kinds of states, especially some highly non-classical states such as coherent superpositions of coherent\nstates.", "reviewer_3": "The authors have identified a class of important quantum states known as the Gaussian-process-entangleable (GE) states, which include important resource states like, ordinary squeezed Gaussian and GKP states. With these states, they proved that learning them is efficient in the sense that the resource complexity is polynomial in the number of modes. A key result is the epsilon-covering bounds in Theorem 2. For states that are not GE, the authors gave an optimal mode-added transformation scheme that turns them into GE ones with more modes interacting with the Gaussian transformation. In doing so, as a corollary, they showed a no-go result of generating a NOON with just a simple beam splitter for N>2.\nOverall, this paper is relevant to the quantum learning community that focuses on learning the parameters of a fixed model (such as that of a state, Hamiltonian, etc.) Together with the consideration of GE states, which are also important in the current quantum computation landscape, I believe that this work will add value to the literature on learning such states efficiently.\nThe results support the message, insofar, so long as pure GE states are concerned. The authors are made note that the learning complexity of mixed GE states shall be relegated to future work. While the methodology is sound and the article indeed carries scientific weight to a degree that could be considered publishable in Nature Communications, I have two concerns which the authors should address before I confirm my aforementioned decision.\nTechnical question: In complexity theory, it is customary to ignore constant prefactors in upper bound big O expressions (by definition). Mathematically, the polynomial scaling behavior in the number of modes m is what matters for astronomically large m. However in practice, the prefactor is actually in important for small to moderate m. It would be instructive if the authors state the prefactor either as a precise expression or an estimate. This would offer the more operationally correct guide for GE-state learning.\nGeneral question: When the pure GE state is corrupted by losses, the resulting state is a noisy mixed non-GE state. Unless I missed something, it appears that the article did not address this important practical issue, which determines the robustness of GE-state learning in the presence of losses. I believe it is rather crucial for authors to address this problem with at least a simulation in order to demonstrate loss robustness."}, "author_response": "RESPONSE TO REVIEWS\n“Complexity of quantum tomography from genuine non-Gaussian entanglement”\nsubmitted to Nature Communications\nby Xiaobin Zhao, Pengcheng Liao, Francesco Anna Mele, Ulysse Chabaud, and Quntao Zhuang\nWe thank the Editor very much for handling our manuscript, as well as the reviewers for their constructive\nfeedback. Wearepleasedthatallthreereviewersrespondedpositivelytoourwork. Inparticular,weareencouraged\nthat Reviewer 1 recommends our study for publication, Reviewer 2 highlights its timeliness, broad impact, and\ninterest, and Reviewer 3 acknowledges its value. We are delighted by this overall positive assessment.\nReviewer 3 suggested that the paper is suitable for publication pending two concerns, which we have fully\naddressed in the revised manuscript. In addition, Reviewer 2 provided several valuable comments, all of which we\nhave carefully considered and incorporated into the revised version of the manuscript. In the following, we provide\na detailed, point-by-point response to each reviewer’s comments and summarize the corresponding changes made\nto the manuscript.\nWe sincerely thank the reviewers for their thoughtful and constructive reports. Having carefully addressed\nall comments, we hope that the revised manuscript will now be considered suitable for publication in Nature\nCommunications.\nReviewer 1\nComment: The manuscript “Complexity of quantum tomography from genuine non-Gaussian entanglement”\nbyX.Zhaoetal. dealswiththerelevanttopicofquantumstatetomographyfocusingonrelationbetweenthequantum\ncorrelation of the state under investigation and the scaling of the number of copies required by the tomographic\nprocess. Theauthorsprovethatthetomographyofamultimodepure“Gaussian-entanglable”staterequiresanumber\nof copies that scales polynomially with respect to the number of modes. To address the non Gaussian-entanglable\nstates they introduce an operation-driven monotone corresponding to the number of ancillary modes that one should\nadd to reduce the states to Gaussian-entanglable ones.\nOverall, the manuscript is indeed technical and it is of interest for the specific field of research, namely, quantum\nstate tomography, but also for the quantum information community, since it points out the relation between the\nquantum correlations of the considered systems and the complexity of quantum tomography.\nThe manuscript is clear, well-written and well-organized. The reference list is rather complete and also the\nmethod section is well-supporting the main text.\nAs a minor comment, at pag. 3, in the definition of the two-mode squeezed vacuum state, there should be “-1/2”\ninstead of “-1”.\nI recommend the publication of the manuscript in its present form.\nReply: We thank Reviewer 1 for the positive comments and for recognizing the relevance to quantum state\ntomography, quantum information, and quantum correlations communities. We also appreciate his/her recommen-\ndation to publish the manuscript in its present form.\nWe thank Reviewer 1 for the comment on the definition of the two-mode squeezed vacuum at page 3. We\nfollow Eq. (5.22) in the book [Serafini, Alessio. Quantum continuous variables: a primer of theoretical methods.\nCRC press, 2023.], under which the “-1” factor is correct. We therefore leave the expression unchanged. The\nnormalization is based on the equality ∞ tanh2kr =cosh2(r).\nk=0\n(cid:80)\nReviewer 2\nComment:\nThe authors theoretically investigate quantum tomography of pure non-Gaussian quantum states. The authors\nshow that efficient reconstruction can be performed if the state is so-called Gaussian-entanglable (GE), i.e. if it can\nbe generated from a product state by Gaussian operation. The authors derive bounds on sample complexity and they\nalso obtain and present several other interesting results concerning GE states and their transformations, and the\nnon-GE states. Intuitively, efficient tomography of the GE states is possible because for m-mode (pure) GE states\nthe number of parameters that characterize the state is given by the number of parameters that characterize the dis-\nentangling Gaussian operation and the number of parameters that characterize m-mode product state, which results\nin polynomial scaling. The concept of GE states appears to be very useful and represents an important addition to\ntheory of quantum resources for bosonic quantum systems. Given the current wide interest in continuous-variable\nquantum information processing at various platforms (light, microwaves combined with superconducting systems,\noptomechanics, etc.) which typically requires engineering of various highly nonclassical and entangled states as a\ncrucial resource, the present work is timely, and potentially of broad impact and interest.\nReply: WethankReviewer2fortheprecisesummaryofourresultsandforrecognizingourcontributiontothe\nstudy of continuous-variable systems and its potential impact.\nComment: The authors should consider and address the following points:\n1) The authors mostly focus on reconstruction of pure states. However, in practice one cannot assume that\nthe experimentally generated state is exactly pure. It would be useful to provide more details on the possibility to\nextend the presented results to mixed states, considering especially the cases when the state is close to a pure state,\nbutitcannotbeconvertedtoaproductof(generallymixed)single-modestatesbyasingleunitaryGaussianoperation.\nReply: We are grateful to Reviewer 2 for bringing up the practical issue that experimentally prepared states\nmay not be perfectly pure and for the helpful suggestion to extend our results to mixed states. In response, we\nadded a new section “Towards efficient learning of mixed GE states” starting on page 6 and a new Supplementary\nNote VIII “Efficient learning of mixed GE states”. In particular, we address two complementary scenarios.\n• First,wedemonstratethatourGDElearningalgorithmisrobust againstsmallperturbationscausedbygeneral\n(potentially non-Gaussian) noise. Even if the unknown state is a slightly perturbed GE state—technically,\na state constructed by sending an input state, which is (cid:15) -close to a pure product state, into a channel,\ninput\nwhich is (cid:15) -close to a Gaussian unitary—our GDE learning algorithm remains effective and tomography\nchannel\nremains efficient. As we now added in the Method section “A useful lemma for mixed GE-state learning”:\nFirst, the state of interest can be defined as follows:\nm\nρˆ = ρˆ (1)\nmix A j\nC  \nj=1\n(cid:79)\n(cid:40) 1 1 2 2 min U m j= A g 1 ∈ ρ G ˆ j U( − HA  )(cid:107) m j= C 1 A | ψ −  j (cid:105) U (cid:104) A g ψ (cid:107) j | E (cid:5) , 1 2 ≤ ≤ (cid:15) (cid:15) i c n h p a u n t nel ,\n(cid:13) (cid:13)\nwhere\nC A\nisanarbitrarym-mo (cid:13) (cid:13)d (cid:78) equantum (cid:78) channelacting (cid:13) (cid:13)onthesubsystems\n{\nA\n1\n,...,A\nm }\n,\nUA\ng(\n·\n)=\nUˆg Uˆg refers to an arbitrary Gaussian unitary channel,\n†\n·\nE,2 = sup ( I )(Ψ Ψ )\n(cid:107)C A (cid:107) (cid:107) C A ⊗ B | AB (cid:105)(cid:104) AB | (cid:107)1\n(cid:5) TrB[ | Ψ (cid:105)(cid:104) Ψ | | Ψ A (cid:105) B A ] B ∈ St( H 2 m,E)\nis the diamond norm with bounded second moment of input energy, ψˆ =Uˆg( m ψ ψ )Uˆg\npure ⊗j=1| j (cid:105)(cid:104) j | †\nis the ideal state, and we have ρˆ ,ψˆ , m ρˆ , m ψ ψ St( 2 ). In the following, we\nmix pure ⊗j=1 j ⊗j=1| j (cid:105)(cid:104) j | ∈ Hm,E\ndenote (cid:15) =(cid:15) +(cid:15) as the total error.\nmix channel input\nThen, we have a lemma concerning all errors.\nLemma 11 Let ρˆ be an approximate GE state described by Eq. (1). The estimate errors of\nmix\nGDE proto- col, designed for pure states, obey\n(cid:15) f( ξ ξ˜ , V V 1/4)+(cid:15) (2)\n≤ (cid:107) − (cid:107)\n2\n(cid:107)\npu\n(cid:15)\nre\n− (cid:107)∞\n(cid:48)ps\nξ ξ˜ 8 E(cid:15) + v , (3)\n(cid:107) − (cid:107) 2 ≤ mix 20√mE (cid:101)\n(cid:112)\nV V 40√2E√(cid:15) +(cid:15) (4)\nmix v\n(cid:107) − (cid:107)∞ ≤\n(cid:15) (cid:15) +f( ξ ξ˜ , V V 1/4)+(cid:15) (5)\n(cid:101)\n(cid:48)ps\n≤\nchannel\n(cid:107) − (cid:107)\n2\n(cid:107)\npure\n− (cid:107)∞\nps\nwhere (cid:15) denotes the reconstruction error in tomography, ξ and V p(cid:101)ure denote the displacement and\ncovariance matrix of the state ψ , f( ξ ξ˜ , V V 1/4)= ( ξ ξ˜ )+ ( V V 1/4)\npure 2 pure 2 pure\n(cid:107) − (cid:107) (cid:107) − (cid:107)∞ O (cid:107) − (cid:107) O (cid:107) − (cid:107)∞\nis the error of estimating moments, the first terms in Eqs. (3)–(5) constitute a fixed error floor, the\nsecond termsofEqs.(3)and(4)givetheerrorsforanidea(cid:101)lpureinput, (cid:15)\n(cid:48)ps\nistheerrorofestim(cid:101)ating\nthe approximated passive separable state after counter-rotation based on estimating moments, (cid:15)\nps\ndenotes the corresponding error where we have an ideal pure GE state.\nWe then specialize to the zero input-error case and summarize it in Theorem 3 of the revised manuscript,\nappended below.\nTheorem 3 Consider approximated pure GE states ρˆ St( 2 ) whose entangling channel is\n∈ Hm,E\n(cid:15) -close to an arbitrary Gaussian unitary. The sample complexity of learning this class of\nchannel\nstate is:\n1 1\nM poly m,E, ,log , (6)\n(cid:39)O (cid:15) (cid:15) δ\n(cid:18) (cid:18) − channel (cid:18) (cid:19)(cid:19)(cid:19)\nfor an error (cid:15)>(cid:15) in trace distance and failure probability δ.\nchannel\nThe above continuity proof fully resolves reviewer 2’s suggestion of states close to pure GE states. When the\nentangling process is a bit noisy and cannot be reversed by a Gaussian unitary, the above results guarantee\nthat we can still estimate the state pretty well, up to the noise in the process.\n• Second, for states obtained by applying a Gaussian unitary to a product of generally mixed input states, we\nhave proven a new Corollary 4 in main text, which we append here for convenience:\nCorollary 4 (Sample complexity of states entangled by a Gaussian unitary) Considerthe\nsetofm-modestates ρˆ obtainedbyapplyingaGaussianunitarytoatensorproductofsingle-mode\nstates m ρˆ , subje { ct } to an energy-moment constraint ρˆ St( 2 ). Meanwhile, the states’\nj=1 j ∈ Hm,E\nsymplectic eigenvalues are nondegenerate, up to a finite eigen-gap ∆. The GDE protocol needs\n(cid:78)\nM poly m,E,∆ 1,(cid:15) 1,log δ 1 copies for this state family, for an error (cid:15) in trace dis-\n− − −\n(cid:39) O\ntance and failure probability δ.\n(cid:0) (cid:0) (cid:0) (cid:1)(cid:1)(cid:1)\nInSupplementalNoteVIII,weshowthatGDEprotocolisabletoachievethesamplecomplexityinCorollary\n4, with local state tomography switched from fidelity witness in the pure-state case to heteredyne local\ntomography. The case of degenerate symplectic eigenvalues is not covered by the above Corollary. In the\nnewly added section “Learning degenerate mixed GE states” in Method, we show that there exists protocols\nrelying on general measurements that in-principle guarantees sample-efficient learning of such degenerate\nstates, appended below:\nTheorem 12 (Application of Theorem 1.5 in [1]) Consider the set of m-mode states ρˆ ob-\ntained by applying a Gaussian unitary to a tensor product of single-mode states m ρˆ , su { bje } ct to\nj=1 j\nan energy-moment constraint ρˆ St( 2 ). There exists an tomography algorithm for this class\n∈ Hm,E (cid:78)\nof states that achieves the sample complexity M = poly m,E ,1,log 1 .\nO II (cid:15) δ\n(cid:0) (cid:0) (cid:0) (cid:1)(cid:1)(cid:1)\nAdditionally, we provide a generalized GDE protocol with the hope of efficient learning such states with\npractical Gaussian operations, however, subject to an unproven conjecture.\nIn this regard, we propose in Supplemental Note VIII a generalized-GDE tomography protocol for\narbitrary mixed states entangled by a Gaussian unitary. The protocol reduces back to the GDE\nprotocol for the non-degenerate case (and thus efficient); for the degenerate case, the protocol is\nefficient up to an unproven conjecture. In the following, we briefly summarize the generalized GDE\nprotocol.\nThe first few steps of generalized-GDE are identical to the GDE protocol, where we estimate the\ncovariance matrix and then perform the disentangling Gaussian unitary. When the symplectic\neigenvalues are non-degenerate, we would conclude with the local tomography as in Corollary 4. If\nthesymplecticeigenvaluesaredegenerate,thenwecontinuewiththenextstepofgeneralized-GDE.\nWeclassicallybuildan(cid:15)-coveringnetofall2m 2msymplecticorthogonalmatrices,sothatforany\n×\nO there exists O in the net with entrywise deviation at most (cid:15) (see Methods). For each O in the\n(cid:48) 2 (cid:48)\nnet S net,O(cid:48) , wereconstructρ ps,O(cid:48) = Uˆ O(cid:48) m j=1 ρˆ j,O(cid:48) Uˆ O† (cid:48) , whereρˆ j,O(cid:48) isthereconstructedreduced\nstatefromheterodynemeasurementinEq(cid:104).(4). Topr(cid:105)eciselyobtainallofthesereconstructedstates,\n(cid:78)\nwe can build up L m(1/(cid:15)(cid:101))4m2 local observables and rely on shadow tomography [2] to estimate\n2\n∝\nthem from heterodyne data on (logL) (poly(m)) copies.\nO ∼O\nFor each of the ρ , we perform a reduced-state test, which is the replacement of fidelity witness\nps,O(cid:48)\nthat only works for pure states. In this test, we apply virtual passive rotations Uˆ O†\nf\nclassically,\ncompute reduce(cid:101)d states of Uˆ O†\nf\nρ ps,O(cid:48) Uˆ Of , and compare them with the locally learned ρˆ j,Of in\nEq. (4). In principle, scanning O U(m) determines whether O is correct. In this regard, we\nf (cid:48)\n∈\nhave the following Lemma 13. (cid:101)\nLemma 13 Given two arbitrary m-mode states ρˆ and ρˆ that satisfy the following condition:\n1 2\nTr j(cid:48)= j 1 (cid:48) ,(cid:54) =\n··\nj\n·\n,m\n(cid:104)\nUˆ O† f ρˆ 1 Uˆ Of\n(cid:105)\n≡ Tr j(cid:48)= j 1 (cid:48) ,(cid:54) =\n··\nj\n·\n,m\n(cid:104)\nUˆ O† f ρˆ 2 Uˆ Of\n(cid:105)\n,\n( j =1, ,m; O U(m)) (7)\nf\n∀ ··· ∀ ∈\nwhereU(m)denotesthesetofallsymplecticorthogonalmatricescorrespondingtoarbitrarym-mode\npassive Gaussian operations. Then, we have:\nρˆ ρˆ . (8)\n1 2\n≡\nThe concrete proof of Lemma 13 is shown in Supplemental Note VIII.\nNote that Lemma 13 requires checking all passive transforms with O U(m); however, we only\nf\n∈\nhave access to the reconstructed reduced states in the r.h.s. of Eq. (7) following a search over all\npossibleO fromtheearlystepsofthegeneralizedGDEprotocol. Nevertheless, asO has\nf\n∈S\nnet,O(cid:48) (cid:48)\nonly m2 independent parameters, we conjecture that going over the O ’s among in Eq. (14)\nf\nS\nnet,O(cid:48)\nis sufficient for guaranteeing small errors in the local reduced-state test. This conjecture is subject\nto future study.\n• Finally, to show case the robustness of our approach to loss. We have added simulations for various exotic\nquantum states subject to uniform loss, as we depict in Fig. 4 of the main text. The figure is also appended\nas Fig. R1 of the reply.\nEntangled number state Entangled cat state Correlated thermal state\nFigure R1: Performance of Gaussian-disentangling protocol (GDE) versus direct tomography (RρR). Here, we\nnumerically simulate the quantum state tomography process on a two-mode GE state with the true state being\nnumber state 0 1 , cat state cat cat and thermal state ρ ρ correlated by a beam-\n| (cid:105)⊗| (cid:105) | (cid:105)0.1 ⊗| (cid:105)1.1 th,n=0.2 ⊗ th,n=0.3\nsplitter (with a rotation angle θ = π/4). The error, defined as the trace distance between the true state and\nthe reconstructed state obtained through direct tomography using the standard RρR algorithm [3], is shown as a\nfunction of the number of copies M used in tomography (black squares for ideal case and grey stars for the lossy\ncase). In the lossy case, pure loss channel (with a transmissivity ratio η = 0.9) is applied to both modes after the\nbeamsplitter. Incontrast, thebottomenvelopeofachievableerrorsintheproposedGaussiandisentangling(GDE)\nalgorithm, after minimising over possible sample number in estimating displacements and the covariance matrix,\nis\ngiven by red circles and pink triangles (lossy case). Their linear fit shows errors in GDE protocol decrease faster\n(with slope k) compared to RρR method in all the three different states, demonstrating the advantage of GDE\nmethod.\nComment: 2) In the text in the Supplementary information below Eq. (18) the symbol omega is used, while it\nappears that a symbol w is used instead before, c.f. e.g. Eq. (17). Please check the notation.\nReply: We thank the referee for pointing out this typo. We have corrected the notation accordingly.\nComment: 3) Proof that the NOON states with N >2 are NGE: The authors first prove that if this state were\nGE it would have to be possible to disentangle it by a passive Gaussian unitary operation. Since passive Gaussian\nunitary operations preserve the total number of photons, this implies that it suffices to prove that a NOON with\nN > 2 state cannot be generated from any two-mode Fock state k > N k > by a passive unitary Gaussian\n| | −\noperation, i.e. by interference at a beam splitter. This is what the authors essentially do in their calculations on\npage 11 of the Supplementary Information but it might be useful to emphasize the above described simple physical\npicture there.\nReply: We appreciate the referee for suggesting presenting the connection to previous no-go results more\nclearly. We have added a comment on page 6:\nNote that, for this specific problem, previous approaches have only established such a no-go result for\nnumberstateinputsinalinearopticalsetting[4–6],wherethephotonnumberispreservedandtherefore\nit suffices to disprove the case with two-mode Fock state. Our approach completes the general case by\nreducing general Gaussian unitaries to passive unitaries.\nWe show that it is sufficient to prove the impossibility of generating NOON states from applying passive operation\nonFockstates. OurapproachcompletesthegeneralcasebyreducinggeneralGaussianunitariestopassiveunitaries.\nComment: 4) Algorithm 1 specifies a procedure for estimation of the covariance matrix and displacement of the\nm-mode state. The authors propose to utilize homodyne and heterodyne measurements and also photon counting.\nIt would be helpful to clarify whether the photon counting is really necessary, since the first and second moments\nof quadrature operators can be estimated from homodyne and/or heterodyne measurements. In fact, by measuring\nsimultaneously each mode with a separate heterodyne detector, one can reconstruct full covariance matrix and vector\nof displacements from this\nmeasurement. Perhaps there is some need to ensure independence of the estimators of each matrix element but additional\nclarification would be very helpful. Note that from the experimental point of view determination of the second quadrature\nmoments by photon counting may be much more challenging than the homodyne or heterodyne measurements.\nReply: We agree that photon counting is more challenging experimentally than homodyne or heterodyne\nmeasurements. We have revised Algorithm 1 in Supplementary Note VII to show how to estimate the number\noperator a†j a\nj\nfrom heterodyne measurement result and to specify the corresponding classical processing to realize\nsuch an observable.\nWe now emphasized that our protocol only relies on Gaussian measurements, as we state in the abstract:\n...by providing an explicit protocol involving only Gaussian unitary, local tomography (e.g. homodyne\nor heterodyne) and classical processing\nComment:5) Another issue related to the previous point: the detectors utilized in experiments need not be per-\nfect and in particular their detection efficiency can be smaller than 1. Can the proposed tomography be performed\nalso with such realistic imperfect detectors?\nReply: We appreciate Referee 2 for raising this point. Building on our changes in the main text in response\nto Comment 1, we added a paragraph on page 6 addressing this issue, appended here for convenience:\nBelow,weapplyCorollary4andTheorem3topracticalexperimentalscenarios. Inthisregard,detector\ninefficiencyisanubiquitousandpracticallyrelevantimperfectionthatleadstomixedNGEstates. Here,\nit can be modeled as ideal detection preceded by pure-loss channels. As uniform loss commutes with\nbeamsplitters, uniform efficiency maps boson-sampling outputs to a passive Gaussian unitary acting\non a product of mixed states. Consequently, Corollary 4 applies in the situation with nondegenerate\nsymplectic eigenvalues, where both the sample and time complexity stay polynomial in the number of\nmodes. For nonuniform loss, the degraded channel remains close to a Gaussian unitary in the small loss\nregime. Thus, Lemma 11 in Methods guarantees a vanishing reconstruction error floor, and Theorem 3\nshows that the sample complexity increases only polynomially with the degradation level.\nComment: 6) In Figure 5 the authors present results of numerical simulations which demonstrate efficiency\nof their tomographic protocol. The results are presented for a thermal state correlated at a beam splitter. It would\nbe interesting to see such results also for other kinds of states, especially some highly non-classical states such as\ncoherent superpositions of coherent states.\nReply: We thank the referee for the suggestion to improve our numerical results. We have ran additional\nsimulations and updated Fig. 5 (Fig. 4 in the revised manuscript) to three panels showing the performance of\nlearning GE states whose input are number states, coherent superpositions of coherent states (i.e., cat states), and\nthermal states. The figure is also appended in Fig. R1 for reference. Here, each panel also includes the cases\nwith detector inefficiency, modeled as a pure-loss channel acting on each mode. We show that, in each case,\nthe\nGaussian-disentangling method achieves smaller error and markedly better scaling than direct tomography.\nReviewer 3\nComment: The authors have identified a class of important quantum states known as the Gaussian-process-\nentangleable (GE) states, which include important resource states like, ordinary squeezed Gaussian and GKP states.\nWith these states, they proved that learning them is efficient in the sense that the resource complexity is polynomial\nin the number of modes. A key result is the epsilon-covering bounds in Theorem 2. For states that are not GE,\nthe authors gave an optimal mode-added transformation scheme that turns them into GE ones with more modes\ninteracting with the Gaussian transformation. In doing so, as a corollary, they showed a no-go result of generating\na NOON with just a simple beam splitter for N >2.\nOverall, this paper is relevant to the quantum learning community that focuses on learning the parameters of a\nfixed model (such as that of a state, Hamiltonian, etc.) Together with the consideration of GE states, which are also\nimportant in the current quantum computation landscape, I believe that this work will add value to the literature on\nlearning such states efficiently.\nThe results support the message, insofar, so long as pure GE states are concerned. The authors are made note\nthat the learning complexity of mixed GE states shall be relegated to future work. While the methodology is sound\nand the article indeed carries scientific weight to a degree that could be considered publishable in Nature Communi-\ncations, I have two concerns which the authors should address before I confirm my aforementioned decision.\nReply: We are gratefull for the reviewer’s detailed summary of our results and his/her recongnition of our\ncontribution\nComment: Technical question: In complexity theory, it is customary to ignore constant prefactors in upper\nbound big O expressions (by definition). Mathematically, the polynomial scaling behavior in the number of modes m\nis what matters for astronomically large m. However in practice, the prefactor is actually in important for small to\nmoderate m. It would be instructive if the authors state the prefactor either as a precise expression or an estimate.\nThis would offer the more operationally correct guide for GE-state learning.\nReply: Weappreciatethereferee’scommentandagreethattheprefactorofsamplecomplexityiscrucialinthe\ncase with a moderate number of modes. We added a new Supplemental Note VIIIF entitled “Prefactor in sample\ncomplexity with moderate mode number”, in which we present the explicit expression of the sample complexity\nupper bound for pure GE state tomography. In the main text accordingly, we add a comment on page 4 pointing\nto Supplemental Note VIIF:\nFor the case with a moderate number of modes in which constant prefactors in the sample complexity\nbecomeimportant,weprovideanexplicitupperboundinSupplementalNoteVIIFtooffermorepractical\nguidance for GE-state learning.\nComment: General question: When the pure GE state is corrupted by losses, the resulting state is a noisy\nmixed non-GE state. Unless I missed something, it appears that the article did not address this important practical\nissue, which determines the robustness of GE-state learning in the presence of losses. I believe it is rather crucial\nfor authors to address this problem with at least a simulation in order to demonstrate loss robustness.\nReply: We thank Reviewer 3 for raising the critical practical issue that, under loss, a pure GE state is not\nexperimentallyfeasibleandbecomesanon-GEstate. Weappreciatethesuggestiontoincludeatleastonesimulation\nto demonstrate robustness to loss. In the revised manuscript, we analyze this issue analytically and numerically.\nIn particular, we added a new section “Towards efficient learning of mixed GE states” beginning on page 6 and a\nnew Supplementary Note VIII “Efficient learning of mixed GE states.” These additions treat two complementary\nscenarios.\n• First,wedemonstratethatourGDElearningalgorithmisrobust againstsmallperturbationscausedbygeneral\n(potentially non-Gaussian) noise. Even if the unknown state is a slightly perturbed GE state—technically, a\nstate constructed by sending an input state, which is (cid:15)-close to a pure product state, into a channel, which\nis (cid:15)-close to a Gaussian unitary—our GDE learning algorithm remains effective and tomography remains\nefficient. In this way, we can include the case where loss is low. We append the formal definition and the\nassociated lemma (which are presented in the Method) below for convenience:\nFirst, the state of interest can be defined as follows:\nm\nρˆ = ρˆ (9)\nmix A j\nC  \nj=1\n(cid:79)\n(cid:40) 1 1 2 2 min U m j= A g 1 ∈ ρ G ˆ j U( − HA  )(cid:107) m j= C 1 A | ψ −  j (cid:105) U (cid:104) A g ψ (cid:107) j | E (cid:5) , 1 2 ≤ ≤ (cid:15) (cid:15) i c n h p a u n t nel ,\n(cid:13) (cid:13)\nwhere isanarbitrarym-mo(cid:13)d (cid:78) equantum (cid:78) channelacting(cid:13)onthesubsystems A ,...,A , g()=\nC A (cid:13) (cid:13) { 1 m } UA ·\nUˆg Uˆg refers to an arbitrary Gaussian unitary channel,\n†\n·\nE,2 = sup ( I )(Ψ Ψ )\n(cid:107)C A (cid:107) (cid:107) C A ⊗ B | AB (cid:105)(cid:104) AB | (cid:107)1\n(cid:5) TrB[ | Ψ (cid:105)(cid:104) Ψ | | Ψ A (cid:105) B A ] B ∈ St( H 2 m,E)\nis the diamond norm with bounded second moment of input energy, ψˆ =Uˆg( m ψ ψ )Uˆg\npure ⊗j=1| j (cid:105)(cid:104) j | †\nis the ideal state, and we have ρˆ ,ψˆ , m ρˆ , m ψ ψ St( 2 ). In the following, we\nmix pure ⊗j=1 j ⊗j=1| j (cid:105)(cid:104) j | ∈ Hm,E\ndenote (cid:15) =(cid:15) +(cid:15) as the total error.\nmix channel input\nThen, we have a lemma concerning all errors.\nLemma 11 Let ρˆ be an approximate GE state described by Eq. (9). The estimate errors of\nmix\nGDE proto- col, designed for pure states, obey\n(cid:15) f( ξ ξ˜ , V V 1/4)+(cid:15) (10)\n≤ (cid:107) − (cid:107)\n2\n(cid:107)\npu\n(cid:15)\nre\n− (cid:107)∞\n(cid:48)ps\nξ ξ˜ 8 E(cid:15) + v , (11)\n(cid:107) − (cid:107) 2 ≤ mix 20√mE (cid:101)\n(cid:112)\nV V 40√2E√(cid:15) +(cid:15) (12)\nmix v\n(cid:107) − (cid:107)∞ ≤\n(cid:15) (cid:15) +f( ξ ξ˜ , V V 1/4)+(cid:15) (13)\n(cid:101)\n(cid:48)ps\n≤\nchannel\n(cid:107) − (cid:107)\n2\n(cid:107)\npure\n− (cid:107)∞\nps\nwhere (cid:15) denotes the reconstruction error in tomography, ξ and V p(cid:101)ure denote the displacement and\ncovariance matrix of the state ψ , f( ξ ξ˜ , V V 1/4)= ( ξ ξ˜ )+ ( V V 1/4)\npure 2 pure 2 pure\n(cid:107) − (cid:107) (cid:107) − (cid:107)∞ O (cid:107) − (cid:107) O (cid:107) − (cid:107)∞\nis the error of estimating moments, the first terms in Eqs. (11)–(13) constitute a fixed error floor,\nthe second terms of Eqs. (11) and (12) give the errors(cid:101)for an ideal pure input, (cid:15)\n(cid:48)ps\nis the(cid:101)error\nof estimating the approximated passive separable state after counter-rotation based on estimating\nmoments, (cid:15) denotes the corresponding error where we have an ideal pure GE state.\nps\nWenotethatapureGEstatecorruptedbyapure-losschannelisthespecialcasewithzeroinputerror. Thus,\nwe added Theorem 3 in the revised manuscript, appended below.\nTheorem 3 Consider approximated pure GE states ρˆ St( 2 ) whose entangling channel is\n∈ Hm,E\n(cid:15) -close to an arbitrary Gaussian unitary. The sample complexity of learning this class of\nchannel\nstate is:\n1 1\nM poly m,E, ,log , (14)\n(cid:39)O (cid:15) (cid:15) δ\n(cid:18) (cid:18) − channel (cid:18) (cid:19)(cid:19)(cid:19)\nfor an error (cid:15)>(cid:15) in trace distance and failure probability δ.\nchannel\n• Second, we observe that the output states of Boson sampling under uniform loss equals a passive Gaus-\nsian unitary applied to a product of mixed states. Motivated by this, Supplementary Note VIII presents a\nlearning protocol for this class, which adapts the reduced-state witness to mixed inputs. For scenarios with\nnondegenerate symplectic eigenvalues, we add Corollary 4 in the main text, appended below for convenience:\nCorollary 4 (Sample complexity of states entangled by a Gaussian unitary) Considerthe\nsetofm-modestates ρˆ obtainedbyapplyingaGaussianunitarytoatensorproductofsingle-mode\nstates m ρˆ , subje { ct } to an energy-moment constraint ρˆ St( 2 ). Meanwhile, the states’\nj=1 j ∈ Hm,E\nsymplectic eigenvalues are nondegenerate, up to a finite eigen-gap ∆. The GDE protocol needs\n(cid:78)\nM poly m,E,∆ 1,(cid:15) 1,log δ 1 copies for this state family, for an error (cid:15) in trace dis-\n− − −\n(cid:39) O\ntance and failure probability δ.\n(cid:0) (cid:0) (cid:0) (cid:1)(cid:1)(cid:1)\nForthedegeneratesymplectic-eigenvaluecase, weaddacorrespondingdiscussioninthenewlyaddedMethod\nsection “Learning degenerate mixed GE states”.\nFinally, we added a paragraph clarifying the relation between photon-loss–corrupted states and our previous\ndiscussion, appended below.\nBelow, we apply Corollary 4 and Theorem 3 to practical experimental scenarios. In this regard,\ndetectorinefficiencyisonetypeubiquitousandpracticallyrelevantimperfectionthatleadstomixed\nNGE states. Here, it can be modeled as ideal detection preceded by pure-loss channels. So uniform\nefficiencymapsboson-samplingoutputstoapassiveGaussianunitaryactingonaproductofmixed\nstates. Consequently,Corollary4appliesinthesituationwithnondegeneratesymplecticeigenvalues,\nwhere both the sample and time complexity stay polynomial in the number of modes. For non-\nuniform loss, the degraded channel remains close to a Gaussian unitary in the small-loss regime.\nThus, Lemma 11 in Methods guarantees a vanishing reconstruction error floor, and Theorem 3\nshows that the sample complexity increases only polynomially with the degradation level.\n• Moreover, to show case the robustness of our approach to loss. We have added simulations for various exotic\nquantum states subject to uniform loss, as we depict in Fig. 4 of the main text. The figure is also appended\nas Fig. R2 of the reply.\nEntangled number state Entangled cat state Correlated thermal state\nFigure R2: Performance of Gaussian-disentangling protocol (GDE) versus direct tomography (RρR). Here, we\nnumerically simulate the quantum state tomography process on a two-mode GE state with the true state being\nnumber state 0 1 , cat state cat cat and thermal state ρ ρ correlated by a beam-\n| (cid:105)⊗| (cid:105) | (cid:105)0.1 ⊗| (cid:105)1.1 th,n=0.2 ⊗ th,n=0.3\nsplitter (with a rotation angle θ = π/4). The error, defined as the trace distance between the true state and\nthe reconstructed state obtained through direct tomography using the standard RρR algorithm [3], is shown as a\nfunction of the number of copies M used in tomography (black squares for ideal case and grey stars for the lossy\ncase). In the lossy case, pure loss channel (with a transmissivity ratio η = 0.9) is applied to both modes after the\nbeamsplitter. Incontrast, thebottomenvelopeofachievableerrorsintheproposedGaussiandisentangling(GDE)\nalgorithm, after minimising over possible sample number in estimating displacements and the covariance matrix,\nis\ngiven by red circles and pink triangles (lossy case). Their linear fit shows errors in GDE protocol decrease faster\n(with slope k) compared to RρR method in all the three different states, demonstrating the advantage of GDE\nmethod.\nList of Changes\nWe attached the marked version with red highlighting all changes. Below is a succinct summary of the changes.\n• a sentence revised in abstract\n• bottom of page 1 introduction revised\n• Table I updated\n• page 3 sentences revised\n• page 4 sentences revised\n• Fig. 2 added, with the paragraph below revised.\n• A new section “Towards efficient learning of mixed GE states” is added, with Fig. 4 revised and moved here.\n• Two sections in Method: “A useful lemmas on mixed GE-state learning” and “Learning degenerate mixed\nGE states” are added, while the section “Additional example of NGE state” is removed.\n• sentence revised in page 15 of the Supplemental Notes\n• Supplemental Notes VIIF and VIII added\nReferences\n[1] CostinB˘adescuandRyanO’Donnell. Improvedquantumdataanalysis. InProceedings of the 53rd Annual ACM\nSIGACT Symposium on Theory of Computing, pages 1398–1411, 2021.\n[2] Simon Becker, Nilanjana Datta, Ludovico Lami, and Cambyse Rouz´e. Classical shadow tomography for contin-\nuous variables quantum systems. IEEE Transactions on Information Theory, 2024.\n[3] Alexander I Lvovsky and Michael G Raymer. Continuous-variable optical quantum-state tomography. Reviews\nof modern physics, 81(1):299–332, 2009.\n[4] Pablo V Parellada, Vicent Gimeno i Garcia, Julio Jos´e Moyano-Fern´andez, and Juan Carlos Garcia-Escartin.\nNo-go theorems forphoton state transformationsin quantum linear optics. Results in Physics, 54:107108, 2023.\n[5] Pablo V Parellada, Julio Jos´e Moyano-Fern´andez, Juan Carlos Garcia-Escartin, et al. Lie algebraic invariants\nin quantum linear optics. arXiv:2409.12223, 2024.\n[6] NM VanMeter, P Lougovski, DB Uskov, K Kieling, J Eisert, and Jonathan P Dowling. General linear-optical\nquantumstategenerationscheme: applicationstomaximallypath-entangledstates.Phys.Rev.As,76(6):063808,\n2007.\n10"}, {"round": 2, "reviewer_comments": {"reviewer_2": "The authors have comprehensively and convincingly addressed all the comments and questions raised in my previous\nreport. The revised manuscript is suitable for publication and I recommend its acceptance for publication in Nature\nCommunications.", "reviewer_3": "The authors have addressed the issues raised adequately. As such, this manuscript is suitable for publication."}, "author_response": ""}]

Decision letter