EUROPEAN ORGANIZATION FOR NUCLEAR RESEARCH (CERN)

 LHCb-DP-2018-001

Selection and processing of calibration samples to measure the particle
identification performance of the LHCb experiment in Run 2

R. Aaij8, L. Anderlini5, S. Benson9, M. Cattaneo8, P. Charpentier8, M.
Clemencic8, A. Falabella3, F. Ferrari,3, M. Fontana4,8, V. V. Gligorov2, D.
Hill13, T. Humair12, C. R. Jones10, O. Lupton8, S. Malde13, C. Marin Benito7,
R. Matev8, A. Pearce8, A. Poluektov11, B. Sciascia6, F. Stagni8, R. Vazquez
Gomez8, Y.-X. Zhang1

1Laboratoire de l’Accelerateur Lineaire, Paris, France
2LPNHE, Université Pierre et Marie Curie, Université Paris Diderot,
CNRS/IN2P3, Paris, France
3Istituto Nazionale di Fisica Nucleare, Sezione di Bologna, Italy
4Istituto Nazionale di Fisica Nucleare, Sezione di Cagliari, Italy
5Istituto Nazionale di Fisica Nucleare, Sezione di Firenze, Italy
6Istituto Nazionale di Fisica Nucleare, Laboratori Nazionali di Frascati,
Italy
7ICCUB, Universitat de Barcelona, Barcelona, Spain
8European Organization for Nuclear Research (CERN), Meyrin, Switzerland
9Nikhef, Amsterdam, Netherlands
10University of Cambridge, Cambridge, United Kingdom
11Department of Physics, University of Warwick, Coventry, United Kingdom
12Imperial College London, London, United Kingdom
13University of Oxford, Oxford, United Kingdom

Since 2015, with the restart of the LHC for its second run of data taking, the
LHCb experiment has been empowered with a dedicated computing model to select
and analyse calibration samples to measure the performance of the particle
identification (PID) detectors and algorithms. The novel technique was
developed within the framework of the innovative trigger model of the LHCb
experiment, which relies on online event reconstruction for most of the
datasets, reserving offline reconstruction to special physics cases. The
strategy to select and process the calibration samples, which includes a
dedicated data-processing scheme combining online and offline reconstruction,
is discussed. The use of the calibration samples to measure the detector PID
performance, and the efficiency of PID requirements across a large range of
decay channels, is described. Applications of the calibration samples in data-
quality monitoring and validation procedures are also detailed.

Submitted to EPJ TI

© CERN on behalf of the LHCb collaboration, licence CC-BY-4.0.

##  1 Introduction

LHCb is a dedicated heavy flavour physics experiment at the LHC. Its main goal
is to search for indirect evidence of new physics in CP-violating processes
and rare decays of beauty and charm hadrons. Among other performance metrics,
like excellent vertex resolution and good momentum and invariant-mass
resolution, charged particle identification (PID) distinguishing electrons,
muons, pions, kaons and protons traversing the detector is essential in the
LHCb physics programme. The required performances range from the per mille
misidentification probability of hadrons as muons in the study of the rare
$B_{(d,s)}\rightarrow\mu^{+}\mu^{-}$ decays [1, 2, 3, 4, 5, 6, 7], to the sub
percent precision, over a wide kinematic range accurate, on the detector
induced asymmetries for the ambitious programme of CP asymmetry measurements
[8, 9, 10, 11].

PID information is extensively used both in the trigger selection and in
offline data analysis. This required the development of a dedicated computing
model and a strategy to select suitable calibration samples, in order to
measure the PID performance and assess systematic effects.

In Section 2, an overview of the LHCb detector is given, together with a
summary of the PID calibration samples required in order to accomplish the
physics goals of LHCb with Run 2 data (2015$-$2018). The article then focuses
on the strategy to select and process PID calibration samples, including a
description of the multivariate classifiers used to combine the response of
calorimeters, RICH and muon system (Section 3); the procedure to measure the
PID performance using dedicated calibration samples, together with the
techniques to determine the selection efficiency on hundreds of different
decay channels, relying on a small number of calibration samples (Section 4);
the dedicated data-processing scheme combining online and offline
reconstruction (Section 5); and the applications of the calibration samples to
data-quality monitoring and validation (Section 6). A brief summary and
outlook are given in Section 7. While this article discusses the calibration
samples specifically for charged particle identification, the general
computing model and selection strategy is also being applied to other
calibration samples in Run 2, such as those for tracking calibration and
neutral pion and photon PID.

##  2 Detector

The LHCb detector is a single-arm forward spectrometer covering the
pseudorapidity range $2<\eta<5$, designed for the study of particlses
containing $b$ or $c$ quarks [12, 13]. The detector includes a high-precision
tracking system consisting of a silicon-strip vertex detector surrounding the
$pp$ interaction region [14], a large-area silicon-strip detector located
upstream of a dipole magnet with a bending power of about $4{\mathrm{\,Tm}}$,
and three stations of silicon-strip detectors and straw drift tubes [15]
placed downstream of the magnet. The tracking system provides a measurement of
momentum, $p$, of charged particles with a relative uncertainty that varies
from 0.5% at low momentum to 1.0% at
200${\mathrm{\,Ge\kern-1.00006ptV\\!/}c}$. The minimum distance of a track to
a primary vertex, the impact parameter, is measured with a resolution of
$(15+29/\mbox{$p_{\mathrm{T}}$}){\,\upmu\mathrm{m}}$, where $p_{\mathrm{T}}$
is the component of the momentum transverse to the beam, in
${\mathrm{\,Ge\kern-1.00006ptV\\!/}c}$. Photons, electrons and hadrons are
identified by a calorimeter system consisting of scintillating-pad and
preshower detectors, an electromagnetic calorimeter and a hadronic
calorimeter. Different types of charged hadrons are distinguished using
information from two ring-imaging Cherenkov (RICH) detectors [16]. Muons are
identified by a system composed of alternating layers of iron and multiwire
proportional chambers [17].

The online event selection is performed by a trigger [18], which consists of a
hardware stage, based on information from the calorimeter and muon systems,
followed by a software stage, which applies a full event reconstruction. Since
2015, in between the hardware and software stages, a real-time procedure
aiming at the alignment and calibration of the detector is performed [19],
making use of a disk buffer [20]. Updated calibration parameters are made
available for the online reconstruction, used in the trigger selection. Online
calibration is of such high quality that it is also used for offline
reconstruction, ensuring consistency between online and offline.

The responses of the calorimeter, RICH, and muon systems, or their
combinations, associated to each track in the reconstruction process are named
for brevity PID variables. They can be used in selections to increase the
signal purity of a sample, reducing the processing time devoted to the
reconstruction of background events, often characterized by high multiplicity,
and helping in fitting into the data storage constraints. Moreover they allow
selections to avoid an explicit bias on quantities of physical interest, such
as decay time [21, 22].

The many contexts in which particle identification is exploited within the
experiment and the difficulties in obtaining a perfect simulation for the PID
detectors, motivate the development of techniques for measuring the PID
performance in suitable PID calibration samples. These samples are datasets
collected by LHCb where decay candidates have a kinematic structure that
allows unambiguous identification of one of the daughters, without the use of
any PID information from the calorimeter, RICH, or muon systems, so that they
are unbiased from the particle identification point of view. Today, most LHCb
physics analyses rely on calibration samples for the determination of PID
efficiencies. In addition, these samples can be used to monitor time
variations in performance, and to test new reconstruction algorithms.

The majority of physics analyses using data collected with the LHCb experiment
rely on the physics quantities as reconstructed in the online trigger
reconstruction. Still, physics analyses with special needs in terms of event
reconstruction, searching for example for interactions of light nuclei or
particles beyond the Standard Model with the detector [23], are able to
reprocess offline the collected calibration datasets with dedicated
reconstruction algorithms.

In order to enable the measurement of selection efficiencies that combine
trigger requirements on the online-computed PID variables and offline
requirements of PID variables obtained through dedicated reconstruction
algorithms, an innovative dedicated data-processing strategy has been
designed. Calibration data are obtained through a real-time selection based on
the online reconstruction without any requirement on PID variables. Each event
belonging to the calibration samples is fully reconstructed independently both
online and offline. The resulting reconstructed particles are then matched,
allowing a measurement of the efficiency of requirements that combine the two
reconstruction types as described in Section 4.

##  3 Global particle identification

The reconstruction algorithms of each of the PID detectors of the LHCb
experiment are very different, but each of them allows the computation of a
likelihood ratio between particle hypotheses for each reconstructed track. The
reconstruction algorithm of the RICH detectors provides the likelihood of the
electron, muon, kaon, proton and deuteron hypotheses relative to the pion
hypothesis. The calorimeter system provides the likelihood of electrons
relative to the pion hypothesis. Finally, the muon system provides the
likelihoods of the muon and non-muon hypotheses. The likelihood ratios of the
three detector systems are combined into _Combined Differential Log-
Likelihoods_ (CombDLL) [13], which are used to define the selection criteria
for the data analyses. Selection strategies based on CombDLL and
$\mathtt{isMuon}$ [24], a binary variable loosely identifying muons, are
widely employed already at the trigger level [25].

Following recent developments in machine learning, more advanced classifiers
have also been designed to combine the likelihoods ratios defined above with
the informations from the tracking system, including the kinematic variables
of the particle, and additional information from the PID detectors not
entering the likelihood computation (_e.g._ the number of hits in the muon
system shared among reconstructed tracks). The classifier with the widest
application in this category, named _ANNPID_ , was developed using _Forward
Feeding Artificial Neural Networks_ [26], structured as a _Multi-Layer
Perceptron_ (MLP) with a single hidden layer composed of roughly 20% more
nodes than the input layer activated, through a sigmoid function. The training
is performed minimising the likelihood of the Bernoulli Cross-Entropy with
Stochastic Gradient Decent as implemented in the TMVA package. The training
sample is obtained from abundant simulated decays of heavy hadrons that
emulate the kinematic distributions of signal samples studied in several
analyses. Depending on the arrangement of the input samples, on the quality of
the simulation, and on the available number of simulated events, the response
of the ANNPID algorithm can vary. As a consequence, the response of the ANNPID
algorithms is provided in several tunings, some for general purpose, and
others specialised for a particular analysis or kinematic range. The variables
combined using the ANNPID classifiers are listed in Table 1.

Table 1:  Input variables of the ANNPID classifiers for the various subsystems
of the LHCb detector.  _Tracking_
---
| Total momentum
| Transverse momentum
| Quality of the track fit
| Number of clusters associated to the track
|  ANN response trained to reject ghost tracks [27]
| Quality of the fit matching track segments upstream and downstream of the
magnet
_RICH detectors_
| Geometrical acceptance of the three radiators, depending on the direction of
the track
| Kinematical acceptance due to Cherenkov threshold for muons and kaons
| Likelihood of the electron, muon, kaon, and proton hypotheses relative to
the pion
| Likelihood ratio of the below-threshold and pion hypotheses
_Electromagnetic calorimeter_
| Likelihood ratio of the electron and hadron hypotheses
| Likelihood ratio of the muon and hadron hypotheses
|  Matching of the track with the clusters in the _preshower_ detector
| Likelihood ratio of the electron and pion hypotheses,
|  after recovery of the Bremsstrahlung photons
_Hadronic calorimeter_
| Likelihood ratio of the electron and hadron hypotheses
| Likelihood ratio of the muon and hadron hypotheses
_Muon system_
| Geometrical acceptance
| Loose binary requirement already available in the hardware trigger
| Likelihood of the muon hypothesis
| Likelihood of the non-muon hypothesis
| Number of clusters associated to at least another tracks

All of the input variables for the ANNPID classifiers are made immediately
available to physics analyses, easing the development of new tunings and
classification algorithms dedicated to single analyses. The many output
variables of the detector reconstruction which are not used as input to ANNPID
can be accessed or even regenerated, relying on the raw detector data stored
on tape.

##  4 Measuring PID performance

More than twenty exclusive trigger selections are designed to select pure
samples of the five most common charged particle species that interact with
the LHCb detector: protons, kaons, pions, muons and electrons [28]. Generally,
low-multiplicity decay modes with large branching fractions are chosen in
order to enhance the statistics and the purity and populate the tails in the
distributions of the PID variables, which are of great relevance when
computing misidentification probabilities. Completely reconstructed final
states composed of charged particles only are preferred, as they are selected
with high purity at LHCb. An overview of the modes utilised is given in Table
2.

The assumption underlying the usage of the calibration samples is that the
distribution of the particle identification variables is independent of the
selection strategy. Simply avoiding explicit requirements on the PID variables
is not sufficient to ensure this. In fact, the hardware trigger relies on
information from the CALO and muon systems to reduce the rate at which the
full detector is read out to around $1{\mathrm{\,MHz}}$, while a first layer
of the software trigger, running before the full event reconstruction,
includes dedicated selection algorithms to identify high $p_{\mathrm{T}}$
muons and muon pairs.

In order to avoid a pre-selection introducing biases in the PID variables, the
selection strategy of the calibration samples imposes requirements on the
algorithms selecting the event in the previous trigger layers. Either the
trigger algorithms do not rely on PID information, or the PID selection in the
trigger is applied to one of the particles not used to measure the
performance.

Table 2: Overview of decay modes that are used to select calibration samples. The high momentum samples are primarily selected, while the low momentum samples are included to maximise the kinematic coverage as much as possible.”  Species | Low momentum | High momentum
---|---|---
$e^{\pm}$ |  ${{{B}^{+}}}\\!\rightarrow{{J\mskip-3.0mu/\mskip-2.0mu\psi\mskip 2.0mu}}{{K}^{+}}$ with ${{J\mskip-3.0mu/\mskip-2.0mu\psi\mskip 2.0mu}}\\!\rightarrow{e^{+}e^{-}}$
${\mu}^{\pm}$ |  ${{{B}^{+}}}\\!\rightarrow{{J\mskip-3.0mu/\mskip-2.0mu\psi\mskip 2.0mu}}{{K}^{+}}$ with ${{J\mskip-3.0mu/\mskip-2.0mu\psi\mskip 2.0mu}}\\!\rightarrow{\mu^{+}\mu^{-}}$ | ${{J\mskip-3.0mu/\mskip-2.0mu\psi\mskip 2.0mu}}\\!\rightarrow{\mu^{+}\mu^{-}}$
${\pi}^{\pm}$ | ${{K}^{0}_{\mathrm{\scriptscriptstyle S}}}\\!\rightarrow{{\pi}^{+}}{{\pi}^{-}}$ |  ${{D}^{*+}}\\!\rightarrow{{D}^{0}}{{\pi}^{+}}$ with ${{D}^{0}}\\!\rightarrow{{K}^{-}}{{\pi}^{+}}$
${K}^{\pm}$ |  ${{D}^{+}_{s}}\\!\rightarrow\phi{{\pi}^{+}}$ with $\phi\\!\rightarrow{{K}^{+}}{{K}^{-}}$ |  ${{D}^{*+}}\\!\rightarrow{{D}^{0}}{{\pi}^{+}}$ with ${{D}^{0}}\\!\rightarrow{{K}^{-}}{{\pi}^{+}}$
$p$, $\overline{p}$ | ${{\mathchar 28931\relax}^{0}}\\!\rightarrow{{p}}{{\pi}^{-}}$ |  ${{\mathchar 28931\relax}^{0}}\\!\rightarrow{{p}}{{\pi}^{-}}$ ; ${{\mathchar 28931\relax}^{+}_{c}}\\!\rightarrow{{p}}{{K}^{-}}{{\pi}^{+}}$

Several of the selection strategies are implemented according to the so-called
_tag-and-probe_ model. Taking the ${{J\mskip-3.0mu/\mskip-2.0mu\psi\mskip
2.0mu}}\\!\rightarrow{\mu^{+}\mu^{-}}$ decay as an example, the tag-and-probe
selection strategy relies on a list of well-identified tag muons of a certain
charge and a list of probe tracks with opposite charge, selected avoiding any
PID requirement. These are combined to form muon pairs with invariant-mass
consistent with the ${J\mskip-3.0mu/\mskip-2.0mu\psi\mskip 2.0mu}$ mass, and
are then filtered further on the basis of the quality of the fit of the decay
vertex, to form the final sample. To extend the $p_{\mathrm{T}}$ range of the
muons in the calibration samples to lower values, where the background from
low momentum pions is difficult to reduce, the
${J\mskip-3.0mu/\mskip-2.0mu\psi\mskip 2.0mu}$ candidates can be combined with
charged kaons to form
${{B}^{+}}\rightarrow{{J\mskip-3.0mu/\mskip-2.0mu\psi\mskip 2.0mu}}{{K}^{+}}$
candidates111Charged-conjugated candidates are implicitly considered here and
throughout the paper. , adding further kinematic constraints related to the
$B$ decay to the final filtering.

Proton calibration samples are obtained from two different decay modes:
${{\mathchar 28931\relax}^{0}}\\!\rightarrow{p}{{\pi}^{-}}$ and ${{\mathchar
28931\relax}^{+}_{c}}\\!\rightarrow{p}{{K}^{-}}{{\pi}^{+}}$. Since the visible
${\mathchar 28931\relax}^{0}$ production cross section in LHCb is very high,
the yield collected at the trigger level exceeds the needs in terms of
statistical precision on the particle identification. This would pose severe
challenges for data storage. Therefore, a large fraction of these signal
candidates is discarded by running the selection only on a randomly selected
fraction of the events. In order to improve the kinematic coverage of the
sample, the fraction of discarded events is defined differently in four bins
of the proton transverse momentum ($p_{\mathrm{T}}$), resulting in a higher
retention rate in the less-populated high-$p_{\mathrm{T}}$ region. The sample
of ${\mathchar 28931\relax}^{+}_{c}$ decays is included to extend the
$p_{\mathrm{T}}$ coverage of the ${\mathchar 28931\relax}^{0}$ samples.

An abundant calibration sample for pions is provided by the decay
${{K}^{0}_{\mathrm{\scriptscriptstyle S}}}\rightarrow{{\pi}^{+}}{{\pi}^{-}}$,
but the spectrum of the probe particles is much softer than what is typical
for hadrons produced in heavy hadron decays. Charm hadron decays allow the
kinematic range to be extended to higher transverse momenta, but the lower
purity of the samples, due to the smaller production cross-section, requires
additional care in the selection and background subtraction strategies.The
decay ${{D}^{*+}}\\!\rightarrow{{D}^{0}}{{\pi}^{+}}$ with
${{D}^{0}}\\!\rightarrow{{K}^{-}}{{\pi}^{+}}$ represents the primary source of
${\pi}^{\pm}$ and ${K}^{\pm}$ calibration samples. The soft pion produced in
the strong decay of the ${D}^{*+}$ hadron allows to tag the flavour of the
${D}^{0}$ and therefore to distinguish the kaon and the pion produced in its
decay without PID requirements on either of the two probe particles. Applying
a requirement on the energy release in the $D^{*}$ decay, which is expected to
be small, enables the rejection of combinatorial background due to the
erroneous combination of ${D}^{0}$ hadrons and pions produced in unrelated
processes. Finally, the ${{D}^{+}_{s}}\\!\rightarrow\phi{{\pi}^{+}}$ decay
with $\phi\\!\rightarrow{{K}^{+}}{{K}^{-}}$ is a further source of kaons. This
sample allows the kinematic range for kaons to be extended to lower momenta,
as the $\phi$ constraint enables the kinematic requirements on the kaons to be
loosened while retaining the purity.

Figure 1:  On the left, mass distributions of the decaying particles with the
results of the fit superimposed; signal contributions are shown by the red
dashed curves, and the total fit functions including background contributions
are shown by the blue solid curves. On the right, the background-subtracted
distributions of the calibration samples for electrons, muons, pions, kaons
and protons as a function of the track pseudorapidity, $\eta$, and momentum
$p$ are shown. The colour scale units are arbitrary.

The residual background that cannot be rejected with an efficient selection
strategy is statistically subtracted using the ${}_{s}\mathcal{P}\mathrm{lot}$
technique [29]. A fit to the invariant-mass of the decaying particle is
performed for each calibration sample, defining a signal component for which
the sample of probe tracks is known to be pure, and one or more background
components of different nature. In several cases, two-dimensional fits are
performed to account for additional background sources. The variables used in
the two-dimensional fits are: the ${{D}^{0}}$ mass and ${{D}^{*+}}-{{D}^{0}}$
mass difference for ${{D}^{*+}}\rightarrow{{D}^{0}}{{\pi}^{+}}$; the
${{{B}^{+}}}$ and $J/\psi$ masses for
${{{B}^{+}}}\rightarrow(J/\psi\rightarrow{\mu^{+}}\mu^{-}){{K}^{+}}$; the
$\phi$ and ${{D}^{+}_{s}}$ masses for
${{D}^{+}_{s}}\rightarrow\phi{{\pi}^{+}}$.

The fit to each calibration sample is repeated (at least) twice; the first
iterations have a large number of free parameters including the means and
widths of the signal components and shape parameters of the background
components, whereas the final iteration fixes all of the parameters apart from
the normalisation of each component (such as signal, misidentified background
and combinatorial background). The covariance matrix produced in the final fit
is used to define a relation between the discriminant variables and a signal
${}_{s}\mathcal{W}\mathrm{eight}$ to be assigned to the daughter candidate.
Correlations between the chosen discriminating variables and the PID variables
do not play a significant role. Figure 1 illustrates the invariant-mass
distributions for some of the calibration samples, as obtained from proton-
proton collision data collected in 2015 corresponding to an integrated
luminosity of 0.17$\mbox{\,fb}^{-1}$. The corresponding kinematic
distributions for the different species of probe particles are also shown.

The response of the PID detectors to a traversing particle depends on the
kinematics of the particle, the occupancy of the detectors (which may be
different event-to-event and for different particle production mechanisms),
and experimental conditions such as alignments, temperature, and gas pressure
(which may modify the response of detectors across runs).

One may assume that the response of a PID variable is fully parameterised by
some known set of variables, such as the track momentum $p$ (which is related
to the Cherenkov angle in the RICH and to the energy deposited in the
calorimeter ) and the track multiplicity, the latter being given by the number
of reconstructed tracks traversing the whole detector. By partitioning the
sample with sufficient granularity in these parameterising variables, the PDF
of the PID variable distribution does not vary significantly within each
subset, such that the efficiency of a selection requirement on that variable
is constant within each subset.

In the trivial case of events that come from the calibration sample, there is
no need to compute per-subset efficiencies, and the average efficiency is
simply given by the fraction of background subtracted events passing the PID
requirement. To compute the PID efficiency on a sample other than the
calibration sample, denoted hereafter as the _reference sample_ , the
parameterising variables in the calibration sample can be weighted to match
those in the reference sample. The PID efficiency can then be computed using
the per-subset weights. The weights are defined as the normalised ratio of
reference to calibration tracks

$w_{i}=\frac{R_{i}}{C_{i}}\times\frac{C}{R}\;,$ (1)

where $R_{i}$ ($C_{i}$) is the number of reference (calibration) tracks in the
$i\text{th}$ subset, and $R$ ($C$) is the total number of reference
(calibration) tracks in the sample.

After applying the PID cut to the weighted calibration sample, the average
efficiency of the PID requirement on the weighted calibration sample is

$\bar{\varepsilon}=\frac{\sum_{i}\varepsilon_{i}w_{i}C_{i}}{\sum_{i}w_{i}C_{i}}\;.$
(2)

where $w_{i}$ is the per-subset weight, $\epsilon_{i}$ is the per-subset
efficiency and $C_{i}$ is the number of calibration tracks in the $i$-th
subset.

The computation of the PID efficiency can be thought of as the reweighting of
the calibration sample to match the reference, or as the assignment of
efficiencies to reference tracks based on the subset they belong to. This can
also be extended to reference samples where PID requirements have been imposed
on multiple tracks, where the efficiency of an ensemble of cuts is required
taking into account the kinematic correlation between tracks.

There are a number of ways in which the calibration samples can be used to
determine PID efficiencies. Three broad strategies have been commonly
implemented by LHCb in the past. The first uses a simulated reference sample
to provide the kinematics of the signal tracks under consideration. This is an
ideal approach to use when the kinematics of the signal tracks are known to be
well modelled in the simulation. If the signal in data can be reliably
separated from the other species in the sample, such that some background
subtraction can be used to extract the signal kinematics, a second approach to
creating the reference sample can be used.

Lastly, the PID response of MC signal samples can be corrected using the PID
calibration data samples. Two options are provided:

  * •

Resampling of PID variables, where the PID response is completely replaced by
the one statistically generated from calibration PDFs.

  * •

Transformation of PID variables, where the PID variables from the simulation
are transformed such that they are distributed as in data.

The PID correction is still considered as a function of track kinematics
($p_{\mathrm{T}}$ and $\eta$) and event multiplicity $N_{\rm evt}$ (such as
the number of tracks in the event). However, unlike in the first two
strategies detailed above, the correction is performed using an unbinned
approach, where the calibration PDFs in four dimensions, the PID variable,
$p_{\mathrm{T}}$, $\eta$, and a measure of $N_{\rm evt}$, are described by a
kernel density estimation procedure using the Meerkat library [30]. The
advantage of resampling and variable transformation is that the corrected PID
response can be used as an input to a multivariate classifier.

However, a limitation of the PID resampling approach is that the PID variables
for the same track are generated independently, and thus no correlations
between them are reproduced. Therefore, only one PID variable per track can be
used in the selection. Correlations between variables for different tracks are
preserved via correlations with the kinematics of tracks, assuming the PID
response is fully parameterised by $p_{\mathrm{T}}$, $\eta$, and $N_{\rm
evt}$.

The PID variable transformation approach aims to remove this limitation [31].
The corrected PID variable ${\rm PID}_{\rm corr}$ is obtained as

${\rm PID}_{\rm corr}=P^{-1}_{\rm exp}(P_{\rm MC}({\rm PID}_{\rm
MC}|\mbox{$p_{\mathrm{T}}$},\eta,N_{\rm
evt})|\mbox{$p_{\mathrm{T}}$},\eta,N_{\rm evt}),$ (3)

where $P_{\rm MC}({\rm PID}_{\rm MC}|\mbox{$p_{\mathrm{T}}$},\eta,N_{\rm
evt})$ is the cumulative distribution function of the simulated PID variable
${\rm PID}_{\rm MC}$, and $P^{-1}_{\rm
exp}(x|\mbox{$p_{\mathrm{T}}$},\eta,N_{\rm evt})$ (where $0<x<1$) is the
inverse cumulative distribution function for the PID variable from the
calibration sample (i.e. for fixed $p_{\mathrm{T}}$, $\eta$ and $N_{\rm evt}$
it returns the PID variable that corresponds to a cumulative probability $x$).
Both functions are obtained from the results of kernel density estimations of
the simulation and calibration PID responses, respectively. The corrected PID
variables obtained in this way follow the PDF of the calibration sample, but
preserve strong correlations with the output of simulation. Through these
correlations in simulation, the ones between PID variables for the same track
are reproduced to first order. The drawback of this approach is that it also
relies on the parametrisation of PID PDFs in simulation, which are extracted
from samples that are typically much smaller than the calibration data.
Although one naively expects this method to perform better due to taking
correlations into account, studies are ongoing to quantify the degree of
agreement between the correlations found in simulation and data. The PID
resampling and variables transformation techniques are schematically
represented in Figure 2.

Figure 2:  Schematic representation of the PID resampling and variable
transformation techniques.

There are a number of sources of uncertainty that affect the measurement of
PID efficiencies. The statistical uncertainty arises from finite statistics in
the input samples used in the calibration procedure, namely the calibration
and reference samples. Due to the large calibration sample sizes, this
uncertainty is usually dominated by the size of the signal reference sample.

Several sources of systematic uncertainty related to the procedure must also
be accounted for, arising from differences between the reference and signal
samples, the specific choice of binning used, and the
${}_{s}\mathcal{W}\mathrm{eight}$ procedure used in the calibration sample
production. The degree to which these uncertainties affect the PID efficiency
precision is analysis dependent, and require specific studies to be carried
out on a case-by-case basis. Moreover the availability of primary and
secondary calibration samples allows to study possible biases coming from
single decay modes.

##  5 Computing model for the calibration samples

Figure 3:  Schematic representation of the computing model for the PID
calibration samples.

In order to face the new challenges of the second run of the LHC, the LHCb
trigger [32, 33] has evolved into a heterogeneous configuration with different
output data formats for different groups of trigger selections. Figure 3 shows
a schematic representation of the computing model that is described in the
following.

Two alternative data formats for physics analyses, named Turbo stream [34] and
Full stream, have been developed. Trigger selections writing to the Turbo
stream are intended for analyses of samples where only the information related
to the candidates and associated reconstructed objects is needed. Trigger
selections that are part of the Turbo stream produce a decay candidate which
is stored for offline analysis, along with a large number of detector-related
variables, while the raw detector data is not kept [35, 36]. When considering
analyses based on the Turbo stream, it is therefore evident that the
calibration samples must provide the PID information as computed online in
order to assess the efficiency of selection requirements applied either in the
trigger selection, or offline on the PID variables retrieved from the trigger
candidate.

Trigger selections for events to be stored in the Full stream are intended for
those measurements and searches for which the Turbo approach is not applied.
While the software trigger fully reconstructs candidates, those are not saved.
If the trigger decision is affirmative, the raw detector data is saved
together with summary information on the trigger decision, including the
CombDLL and isMuon variables, for each particle involved in the trigger
decision. The track and decay candidates are reproduced in a further offline
reconstruction step that accesses the raw detector data. Indeed, some physics
data analyses present special needs in terms of particle identification
algorithms, for example because they explore kinematic regions at the
boundaries of the kinematic acceptance, or because of exceptional requirements
in terms of the accuracy of the efficiency determination. To respond to such
special requirements, dedicated algorithms accessing the raw detector data can
be developed and included in the offline event reconstruction. Hence, the
events selected as part of the calibration samples must include the raw data,
allowing the performance of future algorithms to be measured on data.

An interesting case is presented when a trigger selection targeting the Full
stream includes PID requirements that are then intended to be refined offline.
Potentially, the PID variables computed online can differ from those obtained
from the full event reconstruction performed offline. While accidental
differences in the online and offline algorithms are unlikely thanks to
dedicated checks in the data quality validation procedure, the offline
reconstruction is subject to improvements that provide a slightly different
value for the PID variables. The determination of the efficiency of combined
requirements on online and offline versions of the PID variables, or of
different tunings of the multivariate classifiers adopted in the trigger and
in the statistical data analysis, require the use of calibration samples
combining the information from the online and offline reconstruction, allowing
full offline reprocessing if needed.

A dedicated data format, named TurboCalib, was developed to satisfy the
requirements on PID calibration samples described above [37]. After the online
full event reconstruction, events in which decay candidates useful for
calibration are identified and selected in real-time are stored including both
the trigger candidates themselves and raw detector data. The two output
formats are processed independently for each event, to obtain both decay
candidates propagated from the trigger and decay candidates reconstructed
offline from the raw detector data. The two reconstructions are fully
independent, so that the tracks identified in the two processes must be
matched. This is done according to the fraction of shared clusters in the
detector, or exploiting the TisTos algorithm described in Ref. [38], or with a
combination of the two techniques.

The offline versions of the PID variables can be easily replaced with other
tunings of the multivariate classifiers, or through the output of dedicated
reconstruction sequences. As a result of the matching procedure, each
reconstructed track is associated to two sets of PID variables, obtained
through the online and offline versions of the reconstruction, respectively.
The two sets are available to the analysts to measure the efficiency of
selection requirements that possibly combine the two versions.

As described in Section 4, the measurement of the selection efficiencies from
the selected calibration samples is enabled through the subtraction of the
residual background by means of the ${}_{s}\mathcal{P}\mathrm{lot}$ technique.
In order to overcome to the scalability challenges set by the increasing needs
for precision in many LHCb measurements, resulting in huge calibration samples
to control the statistical uncertainty, the background subtraction is
performed through a dedicated, distributed implementation of the
${}_{s}\mathcal{P}\mathrm{lot}$ technique. Finely binned histograms of the
invariant-mass distributions of the trigger candidates are filled in parallel
on thousands of computing nodes. They are then merged and modeled through a
maximum likelihood fit as the combination of signal and background components.
The relations between the discriminating variables and the
${}_{s}\mathcal{W}\mathrm{eights}$ to be assigned to each candidate are
sampled in fine grids and made available through a distributed file system to
the computing nodes of the LHCb grid [39], where jobs to assign the weights
are run as a final processing step in the calibration sample production
workflow. Such a distributed implementation of the
${}_{s}\mathcal{P}\mathrm{lot}$ technique avoids the storage of the entire
dataset on a single computing node, hence scaling better with the size of the
calibration samples.

The real-time selection strategy, the double-processing scheme combining
event-by-event the online and offline reconstructed variables, and the
distributed approach to background subtraction constitute the main novelties
in the data processing for the calibration samples, overcoming most
scalability issues and making the limited cross-section and the available data
storage resources the only limitations to the statistical precision in the
determination of PID selection efficiencies.

Finally, the _Particle IDentification Calibration_ (PIDCalib) package [40] is
a user interface written in python aiming at a standardization of the
techniques described in Section 4 to transfer the information on PID of the
calibration samples to the reference sample of interest for the many physics
analyses. It includes several reweighing approaches, PID resampling and PID
variable transformation. The set of variables identifying the kinematics of
the tracks, the event multiplicity and the PID response can be chosen case by
case, while the access to calibration samples and the implementation of the
algorithms are maintained centrally.

##  6 Data quality, monitoring and validation

As discussed in Section 4, calibration samples are abundant decays selected at
the trigger level with high purity, representative of all the families of
long-lived particles interacting with the LHCb detector, except deuterons.
Their immediate availability during the data taking and their high statistics
are key ingredients for data-quality monitoring and validation. Since the
reconstruction involves different systems of the LHCb detectors depending on
the nature of the particle, the various samples are used to monitor and
validate different aspects of the reconstruction. For example, the recovery of
Bremsstrahlung photons to improve the momentum resolution of the electrons can
only be monitored and validated using an electron sample. Similarly, the
efficiency of muon identification can be better monitored and validated using
a sample of tagged muons.

In order to add redundancy to the validation procedure, a small fraction of
the calibration samples are reconstructed with the offline procedure in real
time, in order to trigger alarms in a timely manner in case of misalignments
between the online and offline reconstruction, due for example to errors in
the database handling the alignment and calibration constants.

Finally, several checks on the reconstructed quantities of the calibration
samples have been included in the automated validation procedure performed
during data taking, aiming at an early identification of deviations from
standard running conditions, and to check possible temporal variations in
performance due to unstable environmental conditions, or ageing of the
detector [41].

Real-time monitoring on pure decay samples representative of the needs of a
wide physics programme will be of critical importance during Run 3 of the LHC,
when, after a major upgrade of the LHCb experiment, most datasets to perform
physics data analyses will be selected in the trigger and stored as decay
candidates, with no support for raw detector data [36]. Since no further
reprocessing of the reconstruction will be possible, any loss in performance
will unavoidably result in a loss of effectiveness for the resulting physics
measurements.

##  7 Conclusions

The strategy to select and process the calibration samples used to measure the
PID performance has seen several improvements to face the challenges set by
Run 2 of the LHC. The samples are now selected directly in real-time at the
highest level of the software trigger, introducing an important benefit in
terms of statistics and absence of selection bias with respect to the offline
selection strategy adopted in Run 1. The calibration samples are used to
measure the PID performance, to correct the simulated samples, and to monitor
the detector performance during the data-taking.

The computing model to manage and process the calibration samples has been
redesigned in order to overcome the scalability challenges set by the larger
statistics needed to investigate the PID performance for the LHC Run 2.

The new scheme has been later adopted to provide the tracking and the photon
reconstruction performance, paving the way for Run 3.

## References

  * [1] LHCb collaboration, R. Aaij et al., _Search for the rare decays ${{B}^{0}_{s}}\rightarrow{\mu^{+}}{\mu^{-}}$ and ${{B}^{0}}\rightarrow{\mu^{+}}{\mu^{-}}$_, Phys. Lett. B699 (2011) 330, arXiv:1103.2465
  * [2] LHCb collaboration, R. Aaij et al., _Search for the rare decays ${{B}^{0}_{s}}\rightarrow{\mu^{+}}{\mu^{-}}$ and ${{B}^{0}}\rightarrow{\mu^{+}}{\mu^{-}}$_, Phys. Lett. B708 (2012) 55, arXiv:1112.1600
  * [3] LHCb collaboration, R. Aaij et al., _Strong constraints on the rare decays ${{B}^{0}_{s}}\rightarrow{\mu^{+}}{\mu^{-}}$ and ${{B}^{0}}\rightarrow{\mu^{+}}{\mu^{-}}$_, Phys. Rev. Lett. 108 (2012) 231801, arXiv:1203.4493
  * [4] LHCb collaboration, R. Aaij et al., _First evidence for the decay ${{B}^{0}_{s}}\rightarrow{\mu^{+}}{\mu^{-}}$_, Phys. Rev. Lett. 110 (2013) 021801, arXiv:1211.2674
  * [5] LHCb collaboration, R. Aaij et al., _Measurement of the ${{B}^{0}_{s}}\rightarrow{\mu^{+}}{\mu^{-}}$ branching fraction and search for ${{B}^{0}}\rightarrow{\mu^{+}}{\mu^{-}}$ decays at the LHCb experiment_, Phys. Rev. Lett. 111 (2013) 101805, arXiv:1307.5024
  * [6] CMS and LHCb collaborations, V. Khachatryan et al., _Observation of the rare ${{B}^{0}_{s}}\rightarrow{\mu^{+}}{\mu^{-}}$ decay from the combined analysis of CMS and LHCb data_, Nature 522 (2015) 68, arXiv:1411.4413
  * [7] LHCb collaboration, R. Aaij et al., _Measurement of the ${{B}^{0}_{s}}\rightarrow{\mu^{+}}{\mu^{-}}$ branching fraction and effective lifetime and search for ${{B}^{0}}\rightarrow{\mu^{+}}{\mu^{-}}$ decays_, Phys. Rev. Lett. 118 (2017) 191801, arXiv:1703.05747
  * [8] LHCb collaboration, R. Aaij et al., _Measurement of the CKM angle $\gamma$ from a combination of ${{B}^{\pm}}\rightarrow{D}h^{\pm}$ analyses_, Phys. Lett. B726 (2013) 151, arXiv:1305.2050
  * [9] LHCb collaboration, R. Aaij et al., _Measurement of ${C\\!P}$ violation and constraints on the CKM angle $\gamma$ in ${{B}^{\pm}}\rightarrow{D}{{K}^{\pm}}$ with ${D}\rightarrow{{K}^{0}_{\mathrm{\scriptscriptstyle S}}}{{\pi}^{+}}{{\pi}^{-}}$ decays_, Nucl. Phys. B888 (2014) 169, arXiv:1407.6211
  * [10] LHCb collaboration, R. Aaij et al., _Measurement of the CKM angle $\gamma$ using ${{B}^{\pm}}\rightarrow{D}{{K}^{\pm}}$ with ${D}\rightarrow{{K}^{0}_{\mathrm{\scriptscriptstyle S}}}{{\pi}^{+}}{{\pi}^{-}}$, ${{K}^{0}_{\mathrm{\scriptscriptstyle S}}}{{K}^{+}}{{K}^{-}}$ decays_, JHEP 10 (2014) 097, arXiv:1408.2748
  * [11] LHCb collaboration, R. Aaij et al., _Measurement of $C\\!P$ observables in $B^{\pm}\rightarrow D^{(\ast)}K^{\pm}$ and $B^{\pm}\rightarrow D^{(\ast)}\pi^{\pm}$ decays_, Phys. Lett. B777 (2017) 16, arXiv:1708.06370
  * [12] LHCb collaboration, A. A. Alves Jr. et al., _The LHCb detector at the LHC_, JINST 3 (2008) S08005
  * [13] LHCb collaboration, R. Aaij et al., _LHCb detector performance_ , Int. J. Mod. Phys. A30 (2015) 1530022, arXiv:1412.6352
  * [14] R. Aaij et al., _Performance of the LHCb Vertex Locator_ , JINST 9 (2014) P09007, arXiv:1405.7808
  * [15] R. Arink et al., _Performance of the LHCb Outer Tracker_ , JINST 9 (2014) P01002, arXiv:1311.3893
  * [16] M. Adinolfi et al., _Performance of the LHCb RICH detector at the LHC_, Eur. Phys. J. C73 (2013) 2431, arXiv:1211.6759
  * [17] A. A. Alves Jr. et al., _Performance of the LHCb muon system_ , JINST 8 (2013) P02022, arXiv:1211.1346
  * [18] R. Aaij et al., _The LHCb trigger and its performance in 2011_, JINST 8 (2013) P04022, arXiv:1211.3055
  * [19] G. Dujany and B. Storaci, _Real-time alignment and calibration of the LHCb Detector in Run II_ , J. Phys. Conf. Ser. 664 (2015) 082010
  * [20] M. Frank, C. Gaspar, B. Jost, and N. Neufeld, _The LHCb Data Acquisition and High Level Trigger Processing Architecture_ , Journal of Physics: Conference Series 664 (2015), no. 8 082011
  * [21] LHCb collaboration, R. Aaij et al., _Measurement of the ${{B}_{c}^{+}}$ meson lifetime using ${{B}_{c}^{+}}\rightarrow{{J\mskip-3.0mu/\mskip-2.0mu\psi\mskip 2.0mu}}{\mu^{+}}{{\nu}_{\mu}}X$ decays_, Eur. Phys. J. C74 (2014) 2839, arXiv:1401.6932
  * [22] LHCb collaboration, R. Aaij et al., _Measurements of the ${{{B}^{+}}}$, ${{B}^{0}}$, ${{B}^{0}_{s}}$ meson and ${{\mathchar 28931\relax}^{0}_{b}}$ baryon lifetimes_, JHEP 04 (2014) 114, arXiv:1402.2554
  * [23] LHCb collaboration, R. Aaij et al., _Search for long-lived heavy charged particles using a ring-imaging Cherenkov technique at LHCb_ , Eur. Phys. J. C75 (2015) 595, arXiv:1506.09173
  * [24] F. Archilli et al., _Performance of the muon identification at LHCb_ , JINST 8 (2013) P10020, arXiv:1306.0249
  * [25] F. Archilli et al., _Performance of the Muon Identification at LHCb_ , JINST 8 (2013) P10020, arXiv:1306.0249
  * [26] A. Hoecker et al., _TMVA: Toolkit for Multivariate Data Analysis_ , PoS ACAT (2007) 040, arXiv:physics/0703039
  * [27] M. De Cian, S. Farry, P. Seyfert, and S. Stahl, _Fast neural-net based fake track rejection_ , Tech. Rep. LHCb-PUB-2017-011. CERN-LHCb-PUB-2017-011, CERN, Geneva, Dec, 2017
  * [28] O. Lupton, L. Anderlini, B. Sciascia, and V. Gligorov, _Calibration samples for particle identification at LHCb in Run 2_ , Tech. Rep. LHCb-PUB-2016-005. CERN-LHCb-PUB-2016-005, CERN, Geneva, Mar, 2016
  * [29] M. Pivk and F. R. Le Diberder, _sPlot: A statistical tool to unfold data distributions_ , Nucl. Instrum. Meth. A555 (2005) 356, arXiv:physics/0402083
  * [30] A. Poluektov, _Kernel density estimation of a multidimensional efficiency profile_ , JINST 10 (2015), no. 02 P02011, arXiv:1411.5528
  * [31] Particle Data Group, C. Patrignani et al., _Review of Particle Physics_ , 2016\.  doi: 10.1088/1674-1137/40/10/100001
  * [32] R. Aaij et al., _The LHCb Trigger and its Performance in 2011_ , JINST 8 (2013) P04022, arXiv:1211.3055
  * [33] LHCb HLT project, J. Albrecht, V. V. Gligorov, G. Raven, and S. Tolk, _Performance of the LHCb High Level Trigger in 2012_ , J. Phys. Conf. Ser. 513 (2014) 012001, arXiv:1310.8544
  * [34] R. Aaij et al., _Tesla : an application for real-time data analysis in High Energy Physics_ , Comput. Phys. Commun. 208 (2016) 35, arXiv:1604.05596
  * [35] R. Aaij et al., _Tesla: an application for real-time data analysis in High Energy Physics_ , Comput. Phys. Commun. 208 (2016) 35, arXiv:1604.05596
  * [36] _LHCb Trigger and Online Upgrade Technical Design Report_ , Tech. Rep. CERN-LHCC-2014-016. LHCB-TDR-016, May, 2014
  * [37] L. Anderlini et al., _Computing strategy for PID calibration samples for LHCb Run 2_ , Tech. Rep. LHCb-PUB-2016-020. CERN-LHCb-PUB-2016-020, CERN, Geneva, Jul, 2016
  * [38] S. Tolk, J. Albrecht, F. Dettori, and A. Pellegrino, _Data driven trigger efficiency determination at LHCb_ , Tech. Rep. LHCb-PUB-2014-039. CERN-LHCb-PUB-2014-039, CERN, Geneva, May, 2014
  * [39] F. Stagni et al., _DIRAC in Large Particle Physics Experiments_ , Journal of Physics: Conference Series 898 (2017), no. 9 092020
  * [40] L. Anderlini et al., _The PIDCalib package_ , Tech. Rep. LHCb-PUB-2016-021. CERN-LHCb-PUB-2016-021, CERN, Geneva, Jul, 2016
  * [41] M. Adinolfi et al., _Lhcb data quality monitoring_ , Journal of Physics: Conference Series 898 (2017), no. 9 092027


Fine-Gray models specify the subdistribution hazards for one out of multiple competing risks to be proportional.
The estimators of parameters and cumulative incidence functions under Fine-Gray models have a simpler structure when data are censoring-complete than when they are more generally incomplete.
This paper considers the case of incomplete data but it exploits the above-mentioned simpler estimator structure for which there exists a wild bootstrap approach for inferential purposes.
The present idea is to link the methodology under censoring-completeness with the more general right-censoring regime with the help of multiple imputation.
In a simulation study, this approach is compared to the estimation procedure proposed in the original paper by Fine and Gray when it is combined with a bootstrap approach.
An application to a data set about hospital-acquired infections illustrates the method.

Keywords: censored data, confidence regions, multiple imputation, inference, resampling, survival analysis

left=3cm, right=2cm,bottom=3.5cm

[Allignol, 2019]
Arthur Allignol.
kmi: Kaplan-Meier Multiple Imputation for the Analysis of
Cumulative Incidence Functions in the Competing Risks Setting, 2019.
URL <https://CRAN.R-project.org/package=kmi>.
Version 0.5.5.

[Allignol and Beyersmann, 2010]
Arthur Allignol and Jan Beyersmann.
Software for fitting nonstandard proportional subdistribution hazards
Biostatistics, 110 (4):0 674–675, 2010.

[Beyersmann et al., 2006]
Jan Beyersmann, Petra Gastmeier, Hajo Grundmann, Sina Bärwolff, Christine
Geffers, Michael Behnke, Henning Rüden, and Martin Schumacher.
Use of multistate models to assess prolongation of intensive care
unit stay due to nosocomial infection.
Infection Control and Hospital Epidemiology, 270
(5):0 493–499, 2006.

[Beyersmann et al., 2009]
Jan Beyersmann, Aurélien Latouche, Anika Buchholz, and Martin Schumacher.
Simulating competing risks data in survival analysis.
Statistics in Medicine, 280 (6):0 956–971,

[Bärwolff et al., 2005]
Sina Bärwolff, Hajo Grundmann, Frank Schwab, Adriana Tami, Michael Behnke,
Christine Geffers, Elke Halle, Ulf B. Göbel, Reinhold Schiller, Daniel
Jonas, Ingo Klare, Klaus Weist, Wolfgang Witte, E. Dinger, Kathrin Beilecke,
Henning Rüden, and Petra Gastmeier.
Incidence of transmission of pathogens in intensive care units.
results of the SIR 3 study.
Der Anaesthesist, 540 (6):0 560–566, 2005.

[Dietrich et al., 2023]
Marina T. Dietrich, Dennis Dobler, and Mathisca C. M. de Gunst.
Wild bootstrap for counting process-based statistics, 2023.

[Dobler and Pauly, 2014]
Dennis Dobler and Markus Pauly.
Bootstrapping Aalen-Johansen processes for competing risks:
Handicaps, solutions, and limitations.
Electronic Journal of Statistics, 80 (2):0
2779–2803, 2014.

[Fine and Gray, 1999]
Jason P. Fine and Robert J. Gray.
A proportional hazards model for the subdistribution of a competing
Journal of the American Statistical Association, 940
(446):0 496–509, 1999.

[Grundmann et al., 2005]
Hajo Grundmann, Sina Bärwolff, Adriana Tami, Michael Behnke, Frank Schwab,
Christine Geffers, Elke Halle, Ulf B. Göbel, Reinhold Schiller, Daniel
Jonas, Ingo Klare, Klaus Weist, Wolfgang Witte, Kathrin Beck-Beilecke, Martin
Schumacher, Henning Rüden, and Petra Gastmeier.
How many infections are caused by patient-to-patient transmission in
intensive care units?
Critical care medicine, 330 (5):0 946–951,

[Kawaguchi, 2019]
Eric Kawaguchi.
Fine-Gray Regression via Forward-Backward Scan, 2019.
URL <https://CRAN.R-project.org/package=fastcmprsk>.
Version 1.1.1.

[Kawaguchi et al., 2020]
Eric S. Kawaguchi, Jenny I. Shen, Gang Li, and Marc A. Suchard.
A fast and scalable implementation method for competing risks data
with the r package fastcmprsk.
The R Journal, 120 (2):0 163–172, 2020.

[Lin, 1997]
Danyu Lin.
Non-parametric inference for cumulative incidence functions in
competing risks studies.
Statistics in Medicine, 160 (8):0 901–910,

[Little and Rubin, 2020]
Roderick J. A. Little and Donald B. Rubin.
Statistical Analysis with Missing Data.
John Wiley & Sons, Inc., 2020.

[Murray, 2018]
Jared S. Murray.
Multiple Imputation: A Review of Practical and Theoretical
Statistical Science, 330 (2):0 142–159, 2018.

[R Core Team, 2016]
R Core Team.
R: A language and environment for statistical
R Foundation for Statistical Computing, Vienna, Austria,
URL <https://www.R-project.org/>.

[Robins and Rotnitzky, 1992]
James M. Robins and Andrea Rotnitzky.
Recovery of Information and Adjustment for Dependent Censoring
Using Surrogate Markers, pages 297–331.
Birkhäuser Boston, 1992.

[Ruan and Gray, 2008]
Ping K. Ruan and Robert J. Gray.
Analyses of cumulative incidence functions via non-parametric
multiple imputation.
Statistics in Medicine, 270 (27):0 5709–5724,

[Rubin, 1978]
Donald B. Rubin.
Multiple imputations in sample surveys—a phenomenological
Bayesian approach to nonresponse.
Proceedings of the Survey Research Methods Section of the
American Statistical Association, pages 20–34, 1978.

[Rubin, 1987]
Donald B. Rubin.
Multiple Imputation for Nonresponse in Surveys.
John Wiley & Sons, Inc., 1987.

[Schomaker and Heumann, 2018]
Michael Schomaker and Christian Heumann.
Bootstrap inference when using multiple imputation.
Statistics in Medicine, 370 (14):0 2252–2266,

[Taylor et al., 2002]
Jeremy M. G. Taylor, Susan Murray, and Chiu-Hsieh Hsu.
Survival estimation and testing via multiple imputation.
Statistics & Probability Letters, 580 (3):0
221–232, 2002.

[Wolkewitz et al., 2008]
Martin Wolkewitz, Ralf-Peter Vonberg, Hajo Grundmann, Jan Beyersmann, Petra
Gastmeier, Sina Bärwolff, Christine Geffers, Michael Behnke, Henning Rüden,
and Martin Schumacher.
Risk factors for the development of nosocomial pneumonia and
mortality on intensive care units: Application of competing risks models.
Critical care (London, England), 120 (2), 2008.