Quantum circuits for exact unitary $t$-designs and
applications to higher-order randomized benchmarking

Unitary designs have been studied in terms of representation theory RS2009 ; RS2011 since the operator $U^{\otimes t}\otimes U^{\dagger\otimes t}$ in the definition can be regarded as a representation of the unitary group. Our construction is based on representation theory, where irreducible decomposition of the operator plays an important role. A brief introduction of irreducible representations (irreps) will be provided in Section VII.1. Here, we mention a couple of well-known facts that are necessary to state our main result.

Any irrep of the unitary group can be indexed by a non-increasing integer sequence $\lambda:=(\lambda_{1},\lambda_{2},\dots,\lambda_{d})$ of length $d$, i.e., $\lambda_{i}\in\mathbb{Z}$ for $i=1,\dots,d$, and $\lambda_{1}\geq\lambda_{2}\geq\dots\geq\lambda_{d}$ B2004 . In particular, spherical representations of ${\sf U}(d)$ with respect to ${\sf K}:={\sf U}(d_{1})\times{\sf U}(d-d_{1})$ are of great importance in the construction. Let $\Lambda_{\rm sph}(d_{1},d,t)$ be a set of all non-increasing integer sequences $\lambda$ in the form of

where $d_{1}\leq d/2$ and $t\geq\lambda_{1}\geq\dots\geq\lambda_{d_{1}}\geq 0$. The spherical representation is the irrep indexed by $\lambda\in\Lambda_{\rm sph}(d_{1},d,t)$ GW2009 . For a spherical representation $\lambda\in\Lambda_{\rm sph}(d_{1},d,t)$, a zonal spherical function $Z^{(d_{1})}_{\lambda}(z_{1},\dots,z_{d_{1}})$ is defined by the unique bi-${\sf K}$-invariant function JC1974 ; R2010 ; BNOZ2020 . The zonal spherical functions are a certain type of symmetric polynomials, and can be explicitly written down (see, e.g., Appendix A of Ref. BNOZ2020 ).

Our main technical result is to construct a strong unitary $t$-design on ${\sf U}(d)$ from those on ${\sf U}(d_{1})$ and on ${\sf U}(d-d_{1})$.

Theorem 3 follows from a more general result BNOZ2020 shown by some of the authors, which works not only for the unitary group but also for a broader class of compact groups. For the sake of completeness, we provide a direct proof of Theorem 3 in Sec. VII.

We then claim that

where $\omega=\exp[\frac{2\pi}{t+1}]$ is the $(t+1)$-th root of unity, is a strong unitary $t$-design on ${\sf U}(1)$ for any $t$. This is easily checked by direct calculations:

where $\delta_{rs}$ is the Kronecker delta. Hence, we have $\mathbb{E}_{U\sim{\sf W}_{1}}[U^{\otimes r}\otimes U^{\dagger\otimes s}]=\mathbb{E}_{U\sim{\sf H}(1)}[U^{\otimes r}\otimes U^{\dagger\otimes s}]$ for $0\leq\forall s,r\leq t$, implying that ${\sf W}_{1}$ is a strong unitary $t$-design.

From Theorem 3 and ${\sf W}_{1}$, a strong unitary $t$-design on a qudit can be inductively constructed.

In this construction, it is important to obtain zeros for the zonal spherical functions $Z_{\lambda}^{(1)}$. This is computationally feasible since they are polynomials of one variable and are explicitly given (see Appendix A of Ref. BNOZ2020 ). Furthermore, $\Lambda_{\rm sph}(1,d,t)$ contains only $t$ elements. Hence, we need to solve $t$ polynomials with one variable, which is tractable as far as $t$ is not too large.

We now consider a strong unitary $t$-design on $N$ qubits. Again using Theorem 3, we obtain the quantum circuit on $(N+1)$ qubits based on that on $N$ qubits. See also Fig. 1.

Corollary 5 implies that a quantum circuit for an exact unitary $t$-design can be inductively constructed from a strong unitary $t$-design on one qubit, i.e., that on ${\sf U}(2)$. Furthermore, a strong unitary $t$-design on ${\sf U}(2)$ can be constructed using Corollary 4. Thus, combining Corollaries 4 and 5, we obtain a quantum circuit for an exact unitary $t$-design for any $t$ and on an arbitrary number of qubits.

Note that the circuit, constructed in this way, can be explicitly decomposed into two-qubit gates. The controlled unitary ${\rm Ctrl\mathchar 45\relax}{\sf Q}_{N}$ part contains up to three-qubit gates, if the circuit ${\sf Q}_{N}$ on $N$ qubit is already decomposed into two-qubit gates. The three-qubit gates can be easily rewritten as a series of two-qubit gates. Also, the $X$-rotation controlled by $N$ qubits, $R_{X}(\bm{\theta}_{\lambda})$, can be decomposed into a sequence of two-qubit gates of polynomial length using sufficiently many number of ancillary qubits, which is based on a classical oracle that computes the angle $\theta_{\lambda}^{(\bm{j})}$ from $\bm{j}$ (see Appendix A).

In special cases, we can find a much more concise construction based on a similar technique.

Analytically, we can prove that there exist unitaries $V_{1}$ and $V_{2}$ such that ${\sf C}(4)V_{1}{\sf C}(4)V_{2}{\sf C}(4)$ is an exact unitary $4$-design on $2$ qubits BNOZ2020 . Also, an algorithm for computing the unitaries $V_{1}$ and $V_{2}$ is given. It, however, turns out from numerics that it is not necessary to apply two extra unitaries if we choose a proper unitary $U_{c}$, which leads to Proposition 6. An explicit form of the unitary $U_{c}$ is numerically obtained and is provided in Appendix B. Note that the existence of $U_{c}$ is confirmed numerically, so the statement holds up to the numerical precision.

This construction is only for a $4$-design on $2$ qubits, but the number of unitaries in the $4$-design is much smaller than that of Corollary 5. It is an open problem whether a similar construction works for higher-designs on a larger number of qubits.

To quantitatively evaluate the complexity of the quantum circuit for an exact unitary $t$-design obtained in Corollary 5, we provide an order estimate of the number $G(N)$ of two-qubit gates in the circuit. Assuming $2^{N}\gg t$ and using the fact that $|\Lambda_{\rm sph}(D,2D,t)|=O(e^{\pi\sqrt{2t/3}})$ due to the Hardy and Ramanujan formula for the asymptotics of the number of partitions, we obtain

to the leading order of $N$. Hence, it is necessary to use exponentially many two-qubit gates as the number of qubits increases. This inefficiency of the quantum circuit may be intrinsic since the construction is inductive.

There is another source of inefficiency. In Corollary 5, it is necessary to find zeros of zonal spherical functions (see Eq. (19)) for all $\lambda\in\Lambda_{\rm sph}(D,2D,t)$. The zonal spherical function is given in terms of the summation of the (normalized) Schur polynomials (see Appendix A in Ref. BNOZ2020 ). It is unlikely that the Schur polynomials have polynomial size algebraic formulas in general SchurPoly2019 . Moreover, the number of variables for each $Z_{\lambda}^{(D)}$ is $D=2^{N}$. Hence, finding zeros of zonal spherical functions is computationally intractable.

In total, the construction for an exact unitary $t$-design on a large number of qubits is inefficient from both quantum-circuit and classical-computation viewpoints. We, however, think that our construction and our proof technique will form a solid basis of searching more efficient constructions of exact, as well as approximate, unitary designs. We also emphasize that, despite its inefficiency, our construction is of practical use on a few-qubit system, as we seek in the following sections.

Before we conclude this section, we comment on advantages of our construction of an exact unitary $t$-design over a direct implementation of a Haar random unitary. A naive way of implementing a Haar random unitary by a quantum circuit consists of three steps. First, we sample a Haar random unitary as a matrix by a classical computer. A classical algorithm for this is known M2007 , but it is trivially inefficient since the size of the matrix is exponentially large. We then classically compute a decomposition of the unitary matrix into a sequence of two-qubit unitaries, providing a classical description of a quantum circuit for the unitary. This step is also inefficient, and the resulting quantum circuit is almost surely composed of the exponentially many number of two-qubit gates. Finally, we implement the circuit in practice.

This quantum circuit for a Haar random unitary is inefficient in terms of the number of qubits, and so, cannot be of practical use in a large system. Even in a small system, this naive implementation has a crucial difficulty that, every time a unitary is sampled, the above protocol outputs a quantum circuit with a rather different sequence of various two-qubit gates. This implies that, in each sampling, one needs to significantly modify the quantum circuit. This is in a sharp contrast to our quantum circuit for a unitary $t$-design based on Corollary 5 since it has a fixed structure. In each sampling, only what one needs to do is to randomly choose single-qubit gates, or more precisely elements of ${\sf U}(1)$ from ${\sf W}_{1}$ (see Eq. (12)), and to plug them into the quantum circuit with a fixed structure. This will help practical implementations of the circuit in small systems.

It should be also noted that the single-qubit gates in our construction can be sampled from a discrete set, though sampling from a continuous gate set is necessary in the direct implementation of a Haar random unitary. This is another advantage of our construction.

A noise $\mathcal{E}$ acting on a $q$-qubit system is formulated by a completely-positive and trace-preserving (CPTP) map. Let $d$ be defined as $d:=2^{q}$. The average fidelity and the unitarity are defined by

respectively, where $\mathcal{E}^{\prime}(\rho):=\mathcal{E}(\rho-I/d)$, and $|\!|A|\!|_{2}=(\operatorname{Tr}[A^{\dagger}A])^{1/2}$ is the Schatten 2-norm. The average fidelity satisfies $1/(d+1)\leq F(\mathcal{E})\leq 1$, and $F(\mathcal{E})=1$ if and only if the system is noiseless, i.e., $\mathcal{E}$ is the identity channel, while the unitarity satisfies $0<u(\mathcal{E})\leq 1$, and $u(\mathcal{E})=1$ if and only if the noise is coherent, i.e., $\mathcal{E}$ is a unitary channel. The unitarity is an important metric in the context of QEC since coherent noise is known to be hard to correct in general KLDF2016 ; SWS2015 ; SFK2017 .

In the RB-type protocols, it is more natural to use a fidelity parameter $f(\mathcal{E})$ rather than the average fidelity itself. It is defined by

and satisfies $-1/(d^{2}-1)\leq f(\mathcal{E})\leq 1$.

We next introduce a self-adjointness of the noise. For any linear map $\mathcal{E}$, an adjoint map $\mathcal{E}^{\dagger}$ is defined by $\operatorname{Tr}[A\mathcal{E}(B)]=\operatorname{Tr}[\mathcal{E}^{\dagger}(A)B]$. A noise $\mathcal{E}$ is called self-adjoint if $\mathcal{E}=\mathcal{E}^{\dagger}$, which is equivalent to that all the Kraus operators of $\mathcal{E}$ are self-adjoint.

The self-adjointness $H(\mathcal{E})$ of the noise $\mathcal{E}$ is defined by

The normalization constant is chosen such that $0\leq H(\mathcal{E})\leq 1$. Obviously, $H(\mathcal{E})=1$ if and only if $\mathcal{E}$ is self-adjoint, i.e. $\mathcal{E}=\mathcal{E}^{\dagger}$. Note that the self-adjointness has two contributions from the noisy map $\mathcal{E}$, one is from the unital part and the other from the non-unital part. The non-unital part of the noise makes the self-adjointness less than one since, if $\mathcal{E}$ is not unital, then $\mathcal{E}^{\dagger}$ is not trace-preserving, which implies that $\mathcal{E}\neq\mathcal{E}^{\dagger}$.

To clearly separate the two contributions, we introduce a self-adjointness parameter $h(\mathcal{E})$. Using $\mathcal{E}^{\prime}(\rho)=\mathcal{E}(\rho-I/d)$, we defined it by

The self-adjointness parameter $h(\mathcal{E})$ is related to the self-adjointness $H(\mathcal{E})$ and the unitarity $u(\mathcal{E})$ by

where $|\alpha_{\mathcal{E}}|$ is a measure of the non-unital part of the noise (see Subsec. VIII.1 for the definition). We can clearly observe that $H(\mathcal{E})$ consists of two factors, the unital part $h(\mathcal{E})$ and the non-unital part $|\alpha_{\mathcal{E}}|$.

The three metrics of noises, namely, fidelity, unitarity, and self-adjointness, all capture different properties of the noises. The fidelity reveals the first order property of the noises, while the unitarity and the self-adjointness, which are independent to each other, reveal the second-order. In order to improve noisy quantum devices, it is of crucial importance to obtain the information of noise as much as possible. Hence, it is certainly of practical use to introduce the self-adjointness as a new metric of noise. In addition, we argue in the next subsection that the self-adjointness has important implications for QEC.

The most important family of self-adjoint noises is stochastic Pauli noises, whose Kraus operators are all proportional to Pauli matrices. In QEC, Pauli noises are the standard yet most important class of noises both in theory and in practice. From a theoretical perspective, Pauli noises are easy to numerically handle. Hence, most numerical calculations have been carried out by assuming Pauli noises, and it has been confirmed that QEC has preferable features, such as exponential decreases and threshold behaviors of logical error rates, if the noise is Pauli.

The noise being Pauli is also practically preferable in experimental realizations of QEC since it typically simplifies the decoding tasks. This is especially the case for stabilizer codes, such as surface and color codes, whose standard decoders are to estimate what types of Pauli operators should be applied on which physical qubits during recovery operations. For stochastic Pauli noises, if the estimation goes well, the state is fully retrieved with high probability by applying Pauli operators to the suitable physical qubits. In contrast, it is not possible to fully correct non-Pauli noises by applying Pauli operators since they generate undesired coherence between different code spaces. Thus, QEC of non-Pauli noises generally suffers from degradation of logical error rates when the standard decoders are used SFK2017 ; BEKP2018 or requires more complicated algorithms for retrieving the performance of QEC. Neither of them is preferable in practice since it induces additional experimental difficulties.

For these reasons, it is desirable to check that the noise on an experimental system is stochastic Pauli. To this end, the self-adjointness provides useful information since, if $H(\mathcal{E})\ll 1$, then the noise is far from self-adjoint and cannot be approximated by Pauli noises. This implies that the practical situation differs from the standard assumption in theoretical studies of QEC and incurs additional difficulties on decoding procedure. Thus, the self-adjointness provides practical information about the feasibility of QEC using Pauli-based decoders.

Note that the difficulty of QEC for non-Pauli noises, captured by the self-adjointness, highly depends on the assumptions in quantum error correction schemes. When any decoding procedure is available, it would not be so important whether the noise is Pauli or non-Pauli. When this is the case, the unitarity will be a more suitable metric of noise relevant to the feasibility of QEC KLDF2016 ; SWS2015 ; SFK2017 . Note also that non-Pauli noises can be always transformed to a Pauli noise by Pauli-twirling. However, Pauli-twirling induces additional noise onto the system and, as a result, the performance of QEC will degrade. Thus, it is practically desirable to manufacture the system so that the noise is stochastic Pauli.

We also provide a pedagogical example of noise, where performance of QEC can be directly captured by the self-adjointness but not by fidelity nor unitarity. Consider a $\theta$-rotation error around the $X$-axis on one qubit, i.e., $\exp[i\theta X/2]$, where $X$ is the Pauli-$X$ operator. The average fidelity $F_{\theta}$ and the self-adjointness $H_{\theta}$ can be obtained as

The unitarity is $1$ for any $\theta$.

One may expect that the $\pi/2$-rotation error is easier to correct than the $\pi$-rotation since the former has higher fidelity than the latter. However, this is not the case since $\pi$-rotation is simply a perfect bit-flip that can be trivially corrected, while the $\pi/2$-rotation error is known to be particularly hard to correct DP2017 . Thus, neither the average fidelity nor the unitarity, which is $1$ for both errors, is a good metric of the error correctability. In contrast, the self-adjointness clearly captures whether the error can be corrected, at least in this case, since $H_{\pi/2}$ and $H_{\pi}$ are the minimum and the maximum values of the self-adjointness, respectively.

We now introduce the $t$-RB using an exact unitary $2t$-design ${\sf U}_{2t}:=\{U_{i}\}_{i}$. As is the case for the standard RB, we assume that the noise is gate- and time-independent, so that the noisy implementation of ${\sf U}_{2t}$ is given by $\{{\cal G}_{i}:=\mathcal{E}\circ\mathcal{U}_{i}\}_{i}$, where $\mathcal{E}$ is the CPTP map that represents the noise, and we used the notation that $\mathcal{U}(\rho):=U\rho U^{\dagger}$.

Let $O_{\rm ini}$ and $O_{\rm meas}$ be the initial and measurement operators, respectively, which we assume to be Hermitian. We first apply a sequence of unitaries $U_{\bm{i}}=U_{i_{m}}\dots U_{i_{1}}$ onto the initial operator $O_{\rm ini}$. Each $U_{i_{n}}$ is chosen uniformly at random from ${\sf U}_{2t}$, which we denote by $U_{\bm{i}}\sim{\sf U}_{2t}^{\times m}$. We then apply its inverse $U_{i_{m+1}}:=U_{\bm{i}}^{-1}$, and measure $O_{\rm meas}$.

If the system is noiseless, $\mathcal{E}={\rm id}$, this protocol results in a trivial expectation value that

due to the inverse unitary $U_{i_{m+1}}$. However, when the system is noisy, the expectation value becomes

which in general differs from $\operatorname{Tr}[O_{\rm meas}O_{\rm ini}]$. The basic idea of the RB-type protocol is to extract some information about the noise $\mathcal{E}$ from the difference.

In the $t$-RB, we especially focus on the average of the $t$-th power of the expectation value over all choices of the unitary sequence. That is,

Using the representation-theoretic technique, it can be shown that $V^{(t)}(m,\mathcal{E}|O_{\rm ini},O_{\rm meas})$ is generally given in the following form:

where $\lambda$ labels the irreps of the unitary group, $\hat{A}^{(t)}_{\lambda}$ and $\hat{C}^{(t)}_{\lambda}(\mathcal{E})$ are $m_{\lambda}\times m_{\lambda}$ matrices with $m_{\lambda}$ being the multiplicity of the irrep $\lambda$. This is well-known in the literature of RB-type protocols, but we provide a proof in Sec. VIII.3 for completeness.

Despite its abstract expression, Eq. (34) has an important implication that the matrix $\hat{C}^{(t)}_{\lambda}(\mathcal{E})^{m}$ depends only on $\mathcal{E}$ and $m$, but not on $O_{\rm ini}$ and $O_{\rm meas}$. Hence, from the experimental data of $V^{(t)}(m,\mathcal{E}|O_{\rm ini},O_{\rm meas})$ for various $m$, it is in principle possible to estimate the matrix $\hat{C}^{(t)}_{\lambda}(\mathcal{E})$, which contains certain information of the noise $\mathcal{E}$, in the way independent of $O_{\rm ini}$ and $O_{\rm meas}$.

In practice, the most important situation is when the representation is multiplicity-free, i.e., $m_{\lambda}=1$ for any $\lambda$. In this case, $V^{(t)}$ reduces to a much simpler form:

where $A^{(t)},C^{(t)}_{\lambda}(\mathcal{E})\in\mathbb{R}$. Note that $|C^{(t)}_{\lambda}(\mathcal{E})|\leq 1$ since $V^{(t)}$ is a bounded function. Hence, in this case, $V^{(t)}$ becomes a sum of some exponentially decreasing functions with respect to $m$.

To be more concrete, let us recall the standard RB, corresponding to the $1$-RB. As shown in Ref. EAZ2005 , $V^{(1)}$ is given by

where $A_{0}^{(1)}$ and $A_{1}^{(1)}$ depend only on $O_{\rm ini}$ and $\mathcal{E}(O_{\rm meas})$, and $f(\mathcal{E})$ is the fidelity parameter of the noise $\mathcal{E}$. Thus, by fitting experimentally obtained data of $V^{(1)}$ for different $m$ with the fitting function $F(m)=A+B\alpha^{m}$, we can estimate the fidelity parameter $f(\mathcal{E})$.

In the RB protocol, it is important to use exact unitary designs because designs are used many times, sometimes a few hundreds to a thousand, in a single run of the protocol. To illustrate this, let us consider the $1$-RB when the unitary $2$-design in the protocol is $\epsilon$-approximate.

Let $m$ be the length of the unitary sequence as above. It is straightforward to show that,

where $E_{i}$’s are some constants that depend on $O_{\rm ini}$, $\mathcal{E}(O_{\rm meas})$, $f$, and how the design differs from the exact one. See Subsec. VIII.4 for the derivation. Compared to the $1$-RB with exact ones, i.e., Eq. (36), fitting this function with respect to $m$ is much harder since it is not a simple exponential decay.

The fitting may go well if $\epsilon\ll f^{m}/m$. This requires a very high precision of the design since $m$ can be a few hundreds in actual experiments. For instance, when $f=0.95$, the degree $\epsilon$ of approximation of the unitary design should be order $10^{-5}$ or so. Although it is possible to achieve this degree of approximation by a sufficiently long quantum circuit BHH2016 ; NHKW2017 ; HMHEGR2020 , the RB becomes unpractical if we use such a long circuit at every use of a unitary design in the protocol and repeat it a few hundreds times.

There might be a possibility to improve Eq. (37) by using different constructions of approximate unitary designs at every step, by which the differences from the exact design may become random so that they cancel out in total. This will be an interesting question, but at this point, it is not clear if such a technique works. Also, even if it works, we need to assume additional structures of approximate constructions.

The higher-order RB with approximate designs will incur more difficulty in practice. Since it uses higher moment of the outcomes, the fitting function becomes more complicated than Eq. (37) when one uses approximate designs. Similarly to the $1$-RB with approximate $2$-designs, much better degree of approximation, that is, longer quantum circuits, will be needed, which is not practical. Thus, we conclude that exact unitary designs are of crucial importance in a practical implementation of the $t$-RB.

We next focus on the $2$-RB using exact unitary $4$-designs, and show that the $2$-RB reveals the self-adjointness of the noise. To this end, we set the initial operator $O_{\rm ini}$ to a traceless one, i.e., $\operatorname{Tr}[O_{\rm ini}]=0$. This setting, together with the fact that the noise is trace-preserving, makes the representation multiplicity-free (see Appendix C). Hence, the expectation value $V^{(2)}(m,\mathcal{E}|O_{\rm ini},O_{\rm meas})$ for the $2$-RB is given by a sum of exponentially decaying functions as shown in Eq. (35).

Note that the expectation value for a traceless initial operator can be obtained by performing the same experiment for two different quantum states $\rho$ and $\rho^{\prime}$, and by taking the difference of the expectation values before they are squared. That is,

where $\Delta=\rho-\rho^{\prime}$ is a traceless operator.

Our second main result in this paper is about $V^{(2)}(m,\mathcal{E}|\Delta,O_{\rm meas})$ as summarized in Theorem 7.

See Subsec. VIII.5 for the proof.

In the single-qubit case, $V^{(2)}$ is a sum of two exponentially decaying functions with respect to $m$. Hence, from the double-exponential fitting of the experimental data of $V^{(2)}$, we can simultaneously estimate $u(\mathcal{E})$ and $\frac{9}{10}f(\mathcal{E})^{2}-\frac{1}{5}u(\mathcal{E})+\frac{3}{10}h(\mathcal{E})$. Since it can be shown that the former is not less than the latter, we can estimate which of the two decaying rates corresponds to which quantity without any ambiguity. It is also possible to estimate the fidelity parameter $f(\mathcal{E})$ from the same data set by computing $V^{(1)}(m,\mathcal{E}|\Delta,O_{\rm meas})$ because a unitary $2$-design is also a unitary $1$-design. Thus, from the experiment of the $2$-RB on a single qubit, all of $f(\mathcal{E}),$ $u(\mathcal{E})$, and $h(\mathcal{E})$ can be estimated simultaneously.

In multi-qubit systems, $V^{(2)}$ has a little more complicated form and consists of four exponentially decaying functions. Also, the decaying rates do not directly correspond to neither the unitarity nor the self-adjointness parameter. We observe from Eq. (41) that $h(\mathcal{E})$ can be obtained from a linear combination of the decaying rates $C_{\lambda}(\mathcal{E})$, the value of $u(\mathcal{E})$, and $f(\mathcal{E})$.

One may think that, in the case of multiple qubits, it is practically intractable to accurately fit four exponentially decaying functions from experimental data because each data point has an error. This difficulty can be circumvented by choosing appropriate initial and measurement operators. By doing so, we can set some of $A_{\lambda}$ zero in the ideal situation (see Tab. 1). This allows us to estimate the decaying rates one by one. Note that the initial and measurement operators in Tab. 1 are all diagonal in the computational basis. Hence, it suffices to perform the experiments for the four initial operators $|00\rangle,|01\rangle,|10\rangle$, and $|11\rangle$, with the measurement in the computational basis. From the data of these experiments, it is possible to reproduce all cases listed in Tab. 1 by post-processing.

In the multi-qubit case, the ambiguity remains to decide which of the decaying rates corresponds to which quantity. This is the case even when we use the above step-by-step estimation of the rates since, for instance, it is not clear if the unitarity $u(\mathcal{E})$ is larger or smaller than $C_{\rm I}$. In this case, we need to additionally perform the unitarity benchmarking WGHF2015 ; HHFFW2019 to separately estimate $u(\mathcal{E})$. If we have an estimated value of $u(\mathcal{E})$, the step-by-step estimation allows us to decide all decaying rates without any ambiguity.

See Sec. V.1 and Sec. IX.1 for the performance of 2-RB in concrete cases.

The $t$-RB for $t\geq 2$ inherits most of the desired properties of the RB-type protocols. For instance, it is experimentally-friendly since, apart from using higher-designs, the difference of the $t$-RB from the standard RB ($1$-RB) is only taking the $t$-th power of the expectation value before the average. It is also true that the $t$-RB is free from SPAM errors (see Eqs. (34) and (35)).

The property that the standard RB does have and the $t$-RB does not in general is the scarability. This is for two reasons. First, no efficient construction of exact unitary $2t$-designs is known for $t\geq 2$ so far. Second, in the $t$-RB protocol, it is necessary to apply the inverse unitary at the end of the unitary sequence. Hence, we need to beforehand compute the inverse of each sequence. When the system is large, the task is intractable in general. This difficulty is avoided in the standard RB by using the Clifford group, which is an exact unitary $2$-design. Since the inverse is contained in the group, we can find the inverse relatively easily. One may expect that the difficulty of finding the inverse could be also avoided in the $t$-RB by using the $2t$-design that is also a group, which is called a unitary $2t$-group BNRT2020 . However, it is known that unitary $2t$-groups do not exist for $t\geq 2$ if the number of qubits $\geq 3$. Thus, in the $t$-RB for $t\geq 2$, the hardness of finding the inverse in a large system is inevitable.

Nonetheless, we emphasize that, in the current experimental situations, the RB-type protocols for more than three qubits are practically intractable due to the limitation of the coherent time. Thus, the experimental use of the RB-type protocols is currently aiming to characterize the noise on one- or two-qubit systems in a concise manner. Considering this fact, even if the $t$-RB is not scalable, it is practically useful and beneficial: it is as concise as the standard RB and provides more information about the noise, such as self-adjointness.

When the $2$-RB is practically implemented, there are two additional concerns. One is originated from the fact that the expectation value $\langle O_{\rm meas}\rangle_{O_{\rm ini},\bm{i}}$ is obtained from a limited number of measurements in the basis of $O_{\rm meas}$, resulting in an error due to a finite number of measurements. The other originates from the evaluation of the average $\mathbb{E}_{U_{\bm{i}}\sim{\sf U}_{4}^{\times m}}$ over the sequence of unitaries in the $4$-design. Ideally, all sequences in ${\sf U}_{4}^{\times m}$ should be taken, but practically, the average is often evaluated from a small subset in ${\sf U}_{4}^{\times m}$ of randomly chosen sequences, leading to an additional error of estimation.

Taking sufficiently many measurements and samplings of unitary sequences will reproduce the analytical results with high accuracy. However, it is complicated to analytically derive the numbers sufficient for achieving a desired accuracy. We hence perform numerical experiments and show that experimentally-tractable numbers of samplings are sufficient for a reliable $2$-RB.