$B$-anomalies in a twin Pati-Salam theory of flavour including the 2022 LHCb $R_{K^{(*)}}$ analysis

We perform a comprehensive phenomenological analysis of the twin Pati-Salam theory of flavour, focussing on the parameter space relevant for interpreting the $B$-anomalies via vector leptoquark $U_1$ exchange. This model provides a very predictive framework in which the $U_1$ couplings and the Yukawa couplings find a common origin via mixing of chiral quarks and leptons with vector-like fermions, providing a direct link between the $B$-anomalies and the fermion masses and mixing. We propose and study a simplified model with three vector-like fermion families, in the massless first family approximation, and show that the second and third family charged fermion masses and mixings and the $B$-anomalies can be simultaneously explained and related. The model has the proper flavour structure to be compatible with all low-energy observables, and leads to predictions in promising observables such as $\tau\rightarrow3\mu$, $\tau\rightarrow\mu\gamma$ and $B\rightarrow K^{(*)}\nu\bar{\nu}$ at Belle II and LHCb. The model also predicts a rich spectrum of TeV scale gauge bosons and vector-like fermions, all accessible to the LHC. In this updated version we have included an extended analysis considering the new 2022 LHCb data on $R_{K^{(*)}}$, which has slightly shifted the preferred parameter space with respect to the 2021 case. The model can still explain the $R_{D^{(*)}}$ anomalies at 1$\sigma$ in a narrow window, however we expect small deviations from the SM on the $R_{K^{(*)}}$ ratios, to be tested in the future via more precise measurements by the LHCb collaboration. We also predict $\Delta R_{D}= \Delta R_{D^{*}}$, with future measurements shifting the world averages to slightly smaller central values.


Introduction
Fundamental fermions in the Standard Model (SM) come in three copies, denoted as "flavours", which share universal gauge interactions but have different masses and mixings, also known as flavour parameters. The origin of flavour in the SM remains as a complete mystery, as it lacks of any dynamical explanation to the high number of flavour parameters and their hierarchical patterns. A further theory of flavour beyond the SM should provide a solution to the long-lasting "flavour puzzle".
Simultaneously, the non-universal structure of such a theory of flavour could leave its imprints in flavour physics observables, which are becoming accessible up to a high precision level in the current generation of colliders and meson factories. Given the prolific history of flavour physics anticipating the discovery of new physics, searching for the origin of flavour in flavour physics is well motivated. In this direction, a conspicuous series of anomalies in flavour observables emerged in the last years.
Back in 2021, when this project was started, the R K ( * ) ratios had been measured by LHCb to be smaller than 1 [1,2], in good agreement with other anomalies in b → sµµ data which were hinting for flavourful new physics (NP) affecting muons rather than electrons. In particular, R [1.1,6] K was alone in more than 3σ tension with the SM prediction. The breaking of SM lepton flavour universality (LFU) was not only suggested by R K ( * ) , but also the ratios R D ( * ) had been measured to show discrepancy with the SM (see the world averages in [3]), hinting for flavourful NP affecting tau leptons. Although no single measurement of R D ( * ) is very significant, the combination of all of them hints for a consistent deviation from the SM prediction with more than 3σ significance. Both LFU ratios together gave rise to a very consistent picture of hierarchical anomalies, where strong NP mainly coupled to the third family interfere with a SM charged current tree-level effect, while weaker NP couple to the much lighter muons, interfering with 1-loop and GIM-suppressed SM neutral currents.
This picture of "B-anomalies" led to important model building efforts by the community during the last 8 years, in order to interpret these anomalies as a low energy signal of a consistent NP model. A massive, electrically neutral Z ′ vector was identified as a possible explanation of the R K ( * ) anomalies (see e.g. [4][5][6][7][8][9]), while different leptoquarks were proposed to address either R K ( * ) or R D ( * ) separately (see e.g. [10][11][12][13]). Interestingly, the vector leptoquark U 1 (3, 1, 2/3) was identified as the only single mediator capable of addressing both B-anomalies simultaneously [12]. However, the gauge nature of U 1 requires to specify a clear ultra-violet (UV) completion that explains its origin. The original ideas by Pati and Salam (PS) [14], led to tensions with unobserved processes such as K L → µe. Instead, an interesting proposal was firstly laid out in the Appendix of [15], and more formally later in [16], following the idea introduced in [17] that color could appear as a diagonal subgroup of a larger SU (3 + N ) × SU (3) ′ local symmetry valid at high energies. The particular choice N = 1 leads to the so-called "4321" gauge symmetry, with a massive colour octet g ′ (8, 1, 0) and a massive Z ′ (1, 1, 0) with suppressed couplings to light SM fermions. Vector-like (VL) fermions need to be introduced in order to obtain effective couplings of (at least) second family fermions to U 1 . The model, even if not minimal, is very predictive and leads to a rich phenomenology in both low-energy and high-p T searches. However, the flavour structure of the model was rather ad-hoc, and it was hinted that the 4321 gauge group could be the TeV scale effective field theory of a complete model addressing more open questions of the SM. In particular, the 4321 model seemed to be a nice playground to connect the picture of B-anomalies with the flavour puzzle of the SM.
Motivated by the desire to link the origin of the B-anomalies with the origin of Yukawa couplings in the SM, one of us proposed a theory of flavour involving a twin Pati-Salam group [21]. Unlike the other models already present in the market [22][23][24], the twin PS treats all three fermion families in the same way. The basic idea is that all three families of SM chiral fermions transform under one PS group, while families of vector-like fermions transform under the other one. The first PS group, broken at a high scale, provides Pati-Salam unification of all SM quarks and leptons, while a fourth family of vector-like fermions transforms under a second PS group, broken at the TeV scale to the SM, as in Fig. 1.1. The full twin Pati-Salam symmetry, together with the absence of a standard Higgs electroweak (EW) doublet, forbids the usual Yukawa couplings for the SM fermions. Instead, effective Yukawa couplings arise through the mixing between SM fermions and vector-like partners. The same mixing leads to U 1 couplings for SM fermions which could address the B-anomalies. This way, B-anomalies and the flavour puzzle are dynamically and parametrically connected. Furthermore, the twin PS model predicts dominantly left-handed (LH) U 1 currents that were preferred by the 2021 picture of B-anomalies [12,25,26], while the other proposals [22][23][24] predict large couplings for right-handed (RH) third family fermion, which lead to tight constraints from high-p T searches.
In this paper, we studied the phenomenology of the simplified twin PS model presented in [21], which turned out to be incompatible with low-energy data. Afterwards, we performed further model building and presented an extended version of the model that can explain the 2021 picture of B-anomalies and address charged fermion masses and mixings, while being compatible with all existing data. However, during the peer-review process of this paper, the LHCb collaboration presented a reanalysis of the LFU ratios R K ( * ) , with the new measurements in the central q 2 shown below (where q 2 denotes the dilepton invariant-mass squared) [27] R [ 2) with correlation ρ = −0.017. Unexpectedly, the updated results are in good agreement with the SM predictions of R [1.1,6] K ( * ) = 1 ± 0.01 [28], as a result of backgrounds in the electron channel which were misidentified in all the previous analyses. Although our model was originally built to explain large deviations from the SM in both R K ( * ) and R D ( * ) , the new experimental data offers the opportunity to further test the model as a legitimate theory of flavour addressing the origin of quark and lepton masses and mixings. Therefore, we present here an updated analysis which includes the new 2022 LHCb data on R K ( * ) , and we confront the new results versus the previous 2021 picture for which the model was intended. Beyond the R K ( * ) ratios, a new combined measurement of R D ( * ) was presented by LHCb in late 2022 [29]. This measurement is in line with the previous experimental data on R D ( * ) , and does not significantly modify the HFLAV average [3], (1.3) which remains at roughly 3σ discrepancy with the SM predictions.
The layout of the remainder of the paper is as follows. In Section 2 we introduce the simplified twin Pati-Salam model as presented in [30], featuring only one vector-like family, and show that it is unable to explain R D ( * ) in a natural way, while being compatible with the stringent constraints from B s −B s mixing. Instead, in Section 3 we present a new, extended version of the twin Pati-Salam model including three vector-like families and a discrete flavour symmetry, which is compatible with low-energy data and high-p T searches. Section 3.2 shows how effective Yukawa couplings for the SM fermions arise in the model, addressing charged fermions masses and mixings. Similarly, Section 3.3 shows the origin of effective couplings between SM fermions and the exotic gauge bosons. Section 3.4 shows the phenomenological analysis and the discussion of the results, including promising signals to test the model in low-energy observables and high-p T searches, along with a study of the perturbativity of the model. Section 4 includes a comparison of our predictions with other models in the market. Finally, we conclude the paper in Section 5.
On the other hand, according to the matter content in Table 1, there are no standard Higgs fields which transform as (1,2,2) under G II 422 , hence the standard Yukawa couplings involving the chiral fermions are forbidden by the twin PS symmetry. These will be generated effectively via mixing with the fourth family of vector-like fermions which only have quantum numbers under the first PS group, G I 422 . This mixing is facilitated by the non-standard Higgs scalar doublets contained in ϕ, ϕ, H, H in Table 1, via the couplings, are the VL mass terms. These couplings mix the chiral fermions with the VL fermions, and will be responsible for generating effective Yukawa couplings for the second and third families. Moreover, the same mixing leads to effective couplings to TeV scale SU (4) I gauge bosons which violate lepton universality between the second and third families, as we shall see.

High scale symmetry breaking
The twin Pati-Salam symmetry displayed in Eq. (2.2) is spontaneously broken to the "4321" symmetry at the high scale M High 1 PeV (the latter bound due to the non-observation of K L → µe [31]), We can think of this as a two part symmetry breaking: (i) The two pairs of left-right groups break down to their diagonal left-right subgroup, via the VEVs ⟨Φ⟩ ∼ v Φ and ⟨Φ⟩ ∼ v Φ , leading to the symmetry breaking, Since the two SU (4) P S groups remain intact, the above symmetry breaking corresponds to (2.8) (ii) Then we assume the second PS group is broken at a high scale via the Higgs H ′ in Table 1, which under G 4422 transform as and develops VEV in its right-handed neutrino component 1 ,⟨ν H ′ ⟩ 1 PeV, leading to the symmetry breaking where SU (4) II P S is broken to SU (3) II c × U (1) II B−L (at the level of fermion representations, chiral quarks and leptons are split 4 II → (3, 1/6) II ⊕ (1, −1/2) II ), while SU (2) I+II R is broken to U (1) I+II T 3R and the abelian generators are broken to The broken generators of SU (4) II P S are associated with PeV-scale gauge bosons that will mediate processes at acceptable rates, beyond the sensitivity of current experiments and colliders. Instead, the further symmetry breaking of G 4321 will lead to a rich phenomenology at the TeV scale, as we shall see. We anticipate that SU (2) I+II L is already the SU (2) L of the SM gauge group, while SM color and hypercharge are embedded in SU (4) On the other hand, the Yukon scalars ϕ and ϕ in Table 1 decompose under G I 422 × G II 422 → G 4422 → G 4321 as ϕ(4, 2, 1; 4, 2, 1) → ϕ(4, 4, 1 ⊕ 3, 1) → ϕ 3 (4, 3, 1 ⊕ 3, −1/6) ⊕ ϕ 1 (4, 1, 1 ⊕ 3, 1/2) , ϕ(4, 1, 2; 4, 1, 2) → ϕ(4, 4, 1, 1 ⊕ 3) → ϕ 3 (4, 3, 1, 1/6) ⊕ ϕ 1 (4, 1, 1, −1/2) , (2.11) plus extra ϕ 3 and ϕ 1 with different values of Y ′ associated to the breaking of the SU (2) I+II R triplet, that we ignore because they do not couple to fermions. The decomposition above is of phenomenological interest, as the Yukons ϕ 3 , ϕ 3 will couple to quarks while ϕ 1 , ϕ 1 will couple to leptons, allowing non-trivial mixing between SM fermions and VL fermions. They will also lead to a nontrivial breaking of G 4321 down to the SM.
The Higgs scalars H and H in Table 1 where α u,d , β u,d , γ u,d , δ u,d are complex elements of two unitary Higgs mixing matrices. The orthogonal linear combinations are assumed to be very heavy, well above the TeV scale in order to sufficiently suppress the FCNCs. We will further assume that only the light Higgs doublet states get VEVs in order to perform EW symmetry breaking, while the heavy linear combinations do not, i.e. we assume that in the Higgs basis the linear combinations which do not get VEVs are very heavy. The discussion of such Higgs potential is beyond the scope of this paper, for the interested reader a deeper discussion was made in Section 3.4 of [21]. In any case, we shall invert the unitary transformations in Eq. (2.14) to express each of the personal Higgs doublets in terms of the light doublets H u , H d , .. This approach will be used in the next section, when constructing the lowenergy quark and lepton mass matrices.

Effective Yukawa couplings and fermion masses
We have already remarked that the usual Yukawa couplings involving purely chiral fermions are absent in the twin PS model. In this subsection, we show how they may be generated effectively via mixing with the vector-like fermions.
We may write the mass terms and couplings in Eq. (2.5) as a 5 × 5 matrix in flavour space (we also define 5-dimensional vectors as ψ T α and ψ c β ), where extra zeroes had been achieved via suitable rotations that leave unchanged the upper 3 × 3 blocks. There are several distinct mass scales in this matrix: the Higgs VEVs ⟨H⟩ and ⟨H⟩, the Yukon VEVs ⟨ϕ⟩ and ⟨ϕ⟩ and the VL fourth family masses M ψ 4 , M ψ c 4 . Assuming the latter are heavier than all the scalars VEVs, we may integrate out the fourth family, to generate effective Yukawa couplings for chiral quarks and leptons which originate from the diagrams in Fig. 2.1. This is denoted as the mass insertion approximation.
As anticipated in [21], the heavy top mass requires ⟨ϕ⟩ /M ψ 4 ∼ 1 and thus it is necessary to go beyond the mass insertion approximation, where the large mixing angle formalism introduced in Appendix A applies. We shall block-diagonalise the mass matrix in Eq. (2.19) in order to obtain the SM Yukawa couplings for the chiral families,  where the mixing angles are given in Appendix A, we define where i, j = 1, 2, 3. We obtain Until the breaking of the twin PS symmetry, the matrix above is Pati-Salam universal, so all fermions of the same flavour share the same effective Yukawa y ψ eff . If we assume a hierarchy of scales for the VL masses On the other hand, since the sum of the two matrices in Eq. (2.20) has rank 1, the first family will be massless. The masses of first family fermions can arise via the mechanism presented in [21], however it leads to no connections with B-physics and the relevant phenomenology discussed here. Therefore, for the phenomenological purposes of this manuscript, we can safely assume the first family to remain massless.
After the symmetry breaking of the twin PS group to G 4321 , the Yukawa couplings x ψ 34 , x ψ c 42,43 and VL masses M ψ 4 , M ψ c 4 remain universal up to small renormalisation group evolution (RGE) effects, however the Yukons decompose in a different way for lepton and quarks as per Eq. (2.11). Due to this decomposition, the mixing angles in Eq. (2.27) are now different for quark and leptons. The VEVs of the Yukons break the SU (4) symmetry relating quarks and leptons, but an accidental SU (2) q c symmetry relating ψ c quarks remains. Hence, the mixing angles (s u c i4 = s d c i4 ) are the same for up and down quarks, and we define q c = u c , d c . On the other hand, the Higgs fields H, H decompose as personal Higgs doublets for the second and third fermion families as per (2.12) and (2.13). The personal Higgses are introduced in order to break the accidental symmetry SU (2) qc , otherwise the mass matrices in the up and down sector would remain identical. A similar discussion applies to charged leptons and neutrinos, and personal Higgses apply in the same way. Mass terms for second and third family fermions will be obtained after the personal Higgses develop a VEV, see Section 2.2. This way, Eq. (2.27) decomposes for each charged sector as the following effective mass matrices, where the Yukawas y ψ 43 and y ψ 24,34 are Pati-Salam universal, and we have approximated all cosines related to ψ c fields to be 1 due to the hierarchy of VL masses in Eq. (2.28). We obtain a similar Dirac-like matrix for neutrinos. In the complete version of the model presented in [21], a further Majorana matrix for the singlet neutrinos ν c is obtained, and all neutrino masses and mixings are accommodated via a type I seesaw mechanism (see full discussion in Section 4.2 of [21]). However, for the sake of simplicity, we will consider massless neutrinos in this simplified framework, as they are of subleading importance for the B-anomalies and for the phenomenological analysis intended for this article.
Due to the fact that VL fermions are much heavier than SM fermions, the fourth row and column, that we have intentionally ignored when writing Eqs. ( (2.32) The CKM matrix is then predicted as We do not address the mixing involving the first family since we are assuming massless first family fermions, as previously discussed. We are however required to preserve V cb as [32] positive in our parameterisation, where in the last step we have approximated the cosines to be 1. We will not fit V ts , V tb , V cs up to the experimental precision, as corrections related to the first family mixing (and CPV phase) are required. In the following we explore the parameters in the mass matrices of Eqs. (2.29), (2.30), (2.31), and its impact over the diagonalisation of the mass matrices: • In very good approximation, the mass of the top quark is given by the (3,3) entry in the first matrix of Eq. (2.29), i.e.
where v SM = 246 GeV and we have applied as in usual 2HDM. If we consider tan β ≈ 10 and α u ≈ 1, then we obtain (2.37) From the expression above, it is clear that very large or maximal s Q 34 ≈ 1 is required in order to preserve a natural y ψ 43 , and to avoid perturbativity issues. Moreover, we will see that maximal values for s Q 34 are also well motivated by the R D ( * ) anomaly, leading to a clear connection between the B-physics and the flavour puzzle only present in this model.
• In the bullet point above, the effective top Yukawa coupling in the Higgs basis has been estimated as y t ≈ 1. By following the same procedure, we can see that all fermion masses can be accommodated with natural parameters. Remarkably, we obtain that all the effective Yukawa couplings are SM-like in the Higgs basis, explaining the observed pattern of SM Yukawa couplings at low-energy.
• The mixing between left-handed quark fields arise mainly from the off-diagonal (2,3) entry in the quark mass matrices, which is controlled by s q c 34 . This mixing can be estimated for each sector by the ratio of the (2,3) entry over the (3,3) entry, i.e.
obtained under the approximation s q c 34 ≈ s q c 24 . Therefore, the model predicts that V cb originates mainly from the down sector, while the mixing in the up sector is small, suppressed by the heavy top mass. The specific values of the mixing angles can be different if we relax s q c 24 ≈ s q c 34 , but the CKM remains down-dominated in any case.
• The lepton sector follows a similar discussion as that of the quark sector. However, the phenomenological relation ⟨ϕ 3 ⟩ ≫ ⟨ϕ 1 ⟩ will lead to smaller angles than those of quarks  24 , where a larger θ e 23 contributing to large atmospheric neutrino mixing is obtained. In this scenario, interesting signals in lepton flavour-violating (LFV) processes such as τ → 3µ or τ → µγ arise, mediated at tree-level by SU (4) I gauge bosons. This is obtained if x ψ c 43 > x ψ c 42 , without the need of any tuning.
• Unlike private Higgs models [33][34][35][36], the personal Higgs VEVs are not hierarchical, all of order 1-10 GeV, with the exception of the top one whose VEV is approximately that of the SM Higgs doublet, as discussed above. The reason is that the fermion mass hierarchies arise from the hierarchies s ψ 34 ≫ s ψ c 24 , s ψ c 34 , which find their natural origin in the hierarchy of VL masses M ψ 4 ≪ M ψ c 4 in Eq. (2.28). The latter simultaneously leads to dominantly left-handed leptoquark currents, as mentioned before.

The low-energy theory G 4321
In this section we shall discuss the G 4321 theory that breaks to the SM symmetry group at low energies G 4321 → G SM , which is achieved via the scalars ϕ 3 (4,3, 1 + 3, −1/6) and ϕ 1 (4, 1, 1 + 3, 1/2) developing the VEVs where v 1 , v 3 ≲ 1 TeV, and analogously for ϕ 3 and ϕ 1 developing VEVs v 3 and v 1 , leading to the symmetry breaking of G 4321 down to the SM gauge group, Here the SU (4) I P S is broken to The physical massive scalar spectrum includes a real colour octet, three SM singlets and a complex scalar transforming as (3,1,2/3). The heavy gauge boson spectrum includes a vector leptoquark U µ 1 = (3, 1, 2/3), a colour octet g ′ µ = (8, 1, 0) also identified as coloron, and Z ′ µ = (1, 1, 0). The heavy gauge bosons arise from the different steps of the symmetry breaking, The gauge boson masses resulting from the symmetry breaking in Eq. (2.28) are a generalisation of the results in [15,16], where we have assumed v 3 ≈ v 3 and v 1 ≈ v 1 for simplicity. The mass of the coloron depends only on v 3 , and the scenario v 3 ≫ v 1 leads to the approximated relation M g ′ ≈ √ 2M U 1 . This way, the coloron can be slightly heavier than the vector leptoquark, which can help to pass the stringent bounds from high-p T searches.
In the original gauge basis, the heavy gauge bosons couple to the EW doublets (including the EW doublets formed by fourth family VL fermions) via the left-handed interactions and also to the EW singlets, although these couplings are suppressed by small mixing angles connected to the origin of second family fermion masses. Therefore, they can be safely neglected 2 . This way, the U 1 couplings will be purely left-handed, which can alleviate the stringent bounds from highp T . Similar couplings are obtained for the VL partners in the conjugated representations, however those couplings are irrelevant for the phenomenology since the conjugated partners do not mix with the SM fermions. The expressions above can be readily written from CP-conjugated notation to L, R notation via the formulae of Appendix E.
The gauge couplings of SU (3) c and U (1) Y are given by where g 4,3,2,1 are the gauge couplings of G 4321 . The scenario g 4 ≫ g 3,2,1 is well motivated from the phenomenological point of view, since here the flavour-universal couplings of light fermions to the heavy Z ′ and g ′ are suppressed by the ratios g 1 /g 4 and g 3 /g 4 , which will inhibit the direct production of these states at the LHC. In this scenario, the relations above yield the simple expressions g s ≈ g 3 and g Y ≈ g 1 for the SM gauge couplings.
A key feature of the gauge boson couplings in Eqs. (2.46-2.48) is that, while the coloron g ′ µ and the Z ′ µ couple to all chiral and VL quarks and leptons, the vector leptoquark U µ 1 only couples to the fourth family VL fermions. However, the couplings in Eqs. (2.46-2.48) are written in the original gauge basis. We shall perform the transformation to the decoupling basis (primed) as per Eq. (2.21), L gauge where α, β = 1, .., 4 and the indexes of the matrices are implicit. We obtain an effective coupling for the third family due to mixing with the fourth family, where we have omitted the fourth column and row for simplicity. The diagrams in Fig. 2.2 are illustrative, however it must be remembered that the mass insertion approximation is not accurate here due to the heavy top mass, instead we have to work in the large mixing angle formalism. In principle, the couplings in (2.51) can simultaneously contribute to both LFU ratios R K ( * ) and R D ( * ) once the further 2-3 transformations required to diagonalise the quark and lepton mass matrices are taken into account. Such transformations split the SU (2) L doublets and lead to different couplings for the different chiral fermions, included in Appendix B.
In Eq. (B.1) it is shown how the leptoquark couplings that contribute to LFU ratios arise due to the same mixing effects which diagonalise the mass matrices of the model, yielding mass terms for the SM fermions. Therefore, the flavour puzzle and the B-physics anomalies are dynamically and parametrically connected in this model, leading to a predictive framework.
Following the same methodology, we obtain the coloron and Z ′ couplings in the basis of mass eigenstates, which can be found in Appendix B.
The flavour-violating couplings of U 1 in Eq. (B.1) are all proportional to mixing between chiral fermions. In principle, such mixing is of order V cb in the down sector, and of order 0.1V cb in the up sector (see discussion in Section 2.3). The small mixing in the up sector leads to a small U 1 2-3 coupling, possibly too small for R D ( * ) , however a deeper analysis was required and we will perform such analysis in the next section. Moreover, flavour-violating couplings involving the coloron and Z ′ could be sizable in the down 2-3 sector, since the CKM is predicted to be originated from the down sector in this model. We shall study whether this is compatible or not with the stringent constraints coming from B s − B s meson mixing.

R K ( * ) and R D ( * )
New contributions to the R D ( * ) and R K ( * ) ratios arise in our model via tree-level contributions mediated by the U 1 vector leptoquark, see the formulae in Appendix D.2 and D.3. After integrating out U 1 , we obtain the following scaling From Eq. (2.52) it can be seen that our contribution to R D ( * ) is proportional to the mixing angle θ u 23 . Such angle is naturally small in this model, roughly O(0.1V cb ) as per Eq. (2.38), due to the fact that the CKM mixing is originated from the down sector. As a consequence, the contribution to R D ( * ) is suppressed. On the other hand, the contribution of U 1 to R K ( * ) is further suppressed by the O(V cb ) mixing angles θ d 23 and θ e 23 , for a total suppression of O(V 3 cb ).

B s −B s mixing
Flavour-violating couplings involving the coloron and Z ′ could be sizable in the 2-3 down sector, since the CKM is predicted to be originated from the down sector in this model. The formulae, the treatment and the bounds obtained from B s −B s mixing are derived in Appendix D.4. The bounds are highly constraining over this model because both the coloron and Z ′ mediate tree-level contributions to ∆M s , which interfere positively with the SM prediction, while the latter are already larger than the experimental result. We estimate that, in order to satisfy the bound ∆M NP s /∆M SM s < 0.11, the 2-3 down-quark mixing needs to satisfy s d

Results in the simplified model
As anticipated in the previous sections, the contribution of the leptoquark to the R D ( * ) anomaly is strongly suppressed by a naturally small mixing angle θ u 23 ≈ m c /m t , leading to a suppression of O(0.1V cb ). In Fig. 2.3a it can be seen that for a typical benchmark mass M U 1 = 3 TeV, a larger s u

23
4V cb is needed in order to address the R D ( * ) anomaly, provided that the 3-4 mixing is maximal. The contribution to R K ( * ) also suffers from an overall suppression of O(V 3 cb ). We can go beyond the natural value of  x ψ c 34 , or reducing the VL mass M ψ c 4 ), which controls the overall size of the off-diagonal (2,3) entry in the effective mass matrices of Eqs. (2.29) and (2.30). This way, we can explore the parameter space of larger 2-3 mixing angles, provided that the experimental value of V cb is preserved through Eq. (2.34), which entangles both quark mixings θ u 23 and θ d 23 . We further assume s Q 34 = s L 34 and s d 23 = s e 23 to simplify the parameter space. Both assumptions are well motivated, the former due to universality of the Yukawa x ψ 34 and VL mass M ψ 4 , the latter due to both mixing angles being proportional to similar parameters, with the mass matrices having the same mass scale.
Our results are depicted in Fig. 2.3a for a spectrum of heavy gauge boson masses compatible with high-p T searches (see Section 3.4.9). We find that for the given benchmark, a small region of the parameter space is compatible with the 2022 data of R K ( * ) , however the 1σ region of R D ( * ) is not compatible with ∆M s . This version of the model was already unable to explain the 2021 data of both LFU anomalies due to the large constraints from tree-level Z ′ and coloron contributions to ∆M s .
In Fig. 2.3b we have varied the VEVs of ϕ 3 and ϕ 1 , effectively exploring the parameter space of gauge boson masses in line with Eq. (2.45). However, we find that the stringent constraints from ∆M s are only alleviated when ⟨ϕ 3 ⟩ 8 TeV, which corresponds to a coloron with mass M g ′ ≳ 50 TeV and a vector leptoquark with mass M U 1 34 TeV, too heavy to address R D ( * ) .
We conclude that the model in this simplified version is over-constrained by large tree-level contributions to ∆M s mediated by the coloron and Z ′ . Such FCNCs arise due to the 2-3 CKM mixing having its origin in the down sector. Moreover, the same small 2-3 mixing angles suppress the contribution of the model to R D ( * ) . However, we shall show that the proper flavour structure to be compatible with all data is achieved in the extended version of the model presented in Section 3.

Extending the simplified twin Pati-Salam theory of flavour
In this section we present an extended version of the simplified twin Pati-Salam model, featuring extra matter content and a discrete flavour symmetry. This new version can achieve the proper flavour structure required to be compatible with all data, solving the problems of the simplified twin Pati-Salam model discussed in Section 2. Firstly, we will introduce the extended version of the model. Secondly, we will revisit the diagonalisation of the mass matrix, leading to the fermion masses and to the new couplings with the heavy gauge bosons. Finally, we will study the phenomenology, showing that the model is compatible with all data while predicting promising signals in flavour-violating observables, rare B-decays and high-p T searches.

New matter content and discrete flavour symmetry
As identified in [18], when one considers a 4321 model with all chiral fermions transforming as SU (4) singlets (fermiophobic framework), three vector-like fermion families can achieve the proper flavour structure to explain the B-anomalies. Such flavour structure can provide a GIM-like suppression of FCNCs, along with large leptoquark couplings that can contribute to the LFU ratios. Hence, as depicted in Table 2, we extend now the simplified model by two extra vector-like families, to a total of three, ψ 4,5,6 (4, 2, 1; 1, 1, 1) (1,1,α) , ψ 4,5,6 (4, 2, 1; 1, 1, 1) (1,1,α 3 ) , ψ c 4,5,6 (4, 1, 2; 1, 1, 1) (1,1,α) , ψ c 4,5,6 (4, 1, 2; 1, 1, 1) (1,1,α 3 ) , where it can be seen that all VL families originate from the first Pati-Salam group, being singlets under the second. They are indistinguishable under the twin Pati-Salam symmetry in Eq. (2.2), however a newly introduced Z 4 flavour symmetry discriminates the sixth family from the fourth and fifth, via different powers of the Z 4 charge α = e iπ/2 . This way, the total symmetry group of the high energy model is extended to The new Z 4 discrete symmetry is introduced for phenomenological purposes, as it will prevent finetuning, reduce the total number of parameters of the model and protect from FCNCs involving the first family of SM-like chiral fermions. Moreover, Z 4 will simplify the diagonalisation of the full mass matrices and preserve the effective Yukawa couplings for SM fermions given in Section 2.3, with specific modifications. The origin of the chiral fermion families is still the second Pati-Salam group, however now they transform in a non-trivial way under Z 4 , Finally, the scalar content is extended by an additional scalar Ω 15 which transforms in the adjoint representation of SU (4) I , whose VEV ⟨Ω 15 ⟩ = T 15 v 15 splits the vector-like masses, and We also include an additional copy of the Yukon ϕ, denoted as ϕ ′ , featuring α 2 charge under Z 4 . The simplified Lagrangian in Eq. (2.5) is extended by the new matter content to 66 15 Ω 15 ψ 6 ψ 6 +λ aa 15 Ω 15 ψ c a ψ c a +λ 66 15 Ω 15 ψ c 6 ψ c 6 + h.c. , (3.5) where i = 2, 3 and a, b = 4, 5 (terms i = 1 and a, b = 6 forbidden by Z 4 ). The symmetry breaking and the decomposition of the different fields proceeds just like in the simplified model, see Section 2.2, however the VEVs of the additional scalars ϕ ′ and Ω 15 play a role in the spontaneous breaking of the 4321 symmetry, and the corresponding gauge boson masses become (assuming v 1,

Effective Yukawa couplings revisited
In this section, we diagonalise the full mass matrix of the extended model, following the same procedure as in Section 2.3, but including the extra matter content of the extended model. We may write the mass terms and couplings in Eq. (3.5) as a 9 × 9 matrix in flavour space (we also define 9-dimensional vectors as ψ α and ψ c β below), where the diagonal mass parameters M Q,L 44,55,66 are splitted for quarks and leptons due to the VEV of Ω 15 , where a = 4, 5, 6. Similar equations are obtained for the ψ c sector, however in the ψ c sector the mass splitting is minimal due to ⟨Ω 15 ⟩ being of order a few hundreds GeV while M ψ c aa are much heavier due to a generalisation of the hierarchy in Eq. (3.14). In Eq. (3.9) we have achieved an extra zero in the (2,7) entry by rotating ψ 2 and ψ 3 , without loss of generality thanks to the zeros in the upper 3 × 3 block (see Section 2.3).
The matrix in Eq. (3.9) features three different mass scales, the Higgs VEVs ⟨H⟩ and ⟨H⟩, the Yukon VEVs ⟨ϕ⟩, ⟨ϕ⟩, ⟨ϕ ′ ⟩ and the VL mass terms M ψ ab and M ψ c ab . We can block diagonalise the matrix above by taking advantage of the different mass scales. Firstly, we diagonalise the 2 × 2 sub-blocks containing the heavy masses M ψ ab and M ψ c ab , 11) and similarly in the ψ c sector. The 4-5 rotations above just redefine the elements in the 4th, 5th, 7th and 8th rows and columns of the full mass matrix, leaving the upper 3 × 3 blocks unchanged (plus we reintroduce the zero in the (2,7) entry by another rotation of ψ 2 and ψ 3 ). Then we perform a further sequence of rotations to go to the decoupling basis, where no large elements appear apart from the diagonal heavy masses (i.e. those terms in the seventh, eighth and ninth rows and columns involving the fields ϕ and ϕ are all absorbed into a redefinition of the heavy masses), and we obtain a block-diagonal matrix similar to that of Eq. (2.20) but enlarged with the fifth and sixth VL families. The total set of unitary transformations is given by The mixing angles controlling the unitary transformations in Eq. (3.12) are given in Appendix A.
The transformations in the ψ c sector of (3.13) can be described by V ψ c 34 V ψ c 24 in good approximation, whose mixing angles are given by Eqs. (A.10) and (A.11). This approximation is accurate as far as the mixing involving the 5th and 6th ψ c fields is further suppressed by a generalisation of the hierarchy in Eq. (2.28) to three vector-like families, namely (3.14) The hierarchy above will preserve most features of the basic simplified model, such as large third family Yukawa couplings arising from mixing with ψ 4 fermions, and small second family Yukawa couplings arising from mixing with ψ c 4 . The couplings of U 1 to chiral fermions will remain dominantly left-handed, since the couplings to ψ c chiral fermions (or equivalently right-handed fermions) will remain suppressed by small mixing angles. On the other hand, the hierarchy M Q,L 44 ≪ M Q,L 55 will provide hierarchical couplings of U 1 to third family and second family fermions, so we anticipate a small contribution to R K ( * ) and a large contribution to R D ( * ) .
We obtain the effective Yukawa couplings for SM fermions by applying the set of unitary transformations in Eqs. (3.12) and (3.13) to the upper 6 × 6 block of (2.19), in the same way as in Eq. (2.25). In this basis (primed), the mass matrix for each charged sector reads (assuming a small x ψ 35 , see Section 3.4.4 for the motivation, and approximating cosines in the ψ c sector to be 1), • Let us impose that the total (2,3) entry in the down quark mass matrix is small, i.e.
Following the discussion of Section 2.3, a natural benchmark is On the other hand, the mixing angle s Q 25 is very relevant for the B-decays and related phenomenology, and we obtain the typical value s Q 25 ≈ 0.2 in Section 3.4, featuring another connection between the flavour puzzle and B-physics in our model. With this input, we obtain In particular, the benchmark in Table 5 where we have considered y ψ 43 = 1, s Q 34 ≈ 1, as required to explain the top mass (see the discussion in the first bullet point of Section 2.3) and we have neglected the (2,3) term proportional to the smaller energy scale ⟨H c ⟩ when compared with the heavier ⟨H t ⟩. This way, we have taken advantage of the new contribution via the 5th family (and of the different hierarchies m c /m t and m s /m b ) to cancel the dangerous down mixing while preserving the CKM via up mixing.
• The situation in the lepton sector is similar due to Pati-Salam universality of the parameters, i.e.
However, the leptonic mixing angles s e c 24 and s e c 34 are smaller than the quark ones due to the phenomenological relation ⟨ϕ 3 ⟩ ≫ ⟨ϕ 1 ⟩. This leads to ⟨H µ ⟩ being above the scale of the muon mass, which predicts a quick growth of lepton mixing in the scenario s e c 34 > s e c 24 . This can be easily achieved in realistic benchmarks. In this scenario, interesting signals arise in LFV processes such as τ → 3µ or τ → µγ, mediated at tree-level by the Z ′ boson, see Section 3.4.5.
Other than the bullet points above, the mass matrices in Eqs. (3.15), (3.16), (3.17) lead to similar predictions as those of the simplified model in Section 2.3.

Vector-fermion interactions in the extended model
In this section we shall compute the vector-fermion couplings involving the heavy gauge bosons U 1 , g ′ , Z ′ . The complete formulae can be found in Appendix B. We omit the couplings of the vector-like partners in the conjugate representations ψ α and ψ c α , since they do not mix with SM fermions.

U 1 couplings
In the original gauge basis, the vector leptoquark couples to the heavy EW doublets via the lefthanded interactions, where similar couplings to the heavy EW singlets ψ c are also present, however they lead to suppressed couplings to SM fermions due to the hierarchy in Eq. (3.14). This way, we obtain purely left-handed U 1 couplings in good approximation. Now we shall apply the unitary transformations in Eq. (3.12) to rotate the fields from the original gauge basis to the decoupling basis (primed), where The 4-5 rotations are different for quarks and leptons due to ⟨Ω 15 ⟩ splitting the mass terms of the VL fermions. They lead to a non-trivial CKM-like matrix for the U 1 couplings, where s θ LQ depends on the angles s Q 45 and s L 45 , obtained from the diagonalisation in Eq. (3.11). The unitary matrix W LQ can be regarded as a generalisation of the CKM matrix to SU (4) or quarklepton space. Similarly to the CKM case, the W LQ matrix is the only source of flavour-changing transitions among SU (4) I states, and it appears only in interactions mediated by U 1 . In this sense, the vector leptoquark, U 1 , is analogous to the SM W ± bosons. Similarly, the Z ′ , g ′ are analogous to the SM Z boson, and we will show that their interactions are SU (4) I flavour-conserving at treelevel. In analogy to the SM, we will denote U 1 transitions as charged currents and Z ′ , g ′ transitions as neutral currents. As in the SM, flavour-changing neutral currents (FCNCs) proportional to the W LQ matrix are generated at loop level. This mechanism was firstly identified in [18] for a similar 4321 framework.
The same mixing that leads to the SM fermion masses and mixings, see Eq. (3.25), also leads to effective U 1 couplings to SM fermions which can contribute to the LFU ratios, where we have considered that s Q,L 35 are small, see Sections 3.2 and 3.4.4. The first family coupling can be diluted via mixing with vector-like fermions, which is parameterised via the effective parameter ϵ (see Appendix F for more details). The couplings above receive small corrections due to 2-3 fermion mixing arising after diagonalising the effective mass matrices in Eqs. (3.15), (3.16), (3.17). It can be seen from Eq. (3.27) that a large coupling (2,3) coupling β cντ arises now, proportional to the large sines s θ LQ , s L 34 and s Q 25 . This solves one important issue of the simplified model, where the flavour-violating couplings β cντ and β bµ where connected to small 2-3 mixing angles, suppressing the contributions of U 1 to the LFU ratios. In any case, the leptoquark couplings that contribute to B-decays arise due to the same mixing effects which diagonalise the mass matrices of the model, yielding mass terms for the SM fermions. This way, the flavour puzzle and the B-anomalies are dynamically and parametrically connected in this model, leading to a predictive framework.

Coloron couplings and GIM-like mechanism
In the original gauge basis, the coloron couplings are flavour diagonal, featuring the following couplings to EW doublets, where i = 1, 2, 3. Now we rotate to the decoupling basis by applying the transformations in Eq. (3.25), (assuming small x ψ 35 as discussed in Section 3.4.4) obtaining Here V Q 45 cancels due to unitarity and due to the g ′ couplings between VL quarks being flavouruniversal in the original basis of (3.28). Therefore, as anticipated before, the CKM-like matrix W LQ does not affect the neutral currents mediated by g ′ (and similarly by Z ′ ). The coloron couplings in 3.29 receive small corrections due to 2-3 mixing arising after diagonalising the effective mass matrices in Eqs. (3.15), (3.16), (3.17), predominantly in the up sector, due to the down-aligned flavour structure achieved in Section 3.2. We obtain similar couplings for EW singlets, however their mixing angles are suppressed by the hierarchy in Eq. (3.14), and so they remain like in the original gauge basis.
The coloron couplings of Eq. (3.29) are flavour-universal if leading to a GIM-like protection from tree-level FCNCs mediated by the coloron. The condition above was already identified in [18], denoted as full alignment limit. However, we have seen that maximal s Q 34 ≈ 1 is well motivated in our model to protect the perturbativity of the top Yukawa, by the fit of the R D ( * ) anomaly, and furthermore it naturally suppresses s Q 35 via a small c Q 34 . The caveat is that if the condition in Eq. (3.30) is implemented, then s Q 16 and s Q 25 would also be maximal, leading to large couplings to valence quarks which would blow up the production of the coloron at the LHC. This fact was already identified in [18], where large s Q 34 was also suggested by the B-anomalies, and a partial alignment limit was implemented, which suppresses FCNCs between the first and second quark families, proportional to the largest off-diagonal elements of the CKM matrix. FCNCs between the second and third families still arise, however we are protected from the stringent constraints of B s −B s meson mixing due to the downaligned flavour structure achieved in Section 3.2. Finally, FCNCs between the first and third families are also under control, as they are proportional to the smaller elements of the CKM matrix.
The GIM-like condition of Eq. (3.31) translates, in terms of fundamental parameters of our model, into which could be naively achieved with natural couplings and M Q 5 , M Q 6 being of the same order, as allowed by the messenger dominance in Eq. (3.14). The couplings and vector-like mass terms can also be chosen differently, as far as Eq. (3.32) is preserved. At the moment, the GIM-like mechanism is accidental. However, Eq. (3.32) suggests that the sixth and fifth family, and also the first and second families, might transform as doublets under a global SU (2) symmetry, enforcing the parametric relations of Eq. (3.32).
A similar treatment of Z ′ couplings can be found in Appendix B, and a similar condition is obtained to suppress LFV between the first and second lepton families, Remarkably, if the condition of Eq. (3.31) is fulfilled, then Eq. (3.33) would also be fulfilled in good approximation thanks to the underlying twin Pati-Salam symmetry, the small breaking effects given by the splitting of VL masses via ⟨Ω 15 ⟩.

Low-energy phenomenology
The twin PS model features a fermiophobic low-energy 4321 theory with a rich phenomenology. Although extensive analyses of general 4321 models have been performed during the last few years, the vast majority of them have been performed in the framework of non-fermiophobic 4321 models [19,20,22,23,37]. Instead, the twin PS model offers a fermiophobic scenario with a different phenomenology. Being a theory of flavour, extra constraints and correlations arise via the generation of the SM Yukawa couplings and the prediction of fermion masses and mixing, including striking signals in LFV processes. Moreover, the underlying twin Pati-Salam symmetry introduces Observable Experiment/constraint Theory expr.
universality (and perturbativity) constraints over several parameters, which are not present in other models. These features motivate a dedicated analysis. We will highlight key observables for which the intrinsic nature of the model can be disentangled from all alternative proposals. All low-energy observables considered are listed in Table 3, with references to current experimental bounds and links to theory expressions.
The benchmark points BP1 and BP2 in Table 5 address the R D ( * ) anomalies and are compatible with the 2021 and 2022 data on R K ( * ) , respectively, plus all the considered low-energy observables and high-p T searches. They provide a good starting point to study the relevant phenomenology, featuring typical configurations of the model, and allow us to confront the 2021 picture of the model versus the new situation with LFU preserved in µ/e ratios. Moreover, they fit second and third family charged fermion masses and mixings, featuring a down-aligned flavour structure with O(0.1) µ − τ lepton mixing. The latter is more benchmark dependent, with the common range being In the forthcoming sections we will assume the couplings of the fundamental Lagrangian to be universal, such as x ψ 34 and x ψ 25 , however their universality is broken by small RGE effects which we estimate in Section 3.4.8 to be below 8%. We neglect the small RGE effects and preserve universal parameters for the phenomenological analysis, in order to simplify the exploration of the parameter space and highlight the underlying twin Pati-Salam symmetry. vs C µµ 10 to 1σ accuracy derived by using 2022 data on R K (blue region), R K * (green region) and B(B s → µ + µ − ) (yellow region). The red contours denote the 1σ, 2σ and 3σ regions of the χ 2 fit. The black dot denotes the best fit point with ∆χ 2 /dof ≈ 0.56.

Model independent analysis of 2022 clean b → sµµ data
This model was originally built to address and relate the 2021 R K ( * ) and R D ( * ) LFU anomalies, while connecting their origin to the origin of Yukawa couplings in the SM. This picture changed completely after the 2022 LHCb update of the R K ( * ) ratios [27], which are now broadly compatible with the SM predictions (see Eq. (1.2)). New 2022 data of B(B s → µ + µ − ) by CMS [48] is also compatible with the SM, while the previous measurements were hinting for values smaller than the SM prediction, including the 2021 measurement by LHCb [49]. In our analysis, we will consider the global average of B(B s → µ + µ − ) experimental data by Allanach and Davighi [50], which is roughly 1σ below the SM prediction B(B s → µ + µ − ) SM = (3.67±0.15)×10 −9 . Instead, other observables such as the B → Kµµ [51] branching fraction and the angular observable P ′ 5 [52,53] show important tensions with the SM, however such anomalies rely on important assumptions about the hadronic uncertainties, and their study is beyond the scope of this manuscript. Nonetheless, an extended discussion can be found in Section 3.4.3.
Following common practice, we describe b → sµµ transitions in terms of the low-energy Lagrangian containing the usual O 23µµ 9 and O 23µµ 10 operators, defined in Eq. (D.10). More details and all the formulae are included in Appendix D.3. For the sake of clarity, we further simplify the notation by removing quark indexes and denote the corresponding NP Wilson coefficients (WCs) as C µµ 9 and C µµ 10 . To the best of our knowledge, no explicit data for the theoretically clean fit of the WCs considering the new SM-like R K ( * ) is given in the literature, motivating our own model independent analysis. The fits presented in the recent analyses [54,55] include observables which are not theoretically clean, and assumptions about the hadronic uncertainties need to be made. Instead, we need to explore how the theoretically clean, SM-like observables constrain the LFUV C µµ 9 = −C µµ 10 which unavoidably arises in our model. In Fig. 3.1 we show the parameter space in the plane (C µµ 9 , C µµ 10 ) preferred by the 2022 R K ( * ) ratios (in the central q 2 ) and the average of B(B s → µ + µ − ). We also display the result of a combined χ 2 fit to all three observables as the red ellipses, denoting 1σ, 2σ and 3σ intervals. Our results show that a small but non-zero value of C µµ 10 is still preferred by B(B s → µ + µ − ). On the other hand, C µµ 9 is compatible with zero, but small positive and negative values are still allowed by the new R K ( * ) ratios at 1σ. In particular, left-handed NP C µµ 9 = −C µµ 10 are not far away from the 1σ region, and our 1-dimensional fit for the latter is with a best fit value of C µµ 9 = −C µµ 10 = −0.0725 with ∆χ 2 /dof ≈ 0.58. Although left-handed NP are still allowed by the new data, the WCs are much smaller than those preferred by 2021 data (see Table 3).
In our model, the left-handed WC in Eq. (3.35) is obtained after integrating out the heavy gauge bosons, with the overall contribution being dominated by U 1 tree-level exchange. We shall constrain such contribution to the 1σ region in Eq. (3.35), and confront the new results against with the previous picture of 2021 data for which the model was developed.

R D ( * ) and R K ( * )
Beyond the contribution to C µµ 9 = −C µµ 10 (see the EFT of the model in Appendix D. 3), our model also generates a contribution to the WC [C * νedu ] 3332 as defined in Appendix D.2. The latter contribution can accommodate existing tensions between the R D ( * ) ratios and the SM. Namely, in terms of fundamental parameters of the model, the deviations from the SM of the LFU ratios scale as follows, where we have fixed the VL masses and the 4321-breaking VEVs to the values of our benchmark (Table 5). This way, the Yukawa couplings above control the contributions to most of the relevant phenomenology, including the LFU ratios. The Pati-Salam universality of x ψ 34 and x ψ 25 provides here a welcome constraint, not present in other 4321 models. In particular, one can see that both R D ( * ) and R K ( * ) are connected via the same parameters and deviations in both are expected, while in other 4321 models the equivalent of x ψ iα decompose in different parameters for quarks and leptons, which decouple R K ( * ) from R D ( * ) .
Following from Eqs. (3.36) and (3.37), the cubic dependence of R K ( * ) on x ψ 25 anticipates that we can suppress the contribution to R K ( * ) , while preserving a large contribution to R D ( * ) thanks to its linear dependence on x ψ 25 . As a consequence, the yellow band of parameter space preferred by 2022 R K ( * ) is just shifted below the orange band of 2021 R K ( * ) in Fig. 3.3b. The 2022  Table 5. The gray (light gray) region denotes the 1σ (2σ) contour of C U 9 as preferred by a global fit to b → sℓℓ data taken from [56]. The yellow (orange) band denotes the 1σ region preferred by R 2022 K ( * ) (R 2021 K ( * ) ). The blue and red stars denote BP1 and BP2 respectively.
is compatible with R D ( * ) at 1σ only in a narrow region of the parameter space. This is encouraging, given the fact that the model was built to address the 2021 tensions in both LFU ratios. However, in order to explain R D ( * ) , small deviations from the SM in the R K ( * ) ratios are unavoidable, to be tested in the future via more precise measurements of LFU by the LHCb collaboration. Moreover, lower central values for R D ( * ) are also expected.
Remarkably, the fact that the twin PS model only generates the effective operator (c L γ µ b L ) (τ L γ µ ν τ L ) implies that both R D and R D * are corrected in the same direction and with the same size. Instead, non-fermiophobic 4321 models also predict the scalar operator (c L b R )(τ R ν τ L ), which leads to a larger correction of R D than that of R D * (about 5/2 larger for the PS 3 model, see Eq. (27) in [22]).

Off-shell photon penguin with tau leptons
The explanation of R D ( * ) in our model is correlated to new contributions to b → sτ τ due to SU (2) L invariance of the U 1 couplings, to be explored in detail in Section 3.4.7. Interestingly, the same couplings that contribute to b → sτ τ also lead to an off-shell photon penguin diagram with tau leptons running in the loop, which generates a lepton universal contribution to the operator O 23ℓℓ 9 = (s L γ µ b L )(lγ µ ℓ) entering in b → sℓℓ transitions, namely which is explicitly correlated to b → sτ τ , as well as to R D ( * ) since SU (2) L invariance implies β sτ ≈ β cντ for the U 1 couplings. Therefore, the scaling is |C U 9 | ∝ (x ψ 34 ) 3 x ψ 25 , just like R D ( * ) . A similar contribution has been studied in the literature in a model independent framework [57][58][59][60], however in our model we need to add the contributions via the extra VL charged leptons E 4,5 , see Fig. 3.2a. Unfortunately, due to the flavour structure of our model, the contributions via VL leptons interfere negatively with the leading contribution via the tau loop, and hence our overall contribution to C U 9 is smaller than in other models. The contribution from E 4 is negligible, but the contribution from E 5 reduces C U 9 by a 20% factor of the tau loop contribution. In our model, R D ( * ) and R K ( * ) are correlated, as can be seen from Eqs. (3.36) and (3.37). Therefore, C U 9 is not only correlated with R D ( * ) but also with R K ( * ) . Given that deviations from 1 in R K ( * ) are now constrained by the new LHCb measurements, our final contribution to C U 9 is constrained to be C U 9 ≈ −0.4, as can be seen in Fig. 3.2b. However, global fits of b → sℓℓ data (see e.g. [54][55][56]), mostly driven by anomalies in Br(B → Kµµ), Br(B s → ϕµµ) and P ′ 5 (B → K * µµ) (see e.g. [51]), prefer a larger value C U 9 ≈ −0.8. Therefore, we conclude that our model is not able to fully address the anomalies in b → sℓℓ via the off-shell photon penguin, although our contribution to C U 9 ameliorates the tensions. Performing a more ambitious analysis would require to make assumptions about the hadronic uncertainties afflicting Br(B → Kµµ), Br(B s → ϕµµ) and P ′ 5 (B → K * µµ), which is beyond the scope of this paper.

B s −B s mixing
In the twin PS model as presented in Section 3, tree-level contributions to B s −B s mixing via 2-3 quark mixing are suppressed due to the down-aligned flavour structure achieved in Section 3.2. A further 1-loop contribution mediated by U 1 has been studied in the literature [18,20,61] for other 4321 models, and vector-like charged leptons are known to play a crucial role. In [18] a framework with three VL charged leptons was considered, however the loop function was generalised from the SM W box, so the bounds where expected to be slightly overestimated. Instead, in [61] the proper loop function was derived, but a framework with only one VL charged lepton was considered. For this work, we have generalised the loop function of [61] to the case of three VL leptons. The 1-loop contribution mediated by U 1 reads,  Table 5 for both panels. The dashed lines show contours of constant C U 9 . The blue region is excluded by the ∆M s bound, the region excluded only due to the contribution via the 5th lepton is also shown in lighter blue for comparison. The blue and red stars denote BP1 and BP2 respectively. The loop function is dominated by the heavy vector-like partners. In particular, in the motivated scenario with maximal s L 34 , the couplings with the fourth family β * sE 4 β bE 4 are suppressed by the small cosine c L 34 . This way, the loop is dominated by E 5 in good approximation, and we can apply the property (3.40) to obtain The loop function grows with x E 5 (see Appendix D.4.1). However, in the limit of large bare mass term M L 5 the effective coupling β * sE 5 ∝ s Q 25 vanishes (since large M L 5 also implies large M Q 5 due to the Pati-Salam symmetry), hence both the contribution to ∆M s and R D ( * ) go away. In Fig. 3.3a we plot δ(∆M s ) defined in Eq. (D.25) in terms of M L 5 , and we vary x ψ 25 in the ranges compatible with R D ( * ) and R 2022 K ( * ) (R 2021 K ( * ) ). We can see that the ∆M s bound requires a vector-like lepton around 1.5-2 TeV in the 2022 case, while 2021 data was pointing to a VL lepton with a mass around 1 TeV.
In Fig. 3.3b we show that Eq. (3.41) is indeed a good approximation, up to small interference effects between the 4th and 5th family contributions in the small x ψ 34 region, where the fourth lepton is lighter. We also show the parameter space compatible with ∆M s and the LFU ratios in our benchmark scenario. In particular, ∆M s turns out to be the largest constraint over the parameter space other than R 2022 K ( * ) .

LFV processes τ → 3µ
The partial alignment condition of Eq. (3.33) allows for Z ′ -mediated FCNCs in τ µ processes, due to the fact that the model predicts significant mixing between the muon and tau leptons. This is a crucial prediction of the twin PS theory of flavour, not present in general 4321 models. Of particular interest is the process τ → 3µ, which receives a tree-level Z ′ contribution that grows with the τ µ mixing angle s e 23 . Beyond the latter, τ → 3µ also receives a U 1 -mediated 1-loop contribution The effective coupling β D 5 µ is proportional to s L 25 ≈ 0.1, which provides a further suppression of O((s L 25 ) 3 ) that renders the loop negligible against the much larger tree-level Z ′ -mediated contribution. The typical benchmark s L 25 ≈ 0.1 naturally suppresses the µµZ ′ coupling, keeping the Z ′ contribution to τ → 3µ under control, and simultaneously protects from Z ′ → µµ at LHC (see Section 3.4.9).
As depicted in Fig. 3.4a, the Z ′ contribution dominates over the U 1 contribution, and the regions of the parameter space with very large s e 23 are already excluded by the experiment. We have chosen to plot the results of the 2022 case only, since this observable depends mostly on s e 23 and there is little variation with 2021 data. The Belle II collaboration will test a further region of the parameter space [62], setting the bound s e 23 < 2.8V cb if no signal is detected. In all UV incomplete 4321 models (such as [16,18]) the µ − τ mixing is unspecified, so only the small U 1 signal is predicted. Therefore, the large Z ′ signal offers the opportunity to disentangle the twin Pati-Salam model from other 4321  Table 5. The dashed lines show contours of constant C U 9 . The regions excluded by LFV violating processes are displayed. The blue (red) star shows BP1 (BP2). proposals.
As depicted in Fig. 3.5b, τ → 3µ is the most constraining signal over the parameter space out of all the LFV processes, provided that the 2-3 charged lepton mixing is O(0.1).

τ → µγ
The process τ → µγ receives 1-loop contributions via both Z ′ and U 1 , as depicted in Fig. D.2, with formulae reported in Appendix D.5.1. Provided that the 3-4 mixing is maximal, the U 1 loop is dominated by the 5th vector-like quark, and in this situation the couplings β * D 5 µ β D 5 τ are controlled by x ψ 25 . The Z ′ loop is dominated by chiral leptons, in particular by the τ lepton, since the coupling ξ τ τ is maximal while ξ µµ is suppressed. In this scenario, the overall Z ′ contribution is controlled by ξ τ µ which grows with the µ − τ mixing angle s e 23 , and the variation via x ψ 25 is minimal. In Fig. 3.4b we can see that the Z ′ contribution dominates the branching fraction in the range of large x ψ 34 motivated by R D ( * ) , leading to the predictions for B (τ → µγ) being one/three orders of magnitude below the current experimental limit depending on the value of s e 23 . We have also included the projected bound by Belle II (50 ab −1 ) [62], which will partially test the parameter space. In the 4321 models of [16,18] the µ − τ mixing is unspecified, so only the blue U 1 signal is predicted. For non-fermiophobic models, this signal is largely enhanced via a chirality flip with the bottom quark running in the loop [19,20,37,61,63], predicting a larger signal B (τ → µγ) ≈ 10 −8 . Instead, our Z ′ signal lies below, offering the opportunity to disentangle the twin Pati-Salam model from all other proposals.

B s → τ µ, B → Kτ µ and τ → µϕ
The vector leptoquark U 1 mediates tree-level contributions to flavour-violating (semi)leptonic Bdecays to (kaons), taus and muons. The experimental bound for B s → τ µ was obtained by LHCb   Table 5 for both panels, and current exclusion limits are shown. [41], while for B → Kτ µ experimental bounds are only available for the decays B + → K + τ µ [42]. The process τ → µϕ receives tree-level contributions from both U 1 and Z ′ , see Appendix D.3 and D.5. However, τ → µϕ turns out to be suppressed by the small effective couplings β sµ ∝ s Q 25 s L 25 and ξ ss ∝ s Q 25 2 and we find B (τ → µϕ) ≈ 10 −9 , roughly two orders of magnitude below the current experimental bounds, and just below the future sensitivity of Belle II.
As can be seen in Fig. 3.5b, B s → τ + µ − implies the largest constraint over the parameter space out of all semileptonic LFV processes involving τ leptons, followed by B + → K + τ + µ − and τ → µϕ. The present experimental bounds lead to mild constraints over the parameter space compatible with R D ( * ) . As depicted in Fig. 3.5a, the 2021 region for B s → τ + µ − was partially within LHCb projected sensitivity, but the 2022 region will mostly remain untested.

K L → µe
The LFV process K L → µe sets a strong constraint over all models featuring a vector leptoquark U 1 with first and second family couplings [64], (3.43) The first family coupling β de can be diluted via mixing with vector-like fermions, which we parameterised via the effective parameter ϵ in Eq. (3.27), so that β se ≈ s Q 16 s L 16 ϵ. The mechanism to perform this and the definition of ϵ in terms of fundamental parameters of the model is included in Appendix F. In Fig. 3.6a we can see that for the 2022 case, some region of the parameter space is compatible with K L → µe without the need of diluting the coupling. Instead, for the benchmark values BP1 and BP2, a mild suppression is required. In Appendix F, a benchmark with the fundamental parameters   Table 5. Current exclusion limits are displayed, along with their future projections. The blue (red) star shows BP1 (BP2).
of the model that provide such a suppression is included. This signal is a direct consequence of the underlying twin PS symmetry and the GIM-like mechanism, which lead to quasi-degenerate mixing angles s Q 16 ≈ s L 16 that are equal to their 25 counterparts, and as a consequence β de ̸ = 0. Therefore, it is not present in other 4321 models [16,[18][19][20]37].

Tests of universality in leptonic τ decays
NP contributions to R D ( * ) commonly involve large couplings to τ leptons, which can have an important effect over LFU ratios originated from τ decays. Such tests are constructed by performing ratios of the partial widths of a lepton decaying to lighter leptons and/or hadrons. We find all ratios in our model to be well approximated by (see Appendix D.6), where β bτ ≈ cos θ LQ assuming maximal 3-4 mixing. Therefore, it can be seen as a constraint over the β bτ coupling, and hence is not directly related to R K ( * ) so we do not plot two bands here. The high-precision measurements of these effective ratios only allow for per mille modifications, see the HFLAV average [3] in Table 3. As depicted in Fig. 3.6b, this constraint sets the lower bound M U 1 ≳ 2.2 TeV for sin θ LQ = 1/ √ 2 and g 4 = 3.5. This bound becomes more restrictive for cos θ LQ ≈ 1, or equivalently β bτ ≈ 1, for which we find M U 1 ≳ 3.3 TeV if g 4 = 3.5 and M U 1 ≳ 2.9 TeV if g 4 = 3.

B s → τ τ and B → Kτ τ
The explanation of R D ( * ) in our model requires a large U 1 contribution to the 4-fermion operator (c L γ µ b L )(τ L γ µ ν τ L ). Therefore, large contributions to the rare decays B s → τ τ and B → Kτ τ arise via SU (2) L invariance of the U 1 couplings. The respective branching fractions are of order 10 −7 in the SM and mild upper bounds have been obtained by LHCb [44] and BaBar [45], respectively.
In Fig. 3.7, we plot the branching fractions as a function of x ψ 34 , while x ψ 25 is varied in the ranges compatible with 2021 and 2022 R K ( * ) , respectively. We find that the predictions are far below the current bounds, however they lie closer to the expected future bounds from LHCb and Belle II data [62,65]. This prediction is different in non-fermiophobic 4321 models [19,20,37], where these contributions are enhanced and all the parameter space will be tested in B + → K + τ + τ − by Belle II.

B → Kνν
The U 1 leptoquark does not contribute at tree-level to b → sνν transitions, and the tree-level exchange of the Z ′ is suppressed due to the down-aligned flavour structure. However, loop-level corrections can lead to an important enhancement of the channel B → Kν τντ [20]. We parameterise corrections to the SM branching fraction as where the EFT and the Wilson coefficients are defined in Appendix D.7. As depicted in Fig. D.3, the main contributions are a semileptonic box diagram mediated by U 1 and a triangle diagram correction to the flavour-violating bsZ ′ vertex, plus the RGE-induced contribution via the U 1 -mediated operator (s L γ µ b L ) (τ L γ µ τ L ), denoted as C RGE ν,U . The former two 1-loop contributions are dominated by the fifth VL charged lepton and grow with its bare mass, M L 5 . This way, the overall contribution to B → Kνν can be sizable, yielding up to O(1) corrections with respect to the SM value, as depicted in Fig. 3.8. The details of the calculation are found in Appendix D.7. For low M L 5 , the value of δB(B → K ( * ) νν) corresponds to C RGE ν,U . For large M L 5 , however, we have seen that stringent constraints from B s −B s meson mixing play an important role, see Section 3.4.4. This constraint is depicted as the hatched region in Fig. 3.8, correlating B → Kνν and ∆M s , a feature which has not been highlighted in other analyses. In particular, ∆M s rules out the region where δB(B → K ( * ) νν) can reach values close to current experimental limits. Nevertheless, the Belle II collaboration will measure B(B → K ( * ) νν) up to 10% of the SM value [62], hence testing all the parameter space of the model.
Our signal of B → K ( * ) νν also offers a great opportunity to disentangle our twin PS framework from non-fermiophobic 4321 models and the PS 3 model [19,20,22,37], as they predict a much smaller signal (see Fig. 4.4 of [20] and compare their purple region with our Fig. 3.8).

Perturbativity
The explanation of the R D ( * ) anomaly requires large mixing angles s Q 34 and s L 34 , which translate into a sizeable Yukawa coupling x ψ 34 , thus pushing the model close to the boundary of the perturbative domain. Perturbativity is a serious constraint over our model, since we need the low-energy 4321 theory to remain perturbative until the high scale of the twin Pati-Salam symmetry. When assessing the issue of perturbativity, two conditions must be satisfied: • Firstly, the low-energy observables must be calculable in perturbation theory. For Yukawa couplings, we consider the typical bound x ψ 34 < √ 4π. Regarding the gauge coupling g 4 , standard perturbativity criteria imposes the beta function criterium [66] |β g 4 /g 4 | < 1, which yields g 4 < 4π √ 3/ √ 28 ≈ 4.11.
• Secondly, the couplings must remain perturbative up to the energy scale of the UV completion, i.e. we have to check that the couplings of the model do not face a Landau pole below the energy scale of the twin PS symmetry, namely µ ≈ 1 PeV.
The phenomenologically convenient choice of large g 4 is not a problem for the extrapolation in the UV, thanks to the asymptotic freedom of the SU (4) gauge factor (see The Yukawas x ψ 25 and x ψ 16 are not dangerous as they are order 1 of smaller. The problematic Yukawa is x ψ 34 , which is required to be large in order to both the LFU ratios, and also it is connected with the physical mass of the fourth lepton as per Eq. (A.12). Large λ 5 15 is also required to obtain a large splitting of VL masses, which leads to a large θ LQ as required by R D ( * ) . Fig. 3.9 shows that the Yukawas of our benchmark scenario remain perturbative up to the high energy scale µ ≈ 1 PeV, thanks to the choice of a large g 4 = 3.5. However, we have checked that the Landau pole is hit when x ψ 34 > 2.5, hence this region should be considered as disfavoured by perturbativity.
The small RGE effects that break the PS universality of the Yukawa couplings are below 8% in any case, hence the universality of the couplings is preserved at the TeV scale in good approximation.

High-p T signatures
General 4321 models predict a plethora of high-p T signatures involving the heavy gauge bosons and at least one family of vector-like fermions, requiring dedicated analyses such as those in [18,20,67]. In particular, our model predicts a similar high-p T phenomenology as that of [18], which also considers effective U 1 couplings via mixing with three families of vector-like fermions. However, the bounds obtained in the high-p T analysis of [18] are outdated. Moreover, certain differences arise due to the underlying twin Pati-Salam symmetry in our model, plus the different implementation of the scalar sector and VEV structure. Furthermore, we anticipate that some bounds obtained in [20,67] might be overestimated for our model, as they usually consider large couplings to right-handed third family fermions, motivating a dedicated analysis.
We have included the particle spectrum of our benchmark scenario in Fig. 3.10a, as a typical configuration for the masses of the new vectors and fermions. Table 3.10b shows the main decay channels of the new vector bosons, which feature large decay widths Γ/M due to all the available decay channels to vector-like fermions, plus the choice of large g 4 = 3.5 close to perturbativity

Particle
Decay mode bounds.
In this section, we revisit some of the most simple collider signals, such as coloron dijet searches and Z ′ dilepton searches. We will also comment on U 1 searches, coloron ditop searches and vectorlike fermions. We will point out the differences between our framework and general 4321 models, motivating a future manuscript dedicated to specific high-p T signals of the twin PS model.

Coloron signals
The heavy colour octet has a large impact over collider searches for 4321 models, and its production usually sets the lower bound on the scale of the model. In our case, the heavy coloron has a gauge origin, hence the coloron couplings to two gluons are absent at tree-level, reducing the coloron production at the LHC. Moreover, in the motivated scenario ⟨ϕ 3 ⟩ ≫ ⟨ϕ 1 ⟩, the coloron is slightly heavier than the vector leptoquark at roughly M g ′ ≈ √ 2M U 1 , helping to suppress the impact of the coloron over collider searches while preserving a slightly lighter U 1 for R D ( * ) . In the scenario g 4 ≫ g 3,1 , the coupling strength of the coloron is roughly g 4 , which receives NLO corrections via the K-factor [68,69] K NLO ≈ 1 + 2.65 g 2 4 16π 2 + 8.92 We have computed the coloron production cross section from 13 TeV pp collisions with Madgraph5   destructively with the flavour-universal term, allowing for a certain cancellation of the left-handed couplings to light quarks. However, the coloron is still produced via the flavour-universal couplings to right-handed quarks. We estimate analytically the branching fraction to all SM quarks excluding tops, and then we compute the total cross section via the narrow width approximation. Finally, we confront our results with the limits for a qq-initiated spin-1 resonance provided by CMS in Fig. 10 of [72]. The results are displayed in Fig. 3.12a, where we have varied the coupling to light LH quarks κ qq and fixed the rest of parameters as in Table 5. We find bounds ranging from M g ′ 2.5 TeV when κ qq ≈ 0 and M g ′ 3 TeV when κ qq ≈ g 2 s /g 2 4 . These bounds are slightly milder than those obtained in [20], the reason being that in [20] right-handed bottom quarks are assumed to couple maximally to the coloron, while in our model this coupling is suppressed.
We expect to find more stringent bounds in resonant coloron production with tt final states, due to the maximal couplings of the coloron to the third generation EW quark doublet. According to the recent analysis in [20], our benchmark scenario would lie below current bounds, due to the large decay width Γ g ′ /M g ′ ≈ 0.5 provided by extra decay channels to TeV scale vector-like quarks. The limit over the coloron mass is roughly 3.5 TeV, however this bound might be overestimated again for our model due to the different description of RH third family quarks. Reconstructing the tt channel requires a dedicated analysis and a different methodology, which is beyond the scope of this article, and we leave it for a future work.

Z ′ signals
For the Z ′ boson, the flavour-universal couplings to valence quarks are more heavily suppressed than those of the coloron, via the small ratio g 2 Y /g 2 4 . Therefore, cancellation between the lefthanded mixing term proportional to s Q 25 ≈ s Q 16 and the flavour-universal one is not possible here. In contrast with the coloron, the large LH couplings to bottoms can play a role in Z ′ production. The production cross section is estimated via the same methodology as for the coloron above. We do not consider any NLO corrections in this case, following the methodology of [67]. In Fig. 3.11b we , where q = Q 1 , Q 2 . The rest of parameters are fixed as in Table 5 for both panels. The exclusion bound from CMS is shown in green. show that the production via bottoms is larger than the production via valence quarks for a light Z ′ , however the production via valence quarks is bigger for M Z ′ 2 TeV, and shall not be neglected as it commonly happens in the literature (see e.g. [67]).
We estimate the branching fraction to muons and taus, and we compute the total decay width via the narrow width approximation. We confront our results with the limits from the dilepton resonance searches by ATLAS, Fig. 4 Fig. 7 (c) of [74] for taus. We display the results in Figs. 3.13a and 3.13b. In Fig. 3.12b we see that these processes, along with coloron dijet searches, mildly constrain the region of large x ψ 25 . Ditau searches are more competitive than dimuon searches or coloron dijet searches, due to the branching fractions to muons and light quarks being suppressed by mixing angles s Q,L 25 ∼ O(0.1). Instead, the ditau channel is enhanced by maximal 3-4 mixing, and sets bounds of roughly M Z ′ > 1.5 TeV, see Fig. 3.13b.

U 1 signals
Leptoquark pair-production cross sections at the LHC are dominated by QCD dynamics, and thus are largely independent of the leptoquark couplings to fermions. Therefore, we are able to safely compare with the analyses of Refs. [20,67]. A certain model dependence is present in the form of non-minimal couplings to gluons, however these couplings are absent in models where U 1 has a gauge origin. According to Fig. 3.3 of [20], current bounds over direct production exclude M U 1 < 1.7 TeV, and the future bound is expected to exclude M U 1 < 2.1 TeV if no NP signal is found during the high-luminosity phase of the LHC.
An important constraint over U 1 arises from modifications of the high-p T tail in the dilepton invariance mass distribution of the Drell-Yan process pp → τ + τ − + X, induced by t-channel U 1 exchange [12,15,20,67]. This channel is well motivated by the U 1 explanation of R D ( * ) , which unavoidably predicts a large bτ U 1 coupling. The scenario β R bτ = 0 considered in the study of [20, 67] fits well the twin Pati-Salam framework, up to a re-scaling of the U 1 coupling strength as g U → g U β L bτ , in order to account for the fact that our β L bτ is not maximal but β L bτ ≈ c θ LQ ≈ 0.67 in our benchmark scenario, where we obtain g U ≈ 2.3. According to the left panel of Fig. 3.3 in [20], the 3 TeV leptoquark of our benchmark lies well below the current bounds, but within projected limits for the high luminosity phase of LHC. Finding U 1 much below 3 TeV enters in tension with pp → g ′ → tt as explained before, due to the approximate relation M g ′ ≈ √ 2M U 1 that entangles the masses of U 1 and the coloron (although the bound is probably overestimated for our dominantly left-handed model).
The twin Pati-Salam model could provide a good U 1 candidate for the 3σ excess at CMS [75] pointing to a 2 TeV U 1 leptoquark in the well motivated channel pp → U 1 → τ τ , once the extra decay channels to vector-like fermions are considered, assuming that the bound from pp → g ′ → tt is indeed overestimated.

Vector-like fermions
The presence of vector-like fermions is of fundamental importance to discriminate between the different implementations of the 4321 model addressing the B-anomalies. A common constraint arises from ∆F = 2 transitions at low energies, which require that the vector-like charged lepton that mixes with muons is light (see Fig. 3.3a). The natural mass of the quark partner of L 5 should not lie far away due to the approximate Pati-Salam universality, the small breaking effects given by the VEV ⟨Ω 15 ⟩. In particular, in our benchmark scenario we obtained M L 5 ≈ 0.8 TeV and M Q 5 ≈ 1.2 TeV. The flavour structure of the model naturally predicts that both Q 5 and L 5 are strongly coupled to the third generation of SM fermions.
The twin Pati-Salam model features also L 4 and Q 4 as a relevant pair of vector-like fermions, which mix maximally with the third generation in order to obtain the large couplings required for R D ( * ) , and also to fit the top mass without perturbativity issues. This implies that their bare mass terms in the original Lagrangian are small, therefore their physical masses are dominated by x ψ 34 ⟨ϕ 3,1 ⟩, see Eq. (A.6) and (A.12). In the motivated scenario ⟨ϕ 3 ⟩ ≫ ⟨ϕ 1 ⟩ which slightly suppresses the production of the coloron, we found L 4 to be very light, featuring roughly 600 GeV in our benchmark. Instead, Q 4 can lie above 1 TeV, featuring roughly 1.2 TeV in our benchmark. The couplings of L 4 to SM fermions are smaller than those of L 5 , but it is dominantly coupled to third generation fermions.
Interestingly, CMS recently performed a search for the vector-like leptons of the 4321 model [76], finding 2.8σ preference for a vector-like lepton with 600 GeV, however the analysis assumes EW production only and maximal couplings to the third generation. If Z ′ -assisted production is included, L 5 a with mass at 800 GeV could be a good candidate for the anomaly. Furthermore, L 4 at 600 GeV could also provide a good fit once not maximal couplings are considered, however this requires verification in a future dedicated analysis. Non-fermiophobic 4321 models, such as [19,20], predict a heavier vector-like lepton, while [18] also predicts L 5 at around 800 GeV but a heavier L 4 . Regarding the sixth vector-like fermions L 6 and Q 6 , we expect them to have similar masses as L 5 and Q 5 in order to preserve the GIM-like suppression of 1-2 FCNCs, they are feebly coupled to first generation fermions but not to the second or third generation.
Current bounds on vector-like quarks lie around 1 TeV, however they are usually very model dependent. Our vector-like quarks are pair produced through gluon fusion and through the decay of the coloron, which is very likely to be kinematically allowed. Their decays leave a large amount of third generation fermions in the final state, following a similar mechanism as the one discussed in [18]. The twin Pati-Salam model naturally predicts light vector-like quarks with masses around 1 TeV, which is a feature not present in all other 4321 models and could motivate specific searches. Table 4 includes a simplified set of observables that allows to disentangle the twin PS model from the ones that are already in the market. A further discussion can be found in the two following sections.

Non-fermiophobic 4321 models (including the PS 3 model)
The twin PS model is built as a fermiophobic framework, where all chiral fermions are singlets under the TeV scale SU (4). This is a crucial difference between our model and the non-fermiophobic 4321 models [19,20,37,61,63] and its UV completions (including the PS 3 model [22]), where the third family of chiral fermions transforms under the TeV scale SU (4). This implies large left-and right-handed third family couplings to SU (4) gauge bosons. By contrast, in our theory of flavour, the right-handed couplings of SM fermions to U 1 (and also to Z ′ and g ′ ) arise via small mixing angles connected to the origin of second family fermion masses, hence the twin PS model predicts dominantly left-handed U 1 couplings. The low-energy phenomenology between both approaches is radically different.
In terms of the charged current anomalies R D ( * ) , the twin PS model only predicts the effective operator (c L γ µ b L )(τ L γ µ ν τ L ), and hence both R D and R D * are corrected in the same direction and with the same size. Instead, non-fermiophobic 4321 models also predict the scalar operator (c L b R )(τ R ν τ L ). Due to the presence of this operator, the effect on ∆R D is larger than on ∆R D * (about 5/2 larger for the PS 3 model, see Eq. (27) in [22]   The dash (-) indicates that the observable was not considered or numbers were not given in the corresponding references. In the high-p T row we broadly refer to how constrained is the model by high-p T searches.
Another key observable is B → Kνν, for which the twin PS model predicts a larger branching fraction that will be fully tested by Belle II. Instead, non-fermiophobic 4321 models predict a smaller branching fraction, see Fig. 4.4 of [20] and compare their purple region with our Fig. 3.8. Moreover, in our analysis of B → Kνν we have highlighted correlations with B s −B s mixing due to the loops being dominated by the same VL charged lepton, a feature which is missing in the analysis of [20].
Regarding the rest of the observables, broadly speaking the twin PS model predicts larger branching fractions for LFV processes. The exception is τ → µγ, which is enhanced in nonfermiophobic models via a chirality flip with the bottom quark running in the loop. The rare decays B s → τ τ and B → Kτ τ are larger in non-fermiophobic models due to the presence of scalar operators connected to the third family RH couplings. The LHCb and Belle II collaborations will test regions of the parameter space and allow to disentangle between the different 4321 approaches.
High-p T searches also offer a window to disentangle both approaches, since most of the constraints afflicting non-fermiophobic 4321 scenarios are relaxed in the dominantly left-handed scenario of the twin PS model. Particularly relevant are also the different implementations of vector-like fermions.
Regarding the PS 3 model, it predicts the same signals as non-fermiophobic 4321 models plus enhanced LFV. It can be disentangled from the twin PS model via the set of observables already discussed.

Fermiophobic 4321 models
Up to our knowledge, the only fermiophobic 4321 model proposed in the literature is that of Ref. [16], whose phenomenology was studied in detail in [18]. This model presents a simplified fermiophobic scenario with a rather assumed flavour structure motivated by the phenomenology, including an ad-hoc alignment of SM-like Yukawas and VL-chiral fermion mixing. Furthermore, [16] does not address the question of quark and lepton masses (is not a theory of flavour), unlike the model proposed here. It lacks from quark-lepton unification of SM fermions and leads to a less predictive framework than the twin PS model.
The twin Pati-Salam model leads to an effective fermiophobic 4321 model. However, the underlying twin PS symmetry implies universality of key parameters, leading to extra constraints and correlations between observables, which are not present in the analyses of [16,18]. For example, R D ( * ) and R K ( * ) are correlated here due to the universality of x ψ 25 and x ψ 34 , leading to quasi-universal mixing angles s Q,L 25 and s Q,L 34 . By contrast, such mixing angles are free parameters in [16,18] and one can explain R D ( * ) without giving any contribution to R K ( * ) . As a consequence, it was required a dedicated analysis to show that R D ( * ) can be explained in the twin PS model while being compatible with the recent data on R K ( * ) , as we did in this paper.
Since the twin PS model is a theory of flavour, while [16,18] is not, new signals are predicted in LFV processes connected to the origin of fermion masses and mixings. The twin PS model predicts non-vanishing µ − τ mixing, leading to striking signals in τ → 3µ and τ → µγ close to current experimental bounds, summarised in our Figs. 3.4a and 3.4b. The large contributions to the branching fractions of τ → 3µ and τ → µγ are mediated by the Z ′ boson, see the purple region in Figs. 3.4a and 3.4b. On the other hand, in general fermiophobic 4321 models [16,18] only a much smaller 1-loop U 1 mediated signal is predicted. This signal was not computed in Refs. [16,18] as it is very small compared to the experimental bounds, but we have computed it here for the sake of comparison, and it is depicted as the blue region in Figs. 3.4a and 3.4b. In a similar way, we obtain B(τ → µϕ) two orders of magnitude larger than in [16,18].
Finally, the fermion mixing predicted by the twin PS model avoids current constraints coming from CKM unitarity, ∆F = 2 and EW precision observables presented in [78]. The reasons are the absence of SM-like Yukawa couplings for chiral fermions in the original basis (as they will be generated indeed via this mixing), along with the fact that VL quark EW doublets and SM quark EW singlets do not mix, hence the VL quark doublet remains unsplitted. Remarkably, this is different from [16,18], where mixing between the SM quark singlets and the VL (right-handed) quark doublet was induced due to the presence of the SM-like Yukawa couplings for chiral fermions, leading to possible splitting of the VL quark doublet, which constrains the mixing angles for third generation quarks according to the analysis in [78].

Conclusions
We have performed a comprehensive phenomenological analysis of the twin Pati-Salam model, which is capable of explaining anomalies in LFU ratios of B-decays, while simultaneously accounting for the fermion masses and mixings of the SM. The basic idea of this model is that all three families of SM chiral fermions transform under one PS group, while families of vector-like fermions transform under the other one. Vector leptoquark couplings and SM Yukawa couplings emerge together after mixing of the chiral fermions with the vector-like fermions, thereby providing a direct link between B-physics and fermion masses and mixings. The model was originally built to explain the 2021 picture of LFU anomalies in the R K ( * ) and R D ( * ) ratios. In this updated version we have included an extended analysis considering the new 2022 LHCb data on R K ( * ) , which has shifted the preferred parameter space with respect to the 2021 case. The model can still explain the R D ( * ) anomalies at 1σ in a narrow window while being compatible with all data, however we expect small deviations from the SM on the R K ( * ) ratios, to be tested in the future via more precise measurements of LFU by the LHCb collaboration. We also predict ∆R D = ∆R D * , with future measurements shifting the world averages to slightly smaller central values.
Firstly, we presented a simplified version of the model that is able to explain second and third family charged fermion masses and mixings via effective Yukawa couplings, which arise naturally from mixing effects with a fourth vector-like family of fermions. However, with only a single vectorlike family the model is unable to explain the B-anomalies in a natural way, as it does not achieve the flavour structure required by 4321 models, and hence is over constrained by flavour-violating processes such as B s −B s meson mixing. The latter are mediated by a heavy colour octet and a Z ′ that also acquire flavour-violating couplings with chiral fermions.
We then extended the simplified model to include three vector-like families, together with a Z 4 discrete symmetry to control the flavour structure. This version of the model allows for larger flavour-violating and dominantly left-handed U 1 couplings as required to address R D ( * ) , thanks to mixing between a fourth and fifth vector-like families which also mix with the second and third generations of SM fermions. A sixth vector-like family is included to mix with the first SM generation, for the sake of suppressing any FCNCs involving first and second generation fermions. The mechanism resembles the GIM suppression of FCNCs in the SM, featuring a similar Cabbibo-like matrix which is present in leptoquark currents, but not in neutral currents mediated by the coloron and Z ′ .
As in the simplified twin Pati-Salam model, the origin of second and third generation charged fermion masses and mixings remains addressed via effective Yukawa couplings, featuring now a down-aligned flavour structure in the 2-3 sector (requiring a mild tuning as described in the main text) that protects from the dangerous tree-level contributions to B s −B s meson mixing. Non-zero 2-3 mixing in the charged lepton sector is also predicted, leading to interesting signals in τ → 3µ and τ → µγ, mostly due to Z ′ exchange, which are close to present experimental bounds in some region of the parameter space. Signals in LFV semileptonic processes mediated by U 1 at tree-level are found to lie well below current experimental limits, with the exception of K L → eµ which constrains a small region of the parameter space. However, this tension can be alleviated if the first family U 1 coupling is diluted via mixing with vector-like fermions. Tests of LFU in tau decays set important bounds over the mass of U 1 depending on its coupling to third generation fermions. Contributions of U 1 to the rare decays B s → τ τ and B → Kτ τ are broadly below current and projected experimental sensitivity. Instead, the rare decay B → K ( * ) νν offers the opportunity to fully test the model in the near future, since Belle II is expected to cover all the parameter space compatible with the LFU ratios. Remarkably, the model can be easily disentangled from all other proposals via the previous set of observables, as discussed in Section 4.
Apart from the above low-energy predictions at LHCb and Belle II, the model is also testable via high-p T searches at the LHC. The study of the 1-loop contribution of vector leptoquark U 1 exchange to B s −B s mixing revealed that the fifth vector-like lepton has to be light, around 1-2 TeV, to be compatible with the stringent bound from ∆M s . This is easily achieved in the twin Pati-Salam model, where light vector-like fermions are well motivated in order to naturally obtain the large mixing to fit the R D ( * ) anomaly, and also to fit the heavy top mass without perturbativity issues.
In particular, the fourth and fifth charged leptons are suggested as good candidates to explain the CMS excess [76], but further study is required in this direction. Vector-like quarks are found to lie not far above 1 TeV in the suggested benchmark, hence motivating specific searches at LHC to be performed. Regarding the heavy vectors, dijet searches and dilepton searches set mild bounds over the mass of the coloron and Z ′ , respectively. The more stringent bound over the scale of the model arises from the ditop searches in [20,67], which push the mass of the coloron to lie above 3.5 TeV, however those bounds could be slightly overestimated for our model as they only consider non-fermiophobic 4321 scenarios. Finally, the mass range for U 1 is compatible with current bounds, and mostly lie within the projected sensitivity of the high luminosity phase of LHC. A good fit for the 3σ CMS excess in U 1 searches [75] could be provided if the extra decay channels to vector-like fermions, assuming that the bound from ditop searches is indeed overestimated, motivating a future dedicated collider analysis.
The model proposed here features clear connections between the SM fermion masses and the leptoquark couplings which can address anomalies in LFU ratios, along with Pati-Salam universality of most of the parameters, leading to a very predictive and testable framework. The masses and mixings of first family fermion can arise along the lines of the original twin Pati-Salam model [21], involving a new family of vector-like fermions which does not couple to U 1 . Thanks to the GIMlike mechanism implemented, the phenomenology discussed in this paper would remain unaltered. However, the GIM-like mechanism introduced is accidental, and enforcing this mechanism via extra symmetries would involve the first family. Therefore, we leave the origin of first family masses, neutrino masses and the symmetry behind the GIM-like mechanism for a future manuscript, as they are uncorrelated to B-physics and the phenomenology discussed here.

B.2 Extended model
For the couplings to heavy gauge bosons we obtain,    where the up-quark couplings above receive small corrections due to 2-3 mixing arising after diagonalising the effective mass matrices in Eqs. (3.15), (3.16). However, larger 2-3 charged lepton mixing is possible (see Section 3.2), obtaining for charged leptons: at first order in s e 23 and taking c e 23 ≈ 1. The flavour-violating couplings above can lead to interesting signals in LFV processes such as τ → 3µ and τ → µγ, see more in Section 3.4.5. ..

C Benchmarks
In Table 5 we include the benchmark points considered in Section 3.4.

D.1 4-fermion operators in the SMEFT
Here we include the set of 4-fermion operators obtained at tree-level after integrating out the heavy U 1 , Z ′ and g ′ in the extended model of Section 3, and Q ll are defined as in the so-called Warsaw basis [80] of dim-6 operators built out of SM fields. In our model, the Wilson coefficients are given by where we have defined We consider all the fields in (D.1) to be mass eigenstates, as the effects of fermion mixing are encoded into the U 1 (β iα ), g ′ (κ ij ) and Z ′ (ξ ij ,ξ αβ ) couplings given in Eqs. (B.5), (B.6) (B.7) and (B.8).

D.2 b → cτ ν
The charged-current transition b → cτ ν τ is described in our model by the effective Lagrangian, where we have omitted all operators including right-handed fermions, as they receive zero or negligible contributions in our model. The matching with the SMEFT is where the negligible RGE effect is encoded as η ντ V ≈ 1.00144 and has been computed with DsixTools 2.1 [81] for Λ = 1 TeV. The Wilson coefficient [C * νedu ] 3332 (m b ) can provide a very good fit to R D ( * ) , here we take the 1σ interval from [12] (where B(B c → τν) < 30% was also imposed), The effective Lagrangian describing a generic b → sℓℓ transition reads 10) we are interested in the matching to the SMEFT of the left-handed operator, where the RGE is encoded as η ℓℓ V ≈ 0.974 and has been computed with DsixTools 2.1 [81] for Λ = 1 TeV. The expressions of the LFU ratios R K and R K * in terms of the Wilson coefficients (D.14)  where we use the numerical input C SM 10 = −4.17 [83], B B s → τ + τ − SM = (7.73 ± 0.49) · 10 −7 [84], f Bs = 230.3 ± 1.3 MeV [85], m Bs = 5366.92 ± 0.10 MeV [32], τ Bs = 1.515 ± 0.005 ps [32], α EM = 1/137.036 and G F = 1.166 · 10 −5 GeV −2 .

D.4 Bounds from B s −B s mixing
We describe B s −B s mixing with the effective Lagrangian where in the simplified model the Wilson coefficient receives tree-level NP contributions from the Z ′ and g ′ gauge bosons, written in the phenomenological limit of interest g 4 ≫ g 3,1 . Here it is clear that the coloron contribution dominates over the Z ′ one. Even in the motivated scenario ⟨ϕ 3 ⟩ ≫ ⟨ϕ 1 ⟩, where the coloron can be twice heavier than the Z ′ , the coloron contribution is at least four times larger than the Z ′ one. Such a NP contribution is constrained by the results of the mass difference ∆M s of neutral B s mesons. The experimental value is known very precisely, see for example the most recent HFLAV average [3], which is dominated by the updated measurement by LHCb [86]. However, the SM prediction historically suffered from larger uncertainties, and we need a precise knowledge of the SM contribution in order to quantify the impact of possible contributions from new physics. The theoretical determination of ∆M s is limited by our understanding of non-perturbative matrix elements of dimension six operators. The matrix elements can be determined with lattice simulations or sum rules. As discussed in Ref. [38], the 2019 FLAG average [87] is dominated by the lattice results [88][89][90], and suffers from an uncertainty just below 10% with the central value being 1 was computed in [38] as a weighted average of both the FLAG'19 average [87] and sum rule results [91][92][93]. The weighted average shows better agreement with the experiment, and a reduction of the total uncertainty (see the further discussion in [38]). The Average'19 result for ∆M s leaves where in the second step we have assumed real and positive Wilson coefficients. The SM contribution reads with S 0 (x t ) = 2.37 [94]. This way, we obtain the numerical bound where (D.29) The loop function is dominated by the heavy vector-like partners. In particular, in the motivated scenario with maximal s L 34 , the couplings with the fourth family β * sE 4 β bE 4 are suppressed by the small cosine c L 34 . This way, the loop is dominated by E 5 in good approximation. We obtain the effective loop function in this scenario by removing all constants in x α,β , which vanish due to the property (3.40), (D.30)

D.5.1 Dipole operators for τ → µγ
The process τ → µγ is described by the dipole operator, which receives contributions via both U 1 and Z ′ , The running of the dipole operator from Λ = 2 TeV to the scale µ ∼ m τ is given by C NP µτ (m τ ) ≈ 0.92C NP µτ (Λ), as estimated with DsixTools 2.1 [81]. Neglecting the muon mass, the branching ratio is given by (D.42)

D.7 b → sνν
We define the relevant Lagrangian to describe b → sνν transitions as The universal SM contribution reads where X t = 1.48 ± 0.01 [96], and α W = g 2 2 /(4π) with g 2 ≃ 0.65 being the SU (2) L coupling. We further split the NP effects into Z ′ -mediated and U 1 -mediated contributions as follows, and we only obtain sizable contributions in the τ τ channel, The U 1 contribution at NLO accuracy reads [61] C τ τ ν,U ≈ C RGE ν,U + where the second term arises from the semileptonic box diagram in Fig. D.3a, and the first term encodes the RGE-induced contribution from the tree-level leptoquark-mediated operator (s L γ µ b L ) (τ L γ µ τ L ), computed in DsixTools 2.1 [81] as C RGE ν,U = 0.047 (D. 52) In the twin PS framework, we expect extra radial modes associated toφ 3,1 andφ ′ 3,1 , however they only couple to right-handed SM fermions and hence they cannot contribute to the effective operator in Eq. (D. 46).

E From CP-conjugated notation to left-right notation
Under the SM symmetry gauge group, the VL families decompose into fermions with the usual SM quantum numbers of the chiral quarks and leptons, but including partners in conjugate representations, ψ a → (Q a , L a ) ≡ (Q La , L La ) , ψ a → Q a , L where a = 4, 5, 6. In the equations above we show the equivalence between the CP-conjugated notation and the left (L) and right (R) convention, by using a CP transformation where applicable. We use the tilde notation to highlight the partners in conjugate representations. Similarly, we can write the three chiral families of quarks and leptons in L, R convention as, F ϵ dilution of the first family U 1 coupling In Eq. (3.27) we introduced a parameter ϵ which parameterises a possible suppression of the first family U 1 coupling via mixing with vector-like fermions. This idea of suppressing leptoquark couplings via mixing with VL fermions is common in the bibliography, as similar ideas are applied in [97,98] for the same purpose, and also to suppress right-handed couplings in models were the third family is charged under the low-scale SU (4) [19,20,37]. The origin of the first family U 1 coupling β de is mixing between the sixth VL fermion family and the first chiral family, i.e.
However, prior to this mixing, the sixth VL fermion is allowed to mix with another VL fermion family, following a mechanism similar to the one that originated the Cabbibo-like matrix W LQ . Let us assume an extra sixth-primed VL family transforming in the same way as the sixth family under the twin Pati-Salam symmetry, but discriminated by a flavour symmetry which we assume as Z 2 for simplicity, which forbids mixing between the sixth-primed family and any chiral family. Instead, mixing between the sixth and sixth-primed fermion families is allowed via a twin Pati 1.5, 2.5 cos θ 6 0.045 Table 6: (Left) Charge assignments under Z 2 that allow the desired mixing. (Right) Benchmark parameters which lead to a dilution ϵ < 0.1.
If we define We can achieve values of cos θ 6 smaller than 0.1 without any aggressive tuning of the parameters, obtaining the mild suppression desired for K L → µe as per Fig. 3.6a. A suitable benchmark can be found in Table 6. Interestingly, this mechanism does not affect the Z ′ and g ′ interactions, as the unitary matrices V Q 66 ′ and V L 66 ′ cancel in neutral currents. This allows the GIM-like suppression of 1-2 FCNCs to remain in place for both the quark and lepton sector via s Q 16 = s Q 25 and s L 16 = s L 25 , without entering in conflict with K L → µe nor with B-physics.