Mitigation of readout noise in near-term quantum devices by classical post-processing based on detector tomography

Filip B. Maciejewski,1, ∗ Zoltán Zimborás,2, 3, 4, † and Michał Oszmaniec5, ‡ 1University of Warsaw, Faculty of Physics, Ludwika Pasteura 5, 02-093 Warszawa, Poland 2Wigner Research Centre for Physics of the Hungarian Academy of Sciences, H-1525 Budapest, P.O.Box 49, Hungary 3BME-MTA Lendület Quantum Information Theory Research Group, Budapest, Hungary 4Mathematical Institute, Budapest University of Technology and Economics, P.O.Box 91, H-1111, Budapest, Hungary 5International Centre for Theory of Quantum Technologies, University of Gdansk, Wita Stwosza 63, 80-308 Gdansk, Poland


I. INTRODUCTION
In recent years, quantum technologies have been rapidly developing. Scientists and engineers around the world share the hope and fascination caused by the possibility of creating devices that would allow for the manipulation of delicate quantum states with unprecedented precision [1]. Due to the advent of quantum cloud services (IBM [2, 3], Rigetti [4], DWave [5]), any researcher has a possibility to perform experiments on actual quantum devices. However, if one really hopes for utilizing such near-term devices for real-life applications like quantum computation [6] or quantum simulations [7], experimental imperfections must be taken into account. Hence, to properly characterize noise occurring in the devices and to develop error correction and mitigation schemes that may help to fight it have become tasks of fundamental importance [8][9][10][11][12][13]. In the present work, we address this problem for the noise affecting quantum measurements.
In theoretical considerations about quantum information tasks, a quantum device is often assumed to perform perfect measurements. In practice this assumption may be significantly violated due to experimental imperfections and noise. Specifically, broad class of errors occurring in quantum devices may be described as State Preparation And Measurement (SPAM) errors [14,15]. It follows that if FIG. 1. Pictorial representation of our error-mitigation procedure. (i) In the first stage, one performs the tomography of a noisy detector M noisy (red semicircle). (ii) In the next stage, when measuring an arbitrary quantum state ρ, one employs a post-processing procedure on the measured statistics through the application of Λ −1 , the inverse of a stochastic noise map obtained in the QDT. This gives access to the statistics that would have been obtained in an ideal detector M ideal (green semicircle). state preparation errors are negligible, one is able to reconstruct a POVM (Positive-Operator Valued Measure) associated with a given measurement device, and therefore to infer measurement errors. A standard method of such reconstruction is to perform Quantum Detector Tomography (QDT) [16]. Results of the QDT may be used to compare the reconstructed POVM with the ideal one [16] and to analyze the nature of measurement errors [17][18][19][20]. In this work, we propose to do more -we show how to use the knowledge about the POVM associated with a mea-surement device to correct the results obtained in further experiments performed on this device. Such a correction is possible provided the noise which affects the measurement is of a classical stochastic type and non-degenerate (invertible). In other words, such type of noise is equivalent to an invertible classical post-processing of experimental statistics, and it may also be inverted solely by classical post-processing. The general idea of mitigating the effects of such a noise is the following. Assuming that one has access to the aforementioned invertible stochastic map describing a noise, it is straightforward to use the (generally non-stochastic) inversion of this map to classically reverse effects of a noise simply by multiplying the vector of experimental statistics (obtained in a further experiment) by this inverted matrix, see Fig. 1. The main aim of this work is to present such a classical correction scheme together with the analysis of its accuracy. Specifically, we show how the deviations from classical noise model and finite-size statistics affect our procedure, by providing upper bounds for the distance of the corrected probability distribution from the ideal noise-reversed scenario.
Beside the theoretical description, we also test our procedure experimentally. First, we present data that suggests that IBM's and Rigetti's quantum processors are affected by a significant readout noise that can be described (to a good approximation) by the classical model characterized above. Both architectures consist of superconducting transmon qubits [21], which may suggest that the classical noise may be the dominant form of measurement noise in such devices. Second, we test the correction scheme on the IBM's five-qubit device ibmqx4 and conclude that indeed it significantly compensates for the effects of this noise for several one-and two-qubit experiments, including Quantum State Tomography (QST) [6], Quantum Process Tomography (QPT) [6], the implementation of non-projective measurements [22], and two quantum algorithms (Grover's search [23] and the Bernstein-Vazirani algorithm [24]). Furthermore, we test our scheme for the 5-qubit experiments concerning the implementation of certain probability distributions.
We will now comment on related works. There is a variety of recent research concerning mitigation of different types of errors in the contemporary quantum devices [8][9][10][11][12][13]. In particular, in the works [9,10], the classical postprocessing of the data was used to mitigate 'assignment errors', which may refer to the related procedure. Quantum Detector Tomography for superconducting qubits, along with other characterizations of measurement errors, was presented in the Ref. [20] for Rigetti devices. Characterization of SPAM errors for IBM and Rigetti devices has been presented in Ref. [15].
The paper is organized as follows. In Section II, we present the necessary theoretical concepts, including the POVMs formalism, the QDT scheme and the description of classical noise. Section III is devoted to the main idea of this work; it consists of a detailed description of the statistics correction procedure, which is preceded by stating the necessary assumptions. In Section IV, we analyze how violations of the assumptions affect our correction procedure. Section V contains a summary of the correction procedure in the form of pictorial representation. The scheme presents practical steps that need to be done in the case of noisy projective measurement in computational basis. In Section VI, we present experimental results from IBM's and Rigett's devices that provide insight into the magnitude of SPAM errors. Section VII consist of experimental results from applications of our correction scheme for exemplary quantum information tasks in IBM's five-qubit device. Finally, we present the conclusions and some possible future research directions in Section VIII.

II. THEORETICAL BACKGROUND
In this section, the necessary theoretical background and mathematical tools are shortly reviewed. First, we define the notion of Positive-Operator Valued Measures, which are useful tools for modeling measurement noise in quantum devices. Second, we describe a Quantum Detector Tomography procedure, which allows to reconstruct a POVM associated to a given device. Then we discuss a measure of distance between quantum measurements known as the operational distance. Finally, we use the formalism of POVMs to precisely define classical noise affecting measurements, which will be used in Section III to derive our correction procedure.

A. Mathematical description of quantum measurements
We start by a mathematical description of a generalized quantum measurement modeled by a Positive-Operator Valued Measure [25]. A POVM with n outcomes on a d-dimensional Hilbert space may be described by an nelement vector M of operators M i (called effects), such that where M i 's are represented by d × d semi-definite positive matrices, 1 is the operator identity and T denotes transposition which we use because the column form of POVMs will be useful later. If a quantum system was initially described by a state ρ, the probability of obtaining the outcome associated with the effect M i is given by Born's rule In other words, probabilities of obtaining particular outcomes are equal to the expectation values of the associated effects. In quantum information protocols, the considered measurements are often projective [6], which means that the effects M i fulfill the additional requirements In the experimental part of this work, we will focus on projective measurements and their noisy versions, since this fits the set-up of the IBM Q devices. Nevertheless, all our theoretical considerations, including the correction procedure, are formulated for generalized measurements given by Eq. (1).

B. Quantum Detector Tomography
Characterization of quantum devices requires knowledge of the POVM associated with a given measurement performed by the device. The procedure for obtaining this is known as Quantum Detector Tomography (QDT) [16]. The general idea of QDT is as follows. Recall that Born's rule associates the probabilities of particular outcomes with the expectation values of effects. If one performs multiple experiments on the set of quantum states which form a basis for Hermitian operators, then one may use the obtained statistics to decompose all effects in that basis via Born's rule. In the case of a d-dimensional Hilbert space, such basis is d 2 dimensional, and this is the minimal number of different experiments (different quantum states) that must be performed to reconstruct a POVM.
Naturally, in practice one will not have access to perfect statistics, nor to the perfectly prepared quantum states. This may cause that the simple tomographic method returns nonphysical, non-positive POVM elements. To cope with this issue, the method of Maximum Likelihood Estimation (MLE) may be used, which asserts positivity of the reconstructed quantum measurements [26,27]. In order to reconstruct POVMs from experimental data, we have thus implemented an iterative algorithm from [26] which converges to the MLE estimator.

C. Distances
After reconstructing the POVMs, it is useful to have a distance measure between them in order to compare the reconstructed quantum measurements with the ideal detectors and, naturally, to compare detectors with each other. Distances between quantum measurements are usually closely related to the probability distributions that they generate via Born's rule. Therefore, first we need a measure of distance between classical probability distributions. In this work, we will be using so called Total-Variation (TV) distance, which for two arbitrary probability vectors p = (p 1 , . . . , p n ) T and q = (q 1 , . . . , q n ) T , TV distance between them is defined as where ||.|| 1 denotes l 1 norm. Total-variation is a very stringent measure of distance between probability distributions. It plays important role in proposals for attaining quantum computational advantage [28,29] and can be also easily utilized to control expectation values of (classical) observables defined on the sample space of interest.
TV-distance is related to the operational distance [30][31][32] between arbitrary POVMs M and N via the equality where maximization goes over all quantum states ρ, and p M/N denote probability distributions generated by measurement of quantum state via M/N, i.e., Hence, the operational distance is in fact the worst-case scenario of the distance between probability distributions generated by performing those measurements. It can be shown [30] that operational distance may be calculated as where the maximization is over all subsets of indices enumerating effects (i.e., all possible sets of outcomes) and ||.|| ∞ denotes operator norm [33]. Operational distance between POVMs [34] has an interesting operational interpretation through the formula D op (M, N) = 2p disc (M, N) − 1, where p disc (M, N) is the optimal probability of distinguishing between measurements M and N (without using entanglement) [31].

D. Classical noise affecting measurements
Now we will describe model of classical noise that we are considering in this work. Let us denote by M ideal an noutcome quantum measurement that in theory should be associated with our measurement device. In practice, due to the presence of noise, real measurement describing our device is some POVM M exp . In this model, we assume that the relation between M ideal and M exp is given by a (left) stochastic, invertible map Λ, whose element Λ i,j := p (i|j) ∈ [0, 1] are defined by equation The left-stochastic property of Λ means that i Λ i,j = i p (i|j) = 1. Note that this property asserts that M exp is a proper POVM. The above equation is somewhat abstract, since it provides a description of the noise on the measurement operators level. In order to find its operational interpretation, let us use the fact that since the noise affects only a measurement process, Eq. (7) is fulfilled for arbitrary quantum state which is measured on the noisy device. For an arbitrary quantum state ρ, let us denote the ideal vector of probabilities that one would've obtained in the ideal device as p ideal , which should be understood as in Eq (5). Due to the linearity of Born's rule, from Eq. (7) it follows that the vector of probabilities p exp obtained in an experiment on a noisy device is given by This allows us to give an intuitive interpretation of classical noise. Namely, the application of such noise is equivalent to the classical post-processing of the statistics that one would have obtained in an ideal scenario [35]. Such interpretation has been used previously in the context of implementating of general POVMs by projective measurements (see works [22,36,37]).

III. SCHEME OF MITIGATION OF READOUT ERRORS
In this part we lay out the idea of our procedure to mitigate readout errors. We first state the assumptions under which the method works perfectly. Then, we formally formulate the method itself. Finally, we compute operational distances (see Eq. (6)) between ideal and noisy projective measurements for (i) single-qubit classical readout errors and (ii) uncorrelated classical errors affecting multi-qubit projective measurements. These numbers give us the indication on the magnitude errors that can be corrected by our error mitigation technique in realistic devices.

A. Assumptions
Our error-mitigation scheme relies on the following assumptions: 1. (Infinite statistics) We have access to the statistics given by Born's rule (Eq. (2)).

(Classical noise)
The measurement noise occurring in the device has the form of a classical noise described in Section II (Eq. (8)).

(Characterized detector)
We have at our disposal a perfect description of that noise.
The approximate validity of these conditions can be motivated by the following arguments. Assumption of infinite statistics may be fulfilled to a good extent simply by increasing the number of experiments one gathers statistics from. Note that violations of this assumption introduces statistical errors that can also be taken into account.
Moreover, in Section VI (where we present results of QDT for IBM and Rigetti devices) we show that the classical noise model stated in second condition turns out to be a dominant type of noise in the IBM quantum devices, which uses superconducting transmon qubits. This suggests that classical noise may be a good measurement noise model for devices relying on similar architectures.
Finally, the characterization of noise occurring in a detector may in practice be obtained via Quantum Detector Tomography (see Section II). Such reconstruction is of good quality, provided one has access to high fidelity preparations of single-qubit quantum states. For example, the fidelities [38] of single qubit gates in IBM quantum devices are typically high (of the order of 99.6% -see Section VI).
We conclude by noting that the violations of the assumptions affect the applicability of our error-mitigation procedure. The effects of such errors are analyzed in the Section IV.

B. Correction of the statistics
Let us now describe the error-mitigation procedure. In fact, all we need is to notice that since the matrix Λ representing the noise is invertible, we may simply left-multiply the Eq. (8) by its inverse Λ −1 , obtaining What Eq. (9) tells us, is that provided the knowledge of Λ, one can simply left-mulitply p exp by its inverse, in order to obtain statistics p ideal that one would have obtained in the experiments performed on the ideal, non-noisy device.
This result indeed appears very natural if one recalls the fact that such classical noise is equivalent to the classical post-processing of ideal statistics (see Section II).

C. Correction of statistics based on Quantum Detector Tomography in (multi-)qubit devices
So far, we have not assumed any particular structure of the measurement M ideal . In this subsection (and in fact throughout the most part of the paper), we focus on projective measurements in the computational basis in singlequbit and multi-qubit systems. Specifically, we describe in detail how our mitigation method works in these systems. Moreover, we analytically quantify the magnitude of statistical errors that our scheme is capable of correcting under the assumption of independent classical errors affecting the measurements.

Single qubit
For the case of single-qubit projective measurement, without looss of generality we may write effects of a perfect detector POVM M ideal ≡ P in the computational basis as Let us assume that the measurement that is actually implemented is of the form M exp ≡ M = ΛP, for some invertible stochastic transformation Λ. This means that effects of M can be written as where p, q ∈ [0, 1] are probabilities of erroneous detection for outcomes 1 and 2, respectively. It is now straightforward to verify that The operational distance Dop (M, P) between the ideal single-qubit projective measurement P and the measurement with classical noise M plotted against the probability q of erroneous detection. The figure presents the regime in which q > p (see Eq. (11)), which is the typical situation in experiments. The stars (triangles) correspond to the experimental values of q obtained in the detector tomography experiments on IBM's (Rigetti's) device.

From Λ we can construct the correction matrix
which can be used in the future to correct any statistics obtained in one-qubit experiments performed on this device. The adjustment simply amounts to multiplying column vector of statistics from the left by Λ −1 . Finally, to illustrate how error probabilities p and q affect the distance of M from the ideal detector P, we compute the operational distance The plot of this dependence is presented in Fig. 2.

Multiple qubits
Considerations from previous subsection may be easily generalized for multiple-qubit systems. In general, detector errors may exhibit correlations between qubits, therefore QDT should be performed via simultaneous measurements on multiple qubits. The problem of tomographic reconstruction is then exponential in the number of qubits K, since it requires preparation of exponentially many quantum states. However, QDT can be performed efficiently if the readout erros affecting the qubits are uncorrelated, i.e., where Λ i is the stochastic matrix describing the readout noise on the ith qubit. If that is the case, it suffices to perform multiple single-qubit QDTs, which makes the problem linear in the number of qubits. Interestingly, the operational distance between the noisy detector M (K) = Λ (K) P (K) and a multi-qubit projective measurement P (K) = ⊗ K i=1 P i can be computed analytically (see Appendix A for a detailed derivation) is the operational distance between ideal and noisy detector on the ith qubit.
Importantly for small individual errors (i.e.
This motivates the usage of our error-mitigation procedure in this regime.

IV. ERROR ANALYSIS
So far we have considered idealized scenarios in which noise affecting the measurements was perfectly classical and we could repeat experiments infinitely many times. Now we will provide an analysis of the possible deviations from the ideal model for the problem of reconstruction of measurement statistics. Specifically, we will describe three sources of errors: 1. Not entirely classical noise -it might happen that the noise is not of purely classical nature (as described in Sec. II), but also have some non-classical unitary rotation part.
2. We have only access to a finite number of experiments, which introduces the statistical noise when we reconstruct probability distributions.
3. Finding the closest probability vector in the case when the correction yields non-physical probability distribution.
We note that, we omit here errors resulting from the imperfect tomography of measurements, which we shall address in the future work [39].

A. Non-classical noise
Typically noise affecting quantum measurements cannot be described by Eq. (7). Let us consider the situation, in which the measurement that is actually being implemented is of the form where we have decomposed a POVM into a part ΛM ideal that represents ideal POVM affected by the classical noise (as in Eq. (18)) and ∆ which represents every other (nonclassical) errors in reconstructed POVM. Of course, having access to M exp (e.g., from QDT) and M ideal (from the theoretical model), decomposition presented in Eq. (18) may be done in arbitrary way. However, in the case of projective d-outcome measurement in the computational basis, the following anzatz seems natural. Namely, we propose to consider a diagonal part of the POVM as consisting information about Λ, while all off-diagonal terms should be regarded as non-classical part of the noise ∆ (see also Remark 1). Clearly, non-zero ∆ affects our error-mitigation procedure. In the case of a noisy POVM of the form given in Eq. (18), it is impossible to reverse the effects of such noise solely by classical post-processing. Computing the vector of statistics for arbitrary quantum state ρ with the help of Eq. (18) yields (in analogy to Eq. (8)) where∆ denotes a generic disturbance of experimental statistics which arises due to the presence of non-classical part of the noise ∆. Its elements may be defined as ∆ i := Tr (ρ∆ i ), where ρ is an unknown quantum state and ∆ i corresponds to the non-classical part of ith effect. Applying Λ −1 to Eq. (19) gives which clearly does not leave us with the ideal statistics, but consists also some additional term Λ −1∆ . However, if non-classical part ∆ is small compared to the term ΛM (ideal) , we propose to neglect the non-classical part ∆ and perform our error-mitigation procedure as if there were only classical noise. Such action will clearly introduce some error into resulting estimated probability vectors. We quantitatively characterize the error introduced by neglecting non-classical part of the noise by finding an upper bound on D T V Λ −1 p exp , p ideal in terms of the operational distance between the ideal measurement M exp and ΛM ideal . Using Eq. (19) we obtain where in the inequality in the second line we have used the standard inequality ||Ax|| 1 ≤ ||A|| 1→1 ||x|| 1 valid for all vectors in R n and linear transformations [40] A : R n → R n .
Performing the optimization over all quantum states (note that∆ implicitly depends on a quantum state ρ), and using Eq.(4) finally yields the bound We see that upper bound for the error that is introduced by our error-mitigation procedure, due to neglecting the nonclassical part of the noise, can be expressed by 1 → 1 norm of Λ −1 (i.e. maximal l 1 norm of its columns) and by the operational distance between the actually implemented measurement M exp (given by Eq. (18)), and the part of POVM which is affected only by classical noise , i.e., ΛM ideal . Remark 1. The above considerations were based on the decomposition of the imperfect measurements on the 'classical' (ΛM) and the 'non-classical' (∆) part. In general, such decomposition may be done in somewhat arbitrary way. In principle, one could define an optimization problem, which should minimize, e.g., some suitably-chosen cost function depending on ||Λ −1 || 1→1 and ||∆|| 1 so aiming to minimize the upper bound on the error of the error-mitigation procedure.

B. Statistical errors
We will now study how finite-statistics errors affect the performance of our method. First of all, let us note that in experiments we do not have access to the probability vector p exp which we used previously, but only to its estimator. In what follows we will focus on the natural Maximum-Likelihood (ML) estimator for which estimates of individual probabilities p i are given by relative frequencies n i /N of events observed in the experiment repeated multiple times. Therefore, we shall now call the experimental statistics vector p est exp , putting emphasis on the fact that it is an estimator.
The quality of the estimation procedure can be quantified by Pr err ( ) := Pr D T V p est exp , p exp ≥ , i.e., the probability that the TV distance between the estimated statistics vector p est exp and the true vector p exp is greater than some threshold . Hence, the value of 1 − Pr err ( ) may be interpreted as the acceptable confidence level for our estimation. Let us denote by N the number of samples (experimental runs) and by n the number of outcomes of the considered measurement (since we mainly focus on the standard projective measurements, n coincides with the dimension of the considered Hilbert space). In work [41] authors proved that in this scenario with probability at least 1 − Pr err The above inequality gives us the upper bound for the TVdistance between estimated statistics (frequencies) vector p est exp and actual probability vector p exp as a function of the number of experimental runs and the accepted error probability.
Let us finally stress that we are not interested in the statistical error itself as it is inherently present in any estimation scheme, whenever we have on our disposal only a finite number of samples. Instead, we want to understand how finite-sample statistics affects our procedure. Specifically, we want to bound the TV distance of our corrected estimated statistics Λ −1 p est exp to the ideal probability vector p ideal . This can be easily done in a manner analogous to the derivation leading to Eq. (22). The final bound is the following (24) where the first inequality is a consequence of the triangle inequality, and the derivation of the second inequality is fully analogous to that of Eq. (22). In the upper bound above, the first term in the sum accounts for the statistical errors while the second term accounts for the non-classical part of the noise.

C. Non-physical probability vectors
Due to the coherent and statistical errors described above, it may happen that the corrected vector of the estimated statistics Λ −1 p est exp will not be a proper probability vector. Generally, the vector obtained in our procedure may contain negative elements that will however sum up to 1 (this is because the inverse of left-stochastic matrix still has columns of which elements sum up to 1). When this is the case, we propose to find the corrected vector p * exp which is the closest probability vector to Λ −1 p est exp in Euclidean norm [42]. In other words, if obtained vector Λ −1 p exp is not a proper probability vector, we use the corrected vector p * exp that is the solution to the following problem where || · || 2 denotes the Eucledian norm. The above problem can be easily solved by convex optimization solvers like cvxopt [43]. Furthermore, the error α introduced by this method may be quantified in terms of TV-distance between new corrected vector p * exp and nonphysical Λ −1 p est In order to account for the overall error of our procedure, we modify the Eq. (24) by putting p * exp in place of Λ −1 p est exp . The usage of the triangle inequality gives the following upper bound for the overall error in our procedure where δ and α are defined in Eq. (24) and Eq. (26), respectively. The quantity δ + α can be considered as the overall error that our procedure yields for the problem of estimation of the probability vector p ideal . This upper bound on D T V p * exp , p ideal is a function of four variables: the number of experimental runs N used in the procedure of probability estimation, accepted error probability Pr err , the reconstructed measurement M exp , and the vector of statistics obtained in a given experiment p est exp .

D. When is the mitigation successful?
We would like to address now a crucial question. Namely, when can we consider our error-mitigation procedure to be successful? The answer to this can be given by analyzing the the upper bound δ + α on the error of estimation of probability vector that can be potentially obtained by the error-mitigation scheme (resulting from the presence of coherent errors and finite-statistics errors).
Imagine that we have on our disposal only the noisy quantum measurement M exp that produces N samples from the distribution p exp . Let p est exp denote the empirical estimator of this probability vector. By the virtue of considerations given in earlier subsections, the error resulting from imperfect measurement and the finite number of samples (for the assumed error probability Pr err ) can be upper bounded as follows We then propose to consider our mitigation successful if the upper bound in Eq. (27), on the error D T V p * exp , p ideal that can be introduced by the error-mitigation procedure, is smaller then the upper bound on D T V p est exp , p ideal discussed above (for the same number of samples N and error probability Pr err ). In other words we propose the following rule (29) Remark 2. Of course, the actual distance from the perfect probability distribution after post-processing, is highly dependent on the quantum state that is measured. Importantly, it might happen that though inequality in (29) holds, the particular quantum state will give statistics that after post-processing are worse than without it in terms of the distance from the perfect probability distribution. We analyze this problem in some detail in Appendix B where we provide more support of the heuristic procedure presented above.
Let us illustrate the above ideas on the simple example of single-qubit projective measurement. We assume that in detector tomography, noisy POVM describing our device is given by Unlike the case of purely classical noise discussed in Section III, there are also off-diagonal terms which represent non-classical part of the noise. If the their magnitude is much smaller than that of the diagonal elements, we may approximate our POVM by ΛP (where the appropiate stochastic matrix is given in Eq. (12)) at the cost of introducing some additional errors. According to the analysis presented earlier, in order to bound this error we need to compute ||Λ −1 || 1→1 and D op (M, ΛP). Straightforward computation give us the following expressions Inserting these expressions into Eq. (24) we obtain the following upper bound on the error of estimation of the probability distribution with the usage of our error mitigation scheme where statistical error is determined by Eq. (23) and depends on the assumed level of confidence and the number of samples N used in the estimation procedure. For example, if we set N = 8192, which is maximal number of experimental runs for a single quantum circuit in IBM Q devices, and Pr err = 0.01, we obtain ≈ 0.018. Importantly in Eq. (32) we have assumed that the estimated probability vector Λ −1 p est exp is a probability distribution add hence we could set α = 0 in Eq. (26). The visualization of overall error quantifier δ is shown in Fig. 3.

V. SUMMARY OF THE METHOD
In the following, we aim to correct the output statistics obtained from a noisy measurement. The error-mitigation procedure applied to projective measurements is realized by performing the following sequence of practical steps, the graphical illustration of which is given in Fig 4. 1. Perform Quantum Detector Tomography (QDT) experiments. The estimation of the detector's POVM M can be done using, e.g., an algorithm which converges to the Maximum Likelihood Estimation, like the one of [27].
2. Neglect all the off-diagonal terms of the elements of M, obtaining ΛP, where P is the ideal projective measurement we want to implement.
4. After performing arbitrary experiment on the device characterized in the above way, multiply the vector of obtained frequencies p est exp by Λ −1 .
5. Check if the Λ −1 p exp is a proper probability vector, i.e., its elements are positive and sum up to 1.
• If yes, set p * exp = Λ −1 p est exp and proceed to the next step.
• If no, solve the problem formulated in Eq (25), obtaining a vector of corrected statistics p * exp . 6. Calculate δ + α, i.e., the upper bound on the total error magnitude D T V p * exp , p ideal ). Compare this with D op (P, M) + , i.e., with the upper bound on the error occurring without correction (Eq. 28).
Remark 3. We note that if we wanted to mitigate the readout errors in the case of arbitrary ideal measurement M ideal , we would have to modify the second step of our method. Specifically, we would have to propose a decomposition of the noise on the classical (ΛM ideal ) and nonclassical (∆) part (see discussion in Remark 1).

VI. DEVICE CHARACTERIZATION
In this section, we present experimental results [44] that confirm the (approximate) validity of the physical assumptions stated in Section III in the IBM's device. These assumptions have to be satisfied for our error-mitigation procedure to work. We first present the randomized benchmarking data for IBM's devices, which shows that singlequbit gates errors are small and makes Quantum Detector Tomography feasible. Then the results of QDTs performed on IBM devices are provided. This shows that readout errors in these systems are indeed mostly of classical nature (see Section II). We observe that the readout errors between physically connected pairs of qubits are mostly uncorrelated (see Fig. 5 for visual presentation of the connectivity of ibmqx4 device). Furthermore, through out the section, we present analogous data for exemplary five qubits in Rigetti's 16-qubit device Aspen-4-16Q-A. However, in the case of this device, the assumption of perfect detector reconstruction is significantly violated due to the high infidelities of the single-qubit quantum gates. Remark 4. The experimental results presented in this work were obtained on publicly available devices. For this reason, it was not possible to perform all the experiments during a single calibration period. Consequently the data presented throughout the paper may consist of values obtained in experiments performed on different days. Importantly, to correct statistics in all experiment performed in a given calibration period, we have used the data from QDT specific to that period. We note that results of QDT were varying across different calibration periods. However those fluctuations do not change qualitative, nor quantitative conclusions derived throughout our work.

A. Perfect state preparation
In order to be able to perform reliable Quantum Detector Tomography, one needs to assure accurate state preparation. Since tensor product of single-qubit quantum gates suffice to prepare quantum states spanning the whole operator space, it is enough to analyze only single-qubit gate errors. Randomized benchmarking experiments [46,47] allow for inferring such errors in a manner approximately independent on state-preparation and measurement errors. In Table I we list average gate fidelities F RB for ibmqx4 device. The errors are relatively small (i.e., 1 − F RB is of order 0.1%) which indicates that the assumption of perfect state preparation is reasonable. The average gate fidelities in Rigetti's device are, however, worse, which makes it difficult to precisely gauge the errors introduced by our error-mitigation procedure.
We would like to remark that in general the problem of rigorously separating measurement, gate, and statepreparation errors is a difficult task (for recent progress in this direction, see [48][49][50]). Without performing this division, we cannot precisely estimate the accuracy of our error-

FIG. 7. Visualization of single-qubit detectors for IBM's ibmqx4 (stars) and Rigetti's Aspen-4-16Q-A (triangles). We parametrize
the first effect of a measurement as M1 = k∈{0,x,y,z} n k σ k , for σ0 := 1. We set |z| := n 2 x + n 2 y to be the magnitude of coherent errors. In the figures, we plot the two-dimensional projections of the four-dimensional parameter space of M1: part a) depicts the identity (n0) and the σz (nz) components, whereas part b) shows the identity (n0) and the coherent (|z|) components. mitigation procedure. Let us emphasize, however, that for many practical applications (like variational quantum algorithms for example [51]) one can simply implement our scheme and test if it works in practice.

B. Quantum Detector Tomography
Here we present results of detector tomographies performed on all 5 qubits of IBM's five-qubit device ibmqx4 and on first 5 qubits of Rigetti's 16-qubit device Aspen-4-16Q-A.

Single qubit
On Fig. 6 we present the results of QDT performed on individual qubits in IBM's and Rigetti's devices. Three observations can be made. First, readout noise is significant for both platforms. We note that although it cannot be directly compared to the single-qubit gate errors from RB experiments, in the case of IBM ibmqx4 average gate fidelities are so large that the readout-noise can be considered as the predominant type of noise on the level of individual qubits (at least for short quantum circuits). Second, the distance between actual POVM and its diagonal form is in all cases small, which confirms the assertion that clas-    (22)) is always significantly smaller than the distance of the noisy POVM from the ideal detector, which suggests that the mitigation scheme should be beneficial [52]. We also visualize the results of QDT of single-qubit projective measurements on Fig 7. To this end , we parametrize the first effect of a measurement as M 1 = k∈{0,x,y,z} n k σ k , and σ 0 := 1. Then, we set |z| := n 2 x + n 2 y to be the magnitude of coherent errors. In order for M 1 to remain element of POVM, its coefficients must satisfy n 2 0 − 1 ≤ |z| 2 + n 2 z ≤ n 2 0 .

Two qubits
We have also performed QDT for two qubit projective measurements on Rigetti's and IBM's devices. In Figs. 8a, c we present results of detector characterization performed on the basis of joint two-qubit measurement tomography. We also give results of the alternative approach in which two qubit POVMs are reconstructed from tensor products of single-qubit POVMs (see Figs. 8b, d) . In both cases, the data resembles single-qubit scenario -coherent errors are relatively small and correction er- ror ||Λ −1 || 1→1 D op M exp , ΛM ideal (relevant for the infinite statistics scenario) is always much smaller than operational distance of the POVM from ideal detector.
We have also investigated the correlations between measurement errors. Tab. II shows operational distances between POVMs reconstructed via two-qubit tomography and these obtained via tensor product of single-qubit measurements. We observe that for most pairs in IBM's device, the correlations between readout errors are small. However, with two exceptions of the pairs q 3 q 1 and q 2 q 1 . This suggests that mitigation for those pairs should be performed based on joint two-qubit QDT. In the case of Rigetti's device, examined pairs do not show significant level of correlations.

VII. APPLICATIONS ON IBM DEVICES
We have performed numerous experiments on IBM ibqmx4 quantum device to demonstrate the usefulness of our error mitigation method for a number of quantum information tasks. In what follows we first describe briefly the theoretical basis of each of the tasks performed. Then, we proceed to the detailed description of experimental results.

Quantum State Tomography and Quantum Process Tomography
Both Quantum State Tomography (QST) and Quantum Process Tomography (QPT) are based on the same idea as Quantum Detector Tomography, with the clear difference that in the case of QST a quantum state is to be reconstructed, while a QPT provides the characterization of a quantum process, i.e., a channels. In our experiments, the tomographic reconstruction was done using the algorithm from [53], which is available in qiskit.
To compare the target pure state |ψ ψ| with the tomographically reconstructed ρ we used state fidelity [6], F state (|ψ ψ|, ρ) = Tr (|ψ ψ|ρ) . (33) Likewise, in compare the target unitary channel U to-beideally-implemented with the actual quantum channel Λ obtained in the process tomography, we used the entanglement fidelity [54,55], which may be calculated as where Φ U is the Choi-matrix representation of the unitary channel and ρ Λ is the Choi-matrix representation of Λ.

Implementation of non-projective measurements
In our work, we have used generalized quantum measurements as models for noisy projective measurements. However, it is certainly not the only application of POVMs in quantum information. Indeed, POVMs can outperform projective measurements in some tasks, such as unambiguous state discrimination [56,57], minimal error state discrimination [56,57], quantum tomography [58], port-based teleportation of quantum states [59] or quantum computing [60]. Therefore, it may happen that one actually wants to implement a certain generalized measurement.
The standard method to do so is via the Naimark's extension [25,61] and requires adding ancilla system and performing a projective measurement on a whole system. We have implemented in this way three different single-qubit POVMs. We have assessed the quality of implementation by performing quantum measurement tomography and computing the operational distance, see Eq. (6).

Quantum algorithms -Grover's search and the Bernstein-Vazirani algorithm
Grover's [23] and Bernstein-Vazirani [24] algorithms are the canonical examples of quantum algorithms operating withing a oracular (black box) model of computation. Grover's search aims to find unknown element y encoded in the application of the unitary gate U y defined via U y |x = (1 − 2δ x,y )|x , where δ x,y is the Kronecker's delta. BV algorithm algorithm uses Hadamard transform to find (in a single query) a hidden string s ∈ Z n 2 encoded in a n-qubit unitary transformation V s defined via We have based our implementation of these routines on the expository work [62]. Both algorithms were implemented on three qubits, one of which was an ancilla. In that case, Grover's algorithm required only a single query to the oracle, hence it could have been realized in a single quantum circuit. The figure of merit for the quality of implementation in the case of both algorithms is a single number -the probability of obtaining a particular outcome. We have used this number to benchmark our error mitigation scheme.

Probability distributions
We have implemented various five qubit circuits in IBM device aiming to generate specific probability distributions upon measuring all the qubits. As a figure of mearit we have used Total Variation distance (see Eq. (3)) between target probability distribution, and the one estimated from relative frequencies. In order to test our error mitigation method we have implemented these probability distributions with and without the error-mitigation procedure.

B. Experimental results
Now we present experimental results which demonstrate the practical effectiveness of our error-mitigation procedure. For each presented experiment, QDT from which we inferred correction matrix Λ −1 , has been performed in the same calibration period as the corrected experiments. Moreover, for tasks involving measuring multiple-qubits we compare mitigation which assumes lack of correlation between qubit readout errors, with the one which accounts for such correlations.

Quantum State Tomography and Quantum Process Tomography
In what follows we present results of quantum state and quantum process tomographies performed on single-qubit systems and a quantum state tomographies performed on two qubit systems. The procedure of tomography involves performing multiple experiments that are then used to reconstruct the objects in question. We use error-mitigation procedure to correct the statistics (probability vectors) in these experiments in order to enhance the quality of the tomographic reconstruction.
On Fig. 9 we present results for different quantum state and process tomographies on individual qubits in IBM's ibmqx4 device. It is worth noting that although results highly depend on the input state, our correction reduces infidelities of reconstructed states in every tested case.
We observe the improvement also for two-qubit state tomographies, as shown in Fig. 10. We compare the performance of the error-mitigation procedures based on the measurement reconstruction obtained from two-qubit QDTs and the one obtained from the tensor product of singlequbit QDTs. Interestingly, although one would expect that for highly correlated readout errors (in particular the pair q 2 q 1 ), accounting for those correlations in readout errors should provide improvement over the mitigation based on Infidelities obtained from the reconstruction of a) single-qubit states and b) single-qubit channels. Different colors refer to different qubits. The first bar in each pair corresponds to uncorrected statistics, while the second bar corresponds to the corrected one. The three regions separated by dashed lines refer to different Haar-random [63] quantum states (processes) to be reconstructed. Each tomography experiment consisted of several quantum circuits (3 for QST, 12 for QPT), and each quantum circuit was implemented 8192 times. Furthermore, each tomography experiment was repeated 5 times in order to estimate statistical errors, the corresponding 3σ intervals are shown on the barplot as black lines. The correction was based on a QDT which used an over-complete set of states (all Pauli eigenstates), with each circuit implemented 32768 times.
not-correlated tomography, we observed such improvement in the case of only one out of three tested quantum states.
We remark on the possible causes of the result reported above. First, in the course of correction of statistics we sometimes obtain nonphysical probability vectors. As already pointed out in Section IV C dealing with such situations potentially introduce additional errors. These errors can be different for mitigation schemes based on different tomographic reconstructions. The second probable cause is the insufficient number of experiments necessary to perform a reliable QDT of a two-qubit measurement (the sample complexity of this problem is definitely higher as there are are much more parameters needed to describe the two-qubit measurement compared to two single-qubit measurements).
Infidelities obtained from the tomographic reconstruction of Haar-random [63] two-qubit quantum states using error-mitigation procedures based on a) uncorrelated singlequbit QDTs and b) correlated two-qubit QDT. The convention for presenting data is analogous to that of Fig. 9. Each tomography experiment consisted of 9 quantum circuits, corresponding to local measurements performed in different Pauli bases. Each quantum circuit was implemented 8192 times. Furthermore, each tomography experiment was repeated 5 times in order to estimate statistical errors, the corresponding 3σ intervals are shown on the barplot as black lines. The QDT used in the error-mitigation procedure was implemented using an overcomplete set of states (all Pauli eigenstates), with each circuits implemented 32768 times.

Implementation of non-projective measurements
For implementation of non-projective measurements via Naimark's extension, we observe the general improvement of the quality of reconstructed POVMs when errormitigation scheme is used (see Fig. 11). Moreover, the difference between error-mitigation procedure based on the non-correlated and on correlated detector tomography are negligible this time. However, we can offer a simpler explanation of than in the case of two-qubit QSTs described above. Namely, at the time when those particular experiments were performed, correlations between errors for pair q 2 q 1 were relatively small (compared to data regarding quantum process tomography and quantum detector tomography). This points at the importance of performing systematic device characterization and calibration. ; while the second bar shows the distance between the ideal POVM and the one reconstructed from a tomography based on statistics corrected by our scheme (D corr op ). The three regions separated by dashed lines correspond to different to-be-implemented POVMs: the Haar-random 4-outcome POVM, the tetrahedral measurement, and the trine measurement [63]. Each POVM reconstruction required 4 quantum circuits (Pauli 'x+','y+','z+' and 'z-' states) and each quantum circuit was implemented 8192 times. Reconstructions were repeated 5 times in order to estimate standard deviations, the corresponding 3σ intervals are shown here as black bars. The mitigation or measurement errors by our method was done based on a) separate single-qubit QDTs, and b) joint two-qubit QDT. The QDTs were based on the preparation of Pauli 'x+','y+','z+' and 'z-' states for single qubits, and all of their tensor-product combinations for pairs of qubits. Each circuit for the QDT was implemented 32768 times.

Quantum algorithms -Grover's search and the Bernstein-Vazirani algorithm
From the data presented in Table III we observe that our error mitigation method is favorable for implementation of both Grover and BV algorithms. In the case of correction performed on the highly correlated pair q 2 q 1 (Grover's search), the mitigation scheme accounting for the correlations of readout errors provided the significant advantage TABLE III. The guess probabilities for the three-qubit Grover's search and BV algorithms implemented on ibmqx4. The third column shows the mitigation based on uncorrelated single-qubit QDTs, and the fourth column shows a mitigation where the correlations between readout errors were accounted for. For both algorithms, we used the qubits q0, q1, and q2. The measurement was done on the q2q1 pair for Grover's search and on the pair q1q0 for the BV algorithm. For Grover's search the hidden element was chosen to be '11', while for BV it was '01'. Each experiment was repeated 5 times in order to estimate standard deviations, the corresponding 3σ intervals are given in the table. The mitigation was based on a minimal QDT (Pauli 'x+','y+','z+' and 'z-' states), with each circuit implemented 32768 times. TABLE IV. Errors for five-qubit measurements in ibmqx4. The first row presents the lower bound [64] for the operational distances between noisy and ideal measurements including statistical errors when the number of repetitionsis N = 72 × 8192 = 589824 and accepted error probability Prerr = 0.01 (see Eq. (23)). The second row gives the upper bound on δ defined in Eq. (24), which bounds the errors in the observed statistics. The second column corresponds to a tomographic reconstruction based on the tensor product of five single-qubit QDTs, while in the third column the joint QDT for a pair q2q1 was included. The presented data is for POVMs which were used for the correction procedure, hence we do not provide standard deviations. over over the scheme using non-correlated tomography.
We remark that that relative improvement is much higher in the case of Grover's algorithm than in BV's. This may result from two factors. First, the circuit for Grover's algorithm has lower depth than the circuit for the BV algorithm, therefore in this case gate-errors might have smaller impact on the final result. Second, for the implemented instance of Grover's algorithm, the expected result was '11' (which is the state most prone to the noise observed in IBM's device), while for the BV algorithm the theoretically correct result was '01'.
Remark 5. We note that in Ref. [62], authors obtained values ≈ 0.65 and ≈ 0.386 for the correct guessing probabilities in the same instances of Grover's and BV algorithms. Both of these values differ significantly from the ones observed by us without implementing the correction procedure. Such discrepancy is another confirmation that indeed device's performance varies over long periods of time.
TABLE V. Results of experiments for implementing specific probability distributions on the ibmqx4 processor. We use the total variation distance between the target and the reconstructed probability distribution in order to asses the quality of the implementation. The upper table presents results of for the uncorrected case and the case when applying a correction procedure that does not take into account the correlations in the readout errors, while the lower table provides the results of a correction procedure where the correlations between a single qubit pair (q2q1) measurements were taken into account. Each experiment consisted of a single quantum circuit, which was run 72 * 8192 = 589824 times. This high number of repetitions was necessary in order to minimize the statistical errors. Furthermore, each experiment was repeated four times in order to estimate the statistical errors. The errors given in the tables have magnitude 3σ, where σ is the estimated standard deviation. The correction procedure used in the experiments was based on over-complete (all Pauli eigenstates) single-qubit QDTs with each circuit implemented 32768 times. The twoqubit QDT for pair q2q1 was done using all tensor-product combinations of Pauli eigenstates on both qubits, with each circuit implemented the same number of times.

Probability distributions
Finally, we present results on the implementation of probability distributions on five qubits. We start by giving the operational distances and bounds on the correction errors for the five-qubit detector in Table IV. From the data presented there, it is clear that even for 5 qubits we are still in the regime in which correction are beneficial (see the criterion given in Eq. (29)).
The probability distributions tested here were generated by performing a simple circuit (acting on the standard input state |0 ⊗5 ) followed by a measurement in the computational basis. We have chosen to implement the following three probability distributions: The first distribution was the uniform one and was implemented by the tensor product of five Hadamard gates. The second probability distribution, named 'NOT ', was implemented by the product of five X gates. For the last probability distribution, named by us 'Mixed', the circuit consisted of two X gates on q 0 and q 2 , two Hadamards on q 3 and q 4 , and the identity gate on q 1 .
Results of the correction of the probability distributions are presented in Tab. V. There are a few interesting features to be observed. First, let us note that in the case of the uniform probability distribution, our correction procedure is ineffective if we do not account for correlations be-tween qubit readout errors. However, if we do account for such correlations in only a single pair of qubits (the highly correlated pair q 2 q 1 ), it is clear that correction procedure provides a high improvement [65]. Second, in the case of the 'NOT' probability distribution, the correlation-ignorant correction has attained the perfect estimation -yet with remarkably high distance α between first-guess non-physical probability vector and its closest physical neighbor. This moves us very close to the regime in which the correction is unreliable according to the criterion given in Eq. (29). Accounting for correlations this time proves to improve results also remarkably highly, this time with only slightly distant first-guess non-physical distribution. Finally, the method provides advantage for the 'Mixed' distribution in both cases. However, interestingly, this time the difference of improvement after accounting for correlations is negligible.

A. Conclusions
In this work, we have presented a scheme for the mitigation of readout errors which is suitable for noisy and imperfect quantum devices. Our method uses solely classical post-processing that depends on the results obtained in the course of Quantum Detector Tomography. In principle, it can be applied on any quantum hardware provided the readout errors are classical and one has access to approximately perfect detector tomography. We have analyzed how our error-mitigation procedure is affected by nonclassical noise and finite statistics (at the stage of estimation of probability distributions). We have also presented a comprehensive study confirming the (approximate) validity of the classicality of errors in the publicly available prototypes of quantum chips based on architectures with superconducting transmon qubits.
Our method was tested on a variety of quantum tasks and algorithms implemented on IBM's five-qubit quantum processor. We observed that our scheme, despite relying on a purely classical processing of the obtained results, yields substantial improvements in the performance for a number of in single-, two-and five-qubit experiments on the ibmqx4 device. We have compared the error mitigation schemes which accounts for correlations in readout errors with the one that does not account for them. Taking into account correlations proved to be crucial in the case of experiments with single probability vectors (Grover's algorithm for two qubits and the probability distributions in the five-qubit case). We point at the study of readout error correlations and accounting for them as a future research direction.
We are convinced that our method and its future extensions will find applications to mitigate readout errors in realistic near-term quantum devices that, despite being inherently noisy and imperfect, are expected to be capable of effectively solving problems of practical relevance. One can argue that performing quantum detector tomogra-phy and implementing our error-mitigation procedure will not be feasible for systems involving many qubits (due to the exponential growth of the size of this problem). However, many quantum algorithms tailored to near-term applications (such as quantum approximate optimization algorithm [51] or variational quantum eigensolvers [66]) can be performed by implementing just a polynomial number of a few-qubit quantum measurements. We expect that our scheme can be particularly useful in such scenarios.

B. Further research
Finally, we state a number of possible future research directions.
The first question concerns the possibility of physical correction of coherent readout errors. For the case of a singlequbit measurements there exists a natural possibility of correction of the coherent measurement errors. Namely, if the readout errors are uncorrelated one could obtain singlequbit POVMs with diagonal effects by implementing suitably chosen unitaries at the end of every quantum circuit. Moreover, in the case of uncorrelated errors, such physical correction is in principle possible for an arbitrary number of qubits (provided that readout errors dominate over gate errors). Determining whether this procedure will work for realistic near-term devices is an important research question.
Another important problem is to give confidence intervals (with respect to the operational distance defined in Eq. (4)) for the problem of quantum detector tomography. This kind of results are required to give estimates for the sample complexity of QDT and hence are crucial in assessing the feasibility of QDT on realistic near-term quantum devices. So far, such results have been obtained only for the case of quantum state tomography and the trace distance as the figure of merit [67,68].
Last but not least, it is desirable to extend our errormitigation procedure to larger systems with complicated geometry of the physical connections between qubits. To realize this goal it, is necessary to efficiently perform detector tomography of multi-qubit devices and understand if the readout correction is possible without the necessity to perform estimation of the full probability distribution.
the Hungarian National Research, Development and Innovation Office (NKFIH) through Grants No. K124351, K124152, K124176 KH129601, and the Hungarian Quantum Technology National Excellence Program (Project No. 2017-1.2.1-NKP-2017-00001); and he was also partially funded by the János Bolyai and the Bolyai+ Scholarships. We acknowledge use of the IBM Quantum Experience and Rigetti Forest SDK for this work. The views expressed are those of the authors and do not reflect the official policy or position of IBM or the IBM Quantum Experience team, or the Rigetti team.

NOTE ADDED
Upon completion of this manuscript, we became aware of a recent preprint [69] that provided a comprehensive characterization of detectors in the IBM devices and proposed analogous scheme for the mitigation of readout errors. Furthermore, we note that IBM's qiskit package has been recently updated to include the readout error correction scheme, which seems to rely on an analogous procedure. We have studied how fraction f changes with growing ratio δ+α Dop(M,P)+ in the case of single-qubit detectors. To this aim, we have implemented the above algorithm for the range of magnitudes of off-diagonal terms of POVM's elements, while keeping diagonal terms from actual IBM's data. Results are shown in Fig. 12. For all qubits the fraction of corrected statistics crosses 50% in the region in which δ + α ≈ D op (M, P) + , which backs up rule give in Eq. (B1). Furthermore, the fraction of successful error mitigations for actual experimental data lies between 88% and 99%, which suggests that for single-qubit experiments our protocol should be helpful in most cases.
For two-qubit case we have tested POVMs created via tensor product of single-qubit POVMs with increasing off-diagonal terms, again leaving diagonal terms from experimental data. Results of simulation are shown on Figure 13. Similarly to single-qubit case, fraction of successful mitigations crosses 50% around the regime where Eq. (B1) becomes unsatisfied. Furthermore, for actual experimental data f takes values between 98.86% and 99.99%, which suggests that for uncorrelated pairs of qubits our mitigation procedure should be successful for generic quantum states.   Fig. 6 April 28 Fig. 7 April 28 Fig. 8 April 28 Fig. 9 April 28 Fig. 10 April 28 Fig. 11 February 21