Measurement of H → μ+ μ− production in association with a Z boson at the CEPC

The Circular Electron-Positron Collider (CEPC) is a future Higgs factory proposed by the Chinese high energy physics community. It is planned to operate at a center-of-mass energy of 240–250 GeV and is expected to accumulate an integrated luminosity of 5 ab−1 over ten years of operation. At the CEPC, Higgs bosons will be dominantly produced from the ZH associated process. The vast number of Higgs events collected will enable precise studies of its properties, including Yukawa couplings to massive particles. With GEANT4-based simulation of detector effects, we study the feasibility of measuring the Higgs boson decaying into a pair of muons at the CEPC. The results with and without information from the Z boson decay products are provided, showing that a signal significance of over 10 standard deviations can be achieved and the H-μ-μ coupling can be measured within 10% accuracy.


Introduction
The discovery of the Higgs-like boson completes the particle table of the Standard Model (SM) of particle physics. Up-to-date LHC measurements all indicate that the Higgs boson is indeed highly SM-like [1][2][3][4][5][6]. In the SM, Higgs couplings to massive particles are proportional to their mass (squared). Hence, the event rate of Higgs couplings to the first and second generation of massive fermions can be very small, making them difficult to measure at the LHC. The Circular Electron-Positron Collider (CEPC) [7], however, is designed to run at around 240∼250 GeV with a combined instantaneous luminosity of 2 × 10 34 cm −2 s −1 for the two planned experiments, and will deliver 5 ab −1 of integrated luminosity over ten years of running. The huge amount of data will enable precise measurement of the branching ratios of the Higgs to light fermions and determine the associated Yukawa couplings, including H-µ-µ, which is crucial to validate the consistency of the SM Higgs mechanism, since any deviation would indicate the existence of new physics.
Searches for H → µ + µ − production have been performed at the ATLAS and CMS experiments with Run-I and Run-II data [8][9][10]. The most stringent observed (expected) upper limit on the cross-section times branching ratio is found to be 2.8 (2.9) times the SM prediction [10]. Projections have also been made for the High Luminosity-LHC assuming an integrated luminosity of 3000 fb −1 collected by the ATLAS or CMS detector, which can lead to a signal significance of about 7 σ [11] with an accuracy of around 20% [12]. Studies have also been performed for the International Linear Collider (ILC). Considering a center mass energy of 250 GeV and an integrated luminosity of 250 fb −1 , the signal is dominated by the Higgs-strahlung from a Z boson and the signal significances for the sub-processes with a Z boson decaying into νν and qq are found to be 1.8 and 1.1 σ, respectively [13]. Then, accumulating 2000 fb −1 at 250 GeV and 4000 fb −1 at 500 GeV with actual beam polarization, the precision of signal strength becomes 20.5% and 15.4% respectively [14]. At a center-of-mass energy of 1 TeV with an integrated luminosity of 500 fb −1 , the signal is dominated by the WW-fusion process and a sensitivity of 2.75 σ can be achieved [15].
At the CEPC, the signal H → µ + µ − production is dominated by the Higgs-strahlung from a Z boson. We perform a feasibility study based on events generated at leading order accuracy with initial state radiation (ISR), parton shower, hadronization and detector effects simulated.
Considering that 70% of the Z bosons decay hadronically and 20% decay invisibly, we focus on two scenarios, one for Z boson inclusive decay and the other for hadronic decay. The first case maximally exploits the statistics of the produced H→µ + µ − events and the second takes advantage of the major part of the decay kinematics. For both cases, we first perform a cut-based analysis and then improve the measurement using a Boosted Decision Tree (BDT) technique.
This paper is organized as follows. Section 2 describes event generation and simulation. Section 3 presents results for the inclusive measurement. Section 4 presents results for the Z→qq decay channel. Section 5 summarizes the paper.

Monte Carlo simulation
At a 250 GeV CEPC, Higgs bosons will mainly be produced through Higgs-strahlung, i.e. e + e − → ZH. With an integrated luminosity of 5000 fb −1 , about 230 of our signal events H → µ + µ − can be produced. The expected background to the signal production includes 2-fermion processes e + e − → ff, where f can be any SM fermion other than the top quark, and 4-fermion processes, which can be mediated through associated ZZ, WW, ZZ, WW production and a single Z boson production. All Monte Carlo (MC) events are generated with the Whizard V1.9.5 [16] event generator at parton level with ISR and interference effects included. The generated events are interfaced to Pythia 6 [17] for parton shower and hadronization simulation. Detector effects are simulated with the CEPC detector implemented with Mokka/GEANT4 [7,18,19]. The detector is assumed to have a similar structure to the International Large Detector (ILD) [20,21] at the ILC [22]. At the CEPC, the muon identification efficiency is expected to be over 99.5% for P T larger than 10 GeV, and with excellent P T resolution of σ 1/P T = 2×10 −5 ⊕1×10 −3 /(P T sinθ). The fully simulated events are reconstructed with a particleflow algorithm ArborPFA [23]. More details about the CEPC sample set can be found in Ref. [24].
The major SM backgrounds, including all the 2fermion processes (e + e − →ff, where ff refers to all lepton and quark pairs except tt) and 4-fermion processes (ZZ, WW, ZZ or WW, single Z). The initial state radiation (ISR) and all possible interference effects are automatically taken into account in the generation. The classification for four fermion production follows that of LEP [25], depending crucially on the final state. For example, if the final states consist of two mutually charge conjugated fermion pairs that could decay from both WW and ZZ intermediate states, such as e + e − ν eνe , this is classified as a "ZZ or WW" process. If there is e ± together with its partner neutrino and an on-shell W boson in the final state, this is called a "single W" process. Meanwhile, if there is an electron-positron pair and an on-shell Z boson in the final state, it is called a "single Z" process. Detailed information on the 2-fermion and 4-fermion samples used in our analyses are listed in Tables A1 and A2.

Inclusive analysis
A recoil mass method enables a measurement of the H →µ + µ − production without measuring the associated Z boson decay. We define the recoil mass as where √ s is the center-of-mass energy, and M H and E H correspond to the reconstructed mass and energy of the Higgs boson. The ZH (H→µ + µ − ) events form a peak in the M recoil distribution at the Z boson mass window. We select two muons with the largest transverse momenta and consider selections on the following kinematic variables: invariant mass of the di-muon system M µ + µ − , recoil mass of the di-muon system M µ + µ − recoil , transverse momentum of the di-muon system P T µ + µ − , third component of the di-muon momentum P Z µ + µ − , energy of dimuon system E µ + µ − , and angular variables in the laboratory system frame cosθ µ − , cosθ µ + , cosθ µ + µ − , cosθ Zµ − , and cosθ Zµ + , where θ µ − ,θ µ + means the polar angle of µ − , µ + ; θ µ + µ − means the angle between µ − and µ + ; and θ Zµ ± represents the angle between the Z boson and muon leptons.

Cut-count analysis
The event numbers under selection flow, which are determined by maximizing s/ √ s+b, with s and b represent signal and background yields, are summarized in Table 1. The two mass windows M µ + µ − and M µ + µ − recoil are set in accordance with the signal signature. P T µ + µ − and P Z µ + µ − are set to reduce the ZZ events, where one of the Z bosons decays to µ + µ − , and Drell-Yan Z→ µ + µ − background. The Higgs and Z boson decays can lead to different cosθ µ + and cosθ µ − distributions due to the spin-dependence of the couplings and the parity violation of the weak interaction. cosθ µ + µ − selection is chosen to supress the 2f background.
cosθ µ − , and cosθ µ + µ − in the inclusive analysis, after the preselection (2 well identified muons) and 120<M µ + µ − <130 GeV requirements. All the distributions are normalized to 10. An unbinned maximum likelihood fit is performed on the M µ + µ − distribution. The signal is parameterized by a crystal ball function, with parameters fixed by simulated events. The background is parametrized by a second or-der Chebychev function chosen by F-test [28].

BDT optimization
We have also exploited the Toolkit for Multivariate Analysis (TMVA) [26] for further background rejection, using the method of Gradient Boosted Decision Trees (BDTG). After fixing the range of the invariant mass and the recoil mass as mentioned above, 5 variables are taken as inputs to TMVA, including cosθ µ ± Z , cosθ µ ± and P Z µ + µ − . The choice of these variables is based on many tests and importance ranking. The resulted BDT response distribution can be seen in Figure 3, where the agreement between training and testing samples shows no obvious overtraining. We then take the final event selections as: BDTG response > 0.369, 20<P T µ + µ − <64 GeV and cosθ µ + µ − > −0.996. A maximum likelihood fit is performed on the resulting invariant mass of the di-muon system. The signal and background probability functions are parametrized in the same form as in the previous cut-count study. Figure 3 shows the BDT response distribution and the post-fit result of M µ + µ − . The fitted number of signal events is 62.3±10.9. At 68% confidence level, an accuracy from -16% to 17% can be achieved for the signal strength based on a likelihood scan. The signal under the peak at 124.4-125.2 GeV leads to a significance of 10.9 σ.

Z(qq)H(µµ) analysis
Of all the Z boson decay modes, the hadronic channel is most promising, due to its large branching fraction  (∼70%). The exclusive method of the k t algorithm for e + e − collisions in Fastjet [27] is used to reconstruct two jets with the particles expect the chosen µ − and µ + , and the jets are sorted by energy. We perform an analysis on the Z(qq)H(µµ) production. Apart from the previously mentioned variables related to the H(µµ) system, we further exploit the following selections on jets: the third component of di-jet system momentum P Z jj , the recoil mass of the di-jet system M jj recoil mass of jets M j1,2 , and the invariant mass of the di-jet system M jj .

Cut-count analysis
A cut-count analysis is performed for the exclusive analysis. The event flow under selections are summarized in Table 2. Selections on single and di-jet masses eliminate most background without hard jets. The recoil mass cut further reduces the Z(ll)Z(qq) background.
As in the inclusive channel, we perform a likelihood fit to extract the signal yield and strength parameter. The quality of the fit is demonstrated in Fig. 5. The signal yield from the fit is 75.5±12.5. The signal strength can be determined with an uncertainty from -16% to 17%, at 68% confidence level. The signal significance under the peak at 124.3-125.2 GeV is found to be 10.8σ.

BDT improvement
In order to achieve the highest significance, we perform a two-step multivariate analysis. The first step exploits a MLP (multilayer perceptron) [26] method to suppress the fully leptonic WW and ZZ backgrounds. After applying M µ + µ − recoil > 90 GeV, 4 variables including M j1,2 , M jj and M jj recoil are considered as inputs for the MLP. The effectiveness of this MLP is shown in Fig. 6. After requiring the MLP response to be greater than 0.71, we exploit a BDTG to further reduce the backgrounds from semileptonic ZZ and WW. In this second step, the variables cosθ µ ± , cosθ µ ± Z , P Z µ + µ − , P Z jet12 , cosθ j1/j2,H , cosθ j1,2 , and M jj are taken as inputs. After the two-step multivariate analysis, we require a BDTG response>-0.13, 90.4 < M µ + µ − recoil < 93 GeV, and 28<P T µ + µ − <64 GeV. Finally, we perform a likelihood fit to extract the signal yield and strength parameter, as  shown in Fig. 7. The signal yield from the fit is 73.4±12.4. Based on a likelihood scan, the signal strength can be determined with an uncertainty from -16% to 17%, at  68% confidence level. The significance of the signal in the peak region 124.3-125.2 GeV is found to be 10.8σ, which means the best boundary to distinguish signal and background is nearly an N-dimensional rectangle in the parameter phase space.

Summary
The feasibility of measuring H→µ + µ − at the CEPC has been studied considering center-of-mass energy 250 GeV collisions and 5000 fb −1 integrated luminosity. The measurement was performed in two complementary channels: ZH production without measuring the Z boson decay, and ZH production with the Z boson decaying hadronically. For each decay channel, a cut-count analysis was tested and followed with an improvement using multivariate techniques. Similar results are obtained from the two channels.
Finally, we want to mention in brief possible systematic factors affecting this future analysis. By referring mostly to Ref. [19], which shares most of the common factors, the systematics uncertainty should be under control, while the statistical uncertainty will dominate for this analysis at the CEPC.
Over 10 σ significance can be reached for the signal H → µ + µ − process. The accuracy of the signal strength can be measured with ±17% uncertainty and the associated H-µ-µ coupling can be restricted to 10% level. The results are comparable to the High-Luminosity LHC.