A Three-loop Neutrino Model with Leptoquark Triplet Scalars

We propose a three-loop neutrino mass model with a few leptoquark scalars in $SU(2)_L$-triplet form, through which we can explain the anomaly of $B \to K^{(*)} \mu^+ \mu^-$, a sizable muon $g-2$ and a bosonic dark matter candidate, and at the same time satisfying all the constraints from lepton flavor violations. We perform global numerical analyses and show the allowed regions, in which we find somewhat restricted parameter space, such as the mass of dark matter candidate and various components of the Yukawa couplings in the model.


I. INTRODUCTION
Recently, there was an 2.6σ anomaly in lepton-universality violation in the ratio R K ≡ B(B → Kµµ)/B(B → Kee) = 0.745 +0.090 −0.074 ± 0.036 by the LHCb Collaboration [1]. In addition, sizable deviations were observed in angular distributions of B → K * µµ [2]. The results can be interpreted by a large negative contribution to the Wilson coefficient C 9 of the semileptonic operator O 9 , and also contributions to other Wilson coefficients, in particular to C ′ 9 [3][4][5][6]. The discrepancy between the theoretical prediction and experimental value on the muon anomalous magnetic dipole moment has been a long-standing problem, which stands at 3.6σ level with the deviation from the SM prediction at [7]. ∆a µ = a exp µ − a SM µ = 288(63)(49) × 10 −11 .
If one insists on fulfilling the muon g − 2 within 1σ − 2σ of the experimental value in any models, it puts a strong constraint on the parameter space. For example, it requires a relatively light spectrum in the supersymmetric particles in the MSSM in order to bring the prediction to be within 1σ − 2σ of the experimental value. A number of leptoquark models have been proposed to solve the B → K ( * ) µµ anomaly, but however it is very hard to satisfy simultaneously the muon g − 2: see for example Ref. [8].
In this work, we propose a three-loop neutrino mass model with a few leptoquark scalars in SU(2) L -triplet form. We attempt to use the model to explain the anomaly of B ( * ) → Kµ + µ − , to achieve a sizable muon g − 2, and to provide a bosonic dark matter candidate, and at the same time satisfying all the constraints from lepton flavor violations. The concrete model is based on the SM symmetry and a Z 2 symmetry as SU(3) C × SU(2) L × U(1) Y × Z 2 . The model consists of the SM fields, 3 additional leptoquark triplet fields ∆ a 1,2,3 , and one colorless doublet scalar field η. These fields are assigned different Z 2 parities and hypercharges in such a way that each of the Yukawa-type couplings contributes to either neutrino mass, B → K ( * ) µµ anomaly, muon g − 2, or the dark matter interactions. In this way, although the model contains more parameter, it can however explain all the above anomalies. The achievements of the model can be summarized in the following.
coupling terms f, g, h in three-loop diagrams 1 .

The
Yukawa coupling term f can give useful contributions to the Wilson coefficients C 9,10 in such a way that it can explain successfully the B → K ( * ) µµ anomaly.
3. The muon g − 2 receives a large contribution from the Yukawa coupling term r. With some adjustment of the parameters a level of 10 −9 is possible.
4. It provides a dark matter (DM) candidate η R , the real part of the neutral component of the η field with correct relic density. 5. The model can satisfy all the existing constraints from the lepton-flavor violations (LFVs), meson mixings, and rare B decays.
This paper is organized as follows. In Sec. II, we describe the neutrino mass matrix and the solution to the anomaly in b → sµμ. In Sec. III, we discuss various constraints of the model, including lepton-flavor violations, FCNC's, oblique parameters, and dark matter. In Sec. IV, we present the numerical analysis and allowed parameter space, followed by the discussion on collider phenomenology. Sec. IV is devoted for conclusions and discussion.

II. THE MODEL
In this section, we describe the model setup, derive the formulas for the active neutrino mass matrix, and calculate the contributions to b → sµμ.
where the superscript index a = 1 − 3 represents the color.

A. Model setup
We show all the field contents and their charge assignments in Table I for the fermionic sector and in Table II for the bosonic sector. 2 Under this framework, the relevant part of the renormalizable Lagrangian and Higgs potential related to the neutrino masses are given where we have defined L ′ ≡ [N, E] T , σ 2 is the second Pauli matrix and we have abbreviated the trivial terms for the Higgs potential. The scalar fields can be parameterized as where the subscript next to the each field represents the electric charge of the field, v = 246 GeV, and Φ is written in the form after the Goldstone fields are aboserbed as the longitudinal components of W and Z bosons. Notice here that each of the components of ∆ 3 and η is 2 The same contents of the field are found in the systematic analysis in the last part of Table 3 of ref. [12].
in mass eigenstate, since there are no mixing terms that are assured by the Z 2 and U(1) Y symmetries. On the other hand, components of ∆ 1 and ∆ 2 can mix via Φ * Φ * ∆ 1 ∆ 2 term.
In the following analysis, we ignore such mixing effects assuming the relevant coupling is small.
Oblique parameters: Each of the mass components among ∆ i is strongly restricted by the oblique parameters. In order to evade such a strong constraint, we simply assume that each of the components should be of the same mass [13]. Thus, we define m ∆ i as the mass for the components of ∆ i . On the other hand, each component of η cannot have the same mass, because the neutrino mass is proportional to the mass difference between the components of η, as you shall see later. Hence, we consider the oblique parameter constraints on η, which are characterized by ∆T and ∆S. Their formulae are given by [14] ∆T where α em ≈ 1/137 is the fine structure constant, and The experimental bounds are given by [7] (0.05 − 0.09) ≤ ∆S ≤ (0.05 + 0.09), (0.08 − 0.07) ≤ ∆T ≤ (0.08 + 0.07). (II.5) We consider these constraints in the numerical analysis.
Active neutrino mass matrix: The neutrino mass matrix is induced at three-loop level as shown in Fig. 1, and its formula is generally given by where we used the shorthand notation m R/I ≡ m η R/I , and define M Max ≡ , and the three-loop function F III is given in the Appendix. Here we adopt an assumption M max = m ∆ 3 , and require 1TeV m ∆ i (which suggests x d , x u ≈ 0), which is required by the direct bound on leptoquarks [13]. In this case, the neutrino mass matrix can be simplified as where we have abbreviated the symbol of summation and the argument of F III . Then we derive the Yukawa coupling in terms of the experimental values and the parameters by introducing an arbitrary anti-symmetric matrix with complex values A [15], that is A T +A = 0, as follows: where we shall adopt the former formula in the numerical analysis below, and we define Here we assume one massless neutrino (with normal ordering) in the numerical analysis below.
On the term f : The new physics contributions to account for the B → K ( * ) µµ anomaly [2] can be interpreted as the shifts in the Wilson coefficients C 9,10 . In our model, the relevant Wilson coefficients can be calculated as follows [13]: We can then compare them to the best-fit values of C 9,10 from a global analysis based on the LHCb data in Ref. [3] as Here we also have to work within the −0.75 C 9 −0.35 in order to satisfy the the LHCb .036, which shows a 2.6σ deviation from the SM prediction [13]. Notice here that various constraints arising from the f term include B d/s → ℓ + ℓ − , ℓ i → ℓ j γ, which were given in Refs. [13] and [8], and we consider these constraints in the current numerical analysis. Although the muon g − 2 is also induced from this term, the typical order is 10 −12 ∼ 10 −13 with a negative sign [13]. Thus, we neglect this contribution to the muon g − 2.
On the terms g and h: The main constraint on g and h comes from B(b → sγ). The partial decay width for b → sγ is given by 1/3 , a] 2 , (II.14) then the branching ratio and its experimental bound [16] are given by where Γ tot. ≈ 4.02 × 10 −13 GeV is the total decay width of the bottom quark. In our numerical analysis, we consider this constraint only for the g and h terms.
On the term r: This term is very important in our model because it can induce a large contribution to the muon g − 2 and explain the relic density of dark matter (DM) if we assume the η R to be the DM candidate. First of all, let us consider the LFVs processes, ℓ a → ℓ b γ, via one-loop diagrams. The branching ratio is given by where m a(b) is the mass for the charged-lepton eigenstate, C ab ≈ (1, 0.1784, 0.1736) for (a, b) = (2, 1), (3,1), (3,2), and a L(R) is simply given by Q −Q mixing: The forms of K 0 −K 0 , B 0 d −B 0 d , and D 0 −D 0 mixings are, respectively, given by    Dark Matter: Here we identify η R as the DM candidate, and denote its mass by m R ≡ M X .

Direct detection:
We have a Higgs portal contribution to the DM-nucleon scattering process, 3 Since we assume that one of the neutrino masses to be zero with normal ordering that leads to the first column in g to be almost zero, i.e., (g) 11,12,13 ≈ 0, and so these constraints can easily be evaded.
Once we apply this bound on our model, we obtain λ 3 + λ 4 + 2λ 5 2 × 10 −3 . Hence we assume that all the Higgs couplings are small enough to satisfy the constraint, and we neglect DM annihilation modes via Higgs portal in estimating the relic density below. Notice here that photon and Z boson fields transform as V µ → −V µ under charge-parity (CP ) conjugation, while X is CP -even. Thus X − X − γ(Z) couplings are not allowed because they violates the CP invariance.
Relic density: We consider parameter region in which the DM annihilation cross section is d-wave dominant and the dark matter particles annihilate into a pair of charged-leptons, via the process η R η R → ℓ ilj with an E a exchange. Notice that there exist annihilation modes such as η R η R → W + W − /2Z arising from the kinetic term. These modes require a DM mass heavier than at least 500 GeV [22] in order to obtain the correct relic density where coannihilation processes should be included. This case is, however, not in favor of explaining the muon g − 2 anomaly. Thus we assume that M X 80 GeV and η R η R → W + W − /2Z processes are not kinematically allowed. 4 Then the relic density is simply given by , d eff ≈

III. NUMERICAL ANALYSIS
As a first step, we perform the numerical analysis on the r term since this term is independent of the other parameters. We prepare 25 million random sampling points for the relevant input parameters as follows: , r ij ≡ (±1) × 10 r ′ ij and the lower mass bound 1.1 × M X is expected to forbid the coannihilation processes. Under such parameter ranges, we have found 246 allowed points, which are shown in Fig. 2 satisfying all the constraints including LFVs, oblique parameters, and invisible decay of Z boson, as discussed before.
The left panel shows the allowed region to be where the lower bound of η ± comes from the LEP experiment [24,25]. The right panel shows the correlation of r 12 versus ∆a µ , where r 12 is the most relevant parameter to obtain the sizable muon g − 2 and correct relic density of the DM. It suggests that a rather small r 12 is possible to achieve the range 10 −9 ≤ ∆a µ ≤ 2 × 10 −9 , but we need a large r 12 for the range 2 × 10 −9 < ∆a µ ≤ 2 − 4 × 10 −9 . 5 In the next step, we attempt to find the parameter space points that can also solve the B → K ( * ) µµ anomaly by contributing to the C 9 (10) , and at the same time satisfy all the constraints of the LFVs and FCNCs. Since the number of parameters is getting more and more, we are content with a benchmark point that we obtained in the first step. We prepare the benchmark point for the masses to fix the three-loop neutrino function F III 6 as follows:  (III.6) Also, we fix −λ 0 λ ′ 0 v 2 /2 ≈ 1[GeV 2 ] for simplicity. Then, we prepare 0.1 billion random sampling points for the following relevant input parameters: where we define f (h) ij ≡ (±1) × 10 f ′ (h ′ ) ij . In these parameter ranges, we have found 870 allowed points shown in Fig. 3, which satisfy all the constraints as discussed before. If we focus on the best-fit value of C 9 , these figures show that the Yukawa couplings f 1j (j = 1 − 3) are very restricted as |f 11 | 0.1, |f 12 | 0.02, |f 13 | 0.2, (III. 8) where these coupling may affect the collider physics. is the one at 1σ range. Notice here that we have taken the range C 9 ∈ [−0.75, −0.50] to satisfy R K at 1σ confidential level [3].
Collider Phenomenology: There are two types of new particles in this model other than the SM particles: leptoquarks ∆ a 1,2,3 and a set of scalar bosons η ±,0 resulting from an isodoublet scalar field.
Note that ∆ a 1,2 are assigned with Z 2 = −1 while ∆ a 3 are assigned with Z 2 = +1. All ∆ 1,2,3 can be pair produced at hadronic colliders and being searched at the LHC. The direct search bound is roughly 1 TeV [26]. On the other hand, only ∆ a 3 can participate in the 4fermion contact interaction, because of the Z 2 parity. The bound from the 4-fermion contact interaction was worked out in Ref. [8,13] that the bound is currently inferior to the direct search bound of about 1 TeV. Therefore, we shall use 1 TeV as the current bound on the ∆ 1,2,3 bosons.
The isodoublet field η gives rise to a pair of charged bosons η ± , a scalar η R , and a pseudoscalar η I . The charged bosons η ± can run in the triangular loop vertex of Hγγ.
Nevertheless, we can suppress such effects by choosing the λ 5 term in Eq. (II.1) very small.
As we have explained above such a term is small to avoid the conflict of the direct detection bound. Although the interaction between the η field and the SM Higgs field is suppressed for the above reasons, the η can still interact through the kinetic term as it has SU(2) L and U(1) Y interactions. We expect some typical interactions with the gauge bosons: An interesting signature would be Drell-Yan type production of η ± η R via a virtual W . The η ± decays into multi-leptons and η R via the virtual L ′ fields. Therefore, the final state consists of multi-charged-leptons and missing energies. Similarly, in the process pp → Z * → η R η I , the η I would decay into η R eventually with a number of very soft leptons, which may not be detectable. Therefore, the best would be the one produced via virtual W .

IV. CONCLUSIONS
We have investigated a three-loop neutrino mass model with some leptoquark scalars with SU(2) L -triplet, in which we have explained the anomaly in B → K ( * ) µ + µ − , sizable muon g − 2, bosonic dark matter, satisfying all the constraints such as LFVs, FCNCs, invisible decay, and so on. Then we have performed the global numerical analysis and shown the allowed region, in which we have found restricted parameter space, e.g., We find that ∼ O(100) GeV inert doublet scalar is preferred to obtain sizable muon g − 2.
Thus these light inert scalars could be tested by collider experiments such as LHC in which these scalars are produced via electroweak processes. The promising signature of our model comes from the process pp → W ± → η ± η R which provides signals of multi-leptons plus missing transverse momentum.

Appendix A: Loop function
Here we show the explicit form of three-loop function F III , which is given by