New physics models facing lepton flavor violating Higgs decays at the percent level

We speculate about the possible interpretations of the recently observed excess in the h → τμ decay. We derive a robust lower bound on the Higgs boson coupling strength to a tau and a muon, even in presence of the most general new physics affecting other Higgs properties. Then we reevaluate complementary indirect constraints coming from low energy observables as well as from theoretical considerations. In particular, the tentative signal should lead to τ → μγ at rates which could be observed at Belle II. In turn we show that, barring fine-tuned cancellations, the effect can only be accommodated within models with an extended scalar sector. These general conclusions are demonstrated using a number of explicit new physics models. Finally we show how, given the h → τμ signal, the current and future searches for μ → eγ and μ → e nuclear conversions unambiguously constrain the allowed rates for h → τe.


Introduction
The discovery of the Higgs boson at the LHC [1,2] imbues the standard model (SM) of particle physics with completeness and self-consistency. Nonetheless, its failure to account for non-vanishing neutrino masses is one of the main motivations for considering physics beyond the SM. Incidentally, the accidental SM symmetries that prevent neutrinos from acquiring mass also completely suppress lepton flavor violating (LFV) processes. The observation of the former thus provides ample motivation for a rich experimental program to search for the latter. The CMS collaboration has recently reported a slight excess with a significance of 2.4 σ in the search for LFV decay h → τ µ [3]. The best fit for the branching ratio of the Higgs boson to τ µ final state (summed over τ − µ + and τ + µ − ), assuming SM Higgs production, is found to be B(h → τ µ) = 0.84 +0. 39 −0.37 % . (1.1)

JHEP06(2015)108
This recent hint has expectedly received significant amount of attention in the literature [4][5][6][7][8][9][10]. It is thus an imperative to either confirm or reject validity of this tantalizing hint with more data by both ATLAS and CMS experiments. At the same time, it is instructive to revisit expectations for this observable within various new physics (NP) scenarios and in particular re-evaluate the feasibility of obtaining such a large signal in light of severe indirect constraints on LFV Higgs interactions coming from low energy probes.
Without loss of generality, one can parameterize the mass terms and Higgs boson couplings of charged leptons after electroweak symmetry breaking (EWSB) as where the ellipsis denotes non-renormalizable interactions involving more than one Higgs boson and i = e, µ, τ . In the SM, the Higgs couplings are diagonal and given by y ij = (m i /v)δ ij , where v = 246 GeV . On the other hand, non-zero y τ µ and/or y µτ induce h → τ µ decays with a branching ratio of Assuming that the total Higgs boson decay width (Γ h ) is given by its SM value enlarged only by the contribution from h → τ µ itself, i.e., Γ h = Γ SM h /[1 − B(h → τ µ)], where Γ SM h (m h = 125 GeV) = 4.07 MeV [11], the measurement in eq. (1.1) can be interpreted as a two-sided bound on the |y τ µ | 2 + |y µτ | 2 combination of couplings. These limits read (see also the left-hand side panel of figure 1) 0.0019(0.0008) < |y τ µ | 2 + |y µτ | 2 < 0.0032(0.0036) at 68% (95%) C.L. . (1.4) In general, the experimentally measured h → τ µ event yield depends not only on the values of y τ µ and y µτ , but also on other Higgs couplings contributing both to its total decay width Γ h as well as its production cross-section (σ h ). In particular, a given signal can be reproduced for larger (smaller) values of |y τ µ | and |y µτ | by enhancing (suppressing) Γ h and/or suppressing (enhancing) σ h . Since both Γ h and σ h affect other currently measured Higgs observables, their individual effects can be disentangled by performing a global fit to all Higgs production and decay event yields at the LHC. The details of this procedure can be found in appendix A, while the resulting best fit χ 2 values as functions of |y τ µ | 2 + |y µτ | 2 are shown in the right-hand side panel of figure 1. In particular, σ h is well determined by the measurements of both inclusive and separate exclusive Higgs production channels in several different decay modes. Similarly, Γ h is bounded from below by the observations (or at least indications [12]) of the dominant SM Higgs decay modes (h → bb, h → W W * , h → τ + τ − , etc.). Consequently the lower range of allowed |y τ µ | and |y µτ | values does not change much compared to eq. (1.4) when considering the global fit. On the other hand, the fact that the total Higgs decay width is currently only weakly bounded from above [13,14] allows significantly larger |y τ µ | and |y µτ | couplings to reproduce the same observed signal in the general case. 1 Numerically, we find 0.0017(0.0007) < |y τ µ | 2 + |y µτ | 2 < 0.0036(0.0047) at 68% (95%) C.L. . (1.5) The rest of the paper is devoted to interpreting these ranges in terms of hypothetical new physics (NP) effects. In section 2 we review the model independent considerations of flavor violating Higgs couplings using an effective field theory approach. Section 3 is devoted to a case study of the two Higgs doublet model (THDM) with a generic Yukawa structure (the so-called type-III THDM) that can account for the observed anomaly at tree level and be in agreement with present low energy flavor constraints. Scenario that relies on vector-like fermions to explain the data is presented in section 4. Possibility to explain the anomaly via one-loop physics effects is discussed in section 5, where we investigate phenomenology of scalar leptoquarks to demonstrate difficulties with this type of approach. In particular, section 5.4 gives a very special solution that can accommodate the anomaly through an interplay between two different one-loop effects, and satisfy constraints imposed by low energy flavor physics experiments at the price of fine-tuning. We briefly conclude in section 6.

Model independent considerations
In the following we assume the SM contains all the relevant degrees of freedom at energies O(few 100) GeV. 2 In particular, any additional heavy field can be integrated out, so that 1 It has been argued recently [15,16], that the current total Higgs decay width measurements using offpeak data [13,14] introduce some model-specific assumptions into the fit. We have checked that removing the total Higgs width constraint from the fit, no upper bound on |yτµ| 2 + |yµτ | 2 can be set. 2 In the minimal extension of the SM with Dirac neutrinos accounting for the observed neutrino masses, a non-vanishing h → τ µ rate is induced at the one-loop level. However, this contribution is negligibly small since the relevant decay amplitude is suppressed by two powers of neutrino Yukawa couplings in addition to the tau Yukawa.

JHEP06(2015)108
we can apply effective field theory (EFT) methods. If the electroweak symmetry is realized linearly, then current data already strongly suggest that h for the most part should be a neutral CP even component of a linear combination of Y = 1/2 hypercharge weak doublets H α = (h + α , v α + x α h + . . .) T , where the ellipsis denotes possible projections to other heavier scalar mass eigenstates and α |x α | 2 1/2. In general, EWSB can receive contributions from several sources (for example several Higgs multiplet condensates). However, electroweak precision tests severely constrain contributions from higher representations, so that also α v 2 α v 2 /2. In the decoupling limit we are considering, the contributions from the heavier Higgses are power suppressed and can be completely described in terms of higher dimensional operators. The leptonic Yukawa sector of such a theory, including the leading (dimension six) non-renormalizable operators suppressed by the EFT cut-off scale Λ, can be written as where L i and E i correspond to the leptonic weak doublets (L i = ( i L , ν i ) T ) and singlets (E i = i R ), respectively, and we explicitly show the flavor indices i, j = 1, 2, 3. Any additional dimension six operators coupling the Higgses to the leptons can be shown to be either redundant or not to contribute to the effective Lagrangian (1.2) after EWSB [17]. We can thus identify withv α = v α /v and the unitary matrices V L,R diagonalize the leptonic mass term as We can note immediately the two possible sources of non-vanishing y τ µ and/or y µτ . In a theory with multiple Higgs doublets, the first term in the square brackets in eq. (2.3) can be non-vanishing if for at least two terms in the α-sum (with non-commuting λ α ) x α =v α . Thus (lepton) flavor violating Higgs boson interactions are possible even without additional NP contributions above the EFT cut-off. The simplest example of such a theory, the two Higgs doublet model (THDM) of type-III has recently been analyzed in this context in refs. [6,10] and we consider it in more detail in section 3. On the other hand, in a single Higgs theory only v 1 = v/ √ 2 and x 1 = 1/ √ 2 are non-zero so that the first term in the square brackets in eq. (2.3) vanishes identically. Thus, the observation of h → τ µ in this context necessarily implies the presence of the second term in eq. (2.1). In particular, the CMS measurement can be interpreted directly in terms of the effective NP scale of We observe that percent level B(h → τ µ) can be accommodated within the EFT approach (containing a single light Higgs doublet) with a sizable mass gap. At the same time, new degrees of freedom generating the observed h → τ µ excess should be within the kinematic reach of the second LHC run.
Since λ and λ contribute to both as well as m the observed hierarchical structure of the charged lepton masses can be used to define the natural ranges of ij . In particular, without delicate cancellations between the various contributions in eq. (2.4), the observed hierarchy between the muon and tau lepton masses implies [18,19] This constraint is compared to the CMS preferred range for |y τ µ | and |y µτ | from eq. (1.4) in figure 2. We observe that the two indications are not in sharp conflict and thus the observation of B(h → τ µ) at the percent level can in principle be explained by natural NP although a mild hierarchy between y τ µ and y µτ is preferred in this case. Another class of model-independent constraints comes from the fact that the dimension six operators contained in L Y in eq. (2.1) will in general mix with other operators contained in the effective SM. In particular, restricting the discussion to a single Higgs setup with H ≡ H 1 , λ will mix under charged lepton Yukawa renormalization into (H † H) 3 affecting, for example, Higgs boson pair production at the LHC. However, as shown in refs. [20,21], such mixing is suppressed by three powers of the charged lepton Yukawa matrix and is thus (baring fine-tuned cancellations) completely negligible in practice. More phenomenologically relevant are the UV finite one-and two-loop contributions to the operatorsLH(σ·B)E andLτ a H(σ · W a )E, where σ µν = i[γ µ , γ ν ]/2, B µν and W a µν are the hypercharge and weak isospin field strengths, respectively, and τ a are the Pauli matrices. These operators can mediate radiative flavor violating lepton decays as well as contribute to leptonic anomalous electric and magnetic dipole moments, all of which are already tightly constrained by JHEP06(2015)108 experiment. The most stringent constraint comes from the τ → µγ decay [17], mediated by the effective Lagrangian where Q L,Rγ = (e/8π 2 )m τ (μσ αβ P L,R τ )F αβ , P L,R = (1 ∓ γ 5 )/2 and F αβ is the electromagnetic field strength tensor. The Wilson coefficients c L,R receive comparable one-and two-loop contributions. In the experimentally justified approximation y µµ y τ τ (both assumed real) and m µ m τ m h , they are given by [17,22,23] c where y tt is the top quark Yukawa with the SM value of y tt =m t /v = 0.67 for a MS top mass ofm t 164 GeV. The resulting EFT correlation between B(h → τ µ) in eq. (1.3), and B(τ → µγ) given by is shown in figure 3 (diagonal dashed orange line), assuming SM values of all Higgs boson couplings except y τ µ and y µτ . In the same plot, the CMS preferred range of B(h → τ µ) in eq. (1.1) is displayed by the horizontal blue band, while the current (B(τ → µγ) < 4.4 × 10 −8 @ 90% C.L.) [24] and projected future (B(τ → µγ) < 3 × 10 −9 @ 90% C.L.) [25] indirect constraints are shaded in light and dark gray vertical bands, respectively. We observe that within the EFT approach, the CMS signal is well compatible with the nonobservation of τ → µγ at the B factories and will marginally remain so even at Belle II. This is in contrast with the situation in most explicit NP models generating non-zero y τ µ or y µτ , as we demonstrate in sections 3-5 . Before going through explicit examples, let us consider the appearance of (lepton) flavor violating Higgs interactions from a symmetry point of view. In absence of Higgs Yukawa interactions, the SM possesses a large global flavor symmetry G F = SU(3) 5 (where we have omitted the U(1) factors). In particular, the leptons transform as L ∼ (3, 1), In the SM (without neutrino masses), the charged lepton Yukawa matrix λ ∼ (3,3) is the only source of G breaking. Consequently all lepton interactions are flavor conserving in the charged lepton mass basis. Conversely, as also demonstrated explicitly in eq. (2.3), the generation of lepton flavor violating Higgs interactions requires at least two non-aligned sources of lepton flavor symmetry breaking. At the tree level, there are only two possibilities: (1) one can enlarge the SM scalar sector, such that more than one Higgs doublet couples to the leptons (corresponding to the first term in eq. (2.3)); (2) one can extend the leptonic sector by vector-like fermions, whose JHEP06(2015)108  [3]. Current and future projections for B(τ → µγ) experimental sensitivity are represented with vertical light [24] and dark [25] gray bands, respectively. Superimposed are the predictions within the EFT approach (diagonal dashed orange line), in the type-III THDM (green and black bands), in models with vector-like leptons (diagonal dotted purple line) and in models with scalar leptoquarks (diagonal red and orange shaded band). See text for details.
Dirac masses and mixing terms with SM chiral fields can provide additional sources of G breaking. This leads to the appearance of the λ contributions after integrating out the new heavy fermionic states. Both possibilities are explored in the following sections. Example of an enlarged Higgs sector is given in section 3 whereas the vector-like fermion case is discussed in section 4.
Of course many more possibilities exist at the loop level (e.g. [26,27]). In practice however, all such models predicting sizable LFV Higgs interactions necessarily suffer from a severe fine-tuning problem, if they are to simultaneously avoid the stringent τ → µγ constraint. We demonstrate this on an explicit example in section 5.4. A general heuristic argument however goes as follows. The operatorLHE(H † H) has exactly the same transformation properties under both G and chiral transformations of leptons as opera-torsLH(σ · B)E andLτ a H(σ · W a )E. Thus, the simultaneous presence of the former and absence of the latter cannot be protected by a symmetry valid in the low energy theory. If LHE(H † H) is generated at the loop level, it will necessarily involve charged states propagating in the loop. These charged states will couple to photons generating in turn c L,R in eq. (2.7) at the same loop level as y τ µ and y µτ , respectively. The two sets of matching amplitudes will also exhibit identical dependence on the fundamental flavor parameters of the model. The resulting expectation c L,R ∼ 8πy τ µ,µτ /vm τ is clearly in conflict with JHEP06(2015)108 existing experimental results, since a percent level B(h → τ µ) would imply even an order of magnitude bigger B(τ → µγ)! The accidental cancelation required in the matching procedure for the radiative decay in order to accommodate both results is thus expected to be approximately one part in 10 3 at the amplitude level. As we demonstrate in sections 3, 4 and 5, these expectations are confirmed in explicit models.
A positive experimental indication for the h → τ µ decay when combined with existing stringent experimental limits on µ − e LFV processes can be used to constrain LFV τ − e processes. In particular, the product of the h → τ e and h → τ µ branching fractions is constrained from above by the rates of µ → eγ, and µ − e conversion on nuclei. This is due to the fact that tree level Higgs decays to τ µ (τ e) depend on y µτ,τ µ (y eτ,τ e ) while the same sets of couplings contribute at the loop level to µ → eγ and µ − e conversion via diagrams with a virtual τ .
The contributions to the µ → eγ process stemming from a virtual τ are m τ enhanced with respect to diagrams with intermediate µ or e states [17] leading to where we have neglected the effects of the light lepton masses while m τ dependence is retained in its exact form. The coefficient c τ R is obtained from eq. (2.11) by replacing y ij → y * ji . The µ → eγ branching fraction is thus sensitive to a distinct combination of the LFV Yukawas: On the other hand, µ − e conversion on nuclei is most sensitive to vector current effective operators (ēγ ν P L,R µ) (qγ ν q) which probe different combinations of LFV Yukawas. The largest contribution to the vector coupling of the proton follows from eq. (A17) of ref. [17] and reads, in the limit of massless light leptons, The above effective coupling is evaluated in the q 2 → 0 limit, where q is the momentum exchanged in the conversion process. As before the right-handed current coefficient g is obtained from the above result by replacing y ij with y * ji . Note that due to strong experimental limits on B(µ → eγ) we can at present safely neglect the dipole operator Q L,R contributions to µ − e conversion [17]. The relative conversion rate B(µ → e) Au ≡ Γ(µ → e) Au /Γ capture Au is then given by (2.14)
The complementary information on the LFV couplings extracted from the rates of µ → eγ and µ − e conversion is sufficient to completely correlate the Higgs LFV h → τ µ and h → τ e decays: where in the second line we have approximated Γ h Γ SM h . The best experimental limit on B(µ → e) Au < 7×10 −13 (at 90% C.L.) was achieved by the SINDRUM II Collaboration [29] while the best upper bound on B(µ → eγ) < 5.7 × 10 −13 (at 90% C.L.) was recently determined by the MEG Collaboration [30]. Note that with the current experimental data the sum on the right-hand side of eq. (2.15) is completely saturated by the µ − e conversion contribution. Combining the two bounds with the central value for the h → τ µ branching fraction in eq. (1.1) leads to an upper bound B(h → τ e) < 0.26 . This is above the current indirect constraint coming from searches for τ → eγ [24] which reads B(h → τ e) < 0.19 . In the near future the experimental sensitivity to both µ → eγ and especially µ − e conversion in nuclei is expected to improve significantly. In particular the DeeMe [31] experiment aims at a sensitivity of the order of 10 −14 while the Mu2e [32] experiment could improve existing SINDRUM II bounds by four orders of magnitude. Given B(h → τ µ) at the percent level, such a measurement will indirectly probe B(h → τ e) at the order of magnitude 10 −5 .

Type-III two Higgs doublet model
Several well motivated scenarios of physics beyond the SM (e.g. supersymmetric extensions or models addressing the strong CP problem via a PQ symmetry) require the addition of another SU(2) L doublet scalar to the minimal SM Higgs sector (cf. [19] for a recent general review). The most extensively studied version of the THDM is the type-II model in which the absence of tree level flavor changing neutral currents (FCNCs) is imposed by the requirement that only one of the doublets couples to up-type quarks, while the other one couples to down-type quarks and charged leptons. FCNCs in the quark sector are relegated to the one-loop level via the interactions involving charged scalars, whose couplings to quarks follow the CKM structure of the SM. In this model LFV is absent as in the SM due to massless neutrinos and because only a single λ α in eq. (2.1) is non-vanishing.
In order to accommodate LFV Higgs interactions one thus has to depart from this version and consider the THDM with a generic Yukawa structure (type-III THDM). In the following we restrict our discussion to the parameter space with a MSSM-like Higgs potential, for which the following relations hold [33] tan  where 0 < β < π/2 is the angle of the rotation matrix that diagonalizes the mass squared matrices of charged scalars and pseudoscalars and −π/2 < α < 0 is an analogous angle for the neutral scalars. The vacuum expectation values of the doublets are denoted by v u and GeV. In our input we set the mass of the Higgs boson to the observed value m h 125 GeV and choose β and m A to be free parameters of the scalar potential. Such a THDM can arise as a low energy effective theory of the MSSM in the limit in which all the heavy superpartners are integrated out. The LFV interactions of the Higgs boson then originate from diagrams involving heavy sleptons in the loop, cf. [34,35]. Note however that generically in the MSSM, LFV radiative lepton decay contributions will be generated at the same order as discussed in section 2, leading to tighter constraints, cf. [36]. We thus do not attempt to interpret the type-III THDM model parameters in terms of the underlying MSSM parameters, but consider them as essentially uncorrelated, only subject to experimental bounds and theoretical constraints such as perturbativity conditions. The relevant part of the Yukawa Lagrangian in the charged lepton mass eigenstate basis is [33] 3)

JHEP06(2015)108
where the Yukawa couplings can be written as The terms f i parametrize the off-diagonal charged lepton Yukawa couplings and are free parameters of the model. Finally, the coefficients x k q for H k = (H, h, A) are given by [33] x k u = (− sin α, − cos α, i cos β) , x k d = (− cos α, sin α, i sin β) . (3.5) In addition to tree level exchanges of neutral Higgses we have loop contributions to FCNCs due to charged Higgs, whose interactions with leptons are parameterized by the following Yukawa couplings: Using preceding relations, at the tree-level the type-III THDM contributes to the effective Higgs couplings as and we can translate eq. (1.4) into a 1 σ two-sided bound on the parameters For reference we also give the expression for the branching fraction We show in figure 5 the points in the plane spanned by (| τ µ | 2 +| µτ | 2 ) 1/2 and m A for which B(h → τ µ) = 0.84% when all other Higgs decay rates are held SM-like. Interestingly, the decay exhibits small dependence on β for large values of tan β. Indeed, an expansion around the decoupling limit in small ξ = m 2 Z /m 2 A gives sin α tan β + cos α − 2 tan β(tan 2 β − 1) (3.10) For light pseudo-scalar masses (m A ), mixing with the heavy scalar Higgs could affect the couplings of the light Higgs to the gauge bosons. However, the overall modification factor of hW W and hZZ couplings relative to the SM values is given by sin(β − α) and starts to differ from 1 only at order O(ξ 2 ): (3.11) We have checked that even for m A = 150 GeV the hW W coupling deviates less that 2 % from its SM value, for any value of tan β. We conclude that the effect of modified light Higgs couplings to gauge bosons can be safely neglected at the current precision of Higgs properties measurements. The more important effects of possible modifications of the light Higgs (tau) Yukawa couplings are discussed in the next subsection.

Constraints from τ → µγ
The coefficients c L,R in eq. (2.7) are generated via one-loop diagrams with charged or neutral Higgses. In the limit m µ = 0 (for the sake of consistency also y H k µµ = µµ = 0), these coefficients are [33] c The suppression of the above result by two small Yukawa couplings is relaxed at the twoloop level, where Barr-Zee type diagrams require only a single LFV Yukawa coupling along with y tt . Therefore the contributions of the neutral Higgses at the two-loop level can dominate over the one-loop result. Here we employ the Barr-Zee contributions due to neutral Higgses calculated for µ → eγ in ref. [37] and adapted to τ → µγ in appendix B.1. The 2-loop contributions to B(τ → µγ) involve exactly the same combination of elements as B(h → τ µ) in eq. (3.9). For the charged Higgs contributions to the best of our knowledge no complete calculation exists in the literature. A partial result with the charged Higgs contribution inducing the H 0 k γγ * or H 0 k γZ * vertices has been presented in ref. [35]. However, the diagrams with an effective H ± W ∓ γ coupling via loop of t and b are still missing. In the following we use the known partial results to estimate the order of magnitude of the constraint. The missing contributions could potentially further strengthen the τ → µγ bound but are not expected to substantially weaken it.

JHEP06(2015)108
In order to explore the nontrivial correlation between the B(h → τ µ) and B(τ → µγ) we sample the parameter space of τ µ , µτ , τ τ for a specific choice of tan β and m A . The ranges allowed for τ µ,µτ are required to fulfill the naturalness criterium of eq. (2.6) 13) and that each of the effective couplings is within the perturbative range (3.14) The 1-loop amplitude of the τ → µγ process further depends on the diagonal y h τ τ Yukawa coupling which is experimentally constrained by the searches for h → τ τ decays within ATLAS and CMS experiments. In fact both collaborations report strong evidence for the existence of this mode with the signal strengths (normalized to SM expectations) µ τ τ = 1.43 +0.43 −0.37 and µ τ τ = 0.78±0.27, observed by ATLAS [38] and CMS [39], respectively. In the following we use a naïve average of the two experimental results: In the decoupling limit (cf. eq. (3.11)) and with y h tt SM-like, we can assume SM-like Higgs production cross sections and model the departure of µ τ τ from 1 by an appropriate shift of y h τ τ away from − √ 2m τ /v (i.e. non-vanishing τ τ ). A scenario with SM-like y h τ τ coupling corresponding to µ τ τ = 1 for fixed tan β = 10 and m A = 0.3 TeV is presented in figure 3 by a black narrow stripe. This scenario easily passes both experimental constraints. A perfect one-to-one correspondence between the two observables is spoiled by the asymmetry between the left-and right-handed τ → µγ amplitudes at one-loop, see eq. (3.12). On the other hand, for masses m A significantly larger than 500 GeV it is not possible to reconcile both predictions with the corresponding experimental values.
Allowing τ τ to vary within the experimental bounds imposed by the h → τ τ decays, i.e. eq. (3.15), one obtains predictions for h → τ µ and τ → µγ as represented by the green band in figure 3. In particular, this additional freedom in τ → µγ breaks, to some extent, the strict correlation with the h → τ µ rate. The remaining correlation signals the dominance of the 2-loop contributions in B(τ → µγ) which are independent of τ τ .
The interaction vertices that induce the LFV decay of Higgs boson to τ µ pairs also contribute to the LFV decay of tau lepton to three muons, τ − → µ − µ + µ − . The current best experimental upper limit on its branching fraction has been reported by the Belle Collaboration B(τ − → µ − µ + µ − ) ≤ 2.1 × 10 −8 at the 90% C.L. [40]. This process proceeds through the tree level exchange of the neutral scalars h, H and the pseudoscalar A. Since all the relevant amplitudes are suppressed by flavor-diagonal muon Yukawa, this constraint is not competitive with τ → µγ in terms of sensitivity to τ µ and µτ .

h → τ µ vs. h → τ e
Here we repeat the derivation of the bound on h → τ e described in section 2.1 in the context of the type-III THDM model. In this case, both µ → eγ and µ − e conversion receive

JHEP06(2015)108
contributions not from a single Higgs but from three neutral Higgses while the charged Higgs contributions are relatively suppressed by light lepton masses. The effective operator coefficients presented in section 2.1 now have to be summed over the three Higgses with proper mass and effective coupling replacements. Ultimately, the observables depend on t β ≡ tan β, m A , and the LFV parameters ij , where dependence on the first two parameters can be absorbed into an overall coefficient:  Tree level LFV Higgs decays h → τ µ and h → τ e are still properly described in the effective framework of eq. (3.7). In contrast with what has been found in the effective framework (cf. eq. (2.14)), the sensitivity to ij is now similar in both low energy observables. Also, the coefficients B µ→eγ,µe 0 now have a strong dependence on tan β and m A . Finally, Higgs decays are suppressed by a projection factor so that . The above bound is monotonically decreasing with tan β and m A as can be seen on figure 4. Note that at the projected sensitivity of the Mu2e experiment, the right hand side of eq. (3.17) would start to be dominated by the µ → eγ constraint, so an improvement expected from the planned MEG upgrade (MEGII) [41] will be essential in this context, as can be seen by comparing the two bottom-most dashed contours in the plot.
where g 0.65 and c W ≡ cos θ W 0.88 are the weak isospin gauge coupling and the cosine of the weak mixing angle, respectively. In particular one obtains a relation y ij = (m j /v)(X ij − Y ij ) . Severe constraints on X τ µ,µτ , Y τ µ,µτ 10 −3 from searches for τ → 3µ [52] then preclude any sizable contributions to B(h → τ µ). 3 De-correlating LFV Higgs and Z boson couplings requires the introduction of VLLs in both SM gauge representations (Ψ E , Ψ L ) mixing with chiral leptons Integrating out the vector-like states one finds the LFV Higgs and Z boson couplings can now be completely de-correlated in the limit where the direct couplings of chiral fermions to the Higgs vanish, and the SM leptons obtain their masses exclusively through mixing with VLLs in eq. (4.2). In particular the lepton flavor off-diagonal Higgs Yukawas in eq. (2.2) are then of the form [54] All the LFV phenomena are driven by this contribution and one can derive a one-toone correlation between tree-level Higgs decay h → τ µ and the 1-loop radiative decay The 1-loop suppression factor of τ → µγ is not small enough to evade the experimental upper bound and at the same time accommodate the h → τ µ at the percent level, as illustrated by a pink-dotted line in figure 3.

Scalar leptoquarks
A scalar leptoquark state (LQ) can induce h → τ µ decay via quark-LQ penguin diagrams. This type of diagrams requires a helicity flip on one of the fermion lines. One therefore expects suppressed amplitude when all the fermions are substantially lighter than v. It turns out that size of LQ Yukawa couplings must then be well beyond the perturbative regime to result in B(h → τ µ) at the percent level. In this respect, only LQ states that JHEP06(2015)108 couple to charged leptons and top quark are suitable candidates that will be studied in what follows. Inspection of helicity structure of the diagrams reveals that both left and right chiralities of leptons and top quark have to couple to the LQ state. The physical Higgs can couple to the scalar ∆ or to the top quark as shown in the first two diagrams in figure 6. While the strength of the coupling relevant for the latter process is fixed by the top Yukawa, the former process depends on an unknown hLQLQ coupling, λv, that originates from the marginal "Higgs portal" operator, Here ∆ is the scalar LQ and H is the SM Higgs doublet.
The Yukawa couplings of ∆ 1 are given by the following Lagrangian where Q i = (u i L , d i L ) T and U i = u i R are the quark weak doublets and up-type singlets, respectively. We explicitly show flavor indices i, j = 1, 2, 3, and SU(2) indices a, b = 1, 2, with 12 = 1. Also, here y L ij and y R ij are elements of arbitrary complex 3×3 Yukawa coupling matrices. After expanding the SU(2) indices, we obtain where V CKM and V PMNS represent Cabibbo-Kobayashi-Maskawa and Pontecorvo-Maki-Nakagawa-Sakata mixing matrices, respectively. All fields in eq.  4) where N c = 3 is the number of colors. The relevant loop function further depends on the portal coupling λ, given in terms of Passarino-Veltman functions B 0 and C 0 . The h → τ µ decay width in the ∆ 1 LQ scenario is then given by The state ∆ 1 also contributes to the τ → µγ through leptoquark Yukawas present in h → τ µ. The Wilson coefficients of effective Lagrangian introduced in eq. (2.7) are JHEP06 (2015)108 [55]. The τ → µγ branching fraction in presence of ∆ 1 depends on identical Yukawa combination as the h → τ µ decay rate: Additional constraints on the model parameters can also come from the measurements of (g − 2) µ and the electroweak precision observables, i.e. B(Z → bb). The ∆ 1 contribution to the former is dominated by the coupling combination y L tµ y R tµ which can always be suppressed. For the Z → bb the presence of ∆ 1 will modify the effective coupling δg b L , defined as: Here g is the SU(2) L gauge coupling, θ W is the Weinberg angle. g b L,R are the SM effective couplings including the electroweak radiative corrections, while δg b L,R parameterize contributions that originate from beyond the SM. The LQ state ∆ 1 contributes to δg b L through loops with neutrinos, The constraint on LQ Yukawas follows from combining eq. (5.10) with the results of the recent fit to the Z-pole observables [51]: The Z-pole measurements are sensitive to interference terms between g b L,R and δg b L,R , where g b L,R are approximately real numbers, and therefore only real parts of δg b L,R give linear shift in g b L,R . According to ref. [51] the SM point corresponds to χ 2 SM = 7.5, where 0 is defined as the fit with free δg b L,R . The presence of LQ will push δg b L in the negative direction. At 1σ we insist on χ 2 < χ 2 SM + 1 that corresponds to a lower limit ∆g b L > −1.2 × 10 −4 . An approximate constraint on the Yukawas is finally obtained, Clearly, this constraint does not preclude the Yukawa couplings combination entering B(h → τ µ) from being large.
In the left-hand side panel of figure 7 we show the dependence of B(h → τ µ) on the leptoquark Yukawa couplings assuming the total Higgs decay width to be SM-like. Here, we turn off the Higgs portal coupling. We set the mass of ∆ 1 LQ to 650 GeV, consistent with current direct search limits at the LHC [56,57]. In the same panel of figure 7, we show the CMS measurement reported in [3]. The B(τ → µγ) dependence on the leptoquark Yukawa couplings for ∆ 1 case is shown in the right-hand side panel of figure 7 (dashed line), again for m ∆ 1 = 650 GeV. The following transformation needs to be applied when going from the ∆ 1 to the ∆ 2 case: y L ij → y L ji .

Impact of the Higgs portal and its constraints
Next we study the impact of a non-zero Higgs portal coupling in the scenario with ∆ 1 . As an example, for m ∆ 1 = 650 GeV the loop function dependence on the portal coupling is Thus, a positive large λ could in principle relax the leptoquark Yukawa couplings and yield sizable h → τ µ rates without violating the τ → µγ constraint. However, the Higgs portal coupling also induces corrections to the h → γγ decay and to gluon-gluon fusion (ggF) induced Higgs production with the leptoquark running in the triangular loop. The modified ggF production cross section normalized to its SM value is Here, N ∆ i is the number of ∆ i components in the weak multiplet ∆. C(r ∆ ) is the index of color representation r ∆ of ∆ and for the triplet (C(3) = 1/2). We consider heavy enough colored scalars such that the loop function is in the decoupling limit. Similarly, the modified h → γγ decay width, normalized to its SM value, is given by The sum in eq. (5.15) runs over all weak components of the SU(2) L multiplet. d(r ∆ ) and Q ∆ i are the dimension of the color representation of ∆ and the electric charges of weak ∆ i components, respectively. We fit the latest LHC Higgs data (including the CMS signal of h → τ µ) taking λ and leptoquark Yukawa couplings as free parameters. The result for the ∆ 1 LQ model is shown in figure 8. The preferred regions at 1 σ and 2 σ are rendered in solid JHEP06(2015)108 pink and dashed pink, respectively. If besides the λ and leptoquark Yukawa couplings, we also allow other Higgs couplings to vary (see appendix A) we instead get the dark and light grey regions at 1 σ and 2 σ levels, respectively. The correlation between the h → τ µ and τ → µγ branching ratios in presence of the (3, 1, −1/3) leptoquark state for |λ| < 1 and m ∆ 1 > 600 GeV is depicted in figure 3 with a pink stripe.

The ∆ 2 = (3, 2, 7/6) case
The Yukawa couplings of the ∆ 2 leptoquark to SM fermions are where we explicitly show the flavor indices i, j = 1, 2, 3, and SU(2) indices a, b = 1, 2. y L and y R in eq. (5.16) are arbitrary complex 3 × 3 Yukawa matrices. In the mass eigenstate basis we have where the loop function g 1 has been introduced in the previous section. The h → τ µ decay rate is then Allowed values for the Higgs portal coupling λ can be inferred from a global fit to the Higgs data as has been done for the portal coupling of the (3, 1, −1/3) state. (See figure 8.) We do not attempt to repeat the same procedure for the state (3, 2, 7/6) since we do not expect a drastic change in the allowed range of λ. The Wilson coefficients for τ → µγ are again proportional to the couplings responsible for h → τ µ: t /m 2 ∆ 2 and agree with the formulas presented in [55]. Finally, the corresponding branching fraction is given by As in the case of the ∆ 1 leptoquark, contributions to the muon anomalous magnetic moment (g − 2) µ are proportional to y R tµ y L µt and can be always suppressed. On the other hand, contributions of ∆ 2 enter the Z → bb via positive contributions to δg b L that can actually help in easing the tension with the experimental data [58]. Most importantly, these observables do not constrain large possible ∆ 1 effects in h → τ µ .
In the left-hand (right-hand) side panel of figure 7 we show the B(h → τ µ) (B(τ → µγ)) dependence on the leptoquark Yukawa couplings for the ∆ 2 LQ case with the Higgs portal coupling turned off and taking m ∆ 2 = 650 GeV. Also in this leptoquark scenario the bound on B(τ → µγ) excludes sizable B(h → τ µ) due to the strict correlation between the two observables. See the orange stripe in figure 3, where the portal coupling is restricted to |λ| < 1.

Fine-tuning solution
In this section we give an example of a phenomenologically viable model with fine-tuning. We add to the SM a scalar leptoquark ∆ and a vector-like top quark partner T and turn on the following couplings Here, the top partner T has a Dirac mass term m TT T and we assume no mixing with the top quark. Thus, T will have no impact on h → τ µ if the Higgs portal coupling is turned off. The h → τ µ excess is explained by the top quark contribution only and the left-hand side panel of figure 7 is applicable in this case.
On the other hand, T can cancel the top quark contribution to τ → µγ. We calculate the decay rate where h 2 (x) is defined in eq. (5.20).
JHEP06(2015)108 Figure 9. Fine-tuning τ → µγ to zero while keeping h → τ µ at best fit point within the LQ model with a vector-like top quark partner.
As an example we take m ∆ = 650 GeV and y L 33 y R * 32 = 0.8, which gives the best fit value for B(h → τ µ). In figure 9 we show how to avoid the τ → µγ constraint with an appropriate choice of x R * 32 x L 33 and m T . The scenario can be realized in both LQ models. In passing we note that if m T > m ∆ , the T → ∆ decays can produce spectacular signatures at the LHC, while at the same time avoiding traditional searches for vector-like top quark partners, especially if T does not mix with the top.

Conclusions
Prompted by the recent experimental hint of h → τ µ events by the CMS Collaboration, we have carefully examined the implications of LFV Higgs decays at the percent level on possible extensions of the SM. In particular we have shown how a tentative B(h → τ µ) signal can be combined with other Higgs measurements to yield a robust lower bound on the effective LFV Higgs Yukawa couplings to taus and muons. Then we have reexamined the connection between LFV Higgs decays and LFV radiative decays of charged leptons, and demonstrated using EFT methods that the current CMS hint implies τ → µγ at rates, which could be observable at the Belle II experiment. In explicit models, the τ → µγ constraint is generically much more severe. In fact, an eventual observation of h → τ µ at the LHC together with existing indirect constraints would already single out an extended SM scalar sector as a required ingredient in any natural explanation, the minimal example being THDM of type III. We have also examined purely fermionic SM extensions and models where h → τ µ is generated at loop level, only to show that without the introduction of extra Higgs doublets, reconciling all existing indirect constraints with percent level B(h → τ µ), when at all possible, requires a high degree of fine-tuning. Finally, we have shown how a positive signal of h → τ µ can be combined with experimental searches for µ → eγ decays and µ − e conversions in nuclei to yield robust bounds on B(h → τ e). In particular, considering only low energy Higgs EFT effects, the two LFV Higgs decay rates could still JHEP06(2015)108 be comparable. On the other hand, the THDM III cannot accommodate both respective branching ratios at the percent level. Currently planned improvements in experimental searches for µ − e LFV processes will be able to probe the product B(h → τ µ)B(h → τ e) at the 10 −7 level in generic EFT and to order 10 −12 or better within the THDM III.

Acknowledgments
We acknowledge insightful discussions with Vincenzo Chiochia. This work was supported in part by the Slovenian Research Agency. I.D. acknowledges the SNSF support through the SCOPES project No. IZ74Z0 137346 and Croatian Science Foundation support under the project 7118. I.N. acknowledges support in part by the Bundesministerium für Bildung und Forschung.

A Higgs fit
In the following we explain the method used to fit the LHC Higgs data. Ideally, we want the measurements to be reported in terms of the signal strengths normalized to the SM predictions where index i represents a Higgs decay mode, while k denotes a production channel. Such observables could then easily be expressed in terms of new physics parameters. However, the actual experimental categories are never pure and contain events from different production mechanisms. Furthermore, in order to fully reconstruct the total likelihood function, it is necessary to know all the correlations among the different categories, which are available only to the experimental collaborations. It has been argued recently in ref. [59] that the existing measurements for a given decay channel should be presented in terms of two-dimensional likelihoods, in which gluon-gluon fusion (ggF) and associated production with a top pair (ttH) are combined as one signal (µ (ggF +ttH) ), while vector boson fusion (VBF) and associated production with a gauge boson (VH) as another, (µ (V BF +V H) ). This decomposition is particularly useful since, on one hand, ttH is sub-dominant with respect to ggF and poorly explored in the present data set, while VBF and VH receive common corrections in wide class of NP models obeying the custodial invariance. Furthermore, the correlations are automatically provided by the experimental collaborations for a given decay channel.
Following the recommendation, ATLAS and CMS have combined different search categories for a given decay mode to provide separation into production mechanisms. The results are usually presented in 2D plots in (µ (ggF +ttH) , µ (V BF +V H) ) plane. In this case, we parametrize the likelihood with where the sum goes over the decay channels and the covariance matrices are given by We obtain the best-fit values (μ), variances (σ) and correlations (ρ) from the plots provided by the experimental collaborations. These are presented in table 1 and table 2.
If the separation into production modes is not provided, we use the signal strength measurements from existing search categories. These in general target certain production mechanism, which, however, does not imply 100% purity. Inclusive categories are dominated by ggF (90%), while VBF-tagged categories can have 20% to 50% contamination from ggF. VH-and ttH-tagged categories are assumed to be pure. In this case, we write σ A→h σ SM where ξ i represent contributions of the specified production mechanisms to the given category. We do not consider correlations here and add each measurement to total χ 2 separately, (A.5) We approximate the likelihood for the total Higgs decay width measurement [74] by h → invisible VBF+VH 0.14 ± 0.22 [71] h → Zγ inclusive 0.0 ± 4.8 [72] h → µµ inclusive 2.9 ± 2.8 [73] Table 2. CMS Higgs data as used in our fitting procedures. The separation into production mechanisms for a given decay channel is used if provided.
In table 1 and table 2 we finally summarize all the available LHC Higgs data used in our analyses.
The total χ 2 function is given by the sum of all the contributions, namely, χ 2 = χ 2 1 + χ 2 2 + χ 2 Γ . In order to confront the new physics model to data, we express all signal strengths (µ) and total Higgs decay width in terms of model parameters and minimize the χ 2 to find the best fit point. The best fit regions are defined by the appropriate cumulative distribution functions, namely, 68.2% (1 σ) best-fit region for one-(two-) parameter fit satisfies χ 2 − χ 2 min < 1.0 (2.3), while 95.5% (2 σ) best-fit region satisfies 1.0 (2.3) < χ 2 − χ 2 min < 4.0 (6.2). The other parameters in the likelihood are treated as nuisance parameters. The limits quoted in eq. (1.4) are obtained after allowing only y τ µ , y µτ to be free parameter while fixing other Higgs couplings to their SM values. On the other hand, in the limits shown in eq. (1.5) we further allow for arbitrary values of the usual SM tree level Higgs boson couplings; κ t , κ τ , κ b , κ V and κ c where κ V = κ W = κ Z as well as arbitrary new contributions to loop induced Higgs couplings; κ g and κ γ .
Here C H 0 k W W measure the coupling strength of neutral Higgses to W W relative to the SM coupling hW W , i.e., deviation from the vertex igm W g µν . Their values are C HW W = cos(β − α), C hW W = sin(β − α). The coefficients shown above correspond to topologies that contain H k γγ coupling, mediated either by t or W loop. The virtual photon can be replaced by a Z boson, however such contributions are subdominant. Loop integrals f , g, and h were calculated in ref. [37], whereas their arguments are defined as z tk = m 2 t /m 2