Comparatively light extra Higgs states as signature of SUSY SO(10) GUTs with 3rd family Yukawa unification

We study 3rd family Yukawa unification in the context of supersymmetric (SUSY) SO (10) GUTs and SO(10)-motivated boundary conditions for the SUSY-breaking soft terms. We consider μ < 0 such that the SUSY loop-threshold effects enable a good fit to all third family masses of the charged Standard Model (SM) fermions. We find that fitting the third family masses together with the mass of the SM-like Higgs particle, the scenario predicts the masses of the superpartner particles and of the extra Higgs states of the MSSM: while the sparticles are predicted to be comparatively heavy (above the present LHC bound but within reach of future colliders), the spectrum has the characteristic feature that the lightest new particles are the extra MSSM Higgses. We show that this effect is rather robust with respect to many deformations of the GUT boundary conditions, but turns out to be sensitive to the exactness of top-bottom Yukawa unification. Nevertheless, with moderate deviations of a few percent from exact top-bottom Yukawa unification (stemming e.g. from GUT-threshold corrections or higher-dimensional operators), the scenario still predicts extra MSSM Higgs particles with masses not much above 1.5 TeV, which could be tested e.g. by future LHC searches for ditau decays H0/A0→ ττ . Finding the extra MSSM Higges before the other new MSSM particles could thus be a smoking gun for a Yukawa unified SO(10) GUT.


Introduction
Grand Unified Theories (GUTs) [1][2][3] present an attractive setup for Physics Beyond the Standard Model (BSM). While gauge coupling unification in GUT is necessary for consistency, the unification of Yukawa couplings is optional, depending on the GUT operators generating the Yukawa interactions. Conversely, barring a numerical accident, Yukawa unification at high energies might indicate a bigger gauge symmetry.
The most convenient setup for Yukawa unification are supersymmetric (SUSY) GUT models; while supersymmetry helps with gauge coupling unification by modifying the renormalization group (RG) slopes, it can also help with Yukawa unification indirectly via loopthreshold corrections at the SUSY scale M SUSY [4][5][6][7].
The simplest example of some Yukawa couplings unifying would be b-τ unification in the 3rd family within the context of SU(5) GUTs [8]. An even more restrictive and predictive setup is that of t-b-τ (-ν) unification, which is most straightforwardly achieved in SO (10), where all SM fermions of one family, with an addition of a right-handed neutrino, constitute a single irreducible representation 16 of SO (10). In such a setup, the neutrino 3rd family coupling also has the same value as the top, bottom and tau Yukawa coupling, coming from the operator 16 3 ·16 3 ·10, where 16 3 contains the entire Standard Model (SM) 3rd family and the Minimal Supersymmetric SM (MSSM) Higgs doublets are contained in the representation 10. Henceforth, we shall refer to this scenario simply as t-b-τ unification and omit the ν, despite its coupling also unifying.

JHEP06(2020)014
In this work we study t-b-τ unification and assume its origin to be in a SUSY SO (10) GUT. Below the GUT scale, we take the effective theory to be a softly broken MSSM. In such a framework, GUT symmetry would impose relations between the soft breaking terms of the MSSM at the GUT scale. The attractive phenomenological feature of such a setup is that Yukawa unification with GUT-like boundary conditions for the soft terms results potentially in a predictive sparticle spectrum.
In the most direct "vanilla" approach, SO(10) symmetry would result in all the sfermion mass parameters to unify in a single value m 16 , the soft Higgs masses to unify in m 10 , universal gaugino masses M 1/2 , and a universal factor a 0 for the proportionality between the Yukawa and A-matrices. The only other SUSY parameters in the theory would then be the ratio of the MSSM Higgs vacuum expectation values (VEVs) tan β, and the sign of the coupling µ of the term H u · H d in the superpotential. It is known that for t-b-τ unification tan β has to be large (∼ 50) due to the top-bottom mass hierarchy m t m b . Recall that with no SUSY threshold corrections, the Yukawa coupling ratio y τ /y b tends to run via renormalization group equations (RGEs) to a GUT value of 1.3 (see e.g. [9]), and µ < 0 gives the correct sign in the threshold correction of y b to help lower this ratio to 1, see e.g. [10,11]. For this reason we consider µ < 0 to be the better motivated setup for t-b-τ unification. Interestingly enough, fits to low energy data within this specific setup, at least to our knowledge, have not really been attempted, mostly due to the region being disfavored by RGE estimates showing no electroweak symmetry breaking (EWSB), to be discussed later. In this paper, we investigate this "vanilla" region and find it viable from the point of view of EWSB. Furthermore, we obtain good fits to the low energy Yukawa data and the SM Higgs mass, resulting in a predictive sparticle spectrum. The most striking feature of the entire setup is the prediction of a typically ∼ TeV mass for the additional neutral and charged Higgses in the MSSM, a prediction which is now being tested by the LHC. The extra Higgs prediction is very sensitive especially to top-bottom unification, and is very hard to observe with a bottom-up approach, especially if the mass m A 0 is assumed a priori as in some studies, e.g. [12].
To be more specific in what our setup achieves, and to put our results in context, it is necessary to survey the existing extensive literature on the topic of t-b-τ unification. Many early studies [13][14][15][16][17][18][19][20][21] predate the Higgs mass measurement in 2012, or even the top quark mass measurement in 1995. Beside considering the viability of Yukawa unification, they also had to contend with predicting the top quark or Higgs mass, e.g. [6,21,22], or were considering naturalness based criteria [23].
In the literature, a number of important issues have been identified: 1. The µ term: µ > 0 or µ < 0?
The Higgs connecting coupling µ from the superpotential is present in the potential V of the Lagrangian only via |µ|. Assuming no additional CP violation, µ ∈ R, so the choice of the sign of µ is free.

JHEP06(2020)014
The main preference for µ > 0 in the literature stems from considerations of the anomalous magnetic moment of the muon g µ − 2, see e.g. motivation in [31]. This was measured to be above what the Standard Model predicts (see e.g. PDG [48]), and µ > 0 would provide a SUSY contribution in the positive direction, potentially explaining the discrepancy. Despite this there are indications that a fit of g µ − 2 for µ > 0 with universal gaugino masses is difficult to achieve in SO(10) [30].
The study of µ > 0 scenarios, typically within parametrization as close as possible to the constrained MSSM (CMSSM, a.k.a. mSUGRA) with universal gaugino masses, furthermore showed that there is a preferred "funnel" region for the soft MSSM parameters [28], and that the universal gaugino mass parameter should be quite small: M 1/2 500 GeV [27]. Consequently, these scenarios prefer a light gluino mg 450 TeV [34] and suggest an upper bound on the attainable gluino mass of around mg < 2 TeV [30], a constraint coming from fitting the SM Higgs mass. Due to the non-observation of such low gluino mass scenarios at the LHC, the possibility of increasing its mass was investigated in subsequent works: it was found in [26] that the gluino mass can be raised to 2-3 TeV by relaxing the Yukawa unification to be approximate at a few % level, or to introduce a split in the squark mass parameters [33]. Note that all these results are specific to the preferred soft parameter region for µ > 0.
From the point of view of a fit to the data, however, it was already realized a long time ago that µ < 0 is preferred, see e.g. [38,45], since it gives the correct sign to the threshold corrections to the y b Yukawa coupling. Since the sign of the contribution to g µ − 2 depends on Sign(µM 2 ), see e.g. [49], this prompted a consideration of nonuniversal gaugino masses, see [37][38][39][41][42][43], with M 2 < 0. Such boundary conditions can most conveniently be achieved by considering Yukawa unification within the context of the Pati-Salam symmetry instead of fully unified SO (10), see [13,[35][36][37][38]50] for various Pati-Salam setups and studies of Yukawa unification. Another possible approach to g µ − 2 with µ < 0 is to only demand that the g µ − 2 prediction is no worse than in the Standard Model, see [40]. This last case still considered non-universal gaugino masses due to EWSB considerations, see next point.

EWSB and the split between m 2
H d and m 2 Hu at M GUT . Another issue in Yukawa unification models important for their consistency turns out to be electroweak symmetry breaking. In a softly broken MSSM, a necessary condition for EWSB is to obtain m 2 Hu < 0 at the SUSY scale. This is typically automatically achieved by RGE running from M GUT , where this parameter value is positive; the scenario where RG running triggers EWSB is referred to as radiative EWSB (REWSB). Another necessary non-tachyonicity condition, however, also requires m 2 H d > m 2 Hu at the SUSY scale. Assuming the equality m 2 H d = m 2 Hu at the GUT scale, the m 2 H d is driven down faster than the m 2 Hu essentially due to the former having positive contributions to its beta function from both y b and y τ , while the latter has only contributions from y t (and potentially from y ν ), cf. [47].

JHEP06(2020)014
For this reason, most models in the literature introduce a split m 2 H d > m 2 H d already at the GUT scale [27-33, 37-45, 51]. The simplest way to achieve this is by imposing the split ad hoc, which is called "just so" Higgs splitting and assumes m 2 H d,u = m 2 0 ± ∆ at the GUT scale, e.g. [27,29], with the relative split amounting to ∼ 13 %. An alternative mechanism to generate this split is by D-term splitting [17,28,[41][42][43][44][45], which also splits up the other soft scalar masses in a particular way due to D-term contributions to the masses. Attempts to avoid m 2 H d slipping below the value of m 2 Hu have also been studied in the context of adding right-handed neutrinos or introducing a first/third scalar mass split in the GUT boundary conditions, see [51], both options essentially modifying the RGE beta functions for m 2 H d and m 2 Hu . The well known issue regarding REWSB with m 2 H d = m 2 Hu at the GUT scale has been studied in [6,17,52,53], and reiterated later in e.g. [31] based on an approximate expression for m 2 H d − m 2 Hu at low scales taken from [54]. It should be noted, however, that these papers use semi-analytic formulas for RGEs running from the GUT scale to the SUSY scale, which hold only approximately. In the context of the GUT boundary condition m 2 H d = m 2 Hu , successful REWSB was achieved for the case of non-universal gaugino masses [47,55], while the old arguments for the universal gaugino mass case are reiterated. On the other hand, successful REWSB was found for the case of CMSSM with µ < 0 in [46], albeit with only approximate Yukawa unification due to their bottom-up approach of running Yukawa parameters.
In contrast to most considerations in past works presented above, we find that exact Yukawa unification with universal gaugino mass terms and m 2 Hu = m 2 H d is in fact possible. We show this explicitly by performing the RGE running numerically; although we use 2-loop RGEs for the MSSM + soft terms for (most) results, the 1-loop RGE solutions already confirm this qualitative picture. While we agree with prior analyses that RG running just below the GUT scale causes m 2 Hu > m 2 H d in the running parameters, this relation reverses later by RG running a few orders of magnitude above the SUSY scale, thus achieving successful REWSB. This holds true at least in a large part of the soft parameter space. Crucially, however, the running value of m 2 H d − m 2 Hu is typically below (1 TeV) 2 at the SUSY scale, causing the extra MSSM Higgs bosons to be the lightest part of the sparticle spectrum.

Experimental constraints and considerations.
The most obvious type of prediction studied in Yukawa unification models is the MSSM spectroscopy, see [24,36,44,56,57] for studies which focus on this.
Studies which fit GUT models to the experimental data usually consider some or all of these constraints. It was found in many specific realizations of Yukawa unification, however, that potential experimental tensions can usually be relieved by relaxing the demand for exact Yukawa unification and impose it only at a level of some %. This essentially works due to relaxing constraints on the superpartner masses. Such scenarios have been dubbed "quasi-unification", see e.g. [25, 35-37, 50, 55, 61, 62]. Alternative setups to improve fits have also been tried, such as splitting the Aterms [63], considering 4 Higgs doublets instead of 2 [64], introducing certain extra vector-like fermions motivated by an E 6 GUT context [65], or introducing an entire vector-like family of SM fermions [66].
In this paper, as motivated earlier, we consider µ < 0 and numerically find a good solution for REWSB despite the relation m 2 H d = m 2 Hu and universal gaugino masses. In the literature, as far as we are aware, the only case directly comparable with ours is in [46], with the limitation that the SM Higgs mass was not yet measured at the time. One of the scenarios they consider successful (including EWSB) is the CMSSM (implying universal gaugino masses and no GUT split between m 2 H d and m 2 Hu ) with µ < 0. They use, however, a bottom-up approach for Yukawa RGE, and therefore consider only the quasi-unification scenario with a parameters scan. They consequently do not find the low MSSM Higgs mass effect, since it is very sensitive to exact unification, as we show in this paper.
Given the effect of the low extra Higgses we study in this paper, the most acute experimental constraints would come from two possible sources. The first is the B s → µ + µ − decay, with the extra Higgs contribution estimated as, see e.g. [67], compared to the PDG measured value of (3.2 ± 0.7) · 10 −9 [48]. The second constraint is the increasingly competitive LHC searches for ditau decays H 0 /A 0 → τ + τ − of the neutral MSSM Higgses, see [68,69], with current bounds implying m A 1.5 TeV (for tan β = 50). Given this most recent estimate and future trends of bounds, we find the ditau search to be comparable or more stringent than the B s → µ + µ − process; we thus focus only on the ditau decay in this paper for simplicity. The other parts of the SUSY spectrum in our setup are heavy, larger than 4 TeV for gluinos and squarks, far above the present ATLAS and CMS bounds but within reach of future colliders such as the FCC-hh or SppC. The organization of the paper is as follows: in section 2 we introduce our notation and conventions, and analyze the salient points regarding EWSB and the masses of the extra Higgs bosons in the MSSM. In section 3, we perform an RGE analysis of the quantity m 2 H d − m 2 Hu relevant for both those aspects and perform a sensitivity analysis to deformations of various parameter relations around an example point. In section 4, we perform a more general investigation of the CMSSM parameter space and show that the masses of the extra Higgses are predicted to be low in general. Finally, in section 5, we analyze how JHEP06(2020)014 constraints from the LHC challenge exact Yukawa unification and how a quasi-unification scenario helps in this regard. Then we conclude. For completeness, we also include two appendices. In appendix A the general 1-loop RGEs for a softly broken MSSM with righthanded neutrinos are presented. In appendix B a simplified version of the RGEs neglecting the Yukawa couplings of the first 2 families is given.
The indices i and j are family indices, the SU(2) contractions between doublets are denoted by a dot and defined by Φ · Ψ ≡ ab Φ a Ψ b with 12 = − 21 = 1, while the SU(3) indices are suppressed. Also note that a left-chiral superfield Φ c contains the charge conjugated fermion field ψ † , as well as the conjugated complex scalar fieldφ * R . The soft-breaking terms consist of gaugino mass terms, the scalar trilinear A-terms, the scalar soft-mass terms, and the b-term:

JHEP06(2020)014
We labeled the SU(3) C , SU(2) L and U(1) Y gauginos by λ a 3 , λ b 2 and λ 1 , respectively. The tildes above the fields indicate the scalar component of the superfield, with the exception of H u and H d , which also indicate scalar parts.
The neutral components of H u and H d each acquire an EW breaking VEV: which -motivated by EW symmetry breaking in the SM -are parametrized by This leaves tan β as the only free parameter, and v u , v d ∈ R.
Minimization of the potential with respect to the electrically neutral components H 0 u and H 0 d of SU(2) doublets leads to a (tree-level) vaccum solution Note that we have solved the vacuum equations for the superpotential parameter |µ| 2 and the soft parameter b, while treating the unknown VEVs v u and v d as independent variables, appearing implicitly via v u /v d = tan β. In the large tan β regime, we can make the approximation implying that a solution to EWSB (at tree level) is possible only if the soft mass parameter is negative at the energy scale of computation, i.e. m 2 Hu < 0 at the SUSY scale. After EW symmetry breaking, 3 real scalar degrees of freedom in H u and H d become part of the longitudinal components of the massive gauge bosons W ± and Z 0 via the Higgs mechanism, leaving 5 real degrees of freedom to be physical. We label them in the standard way by h 0 , H 0 , A 0 , H + and H − , where their superscripts denote their EM charge. The low mass Higgs at 125 GeV is denoted by h 0 , while H 0 and A 0 denote heavier neutral scalars with even and odd parity P , respectively. We get the following well-known expressions for their tree-level masses: (2.14) Considering the regime m 2 A 0 m 2 Z , m 2 W leads in leading order to showing that all extra Higgs particles H 0 , A 0 and H ± are near the scale m 2 A 0 . The scale of m 2 A 0 in turn depends on the vacuum solution for |µ 2 |; combining eq. (2.11) and (2.8) gives the tree level value We see that, crucially, the scale m 2 A 0 depends on the difference m 2 H d − m 2 Hu of the masssquare soft parameters. In the large tan β regime, this approximates to so that a non-tachyonic tree-level mass for A 0 requires m 2 Hu as a necessary condition. We now briefly turn to a discussion of the scale of masses at 1-loop level. The vacuum solutions at 1-loop become (see [71,72]) The hatted quantities, includingm 2 W for later convenience, are defined bŷ where t u and t d are 1-loop tadpole expressions, and Π T ZZ and Π T W W are the transverse Z and W -boson 1-loop self-energies. The hatted massesm 2 Z andm 2 W are the 1-loop masses computed in the DR renormalization scheme. Their explicit expressions can be found in [71] and will not be reproduced here. For a consistent loop calculation, the quantities in the expressions for 1-loop corrections can be taken to be the parameters at tree-level.
When the quantities in the superpotential of eq. (2.3) are complex, the neutral states h 0 , H 0 and A 0 mix: with the 1-loop correction, the masses may no longer be CP eigenstates.

JHEP06(2020)014
We shall not be considering complex phases in the SUSY parameters, so this complication need not be considered.
Due to the breaking of CP symmetry at next to leading order in the general case, rather than the mass m 2 A 0 ,tree from (2.16), a more convenient quantity to consider is the mass of the charged Higgses H ± , since the charged Higgses H ± have no other states to mix with. The expression at 1-loop order for the mass of H + is known to be with Π H + H − denoting the self-energy of H ± , see [71].

RGE analysis of m
As a first step in assessing models with Yukawa unification and SO(10) boundary conditions for soft parameters, we study the RG running of the quantity m 2 H d − m 2 Hu . This quantity must be positive at the SUSY scale, a feature crucial for EWSB, and its magnitude sets the mass scale of the extra MSSM Higgs states H 0 , A 0 and H ± , as was discussed in section 2. An often cited requirement in the literature for REWSB to occur is a split in the GUT scale boundary conditions for m 2 H d and m 2 Hu , see section 1 and references therein. We show here, however, that such a split is not necessary, since we obtain m 2 H d − m 2 Hu > 0 at the SUSY scale regardless. The value of this difference, however, is small compared to the magnitude of each term, implying low lying extra Higgs states in the MSSM, an effect that we show to be especially sensitive to t-b unification.
To facilitate the RGE analysis, we make use of simplified RGEs at 1-loop and CMSSM boundary conditions, as explained in separate subsections below. Note that these simplifications are specific to this section of the paper and do not change the general conclusions, confirmed by comprehensive analyses in later sections by use of 2-loop RGEs and SO(10) motivated boundary conditions. The analysis of the simplified case nevertheless gives valuable insights into EWSB and the low spectrum of the extra MSSM Higgses, confirming that this striking feature can be understood as an RGE effect, and is seen already at 1-loop order.

The simplified boundary conditions -CMSSM
In this section we make a slight simplification and consider the CMSSM boundary conditions (see e.g. [73]) as the default scenario, instead of the SO(10) motivated split in the sfermion and Higgs soft masses to be studied later. We also study how RG running changes under various deformations of the default CMSSM boundary conditions, obtaining a number of important conclusions applicable to the more general scenario beyond CMSSM.
More explicitly, we assume the following for the RGE analysis in this section: • The boundary conditions are set at a high energy: M GUT = 2 · 10 16 GeV.
• The MSSM is extended by right-handed neutrinos at a scale M R , with M R ≤ M GUT , below which they are integrated out.
• The boundary conditions of the soft parameters are those of CMSSM: The RGE boundary conditions for the soft parameters are thus parametrized by the 3 CMSSM parameters m 2 0 , M 1/2 and a 0 .
• Unification of 3rd family Yukawa couplings at the scale M GUT : The above assumptions are a simplified version of the "SO(10) boundary conditions" with only one soft scalar mass parameter m 0 and with universal sfermion soft matrices (typical leading order pattern in "flavored GUTs" with family symmetry): the constraints are implied in the unification of all fermion sectors, and t-b-τ unification arises in the simple case when the Yukawa contribution to the 3rd family of 16 F comes from the 16 F 3 ·16 F 3 ·10 H operator in SO (10). We note that although the stated class of SO(10) models gives rise to the MSSM setup described below M GUT , we do not necessarily commit to a particular SO(10) UV completion. In this context, we would also like to remark that the exact Yukawa unification will be subject to model-dependent corrections such as e.g. GUT threshold corrections, which however depend on the details of the UV completion. We will study the effects of such perturbations of the scenario later in the paper.

The simplified 1-loop RGE
The complete set of RGEs for the neutrino-extended and softly-broken MSSM are given in appendix A (also cf. [71]). The full RGEs can be simplified by eliminating some degrees of freedom which are either numerically irrelevant or unnecessary for our considerations. In the quark sector, for example, there is little mixing, and the Yukawa matrices in both quark JHEP06(2020)014 sectors as well as the charged lepton sector have hierarchical masses. A good approximation is therefore to consider only the 3rd family of fermions. Also, we assume family universality in all sfermion mass matrices at the GUT scale.
To simplify the RGE, we consider the minimal amount of variables consistent with the above assumptions. It turns out that the following 28 variables in the RGEs are required: • The 3 gauge couplings g 1 , g 2 and g 3 .
• 4 Yukawa couplings of the 3rd family y t , y b , y τ , y ν .
• The 6 × 2 + 2 soft mass parameters: m 2 x i , where x ∈ {Q, L, u, d, e, ν} and i ∈ {1, 3} are independent, and the Higgs mass parameters m 2 H d and m 2 Hu . The case i = 2 does not have to be studied separately since, in our setup, the i = 2 quantities have exactly the same running and boundary conditions as those for i = 1.
The resulting simplified 1-loop RGE are presented in appendix B, which contains also more details on the above variables, cf. eq. (B.1)-(B.5). Making use of the RGEs from appendix B, the running of the expression m 2 Hu is then determined to be where c 1 is the loop factor and S is a linear combination of soft masses: We see that the first 4 terms of the result in eq. (3.6) are analogous to each other, the quantities in the terms correspond respectively to the particles b, t, τ and ν τ (and their superpartners). Each term contains the modulus-squared of its Yukawa coupling, and the factor next to it contains a modulus-squared of the appropriate A-term factor, as well 3 more terms with the soft masses of particles present in the corresponding superpotential Yukawa term. The b and t terms have an additional numerical factor 3 compared to τ and ν due to the 3 possible SU(3) colors they can take. Crucially, the terms also come into the RG beta function with different signs, so it may happen that they cancel. Below the right-handed neutrino mass scale M R , the ν term vanishes. The boundary conditions JHEP06(2020)014 imply that at exactly M GUT , the last term vanishes due to S = 0, and the b and t terms cancel each other, and as well as the τ and ν terms, such that we have As already stated, the scale of the masses of the extra MSSM Higgs bosons will be determined by This same quantity must be positive at low energies also for successful EWSB. It is computed numerically by solving the RGE differential equations of appendix B. We shall often allude to eq. (3.6) for a better understanding of the numerical results, which we now consider.

Numerical RGE results
We now investigate the RGE properties of the system numerically. To do this as explicitly as possible, we take an example parameter point, whose neighborhood we study. We stress that the conclusions of the RGE behavior in this section nevertheless hold generally, i.e. different example points of Yukawa unification at high energies and consistent with experimental data at low energies yield the same qualitative conclusions, which we checked explicitly by considering different parameter points. Furthermore, we identify the underlying reasons for certain RG behaviors throughout this section, and the generality (where applicable) is also confirmed by results in later sections. We take the following boundary values for the parameters at the scale M GUT = 2.0 · 10 16 GeV: The gauge coupling g 1 is given in the GUT normalization, and M R is the mass of the added right-handed neutrino. The above values are to be understood as boundary conditions for the RGE in appendix B. At the scale M R , the right-handed neutrino is integrated out;

JHEP06(2020)014
below this threshold, the RGE are corrected by removing all terms containing y ν . For the example point under consideration, we have taken M R = M GUT so that by default no effects arise due to the right-handed neutrinos, since the y ν term with the large 3rd family neutrino Yukawa coupling is removed already at the GUT scale; its effect is studied separately below. The values of the gauge couplings at the GUT scale are taken from the high-energy data provided by [9], which uses 2-loop RGEs and takes the SUSY scale at 3 TeV; note that their values are consistent with a typical unified gauge coupling value of ≈ 0.7.
The overall scale of the soft parameters m 0 , M 1/2 and a 0 has been taken at the order of a few TeV, which tends to be the preferred scale for the fits to low energy data, as will be seen in the next sections. Also, the main effect we are after in this paper is that the extra MSSM Higgs particles are unexpectedly light compared to the SUSY scale, for example 1 TeV; this effect will be obscured if the SUSY scale is also taken to be lighter than 1 TeV, as used to be popular in past SUSY studies. The few TeV scale for sparticles is compatible with (as of yet) non-observation of SUSY particles at the LHC.
Note that the chosen point is such that it gives the correct 3rd generation Yukawa couplings y t , y b and y τ at the scale M Z in the MS scheme, based on the data from [9]. An intuitive qualitative description of how the GUT scale parameters control the fit of the 3rd family Yukawa parameters is the following: • The value y 0 controls the overall scale of the 3 Yukawa couplings, and needs to have the value y 0 ≈ 0.5.
• The effect of the soft parameters m 2 0 , M 1/2 and a 0 is to control the SUSY spectrum, through which SUSY threshold effects give the correct ratio y τ /y b .
• The quantity tan β controls for the ratio y t /y b (alongside SUSY threshold corrections).
Low energy data demands a large value of tan β ≈ 50, a well-known feature of MSSM based t-b-τ unification models.
We plot the running under 1-loop RGE from appendix B for the various quantities of the MSSM, with the boundary conditions at M GUT given by the example parameter point in eq. (3.11)- (3.20). We shall also investigate the effect of changing one feature of the boundary conditions at a time, understanding its impact; note that we do not evaluate the worsening of the fit to low energy data under such a deformation, since we are for now interested only in the (numerical) effect on the RGE running. We plot quantities in the range [M SUSY , M GUT ]; note that the lower scale is the SUSY scale, since that is the scale where the sparticle spectrum is computed. This scale is also where a match between the SM and MSSM theories is performed, and it is taken to be the geometric mean of the masses of the two stops (computed for our example point using SusyTC [71] to be JHEP06(2020)014 M SUSY = 5901 GeV). While we used a custom computer code for RGE running based on appendix B for greater control, the results were compared and confirmed with SusyTC when applicable.
The RGE running of the system, based on the results of the example point, turns out to have the following properties: 1. Running of gauge and Yukawa couplings, gaugino masses and the A-terms.
The RGE running of the gauge couplings, Yukawa couplings, gaugino mass parameters, as well as the the A-term factors a x from eq. (B.3) is shown in figure 1.
As always in the MSSM, each of the gauge couplings evolves independently from other quantities (at 1-loop level); the couplings approximately meet at ∼ 0.7, and their running values are determined; when the renormalization scale µ r decreases to low energies, g 3 runs upwards and g 1 and g 2 run downwards, see eq. (B.6), due to the signs of MSSM beta coefficients β 3 < 0 and β 1 , β 2 > 0 from eq. (A.24).
The running of gaugino mass parameters, according to eq. (B.7), is influenced by the gauge couplings. It is the differences in gauge couplings which drive the gaugino mass-parameter differences from a common boundary point M 1/2 at M GUT . This explains why the gluino mass parameter M 3 increases when approaching M SUSY , while M 1 and M 2 decrease, but all are at a scale of 2 TeV or higher.
The RGEs of the Yukawas have two competing contributions to the beta functions, cf. (B.9)-(B.11): a positive contribution from the Yukawas themselves, and a negative contribution from gauge bosons (terms proportional to g 2 i ). The Yukawa couplings can then rise or fall with smaller µ r , depending on whether the gauge or Yukawa contributions to the beta function are dominant, respectively.
The 3rd family Yukawa parameters y t and y b rise with lower scale µ r essentially due to the relatively large negative g 2 3 term from the gluons, while y τ stays mostly flat, since realistic unified values of the gauge couplings of ≈ 0.7 and Yukawa couplings of ≈ 0.5 give the Yukawa and gauge contributions approximately equal. The difference between the top and bottom Yukawa, on the other hand, is small and is essentially driven by the |y τ | 2 term in β(y b ) and the difference in the g 2 1 terms in β(y t ) and β(y b ), see eq. (B.9) and (B.10). This ensures a small relative difference y t − y b , with y t > y b at all energies; the very different values of y t and y b at M Z , as implied by the different masses of the t and b quarks, must thus come from the MSSM to SM matching at M SUSY , implying a large tan β of around 50.
The RGEs for the A-term factors are given in eq. (B.16)-(B.18). We can see that the difference between a u and a d is essentially driven by the difference between y t and y b , as well as the |y τ | 2 and g 2 1 terms, which essentially already drive the y t and y b difference, as discussed earlier. For this reason, there is again only a small deviation between a u and a d . The slope of a e in absolute terms is smaller due to no gluino related terms, and because of smaller numerical factors in front of the Yukawa terms. JHEP06(2020)014  For m 2 Hu and m 2 H d , the positive Yukawa term contributions to the β functions dominate, leading to a positive slope and thus the parameters becoming smaller and eventually negative with smaller µ r . The drive to m 2 Hu < 0 at low µ r confirms that the EWSB is radiative. Crucially, the necessary condition for EWSB m 2 H d > m 2 Hu is also satisfied at low scales, as will be discussed in more detail later.
The soft mass parameters related to the squarks grow fast with smaller µ r due to the large negative contribution of the gluino related terms g 2 3 |M 3 | 2 . These terms are not present in the β function for soft-mass parameters of leptons, so the slepton masses stay almost flat.
Another general feature of the soft-mass parameter running is that the masses of the 1st and 2nd family of squarks and sleptons (index 1) become larger than those of the JHEP06(2020)014 3rd family (index 3); we are comparing here the soft-mass parameters of particles of the same flavor, but from different families. The simple reason is the additional positive terms proportional to squares of Yukawa couplings, which appear only for 3rd family squarks and sleptons (since the 1st and 2nd family Yukawa coupling are negligible compared to the 3rd family, and they are set to zero in our simple scenario). We thus have the usual inverted hierarchy in the squark and slepton masses.
We now discuss how the scenario of t-b unification and y t = 1.1 y b compare. We see that there is little qualitative difference for the values of any one soft parameter taken on its own. Visually though, major quantitative changes in relative terms can be spotted when comparing the quantity m 2 H d − m 2 Hu in the two scenarios, as well as changes in the quantity m 2 u 3 − m 2 d 3 . These changes might be deemed to have an insignificant effect on the low energy observables. But as shown in the previous section, the difference m 2 H d − m 2 Hu turns out to determine the mass scale of the extra MSSM Higgs bosons. That means that the exactness of t-b unification at the GUT scale, as demonstrated by the two scenarios in figure 2, has a big impact on the sparticle spectrum, i.e. on the extra Higgs sector to be precise. This is the major effect that this paper investigates.

Effect of t-b unification on m 2
H d − m 2 Hu . We have seen from the RGE of the soft masses in the previous step that t-b unification 1 has little qualitative effect on the running of these parameters taken in isolation, but has a crucial effect on m 2 H d − m 2 Hu . Figure 3 shows RGE trajectories for m 2 Hu under different y t /y b ratio boundary conditions at M GUT , essentially demonstrating the sensitivity of this quantity to t-b unification. We see that for our example point, the running expression m 2 H d − m 2 Hu increases essentially linearly with the y t − y b difference (at least when relative differences are small), and with a substantial increase already when y t and y b differ at the percent level. The impact is even more dramatic when considered in terms of relative increases of m 2 H d − m 2 Hu : a deviation of a mere 10 % from t-b unification raises the value by a factor 4, and consequently the masses of the extra MSSM Higgs particles by a factor of 2. Looking at this from a reverse perspective, when approaching t-b-τ unification from a t-b deformation direction, the predicted masses of the extra Higgses drop very quickly, typically below 1 TeV.
Hu > 0 is necessary for (tree-level) EWSB. We can see in figure 3 that this condition is fulfilled even for exact Yukawa unification (the y t = y b curve), at least for this particular example point. This shows that there exist parameter points with exact Yukawa unification and successful EWSB. It is important to note that a successful EWSB with the m 2 H d = m 2 Hu GUT boundary condition (and universal gaugino masses) was not found in some of the prior literature [6,17,31,52,53] due to extensive use of semi-analytic approximate formulas from e.g. [54], as was of the stop masses. Strictly speaking, the scale M SUSY shifts slightly with different ratios y t /y b , so that comparing the running quantity m 2 H d − m 2 Hu at a fixed scale is not exactly the same as comparing the mass scales of the extra Higgses. This shift, however, is negligible, since the quantities determining the stop masses run logarithmically with µ r and change only slightly with the ratio y t /y b , as argued in the previous analysis step. It is thus justified to compare the running expression for different curves at a fixed scale M SUSY for qualitative considerations.

Contributions to β(m 2
H d − m 2 Hu ). To understand the effect that t-b-τ unification has on the RGE running, we consider the various contributions to β(m 2 H d − m 2 Hu ). One can combine the separate RGE in eq. (B.21) and (B.22) into the β function of eq. (3.6). For our example point, where the right-handed neutrinos are already integrated out at M GUT , there are 4 terms: terms proportional to |y t | 2 , |y b | 2 and |y τ | 2 , as well as a term proportional to S, which is a linear combination of scalar soft masses, see eq. (A.25). We plot these contributions for the cases y t = y b and y t = 1.1 y b in figure 4.
The results show that in absolute terms the |y t | 2 and |y b | 2 contributions dominate over the |y τ | 2 one at M GUT , an effect which only increases when running to lower µ r , while the contribution from the S term stays numerically negligible throughout and will thus be ignored in the following discussion. The larger contributions of the t and b terms start out due to larger numeric prefactors (due to color) compared to the τ term. Furthermore, at lower energies the Yukawa couplings y t and y b rise with smaller scale, while y τ falls, see figure 1. In addition, also the soft masses show the same trend, see figure 2. Note, however, that these terms in β(m 2 H d − m 2 Hu ) come with different signs; in particular, the t and b contributions have opposite signs.

JHEP06(2020)014
It is thus convenient to compare the difference of the t and b terms (red curve) with the τ contribution (green curve), see right panels of figure 4. We shall refer to these two contributions as the t-b and τ contributions, respectively. The t-b contribution comes into the β function with a negative sign, so whenever the red curve dominates over the green curve, the beta function value becomes negative, i.e. the RGE running of m 2 H d − m 2 Hu has a negative slope. Conversely, when the τ contribution dominates and the green curve is above the red, the slope is positive. As the figure shows, the slope is positive at large µ r and negative at small µ r , which is consistent with figure 3.
At low enough µ r the t-b contribution is expected to dominate over the τ contribution regardless of the starting y t /y b ratio simply due to Yukawa coupling values at those energies, and that typically the squark soft masses are larger than the corresponding lepton ones. The ratio y t /y b is crucial, however, for the t-b contribution at energies near the GUT scale: when y t = y b the t-b contribution starts at zero, while y t /y b > 1 implies a non-vanishing starting value for the RGE. 2 This crucially impacts the scale at which the t-b contribution becomes bigger than the τ one, i.e. when the red and green curves on the right panels of figure 4 cross. We see that for y t = 1.1y b the t-b contributions already starts out almost as big as the τ contribution at M GUT , so the curves intersect above 10 14 GeV, while t-b unification delays this until below 10 11 GeV. Consequently, with t-b unification the value of m 2 H d − m 2 Hu will be much lower, since the rise in its running value is delayed by several orders of magnitude in the energy scale µ r .
This completes our understanding of the effect of t-b unification on m 2 H d − m 2 Hu . Yukawa unification delays when the t-b contribution in the beta function rises enough to dominate over the τ contribution, allowing for the running value of m 2 H d − m 2 Hu to rise much less by the scale µ r = M SUSY . We emphasize that this effect is an indirect consequence of RG running of all parameters, and can thus be seen only when solving for the entire system of RGE numerically and evolving it over multiple orders of magnitude of µ r . In simplified analyses, such as studying the local RG behavior at M GUT by Taylor expansion or taking some running quantities in the beta function as constant to derive a linear-log semi-analytic approximation [54], not even the m 2 H d > m 2 Hu property at low µ r is reproduced, let alone the more subtle effect of the t-b deformation.

Effect of
Hu . An interesting question now is what impact b-τ unification of couplings has on lowering the value m 2 H d − m 2 Hu . It turns out that while t-b unification is crucial for this effect, b-τ unification is not.
We plot the RGE flow of m 2 H d − m 2 Hu for different ∆ τ := y b − y τ = y 0 − y τ in figure 6. The results clearly show that b-τ unification has minimal effect on that quantity at the SUSY scale. The two sets of trajectories on the plot correspond to the y t = y 0 JHEP06(2020)014 2 |yτ 2 (a e 2 + m Hd 2 + m L3 2 + m e3 2 )  case (red-blue) and the y t = 1.1y 0 case (green-cyan); trajectories in the same set differ in ∆ τ from 0 to 0.2, which presents a relative drop in y τ compared to b-τ unification of more than 40 %, but trajectories in the same set nevertheless cluster together at M SUSY , despite diverging at first at intermediate energies.
Hu . We see from figure 5 that the scale of right-handed neutrino M R , associated with the large 3rd family neutrino Yukawa coupling y ν , has a comparatively small effect on the value of m 2 H d − m 2 Hu at M SUSY , relative to effect of the t-b deformation. The discontinuous changes in the slope happen at scales when the right-handed neutrino is integrated out, i.e. at the scale M R . We conclude that the right-handed neutrinos do not have a large direct effect on the mass scale of the extra Higgs particles, and we therefore do not include them in the analyses of sections 3 and 4. It should be noted though that an indirect effect turns out to be possible, since their presence shifts the region of parameter space where good fits to low energy data are obtained, see section 5. We investigate whether having a simplified set of CMSSM parameters for the soft term boundary conditions is crucial for having light extra Higges. A more realistic, yet JHEP06(2020)014 We therefore consider a slightly more general case of parametrization for the soft terms, which we refer to as "SO(10) boundary conditions". We keep the M 1/2 and a 0 parameters, but have two different soft mass parameters m 16 and m 10 for the sfermions and Higgses, respectively: (3.25) The notation for m 16 and m 10 signifies which SO(10) representation the scalars of the soft term are part of. It is presumed here that H u and H d come from a 10 of SO(10), which allows for t-b-τ unification with the simple renormalizable 3rd family Yukawa operator 16 3 · 16 3 · 10.
We investigate the effect of such an SO (10)  Hu at the GUT scale, as is common in the literature [27-33, 37-40, 44, 45, 51], can erase the effect of low extra Higgs masses.

JHEP06(2020)014
The above discussion also makes it clear that unless one studies scenarios of t-b unification, usually in the context of t-b-τ unification, this effect will be missed. In particular, this effect will not be present in any kind of SU(5) SUSY GUT model attempting merely b-τ unification.

The typical mass scales of the extra Higgs particles
In this section, we turn to the broader question of the predicted mass range of the extra MSSM Higgses when considering the entire region of parameter space that yields good fits to low energy data.
We The next step is a more precise calculation going beyond the proxy quantity m 2 H d −m 2 Hu , instead considering the masses of the extra Higgs particles directly. We make the following improvements in the analysis for estimating the Higgs masses as accurately as possible: 1. The RGE running of the softly broken MSSM is performed at 2-loop level.

The masses of the extra Higgses are computed at 1-loop instead of tree level.
To perform such improved calculations, we make use of the following tools:

JHEP06(2020)014
• The MSSM Higgs sector is computed to higher loop order by the program Feyn-Higgs [74][75][76][77][78][79][80], version 2.13.0. The output of SusyTC gives the Higgs masses at tree level, with the exception of m 2 H ± given at 1-loop by using eq. (2.21). Using the output values of SusyTC as input for FeynHiggs, the SM Higgs mass is computed to 2-loop and the extra Higgs particles' masses are computed to 1-loop.
• For the computation of EW vacuum stability we make use of Vevacious [81]. We use SusyTC to produce an SLHA file, amended with values of the MSSM µ and b terms at tree and loop level, computed from the VIN file of the tree and 1-loop potential for EW breaking produced by SARAH 4.14.1 [82,83]. We use the SARAH predefined model with possible charge breaking via stau VEVs.
We use these tools for improved computations of the t-b-τ unification model, where we still consider only the 3rd family Yukawa couplings to be non-vanishing as in section 3, and assume the right-handed neutrinos are integrated out at the GUT scale. The GUT scale values of the gauge couplings are taken to be those from eq. (3.11)-(3.13). We shall consider two scenarios of boundary conditions: the CMSSM scenario (5 parameters) and the SO (10)  We take µ < 0 in all cases. The standard notation of CMSSM parameters applies, the parameter y 0 is the t-b-τ unified Yukawa coupling, while m 16 and m 10 are defined according to eq. (3.24) and (3.25). Each parameter point in a scenario allows the computation of the Yukawa couplings at M Z , the Higgs mass, as well as the SUSY spectrum. The part of the SUSY spectrum that is of greatest interest to us is the one of the masses of the extra MSSM Higgs particles; we would like to confirm that due to t-b-τ unification they should indeed be comparatively low.
As a first check, we recompute the example point from eq. (3.11)-(3.20) with improvements of higher loop order. The results for the mass prediction of the CP-odd Higgs A 0 are the following: (4. 3) The result I corresponds to the tree level mass from eq. (2.11) and 1-loop RGE, the result II corresponds to tree level mass and 2-loop RGE, while result III is the most accurate with the 2-loop RGE and 1-loop mass from FeynHiggs. We see that the predicted mass reduced after every improvement, which we find happens generically. This confirms that the low MSSM Higgs mass phenomenon persists (and may be further enhanced) even with the improved loop order in the calculation. We now turn to a more general study of the parameter space beyond just the example point. In the subsequent analysis, the 3rd family Yukawa couplings and the SM Higgs mass JHEP06(2020)014 are considered to be observables: As a measure of goodness of fit we make use of the χ 2 function: where the vector x represents the input parameters of the model from either eq.  1 in [9], with relative errors adjusted upwards to 1 % due to limited precision of our RGE procedure from M GU T to M Z . The SM Higgs mass central value was taken to be m h = 125.09 GeV [84], with a 3 GeV error due to theoretical uncertainties in the computation.
We show that the prediction of a low extra Higgs mass is a generic feature of t-b-τ unification rather than of just the example point from the previous section. For this reason we search for a number of other points in the parameter space of CMSSM, which provide a good fit of the observables. We do this by a systematic search in the m 0 -a 0 plane of parameters. For a fixed m 0 and a 0 , we perform a minimization of the χ 2 for the other 3 input parameters M 1/2 , y 0 and tan β in eq. (4.1). Remember that these 3 free parameters are used to fit 4 observables of eq. (4.4), which may not necessarily be possible for an arbitrary point in the m 0 -a 0 plane. The computation involves a minimization of χ 2 for each point in a 25 × 37 grid and subsequent interpolation between grid points; the points were taken equidistant and in the range 100 GeV ≤ m 0 ≤ 5500 GeV, − 12000 GeV ≤ a 0 ≤ 6000 GeV, (4.6) and include the edge points of these intervals. As we shall see, this range includes the entire region of admissibly low χ 2 , at least in the CMSSM context. The relevant results of this fit are summarized in figures 9, 10 and 11. We analyze them below: • Figure 9 shows the contours of the minimal attainable χ 2 for a point in the m 0 -a 0 plane, with the shaded region excluding points due to vacuum stability, to be discussed below. Contour regions from blue to white represent points where a reasonable fit can be obtained: the darkest shade of blue represents almost perfect fits of χ 2 < 1, while the white region represents the edge points where χ 2 < 9, such that the deviation in any one observable cannot be more than 3σ. We see that the allowed region in the m 0 -a 0 plane is compact: the ranges are roughly m 0 < 4 TeV, − 12 TeV < a 0 < 5 TeV, (4.7)

JHEP06(2020)014
i.e. the regions involve scales of a few TeV. The numeric values of µ at the SUSY scale are computed to be in the interval (−12 TeV, −2 TeV) for all points.
• The darkly shaded region in figure 9 corresponds to points in the m 0 -a 0 plane for which χ 2 has been minimized, but the vacuum is not sufficiently stable. The threshold is taken to be at 10× the current age of the universe, but the exponential sensitivity of the lifetime to the bounce action (see [85][86][87]) means that one order of magnitude difference in the threshold does not appreciably change the excluded area. The unshaded region thus represents points with the EW vacuum either being metastable with a sufficiently long lifetime or stable. Note that the instability in the shaded region does not necessarily exclude all possible points with a given m 0 and a 0 , but only the one minimizing χ 2 . Although an improved approach would be to include a sufficiently long vacuum lifetime as a necessary condition in the minimization of χ 2 , this would be much more demanding computationally. Ultimately, the vacuum computation performed here is sufficient to show that most of the low χ 2 region consists of allowed points.
• The minimization of χ 2 gives the following ranges for tan β and y 0 for all best-fit points: 48 < tan β < 55, 0.44 < y 0 < 0.50. (4.8) These two parameters thus have small relative changes for best-fit points with different CMSSM soft parameters. The results are compatible with the well-known fact that t-b-τ unification requires tan β ≈ 50, while the unified coupling is approximately y 0 ≈ 0.5. A more interesting input parameter to track for different best-fit points in the m 0 -a 0 plane, however, is the gaugino mass parameter M 1/2 , since this provides the information for all CMSSM soft parameters of the well-fit points. A contour plot of the M 1/2 values is presented in figure 10; this data represents a 2D surface of best (3rd family) Yukawa fits in the CMSSM soft-parameter space of m 0 , M 1/2 and a 0 . Any good fit of t-b-τ unification in the CMSSM would thus be expected to always lie in a compact region around the hypersurface: the m 0 and a 0 values would need to lie in the region of low χ 2 , while the M 1/2 value would need to lie near the one for the best-fit point. Results show that M 1/2 values of most best-fit points with χ 2 < 9 lie in the range between 2.5 TeV and 6 TeV, with the value increasing with increasing m 0 and |a 0 |.
• Figure 11 shows the predicted mass m A 0 (at 1-loop) of the neutral CP-odd MSSM Higgs A 0 , which is the main result of interest. Note that CP is not broken at 1-loop, because our parameters do not have complex phases. We see that all best-fit points in the allowed region of the m 0 -a 0 plane give a relatively low mass m A 0 , roughly in the range between 150 GeV and 1200 GeV. Important note: the m A 0 values are given only for the best-fit points, so one should be careful not to interpret the figure as a precise prediction of the CP-odd Higgs mass as a function of only a 0 and m 0 . Note the following important reservation about the results: they merely show the "naive" predicted mass of the extra Higgs particles in the CMSSM model. Potential experimental constraints have not been considered in this plot. In fact, as shall be discussed in the next section, practically the entire region predicted here (assuming exact t-b-τ unification) is under severe stress from ATLAS and CMS searches of H 0 → τ τ .

Challenges to t-b-τ unification
We have seen in section 4 that the scale of the extra MSSM Higgses is generically expected to be low in t-b-τ unification. The ultimate reason lies in the RG flow of the quantity m 2 H d − m 2 Hu , which was analyzed in section 3, and found to have a relatively small yet positive value, the latter being important for consistent EWSB. In this section, we analyze the predictions of t-b-τ further and confront them with experimental data from the LHC. SUSY spectrum with SO(10) boundary conditions Figure 12. As a first step, we extend the CMSSM scenario to the more general one with SO(10) boundary conditions, where the parameters consist of those in eq. (4.2), while the χ 2 is again defined with the observables of eq. (4.4). The standard deviations are taken as follows: the relative errors of the 3rd family Yukawa couplings are taken to be 1 %, while the error of for the SM Higgs mass is taken to be 2 GeV due to theoretical uncertainties in the computation.
This time we compute the overall expectations from this setup (with no fixed parameter values) by computing posterior probability densities of quantities of interest in a Bayesian approach by use of the Markov Chain Monte Carlo algorithm.
This paragraph contains some technical details of the computation. The MCMC algorithm was performed with 12 parallel chains, each yielding 1.3 · 10 5 points after discarding the initial bunch of 10 4 in the burn-in period. The total number of used data points is thus 1.56 million. Vacuum existence at 1-loop was checked, but not vacuum stability under EM charge breaking.
The result of interest from the MCMC computation is the SUSY sparticle spectrum, which turns out to be quite predictive, due to good fits obtained only in a compact region of parameter space, analogously to section 4. The results are presented in figure 12, where we draw the 1-σ and 2-σ highest posterior density (HPD) intervals for the masses of the sparticles. We use the labelsg for gluinos,χ 0 i for neutralinos,χ ± i for charginos,ũ i for up-type squarks,d i for down-type squarks,ẽ i for charged sleptons andν i for sneutrinos, where the index i goes over different ranges for different types of superpartners, but always corresponds to increasing mass (these are mass eigenstates, so the index i is not directly related to flavor).

JHEP06(2020)014
We make the following comments on the sparticle spectrum results: • The lowest part of the SUSY spectrum are the extra Higgs particles H 0 , A 0 and H + . They are expected in the rough range between 500 GeV and 1000 GeV. This reproduces the results for the case of CMSSM from section 4.
• The next lightest states are the lightest neutralinoχ 0 1 and the lightest charged sleptoñ e 1 . We see from the expected ranges that the lightest supersymmetric particle (LSP) for some points must be the lightest charged slepton (i.e. the stau) instead of the neutralino. Such points are experimentally problematic, since they would predict a charged LSP as a dark matter candidate. We performed a second MCMC analysis with the added constraint that the LSP must be the neutralino; this addition only minimally changes the quantitative predictions for HPD intervals of the other parts of the spectrum, so we choose not to include a separate plot.
• The rest of the spectrum is higher than 2 TeV, with gluinos typically at > 5 TeV.
An interesting feature is that the sleptons are expected to have lower masses than squarks.
The predicted sparticle spectrum is mostly compatible with the LHC data and searches for these particles, with one notable exception: the extra MSSM Higgs particles. The most stringent constraint comes from the possible ditau decay of neutral Higgses H 0 /A 0 → τ τ . The general scenario relevant in our case is the so called hMSSM [88], which assumes for all SUSY particles other than Higgses to be above 1 TeV. It was shown that specifying only two parameters, tan β and m A 0 , is sufficient to uniquely predict other tree-level quantities. The observed ditau rate is consistent with the SM background, so the non-observation of H 0 or A 0 is summarized by upper bounds on tan β for a given m A 0 in the m A 0 -tan β plane. The latest ATLAS [68] and CMS [69] results on this, using the dataset with 36 fb −1 of integrated luminosity at √ s = 14 TeV, suggest a bound of m A 0 1.5 TeV at tan β ≈ 50. Based on figure 12, the t-b-τ model prediction for the mass of H 0 and A 0 is clearly in tension with the experimental bounds, at least for most of the otherwise available parameter space. In fact, a search among computed MCMC points showed that the extra Higgs masses in the scenario of SO(10) boundary conditions cannot go much higher than 1200 GeV (since that would incur a severe χ 2 penalty). Comparing the various contributions to χ 2 shows that the tension comes from the SM Higgs mass, which tends to be dragged too high for high values of the extra Higgses.
This result is consistent with the upper limit for the best fit points in the more constrained CMSSM scenario, see figure 11; the additional parameter gained by the split of m 0 to m 16 and m 10 in the SO(10) boundary conditions thus does not appear to gain much maneuvering space over CMSSM for increasing the masses of the extra Higgs states. The CMSSM region in figure 11 with high extra Higgs masses is located at small m 0 , i.e. m 0 500 GeV, while a 0 ∼ −5 TeV.
This result indicates that exact t-b-τ unification, at least within the SO(10) boundary conditions scenario, is under strain exactly because of the low masses of the extra MSSM Higgses, the very feature pointed out and studied in this paper. Since we are now interested also in the masses of the extra Higgses, we perform the minimization with more observables in the χ 2 . For the input we have the SO(10) boundary condition parameters, now also assuming a possible split in t-b and one right-handed neutrino (the one with the largest Yukawa coupling, i.e. the unified coupling, in the Dirac mass term) at the scale M R , which may now be below M GUT . The other two Majorana type masses of the right-handed neutrinos are again set at the GUT scale. The input parameters are now Deformation scenario parameters: tan β, y 0 , y t , where the unified Yukawa coupling now excludes the top coupling y t : As for the χ 2 , we consider the observables from eq. (4.4), with two additional penalty terms. The first penalty term is associated to the non-observation of H 0 /A 0 → τ τ at the LHC, and is present only if tan β is too high given the value of m A 0 . The expected values of the tan β upper bound and 1-σ upper error of the constraint (extended to bigger errors assuming a Gaussian profile) are taken from figure 10b from the ATLAS analysis [68]. The other penalty basically enforces the neutralino to be the LSP, which turns out to be easily possible.

JHEP06(2020)014
We now fix the t-b deformation quantity y t /y 0 −1 and M R , and perform a minimization in the other parameters. We do so for each point in a 7 × 7 grid of equidistant points in the "deformation plane" of y t /y 0 − 1 and M R . The results of the minimized χ 2 (using interpolation of the grid results to show contours) is shown in figure 13. The range of t-b deformations is taken from 0 to 6 %, while the right-handed neutrino scale M R is considered on a logarithmic axis in the range between 10 13 GeV and 10 16 GeV. Note: the points were checked for the existence of the EW vacuum at 1-loop, but not explicitly for vacuum stability due to too excessive computation time. On the other hand, the points are close to points which have been checked with Vevacious, and overall in an unproblematic region with respect to vacuum stability. The numeric value of the µ parameter at the SUSY scale is in the interval (−10 TeV, −6.7 TeV) for all points. All points in the figure have the extra MSSM Higgs particles as the lowest lying states at around 1.3-1.5 TeV in the sparticle spectrum, followed by the neutralino with a mass > 2 TeV.
As stated earlier, the main difficulty is the reconciliation of the SM Higgs mass with the H 0 /A 0 → τ τ constraint on extra Higgs masses. The best fit points all have small m 10 , i.e. m 10 < 500 GeV, as in the CMSSM case, but the m 16 -m 10 split now allows for a bit bigger a 0 in magnitude without compromising χ 2 : a 0 ∼ −10 TeV.
The results clearly show that the t-b deformation at a few percent level can indeed greatly reduce the tension (for example the blue region in the plot corresponding to χ 2 < 6). This actually happens in two ways: first, it increases the masses of the extra Higgs particles and thus m A 0 (RGE effect), and second, it allows for a smaller tan β of around 46, which also relaxes tension, since H 0 /A 0 → τ τ constraints are in the form of an upper bound on tan β. In addition, figure 13 also shows that the fit is improved by a lower right-handed neutrino scale, but the effect is sub-dominant compared to the t-b deformation.
Another important result of the minimization in the grid worth stating is also the following: the best fit points still tend to have the extra Higgs masses at the lower end of the allowed range. The non-deformed points under tension have the Higgs just above 1300 GeV, while the deformed points not-under tension have those masses up to 1500 GeV. Though the ditau constraint did not require them to be higher than around 1500 GeV, this still shows that the deformed points have a preference for lower rather than higher masses of m A 0 . A continuing non-observation of the ditau decay coming from H 0 /A 0 neutral MSSM Higgses at the LHC would thus put the other points under increasing strain as well, requiring an ever larger t-b deformation.

Conclusions
We considered in this paper t-b-τ Yukawa unification in the context of SO(10) SUSY GUTs with µ < 0. The µ < 0 is the preferred sign for Yukawa unification, since it provides the SUSY threshold corrections to the b quark in the correct direction. Below the GUT scale, a good effective description is a softly broken MSSM possibly extended by right-handed neutrinos (if they are not yet integrated out). The boundary condition for the soft parameters at the GUT scale are assumed to be CMSSM-like, except for an additional split of the scalar soft mass parameter m 0 into sfermion masses m 16 and the mass parameter m 10 of JHEP06(2020)014 the Higgs doublets H u and H d , since these two soft mass parameters involve particles from different SO (10) representations. In particular, the features most important for comparison with the existing literature are exact Yukawa unification as opposed to quasi-unification, Hu at the GUT scale, µ < 0, and universal gaugino masses. We consider the above scenario to be the vanilla setup for Yukawa unification in SO(10), yet this has remained a largely unexplored possibility in the literature, where one or more of our stated assumptions are violated in an important way. The reason for that was a pessimistic outlook on the possibility of REWSB, based on approximate semi-analytic solutions of RGEs. In contrast, we show in this paper that REWSB is in fact possible to achieve by solving the full set of RGEs numerically.
The quantity of interest for successful EWSB is m 2 H d − m 2 Hu , which must be positive at the SUSY scale. In the large tan β regime needed for Yukawa unification, this same quantity determines also the mass scale of the extra MSSM Higgs particles H 0 , A 0 and H ± (cf. section 2). We find that the running quantity m 2 H d − m 2 Hu vanishes at the GUT scale due to the boundary conditions, first runs to negative values at lower scales, but the trend then reverses and it results in a positive value at M SUSY . Crucially, this positive value is smaller than might be expected based on the scale of the soft parameters, typically below TeV (when assuming exact t-b-τ Yukawa unification at the GUT scale). This yields a SUSY mass spectrum with the characteristic feature that the extra Higgs states are the lowest lying sparticle states, a feature that we focused on in this paper.
We study in detail the 1-loop RGE running of the quantity m 2 H d − m 2 Hu in section 3; we analyze the various contributions to its beta function, as well as determine the sensitivity to various deformations of boundary conditions. We find that the low mass feature for the extra MSSM Higgs particles is very sensitive to the exactness of t-b unification, with a 10 % percent deformation easily raising the scale by a factor of 2. The b-τ unification, presence of right-handed neutrinos, or a split of a universal scalar soft mass m 0 into the sfermion and Higgs parameters m 16 and m 10 , on the other hand, produce numerically a far more modest effect. Given the large sensitivity to t-b deformations, we conclude that a top-down RGE calculation is more suitable to accurately model the extra Higgs masses in exact t-b-τ unification.
This effect of low extra Higgs masses is ubiquitous in the entire parameter space, at least where t-b-τ unification leads to realistic Yukawa values at low energies. Most of the parameter space, both in the CMSSM and in the SO(10) boundary condition scenario, where a good fit to the 3rd family Yukawa couplings and the SM Higgs mass can be obtained, favors the extra Higgs masses at less than 1 TeV (for the case of exact t-b-τ unification), as presented in sections 4 and 5.
These model predictions, however, are in tension with ATLAS and CMS searches of ditau decays of neutral extra Higgses, i.e. H 0 /A 0 → τ τ . The experimental searches result in upper bounds on tan β as a function of m A 0 . Since t-b-τ unification requires a large tan β ≈ 50, this suggests the extra Higgses to be above roughly 1.5 TeV. In exact t-bτ unification with correct Yukawa predictions at low scales, it is hard to achieve masses above ∼ 1.3 TeV; the main obstacle turns out to simultaneously obtain heavy extra Higgses alongside a sufficiently low SM Higgs mass near 125 GeV.

JHEP06(2020)014
The tension with experiment can be reduced by relaxing exact t-b-τ unification. As shown in section 5, a deformation of t-b unification at a level of a few percent can completely relieve the tension with experiment, both by raising the masses of the extra Higgs particles and lowering the required tan β. Such a deformation of a few percent could come about from GUT threshold corrections, especially given the large numbers of particles in the SO (10) representations in the Higgs sector (which are of course model dependent), or Planck scale suppressed operators. It should be noted, however, that even deformed t-b-τ unification prefers lower rather than higher extra Higgs masses.
In summary, we have shown that t-b-τ (quasi-)unification in SO(10) SUSY GUTs with µ < 0 generically features comparably light extra MSSM Higgs particles. For exact t-b-τ unification we find a tension with LHC constraints from H 0 /A 0 → τ τ , due to predicting too light masses of the extra MSSM Higgses. The tension can be successfully alleviated by relaxing the scenario to quasi-unification of Yukawa couplings: a few percent split of the top Yukawa from the unified value (most importantly from the bottom Yukawa) can bring the extra Higgs states to sufficiently high values to avoid the present experimental constraints. Nevertheless, masses of these states close to the present bounds are still preferred. This implies that a continuing non-observation of the extra MSSM Higgses would require ever bigger deformation of t-b-τ unification, finally disfavoring the scenario. Conversely, an observation of an extra Higgs state in the ditau decay channel could be the first sparticle observation of the t-b-τ unified SO(10) SUSY GUT model, and measuring a sparticle spectrum with extra Higgses having the lowest masses could be a hint for the realization of this scenario in nature.

JHEP06(2020)014
(A.10) JHEP06(2020)014 The loop factor c 1 is defined as  For the Majorana neutrino mass associated to the large 3rd family neutrino Yukawa coupling, we assume the value M ν3 = M R at the scale M R , implying that this heavy neutrino is integrated out at the scale M R . The M ν3 does not appear in the RGE of any other quantity.
Open Access. This article is distributed under the terms of the Creative Commons Attribution License (CC-BY 4.0), which permits any use, distribution and reproduction in any medium, provided the original author(s) and source are credited.