Anomalous Higgs boson couplings in weak boson fusion production at NNLO in QCD

The production of Higgs bosons in weak boson fusion has the second largest cross section among Higgs-production processes at the LHC. As such, this process plays an important role in detailed studies of Higgs interactions with vector bosons. In this paper we extend the available description of Higgs boson production in weak boson fusion by considering anomalous $HVV$ interactions and NNLO QCD radiative corrections at the same time. We find that, while leading order QCD predictions are too uncertain to allow for detailed studies of the anomalous couplings, NLO QCD results are sufficiently precise, most of the time. The NNLO QCD corrections alter the NLO QCD predictions only marginally, but their availability enhances the credibility of conclusions based on NLO QCD computations.


I. INTRODUCTION
The discovery of the Higgs boson nearly ten years ago completed the Standard Model (SM) of particle physics, providing, for the first time, experimental support for the hypothesis that electroweak symmetry is broken by a scalar field. By now quantum numbers of the discovered Higgs boson, as well as its couplings to gauge and matter fields, have been studied in great detail and it appears that the properties of the Higgs boson are closely aligned with the SM expectations. For example, the Higgs couplings to massive electroweak vector bosons have been measured to within O(30) percent of their Standard Model values (see e.g. Refs [1,2]).
Verifying the Higgs boson couplings to a better, perhaps a few percent, precision is the goal of the LHC Run III and, especially, of the high-luminosity LHC. Reaching this goal will not be easy as it will require combining precise measurements with very detailed theoretical predictions. In addition, assuming that deviations in Higgs couplings are found, understanding their origin and what they imply becomes important. A convenient way to investigate this in a relatively model-independent way is to use an Effective Field Theory (EFT) framework and parameterize deviations from the Standard Model by higher-dimensional operators that depend on the Standard Model fields. These operators modify interactions * Electronic address: kasteriad@bnl.gov † Electronic address: fabrizio.caola@physics.ox.ac.uk ‡ Electronic address: kirill.melnikov@kit.edu § Electronic address: raoul.rontsch@unimi.it between Standard Model particles and affect production cross sections and kinematic distributions that are observed at the LHC. An EFT that extends the SM is known as SMEFT [3][4][5].
Assuming that SMEFT provides a faithful description of beyond the Standard Model (BSM) physics, it becomes important to extract Wilson coefficients of higherdimensional operators from the experimental data. Similar to measurements of Standard Model parameters, such an extraction of Wilson coefficients is affected by radiative corrections, of which QCD corrections are especially important. In this paper we investigate the interplay between anomalous couplings and radiative corrections for Higgs production in weak boson fusion (WBF). Such an interplay can be quite subtle for this process. Indeed, corrections to inclusive Higgs production in WBF are known to be small, but their impact on kinematic distributions can be larger. Anomalous couplings also distort kinematic distributions and, in fact, it is known that shapes of various observables often provide the best means to distinguish between different anomalous couplings (see e.g. Refs [6][7][8]). Hence, for a reliable EFT analysis it becomes important to understand to what extent the effects of anomalous couplings and QCD corrections can be disentangled.
The theoretical description of Higgs production in weak boson fusion is very advanced. At the inclusive level, N 3 LO QCD corrections are known in the so-called factorized approximation [9]. At the differential level, factorized NNLO QCD corrections [10][11][12] as well as NLO EW corrections [13][14][15]  Theoretical predictions that include the anomalous couplings together with the NLO corrections to Higgs boson production in weak boson fusion can be obtained using such programs as HAWK [20], VBFNLO [8,21] and Mad-graph5 [22].
In this paper we extend these results by computing NNLO QCD corrections to weak boson fusion in the factorized approximation in the presence of anomalous HV V couplings. Although the HV V vertex is not affected by QCD effects, such a computation is non-trivial since obtaining fully-differential predictions for the WBF process with NNLO QCD accuracy [10][11][12] is quite demanding. The reason behind this is the relative smallness of radiative corrections and a very large phase space (typical for a 2 → 5 process) whose efficient generation is challenging.
In this paper, we employ the calculation of NNLO QCD corrections reported in Ref. [12]. It is performed using the nested soft-collinear subtraction scheme [23] which was adapted to WBF kinematics in Ref. [24]. It provides a sufficiently efficient implementation and phasespace sampling that allowed us to incorporate decays of the Higgs boson [12] into the calculation. This efficiency is quite important for studying the anomalous couplings, especially for performing scans in the parameter space.
The rest of the paper is organized as follows. In Section II we describe the EFT framework which we employed in the calculations reported in this paper, and discuss the parametrization of the anomalous couplings adopted for this analysis. In Section III we present phenomenological results for the 13 TeV LHC. After describing our setup, we show predictions for cross sections and scans in the space of the anomalous couplings. We then focus on two scenarios where anomalous interactions are present but fiducial cross sections are indistinguishable from their Standard Model values. We then consider kinematic distributions which are sensitive to the anomalous couplings of the Higgs boson, and discuss the impact of higher order QCD corrections to them. We conclude in Section IV.

II. ANOMALOUS HVV INTERACTIONS
For our analysis, we focus on anomalous HV V interactions, where V is a massive vector boson. 1 Lorentz invariance and Bose symmetry imply that the most generic 1 For definiteness, we do not consider photon-mediated WBF in our study. We stress however that its inclusion in our framework is straightforward.
where A,B and C are arbitrary functions which, in the context of an effective theory, are Taylor-expandable in p 2 i and p i · p j (i, j = 1, 2). It turns out that the dependence on these kinematic invariants can be further restricted. This point has been studied quite recently through an all-order EFT expansion in the context of the Standard Model Effective Theory [25]. In case of the dominant dimension-six operators, it has been shown [25] that the functions A,B and C can be written as follows where Λ is the EFT scale and g where g is the weak coupling and θ W is the weak mixing angle. 3 We find it convenient to parametrize the anomalous couplings as follows The new constants c (i) HV V are dimensionless; in terms of these constants, the coupling of the Higgs boson and the 2 We note that, in principle, tensor couplings of the form p µ 1 p ν 2 are allowed in Eq. (1). However, neglecting contributions of heavy quarks in q → qV and similar transitions, the vector bosons couple to conserved currents in weak boson fusion, so that such couplings do not contribute to the final result. 3 We note that the g (i) HV V couplings introduced in Eq. (2) are defined exactly as in Ref. [26]. In this reference, an additional coefficient g (2) HV V is also present. However, this Wilson coefficient is not generated at dimension six [25]. For this reason we do not consider it in our analysis, although it is straightforward to include it if needed. electroweak vector bosons reads, c.f. Eq. (1), The factor (6π) introduced in the CP -odd term is chosen for convenience; indeed, as we will see later, with this normalization and for Λ = 1 TeV,c HV V ∼ 1 induces O(1%) corrections to the fiducial cross section.
In general, the HZZ and HW W anomalous couplings do not have to be the same. However, for simplicity, in the analysis below we will assume that, by factorizing the SM couplings, we already account for main differences between them. Hence, in what follows, we restrict ourselves to the case where c HV V and c HZZ =c HW W =c HV V .

Setup
For the phenomenological results reported in this paper, we use the same parameters and kinematic selection criteria as in Ref. [24]. We collect them here for completeness. We consider proton-proton collisions with the center-ofmass energy 13 TeV and treat the Higgs boson as stable. We use the Higgs mass m H = 125 GeV, the vector boson masses M W = 80.398 GeV and M Z = 91.1876 GeV, and their widths Γ W = 2.105 GeV and Γ Z = 2.4952 GeV. We use the Fermi constant G F = 1.16639×10 −5 GeV −2 to derive the weak couplings and we set the CKM matrix to an identity matrix. We employ the NNPDF31-nnlo-as-118 parton distribution functions [27] and use them for all calculations reported in this paper irrespective of the nominal perturbative order. We also use α s (M Z ) = 0.118. The evolution of the parton distribution functions and the strong coupling constant is obtained from LHAPDF [28]. Finally, we employ dynamical renormalization and factorization scales µ R = µ F = µ using a central value [10] To define the weak boson fusion fiducial volume we reconstruct jets using the inclusive anti-k ⊥ algorithm [29] Table I. Fiducial cross sections, in fb, at various orders of perturbative QCD using µ = µ0. The sub-and super-scripts indicate the results computed with µ = µ0/2 and µ = 2µ0, respectively. Numerical uncertainties are much smaller than scale uncertainties and are not shown. See text for further details.
with R = 0.4. We require events to contain at least two jets with transverse momenta p ⊥,j > 25 GeV and rapidities |y j | < 4.5. Also, the two leading-p ⊥ jets should be separated by a large rapidity interval |y j1 − y j2 | > 4. 5 and their invariant mass should be larger than 600 GeV.
In addition, the two leading jets should be in different hemispheres in the laboratory frame; to enforce this, we require that the product of their rapidities in the laboratory frame is negative, y j1 y j2 < 0. Finally, for definiteness we set Λ = 1 TeV in all computations that we report below.

Fiducial cross sections
We now present our results for the cross sections for Higgs production in weak boson fusion in the fiducial region described above, through NNLO QCD. Since three independent couplings appear in the HV V vertex, cf. Eq. (4), the cross section naturally separates into six terms We note that X 1 is the Standard Model cross section. Two of the other terms, X 5,6 describe the interference of CP -even and CP -odd couplings and hence vanish for the azimuthally-symmetric cuts that we employ for computing fiducial cross sections. 4 The results for the non-vanishing coefficients at leading, next-to-leading and next-to-next-to-leading orders are shown in Table I.
We note that the radiative corrections to the various cross sections X 1,2,3,4 are similar but not identical. In general, when moving from LO to NLO we also observe a significant reduction in the dependence of the cross sections X 1,2,3,4 on the renormalization and factorization scales, which we estimate by changing the central scale µ 0 in Eq. (5) by a factor two in either direction. Instead, there is no substantial scale-uncertainty reduction when moving from NLO to NNLO, and the NLO scale variation bands do not contain the NNLO result. This feature is well-known for the SM case, where it is also known that the NNLO and N 3 LO scale variation bands of the inclusive cross section do overlap [9]. Because of this, it is understood in the context of SM studies that drawing conclusions based only on NLO QCD predictions and their scale uncertainty is delicate, while NNLO analyses should be more robust. Table I suggests that this also holds in the presence of anomalous HV V interactions. These features will play an important role in the discussion below.
To illustrate the importance of higher order corrections, we now study the dependence of fiducial cross sections on the anomalous coupling constants at different orders of QCD perturbation theory. To this end, we take the Standard Model cross section as a reference point. We then consider a hypothetical BSM scenario with nonzero We show the results of this analysis in Figs 1,2 where different pairs of the anomalous couplings are chosen for two-dimensional projections. In those figures, the colored regions represent the allowed values for anomalous couplings according to the criterion just explained. We note that for this analysis we take Eq. (6) exactly as written, i.e. we include the contribution of dimension-six operators squared. One may worry about the consistency of this approach since we are neglecting the contribution of dimension-eight operators when defining the anomalous couplings. While this worry is clearly justified, the goal of our analysis is to study the possible interplay between higher order QCD corrections and EFT effects, and for this dealing with a representative set of anomalous couplings is sufficient. We stress however that the results of Table I make it straightforward for anyone to repeat the analysis including only the interference of dimension-six and SM operators. Similarly, it is immediate to change the definition of compatibility between SM and BSM predictions.
Having clarified our procedure, we now comment on the results, c.f. Figs 1,2. We note that since the cross sections X 2,3,4 describing the contributions of the anomalous couplings change the Standard Model cross section by about 2 − 3 percent and the scale uncertainty of the leading order Standard Model cross section is about 10 percent, only relatively large anomalous couplings can be excluded at leading order. The situation changes dramatically at NLO since theoretical uncertainties are reduced to about 1 − 2 percent. There are further improvements at NNLO but they are minor; such improvements would imply changes in the anomalous couplings at the level of δc (1,2) HV V ∼ 0.1 which, depending on the value of c HV V can be about ten percent. However, as we noted above, the NNLO shifts in the cross sections shown in Table I are often outside the NLO scale uncertainty bands. We then take the consistency of the NLO and NNLO results as an indicator of the reliability of the NLO analysis.
It is clear that the above study can be made more complex by, for example, considering three couplings at the same time, or by comparing results that include the contribution of dimension-six operators squared and results that do not, or by defining a more refined estimate of the compatibility of SM and BSM results. We do not pursue these investigations in this paper but we point out that once the fiducial cross sections X 1,2,3,4 are known, defining selection criteria and scanning the parameter space become straightforward. Obviously, if the fiducial cuts are to be changed, then all the cross sections X 1,2,3,4 have to be recalculated. We note in this respect that it takes about 25 000 CPU hours to compute each X i with good precision. Hence, if O(1000) CPU's are available, about a few days of running is required.

Kinematic distributions
In the previous section, we discussed fiducial cross sections. This, of course, does not exhaust all opportunities for the analysis as kinematic distributions can also help to distinguish between the different couplings. To illustrate this point and to show how higher-order QCD corrections impact such an analysis, we consider two sets of  anomalous couplings which lead to nearly identical cross sections. We choose the following scenarios The results for fiducial cross sections at various orders of QCD perturbation theory are shown in Table II. It is clear that, at leading order, the two scenarios cannot be distinguished from each other and from the Standard Model. Even with the significant reduction of scale variation uncertainties at next-to-leading and next-to-nextto-leading orders, the fiducial cross sections remain compatible.
To understand if the two scenarios and the Standard Model can still be distinguished from each other, we consider the kinematic distributions. For most observables, however, there is very little difference between the various scenarios even at higher orders. To give an example, in Fig. 3 we consider the Higgs transverse momentum distribution. We report predictions at LO (left), NLO (middle) and NNLO (right), for the scenarios A and B and the SM. It is clear from this figure that one cannot disentangle the different models using LO predictions. At NLO and NNLO one starts seeing hints of slightly different shapes, but this is a mild effect and the residual scale variation uncertainty is still too large to allow any definite conclusion to be drawn. We note that the predictions for the three models start deviating significantly at large transverse momenta p ⊥,H 250 GeV. This is expected as the EFT effects are enhanced at higher energies. However, in this region the cross section is already quite small, having dropped off by about one order of magnitude from its peak value.
On the other hand, it is well-known that angular variables are especially well-suited for discriminating between the SM and the various EFT scenarios [6][7][8]. To this end, we study the azimuthal separation between the two hardest jets ∆φ j1j2 = φ j1 −φ j2 , where j 1 is the forward jet and j 2  the backward jet, i.e. y j1 > 0 and y j2 < 0. We show this distribution in Fig. 4. It follows from this figure that, at LO, while the differences between the SM and scenario A is just covered by scale variation uncertainties, scenario B is clearly distinguishable from the SM as well as from scenario A. The situation improves at NLO. Indeed, in this case the scale uncertainty is significantly reduced and the three scenarios become clearly distinguishable. The situation remains the same at NNLO. We note that similar to the situation with the fiducial cross section, there is no significant reduction in NNLO scale variation uncertainties with respect to NLO. It can also be seen in Fig. 4 that the NLO shapes are quite stable under radiative corrections. This is a welcome feature, as it implies that the perturbative expansion for this observable is under very good control and NLO results are reliable. To quantify this statement, in Fig. 5 we plot the ratio of LO and NNLO to NLO predictions, for both the Standard Model (left), scenario A (middle) and scenario B (right). We see that in all the three cases there is no overlap between the scale uncertainty bands, but the corrections are rather flat and can be captured by a global K-factor. Interestingly, both the very mild shape distortion and the K-factor seem independent on the value of the anomalous couplings, at least for the reference values chosen here. These results imply that NLO predictions, possibly augmented by a global K-factor rescaling, seem to provide a robust enough framework for performing anomalous coupling studies with this observable.
Differences between scenario A, scenario B and the SM in the distributions shown in Fig. 4 are primarily of an anti-symmetric nature since they are dominated by the interference of CP -odd and CP -even couplings, which are absent in azimuthally-averaged observables. Hence, in the presence of non-vanishing CP -odd EFT operators, this observable may not be optimal to study effects of CP -even couplings. To this end, it is useful to consider the absolute value of the azimuthal separation, in which the interference between CP -odd and CP -even operators drops out [6][7][8]. We show this distribution in Fig. 6.
As expected, the deviations of scenarios A and B from the SM for |∆φ j1j2 | in Fig. 6 are much smaller than those observed in Fig. 4. In fact, at LO, the differences between the various scenarios are entirely swamped by scale variation uncertainties. However, despite being less sensitive to non-SM couplings, already at NLO scale uncertainties are sufficiently reduced and the three scenarios become distinguishable. Different shapes of scenarios A and B, due to different CP -even contributions, can now be clearly observed. Ratios in Fig. 7 imply that these shapes are also relatively stable under radiative corrections and might, therefore, be captured by global K-factors (as was the case for the ∆φ j1j2 observable).

IV. CONCLUSIONS
We have studied the combined impact of anomalous Higgs interactions and NNLO QCD corrections on Higgs production in weak boson fusion. For simplicity, we have only considered modifications to the HV V vertex and have performed phenomenological analyses assuming that the modification of the HW W and HZZ couplings are correlated. We have found that the patterns of radiative corrections for the WBF cross section in a typical experimental fiducial region are similar in the Standard Model and when anomalous couplings are present. Namely, in all cases scale variation uncertainties underestimate the size of the correction at the next perturbative order, but predictions seem to be reasonably stable when moving from NLO to NNLO. Keeping in mind that for inclusive SM cross sections the NNLO scale uncer-tainty band does overlap with N 3 LO QCD predictions, we have investigated the effect of a reduced theoretical uncertainty on extractions of anomalous couplings. We have studied the constraining power of the fiducial cross section, and found that NLO and NNLO results lead to a similar discriminating power, which is significantly better than the LO one. The relative stability of this picture, when going from NLO to NNLO, gives us confidence that for the analysis of the anomalous couplings in Higgs production through WBF, QCD radiative corrections are sufficiently well understood.
We have also investigated the constraining power of differential distributions, focusing on scenarios where different choices of anomalous couplings led to the same fiducial cross section. It is known that angular distributions have a strong constraining power. However, including NLO QCD corrections is mandatory for differences in shapes induced by small anomalous couplings not to be swamped by theoretical uncertainties. We have shown that theoretical predictions for differential distributions are reasonably stable when moving from NLO to NNLO. In particular, we have found that shape distortions are quite small, and the bulk of the effect of NNLO QCD corrections seems to be captured by a global K-factor. We have found that this holds irrespective of the presence or absence of anomalous couplings, at least for the scenarios studied in this paper.
There are several ways in which the analysis presented in this paper can be improved. First, one could study more general scenarios where modifications of the HZZ and HW W couplings are different. Also, one may study the impact of dimension-eight operators in HV V couplings on Higgs signal in WBF. Implementing both of these extensions in our framework is conceptually straightforward. A more involved improvement would be to con-sider a wider class of higher dimensional operators, especially those that are affected by higher-order QCD corrections. Also, one may consider more realistic setups, where e.g. the decay of the Higgs boson is also taken into account. We look forward to investigating these and other issues relevant for EFT studies in weak boson fusion in the future.