Search for $\mathrm{t}\overline{\mathrm{t}}$H production in the all-jet final state in proton-proton collisions at $\sqrt{s}=$ 13 TeV

A search is presented for the associated production of a Higgs boson with a top quark pair in the all-jet final state. Events containing seven or more jets are selected from a sample of proton-proton collisions at $\sqrt{s} =$ 13 TeV collected with the CMS detector at the LHC in 2016, corresponding to an integrated luminosity of 35.9 fb$^{-1}$. To separate the $\mathrm{t}\overline{\mathrm{t}}$H signal from the irreducible $\mathrm{t}\overline{\mathrm{t}}+\mathrm{b}\overline{\mathrm{b}}$ background, the analysis assigns leading order matrix element signal and background probability densities to each event. A likelihood-ratio statistic based on these probability densities is used to extract the signal. The results are provided in terms of an observed $\mathrm{t}\overline{\mathrm{t}}$H signal strength relative to the standard model production cross section $\mu=\sigma/\sigma_\mathrm{SM}$, assuming a Higgs boson mass of 125 GeV. The best fit value is $\hat{\mu} =$ 0.9 $\pm$ 0.7 (stat) $\pm$ 1.3 (syst) = 0.9 $\pm$ 1.5 (tot), and the observed and expected upper limits are, respectively, $\mu<$ 3.8 and $<$ 3.1 at 95% confidence levels.


Introduction
Since the observation of a new boson with a mass of approximately 125 GeV [1][2][3], a major focus of the program of the ATLAS and CMS experiments at the CERN LHC has been to measure the properties of this particle.To date, all measurements remain consistent with the standard model (SM) Higgs boson (H) hypothesis.The new particle has been observed in γγ, ZZ, WW, and ττ decays [4][5][6][7][8][9][10], and evidence has been reported in bb final states [11,12].The observed production and decay rates and spin-parity properties are compatible with SM expectations [13][14][15][16][17] and with the measured mass of m H = 125.26± 0.21 GeV [18].This paper describes a search for the SM Higgs boson produced in association with a pair of top quarks (ttH production) in all-jet events, where the Higgs boson decays exclusively to bb, and each top quark decays to a bottom quark and a W boson, which in turn decays to two quarks.The analysis uses √ s = 13 TeV proton-proton (pp) collision data collected by the CMS detector in 2016, corresponding to an integrated luminosity of 35.9 fb −1 .
In the context of probing the properties of the Higgs boson, the ttH process is interesting for several reasons.Because of the large top quark mass of m t = 172.44 ± 0.48 GeV [19], the top quark Yukawa coupling (y t ) is close to unity, and it therefore plays a major role in many beyond the SM theories.Although model-dependent indirect measurements of y t have been made in processes involving top quark loops, such as in Higgs boson production through gluongluon fusion and in Higgs boson decays to photons, these loop processes can be influenced by new contributions.A direct measurement of y t is therefore crucial to constrain many beyond the SM models.As the direct measurements of y t have not yet reached suitable precision in specific final-state searches nor in their combinations, the addition of an independent search will enhance the overall sensitivity of such measurements.
The ATLAS and CMS experiments have previously searched for ttH production in data corresponding to ≈5, 20, and 36 fb −1 at centre-of-mass energies of 7, 8, and 13 TeV, respectively [6,[20][21][22].The most recent searches in the H → bb decay mode have been performed in the leptonic channels by CMS [23] and ATLAS [24].In addition, ATLAS dedicated a search for ttH in the all-jet final state at √ s = 8 TeV, in which the observed and expected upper limits were set at multiples of 6.4 and 5.4 of the respective SM cross sections (σ SM ) at 95% confidence levels (CL), and a best fit value was obtained for the signal strength relative to the SM value of μ = σ/σ SM = 1.6 ± 2.6 [25].
The leading order (LO) Feynman diagram for ttH production in the all-jet final state is shown in Fig. 1.Ideally, signal events contain eight jets, of which four are tagged as originating from b quarks via the algorithm described in Section 4. To accommodate jets lost due to geometrical acceptance, merging of jets, and b tagging inefficiency, events with seven or more jets and three or more b-tagged jets are included in the defined signal region (SR).To account for extra jets from initial-or final-state radiation, up to nine jets are considered per event.
Despite that the signature involves a large number of high transverse momentum (p T ) finalstate jets, the analysis is constrained by background, mainly from jets produced through the strong interaction, referred to as quantum chromodynamic (QCD) multijet production.The next largest background contribution is from tt+jets, which includes tt+light-flavour jets, where one or more of the jets are incorrectly identified as b jets, tt+cc, and the topologically irreducible tt+bb background.Smaller background contributions arise from single top quark and vector boson production in association with jets.Although the all-jet final state has significantly more background than, for example the lepton+jets final state, it offers a fully reconstructed event, and a larger contribution from the signal.Because background rates are much larger than the ttH signal, it is important that all available information on the differences between the signal and backgrounds is incorporated into the analysis.A quark-gluon discriminant is used to differentiate between events containing jets originating from light-flavour quarks and events containing jets from gluons.Given the many combinations of jet-quark associations, it is difficult to resolve a clear, resonant Higgs boson mass peak.Nevertheless, there are underlying kinematic differences between the ttH signal and the multijet background and, to a lesser extent, the tt+jets background.These differences are exploited through the matrix element method (MEM) [23,[26][27][28] to distinguish the signal from background.Specifically, events are assigned a probability density according to how compatible they are with the ttH process.Although this probability density suffices to separate the signal from most backgrounds, a second probability density is assigned to each event according to its compatibility with the tt+bb process, which provides additional discrimination against the topologically irreducible tt+bb background.The two probability densities are combined in a likelihood ratio to form the final discriminant between signal and background.
The body of this paper is devoted to the details of the analysis.The CMS detector is described in Section 2, while the description of the Monte Carlo (MC) simulation follows in Section 3. The reconstruction of physical objects and event selections are discussed in Section 4. The QCD multijet background, estimated from control regions in data, is described in Section 5, followed by the description of the analysis strategy and the MEM in Section 6.The treatment of the systematic uncertainties is discussed in Section 7, the results are presented in Section 8, and a brief summary is provided in Section 9.

The CMS detector
The central feature of the CMS apparatus is a superconducting solenoid of 6 m internal diameter, providing a magnetic field of 3.8 T. 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 end sections, reside within the field volume.Forward calorimeters extend the pseudorapidity (η) coverage beyond these barrel and end sections.Muons are measured in gas-ionization detectors embedded in the steel flux-return yoke outside the solenoid.Events of interest are selected using a two-tiered trigger system [29].The first level, composed of specialized hardware processors, uses information from the calorime-ters and muon detectors to select events at a rate close to 100 kHz, ascertained within a time interval of less than 4 µs.The second level, known as the high-level trigger, consists of a farm of processors running a version of the full event reconstruction software optimized for fast processing, and reduces the event rate to ≈1 kHz before data storage.A more detailed description of the CMS detector, together with a definition of the coordinate system and the kinematic variables, can be found in Ref. [30].

Simulated samples
This analysis uses MC event generators and a detailed detector simulation based on GEANT4 (v9.4) [31] to model experimental effects, such as reconstruction, selection efficiencies, and resolutions in the CMS detector.
The ttH signal is modelled at next-to-LO (NLO) perturbative QCD via the POWHEG (v2) event generator [32][33][34][35].The value of the Higgs boson mass is set to m H = 125 GeV, the value of the top quark mass to m t = 172.5 GeV, and the parton distribution functions (PDFs) of the proton are described using NNPDF3.0[36], at the same level of accuracy in perturbative QCD as used to generate the event samples.These events are subsequently processed through PYTHIA (v8.200) [37] for parton showering and hadronization.
The POWHEG generator [38][39][40] is used to simulate tt and single top quark (single t) production at NLO.The associated production of tt with a vector boson (referred to as tt+W and tt+Z, or commonly as tt+V in the following) is simulated at NLO using MADGRAPH5 aMC@NLO (v2.2.2) [41] with FxFx jet matching and merging [42].The MADGRAPH5 aMC@NLO generator is also used at LO accuracy with the MLM matching scheme [43] to simulate the production of W and Z bosons with additional jets (V+jets) as well as QCD multijet events.The multijet background is estimated from control regions in data, as described in Section 5, while simulated events are used to ensure the self-consistency of the method.Diboson (WW, WZ, and ZZ) events are simulated at LO using PYTHIA.All background events also use PYTHIA for parton showering and hadronization.The underlying event in the simulation of all backgrounds, except for tt, is characterized through the PYTHIA CUETP8M1 event tune [44,45].For ttH and tt production, we use a custom tune CUETP8M2, developed by CMS with an updated strong coupling α s for initial-state radiation, to better model the jet multiplicity spectrum.
For comparing with measured distributions, the events in the simulated samples are normalized to the integrated luminosity of the data according to their predicted cross sections, taken from calculations at next-to-NLO (NNLO) for V+jets production, at approximate NNLO for single t production [46], and at NLO for diboson [47] and tt+V production [48].The ttH cross section and Higgs boson branching fractions are also generated at NLO accuracy [49].The simulated tt events are normalized to the full NNLO calculation, including soft-gluon resummation at next-to-next-to-leading logarithmic accuracy [50], assuming m t = 172.5 GeV.The dependence on choice of the PDFs was also studied and found to have negligible impact on all observed distributions; however, their impact on systematic uncertainties enters through their normalization.As the differential measurements of top quark pair production indicate that the measured p T spectrum of top quarks is softer than in simulation [51], the top quark p T distribution in POWHEG-generated tt events is reweighted to the data using a control region enriched in lepton+jets tt events.The events in the tt sample are also separated according to the flavour of additional jets that do not originate from top quark decays, namely: tt+bb, where two additional b jets are within the acceptance defined in Section 4 and originate from one or more b quarks; tt+b, where only one b jet is within the acceptance and originates from a single b quark; tt+2b, where a single b jet is within the acceptance and originates from two spatially close b quarks; tt+cc, where there is at least one c jet and no additional b jets within the acceptance; and tt+light-flavour partons (tt+lf).The separation is motivated by the fact that the subsamples originate from different processes and therefore have different associated systematic uncertainties.
Effects from additional pp interactions in the same or neighbouring bunch crossings (referred to as pileup) are modelled by adding simulated minimum-bias events generated in PYTHIA to all simulated processes.The pileup multiplicity distribution in simulation is reweighted to reflect the one observed in pp collision data, which is distributed between about 10 and 40 and peaks at ≈20.Additional rescaling factors are applied where necessary to improve the description of data in simulation, as discussed in the next section.

Event reconstruction and selection
The particle-flow event algorithm [52] is used to reconstruct and identify physical objects from an optimized combination of information from the various components of the CMS detector.The reconstructed vertex with the largest value of summed object p 2 T is taken to be the primary pp interaction vertex.These objects correspond to jets, clustered with the jet finding algorithm [53,54] using the tracks assigned to the vertex as inputs, and to the associated missing transverse momentum ( p miss T ), defined as the negative vector sum of the p T of those jets.The energy of charged hadrons is determined from a combination of their momenta measured in the tracker and the corresponding ECAL and HCAL energy depositions, corrected for the response function of the calorimeters to hadronic showers, while the energy of neutral hadrons is obtained from the corresponding corrected ECAL and HCAL energies.Charged hadrons from pileup events are omitted in the subsequent event reconstruction.The energy of electrons is determined from a combination of the electron momentum at the primary interaction vertex obtained from the tracker, the energy in the corresponding ECAL clusters, and the energy sum of all bremsstrahlung photons spatially compatible with the electron track.The p T of muons is based on the curvature of the corresponding track.
The selection criteria reflect the all-jet final state in ttH production where the Higgs boson decays into a bb pair.The analysis requires the reconstruction of all objects originating from a common vertex, i.e. jets from quark hadronization, including b-tagged jets, and electrons and muons to veto any events with leptons.The lepton veto ensures that leptonic final states in ttH production are not considered, as these are covered by separate searches at CMS. Jets are reconstructed from particle-flow candidates using the anti-k T clustering algorithm [53], implemented in the FASTJET [54] package, using a distance parameter of 0.4.In the offline reconstruction, each jet is required to have p T > 30 GeV and |η| < 2.4, and ECAL and HCAL energy fractions of at least 1%, which removes jets originating from instrumental effects.In addition, jet energy corrections are applied that depend on p T , η, pileup, and residual differences between data and MC simulations [55,56].
Jets originating from decays of b quarks are identified using the combined secondary vertex (CSVv2) algorithm [57,58], which uses information from the tracker.The algorithm combines information on the impact parameters at the primary vertex with features of secondary vertices within the jet, via a neural network discriminant.The least restrictive ("loose") and "medium" working points of the CSVv2 algorithm have respective efficiencies of ≈80 and 65% for tagging jets originating from b quarks, ≈35 and 10% for jets originating from c quarks, and probabilities of ≈10 and 1% for jets from light-flavour quarks or gluons to be misidentified as b jets.In the remainder of this paper, jets passing the loose and medium b tagging working points are referred to as CSVL and CSVM jets, respectively.The distribution in the output of the CSVv2 discriminant in simulation is corrected using scale factors to provide a better description of the observed features in data.This correction is obtained separately for light-flavour and b jets via a "tag-and-probe" approach [59] in dilepton control samples enriched in events with a Z boson and two jets, and in tt events with no additional jets, respectively.No correction is applied for c jets, although a larger systematic uncertainty is assigned to such jets, as described in Section 7.
Two dedicated triggers have been developed to efficiently select signal events, both requiring at least six jets with |η| < 2.6.The first trigger has kinematic requirements of jet p T > 40 GeV and H T > 450 GeV, and requires at least one of the jets to be b tagged; the second trigger is complementary in that it has less stringent kinematic conditions of jet p T > 30 GeV and H T > 400 GeV, but requires at least two b-tagged jets; the H T being defined by the scalar sum of jet p T of all jets in the event satisfying the given thresholds.The b jets are tagged online by the triggers at an efficiency of ≈70-80%, with a misidentification rate of ≈6% for light-flavour quark and gluon jets.Events in data and simulation are selected if they pass either of these two triggers.The efficiency in data and in simulation is measured in bins of the number of CSVM jets, the p T of the jet with the sixth-highest p T , and the H T in control samples collected using single-muon triggers.A bin-by-bin scale factor, based on the ratio of this efficiency in data to that in simulation, is applied to simulated events to correct for any remaining differences.The overall trigger efficiency for signal events that pass the offline event selection is 99%.
To ensure that the trigger selection is close to full efficiency relative to the offline selection, thereby reducing the uncertainty in any efficiency differences between data and simulation, the offline analysis selects events that contain at least six jets with p T > 40 GeV, requiring H T > 500 GeV and at least 2 b-tagged CSVM jets.Events with reconstructed electrons or muons in the final state are vetoed using loose lepton reconstruction and identification criteria, in particular, electrons and muons are required to have p T > 15 GeV, |η| < 2.4, and to be isolated from jet activity.If at least one reconstructed electron or muon is found in an event, the event is rejected.The distributions in H T and jet multiplicity after the preselection are shown in Fig. 2 for the data and simulation, with the signal contribution scaled to the total background yield for illustration.
To reject events that are unlikely to include a W boson from top quark decay, a cutoff is placed on the dijet invariant mass.Any jets that are not b tagged are considered in this calculation, and the pair with invariant mass m qq closest to the W boson mass is chosen as the W boson candidate in the event.Requirements of 60 < m qq < 100, and 70 < m qq < 92 GeV are applied in the 7-or 8-jet, and 9-or-more-jet categories, respectively.
Jets are also classified according to a quark-gluon likelihood (QGL) [60,61].The QGL discriminates light-flavour quark jets from gluon jets, exploiting the higher track multiplicity of gluon jets, the lower p T of their constituents, as well as their less collimated spatial profile.A jet is assigned a QGL discriminant value near 1 if it is quark-like, while gluon-like jets receive a value near 0. This is used to separate events containing light-flavour jets from qq decays of W bosons from events containing gluon jets produced in multijet interactions or via initial-state radiation.The QGL discriminant is used only for jets not identified as originating from b quarks, and is combined in a likelihood ratio (LR) defined as: where N 1 is the number of jets as well as the number of quarks in the all-quark hypothesis, and N 2 is the number of quarks in the quark+gluon hypothesis.The individual likelihoods are defined by: where ζ i is the QGL discriminant for the i th jet, and f q and f g are the probability density functions for ζ i when the i th jet originates from a quark or gluon, respectively.The former include u, d, s, and c quarks.The sum in Eq. ( 2) runs over all unique permutations in assigning N q jets to quarks and N g jets to gluons.We use q LR(N,0) to compare the likelihood of N reconstructed light-flavour jets that arise from N quarks, to the likelihood of N reconstructed light-flavour jets that arise from N gluons.After excluding either the first 3 or 4 CSVM jets in an event (based on the highest CSVv2 value), up to N light-flavour jets, with N = 3, 4, or 5 are used to calculate the LR, with a requirement of q LR(N,0) > 0.5 set for the SR.To correct the modelling of the QGL distribution in simulation, a reweighting based on a control sample of µ + µ − + jet and dijet events is applied to each event according to the type (quark or gluon) and QGL value of all jets in the event.The distribution in the quark gluon LR (QGLR), calculated excluding the first 3 b-tagged jets, is shown in Fig. 3 for data and simulation for events passing the preselection.Good agreement is observed between prediction and data.
Events / 100 GeV      The simulated backgrounds are first scaled to the integrated luminosity of the data, and then the simulated multijet background is rescaled to match the yield in data.The contribution from signal (blue line) is scaled to the total background yield (equivalent to the yield in data) to enhance readability.The hatched bands reflect the total statistical and systematic uncertainties in the background prediction, prior to the fit to data, which are dominated by the systematic uncertainties in the simulated multijet background.The ratio of data to background is given below the main panel.

Background estimation
The main backgrounds stem from multijet and tt production associated with additional gluons, light-flavour, charm, or bottom quarks (tt+jets).The background from tt+jets as well as other minor backgrounds (single t, V+jets, tt+V, and diboson events) are estimated from MC simulation, while the background from multijet events is obtained from control regions in data, as described below.The approach uses a control region (CR) with low b tag multiplicity to estimate the contribution from multijet events in the SR.The CR is enriched in multijet events, and the remaining contribution from other backgrounds (mainly tt+jets) is subtracted using simulation.
The CR is defined by events with two CSVM jets and at least three CSVL jets.We define a validation region (VR) using events with QGLR < 0.5.This definition provides four orthogonal regions, summarized in Table 1, from which we can obtain and check the multijet background estimate.The use of the VR relies on the fact that the QGLR and the number of additional CSVL jets are uncorrelated by construction, since only the non-CSVL jets are used in the calculation of the QGLR, except for the rare case of events with 5 or more CSVL jets.
Table 1: Definition and description of the four orthogonal regions in the analysis.The four orthogonal regions are used independently in each of the six analysis categories described in Section 4. For a given variable, the distribution of the multijet events in the SR of each category is estimated from the data in the CR, after subtracting tt+jets and other minor backgrounds.Since the kinematic properties of b-tagged and untagged jets differ in the CR and SR because of different heavy-flavour composition, corrections as a function of jet p T , η, and the minimum distance between the jet and the first two b-tagged jets (∆R min ) are applied to the CSVL jets in the CR.The correction function is obtained from jets in events passing the preselection, excluding the first two jets, ordered according to their CSVv2 discriminants.This is used to reweight the kinematic distributions of CSVL jets to match those of CSVM jets.The corrected multijet distribution in the CR is scaled to provide an estimate of the distribution in the SR.Specifically, the multijet yield in each category is left floating in the final fit, as discussed in Section 8.
A consistency check of the procedure used to estimate the multijet background is performed through simulation.Since the power of this test is limited by the statistical uncertainty in simulated QCD multijet events, the method is validated in data using events with QGLR < 0.5, by applying the same procedure used to estimate the multijet background in the SR.The distributions in the MEM discriminant, described in Section 6, as well as several kinematic variables in the VR in data, together with the multijet background estimate obtained from data, and other simulated contributions, are shown for the most sensitive event category, (≥9j, ≥4b), in Fig. 4. All discrepancies are taken into account via the systematic uncertainties in the multijet background, particularly in events with low H T in the (7j, 3b) category and in all ≥4b categories, in which two uncorrelated uncertainties are applied, as discussed in Section 7.
To verify that the good performance demonstrated in the gluon jet enriched VR (QGLR < 0.5) also holds in the quark jet enriched SR (QGLR > 0.5), we investigate another CR.Specifically, the b tagging criteria are changed by selecting jets that fulfil some intermediate but not the formal medium b tagging requirement, which corresponds to an efficiency for tagging b jets of ≈75% with a misidentification rate of ≈4%.The classified jets are then used to form jet and b jet multiplicity categories in analogy with the SR.Since these categories are orthogonal to the categories in the SR, they are used to verify that the background estimation is valid for QGLR > 0.5.In all validation regions, this background estimation reproduces the kinematic distributions measured in data.In the remainder of this paper, all multijet estimates are based on data, unless stated to the contrary.

Signal extraction
A likelihood technique based on the LO matrix elements for the ttH signal and the tt+bb background processes is used to extract the signal.This method utilizes the full kinematic properties of each event to provide a powerful MEM discriminant between the signal and background.Although the discriminant is constructed to discriminate against the tt+bb background, it performs well against tt+light-flavour jets and against multijet events.Consequently, it is used as the single discriminant against all background sources.The MEM algorithm is similar to the one documented in Ref. [23], but adapted to the all-jet final state.
The MEM probability density functions for signal and background are constructed at LO assuming that both reactions proceed via gluon-gluon fusion, which is valid, as this production mechanism represents most of the event yield.For example, at √ s = 14 TeV, the fraction of gluon-gluon initiated ttH subprocesses corresponds to ≈80% of the inclusive NLO cross section [62].Example Feynman diagrams for these processes are shown in Fig. 5.All inequivalent jet-quark associations in the reconstruction of the final state are considered in the analysis.The level of agreement between data and estimation is expressed in terms of a χ 2 divided by the number of degrees of freedom (dof), which are given along with their corresponding statistical p-values.
The differences between data and the total estimates divided by the quadratic sum of the statistical and systematic uncertainties in the data and in the estimates (pulls) are given below the main panels.The numbers in parenthesis in the legends represent the total yields for the corresponding entries.The last bin includes event overflows.
Hadronization and detector effects are taken into account via transfer functions obtained from simulation that map the measured jet four-momenta to the final-state particles.Each event contains three or four jets that most likely originate from b quarks (according to their CSVv2 discriminant values), and are considered as candidates for b quarks from the H → bb and t → Wb decays, whereas untagged jets are considered as candidates for the W → ud or cs quarks (and their charge conjugates).The accepted jet permutations must retain the assigned "b" status of their quarks.To account for the loss of jets because of limited detector acceptance, as well as the presence of additional jets from gluon radiation, the method integrates over certain final-state variables in many categories, specifically over quark directions.In events with eight jets and four b jets, one light-flavour jet is excluded from our MEM calculation in turn, and a sum is taken over the additional permutations.This approach has been checked and shown to provide improved performance relative to the fully reconstructed hypothesis in the (8j, ≥4b) category.Events with only three b jets are assumed to have lost a bottom quark from the decay of a top quark.Up to five untagged jets are considered as other than b quark candidates, while additional untagged jets are ignored.For five such quark candidates, one is excluded in turn and the number of permutations is increased by a factor of five.The final choice of hypothesis for each category, optimized according to discrimination power and computing performance, is given in Table 2.
where y represents the set of measured observables, i.e. jet momenta, and s is the square of the pp centre-of-mass energy.The sum is taken over the possible permutations of matching observed jets to final-state quarks.The integration is performed over the phase space volume, x, of the final-state quarks, with momenta p = (p 1 , . . ., p 8 ), and over the incident gluon momentum fractions, x a and x b , using the VEGAS [63] algorithm.The gluon PDF g is calculated in LHAPDF [64] at the factorization scale Q, which is the momentum transfer in the event, using the CTEQ6.6 set of PDFs [65].The four-dimensional delta function is reduced to two dimensions, which ensures conservation of energy and longitudinal momentum, while the total p T of the final-state quarks is weakly constrained via a resolution function, R( ρ T ), where ρ T is the measured transverse recoil, equal to the difference between p miss T and the vector p T sum of the specifically selected jets.The scattering amplitude for the given hypothesis, M H , is calculated at LO using OPENLOOPS (v1.3.1)[66], and all resonances are treated using the narrow-width approximation [67].The masses of the Higgs boson, the W boson, and the top quark are fixed, thereby constraining the energy of certain final-state quarks.The final integrations of the fully reconstructed ttH and tt+bb hypotheses are therefore performed over the incoming gluon momentum fractions and three or four respective final-state quark energies, with an integration range of a factor of three that of the associated jet energy resolution (JER).
The transfer function W(y|p) corresponds to the likelihood of measuring the set of variables y given the final-state quark momenta p.Given the excellent angular resolution of reconstructed jets, the direction of the quarks is assumed to be precisely measured, which reduces the total transfer function to a product of quark energy transfer functions.Each quark energy transfer function is modelled via a double Gaussian with different parameterizations for jets associated with b and other quarks, all obtained from MC simulation.Using Eq. ( 3), the MEM event discriminant is defined as: where κ S/B is a positive constant that adjusts the relative normalization between w(y|S) and w(y|B), given that the individual probability densities are not normalized to unity.The κ S/B constant is optimized category-by-category by minimizing the expected exclusion limit on the signal strength discussed in Section 8. Starting with a value that provides good visual discrimination between signal and background, the optimization is performed by testing different values of κ S/B above and below its original value until a minimum in the expected limit is found in each category as well as in combinations of categories.This optimization was proven to be equivalent to optimizing the significance of discovery.By construction, the MEM discriminant satisfies the condition 0 ≤ P S/B ≤ 1.
The event yields expected for the signal and the different background processes, and the yields observed in data, are listed in Table 3 for each analysis category.The proportion of expected ttH events increases with the number of jets and of b jets.To better represent the signal strength in each category, the expected yields for signal and background in the last three bins of the MEM discriminant are listed in Table 4.

Systematic uncertainties
Each source of systematic uncertainty is associated with a "nuisance" parameter that modifies the likelihood in the final fit discussed in Section 8, and can either affect the yield from a process (rate uncertainty), the distribution in the MEM discriminant, or both.In the latter case, the effects from the rate and distribution are treated simultaneously, and can be considered completely correlated.Each individual source of systematic uncertainty is independent of other Table 3: Expected numbers of ttH signal and background events, and the observed event yields for the six analysis categories, following the fit to data.The signal contributions are given at the best fit value.The quoted uncertainties contain all the contributions described in Section 7 added in quadrature, considering all correlations among the processes.Table 4: Expected numbers of ttH signal and background events in the last three bins of the MEM discriminant for the six analysis categories, following the fit to data.The signal contributions are given at the best fit value.The quoted uncertainties contain all the contributions described in Section 7 added in quadrature, considering all correlations among the processes.Also given are the signal (S) to total background (B) ratios for the SM ttH expectation (µ = 1).sources, and its effect on signal and background is completely correlated across the relevant processes.
The uncertainty in the integrated luminosity is estimated to be 2.5% [68].The uncertainties in trigger scale factors are determined from the bin-by-bin uncertainties in the ratio of efficiency in data relative to simulation, and are ≈1.0%on average, with some being as large as 15%.The uncertainty in the distribution in the number of pileup interactions is evaluated by changing the total inelastic cross section in the MC simulation by its uncertainty of 4.6%.The changes in the weight factors are propagated to the MEM discriminant distributions and treated as fully correlated among all processes.The impact of the uncertainty in the jet energy scale (JES) [55] is evaluated for each jet in the simulated events by changing the correction factors by their uncertainties, and propagating the effect to the MEM discriminant by recalculating all kinematic quantities.Various independent sources contribute to the overall JES uncertainty, and their impact is therefore evaluated separately and they are treated as uncorrelated in the final fit.A similar procedure is applied to account for the uncertainty related to the JER, which ranges between about 1 and 5% of the expected energy resolution, depending on the jet direction.Since the analysis categories are defined in terms of jet multiplicity and kinematics, a change in either JES or JER can induce a migration of events between categories, as well as in the SR.The fractional changes in event yields induced by one standard deviation shifts in JES or JER range, respectively, from 3-11% and 2-11%, depending on the process and the category.
The scale factors applied to correct the CSVv2 discriminant, described in Section 4, are affected by several components of systematic uncertainty, which are attributed to three main sources: JES, purity of heavy-or light-flavour jets in the control sample used to obtain the scale factors, and the statistical uncertainty in the event sample used in their extraction [58].A separate, large uncertainty is applied to charm jets, owing to the lack of a reliable calibration based on data.Each component of these systematic b tagging uncertainties is considered as a separate nuisance parameter in the final fit.
The systematic uncertainties from reweighting the simulated QGL and top quark p T distributions affect both the rate and the distribution in the MEM discriminant, and are taken into account by changing the functions used in the reweighting procedure by their uncertainties.Many uncertainties that would normally be related to the MC simulation of the multijet background are avoided by estimating this contribution from data.Nevertheless, a few small systematic uncertainties remain: an uncertainty in the correction applied to CSVL jets that is estimated via an alternative correction obtained by applying a second η correction; a 2.5 or 2.0% uncertainty in the first bin of the distribution in the MEM discriminant across ≥4b and 3b categories, applied to account for respective systematic over-or underestimation of this bin in the VR; a reweighting based on the H T distribution (considering the p T of just the first six jets) in the (7j, 3b) and all ≥4b categories of the VR, to account for mismodelling at low H T ; and the total normalization in each category that is left unconstrained in the final fit.The uncertainties in the multijet background normalization have the largest impact on the sensitivity of the analysis, and setting the normalization to a fixed value in each category improves the expected limit by 20 to 30%.
Theoretical uncertainties are treated as process-specific if they only impact a single simulated event sample, otherwise they are treated as correlated across all samples to which they apply.For example, the PDF uncertainty for quark-quark initiated processes, such as V+jets, tt+W, and diboson events, are treated as 100% correlated.The tt+heavy-flavour processes represent important sources of irreducible background, which have not yet been reliably measured, nor subjected to higher-order calculations that constrain these contributions.In fact, the most recent direct measurement of the tt+bb cross section in the dilepton final state has an accuracy of ≈35% [69].A 50% uncertainty in the production rate is therefore assigned separately to the tt+bb, tt+2b, tt+b, and tt+cc processes.These uncertainties have a significant impact on the expected limits, and ignoring these improves the overall limit by about 5%.
Uncertainties from factorization (µ F ) and renormalization (µ R ) scales in the inclusive cross sections range from 1 to 13%, depending on the process.The PDF uncertainties [70] range between 2 and 4%, and are treated as fully correlated for all processes that share the same dominant initial state (i.e.gg, gq, or qq), except for ttH, which is considered separately.Finally, the limited number of simulated background and signal events leads to statistical fluctuations in the nominal predictions.This is taken into account by assigning a nuisance parameter for each bin of each sample that can be changed by its uncertainty as described in Refs.[71,72].Table 5 summarizes the various sources of systematic uncertainty and their impact on yields.

Results
The MEM discriminant is used to extract the signal.The results are interpreted in terms of the signal strength modifier, defined as the ratio of the measured ttH production cross section to the SM expectation for m H = 125 GeV (µ = σ/σ SM ).The statistical method used to determine the signal strength is documented in Ref. [73].In this method, a likelihood function L(µ, θ) is constructed from the product of Poisson likelihoods for each bin of the MEM discriminant distribution for all six event categories, multiplied by a probability density function for the nuisance parameters θ, discussed in Section 7. The following test statistic, based on the LR, is used to compare competing hypotheses for the value of µ: where θ and θµ represent the best fit estimate of θ for a floating or fixed µ, respectively, following the procedure described in Section 3.2 of Ref. [14].
The final MEM discriminants in each category, following the combined fit to data, are displayed in Fig. 6.Each background contribution is initially normalized to an integrated luminosity of 35.9 fb −1 , while the multijet contribution is free to float in the fit.The fit to the nuisance parameters constrains many uncertainties including the multijet normalization.The total uncertainty contains the correlations among the fitted parameters that are sampled through the covariance matrix.
The observed distributions are compared to the signal+background expectation to extract a best fit value of μ = 0.9 ± 0.7 (stat) ± 1.3 (syst) = 0.9 ± 1.5 (tot), which is consistent with the SM expectation of µ = 1.The statistical component of the total uncertainty is estimated by including in the fit only a limited set of nuisance parameters that are statistical in nature, namely the per-category multijet normalizations and the multijet bin-by-bin uncertainties, which are dominated by the limited event count in the CR.The systematic component is then calculated as the difference in quadrature between the total uncertainty and the statistical component.
The best fit values of µ in each analysis category and in the combination are listed in Table 6, and displayed in Fig. 7 (left).Since the value is also compatible with the background-only hypothesis, an exclusion limit at the 95% CL can be set using the asymptotic approximation [74] of the modified frequentist CL s method [75,76].Combining all categories, we obtain the observed and expected upper limits of µ < 3.8 and < 3.1, respectively.These limits in each category, as well as the combined limit, are listed in Table 6 and displayed in Fig. 7 (right).
The signal is displayed in Fig. 8 in terms of the distribution in log 10 (S/B), where S/B is the ratio of the signal to background yields in each bin of the six MEM discriminant histograms, obtained from a combined fit constrained to the SM cross section of µ = 1.Good agreement is found between the data and the sum of the SM signal and background over the whole range of this variable.

CMS
Figure 7: Best fit values in the signal strength modifiers μ, and their 68% CL intervals as split into the statistical and systematic components (left), and median expected and observed 95% CL upper limits on µ (right).The expected limits are displayed with their 68% and 95% CL intervals, as well as with the expectation for an injected SM signal of µ = 1.

Summary
A search for the associated production of a Higgs boson with a top quark pair is performed using proton-proton collision data collected by the CMS experiment at a centre-of-mass energy of 13 TeV, corresponding to an integrated luminosity of 35.9 fb −1 .Events are selected to be compatible with the H → bb decay and the all-jet final state of the tt pair, and are divided into six categories according to their reconstructed jet and b jet multiplicities.
The result of the search is presented in terms of the signal strength modifier µ for ttH production, defined as the ratio of the measured ttH production cross section to the one expected for a standard model Higgs boson with a mass of 125 GeV.From a combined fit of signal and background templates to the data in all event categories, observed and expected upper limits of µ < 3.8 and < 3.1, respectively, are obtained at 95% confidence levels.These limits correspond to a best fit value of μ = 0.9 ± 0.7 (stat) ± 1.3 (syst) = 0.9 ± 1.5 (tot), which is compatible with the standard model expectation.

Figure 1 :
Figure 1: An example of an LO Feynman diagram for ttH production, including the subsequent decay of the top quark-antiquark pair, as well as that of the Higgs boson into a bottom quarkantiquark pair.

Figure 2 :
Figure2: Distribution in H T (left) and jet multiplicity (right) in data (black points) and in simulation (stacked histograms), after implementing the preselection.The simulated backgrounds are first scaled to the integrated luminosity of the data, and then the simulated QCD multijet background is rescaled to match the yield in data.The contribution from ttH signal (blue line) is scaled to the total background yield (equivalent to the yield in data) to enhance readability.The hatched bands reflect the total statistical and systematic uncertainties in the background prediction, prior to the fit to data, which are dominated by the systematic uncertainties in the simulated multijet background.The last bin includes event overflows.The ratios of data to background are given below the main panels.

Figure 3 :
Figure3: Comparison of the distributions in the quark-gluon likelihood ratio in data (black points) and in simulation (stacked histograms), after preselection, excluding the first 3 b-tagged jets.The simulated backgrounds are first scaled to the integrated luminosity of the data, and then the simulated multijet background is rescaled to match the yield in data.The contribution from signal (blue line) is scaled to the total background yield (equivalent to the yield in data) to enhance readability.The hatched bands reflect the total statistical and systematic uncertainties in the background prediction, prior to the fit to data, which are dominated by the systematic uncertainties in the simulated multijet background.The ratio of data to background is given below the main panel.

Figure 4 :
Figure 4: Distributions in data, in simulated backgrounds, and in the estimated multijet background for the (≥9j, ≥4b) VR category.The MEM discriminant (upper left), H T (upper right), p T of the leading jet (middle left), minimum ∆R between b jets (middle right), invariant mass of the closest b jet pair (lower left), and minimum mass of all jet pairs (lower right).The level of agreement between data and estimation is expressed in terms of a χ 2 divided by the number of degrees of freedom (dof), which are given along with their corresponding statistical p-values.The differences between data and the total estimates divided by the quadratic sum of the statistical and systematic uncertainties in the data and in the estimates (pulls) are given below the main panels.The numbers in parenthesis in the legends represent the total yields for the corresponding entries.The last bin includes event overflows.

Figure 5 :
Figure 5: Examples of LO Feynman diagrams for the partonic processes of gg → ttH and gg → tt+bb.

Figure 6 :
Figure 6: Distributions in the fitted MEM discriminant for each analysis category.The contributions expected from signal and background processes (filled histograms) are shown stacked.The tt+b and tt+2b background process are shown combined with tt+bb, the single t and tt+V processes are shown combined as "Other t", and the V+jets and diboson processes are shown combined as "Electroweak".The signal distributions (lines) for a Higgs boson mass of m H = 125 GeV are multiplied by a factor of 500 and superimposed on the data.The distributions in data are indicated by the black points.The ratios of data to background appear below the main panels.

Figure 8 :
Figure 8: Distribution in log 10 (S/B), where S and B indicate the respective bin-by-bin yields of the signal and background expected in the MEM discriminant distributions, obtained from a combined fit with the constraint in the cross section of µ = 1.

Table 2 :
Selected MEM hypotheses for each event topology.The 4W2H1T hypothesis assumes 1 b quark from a top quark is lost, 3W2H2T assumes that 1 quark from a W boson is lost, and 4W2H2T represents the fully reconstructed hypothesis requiring at least 8 jets.

Table 5 :
Summary of the systematic uncertainties affecting the signal and background expectations.The second column indicates the range in yield of the affected processes, caused by changing the nuisance parameters by their uncertainties.The third column indicates whether the uncertainties impact the distribution in the MEM discriminant.A check mark ( ) indicates that the uncertainty applies to the stated processes.An asterisk (*) indicates that the uncertainty affects the multijet distribution indirectly, because of the subtraction of directly affected backgrounds in the CR in data.

Table 6 :
Best fit value of the signal-strength modifier μ and the median expected and observed 95% CL upper limits (UL) on µ in each of the six analysis categories, as well as the combined results.The best fit values are shown with their total uncertainties and the breakdown into the statistical and systematic components.The expected limits are given together with their 68% CL intervals.