---
title: 'Langevin Sampler: Scalable Quantum Tomography'
url: https://www.emergentmind.com/topics/langevin-sampler-for-quantum-tomography
type: topic
---

# Langevin Sampler: Scalable Quantum Tomography

A Langevin sampler for quantum tomography is a Bayesian computational approach that leverages the Burer–Monteiro factorization to efficiently estimate quantum states from measurement data, with explicit low-rank structure and scalability guarantees. The method operates by parameterizing a \(d \times d\) Hermitian positive semidefinite (PSD) density matrix via a product of a complex matrix and its conjugate transpose, imposing physical constraints directly in the parameter space. This enables the construction of a posterior distribution restricted to matrices of known or bounded rank and, through the use of a spectral Student–$t$ prior, promotes solutions of even lower rank when the true rank is unknown. The posterior is explored via an unadjusted Langevin algorithm (ULA), with rigorous PAC–Bayesian risk bounds that match state-of-the-art rates, and the algorithm achieves substantial computational savings compared to conventional Markov chain Monte Carlo (MCMC) techniques when the target density is low rank [2601.08775].

## 1. Parameterization via Burer–Monteiro Factorization

Quantum tomographic inference seeks a density matrix $\rho \in \C^{d \times d}$ satisfying $\rho = \rho^\dagger$, $\rho \succeq 0$, and $\tr(\rho) = 1$. If $\rho$ is known (or assumed) to have rank $r \leq d$, it is parameterized as $\rho = Z Z^\dagger$ with $Z \in \C^{d \times r}$. The unit-trace condition translates to $\|Z\|_F^2 = 1$, so $Z$ resides on the complex hypersphere.

Measurement data is modeled by $A$ observables, each with $S$ possible outcomes. The empirical frequencies $\hat p_{a,s}$ for observable $a$ and outcome $s$ are related to the Born prediction $p_{a,s}(Z) = \tr(\mathcal P_s^a Z Z^\dagger)$, where $\mathcal P_s^a$ denote the POVM elements. A pseudo-likelihood corresponds to a sum-of-squares loss:
$$
L(Z) = \sum_{a=1}^A \sum_{s=1}^S [\hat p_{a,s} - \tr(\mathcal P_s^a Z Z^\dagger)]^2
$$
Bayesian inference proceeds with the posterior $p(Z | \text{data}) \propto \exp(-\lambda L(Z)) p(Z)$, leading to potential $U(Z) = \lambda L(Z) - \log p(Z)$.

This factorization reduces the ambient parameter space from $d^2$ to $d r$ dimensions and automatically maintains Hermitian PSD structure.

## 2. Low-Rank–Promoting Spectral Prior

When the rank of $\rho$ is unknown but an upper bound $r$ is available, a spectral Student–$t$ prior is used to promote low-rank solutions:
$$
p(Z) \propto \det(\theta^2 I_d + Z Z^\dagger)^{-\alpha},\quad \alpha = \tfrac{2d + r + 2}{2},\;\; \theta > 0
$$
The prior decomposes as a product over the singular values $\sigma_i$ of $Z$: $p(Z) \propto \prod_{i=1}^r (\theta^2 + \sigma_i^2)^{-\alpha}$, heavily penalizing small $\sigma_i$ and thus favoring low-effective-rank $\rho$. The gradient of $\log p(Z)$ is available in closed form, facilitating efficient implementation:
$$
\nabla_Z \log p(Z) = - (2d + r + 2)\, (\theta^2 I_d + Z Z^\dagger)^{-1} Z
$$

This prior is a complex generalization of that studied by Dalalyan (2020) for promoting low-rank matrix estimation.

## 3. Langevin Sampler: Stochastic Dynamics and Discretization

The posterior on $Z$ is sampled via complex Langevin dynamics governed by the SDE:
$$
\mathrm{d}Z_t = - \nabla_Z U(Z_t) \, \mathrm{d}t + \sqrt{2}\, \mathrm{d}W_t
$$
where $W_t$ is Brownian motion in $\C^{d \times r}$.

Discretization through the unadjusted Langevin algorithm (ULA) with step size $\eta$ yields
$$
Z_{k+1} = Z_k - \eta \nabla_Z U(Z_k) + \sqrt{2\eta}\, \Xi_k,\quad \Xi_k \sim \mathcal N(0, I_{dr})
$$
To maintain the trace constraint $\|Z\|_F = 1$, two strategies are used: (1) projected Langevin, normalizing $Z_{k+1}$ after every step, or (2) unconstrained iteration with trace normalization applied only to the final estimator. Empirically, per-step drift in $\|Z\|_F$ is negligible, so the latter is often preferable for computational simplicity.

## 4. Implementation Workflow

The following high-level pseudocode summarizes the Langevin sampling scheme for quantum tomography:

```python
Input: r, η, λ, θ, N, B, initial Z₀ ∈ ℂ^{d×r}, ‖Z₀‖_F=1
for k = 0 to N−1:
    G = ∇_Z [ λ L(Z_k) − log p(Z_k) ]
    Ξ ∼ 𝒩(0, I_{dr})     # i.i.d. complex Gaussian noise
    Z_{k+1} = Z_k − η G + √(2η) Ξ
    # Optionally: Z_{k+1} ← Z_{k+1} / ‖Z_{k+1}‖_F
end for
ρ̄ = (1/(N−B)) ∑_{k=B}^{N−1} Z_k Z_k^†
return ρ̂ = ρ̄ / tr(ρ̄)
```

Key hyperparameter choices:
- λ controls the data-prior trade-off; for complete measurements, λ ≈ m/2 or 3m/8 is recommended.
- θ tunes the rank penalty; small θ strongly penalizes rank, large θ recovers a nearly uniform prior.
- η is chosen empirically for stability, typically $10^{-5}$–$10^{-4}$.
- N (number of iterations) and B (burn-in) are set to ensure convergence and posterior mixing.

## 5. PAC–Bayesian Risk Guarantees

In the full Pauli measurement setting with $N_{\text{tot}} = 3^n m$ samples and true rank-$r$ density $\rho^0 = Y^0 Y^{0†}$ with $Y^0 \in \C^{d \times r}$, a PAC–Bayesian Frobenius risk bound holds. For any $\epsilon \in (0,1)$ and comparison factor $\bar Y$ of rank $p$ such that
$$
\|\bar Y - Y^0\|_F\,\|\bar Y + Y^0\|_F \leq 3^{-3n/2}2^{-n/2} m
$$
it holds with probability at least $1-\epsilon$ that
\begin{align*}
\|\hat\rho-\rho^0\|_F^2
&\leq
\frac{3}{N_{\text{tot}}\bigl(3^{3n/4}2^{(n+6)/4}(r+\sqrt{r}\,\|\bar Y\|_F)+2r/m+1\bigr)} \\
&\qquad +
\frac{8\cdot3^n}{2^n N_{\text{tot}}}
\frac{\log(2/\epsilon)+2p(2^{n+1}+r+2) \log(1+\|\bar Y\|_2/\theta)}{}
\end{align*}
The analysis leverages exponential moment inequalities for the empirical squared error, KL-divergence bounds for shifted priors, and spectral properties of measurement operators. The leading term (up to log factors) matches the minimax rate $\frac{3^n\,\rank(\rho^0)}{N_{\text{tot}}}$ known from the literature (Mai & Alquier 2017, Mai 2021).

## 6. Computational Complexity and Empirical Performance

Each Langevin iteration involves:
- **Data drift term**: $\nabla L(Z)$ costs $O(A d r)$ operations with $A$ measurements (for complete Pauli, $A = 3^n$).
- **Prior gradient**: solving $(\theta^2 I_d + Z Z^\dagger)^{-1} Z$ is reduced via Sherman–Morrison–Woodbury to $O(d r^2 + r^3)$ for $r \ll d$.
- **Noise sampling**: $O(d r)$.

Total per-step cost: $O(A d r + d r^2 + r^3)$, dominated by $O(d r^2)$ for low-rank setting and moderate $A$.

Overall, runtime scales with the number of iterations times this per-step cost. Empirical benchmarks indicate:
- **Scalability**: With small $r$, each Langevin update is vastly faster than in full-rank MCMC schemes.
- **Mixing**: The chain mixes in fewer than $10^4$ steps, substantially fewer than required by Metropolis–Hastings methods for comparable accuracy.
- **Estimation accuracy**: The final estimator achieves Frobenius norm errors competitive with, or superior to, existing Bayesian quantum tomography algorithms.

## 7. Extensions and Practice Considerations

Potential refinements include step-size (η) annealing or adaptive temperature control for improved mixing. Metropolis-adjusted Langevin (MALA) or Riemannian variants may strengthen theoretical convergence guarantees. The method applies directly to process tomography (Choi matrix estimation) and can address incomplete measurement regimes. The reduction from $d^2$ to $d r$ dimensional parameter space, inherent PSD constraint enforcement, and low-rank–favoring prior make the Langevin sampler a scalable and theoretically sound approach for large-scale quantum tomography with explicit risk guarantees [2601.08775].

Source: https://www.emergentmind.com/topics/langevin-sampler-for-quantum-tomography