On the stability of the electroweak vacuum in the presence of low-scale seesaw models

The scale of neutrino masses and the Planck scale are separated by more than twenty-seven order of magnitudes. However, they can be linked by imposing the stability of the electroweak (EW) vacuum. The crucial ingredient is provided by the generation of neutrino masses via a seesaw mechanism triggered by Yukawa interactions between the standard model (SM) Higgs and lepton doublets and additional heavy right-handed neutrinos. These neutrinos participate to the renormalization group (RG) running of the dimensionless SM couplings, affecting their high-energy behavior. The Higgs quartic coupling is dragged towards negative values, thus altering the stability of the EW vacuum. In the usual type-I seesaw model, this effect is too small to be a threat since, in order to comply with low-energy neutrino data, one is forced to consider either too small Yukawa couplings or too heavy right-handed neutrinos. In this paper we explore this general idea in the context of low-scale seesaw models. These models are characterized by sizable Yukawa couplings and right-handed neutrinos with mass of the order of the EW scale, thus maximizing their impact on the RG flow. As a general result, we find that Yukawa couplings such that ${\rm Tr}(Y_{\nu}^{\dag}Y_{\nu}) \gtrsim 0.4$ are excluded. We discuss the impact of this bound on several observables, with a special focus on the lepton flavor violating process $\mu\to e\gamma$ and the neutrino-less double beta decay.


Introduction
The discovery of the Higgs boson at the LHC [1, 2] represents a cornerstone in particle physics, finally achieved after a nearly half-century quest. Paradoxically though it may seem, this astonishing discovery did not clarify the true origin of the electroweak (EW) symmetry breaking: after LHC run I, any clue of new physics beyond the Standard Model (SM) is still missing.
To crown it all, the measured value of the Higgs boson mass, M h = 125.09 ± 0.21 stat. ± 0.11 syst. GeV according to the latest combination of both ATLAS and CMS experiments [3], is consistent, in spite of any naturalness argument, even if the SM is extrapolated up to the Planck scale. This remarkable conclusion is based on the fact that the value of the Higgs boson mass is caught in the vice of two solid theoretical constraints if one tries to extend the SM up to arbitrarily high energies by means of its Renormalization Group Equations (RGEs). On the one hand, for large values of M h , the Higgs quartic coupling λ is driven towards greater and greater values, eventually becoming non-perturbative. On the other one, for small values of M h , λ is driven towards smaller values, eventually becoming negative. The occurrence of the latter condition may have dramatic consequences since it implies that the EW vacuum is only a local minimum, and there exists a global minimum at higher energy in the Higgs potential [4][5][6][7][8][9][10][11][12][13][14][15][16][17]. For a sufficiently old Universe, nothing would prevent the vacuum expectation value (vev) of the Higgs field to be at the global minimum, as a mere consequence of a tunneling transition between the two vacua [18][19][20]. Fortunately, Nature spared us from this unspeakable catastrophe: in the SM with M h ≃ 125 GeV λ does become negative at high energies, but its absolute value remains far too small to trigger the tunneling. In other words, in the SM the EW vacuum is metastable [21][22][23].
The metastability of the EW vacuum holds under the assumption that the SM is the full theory all the way up to the Planck scale. This argument poses a legitimate question: what is the role of new physics beyond the SM with respect to the stability of the EW vacuum [24][25][26][27][28][29][30][31][32][33][34][35]? First and foremost, to address this question we have to carefully specify which kind of new physics we are referring to. Undeniably, the most compelling case is provided by the observation of non-zero neutrino masses and oscillations. This fact cries out for the inclusion in the SM of a mechanism responsible for the generation of neutrino masses [36]. A very economical way is provided by the so-called type-I seesaw mechanism [37][38][39][40]. In a nutshell, the SM is enlarged with the addition of three heavy right-handed neutrinos 1 with a Majorana mass M R breaking lepton number, coupled to the SM Higgs doublet H and leptonic doublets L via a Yukawa coupling Y ν . Integrating out the heavy fields, the effective light neutrino mass m ν ≈ v 2 Y ν M R is generated after EW symmetry breaking, where v ≃ 246 GeV is the EW vev.
In principle, this additional Yukawa interaction pushes the SM towards the instability region, since, at least at one-loop, it exactly mimics the effect of the top Yukawa coupling on the running of λ: it introduces a negative contribution which drags the Higgs quartic coupling towards negative values [22,[42][43][44]. However, in order to have a sizable effect, the Yukawa coupling Y ν must be of order one. This is a necessary but not a sufficient requirement. Indeed, it implies M R ∼ O (10 15 ) GeV if we assume m ν ≃ 0.1 eV for the light neutrino mass scale. This means that the right-handed neutrinos actively participate to the running of λ only for values of the RG scale larger than M R . Or, to put it another way, there is not enough time, in terms of RG evolution, to sizably alter the SM picture. This result does not change trying to lower the mass scale M R since in this case, in order to reproduce the correct order of magnitude for the mass scale of light neutrinos, one is forced to consider Y ν ∼ O(10 −5 ). Of course, this value is too small to modify the running of λ.
However, this is anything but the end of the story. Despite its simplicity, the minimal type-I seesaw suffers from a penalizing phenomenological issue: it is not testable. As discussed before, it introduces either unaccessible large mass scale or extremely tiny Yukawa couplings. The necessity to overcome this unpleasant problem led to the introduction of extended seesaw models with both TeV-scale right-handed neutrinos and sizable Yukawa couplings [45][46][47][48][49][50]. The key feature of these models is that lepton number is softly broken via the introduction of extra singlet fermions in addition to the usual right-handed neutrinos. The mass M R can be brought down to the EW scale without neither causing problem with low-energy neutrino phenomenology nor lowering the Yukawa coupling Y ν . As a byproduct of this construction, a very rich low-energy phenomenology emerges. Potentially interesting signals include, for instance, lepton flavor violating radiative decays, deviations from EW precision observables, and production at colliders of heavy Majorana fermions.
As already clear from these introductory comments, low-scale seesaw models may alter the metastability of the EW vacuum since they feature, at the same time, sizable Yukawa couplings and relatively low mass thresholds [51]. In this work we study in detail the stability of the EW vacuum considering the so-called inverse seesaw model (ISS) [45,46] as a reference model with low-scale seesaw. We include in our analysis the constraints coming from low-energy neutrino phenomenology with the aim to provide a complete and realistic description of the physics involved.
Schematically, this paper is organized as follows. In section 2 we discuss the most relevant features of the ISS model, with a special focus on the phenomenological constraints included in our analysis. In section 3 we study the stability of the EW vacuum in the ISS model, using the results of section 2. In section 4 we present our results, and finally we conclude in section 5. In appendix A we generalize our results to the linear and double seesaw models.

The inverse seesaw model
In this section we start our analysis describing the most relevant prerogatives of the ISS model. In section 2.1 we introduce the Lagrangian, the mass matrix and the corresponding diagonalization. In section 2.2 we discuss the constraints coming from low-energy neutrino phenomenology included in our analysis while in section 2.3 we present our strategy for the numerical analysis.

The model
In the ISS model the SM field content is extended to incorporate n R right-handed neutrinos N i R and n S singlet fermionic fields S i . The Lagrangian in the ISS model is given by where L ≡ (L e , L µ , L τ ) T represents the left-handed lepton doublets with the usual contents L l=e,µ,τ = (ν lL , l L ) T while H is the Higgs field, withH ≡ iσ 2 H * . Y ν is the n R × 3 Yukawa matrix mediating the interactions between the SM leptons and the right-handed neutrinos while M R and µ S are, respectively, n R × n S and n S × n S mass matrices. Both right-handed neutrinos and singlet fermions have lepton number L = 1; consequently, the mass term S C µ S S violates lepton number for two units.
Introducing the left-handed basis N L ≡ (ν L , N C R , S) T we have, after EW symmetry breaking, the following mass matrix The mass matrix in eq. (2.2) can be diagonalized by means of the following unitary transformation where M ∆ is the diagonal matrix referred to the mass eigenstates N ′ L . The first three eigenstates correspond to the standard light active neutrinos while the remaining n R + n S states are additional heavy sterile neutrinos.
Following the standard seesaw calculation [52,53] and assuming the hierarchy M R ≫ m D ≫ µ S , it is possible to extract the effective light neutrino mass matrix In the following we consider the case n R = n S = 3. Moreover, without loss of generality, we take M R to be real and diagonal. We also work in a basis in which the mass matrix of charged SM leptons is diagonal. Within this framework, in the µ S → 0 limit the three light neutrinos are massless while the six heavy neutrinos can be recast into three pairs of Majorana particles with three (double degenerate) masses M Ri , i = 1, 2, 3. On a general ground, from eq. (2.4) it follows that the order of magnitude of the light neutrino mass is . Assuming µ S ∼ O(1) keV, the model can accommodate sub-eV light neutrino masses with Y ν ∼ O(1) couplings and M R ∼ O(1−10) TeV seesaw scale. These order of magnitude estimates lie at the hearth of the ISS scenario. Small values of µ S are expected by virtue of the 't Hooft naturalness criterium [54], since the limit µ S → 0 increases the symmetry of the theory. Interestingly, the keV scale nicely fits the typical mass scale characterizing warm dark matter; in the ISS models with n R ≠ n S the spectrum -in addition to light and heavy neutrinos -also contains intermediate states with keV mass that are valuable warm dark matter candidates [55]. The estimates Y ν ∼ O (1) and M R ∼ O(1 − 10) TeV represent the most relevant phenomenological properties of the ISS model since they allow for sizable (and, in principle, measurable) mixing effects between light active and heavy sterile neutrinos. This issue is particularly striking if compared with the typical high-scale characterizing the minimal type-I seesaw [37][38][39][40] -that is M R ∼ O(10 15 ) GeV for order one Yukawa couplingsin which mixing effects, typically of order O(m 2 D M 2 R ), are negligible. As we shall explain in detail in section 3, the occurrence of both the peculiar conditions Y ν ∼ O(1) and M R ∼ O(1 − 10) TeV is of fundamental importance to determine the stability of the EW vacuum in the context of the ISS model.
Of particular relevance for many phenomenological applications is the generalized Casas-Ibarra parametrization [56]  where c ij ≡ cos θ ij and s ij ≡ sin θ ij . In addition to the Dirac CP violation phase δ CP there are also two Majorana CP violation phases (not shown in eq. (2.6)). The latter are physical only if light neutrinos are Majorana particles, otherwise they can be always rotated away from the Lagrangian in the mass basis. In the following we omit the three phases of the PMNS matrix, since the Majorana phases are completely unknown and there are only preliminary hints about a non-zero Dirac phase.
In this work, we use the latest results of the νfit group [64]. As customary, we define ∆m 2 ij ≡ m 2 νi − m 2 νj , and we adopt the convention that results for the mass squared differences are reported with respect to the one with the largest absolute value.
○ Neutrino mass squared differences 3-σ C.L. ranges on the mass squared differences

○ Unitarity
In the basis where the charged lepton mass matrix is diagonal, the leptonic mixing matrix in the ISS model is given by the rectangular 3 × 9 sub-matrix corresponding to the first three rows of the matrix U defined in eq. (2.3), with the 3 × 3 block corresponding to the (non-unitary)Ũ PMNS . Bounds on the nonunitarity of the matrixŨ PMNS were derived in [65][66][67] using an effective field theory approach. These bounds can be recast as follows 2

○ Additional constraints
The absolute values of neutrino masses m νi are unknown. Cosmology sets the most stringent upper bounds using data from the Cosmic Microwave Background (CMB) radiation, supernovae and galaxy clustering. Assuming the validity of the ΛCDM model [68], the Planck collaboration placed the upper bound ∑ i m νi < 0.66 eV at 95% C.L. [68]; this bound becomes even more stringent adding data on the Baryon Acoustic Oscillation, ∑ i m νi < 0.23 eV at 95% C.L. [68]. In our analysis we scan over the interval 10 −4 eV ⩽ m ν1 ⩽ 10 −1 eV for the mass of the lightest neutrino.
In this paper, in order to investigate the stability of the EW vacuum in a region of the parameter space of particular interest for present and future experimental prospects, we focus on the lepton flavor violating process µ → eγ and the neutrinoless double beta decay (0ν2β hereafter).
As far as the radiative µ → eγ decay is concerned, the rate induced by the presence of sterile neutrinos is given by [85,86] where α is the electromagnetic fine structure constant and M W the W mass. The loop function is given by The present experimental upper bound, reported by the MEG collaboration, is Br(µ → eγ) < 5.7 × 10 −13 at 90% C.L. [87]. The amplitude of the 0ν2β process is proportional to the so-called effective neutrino mass, m νe eff . Current experiments (among others, GERDA [88], EXO-200 [89,90], and KamLAND-ZEN [91]), put an upper limit in the range m νe eff ≲ 140 − 700 MeV. In the presence of sterile state the effective neutrino mass is given by [92] where p 2 ≃ −(125 MeV) 2 is the momentum of the virtual neutrino. Notice that, since we are considering the regime m 2 νi=4,...,9 ≫ p 2 , heavy neutrinos decouple in eq. (2.13) and the dominant contribution to m νe eff comes from the light active neutrinos.

Strategy and first numerical results
We perform a scan over the parameter space of the model, and our procedure goes as follows. First, we randomly generate i) the light neutrino masses m νi=1,2,3 and the leptonic mixing angles θ 12 , θ 23 , θ 13 according to the corresponding 1-σ intervals allowed by the analysis of present experimental data, ii) the entries of the matrices M R and µ S in the intervals 10 2 GeV ⩽ M Ri ⩽ 10 2 TeV, 10 −1 keV ⩽ (µ S ) ij ⩽ 10 2 keV and iii) the complex angles defining the arbitrary matrix R in the interval [0, 2π]. Second, we reconstruct the full Yukawa matrix Y ν using the generalized Casas-Ibarra parametrization in eq. (2.5). Finally, plugging back M R , µ S and Y ν into eq. (2.2), we diagonalize the mass matrix M in order to find the full 9 × 9 mixing matrix U . The phases of the mixing matrix are fixed using eq. (2.3), by means of the condition m νi ⩾ 0 for all i. As a consistency test, for each point of the scan we check that the mass matrix M, randomly generated as discussed above, correctly reproduces after diagonalization light neutrino masses and mixing angles in agreement with the bounds discussed in section 2.2. Equipped by these results we can easily compute the branching ratio Br(µ → eγ) in eq. (2.12) and the effective neutrino mass in eq. (2.13). We show our results in fig. 1 for the normal ordering and in fig. 2 for the inverted ordering. In the upper (lower) panels of both figures we show the branching ratio Br(µ → eγ) (the effective neutrino mass m νe eff ) as a function of the lightest heavy neutrino mass m ν4 (plot on the left) and the trace of the Yukawa couplings Tr(Y † ν Y ν ) (plot on the right). Few comments are in order.

Normal and inverted ordering produce very similar distributions considering
the radiative decay Br(µ → eγ). This is caused by the well-known fact that the contribution of light active neutrinos is strongly suppressed by the extremely small value of light neutrino masses. In the ISS model a non-zero contribution to Br(µ → eγ) is entirely generated by the additional heavy neutrinos, and controlled by the mixings U µi , U ei ∼ m D M R (see eq. (3.11) in section 3.2). We notice that in our scan we can obtain a signal close to the present experimental bound even considering m ν4 as large as 10 TeV and Tr(Y † ν Y ν ) as small as  [90]. Details about the numerical scan are given in the text. All points comply with the bounds discussed in section 2.2. Br(µ → eγ) ⩾ 10 −13 . Points with large Yukawa couplings (e.g. Tr(Y † ν Y ν ) ≳ 0.5) and sizable µ → eγ rate are generated in the whole interval of analyzed masses for the right-handed neutrinos.
2. Normal and inverted ordering produce completely different distributions considering the effective neutrino mass m νe eff . In this case the contributions of additional heavy neutrinos decouple since their masses are much larger than the typical momentum scale p 2 ≃ −(125 MeV) 2 . Therefore, in our numerical scan the ISS model resembles the typical scenario with only three light active neutrinos. The situation is well represented by the right panel of fig. 3 where we show the effective neutrino mass for the normal and inverted ordering as a function of minimal neutrino mass. The normal ordering is suppressed since the largest neutrino mass is multiplied by the small value of s 13 . However, in both cases the effective neutrino mass is close to the future sensitivity of the EXO-200 experiment. Distribution of our numerical scan in the plane [Tr(Y † ν Y ν ), m ν4 ] considering the normal ordering (the inverted ordering gives an analogous result). All points comply with the bounds discussed in section 2.2. Points with Br(µ → eγ) ⩾ 10 −13 are marked in dark cyan. On the qualitative level, the plot shows that Yukawas Tr(Y † ν Y ν ) ≳ 0.5 arise in the whole range of analyzed masses for the right-handed neutrinos. Right panel. Effective neutrino mass as a function of the lightest neutrino mass considering both normal and inverted ordering. Our scan correctly reproduces, as expected if the heavy neutrinos decouple if compared with the typical virtual momentum p 2 ≃ −(125 MeV) 2 , the well-known result characterizing the presence of only three light active neutrinos.

Additional remarks
Let us close this section summarizing further predictions and constraints on the inverse seesaw scenario. The aim of the following discussion is to provide additional motivations enforcing the phenomenological relevance of EW-scale right-handed neutrinos with O(1) Yukawa couplings.
Collider searches at the LHC. At the LHC right-handed neutrinos with a mass not far above or below the Higgs mass and with a sizable Yukawa coupling Y ν ≳ O(10 −2 ) affect the Higgs decay h → llνν (see [80] for a recent analysis). Present bounds hold in the range 60 ⩽ M R ⩽ 200 GeV with 10 −2 ⩽ Y ν ⩽ 2. 3 Furthermore, the CMS collaboration placed upper limits on the active-sterile neutrino mixings in the same mass range for M R considering direct production of heavy neutrinos [93,94]. 4 As stated in the introduction, these values of masses and couplings may have an impact on the stability of the EW vacuum, thus providing an additional motivation for the analysis that we shall perform in the next section.
Fit of LEP data via oblique parameters. The fit of LEP data still provides today an important constraint on beyond the SM physics. The presence of additional sterile neutrinos modifies the oblique radiative corrections [96][97][98][99][100][101][102][103]. In [79] it was shown that right-handed neutrino masses of the order M Ri ∼ O(10) TeV, together with violation of unitarity of the order αβ ∼ 10 −3 ÷ 10 −6 , can improve the fit of LEP data with respect to the SM.
Leptogenesis. In the inverse seesaw scenario the decay of (nearly degenerate) heavy Majorana neutrinos can realize the so-called resonant leptogenesis [104]. Remarkably, resonant leptogenesis can be realized with heavy Majorana neutrinos as light as 1 TeV [117] (in contrast with the usual thermal leptogenesis, realized in the type-I seesaw, in which M R ≳ 10 9 GeV [106]).
Naturalness. On a general ground, whenever a threshold with particles of mass M coupled to the Higgs with strength ξ is present, quantum corrections generate a contribution δm 2 H ≈ ξ 2 M 2 16π 2 to the renormalized Higgs boson mass. If δm 2 H ≫ v 2 , an unnatural cancellation between δm 2 H and the bare Higgs mass is required in order to reproduce the observed value of the Higgs boson mass. The condition δm 2 H ≲ v 2 can be used as a criterium to construct natural model of new physics [107,108]. In the context of inverse seesaw models, the scale M is the mass of right-handed heavy neutrinos M R , and the coupling ξ is the Yukawa coupling Y ν . TeV-scale values of M R and O(1) Yukawa couplings satisfy the naturalness condition.
To sum up, the presence of heavy neutrinos with a mass not far above the EW scale and a sizable Yukawa coupling with the Higgs boson has extremely rich phenomenological consequences. Motivated by the results achieved in this section, we are now in the position to tackle the second part of our analysis in which we aim to investigate the impact of heavy neutrinos on the stability of the EW vacuum.

Stability of the electroweak vacuum and low-scale seesaw
The key ingredient in the study of the stability of the EW potential is represented by the quantum effective potential V eff . In particular, we search for the instability scale Λ corresponding to the Higgs field value for which the potential becomes smaller than its value at the EW minimum. If such scale does not exist the vacuum is called absolute stable, otherwise the Higgs vacuum is not the global minimum, and a quantum tunnelling to the true vacuum may occur. The decay probability can be computed from the bounce solution of the euclidean equations of motion of the Higgs field [18,19]. If the corresponding lifetime is much bigger that the age of the Universe the EW vacuum is called metastable. This can be translated into a lower bound on the Higgs effective quartic coupling λ eff , driven to negative values, that at leading order reads as where τ is the age of the Universe τ = 4.35×10 17 sec [109] and µ is the renormalization scale of the RG running. The Higgs effective quartic coupling will be defined in the next section. The SM quantum corrections at zero temperature to the metastable condition have been computed in [110] but are not considered in this work. The measured values of the Higgs and the top masses place the SM in a metastable position in the phase space diagram [111]. The inclusion of extra fermionic degrees of freedom, as right-handed neutrinos in seesaw extensions, can clearly change this picture possibly driving the model to an instability region. Therefore, by requiring the metastability of the EW vacuum, together with the perturbativity of the gauge, scalar and Yukawa couplings up to the Planck mass, we can constrain the parameter space of low-scale seesaw models.
In the following sections we will describe the theoretical tools needed in our stability analysis: the Higgs effective quartic couplings which determines the instability scale, the RGEs in the MS describing the running of the couplings from the EW scale up to the Planck scale and, finally, the matching conditions which provide the initial values of the MS parameters from the EW physical observables.

The Higgs effective quartic coupling
From a preliminary analysis one can show that the vacuum instability appears at a scale much bigger than the EW minimum. This allows us to consider the effective potential at large field values, φ ≫ v, and to use the following approximation where the Higgs quadratic mass term has been neglected and t is the logarithm of the renormalization running scale µ. The effective quartic coupling λ eff is extracted from the RG-improved effective potential at two-loop order in the SM [112] and at oneloop for the right-handed neutrino corrections, computed in the MS renormalization scheme and in the Landau gauge. The complete two-loop expression is too lengthy to be given here. We show, instead, the effective quartic coupling at one-loop order in the SM [113] where the p-coefficients are summarized in table 1 and Γ(t) is defined as with γ the anomalous dimension of the Higgs field. The contribution of the heavy neutrinos to the RG-improved effective potential is [114] where θ νi (t) = θ(t − ln M Ri µ 0 ) and M νi (φ, t) are the three (double degenerate) nonzero eigenvalues of the 9 × 9 mass matrix withM Ri ≡ diag(M R1 , M R2 , M R3 ) and φ(t) = e Γ(t) φ. Compared with eq. (2.2), we are considering the µ S → 0 limit in which the three light neutrinos are massless. The factor of two in eq. (3.5) comes from the fact that each non-zero eigenvalue is double degenerate.
The RG running parameter t can be chosen in such a way that the convergence of perturbation theory -otherwise spoiled by the presence of large logs -is improved. We follow the prescription µ(t) = φ widely used in vacuum stability analyses.
The contribution of heavy neutrinos to the effective potential produces two distinctive effects on the RG-evolution of the effective quartic coupling in eq. (3.3).

○ Above each threshold M Ri and at high values of renormalization scale µ(t) ≫
M Ri it contributes explicitly to λ eff introducing the correction where κ νi (t) is the coefficient of the φ-dependent part of M νi (φ, t).
○ Below each threshold M Ri the corresponding heavy neutrinos are integrated out. The matching produces the threshold correction to the effective potential which translates into a threshold correction ∆λ νi th to the Higgs quartic coupling at the M Ri mass scale which can be extracted from the φ 4 term of eq. 3.8 and it is explicitly given by (3.9)

The matching conditions
The RG equations employed in this work are computed in the MS renormalization scheme and must be equipped with suitable initial conditions for the running parameters evaluated in the same scheme. These parameters are expressed in terms of physical observables, defined in the on-shell OS scheme, through appropriate relations called matching conditions. For a generic parameter α the matching at the scale µ is given by where α(µ) and α OS are respectively the MS and the OS parameters, while δα MS and δα OS are the corresponding counterterms in the two schemes. The difference between the two of them is an ultraviolet-finite correction. The SM parameters obtained from the matching procedure are (λ, y t , g 1 , g 2 ), while the QCD coupling constant g 3 is directly extracted from the α 3 (M Z ) which is already defined in the MS scheme. The physical observables used to compute the input values of the MS parameters are given in table 2. The details of the strategy can be found in [111] where the SM two-loop (NNLO) corrections to the matching conditions have been discussed. In particular, in [111] a complete two-loop analysis has been performed in the EW sector and the N3LO (three-loop) pure QCD effect has been included in the matching of the top Yukawa coupling and the strong coupling constant. The running of the latter from M Z to   [3], for all the other parameters we refer to [111].
M t , which is the starting scale of our stability analysis, has been performed including the QCD four-loop β function.
In low-scale seesaw extensions, the right-handed neutrino corrections to the matching conditions can be important and must be taken into account. In our analysis we have considered all the SM results given in [111], supplemented by new physics contributions computed at one-loop order from the Lagrangian in eq. (2.1). In particular, the additional neutrinos introduce corrections to the masses of the gauge bosons M Z , M W , the Higgs M h , the quark top M t and to the muon decay, from which the Fermi constant G µ is extracted. These corrections depend on the masses of the heavy neutrinos, their interactions with the SM fields mediated by the Yukawa couplings Y ν (obtained from the parameterization in eq. (2.5)), and the full 9 × 9 mixing matrix U . Due to the mass hierarchy M R ≫ m D , and to the smallness of the light neutrino masses -which can be safely neglected in the computation of the matching conditions at the EW scale -the mixing matrix U can be expanded in the ratio m D M R , with µ S → 0, as Notice that the PMNS block has been set to the unit matrix, consistently with the approximation m νi ≃ 0. Moreover, we have verified the Appelquist-Carazzone theorem in the new physics sector. In particular, we have checked the decoupling of right-handed neutrino contributions from the matching conditions in the limit of large Majorana masses. This has to be expected since these masses are not generated by spontaneous symmetry breaking of the Higgs field.
Concerning the matching conditions of the Yukawa matrix Y ν , computing perturbative corrections to the matching conditions at the EW scale will not considerably improve the precision on the determination of Y ν in the MS scheme (since in any case Y ν turns out to be related to unknown parameters by means of eq. (2.5)). Therefore, for the sake of simplicity, we decided to match the MS Yukawa matrix Y ν to its OS version at the tree-level, namely Y ν (M t ) ≃ Y ν .

The renormalization group equations
All the dimensionless couplings (λ, y t , g i , Y ν ) are evolved from the top-mass scale, M t , up to the Planck scale using the three-loop (NNLO) RG equations for the SM parameters and the two-loop (NLO) β functions for the Yukawa matrix Y ν . Here g i stands for the three gauge coupling constants and we have retained only the topquark contribution in the SM Yukawa sector. The system of coupled RG equations is then given by where κ = 1 (16π 2 ) and the abelian gauge coupling is given in the GUT normalization g 1 = 5 3g Y . Notice that the RG equations given above are defined for a renormalization scale µ bigger than any particle mass of the model. For lower scales, heavy degrees of freedom are integrated out and the corresponding coupling must be removed by hand from the β functions. Indeed, in the MS renormalization scheme the decoupling of heavy degrees of freedom is not automatic and has to be explicitly implemented at the different particle thresholds. As far as the vacuum stability is concerned, the right-handed neutrinos behave like the top-quark and drive λ to negative values faster than the SM case. This is clear from eq. (3.16). On the other hand, their impact on the top Yukawa couplingsee eq. (3.17) -is to increase y t all along the RG evolution. Actually, this is another source of vacuum destabilization with respect to the SM picture, because a bigger value of y t has a bigger decreasing effect on the Higgs quartic coupling. Nevertheless, the overall behavior of the top-quark Yukawa coupling is dictated by the large and negative QCD corrections which lead to a decreasing y t . These features can be easily deduced from fig. 4, left panel, where the running of the SM couplings is depicted. Here solid lines represent the evolution in the seesaw extended scenario, while dashed curves correspond to the pure SM. Contrary to the top Yukawa case, the QCD contributions are obviously absentat least in the leading one-loop approximation -from the evolution of Y ν , and the Tr(Y † ν Y ν ) term, which affects β λ e β yt , is always increasing. This feature, shown in the right panel of fig. 4, has a negative impact both on the vacuum stability and on the perturbativity of Y ν which can be violated, if Y ν ij > √ 4π, during the RG evolution for sufficiently big values of the Yukawa coupling at the EW scale. In the same figure the decoupling of the heavy right-handed neutrinos below their mass thresholds is also manifest. Indeed, for µ ≪ M Ri , the N i R neutrino is integrated out and does not contribute to the RG running: the corresponding row in the Yukawa matrix Y ν is frozen, and enters in the β functions only above the threshold scale M Ri as shown in eq. (3.20) where we mark with a generic × non-zero entries for the Yukawa matrix Y ν .
Finally, we show in fig. 5 the evolution of the effective quartic coupling λ eff in two different scenarios. In the left panel the inverse seesaw is realized with Tr(Y † ν Y ν ) ≃ 0.36, while in the right panel the Yukawa matrix is such that Tr(Y † ν Y ν ) ≃ 0.6. In the latter case the effects of Y ν , which affects the RG running above the threshold scales, are quite large and λ eff is driven outside the metastability region below the Planck scale.

Results and discussion
In this section we discuss the results of our analysis. First, let us briefly remind our strategy. The Yukawa matrices generated following the prescription outlined in section 2.3 enter as initial conditions, together with all the other SM external parameters, in the solution of the RGEs in eqs. (3.12-3.15). From the effective potential, improved by the running couplings previously computed, we extract the Higgs effective quartic coupling λ eff in eq. (3.3). Finally, we use eq. (3.1) to check whether the analyzed points violate the metastability bound on the EW vacuum.
We show our results in fig. 6 for the normal ordering and in fig. 7 for the inverted ordering. In order to make contact with phenomenology, in both cases we present the impact of the metastability bound with respect to the observables targeted in section 2.3.1, namely the branching ratio Br(µ → eγ) (left panel) and the effective neutrino mass (right panel). The red points are excluded by the metastability bound: the EW vacuum would decay too fast in the true vacuum of the EW potential.
Inverse seesaw: Normal Ordering Our numerical analysis clearly indicates that points with Yukawa couplings such that Tr(Y † ν Y ν ) ≳ 0.4 are excluded. This bound does not depend on the assumed hierarchy, since light neutrino masses are irrelevant. Most importantly, the excluded points lie in a region of the parameter space that is close to the present bounds and future experimental sensitivities for both the analyzed observables. Moreover, as clear from the left panel of fig. 3, the metastability bound applies in the whole range of analyzed masses for the right-handed neutrinos. This is an interesting piece of information since, for instance, searches for heavy neutrinos with m ν4 ∼ O(100) GeV and sizable Yukawa couplings are currently ongoing at the LHC (see section 2.3.3).

EXO 200
Figure 7: Same as in fig. 6, but considering the inverted ordering.
Before moving to conclusions, we refer the reader to the appendix A for a generalization of these results to the linear and double seesaw models. In fact, one might be concerned about the model-dependence of our bound, since we worked in a well-defined model. By extending our analysis to the two aforementioned models, we intend to prove that the significance of the metastability bound is a general aspect of low-scale seesaw models.

Conclusions and outlook
In this paper we analyzed the stability of the EW vacuum in the presence of low-scale seesaw models. For definiteness, in the main part of this manuscript we focused on the ISS model, but in appendix A we generalized our results to other popular lowscale scenarios. In the ISS model lepton number is broken at the scale µ S , and the scale of light neutrino masses is given by m ν ≈ µ S (vY ν M R ) 2 where M R is the mass of the heavy neutrinos. Our analysis can be divided in two parts.
First, in section 2, we carefully chose the parameter space of our numerical scan in order to comply with low-energy neutrino data. We performed this selection by means of the Casas-Ibarra parametrization, and we considered the intervals 10 2 GeV ⩽ M R ⩽ 10 2 TeV and 10 −1 keV ⩽ µ S ⩽ 10 2 keV in order to generate order one Yukawa couplings. As expected, the output of our numerical scan shows that the ISS model presents an extremely rich low-energy phenomenology, and we decided to focus our attention on the lepton flavor radiative decay µ → eγ and on the 0ν2β.
Second, in section 3, we computed the contribution of the heavy neutrinos to the RG evolution of the dimensionless SM couplings and to the RG-improved effective potential. On the technical level, we worked including the most updated results for the computation of β functions and matching conditions in the SM, and we added the contributions of the right-handed neutrinos at two-loops in the β functions and one-loop in the matching.
With such heavy artillery, we analyzed the Yukawa couplings generated in section 2 by imposing the metastability condition on the lifetime of the EW vacuum. Our conclusions are summarized in fig. 6 and fig. 7: Yukawa couplings such that Tr(Y † ν Y ν ) ≳ 0.4 are excluded by the metastability condition. The bound applies in the whole range of analyzed right neutrino masses, and affects a region of the parameter space close to present or future sensitivities for both µ → eγ rate and 0ν2β.
To sum up, we argue that the metastability bound represents an important consistency condition that should be included in all the phenomenological analysis of low-scale seesaw models featuring EW-scale right-handed neutrinos with O(1) Yukawa couplings.

A Linear and double seesaw
In full generality, we introduce a Majorana mass term for the right-handed neutrinos and a Yukawa coupling for the fermionic singlets In the left-handed basis N L ≡ (ν L , N C R , S) T we have, after EW symmetry breaking, the following mass matrix ○ The linear seesaw [115] corresponds to µ S = 0, M N = 0 In this case the only source of lepton number violation comes from m S . Following the standard seesaw approximation we find Generalizing the Casas-Ibarra parametrization [56], we obtain [53] where A is a 3 × 3 matrix satisfying the equation A + A T = 1. Consequently, A ii = 1 2 and A ij = −A ji .
○ The double seesaw corresponds to m S = 0 The Majorana mass term M N , in addition to µ S , violates lepton number for two units. We assume the hierarchy M N ≫ M R ≫ m D ≫ µ S . Integrating out the heavy fields N R , we find the effective Lagrangian In the left-handed basis N L = (ν L , S) T eq. (A.7) corresponds, after EW symmetry breaking, to the mass matrix After block-diagonalization we find From eq. (A.10) it follows that m ν → 0 if µ S → 0 since the type-I contribution m T D M −1 N m D cancels out between the two remaining terms. Considering a perturbative expansion in µ S , we find As a result, the Casas-Ibarra parametrization is analogous to eq. (2.5). The unitary mixing matrix V in eq. (A.9) has the general structure and, at the lowest order, we find B * = m T D (M T R ) −1 . For the sake of simplicity, we focus on the case with n R = n S = 3, and we assume M R = diag(M Ri ), and M N = diag(M N i ) with i = 1, 2, 3. In this case the heavy block in eq. (A.9) simplifies to In eq. (A.7) the interactions between the lepton doublets L i and the three singlet fermions S i are mediated by an effective

A.1 Linear seesaw: results
For definiteness, we focus on the case with n R = n S = 3. Considering the low-energy neutrino data, we follow the same strategy outlined in section 2 (see in particular section 2.3). We make use of the Casas-Ibarra parametrization in eq. (A.5) to numerically reconstruct the Yukawa matrix Y ν , and we randomly scan over the intervals 10 2 GeV ⩽ M Ri ⩽ 10 2 TeV, 10 −2 keV ⩽ (m S ) ij ⩽ 10 2 keV, and 10 −1 ⩽ A ij ⩽ 10 2 . As done for the inverse seesaw, we discard points unable to comply with the bounds discussed in section 2.2. As far as the stability of the EW vacuum is concerned, in the m S → 0 limit the mass matrix M LSS in eq. (A.3) reduces to the same structure already studied in the inverse seesaw case (see eq. (2.2)). Consequently, in the definition of the Higgs effective quartic coupling we employ the same RG equations, matching conditions and effective potential used in the inverse seesaw analysis. We show the final result of our analysis in fig. 8 considering a normal ordering of light neutrino masses. As expected, we find the same quantitative conclusion if compared with the inverse seesaw case. Yukawa couplings such that Tr(Y † ν Y ν ) ≳ 0.4 are excluded by the metastability bound. Remarkably, considering both the lepton flavor violating process µ → eγ and the 0ν2β, these points lie in a region of the parameter space close to present or future sensitivities. Therefore, we conclude that the metastability bound represents an important consistency condition that should be included in all the phenomenological analysis of the linear seesaw model featuring EW-scale righthanded neutrinos with O(1) Yukawa couplings.

A.2 Double seesaw: results
We follow the same approach already exploited for the inverse and linear seesaw models. However, in the double seesaw case there are few remarkable differences. We use the Casas-Ibarra parametrization in eq. (2.5) to sample the Yukawa matrix Y ν , and we randomly scan over the intervals 10 TeV ⩽ M Ri ⩽ 10 3 TeV, 10 9 GeV ⩽ M N i ⩽ 10 11 GeV, and 10 keV ⩽ (µ S ) ij ⩽ 10 4 keV. Few comments are in order. First, notice that this choice of parameters -optimized in order to obtain O(1) Yukawa couplings -respects the hierarchy M N ≫ M R ≫ m D ≫ µ S assumed above (see discussion below eq. (A.6)). Second, we expect the following order of magnitude estimates: for the mixing parameter in eq. (A.13), B ∼ O(10 −3 ); for the mass of the heavy neutrinos in eq. (A.11), M heavy ∼ O(1) GeV; for the effective Yukawa coupling in eq. (A.7),Ỹ ν ∼ 10 −5 Y ν . Armed with these numbers, we can outline as follows. At large renormalization scale values, µ ≫ M N i , the model is described by the full Lagrangian in eq. (A.1) (with Y S = 0). In terms of the RG running, the only relevant parameter is the Yukawa matrix Y ν describing the interactions between the Higgs doublet, the lepton doublets and the three right-handed heavy neutrinos. At this stage, the situation is formally equivalent to the familiar type-I seesaw. 5 There is, however, one remarkable difference. In the type-I seesaw righthanded neutrino masses M N i ∼ O(10 9 ) GeV require, in order to reproduce lowenergy neutrino phenomenology, small Yukawa couplings (typically Y ν ∼ 10 −5 ). As a consequence, the impact of the interactions N R Y νH † L + h.c. on the running of the Higgs quartic couplings is negligible. In the double seesaw case, on the contrary, we are allowed to consider O(1) Yukawa couplings since the mass of light neutrinos is set by µ S (see eq. (A.12)). Below the thresholds M N i the heavy right-handed neutrinos are integrated out, and eventually the model is described by the effective Lagrangian in eq. (A.7). Given our choice of parameters, in this region the running is approximately equivalent to the pure SM since for the effective Yukawa interactions L CH * Ỹ ν S we expectỸ ν ∼ 10 −5 Y ν .
We summarize our results in fig. 9. In the left panel, we show the running of the effective Higgs quartic coupling for a specific realization of the double seesaw model with Tr(Y † ν Y ν ) ≃ 0.62. Above the thresholds M N i the Yukawa couplings Y ν sizably affect the running of λ eff eventually violating the metastability bound before the Planck scale. In the right panel we show the result of our numerical scan focusing on the effective neutrino mass relevant for the 0ν2β. The most striking difference with respect to the inverse and linear seesaw models (see, respectively, figs. 6,8) is that the presence of additional neutrinos with mass M heavy ∼ O(1) GeV gives a sizable contribution to m νe eff . As a result, numerous points in our numerical analysis are close to (or even exceed) the present experimental bound. We find that Yukawa couplings such that Tr(Y † ν Y ν ) ≳ 0.6 are excluded by the metastability bound. Before concluding, let us stress that the aim of this discussion is not to be exhaustive since a careful phenomenological analysis is clearly far beyond the scope of this appendix. The most important message of our discussion is that even in the presence of very heavy right-handed neutrinos (with mass of the same size as the one expected in standard type-I scenarii) large modifications of the effective Higgs quartic couplings are a concrete and realistic possibility. As a final remark, notice that double seesaw models represent an appealing setup to realize leptogenesis [116,117]. With this respect, a more detailed study of the parameter space will be addressed in a forthcoming publication.