Combination of CMS searches for heavy resonances decaying to pairs of bosons or leptons

A statistical combination of searches for heavy resonances decaying to pairs of bosons or leptons is presented. The data correspond to an integrated luminosity of 35.9 fb − 1 collected during 2016 by the CMS experiment at the LHC in proton-proton collisions at a center-of-mass energy of 13 TeV. The data are found to be consistent with expectations from the standard model background. Exclusion limits are set in the context of models of spin-1 heavy vector triplets and of spin-2 bulk gravitons. For mass- degenerate W (cid:2) and Z (cid:2) resonances that predominantly couple to the standard model gauge bosons, the mass exclusion at 95% conﬁdence level of heavy vector bosons is extended to 4.5 TeV as compared to 3.8 TeV determined from the best individual channel. This excluded mass increases to 5.0 TeV if the resonances couple predominantly to fermions. © 2019 The Author(s). Published by Elsevier B.V. This is an open access article under the CC BY license (http://creativecommons.org/licenses/by/4.0/). Funded by SCOAP 3 .

This Letter describes a statistical combination of CMS searches for heavy resonances that decay to pairs of bosons or leptons [1- 10,19,20]. The event selection, the simulated samples, the background estimation, and the systematic uncertainties of the individual analyses are unchanged. The results are used to set exclusion E-mail address: cms -publication -committee -chair @cern .ch. limits on models that invoke heavy vector resonances and on a model with a bulk graviton.
The SM bosons produced in the decay of resonances X with masses m X that exceed 1 TeV are expected to have a large Lorentz boost [44]. The decay products of the SM bosons are therefore highly collimated, requiring dedicated techniques for their identification and reconstruction. In the case of hadronic decays, the pair of quarks produce a single large-cone jet with a two-prong structure. In addition, Higgs bosons can also be identified by tagging the b quarks originating from their decays. In models where heavy vector bosons couple predominantly to fermions, the combined contribution from the leptonic decays W → ν and Z → dominates, and these channels [19,20] are included in the combination, providing complementary sensitivity to the searches in diboson channels. The analyses considered in this Letter are based on a data sample of proton-proton (pp) collisions at a center-ofmass energy of 1 TeV collected by the CMS experiment in 2016, and corresponding to an integrated luminosity of 35.9 fb −1 . A similar combination performed on a comparable set of data has been recently published by ATLAS [45].

The CMS detector and event reconstruction
The CMS detector features a silicon pixel and strip tracker, a lead tungstate crystal electromagnetic calorimeter (ECAL), and a brass and scintillator hadron calorimeter (HCAL), each composed of a barrel and two endcap sections. These detectors reside within a superconducting solenoid, which provides a magnetic field of 3.8 T. Forward calorimeters extend the pseudorapidity coverage up to |η| < 5.2. Muons are measured in gas-ionization detectors embedded in the steel flux-return yoke outside the solenoid. A detailed description of the CMS detector, together with a definition of the coordinate system and the kinematic variables, can be found in Ref. [46].
The information from various elements of the CMS detector is used by a particle-flow (PF) algorithm [47] to identify stable particles reconstructed in the detector as electrons, muons, photons, and charged or neutral hadrons. The energy of electrons is determined from a combination of the electron momentum, as determined in the tracker, the energy of the corresponding ECAL cluster, and the energy sum of all bremsstrahlung photons spatially compatible with originating from the electron track. The dielectron mass resolution for Z → ee decays ranges from 1.9% when both electrons are in the ECAL barrel, to 2.9% when both electrons are in the endcaps [48]. The momentum of muons is obtained from the curvature of the corresponding track. The p T resolution in the barrel is better than 7% for muons with p T up to 1 TeV [49]. The τ leptons that decay to hadrons and a neutrino (labeled τ h ) are reconstructed by combining one or three charged PF candidates with up to two neutral pion candidates [50]. The e, μ, and τ h with large momenta are required to be isolated from other hadronic activity in the event. The energy of charged hadrons is determined from a combination of their momenta measured in the tracker and the matching ECAL and HCAL energy deposits. Finally, the energy of neutral hadrons is obtained from the corresponding corrected ECAL and HCAL energies.
Jets are reconstructed from PF candidates clustered with the anti-k T algorithm [51], with a distance parameter 0.4 (AK4 jets) or 0.8 (AK8 jets), using the FastJet 3.0 package [52,53]. The geometrical overlaps between jets and isolated leptons, and between AK4 jets and the AK8 jet, are handled by assigning the overlapping four-momentum to the lepton, and the AK8 jet, respectively. The four-momenta of the AK4 and AK8 jets are obtained by clustering candidates passing the charged-hadron subtraction algorithm [54]. The contribution of neutral particles originating from additional pp interactions within the same or neighboring bunch crossings (pileup) is proportional to the jet area, and is estimated using the median area method implemented in FastJet, and then subtracted from the jet energy. The jet energy resolution, after the application of corrections to the jet energy, is 4% at 1 TeV [55]. The missing transverse momentum vector p miss T is computed as the negative vector sum of the transverse momenta of all the PF candidates in an event, and its magnitude is denoted as p miss T [56]. The p miss T is corrected for adjustments to the energy scale of the reconstructed AK4 jets in the event.
While AK4 jets are used for single quarks, AK8 jets are adopted to reconstruct large momentum SM bosons that decay to quarks. Reconstructing the AK8 jet mass (m j ) and substructure relies on the pileup-per-particle identification algorithm [54,57]. The contributions from soft radiation and additional interactions are removed using the soft-drop algorithm [58,59], with parameters β = 0 and z cut = 0.1. Dedicated mass corrections, obtained from simulation and data in a region enriched with tt events with merged W(qq ) decays, are applied to the jet mass to reduce any residual jet p T dependence [6,60], and to match the jet mass scale and resolution observed in data. The measured soft-drop jet mass resolution is approximately constant at 10% in the considered range of jet p T [60]. Exclusive m j intervals labeled m W , m Z , and m H , which range from 65 to 85, 85 to 105, and 105 to 135 GeV, respectively, are defined according to the SM boson masses and the jet mass resolution.
The variables m j and τ 21 are calibrated using a top quark-antiquark sample enriched in hadronically decaying W bosons [62]. The decay of a Higgs boson to a pair of b quarks is identified using one of two b tagging algorithms, depending on the background composition. The first consists of a dedicated b tagging discriminator, specifically designed to identify a pair of b quarks clustered in a single jet [63]. The second relies on splitting the AK8 jet into two subjets, then applying the combined secondary vertex algorithm [63] to the subjets. The latter is also applied to AK4 jets to identify isolated b quarks in the event.

Signal modeling
The response of the CMS detector to the production and decay of heavy resonances is evaluated through simulated events, which are reconstructed using the same algorithms as used in collision data. The spin-1 gauge bosons, W and Z , are simulated at leading order (LO) using the MadGraph5_amc@nlo 2.4.2 matrix element generator [64], in the heavy vector triplet (HVT) framework [65,66], which introduces a triplet of heavy vector bosons with similar mass, of which one is neutral (Z ) and two are electrically charged (W ± ). The coupling strength of the heavy vector bosons to SM bosons and fermions is determined from the combinations g H = g V c H and g f = g 2 c F /g V , respectively. The parameter g V is the strength of the new interaction; c H characterizes the interaction between the HVT bosons, the Higgs boson, and longitudinally polarized SM vector bosons; c F represents the direct interaction between the V bosons and the SM fermions; g is the SM SU(2) L gauge coupling constant. The HVT framework is presented in two scenarios, henceforth referred to as model A and model B, depending on the couplings to the SM particles [66]. In the former, the coupling strengths to the SM bosons and fermions are comparable, and the new particles decay primarily to fermions. In the latter, the couplings to the SM fermions are small, and the branching fraction to the SM bosons is nearly 100%. Events are simulated with different m X hypotheses from 800 to 4500 GeV for X decays to bosons, and 800 to 5500 GeV for X decays to leptons, assuming a negligible resonance width compared to the experimental resolution (3-20%). The kinematic distributions of the signal do not depend on the choice of benchmark scenario, and the samples are reweighted according only to the cross sections and branching fractions.
The bulk graviton events are also simulated at LO using the same MadGraph5_amc@nlo generator. The cross section and the width of the bulk graviton mainly depend on its mass and the ratio κ ≡ κ/M Pl , where κ is a curvature factor of the model and M Pl is the reduced Planck mass [43,67]. The graviton signals are generated assuming κ = 0.5, which guarantees that the width of the graviton is smaller than the experimental resolution.
The signal events are generated using the NNPDF 3.0 [68] parton distribution functions (PDFs), and are interfaced to pythia 8.205 [69] for parton showering and hadronization, adopting the MLM matching scheme [70]. Simulated pileup interactions are superimposed on the hard process, and their frequency distribution is weighted to match the number of interactions per bunch crossing observed in data. Generated events are processed through the CMS detector simulation, based on Geant4 [71].

Search channels
The search strategies for the analyses contributing to the combination are summarized in this Section. More details are provided in the relevant publications [1-10,19,20].

Fully hadronic diboson channels
Diboson resonances have been searched for in several final states, depending on the decay modes of the bosons. The final states targeted were VV → qqqq [1], VH → qqbb [6], and HH → bbbb final states [9,10]. Each boson is reconstructed as a twoprong, large-cone jet, so that for diboson resonances, the two jets would recoil against each other. The presence of a diboson resonance would be observed in the dijet invariant mass spectrum m jj .
The W and Z bosons are identified (V tagged) via the m j and τ 21 variables, and b tagging is used to identify b quarks from Higgs bosons in addition to m j . Although the signal yield is large, because of the large fraction of Higgs boson decays to b quarks, these channels are subject to an overwhelming background from quantum chromodynamics (QCD) multijet production.
In the VV and VH analyses, the background is estimated directly from data, assuming that the invariant mass distribution of the background can be described by a smooth, parameterizable, monotonically decreasing function of m jj . The signal template, based on a Gaussian core, is fitted to the data simultaneously with the background function. The HH analyses also use an additional region in the fit, obtained by inverting the b tagging selection on the H candidates, which constrains the parameters of the background function.

Semi-leptonic diboson channels
Searches for VV, VH, and HH resonances have been performed in channels where one of the SM bosons decays to leptons and the other to quarks (ZV → ννqq [2], WV → νqq [3], ZV → qq [4], ZH → ννbb [7], WH → νbb [7], ZH → bb [7], VH → qqτ τ [8] and HH → bbτ τ [8]). These final states represent an attractive alternative to all-jet final states, thanks to the large selection efficiencies and natural discrimination against multijet background stemming from the presence in the signal of energetic and isolated leptons or neutrinos.
The decay of a Z boson to neutrinos can be identified through its large p miss T , and the resonance mass can be inferred from the transverse mass between the p miss T and the jet originating from the hadronic decay of the other boson. For the W → ν decay, there is a single, isolated lepton associated with a moderate p miss T , and the vector boson can therefore be reconstructed by imposing a constraint from the W boson mass to recover the longitudinal momentum of the neutrino. In Z → decays, two opposite-sign, same-flavor leptons with a combined invariant mass compatible with the Z boson mass are used to determine accurately the Z boson four-momentum. A Higgs boson decaying to τ leptons is identified from dedicated τ h decay reconstruction and isolation techniques [8], and its mass is estimated through the measured momenta of the visible decay products of the τ leptons and the p miss T , with the SVfit algorithm [72]. The boson that decays to a pair of quarks is reconstructed as an AK8 jet. The AK8 jet mass is required to be compatible with either the W and Z boson masses, denoted as m V (65 < m j < 105 GeV), or the Higgs boson mass (105 < m j < 135 GeV).
In the semi-leptonic analyses, the main V + jets background is estimated from a fit to data in the m j sidebands of the hadronic jet, and extrapolated to the signal region using a transfer function ("α function") obtained from simulation. The top quark pair production is estimated from simulation, but its normalization is rescaled to match the data in control regions obtained by requiring an additional b-tagged AK4 jet in the event [2,7,8].
In addition, the WV → νqq analysis introduces a novel signal extraction method based on a two-dimensional (2D) fit to data [3]. The backgrounds are separated into non-resonant and resonant categories depending on the presence or absence of W bosons and top quarks in the jet mass spectrum, and are fitted simultaneously in the space of m j and m WV , accounting for the correlation between the two variables.

Fully leptonic diboson channels
Searches for diboson resonances decaying to a pair of Z bosons have been performed in fully leptonic final states, with one boson undergoing the decay Z → and the other Z → νν [5]. The presence of the leptons and neutrinos defines a very clean final state with reduced backgrounds, but the small branching fraction makes this channel competitive only for small resonant mass values.

Decays to a pair of fermions
The decay of a heavy resonance to a pair of fermions can be sizable when the couplings to SM fermions are large. If the resonance is electrically charged, as in the case of W ± , the decay to a neutrino and an electron or muon yields a broad excess in the ν transverse mass spectrum [19]. If the new state is neutral, as in the case of Z , a narrow resonance would emerge from the dielectron or dimuon ( ) invariant mass spectra [20]. The analyses of these fermionic decays extend to masses above 5 TeV, and employ selection techniques optimized to identify and measure very energetic electrons and muons in the detector [19,20].

Event selection
The search regions of the analyses entering in the combination are statistically independent because of mutually exclusive selections on the number of leptons and their flavor, number of AK8 jets, and jet mass intervals. Analyses with hadronic final states reject events with isolated leptons or with large p miss T .
Overlaps between channels that share the same lepton multiplicity are avoided by selecting different jet mass ranges. The W → ν search does not share any events with the W → VW and W → WH analyses because of the requirements on the angular separation φ( , p miss T ) between the p miss T and the lepton direction. The Z → analysis includes events with a dilepton invariant mass m > 120 GeV, which is incompatible with the 70 < m < 110 GeV selection used in the diboson channels, in which the Z boson is on-shell. The two searches for resonant HH bosons that decay to b quarks have common events explicitly removed [10].
In the WV → νqq channel [3] the background is estimated using a 2D fitting technique that scans the full jet mass range, and is therefore not independent of the WH → νbb channel. For this reason, in the W , Z , and V interpretations, where the two signals might be present simultaneously, the "α function" is used to estimate the background instead. This method considers only events in the m j regions of the W and Z bosons, thereby preventing double counting of events in the Higgs boson mass region. The results of the alternative background estimation method are consistent with those obtained in the 2D fit, but the method is about 10% less sensitive. The main selections that define the exclusivity of the analyses are summarized in Table 1.

Systematic uncertainties
The systematic uncertainties originating from the background estimation in the individual channels are considered uncorrelated, because the backgrounds are determined from statistically independent control regions. The uncertainties arising from the reconstruction and calibration are instead correlated among the different Table 1 Summary of the main selections that guarantee the exclusivity between individual final states. The symbol represents an electron or a muon; τ leptons are considered separately. If an entry is not preceded by "≥", the exact value is required. The AK4 b jets are additional b tagged AK4 jets that do not geometrically overlap with AK8 jets. The symbol "-" implies that no selection is applied. Table 2 Summary of the main systematic uncertainties. The second column reports whether a systematic uncertainty is considered fully correlated or not across different channels. The third column indicates whether the uncertainty affects the yield, the shape of the distributions, or both, or if it induces migration (migr.) effects across search regions. The fourth column reports the smallest and largest effect of the uncertainty, among the yield, the migration, and the signal shape parameters. The symbols "s", "b" indicate that the uncertainty affects the signal, the main backgrounds of the analysis, respectively. The treatment of non-dominant backgrounds is often different and not reported here. The symbol "f" indicates that the parameters are not constrained, or associated with large uncertainties as in the case of multi-dimensional fits. The entries labeled with "t" are treated differently depending on the interpretation of the exclusion limit, as discussed in Section 7. Uncertainties marked with "-" are not applicable or are negligible.

Correlation
Type Variation ννqq [2] νqq νν [5] qqbb [6] (νν, bbbb [9,10] ν [19] [20] Bkg. modeling no shape - channels. These include the uncertainties in jet energy and resolution, and in the e, μ, and τ h lepton energy, reconstruction, and identification. The uncertainties in the identification of the SM bosons that decay to quarks are dominant in final states containing at least one such decay, and originate from the m j scale and resolution, the τ 21 selection, and the b tagging. The uncertainties in the V tagging extrapolation to large jet p T and the selection of events in the jet mass window of the Higgs boson are estimated by comparing the results with those obtained using the alternative herwig++ [73] shower model. These uncertainties affect the signal distribution, selection efficiency, and induce event migration effects between search regions. The uncertainties in the proton-proton inelastic cross section [74], the integrated luminosity during the 2016 data-taking [75] and the kinematic acceptance of final-state particles, which affect the signal normalization, are also considered as correlated among channels. Theoretical uncertainties in the cross section and in the signal geometric acceptance related to the choice of PDFs used in the event generators [76], and the uncertainties in the factorization and renormalization scales, are evaluated according to the PDF4LHC recommendations [76]. The impact of these uncertainties on the signal cross section can be as large as 78%, depending on the signal mass and the initial state (qq or gg). A summary of the main systematic uncertainties is given in Table 2.

Statistical combination
No significant excess is observed above the SM background expectations. Upper limits are set at 95% confidence level (CL) [77,78] on the cross section of a heavy resonance, which is rescaled by a signal strength modifier parameter μ with uniform prior. Systematic uncertainties are represented by nuisance parameters θ , and affect both signal and background expectations [79]. Systematic uncertainties are considered fully (anti-)correlated when related to a common nuisance parameter, or uncorrelated when different nuisance parameters are used. The prior on these parameters is either flat (represented by the symbol "f" in Table 2) or log-normal distributed (identified with "s", "b", "t" in Table 2). The statistical procedure is based on a likelihood constructed as: where P represents the Poisson probability, and p j (θ j |θ j ) is the frequentist pdf of the nuisance parameter θ j and its default value θ j associated with the jth uncertainty. The values s c,i (θ) and b c,i (θ) represent the number of signal and background events in the channel c and bin i, respectively. The test statistic is based on the profile likelihood ratio q: and the quantities μ, θ μ , θμ are fixed to their best-fit value, and the range of μ is limited in the 0 ≤μ ≤ μ interval. The uncertainties that affect the signal normalization (PDFs and factorization and renormalization scales, marked with "t" in Table 2) are treated differently depending on how the exclusion is presented. When deriving upper limits on the cross section, these uncertainties are not varied in the fit, but are reported separately as the uncertainty of the theoretical cross sections from the model. When placing limits on the model parameters, these nuisance parameters are fixed at the best-fit values, in the same manner as to the other systematic uncertainties.
The upper limit on μ is derived from the 95% CL s criterion, defined as CL s (μ) = p(μ)/(1 − p(0)), such that CL s (μ) = 0.05. The quantities p(μ) and 1 − p(0) represent the probabilities to have a value of q equal to, or larger than the observed value in the signal or background (μ = 0) hypotheses, respectively, and are derived from analytical functions using the asymptotic approximation [80]. This approximation leads to limits that are up to 30% stronger in regions with a small number of data events, compared to those obtained from the Monte Carlo generation of pseudoexperiments [79].

Results and interpretation
The data are in agreement with the expected background from SM processes. The largest deviation from the expected limit is observed in the V model A at a mass of 1.3 TeV, with a local significance of 2.7 standard deviations, corresponding to a global significance of 1.6 standard deviations. The local significances are obtained in the asymptotic approximation [80], and the global significances are evaluated using the trial factors method [81].
The exclusion limit on the cross section of each diboson channel (WW, WZ, ZZ, WH, ZH, HH) is depicted in Fig. 1 for the combination of all contributing channels, reported in Table 1, according to the spin of the new resonance. The generated signal can be either a spin-1 heavy vector (W or Z , as in the HVT model) or a spin-2 boson (as in the bulk graviton model). In fact, the spin and polarization of the heavy resonance does affect the final state, signal acceptance, and selection efficiencies. The exclusion limits are not presented above 4.5 TeV, since at larger masses the background estimation procedure used in diboson analyses becomes less reliable because of the lack of events in data.
The combined exclusion limits for the spin-1 singlet hypotheses (W or Z ) in the HVT model B framework, where the branching fractions to SM bosons are dominant, are shown in Fig. 2. In this scenario, the contribution of the dilepton channels is negligible because their branching fraction is of the order of a few permil.
In the VV analyses, the V → VH signal is added to the V → VV signal because the fraction of H jets passing the m W and m Z selections is not negligible and can contribute up to 30% to the total signal yield. The predictions of the HVT model B are superimposed on the exclusion limits, showing that a W boson of mass below 4.3 TeV, and a Z boson with mass below 3.7 TeV are excluded at 95% CL.
The HVT hypothesis is tested in Fig. 3 Fig. 4. The dilepton and diboson searches constrain different regions of the parameter space, as the dilepton searches can probe the region where the coupling to the SM bosons approaches zero. In the triplet interpretation, the ratio of the W to Z cross sections is assumed to be determined by the ratio of the partonic luminosities, and to depend only weakly on the model parameters.
The fraction of the parameter space where the natural width of the resonances is larger than the average experimental resolution of 5%, and the narrow-width approximation is thus invalid, is also indicated in Fig. 4. In the spin-2 bulk graviton model, the WW, ZZ, and HH channels are combined, setting upper limits of up to 1.1 fb on the cross section of a graviton with mass up to 4.5 TeV. In the κ = 0.5 sce- nario, a graviton with a mass smaller than 850 GeV is excluded at 95% CL, as shown in Fig. 5. Larger κ values increase the production cross sections, but also the graviton natural width, which can be comparable or larger than the experimental resolution. In these cases, the narrow-width approximation is no longer valid.
These results represent the most stringent limits on heavy vector models set by CMS and are comparable at large V mass to the limits obtained the ATLAS in the combination of similar channels [45]. At lower masses, the ATLAS combination excludes smaller cross sections because of the inclusion of final states with three or more leptons. The exclusion limits in the bulk graviton model are not directly comparable because of the large κ = 1.0 parameter adopted by ATLAS, which increases the cross section but also implies a non-negligible natural width.   limits of up to 1.1 fb are set at 95% confidence level on the production cross section of the spin-2 bulk graviton. For models with a triplet of narrow spin-1 resonances, heavy vector bosons with masses below 5.0 and 4.5 TeV are excluded at 95% confidence level in models where the W and Z bosons couple predominantly to fermions and bosons, respectively. In the latter, the statistical combination extends the exclusion limit by 700 GeV as compared to the best individual channel.

Summary
In addition, we gratefully acknowledge the computing centers and personnel of the Worldwide LHC Computing Grid for delivering so effectively the computing infrastructure essential to our analyses. Finally, we acknowledge the enduring support for the construction and operation of the LHC and the CMS detector provided by the following funding agencies: