Tight finite-key security for twin-field quantum key distribution

Quantum key distribution (QKD) offers a reliable solution to communication problems that require long-term data security. For its widespread use, however, the rate and reach of QKD systems must be improved. Twin-field (TF) QKD is a step forward toward this direction, with early demonstrations suggesting it can beat the current rate-versus-distance records. A recently introduced variant of TF-QKD is particularly suited for experimental implementation, and has been shown to offer a higher key rate than other variants in the asymptotic regime where users exchange an infinite number of signals. Here, we extend the security of this protocol to the finite-key regime, showing that it can overcome the fundamental bounds on point-to-point QKD with around $10^{10}$ transmitted signals. Within distance regimes of interest, our analysis offers higher key rates than those of alternative variants. Moreover, some of the techniques we develop are applicable to the finite-key analysis of other QKD protocols.

Quantum key distribution (QKD) enables two remote parties, Alice and Bob, to generate a shared secret key in the presence of an eavesdropper, Eve, who may have unbounded computational power at her disposal [1][2][3]. While, ideally, the two parties can be at any distance, in practice, due to the loss and noise in the channel, point-to-point QKD is limited to a certain maximum distance at which secret key bits can securely be exchanged. In fact, the longest distance achieved to date in a terrestrial QKD experiment is about 400 km [4,5]. The main limitation is the exponential decrease of the transmittance, η, with the channel length in optical fibres. Even with a high repetition rate of 10 GHz, it would take an average of about two minutes to send a single photon over a distance of 600 km of standard optical fibres, and about 300 years to send it over 1000 km [6]. Indeed, fundamental bounds [7][8][9][10][11] on the private capacity of repeaterless point-to-point QKD protocols show that their secret-key rate scales at best approximately linearly with η. A protocol that aims to overcome this linear scaling must then include at least one middle node. Interestingly, this is not a sufficient condition. A well-known counterexample is the so-called measurement-device independent QKD (MDI-QKD) [12], which uses the middle node for an untrusted Bell-state measurement operation. There are, however, extensions of MDI-QKD that can improve its rate scaling from η to √ η by either using quantum memories [13,14] or quantum non-demolition measurements [15]. Such setups can, in fact, be considered to be the simplest examples of quantum repeaters [6,16], which are the ultimate solution to trust-free long-distance quantum communications [17]. However, even these simple versions may need more time to be efficiently implemented in practice [18,19].
Remarkably, the recently proposed twin-field QKD (TF-QKD) [20] can also overcome this linear scaling while using a relatively simple setup. TF-QKD is related to MDI-QKD, and it inherits its immunity to detector side-channels. However, it relies on single-photon, rather than two-photon, interference for its entanglement swapping operation. The secret-key rate of this protocol was first conjectured [20] and then proven [21,22] to scale with √ η too, making this approach a strong candidate to beat the current QKD records [23][24][25][26] with today's technology. The main experimental challenge is that single-photon interference needs very precise phase stability, which makes it more demanding than two-photon interference. Also, some of its current security proofs [21,22] need Alice and Bob to randomly choose a global phase, and then post-select only those rounds in which their choices match, which causes a drop in the secret key rate. Since the original proposal, several variants of TF-QKD have been developed [27][28][29][30], sharing the single-photon interference idea and its consequent √ η scaling, but differing in their experimental setups and security proofs. Moreover, some of these variants have been shown to be robust against phase misalignment [28][29][30], which simplifies their experimental implementation.
In this paper, we focus on the TF-QKD variant introduced in [28], which has two key features: (i) it does not Figure 1. Setup for the TF-QKD protocol [28] described in Box 1. Alice and Bob generate their raw key from the rounds in which they both select the X basis and Charlie declares that a single detector has clicked. The key bit is encoded in the phase of their coherent state. When the users select the same (a different) bit, the constructive (destructive) at Charlie's 50:50 beamsplitter interference should cause a click in detector Dc (D d ). The Z-basis PRCSs are only used to estimate the phase-error rate of the X-basis emissions.
Alice (Bob) chooses the key-generation basis X with probability pX or the parameter-estimation basis Z with probability pZ = 1 − pX , and (1.1) If she (he) chooses the X basis, she (he) generates a random bit bA (bB), prepares an optical pulse in the coherent state (−1) b A α ( (−1) b B α ), and sends it to Charlie.
(1.2) If she (he) chooses the Z basis, she (he) sends an optical pulse in a PRCS of intensity µ, selected from the set µ = {µ0, µ1, . . . , µ d−1 } with probability pµ, where d is the number of decoy intensities used.
They repeat step (1) for N rounds.
(2) Detection An honest Charlie measures each round separately by interfering Alice and Bob's signals at a 50:50 beamsplitter, followed by threshold detectors Dc and D d placed at the output ports corresponding to constructive and destructive interference, respectively. After the measurement, Charlie reports the pair (kc, k d ), where kc = 1 (k d = 1) if detector Dc (D d ) clicks and kc = 0 (k d = 0) otherwise. If he is dishonest, Charlie can measure all rounds coherently using an arbitrary quantum measurement, and report N pairs (kc, k d ) depending on the result. A round is considered successful (unsuccessful) if kc = k d (kc = k d ).
(3) Sifting For all successful rounds, Alice and Bob disclose their basis choices, keeping only those in which they have used the same basis. Let MX (MZ ) be the set of successful rounds in which both users employed the X (Z) basis, and let MX = |MX | (MZ = |MZ |) be the size of this set. Alice and Bob disclose their intensity choices for the rounds in MZ and learn the number of rounds M µν in MZ in which they selected intensities µ ∈ µ and ν ∈ µ, respectively. Also, they generate their sifted keys from the values of bA and bB corresponding to the rounds in MX . For those rounds in which kc = 0 and k d = 1, Bob flips his sifted key bit.
(4) Parameter estimation Alice and Bob apply the decoy-state method to M µν , for µ, ν ∈ µ, obtaining upper-bounds M U nm on the number of rounds Mnm in MZ in which they sent n and m photons, respectively. They do this for all n, m ≥ 0 such that n + m is even and n + m ≤ Scut for a prefixed parameter Scut. Then, they use this data to obtain an upper bound N U ph on the number of phase errors, N ph , in their sifted keys, and check if the upper bound e U ph = N U ph /MX is lower than a predetermined threshold value. If so, they continue to the last step; otherwise they abort the protocol.  5.1) Error correction: Alice sends Bob a pre-fixed amount λEC of syndrome information bits through an authenticated public channel, which Bob uses to correct errors in his sifted key.
(5.2) Error verification: Alice and Bob compute a hash of their error-corrected keys using a random universal hash function, and check whether they are equal. If so, they continue to the next step; otherwise, they abort the protocol.
(5.3) Privacy amplification: Alice and Bob extract a secret key pair (SA, SB) of length |SA| = |SB| = from their error-corrected keys using a random two-universal hash function.
Parameter estimation and Secret-key rate analysis The main contribution of this work-see Methods for the details-is a procedure to obtain a tight upper-bound N U ph on the total number of phase errors N ph in the finite-key regime for the protocol described in Box 1. Namely, we find that except for an arbitrarily small failure probability ε, it holds that where p nm|X (p nm|Z ) is the probability that Alice and Bob's joint X (Z) basis pulses contain n and m photons, respectively, given by with p n|µ = µ n exp(−µ)/n! being the Poisson probability that a PRCS pulse of intensity µ will contain n photons; ∆ and ∆ nm are statistical fluctuation terms defined in step 4 of Box 2; N 0 (N 1 ) is the set of non-negative even (odd) integers; and the rest of the parameters have been introduced in Box 1. The phase-error rate is then simply upperbounded by e U ph := N U ph /M X . Box 2 provides a step-by-step instruction list to apply our results to the measurement data obtained in an experimental setup.
When it comes to finite-key analysis, there is one key difference between the protocol in Box 1 and several other protocols, such as, for example, decoy-state BB84 [44], decoy-state MDI-QKD [35], and sending-or-not-sending TF-QKD [41]. In all the latter setups, when there are no state-preparation flaws, the single-photon components of the two encoding bases are mutually unbiased; in other words, they look identical to Eve once averaged by the bit selection probabilities. This implies that such states could have been generated from a maximally entangled bipartite state, where one of its components is measured in one of the two orthogonal bases, and the other half represents an encoded key bit. In fact, the user(s) could even wait until they learn which rounds have been successfully detected to decide their measurement basis, effectively delaying their choice of encoding basis. This possibility allows the application of a random sampling argument: since the choice of the encoding basis is independent of Eve's attack, the bit error rate of the successful X-basis emissions provides a random sample of the phase-error rate of the successful Z-basis emissions, and vice-versa. Then, one can apply tight statistical results such as the Serfling inequality [45] to bound the phase-error rate in one basis using the measured bit-error rate in the other basis. This approach, however, is not directly applicable to the protocol in Box 1, in which the secret key is extracted from all successfully detected X-basis signals, not just from their single-photon components. Moreover, the encoding bases are not mutually unbiased: the Z-basis states are diagonal in the Fock basis, while the X-basis states are not. This will require a different, perhaps more cumbersome, analysis as we highlight below.
To estimate the X-basis phase-error rate from the Z-basis measurement data, we construct a virtual protocol (see Box 3) in which the users learn their basis choice by measuring a quantum coin after Charlie/Eve reveals which rounds were successful. Note that, because of the biased basis feature of the protocol, the statistics of the quantum coins associated to the successful rounds could depend on Eve's attack. This means that the users cannot delay their choice of basis, which prevents us from applying the random sampling argument. Still, it turns out that the quantum coin technique now allows us to upper-bound the average number of successful rounds in which the users had selected the X basis and undergone a phase error. This bound is a non-linear function of the average number of successful rounds in which they had selected the Z basis and respectively sent n and m photons, with n + m even. More details can be found in the Methods Section; see Eq. (20).
The main tool we use to relate each of the above average terms to their actual occurrences, N ph and M nm , is Azuma's inequality [39], which is widely used in security analyses of QKD to bound sums of observables over a set of rounds of the protocol (in our case, the set of successful rounds after sifting), when the independence between the observables corresponding to different rounds cannot be guaranteed. When using Azuma's inequality, the deviation term ∆ scales with the square root of the number of terms in the sum. In our case, ∆ scales with √ M s , where M s is the number of successful rounds after sifting. For parameters of comparable magnitude to M s , this provides us with a reasonably tight bound. Whenever the parameter of interest is small, however, the provided bound could instead be loose. This is the case for the crucial term M U 00 in Eq. (1), as vacuum states are unlikely to result in successful detection events, thus the bound obtained with Azuma's inequality can be loose. This is important because, in Eq. (1), the coefficient associated to the vacuum term is typically the largest. It is then essential to find a tighter bound for this term.
For this, we employ a remarkable recent technique to bound the deviation between a sum of dependent random variables and its expected value [38]. This technique provides a much tighter bound than Azuma's inequality when the value of the sum is much lower than the number of terms in the sum. In particular, it provides a tight upper-bound for the vacuum component M 00 . In Methods, we provide a statement of the result and we explain how we apply it to our protocol.
Having obtained e U ph , we show in Supplementary Note A that the secret-key length, , can be lower bounded by while guaranteeing that the protocol is c -correct and s -secret, with s = 2ε is the Shannon binary entropy function, λ EC is number of bits that are spent in the error-correction procedure, PA is the failure probability of the privacy amplification scheme, and ε is the failure probability associated to the estimation of the phase-error rate. Here, our security analysis follows the universal composable security framework [46,47], according to which a protocol is sec -secure if it is both c -correct and s -secret, with sec ≥ c + s . The correctness criterion is met when Alice and Bob's secret keys S A and S B are identical, and the protocol is ccorrect when Pr[S A = S B ] ≤ c . The secrecy criterion is met when the classical-quantum state ρ AE describing Alice's secret key and Eve's side information is of the form ρ AE = U A ⊗ ρ E , where U A is the uniform distribution over all bit strings, and ρ E is an arbitrary quantum state. The protocol is s -secret if where · is the trace norm.
Box 2: Instructions for experimentalists 1. Set the security parameters c and PA, as well as the failure probabilities εc and εa for the inverse multiplicative Chernoff bound and the concentration bounds for sums of dependent random variables, respectively. Calculate the overall failure probability ε of the parameter estimation process, which depends on the number of times that the previous two inequalities are applied. In general, ε = d 2 εc + S cut 2 + 1 2 εa + εa, where d is the number of decoy intensities employed by each user. For Scut = 4 and three decoy intensities, we have that ε = 9εc + 10εa.

2.
Use prior information about the channel to obtain a predictionM U 00 on M U 00 , the upper bound on the number of Z-basis vacuum events that will be obtained after applying the decoy-state method.
6. Use Eq. (1) to find N U ph and set e U ph = N U ph /MX . 7. Use Eq. (4) to specify the required amount of privacy amplification and to find the corresponding length of the secret key that can be extracted. The key obtained is sec-secure, with sec ≥ c + s and s = 2ε + PA.

DISCUSSION
In this section, we analyse the behaviour of the secret-key rate as a function of the total loss. We simulate the nominal no-Eve scenario with an honest Charlie. In this case, the total Alice-Bob loss includes the loss in the quantum channels as well as the inefficiency of Charlie's detectors. We compare the key rate for the protocol in Box 1, using the finite-key security analysis introduced in the previous section, with that of the sending-or-not-sending TF-QKD protocol [30,41], as well as with the finite-key analysis presented in Ref. [40]. We also include the asymptotic secret key capacity for repeaterless QKD systems over lossy channels, known as the PLOB bound [9], for comparison. While specific bounds for the finite-key setting have recently been studied [10,48], in the practical regimes of interest to this work, they numerically offer a negligible difference to the PLOB bound. The latter has then been used in all relevant graphs for consistency. To simulate the data that would be obtained in all protocols, we use the simple channel model described in Supplementary Note C, which accounts for phase and polarisation misalignments. Also, we assume that both users employ three decoy-state intensities µ 0 > µ 1 > µ 2 . Since the optimal value µ 2 = 0 is typically difficult to achieve in practice, we set µ 2 = 10 −4 and optimise the secret-key rate over the value of µ 0 and µ 1 . We also optimise it over the selection probabilities, as well as over p X and α.  Table I. List of parameters used in the numerical simulations. Here, p d is the dark count probability, per pulse, of the detectors, δ ph is the phase misalignment of the system, δ pol is the polarisation misalignment of the system, and f is the error-correction inefficiency. In our numerical simulations, we set ε = PA = s/3.
The nominal values for system parameters are summarised in Table I. We assume a phase mismatch of 9.1% between Alice and Bob's signals, corresponding to a QBER of around 2% for most attenuations, matching the experimental results in [23]. For brevity, we do not consider the effect of polarisation misalignment in our numerical results, but one can use the provided analytical model to study different scenarios of interest. In principle, even if the mechanism used for polarisation stability is not perfect, one can use polarisation filters to ensure that the same polarisation modes are being coupled at the 50:50 beamsplitter, at the cost of introducing additional loss. We assume an error correction where e X is the bit error rate of the sifted key, and f is the error correction inefficiency. For the security bounds, we set c = s = 10 −10 , and for simplicity we set ε = PA = s /3.   Table I. In Fig. 2, we display the secret key rate per pulse for the protocol in Box 1 for different values of the block size, N , of transmitted signals. It can be seen that the protocol can outperform the repeaterless bound for a block size of just 10 10 transmitted signals per user, at an approximate total loss of 50 dB. For standard optical fibres, this corresponds to a total distance of 250 km, if we neglect the loss in the photodetectors. At a 1 GHz clock rate, it takes only ten seconds to collect the required data. For a block size of 10 11 transmitted signals, the protocol can already outperform Secret key rate (bits per pulse)  Table I. the repeaterless bound for a total loss ranging from 45 dB to over 80 dB. By increasing N , we approach the asymptotic performance of the protocol.
The dependence of the secret key rate on the block size N has been shown in Fig. 3, at a fixed total loss of 50 dB and for several values of phase misalignment. In all cases, there is a minimum required block size to obtain a positive key rate. This minimum block size can be even lower than 10 9 in the ideal case of no phase misalignment, and it goes up to around 10 10 at δ ph = 20%. There is a sharp increase in the secret key rate once one goes over this minimum required block size after which one slowly approaches the key rate in the asymptotic limit. The latter behaviour is likely due to the use of Azuma's inequality. One can, nevertheless, overcome the repeaterless bound at a reasonable block size in a practical regime where δ ph ≤ 15%. At higher values of total loss this crossover happens at even larger values of δ ph .
In Fig. 4, we compare the performance of the protocol in Box 1 with that of the sending-or-not-sending TF-QKD protocol presented in [30,41]. In the asymptotic regime, the former protocol outperforms the latter at all values of total loss. For a block size of 10 12 transmitted signals, this is still the case up to 80 dB of total loss, after which the key rate is perhaps too low to be of practical relevance. For a block size of 10 10 transmitted signals, however, the curves for the two protocols cross at around 55 dB, where they happen to cross the repeaterless bound as well. In this case, the sending-or-not-sending protocol offers a better performance after this point. This behaviour is due to the different statistical fluctuation analyses applied to the two protocols. As explained in the Result section, the single-photon components in the sending-or-not-sending protocol are mutually unbiased, allowing for a simpler and tighter estimation of the phase-error rate. This is not the case for the TF-QKD protocol in Box 1, for which this estimation involves the application of somewhat looser bounds for several terms in Eq. (1). We conclude that for sufficiently large block sizes, the protocol in Box 1 maintains its better performance over the sending-or-not sending variant.
Finally, in Fig. 5, we compare our results with those of the alternative analysis in [40]. To compute the secret-key rate of the latter, we use the code provided by the authors, except for the adjustments needed to match it to the channel model described in Supplementary Note C. It can be seen that, at all cases considered and for most practical regimes of interest, the analysis introduced in this paper provides a higher key rate than that of [40]. Moreover, we remark that the security proof presented in [40], in its current form, is only applicable when the state generated by the weakest decoy intensity µ 2 is a perfect vacuum state of intensity µ 2 = 0. The security analysis presented in this work, however, can be applied to any experimental value of µ 2 , and we assume a value of µ 2 = 10 −4 , which may be easier to achieve in practice. That said, the security proof in [40] adopts an interesting approach that results in a somehow simpler statistical analysis. In particular, unlike in the analysis presented in this paper, the authors in [40] Table I. do not estimate the detection statistics of photon-number states as an intermediate step to bounding the phase-error rate. Instead, they show that the operator corresponding to a phase-error can be bounded by a linear combination of the Z-basis decoy states. While this linear bound is asymptotically looser than the non-linear formula in Eq. (1), it allows the application of a simpler statistical analysis based on a double use of Bernoulli sampling. Given that the finite-key analysis of a protocol could be part of the software package of a product, we believe that the additional key rate achievable by our analysis justifies its slightly more complex approach.
In conclusion, we have proven the security of the protocol in Box 1 in the finite-key regime scenario against coherent attacks. Our results show that, under nominal working conditions experimentally achievable by today's technology, this scheme could outperform the repeaterless secret-key rate bound in a key exchange run of only ten seconds, assuming a 1 GHz clock rate. It would also outperform other TF-QKD variants, as well as alternative security proofs in practical regimes of interest.

METHODS
In this section, we introduce the procedure that we use to bound the phase-error rate of the protocol in Box 1, referring to the Supplementary Notes when appropriate. For notation clarity, we assume the symmetric scenario in which Alice and Bob employ the same X-basis amplitude α and the same Z-basis intensities µ ∈ µ. However, the analysis can be applied as well to the asymmetric scenario [42,43] by appropriately redefining the parameters p nm|X and p nm|Z .

Virtual protocol
To bound the phase-error rate, we construct a virtual protocol in which Alice and Bob measure an observable that is conjugate to that used to generate the key. By the complementarity argument [31], the bit-error rate of this virtual protocol is identical to the phase-error rate of the actual protocol, provided that the two protocols are equivalent. The equivalence condition is satisfied if the two protocols send the same quantum and classical information, thus Eve cannot tell which of the two is being performed. More concretely, our virtual protocol replaces Alice's X-basis Analysis in [40], N = 10 10 This work, N = 10 12 Analysis in [40], N = 10 12 This work, N → ∞ Analysis in [40], N → ∞ Repeaterless bound [9] Figure 5. Comparison between this work and the alternative analysis in [40], for different block size values N of transmitted signals. The simulation parameters are listed in Table I. emissions by the preparation of the state where A is an ancilla system at Alice's lab, a is the photonic system sent to Eve, and |± = 1 √ 2 (|0 ± |1 ); while Bob's X basis emissions are replaced by a similarly defined |ψ x Bb . After Eve's attack, Alice and Bob measure systems A and B in the Z basis {|0 , |1 }, which is conjugate to the X basis {|+ , |− } that they would use to generate the key. It is useful to write the state in Eq. (7) as where |C 0 and |C 1 are the (unnormalised) cat states Alice's Z-basis emissions are diagonal in the Fock basis, and the virtual protocol replaces them by their purification where p n|Z = µ∈µ p µ p n|µ is the probability that Alice's Z basis pulse contains n photons, averaged over the selection of µ. Unlike in the actual protocol, in the virtual protocol Alice and Bob learn the photon number of their signals by measuring systems A and B after Eve's attack. Lastly, Alice's emission of |ψ x Aa with probability p X and |ψ z Aa with probability p Z is replaced by the generation of the state where A c is a quantum coin ancilla at Alice's lab; while Bob's is replaced by an equally defined |ψ BcBb . Alice and Bob measure systems A c and B c after Eve's attack, delaying the reveal of their basis choice. The different steps of the virtual protocol are described in Box 3.

Box 3: Virtual protocol
(1) Preparation Alice and Bob prepare N copies of the state |φ = |ψ AcAa ⊗ |ψ BcBb and send all systems a and b to Eve over the quantum channel.
(2) Detection Eve performs an arbitrary general measurement on all the subsystems a and b of |φ ⊗N and publicly announces N bit pairs (kc, k d ). Without loss of generality, we assume that there is a one-to-one correspondence between her measurement outcome and her set of announcements. A round is considered successful (unsuccessful) if kc = k d (kc = k d ). Let M (M) represent the set of successful (unsuccessful) rounds.
(3) Virtual sifting For all rounds, Alice and Bob jointly measure the systems Ac and Bc, learning whether they used the same or different bases, but not the specific basis they used. Let Ms (M d ) denote the set of successful rounds in which they used the same (different) bases.
(4) Ancilla measurement  Two points from the virtual protocol in Box 3 require further explanation. The first is that, in the real protocol, Bob flips his key bit when Eve reports k c = 0 and k d = 1. This step is omitted from the virtual protocol, since the X-basis bit flip gate σ z has no effect on Bob's Z-basis measurement result. The second point concerns step 5, which may appear to serve no purpose, but it is needed to ensure that the classical information exchanged between Alice and Bob is equivalent to that of the real protocol. The term p µ|n is the probability that Alice's (Bob's) Z-basis n-photon pulse originated from intensity µ, and it is given by Phase-error rate estimation We now turn our attention to Alice and Bob's measurements in step (4.1) in Box 3. Let u ∈ {1, 2, ..., M s } index the rounds in M s , and let ξ u denote the measurement outcome of the u-th round. The possible outcomes are ξ u = X ij , corresponding to |00 AcBc |ij AB , where i, j ∈ {0, 1}; and ξ u = Z nm , corresponding to |11 AcBc |n, m AB , where n and m are any positive integers. Note that the outcomes |10 AcBc and |01 AcBc are not possible due to the previous virtual sifting step. A phase error occurs when ξ u ∈ {X 00 , X 11 }. In Supplementary Note D, we prove that the probability to obtain a phase error in the u-th round, conditioned on all previous measurement outcomes in the protocol, is upper-bounded by where F u−1 is the σ-algebra generated by random variables ξ 1 , ..., ξ u−1 , N 0 (N 1 ) is the set of non-negative even (odd) numbers, and the probability terms p nm|X and p nm|Z have been defined in Eqs. (2) and (3). In Eq. (13), for notation clarity, we have omitted the dependence of all probability terms on the outcomes of the measurements performed in steps (2) and (3) in Box 3.
From Eq. (31), we have that, except with probability ε a , Pr (ξ u ∈ {X 00 , X 11 }|F u−1 ) + ∆, (14) where N ph is the number of events of the form ξ u ∈ {X 00 , X 11 } in M s , and ∆ = 1 2 M s ln ε −1 a is a deviation term. Similarly, from Eq. (31), we have that, except with probability ε a ,

Ms u=1
Pr where M nm is the number of events of the form ξ u = Z nm in M s . As we will explain later, this bound is not tight when applied to the vacuum counts M 00 . For this term, we use Eq. (34), according to which, except with probability ε a ,
In this case, the deviation term is given by where a and b can be found by substitutingΛ n byM U 00 in Eq. (32). Now we will transform Eq. (13) to apply Eqs. (14) to (16). Let us denote the right-hand side of Eq.
We are now ready to apply Eqs. (14) to (16) to substitute the sums of probabilities by N ph and M nm . However, note that Alice and Bob only estimate the value of M nm for terms of the form n + m ≤ S cut and it is only useful to substitute Eq. (15) for these terms. With this in mind, we obtain where ∆ nm = ∆ except for ∆ 00 .
We still need to deal with the sum over the infinitely many remaining terms of the form n + m > S cut . For them, we apply the following upper bound where ξ u = Z denotes that Alice and Bob learn that they have used the Z basis in the uth round in M s ; and M Z is the number of events of the form ξ u = Z obtained by Alice and Bob. In the last step, we have used Eq. (31), using an identical argument as in Eq. (14). When we apply Eq. (21) Substituting Eq. (21) into Eq. (20), and isolating N ph , we obtain Note that the right hand side of Eq. (24) is a function of the measurement counts M nm , which are unknown to the users. They must be substituted by the upper bounds M U nm obtained via the decoy-state analysis. After doing so, we obtain Eq. (1). The failure probability ε associated to the estimation of N ph is upper-bounded by summing the failure probabilities of all concentration inequalities used. That includes each application of Eqs. (31) and (34), which fail with probability ε a ; and each application of the multiplicative Chernoff bound, which fails with probability ε c . In the case of three decoy intensities and S cut = 4, we have ε = 9ε c + 10ε a . In our simulations, we set ε c = ε a for simplicity.

Decoy-state analysis
Since Alice and Bob's Z-basis emissions are a mixture of Fock states, the measurement counts M nm have a fixed value, which is nevertheless unknown to them. Instead, the users have access to the measurement counts M µν , the number of rounds in M Z in which they selected intensities µ and ν, respectively. To bound M nm , we use the decoy-state method [32][33][34]. This technique exploits the fact that Alice and Bob could have run an equivalent virtual scenario in which they directly send Fock states |n, m with probability p nm|Z . Then, after Eve's attack, they know the number M nm of successful events in which they respectively sent n and m photons, and they assign each of them to intensities µ and ν with probability where p µν = p µ p ν and p nm|µν = p n|µ p m|ν . Each of these assignments can be regarded as an independent Bernoulli random variable, and the expected number of events M µν assigned to intensities µ and ν is In the actual protocol, Alice and Bob know the realisations M µν of these random variables. By using the inverse multiplicative Chernoff bound [50,51] This problem can be solved numerically using linear programming techniques, as described in the Supplementary Note 2 of [35]. While accurate, this method can be computationally demanding. For this reason, we have instead adapted the asymptotic analytical bounds of [42,52] to the finite-key scenario and used them in our simulations. The results obtained using these analytical bounds are very close to those achieved by numerically solving Eq. (27). This analytical method is described in Supplementary Note B.

Concentration inequality for sums of dependent random variables
A crucial step in our analysis is the substitution of the sums of probabilities in Eq. (19) by their actual realisation in the protocol. Typically, this is done by applying the well-known Azuma's inequality [39]. Instead, we use the following novel result [38]: Let ξ 1 , ..., ξ n be a sequence of random variables satisfying 0 ≤ ξ l ≤ 1, and let Λ l = l u=1 ξ u . Let F l be its natural filtration, i.e. the σ-algebra generated by {ξ 1 , ..., ξ l }. For any n, and any a, b such that b ≥ |a|, By replacing ξ l → 1 − ξ l and a → −a, we also derive In our analysis, we apply Eqs. (28) and (29) to sequences ξ 1 , ..., ξ n of Bernoulli random variables, for which E(ξ u |F u−1 ) = Pr(ξ u = 1|F u−1 ). Now, if we set a = 0 on Eqs. (28) and (29), we obtain This is a slightly improved version of the original Azuma's inequality, whose the right-hand side is exp − 1 2 b 2 . Equating the right hand sides of Eq. (30) to ε a , and solving for b, we have that n u=1 Pr (ξ u = 1|ξ 1 , ..., ξ u−1 ) ≤ Λ n + ∆, Pr (ξ u = 1|ξ 1 , ..., ξ u−1 ) + ∆, (31) with ∆ = 1 2 n ln ε −1 a , and where each of the bounds in Eq. (31) fail with probability at most ε a . The bound in Eq. (31) scales with √ n, and it is only tight when Λ n is of comparable magnitude to n. When Λ n n, one can set a and b in Eq. (28) appropriately to obtain a much tighter bound. To do so, one can use previous knowledge about the channel to come up with a predictionΛ n of Λ n before running the experiment. Then, one obtains the values of a and b that would minimise the deviation term if the realisation of Λ n equalledΛ n , by solving the optimisation problem The solution to Eq. (32) is After fixing a and b, we have that except with probability ε a , where In our numerical simulations, we have found the simple bound in Eq. (31) to be sufficiently tight for all components except the vacuum contribution M 00 . For this latter component, we use Eq. (34) instead. However, note that the users do not know the true value of M 00 , even after running the experiment. Instead, they will obtain an upper-bound M U 00 on M 00 via the decoy-state method, and they will apply Eq. (34) to this upper bound. Therefore, to optimise the bound, the users should come up with a predictionM U 00 on the value of M U 00 that they expect to obtain after running the experiment and performing the decoy-state analysis, and then substituteΛ n →M U 00 in Eq. (32) to obtain the optimal values of a and b. To findM U 00 , one can simply use their previous knowledge of the channel to come up with predictionsM µν of M µν , and run the decoy-state analysis using these values to obtainM U 00 .

Supplementary Note A: Security bounds
Let X (X ) denote Alice's (Bob's) sifted key of length M X before the post-processing step of the protocol. After the error correction and verification steps, Bob should have turned X into a copy of X. Then, Alice and Bob apply a privacy amplification scheme based on two-universal hashing to obtain a shorter secret key of length . The protocol is s -secret if [53] s ≤ 2ε + where E represents Eve's total side information about X, and H ε min (X|E ) is the ε-smooth min entropy of X conditioned on E . Let E denote Eve's side information before the error correction step. By the chain rule for smooth min-entropies [53], where λ EC (log 2 2 c ) is the number of bits revealed in the error correction (verification) step of the protocol. We now make use of the following theorem, introduced in [54], which we reproduce here for completeness.
Theorem [54]: Let ε > 0, ρ AEB be a tripartite quantum state, X = {M x } and Z = {N z } be two POVMs on A, and X (Z) be the result of the measurement of X (Z). Then, To apply this theorem, we consider a slight modification to our virtual protocol in Box 3. In step (4.1), Alice and Bob now first measure all the basis ancillas A c and B c in M s , keeping only the M X successful rounds in which they used the X basis, which we denote as the set M X . Let ρ AEB be the quantum state that describes all systems A and B in M X , as well as Eve's side information E on them. Note that if Alice measures all her systems A in the X basis, she will obtain a raw key X that is identical to the one she would have obtained in the real protocol; while if she measures them in the Z basis, she will obtain a raw key Z that is identical to that of the virtual protocol.
Let X = {M x } (Z = {N z }) denote Alice's overall POVM if she chooses to measure all her systems A in the X (Z) basis. The elements M x of X are of the form |x 1 x 2 x 3 ... x 1 x 2 x 3 ...|, where x n ∈ {+, −} is the result of the measurement of round n ∈ {1, ..., M X }. Conversely, the elements N z of Z are of the form |z 1 z 2 z 3 ... z 1 z 2 z 3 ...|, with z n ∈ {0, 1}. Since M x and N z are rank 1 projective measurements, we have that where in the last step we have used the fact that x n |z n 2 = 1/2, independently of the value of x n and z n . From Eq. (A3), it follows that Now, let us assume that Bob measures his systems B using POVM Z, obtaining a string Z that is identical to the one that he would obtain in the virtual protocol. Clearly, the result of a measurement of B cannot contain more information about Z than system B itself, and therefore from which we finally obtain In the Methods section, we show that the error rate between Z and Z is bounded by e U ph , except with a failure probability ε. Therefore [53], and by combining Eqs. (A2), (A8) and (A9), we have that By substituting Eq. (A10) and the secret key length given in Eq. (4) of the main text into Eq. (A1), we finally obtain that the protocol is s -secret if we choose s such that s ≤ 2ε + PA . (A11)

Supplementary Note B: Analytical estimation method
In this Note, we present an analytical method to obtain the upper bounds M U nm in Eq. (1), using the observed quantities M µν . First, we explain the general idea behind the procedure, and then we obtain specific analytical bounds for the case of three decoy intensities and S cut = 4, which we use in our simulations. We have numerically verified that the choice of three decoy intensities is optimal for reasonable block size values below 10 12 transmitted signals.
Our starting point is Eq. (26), which we rewrite aŝ where we have singled out the index (i, j), ensured that c ij > 0, and defined S + (S − ) as the set of pairs (n, m) = (i, j) such that c nm is a positive (negative) number. From Eq. (B4), one can obtain the following upper bound on M ij where we have chosen c max such that c max ≥ |c nm | for all the pairs (n, m) ∈ S − , and the last lower bound depends on the particular M ij that we are trying to estimate, as we will show later. Substituting the three bounds in Eq. (B5) and isolating M ij , we obtain (B7)

Bounds for three decoy intensities
Now, we obtain explicit lower bounds for the case in which S cut = 4 and each of Alice and Bob use three different intensity settings, satisfying µ 0 > µ 1 > µ 2 and ν 0 > ν 1 > ν 2 , respectively. For this, we take inspiration from the asymptotic analytical bounds derived in [42]. First, we define which is a function of some intensities that satisfy a S > a I and b S > b I , with a S , a I ∈ {µ 0 , µ 1 , µ 2 } and b S , b I ∈ {ν 0 , ν 1 , ν 2 }. The coefficients κ µ A and κ ν B depend on the specific M ij that is to be estimated, but we have omitted this dependence from the notation for simplicity. Using the previous equation, we now define where the coefficients w µν A and w µν B also depend on the particular M ij that we want to estimate. If we rewrite Ω as Ω = 2 k,l=0ĉ µ k ν lM µ k ν l , it is easy to prove that if the coefficients w µν A , w µν B , κ µ A and κ µ B are all positive, the coefficientŝ c µ k ν l are always positive (negative) when k + l is even (odd). Thus, one can find upper (Ω U ) and lower (Ω L ) bounds on Ω by properly replacing eachM µ k ν l by either its upper or lower bound, as explained in the introduction of this Note.

a. Upper bound on M00
By substituting κ µ A = κ µ B = µ and w µν A = w µν B = (µ 2 ν − ν 2 µ) in Ω, we obtain the following function Ω 00 : where the coefficients can be shown to be non-negative for all n, m [42]. Then, an upper bound on M 00 is straightforwardly given by where we have lower bounded the term c nm M nm by zero since all the coefficients satisfy c nm ≥ 0.

b. Upper bound on M11
By substituting κ µ A = κ µ B = 1 and w µν A = w µν B = (µ 2 − ν 2 ) in Ω, we obtain the following function Ω 11 : where the coefficients where c max ≥ |c nm | for all the pairs (n, m) ∈ S − and S o = {(n, 0)|n ≥ 0} ∪ {(0, m)|m ≥ 0}. In Eq. (B15), the second inequality comes from the fact that we have set to zero all those terms M nm , with (m, n) = (1, 1), which do not belong to S − nor to S o because M nm ≥ 0, ∀n, m. A valid c max can be obtained by noticing that, for n > s, This means that, from Eqs. (B14) and (B16), c max is given by Finally from Eqs. (B5) and (B15), an upper bound on M 11 is given by where M L 0A and M L 0B are lower bounds on the quantities M 0A = ∞ m=0 M 0m and M 0B = ∞ n=0 M n0 , respectively, and we have lower bounded the term ∞ n,m=3 c nm M nm by zero. Since M 0A and M 0B depend only on a single emitter, we can estimate them using the same method as for the vacuum component in BB84. Using the results of [44], we have that pν , M µ = ν M µν and M ν = µ M µν , with µ ∈ {µ 0 , µ 1 , µ 2 } and ν ∈ {ν 0 , ν 1 , ν 2 }; and the upper and lower bounds included in Eqs. (B19) and (B20) are obtained accordingly to Eq. (E1).

Supplementary Note C: Channel model
For our simulations, we use the channel model of [28], which we summarize here. We model the overall loss between Alice (Bob) and Charlie by a beamsplitter of transmittance √ η, which includes the channel transmissivity and the quantum efficiency of Charlie's detectors. We consider that the quantum channels connecting Alice and Bob with Charlie introduce both phase and polarisation misalignments. We model the phase mismatch between Alice and Bob's pulses by shifting Bob's signals by an angle φ = δ ph π. We model polarisation misalignment as a unitary operation that transforms Alice's (Bob's) polarisation input mode a † in (b † in ) into the orthogonal polarisation output modes a † out and a † out⊥ (b † out and b † out⊥ ) as follows: . The rotation angles are assumed to be θ A = −θ B = arcsin δ pol . With this channel model, it can be shown [28] that the probability that Charlie reports a successful detection, given that both users employ the X basis, is given by where γ = √ ηα 2 , θ = θ A − θ B , and Ω(φ, θ) = cos φ cos θ. The probability that Alice and Bob end up with different key bits is given by e X = e −γΩ(φ,θ) − (1 − p d )e −γ e −γΩ(φ,θ) + e γΩ(φ,θ) − 2(1 − p d )e −γ , while the probability that Charlie reports a successful detection, given that both users employ the Z basis and select the intensities µ and ν, respectively, is where I 0 (z) = 1 2πi e (z/2)(t+1/t) t −1 dt is the modified Bessel function of the first kind. In our simulations, we assume that the observed measurement counts equal their expected value, that is, we set M X = N p 2 X Q X and M µν = N p 2 Z p µ p ν Q µν , where M µν denotes the number of successful rounds in which Alice and Bob select the Z basis and the intensities µ and ν, respectively. Also, we assume that the bit-error rate of the sifted-key equals the probability given by Eq. (C2).
Supplementary Note D: Proof of Equation (13) Let us consider the evolution that the initial quantum state |Φ = |φ ⊗N , where |φ = |ψ AcAa ⊗ |ψ BcBb and |ψ is given by Eq. (11), experiences before step (4.1) in Box 3, by taking into account all operations applied to it. After Eve's measurement in step 2 of Box 3, it is transformed toM eve |Φ , whereM eve is the operator associated with her outcome. Let us reorder |Φ as |Φ = |φ ⊗M |φ ⊗M , writing first (last) the M (M ) successful (unsuccessful) rounds. In the virtual sifting step, Alice and Bob measure all subsystems A c and B c , using measurement operators Using an identical approach, we can show that the probability that Alice and Bob will learn that they used the X basis and sent cat states |C i |C j in the u-th successful round is Pr (ξ u = X ij |ξ, ξ 1 , ..., ξ u−1 ) = p 2 X Ê u |C i a |C j b 2 Pr(ξ, ξ 1 , ..., ξ u−1 ) .
Now, we want to relate the probability terms on the left hand side of Eqs. (D6) and (D7). For this, we use the approach of [28] and apply the Cauchy-Schwartz inequality to show that (D8) Combining the three previous equations, we obtain except with a probability ε c , where δ L and δ U are the solutions of the following equations The solutions to Eq. (E2) satisfy [51] 1 1 + δ L = W 0 (−e ln(εc/2−χ)/χ ), where W 0 and W −1 are branches of the Lambert W function, which is the inverse of the function f (z) = ze z .