Inverse counting statistics based on generalized factorial cumulants

We propose a procedure to reconstruct characteristic features of an unknown stochastic system from the long-time full counting statistics of some of the system's transitions that are monitored by a detector. The full counting statistics is conveniently parametrized by so-called generalized factorial cumulants. Taking only a few of them as input information is sufficient to reconstruct important features such as the lower bound of the system dimension and the full spectrum of relaxation rates. The use of generalized factorial cumulants reveals system dimensions and rates that are hidden for ordinary cumulants. We illustrate the inverse counting-statistics procedure for two model systems: a single-level quantum dot in a Zeeman field and a single-electron box subjected to sequential and Andreev tunneling.

Counting the number N of transitions within a time interval [0; t] repeatedly many times yields the probability distribution P N (t), referred to as full counting statistics. For a well-characterized system, the states, the rates between them, and the coupling to the detector are known. It is, then, straightforward to compute the full counting statistics and compare it with experimentally measured data. Suppose, however, that the underlying model for a stochastic system is unclear and the only information available is the counting statistics measured by the detector. It is, then, desirable to have a systematic approach to distill out of the measured counting statistics the relevant information for reconstructing properties of the underlying model, including basic properties such as the number of states of the stochastic system. Such an approach can be dubbed inverse counting statistics [43].
What are the properties of the stochastic system that one may hope to reconstruct by inverse counting statistics? First of all, there is the number M of possible system states. Second, the stochastic system is characterized by the spectrum of relaxation rates with which it relaxes back to its (equilibrium or nonequilibrium) steady state after being disturbed externally [44][45][46][47][48][49]. Of course, inverse counting statistics cannot distinguish between different stochastic systems that are equivalent in the sense that they produce the same counting statistics, even if the stochastic systems possess a different numbers of states. Therefore, inverse counting statistics can at most deliver the minimal number of system states necessary that is compatible with the observed counting statistics.
To make the inverse counting statistics a powerful and practical tool, one should keep the circumstances for the acquisition of the input data as transparent and as simple as possible. Therefore, we restrict ourselves to the following scenario. First, only steady-state counting statistics is considered, i.e., we assume that the system has already relaxed before counting starts. This excludes studying transient behavior after a perturbation of the system. The latter would, on the one hand, offer a direct access to some relaxation rate of the system [50][51][52][53][54]. On the other hand, the determination of the full spectrum of relaxation rates would require the knowledge of how to perturb the system in order to probe a specific relaxation rate.
Second, we concentrate on the limit of long measuring-time intervals [0; t], for which the system dynamics is dominated by the slowest relaxation rate only. Nevertheless, the procedure of inverse counting statistics yields the full spectrum of all relaxation rates, as explained below.
Given a measured distribution P N (t), what are the input data for the inverse counting statistics? In [43], it was suggested to use the (long-time) cumulants of the distribution function. Here, we propose to use generalized factorial cumulants [55,56] instead. The advantage of the latter is that they depend on an arbitrarily chosen parameter s. The outcome of the inverse counting statistics (such as the number of system states or the spectrum of relaxation rates) should, however, not depend on this parameter s. Therefore, the s-independence of the results defines a powerful consistency criterion. Furthermore, as we will see in section 5, there are special cases in which part of the relaxation-rate spectrum is not accessible by inverse counting statistics with ordinary cumulants but is detectable by using generalized factorial cumulants with properly chosen parameter s.
As another difference to [43], we allow for a more general system-detector coupling by introducing the counting power m. In [43], the detector is assumed to be sensitive to only a single transition between two specific states increasing the detector counter just by one, the counting power is m = 1. If this transition increases the detector counter by k, the counting power is m = k. If several transitions are counted by the detector, the counting power can be even larger. We allow for detectors counting arbitrarily many transitions between arbitrarily many states increasing the detector counter by an arbitrary amount. Therefore, our inverse counting procedure does not only test compatibility with the number M of system states but also with the counting power m.
The paper is organized as follows. In section 2, we give a short introduction how full counting statistics, especially generalized factorial cumulants, are calculated by means of a Markovian master equation. Subsequently, in section 3, we explain the general procedure of inverse counting statistics. This procedure is, then, illustrated in sections 4 and 5 for two model systems: a single-level quantum dot in a Zeeman field and a singleelectron box subjected to sequential and Andreev tunneling. In section 6, we give a short summary of the inverse counting procedure introduced in this paper.

Full counting statistics for stochastic systems
In this paper, we consider stochastic systems as represented in figure 1, i.e., systems whose time evolution is completely governed by a Markovian master equation. There is a finite number M of system states, labeled by χ. Systems with coherent superposition of different χ are not taken into account. Arrows indicate transitions from state χ to state χ with a transition rate Γ χ χ . Non-Markovian effects [57][58][59] are not taken into account. The states of the master equation can be, e.g., three charge states of a metallic island coupled to one superconducting lead as discussed in the section 5, while the continuum of electronic states on the metallic island and lead enters only effectively via the transition rates.
The counting factors z k with k = 1, 2, . . . next to the dashed arrows specify the coupling of the detector to the system: if the system undergoes such a transition, the number N counted by the detector increases by k. There may be also processes with rate Γ χχ monitored by the detector in which the system ends up in the same state χ as it started, which is, e.g., the case for cotunneling through a magnetic atom [60].
Without the detector, the master equation for the probabilities p χ (t) to find the system in state χ takes the forṁ The first term of the sum describes transitions into state χ, the second term transitions out of state χ. To account for the detector, we need to replace this master equation for p χ (t) by an N -resolved one for p χ N (t), where N is the number of detector counts in the time interval [0; t]. For this, we introduce the coupling constants d k χχ which is 1 if the detector count is increased by k for the transition from χ to χ, and 0 otherwise. We obtainṗ Solving for the N -resolved probabilities yields the counting statistics of the detector via P N (t) = χ p χ N (t). We now turn to the question of how to extract from a given distribution P N (t) the information that serves as an input for the inverse counting statistics. One possibility would be to use the ordinary (long-time) cumulants of the distribution function, as suggested in [43]. The cumulants are obtained as derivatives C k (t) = ∂ k z ln M(z, t)| z=0 of the generating function M(z, t) = ∞ N =0 e N z P N (t). Here, we propose, as an alternative, to employ generalized factorial cumulants [55] in the long-time limit. Generalized factorial cumulants are derived from the generating function by evaluating the derivatives Counting factors enter the generating function (3) as powers of z and not e z as in the case of ordinary cumulants, which simplifies analytic expressions in the following. The special case s = 1 recovers the factorial cumulants recently introduced in the context of mesoscopic transport [61,62]. They are called factorial because the corresponding moments ∂ k z M 1 (z, t)| z=0 = N (k) are expectation values of the factorial power N (k) := N (N − 1) . . . (N − k + 1) instead of the ordinary power N k which defines ordinary moments. Generalized factorial cumulants depend on an extra parameter s that describes a shift of the complex variable z by the amount s − 1 along the real axis in the complex plane. In contrast to previous works [55,56], we define the generating function without a normalization factor 1/M s (0, t). Such a z-independent normalization factor would not influence the cumulants of order k > 0, but it would set, per definition, C s,0 = 0. Without this normalization factor, C s,0 (t) = ln N s N P N (t) contains nontrivial information that can be used for the inverse counting statistics. Both factorial (s = 1) and generalized factorial cumulants (s = 1) have been utilized to detect correlations in charge-transfer statistics [55,56,[61][62][63].
In the context of inverse counting statistics, the parameter s will play a very important role in two respects. First, since the underlying stochastic model for a measured distribution P N (t) is independent of s, the outcome of the inverse counting statistics must also be independent of the parameter s. Therefore, the required sindependence of the obtained results defines a criterion for the compatibility of the measured data with the assumed underlying model. Second, as we will see below, there are special cases in which the inverse counting statistics with factorial cumulants (s = 1) would indicate compatibility with a too small stochastic system and only generalized factorial cumulants (s = 1) reveal the higher dimension of the underlying stochastic model.
To calculate the generalized factorial cumulants for a given stochastic system, it is convenient to first perform a z-transform of the N -resolved master equation (2), i.e., multiply with z N and then sum over N . If we combine the z-transformed N -resolved probabilities p χ z = N z N p χ N of the different states χ in a vector p z (t), the z-transformed master equation can be written in the forṁ with matrix elements Two examples for W z are given in equations (14) and (18). For z = 1, equation (5) is nothing but the master equation (1).
The solution of equation (5) is p z (t) = exp (W z t) p z (0). Since p χ N (0) ∼ δ N,0 , the initial vector p z (0) is independent of z and describes the initial probability distribution. The matrix exponential exp (W z t) = M j=1 exp [λ j (z)t] r j,z ⊗ l T j,z can be expressed in terms of the eigenvalue spectrum {λ j (z)} of W z by making use of the decomposition into the left and right eigenvector l j,z and r j,z with normalization l T j,z · r j ,z = δ jj . For an arbitrary initial distribution p z (0) of the system, the systems relaxes exponentially in time to its steady state, governed by the eigenvalues λ j (z). The eigenvalues λ j (z) at z = 1 are the system's relaxation rates mentioned in the introduction. They are either real or they appear as complex-conjugated pairs. In the following, we assume that counting starts only after the system has reached its steady state, i.e., p z (0) is the stationary probability distribution, determined by W 1 p z (0) = 0 and e T · p z (0) = 1, where we defined e T = (1, . . . , 1) to sum over all states χ in p z (0).
Finally, taking into account that the generating function can be written as The summation over j complicates the time dependence of the generalized factorial cumulants. In the long-time limit, however, the above expression becomes considerably simpler since the exponential factors suppress all terms of the sum except the ones with the largest real part Re [λ j (z)] for z = s. For z = 1 and systems with a unique stationary state, the dominant eigenvalue is 0, i.e., all other eigenvalues have a negative real part. Around z = 1, the dominant eigenvalue, denoted by λ max (z), remains real and the limit provides well-defined constants, referred to as scaled long-time (generalized factorial) cumulants. These scaled long-time cumulants c s,k define the input information for the inverse counting statistics. Of course, in an experiment, the time t is always finite and, moreover, λ max (z) is not directly accessible. Therefore, one has to use scaled finite-time cumulants C s,k (t)/t, with C s,k (t) obtained via equations (3) and (4) from the measured P N (t). For large t, the scaled cumulants become time independent such that C s,k (t)/t ≈ c s,k .
If s is chosen very negative, it may happen that the dominant eigenvalues are given by a complex-conjugated pair. In this case, the limit is not well defined and the inverse counting-statistics procedure derived below cannot be applied.

Inverse Counting Statistics
In the previous section, we have shown how to calculate for a given stochastic model defined by the matrix W z (referred to as the generator of the stochastic system) the long-time (generalized factorial) cumulants. Inverse counting statistics deals with the opposite problem: how much can we learn about the stochastic system if only a few numbers, namely the experimentally determined values of the scaled long-time cumulants (up to some order), are given? To be more specific, we aim at the following properties of the stochastic system. First, a very important feature of the stochastic system is the dimension M of W z , i.e., the number of participating states in the stochastic process. Furthermore, the coupling to the detector is described by powers of z attached to some matrix elements of the generator. As a consequence, the characteristic polynomial det (λ1 − W z ) is of order m in z. Thus, as a second feature, we identify the counting power m. We will show below that the values of the first (m+1)M scaled longtime cumulants are enough to check compatibility with a stochastic system of dimension M and the counting power m characterizing the coupling to the detector.
But, with inverse counting statistics, we can get much more. From the (m + 1)M input parameters c s,k it is possible to determine the full spectrum of W z , i.e., the full z-dependence of the eigenvalues λ j (z). To appreciate how remarkable this statement is, let us remind that the input parameters are only a finite [(m + 1)M ] number of derivatives of only one eigenvalue λ max at only one value of z, namely the arbitrarily chosen s. From this rather restricted amount of information, we aim at reconstructing also the other eigenvalues different from λ max at all values of z different from s. How is this possible and how does it work in practice?
To answer this question, we observe that the characteristic function of the generator W z , is a polynomial both in λ (of order M ) and in z (of order m). The eigenvalues λ j (z), i.e., the zeros of the characteristic function, χ(z, λ j (z)) = 0, are, in general, nonanalytic functions in z. The characteristic function χ(z, λ) itself, however, is a polynomial in z and can, therefore, be written in the form where s is the arbitrarily chosen parameter of the generalized factorial cumulants. As a consequence, the (s-independent) characteristic function is fully determined by the (m + 1)M real (and s-dependent) coefficients a µν . This fixes all z-dependent (but sindependent) eigenvalues λ j (z) of the generator W z . For this reason, (m + 1)M input parameters are enough to fully determine the spectrum of W z .
Suppose that M and m are already known (we will discuss below how this is done with the help of inverse counting statistics). How do we get the spectrum of W z ? As input parameters we use the scaled generalized factorial cumulants c s,k for k = 0, . . . , (m + 1)M − 1 in the long-time limit. To determine the coefficients a µν , we perform l = 0, . . . , (m + 1)M − 1 times a derivative of χ(z, λ(z)) ≡ 0 with respect to z and set z = s afterwards. For technical reasons, it is convenient to divide the resulting equation by l!. Then, we arrive at the set of linear equations for l = 0, . . . , (m + 1)M − 1. The coefficients A l,µν , defined for nonnegative l, µ, and ν, are given by and otherwise (i.e., l ≥ µ together with ν ≥ 1) by Obviously, A l,µν depends on l and µ only via the difference l − µ (this was the reason of dividing by l!). The multiple sum over the α's is constrained by An alternative expression for the A l,µν can be found in Appendix A.
To be explicit, let us write down all the terms that are relevant for the case M = 3 and m = 2. We need ν up to 2 and l up to 8. For ν = 0 we get A l,µ0 = δ lµ , for ν = 1 we have A l,µ1 = c s,l−µ /(l − µ)!, and for ν = 2 we find To obtain the full spectrum {λ j (z)}, one first solves the set of linear equations (11) for a µν . Second, the result for a µν is inserted into equation (10) to get the characteristic function. Finally, the zeros of the characteristic function are determined. It all works because the characteristic function is fully determined by a finite number of coefficients a µν only.
It is important to remark that, in general, there is no guarantee that the set of linear equations (11) provides a unique solution. There may be special situations (we will discuss such a case below) in which the generator W z is separable in the sense that its characteristic function can be written as a product of two polynomials, the first one of order m and M in z and λ, the second one of order m − m and M − M . Of course, it is trivial that the characteristic function can always be written as a product of two polynomials in λ, but separability requires that, in addition, the two polynomials in λ are polynomials in z as well. For separable generators, inverse counting statistics has only access to the part of the spectrum to which the eigenvalue λ max (z) with the largest real part belongs, i.e., the effective problem has reduced dimensions M and m . The remaining question to be answered is how to determine the dimensions M and m of the stochastic system. In the spirit of [43], one may suggest that after having determined λ max (z) from the first (m + 1)M scaled long-time cumulants c s,k one can calculate the c s,k of higher order and compare them with experimentally measured values C s,k (t)/t in the long-time limit. If M and m were chosen correctly, then one expects a coincidence of calculated and measured values. This has, however, two downsides. First, measuring cumulants of increasingly higher order may become more and more difficult. Second, comparing just numbers to establish a consistency criterion may be of only limited significance in view of experimentally unavoidable noise.
With the use of generalized factorial cumulant, however, we can do much better. Remember that the parameter s was chosen arbitrarily and the result, i.e., the spectrum of the generator, must be independent of this parameter s. Therefore, we can use the very same measured time trace of the detector to determine scaled cumulants C s,k (t)/t in the long-time limit for different values of s and, afterwards, run the inverse countingstatistics procedure for each s. If M and m were chosen correctly, then the full z-dependence of the full spectrum should be independent of the choice of s. This establishes a much stronger consistency check than comparing just a few numbers.
In the following sections 4 and 5, we illustrate the inverse counting-statistics procedure for two model systems: a single-level quantum dot in a Zeeman field and a single-electron box subjected to sequential and Andreev tunneling.

Single-level quantum dot in a Zeeman field
The first model system is depicted in figure 2(a). A single-level quantum dot is weakly tunnel coupled to one normal-state metallic lead N and subjected to a magnetic field. The orbital level ε measured relative to lead's electrochemical potential is splitted by the Zeeman energy ∆ into ε σ = ε ± ∆/2. The positive (negative) sign applies to a spin σ =↑ (↓) electron on the quantum dot. The current through a nearby quantum point contact is sensitive to the dot charge which allows to monitor the charge transfer between dot and lead as function of time. The empty dot 0 can be occupied by a spin σ electron with the sequential tunneling rate Γ σ0 = Γf (ε σ ) given by Fermi's golden rule. The reverse transition occurs with the rate The temperature T is as large that both Γ σ0 and Γ 0σ are nonvanishing, but transitions to higher charge states are negligible due to charging energy.
The stochastic system is depicted in figure 2(b). Its generator is given by Each counting factor z in the upper-right off-diagonal matrix elements correspond to counting an electron leaving the dot. The characteristic function is a polynomial of order M = 3 in λ and of order m = 1 in z.
The case of vanishing magnetic field ∆ = 0, however, is special because of spin degeneracy ε ↑ = ε ↓ , and only two different rates Γ ↑0 = Γ ↓0 = Γ 10 and Γ 0↑ = Γ 0↓ = Γ 01 appear. As a consequence, the characteristic function becomes separable, The first factor is a polynomial of order M = 2 in λ and of order m = 1 in z, while the second factor is of order M = 1 and independent of z. Due to the z-independency, the second factor does not influence counting statistics and thus can not be detected anymore. The electron transfer can be completely described by a spinless orbital which is occupied with rate 2Γ 10 and emptied with rate Γ 01 . The factor 2 reflects the spin degeneracy [52]. Note that a separable characteristic function does not always separate in a z-dependent and z-independent factor, an example will be given in section 5.

Nonvanishing magnetic field
We start with discussing the case of nonvanishing magnetic field for which we choose ε = −k B T and ∆ = k B T /2. As unit of time we chose the inverse of Γ = Γ σ0 + Γ 0σ . The input information for the inverse counting statistics (for M = 3 and m = 1) is given Figure 3. (Color online) Eigenvalue spectrum obtained via inverse counting statistics for a single-level quantum dot ε = −k B T (a) inside a magnetic field ∆ = k B T /2 and (b) without magnetic field. The input data are scaled cumulants C s,k (t)/t for finite times t. Only in the long-time limit Γt 100, the output {λ j } is indeed the eigenvalue spectrum of the generator (14). Values λ j with a finite imaginary part (occurring for very negative z) are not depicted.
by the scaled long-time generalized factorial cumulants c s,k from order 0 up to order (m + 1)M − 1 = 5. Since in experiments, the measurement time is always finite, we do not take as input parameters the exact scaled long-time cumulants of the defined model but calculate, instead, the scaled finite-time cumulants C s,k (t)/t at some large but finite time t. Hence, the scaled cumulants are close but not identical to the exact values in the long-time limit. For the exact long-time scaled cumulants, a value for the dimension M that was assumed to be too large, can be immediately identified by trying to solve the system of linear equations (11). Then, no unique solution for every a µν is obtained, i.e., the linear equations are not independent from each other. However, in the following, we stick to the experimental situation that, due to the finite measuring time or due to experimental noise, a unique solution for the system of linear equations (11) is even obtained if M is assumed too large.
First, we determine the eigenvalue spectrum from the inverse counting statistics performed at s = 1 and assuming the correct values M = 3 and m = 1. As input parameter, we take the calculated scaled generalized factorial cumulants C s,k (t)/t at times Γt = 10, 20, 100. The result is shown in figure 3(a). For Γt = 100, the long-time limit is reached and no difference to the exact eigenvalues can be recognized anymore.
Next, we demonstrate the consistency check for the dimensions M and m. For this, we check the required s-independence of the eigenvalue spectrum. In the following, we always use as input information the calculated scaled cumulants at Γt = 6000 to ensure convergence to the asymptotic long-time behavior not only for s = 1 but also for other s.

Vanishing magnetic field
We now turn to the case of vanishing magnetic field ∆ = 0, for which the characteristic function is separable, i.e., the characteristic function is a product of two polynomials, one of order M = 2 and m = 1 in λ and z, and the other one is of order M = 1 and independent of z. The latter polynomial, i.e., the third dimension, does not influence transport anymore because the system can be described in a spineless model with M = 2 and m = 1.
Performing the inverse counting statistics at s = 1 with the correct values M = 2 and m = 1 gives the spectrum depicted in figure 3(b). The input information is given by the scaled finite-time cumulants C s,k (t)/t from order 0 to 3. For Γt = 100, the long-time limit is reached and no difference to the exact eigenvalues can be recognized anymore.
For the consistency check of the dimensions M and m, we use as input information the calculated scaled cumulants at Γt = 100. Horizontal contour lines in figure 5(a), (b) indicate that the eigenvalues are independent of s, i.e., the assumed dimensions M = 2 and m = 1 are compatible with the input data. In contrast for the choice M = 3 and m = 1, we obtain no s-independent spectrum of eigenvalues. The dimension M = 3 is too large.
For completeness, we discuss in Appendix B how close the magnetic field has to be tuned to ∆ = 0 in order to observe the discussed behavior.

Sequential and Andreev tunneling in a single-electron box
The second example illustrating the inverse counting-statistics procedure is a model system that has been already experimentally realized in [64][65][66]. The set up is depicted in figure 6(a). A single-electron box (SEB) is formed by one superconducting lead weakly coupled (characterized by the tunnel resistance R T ) to a normal-state metallic island. The energy required to bring n excess electrons on the island is E C (n − n G ) 2 . By applying a voltage V G to a gate electrode, the gate charge is tuned near n G = 0. The charge n on the island is monitored by an electrostatically coupled single-electron transistor (SET): each value of n results in a characteristic value of the current through the SET. Due to a finite temperature, transitions n = 0 → ±1 with the rate Γ ± u are possible (but temperature is small enough so that transitions to further charge states −2 and 2 are negligible). If the island is in one of these excited states, Andreev tunneling n = ±1 → ∓1 with the rate Γ ∓ A is possible until the island relaxes back to the ground state n = ±1 → 0 with rate Γ ∓ d . The stochastic system is depicted in figure 6(b). Its generator is given by Each counting factor z in the upper-right off-diagonal matrix elements correspond to counting an electron leaving the island. Since an Andreev-tunneling process transfers two electrons, Γ − A is multiplied with the square z 2 instead of a single z. This counting procedure is, especially, reasonable if the detector can not resolve two separate transitions −1 → 0 followed directly by 0 → +1 from the Andreev-tunneling process −1 → +1. For the chosen counting procedure, theses transitions need not to be distinguished.
If transitions are counted in which the number of electrons on the island is increased, the counting factors are removed from the upper-right off-diagonal matrix elements and put, instead, to the lower-left terms. The cumulants are identical for both counting procedures in the long-time limit. However, for finite times and asymmetric transition rates, the cumulants can differ, which indicates a violation of detailed balance [67].
The characteristic function is a polynomial of order M = 3 in λ and of order m = 2 in z. In general, the inverse counting statistics should deliver the z-dependence of all three eigenvalues.
The symmetry point n G = 0, however, is special. At this point, there are only three different rates As a consequence, the characteristic function becomes separable, The first factor is a polynomial of order M = 2 in λ and of order m = 1 in z, while the second one is of order M = 1 and m = 1. Depending on the choice of s, the inverse counting statistics will only provide the two eigenvalues of the first or the single eigenvalue of the second factor. Away from the symmetry point, n G = 0, the generator is not separable. Figure 7. (Color online) Eigenvalue spectrum obtained via inverse counting statistics for a normal-state metallic island weakly tunnel coupled to one superconducting lead. The gate charge is (a) n G = 0.09 and (b) n G = 0.00. The input data are scaled cumulants C s,k (t)/t for finite times t. Only in the long-time limit, (a) t 10 0 s and (b) t 10 −1 s, the output {λ j } is indeed the eigenvalue spectrum of the generator (18). Values λ j with a finite imaginary part (occurring for very negative z) are not depicted.
Instead of calculating the tunneling rates Γ ± u , Γ ± d , and Γ ± A in the presence of an electromagnetic environment [64,68], we rely on experimentally measured rates for n G = 0.09 in [65] and n G = 0.00 in [66]. In the former case, the experimental parameters are n G = 0.09, E C = 43µeV, ∆ = 216µeV, and R T = 2000 kΩ at 60 mK temperature. The measured rates are Γ + u = 10.

Non-symmetric case
We start with discussing the generic, non-symmetric case, for which we choose the n G = 0.09. The input information for the inverse counting statistics (for M = 3 and m = 2) is given by the scaled long-time generalized factorial cumulants c s,k from order 0 up to order (m + 1)M − 1 = 8. Similar to the case of the single-level quantum dot discussed in section 4, we take into account that the measurement time is always finite in an experiment and do not take as input parameters the exact scaled long-time cumulants of the defined model but calculate, instead, the scaled cumulants C s,k (t)/t at some large but finite time.
First, we determine the eigenvalue spectrum from the inverse counting statistics performed at s = 1 and assuming the correct values M = 3 and m = 2. As input parameter we take the calculated scaled generalized factorial cumulants at times t = 0.01, 0.1, 1 s. The result is shown in figure 7(a). For t = 1 s no difference to the exact eigenvalues can be recognized anymore.
Next, we demonstrate the consistency check for the dimensions M and m. For the

Symmetric case
We now turn to the symmetric case, n G = 0, for which the characteristic function is separable, i.e., the characteristic function is a product of two polynomials, one of order M = 2 and m = 1 in λ and z, and the other one of order M = 1 and m = 1. The choice of s determines whether the eigenvalue with the largest real part is a zero of the first or the second factor and, therefore, which and how many of the eigenvalues are accessible via the inverse counting statistics. If we choose s = 1 (corresponding to factorial cumulants), then we obtain only two of the three eigenvalues. For small s < −1.36 (for figure 7(b), we choose s = −2), we get only the third one. The number of required cumulants is 4 in the former and 2 in the latter case. The resulting spectrum of all three eigenvalues is depicted in figure 7(b) for the times t = 0.001, 0.01, 0.1 s. For t = 0.1 s, no difference to the exact eigenvalues can be recognized anymore. The consistency check of the dimensions M and m (for t = 20 s) is depicted in figure 9(a) for M = 1 and m = 1 and in figure 9(b), (c) for M = 2 and m = 1. If we perform the inverse counting statistics around s = 1 (or for any s > −1.36), then we conclude that the dimension M = 1 is too small (no horizontal contour lines in figure 9(a) for s > −1.36), but M = 2 seems to be sufficient (horizontal contour lines in Figs. 9(b) and (c) for s > −1.36). Thus, by employing only factorial cumulants (s = 1) one may be tempted to conclude that the dimension of the stochastic system is M = 2 only. If, on the other hand, inverse counting statistics is also done for s < −1.36, then the horizontal contour lines in figure 9(a) indicate that there is another eigenvalue. Since the obtained z-dependent eigenvalues are different from each other, we conclude that there must be, in total, three eigenvalues.
For completeness, we estimate in Appendix C how close n G has to be tuned to zero in order to observe the discussed behavior.

Inverse-counting-statistics manual
For practical use, we summarize the inverse-counting-statistics procedure introduced in this paper in the following step-by-step manual: (i) Some positive integer values for the number of system states M and the counting power m are assumed.
(ii) From the measured probability distribution P N (t), the scaled finite time cumulants C s,k (t)/t of order k = 0, . . . , (m + 1)M − 1 for some real s are calculated, see equations (3) and (4). The time t is chosen large enough for C s,k (t)/t to become time independent such that C s,k (t)/t ≈ c s,k .
(iii) From the c s,k the set of linear equations (11) is derived.
(iv) Solving this set for the coefficients a µν yields the characteristic function χ(z, λ) defined in equation (9).
If the z-dependence of each zero λ j (z) is independent of the parameter s chosen in step (ii), the values assumed in step (i) are lower bounds for M and m. Moreover, the λ j (z) are, indeed, eigenvalues of the system's generator W z and λ − λ j (z) linear factors of the generator's characteristic function. The λ j (1) are the system's relaxation rates. If the z-dependence is not independent of the parameter s, the assumption in (i) is falsified. One needs to start again at step (i) increasing, successively, from small values for M and/or m to larger values.
There are special cases of separable characteristic functions, see section 5, where different intervals of s-values reveal different dimensions and eigenvalues of the generator. The total dimension and full spectrum is, then, obtained by combining the results from the different s-intervals.

Conclusions
In this paper, we propose inverse counting statistics based on generalized factorial cumulants as a convenient and powerful tool to reconstruct characteristic features of a stochastic system from measured counting statistics of some of the system's transitions. Such a method is particularly useful in cases in which very little is a priori known about the stochastic system under investigation. As the only input information for the inverse counting-statistics procedure, we use a few experimentally determined numbers, namely the scaled generalized factorial cumulants in the long-time limit. Despite the limited amount of input, the inverse counting-statistics procedure yields a remarkable extended amount of output. First, we can determine a lower bound of M , the dimension of the stochastic system. Second, we can find a lower bound of the counting power m, which characterizes the coupling between stochastic system and detector. Third, we can reconstruct the characteristic function χ(z, λ) of the generator W z , which is a polynomial of order M in λ and a polynomial of order m in z. From the zeros of the characteristic function, we can, then, determine the full z-dependence of the full spectrum of eigenvalues of W z . This is quite a remarkable result, since the long-time cumulants used as input depend only on one of the eigenvalues, λ max , determined around one value of z.
The use of generalized factorial cumulants instead of ordinary ones is crucial for our proposal. While the evaluation of generalized factorial cumulants from a measured time trace of the detector signal does not introduce any extra complication as compared to the evaluation of ordinary cumulants, the benefit of having a free parameter s in the definition of the generalized factorial cumulants is immense. First, the outcome of the inverse counting statistics must be s-independent. Therefore, an s-dependent outcome of the inverse counting statistics immediately indicates a wrong choice of M or m. Second, there are special cases of separable characteristic functions, for which the inverse counting statistics with ordinary cumulants would only reveal a part of the spectrum of eigenvalues of the generator, while the variation of s makes it possible to access the full spectrum.
The proposed inverse counting-statistics procedure is quite general and, therefore, applicable to large variety of systems. To illustrate the procedure we choose two examples from electronic transport in nanostructures: a single-level quantum dot in a Zeeman field and a single-electron box subjected to sequential and Andreev tunneling. For the latter case, the full dimension of the system's state space and the full spectrum of eigenvalues can only be revealed by varying the parameter s. Appendix C. Consistency check for n G = 0.001 For the single-electron box, we estimate that at least for |n G | 0.001 the n G = 0 case is already reached in good approximation for s −1.5 or s −0.5 (compare figure C1 to figure 9). For figure C1, the Andreev-tunneling rates are approximated by Γ ± A ≈ (615.11 ± 11.42) Hz [69]. The sequential tunneling rates are estimated via an interpolation between the experimental values of [65]: Γ ± u ≈ (12.00 ± 0.03) Hz and Γ ± d ≈ (252.00 ± 0.83) Hz.