Search for direct pair production of supersymmetric partners to the τ lepton in proton-proton collisions at √ s = 13 TeV

A search is presented for τ slepton pairs produced in proton-proton collisions at a center-of-mass energy of 13 TeV. The search is carried out in events containing two τ leptons in the final state, on the assumption that each τ slepton decays primarily to a τ lepton and a neutralino. Events are considered in which each τ lepton decays to one or more hadrons and a neutrino, or in which one of the τ leptons decays instead to an electron or a muon and two neutrinos. The data, collected with the CMS detector in 2016 and 2017, correspond to an integrated luminosity of 77.2 fb−1. The observed data are consistent with the standard model background expectation. The results are used to set 95% confidence level upper limits on the cross section for τ slepton pair production in various models for τ slepton masses between 90 and 200 GeV and neutralino masses of 1, 10, and 20 GeV. In the case of purely left-handed τ slepton production and decay to a τ lepton and a neutralino with a mass of 1 GeV, the strongest limit is obtained for a τ slepton mass of 125 GeV at a factor of 1.14 larger than the theoretical cross section. Submitted to the European Physical Journal C c © 2019 CERN for the benefit of the CMS Collaboration. CC-BY-4.0 license ∗See Appendix A for the list of collaboration members ar X iv :1 90 7. 13 17 9v 1 [ he pex ] 3 0 Ju l 2 01 9


Introduction
Supersymmetry (SUSY) [1][2][3][4][5][6][7][8] is a possible extension of the standard model (SM) of particle physics, characterized by the presence of superpartners for SM particles. The superpartners have the same quantum numbers as their SM counterparts, except for the spin, which differs by half a unit. One appealing feature of SUSY is that the cancellation of quadratic divergences in quantum corrections to the Higgs boson mass from SM particles and their superpartners could resolve the fine tuning problem [9][10][11][12][13][14]. Another feature is that the lightest supersymmetric particle (LSP) is stable [15,16] in SUSY models with R-parity conservation [17], and could be a dark matter (DM) candidate [18].
The hypothetical superpartner of the τ lepton, the τ slepton ( τ ), is the focus of the search reported in this paper. Supersymmetric models where a light τ is the next-to-lightest supersymmetric particle are well motivated in early universe τ -neutralino coannihilation models that can accommodate the observed DM relic density [19][20][21][22][23][24]. The existence of a light τ would enhance the rate of production of final states with τ leptons in collider experiments [25,26].
In this analysis, we study the simplified model [27][28][29] of direct τ pair production shown in Fig. 1. We assume that the τ decays to a τ lepton and χ 0 1 , the lightest neutralino, which is the LSP in this model. The search is challenging because of the extremely small production cross section expected for this signal, as well as the large backgrounds. The most sensitive previous searches for direct τ pair production were performed at the CERN LEP collider [30][31][32][33], excluding τ masses at 95% confidence level (CL) up to ≈90 GeV for neutralino masses up to 80 GeV in some models. At the LHC, the ATLAS [34, 35] and CMS [36] Collaborations have also performed searches for direct τ pair production using 8 TeV data, and the CMS Collaboration has reported a search for direct τ pair production in an initial sample of 35.9 fb −1 at 13 TeV collected in 2016 [37]. This paper presents a significant improvement in search sensitivity, which was limited by the small signal production rates, through the incorporation of improved analysis techniques and the inclusion of the data collected in 2017. The data used correspond to a total integrated luminosity of 77.2 fb −1 .
Events with two τ leptons are used. We consider both hadronic and leptonic decay modes of the τ lepton, in which it decays to one or more hadrons and a neutrino, or to an electron or muon and two neutrinos, respectively. Independent analyses are carried out in the final states with two hadronically decaying τ leptons (τ h τ h ) and with one τ h and an electron or a muon ( τ h , where = e or µ). The presence of missing transverse momentum, which can originate from stable neutralinos as well as neutrinos from τ lepton decays, provides an important source of discriminating power between signal and background.
We have introduced several improvements with respect to the analysis presented in Ref. [37] that are applied to both 2016 and 2017 data. We make use of dedicated machine learning techniques to enhance the search sensitivity. These include the incorporation of an improved τ h selection method that makes use of a deep neural network (DNN) for the τ h τ h analysis, and of a boosted decision tree (BDT) for event selection in the τ h analyses. Improvements have also been made to the background-estimation techniques and to the search region (SR) definitions. The incorporation of these enhancements is expected to improve the search sensitivity by up to 50%, where the figure of merit considered is the 95% CL upper limit on the cross section for τ pair production obtained with the data collected in 2016. The improvement is less significant than expected, since it is found that the estimated signal acceptance is reduced when the fast detector simulation that was previously used to model signal events is replaced in this search with the more realistic, full GEANT4-based detector simulation [38]. Differences in the signal acceptance for the fast and more accurate full detector simulations are mainly caused by differ-ences in the reconstructed τ h visible transverse momentum (p T ), which is found to have larger values in the case of the fast simulation.
We consider the superpartners of both left-and right-handed τ leptons, τ L and τ R . The cross section for τ L pair production is expected to be about a factor of three larger than for τ R pairs [39]. The experimental acceptance is also expected to be different for left-and righthanded assignments because of the differences in the polarization of the τ leptons produced in τ L and τ R decays. The decay products of hadronically and leptonically decaying τ leptons originating from τ R decays are predicted to have larger and smaller p T , respectively, than those originating from τ L decays. Two simplified models are studied for direct τ pair production. One model involves production of only τ L pairs and the other is for the degenerate case in which both τ L and τ R pairs are produced. No mixing is introduced between left-and right-handed states. We study models with τ masses ranging from 90 to 200 GeV. The LEP limits [30][31][32][33] place strong constraints on the allowed values of the τ mass below this range, while the search sensitivity for τ masses above this range is low as a result of the decrease in production cross section with increased mass. We also consider different assumptions for the χ 0 1 mass, namely 1, 10, and 20 GeV. The search sensitivity decreases when the mass difference between the τ and χ 0 1 becomes small, since the visible decay products in such cases have lower momentum, resulting in a loss of experimental acceptance for such signals.  Figure 1: Diagram for direct τ pair production, followed by decay of each τ to a τ lepton and a χ 0 1 .

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, each composed of a barrel and two endcap sections, reside within the solenoid volume. Forward calorimeters extend the pseudorapidity (η) coverage provided by the barrel and endcap detectors. Muons are detected in gas-ionization chambers embedded in the steel flux-return yoke outside the solenoid. Events of interest are selected using a two-tiered trigger system [40]. The first level, composed of custom hardware processors, uses information from the calorimeters and muon detectors to select events at a rate of around 100 kHz 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, which reduces the event rate to about 1 kHz before data storage. A more detailed description of the CMS detector, together with definitions of the coordinate system and kinematic variables, can be found in Ref. [41].

Event reconstruction and simulation
The event reconstruction uses a particle-flow (PF) algorithm [42] that combines information from the tracker, calorimeter, and muon systems to identify charged and neutral hadrons, photons, electrons, and muons in an event. The missing transverse momentum vector, p miss T , is computed as the negative of the vector sum of the p T of all PF candidates reconstructed in an event, and its magnitude p miss T is used in the search as a discriminator between signal and SM background. Events selected for the search are required to pass filters [43] designed to remove detector-and beam-related backgrounds, and must have at least one reconstructed vertex. Usually, more than one such vertex is reconstructed because of pileup, i.e., multiple proton-proton (pp) collisions within the same or neighboring bunch crossings. The mean number of interactions per bunch crossing was 27 in 2016, and increased to 37 in 2017, assuming a total inelastic pp cross section of 80 mb. The reconstructed vertex with the largest value in summed object p 2 T is selected to be the primary pp interaction vertex (PV). These objects are defined by tracks associated with a given vertex that are clustered using a jet finding algorithm [44,45], and a more restricted form of the vector missing transverse momentum that is calculated from these track-based jets.
Charged particles that originate from the PV, photons, and neutral hadrons are clustered into jets using the anti-k T algorithm [44] with a distance parameter of 0.4, as implemented in the FASTJET package [45]. The jet energies are corrected to account for the contribution from pileup interactions and to compensate for variations in the detector response [45,46]. To mitigate issues related to noise in the ECAL endcaps that led to significantly worse modeling of the p miss T distribution, particularly for events with large values of p miss T in 2017 data, PF candidates that are clustered in jets in 2.65 < |η| < 3.14 with uncorrected p T < 50 GeV are not used in the calculation of p miss T in 2017 data and simulation. This improves the agreement in p miss T between simulated events and data.
Jets in the search are required to have their axes within the tracker volume of |η| < 2.4. For the τ h τ h analysis, we use jets with p T > 30 GeV, while for the τ h analyses, we veto events containing jets with p T > 20 GeV to provide efficient background rejection. Jets are required to be separated in η and azimuthal angle (φ) by ∆R ≡ √ (∆η) 2 + (∆φ) 2 > 0.4 from electron, muon, or τ h candidates in order to minimize double counting of objects. Jets originating from the hadronization of b quarks are "tagged" in the τ h τ h analysis through the DNN-based combined secondary vertex algorithm (DeepCSV) [47] to reject events with b quark jets that are likely to originate from backgrounds with top quarks. The efficiency for tagging b quarks originating from top quark decays is about 84%, while the misidentification rates for jets from charm quarks, and from light quarks or gluons, are about 41 and 11%, respectively. In the τ h analyses, the CSVv2 tagger [47] is used to identify b quark jets for the selection of backgroundenriched control regions (CRs). The working point that is used corresponds to an efficiency of 63% and misidentification rates of 12 and 0.9% for jets from charm quarks and light quarks or gluons, respectively. Electron candidates are reconstructed by first matching reconstructed tracks to clusters of energy deposited in the ECAL. Selections based on the spatial distribution of the shower, trackcluster matching criteria, and consistency between the cluster energy and the track momentum are then used in the identification of electron candidates [48]. Muon candidates are reconstructed by requiring reconstructed tracks in the muon detector to be matched to the tracks found in the inner tracker [49]. We require the origin of electron and muon candidates to be consistent with the PV. Restrictions are imposed on the magnitude of the impact parameters of their tracks relative to the PV in the transverse plane (d xy ), and on the longitudinal dis-placement (d z ) of the point of closest approach. To ensure that electron or muon candidates are isolated from jet activity, we define a relative isolation quantity (I rel ) as the ratio of the scalar p T sum of hadron and photon PF candidates, in an η-φ cone of radius 0.3 or 0.4 around the candidate electron or muon, to the candidate p T , requiring it to be below an upper bound appropriate for the selection. The quantity I rel is adjusted to account for the contributions of particles originating from pileup interactions. The electron and muon selection criteria applied in the analysis are the same as those described in Ref. [37].
The τ h candidates are reconstructed using the CMS hadrons-plus-strips algorithm [50]. The constituents of the reconstructed jets are used to identify individual τ lepton decay modes with one charged hadron and up to two neutral pions, or three charged hadrons. The τ h candidate momentum is determined from the reconstructed visible τ lepton decay products. The presence of extra particles within the jet that are incompatible with the reconstructed decay mode is used as a criterion to discriminate jets from τ h decays. A multivariate-analysis (MVA) based discriminant [50], which contains isolation as well as lifetime information, is used to suppress the rate for quark and gluon jets to be misidentified as τ h candidates. We employ a relaxed ("very loose") working point of this discriminant as a preselection requirement for the τ h candidates selected in the τ h τ h analysis, as well as in the extrapolation used to estimate the contributions of events to the background in which quark or gluon jets are misidentified as τ h candidates. This working point corresponds to an efficiency of ≈70% for a genuine τ h , and a misidentification rate of ≈1% for quark or gluon jets. A DNN is used to improve the discrimination of signal τ h candidates from background, as discussed in more detail below. Two working points are used in the τ h analysis: a "very tight" working point for selecting signal τ h candidates that provides stringent background rejection, and a "loose" working point for the extrapolation procedure to estimate the misidentified τ h background that provides higher efficiency and less background rejection. These working points, respectively, typically have efficiencies close to 45 and 67% for a genuine τ h , with misidentification rates of ≈0.2 and 1% for quark or gluon jets. Electrons and muons misidentified as a τ h are suppressed via criteria specifically developed for this purpose that are based on the consistency of information from the tracker, calorimeters, and muon detectors [50].
The dominant background in the τ h τ h final state originates from misidentification of jets as τ h candidates, mainly in SM events exclusively comprising jets produced through the strong interaction of quantum chromodynamics (QCD). These are referred to as QCD multijet events in what follows. To further improve the suppression of this background while retaining high signal efficiency, we have pursued a new approach for τ h isolation in the τ h τ h analysis that is based upon the application of a DNN that is fed information about the properties of PF candidates within an isolation cone with ∆R < 0.5 around the τ h candidate. We refer to this as "Deep Particle Flow" (DeepPF) isolation. Charged PF candidates consistent with having originated from the PV, photon candidates, and neutral hadron candidates with p T > 0.5, 1, and 1.25 GeV, respectively, provide the inputs to the DeepPF algorithm. The list of observables incorporated for each PF candidate includes its p T relative to the τ h jet, ∆R between the candidate and τ h , particle type, track quality information, and d xy , d z and their uncertainties, σ(d xy ) and σ(d z ). A convolutional DNN [51] is trained with simulated signal and background events. Signal τ h candidates are those that are matched to generator-level τ leptons from a mixture of processes that give rise to genuine τ leptons. Background candidates that fail the matching are taken from simulated W+jets and QCD multijet events. The DeepPF discriminator value is obtained by averaging the DNN output with the nominal MVA-based discriminant described above. The working point for DeepPF isolation is chosen to maintain a constant efficiency of ≈50%, 56%, and 56% as a function of p T for the three respective τ h decay modes: one charged hadron, one charged hadron with neutral pions, and three charged hadrons. Since the τ h candidate p T distribution in signal events depends on the τ and χ 0 1 masses, this choice of discriminator and working points allows us to maintain high efficiency for τ pair production signals under a large range of mass hypotheses. The overall misidentification rate for jets not originating from τ leptons ranges from 0.15% to 0.4% depending on p T and decay mode.
Significant contributions to the SM background originate from Drell-Yan+jets (DY+jets), W+jets, tt, and diboson processes, as well as from QCD multijet events, where DY corresponds to processes such as qq → + − . Smaller contributions arise from single top quark production and rare SM processes, such as triboson and Higgs boson production, and top quark pair production in association with vector bosons. We rely on a combination of measurements in data CRs and Monte Carlo (MC) simulation to estimate contributions of each source of background. The MC simulation is also used to model the signal. package. Finally, uncertainties in renormalization and factorization scale, and parton distribution functions (PDFs) have been obtained using the SYSCALC package [58]. Models of direct τ pair production are generated with MADGRAPH5 aMC@NLO at LO precision up to the production of τ leptons, with their decay modeled by PYTHIA 8.212 and 8.230 for the analysis of 2016 and 2017 data, respectively. The CUETP8M1 [59] (CUETP8M2T4 [60] for tt) and CP5 [61] underlying-event tunes are used with PYTHIA for the 2016 and 2017 analyses, respectively. The 2016 analysis uses the NNPDF3.0LO [62] set of PDFs in generating W+jets, DY+jets, and signal events, while the NNPDF3.0NLO PDFs are used for other processes. The NNPDF3.1NLO PDFs are used for all simulated events in the 2017 analysis.
Simulated events are reweighted to match the pileup profile observed in data. Differences between data and simulation in electron, muon, and τ h identification and isolation efficiencies, jet, electron, muon, and τ h energy scales, and b tagging efficiency are taken into account by applying scale factors to the simulation. We improve the modeling of initial-state radiation (ISR) in simulated signal events by reweighting the p ISR T distribution, where p ISR T corresponds to the total transverse momentum of the system of SUSY particles. This reweighting procedure is based on studies of the p T of Z bosons [63]. The signal production cross sections are calculated at NLO using next-to-leading logarithmic (NLL) soft-gluon resummations [39]. The most precise calculated cross sections available are used to normalize the simulated SM background samples, often corresponding to next-to-next-to-leading order accuracy.

Event selection
The data used in this search are selected through triggers that require the presence of isolated electrons, muons, τ h candidates, or p miss T . The data used for the τ h τ h analysis are collected with two sets of triggers. Events with p miss T < 200 GeV are selected using a trigger that requires the presence of two τ h candidates, each with p T > 35 and >40 GeV in 2016 and 2017 data, re-spectively. We gain up to 7% additional signal efficiency for events with p miss T > 200 GeV with the help of a trigger that requires the presence of substantial p miss T , with a threshold varying between 100 and 140 GeV during the 2016 and 2017 data-taking periods. For the eτ h final state, the trigger relies on the presence of an isolated electron satisfying stringent identification criteria and passing p T > 25 or >35 GeV in 2016 and 2017 data, respectively. For the µτ h final state, the trigger is based on the presence of an isolated muon with p T > 24 and >27 GeV in 2016 and 2017 data, respectively. Trigger efficiencies are measured in data and simulation. In addition to corrections mentioned in Section 3, we apply scale factors to the simulation to account for any discrepancies in trigger efficiency with data. These scale factors are parameterized in the p T and η of the reconstructed electron, muon, or τ h candidates, or the reconstructed p miss T for events selected using p miss T triggers.

Event selection and search regions in the τ h τ h final state
Beyond the trigger selection, the baseline event selection for the τ h τ h analysis requires the presence of exactly two isolated τ h candidates of opposite charge, satisfying the DeepPF selection described in Section 3, with |η| < 2.3 and p T > 40 and >45 GeV in the 2016 and 2017 analysis, respectively, as well as no additional τ h candidates with p T > 30 GeV satisfying the very loose working point of the MVA-based discriminant. We veto events with additional electrons or muons with p T > 20 GeV and |η| < 2.5 or <2.4 for electrons and muons, respectively, and reject any events with a b-tagged jet to suppress top quark backgrounds. A requirement of |∆φ(τ (1) h , τ (2) h )| > 1.5 helps to suppress the DY+jets background, while retaining high signal efficiency. Finally, we require p miss T > 50 GeV to suppress the QCD multijet background.
The removal of low-p T jets in the forward ECAL region from the p miss T calculation in 2017 (see Section 3) causes the background originating from DY+jets and other sources to increase in the SRs, since events with low-p T jet activity in that region are assigned larger values of reconstructed p miss T . We recover some of the corresponding loss in sensitivity in the 2017 analysis by placing an upper bound of 50 GeV on the scalar p T sum of low-p T jets excluded from the p miss T calculation (H low T ). This restriction reduces the impact of background events with significant low-p T jet activity in the forward region, for which the p miss T would be overestimated. To ensure that the efficiency of this requirement is correctly estimated in simulation, a Z → µ + µ − CR is used to extract correction factors for the H low T distribution in simulation that account for discrepancies with the distribution observed in data. In addition, to avoid effects related to jet mismeasurement that can contribute to spurious p miss T , we require the p miss T to have a minimum separation of 0.25 in |∆φ| from jets with p T > 30 GeV and |η| < 2.4, as well as from those with uncorrected p T > 50 GeV in the region 2.4 < |η| < 3.14.
Events satisfying the baseline selection criteria are subdivided into exclusive SRs using several discriminants. To improve the discrimination of signal from SM background, we take advantage of the expected presence of two χ 0 1 in the final state of signal events and their contribution to p miss T . Their presence skews the correlations between p miss T and the reconstructed leptons to be different from background processes, even for those backgrounds with genuine p miss T . These differences can be exploited by mass observables calculated from the reconstructed lepton transverse momenta and p miss T to provide discrimination of signal from background. For a particle decaying to a visible and an invisible particle, the transverse mass (m T ) calculated from the p T of the visible decay products should have a kinematic endpoint at the mass of the parent particle. Assuming that the p miss T corresponds to the p T of the invisible particle, we calculate the m T observable for the visible particle q and the invisible particle as follows: We use as a discriminant the sum of the transverse masses calculated for each τ h with p miss T , Σm T , given by Another variable found to be useful in the discrimination of signal from background is the "stransverse mass" m T2 [64][65][66]. This mass variable is a generalization of m T in the case of multiple invisible particles. It serves as an estimator of the mass of pair-produced particles when both particles decay to a final state containing the same invisible particle. It is given by: where p are the unknown transverse momenta of the two undetected particles, X(1) and X(2), corresponding to the neutralinos in our signal models, and m (i) T are the transverse masses obtained by pairing either of the two invisible particles with one of the two leptons. The minimization (min) is over the possible momenta of the invisible particles, taken to be massless, which are constrained to add up to the p miss T in the event. For direct τ pair production, with each τ decaying to a τ lepton and a χ 0 1 , m T2 should be correlated with the mass difference between the τ and χ 0 1 . A large value of m T2 is thus common in signal events for models with larger τ masses and relatively rare in SM background events.
The SR definitions for the τ h τ h analysis, shown in Table 1, are based on a cut-and-count analysis of the sample satisfying the baseline selections. The regions are defined through criteria imposed on m T2 , Σm T , and the number of reconstructed jets in an event, N j . The Σm T and m T2 distributions of events in the τ h τ h final state surviving the baseline selections are shown in Fig. 2. The distributions obtained for 2016 and 2017 data are combined. Separate sets of simulated events are used to model signal and background events in 2016 and 2017 data using the methods described in Section 3. In all distributions, the last bin includes overflow events. After applying a minimum requirement of m T2 > 25 GeV in all SRs, we subdivide events into low (25-50 GeV) and high (>50 GeV) m T2 regions, to improve the sensitivity to lower and higher τ mass signals, respectively. For each m T2 region, the Σm T distribution is exploited to provide sensitivity for a large range of τ mass signals. We define three bins in Σm T : 200-250, 250-300, and >300 GeV. Finally, we subdivide events in each m T2 and Σm T region into the categories N j = 0 and N j ≥ 1. This binning is beneficial as background events passing the SR kinematic selections are largely characterized by additional jet activity, while signal contains very few additional jets. The 0-jet category therefore provides nearly background-free SRs. However, we retain the SRs with N j ≥ 1 that are also expected to contain signal events with ISR or pileup jets.

Event selection in the τ h final states
The baseline event selections for the τ h analyses require either an electron with p T > 26 (35) GeV and |η| < 2.1 or a muon with p T > 25 (28)   Other SM Bkg. uncertainty Obs. / Pred. have |d z | < 0.2 cm, and electrons and muons are also required to have |d xy | < 0.045 cm. Electrons and muons have to satisfy I rel < 0.15 and <0.1, respectively. Backgrounds from tt and W+jets are greatly reduced by vetoing events that contain jets with p T > 20 GeV. Events from the W+jets background are further reduced by requiring the transverse mass m T ( , p miss T ), calculated using the electron or muon momentum vector and p miss T , to be between 20 and 60 GeV or above 120 GeV. A significant background from DY+jets events is reduced by requiring the invariant mass of the electron or muon and the τ h , m τ h to be above 50 GeV. To reduce background from QCD multijet events, we require 2.0 < ∆R( , τ h ) < 3.5.
With these preselection criteria in place, we train several BDTs corresponding to different signal hypotheses to classify signal and background events. The input variables are the p T of the electron or muon, the p T of the τ h candidate, p miss . We also include m T2 and the contransverse mass (m CT ) [67,68], computed from the visible decay products and defined as For signal events, m CT is expected to have an endpoint near (m with ζ being the bisector of the directions of the transverse momenta of the electron or muon and the τ h candidate [69,70]. The value of 0.85 reflects an optimization to efficiently distinguish DY+jets events from other backgrounds and the signal. Figure 3 shows the distributions of events passing the baseline selections in the µτ h final state in two of the BDT input variables that provide the highest discriminating power, p miss T and m tot T . The distributions observed in the eτ h final state are similar. Since the signal kinematics depend on mass, we train BDTs for signals with τ masses of 100, 150, and 200 GeV. In all cases we use a χ 0 1 mass of 1 GeV. As the results of the training depend critically on the number of input events, we relax the τ h MVA-based isolation criteria and reduce the p T threshold for the τ h to 20 GeV for the training sample in order to increase the number of training and test events. The "very tight" isolation and a p T threshold of 30 GeV for the τ h are applied in the final analysis. For a given signal hypothesis, we choose the BDT trained with the same τ mass for models with τ masses of 100, 150, and 200 GeV, or the one that provides optimal sensitivity for models with other τ mass values. For signal models with τ masses of 90 and 125 GeV, we use the BDT trained for m( τ ) = 100 GeV, while for those with a τ mass of 175 GeV, we use the BDT trained for m( τ ) = 200 GeV. While signal events are largely expected to have high BDT output values, we include the full BDT distribution in a binned fit for the statistical interpretation of the analysis as described in Section 7. The binning is chosen to optimize signal significance. Other SM Bkg. uncertainty Obs. / Pred.

Background estimation
Our most significant backgrounds are from DY+jets, W+jets, QCD multijet, tt, and diboson processes. They have relative contributions that vary with final state. For the τ h τ h final state, the dominant background arises from the misidentification of jets as τ h candidates in QCD multijet and W+jets events, constituting ≈65% of background after the baseline selection. For the τ h final states after the baseline selection, the main backgrounds are from DY+jets (≈50%), W+jets (≈30%), and QCD multijet (≈10%) events. The DY+jets contribution, which is also a major background in the τ h τ h final state (≈20%), usually consists of events with two prompt τ leptons. This background is determined with simulation samples after applying corrections to match the normalization and to be consistent with variable distributions in collider data. The W+jets and QCD multijet backgrounds usually contain one or more jets misidentified as τ h and their contributions are determined via methods that rely on data. Finally, we have smaller contributions from other SM processes such as the production of Higgs bosons, dibosons, and top quark pairs with or without vector bosons. These are estimated via MC simulation with appropriate correction factors applied as described in Section 3. For the τ h analyses, dedi-cated CRs that are each enriched in one of the major background processes are used to validate the modeling of the BDT distribution and to extract uncertainties that are used to account for any potential mismodeling of the distributions in simulation. These CRs are described in the following subsections below.

Misidentified jets in the τ h τ h final state
After requiring two τ h candidates with high p T , events with misidentified τ h candidates are the dominant background in the τ h τ h final state. This background, which originates predominantly from QCD multijet and W+jets production, is predicted by extrapolating the event count in a data sample selected with a relaxed isolation requirement into the SR. The fraction of nonprompt or misidentified τ h candidates selected with the very loose MVA-based isolation working point that also pass the tight DeepPF isolation requirement is measured in a QCD multijetenriched sample of same-charge τ h τ h events. The same-charge τ h τ h events are collected with the same τ h τ h trigger as opposite-charge τ h τ h events to avoid additional trigger-related biases. We also require m T2 to be low (<40 GeV) to reduce potential contributions from signal events. We find that roughly 20% of the same-charge events with misidentified τ h candidates selected with very loose isolation also pass the tight isolation requirement. However, the rate depends on the p T and decay mode (one-or three-prongs) of the τ h candidate, as well as the jet flavor, i.e., whether the misidentified jet originates from the hadronization of light-flavor quarks, heavy-flavor quarks, or gluons. The τ h misidentification rate is therefore measured in bins of p T and decay mode to mitigate the dependence on these factors. The measurement is also binned in the number of primary vertices (N PV ) to capture the effects of pileup. From studies performed with MC simulation samples, a systematic uncertainty of ≈30% is assigned to account for the dependence of the misidentification rate on jet flavor.
Since the isolation efficiency for prompt τ h candidates is only around 70-80%, processes containing genuine τ h candidates can enter the sideband regions in events that are selected with the relaxed isolation requirement. To take this into account when calculating the final background estimate, we define three categories of events with at least two loosely isolated τ h candidates: (i) events in which both τ h candidates pass the tight DeepPF isolation requirement, (ii) events in which one passes and one fails the tight isolation requirement, and (iii) events in which both τ h candidates fail the tight isolation requirement. We then equate the count of events in each of these three event categories to the sum of expected counts for the events with two prompt τ h candidates, two jets misidentified as τ h candidates, or one prompt τ h candidate and one jet misidentified as a τ h candidate, that contribute to each category. The contributions from backgrounds with one or two jets misidentified as τ h candidates in the SRs are then determined analytically by solving a set of linear equations.

Misidentified jets in the eτ h and µτ h final states
The misidentification of jets as τ h candidates also gives rise to a major source of background in the eτ h and µτ h final states that arises mainly from W+jets events with leptonic W boson decays. We estimate this background from a sideband region in data selected using the SR selection criteria, with the exception that the τ h candidates are required to satisfy the loose isolation working point and not the very tight working point. A transfer factor for the extrapolation of event counts from this τ h -isolation range into the tight isolation range of the SR is determined with a W+jets CR selected from events with one muon and at least one τ h candidate that passes the loose isolation requirement. In events with more than one τ h candidate, the candidate with the highest value of the MVA-based isolation discriminant is used. To increase the purity of W+jets events in this region, we reduce the contribution from tt and QCD multijet events by requiring 60 < m T ( , p miss T ) < 120 GeV, p miss T > 40 GeV, no more than two jets, and an azimuthal separation of at least 2.5 radians between any jet and the W boson reconstructed from the muon and p miss T (∆φ(W, jet) > 2.5). We also reject events with additional electrons or muons satisfying looser identification criteria. The remaining sample has an expected purity of ≈85% for W+jets events. The transfer factor, R, is then determined from this control sample after subtracting the remaining non-W+jets background contributions estimated from simulation, as follows: where N CR data corresponds to the number of events in the CR in data. The parenthetical argument VT denotes events in which the τ h candidate satisfies the very tight isolation working point, while LVT denotes those that satisfy the loose, but not the very tight requirement. Transfer factors are determined separately in bins of p T and η of τ h candidates.
The contribution of the background originating from a jet misidentified as a τ h candidate in the SR is then determined from the corresponding sideband in data: where N sideband data is the number of events in the sideband in data, from which N sideband MC,τ , the number with genuine τ leptons as estimated with MC simulation by generator-level matching, is subtracted. We validate the estimation of jets misidentified as τ h in a CR requiring 60 < m T ( , p miss T ) < 120 GeV and ∆φ(W, jet) < 2.5 to ensure that the region is independent of the region described above that is used to estimate the background.

Estimation of background from Drell-Yan+jets
The DY+jets background comes primarily from Z → τ + τ − decays. We estimate this contribution via simulation, after applying corrections based on CRs in data. Mismodeling of the Z boson mass or p T distribution in simulation can lead to significant differences between data and simulation in kinematic discriminant distributions, especially when considering the large values of these variables that are relevant for the τ h τ h SRs. We therefore use a high-purity Z → µ + µ − CR to compare the dimuon mass and p T spectra between data and simulation and use the observed differences to correct the simulation in the SRs with weights parameterized by generator-level Z boson mass and p T . The correction factors range up to 30% for high-mass and high-p T values. Because these factors are intended to compensate for missing higher-order effects in the simulation, we assign the differences between the generator-level Z boson mass and p T distributions in LO and NLO simulated events as systematic uncertainties. The differences between data and simulation are taken into account through the use of scale factors, as described in Section 3. The uncertainties in these corrections are propagated to the final background estimate. The corrected simulation is validated in the τ h τ h final state using a Z → τ + τ − CR selected by inverting the m T2 and Σm T requirements used to define the SRs. In addition, requiring a p T of at least 50 GeV for the τ h τ h system reduces the QCD multijet background and improves the purity of this CR. This choice makes it possible to increase the statistical power of this region by removing the p miss T > 50 GeV requirement. The visible mass distribution of the τ h τ h system shown in Fig. 4 (left) demonstrates that the corrected simulation agrees with the data within experimental uncertainties.
For the analysis in the τ h final states, a normalization scale factor, as well as corrections to the p T distribution of the Z boson in simulation are obtained from a very pure Z → µ + µ − CR in data. These events are selected by requiring two isolated muons and no additional leptons, at most one jet, no b-tagged jets, and a dimuon mass in a window of 75-105 GeV, to increase the probability to >99% that they originate from Z → µ + µ − decays. After subtracting all other contributions estimated from simulation, a normalization scale factor that is compatible with unity is extracted from the ratio of data to simulated events. The uncertainty in the scale factor is determined by varying systematic uncertainties associated with objects such as the muon efficiency and jet energy uncertainties.
To validate the DY+jets background prediction in the τ h analyses, we construct a CR in µτ h events with m T (µ, p miss T ) < 20 GeV, 50 < m(µτ h ) < 80 GeV, and N j = 0. These requirements are chosen to obtain a Z → τ + τ − sample with good purity. The m(µτ h ) range is chosen to select the Z boson peak, low m T (µ, p miss T ) helps to remove W+jets and potential signal contamination while the 0-jet requirement helps remove other backgrounds. The p miss T distribution of these events is shown in Fig. 4 (right). We observe good agreement between data and the predicted background.

Estimation of other backgrounds
Smaller contributions are expected from other SM backgrounds, including diboson, triboson, and Higgs boson production. There are also contributions from tt and single top quark production, or top quark pair production in association with a vector boson. These are estimated via MC simulation after application of efficiency and energy-scale corrections. Experimental and theoretical uncertainties are evaluated as described below in Section 6.
For the τ h analyses, we check the BDT distribution in a tt-enriched CR that is defined by requiring the event selection to be the same as in the SR, except for a requirement of one or two b-tagged jets. To validate the WW background prediction, we construct a CR of events with oppositely charged muon-electron pairs that have m µe > 90 GeV and N j = 0. We obtain systematic uncertainties for the normalization of the corresponding backgrounds and any potential mismodeling of the BDT distribution in these CRs. The latter is done by constructing a χ 2 test for all CRs with the BDT modeling taken into account by including an additional floating uncertainty that is determined by requiring a p-value [71] of at least 68% in all CRs. In this way, the BDT shape uncertainty is estimated to be 9%.

Systematic uncertainties
The dominant uncertainties in this analysis are the statistical uncertainties resulting from limited event counts in data sidebands or in simulated event samples used to obtain background estimates and the systematic uncertainties in the estimated rates for jets to be misidentified as τ h candidates. We rely on an extrapolation in τ h isolation to obtain an estimate of the background originating from jets misidentified as τ h candidates. In the τ h τ h analysis, the uncertainty in this extrapolation is dominated by the dependence of isolation on jet flavor. It also includes the statistical uncertainty associated with the CR samples from which the extrapolation factors are obtained. The uncertainty in the combined identification and isolation efficiency for prompt τ h candidates is also propagated to the final estimated uncertainty. In the τ h analyses, we estimate a transfer factor for the extrapolation in τ h isolation from a W+jets-enriched CR. The purity of W+jets events this region is ≈85% as determined from simulation. We therefore propagate a relative uncertainty of 15% to account for contamination from other sources.
We use simulation to obtain estimates of the yields from other background contributions and to estimate the potential signal contributions. We propagate uncertainties related to the b tagging, trigger, and selection efficiencies, the renormalization and factorization scales, PDFs, jet energy scale and resolution, unclustered energy contributing to p miss T , and the energy scales of electrons, muons, and τ h candidates. For the DY+jets background, we have an additional uncertainty associated with the corrections applied to the mass and p T distributions. We assign a 15% normalization uncertainty in the τ h τ h final state for the cross sections of processes estimated from simulation, namely DY+jets, tt, diboson, and rare SM processes, based on the results of CMS differential cross section measurements [72,73]. For the τ h analyses, we extract normalization uncertainties of 5, 5, and 20% for the DY+jets, tt, and WW backgrounds, respectively, based on the estimated impurity of the corresponding process-enriched CRs. An additional uncertainty of 9% is assigned to cover potential mismodeling of the BDT distribution in simulation that is based on studies in CRs.
The categorization of events in the τ h τ h final state by the number of reconstructed jets induces sensitivity to the modeling of ISR in the signal simulation. The p ISR T distribution of simulated signal events is reweighted to improve the ISR modeling. The reweighting factors are obtained from studies of Z boson events. We take the deviation of the reweighting factors from unity as a systematic uncertainty.
The uncertainty in the integrated luminosity is taken into account in all background estimates for which we do not extract normalization scale factors in dedicated data CRs, as well as for signal estimates. This uncertainty corresponds to 2.5% [74] and 2.3% [75] for the 2016 and 2017 data, respectively. The main systematic uncertainties for signal and background are summarized in Table 2.
In general, we treat all statistical uncertainties as uncorrelated. In addition, all systematic uncertainties arising from statistical limitations in the 2016 and 2017 data are assumed to be uncorrelated while systematic uncertainties from similar sources are treated as correlated or partially correlated across the various background and signal predictions. For the combination of the τ h τ h and τ h analyses, we correlate uncertainties related to object reconstruction, with the exception of the τ h selection efficiency, which is treated as uncorrelated because of the use of different isolation algorithms.

Results and interpretation
The results of the search in the τ h τ h final state are presented in Fig. 5 and summarized in Tables 3 and 4. The background predictions resulting from a maximum likelihood fit to the data under the background-only hypothesis are shown in the lower row of Fig. 5. The BDT distributions corresponding to a training for a τ mass of 100 GeV and a χ 0 1 mass of 1 GeV are shown before and after the maximum-likelihood fit to the data in Figs. 6 and 7 for the µτ h and eτ h final states, respectively. The data are consistent with the prediction for SM background. The predicted and observed event yields in the last, most sensitive BDT bins are summarized in Tables 5 and 6 for τ h final states. For the statistical interpretation of these results, the normalization uncertainties affecting background and signal predictions are generally assumed to be log-normally distributed. For statistical uncertainties limited by small event counts in data or simulation, we use a Γ distribution.
The results are used to set upper limits on the cross section for the production of τ pairs in the context of simplified models [27][28][29]76] using all of the exclusive τ h τ h SRs and the τ h BDT distributions in a full statistical combination. The limits are evaluated using likelihood fits with the signal strength, background event yields, and nuisance parameters corresponding to the uncertainties in the signal and background estimates as fitted parameters. The nuisance parameters are constrained within their uncertainties in the fit. We assume that the τ decays with 100% branching fraction to a τ lepton and a χ 0 1 . The 95% CL upper limits on SUSY production cross sections are calculated using a modified frequentist approach with the CL s criterion [77,78]. An asymptotic approximation is used for the test statistic [79,80], q µ = −2 ln L µ /L max , where L max is the maximum likelihood determined by allowing all fitted parameters, including the signal strength µ, to vary, and L µ is the maximum likelihood for a fixed signal strength. Figure 8 shows the limits obtained for purely left-handed τ pair production, while Fig. 9 shows the limits obtained for the degenerate τ model in which both left-and right-handed τ pairs are produced. The τ h τ h analysis makes the dominant contribution to the search sensitivity. A slight excess of events over the background expectation in the τ h τ h SRs results in an observed limit that is weaker than the expected limit. The strongest limits are observed in the case of a nearly massless χ 0 1 . In general, the constraints are weaker for higher values of the χ 0 1 mass because of smaller experimental acceptances. In the purely left-handed model, the strongest limits are observed for a τ mass of 125 GeV where we exclude a τ pair production cross section of 132 fb. This value is a factor of 1.14 larger than the theoretical cross section. In the degenerate τ model we exclude τ masses between 90 and 150 GeV under the assumption of a nearly massless χ 0 1 .

Summary
A search for direct τ slepton ( τ ) pair production has been performed in proton-proton collisions at a center-of-mass energy of 13 TeV in events with a τ lepton pair and significant missing transverse momentum. Search regions are defined using kinematic observables that exploit For the background estimates with no events in the sideband or in the simulated sample, we calculate the 68% CL upper limit on the yield. The first and second uncertainties given are statistical and systematic, respectively. We also list the predicted signal yields corresponding to the purely left-handed model for a τ mass of 100 GeV and a χ 0 1 mass of 1 GeV.

Acknowledgments
We congratulate our colleagues in the CERN accelerator departments for the excellent performance of the LHC and thank the technical and administrative staffs at CERN and at other CMS institutes for their contributions to the success of the CMS effort. 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.  For the background estimates with no events in the sideband or in the simulated sample, we calculate the 68% CL upper limit on the yield. The first and second uncertainties given are statistical and systematic, respectively. We also list the predicted signal yields corresponding to the purely left-handed model for a τ mass of 100 GeV and a χ 0 1 mass of 1 GeV.     Table 6: Predicted background yields and observed event counts in the most sensitive last bins of the BDT distributions in the eτ h and µτ h final states, in data collected in 2017. The numbers in parentheses in the first row are the τ and χ 0 1 masses corresponding to the signal model for left-handed τ pair production that is used to train the BDT. In the bottom row, we list the corresponding predicted signal yields in the last bin of the BDT distribution. The first and second uncertainties given are statistical and systematic, respectively. BDT           [40] CMS Collaboration, "The CMS trigger system", JINST 12 (2017) P01020, doi:10.1088/1748-0221/12/01/P01020, arXiv:1609.02366.
[52] J. Alwall et al., "The automated computation of tree-level and next-to-leading order differential cross sections, and their matching to parton shower simulations", JHEP 07