NLO QCD bottom corrections to Higgs boson production in the MSSM

We present a calculation of the two-loop bottom-sbottom-gluino contributions to Higgs boson production via gluon fusion in the MSSM. The calculation is based on an asymptotic expansion in the masses of the supersymmetric particles, which are assumed to be much heavier than the bottom quark and the Higgs bosons. We obtain explicit analytic results that allow for a straightforward identification of the dominant contributions in the NLO bottom corrections. We emphasize the interplay between the calculations of the masses and the production cross sections of the Higgs bosons, discussing sensible choices of renormalization scheme for the parameters in the bottom/sbottom sector.


Introduction
With the coming into operation of the Large Hadron Collider (LHC), a new era has begun in the search for the Higgs boson(s). At the LHC the main production mechanism for the Standard Model (SM) Higgs boson, H SM , is the loop-induced gluon fusion mechanism [1], gg → H SM , where the coupling of the gluons to the Higgs is mediated by loops of colored fermions, primarily the top quark. The knowledge of this process in the SM includes the full next-to-leading order (NLO) QCD corrections [2,3], the next-to-next-to-leading order (NNLO) QCD corrections [4] including finite top mass effects [5], softgluon resummation effects [6], an estimate of the next-to-next-to-next-to-leading order (NNNLO) QCD effects [7] and also the first-order electroweak corrections [8,9,10].
The Minimal Supersymmetric extension of the Standard Model, or MSSM, features a richer Higgs spectrum which consists of two neutral CP-even bosons h, H, one neutral CP-odd boson A and two charged scalars H ± . The gluon-fusion process is one of the most important production mechanisms for the neutral Higgs bosons, whose couplings to the gluons are mediated by colored fermions and their supersymmetric partners. The gluon-fusion cross section in the MSSM is known at the NLO. The contributions arising from diagrams with squarks and gluons were first computed under the assumption of vanishing Higgs mass in ref. [11]. The complete top/stop contributions, including stop mixing and gluino effects, were computed under the same assumption in ref. [12], and the result was cast in a compact analytic form in ref. [13]. Later, more refined calculations aimed at the inclusion of the full Higgs-mass dependence. In particular, the full squark-gluon contribution is known in a closed analytic form [14,15,16], while the full quark-squark-gluino contribution has been computed in ref. [17] via a combination of analytic and numerical methods.
It should be stressed that, at least for the case of the light Higgs, the exact two-loop QCD Higgsgluon-gluon amplitude is in general well approximated by the amplitude evaluated in the limit of neglecting the Higgs mass. The latter is much easier to compute and the corresponding result can be straightforwardly implemented in computer codes that aim to evaluate the Higgs boson production cross section in a fast and efficient way. Indeed, it was noticed several years ago for the SM case [18] that the exact K factor, defined as the ratio between the NLO and leading-order (LO) cross sections, is well approximated by the so-called effective K factor that can be obtained via an improved effectivetheory calculation. By the latter we mean a result in which the effective NLO cross section is obtained by multiplying the exact LO partonic cross section by the O(α s ) corrections evaluated in the limit of vanishing Higgs mass. For the SM case this approximation works at the level of few per cent for Higgs mass values below the 2 m t threshold, and up to 10% for any Higgs mass value. The same level of accuracy is reached when the Higgs couples to a generic scalar particle with mass m S , with the exception of a narrow region close to the m H ≃ 2 m S threshold [19].
There is only one case in which the effective approximation does not work sufficiently well, namely when the bottom contribution becomes very relevant. This can happen in the MSSM when tan β, i.e. the ratio of the vacuum expectation values (vev) of the neutral components of the two Higgs doublets, becomes large. In such a situation, in principle, the exact computation of the NLO bottom contribution to the cross section should be employed. However, explicit analytic results for the part of the NLO bottom contribution that is induced by two-loop bottom-sbottom-gluino diagrams have not been made available so far.
In this paper we present an approximate evaluation of the bottom-sbottom-gluino diagrams, based on an asymptotic expansion in the large supersymmetric masses that is valid up to and including terms of where m φ denotes a Higgs-boson mass and M denotes a generic superparticle mass (M = mg, mb 1 , mb 2 ). Together with the known exact results for the (s)bottom-gluon diagrams [3,14,15,16,19], our new result allows us to obtain effective K factors that can be easily implemented in computer codes to provide an accurate and efficient evaluation of the cross section for Higgs boson production in the MSSM. Such K factors are expected to be at the same level of accuracy as in the SM, i.e. within a few per cent of the exact value for the light Higgs and up to ten per cent for the heavy Higgs, with however a much better accuracy expected if the heavy-Higgs mass is below all the heavy-particle thresholds.
The paper is organized as follows: in section 2 we summarize general results on the cross section for Higgs boson production via gluon fusion. In section 3 we outline the way we perform the asymptotic expansion of the bottom-sbottom-gluino diagrams. Section 4 contains the explicit results for the NLO bottom contribution both in DR and in the on-shell (OS) scheme. In section 5 we discuss the numerical relevance of the NLO bottom corrections devoting particular attention to the interplay between the computation of the Higgs mass and that of its production cross section. In the last section we discuss an approximate way to take the NLO bottom contribution into account, based on an improved LO term. Finally, in appendix A we specialize to the MSSM case the general exact results for the real radiation derived in ref. [19], while appendix B contains the explicit expressions for the shifts from the DR to the OS parameters in the sbottom sector.

Higgs boson production via gluon fusion at NLO in the MSSM
In this section we recall for completeness some general results on Higgs boson production via gluon fusion. The hadronic cross section for Higgs boson production at center-of-mass energy √ s can be written as the parton density of the colliding hadron h i for the parton of type a (for a = g, q,q), andσ ab the cross section for the partonic subprocess ab → φ + X at the center-of-mass energyŝ = x 1 x 2 s = m 2 φ /z. The latter can be written in terms of the LO contribution σ (0) and a coefficient function G ab (z) aŝ We consider now the production of the lightest CP-even Higgs boson, h, through gluon fusion 1 . The LO term can be written as where τ k ≡ 4 m 2 k /m 2 h , and the functions G 1ℓ 0 and G 1ℓ 1/2 read The analytic continuations are obtained with the replacement m 2 h → m 2 h + iǫ . For later convenience, we recall the behavior of G 1ℓ 0 and G 1ℓ 1/2 in the limit in which the Higgs boson mass is much smaller or much larger than the mass of the particle running in the loop. In the first case, i.e. τ ≫ 1, which applies to the top and squark contributions for the light-Higgs case, while in the opposite case, i.e. τ ≪ 1, which is relevant for the bottom quark, The coefficient function G ab (z) in eq. (2) can be decomposed, up to NLO terms, as with the LO contribution given only by the gluon-fusion channel: The NLO terms include, besides the gg channel, also the one-loop induced processes gq → qh and qq → gh: where the LO Altarelli-Parisi splitting functions are In the equations above, C A = N c and C F = (N 2 c − 1)/(2 N c ) (N c being the number of colors), β 0 = (11 C A − 2 N f )/6 (N f being the number of active flavors) is the one-loop β-function of the strong coupling in the SM, and The gg-channel contribution, eq. (17), involves two-loop virtual corrections to gg → h and oneloop real corrections from gg → hg. The former, regularized by the infrared-singular part of the real emission cross section, are displayed in the first line of eq. (17). The second line contains the nonsingular contribution from the real gluon emission in the gluon fusion process. The latter contribution as well as the ones due to the qq → hg annihilation channel and the quark-gluon scattering channel, eq. (18), are obtained from one-loop diagrams where only quarks or squarks circulate in the loop. General expressions for the functions R gg , R qq , R qg can be found in ref. [19] (see also refs. [20,21]). In appendix A we provide expressions in which the contribution of the bottom quark is kept exact while those of the top quark and of the squarks are evaluated in the limit of vanishing Higgs mass.
The two-loop top/stop contributions to the form factors H 2ℓ 1,2 entering eq. (17) are fully under control in the light-Higgs case. Typically, the mass ratios between the Higgs and the particles running in the loops allow for the evaluation of the relevant diagrams via a Taylor expansion in the Higgs mass, with the zero-order term in the series already a very good approximation of the full result. The case of the two-loop bottom/sbottom contributions is obviously different. In general, Taylor-expanded evaluations of the relevant diagrams are no longer viable, due to the presence of the light bottom quark in the loops. Thus, the diagrams must be evaluated either exactly or via an asymptotic expansion in a large mass or momentum.
In the following sections we present the result for the NLO bottom contribution, combining earlier results in the literature with our new calculation of the bottom-sbottom-gluino contribution. The latter has been obtained via an asymptotic expansion, retaining terms of

Outline of the calculation
An exact analytic evaluation of the bottom-sbottom-gluino contribution to H 2ℓ 1,2 is, at the moment, beyond our computational ability. However, it is reasonable to assume that all of the supersymmetric particles are much heavier than the lightest Higgs boson and the bottom quark, and look for an approximate evaluation of the diagrams in terms of a small-momentum (large-mass) expansion. We follow this path by performing a large-mass expansion, assuming all the supersymmetric particles to be heavy but without requiring any specific hierarchy among them.
After generating the two-loop diagrams involving bottom, sbottom and gluino that contribute to the process g(q 1 ) + g(q 2 ) → h(q) with the help of FeynArts [22], we separate them in two classes: i) those that can be evaluated via an ordinary Taylor expansion in powers of q 2 /M 2 , of which we keep only the term of order zero; ii) the diagrams that require an asymptotic expansion. We recall that a Taylor expansion of a two-loop diagram in the external momentum q 2 is viable for values of q 2 up to the first physical threshold. In our case, diagrams with a physical threshold at q 2 = 4 m 2 b , when Taylor-expanded in q 2 , exhibit an infrared (IR) divergent behavior as the bottom mass is sent to zero. Thus, these diagrams belong to the class ii.
Class-i diagrams are expressed in terms of two-loop vacuum integrals that can be evaluated using the results of ref. [23]. Concerning the diagrams belonging to class ii, reviews of the method of asymptotic expansions of Feynman diagrams with respect to masses or momenta can be found in ref. [24]. In practice, we generate the expansion of a diagram by adding and subtracting to it the part of the diagram itself that becomes IR-divergent when m b and q 2 are sent to zero. Formally we are adding nothing to the original diagram but, as graphically 2 exemplified in fig. 1, this construction allows us to separate the diagram in two parts: part A in fig. 1 which, being by construction IRsafe, can be evaluated via a Taylor expansion in the same way as class-i diagrams; part B in fig. 1, containing the IR-divergent contribution, which should be evaluated exactly.
The IR-divergent part of a diagram is constructed in the following way. We first note that in all the diagrams entering our calculation one can choose a routing of momenta such that the connecting propagators, i.e. the propagators that contain both integration momenta k 1 and k 2 , are always accompanied by a heavy mass M . Furthermore, only one subintegration, let us assume the one on k 2 , is IR divergent. Then, one can rewrite the connecting propagators using the identity The first term on the r.h.s. of eq. (21) leads to a disconnected integral (product of two one-loop integrals) that contains the IR-divergent contributions present in the original diagram. This term can be evaluated exactly, i.e. for arbitrary q 2 , giving rise to the ln(q 2 /m 2 b ) terms that describe the physical threshold. The second term, instead, leads to a two-loop integral with improved infrared convergence in the k 2 integration and improved ultraviolet convergence in the k 1 integration. Therefore, if, for example, the original integral is logarithmically IR divergent in the k 2 integration when q → 0 and m b → 0, the corresponding two-loop integral associated with the second term in eq. (21) evaluated at q 2 = m b = 0 is no longer IR divergent, but it actually gives a finite result that differs from the result valid for q 2 = 0 and m b = 0 by terms of O(m 2 /M 2 ln(m 2 /M 2 )), where m 2 denotes either q 2 or m 2 b . In general, a repeated application of eq. (21), controlled by the power counting in the IR-divergent terms, allows us to construct the IR-divergent part of any diagram in terms of products of one-loop integrals with numerators that contain terms of the form (k i · q j ) m , (k i · k j ) n (i, j = 1, 2) where m, n are generic powers. The Passarino-Veltman reduction method is then applied to eliminate the numerators and express the result in terms of the known one-loop scalar integrals [26]. A check of the validity of our construction of the IR-divergent part of a diagram is given by the evaluation of its part A. Indeed, one verifies explicitly that the IR-divergent contributions of the original diagram are canceled by the terms constructed via eq. (21), so that the final result for part A is free of any ln(q term.

Two-loop bottom/sbottom contributions
In this section we present the result for the two-loop bottom/sbottom contribution to the form factor for Higgs boson production via gluon fusion. We stress that, in the MSSM, the result for the production cross section of a Higgs boson is strictly linked to the computation of its mass, i.e., both observables should be computed in terms of the same set of SUSY parameters, defined in the same way beyond tree level. There is however an important difference between the two calculations. In the computation of the one-loop corrections to the Higgs masses, the diagrams involving the bottom quark are suppressed by the bottom mass and can be safely neglected, resulting in a one-loop contribution that is actually due to the sbottom diagrams only. This implies that in the two-loop calculation of the Higgs masses the only couplings that require a one-loop renormalization are the trilinear sbottom-Higgs couplings, while the definition of the bottom-Higgs Yukawa coupling beyond tree level is irrelevant. On the other hand, in the one-loop calculation of the amplitude for Higgs production both the bottom-Higgs Yukawa coupling and the trilinear sbottom-Higgs couplings play a role, thus they both require a one-loop definition when the two-loop contributions are computed.   (5) and (6). For both kinds of contributions, we first report the results obtained in the DR renormalization scheme, which is the scheme employed by several public computer codes [27] that return the MSSM mass spectrum starting from a set of high-energy boundary conditions for the SUSY-breaking parameters. We then discuss how to move from the DR scheme to a different renormalization prescription, which generalizes the one introduced for the calculation of the Higgs masses in ref. [28] and implemented in the code FeynHiggs [29]. With a slight abuse of language, since some of the amplitudes involved in the renormalization of the SUSY parameters are in fact evaluated off mass shell, we refer to this prescription as to the "on-shell" scheme.

Contributions controlled by the Higgs-bottom coupling
We start by considering the contributions of the two-loop diagrams involving the Higgs-bottom coupling. In our approximation those are the only diagrams that contribute to the function G b in eq. (5). The two-loop part of the function can be decomposed as Assuming that the one-loop form factor H 1ℓ 1 is expressed in terms of DR-renormalized parameters evaluated at the scale Q 2 , the contribution of the two-loop diagrams with bottom quarks and gluons Exact expressions for the functions F (2ℓ,a) are given in eqs. (2.12), (2.13) and (3.8) of ref. [15], respectively. In the limit τ ≪ 1 they reduce to where 1 is expressed in terms of DR-renormalized parameters evaluated at the scale Q 2 , they read: where x i = m 2 b i /m 2 g , and (δm b ) SU SY denotes the SUSY contribution to the bottom self-energy, in units of C F α s /π and in the limit of vanishing m b : where If the bottom-quark contribution to H 1ℓ 1 is expressed in terms of the pole bottom mass M b , the functions G are shifted with respect to their expressions in eqs. (23) and (27). In particular, the former becomes and the term proportional to (δm b ) SU SY in the first line of eq. (27) is canceled out. Eqs. (27) and (29) This manipulation amounts to differentiating, in the one-loop contribution, between the parameter that describes the mass of the bottom quark running in the loop -which is identified with M b -and the parameter that describes the Yukawa coupling of the bottom quark to the Higgs boson -which is identified with m b . We recall that, in the MSSM, the running bottom mass m b can be related to the corresponding SM parameter m b as [30] where δ b denotes terms that are not enhanced by tan β and, to O(α s ), Since ) , it is easy to see that the terms enhanced by mg/m b in do indeed contain ǫ b tan β. In the effective-theory language of ref. [31] we can argue that, by expressing the bottom Yukawa coupling entering H 1ℓ 1 in terms of m b as defined in eq. (33), we "resum" in the one-loop part of the form factor the tan β-enhanced threshold corrections to the relation between the mass and the Yukawa coupling of the bottom quark. As a result of this special choice of parameters, all terms of the form ǫ b tan β drop out of the two-loop part of the form factor.

Contributions controlled by the Higgs-sbottom coupling
We now turn our attention to the diagrams that involve the Higgs-sbottom coupling. In our approximation those diagrams contribute only to the functions F b , F b and G b in eqs. (5) and (6). In analogy with ref. [13], the two-loop parts of the functions can be written as Furthermore, the various terms in eqs. (35)- (37) can be split in the contributions coming from diagrams with gluons (g, fig. 3a), with strong, D-term-induced quartic bottom couplings (4b, fig. 3b), and with gluinos (g, fig. 3c), The first two terms in eq. (38) can be obtained from the analytic expressions presented for the stop contributions in eqs. (27)-(30) of ref. [13], identifying ∂Z a ∂x with Y a x after making the trivial replacement t →b. The gluino contributions, on the other hand, require a dedicated calculation. Writing and assuming that the parameters in H 1ℓ 1 and H 1ℓ 2 are expressed in the DR scheme at the renormalization scale Q 2 , we find for the functions Yg x : where we retained only terms that induce The expression for Yb 2 can be obtained from the expression for Yb 1 through the replacements x 1 → x 2 and s 2θ b → −s 2θ b . Comparing eqs. (40) and (41) with eq. (42) of ref. [13] we notice that, contrary to what we stated in section 3.3 of that paper, even for θ b = 0 the two-loop bottom-sbottomgluino contribution to F b and G b cannot be obtained by taking the limit m t → 0 in the corresponding top contribution.

On-shell renormalization scheme for the sbottom parameters
We now discuss a suitable OS renormalization scheme for the parameters that determine the sbottom contribution to H 1ℓ 1,2 . We recall that, at the one-loop level, the vev v 1 and v 2 , the Z boson mass, the Weinberg angle and the parameter µ are not renormalized by the strong interactions. Therefore, the only parameters that require a one-loop definition are (h b , A b , s 2θ b , mb 1 , mb 2 ), where by h b we denote the coupling constant entering the cubic and quartic sbottom-Higgs interactions, which at tree level is related to the bottom mass by m b = h b v 1 / √ 2. Indeed, the factor m b that multiplies the function F b in eqs. (5) and (6) has to be interpreted as a bookmark for h b . In fact, only four of those parameters are independent, because of the relation In the analysis of the sbottom corrections to the neutral MSSM Higgs boson masses presented in ref. [28] it was pointed out that, while the sbottom masses can be naturally identified with the pole masses, an OS definition of (h b , A b , s 2θ b ) is less easily singled out. Proceeding in analogy with the OS renormalization of the stop sector (see, e.g., ref. [13]), we might choose as independent parameters a conveniently defined bottom mixing angle, s 2θ b , and the bottom Yukawa coupling h pole b , as defined by the pole bottom mass Then, eq. (44) might be used to establish the one-loop definition of A b in terms of the pole bottom and sbottom masses and the sbottom mixing angle. However, for large values of tan β such definition would produce very large shifts in A b with respect to its DR value 3 , δA b = O(α s µ 2 tan 2 β/mg) [32]. This is related to the fact that, in the large-tan β limit (i.e., v 1 → 0), s 2θ b becomes independent of A b , as can be easily seen from eq. (44). To cure the problem, it was suggested in ref. [28] (see also ref. [33]) to take s 2θ b and A b as independent parameters, while considering h b as a derived quantity via eq. (44). Suitable renormalization conditions were then proposed for s 2θ b and A b .
In the OS analysis of the cross section for Higgs boson production we want to retain the convenient features of the renormalization prescription employed in ref. [28]. However, that prescription needs to be expanded: first of all, the renormalization conditions in ref. [28] were defined in the limit tan β → ∞ (i.e., v 1 = 0), while in the case at hand we do not impose constraints on tan β. Moreover, while in the calculation of the one-loop corrections to the Higgs masses the contributions controlled by the bottom-Higgs Yukawa coupling (which we denote as h Y b to distinguish it from h b ) are suppressed, the one-loop diagram controlled by h Y b gives an important contribution to the production cross section. Therefore, a one-loop definition of h Y b is required. Since it does not seem appropriate to define the bottom-Higgs coupling h Y b in terms of quantities of the sbottom sector, as would happen if we imposed on it the same renormalization condition used for h b , we need to impose different renormalization conditions on h Y b and h b , or, equivalently, on the bottom mass that enters the one-loop bottom contribution and the one that enters the one-loop sbottom contribution. In particular, we identify the former with the pole mass M b (the resulting shifts in the function G 2ℓ b are discussed in section 4.1). Concerning the bottom mass in the one-loop sbottom contribution, we follow ref. [28], extending the prescription presented in that paper to the case of finite tan β.
To obtain definitions for δh b and δA b , we consider two quantities that allow for a natural interpretation: X b , at the classical level, is the off-diagonal term in the sbottom mass matrix, related to the mixing angle s 2θ b via eq. (44); Y b is proportional to the coefficient of the trilinear interaction ( b 1 b * 2 A). A definition of the mixing angle θ b like the one proposed in ref. [34], together with the identification of the sbottom masses as pole masses, can be immediately translated, using eq. (44), into a prescription for X b : In eqs. (46) and (47), Π ij (q 2 )(i, j = 1, 2) denotes the finite part of the (i, j) self-energy of the sbottoms.
Recalling that in the tan the extension to the case of finite tan β of the prescription for A b introduced in eq. (15) of ref. [28] reads: where iΛ 12A (p 2 1 , p 2 2 , p 2 A ) denotes the proper vertexb 1b * 2 A. Finally, the shifts of the parameters h b and A b are related to those of X b and Y b by where v = v 2 1 + v 2 2 . Explicit expressions for δh b and δA b , as well as for δs 2θ b and δm 2 b i , can be found in appendix B.
If the one-loop sbottom contribution to H 1ℓ 1 and H 1ℓ 2 is evaluated in terms of OS quantities, the two-loop functions in eqs. (35)-(37) must be replaced by In addition, the two-loop form factor H 2ℓ 1 receives a contribution originating from the shift in A b :

A numerical example
We will now illustrate the effect of the two-loop bottom/sbottom contributions to the form factors for the production of a Higgs boson in a representative region of the MSSM parameter space. The SM parameters entering our calculation are the Z boson mass m Z = 91.1876 GeV, the Fermi parameter G F = 1.16637 × 10 −5 GeV −2 , the sine of the Weinberg angle s 2 θ W = 0.223 and the strong coupling constant α s (m Z ) = 0.118 [35]. For the pole masses of the top and bottom quarks we take M t = 173.1 GeV [36] and M b = 4.49 GeV, the latter corresponding to the SM running mass (in the MS scheme) m b (m b ) = 4.16 GeV [37]. The tree-level mass matrix for the CP-even Higgs bosons can be expressed in terms of the physical pseudoscalar mass m A and the DR-renormalized parameter tan β, in addition to m Z . In the calculation of the physical Higgs boson masses and of the mixing angle α we include the one-loop O(α t + α b ) and two-loop O(α t α s + α b α s ) corrections as in refs. [28,38].
When computing the two-loop corrections to both mass matrix and production form factors for the CP-even Higgs bosons, the parameters that determine the stop and sbottom masses and mixing angle and are subject to O(α s ) corrections require a one-loop specification. For the stop sector we adopt the OS scheme described e.g. in ref. [13]. In particular, we take as input the pole top mass M t and the soft SUSY-breaking parameters (m Q,t , m U , A t ) that can be derived by rotating the diagonal matrix of the OS stop masses by the angle θ t , defined as in eq. (37) of ref. [13]. Concerning the corresponding parameters of sbottom sector (h b , m Q,b , m D , A b ) additional care is required, because of our non-trivial definition of h b and of the fact that, at O(α s ), the parameter m Q,b entering the sbottom mass matrix differs from the corresponding stop parameter m Q,t by a finite shift [32]. We start by computing the renormalized Higgs-sbottom coupling as given by h b =ĥ b − δh b , whereĥ b is the DR-renormalized running coupling that can be trivially extracted from m b computed via eq. (33), and δh b is defined in eq. (49). Then we compute m Q,b following the prescription of [32]. Finally, we use the parameters h b and m Q,b to compute the actual values of the OS sbottom masses and mixing angle.  To set the stage for further discussion, we show in fig. 4 the masses and mixing angle of the CP-even Higgs bosons as a function of tan β. All the relevant SUSY-breaking parameters, as well as the supersymmetric mass parameter µ, are set to a common value M = 500 GeV, and the physical pseudoscalar mass m A is set to 150 GeV. The left panel of fig. 4 shows the masses of the two Higgs bosons h and H in the range 25 < tan β < 50, while the right panel shows the combinations − sin α/ cos β and cos α/ sin β, which determine the strength of the coupling of h to the bottom and top quarks, respectively, relative to the corresponding SM couplings. For each set of curves, the solid line represents the result obtained in the OS renormalization scheme of ref. [28], described in section 4.3. For comparison, we also show as a dotted line the result that would be obtained if the sbottom parameters h b and A b were renormalized in the same way as the corresponding stop parameters. The left plot shows the well-known fact that, at large tan β, the radiative corrections from sbottom loops tend to reduce m h . The right plot shows that, for the chosen values of m A and tan β, the coupling of h to the bottom quark is still substantially enhanced with respect to its SM value. This has to be contrasted with the couplings of h to the top quark and to the gauge bosons (not shown), which are already very close to the SM values they tend to in the "decoupling" limit m A ≫ m Z .
The comparison between the dotted and solid lines in fig. 4 shows that, if we had adopted for the sbottom parameters h b and A b the renormalization scheme used for the stop parameters, the results for m H and for − sin α/ cos β would differ wildly from the ones obtained with the renormalization scheme discussed in section 4.3 (conversely, we checked that the results obtained in the DR scheme would be in 25   good qualitative agreement with the solid lines). The discrepancy is due to the fact that in the "bad" OS scheme the (1,1) and (1,2) entries of the CP-even Higgs mass matrix are subject to very large two-loop corrections scaling like M 2 tan β 2 , induced by the contribution of the counterterm δA b . It is interesting to note that, since the contribution of δA b to the form factor H 2ℓ 1 in eq. (54) is suppressed by a factor m b , its impact on the Higgs boson production cross section in the "bad" OS scheme is not as extreme as the impact on the Higgs mass. However, we stress that a consistent determination of the properties of the Higgs bosons requires that the same definition of input parameters be used in the calculations of mass and production cross section. Since the naive choice of using the same OS renormalization scheme for the stop and sbottom sectors is not viable in the calculation of the Higgs masses, it should not be applied to the calculation of the cross section either.
We are now ready to discuss the effect of the two-loop bottom/sbottom contributions to the form factor for Higgs-boson production. To this purpose, we define a factor K h that contains the ratio of two-loop to one-loop form factors appearing in eq. (17): In the left panel of fig. 5 we plot K h as a function of tan β, with the same choice of SUSY parameters as in fig. 4, in the OS renormalization scheme described in section 4.3. The one-loop form factors in the denominator of the term between parentheses in eq. (55) contain both the top/stop and bottom/sbottom contributions, computed under the approximations of eqs. (13) and (14). The lines in the plot correspond to different computations of the two-loop form factors in the numerator: the dotted line includes only the contributions of the top/stop sector, as computed in ref. [13]; the dashed line includes also the contribution of two-loop diagrams with bottom quarks and gluons; the solid line includes the full two-loop contribution of the bottom/sbottom sector as computed in section 4; finally, the dot-dashed line is obtained by approximating the bottom/sbottom contribution (with the exception of the bottom-gluon diagrams) with just the terms enhanced by mg/m b in eq. (27). From the comparison between the dotted and dashed lines it can be seen that, in the OS renormalization scheme, the contribution to H 2ℓ 1 of the two-loop diagrams with bottom quarks and gluons is very small. This is due to a partial cancellation between the terms C F F is not enhanced by the potentially large logarithm of the ratio between the bottom mass and the renormalization scale, as can be seen by comparing eqs. (23) and (31). The solid line shows that the effect of the diagrams involving sbottoms can be very sizable at large tan β, more than doubling K h . Indeed, for large tan β the coupling of the light Higgs boson to the (s)bottom is considerably enhanced with respect to its SM value, as can be seen in the right panel of fig. 4. However, the proximity between the solid and dot-dashed lines shows that this sizable effect is almost entirely due to the terms enhanced by mg/m b in the contribution of the two-loop bottom-sbottom-gluino diagrams in which the light Higgs boson couples to the bottom quark.
As discussed in section 4.1, the terms enhanced by mg/m b in the OS result can be canceled out if the Higgs-bottom Yukawa coupling in the one-loop part of the result is identified with the DRrenormalized MSSM bottom mass m b instead of the physical mass M b . To this effect, the factor m 2 b multiplying the function G b in eq. (5) must be expressed as m b M b , and the two-loop part of G b must be shifted as in eq. (32). In the right panel of fig. 5 we present the result of this manipulation, with m b evaluated at the scale Q = m h . The input parameters and the meaning of the different lines are the same as for the plot in the left panel. The proximity between the dashed and solid lines shows that the contribution of the two-loop diagrams involving sbottoms is rather small in this renormalization scheme, at least for our choice of input parameters. However, K h still shows a sizable increase at large tan β. This is due to the fact that the shift in eq. (32) brings back a large logarithm, ln(m 2 b /m 2 h ), in the contribution of the two-loop diagrams with bottom and gluon (this logarithm compensates the scale dependence of the running mass m b ).
To conclude this section, we show in fig. 6 the factor K H for the production of the heavy CP-even Higgs boson, in the range 2 < tan β < 50. The definition of K H can be obtained from the one of K h in eq. (55) via the replacements sin α → − cos α, cos α → sin α. The input parameters are chosen exactly as in figs. 4 and 5, and the meaning of the lines in the left and right panels is the same as in fig. 5. Since in this example the mass of the heavy Higgs boson is of the order of 150 GeV (see fig. 4), i.e. well below any threshold for heavy-particle production, we expect the approximation of vanishing Higgs mass to hold reasonably well even for H. From fig. 6 it appears that the balance of the various contributions to K H in the two different renormalization schemes is qualitatively similar   to the one for K h shown in fig. 5: in the OS scheme the factor K H receives a sizable contribution from the sbottom diagrams, largely dominated by the terms enhanced by mg/m b in the diagrams controlled by the Higgs-bottom coupling; in the "mixed" scheme, on the other hand, the sbottom contribution is rather small, but there is a sizable contribution from the diagrams with bottom and gluon.
This said, the factor K H shows a peculiar dependence on tan β: for sufficiently large values of tan β, it grows linearly in the OS scheme, while it reaches a plateau in the mixed scheme. This can be easily understood by recalling that, for moderate-to-large tan β and for our choice of m A , the Yukawa coupling of the heavy MSSM Higgs to bottom quarks is enhanced by tan β with respect to the SM value, while the coupling to top quarks is suppressed by tan β. Consequently, both the one-loop and the two-loop form factors in K H are dominated by the contribution of the diagrams controlled by the Higgs-bottom coupling, with the result that the coupling itself cancels out in the ratio. However, the dominant contribution from the bottom-sbottom-gluino diagrams in the OS scheme contains an additional tan β-enhancement hidden in the product s 2θ b mg/m b (see the discussion at the end of section 4.1), which explains the linear rise of K H . On the other hand, the dominant contribution of the bottom-gluon diagrams in the mixed scheme possesses no further tan β-enhancement, which explains the plateau.

Conclusions and discussion
In this paper we presented analytic results for the NLO bottom-sbottom-gluino contribution to the cross section for Higgs boson production in gluon fusion, obtained using an asymptotic expansion in the large supersymmetric masses. This approximation is fully valid for the light-Higgs case, while for the heavy Higgs it covers the mass region where m H is below all the heavy-particle thresholds. Together with the previously known results for the NLO corrections in the MSSM, our expressions can be easily implemented in computer codes that aim to provide an accurate and efficient evaluation of the cross section for Higgs boson production in the MSSM.
In our analysis we paid special attention to the consistency between the calculations of the masses and the production cross sections of the MSSM Higgs bosons, i.e. to the fact that the same input parameters, defined at the one-loop level, should be used in both calculations. The OS definition of the parameters of the bottom sector is delicate, as discussed for the case of the Higgs masses in ref. [28]. The choice of treating the top and the bottom sectors on the same footing suffers from the fact that large two-loop corrections proportional to tan 2 β are generated in the contributions controlled by the Higgs-sbottom couplings, affecting both the calculation of the Higgs masses and that of the production cross sections. To avoid such large two-loop effects, a convenient OS renormalization prescription was proposed in ref. [28] for the calculation of the Higgs masses. In the present paper we have extended that prescription to cover also the calculation of the production cross sections.
Our analysis of the NLO bottom contribution to the gluon-fusion production cross section shows that, with our choice of OS renormalization conditions, the bulk of the corrections comes from the two-loop diagrams involving the Higgs-bottom Yukawa coupling, while the diagrams controlled by the Higgs-sbottom coupling play a secondary role. The contribution controlled by the Higgs-bottom Yukawa coupling can be further divided in two parts: diagrams with only bottom and gluons and diagrams involving bottom, sbottom and gluino. By far, the most important pieces of the latter diagrams are the terms of O(m b mg/m 2 h,H ), i.e. the ones in which the helicity flip on the fermion line is achieved via a gluino mass insertion instead of a bottom mass insertion.
It is natural to wonder if it is possible to absorb most of the NLO bottom contribution into the LO term with a suitable choice of the input parameters. In such a situation the factor K h,H , as defined in eq. h,H ) terms in the two-loop contribution cancel out, but the bottom-gluon diagrams give a relevant contribution because of the presence of large logarithms of the ratio between m b and the renormalization scale. However, the explicit knowledge of the NLO bottom contribution allows us to devise a simple recipe to absorb the bulk of the NLO contribution into the LO term. It amounts to writing the LO bottom contribution entirely in terms of the pole bottom mass M b , then rescaling it by a factor 1/(1 + ǫ b tan β). Once this manipulation is implemented, we expect the remaining NLO bottom/sbottom contributions to be quite small -at least in large regions of the parameter space -in which case they can be neglected in the evaluation of the form factors H 2ℓ 1,2 without introducing large errors 4 . We stress that the validity of this simple recipe is strictly linked to the absence of spuriously large corrections to the Higgs-sbottom coupling. This is realized with our choice of OS renormalization conditions for the sbottom sector (and also in the DR scheme) but it is not guaranteed with other renormalization conditions. Finally, the results derived in this paper for the production cross section can be straightforwardly applied to the NLO computation of the gluonic and photonic decay widths of the CP-even Higgs boson in the MSSM, as described in section 5 of ref. [13].