A nonunitary interpretation for a single vector leptoquark combined explanation to the B-decay anomalies

In order to simultaneously account for both RD(∗) and RK (∗) anomalies in B-decays, we consider an extension of the Standard Model by a single vector leptoquark field, and study how one can achieve the required lepton flavour non-universality, starting from a priori universal gauge couplings. While the unitary quark-lepton mixing induced by SU(2)L breaking is insufficient, we find that effectively nonunitary mixings hold the key to simultaneously address the RK (∗) and RD(∗) anomalies. As an intermediate step towards various UV-complete models, we show that the mixings of charged leptons with additional vector-like heavy leptons successfully provide a nonunitary framework to explain RK (∗) and RD(∗) . These realisations have a strong impact for electroweak precision observables and for flavour violating ones: isosinglet heavy lepton realisations are already excluded due to excessive contributions to lepton flavour violating Z -decays. Furthermore, in the near future, the expected progress in the sensitivity of charged lepton flavour violation experiments should allow to fully probe this class of vector leptoquark models.


JHEP12(2019)006
2.1 Accounting for R K ( * ) and R D ( * ) in a minimal V 1 leptoquark framework 5 2.2 Vector-like fermions and "effective" nonunitary mixings in the light sector 5 3 Explaining LFUV data with nonunitary couplings: phenomenological viability 8 3.1 New contributions to R K ( * ) and R D ( * ) 8  In the Standard Model (SM), gauge interactions are strictly flavour universal, as confirmed by precision measurements of several electroweak observables, such as Z → decays [1,2].
Although these anomalies by no means invalidate the SM at this stage, their persistence and relatively coherent pattern inevitably raise the question of which (minimal) new ingredients beyond the SM (BSM) would be required to explain them. A first modelindependent approach [19,20,[26][27][28][29][30][31][32][33][34][35][36][37][38][39] is to introduce higher dimensional effective operators, coupling two quarks with two leptons. Despite the large number of possibilities, it is nevertheless remarkable that only a reduced number of such non-standard couplings significantly eases the tensions with the SM predictions.
In this work, we focus on the exchange of a vector leptoquark V 1 transforming as (3, 1, 2/3) under the SM gauge group, which has been shown to be particularly attractive for its ability to provide a single particle solution simultaneously to both charged and neutral current anomalies [96][97][98][99][100][101][102][103][104][105][106][107][108]. Complying with the experimental measurements suggests that V 1 should have non-universal couplings to quarks and leptons. We assume V 1 to be an elementary spin-1 gauge boson; since it carries charges (as a leptoquark must), the underlying gauge symmetry is necessarily non-abelian, with universal (gauge) couplings as long JHEP12(2019)006 as it remains unbroken. 1 As an example, such a field V 1 is naturally contained within the theoretically well-motivated Pati-Salam model (PS) as an SU(4) gauge boson. However, the current bounds on the charged lepton flavour violating (cLFV) decays K L → µe and K → πµe lead to dramatic (lower) bounds on the mass of such a vector leptoquark (m V ), typically above the 100 TeV scale for O(1) couplings [110][111][112][113][114][115]. In turn, this renders the new state excessively heavy to account for the B-meson decay anomalies.
The cLFV bound on m V turns out to effectively preclude a viable solution to both charged and neutral current anomalies: in the unbroken phase, V 1 has a single universal coupling to matter; SU(2)-breaking introduces a possible misalignment of the quark and lepton mass eigenbases, thus resulting in LFU-violating V 1 couplings, proportional to a 3 × 3 unitary matrix. In order to explain the R D ( * ) anomalies, the bτ and sτ couplings are required to be large, 2 which in turn leads to large couplings between the first two generations of quarks and leptons (as a consequence of the unitarity of the mixing matrix), leading to excessive contributions to cLFV. The question that naturally emerges is whether one can find a minimal embedding of V 1 that successfully allows to overcome the cLFV constraints and address both R K ( * ) and R D ( * ) anomalies. In other words, can one go beyond the tight constraints arising from a (3 × 3) unitary mixing of quarks and leptons? This necessarily requires the addition of new fields, beyond V 1 , and along these lines, one possibility is to add other vector leptoquark fields (thus implying a larger gauge group), whose mixing would allow to overcome the above mentioned constraints, as explored in [99].
In the present study, we avoid this further enlarging of the gauge group, adhering to the single vector leptoquark hypothesis, and pursue a distinct avenue. In particular, and motivated by the phenomenological impact of having nonunitary left-handed leptonic mixings in the presence of (heavy) sterile neutral leptons [117][118][119][120], we consider the possibility of nonunitary V 1 couplings, as arising from the presence of n additional vector-like heavy leptons L (also present in the construction of [99]). In the broken phase, the V 1 couplings are then given by a (3 + n) × (3 + n) mixing matrix, so that the couplings to SM fermions now correspond to a 3 × 3 sub-block, which is no longer unitary. We argue that this departure from unitary mixings might indeed hold the key to simultaneously address R K ( * ) and R D ( * ) data, while satisfying existing cLFV constraints.
The addition of vector-like heavy charged leptons 3 can be seen as an intermediary step towards a full ultraviolet-complete model, providing a better framework for the pecu- 1 In this work we are interested in the minimal (gauge extension) scenario where the vector leptoquark is an elementary gauge boson corresponding to a gauge group under which the SM fermion generations are universally charged and no additional protection or symmetry is introduced to induce non-universality. For models in which the vector leptoquark appears as a composite field, see for instance [109]; for other models where the gauge group is non-minimal and/or the gauge charges of the SM fermion generations are non-universal, see e.g. [70]. 2 To satisfy the constraints from the τ decays, the cν coupling induced by bτ via CKM mixing is in general not sufficient to comfortably explain R D ( * ) [116]; on the other hand the maximum cν coupling induced by diµ and die (for neutrino flavour in cν different from ντ ) are fixed by R

JHEP12(2019)006
liar structure of leptoquark couplings required by the anomalies. In this framework, the nonunitary mixings will also lead to the modification of SM-like charged and neutral lepton currents, establishing an inevitable link to electroweak precision (EWP) observables, such as lepton flavour violating and/or LFUV Z-decays. The latter observables will prove to be extremely constraining, ultimately leading to the exclusion of isosinglet vector-like heavy leptons as a source of non-universality in B-meson decays. These constraints are much milder for isodoublet heavy leptons: after arguing that for a single additional heavy charged lepton, cLFV constraints exclude an explanation of even R D ( * ) alone, we show that the addition of n = 3 vector-like isodoublet leptons allows a simultaneous explanation of both R K ( * ) and R D ( * ) anomalies, while respecting all available constraints.
This work is organised as follows: in section 2 we describe the underlying framework; section 3 is devoted to a comprehensive analysis of the phenomenological implications of the nonunitary framework, regarding the B-meson anomalies, and several flavour and EWP observables. A summary and concluding remarks can be found in section 4.

Towards a nonunitarity interpretation of vector leptoquark couplings
As mentioned in the Introduction, we consider here a SM extension by a single vector leptoquark V 1 , which transforms under the SM gauge group SU(3) c × SU(2) L × U(1) Y as (3, 1, 2/3). Without loss of generality, we assume that V 1 is a gauge boson of an unspecified gauge extension of SU(3) c with a universal (i.e. flavour blind) gauge coupling; without relying on a specific gauge embedding and/or Higgs sector, our only working assumption is that all fermions acquire a mass after electroweak symmetry breaking (EWSB), and that the physical eigenstates are obtained from the diagonalisation of the corresponding (generic) mass matrices. In the weak basis, the interaction of V 1 with the SM matter fields can be written as in which the "0" superscript denotes interaction states, and i = 1 − 3 are family indices. The couplings κ L,R are flavour diagonal, and universal. Since left-handed couplings are the minimal essential ingredient frequently called upon to simultaneously explain the neutral and charged current anomalies [97], for simplicity we will henceforth only consider the latter (i.e., taking κ L = 0 and κ R =κ R = 0). Furthermore, notice that this can be easily realised in chiral PS models [97,104,105], and is moreover phenomenologically well-motivated. 4 In terms of physical fields, the Lagrangian can be written as In the context of PS unification it has been noted in the literature that if the vector leptoquark couples to both left-and right-handed fermion fields with similar gauge strength, then in the absence of some helicity suppression, bounds from various searches for lepton flavour violating mesonic decay modes put a lower limit on the vector leptoquark mass around 100 TeV [111][112][113][114][115].

JHEP12(2019)006
where V is the Cabibbo-Kobayashi-Maskawa (CKM) mixing matrix and U P ≡ U † L U ν L the Pontecorvo-Maki-Nakagawa-Sakata (PMNS) leptonic mixing matrix; we have also introduced K 1 ≡ κ L √ 2 U L to denote the "effective" leptoquark couplings in the physical fermion basis. Being proportional to an arbitrary unitary matrix (which we hereby denote V 0 ), K 1 can be further cast as in which we used the standard parametrisation of a real 3 × 3 unitary matrix in terms of three angles θ 12,23,13 (with c ij and s ij respectively denoting cos θ ij and sin θ ij ), and we restrict ourselves to real parameter space in all our analysis. (We emphasise here that V 0 is not the PMNS matrix, and that the above angles are not those associated with neutrino oscillation data.) 2.1 Accounting for R K ( * ) and R D ( * ) in a minimal V 1 leptoquark framework The presence of a vector leptoquark, whose interactions with quarks and leptons are defined in eqs. (2.1), (2.2), can induce new operators, contributing to b-decays (both neutral and charged currents). In the SM, b → s and b → c ν decays respectively occur at one-loop and at tree-level; on the other hand, the new V 1 -mediated contributions to both decays arise at the tree-level. Thus, the new contributions required to explain R K ( * ) data are comparatively smaller than those needed to account for the discrepancy in R D ( * ) data: in particular, R D ( * ) requires the mass scale of V 1 to be quite low ∼ O(1 TeV), while it is possible to explain R K ( * ) for leptoquark masses m V ∼ O(10 TeV) (taking into account all the constraints from rare transitions and decays). The low mass scale required to explain the R D ( * ) anomaly effectively precludes a simultaneous (combined) explanation for both anomalies, due to the excessive associated contributions to cLFV kaon decays, in particular to K L → e ± µ ∓ (which occurs at the tree level). Consequently, both modes (K L → e + µ − and K L → e − µ + ) have to be suppressed separately. In terms of the parametrisation of eq. (2.3), saturating R D ( * ) requires maximising the 23 and 33 entries of V 0 (thus leading to θ 13 ∼ 0 and θ 23 ∼ π 4 ). This implies that the branching fractions of the tree-level kaon decay modes are proportional to sin 2 θ 12 and cos 2 θ 12 , respectively. A sufficient and simultaneous suppression of contributions to these modes is then clearly not possible.
Likewise, excessively large contributions to µ−e conversion (also occurring at tree-level) further exclude a low scale realisation, with m V ∼ O(1 TeV). The above arguments are illustrated by figure 1, in which we display the predictions for neutrinoless µ − e conversion and K L → e ± µ ∓ associated with having contributions to R D ( * ) within 3σ of the current best fit (for m V ∼ O(1 TeV) and three different values of κ L / √ 2).

Vector-like fermions and "effective" nonunitary mixings in the light sector
The above discussion suggests that the minimal flavour structure encoded in the (unitary) parametrisation of the leptoquark-quark-lepton currents (eq. (2.2)) is insufficient to account  2), and the dotted line a benchmark future sensitivity to CR(µ − e, Al). All mixing angles have been varied randomly between −π and π and the leptoquark mass is set to m V ∼ 1.5TeV. The blue, orange and green points respectively correspond to three benchmark choices, κ L √ 2 = 1, 0.5, 0.1.
for both anomalies. A stronger enhancement of LFUV in the leptoquark couplings can be achieved if one hypothesises that the "effective" leptoquark mixings -i.e. the 3 × 3 matrix V 0 is nonunitary. As we proceed to discuss, in order to explain R K ( * ) and R D ( * ) data simultaneously, and for universal gauge couplings, a highly nonunitary flavour misalignment between quarks and leptons is in fact required.
Such a nonunitary flavour misalignment can be understood in the presence of heavy vector-like fermions, SU(2) L singlets or doublets, which have non-negligible mixings with the SM fermions. This can be encoded by generalising the charged lepton mixing matrix to a 3 × (3 + n) semi-unitary matrix, so that SM interaction fields and physical states are related as 0 L = U L L (for n additional heavy states). The Lagrangian of eq. (2.2) can thus be recast as Notice that in the above equation, the effective leptoquark coupling K L generalises K 1 of eq. (2.3), and now corresponds to a rectangular 3 × (3 + n) matrix, which can be written in terms of U L as K L ≡ κ L √ 2 U L . Finally, K L can be further decomposed as K L = (K 1 , K 2 ), so that K 1 can be now identified with the nonunitary mixings in the light sectors (contrary to the simple limit of eq. (2.3)). K 2 is a 3 × n matrix which corresponds to the n heavy degrees of freedom describing the coupling parameters of the heavy (vector-like) states. Inspired by the approach frequently adopted in the context of neutrino physics, the deviation from unitarity JHEP12(2019)006 in the K 1 block can now be parametrised as [117][118][119][120] with V 0 given in eq. (2.3). The left-triangle matrix A, characterises the deviation from unitarity and encodes the effects of the mixings with the heavy states. As already mentioned in the Introduction, we assume that the vector leptoquark V 1 appears as a gauge boson in an unspecified SU(3) c extension. Since neither the gauge embedding nor the Higgs sector is explicitly specified, our only assumption is that after EWSB all fermions (SM and vector-like) are massive, and that the physical eigenstates are obtained from the diagonalisation of an (effective) generic (3 + n) × (3 + n) lepton mass matrix. For simplicity (see section 3.3), we take n = 3 generations of heavy leptons in what follows; the 6 × 6 charged lepton mass matrix M can be diagonalised by a bi-unitary transformation Being a unitary 6 × 6 matrix, U L can be parametrised by 15 real angles and 10 phases, and cast as a the product of 15 unitary rotations, R ij . By choosing a convenient ordering for the products of the complex rotation matrices, one can establish a parametrisation that allows isolating the information relative to the heavy leptons in a simple and compact form. Schematically, this can be described by the following (2 × 2 block matrix) decomposition [117], to which we adhere for the remainder of our discussion, further defining (2.8) Under the above decomposition, one can still identify the SM-like mixings, given by V 0 (cf. eq. (2.3)); the leptoquark couplings 5 are now parametrised by the 3×6 (rectangular) matrix, The diagonal elements of the triangular matrix A, α ii , can be expressed as

JHEP12(2019)006
in terms of the cosines of the mixing angles, c ij = cos θ ij . (The SM-like limit can be recovered for A → 1.) The off-diagonal elements can be cast as [117] whereŝ ij ≡ e iδ ij sin θ ij , with θ ij and δ ij respectively being the angles and CP phases associated with the R ij rotation. Finally, it is worth emphasising that not only the full 6 × 6 matrix U L is unitary, but its upper 3 × 6 block (AV 0 , R) is also semi-unitary on its own, with 2 This formalism, which can be easily generalised to n extra generations, offers the possibility of successfully separating the information relative to the heavy leptons in a simple and compact form. Although the couplings (in particular the α ij entries) can be in general complex, in what follows we consider a minimal scenario where all couplings are taken to be real.

Explaining LFUV data with nonunitary couplings: phenomenological viability
We recall that, as mentioned in the Introduction, we work under the minimal assumptions that the singlet vector leptoquark V 1 (colour triplet) should correspond to a gauge extension of SU(3) c unifying quarks and leptons with a universal (i.e., flavour independent) gauge coupling. We first describe the effects of the vector leptoquark on the neutral and charged current b decays, and then summarise the most stringent constraints arising from numerous flavour violating and flavour conserving observables (meson oscillations and decays, as well as charged lepton flavour violation processes). We then present our main numerical results.

New contributions to R K ( * ) and R D ( * )
In what follows, we proceed to explore whether the relaxation of the unitarity requirement on the SU(2) L -singlet vector leptoquark V 1 couplings to SM matter does allow addressing R K ( * ) and R D ( * ) data simultaneously.
Anomalies in neutral current b decays: R K ( * ) . As mentioned in the Introduction, several measurements of the ratio of branching ratios of B → K ( * ) ( = e, µ) exhibit tensions when compared to the SM predictions. The most recent averages (and SM estimations) are associated with the following deviations [11][12][13], in which the dilepton invariant mass squared bin (in GeV 2 ) is identified by the subscript:

JHEP12(2019)006
Other anomalies in the neutral current mode of B meson decays have also emerged concerning the observable Φ ≡ dBR(B s → φµµ)/dm 2 µµ in the analogous bin (m 2 µµ ∈ [1, 6] GeV 2 ) [14], with a similar deviation (at a level of approximately 3σ), as well as in the angular observable P 5 in B → K * + − processes. While LHCb's results for P 5 in B → K * µ + µ − decays manifest a slight discrepancy with respect to the SM, the Belle Collaboration [15] reported that, when compared to the muon case, P 5 results for electrons show a better agreement with theoretical SM expectations. Nevertheless, there has been an ongoing discussion about the possibility that incorrectly estimated hadronic uncertainties might be at the origin of the observed anomalies in the mode B → K * µ + µ − , due to power corrections to the form factors, or charm-loop contributions [123][124][125][126].
The effective Hamiltonian describing the neutral current effects at the level of quark transitions d j → d i − + is given in appendix A. A model-independent analysis of the current data at the b-quark mass scale can be made using the 2D hypothesis C µµ 9,NP = − C µµ 10,NP and C ee 9,NP = − C ee 10,NP , allowing for the possibility of (V − A) new physics (NP) effects in both electron and muon channels. Using the available experimental measurements of R K ( * ) in different high and low q 2 bins, including the latest updates from LHCb and Belle collaborations, the available experimental measurements on the angular observables for b → sµµ and b → see, and the latest average for BR(B s → µµ) -which already includes the latest measurement from the ATLAS collaboration -we find the global fit ranges at the 2 σ (1 σ) level. In figure 2 we show the 1σ likelihood contours from the latest experimental measurements of R ( * ) K , b → sµµ and b → see observables, as well as their combined global fit (1σ and 2σ) in the C µµ 9,sb vs. C ee 9,sb plane, assuming the 2D hypothesis C µµ 9,sb = −C µµ 10,sb and C ee 9,sb = −C ee 10,sb . The vector leptoquark V 1 contributes at tree level, yielding the following Wilson coefficients at the leptoquark mass scale for the O ij; 9,10 operators [127] in which G F denotes the Fermi constant, α em is the fine-structure constant, V the CKM matrix, K L the "effective" leptoquark couplings (cf. eq. (2.9)) and m V the vector leptoquark mass. The matching of the model parameters with the Wilson coefficients is performed at the scale of leptoquark mass and the Wilson coefficients are subsequently run down to the b-quark mass scale [116]. In particular, it is interesting to note that due to RG-running a large K iτ L coupling can potentially induce a non-negligible lepton-universal contribution to b → s transitions via a log-enhanced anapole photon penguin contribution, as noted in ref. [128]. The dominant log-enhanced contribution is given by K , b → sµµ, and b → see observables, as well as the combined global fit (1σ and 2σ) in the C µµ 9,sb vs. C ee 9,sb plane, assuming the 2D hypothesis C µµ 9,sb = −C µµ 10,sb and C ee 9,sb = −C ee 10,sb . The global fit is obtained at the b-quark mass scale using the package flavio [122].
In our analysis, the running from the scale of leptoquark mass to the b-quark mass scale, and to any other relevant process (observable) scale, is taken into account using the wilson package [121] in association with the flavio package [122].
Anomalies in charged current b → c ν transitions: R D ( * ) . Important deviations from the SM prediction of lepton flavour universality in B → D ( * ) ν decays have also been reported by several experimental collaborations. The most recent averages for the R D and R D ( * ) ratios by the HFLAV Collaboration [16,17] are We define the effective Hamiltonian for the charged current transitions d k → u jν − as where, in the SM, C f i jk,SM = δ f i . The contribution from the vector leptoquark V 1 is given by One can further construct the double ratios Figure 3. Likelihood 1σ and 2σ contours from R K ( * ) , R D ( * ) , b → sµµ observables and the combined global fit (1σ, 2σ and 3σ) in the plane K 33 L , defined at the leptoquark mass scale cf. eq. (2.2). The first generation lepton and quark couplings are set to zero.
(equal to unity when NP decouples, i.e., κ L → 0 and m V → ∞). After combining current experimental world averages with the SM predictions, the current anomalous data can be summarised as R D /R SM D = 1.14 ± 0.10 , R D * /R SM D * = 1.14 ± 0.06 , in which the statistical and systematical errors have been added in quadrature. In figure 3, we display the 1σ and 2σ likelihood contours from R K ( * ) , R D ( * ) , b → sµµ observables and the combined global fit (1σ, 2σ and 3σ) in the plane K 33 where the couplings (see eq. (2.2)) are varied independently and the others are set to zero.
As can be seen in figure 3, both b → c ν and b → s anomalies can be indeed accommodated simultaneously in a minimal V 1 model. However, we will subsequently show how this is realised in the nonunitary framework taking all couplings into account, as suggested by the b → s data.

Constraints from (rare) flavour processes, EW precision observables and direct searches
The extended framework called upon to address the B meson decay anomalies -not only the additional vector leptoquark, but also the presence of extra vector-like fermions, which are the origin of the nonunitarity of the V 1 effective couplings -opens the door to extensive contributions to numerous observables. While most of the NP contributions occur via higher order (loop) exchanges, it is important to notice that V 1 can also mediate very rare (or even SM forbidden) processes already at the tree level. As we proceed to discuss, the latter observables prove to be particularly constraining, and put stringent bounds on the degrees of freedom of these leptoquark realisations.
Leptoquark SM extensions aiming at addressing the anomalies in R K ( * ) and R D ( * ) data receive strong constraints from d j → d iν ν transitions (in particular s → dνν and b → sνν). However, the vector leptoquark V 1 does not generate contributions at tree level, and the first non-vanishing contribution appears at one loop. Consequently, we find that even with significant uncertainties, the semileptonic decays into charged dileptons d j → d i Observables SM prediction Experimental data Table 1. Relevant observables and current experimental status for leptonic and semi-leptonic meson decays; when appropriate, the associated SM prediction is also included.
violating modes). In the present analysis, we therefore include bounds from K → π and B → K , as well as the stringent limits on NP contributions arising from the observed decay mode B s → µ + µ − . Furthermore, we also take into account the lepton flavour violating leptonic decays B → e ± µ ∓ , B s → e ± µ ∓ and K L → e ± µ ∓ . Leptoquark contributions to neutral meson oscillations and mixings, such as K 0 −K 0 and B 0 s −B 0 s , are also evaluated. Table 1 provides a brief summary of the current experimental status for these mesonic observables (current bounds and future sensitivities). A detailed discussion of the formalism used to evaluate the vector leptoquark contributions is provided in appendix A.
The lepton flavour non-universal couplings of vector leptoquarks (in general nonunitary in the present framework) induce new contributions to cLFV observables: radiative decays JHEP12(2019)006 i → j γ and 3-body decays i → 3 j at loop level, and neutrinoless µ − e conversion in nuclei both at tree and loop level. Further taking into account the impressive associated experimental sensitivity, it is clear that these observables lead to important constraints on the vector leptoquark couplings to SM fermions. It is important to stress that although the radiative decays are generated at higher order, relevant anapole contributions can add to the Wilson coefficients accounting for the tree-level contributions to neutrinoless µ − e conversion and µ → 3e. The higher order anapole contributions can have a magnitude comparable to the tree level ones (or even account for the dominant contribution). In addition, dipole operators also contribute significantly to radiative decays and to neutrinoless µ − e conversion. Although we do take tauonic modes into account, we notice here that due to the associated current experimental sensitivity, (semi)leptonic tau decays in general lead to comparatively looser constraints; likewise, semileptonic meson decays into final states including tau leptons are typically less constraining. However, the expected improvements in sensitivity from dedicated experiments might render the tau modes important probes of SM extensions via vector leptoquarks. 6 As done for the first time in this work, all these contributions must be systematically included to thoroughly constrain the vector leptoquark couplings.
A summary of the current experimental status (current bounds and future sensitivities) is given in table 3; for simplicity, in the numerical analysis we consider a benchmark future sensitivity to µ − e conversion in Aluminium of O(10 −17 ). The relevant details of the computation of the cLFV observables considered in this work can be found in appendix B.
The presence of heavy vector-like fermions (at the source of the nonunitary couplings of the vector leptoquark to the light fermions) can have a non-negligible impact on the couplings of SM fermions to gauge bosons. In turn, this can be manifest in new contributions to several EW precision observables -potentially in conflict with SM expectations and precision data, -which will prove to play a key role in constraining the mixings of the SM charged leptons with the heavy states.
For the Z-couplings, which are modified at the tree level, the most stringent constraints are expected to arise from leptonic Z decays; in our analysis, we take into account the LFU ratios and cLFV decay modes of the Z boson. We summarise in table 2 the EWP observables which are of relevance for our study (experimental measurements and SM predictions). A discussion on how thef f Z couplings are modified as a consequence of nonunitary effective leptoquark couplings, as well as their impact for several observables, is detailed in appendix C.
Finally, it is clear that the (negative) results of direct searches for the exotic states must be taken into account. At the LHC, pairs of vector leptoquarks can be abundantly produced in various processes (via t-channel lepton exchange and direct couplings to one or two gluons). Due to the underlying gauge structure of possible ultraviolet (UV) completions, the production cross section strongly depends on the coupling to gluons which makes the theoretical predictions for the production of vector leptoquarks less robust than those for scalar leptoquarks.   Table 3. Current experimental bounds and future sensitivities of various cLFV processes considered in the analysis.
spontaneously broken non-abelian gauge symmetry, gauge invariance completely fixes the couplings between the vector leptoquark and the gluons, implying lower limits on the vector leptoquark mass (albeit still depending on the branching fractions of the leptoquark) from the negative results in the direct searches for pair production at the LHC, see for instance [103]. As a natural consequence of the favoured structure called upon to maximise the effects on B-meson decay anomalies, V 1 is expected to dominantly decay into either tν τ or bτ . The ATLAS and CMS collaborations have conducted extensive searches, assuming that leptoquarks couple exclusively to third generation quarks and leptons [151][152][153][154][155], which have led to lower limits on the V 1 mass ∼ 1.5 TeV. Further important collider signatures are pp → ττ + X, arising from t-channel leptoquark exchange or from single leptoquark production in association with a charged lepton. As argued in [156,157], the projected vector leptoquark reach of HL-LHC with 3 ab −1 is close to 1.8 TeV. Much higher luminosi-

JHEP12(2019)006
ties and/or more sophisticated search strategies are required to probe the preferred mass scale and couplings of states at the origin of a combined explanation of the anomalies (for instance, as suggested by the study of [97]). Other potentially interesting search modes could include bµbτ and bµbµ final states.
In the present analysis, we select (working) benchmark values for the mass and gauge coupling of the vector leptoquark allowing to comply with the current available limits. In particular, for the following numerical analysis we set κ L √ 2 = 1 as a benchmark choice. Nevertheless, for any other choice consistent with the constraints from direct searches (as discussed above) the qualitative behaviour and the conclusions drawn remain the same. However, for very small values κ L (κ L 0.1 for m V 1.5 TeV) the number of points in the best fit region for R K ( * ) and R D ( * ) anomalies starts to decrease drastically, so that the NP effects become negligible with respect to the SM.

Results and discussion
We now finally address the question of whether a SM extension via vector leptoquarks can simultaneously provide an explanation to both the R K ( * ) and R D ( * ) data, working under the hypothesis of universal gauge couplings for the vector leptoquark V 1 . In this framework, the required flavour non-universality arises from nonunitary mixings among the SM fermions -a consequence of the existence of heavy vector-like fermions which have non-negligible mixings with the light fields.
We first begin by considering the most minimal scenario where the nonunitary flavour misalignment is due to the presence of a single generation of heavy vector-like charged leptons (i.e., n = 1 in eq. (2.4)). Such a minimal field content already leads to a sufficient amount of LFUV to account for both R K ( * ) and R D ( * ) anomalies. Although new contributions to rare meson decays and transitions are still in good agreement with current experimental bounds, this scenario is ruled out due to the stringent constraints on cLFV 8 modes. Excessive contributions to (tree-level) muon-electron conversion in nuclei play a crucial role in ruling out this realisation, as well as the closely related radiative decays. In order to reconcile the model's prediction with the current bounds on CR(µ − e, Au), the photon-penguin contributions must (at least) partially cancel the sizeable tree-level ones; however, such large photon-penguin contributions then translate into unacceptably large µ → eγ decay rates, already in conflict with current bounds.
The required flavour non-universality can be recovered for a less dramatic unitarity violation; this can be achieved by extending the particle content by two or more additional heavy charged lepton states, or formally for n ≥ 2 in eq. (2.4). Although n = 2 provides more freedom to evade the constraints found in the case n = 1, no generic solution was found in this case, which might then require an extreme fine tuning to become viable. In the subsequent discussion, we therefore take n = 3 which is the minimal case and conveniently replicates the number of generations in the SM.  For n = 3, our study suggests that it is in general possible to find regions in the parameter space in which the required non-universal flavour structure to explain both R K ( * ) and R D ( * ) arises in a natural way, while still complying with all constraints from flavour violating processes (meson and lepton sectors). This can be seen from the upper row of figure 4, in which we display the regimes for the entries of the matrix A (cf. eqs. (2.7)-(2.11)) which account for both R K ( * ) and R D ( * ) data, as well as regions respecting the constraints arising from the several flavour violating modes considered in our analysis. Concerning the latter, we find it worth mentioning that the most stringent constraints arise from K L → µ ± e ∓ , µ → eγ, and µ − e conversion in nuclei; B-meson cLFV decays, or (semi)leptonic B and K decays lead to comparatively milder constraints (or are systematically satisfied).
Finally, one should address the compatibility of the considered SM leptoquark extension with the constraints arising from EW precision tests; as discussed in the previous subsection (and in appendix C), nonunitary mixings can modify the couplings of the Z boson. In particular, for isosinglet heavy leptons, the entries of the matrix A (see eqs. (2.10) and (2.11)) are severely constrained by the Z width and by bounds on its cLFV decays (Z → ). This is shown on the lower row of figure 4, which illustrates the tension between LFUV and Z bounds -a tension which ultimately leads to disfavouring this class of extensions as a phenomenologically viable NP model to explain both R K ( * ) and R D ( * ) discrepancies.
If one foregoes a solution to the charged current anomalies (i.e., R D ( * ) ), it is possible to accommodate R K ( * ) in full agreement with constraints from flavour bounds, relying on very JHEP12(2019)006 mild deviations from unitarity, and thus evading constraints from universality violation in Z decays. However, one would be led to regions with considerably heavier leptoquarks, m V 15 TeV. This is depicted on the left panel of figure 5, in which we display regimes complying with R K ( * ) at the 3σ level) in the plane spanned by two particularly constraining observables, BR(K L → µ ± e ∓ ) and CR(µ − e, N), for 15 TeV m V 45 TeV. As can be verified, a small subset of points (consistent with R K ( * ) and respecting universality in Z decays) is compatible with current bounds on the cLFV processes. This is in agreement with the analyses of various UV-complete models, such as [104,105].
For completeness, the right panel of figure 5 shows a similar study for R D ( * ) . In order to accommodate R D ( * ) data, lower leptoquark masses are required (in this case we have taken 1 TeV m V 6 TeV), and it is no longer possible to evade K L → µ ± e ∓ and µ − e conversion bounds while being consistent with leptonic Z-decay universality. The data displayed in the panels of figure 4 was obtained for vector leptoquark masses m V ∼ 1.5 TeV; analogous conclusions can be inferred for m V ∼ (1 − 3) TeV, albeit for different α ij ranges.
Since for the case of isosinglet leptons an explanation of R D ( * ) is excluded by bounds on Z decays, we now consider isodoublet heavy charged leptons. The nonunitarity in the couplings of the vector leptoquark to the SM charged leptons can simultaneously explain R K ( * ) and R D ( * ) data. Moreover, and by construction, in the case of isodoublet heavy charged lepton states Z couplings remain universal (in the absence of mixings between right-handed SM charged leptons and vector-like doublets, ∆g R = 0, see appendix C); nevertheless, flavour observables still play a crucial role, and are (as expected) responsible for severe constraints on the NP degrees of freedom.  Figure 6. On the left, ∆χ 2 distribution for the fit to R K ( * ) and R D ( * ) data (1σ) in the plane of the (K L ) ij couplings. All points comply with the different (flavour) constraints. On the right, regions in the plane spanned by CR(µ − e,N) and BR(K L → µ ± e ∓ ), accommodating both R K ( * ) and R D ( * ) (blue) and those in addition complying with LFV constraints (yellow). Both panels correspond to a heavy sector composed of three isodoublet vector-like charged lepton states, and to having set m V ∼ 1.5 TeV. The ∆χ 2 corresponds to the 1 σ-region around the best fit point.
The left panel of figure 6 offers a global view of this case, showing the ∆χ 2 distribution for the fit to R K ( * ) and R D ( * ) data, in the plane spanned by (K L ) ij "muon and tau couplings" (K 22 K 32 − K 23 K 33 ), marginalising over the other couplings. The leptoquark mass is set to m V ∼ 1.5 TeV. We stress that leading to this plot all couplings were determined by the underlying nonunitarity parametrisation (with all mixing angles randomly sampled); in particular, we have not set the leptoquark couplings to the first generation of quark and leptons to zero. The displayed points comply with all flavour bounds included in our study, as described in section 3.2.
The lowest ∆χ 2 region (dark red ellipsoid) suggests that the best fit scenario corresponds to new physics dominantly coupling to muons and taus. We stress that the patterns emerging from the ∆χ 2 distribution are not an artefact of some particular assumption imposed on the couplings, but rather the result of a very general scan over the full set of (mixing) parameters.
The right panel of figure 6 offers a projection of the viable points (displayed on the left panel) in the plane of the most constraining observables, CR(µ − e,N) and BR(K L → µ ± e ∓ ). It is interesting to notice that, to a very good approximation, most of the currently phenomenologically viable points lie within future reach of the upcoming muon-electron conversion dedicated facilities (COMET and Mu2e).
In the near future, and should the B-meson decay anomalies be confirmed, an explanation in terms of such a minimal leptoquark framework could be probed via its impact for cLFV observables, in particular µ−e conversion in nuclei. Although the cLFV bounds could be evaded by increasing the mass of the vector leptoquark, this would however prevent a viable explanation of the B-meson decay anomalies, especially of R D ( * ) .

JHEP12(2019)006 4 Concluding remarks
In this study we have considered a minimal SM extension via one vector leptoquark V 1 and n generations of heavy vector-like charged leptons, as a candidate framework to explain the current B-meson anomalies, R K ( * ) and R D ( * ) .
Minimal extensions by a single V 1 leptoquark are in general disfavoured due to the strong cLFV constraints on the (unitary) quark-lepton-V 1 couplings. Here we have suggested that the pattern of mixings required to simultaneously address R K ( * ) and R D ( * ) with a single V 1 could be interpreted within a framework of nonunitary V 1 q couplings: the mixings of the SM charged leptons with the additional vector-like heavy leptons can lead to effectively nonunitary V 1 q couplings, offering the required amount of LFUV to account for both anomalies.
As we have argued, the most minimal nonunitary scenario (i.e. n = 1) consistent with both R K ( * ) and R D ( * ) , is ruled out as it leads to excessive contributions to cLFV observables such as muon-electron conversion in nuclei. We have thus considered three families of vector-like heavy leptons, and we have carried out a detailed analysis of the impact for an extensive array of flavour violation and EW precision observables. Our findings revealed that the SU(2) L charges of the heavy charged leptons are of paramount importance for the model's viability: for isosinglet heavy leptons, the mass of the leptoquark must be sufficiently large to avoid excessive contributions to Z decays, which then prevents an explanation of R D ( * ) . This is expected to happen for heavy leptons of any SU(2) L representation except for isodoublets. In this particular case, the Z couplings remain universal, and we have shown that the nonunitarity in the couplings allows to successfully explain both sets of anomalies, while complying with all considered current bounds.
The nonunitary framework here considered in fact provides one of the most economical options to simultaneously address R K ( * ) and R D ( * ) data via a single vector leptoquark. Although our results arise from a phenomenological analysis, our findings can be easily related to well-motivated UV complete models, in which the nonunitarity framework allows to obtain the required flavour structure without further enlarging the gauge group. Interestingly, this is the case of the Pati-Salam model [98][99][100][101][102][103][104][105][106][107], in which the existence of a vector leptoquark with a mass around the TeV scale (as required for a successful simultaneous explanation of both the neutral and charged current anomalies) naturally motivates the existence of such vector-like heavy leptons.
Finally, it is worth stressing that in view of the important progress expected in the near future, cLFV searches can play a crucial role in falsifying these minimal leptoquark frameworks, in view of their strong impact for cLFV observables, in particular µ − e conversion in nuclei. As we have discussed, should the B-meson decay anomalies be verified, and should this be a valid framework to explain them, then one should expect the observation of cLFV signals, in particular from µ − e conversion in Aluminium targets, either at COMET or Mu2e.

JHEP12(2019)006
Acknowledgments We acknowledge support within the framework of the European Union's Horizon 2020 research and innovation programme under the Marie Sklodowska-Curie grant agreements No 690575 and No 674896.

A Constraints from flavour violating rare meson decays and neutral meson mixing
New physics frameworks including leptoquarks give rise to important contributions to observables in the meson sector. These include flavour conserving and flavour violating leptonic and semileptonic decays 9 (including final state neutrinos), as well as neutral meson mixings. In this appendix we discuss the different leptoquark contributions to leptonic and semi-leptonic meson decays (arising at tree-level) and to modes with final state neutrinos (at one-loop level). We recall that the SM predictions and current experimental limits are summarised in table 1.
Note that the matching of the model's parameters with the Wilson coefficients is performed at the scale of the leptoquark mass, and then these are run down to the b-quark mass scale. In our analysis, the running from the scale of leptoquark mass to the b-quark mass scale, and to any other relevant process (observable) scale, is taken into account using the wilson package [121] in association with the flavio package [122].
The effective Hamiltonian for d j → d i − + transitions, including the LFV operators, can be cast as [158][159][160][161][162][163] (A.1) 9 We do not consider the Bs → µµ decay explicitly, since it is already included in deriving the global ranges for the new physics contributions to b → s anomalies.

JHEP12(2019)006
In the above, we recall that V ij denotes the CKM matrix; the operators are defined as with the primed operators O 7,9,10,S,P following from the replacement P L ↔ P R .
Regarding lepton flavour violating processes, only the operators O ( ) 9,10,S,P play a relevant role: for vector leptoquarks, the complete set of associated Wilson coefficients (not present in the SM) is given by [127] C ij; where K i L (K i R ) denoting the left-handed (right-handed) leptoquark couplings. Here, and for completeness, we include the K R mixing matrix, although we do not take it into account in our analysis.

A.2 P → − + decays
Leptonic decays of pseudoscalar mesons lead to important constraints on the (vector) leptoquark couplings. Here, we summarise the computation of the P → − + rates, following the derivation of [164]. With the standard decomposition of the hadronic matrix element 0 |d j γ µ γ 5 d i | P (p) = i p µ f P , (A.4)

JHEP12(2019)006
in which f P denotes the P meson decay constant, the branching fraction 10 is given by where the Källén-function is defined as Since this FCNC transition is generated at tree-level from V 1 interactions, the modes leading to different final state lepton charge assignments must both be included and treated separately. Leptonic pseudoscalar meson decays thus provide extensive (and very tight) constraints on the leptoquark couplings.

A.3 P → P − + decays
Again working in the standard basis, the hadronic matrix elements are parametrised as In the above, the hadronic form factors f +,0,T (q 2 ) depend on the momentum transfer, which lies in the range (m + m ) 2 ≤ q 2 ≤ (m P − m P ) 2 . In our evaluation of the form factors we closely follow the results of [165]; since we focus on decays with heavy-to-light meson transitions, we further assume the scale to be µ = m d j , i.e. µ = m b . With the quantities

A.4 Loop effects in neutrino modes
Both s → dνν and b → sνν transitions are known to provide some of the most important constraints on NP scenarios aiming at addressing the anomalies in R K ( * ) and R D ( * ) data.
Following the convention of [166], at the quark level the |∆S| = 1 rare decays K + (K L ) → π + (π 0 ) ν ν and B → K ( * ) ν ν can be described by the following short-distance effective Hamiltonian for d j → d i ν ν transitions [167,168] in which i, j denote the down-type flavour content of the final and initial state meson, respectively. For vector leptoquarks, the contributions are generated at one loop-level and are a priori divergent; consequently, the calculation should be carried in a given gaugeembedding, including the corresponding would-be Goldstone modes. Here we follow the computation of [128], and the coefficient for d a → d fνi ν j is thus given by where M W is the mass of the W boson and m t the mass of the top quark. The branching fractions for the neutral and charged kaon decay modes are given by [129,169] BR(K ± → π ± νν) = 1 In this convention, λ is one of the Wolfenstein parameters related to the Cabibbo angle, and λ c = V * cs V cd and λ t = V * ts V td . The B → K ( * ) νν decay width has been derived in [166], leading to C SM,f i L,sb ≈ −1.47/s 2 W δ f i . Finally, it is convenient to express the BRs normalised to the SM predictions,

A.5 Loop effects in neutral meson mixing
In this scenario a contribution to |∆F | = 2 amplitudes is generated at one-loop level.
Contributions to neutral meson mixings, P −P with P = B 0 s , B 0 d , K 0 , arise both from SM box diagrams involving top quarks and W 's, and from NP box diagrams involving leptons and vector leptoquarks. These contributions can be described in terms of the following effective Hamiltonian for |∆F | = 2 transitions with {i, j} respectively denoting {b, s}, {b, d} or {d, s} for P = B 0 s , B 0 d or K 0 mesons. The |∆F | = 2 transitions are sensitive to the mass scale of the heavy vector-like fermions, and the widths scale proportionally to m 2 V (similarly to the SM contribution, itself proportional to m 2 t ). A complete evaluation of the contributions must further include the effects of the (physical) Higgs fields; therefore the computation of these observables requires specifying a particular UV completion. Nevertheless, it is possible to draw preliminary conclusions on the mass scales of the vector leptoquarks and heavy leptonic states (here denoted by M ) based on the new physics contribution to the diagrams involving V 1 . For example, taking P = B 0 s , one obtains [97,98] C NP .

JHEP12(2019)006
The |∆F | = 2 transitions thus lead to a (lower) bound on the heavy leptonic mass scales of around 500 GeV, while the vector leptoquark mass should lie above the TeV to keep new physics contributions to ∆M B s,d below O(10%), given the experimental constraints.

B Constraints from cLFV decays
Due to the LFUV couplings of the leptoquarks, sizeable contributions can be generated for cLFV observables, including radiative decays i → j γ, three-body decays i → 3 j , and neutrinoless µ − e conversion in nuclei, which then lead to important constraints. As already mentioned in section 3, although the radiative decays are generated at loop level, anapole operators can induce additional contributions to the Wilson coefficients (comparable to the tree-level contributions to neutrinoless µ − e conversion), and dipole operators can also be responsible to significant contributions to radiative decays and neutrinoless µ − e conversion.
Notice that the one-loop dipole and anapole contributions due to the exchange of vector bosons generically diverge, and a UV completion must be specified to obtain a convergent result in a gauge independent manner. We have computed the anapole and dipole contributions in the Feynman gauge, for which it is necessary to include the relevant contributions from the Goldstone modes. To consistently compute these contributions for vector leptoquarks, we make the minimal working assumption that the new state corresponds to a (non-abelian) SU(3) c gauge extension, whose breaking gives rise to a would-be Goldstone boson degree of freedom, subsequently absorbed by the massive vector leptoquark. We thus include this Goldstone mode (degenerate in mass with V 1 ) to obtain the gauge invariant (finite) form factors for the dipole and anapole contributions.

B.1 Radiative lepton decays
The relevant terms of the effective Lagrangian for radiative lepton decays i → j γ can be written as in which F µν is the electromagnetic field strength tensor. The corresponding Wilson coefficients are related to the form factors σ i j L(R) as where the form factors σ i j L(R) can be computed following the prescription of [170]. At the one-loop level, there are 10 distinct diagrams that must be included. The gauge-invariant amplitude can be decomposed as

JHEP12(2019)006
in which and q denote the photon polarisation and its momentum. The i → j γ decay width is then given by In the physical (mass) basis, the relevant part of the Lagrangian leading to the computation of σ L,R can be written as (cf.eq. (2.4)) with K L (K R ) the left-handed (right-handed) coupling matrix. The Goldstone (ϕ) interaction terms in the Lagrangian are then given by The form factors σ L,R can be cast in the following (compact) form: and N c the number of colours of the internal fermion. The loop functions are given by

B.2 Three body decays →
At the loop level, three body decays can receive contributions from photon penguins (dipole and off-shell "anapole"), Z penguins and box diagrams, arising from flavour violating interactions involving the vector leptoquark V 1 and quarks. The relevant low-energy effective

JHEP12(2019)006
Lagrangian inducing the four-fermion operators responsible for → decays can be written as [171,172] to which the photonic dipole terms entering in L i → j γ eff , cf. eq. (B.1), must be added; the corresponding coefficients parametrised by C i j L(R) have already been discussed in detail in the previous subsection. Neglecting Higgs-mediated exchanges, the off-shell anapole photon penguins, Z penguins and box diagrams will give rise to non-vanishing contributions to the above g 3 , g 4 , g 5 and g 6 coefficients. Note that in the large V 1 mass limit, the offshell anapole photon-penguin diagrams scale proportionally to |K| 2 ln(m 2 q /M 2 )/M 2 , in contrast with the contributions from the Z-penguins and box diagrams, which are (naïvely) proportional to |K| 2 m 2 q /M 4 and |K| 4 m 2 q /M 4 respectively [173]. Therefore we only include in our computation the log-enhanced photonic anapole contributions, in addition to the dipole ones. Neglecting right-handed couplings of the leptoquark as before, the only nonvanishing coefficients (at one-loop) are g 4 = g 6 . The relevant amplitude for the anapole contribution can be written as where M µ anapole can be parametrised in terms of a form factor F γ L as with q the off-shell photon momentum. In this convention the F γ L form factor is independent of q 2 . After performing the calculation in the Feynman gauge, we obtain (in the limit of vanishing external lepton masses) in which x i = m 2 d i /m 2 V and N c is the colour factor (corresponding to the coloured fields entering in the loop). Finally, the loop function f a (x) is given by The leptoquark-induced contributions to the 4-fermion operators are given by (B. 16) In the case of the → 3 decays, Q f = Q denotes the charge of the fermion pair at the end of the off-shell photon (in units of e). As an example, for the case of µ → 3e decays, one obtains the following branching ratio [171,172] BR(µ → eee) = 2 |g 3 | 2 + |g 4 | 2 + |g 5 | 2 + |g 6 | 2 + +8 e Re C µe R (2g * 4 + g * 6 ) + C µe L (2g * 3 + g * 5 ) + Similar expressions can be easily inferred for the other cLFV 3-body decay channels.

B.3 Neutrinoless µ − e conversion
In terms of the relevant effective Wilson coefficients, the general contribution to the neutrinoless µ − e conversion rate is given by [127] Γ µ−e,N = 2 G 2 in which the photonic dipole Wilson coefficients C i j L(R) have been given in eq. (B.1); the other non-vanishing coefficients, induced by the tree-level leptoquark exchange or arising from the photonic anapole contributions, are given by with Q d = − 1 3 and Q u = 2 3 . The values for the overlap integrals (D, V, S) are given in table 4 [174], and the scalar coefficients G can be found in [175]. We again emphasise here that the off-shell anapole contributions, often neglected in the literature, can have a contribution comparable to the tree-level leptoquark exchange, and therefore should be included for a thorough estimation of the rate of µ − e conversion in nuclei.

C Electroweak precision observables
As mentioned in the main body of the paper, the hypothetical heavy vector-like fermions can modify the couplings of SM fermions to gauge bosons 11 Table 4. Overlap integrals D, V, S and Γ capture for Gold and Titanium nuclei, as reported in [174] (tables I and VIII).
the theoretical predictions of the SM, which are mostly in remarkable agreement with EW precision data.

C.1 Couplings of the Z boson and photon
If the heavy vector-like fermion states are SU(2) L singlets, mixings with the light SU(2) L doublets can lead to modified couplings of the latter to the Z boson (f f Z). For the case of charged leptons, the relevant couplings can be obtained from the kinetic terms, where {j, k} denotes the physical fields and {a, b} = 1 . . . 6 the interaction states (with a ∈ {1, 2, 3} corresponding to SM fermions, and a ∈ {4, 5, 6} to the new heavy vector-like states). The covariant derivative associated with the charges of a given state a can thus be written where T 3 and Q respectively denote the weak isospin and the electric charge. Since the electric charge is the same for all lepton states (Q a = −1), the couplings of the photon are not modified. Let us now introduce the "effective" Z boson couplings, where T 3 L(R) is the weak isospin of a left-handed (right-handed) lepton. Should the SM fermions and the heavy states belong to the same SU(2) L representation, universality is trivially restored (by unitarity) for both g Z j k L,R effective couplings, and one recovers the SM universal couplings. For heavy isodoublet vector-like states, one has (g Z j k L ) eff = g cos θ W − 1 2 + sin 2 θ W δ jk , (C.4)

JHEP12(2019)006
However, if the new fields transform differently (have distinct charges) under SU(2) L , the g Z j k L,R couplings are modified. In particular, in the presence of isosinglet heavy states, one finds (g Z j k L ) eff = g cos θ W − 1 2 + sin 2 θ W δ jk + ∆g jk L , with ∆g jk L = 6 a=4 1 2 g cos θ W (U L ) † ja (U L ) ak .
(C.5) Likewise, vector-like doublets also lead to the modification of the g Z j k R couplings: (C.6)

C.2 Couplings of the W boson
The possible mixings with the heavy vector-like leptons can also modify the couplings to the W boson. The charged current interaction terms can be written so that the corresponding charged current couplings are then given by where U P denotes the (generalised) PMNS mixing matrix. A priori, the branching ratios of W → ν can constrain the mixings of the heavy leptons (see, e.g. [176,177]). However, these strongly depend on the neutrino mass generation mechanism (as well as on the structure of the Higgs sector), and in the present analysis we will not take them into account; we nevertheless mention that for a given Higgs sector, in the presence of additional isosinglet heavy neutrinos, the modified charged current vertex can impact several observables. Therefore, in addition to electroweak precision measurements of BR(W → ν) [178], other decays or collider processes with one or two neutrinos in the final state, as for example τ decays, leptonic and semileptonic meson decays [179][180][181], production and decay of W bosons to dilepton and two jets at the LHC [182,183], can also lead to interesting constraints depending on the mass scales of the additional isosinglet neutral states.

C.3 Constraining EWP observables
Due to the tree-level modified Z-couplings (a consequence of the mixing of SM fermions with the heavy vector-like fermions), strong constraints from EWP observables are expected to arise from the observed lepton universality in charged leptonic Z-decays.

JHEP12(2019)006
At tree level, the decay width of a massive vector boson to fermions is given by [178] Γ where g L(R) are the chiral couplings and the Källén function is defined in eq. (A.6). In the case of the Z-boson, the relevant couplings have been introduced in eq. (C.3). From eq. (C.9), and in view of the very good agreement between the SM predictions and experiment (cf. table 2), it is clear that any modification of the tree-level couplings of the Z-boson will be subject to very stringent constraints, which in turn translate into bounds on the mixing parameters responsible for the ∆g L(R) terms (see eqs. (C.5), (C.6)). Using eq. (C.9), one can derive conservative constraints on ∆g L(R) from the requirement of compatibility with the bounds of table 2. As an example, the current experimental data on BR(Z → e ± µ ∓ ) and Γ( + − ) leads to |∆g eµ L | 2 + |∆g eµ R | 2 1.55 × 10 −6 , |∆g L | 5.6 × 10 −4 , |∆g R | 3.5 × 10 −4 . (C.10) The above constraints allow in turn to infer bounds on the elements of the matrix A (see eqs.

D Details of the numerical analysis
In what follows, we detail relevant aspects of the numerical analysis done in sections 2 and 3, in particular concerning the global fits and the scans over the parameter space leading to the different plots. The global fits displayed in figures 2 and 3 are obtained using the "fastfit" method of the flavio package [122]. This is based on the approximation which assumes the likelihood to be of the form L = e −χ 2 ( ξ)/2 where

JHEP12(2019)006
and C( ξ) = C exp + C th ( ξ) being the combined (theoretical and experimental) covariance matrix of the observables x i and the theoretical and experimental uncertainties are approximated as Gaussian. For a more detailed description of the statistical treatment, we refer the reader to ref. [184]. Likelihood contours are obtained by calculating the ∆χ 2 deviation around the best-fit point.
To obtain the various scatter plots, we have varied all mixing parameters of the parametrisation as stated in the captions. The mass of the vector leptoquark is varied in the ranges described in the captions of the figures (or set to a fixed value). The mixing angles were varied between −π and π. To overcome the difficulties associated with presenting a parameter space spanned by 12 mixing angles in the nonunitary parametrisation, we perform a random scan using ∼ 10 10 n-tuples of random numbers and compute the predictions of the model for a large number of flavour violating observables, including those shown in tables 1 and 3. If any of the predictions for a certain parameter n-tuple exceeds the experimental bounds in a given set of constraints, the point is filtered in the appropriate category. Points complying with the imposed constraints are further processed to calculate their global likelihood regarding the b → s and b → c ν observables using flavio [122]. Relevant effects due to renormalisation group running, in all the relevant processes (which we have used as constraints or for fitting procedures in our analysis) were computed with the wilson package [121] in association with the flavio package [122].
Open Access. This article is distributed under the terms of the Creative Commons Attribution License (CC-BY 4.0), which permits any use, distribution and reproduction in any medium, provided the original author(s) and source are credited.