The role of the $S_3$ GUT leptoquark in flavor universality and collider searches

We investigate the ability of the $S_3$ scalar leptoquark to address the recent hints of lepton universality violation in $B$ meson decays. The $S_3$ leptoquark with quantum numbers $(\overline{\mathbf{3}},\mathbf{3},1/3)$ naturally emerges in the context of an $SU(5)$ GUT model without any conflict with the stringent limits from observed nucleon stability. Scalar leptoquark $S_3$ with left-handed couplings to 2nd and 3rd generations of charged leptons and down-type quarks seems well-suited to address both $R_{K^{(*)}}$ and $R_{D^{(*)}}$. We quantify this suitability with numerical fits to a plethora of relevant flavor observables. The proposed $SU(5)$ model calls for a second leptoquark state, i.e., $\tilde{R}_2$ with quantum numbers $(\mathbf{3},\mathbf{2},1/6)$, if one is to generate gauge coupling unification and neutrino mass. We accordingly include it in our study to investigate $\tilde{R}_2$'s ability to offset adverse effects of $S_3$ and thus improve a quality of numerical fits. A global fit of the leptoquark Yukawa couplings shows that large couplings of light $S_3$ to $\tau$ leptons are preferred. We furthermore identify $B \to K^{(*)} \bar\nu\nu$ as the most sensitive channel to probe the preferred region of parameter space. Large couplings of $S_3$ to $\tau$ leptons are finally confronted with the experimental searches for $\tau$ final states at the Large Hadron Collider. These searches comprise a study of decay products of the leptoquark pair production, as well as, and more importantly, an analysis of the high-mass $\tau\tau$ final states.


Introduction
At low energies there are a few experimentally measured observables that exhibit deviation from the Standard Model (SM) predictions. Among them, the three B meson anomalies, indicating possible lepton flavor universality (LFU) violation, particularly stand out. One of these anomalies manifests itself in the ratios according to the experimental results in Refs. [1][2][3][4][5][6][7]. The result for R D appears to be 1.9 σ larger than the SM prediction, i.e., R SM D = 0.286 ± 0.012, that is obtained by relying on the lattice QCD results for both the vector and the scalar form factors [8] (see also Ref. [9]). The experimentally established R D * = 0.304 ± 0.020 has also been confirmed [6,7,10], and it appears to be ∼ 3 σ larger than predicted value R SM D * = 0.252 ± 0.003 [11]. The deviation from the SM prediction in the R D -R D * plane is at 3.9 σ level [10,12,13] and it has accordingly attracted a lot of attention recently [11,[14][15][16][17][18][19][20].
The remaining two B meson anomalies are related to the b → s + − transition. Namely, the LHCb experiment has found that there are slight discrepancies between the SM prediction and the experimental results for the angular observable known as P 5 in B → K * µ + µ − process. In many approaches this disagreement has been attributed to new physics (NP), although the tension might be a result of the SM QCD effects (see e.g. Ref. [21] and references therein). The second of the two b → s + − transition anomalies has been found in the ratio of the branching fractions, R K ≡ B(B + → K + µ + µ − ) q 2 ∈[1,6] GeV 2 B(B + → K + e + e − ) q 2 ∈[1,6] GeV 2 = 0.745 ± 0.090 0.074 ±0.036 [22], (1. 2) The values the LHCb experiment measured for these ratios are consistently lower than the SM prediction, i.e., R SM K = 1.00 ± 0.03, in which the next-to-next-to-leading QCD corrections and soft QED effects have been included [24,25] (for R K * see Table 1 and references in [23]). In other words, the LHCb results point towards a significant effect of the lepton flavor universality violation in this process. Recently, Belle Collaboration [26] found out that the angular observable P 5 agrees with the SM prediction much better for electrons than for muons. This important result suggests that it is much more likely that beyond the SM effects are present in the second generation of leptons, and that there are currently no effects in b → se + e − which would not be accounted for in the SM.
Many scenarios of NP [8,16,17, have been investigated in order to explain either R K ( * ) and P 5 , or R D ( * ) anomalies. An interesting observation was found in Ref. [17] that R K ( * ) and P 5 can be explained if NP couples only to the third generation of quarks and leptons. Furthermore, the authors of Refs. [32,50] noticed that both R D ( * ) and R K ( * ) puzzles can be correlated if the effective four-fermion semileptonic operators consist of lefthanded doublets.
In this work we consider a Grand Unified Theory (GUT) inspired setting with a light scalar S 3 leptoquark (LQ) that transforms as (3, 3, 1/3) under the SM gauge group SU (3) × SU (2) × U (1). The state S 3 is rendered baryon number conserving due to the GUT symmetry, as we discuss in Sec. 7, and generates purely left-handed currentLLQQ operators which seem to be well-suited to explain the LFU puzzles if the S 3 mass is at the TeV scale. The need for the second light LQ state,R 2 in representation (3, 2, 1/6), emerges naturally from the requirement of neutrino masses generation in the advocated GUT model as well as from the gauge coupling unification. We accordingly includeR 2 in our study and investigate whether it could partially compensate for the adverse low-energy effects of S 3 . In Sec. 2 we introduce relevant couplings of these two LQs with the SM fermions. We proceed to show how S 3 could, in principle, address the LFU puzzles in Sec. 3. We then present relevant additional constraints on the LQ parameters in Sec. 4. The low-energy flavor analysis is concluded in Sec. 5 with the global fit of the relevant couplings of the two LQs with quark-lepton pairs for three specific Yukawa structures. Sec. 6 is devoted to collider study of the model signatures in the LQ resonant pair production and in a t-channel LQ exchange contributing to τ τ final states at LHC. We elaborate on the GUT construction behind the two LQ states in Sec. 7. Finally, we conclude in Sec. 8.

Model setup
The LQ multiplet S 3 (3, 3, 1/3) interacts with the SM fermions in accordance with its quantum numbers, given in the brackets. The three charge eigenstate components of S 3 , i.e., where V is the Cabibbo-Kobayashi-Maskawa (CKM) mixing matrix. Note that S 3 has purely left-handed couplings. The diquark interactions with S 3 are not shown in Eq. (2.1) since we assume that S 3 and its interactions originate from the GUT construction presented in Ref. [52] where the baryon number violating diquark couplings are forbidden due to the grand unified symmetry. 1 The main goal of our study is to address the puzzles observed in neutral current LFU tests in the R K ratio (and related anomalies in b → sµ + µ − ) as well as in charged-current LFU ratios R D ( * ) . Thus we have clear target observables that we can affect with a small number of LQ Yukawa couplings.
In the context of SM complemented with effective operators (SM-EFT) it has been shown that NP models contributing to dimension-6 operators made out of left-handed quark and lepton doublets can explain both neutral-and charged-current LFU anomalies [17,32,50,53,54]. However, in an explicit NP model these effective interactions could be correlated, unlike in the effective theory 2 , with other observables. Our intention is to quantify this correlation within this particular model.
The LQ state S 3 can affect all the target LFU observables with a minimal set of parameters, e.g., y sµ , y bµ , and y bτ . In this work, however, we also study the effect of the coupling y sτ which enables a handle on the semitauonic modes entering R D ( * ) . The couplings of S 3 to d L and e L have to be rather small in order to avoid pressing bounds from LFV and kaon physics. We opt to set those couplings to zero to obtain the following flavor structure: Note that the Yukawa couplings of S 3 to up-type quarks are spread over generations due to CKM rotation. In what follows all Yukawa couplings are assumed to be real. The ansatz of Eq. (2.2) summarizes the most general S 3 scenario studied within our work, although we will also comment on more restricted scenarios, where some additional elements of y will be set to zero.
Having only one LQ with mass around the 1 TeV scale would invalidate unification of gauge couplings, thus a second LQ state -R 2 in our case -is needed. The two electric charge eigenstates ofR 2 couple only to down-type quarks: The doubletR 2 can accomodate the measured value of R K , but its right-handed current contributions cause tensions with the reported value for R K * . In the current setting with strictly left-handed neutrinosR 2 does not interact with up-type quarks and thus cannot affect R D ( * ) . In our approach it is S 3 that could, in principle, address both LFU anomalies, whereas its side-effects in other well-constrained observables (e.g. B s -B s mixing and B → K ( * )ν ν) might be, hopefully, cancelled byR 2 . Since S 3 will have largest effects in the τ sector we have to introduce couplings ofR 2 to τ in order to compensate for potentially unwanted effects. In the following analysis we will analyze a light S 3 scenario with the couplings texture (2.2) and along with it test the viability of having lightR 2 with nonzero Yukawas involving the τ lepton. Namely, we takẽ The mass ofR 2 should be at around 1 TeV in order to affect low-energy phenomenology, if required at all. We consistently take this to be the case when we discuss the role ofR 2 in gauge coupling unification and the neutrino mass generation.
For both LQ states the rotations with the CKM matrix V , left over from the transition to the mass basis of fermions, have been assigned to the u L fields. For the study of flavor phenomenology the neutrinos can be safely considered as massless. Thus, Lagrangians in Eqs. (2.1) and (2.3) are written in the fermion mass basis with the exception of ν L whose mass basis is ill-defined. We use flavor basis for the neutrinos, such that the Pontecorvo-Maki-Nakagawa-Sakata (PMNS) matrix becomes unity.

LFU violating contributions
In this section we focus on how the two light LQs would affect the LFU violating anomalies measured in B meson decays. The gross features required of Yukawa matrices will be presented. The detailed discussion of additional observables and their interplay with the LFU anomalies will be presented in the next section.

Charged currents LFU: R D ( * )
The largest LFU violating effect is in the charged current observables R D ( * ) . For a NPinduced effective operator that follows the chirality structure of the SM it has been shown that the dimensionless coupling of ∼ 0.1 is needed, if new particles have mass of Λ = 1 TeV and contribute at tree level [41]. The matched contributions of S 3 generate left-handed current operator, whereasR 2 cannot contribute to charged currents in this setup 3 . In particular in b → c ν transition the S 3 presence leads to the modification of the left-handed charged-current operators: The effect of S 3 may be also understood as nonuniversal CKM elements in semileptonic charged-current processes: (3.3) One also has lepton flavor violating S 3 contributions parameterized by g U D; ν , with their effect being much smaller since they do not interfere with the SM amplitude. They contribute at subleading order, namely at v 4 /m 4 S 3 that we neglect in comparison to the interference terms. Here v = 246 GeV is the electroweak vacuum expectation value. Notice that the form of interaction imposed in Eq. (2.2) implies that both decay modes B → D ( * ) τ ν τ and B → D ( * ) µν µ are affected. From the fit to the measured ratio R D ( * ) , performed in Ref. [41], we learn that at 1 σ we have the following constraint on the S 3 Yukawas: (3.4) The R D ( * ) constraint of Eq. (3.4) includes effects from τν τ and µν µ states. It is important to notice definite signs of the LQ-SM interference contributions which are proportional to V cb . Large y bτ is clearly disfavoured by (3.4) while y bµ results in negative interference term in semi-muonic modes that would be welcome from the R D ( * ) point of view, however this possibility could be in conflict with precise measurements of LFU in R µ/e D ( * ) (studied below in Sec. 4). Out of the remaining two terms y bµ y * sµ is negligible in Eq.(3.4) as required by the b → sµ + µ − . The only numerical scenario with positive interference term for the semi-tauonic mode is the one with large Cabibbo favored contribution, In the next section we will introduce constraints that put important bound on the above product of Yukawas.
3.2 Neutral currents: R K ( * ) , B(B → K ( * ) µ + µ − ) and related observables The R K anomaly can be accounted for by the additional contribution of S 3 state to the effective four-Fermi operators that are a product of left-handed quark and lepton currents [51]. TheR 2 state alone can also explain R K via the right-handed current operators [56], but the recent measurement of R K * being significantly smaller than 1 [23] implies that these operators' contributions must be small [29,56]. If we expand our analysis to a whole family of observables driven by b → sµ + µ − process the scenario with left-handed currents (S 3 state) presents a good fit and prefers the following range at 1 σ [49] (see also [53,57]): The exchange of S 4/3 3 contributes towards the above effective coefficients as For a range (3.6) of Wilson coefficients we find Contrary to S 3 , the right-handed quark currents generated byR 2 do not improve significantly the global agreement between theory predictions and observables related to the b → sµ + µ − . Tree-level matching ofR 2 amplitudes yields Using the result of the global fit from [49] we have checked that including non-zeroỹ sµ and y bµ does not improve the fit considerably.

Constraints on the LQ couplings
The introduction of the two LQ states with sizable couplings to explain the LFU observables, as presented above, inevitably causes side effects in related flavor observables which we will focus on in this section.
where we have averaged over the two Belle results. Due to its smallness the term y sµ y * bµ is irrelevant in the above equation (see Eq. (3.8)), albeit the factor ∼ 20 enhancement due to CKM. After this simplification Eq. (4.1) becomes a rather weak limit, i.e., |y bµ | 1.5(m S 3 /TeV). It is, however, clear that y bµ , in spite of its large value, is not sufficient to explain the R D ( * ) constraint of Eq. (3.4). Since the largest effects are concentrated in the τ flavor, we expect large effect in leptonic decay of B − → τν which is sensitive to |V , as given in Eq. (3.3). The B − → τν rate is thus enhanced by the same combination of Yukawas (and same order of Cabibbo angle) that also drives the B → D ( * ) τν rate. The current experimental average B(B − → τν) = (1.09 ± 0.24) × 10 −4 is indeed slightly higher than the SM prediction B(B − → τν) SM = (0.78 ± 0.07) × 10 −4 . If we assume that LQ Yukawas are real numbers then the leading contribution y * sτ y bτ in both observables leads to correlation where the CKM factor relating the two observables is close to unity.

Semileptonic K and τ decays
On the other hand, LFU in kaon decays has been tested and confirmed with better precision through the following ratios: As pointed out in Ref. [42] these observables enable us to put strong constraints on the corrections arising within models of NP. In the e/µ sector the experimental result [60] agrees well with the SM prediction [61]: (4.5) R K e/µ is most sensitive to |y sµ | since the product y * bµ y sµ must be small as dictated by b → sµµ sector and comes with an additional CKM suppression. The agreement of experiment [60] with the SM prediction [62] in the τ /µ exhibits a ∼ 2 σ tension: where the dominant error of the experimental ratio is due to the τ lifetime uncertainty, whereas on the theory side it is the radiative correction δR τ /K = (0.90 ± 0.22)% [63] which is the source of uncertainty. The constraint is expressed as: (4.7)

Leptonic decays: W → τν, τ → νν
The SM tree-level vertexτ νW is rescaled due to penguin-like contribution of both S 3 and R 2 . As we integrate out S 3 andR 2 at the weak scale the W vertex with τ leptons reads Free color index in the loops graphs results in the N c = 3 factor in front. We have neglected the quark masses in the above calculation and presented only the leading terms in x and x. The contribution of S 3 with mass of 1 TeV shifts the W → τ ν decay width relatively by 4 × 10 −4 (|y bτ | 2 + |y sτ | 2 ) which is well below the current ∼ 2% experimental precision. The W → µν is also rescaled by an analogous δ (µ) W factor. At low energies the effective W → τ ν vertex would, together with direct box contributions with LQs, manifest itself in the τ → ν ν τ decays. Only S 3 may participate in the box diagrams sinceR 2 has no direct couplings to . The effective interaction term of τ → ν τν then reads −g 2 As it has been pointed out recently in the literature [54,55,62] the LFU observable R τ / τ , defined as a ratio B(τ → νν)/B(µ → eνν), and normalized to the SM prediction of this ratio, is very sensitive to models modifying couplings of the τ lepton. Experimentally, R τ /µ τ = 1.0022 ± 0.0030, R τ /e τ = 1.0060 ± 0.0030, while in the present model the leading interference terms shift the ratios as (4.10)

Semileptonic decays of D and t
We have checked the effect of S 3 on the leptonic charm meson decays D s → ν. Using the bounds from kaon LFU observables presented above we find that the S 3 correction to the D → µν width is below 1%, while the experimental uncertainty of D s → τ ν is 4% and can easily accommodate |y sτ | 1.2(m S 3 /TeV) without even taking into account the uncertainty in the decay constant f Ds . For the semileptonic top decay process among the third generation fermions, t → bτ + ν, the correction is also below the current sensitivity [64].

τ → µγ
Current bound B(τ → µγ) ≤ 4.4 × 10 −8 has been determined by the BABAR collaboration [65]. The S 3 LQ contributes to the τ → µγ amplitude via b and s quarks and S 4/3 3 in the loop and also via up quarks u, c, and t mediated by the S −1/3 3 component. Using the loop functions in the small quark mass limit as in Ref. [51] we determine where the effective coupling reads At loop level, S 3 andR 2 modify the Z → f 1f2 decay widths which were precisely measured at LEP-2. The largest effects in presented LQ model are expected for third generation final states both in flavor conserving decays, as in Z → τ + τ − , which has been shown to have only weak constraining power in Ref. [66], as well as in LFV modes Z → τ ± µ ∓ . The latter decay happens due to penguin diagrams with S 3 as well as 1-particle reducible diagrams and is suppressed by a loop factor and small ratio x = m 2 Z /m 2 S 3 , in which we expand to leading order: (4.13) We have checked that B(Z → µτ ) is well below the current experimental bound at 10 −5 . Compared to the closely related τ → µγ decay, this channel is less stringently constrained and thus we do not include it in the fit. On the other hand, the B(τ → 3µ) < 2.1 × 10 −8 at 90% C.L. [60], and can be mediated by the above mentioned LFV Z vertex or via box diagram containing S 3 and quarks. They are both encompassed in the low-energy effective Lagrangian: where, again, The mixed chirality stems from the Z coupling toμ R µ R . In the limit of m µ /m τ → 0 the two terms above do not interfere. We notice that when all couplings are ∼ 1 and m S 3 = 1 TeV then the B(τ → 3µ) is in the ballpark of current experimental upper bound. As will be shown in Sec. 5, realistic values of the Yukawas result in much smaller contribution to this channel, and that is why we omit this channel from the fit.

(g − 2) µ
The difference between the experimental value and the one predicted by the SM is δa µ = a exp µ − a SM µ = (2.8 ± 0.9) × 10 −9 [60]. Following [51] and using the Lagrangian of Eq. (2.1), we derive the contribution of S 3 to the muon anomalous magnetic moment: Since the above contribution has wrong sign with respect to the experimental pull individual Yukawa couplings of S 3 to the µ should be small. Notice that the contribution ofR 2 to (g − 2) µ is greatly suppressed and vanishes at m s,b /mR 2 → 0 [51, 67].

B → Kµτ decays
The lepton flavor violation can be induced by the LQ presence at tree level in B → Kµτ and also in decays of bottomonium to τ µ. As noticed in [42] the latter process has been constrained at the level of 10 −6 however these bounds are not competitive with the bound B(B − → K − µ ± τ ∓ ) < 4.8 × 10 −5 at 90% C.L. [68]. This inclusive mode is sensitive to both y sτ y bµ and y sµ y bτ as B( when the form factors of Ref. [69] are used. The constraint then reads Despite being a loop observable in the LQ scenarios, the B s meson mixing frequency is one of the most important constraints in our particular setup where the product of S 3 Yukawas y bτ y sτ is large. This product alone would lead to uncomfortably large effect in the B s -B s oscillation frequency ∆m s . However, there is an additional box amplitude due toR 2 as well as an amplitude with both S 3 andR 2 propagating in the box, as shown Figure 1. Three types of box-diagrams with S 3 andR 2 contributing to ∆m s .
in Fig. 1. Amplitudes that correspond to the first and second diagram in Fig. 1 can be found in Refs. [8,51] and contribute to operators C 1 andC 1 of the effective Hamiltonian, respectively: The third diagram in Fig. 1 in which both LQs are present, but couple with opposite chirality to the fermions, contributes to the Wilson coefficient C 5 . There the color indices α and β are summed across ∆B = 1 currents. The box diagrams in Fig. 1 are well approximated using a limit of massless virtual leptons and match onto the effective Hamiltonian at scale Evaluation of hadronic matrix elements for B s -B s mixing is performed at the scale µ = m b (m b ) = 4.2 GeV. Utilizing parameterization in terms of bag parameters as in Ref. [70], we find for the oscillation frequency For the SM prediction we use the perturbative QCD renormalization at next-to-leading order whose effect is subsumed in η 2B = 0.55(1) [71]. The non-perturbative parameters and perturbative RG running effects of C 1 are combined into a scale-invariant combination, where the value of renormalization-group invariant bag parameter is taken from the QCD lattice simulation with three dynamical quarks [72]: f 2

BsB
(1) Bs = 0.0754(46)(15) GeV 2 4 . First number in the brackets represents statistical and systematic error, apart from systematic error due to omission of dynamical charm-quark, which is shown in the second bracket. The SM prediction is then ∆m SM s = (19.6 ± 1.6) ps −1 . For the LQ contributions in Eq. (4.19) we use the values of B (i) Bs (µ) from Ref. [72]. For the multiplicative renormalization of coefficients C S 3 1 andCR 2 1 we neglect the running from Λ to m t , such that running effect to low scale is the same as in the SM, whereas for CR 2 S 3 4,5 we use the leading order mixing [74] to find

B → K ( * ) νν
The B → K ( * ) νν decay offers an excellent probe of the lepton flavor conserving as well as lepton flavor violating combination of the LQ couplings. Following [42] and with the help of notation in Refs. [56,75,76], we write the effective Lagrangian: In the SM we have a contribution for each pair of neutrinos and therefore C SM,ij L = C SM L δ ij where C SM L = −6.38 ± 0.06 [75]. The respective contributions of S 3 andR 2 to the left-and right-handed operators are [51]: While the amplitude of B → Kνν depends only on the vectorial part of Wilson coefficients (4.22), the B → K * ν ν amplitude is also sensitive to axial current, and the two decay modes constrain the right-handed Wilson coefficient differently. We follow Ref. [76] and introduce ij = . (4.23) 4 We prefer to use the results of Ref. [72] that include bag parameters for the whole operator basis.
The measurements of µ + µ − spectra at high invariant mass m µµ are sensitive to large couplings y sµ or y bµ . The relevant channel in our case is bb → µ + µ − which directly limits y bµ . If we assume that effective dim-6 operator description is a good approximation to the t-channel S 3 exchange at LHC energy, then we can use a 1 σ bound derived in Ref. [78] (4.25)

D decays
The weak triplet nature of S 3 implies couplings only to the weak doublets of quarks and leptons, and thus corrections to the charged current processes only rescale the SM charged current contributions. The dominant modification of V cs element associated with semimuonic decays follows from Eq. (3.2): (for processes with µν µ ). Assuming that the CKM-suppressed y bµ term can be neglected in Eq. (4.26) and using the fact that current precision on the semileptonically determined V cs reaches 1 per-mille [60], we find y sµ 0.3(m S 3 /TeV). Rare charm decays with two leptons, e.g. D 0 → µ + µ − and D → M µ + µ − , are most constraining at the moment (for dineutrino modes cf. [79]), where M can be a pseudoscalar or a vector meson. The effective Wilson coefficient of the left-handed current, C 9 = −C 10 ≈ (V us πv 2 )/(αV ub V * cb m 2 S 3 )y 2 sµ can be compared to the bounds, |C 9 |, |C 10 | 1.0/|V ub V cb |, obtained in Ref. [80]. We learn that the ensuing bound y sµ 0.5(m S 3 /TeV) from rare decays is weaker than the abovementioned bound from semileptonic decays.

Flavor couplings
In this section we study three scenarios differing in the number of variable Yukawas. For each scenario we report a minimum of χ 2 function, which is a sum of terms corresponding to all observables discussed in the preceding sections. We also report 1 σ regions for the interesting two-dimensional projections of parameter space. While performing these fits we limit all free Yukawa couplings to be smaller than 3. Introduction of this artificial cut-off is guided by the constraints posed by the LHC searches, discussed in Sec. 6. The SM point has χ 2 = 71.6 and serves as a reference value to which χ 2 of the three fits are compared.

S 3 coupled to muons (2 parameters)
In this scenario we consider only the effect of S 3 with non-zero muonic couplings: Left panel in Fig. 2 exposes tension between R D ( * ) (2.8 σ pull) and R µ/e D ( * ) (1.8 σ pull) which is even more exacerbated when we include the direct constraints on y bµ from LHC (right panel of Fig. 2). The latter scenario with all constraints included has χ 2 = 42.4 at point (y sµ , y bµ ) = ±(2 × 10 −3 , 0.46) which corresponds to the 5.0 σ pull of the SM hypothesis. One can observe in the right panel in Fig. 2 that in this case the preferred region is drawn further away from R D ( * ) . The results indicate that R D ( * ) cannot be explained by omitting couplings to τ . Detailed results on the pulls are given in the third column of Tab Table 1. Observables that enter the global fit with their pulls in σ in the SM and S 3 scenarios. Third column represents the case when m S3 = 1 TeV and only y sµ , y bµ are allowed, without/with taking into account bb → µµ constraint. Fourth column represents the fit of the y sµ , y bµ , y sτ , y bτ scenario for m S3 = 1.0/1.5 TeV. The constraints with negligible pulls are not shown in this table.

S 3 coupled to muons and taus (4 parameters)
Since the purely muonic couplings are in conflict with R D ( * ) we allow in addition for tauonic couplings of S 3 : In this case both couplings with the muons tend to be small, below 0.1, and are relevant only in b → sµµ, whereas the couplings to τ are ∼ 1 in order to enhance R D ( * ) . For m S 3 = 1 TeV we find that the minimal χ 2 of this scenario with 4 degrees of freedom is 36.8 reached at (y sµ , y bµ , y sτ , y bτ ) = (0.047, 0.020, 0.87, −0.048) 5 which makes the SM point excluded at 5.0 σ (pull). In Fig. 3 the fit in the tauonic couplings' plane shows how the optimal region is still far from the central value of R D ( * ) , mostly due to R νν and ∆m s , which do not allow for large products of y bτ y sτ . Pulls of individual observables for m S 3 = 1.0/1.5 TeV are presented in the fourth column of Table 1.

S 3 andR 2 (6 parameters)
In order to relax the tension in the y sτ -y bτ plane between large effect in R D ( * ) and well constrained R ( * ) νν and ∆m s , we could invoke a lightR 2 with couplings to τ . We consider a case m S 3 = mR 2 = 1 TeV with six free Yukawa couplings (y ij from the previous subsection and (ỹ sτ ,ỹ bτ ) pair) to find χ 2 = 33.4 at (y sµ , y bµ , y sτ , y bτ ) = (0.051, 0.019, 0.86, −0.069), (ỹ sτ ,ỹ bτ ) = (3, 0.0026) 6 that represents a 4.9 σ pull of the SM. Most importantly, the tension in R D ( * ) is only marginally improved and stands at 3.7 σ. The presence ofR 2 allows for partial cancellation in ∆m s between large tauonic couplings of S 3 andR 2 , which is not the case in both R νν and R * νν , where cancellation in one observable necessary spoils the other (cf. (4.24)). We thus conclude that lightR 2 with relatively large couplings to the SM fermions cannot improve substantially the agreement with data. We accordingly assume that the couplings ofR 2 , i.e.,ỹ ij of Eq. (2.3), are small enough as not to affect flavor observables. With lightR 2 and its Yukawa couplings sufficiently small to avoid flavor constraints we can still aid gauge coupling unification and generate viable neutrino masses in the underlying GUT model as we show in Sec. 7.

Collider constrains
In what follows we confront our model, comprising two light LQs, with collider constraints while taking into account the particularities of the flavor structure derived in the previous section. We demonstrate the viability of the proposed model and present bounds from third generation LQ pair production as well as high-mass τ τ production searches at the 13 TeV LHC for current and projected luminosities. We show that a large portion of the relevant parameter space can be covered by the HL-LHC.

LQ pair production
The current best mass limit for LQs that decay to the third-generation leptons has been recently reported by the CMS Collaboration while searching for a pair of QCD produced LQs decaying into the τ + τ − bb channel [81]. This search excludes an LQ with mass bellow 850 GeV (550 GeV) for a branching ratio (BR) of β = 1 (β = 0.5). This search can set limits on the parameter space of our model via pp → S We focus on the scenario presented in Sec. 5.2 when S 3 is at the TeV scale andR 2 is assumed not to feed significantly in the τ + τ − bb signal. In this case the CMS bound can be applied directly to the S 4/3 3 state, at a benchmark mass of m S 3 = 1 TeV, decaying into a τ b pair with a BR given by β ≈ |y bτ | 2 |y bτ | 2 + |y sτ | 2 . (6.1) Here we neglect the small widths of S 4/3 3 into both muonic channels. Results are given in Fig. 4 (left panel), where the dashed blue contours represent the 95% C.L. exclusion limits for different LHC luminosities and the red region represents the 1 σ region for the mS 3 = 1 TeV   → τ + τ − bs. 7 The search starts losing sensitivity for m S 3 above 1 TeV, while for masses larger than 1.2 TeV the search does not produce any useful limits. In conclusion, a third generation LQ pair production search at the LHC is not a sensitive probe for this particular flavor structure of our LQ model.

LHC 13 TeV
We now turn to the scenario where two generic third generation LQs (like e.g. S 3 and R 2 ) are at the TeV scale and both contribute to LQ pair production. A naive reinterpretation of the CMS limits [81] can be performed when (i) both LQ components are degenerate in mass 8 , (ii) interference terms between the final state τ b pairs at the amplitude level are negligible and, (iii) the LQ decay widths are small enough in order to guarantee the narrow width approximation (NWA) assumed in the experimental search. This scenario would thus correspond to the particular case we investigated in Sec. 5.3.
As shown in Appendix A, the CMS bound for one LQ can be directly mapped into a bound for two degenerate LQs we denote with LQ 1 and LQ 2 , for simplicity, when the associated BRs are β 1 and β 2 , respectively. (See Fig. 2 (right panel) of Ref. [81] for the experimental limit.) The inferred limits, which apply in general to two third-generation LQs with non-interfering final states, are presented in Fig. 4 (right panel) for an integrated luminosity of 12.9 fb −1 (solid green contour). For example, we find that both LQs with equal masses bellow 930 GeV (600 GeV) are excluded at 95% C.L. if β eff = 1 (β eff = 0.5), where β eff = (β 2 1 + β 2 2 )/2. We also include LHC projections for an integrated luminosity of 300 fb −1 (dashed green contour) by rescaling the CMS Collaboration limits with the square root of the luminosity ratio. At this projected luminosity the LQ pair production search is not sensitive for masses above 1.3 TeV. These results, again, can be used for our LQ model for degenerate S 3 andR 2 LQs at the TeV scale. Unfortunately, the LHC bounds on the couplings extracted from LQ pair production search are not strong enough to probe this scenario either.

High-mass τ τ production
In this section we study the implication of light S 3 andR 2 leptoquarks for high-p T τ τ production at the LHC. It was shown in Refs. [32,82] that τ τ resonance searches at the LHC produce stringent constrains on a large class of models explaining the R D ( * ) anomaly. In what follows we give predictions for the deviation from the SM in the invariant mass tails of pp → τ + τ − and derive bounds at different luminosities for the parameter space of the present LQ model from the 13 TeV ATLAS Collaboration resonance search at 3.2 fb −1 [83].
Both LQs contribute to pp → τ + τ − production exclusively through Yukawa interactions by exchanging S This process will produce events in the signal region defined in Ref. [81], which is based on only one b-tagged jet and not two. We have also excluded from our analysis contributions coming from the non-QCD LQ pair production. 8 For the non-degenerate case the results of Ref. [81] are not directly applicable given that each LQ will have different kinematic distributions leading to different selection efficiencies in the signal region. contributions may come from the processes with incoming strange quarks ss → τ + τ − and sb (bs) → τ + τ − , followed by sub-leading contributions from bottom, charm and up quark initiated processes bb (cc) (uū) → τ + τ − . The flavor structure in Eq. (2.2) also allows for S 1/3 3 to couple to u and τ via the CKM mixing. Nevertheless, this coupling is proportional to |V us | y sτ , meaning that τ τ production from incoming up quarks is Cabibbo suppressed leading to negligible cross-sections of order |V us | 2 and |V us | 4 for the processes cū → τ + τ − and uū → τ + τ − , respectively. The Cabibbo suppressed vertices are shown in red in Feynman diagrams of Fig. 5. On the other hand, at high-x the large proton PDF of the valence up quark in the process uc → τ + τ − can marginally compensate for the |V us | suppression in the amplitude giving a contribution comparable to cc → τ + τ − in the total cross-section.
We now focus on the total cross-section σ fid TOT of pp → τ + τ − far from the Z-pole in the high-mass tails of the τ τ invariant mass distribution. We will, for definiteness, study the scenario where only S 3 contributes to τ τ production. The couplings ofR 2 are assumed to be small and can thus be safely neglected for this collider study. This is in accordance with the outcome of the numerical study presented in Sec. 5.2.
At leading-order (LO), τ τ production will receive contributions from the t-channel exchange of S 3 , from the s-channel SM Drell-Yan pp → Z/γ * → τ τ production, and from interference effects between these processes. The high-mass kinematic region is defined by the following fiducial cuts on the final states: p T > 100 GeV (50 GeV) for the leading (subleading) τ and a high invariant mass cut for the τ τ pair of m τ τ > 600 GeV. We define the signal strength µ pp→τ τ as the ratio of σ fid TOT with the the SM Drell-Yan fiducial cross-section σ fid SM : Here the fiducial cross-section σ fid LQ includes all NP contributions from both the LQ squared and LQ-SM interference amplitudes, i.e., σ fid LQ = 2 Re(A * SM A LQ ) + |A LQ | 2 . The ratio σ fid LQ / σ fid SM quantifies the NP deviation of the total fiducial cross-section from the expected SM prediction. The LQ Yukawa couplings enter in σ fid LQ as σ fid LQ (y sτ , y bτ ) = σ ss (y sτ ) + σ sb (y sτ , y bτ ) + σ bb (y bτ ) + σ cc,uū,uc (y sτ ), (6.3) In order to conform with the analysis in Sec. 5 we assume all Yukawa couplings to be real.

4)
σ sb (y sτ , y bτ ) = y 2 sτ y 2 bτ A 2 , (6.5) The polynomial coefficients A i and B i are functions of the mass m S 3 describing the LQ squared amplitudes and LQ-SM interference amplitudes, respectively. In Eqs. (6.4)-(6.7) we define all polynomial coefficients to be positive and include the explicit signs of the LQ-SM interference coefficients B i , indicating the presence of either destructive or constructive interference amplitudes. The origin of the sign of these interference terms can be traced back to the SM amplitude proportional to the T 3 − sin 2 θ W Q coupling of the Z boson with the incoming quarks along with the sign of the LQ Yukawa interactions of S The interference arising between the Z boson and the S

4/3 (1/3) 3
components is driven by the weak isospin of the quark doublets, i.e., positive (negative) for up(down)-type quarks. This translates into a constructive interference for down-type quarks in dd, ss , bb → τ τ processes and destructive interference for up-type quarks in cc , uū → τ τ processes. As a mττ > 600 GeV pT (τ1) > 100 GeV pT (τ2) > 50 GeV  consequence, this relative sign that appears between the LQ-SM interferences in different τ τ production channels leads to a partial cancellation in the cross-section σ fid LQ . For more details see Appendix B.
For the calculations below we choose a benchmark mass value of 1 TeV. In order to extract the values of the polynomial coefficients A i and B i we generate in FeynRules [84] the UFO model file for the Lagrangian for S 3 and simulate in MadGraph5 [85] 13 TeV event samples of pp → τ + τ − subject to the high-mass fiducial cuts discussed above for different values of the couplings. The specific values of the coefficients can be found in Appendix B. In Fig. 6 (left panel) we show results for different contours of constant signal strength µ pp→τ τ in the y sτ -y bτ plane for the 1 TeV mass scenario. This shows that our LQ model equipped with the proposed flavor structure with y sτ,bτ couplings of O(1), predicts at the LHC an enhancement of O(10%) in the τ τ tails when compared to the SM.
In the remaining part of this section we confront our model with existing LHC data and extract 95% C.L. limits for the model parameters at different integrated luminosities. For this we recast a heavy Z search by the ATLAS Collaboration with 3.2 fb −1 of data in the fully hadronic τ had τ had channel [83]. For the recast we generate at LO in MadGraph5 a set of the LQ τ had τ had signal samples followed by parton showering and hadronization in Pythia8 [86]. We have included interference effects between the LQ signals and the SM Drell-Yan background process. Detector effects are simulated in Delphes3 [87] with settings tuned according to the experimental environment of the τ had τ had inclusive category as described in Ref. [83]. For this category, events are selected if the reconstructed objects satisfies the following requirements: • p T > 110 GeV (55 GeV) for the leading (sub-leading) τ had .
• Total transverse mass cut of m tot T > 350 GeV.
Here m tot T is the dynamical variable used to reconstruct the invariant mass of the visible part of τ had τ had defined as (m tot is the total missing energy in the event and m 2 cos ∆φ(A, B)) is the squared transverse mass of objects A and B. In order to extract limits from the tails of the m tot T distributions we use the statistical analysis presented in Ref. [82]. For m LQ > 1 TeV the most sensitive bin corresponds to the last one (m tot T > 684 GeV) where a point in parameter space is excluded at 95% C.L. if the number of events in that bin exceeds 11, 36 and 113 at integrated luminosities of 30, 300, 3000 fb −1 respectively. Here we have applied a naive scaling of the limits for the 3.2 fb −1 search to arbitrary luminosities with L int /3.2 fb −1 . The results are given in Fig. 6 (Right panel) for the benchmark mass of m LQ = 1 TeV. There the regions to the right of the dotted blue boundaries are excluded at 95% C.L. for future LHC luminosities of 30 fb −1 , 300 fb −1 and 3000 fb −1 . In the plot we have also included the 1 σ region in solid red obtained from the low-energy fit to all flavor experiments performed in Sec. 5.2. Notice that these bounds are conservative given that we have not included NLO QCD corrections and have not taken into account other processes such as non-resonant qg → τ τ q that can produce sizeable contributions to the inclusive τ had τ had category 9 . From these results we conclude that the High-Luminosity LHC can probe a large portion of the parameter space for the Yukawa coupling ansatz of Sec. 5.2.

GUT completion
The preceding sections were devoted to the study of the impact of light scalar fields S 3 and R 2 with potentially sizeable couplings to the quark-lepton pairs on the flavor physics processes and LHC observables. Here we want to demonstrate that these fields and associated couplings can originate from a consistent grand unified theory (GUT) model.
To insure that the LHC accessible leptoquarks S 3 andR 2 are not in conflict with stringent limits on matter stability it is necessary that both S 3 andR 2 do not couple to the quark-quark pairs either directly or through the mixing with other scalars in a specific model of unification. It turns out that one can meet this requirement in an SU (5) model that comprises 5-, 15-, 24-, and 45-dimensional scalar representations [52]. The decomposition of the scalar sector of that model is and ∆ 3 are identified withR 2 and S * 3 , respectively. The fermions of the SM, on the other hand, are embedded within the tenplets and fiveplets in the usual manner [88].
We first show that this GUT scenario is compatible with the viable gauge coupling unification. To this end, we take all scalar fields in the model that mediate proton decay at tree-level to reside at or above 10 12 GeV. These fields are Ψ T , ∆ 5 , and ∆ 6 [52]. We furthermore set the masses of both S 3 andR 2 at 1 TeV and constrain all remaining scalar fields to be at or above one scale we simply denote m that is to be determined through the requirement that the gauge coupling unification takes place at the one-loop level. Note that Σ 24 does not affect unification. Also, Σ (3,2) and Σ (3,2) are not physical fields since they provide necessary degrees of freedom for the baryon and lepton number violating gauge bosons X and Y of the SU (5) origin to become massive fields.
The gauge couplings meet at the unification scale m GUT when the following equation is satisfied [89] B 23 B 12 = 5 8 where the right-hand side is evaluated using α S (m Z ) = 0.1193 ± 0.0016, α −1 (m Z ) = 127.906 ± 0.019, and sin 2 θ W = 0.23126 ± 0.00005 [90]. The left-hand side depends on the particle content and the mass spectrum of the model. Namely, coefficients B ij are where b J i are the well-known β-function coefficients of particle J with mass m J and r J = (ln m GUT /m J )/(ln m GUT /m Z ). The sum goes through all particles beside the SM ones that reside between Z boson mass m Z and m GUT . The convention is such that b J 1 , b J 2 , and b J 3 are associated with U (1), SU (2), and SU (3) of the SM, respectively. We identify m GUT not only with the gauge coupling unification scale but with the masses of the proton decay mediating gauge boson fields X and Y .
If and when unification takes place for a given m we evaluate m GUT using equation [89] ln m GUT m Z = 16π 5α to check that m GUT ≥ 5 × 10 15 GeV in order to satisfy stringent bounds on the X and Y gauge boson mediated proton decay. To actually set a lower bound on m we fix m GUT = 5 × 10 15 GeV in our analysis and maximise m. We find that m = 3.1 × 10 10 GeV when the masses of bothR 2 and S 3 are at 1 TeV. The masses of all other scalar particles in the model GeV, m ∆ 5 = 10 12 GeV, and m ∆ 6 = 10 12 GeV, m ∆ 7 = m GUT . Note that the SM Higgs is in principle a mixture of Ψ T and ∆ 7 . We accordingly take one state to be light and treat the mass of the other as a free parameter that is between m and m GUT . The fact that viable unification can take place when S 3 andR 2 are both light does not come as a surprise. Note that the SM field content yields B SM 23 /B SM 12 = 0.53 instead of the experimentally required value given in Eq. (7.1). The nice feature of the set-up with light S 3 andR 2 is that both fields have positive b J 23 and negative b J 12 coefficients. This not only helps in bringing the left-hand side of Eq. (7.1) in agreement with the required experimental value but simultaneously raises the GUT scale m GUT through Eq. (7.2). The relevant coefficients are b S 3 23 = 9/6, b S 3 12 = −27/15, bR 2 23 = 1/6, and bR 2 12 = −7/15. Again, our findings demonstrate that the gauge coupling unification is possible for light S 3 andR 2 in this particular model. We next demonstrate that the explicit forms of the Yukawa couplings of S 3 andR 2 that are used in Sec. 5 to produce numerical fits can originate from the appropriate SU (5) operators. It is also argued that the model can accommodate realistic masses of the SM fermions.
The S 3 ∈ 45 lepton-quark couplings originate from the SU (5) contraction y 45 ij 10 i 5 j 45, where 10 i are the usual fermionic tenplets, i(= 1, 2, 3) is the generation index, and y 45 is a 3 × 3 matrix in flavor space. We can thus identify y of Eq. (2.1) with y 45 / √ 2, where y 45 is related to the difference in masses between charged fermions and down-type quarks [52]. This follows from the fact that there are actually two operators that contribute towards the charged fermion and down-type quark masses in this SU (5) model [91]. One is y 45 ij 10 i 5 j 45 and the other is y 5 ij 10 i 5 j 5, where y 5 is an arbitrary complex matrix. Clearly, y 45 = √ 2y and y 5 together contain enough parameters to easily address observed mismatch between the charged fermion and down-type quark masses. For completeness we specify that the uptype quark masses originate from a single contraction x ij 10 i 10 j 5, where x ij is a symmetric complex 3 × 3 matrix. TheR 2 ∈ 15 lepton-quark couplings are symmetric in flavor space since they originate from y 15 ij 5 i 5 j 15, where 5 i are the usual fermionic fiveplets. We identifyỹ ij of Eq. (2.3) with −(D R y 15 ) ij / √ 2 in the physical basis for the down-type quarks and charged leptons, where D R represents unitary transformation of the right-chiral down-type quarks. If we take that y 15 33 = 0 we obtain the form ofỹ that is used in the fit of Sec. 5 when we consider joint effect of S 3 andR 2 on flavor observables.
It is worth mentioning that it is possible to address neutrino masses within this model. Namely, if one turns on a vacuum expectation value of the electrically neutral field Φ a ∈ 15 one can generate neutrino masses of Majorana nature via type II see-saw mechanism [92,93] through the same operator that yields theR 2 lepton-quark couplings, i.e., y 15 ij 5 i 5 j 15. In this particular instance the entries in y 15 would need to be responsible for the observed masssquared differences and mixing angles in the neutrino sector. That requirement would not be compatible with a simple ansatz for the structure ofỹ given in Eq. (2.4). Again, viable neutrino masses would only be possible if we depart from that ansatz and assume that thẽ y ij entries are sufficiently small to avoid flavor constraints for lightR 2 . This is in agreement with the findings we presented in Sec. 5.3. Note that the neutrino Majorana masses could also receive partial contribution through the one-loop processes, where the particles in the loop are down-quarks and a mixture of S 3 andR 2 [52,94,95].
The preceding discussion demonstrates that the SU (5) GUT model comprising 5-, 15-, 24-, and 45-dimensional scalar representations, with the canonical embedding of the SM fermions, can accommodate light S 3 andR 2 and describe observed fermion masses of the SM without any conflict with relevant experimental constraints.

Conclusion
Our aim, in the present work, is to accommodate the observed lepton non-universality in charged current processes, signalled by R D ( * ) , as well as lepton non-universality and the global tension in the b → sµ + µ − sector through the introduction of light scalar LQ S 3 . This LQ emerges naturally in the context of a specific SU (5) GUT model and has to be accompanied by another light scalar LQR 2 which improves gauge coupling unification and aids neutrino mass generation.
The first state, S 3 , couples left-handed s and b to left-handed µ and τ and is capable of accommodating b → s + − sector. Because of its weak triplet nature it also couples to up-type quarks and neutrinos which are precisely the additional couplings needed to address R D ( * ) . Large couplings needed for R D ( * ) cause the weak triplet S 3 to inevitably contribute to other well constrained flavor observables that agree with the SM predictions. We have analyzed those in detail and demonstrated that the most pressing ones, R ( * ) νν = B(B → K ( * )ν ν)/B(B → K ( * )ν ν) SM and ∆m s , allow only for minor improvement of R D ( * ) puzzle. We furthermore show that the second state,R 2 , cannot significantly improve the agreement with data.
Based on the numerical values of the LQ Yukawa couplings as obtained in the flavor fit we recast two LQ collider searches: (i) search for pair produced LQs decaying to bτ bτ , (ii) search for high-mass τ τ final state which is sensitive to the t-channel LQ exchange.
From the recast of the search for the LQ pair production we find that the proposed scenario with m S 3 = 1 TeV cannot be significantly probed in a large portion of the parameter space even with 300 fb −1 of integrated luminosity at the Large Hadron Collider. Complementary searches for τ τ final states produced via a single LQ exchange, on the other hand, are not as hampered by large LQ masses and are already excluding corners of the parameter space with largest individual Yukawa couplings. Moreover, since the flavor fit requires increasing Yukawa couplings for larger LQ masses the sensitivity of high-mass τ τ final state search does not degrade at higher masses as compared to the pair production mechanism. This method can probe almost entire parameter space of the model at 3000 fb −1 of integrated luminosity.
A natural ultraviolet completion for the two LQ states is an SU (5) GUT. We demonstrate that a particular setting with 5-, 15-, 24-, and 45-dimensional scalar representations is consistent with unification of the gauge couplings, where light S 3 andR 2 leptoquarks reside in 45and 15-dimensional representations, respectively. Furthermore, baryon number violation is sufficiently suppressed by lack of diquark couplings of S 3 and high enough scale of the rest of the GUT spectrum. The model also accommodates the masses of all fermions of the SM.

A LQ pair production recast
In this appendix we give a reinterpretation of the results by the CMS Collaboration [81] for the case of two LQs, denoted LQ 1 and LQ 2 . When addressing the LQ pair production from QCD interactions, a model with two LQs of degenerate mass m LQ ≡ m LQ 1 = m LQ 2 with BRs β 1 and β 2 for LQ 1 → τ b and LQ 2 → τ b, respectively, can be consistently mapped to a model with only one LQ, denoted here as LQ, with an effective mass m eff and an effective BR β eff for LQ → τ b. In this case, assuming the NWA, the total cross-section for pp → LQ * LQ → τ + τ − bb is factorized into production and decay modes as where σ pair is the pp → LQ * LQ pair production cross-section that depends exclusively on the LQ mass when only QCD interactions are taken into account. For the numerical calculations we use the approximate expression from Ref. [96] for the cross-section at NLO σ pair (m) ≈ exp where σ −1 is the inverse function of Eq. (A.2). Here we assume negligible interference effects between the decay products of the LQ 1,2 and simply add two cross-sections together. After calculating σ −1 numerically we can use Eq. (A.3) to map the CMS Collaboration 12.9 fb −1 exclusion limits in the β-m LQ plane as reported in Fig. 9 of Ref. [81] into the exclusion limits for two generic non-interfering third-generation LQs with degenerate mass. These limits are shown in Fig. 4.

B High-mass τ τ production cross-sections
We obtain the following fiducial cross-sections in fb for the process pp → τ τ for m LQ = 1 TeV: Notice that in each individual production channel the interferences can be large. In particular, these dominate in cc (uū)(uc) → τ τ production over the squared LQ terms for Yukawa couplings of order one, as shown in Eq. (B.4). Only after summing across all channels the total interference is found to be sub-leading when compared to the total LQ squared amplitudes in most portions of parameter space. This happens because of an accidental cancellation between the constructive S 3 -Z interference in ss → τ τ given by the second term in Eq. (B.1) and the destructive S 3 -Z interference in cc (uū)(uc) → τ τ given by the second term in Eq. (B.4). The remaining small (constructive) interference after cancellations is mostly given by τ τ production from bottom fusion and is negligible in high-mass τ τ searches for the current level of experimental uncertainties.