Quantum speedup for nonreversible Markov chains

Domain: physics
Published: 2025-11-28
Editorial outcome: major_revision

Abstract
Quantum algorithms can potentially solve a handful of problems more efficiently than their classical counterparts. In that context, it has been discussed that Markov chains problems could be solved significantly faster using quantum computing. Indeed, previous work suggests that quantum computers could accelerate sampling from the stationary distribution of reversible Markov chains. However, in practice, certain physical processes of interest are nonreversible in the probabilistic sense and reversible Markov chains can sometimes be replaced by more efficient nonreversible chains targeting the same stationary distribution. This study constructs Markov chain reversibilizations and develops quantum algorithmic techniques to accelerate nonreversible processes. Such an up-to-exponential quantum speedup goes beyond the predicted quadratic quantum acceleration for reversible chains and is likely to have a decisive impact on many applications ranging from statistics and machine learning to computational modeling in physics, chemistry, biology and finance.

Introduction
Quantum algorithms recently gained attention due to their asymptotic advantage over the best known classical methods for certain problems1,2. Such speedups turned quantum computing into an active field of research. Monte Carlo simulations both offer a wide range of applications and promise to be accelerated by quantum algorithms3.

The speedup is typically provided by the Quantum Amplitude Estimation (QAE) routine4,5. The routine requires access to samples from a target probability distribution π. When dealing with complex, high-dimensional distributions known only up to a normalization constant (as in statistical physics), Markov Chain Monte Carlo (MCMC) is among the most popular methods to provide such classical samples6. If not used solely, it is often part of more elaborate algorithms such as Sequential Monte Carlo7 or Annealed Importance Sampling8. The quantum coherent state |π⟩ associated with the distribution may be prepared using the Quantum Rejection Sampling (QRS) algorithm9, but QRS requires bounds on the normalization constants for efficiency. The present introductory paragraph uses common definitions, which are recalled in the “Methods” section “Markov chains”. In MCMC, we start by designing an ergodic Markov chain that converges to the sought distribution π10. Performing sufficiently many steps of the chain yields an approximate sample from π. The main difficulty is that the mixing time, the minimum number of steps for the sample to be distributed according to π, may be large. It is well-known that, to each reversible chain, we can associate a quantum walk with quadratically larger spectral gap11. It is then possible to construct the reflection through the stationary distribution in the inverse spectral gap of this quantum walk (see, for example12,13). Since the mixing time is of the order of the inverse spectral gap of the chain, the reflection can readily be used to provide a sample from π with quadratic speedup14.

First and foremost, we are often interested in studying physical processes that are nonreversible in the probabilistic sense. Such processes may have a known stationary distribution, as for the underdamped Langevin dynamics, or not, as in the study of out-of-equilibrium systems in statistical physics15,16. In addition, reversible chains are known to show a diffusive behavior17. It is sometimes possible to construct a nonreversible chain that converges faster than the reversible counterpart to the same stationary distribution18,19. One way to do this is to use a so-called lifting procedure. An optimal lift can offer up to a quadratic speedup (like the quantum algorithms)20,21 but typically requires a detailed knowledge of the chain20. Because there is no downside to using lifts, nonreversible processes are generally preferred for sampling purposes22. Nonreversible processes may also be used to encode reversible ones, such as Metropolis-Hastings kernels, more efficiently23. Because of the growing interest in such processes, we ask the following question. Can quantum algorithms accelerate the mixing of nonreversible Markov processes?

In summary, we answer this question by:
• analyzing known quantum singular value transforms in the context of nonreversible Markov chains (which requires simulating the time-reversal), and
• introducing quantum eigenvalue transforms to retrieve an approximation of the stationary distribution without simulating the time-reversal.

More precisely, this work provides two ways to construct an approximate reflection through the target coherent state from the Szegedy quantum walk operators11 (section “Encoding Markov chains in quantum computers” provides further details on the oracle access model). Its main contributions are the algorithmic workflows presented in Fig. 1. In both methods, we perform classical steps of the walk until a certain criterion is met and the chain is “reversible enough.” Then, we perform a suitable quantum singular value or eigenvalue transform, producing the desired reflection with a speedup. Amplitude amplification algorithms use such a reflection to prepare the sought state24. Since QAE relies on such a state preparation unitary4,5, the reflection is readily usable in Quantum Monte Carlo (QMC) routines3. QMC requires quadratically fewer samples from the stationary distribution than classical MCMC, and each sample can be provided more efficiently.

In the first method, the number of steps to be performed classically is called the reversibilization time τrev and has a precise definition. The algorithm constructs the desired reflection in a time of the order of √τrevτ, where τ is the mixing time of the chain (see the section “Singular value transformations of non-normal Markov kernels”). The second method builds on the same intuition to define the number of steps to be performed classically before the quantum eigenvalue transformation is applied (see Fig. 1). We introduce a Markov kernel called the geometric reversibilization of the chain and a condition of reversibility on π-average. This condition guarantees both the efficiency of the method and the quality of the produced reflection. Indeed, if this condition is met, the geometric reversibilization reaches its stationary distribution in a time of the order of the mixing time of the additive reversibilization of the chain, the process that goes forward and backward in time with equal probabilities. This stationary distribution is called the most reversible distribution and is guaranteed to be close to π. Figure 2 illustrates an example where the condition of reversibility on π-average is verified long before the mixing time. Figure 2a displays that the eigenvalues of the geometric reversibilization, whose properties govern the runtime of the quantum algorithm, rapidly approach those of the additive reversibilization. Figure 2b shows that this phenomenon occurs long before the mixing time of the kernel, resulting in a runtime of the square root of the mixing time of the additive reversibilization (see the section “Reversibility on π-average and another notion of reversibilization time”).

By efficiently constructing reflections through the stationary distribution of nonreversible chains, this study expands the set of Markov chain computations that can be achieved faster using a quantum computer. These reflections could ultimately be applied to accelerate molecular dynamics, with implications in drug design25,26, in protein folding studies27–30 and in any subfield using molecular simulations. They could also help simulating the limiting behavior of stochastic differential equations, with implications in financial modeling31,32.

Results
Technical background
In this section, we first summarize Markov chain concepts17. We then describe core quantum algorithmic tools such as projected unitary encodings and qubitized walk operators, more detailed expositions of which can be found in refs. 33,34. We observe that, in order to get a reflection through the stationary distribution, it is sufficient to obtain a projected unitary encoding of the projector on this state. We express the latter as polynomials of discriminant matrices, which are related to Markov kernels and can be encoded through projected unitary encodings. Finally, we recall how the Generalized Quantum Eigenvalue Transform (GQET) and the Generalized Quantum Singular Value Transform (GQSVT) allow us to construct a projected unitary encoding of such polynomials.

Markov chains and mixing time. A Markov kernel on a finite set S is a matrix P of positive numbers such that each row sums to 1. The Markov chain with initial condition X0 ∈ S associated with P is a stochastic process (Xt)t∈N with values in S such that when it is in state x ∈ S, it goes to y ∈ S with probability P(x, y). If the kernel allows reaching any state y ∈ S from any state x ∈ S, we say that it is irreducible. If it has a single eigenvalue on the unit circle, we say that it is aperiodic. If P is both irreducible and aperiodic, it is said to be ergodic. An ergodic Markov kernel is such that, after waiting for a sufficiently long time τ ∈ N, the probability distribution of the system state Xτ is approximately given by a probability distribution π that is independent of the starting state (see the Methods “Markov chains” for a precise definition of the mixing time τ). This distribution is called the stationary distribution. A sufficient condition for an ergodic kernel P to have a stationary distribution π is to be reversible with respect to π: such that π(x)P(x, y) = π(y)P(y, x) for any states x, y ∈ S. Defining the time-reversal P⋆ by P⋆(x, y) = π(y)P(y, x)/π(x) for every x, y ∈ S, the previous condition can be rewritten P = P⋆. We will refer to (P + P⋆)/2 as the additive reversibilization of P, and to PP⋆ as the multiplicative reversibilization of P. MCMC algorithms compute the expectation values of functions of random variables sampled according to π. Therefore, the number of Markov chain steps required to provide a sample from π is a natural performance metric.

Mixing times and notion of spectral gap. The mixing time of ergodic reversible Markov kernels is related to their spectral gaps. If P is an ergodic reversible Markov kernel, define its spectral gap to be γ(P) = 1 − maxλ∈σ(P)n{1} |λ|. Then, the mixing time of P is essentially of the order of 1/γ(P)17. As it turns out, this relationship does not hold true for general Markov kernels. If P is an ergodic kernel, we may define its pseudo-spectral gap γ∞(P) by γ∞(P) = maxk≥1 γ(Pk(P⋆)k)1/35. Then, the previous relationship can be generalized in the sense that the mixing time of P is of the order of 1/γ∞(P).

Projected unitary encodings and qubitized walk operators. A unitary U and two (partial) isometries ΠL, ΠR are said to be a Projected Unitary Encoding (PUE) (U, ΠL, ΠR) of A if ΠL† U ΠR = A (the superscript † indicates the conjugate transpose matrix). If ΠL = ΠR = Π and U is also symmetric, then (U, Π) is said to be a Symmetric Unitary Projected Encoding (SPUE) of A. Given (U, Π) a SPUE of A, define the qubitized walk operator W = (2ΠΠ† − 1)U. This operator is going to be particularly useful when A is the projector on a state of interest. Indeed, if A = |ϕ⟩⟨ϕ| for some state |ϕ⟩, then W2, Π is a SPUE of 2|ϕ⟩⟨ϕ| − 1 (see equation (12) in the “Methods”).

Encoding Markov chains in quantum computers. Let us now define PUEs, namely unitary operators and partial isometries, that encode information about Markov kernels. Consider the Hilbert space spanned by the set S. For every probability distribution η on S, we write |η⟩ = ∑x∈S η(x)|x⟩. If P1, P2 : S2 → [0, 1] are two Markov kernels, then for each x, y ∈ S:
D(x, y) = √P1(x, y)P2(y, x) : (1)
In the case P1 = P, P2 = P⋆, we will refer to the resulting matrix D as the curved discriminant D = DP,P⋆, with PUE (S, Π⋆, Π) = (S, ΠP⋆, ΠP). We will be interested in cases where one of the eigenvectors (resp. singular vector) of D (resp. D) is a good approximation |μ⟩ for |π⟩. Then, we may write our target projector |μ⟩⟨μ| = υ(D), where υ is a polynomial such that υ(μ) = 1 and υ(λ) = 0 for other eigenvalues λ ≠ μ of D. Constructing encodings of such υ(D) is precisely the aim of the GQET and GQSVT presented below. It will make use of the Szegedy quantum walk operators11, R = 2ΠΠ† − 1 and R⋆ = 2Π⋆(Π⋆)† − 1. The operators R and R⋆ are typically constructed with arithmetic operators capable of performing the transformation
|x⟩ → ∑y∈S √P(x, y)|x, y⟩, (2)
for every state x ∈ S (see, for example, ref. 36). Such operators are efficiently implementable whenever the transition probabilities are efficiently computable. If the stationary distribution π can be computed up to a multiplicative constant, ratios of the form π(y)/π(x) can be efficiently evaluated for any pair of states x, y ∈ S. Consequently, the time-reversed transition probabilities P⋆(x, y) = P(y, x)π(y)/π(x) can also be computed efficiently, implying that P⋆ is accessible under the same oracle model. As demonstrated in refs. 37–39 and exemplified in the Supplementary Information, it is often possible to construct these operators much more efficiently. For kinetic processes, where the state space includes both position and velocity components (as in underdamped Langevin dynamics), P⋆ can be efficiently implemented by reversing the initial velocity, applying one step of P, and then reversing the velocity again.

Generalized Quantum Eigenvalue Transform. Let (U, Π) be a SPUE of an operator H. Let υ be a polynomial of degree d with complex coefficients. The GQET is a unitary U built from d uses of controlled-W operations such that (U, |0⟩ ⊗ Π, |0⟩ ⊗ Π) is a PUE of υ(H)/β. W denotes the qubitized walk operator of (U, Π) and 1 ≤ β ≥ O(log(d)) is called the scaling factor of υ34.

Generalized Quantum Singular Value Transform. Let (U, ΠL, ΠR) be a PUE of an operator A and υ be an even polynomial of degree d with complex coefficients. The GQSVT is a unitary U built from d controlled-U, controlled-U†, controlled-(2ΠL ΠL† − 1) and controlled-(2ΠR ΠR† − 1) such that (U, |01⟩ ⊗ ΠR, |01⟩ ⊗ ΠR) a PUE of V†υ(Σ)V/β, where A = W†ΣV is the singular value decomposition of A and β is the same scaling factor as in GQET.

Applying the fast-forwarding polynomial. Let (S, Π⋆, Π) be a PUE of D, the curved discriminant of an ergodic Markov kernel. Without loss of generality, we assume that the singular values of D are in [0, 1 − δ] ∪ {1}, for some δ ∈ ]0, 1[, with unique left and right singular vectors with singular value 1 denoted by |π⟩. We want to apply the GQSVT with a polynomial υ of low degree such that V†υ(Σ)V/β is an approximation of |π⟩⟨π|, where D = WΣV is the singular value decomposition of D. As detailed in the Methods “Quadratic speedup polynomial”, for any ϵ > 0, there exists a polynomial υ of degree O(δ−1/2 log(1/ϵ)) such that υ(1) = 1 and |υ(x)| ≤ ϵ for x ∈ [0, 1 − δ]. An example of such a polynomial of degree 8 is illustrated in Fig. 1, with ϵ = 0.2. The polynomial is more efficient than the monomial x8, corresponding to the execution of 8 steps of the Markov chain. In the rest of the section, we will refer to υ as the fast-forwarding polynomial.

Singular value transformations of non-normal Markov kernels
The GQSVT algorithm allows to apply fixed-parity polynomials to the singular values of P, the square roots of the eigenvalues of PP⋆. Therefore, the performances of the algorithm strongly depend on the mixing time of the multiplicative reversibilization PP⋆ of P. However, up to considering the lazy version of P, PP⋆ may mix quadratically slower than P18. The quantum speedup provided by the GQSVT transformation being quadratic, there are no benefits to directly transforming the singular values of P. In order to optimize the procedure, we suggest to apply the GQSVT algorithm to Pk, for a number of steps k left unspecified, and at the expense of multiplying the cost of the method by k. Is there an integer k such that the resulting construction requires much fewer steps than the mixing time?

As described in the “Technical background” section “Markov chains and mixing time”, the mixing time is related to the inverse of the pseudo-spectral gap. By definition of the pseudo-spectral gap, there exists a smallest integer τrev ≥ 1 such that γ∞(P) = γ(Pτrev(P⋆)τrev) = τrev. Applying the GQSVT transformation to Pτrev to encode an approximation of |π⟩⟨π| with O(1/γ(Pτrev(P⋆)τrev)) uses the PUE of Pτrev. Thus, the overall procedure has complexity:
O(√τrev/γ(Pτrev(P⋆)τrev)) = O(1/γ∞(P)) : (3)
Recalling that 1/γ∞(P) is upper bounded by the mixing time, the overall complexity is in:
O(√τrevτ) : (4)
Proposition 1 summarizes the derivation.

Proposition 1. Let P be an ergodic kernel on finite state space S. A quantum circuit approximating 2|π⟩⟨π| − 1 up to spectral norm error ϵ > 0 can be constructed with O(√τrevτ(ϵ) log(1/ϵ)) uses of the Szegedy quantum walk operators R and R⋆, where τrev is the reversibilization time of P, and τ(ϵ) is its ϵ-mixing time.

To summarize, we may reflect through the stationary measure of a Markov process with speedup whenever the reversibilization time of the process is smaller than its mixing time.

Eigenvalue transformations of the flat discriminant
In practice, we may not be able to access R⋆. Indeed, the time-reversal of the chain may not be easy to sample from37 or the stationary distribution may not be known16. Let us discuss how to implement 2|π⟩⟨π| − 1 using only R. Recall that (S, Π) is a SPUE of the flat discriminant D, a symmetric matrix with eigenvalues in [−1, 1].

Let us now construct a new reversibilization of a Markov kernel. If D is primitive, meaning that Dm > 0 element-wise for some m ∈ N, then the Perron-Frobenius Theorem (Theorem 6 in the Supplementary Information) and the variational principle imply the existence of a strictly positive probability distribution μ on S such that:
μ = argmaxν ⟨ν|D|ν⟩, (5)
where the argmax ranges over all strictly positive probability distributions on S. The geometric reversibilization of P, denoted by Q : S2 → [0, 1], and defined by:
Q(x, y) = √μ(y)/μ(x) D(x, y) / ⟨μ|D|μ⟩ (6)
is well-defined. Moreover, Q is a Markov kernel that is reversible with respect to μ. D being primitive, it is in fact ergodic with a unique stationary distribution μ. Also, the spectra of D and Q are related through the equation σ(Q) = σ(D)/⟨μ|D|μ⟩. According to the variational principle, ⟨π|D|π⟩ is a lower bound for the largest eigenvalue of D. Moreover, if ⟨π|D|π⟩ is much closer to 1 than to 1 − γ(Q), where γ(Q) is the spectral gap of Q, then the overlap ⟨π|μ⟩ between the stationary distribution and the most reversible distribution is large. We refer to this condition as reversibility on π-average. If a kernel is reversible on π-average, 1 − ⟨π|D|π⟩ ≪ γ(Q). In order to implement the reflection 2|μ⟩⟨μ| − 1 up to error ϵ using the GQET formalism, we need a polynomial υ of low degree such that υ(⟨μ|D|μ⟩) ≥ 1 − ϵ and |υ(x)| ≤ ϵ for all x ∈ [−1, ⟨μ|D|μ⟩(1 − γ(Q))]. What is the minimal degree for a polynomial satisfying this property?

Let us construct such a polynomial as the composition of υ, which was applied to the singular values of the curved discriminant, with a function f : [−1, 1] → R such that f(x) = 1 for all x ∈ [⟨μ|D|μ⟩, 1] and f(x) = 0 for all x ∈ [−1, ⟨μ|D|μ⟩[. Applying f ∘ υ to the eigenvalues of the flat discriminant would construct a PUE of |μ⟩⟨μ| from which we can recover 2|μ⟩⟨μ| − 1. f is not a polynomial but can be efficiently approximated by a polynomial, called a leading-eigenvalue-selection polynomial. The resulting composite polynomial is illustrated in Fig. 1. It is of low degree if the gap between the two leading eigenvalues of D does not close more rapidly than the gap between the first eigenvalue and 1, as stated by Proposition 2. A proof is given in the “Methods” section “Quadratic speedup polynomial”.

Proposition 2. Let 0 < δ < 1 and ϵ > 0. Then, there exists a constant c ∈ ]0, 1[ and a real polynomial υ of degree O(δ−1/2 log(1/ϵ)) such that:
• |υ(x)| ≤ 1 for all x ∈ [−1, 1],
• |υ(x)| ≤ ϵ for all x ∈ [−1, 1 − δ] and
• υ(x) ≥ 1 − ϵ for all x ∈ [1 − cδ, 1].

If Q has a mixing time that is small compared to 1/(1 − ⟨π|D|π⟩), Proposition 2 provides a polynomial q such that q(D) is a projection on |μ⟩ and ensures a large overlap ⟨π|μ⟩. The result is summarized in Proposition 3 (see the Methods “Quadratic speedup polynomial” for a proof).

Proposition 3. Let P be an ergodic kernel on finite state space S with primitive flat discriminant, most reversible distribution μ, and geometric reversibilization Q with spectral gap γ(Q). Let c ∈ ]0, 1[ be the constant introduced in Proposition 2. If
1 − c(1 − γ(Q)) ≤ ⟨μ|D|μ⟩, (7)
then 2|μ⟩⟨μ| − 1 can be prepared up to spectral norm ϵ > 0 with O(γ(Q)−1/2 log(1/ϵ)) uses of the Szegedy quantum walk operator R. Moreover,
⟨π|μ⟩2 ≥ 1 − ⟨μ|D|μ⟩ / ⟨π|D|π⟩ (8)

Reversibility on π-average and another notion of reversibilization time
If P is reversible, then D is equal to P up to a change of basis. If P is nonreversible, then D may contain no information about the stationary distribution. For example, consider a Markov chain on the N-point discrete circle that goes clockwise with probability 1/2 and stays where it is otherwise: D is half the identity matrix. Indeed, for every x ≠ y, P(x, y) > 0 implies P(y, x) = 0 and therefore D(x, y) = √P(x, y)P(y, x) = 0. However, if P(x, y) = π(y) for each x, y ∈ S, meaning that it is perfectly mixed, then D = |π⟩⟨π| and R = 2|π⟩⟨π| − 1 = 2(|π⟩⟨π|) − 1. In consequence, it is much more appropriate to study the sequence (Dj)j∈N of the flat discriminants associated to Pj instead of simply D = D1. Indeed, we found that applying GQSVT to Pτrev instead of P improved the efficiency of the reflection construction. We want to build on this intuition and apply a GQET to the flat discriminant of Pj, for some integer j ≥ 1. In light of Proposition 3, we will choose j to be the least integer such that Pj is reversible on π-average: such that 1 − ⟨π|Dj|π⟩ ≪ γ(Qj). Because it plays a role similar to that of τrev, we will denote this integer j by trev. The remainder of the section will describe the construction of an approximation of the desired reflection and presents a complexity that is smaller than the mixing time of the chain.

Figure 3 illustrates the phenomenon for a nonreversible walk on a graph with bottleneck. We will consider two non-symmetric (therefore nonreversible) walks on two circles connected by a bridge that is only taken with small probability. Let N be odd (N = 31 for the numerical experiment) and consider the state space S = {0, ..., N − 1} ∪ {∂} ∪ {N, ..., 2N − 1}. Define the nonreversible Markov chain P by:
P(x, [x + 1]N) = P(N + x, N + [x + 1]N) = 3/4,
P(x, [x − 1]N) = P(N + x, N + [x − 1]N) = 1/4, (9)
for all x ∈ {1, ..., N − 1}, where [N] denotes modulo N and
P(0, ∂) = P(N, ∂) = 1/N3, P(0, 1) = P(N, N + 1) = 3/4,
P(0, N − 1) = P(N, 2N − 1) = (1 − 1/N3)/4,
P(∂, 0) = P(∂, N) = 1/2: (10)
If the process starts in {0, N − 1} and does not exit this set, then it converges to a quasi-stationary distribution in O(N2) steps. The same holds if the starting set is {N, 2N − 1}. Since it takes of the order of N3 steps to exit the set, the overall mixing time is at least of order N3. Figure 3 shows that the fast convergence to a quasi-stationary distribution is sufficient for D to present the desired spectral properties. In particular, the figure shows that 1/γ(Qj) is of the order of the square root of the mixing time of Pj.

Note that ⟨π|Dj|π⟩ is the π-averaged overlap between the laws Pj(x, ⋅) and (P⋆)j(x, ⋅) as x follows π: ⟨π|Dj|π⟩ = Eπ[⟨Pj(x, ⋅)|(P⋆)j(x, ⋅)⟩]. In practice, a chain with local moves tends to stay for a long time near a local maxima of probability (according to the distribution π), regardless of whether it evolves according to P or P⋆. In such regions, the overlap between the laws of P and of P⋆ tends to 1 very fast (1 minus the overlap is in fact the square of a distance). For x ∈ S in a low probability region, Pj(x, ⋅) and (P⋆)j(x, ⋅) can remain different until j equals the mixing time; e.g., kinetic processes and their time-reversal tend to leave local minima of probabilities in opposite directions. Since such initial conditions are of low probability, their contribution to the mean is of low importance and ⟨π|Dj|π⟩ = Eπ[⟨Pj(x, ⋅)|(P⋆)j(x, ⋅)⟩] approaches 1. Such ideas are made rigorous using the concept of quasi-stationary distributions in Proposition 8. A numerical example is given in the Supplementary Information. Many discrete approximations of stochastic differential equations also typically have quasi-stationary distributions and rapidly increasing ⟨π|Dj|π⟩ parameter. More details regarding such processes are given in the Supplementary Information.

Note that the flat discriminant of a Markov kernel is equal to the flat discriminant of its adjoint. Moreover, a kernel and its adjoint do not have the same mixing time in general. Corollary 2, given in the “Methods” section “Mixed kernels” and proven in the Supplementary Information, show that if either a kernel or its adjoint is sufficiently mixed, then the reflection through the most stationary distribution can be constructed with constant cost. In the Supplementary Information, we give an example of a Markov chain such that the reflection construction time is exponentially smaller than the mixing time.

Let us now assume that the Markov chain possesses a group structure, then ∥Dj − (Dj + (D⋆)j)/2∥ = 1 − ⟨π|Dj|π⟩. This implies that the eigenvalues of the flat discriminant are all within distance 1 − ⟨π|Dj|π⟩ of an eigenvalue of (Pj + (P⋆)j)/2. A precise statement is given in Proposition 7. For many applications, it is therefore sufficient to focus on ⟨π|Dj|π⟩.

Methods
Let us start by recalling central definitions. Let μ, ν : S → [0, 1] be such that ∑x∈S μ(x) = ∑x∈S ν(x) = 1. The total variation distance separating them is defined as dTV(μ, ν) = 1/2 ∑x∈S |μ(x) − ν(x)|. The ϵ-mixing time of an ergodic Markov kernel P with stationary distribution π is defined by τ(ϵ) = min{t ∈ N : maxx∈S dTV(Pt(x, ⋅), π) ≤ ϵ}. The default precision is ϵ = 1/4 so that the notation τ means τ = τ(1/4). A kernel P is said to be lazy if P(x, x) ≥ 1/2 for each x ∈ S. The mixing time of P is related to its pseudo-spectral gap γ∞(P) through the inequalities:
1/γ∞(P) ≤ τ(ϵ) ≤ 1/γ∞(P) ln(2/ϵπ*), (11)
where π* = minx∈S π(x).

PUE of projectors and reflections
Let us explain why it is sufficient to construct the SPUE of |π⟩⟨π|. Assume that (U, Π) is an SPUE of A = |ϕ⟩⟨ϕ|, for a normalized state |ϕ⟩. Recall the definition of the qubitized walk operator W = (2ΠΠ† − 1)U. Compute:
Π† W2 Π = Π† U(2ΠΠ† − 1)U Π
= 2Π† UΠΠ† UΠ − Π† U2 Π (12)
= 2|ϕ⟩⟨ϕ|ϕ⟩⟨ϕ| − Π† Π
= 2|ϕ⟩⟨ϕ| − 1:
Since (U, Π) is a SPUE, Π† Π is the projection on the support of Π.

Quantum signal processing
We can now describe precisely the quantum signal processing tools at our disposal. The central result provides PUEs of polynomials of unitaries, and is summarized in Theorem 1 (see ref. 33).

Theorem 1. Let ϒ be a degree-d ∈ N complex polynomial and U be a unitary. Then, using d controlled-U operations, O(d) single-qubit gates and 1 ancilla qubit, one can construct a PUE of ϒ(U).

Theorem 1 can be used to apply polynomials to non-unitary symmetric matrices. Theorem 2 shows how to construct PUEs of a polynomial of a symmetric matrix from one of its SPUEs.

Theorem 2. Let (UP, Π) be a SPUE of A. Consider a degree-d complex polynomial υ(x) = ∑n=0 d an Tn(x), written in the basis of Chebyshev polynomials. Consider its associated signal processing polynomial ϒ(z) = ∑n=0 d an zn. Assume that maxz∈∂B1(0) |ϒ(z)| ≤ 1. Let G be the GQSP (Theorem 1) transformation of the qubitized walk operator of (U, Π) applying the polynomial ϒ. Then, (G, |0⟩ ⊗ Π) is a PUE of υ(A).

Performances of the previous construction depend on the scaling factor of the target polynomial, as defined in Definition 1.

Definition 1. Let υ be a complex-coefficient polynomial and ϒ its signal processing polynomial, as introduced in Theorem 2. Define the scaling factor β of υ by:
β = maxz∈∂B1(0) |ϒ(z)| / maxx∈[−1, 1] |υ(x)| : (13)
Finally, Proposition 4 defines polynomials of non-symmetric matrices and provides a construction of them using their hermitianizations4.

Proposition 4. Let (U, ΠL, ΠR) be a PUE of A. Define Uherm = |0⟩⟨0| ⊗ U + |1⟩⟨1| ⊗ U† + |0⟩⟨1| ⊗ ΠL + |1⟩⟨0| ⊗ ΠR. In particular, Uherm is a SPUE of A. It can be constructed from a controlled-U gate, its adjoint, and an X gate. Uherm is called the hermitianization of (U, ΠL, ΠR).

Proposition 5. Let (U, ΠL, ΠR) be a PUE of an operator A. Let (Uherm, Π) be the hermitianization of (U, ΠL, ΠR) and consider a degree-d complex polynomial υ with even part υe and odd part υo. Let G be the GQET applying the polynomial υ to (Uherm, Π). Then,
⟨0|Π† G |0⟩Π = W† υe(Σ)V / β
V† υo(Σ)V / β (14)
where A = W†ΣV is the singular value decomposition of A and β is the same scaling factor as in GQET.

Quadratic speedup polynomial
Let us now state formally certain properties of the eigenvalues and singular values of the curved discriminant D introduced in Paragraph B. The proof of Proposition 6 is simple and given in the Supplementary Information.

Proposition 6. Let D be the discriminant of an ergodic Markov kernel P. Then, D and P are similar. In particular, they share the same spectrum. Also, D†D is similar to P⋆P.

The quantum algorithm described in Paragraph B applies a real polynomial υϵ,d(x) = ϵTd(xT1/d(1/ϵ)) to the singular values of D, where ϵ ∈ ]0, 1[ and for y ∈ ]0, 1], Ty : [1, ∞[ → R is x 7→ cosh(y cosh−1(x)). The polynomial is such that maxx∈[−1, 1] |υϵ,d(x)| = 1, maxx∈[−1+δ, 1−δ] |υϵ,d(x)| = ϵ, |υϵ,d(±1)| = 1 and of degree d ∈ O(δ−1/2 log(1/ϵ)). Applying such a polynomial using the GQSVT formalism requires the signal processing polynomial ϒϵ,d to verify the scaling condition ∥ϒϵ,d(eiθ)∥ ≤ 1 for all θ ∈ R. If this condition is not met, υϵ,d must be divided by some constant, implying a loss in the success probability of the algorithm. The optimal constant is called the scaling factor βϵ,d and defined by:
βϵ,d = maxz∈∂B1(0) |ϒϵ,d(z)| / maxx∈[−1, 1] |υϵ,d(x)| : (15)
In the Supplementary Information, we show the following Corollary 1 stating that βϵ,d = 1.

Corollary 1. The polynomial υϵ,d(x) = ϵTd(xT1/d(1/ϵ)) has scaling factor βϵ,d = 1 for any d ≥ 1 and 0 < ϵ ≤ 1. Moreover, maxx∈[−1, 1] |υϵ,d(x)| = 1.

Let us now prove Proposition 2. The proof is based on the following theorem, whose proof is given in the Supplementary Information. Theorem 3 is inspired by a result from ref. 41. A more precise version of Proposition 2 can be found in the Supplementary Information.

Theorem 3. Let 0 < ϵ < 1. Let (δk)k∈N be the sequence of positive real numbers defined for each k ∈ N by:
δk = T1/k(1/ϵ)−1 : (16)
Let C(ϵ, k) be the set of defined parity real polynomials υ having all their roots in ] − 1, 1[ and such that:
• υ(1) = 1,
• ∀x ∈ [−1 + δk, 1 − δk], |υ(x)| ≤ ϵ.
Then, there is no polynomial of degree less than k in C(ϵ, k). Moreover, the polynomial υϵ,k(x) = ϵTk(T1/k(1/ϵ)x) is such that max υ(1 − δk) = 1 − δk and there is no polynomial of degree less than k in C(ϵ, k).

Consider the sequence of polynomials (υ1/4,k)k∈N. For each k ∈ N, let yk = max υ−1(3/4). Let m(k) be the smallest integer m such that yk ≤ 1 − δm. Note that υ1/4,k ∈ C(3/4, m(k) − 1) so that k ≥ m(k) − 1 by Theorem 3. Moreover, the estimate k ≈ 1/√δk (see the Supplementary Information) implies the existence of constants c1, c2 > 0 such that k ≥ c1/√δk and k + 1 ≤ c2/√δk for k large enough. By definition of m(k),
yk ≤ 1 − δm(k)−1 ≈ 1 − c2/k2 ) 1 − yk ≤ c2/k2 ≤ m(k) : (17)
Combining the inequalities:
1 − δk yk ≥ c1/c2 > 0 : (18)
Proposition 2 directly implies Proposition 3.

Proof of Proposition 3.
First, note that the spectra σ(D) and σ(Q) of D and Q are related through the equation σ(D) = ⟨μ|D|μ⟩ σ(Q). Therefore, the second largest eigenvalue of D can be written 1 − δ = ⟨μ|D|μ⟩ (1 − γ(Q)). In order to apply Proposition 2 to the eigenvalues of D, it is sufficient to have δ ≤ ⟨μ|D|μ⟩. Reformulating the inequality yields the first claim. The second claim follows from the application of a more general corollary of the variational principle given, stated in the Supplementary Information.

Mixed kernels
Note that the flat discriminant is common to a Markov kernel and its time-reversal. Indeed, for each x, y ∈ S,
√P(x, y)P(y, x) = √π(y)P⋆(y, x)π(x)P⋆(x, y) / √π(x)π(y) = P⋆(x, y)P⋆(y, x) : (19)
Assume, for example, that the kernel we are interested in is not mixed, but its time-reversal P⋆ is. Then, the reflection through the most reversible distribution can be constructed with a constant number of steps. Moreover, this distribution has a large overlap with π. This property is especially interesting if P⋆ is not explicitly known. Corollary 2 states it more formally.

Corollary 2. Let 0 < ϵ < 1/(2√8) and t = min(τ(ϵ), τ⋆(ϵ)) where τ(ϵ) is the Hellinger mixing time of P and τ⋆(ϵ) is the mixing time of P⋆. It is possible to implement the reflection 2|μ(t)⟩⟨μ(t)| − 1 up to spectral norm ϵ > 0 with O(log(1/ϵ)) uses of the Szegedy quantum walk. Moreover, ⟨π|μ(t)⟩2 ≥ 1 − 3√8ϵ + O(ϵ2).

Proof. Without loss of generality, assume that t = τ(ϵ). Following the reasoning from the proof of Lemma 20 in ref. 43, we get that
∥Πt − Π∥ ≤ 2√ϵ : (20)
As a consequence, writing Π for the perfectly mixed kernel,
∥D − |π⟩⟨π|∥ = ∥Πt S Πt − Π S Π∥
≤ ∥Πt S Πt − Π S Πt∥ + ∥Π S Πt − Π S Π∥
≤ 2∥Πt − Π∥ ≤ 4√ϵ : (21)
By the Bauer-Fike Theorem (see Supplementary Information), the eigenvalues of D are at distance at most 4√ϵ from the eigenvalues of |π⟩⟨π|. In particular, the second largest eigenvalue λ2 ≤ 1 + 4√ϵ. As detailed in the Supplementary Information, we arrive at:
⟨π|μ(t)⟩2 ≥ 1 − 3√8ϵ + O(ϵ2) : (22)

Random walks on groups
The presented quantum algorithms have a complexity that is governed by spectral properties of D. Being a symmetric matrix, those eigenvalues are necessarily in [−1, 1], whereas the nontrivial eigenvalues of an arbitrary Markov kernel can lie anywhere in B1(0), the unit disk. As such, the eigenvalues of a kernel and its flat discriminant may be quite different. This paragraph introduces random walks on groups and shows that in the presence of symmetries, the eigenvalues of the flat discriminant are in fact very close to those of the additive reversibilization of the chain.

Definition 2. Consider a finite group (S, ⋆) and a probability distribution ν on S. Let (Zt)t≥1 be independent and identically distributed random variables with law ν, and X0 ∈ S. Define the process X = (Xt)t≥1 with
Xt = Zt ⋆ Zt−1 ⋆ ... ⋆ Z1 ⋆ X0 : (23)
Then X is a Markov chain on S called the random walk on (S, ⋆) with increment law ν. Its transition kernel is given for all x, y ∈ S by P(x, y) = ν(y ⋆ x−1).

We prove the following result in the Supplementary Information.

Proposition 7. Assume that P is a random walk on a group (S, ⋆). Then, ⟨π|μ⟩ = 1 and ∥D − DA∥ = 1 − ⟨π|D|π⟩, where D is the flat discriminant of P and DA is the curved discriminant of PA = (P + P⋆)/2. In particular,
σ(D) ⊂ B1−⟨π|D|π⟩(λ) : (24)
λ∈σ(PA)
Also, if I ⊂ σ(PA) and
∪λ∈I B1−⟨π|D|π⟩(λ) ∩ ∪λ∈σ(PA)nI B1−⟨π|D|π⟩(λ) = ∅, (25)
then there are exactly |I| eigenvalues of D in ∪λ∈I B1−⟨π|D|π⟩(λ).

Quasi-stationary distributions
In practice, Markov kernels often present a long mixing time because of the presence of regions of the state space that are hard to escape. Within such a region, the probability measure of the process conditioned on the fact that it did not exit may however, converge rapidly to some locally supported distribution. Definition 3 defines such distributions properly and Proposition 8 shows that their existence is favorable for ⟨π|Dj|π⟩ to be close to 1 for j much smaller than the mixing time.

Definition 3. Let (Xt)t∈N be a Markov chain with state space E ∪ ∂ (with ∂ ∩ E = ∅) which is absorbed at ∂ (i.e., P∂(X1 = ∂) = 1). A quasi-stationary distribution is a probability measure ν on E such that:
∀t ∈ N, ∀A ⊂ E : ν(Xt ∈ A | t < τ∂) = ν(A), (26)
where τ∂ = inf{t ≥ 0, Xt ∈ ∂} is the absorption time of the absorption time.

Write the state space as S = ∪mi=1 Ei ∪ ∂ for some m ∈ N and assume the Markov kernel has quasi-stationary distributions in each of the Ei (with stationary distribution νi). E∂ = ∪mi=1 Ei ∪ ∂ is a set with absorption time.

Proposition 8. The following inequality holds for any integer j ≥ 1:
⟨π|Dj|π⟩ ≥ π(E) min1≤i≤m minx∈Ei Px(j < τ∂i)
min1≤i≤m (1 − 2 Eπ|Ei [dTV(Px(Xj = ⋅ | j < τ∂i), νi)] − dTV(π|Ei, νi)) : (27)

In the Supplementary Information, we prove Proposition 8. We also prove an often tighter result. However, the probabilistic interpretation of this other result is not as clear because the total variation distance is replaced by another distance.

Peer review rounds
[{"round": 1, "reviewer_comments": {"reviewer_1": "- What are the noteworthy results?\nThe authors propose two quantum algorithms for sampling the stationary distribution of a Markov chain. When the Markov\nchain is irreversible with respect to the target stationary distribution, the second algorithm can provide an up-to-exponential\nspeedup and does not require knowing the stationary distribution up to a normalization constant. The resulting algorithmic\nacceleration goes beyond the quadratic speedup known for reversible chains but hinges on the reversibility on π-average\ncondition.\n- Will the work be of significance to the field and related fields? How does it compare to the established literature? If the work\nis not original, please provide relevant references.\nYes, the work addresses the efficient sampling of stationary distributions of generic Markovian systems, including those\npossibly breaking the detailed balance, which yields important implications in various domains of physical and information\nsciences. The first result summarized in Proposition 1 appears related to previous results on Grover-type quadratic\nspeedups. It would be interesting to compare and distinguish from the established literature.\n- Does the work support the conclusions and claims, or is additional evidence needed?\nAnalytical and numerical examples were demonstrated to showcase the more-than-quadratic speedup of the sampling\nalgorithm. However, the work does not seem to fully support, on a broader scope, the claim of impact across various fields as\nit remains unclear for which applications the assumptions made in the analysis are valid.\n- Are there any flaws in the data analysis, interpretation and conclusions? - Do these prohibit publication or require revision?\nThe paper reads a bit dense and heavy of mathematical notations. It would be helpful to (i) restructure the paper to highlight\nthe main result, (ii) include an overview figure that schematically illustrates the core algorithm, and (iii) summarize and\nclearly distinguish the novel contributions in the paper.\nIt remains a bit unclear what are some sufficient and necessary conditions for reversibility on π-average, since this is a key\ncondition for the efficiency of the sampling. Related to this, the method is illustrated only for a few contrived toy examples\nand the claim of broader applicability to statistics, machine learning, and computational modeling is not immediately\nobvious.\nMinor typos need to be corrected throughout the main text E.g.,\n∂B1(0) as the complex unit circle in the right column of page 2 is undefined.\n1-|P(z)|2 in the paragraph discussing the complementarity condition on page 4 needs to be revised, as P is now defined as\nthe Markov transition kernel.\n- Is the methodology sound? Does the work meet the expected standards in your field?\nThe mathematical (Markov chain theory) and algorithmic (projected unitary encoding and generalized quantum eigenvalue\ntransformation) tools used for implementing the sampling algorithm are well established and developed.\n- Is there enough detail provided in the methods for the work to be reproduced?\nThe work includes extensive mathematical proofs but does not provide source code. However, it offers comprehensive\ndetails presented mostly in a self-consistent manner, likely enabling the exact reproduction of the algorithm and numerical\nresults.", "reviewer_2": "I co-reviewed this manuscript with one of the reviewers who provided the listed reports. This is part of the Nature\nCommunications initiative to facilitate training in peer review and to provide appropriate recognition for Early Career\nResearchers who co-review manuscripts.", "reviewer_3": "This paper introduces quantum algorithms that provide up-to-exponential speedup for sampling from the stationary\ndistribution of nonreversible Markov chains. It extends existing method to non-reversible chains.\nThe core algorithmic components are the Generalized Quantum Eigenvalue Transform (GQET) and the Generalized\nQuantum Singular Value Transform (GQSVT), which facilitate the construction of reflections through the stationary\ndistribution.\nQuantizing nonreversible Markov chains is an important problem since nonreversible chains can be more efficient and may\nexhibit faster mixing. While I appreciate the effort put into this work, I have several concerns about its technical aspects.\n1) I would like to see a stronger motivation for the second method, which addresses cases where the stationary distribution\nis unknown up to a normalizing constant. From my experience, most papers on Markov chain algorithms in physics, machine\nlearning, and statistics assume that the stationary distribution is known up to a normalizing constant, as is standard in\nstatistical physics, Bayesian learning, and many relevant models. This makes me question why so much effort is dedicated\nto the second case.\n2) The technical section is difficult to follow, and the novelty appears somewhat limited. From my reading, it seems (at least\nfor the first method) that the authors simply apply the GQSVT algorithm, and the results follow directly. Therefore I have three\nsuggestions\n2.1) The authors may consider providing a more comprehensive overview of existing work, particularly on GQSVT and\nprojected unitary encoding, as these appear to be central to their approach. Currently, the technical background section is\ntoo terse. In fact, while I am very familiar with about 70% of the material, the remaining 30% is impossible to follow. (see also\ncomment 3)\n2.2) The authors should clarify how they utilize existing algorithms and whether their approach is merely a direct application\nof these methods.\n2.3) The authors could improve clarity by helping readers \"connect the dots.\" Even after reviewing the original GQSVT\npaper, it is not entirely clear to me why Proposition 1 holds. For instance, I assume the authors apply GQSVT to the matrix \\( D = \\sqrt{PP^\\star} \\), but it is unclear why this corresponds to a reflection over \\( \\pi \\) and how it approximates \\( \\pi \\).\n3) I may be overlooking something, but I am unable to verify some of the claims made in the paper. For example, at the end\nof the second paragraph in Section 2, the statement \"... Indeed, ..., .. a comprehensive exposition\" does not seem to hold in a\nstraightforward manner. If I take the naive encoding with □ = I and U = A, the claim suggests that (A^2, I) is a projected\nunitary encoding of 2A-I. However, since A^2 = A, this does not match 2A - I?\n4) I wonder whether GQSVT is necessary? For example, can we directly use the easier operator in e.g., Section 17.2 of\nhttps://www.cs.umd.edu/~amchilds/qa/qa.pdf ? It seems this operator does not explicitly require the chain to be reversible.", "reviewer_4": "This work generalizes the quantum speedup for Markov chain simulations, by proposing a way to apply it not only to\nreversible Markov chains but also to non-reversible ones. This aim and the associated important applications are quite clear.\nThe results themselves though are not so clear to follow.\nA first result, it is said, needs to know the stationary distribution (i.e. the one that we are trying to sample from). Once this is\nknown, there should be much more efficient ways to sample from it than MCMC! E.g. without taking any walk constraints.\nThis is not discussed.\nA second result does not need such knowledge. Still, the formulation of the propositions (for instance, starting by postulating\nthe Markov kernel P) does not make it clear which knowledge of P is required in order to select the polynomial and other\nparameters.\nIn particular, when making statements like \"P and P^* have roughly the same high maxima\", it is easy to imagine\ncounterexamples to this. This is probably why conditions are needed for the scheme to work. But how would one know if the\nconditions are satisfied or not, on a given problem instance? if/when we can trust the results of such scheme?\nIn this sense, the technical properties highlighted in this paper contain interesting points, but their relevance in a full\napplication context is not explained as clearly as it should be. I would advise the authors to revise the writing style in Section\n2, putting emphasis on what must be known or not and how one would use the construction in practice (facing a precise\noracle). Maybe, illustrate on a running example how all things work what is needed or not. Especially for this type of journal,\nthis sort of intermediate yet precise level, between introduction and technical proofs, is the most important part for the\nreaders and it is rather missing here."}, "author_response": "Response to Referee comments\nReferee 1\nComment: - What are the noteworthy results?\nThe authors propose two quantum algorithms for sampling the stationary distribution of a Markov chain.\nWhen the Markov chain is irreversible with respect to the target stationary distribution, the second algorithm\ncan provide an up-to-exponential speedup and does not require knowing the stationary distribution up to a\nnormalization constant. The resulting algorithmic acceleration goes beyond the quadratic speedup known for\nreversible chains but hinges on the reversibility on π-average condition.\nComment: - Will the work be of significance to the field and related fields? How does it compare to the\nestablished literature? If the work is not original, please provide relevant references.\nYes,theworkaddressestheefficientsamplingofstationarydistributionsofgenericMarkoviansystems,including\nthosepossiblybreakingthedetailedbalance,whichyieldsimportantimplicationsinvariousdomainsofphysical\nand information sciences. The first result summarized in Proposition 1 appears related to previous results\non Grover-type quadratic speedups. It would be interesting to compare and distinguish from the established\nliterature.\nResponse: We thank the reviewer for the suggestion. We have significantly revised the para-\ngraph surrounding Proposition 1 to clarify its motivation and contribution. In particular, we now\nemphasize that performance analysis in the nonreversible (and more generally nonnormal) case can\nonly be made rigorous through a newly introduced notion of spectral gap adapted to this broader\nsetting. This allows us to quantify algorithmic performance even when the relationship between the\neigenvaluesofthekernelandthatofthequbitizedwalkoperatorisnotassimpleasinthereversiblecase.\nChanges in the manuscript:\n• Paragraph 3 entirely rewritten.\nComment: - Does the work support the conclusions and claims, or is additional evidence needed?\nAnalytical and numerical examples were demonstrated to showcase the more-than-quadratic speedup of the\nsampling algorithm. However, the work does not seem to fully support, on a broader scope, the claim of\nimpact across various fields as it remains unclear for which applications the assumptions made in the analysis\nare valid.\nResponse: We agree that full theoretical guarantees for end-to-end application scenarios are difficult\nto establish. This limitation, however, is not unique to our work and is a common challenge when\napplying Markov chain methods in practice. To address this, we have expanded the discussion section\nto better highlight open questions and potential avenues for future work, especially concerning the\npractical validity of our assumptions in real-world applications.\nChanges in the manuscript:\n• Added section 6.b. entitled ”Research directions”.\nComment: - Are there any flaws in the data analysis, interpretation and conclusions? - Do these prohibit\npublication or require revision?\nThepaperreadsabitdenseandheavyofmathematicalnotations. Itwouldbehelpfulto(i)restructurethepaper\nto highlight the main result, (ii) include an overview figure that schematically illustrates the core algorithm,\nand (iii) summarize and clearly distinguish the novel contributions in the paper. It remains a bit unclear what\nare some sufficient and necessary conditions for reversibility on π-average, since this is a key condition for the\nefficiency of the sampling. Related to this, the method is illustrated only for a few contrived toy examples and\ntheclaimofbroaderapplicabilitytostatistics,machinelearning,andcomputationalmodelingisnotimmediately\nobvious. MinortyposneedtobecorrectedthroughoutthemaintextE.g.,∂B1(0)asthecomplexunitcirclein\ntherightcolumnofpage2isundefined. 1-—P(z)—2intheparagraphdiscussingthecomplementaritycondition\non page 4 needs to be revised, as P is now defined as the Markov transition kernel.\nResponse: In response to this comment, we have reorganized the manuscript to improve clarity\nand accessibility. (i) We have compiled all relevant material from the existing literature into an\nexpandedtechnicalbackgroundsection. Sections3,4,and5containthecorenovelresultsandarenow\nstructured to directly address the central research question posed in the introduction. (ii) We have\nincluded a new figure to illustrate the algorithmic flows of the two approaches. (iii) A new concluding\nsectionexplicitlysummarizesthenovelcontributionsandplacestheminthebroadercontextofexisting\nwork.\nChanges in the manuscript:\n• Figure to illustrate the algorithmic workflow.\n• Revised Conclusion section.\nComment: - Is the methodology sound? Does the work meet the expected standards in your field?\nThe mathematical (Markov chain theory) and algorithmic (projected unitary encoding and generalized quan-\ntum eigenvalue transformation) tools used for implementing the sampling algorithm are well established and\ndeveloped.\nComment: - Is there enough detail provided in the methods for the work to be reproduced?\nThe work includes extensive mathematical proofs but does not provide source code. However, it offers com-\nprehensive details presented mostly in a self-consistent manner, likely enabling the exact reproduction of the\nalgorithm and numerical results.\nReferee 2\nComment: I co-reviewed this manuscript with one of the reviewers who provided the listed reports. This\nis part of the Nature Communications initiative to facilitate training in peer review and to provide appropriate\nrecognition for Early Career Researchers who co-review manuscripts.\nResponse: We thank the reviewer for evaluating this work.\nReferee 3\nComment: This paper introduces quantum algorithms that provide up-to-exponential speedup for sampling\nfrom the stationary distribution of nonreversible Markov chains. It extends existing method to non-reversible\nchains.\nThe core algorithmic components are the Generalized Quantum Eigenvalue Transform (GQET) and the Gen-\neralized Quantum Singular Value Transform (GQSVT), which facilitate the construction of reflections through\nthe stationary distribution.\nQuantizingnonreversibleMarkovchainsisanimportantproblemsincenonreversiblechainscanbemoreefficient\nand may exhibit faster mixing. While I appreciate the effort put into this work, I have several concerns about\nits technical aspects.\nComment: 1)Iwouldliketoseeastrongermotivationforthesecondmethod,whichaddressescaseswhere\nthe stationary distribution is unknown up to a normalizing constant. From my experience, most papers on\nMarkov chain algorithms in physics, machine learning, and statistics assume that the stationary distribution is\nknown up to a normalizing constant, as is standard in statistical physics, Bayesian learning, and many relevant\nmodels. This makes me question why so much effort is dedicated to the second case.\nResponse: To strengthen the motivation for the second method, which does not assume knowledge\nofthestationarydistributionuptonormalization, wehaveaddedexamplesfrombiologicalsystemsand\nfinancial modeling where such knowledge is typically unavailable. Moreover, we discuss its relevance in\nout-of-equilibrium statistical physics, where the stationary distribution is often complex or analytically\ninaccessible, yet one is still interested in efficient sampling.\nChanges in the manuscript:\n• 6. Discussion b. Research direction.\nComment: 2) The technical section is difficult to follow, and the novelty appears somewhat limited. From\nmy reading, it seems (at least for the first method) that the authors simply apply the GQSVT algorithm, and\nthe results follow directly. Therefore I have three suggestions\nResponse: We have revised the structure of the paper to center it around a guiding question posed\nin the introduction. This makes the logical progression of our results more transparent and improves\nthe overall coherence of the presentation.\nChanges in the manuscript:\n• In the introduction: Can quantum algorithms accelerate the mixing of nonreversible Markov\nprocesses? In summary, we answer this question by:\n– analyzing known quantum singular value transforms in the context of nonreversible Markov\nchains (which requires simulating the time-reversal),\n– introducing quantum eigenvalue transforms to retrieve an approximation of the stationary\ndistribution without simulating the time-reversal.\nComment: 2.1) The authors may consider providing a more comprehensive overview of existing work, par-\nticularlyonGQSVTandprojectedunitaryencoding,astheseappeartobecentraltotheirapproach. Currently,\nthetechnicalbackgroundsectionistooterse. Infact,whileIamveryfamiliarwithabout70%ofthematerial,\nthe remaining 30% is impossible to follow. (see also comment 3)\nResponse: In response to the reviewer’s concern, we have significantly expanded and clarified the\ntechnical background section. The updated version now provides a more comprehensive overview of\nthe tools used, including Generalized Quantum Singular Value Transformation (GQSVT) and projected\nunitary encodings, ensuring that readers unfamiliar with all aspects of these techniques can still follow\nthe technical development.\nChanges in the manuscript:\n• Clarified technical background section.\nComment: 2.2) The authors should clarify how they utilize existing algorithms and whether their approach\nis merely a direct application of these methods.\nResponse: We have now clarified how our approach builds upon and extends existing quantum\nalgorithmic techniques. On the classical side, we study the spectral properties of the flat discriminant\nmatrix, while on the quantum side, we construct specific polynomials to enable spectral gap amplifica-\ntion in a manner reminiscent of Grover’s algorithm. These clarifications are integrated into the revised\nsections and proofs.\nComment: 2.3)Theauthorscouldimproveclaritybyhelpingreaders”connectthedots.”Evenafterreviewing\nthe original GQSVT paper, it is not entirely clear to me why Proposition 1 holds. For instance, I assume the\n√\nauthors apply GQSVT to the matrix D= PP⋆, but it is unclear why this corresponds to a reflection over π\nand how it approximates π.\nResponse: To improve readability and help readers “connect the dots,” we have added a new ex-\nplanatoryparagraphbeforeProposition1. ThisparagraphclarifiestheuseofGQSVT,theconstruction\nof the reflection operator, and how this enables the preparation of the stationary distribution. We have\nalsomovedthedetailedproofsofPropositions2and3totheMethodssectionforbetterstructuralflow.\nChanges in the manuscript:\n• New paragraph before Proposition 1.\n• Proofs of Propositions 2 and 3 moved to the Methods section.\nComment: 3) I may be overlooking something, but I am unable to verify some of the claims made in the\npaper. For example, at the end of the second paragraph in Section 2, the statement ”... Indeed, ..., .. a\ncomprehensiveexposition”doesnotseemtoholdinastraightforwardmanner. IfItakethenaiveencodingwith\n□ = I and U = A, the claim suggests that (A2,I) is a projected unitary encoding of 2A−I. However, since\nA2 =A, this does not match 2A−I?\nResponse: We appreciate the detailed reading of our assumptions. We have clarified in the\nmanuscript that the projected unitary encoding (PUE) must in fact be a symmetric PUE (SPUE).\nIn the example given, the operator A is assumed to be both unitary and symmetric, which ensures\nA2 =1. This corresponds to a reflection through the entire space. For the case where A is a projector\n(e.g., A=|ϕ⟩⟨ϕ|), we now provide an explicit computation.\nChanges in the manuscript:\n• Let us explain why it is sufficient to construct SPUE of |π⟩⟨π|. Assume that (U,□) is a SPUE\nof A = |ϕ⟩⟨ϕ|, for a normalized state |ϕ|. Recall the definition of the qubitized walk operator\nW =(2□□†−1)U. Compute:\n□†W2□=□†U(2□□†−1)U□\n=2□†U□□†U□−□†U2□\n(1)\n=2|ϕ⟩⟨ϕ|ϕ⟩⟨ϕ|−□†□\n=2|ϕ⟩⟨ϕ|−1.\nSince □ is a partial isometry, □†□ is the projection on its support. Thus, □†W2□ acts as\n2|ϕ⟩⟨ϕ|−1 on the support of □.\nComment: 4) I wonder whether GQSVT is necessary? For example, can we directly use the easier operator\nin e.g., Section 17.2 of https://www.cs.umd.edu/~amchilds/qa/qa.pdf ? It seems this operator does\nnot explicitly require the chain to be reversible.\nResponse: We thank the reviewer for raising this important point. We now emphasize in the\nintroduction that most prior quantum algorithms for Markov chains, particularly those based on\ncontinuous-time quantum walks (CTQWs), assume reversibility. Indeed, CTQWs are well-defined\nif equation 17.7 is a Schr¨odinger equation. In particular, the Laplacian should be symmetric. The\nunderlying walk must therefore also be symmetric, namely reversible with respect to the uniform\ndistribution. This structural constraint limits the applicability of such approaches to a relatively narrow\nsubclass of chains. By contrast, our methods are designed to handle more general nonreversible\nkernels,therebyextendingquantumsamplingtechniquesbeyondthescopeofCTQW-basedframeworks.\nReferee 4\nComment: This work generalizes the quantum speedup for Markov chain simulations, by proposing a way\nto apply it not only to reversible Markov chains but also to non-reversible ones. This aim and the associated\nimportant applications are quite clear.\nThe results themselves though are not so clear to follow. A first result, it is said, needs to know the stationary\ndistribution (i.e. the one that we are trying to sample from). Once this is known, there should be much more\nefficient ways to sample from it than MCMC! E.g. without taking any walk constraints. This is not discussed.\nAsecondresultdoesnotneedsuchknowledge. Still, theformulationofthepropositions(forinstance, starting\nby postulating the Markov kernel P) does not make it clear which knowledge of P is required in order to select\nthe polynomial and other parameters. In particular, when making statements like ”P and P∗ have roughly the\nsame high maxima”, it is easy to imagine counterexamples to this. This is probably why conditions are needed\nfor the scheme to work. But how would one know if the conditions are satisfied or not, on a given problem\ninstance? if/when we can trust the results of such scheme?\nResponse: We thank the reviewer for emphasizing the importance of clearly stating what is assumed\ntobeknownaboutthestationarydistribution. Wenowspecifyearlyinthepaperthat,unlessotherwise\nnoted, we understand “knowing the stationary distribution” to mean knowing it up to a multiplicative\nconstant. MCMC are well suited for this setting.\nWe have also expanded the introduction to better explain the high-dimensional motivation behind our\napproach. While rigorous parameter selection is indeed challenging, especially in complex or high-\ndimensional systems, we emphasize that this challenge is not unique to our work and is commonly\naddressed heuristically in practical applications.\nTo support practical implementation, we now include a discussion of how Quantum Phase Estimation\n(QPE)andanalyticaltechniquescanassistinestimatingspectralpropertiesoftheunderlyingoperators.\nAdditionally, we clarify that in the case of kinetic processes, the adjoint kernel P⋆ corresponds to\nthe process where initial and final velocities are flipped. This physical interpretation underpins the\nassumption that P and P⋆ exhibit similar behavior in such systems. A more detailed explanation and\nan illustrative example are now included in the Appendix to help ground this intuition.\nChanges in the manuscript:\n• In the introduction: The first method requires knowledge of π up to a multiplicative constant (in\nthis work, knowledge of π is always understood to be up to a multiplicative constant).\n• Expanded Introduction.\n• Expanded Discussion.\nComment: In this sense, the technical properties highlighted in this paper contain interesting points, but\ntheirrelevanceinafullapplicationcontextisnotexplainedasclearlyasitshouldbe. Iwouldadvisetheauthors\ntorevisethewritingstyleinSection2,puttingemphasisonwhatmustbeknownornotandhowonewoulduse\ntheconstructioninpractice(facingapreciseoracle). Maybe,illustrateonarunningexamplehowallthingswork\nwhat is needed or not. Especially for this type of journal, this sort of intermediate yet precise level, between\nintroduction and technical proofs, is the most important part for the readers and it is rather missing here.\nResponse: We have revised Section 2 to clearly distinguish the roles of R and R⋆, the reflection\noperators used in our construction. Additionally, we now elaborate on how eigenvalue estimates can\nbe obtained in practice using QPE or classical analytic techniques. To build intuition, we include a\nrunningexamplebasedonabottleneckprocess, whichhelpsillustratehowthemethodworksandwhat\ninputs are required from a user-facing oracle perspective.\nChanges in the manuscript:\n• Figure to illustrate the algorithmic workflow.\nReferee 1\nComment: I thank the authors for their efforts in revising the manuscript, ”Quantum Speedup for Nonre-\nversibleMarkovChains.”Theyhavesuccessfullyaddressedseveralkeyissues,particularlyinclarifyingthepaper’s\nmotivations and contributions. While the exposition of the technical background is significantly improved, it\nwould benefit from some further revisions. My specific recommendations are listed below.\nResponse: We thank the reviewer for their constructive feedback and valuable recommendations.\nComment: **Figure 1**\nThisfigurehelpsmakethepaper’sideasmoreaccessible. However"}, {"round": 2, "reviewer_comments": {"reviewer_1": "I thank the authors for their efforts in revising the manuscript and addressing my concerns.", "reviewer_4": "After reading the reply to referees, I must say that the paper seems now much clearer to me; and hopefully this is\nrepresentative of more general readers' opinion.\nMost arguments are clarified and I appreciate the additions about the context.\nThere is just one technical point which I found too sketchy, namely about the oracle access. Indeed, the reference provided\nexplains how to use the oracle of a reversible Markov chain, while here the purpose is precisely to treat non-reversible ones.\nIt would seem necessary to explain how, having access to the \"P\" oracle, one would construct the \"P ^star\" oracle. Indeed,\nin general, it is not obvious (to me) that being able to predict where the walk will go next, implies you could as easily deduce\nwhere the walk is coming from.\nI may just be overlooking some obvious property; or maybe not. If this point was clarified, all my comments would be\nresolved."}, "author_response": ""}, {"round": 3, "reviewer_comments": {"reviewer_4": "I just had one remaining concern: how, in general, an oracle model for the direct chain would imply the oracle model for the\nreverse chain.\nWith respect to the material presented now by the authors, I was worried about two things:\n- how to convert the transition probability, which thus indeed assumes we can evaluate pi(x)/pi(y) for neighboring states;\nwhy\nnot, in many settings.\n- but also: if given x, we can easily find all y for which P(x,y)≠0, this would not necessarily imply that given y, we can\nefficiently find all x for which P(x,y)≠0.\nI therefore found it important to clarify that, with respect to the most general setting, we have to assume these things.\nThe authors now correctly clarify that this is an additional assumption, which often holds in practice. With this I see no further\nissues."}, "author_response": ""}]

Decision letter