<figure><img src="image_1.png" /></figure>

HAL
open science

# False clustering rate control in mixture models

Ariane Marandon, Tabea Rebafka, Etienne Roquain, Nataliya Sokolovska

► To cite this version:

Ariane Marandon, Tabea Rebafka, Etienne Roquain, Nataliya Sokolovska. False clustering rate control in mixture models. 2022. hal-03625469

HAL Id: hal-03625469

https://hal.science/hal-03625469

Preprint submitted on 30 Mar 2022

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.
# False selection rate control in mixture models

Ariane Marandon, Tabea Rebafka, Etienne Roquain, Nataliya Sokolovska

► To cite this version:

Ariane Marandon, Tabea Rebafka, Etienne Roquain, Nataliya Sokolovska. False selection rate control in mixture models. Scandinavian Journal of Statistics, 2025, 10.1111/sjos.70017. hal-03625469v2

HAL Id: hal-03625469

https://hal.science/hal-03625469v2

Submitted on 16 Oct 2025

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# False clustering rate control in mixture models

Ariane Marandon¹, Tabea Rebafka¹, Etienne Roquain¹, Nataliya Sokolovska²

¹ Sorbonne Université, LPSM

² Sorbonne Université, NutriOmics

February 28, 2022

## Abstract

The clustering task consists in delivering labels to the members of a sample. For most data sets, some individuals are ambiguous and intrinsically difficult to attribute to one or another cluster. However, in practical applications, misclassifying individuals is potentially disastrous. To overcome this difficulty, the idea followed here is to classify only a part of the sample in order to obtain a small misclassification rate. This approach is well known in the supervised setting, and referred to as classification with an abstention option. The purpose of this paper is to revisit this approach in an unsupervised mixture-model framework. The problem is formalized in terms of controlling the false clustering rate (FCR) below a prescribed level $\alpha$, while maximizing the number of classified items. New procedures are introduced and their behavior is shown to be close to the optimal one by establishing theoretical results and conducting numerical experiments. An application to breast cancer data illustrates the benefits of the new approach from a practical viewpoint.

# 1 Introduction

## 1.1 Background

Clustering is a standard statistical task that aims at grouping together individuals with sim- ilar features. However, it is common that data sets include ambiguous individuals that are inherently difficult to classify, which makes the clustering result potentially unreliable. To il- lustrate this point, consider a Gaussian mixture model with overlapping mixture components. Then it is difficult, or even impossible, to assign the correct cluster label to data points that fall in the overlap of those clusters, see Figure 1. Hence, when the overlap is large (Figure 1 panel (b)), the misclassification rate of a standard clustering method is inevitably elevated.

This issue is critical in applications where misclassifications come with a high cost for the user and should be avoided. This is for example the case for medical diagnosis, where an error can have severe consequences on the individual's health. When there is too much uncertainty, a solution is to avoid classification for such individuals, and to adopt a wiser "abstention decision", that leaves the door open for further medical exams.

In a supervised setting, classification with a reject (or abstention) option is a long-standing statistical paradigm, that can be traced back to Chow (1970), with more recent works includ- ing Herbei and Wegkamp (2006); Bartlett and Wegkamp (2008); Wegkamp and Yuan (2011), among others. In this line of research, rejection is accounted for by adding a term to the risk that penalizes any rejection (i.e., non classification). However, this approach does not provide a prescribed control of the classification error among the classified items (those that are not rejected).

More confidence in the labelling can be brought via conformal prediction, Vovk et al. (1999, 2005), see also Angelopoulos et al. (2021) for a review. This generic technique can be applied to supervised classification for producing confidence sets for the label vector. However, such inference does not incorporate the abstention rule. Very recently, Denis and Hebiri (2020)
<figure><img src="image_2.png" /></figure>

ORIGINAL ARTICLE

Scandinavian Journal of Statistics

# False selection rate control in mixture models

Ariane Marandon<fnref n="1" /> ID |

Tabea Rebafka<fnref n="2" /> |

Etienne Roquain<fnref n="2" /> ID |

Nataliya Sokolovska<fnref n="3" />

**Correspondence**

Ariane Marandon, The Alan Turing Institute, London, UK.

Email: amarandon-carlhian@turing.ac.uk

**Funding information**

DIM MATH INNOV; Turing-Roche partnership

**Abstract**

The clustering task consists in partitioning elements of a sample into homogeneous groups. Most datasets contain individuals that are ambiguous and intrinsically difficult to attribute to one or another cluster. However, in practical applications, misclassifying individuals is potentially disastrous and should be avoided. To keep the misclassification rate small, one can decide to classify only a part of the sample. In the supervised setting, this approach is well known and referred to as classification with an abstention option. In this paper, the approach is revisited in an unsupervised mixture-model framework. The purpose is to develop a method that guarantees the false selection rate (FSR) does not exceed a predefined level $\alpha$. We propose a plug-in procedure and provide a theoretical analysis, quantifying the deviation of the FSR from the target $\alpha$ with explicit remainder terms. Bootstrap versions of the procedure are shown to improve the performance in numerical experiments.

**KEYWORDS**

abstention option, bootstrap, clustering, false discovery rate, mixture models

<footnote n="1">¹The Alan Turing Institute, London, UK</footnote>

<footnote n="2">²LPSM, Sorbonne Université, Université Paris Cité & CNRS, Paris, France</footnote>

<footnote n="3">³LCQB, Sorbonne Université, CNRS, Paris, France</footnote>

<footnote>This is an open access article under the terms of the Creative Commons Attribution License, which permits use, distribution and reproduction in any medium, provided the original work is properly cited.</footnote>

<footnote>© 2025 The Author(s). *Scandinavian Journal of Statistics* published by John Wiley & Sons Ltd on behalf of The Board of the Foundation of the Scandinavian Journal of Statistics.</footnote>
proposed a conformal approach with a rejection option. The method consists in minimizing
the probability of misclassification conditionally on the fact that it has been classified, while
controlling the probability of being classified. This brings a safety guarantee on the label
assigned to the classified items.

## 1.2 Aim and approach

The goal of the present work is to propose a labelling guarantee on the classified items in the vein of Denis and Hebiri (2020), but in the more challenging unsupervised setting, where no training set is available and data are assumed to be generated from a finite mixture model. This is achieved by the possibility to refuse to cluster ambiguous individuals and by using the false clustering rate (FCR), which is defined as the average proportion of misclassifications among the classified objects. Our procedures are devised to keep the FCR below some nominal level $\alpha$, while classifying a maximum number of items.

It is important to understand the role of the nominal level $\alpha$ in our approach. It is chosen by the user and depends on their acceptance or tolerance for misclassified objects. Since the FCR is the misclassification risk that is allowed on the classified objects, the final interpretation of FCR control at level $\alpha$ is clear: if for instance $\alpha$ is chosen to be 5% and 100 items are finally chosen to be classified by the method, then the number of misclassified items is expected to be at most 5. This high interpretability is similar to the one of the false discovery rate (FDR) in multiple testing, which has known a great success in applications since its introduction by Benjamini and Hochberg (1995). This is a clear advantage of our approach for practical use compared to the methods with a rejection option that are based on a penalized risk.

In our framework, a procedure is composed of two intertwined decisions:

* a selection rule deciding which of the items are labelled;

* a clustering method inferring the labels for the selected items.

Both decisions heavily rely on the appropriate quantification of the uncertainty of the cluster
labels. As such, our approach is model-based, and can be viewed as a method that thresholds
the posterior probabilities of the cluster labels with a data-driven choice of the threshold.
The performance of the method will depend on the quality of the estimates of these posterior
probabilities in the mixture model.

The adaptive character of our method is illustrated in Figure 1: when the clusters are well separated (panel (a)), the new procedure only discards few items and provides a clustering close to the correct one. However, when the clusters are overlapping (panel (b)), to avoid a high misclassification error, the procedure discards most of the items and only provides few labels, for which the uncertainty is low. In both cases, the proportion of misclassified items among the selected ones is small and in particular close to the target level $\alpha$ (here 10%). Hence, by adapting the amount of labeled or discarded items, our method always delivers a reliable clustering result, inspite of the varying intrinsic difficulty of the clustering task.

## 1.3 Presentation of the results

Let us now describe in more details the main contributions of the paper.

* We introduce three new data-driven procedures that perform simultaneously selection and clustering: the plug-in procedure (illustrated in Figure 1) and two bootstrap procedures (parametric and non-parametric), see Section 3.2.
# 1 | INTRODUCTION

## 1.1 | Background

Clustering is the standard statistical task that aims to group together items with similar features. However, it is common that datasets include ambiguous items that are inherently difficult to classify, which makes the clustering result potentially unreliable. To illustrate this point, consider a Gaussian mixture model with overlapping mixture components, see Figure 1. In this setting, it is impossible to construct a clustering method that consistently assigns the correct cluster label to data points that fall in the overlap of those clusters. That is, when the overlap is large, the misclassification rate (the proportion of elements which are erroneously clustered up to label switching, see (1) below for a formal definition) of a standard clustering method is inevitably elevated. This issue can be critical in applications where falsely clustered items induce high costs for the user or make the interpretation of the clusters difficult. In particular in an early-stage analysis of a new dataset, it might be much easier to interpret the clusters when they do not contain spurious observations.

To be more precise about the motivation of our unsupervised approach, let us consider the common situation where a medical doctor considers a clustering of a group of patients based on phenotyping, i.e. clusters are homogeneous groups of individuals in terms of their phenotypic profiles. The doctor may analyze each cluster to decide on a treatment for all patients in the same cluster. For instance, in Eckardt et al. (2023), a clustering method is applied to stratify leukemia patients into groups according to risk; these groups differ in genetics and clinical values, which leads to differences in patient survival and treatment. Another example concerns patients suffering from injuries who are to be clustered according to their physical and psychological profiles to identify optimal recovery for each group (Stoitsas et al., 2022). Our approach is motivated by the preliminary unsupervised stage, for which it is important to have patient clusters that we can clearly identify. Hence, patients with an ambiguous profile should not be considered at this stage and put aside, as taking a decision for these patients might require additional information. This corresponds to following an approach with abstention decision. In a supervised setting, classification with a reject (or abstention) option is a long-standing statistical paradigm, that can be traced back to Chow (1970), with more recent works including Herbei and Wegkamp (2006), Bartlett and

<figure><img src="image_3.png" /><figcaption>**FIGURE 1** Data from Gaussian mixtures with three components (n = 200), in a fairly separated case (panel (a)) and an ambiguous case (panel (b)). In each panel, the left part displays the true clustering (colors/point characters correspond to the different true clusters), while the right part illustrates the new procedure (plug-in procedure at level α = 10%), that does not cluster all items (colors/point characters correspond to the different estimated clusters). The points not classified are depicted by grey crosses. Red circles indicate erroneous labels.</figcaption></figure>
<figure><img src="image_2.png" /><figcaption>Figure 1: Data from Gaussian mixtures with three components ($n = 200$), in a fairly separated case (panel (a)) and an ambiguous case (panel (b)). In each panel, the left part displays the true clustering, while the right part illustrates the new procedure (plug-in procedure at level $\alpha = 10\%$), that does not cluster all items. The points not classified are depicted by grey crosses. Points circled in red correspond to erroneous labels.</figcaption></figure>

* We prove that the plug-in procedure controls the FCR at the desired level $\alpha$, up to a reminder term which becomes small when the sample size grows. In addition, this procedure is shown to satisfy the following optimality property: any other procedure that provides an FCR control necessarily classifies equal or less items than the plug-in procedure, up to a small remainder term. All results are non-asymptotic in the sample size (Theorem 6), and it is shown that the reminder terms tend to zero when the sample size tends to infinity (Corollary 7). We also put forward specific convergence rates depending on the parameter estimation accuracy (Corollary 8).

* Numerical experiments establish that the bootstrap procedures improve the plug-in procedure in terms of non-asymptotic FCR control, and thus are more reliable for practical use, where the sample size may be moderate, see Section 5.1.

* Our analysis also shows that the fixed threshold procedure that only labels items with a maximum posterior probability larger than $1 - \alpha$ is generally suboptimal for an FCR control at level $\alpha$, see Section 5.1. To this extent, our procedures can be seen as more refined algorithms that classify more individuals while maintaining the FCR control.

* The practical impact of our approach is demonstrated on a real data set, see Section 5.2.

## 1.4 Relation to previous work

**Other clustering guarantees in unsupervised learning** While we provide a specific FCR control guarantee on the clustering, other criteria, not particularly linked to a rejection option, have been previously proposed in an unsupervised setting. Considering the cluster labels as fixed parameters, minimizing the misclassification risk has been proposed in Lei and Rinaldo (2015); Lu and Zhou (2016), among others. When clusters are rather considered as latent variables generated concomitantly with the observations (as we do here), other criteria have been investigated: the probability to make a different decision than the Bayes rule (Azizyan et al., 2013), the probability of exact recovery (Arora and Kannan, 2005; Bickel and Chen, 2009), or the fact that all clusters are mostly homogeneous with high probability (Najafi et al., 2020). All these works provide a guarantee only if the setting is favorable enough. By contrast, providing a rejection option is the key to obtain a guarantee in any setting (in the worst situation, the procedure will not classify any item).
Wegkamp (2008), Wegkamp and Yuan (2011), Cortes et al. (Cortes et al., 2016), Ni et al. (2019), Cao et al. (2022), among others. In this line of research, rejection is accounted for by adding a term to the risk that penalizes any rejection (i.e., nonclassification).

Recently, still in the supervised setting, Geifman and El-Yaniv (2017) and Angelopoulos et al. (2025) have considered the problem of having a prescribed control of the classification error conditional on the item being classified (not rejected). In these works, the proposed method consists in thresholding the class probabilities estimated by a pretrained classifier, in a data-driven manner. The computation of the risk-controlling threshold is based on concentration inequalities that bound the true risk with respect to the empirical one. This criterion is also studied in Denis and Hebiri (2020) in a semisupervised setting. However, the control problem is different, because it rather aims at a type-II error control concerning the probability of classifying an item. Finally, a related line of research is conformal prediction (Angelopoulos & Bates, 2021; Lei, 2014; Romano et al., 2020; Sadinle et al., 2019; Vovk et al., 2005). This generic technique builds prediction sets that contain the true label with probability above a fixed confidence level. In that framework, one trades single-valued classification for set-valued predictions to have a prespecified prediction coverage, and finite-sample control is achieved using hold-out labeled data (calibration sample).

The supervised approaches can only be applied when classes are known and already have a specific interpretation, while our unsupervised method is an exploratory approach, where the goal is to identify unknown clusters of individuals that we wish to interpret for a better understanding of the dataset.

## 1.2 | Aim and approach

The goal of the present work is to propose a labeling guarantee on the classified items in the more challenging unsupervised setting, where no labeled training set is available (Eckardt et al., 2023; Stoitsas et al., 2022) and data are assumed to be generated from a finite mixture model (Grün, 2019; Ng, 2019). This is achieved by the possibility to refuse to cluster ambiguous individuals and by using the false selection rate (FSR), which is defined as the average proportion of misclassifications among the classified objects. Our procedures are devised to keep the FSR below some nominal level $\alpha$, while classifying a maximum number of items.

It is important to understand the role of the nominal level $\alpha$ in our approach. It is chosen by the user and depends on their acceptance or tolerance for misclassified objects. Since the FSR is the misclassification risk that is allowed on the classified objects, the final interpretation of an FSR control at level $\alpha$ is clear: if, for instance, $\alpha$ is set to 5% and 100 items are finally chosen to be classified by the method, then the number of misclassified items is expected to be at most 5. This high interpretability is similar to the one of the false discovery rate (FDR) in multiple testing, which has known a great success in applications since its introduction by Benjamini and Hochberg (1995). This is a clear advantage of our approach for practical use compared to the methods with a rejection option that are based on a penalized risk.

In our framework, a procedure is composed of two components:

* a clustering method inferring the labels;

* a selection rule deciding which labels to keep.

Importantly, the selection rule is only applied *after* a clustering method is fitted on the (entire) sample. In other words, the procedure consists of two subsequent steps: a clustering step, that
**Comparison to Denis and Hebiri (2020) and Mary-Huard et al. (2021)** We describe here two recent studies that are related to ours, because they also use a FCR-like criterion. The first one is the work of Denis and Hebiri (2020), which also relies on a thresholding of the (estimated) posterior probabilities. However, the control is different, because it does not provide an FCR control, but rather a type-II error control concerning the probability of classifying an item. Also, the proposed procedure therein uses an additional labeled sample (semi-supervised setting), which is forbidden in our case.

The work of Mary-Huard et al. (2021) also proposes a control of the FCR. However, the analysis therein is solely based on the case where the model parameters are known (thus corresponding to the oracle case developed in Section 3.1 here). Compared to Mary-Huard et al. (2021), the present work provides number of new contributions, which are all given in Section 1.3. Let us also emphasize that we handle the label switching problem in the FCR, which seems to be overlooked in Mary-Huard et al. (2021).

**Relation to false discovery rate** The FCR is closely related to the false discovery rate (FDR), defined as the average proportion of errors among the discoveries, in a multiple testing context. In fact, we can roughly view the problem of designing an abstention rule as testing, for each item $i$, whether the clustering rule correctly classifies item $i$ or not. With this analogy, our selection rule is based on quantities similar to the local FDR values (Efron et al., 2001), a key quantity to build optimal FDR controlling procedure in multiple testing mixture models, see, e.g., Storey (2003); Sun and Cai (2007); Cai et al. (2019); Rebafka et al. (2019). In particular, our final selection procedure shares similarities with the procedure introduced in Sun and Cai (2007), also named cumulative $\ell$-value procedure Abraham et al. (2021). In addition, our theoretical analysis is related to the work of Rebafka et al. (2019), although the nature of the algorithm developed therein is different from here: they use the $q$-value procedure of Storey (2003), while our method rather relies on the cumulative $\ell$-value procedure.

## 1.5 Organization of the paper

The paper is organized as follows: Section 2 introduces the model and relevant notation, namely the FCR criterion, with a particular care of the label switching problem. Section 3 presents the methods: the oracle, plug-in and the bootstrap approaches. Our main theoretical results are provided in Section 4, after introducing the appropriate assumptions. Section 5 presents numerical experiments and an application to a real data set, while a conclusion is given in Section 6. Proofs of the results and other technical details are deferred to appendices.

# 2 Setting

This section presents the notation, model, procedures and criteria that will be used throughout the manuscript.

## 2.1 Model

Let $\mathbf{X} = (X_1, \dots, X_n)$ be an observed random sample of size $n$. Each $X_i$ is an i.i.d. copy of a $d$-dimensional real random vector, which is assumed to follow the standard mixture model:

$$ Z \sim \mathcal{M}(\pi_1, \dots, \pi_Q), \\ X|Z = q \sim F_{\phi_q}, \quad 1 \le q \le Q, $$

where $\mathcal{M}(\pi_1, \dots, \pi_Q)$ denotes the multinomial distribution of parameter $\pi$ (equivalently, $\pi_q = P(Z = q)$ for each $q$). The model parameters are given by
assigns a cluster label to every item and a selection step, that chooses which items to classify—in which case, the label from the previous clustering step is assigned. For the items that are not selected, we discard the cluster label, that is, we effectively abstain to make a classification decision for those items. In particular, we emphasize that the clustering method is not fitted again after selection (which would in general lead to bias).

The quality of the selection heavily relies on the appropriate quantification of the uncertainty of the cluster labels. For this, our approach is model-based and can be viewed as a method that thresholds the posterior probabilities of the cluster labels with a data-driven choice of the threshold. The performance of the method will depend on the quality of the estimates of these posterior probabilities in the mixture model.

The adaptive character of our method is illustrated in Figure 1: when the clusters are well separated (panel (a)), the new procedure only discards few items and provides a clustering close to the correct one. However, when clusters are overlapping (panel (b)), to avoid a high misclassification error, the procedure discards most of the items and only provides few labels, for which the uncertainty is low. In both cases, the proportion of misclassified items among the selected ones is small and in particular close to the target level $\alpha$ (here 10%). Hence, by adapting the amount of labeled or discarded items, our method always delivers a reliable clustering result, inspite of the varying intrinsic difficulty of the clustering task.

## 1.3 | Presentation of the results

Let us now describe in more details the main contributions of the paper.

* We introduce three new data-driven procedures that perform simultaneously selection and clustering: the plug-in procedure (illustrated in Figure 1) and two bootstrap procedures (parametric and nonparametric), see Section 3.2.

* We provide a theoretical analysis of the plug-in procedure, quantifying the FSR deviation with respect to the target level $\alpha$ with explicit remainder terms, which become small when the sample size grows. In addition, this procedure is shown to satisfy the following optimality property: any other procedure that provides an FSR control necessarily classifies as many or less items than the plug-in procedure, up to a small remainder term (Theorem 2).

* Numerical experiments establish that the bootstrap procedures improve the plug-in procedure, and thus are more reliable for practical use, where the sample size may be moderate, see Section 5.1. In particular, the FSR control is shown to be valid in scenarios with various levels of difficulty.

* Our analysis also shows that a fixed threshold procedure that only clusters items with a maximum posterior probability larger than $1 - \alpha$ typically clusters less items than our procedures, see Section 5.1. To this extent, our procedures can be seen as refined algorithms that classify more individuals while maintaining the FSR control.

* The practical impact of our approach is demonstrated on a real dataset, see Section 5.2.

## 1.4 | Relation to previous work

### 1.4.1 | Misclassification criterion

In the unsupervised setting, defining a notion of error for a clustering procedure is not straightforward and there exists multiple criteria (Fahad et al., 2014). We adopt the proportion of
* the probability distribution $\pi$ on $\{1, \dots, Q\}$ that is assumed to satisfy $\pi_q > 0$ for all $q$. Hence, $\pi_q$ corresponds to the probability of being in class $q$;

* the parameter $\phi = (\phi_1, \dots, \phi_Q) \in \mathcal{U}^Q$, where $\{F_u, u \in \mathcal{U}\}$ is a collection of distributions on $\mathbb{R}^d$. We moreover assume that each distribution $F_u$ has a density with respect to the Lebesgue measure on $\mathbb{R}^d$, that we denote by $f_u$. Also, we assume that the $\phi_q$'s are all distinct.

The parameter $Q$, coding for the number of classes, is assumed to be known and fixed throughout the manuscript (see Section 6 for a discussion). The overall parameter is thus $\theta = (\pi, \phi)$, the parameter set is denoted by $\Theta$, and the distribution of $(Z, X)$ is denoted by $P_\theta$. The distribution family $\{P_\theta, \theta \in \Theta\}$ is the considered statistical model. We also assume that $\Theta$ is an open subset of $\mathbb{R}^K$ for some $K \ge 1$ with the corresponding topology.

In this mixture model, the latent vector $\mathbf{Z} = (Z_1, \dots, Z_n)$ encodes a partition of the $n$ observations into $Q$ classes given by $\{1 \le i \le n : Z_i = q\}$, $1 \le q \le Q$. We refer to this model-based, random partition as the true latent clustering in the sequel.

In what follows, the “true” parameter that generates $(Z, X)$ is assumed to be fixed and is denoted by $\theta^* \in \Theta$.

## 2.2 Procedure and criteria

Our approach starts with a given clustering rule, that aims at recovering the true latent clustering. In general, a clustering rule is defined as a (measurable) function of the observation $\mathbf{X}$ returning a vector label $\hat{\mathbf{Z}} = (\hat{Z}_i)_{1 \le i \le n} \in \{1, \dots, Q\}^n$ for which the label $q$ is given to individual $i$ if and only if $\hat{Z}_i = q$.

The classification error of $\hat{Z}$ is thus given by $\varepsilon(\hat{\mathbf{Z}}, \mathbf{Z}) = \sum_{i=1}^n \mathbf{1}\{Z_i \ne \hat{Z}_i\}$. Obviously, since switching the labels in $\mathbf{Z}$ does not change the true latent partition, we should allow for switching the labels when making the decision $\hat{\mathbf{Z}}$ from the data. The clustering risk of $\hat{\mathbf{Z}}$ is thus defined by

$$R(\hat{\mathbf{Z}}) = \mathbb{E}_{\theta^*} \left( \min_{\sigma \in [Q]} \mathbb{E}_{\theta^*} \left( n^{-1} \varepsilon(\sigma(\hat{\mathbf{Z}}), \mathbf{Z}) \mid \mathbf{X} \right) \right), \quad (1)$$

where $[Q]$ denotes the set of all permutations on $\{1, \dots, Q\}$. Above, the minimum over $\sigma$ is added to handle the aforementioned label switching problem.

**Remark 1.** Note that the position of the minimum w.r.t. $\sigma$ matters in the risk (1): $\sigma$ is allowed to depend on $X$ but not on $Z$. Hence, this risk has to be understood as being computed up to a data dependent label switching. This definition coincides with the usual definition of misclassification risk in the situation where the true clustering is deterministic, see Lei and Rinaldo (2015); Lu and Zhou (2016). Hence, it can be seen as a natural extension of the latter to a mixture model where the true clustering is random.

Classically, we aim to find a clustering rule $\hat{\mathbf{Z}}$ such that the clustering risk is “small”. However, as mentioned above, whether this is possible or not depends on the intrinsic difficulty of the clustering problem and thus of the true parameter $\theta^*$ (see Figure 1). Therefore, the idea is to provide a selection rule, that is, a (measurable) function of the observation $\mathbf{X}$ returning a subset $S \subset \{1, \dots, n\}$, such that the clustering risk with restriction to $S$ is small. Throughout the paper, a *procedure* refers to a couple $C = (\hat{\mathbf{Z}}, S)$, where $\hat{\mathbf{Z}}$ is a clustering rule and $S$ is a selection rule.
misclassifications with respect to predefined class labels (computed up to a permutation of the labels). This criterion takes a pointwise point of view where each data point is assigned a particular label independently of other data points, which can be motivated from a mixture model setup. Alternatively, many popular criteria are based on pairwise comparisons, assessing whether pairs of similar data points are placed in the same cluster and dissimilar data points in different clusters. In the binary case, Bao et al. (2020) established that the misclassification error is actually a monotonic function of the pairwise classification error, highlighting the strong relationship between these criteria.

## 1.4.2 | Clustering guarantees in unsupervised learning

While we provide a specific FSR control guarantee on the clustering, other criteria, not linked to a rejection option, have been previously proposed in an unsupervised setting. Previous works provided essentially two types of guarantees: while early works focused on the probability of exact recovery (Abbe, 2018; Arora & Kannan, 2005; Vempala & Wang, 2004), recent contributions rather considered minimizing the misclassification risk (Chretien et al., 2019; Giraud & Verzelen, 2018; Lei & Rinaldo, 2015; Lu & Zhou, 2016). Other criteria include the probability to make a different decision than the Bayes rule (Azizyan et al., 2013), or the fact that all clusters are mostly homogeneous with high probability (Najafi et al., 2020). In these works, the error rates can be small only in regions of the parameter space that make the clusters separated enough. By contrast, FSR control can be provided even when the clusters overlap by selecting only items that are not ambiguous.

## 1.4.3 | Comparison to Mary-Huard et al. (2021)

The recent work of Mary-Huard et al. (2021) also proposes a control of the FSR. However, the analysis therein is solely based on the case where the model parameters are known (thus corresponding to the oracle case developed in Section 3.1 here). Compared to Mary-Huard et al. (2021), the present work provides a number of new contributions, which are all outlined in Section 1.3. Let us also emphasize that we handle the label switching problem in the FSR, which seems to be overlooked in Mary-Huard et al. (2021).

## 1.4.4 | Relation to the false discovery rate

The FSR is closely related to the false discovery rate (FDR) in multiple testing, defined as the average proportion of errors among the discoveries. In fact, we can roughly view the problem of designing an abstention rule as testing, for each item i, whether the clustering rule correctly classifies item i or not. With this analogy, our selection rule is based on quantities similar to the local FDR values (Efron et al., 2001), a key quantity to build optimal FDR controlling procedures in multiple testing mixture models, see, e.g., Storey (2003), Sun and Cai (2007), Cai et al. (2019), Rebafka et al. (2022). In particular, our final selection procedure shares similarities with the procedure introduced in Sun and Cai (2007), also named cumulative *l*-value procedure (Abraham et al., 2022). In addition, our theoretical analysis is related to the work of Rebafka et al. (2022), although the nature of the algorithm developed therein is different from here: they
**Definition 1** (False clustering rate). The false clustering rate (FCR) of a procedure $\mathcal{C} = (\hat{\mathbf{Z}}, S)$ is given by

$$ \mathrm{FCR}_{\theta^*}(\mathcal{C}) = \mathbb{E}_{\theta^*} \left( \min_{\sigma \in [Q]} \mathbb{E}_{\theta^*} \left( \frac{\varepsilon_S(\sigma(\hat{\mathbf{Z}}), \mathbf{Z})}{\max(|S|, 1)} \middle| \mathbf{X} \right) \right), \quad (2) $$

where $\varepsilon_S(\hat{\mathbf{Z}}, \mathbf{Z}) = \sum_{i \in S} \mathbb{1}\{Z_i \neq \hat{Z}_i\}$ denotes the misclassification error restricted to the subset $S$.

In this work, the aim is to find a procedure $\mathcal{C}$ such that the false clustering rate is controlled at a nominal level $\alpha$, that is, $\mathrm{FCR}_{\theta^*}(\mathcal{C}) \le \alpha$. Obviously, choosing $S$ empty implies $\varepsilon_S(\sigma(\hat{\mathbf{Z}}), \mathbf{Z}) = 0$ a.s. for any permutation $\sigma$ and thus satisfies this control. Hence, while maintaining the control $\mathrm{FCR}_{\theta^*}(\mathcal{C}) \le \alpha$, we aim to classify as much individuals as possible, that is, to make $\mathbb{E}_{\theta^*}|S|$ as large as possible.

The FCR (2) involves an expectation of a ratio, which is more difficult to handle that a ratio of expectations. Hence, the following simpler alternative criterion will also be useful in our analysis.

**Definition 2** (Marginal false clustering rate). The marginal false clustering rate (mFCR) of a procedure $\mathcal{C} = (\hat{\mathbf{Z}}, S)$ is given by

$$ \mathrm{mFCR}_{\theta^*}(\mathcal{C}) = \frac{\mathbb{E}_{\theta^*} \left( \min_{\sigma \in [Q]} \mathbb{E}_{\theta^*} \left( \varepsilon_S(\sigma(\hat{\mathbf{Z}}), \mathbf{Z}) \middle| \mathbf{X} \right) \right)}{\mathbb{E}_{\theta^*}(|S|)}, \quad (3) $$

with the convention $0/0 = 0$.

Note that the mFCR is similar to the criterion introduced in Denis and Hebiri (2020) in the supervised setting.

## 2.3 Notation

We will extensively use the following notation: for all $q \in \{1, \dots, Q\}$ and $\theta = (\pi, \phi) \in \Theta$, we let

$$ \ell_q(X, \theta) = \mathbb{P}_{\theta}(Z = q | X) = \frac{\pi_q f_{\phi_q}(X)}{\sum_{\ell=1}^{Q} \pi_{\ell} f_{\phi_{\ell}}(X)}; \quad (4) $$

$$ T(X, \theta) = 1 - \max_{q \in \{1, \dots, Q\}} \ell_q(X, \theta) \in [0, 1 - 1/Q]. \quad (5) $$

Hence, $\ell_q(X, \theta)$ is the posterior probability of being in class $q$ given the measurement $X$ under the distribution $P_{\theta}$. The quantity $T(X, \theta)$ is a measure of the risk when classifying $X$: it is close to 0 when there exists a class $q$ such that $\ell_q(X, \theta)$ is close to 1, that is, when $X$ can be classified with a large confidence.

# 3 Methods

In this section, we introduce new methods for controlling the FCR. We start by identifying an *oracle* method, in the sense that it uses the true value of the parameter $\theta^*$. Substituting in that oracle the unknown parameter $\theta^*$ by an estimator provides our first method, called the *plug-in* procedure. We then define a refined version of the plug-in procedure, that accounts for the variability of the estimator and is based on a *bootstrap* approach.
use the q-value procedure of Storey (2003), while our method rather relies on the cumulative ℓ-value procedure.

## 1.5 | Organization of the paper

The paper is organized as follows: Section 2 introduces the model and relevant notation, namely the FSR criterion, with a particular care for the label switching problem. Section 3 presents the new methods: the oracle, the plug-in and the bootstrap approaches. Our main theoretical results are provided in Section 4, after introducing appropriate assumptions. Section 5 presents numerical experiments and an application to a real dataset, while a conclusion is given in Section 6. Proofs of the results and technical details are deferred to appendices.

## 2 | SETTING

This section presents the notation, model, procedures, and criteria that will be used throughout the manuscript.

### 2.1 | Model

Let $\mathbf{X} = (X_1, \dots, X_n)$ be an observed random sample of size $n$. Each $X_i$ is an i.i.d. copy of a $d$-dimensional real random vector, which is assumed to follow the standard mixture model:

$$ Z \sim \mathcal{M}(\pi_1, \dots, \pi_Q), \\ X|Z = q \sim F_{\phi_q}, \quad 1 \le q \le Q, $$

where $\mathcal{M}(\pi_1, \dots, \pi_Q)$ denotes the multinomial distribution of parameter $\pi$ (equivalently, $\pi_q = \mathbb{P}(Z = q)$ for each $q$). The model parameters are given by

* the probability distribution $\pi$ on $\{1, \dots, Q\}$ that is assumed to satisfy $\pi_q > 0$ for all $q$. Hence, $\pi_q$ corresponds to the probability of being in class $q$;

* the parameter $\phi = (\phi_1, \dots, \phi_Q) \in \mathcal{V}^Q$, where $\{F_u, u \in \mathcal{V}\}$ is a collection of distributions on $\mathbb{R}^d$. Every distribution $F_u$ is assumed to have a density with respect to the Lebesgue measure on $\mathbb{R}^d$, denoted by $f_u$. Moreover, we assume that the $\phi_q$'s are all distinct (and that the $F_u, u \in \mathcal{V}$ are all distinct).

The number $Q \ge 2$ of classes is assumed to be known and fixed throughout the manuscript (see Section 6 for a discussion). Thus, the overall parameter is $\theta = (\pi, \phi)$, the parameter set is denoted by $\Theta$, and the distribution of $(Z, X)$ is denoted by $P_\theta$. The distribution family $\{P_\theta, \theta \in \Theta\}$ is the considered statistical model. We also assume that $\Theta$ is an open subset of $\mathbb{R}^K$ for some $K \ge 1$. In this mixture model, the latent vector $\mathbf{Z} = (Z_1, \dots, Z_n)$ is called the true latent clustering in the sequel. In what follows, the “true” parameter that generates $(Z, X)$ is assumed to be fixed and is denoted by $\theta^* \in \Theta$. Furthermore, we use the common notation $\mathbb{P}_{\theta^*}(A)$ or $\mathbb{P}_{\mathbf{X} \sim P_{\theta^*}}(A)$ to denote the probability of an event $A$ involving the sample $\mathbf{X} = (X_1, \dots, X_n)$ (and possibly also the
## 3.1 Oracle procedures

**MAP clustering** Here, we proceed as if an oracle had given us the true value of $\theta^*$ and we introduce an oracle procedure $C^* = (\hat{\mathbf{Z}}^*, S^*)$ based on this value. As the following lemma shows, the best clustering rule is well-known and given by the Bayes clustering $\hat{\mathbf{Z}}^* = (\hat{Z}_1^*, \dots, \hat{Z}_n^*)$ which can be written as

$$ \hat{Z}_i^* \in \underset{q \in \{1, \dots, Q\}}{\operatorname{argmax}} \ell_q(X_i, \theta^*), \quad i \in \{1, \dots, n\}, \qquad (6) $$

where $\ell_q(\cdot)$ is the posterior probability given by (4).

**Lemma 1.** We have $\min_{\hat{\mathbf{Z}}} R(\hat{\mathbf{Z}}) = R(\hat{\mathbf{Z}}^*) = n^{-1} \sum_{i=1}^n \mathbb{E}_{\theta^*}(T_i^*)$, for the Bayes clustering $\hat{\mathbf{Z}}^*$ defined by (6) and for

$$ T_i^* = T(X_i, \theta^*) = \mathbb{P}_{\theta^*}(Z_i \neq \hat{Z}_i^* | X_i), \quad i \in \{1, \dots, n\}, \qquad (7) $$

where $T(\cdot)$ is given by (5).

In words, Lemma 1 establishes that the oracle statistics $T_i^*$ corresponds to the posterior misclassification probabilities of the Bayes clustering. To decrease the overall misclassification risk, it is natural to avoid classification for points with a high value of the test statistic $T_i^*$.

**Thresholding selection rules** In this section, we introduce the selection rule, which assesses which items have to be classified. From the above paragraph, it is natural to consider a thresholding based selection rule $S = \{i \in \{1, \dots, n\} : T_i^* \le t\}$, for some threshold $t$ to be suitably chosen. The following result gives insights for choosing such a threshold $t$.

**Lemma 2.** For a procedure $C = (\hat{\mathbf{Z}}^*, S)$ with Bayes clustering and an arbitrary selection $S$,

$$ \text{FCR}_{\theta^*}(C) = \mathbb{E}_{\theta^*} \left( \frac{\sum_{i \in S} T_i^*}{\max(|S|, 1)} \right). \qquad (8) $$

As a consequence, a first way to build an (oracle) selection is to set

$$ S = \{i \in \{1, \dots, n\} : T_i^* \le \alpha\}. $$

Since an average of numbers smaller than $\alpha$ is also smaller than $\alpha$, the corresponding procedure controls the FCR at level $\alpha$. This procedure is referred to as the procedure with fixed threshold in the sequel. It corresponds to the following naive approach: to get a clustering with a risk of $\alpha$, we only keep the items that are in their corresponding class with a posterior probability at least $1 - \alpha$. By contrast, the selection rule considered here is rather

$$ S = \{i \in \{1, \dots, n\} : T_i^* \le t(\alpha)\}, $$

for a threshold $t(\alpha) \ge \alpha$ maximizing $|S|$ under the constraint $\sum_{i \in S} T_i^* \le \alpha|S|$. It uniformly improves the procedure with fixed threshold and will in general lead to a (much) broader selection. This gives rise to the oracle procedure, that can be easily implemented by ordering the $T_i^*$'s, see Algorithm 1.
corresponding $\mathbf{Z} = (Z_1, \dots, Z_n)$), when this sample is iid and generated from $P_{\theta^*}$. We also use $\mathbb{E}_{P_{\theta^*}}$ or $\mathbb{E}_{\mathbf{X} \sim P_{\theta^*}}$ in a similar way.

## 2.2 | Procedure and criteria

Our approach starts with a given clustering rule that aims at recovering the true latent clustering for all observed items. In general, a clustering rule $\hat{\mathbf{Z}} = (\hat{Z}_i)_{1 \le i \le n} \in \{1, \dots, Q\}^n$ is defined as a (measurable) function of the observation $\mathbf{X}$ returning a vector of labels $\hat{\mathbf{Z}}(\mathbf{X})$ for which the label $q$ is assigned to individual $i$ if and only if $\hat{Z}_i(\mathbf{X}) = q$. For better readability, we write $\hat{\mathbf{Z}}$ instead of $\hat{\mathbf{Z}}(\mathbf{X})$ in the sequel. The classification error of $\hat{\mathbf{Z}}$, with respect to specific labels, is given by $\varepsilon(\hat{\mathbf{Z}}, \mathbf{Z}) = \sum_{i=1}^n \mathbb{1}\{Z_i \neq \hat{Z}_i\}$.

In the unsupervised setting, only the partition of the observations is of interest, not the labels themselves. Switching the labels of $\hat{\mathbf{Z}}$ does not change the corresponding partition. This is in stark difference to the supervised setting. An invariant label-switching error is the clustering risk of $\hat{\mathbf{Z}}$ defined by

$$R(\hat{\mathbf{Z}}) := \mathbb{E}_{\theta^*} \left( \min_{\sigma \in [Q]} \mathbb{E}_{\theta^*} \left( n^{-1} \varepsilon(\sigma(\hat{\mathbf{Z}}), \mathbf{Z}) \mid \mathbf{X} \right) \right), \quad (1)$$

where $[Q]$ denotes the set of all permutations on $\{1, \dots, Q\}$. The minimum over all permutations $\sigma$ is the way to handle the aforementioned label-switching problem.

**Remark 1.** The position of the minimum w.r.t. $\sigma$ in the expression of the risk (1) matters: the permutation $\sigma$ is allowed to depend on $\mathbf{X}$ but not on $\mathbf{Z}$. Hence, this risk has to be understood as being computed up to a data-dependent label switching. This definition coincides with the usual definition of the misclassification risk in the situation where the true clustering is deterministic, see Lei and Rinaldo (2015), Lu and Zhou (2016). Hence, it can be seen as a natural extension of the latter to a mixture model where the true clustering is random.

Classically, we aim to find a clustering rule $\hat{\mathbf{Z}}$ such that the clustering risk is “small”. However, as mentioned above, whether this is possible or not depends on the intrinsic difficulty of the clustering problem and thus of the true parameter $\theta^*$ (see Figure 1). Therefore, the idea is to provide a selection rule $S$, that is, a (measurable) function of the observation $\mathbf{X}$ returning a subset of indices $S(\mathbf{X}) \subset \{1, \dots, n\}$, such that the clustering risk with restriction to $S(\mathbf{X})$ is small. Again, we write $S$ instead of $S(\mathbf{X})$ for short. Throughout the paper, a procedure refers to a couple $C = (\hat{\mathbf{Z}}, S)$, where $\hat{\mathbf{Z}}$ is a clustering rule and $S$ is a selection rule.

**Definition 1** (False selection rate). The false selection rate (FSR) of a procedure $C = (\hat{\mathbf{Z}}, S)$ is given by

$$\text{FSR}_{\theta^*}(C) := \mathbb{E}_{\theta^*} \left( \min_{\sigma \in [Q]} \mathbb{E}_{\theta^*} \left( \frac{\varepsilon_S(\sigma(\hat{\mathbf{Z}}), \mathbf{Z})}{\max(|S|, 1)} \mid \mathbf{X} \right) \right), \quad (2)$$

where $\varepsilon_S(\hat{\mathbf{Z}}, \mathbf{Z}) = \sum_{i \in S} \mathbb{1}\{Z_i \neq \hat{Z}_i\}$ denotes the misclassification error restricted to subset $S$.
**Algorithm 1 Oracle procedure**

Input: Parameter $\theta^*$, sample $(X_1, \dots, X_n)$, level $\alpha$.

1. Compute the posterior probabilities $P_{\theta^*}(Z_i = q | X_i)$, $1 \le i \le n$, $1 \le q \le Q$;

2. Compute the Bayes clustering $\hat{Z}_i^*$, $1 \le i \le n$, given by (6);

3. Compute the probabilities $T_i^*$, $1 \le i \le n$, given by (7);

4. Order these probabilities in increasing order $T_{(1)}^* \le \dots \le T_{(n)}^*$;

5. Choose $k^*$ the maximum of $k \in \{0, \dots, n\}$ such that $\max(k, 1)^{-1} \sum_{j=1}^k T_{(j)}^*(\mathbf{X}) \le \alpha$;

6. Select $S_\alpha^*$, the index corresponding to the $k^*$ smallest elements among the $T_i^*$'s.

Output: Oracle procedure $C_\alpha = (\hat{\mathbf{Z}}^*, S_\alpha^*)$.

**Algorithm 2 Plug-in procedure**

Input: Sample $(X_1, \dots, X_n)$, level $\alpha$.

1. Compute an estimator $\hat{\theta}$ of $\theta$;

2. Run the oracle procedure given in Algorithm 1 with $\hat{\theta}$ in place of $\theta^*$.

Output: Plug-in procedure $\hat{C}_\alpha^{\text{PI}} = (\hat{\mathbf{Z}}^{\text{PI}}, \hat{S}_\alpha^{\text{PI}})$.

## 3.2 Empirical procedures

**Plug-in procedure** The oracle procedure cannot be used in practice since $\theta^*$ is generally unknown. A natural idea is then to approach $\theta^*$ by an estimator $\hat{\theta}$ and then to plug this estimate into the oracle procedure. The resulting procedure, denoted $\hat{C}_\alpha^{\text{PI}} = (\hat{\mathbf{Z}}^{\text{PI}}, \hat{S}_\alpha^{\text{PI}})$, is called the *plug-in procedure* and is implemented in Algorithm 2.

Obviously, the quality of the plug-in procedure relies on the chosen estimator $\hat{\theta}$. In a mixture model, the maximum likelihood estimator (MLE) is not computable in general. A standard way to approach the MLE is to use an Expectation Maximization (EM) algorithm (Dempster et al., 1977). The property of the EM estimator $\hat{\theta}$ in the current framework has been extensively studied both empirically and theoretically in the literature, see Balakrishnan et al. (2017) and references therein.

In Section 4, we establish that the plug-in procedure has suitable properties: when $n$ tends to infinity, provided that the chosen estimator $\hat{\theta}$ behaves well and under mild regularity assumptions on the model, the FCR of the plug-in procedure is close to the level $\alpha$ while it is nearly optimal in terms of average selection number.

**Bootstrap procedure** Despite the favorable theoretical properties shown in Section 4, the plug-in procedure achieves an FCR that can exceed $\alpha$ in some situations, as we will see in our numerical experiments (Section 5). This is in particular the case when the estimator $\hat{\theta}$ is too rough. Indeed, the uncertainty of $\hat{\theta}$ near $\theta^*$ is ignored by the plug-in procedure.

To take into account this effect, we propose to use a bootstrap approach. It is based on the following result.

**Lemma 3.** For a given level $\alpha \in (0, 1)$, the FCR of the plug-in procedure $\hat{C}_\alpha^{\text{PI}}$ is given by

$$ \text{FCR}(\hat{C}_\alpha^{\text{PI}}) = \mathbb{E}_{\mathbf{X} \sim P_{\theta^*}} \left( \min_{\sigma \in [Q]} \frac{\sum_{i=1}^n \{1 - \ell_{\sigma}(\hat{Z}_i^{\text{PI}}(\mathbf{X})) (X_i, \theta^*)\} 1\{i \in \hat{S}_\alpha^{\text{PI}}(\mathbf{X})\}}{\max(|\hat{S}_\alpha^{\text{PI}}(\mathbf{X})|, 1)} \right). \quad (9) $$

The general idea is as follows: since FCR($\hat{C}_\alpha^{\text{PI}}$) can exceed $\alpha$, we choose $\alpha'$ as large as possible such that $\widehat{\text{FCR}}_{\alpha'} \le \alpha$, for which $\widehat{\text{FCR}}_{\alpha'}$ is a bootstrap approximation of FCR($\hat{C}_{\alpha'}^{\text{PI}}$) based on (9).
In this work, the aim is to find a procedure C such that the false selection rate is controlled at a nominal level $\alpha$, that is, $FSR_{\theta^*}(C) \le \alpha$. Obviously, choosing S empty implies $\epsilon_S(\sigma(\hat{\mathbf{Z}}), \mathbf{Z}) = 0$ for any permutation $\sigma$ and thus satisfies this control. Hence, while maintaining the control $FSR_{\theta^*}(C) \le \alpha$, we aim to classify as much individuals as possible, that is, to make $\mathbb{E}_{\theta^*}|S|$ as large as possible.

The definition of the FSR (2) involves an expectation of a ratio, which is more difficult to handle than a ratio of expectations. Hence, the following simpler alternative criterion will also be useful in our analysis.

**Definition 2** (Marginal false selection rate). The marginal false selection rate (mFSR) of a procedure $C = (\hat{\mathbf{Z}}, S)$ is given by

$$ \text{mFSR}_{\theta^*}(C) := \frac{\mathbb{E}_{\theta^*}\left(\min_{\sigma \in [Q]} \mathbb{E}_{\theta^*}\left(\epsilon_S(\sigma(\hat{\mathbf{Z}}), \mathbf{Z}) \mid \mathbf{X}\right)\right)}{\mathbb{E}_{\theta^*}(|S|)}, \quad (3) $$

with the convention $0/0 = 0$.

Note that the mFSR is similar to the criterion introduced in Denis and Hebiri (2020) in the supervised setting.

## 2.3 | Notation

We will use the following notation: for all $q \in \{1, \dots, Q\}$ and $\theta = (\pi, \phi) \in \Theta$, we let

$$ \ell_q(X, \theta) := \mathbb{P}_{\theta}(Z = q | X) = \frac{\pi_q f_{\phi_q}(X)}{\sum_{\ell=1}^{Q} \pi_{\ell} f_{\phi_{\ell}}(X)}; \quad (4) $$

$$ T(X, \theta) := 1 - \max_{q \in \{1, \dots, Q\}} \ell_q(X, \theta) \in [0, 1 - 1/Q]. \quad (5) $$

Hence, $\ell_q(X, \theta)$ is the posterior probability of belonging to class $q$ given the measurement $X$ under the distribution $P_{\theta}$. The quantity $T(X, \theta)$ is a measure of the risk when classifying $X$: it is close to 0 when there exists a class $q$ such that $\ell_q(X, \theta)$ is close to 1, that is, when $X$ can be classified with large confidence.

# 3 | METHODS

In this section, we introduce new methods for controlling the FSR. We start by identifying an oracle method, that uses the true value of the parameter $\theta^*$. Substituting the unknown parameter $\theta^*$ by an estimator in that oracle provides our first method, called the plug-in procedure. We then define a refined version of the plug-in procedure, that accounts for the variability of the estimator and is based on a bootstrap approach.

## 3.1 | Oracle procedures

Here, we proceed as if an oracle had given us the true value of $\theta^*$ and we introduce an oracle procedure based on this value. The oracle procedure $C_{\alpha}^{\text{or}}$ is defined in Algorithm 1 and described in detail below.
**Algorithm 3** Bootstrap procedure

Input: Sample $(X_1, \dots, X_n)$, level $\alpha$, number $B$ of bootstrap runs.

1. Choose a grid of increasing level $(\alpha(k))_{1 \le k \le K}$;

2. Compute $\widehat{\text{FCR}}_{\alpha(k)}^B$, $1 \le k \le K$, according to (10);

3. Choose $\tilde{k}$ according to (11).

Output: Bootstrap procedure $\hat{C}_{\alpha}^{\text{boot}} = \hat{C}_{\alpha(\tilde{k})}^{\text{PI}}$.

The bootstrap approximation reads as follows: in the RHS of (9), we replace the true parameter $\theta^*$ by $\hat{\theta}$ and $\mathbf{X} \sim P_{\theta^*}$ by $\mathbf{X}' \sim \hat{P}$, where $\hat{P}$ is an empirical substitute of $P_{\theta^*}$. This empirical distribution $\hat{P}$ is $P_{\hat{\theta}}$ for the parametric bootstrap and the uniform distribution over the $X_i$'s for the non-parametric bootstrap. This yields the bootstrap approximation of FCR($\hat{C}_{\alpha}^{\text{PI}}$) given by

$$ \widehat{\text{FCR}}_{\alpha} = \mathbb{E}_{\mathbf{X}' \sim \hat{P}} \left( \min_{\sigma \in [Q]} \frac{\sum_{i=1}^{n} \{1 - \ell_{\sigma}(\hat{Z}_{i}^{\text{PI}}(\mathbf{X}'))(X_i^b, \hat{\theta}(\mathbf{X}))\} \mathbf{1}\{i \in \hat{S}_{\alpha}^{\text{PI}}(\mathbf{X}')\}}{\max(|\hat{S}_{\alpha}^{\text{PI}}(\mathbf{X}')|, 1)} \middle| \mathbf{X} \right). $$

Classically, the latter is itself approximated by a Monte Carlo scheme:

$$ \widehat{\text{FCR}}_{\alpha}^{B} = \frac{1}{B} \sum_{b=1}^{B} \min_{\sigma \in [Q]} \frac{\sum_{i=1}^{n} \{1 - \ell_{\sigma}(\hat{Z}_{i}^{\text{PI}}(\mathbf{X}^{b}))(X_{i}^{b}, \hat{\theta}(\mathbf{X}))\} \mathbf{1}\{i \in \hat{S}_{\alpha}^{\text{PI}}(\mathbf{X}^{b})\}}{\max(|\hat{S}_{\alpha}^{\text{PI}}(\mathbf{X}^{b})|, 1)}, \quad (10) $$

with $\mathbf{X}^1, \dots, \mathbf{X}^B$ i.i.d. $\sim \hat{P}$ corresponding to the bootstrap samples of $\mathbf{X}$.

Let $(\alpha(k))_{1 \le k \le K} \in (0, 1)^K$ a grid of increasing nominal levels (possibly in restriction to values slightly below the targeted level $\alpha$). Then, the bootstrap procedure at level $\alpha$ is defined as $\hat{C}_{\alpha}^{\text{boot}} = \hat{C}_{\alpha(\tilde{k})}^{\text{PI}}$ where

$$ \tilde{k} = \max \left\{ k \in \{1, \dots, K\} : \widehat{\text{FCR}}_{\alpha(k)}^{B} \le \alpha \right\}. \quad (11) $$

This procedure is implemented in Algorithm 3.

**Remark 2** (Parametric versus non parametric bootstrap). The usual difference between parametric and non parametric bootstrap also holds in our context: the parametric bootstrap is fully based on $P_{\hat{\theta}}$, while the non parametric bootstrap builds an artificial sample (with replacement) from the original sample, which does not come from a $P_{\theta}$-type distribution. This gives rise to different behaviors in practice: when $\hat{\theta}$ is too optimistic (which will be typically the case here when the estimation error is large), the correction brought by the parametric bootstrap (based on $P_{\hat{\theta}}$) is often weaker than that of the non parametric one. By contrast, when $\hat{\theta}$ is close to the true parameter, the parametric bootstrap approximation is more faithful because it uses the model, see Section 5.

# 4 Theoretical guarantees for the plug-in procedure

In this section, we derive theoretical properties for the plug-in procedure: we show that its FCR and mFCR is close to $\alpha$ while its expected selection number is close to be optimal under some conditions.
### 3.1.1 | MAP clustering

As the following lemma shows, the best clustering rule is well-known and given by the Bayes clustering $\hat{\mathbf{Z}}^* = (\hat{Z}_1^*,...,\hat{Z}_n^*)$, which can be written as the maximum a posteriori (MAP) rule

$$ \hat{Z}_i^* \in \underset{q \in \{1, \dots, Q\}}{\operatorname{argmax}} \ell_q(X_i, \theta^*), \quad i \in \{1, \dots, n\}, \qquad (6) $$

where $\ell_q(\cdot)$ is the posterior probability given by (4).

**Lemma 1.** We have $\min_{\hat{\mathbf{Z}}} R(\hat{\mathbf{Z}}) = R(\hat{\mathbf{Z}}^*) = n^{-1} \sum_{i=1}^n \mathbb{E}_{\theta^*}(T_i^*) = \mathbb{E}_{\theta^*}(T_1^*)$, for the Bayes clustering $\hat{\mathbf{Z}}^*$ defined by (6) and for

$$ T_i^* := T(X_i, \theta^*) = \mathbb{P}_{\theta^*}(Z_i \neq \hat{Z}_i^* | X_i), \quad i \in \{1, \dots, n\}, \qquad (7) $$

where $T(\cdot)$ is given by (5).

Lemma 1 is proved in Appendix B1. In words, it states that the oracle statistics $T_i^*$ correspond to the posterior misclassification probabilities of the Bayes clustering. To decrease the overall misclassification risk, it is natural to avoid the classification of points with a high value of the test statistic $T_i^*$.

### 3.1.2 | Thresholding selection rules

In this section, we introduce the selection rule that decides which items are to be classified. From the above paragraph, it is natural to consider a thresholding-based selection rule of the form $S_t^* := \{i \in \{1, \dots, n\} : T_i^* \le t\}$, for some threshold $t$ to be chosen suitably. The following result gives insights for the choice of such a threshold $t$.

**Lemma 2.** For a procedure $C = (\hat{\mathbf{Z}}^*, S)$ with Bayes clustering and an arbitrary selection $S$,

$$ \text{FSR}_{\theta^*}(C) = \mathbb{E}_{\theta^*} \left( \frac{\sum_{i \in S} T_i^*}{\max(|S|, 1)} \right). \qquad (8) $$

Lemma 2 is proved in Appendix B2. As a consequence, a first way to build an (oracle) selection with $\text{FSR}_{\theta^*}(C) \le \alpha$ is to set $S_\alpha^* = \{i \in \{1, \dots, n\} : T_i^* \le \alpha\}$ (an average of numbers smaller than $\alpha$ is also smaller than $\alpha$). This procedure is referred to as the procedure with fixed threshold in the sequel. It corresponds to the following naive approach: to get a clustering with a risk of $\alpha$, we only keep the items that are in their corresponding class with a posterior probability of at least $1 - \alpha$.

By contrast, the oracle selection rule considered here is given by $S^* := \{i \in \{1, \dots, n\} : T_i^* \le t(\alpha)\}$, for a threshold $t(\alpha) \ge \alpha$ maximizing $|S|$ under the constraint $\sum_{i \in S} T_i^* \le \alpha|S|$. It has a threshold that adapts to $\theta^*$, uniformly improves the procedure with fixed threshold and will in general lead to a (much) broader selection. This gives rise to the oracle procedure, that can be easily implemented by ordering the $T_i^*$'s, see Algorithm 1.
## 4.1 Additional notation and assumptions

We make use of an optimality theory for mFCR control, that will be developed in detail in Section A.1. This approach extensively relies on the following quantities (recall the definition of $T(X, \theta)$ in (5)):

$$ \mathrm{mFCR}_{t}^{*} = \mathbb{E}_{\theta^{*}} (T(X, \theta^{*}) | T(X, \theta^{*}) < t); \qquad (12) $$

$$ t^{*}(\alpha) = \sup \{t \in [0, 1] : \mathrm{mFCR}_{t}^{*} \le \alpha\} \qquad (13) $$

$$ \alpha_c = \inf\{\mathrm{mFCR}_t^* : t \in (0, 1], \mathrm{mFCR}_t^* > 0\}; \qquad (14) $$

$$ \bar{\alpha} = \mathrm{mFCR}_1^*. \qquad (15) $$

In words, $\mathrm{mFCR}_t^*$ is the mFCR of an oracle procedure that selects the $T_t^*$ smaller than a threshold $t$ (Lemma 11). Then, $t^*(\alpha)$ is the optimal threshold such that this procedure has a mFCR controlled at level $\alpha$. Next, $\alpha_c$ and $\bar{\alpha}$ are the lower and upper bounds for the level $\alpha$, respectively, for which the optimality theory can be applied.

Now, we introduce our main assumption, which will be ubiquitous in our analysis.

**Assumption 1.** For all $\theta \in \Theta$, under $P_{\theta^*}$, the r.v. $T(X, \theta)$ given by (5) is continuous, that is, $t \mapsto P_{\theta^*}(T(X, \theta) < t)$ is continuous on $[0, 1]$. In addition, it is increasing on $(\alpha_c, \bar{\alpha})$.

Assumption 1 is useful to several regards: first, it prohibits ties in the $T(X_i, \theta)$'s, $1 \le i \le m$, so that the selection rule (see Algorithm 1) can be truly formulated as a thresholding rule. Second, it entails interesting properties for function $t \mapsto \mathrm{mFCR}_t^*$, see Lemma 11 (this in particular ensures that the supremum in (13) is a maximum). Also note that the inequality $0 \le \alpha_c < \bar{\alpha} < 1 - 1/Q$ holds under Assumption 1.

The next assumption ensures that the density family $\{f_u, u \in \mathcal{U}\}$ is smooth, and will be useful to establish consistency results.

**Assumption 2.** For $P_{\theta^*}$-almost all $x \in \mathbb{R}^d$, $u \in \mathcal{U} \mapsto f_u(x)$ is continuous.

Let us now introduce some additional quantities measuring the regularity of the model and that will appear in our remainder terms (recall definition (4), (5) and (13) of $l_q(X, \theta)$, $T(X, \theta)$ and $t^*(\alpha)$, respectively): for $\epsilon, \delta, v > 0$,

$$ W_{\ell}(\epsilon) = \sup \left\{ \mathbb{E}_{\theta^*} \left[ \max_{1 \le q \le Q} |l_q(X, \theta^*) - l_q(X, \theta)| \right], \|\theta - \theta^*\|_2 \le \epsilon, \theta \in \Theta \right\}; \qquad (16) $$

$$ W_T(\delta) = \sup\{|P_{\theta^*}(T(X, \theta^*) < t') - P_{\theta^*}(T(X, \theta^*) < t)|, t, t' \in [0, 1], |t' - t| \le \delta\}; \qquad (17) $$

$$ \Psi(\epsilon) = \inf_{\delta \in (0,1)} \{W_T(\delta) + W_{\ell}(\epsilon)/\delta\}; \qquad (18) $$

$$ W_{t^*, \alpha}(v) = \sup \{|t^*(\alpha + \beta) - t^*(\alpha)|, |\beta| \le v\}. \qquad (19) $$

Note that $\Psi(\epsilon) \le W_T(W_{\ell}(\epsilon)^{1/2}) + W_{\ell}(\epsilon)^{1/2}$ for $\epsilon$ small enough. The following lemma is straightforward.

**Lemma 4.** Under Assumption 1, we have $\lim_{\delta \to 0} W_T(\delta) = 0$, $\lim_{v \to 0} W_{t^*, \alpha}(v) = 0$. Under Assumption 2, we have $\lim_{\epsilon \to 0} W_{\ell}(\epsilon) = 0$. Under both assumptions, we have $\lim_{\epsilon \to 0} \Psi(\epsilon) = 0$.

We can derive rates for the convergences in Lemma 4 under the following additional regularity condition.

**Assumption 3.** There exists positive constants $r = r(\theta^*)$, $C_1 = C_1(\theta^*)$, $C_2 = C_2(\theta^*, \alpha)$, $C_3 = C_3(\theta^*, \alpha)$ such that
**Algorithm 1.** Oracle procedure

Input: Parameter $\theta^*$, sample $(X_1, \dots, X_n)$, level $\alpha$.

1. Compute the posterior probabilities $P_{\theta^*}(Z_i = q | X_i)$, $1 \le i \le n$, $1 \le q \le Q$;

2. Compute the Bayes clustering $\hat{Z}_i^*$, $1 \le i \le n$, according to (6);

3. Compute the probabilities $T_i^*$, $1 \le i \le n$, according to (7);

4. Order these probabilities in increasing order $T_{(1)}^* \le \dots \le T_{(n)}^*$;

5. Choose $k^*$ the maximum of $k \in \{0, \dots, n\}$ such that $\max(k, 1)^{-1} \sum_{j=1}^k T_{(j)}^*(\mathbf{X}) \le \alpha$;

6. Set $S^*$ to the set of indices of the $k^*$ smallest elements among the $T_i^*$'s.

Output: Oracle procedure $C_{\alpha}^{or} = (\hat{Z}^*, S^*)$.

**Algorithm 2.** Plug-in procedure

Input: Sample $(X_1, \dots, X_n)$, level $\alpha$.

1. Compute an estimator $\hat{\theta}$ of $\theta$ based on $(X_1, \dots, X_n)$;

2. Run the oracle procedure given in Algorithm 1 with $\hat{\theta}$ in place of $\theta^*$, returning a procedure that is denoted $\hat{C}_{\alpha}^{\mathrm{PI}} = (\hat{Z}^{\mathrm{PI}}, \hat{S}_{\alpha}^{\mathrm{PI}})$.

Output: Plug-in procedure $\hat{C}_{\alpha}^{\mathrm{PI}} = (\hat{Z}^{\mathrm{PI}}, \hat{S}_{\alpha}^{\mathrm{PI}})$.

## 3.2 | Empirical procedures

### 3.2.1 | Plug-in procedure

The oracle procedure cannot be used in practice since $\theta^*$ is generally unknown. A natural idea then is to replace $\theta^*$ by an estimator $\hat{\theta}$ and to plug this estimate into the oracle procedure. The resulting procedure, denoted $\hat{C}^{\mathrm{PI}} = (\hat{Z}^{\mathrm{PI}}, \hat{S}_{\alpha}^{\mathrm{PI}})$, is called the *plug-in procedure* and is implemented in Algorithm 2.

In Section 4, we establish that the plug-in procedure has convenient properties: when $n$ tends to infinity, provided that the chosen estimator $\hat{\theta}$ behaves well and under mild regularity assumptions on the model, the FSR of the plug-in procedure is close to the level $\alpha$, while it is nearly optimal in terms of the average selection number.

### 3.2.2 | Bootstrap procedure

Despite the favorable theoretical properties shown in Section 4, the plug-in procedure achieves an FSR that can exceed $\alpha$ in some situations, as we will see in our numerical experiments (Section 5). This is in particular the case when the estimator $\hat{\theta}$ is too rough. Indeed, the uncertainty of $\hat{\theta}$ is ignored by the plug-in procedure.

To take into account this effect, we propose to use a bootstrap approach. It is based on the following result.

**Lemma 3.** For a given level $\alpha \in (0, 1)$, the FSR of the plug-in procedure $\hat{C}_{\alpha}^{\mathrm{PI}}$ is given by

$$ \mathrm{FSR}_{\theta^*}(\hat{C}_{\alpha}^{\mathrm{PI}}) = \mathbb{E}_{\mathbf{X} \sim P_{\theta^*}} \left( \min_{\sigma \in [Q]} \frac{\sum_{i=1}^n \{1 - \ell_{\sigma}(\hat{Z}_i^{\mathrm{PI}})(\mathbf{X}_i, \theta^*)\} \mathbb{1}\{i \in \hat{S}_{\alpha}^{\mathrm{PI}}(\mathbf{X}_i)\}}{\max(|\hat{S}_{\alpha}^{\mathrm{PI}}(\mathbf{X}_i)|, 1)} \right). \quad (9) $$
(i) for $\mathbb{P}_{\theta^*}$-almost all $x$, $u \in \mathcal{U} \mapsto f_u(x)$ is continuously differentiable, and

$$ \sum_{1 \le q \le Q} \mathbb{E}_{\theta^*} \left( \sup_{\substack{\theta \in \Theta \\ \|\theta - \theta^*\| \le r}} \| \nabla_{\theta} \ell_q(X, \theta) \| \right) \le C_1; $$

(ii) for all $t, t' \in [0, 1]$, $|\mathbb{P}_{\theta^*}(T(X, \theta^*) < t) - \mathbb{P}_{\theta^*}(T(X, \theta^*) < t')| \le C_2|t - t|$;

(iii) for all $\beta \in [(\alpha_c + \alpha)/2, (\alpha + \bar{\alpha})/2]$, $|t^*(\beta) - t^*(\alpha)| \le C_3|\beta - \alpha|$.

**Lemma 5.** Under Assumption 3, we have $W_\ell(\epsilon) \le C_1\epsilon$, $W_T(\delta) \le C_2\delta$, $W_{t^*,\alpha}(v) \le C_3v$ and $\Psi(\epsilon) \le (1+C_2)C_1^{1/2}\epsilon^{1/2}$ for $\epsilon, \delta, v$ small enough.

In Appendix D, it is proved that Assumptions 1, 2 and 3 hold true in the homoscedastic two-component multivariate Gaussian model, see Lemma 18.

Let us now discuss conditions on the estimator $\hat{\theta}$ on which the plug-in procedure is based. We start by introducing the following assumption (used in the concentration part of the proof, see Lemma 14).

**Assumption 4.** The estimator $\hat{\theta}$ is assumed to be valued in a countable subset $\mathcal{D}$ of $\Theta$.

This assumption is a minor restriction, because we can always choose $\mathcal{D} \subset \mathbb{Q}^K$ (recall $\Theta \subset \mathbb{R}^K$). Next, we additionally define a quantity measuring the quality of the estimator: for all $\epsilon > 0$,

$$ \eta(\epsilon, \theta^*) = \mathbb{P}_{\theta^*} \left( \min_{\sigma \in [Q]} \|\hat{\theta}^\sigma - \theta^*\|_2 \ge \epsilon \right). \qquad (20) $$

In the Gaussian case, the quantity $\eta(\epsilon, \theta^*)$ is provided to be small for a constrained MLE, see Proposition 16.

## 4.2 A general result

The following theorem is our main theoretical result: it provides non-asymptotical bounds on the mFCR/FCR of the plug-in procedure and on its average selection number.

**Theorem 6.** Let Assumption 1 be true. For any $\alpha \in (\alpha_c, \bar{\alpha})$, for $s^* = s^*(\alpha, \theta^*) \in (0, 1)$ and $e^* = e(\alpha, \theta^*) > 0$ two constants only depending on $\alpha$ and $\theta^*$, the following holds. Consider the plug-in procedure $\hat{C}_\alpha^{PI} = (\hat{Z}^{PI}, \hat{S}_\alpha^{PI})$ introduced in Algorithm 2 and based on an estimator $\hat{\theta}$ satisfying Assumption 4, with $\eta(\epsilon, \theta^*)$ defined by (20). Then for $\epsilon \le e^*$ and $n \ge 4$, letting

$$ \Delta_n(\epsilon) = 2(W_T(W_{t^*,\alpha}(2\delta_n + 8\Psi(\epsilon)/s^*)) + 4\Psi(\epsilon) + 2\delta_n) + 50/n + 10\eta(\epsilon, \theta^*), $$

for $\delta_n = (2/s^*)\sqrt{(\log n)/n}$ and with the quantities $W_T$, $W_\ell$, $\Psi$, $W_{t^*,\alpha}$ defined by (17), (16), (18), (19), respectively, it holds:

* $\hat{C}_\alpha^{PI}$ controls both the FCR and mFCR at level close to $\alpha$ in the following sense:

$$ \max(\text{mFCR}(\hat{C}_\alpha^{PI}), \text{FCR}(\hat{C}_\alpha^{PI})) \le \alpha + \Delta_n(\epsilon)/s^*; $$
**Algorithm 3.** Bootstrap procedure

Input: Sample $(X_1, \dots, X_n)$, level $\alpha$, number $B$ of bootstrap runs.

1. Choose a grid of increasing levels $(\alpha(k))_{1 \le k \le K}$;

2. Draw $B$ bootstrap samples and apply Algorithm 2 to each of them.

3. Compute $\widehat{\text{FSR}}_{\alpha(k)}^B$, $1 \le k \le K$, according to (10);

4. Choose $\tilde{k}$ according to (11).

Output: Bootstrap procedure $\hat{C}_{\alpha}^{\text{boot}} = \hat{C}_{\alpha(\tilde{k})}^{\text{PI}}$.

Lemma 3 is proved in Appendix B3. The general idea is as follows: since $\text{FSR}_{\theta^*}(\hat{C}_{\alpha}^{\text{PI}})$ can exceed $\alpha$, we choose $\alpha'$ as large as possible such that $\widehat{\text{FSR}}_{\alpha'} \le \alpha$, where $\widehat{\text{FSR}}_{\alpha'}$ is a bootstrap approximation of $\text{FSR}_{\theta^*}(\hat{C}_{\alpha'}^{\text{PI}})$ based on (9).

The bootstrap approximation reads as follows: in the RHS of (9), we replace the true parameter $\theta^*$ by $\hat{\theta}$ and $\mathbf{X} \sim P_{\theta^*}$ by $\mathbf{X}' \sim \hat{P}$, where $\hat{P}$ is an empirical substitute of $P_{\theta^*}$. This empirical distribution $\hat{P}$ is $P_{\hat{\theta}}$ for the parametric bootstrap and the uniform distribution over the $X_i$'s for the nonparametric bootstrap. This yields the bootstrap approximation of $\text{FSR}_{\theta^*}(\hat{C}_{\alpha}^{\text{PI}})$ given by

$$ \widehat{\text{FSR}}_{\alpha} := \mathbb{E}_{\mathbf{X}' \sim \hat{P}} \left( \min_{\sigma \in [Q]} \frac{\sum_{i=1}^n \{1 - \ell_{\sigma}(\hat{Z}_i^{\text{PI}}(\mathbf{X}'))(\mathbf{X}'_i, \hat{\theta}(\mathbf{X}))\} \mathbb{1}\{i \in \hat{S}_{\alpha}^{\text{PI}}(\mathbf{X}')\}}{\max(|\hat{S}_{\alpha}^{\text{PI}}(\mathbf{X}')|, 1)} \middle| \mathbf{X} \right). $$

Above, $\mathbb{E}_{\mathbf{X}' \sim \hat{P}}$ denotes the expectation operator when $\mathbf{X}' = (X'_1, \dots, X'_n)$ is iid $\sim \hat{P}$. Classically, the latter is itself approximated by a Monte-Carlo scheme:

$$ \widehat{\text{FSR}}_{\alpha}^{B} := \frac{1}{B} \sum_{b=1}^{B} \min_{\sigma \in [Q]} \frac{\sum_{i=1}^{n} \{1 - \ell_{\sigma}(\hat{Z}_{i}^{\text{PI}}(\mathbf{X}^{b}))(\mathbf{X}_{i}^{b}, \hat{\theta}(\mathbf{X}))\} \mathbb{1}\{i \in \hat{S}_{\alpha}^{\text{PI}}(\mathbf{X}^{b})\}}{\max(|\hat{S}_{\alpha}^{\text{PI}}(\mathbf{X}^{b})|, 1)}, \quad (10) $$

with bootstrap samples $\mathbf{X}^1, \dots, \mathbf{X}^B$ i.i.d. $\sim \hat{P}$.

Let $(\alpha(k))_{1 \le k \le K} \in (0, 1)^K$ be a grid of increasing nominal levels (possibly with restriction to values slightly below the target level $\alpha$). Then, the bootstrap procedure at level $\alpha$ is defined as $\hat{C}_{\alpha}^{\text{boot}} = \hat{C}_{\alpha(\tilde{k})}^{\text{PI}}$, where

$$ \tilde{k} = \max \left\{ k \in \{1, \dots, K\} : \widehat{\text{FSR}}_{\alpha(k)}^{B} \le \alpha \right\}. \quad (11) $$

This procedure is described in Algorithm 3.

**Remark 2** (Parametric versus nonparametric bootstrap). The usual difference between parametric and nonparametric bootstrap also holds in our context: the parametric bootstrap is fully based on $P_{\hat{\theta}}$, while the nonparametric bootstrap builds an artificial sample (with replacement) from the original sample, which does not come from a $P_{\theta}$-type distribution. This gives rise to different behaviors in practice: when $\hat{\theta}$ is too optimistic (that is, in cases where the estimated posterior probability to be in each class is biased towards 0 or 1, which will be typically the case when the estimation error is large), the correction brought by the parametric bootstrap (based on $P_{\hat{\theta}}$) is often weaker and insufficient compared to that of the nonparametric one. By contrast,
• $\hat{C}_{\alpha}^{PI}$ is nearly optimal in the following sense: for any other procedure $C = (\hat{Z}, S)$ that controls the mFCR at level $\alpha$, that is, with mFCR($C$) $\le \alpha$, we have

$$n^{-1} \mathbb{E}_{\theta^*}(|\hat{S}_{\alpha}^{PI}|) \ge n^{-1} \mathbb{E}_{\theta^*}(|S|) - \Delta_n(\epsilon).$$

The proof employs techniques that share similarities with the work of Rebafka et al. (2019) developed in a different context, see Appendix A. Here, a difficulty is to handle the new statistics $T(X_i, \hat{\theta})$ which is defined as an extremum, see (5).

## 4.3 Corollaries

In this section, we present two straightforward consequences of the main result. By Lemma 4, assuming Assumptions 1 and 2, we have that for all $\epsilon$ small enough, the limit of the remainder term $\limsup_n \Delta_n(\epsilon)$ can be bounded by a function of $\epsilon$ that vanishes when $\epsilon$ tends to 0. This yields the following consistency result.

**Corollary 7 (Asymptotic optimality of the plug-in procedure).** Consider the setting of Theorem 6 where, in addition, Assumption 2 is supposed to hold. Then, we have

$$\limsup_n \text{FCR}(\hat{C}_{\alpha}^{PI}) \le \alpha, \quad \limsup_n \text{mFCR}(\hat{C}_{\alpha}^{PI}) \le \alpha,$$

and for any procedure $C = (\hat{Z}, S)$ that controls the mFCR at level $\alpha$, we have

$$\liminf_n \{n^{-1} \mathbb{E}_{\theta^*}(|\hat{S}_{\alpha}^{PI}|) - n^{-1} \mathbb{E}_{\theta^*}(|S|)\} \ge 0.$$

As shown by Lemma 5, by strengthening the assumptions, Theorem 6 delivers in addition convergence rates.

**Corollary 8 (Optimality of the plug-in procedure with rates).** Consider the setting of Theorem 6, where in addition Assumption 3 holds. Then, with constants $C_1 = C_1(\alpha, \theta^*) > 0$ and $C_2 = C_2(\alpha, \theta^*) > 0$, we have for any sequence $\epsilon_n > 0$ tending to zero, for $n$ larger than a constant only depending on $\alpha$ and $\theta^*$,

$$\max((\text{FCR}(\hat{C}_{\alpha}^{PI}), \text{mFCR}(\hat{C}_{\alpha}^{PI})) \le \alpha + C_1 \left( \epsilon_n^{1/2} + \sqrt{(\log n)/n} + \eta(\epsilon_n, \theta^*) \right) \quad (21)$$

$$n^{-1} \mathbb{E}_{\theta^*}(|\hat{S}_{\alpha}^{PI}|) - n^{-1} \mathbb{E}_{\theta^*}(|S|) \ge -C_2 \left( \epsilon_n^{1/2} + \sqrt{(\log n)/n} + \eta(\epsilon_n, \theta^*) \right), \quad (22)$$

for any procedure $C = (\hat{Z}, S)$ that controls (non asymptotically) the mFCR at level $\alpha$.

The assumptions of Corollary 8 are all satisfied in the homoscedastic two-component multivariate Gaussian model. Applying this result in this case with $\epsilon_n = c\sqrt{(\log n)/n}$ (for some constant $c > 0$ not depending on $n$) and with $\hat{\theta}$ being a constrained MLE, Proposition 16 ensures that $\eta(\epsilon_n, \theta^*)$ is smaller than $1/n$. This indicates that the convergence rates in (21) and (22) are of order $((\log n)/n)^{1/4}$ in this case.

## 5 Experiments

In this section, we evaluate the behavior of the new procedures: plug-in (Algorithm 2), parametric bootstrap and non parametric bootstrap (Algorithm 3). For this, we use both synthetic and real data.
when $\hat{\theta}$ is close to the true parameter, the parametric bootstrap approximation is more faithful because it uses the correct model, see Section 5.

Remark 3 (Weighted likelihood-based nonparametric bootstrap). We explored another bootstrap approach in Appendix E1, which is suitable in the presence of small or overlapping clusters. It is based on weighting the likelihood which is conceptually different, see O'Hagan et al. (2019). However, it seems to provide only a minor improvement of the plug-in procedure compared to the parametric/nonparametric bootstrap approaches mentioned above.

# 4 | THEORETICAL GUARANTEES FOR THE PLUG-IN PROCEDURE

In this section, we derive theoretical properties for the plug-in procedure: we show that its FSR and mFSR are close to $\alpha$, while its expected selection number is close to be optimal under some conditions.

## 4.1 | Additional notation and assumptions

We make use of an optimality theory for mFSR control, that will be developed in detail in Appendix A4. This approach extensively relies on the following quantities (recall the definition of $T(X, \theta)$ in (5)):

$$ \text{mFSR}_{t}^{*} := \mathbb{E}_{\theta^{*}}(T(X, \theta^{*}) | T(X, \theta^{*}) < t); \qquad (12) $$

$$ t^{*}(\alpha) := \sup \{t \in [0, 1] : \text{mFSR}_{t}^{*} \le \alpha\}; \qquad (13) $$

$$ \alpha_c := \inf\{\text{mFSR}_t^* : t \in (0, 1], \text{mFSR}_t^* > 0\}; \qquad (14) $$

$$ \bar{\alpha} := \text{mFSR}_1^*. \qquad (15) $$

In words, $\text{mFSR}_t^*$ is the mFSR of a procedure that selects the items with $T_t^*$ smaller than some threshold $t$ (Lemma 8). Then, $t^*(\alpha)$ is the optimal threshold such that this procedure has an mFSR controlled at level $\alpha$. (This optimal procedure, studied in Appendix A4, is not the same as the oracle procedure $C_\alpha^{or}$ defined in Section 3.1, although these two procedures are expected to behave roughly in the same way for a large $n$.) Next, $\alpha_c$ and $\bar{\alpha}$ are the lower and upper bounds for the nominal level $\alpha$, respectively, for which the optimality theory can be applied.

Now, we introduce our main assumption, which will be ubiquitous in our analysis.

**Assumption 1.** For all $\theta \in \Theta$ and $q \in \{1, \dots, Q\}$, under $P_{\theta^*}$, the r.v. $\ell_q(X, \theta)$ given by (4) is continuous. In addition, the function $t \mapsto \mathbb{P}_{\theta^*}(T(X, \theta) < t)$ is increasing on $(\alpha_c, \bar{\alpha})$, where $T(X, \theta)$ is given by (5).

Note that Assumption 1 implies the continuity of the r.v. $T(X, \theta)$. Indeed, $\mathbb{P}(T(X, \theta) = t) \le \sum_{q=1}^Q \mathbb{P}(\ell_q(X, \theta) = 1 - t)$. Hence, this assumption implies that $t \mapsto \mathbb{P}_{\theta^*}(T(X, \theta) < t)$ is both continuous on $[0, 1]$ and increasing on $(\alpha_c, \bar{\alpha})$. This is useful in several regards: first, it prohibits ties in the $T(X_i, \theta)$'s, $1 \le i \le m$, so that the selection rule (see Algorithm 1) can be truly formulated as a thresholding rule (see Lemma 9). Second, it entails interesting properties for function $t \mapsto \text{mFSR}_t^*$,
<figure><img src="image_3.png" /><figcaption>Figure 2: FCR (left panel) and selection frequency (right panel) as a function of the mean separation $\epsilon$. Known mixture proportions and covariances setting with $Q = 2$, $n = 100$, $d = 2$, $\alpha = 0.1$.</figcaption></figure>

## 5.1 Synthetic data set

The performance of our procedures is studied via simulations in different settings with various difficulties. All of them are gaussian mixture models, with possible restrictions on the parameter space. For parameter estimation, the classical EM algorithm is applied with 100 iterations and 10 starting points chosen with Kmeans++ (Arthur and Vassilvitskii, 2006). In the bootstrap procedures $B = 1000$ bootstrap samples are generated. The performance of all procedures is assessed via the sample FCR and the proportion of classified data points, which is referred to as the selection frequency. For every setting and every set of parameters, depicted results display the mean over 100 simulated datasets. As a baseline, we consider the fixed threshold procedure in which one selects data points that have a maximum posterior group membership probability that exceeds $1 - \alpha$. The oracle procedure (Algorithm 1) is also considered in our experiments for comparison.

Finally note that some figures are postponed in Appendix E for space reasons.

**Known proportions and covariances** In the first setting, the true mixture proportions and covariance matrices are known and used in the EM algorithm. We consider the case $Q = 2$, $\pi_1 = \pi_2 = 1/2$ and $\Sigma_1 = \Sigma_2 = I_d$ with $I_d$ the ($d \times d$)-identity matrix. For the mean vectors, we set $\mu_1 = 0$ and $\mu_2 = (\epsilon/\sqrt{d}, \dots, \epsilon/\sqrt{d})$. The quantity $\epsilon$ corresponds to the mean separation, that is, $\|\mu_1 - \mu_2\|_2 = \epsilon$ and accounts for the difficulty of the clustering problem.

Figure 2 displays the FCR for nominal level $\alpha = 0.1$, sample size $n = 100$, dimension $d = 2$ and varying mean separation $\epsilon \in \{1, \sqrt{2}, 2, 4\}$. Globally, our procedures all have an FCR close to the target level $\alpha$ (excepted for the very well separated case $\epsilon = 4$ for which the FCR is much smaller because a large part of the items can be trivially classified). In addition, the selection rate is always close to the one of the oracle procedure. On the other hand, the baseline procedure is too conservative: its FCR can be well below the nominal level and it selects up to 50% less than the other procedures. This is well expected, because unlike our procedures, the baseline has a fixed threshold and thus does not adapt to the difficulty of the problem.

We also note that the FCR of the plug-in approach is slightly inflated for a weak separation ($\epsilon = 1$). This comes from the parameter estimation, which is difficult in that case. This also illustrates the interest of the bootstrap methods, that allow to recover the correct level in that case, by appropriately correcting the plug-in approach.
see Lemma 8 (this in particular ensures that the supremum in (13) is a maximum). Also note that the inequality $0 \le \alpha_c < \bar{\alpha} < 1 - 1/Q$ holds under Assumption 1.

The next assumption ensures that the density family $\{f_u, u \in \mathcal{U}\}$ is smooth, and will be useful to establish consistency results.

**Assumption 2.** For $P_{\theta^*}$-almost all $x \in \mathbb{R}^d$, $u \in \mathcal{U} \mapsto f_u(x)$ is continuous.

Moreover, we can derive convergence rates under the following additional regularity conditions.

**Assumption 3.** There exist positive constants $r = r(\theta^*), C_1 = C_1(\theta^*), C_2 = C_2(\theta^*, \alpha), C_3 = C_3(\theta^*, \alpha)$ such that

(i) for $P_{\theta^*}$-almost all $x \in \mathbb{R}^d$, $u \in \mathcal{U} \mapsto f_u(x)$ is continuously differentiable, and

$$ \sum_{1 \le q \le Q} \mathbb{E}_{\theta^*} \sup_{\substack{\theta \in \Theta \\ \|\theta - \theta^*\|_2 \le r}} \|\nabla_{\theta} \ell_q(X, \theta)\|_2 \le C_1; $$

(ii) for all $t, t' \in [0, 1]$, $|P_{\theta^*}(T(X, \theta^*) < t) - P_{\theta^*}(T(X, \theta^*) < t')| \le C_2|t - t'|;

(iii) for all $\beta \in [(\alpha_c + \alpha)/2, (\alpha + \bar{\alpha})/2]$, $|t^*(\beta) - t^*(\alpha)| \le C_3|\beta - \alpha|$.

**Example 1.** In Appendix D, it is proved that Assumptions 1, 2 and 3 hold true in the homoscedastic two-component multivariate Gaussian mixture model, see Lemma 17.

Next, we consider the following complexity assumption to ensure concentration of the underlying empirical processes. It is given in terms of the Vapnik-Chervonenkis (VC) dimension of specific function classes involving $\ell_q$. In the sequel, the VC dimension of a function set $\mathcal{F}$ is defined as the VC dimension of the set family $\{\{x \in \mathbb{R}^d : f(x) \ge u\}, f \in \mathcal{F}, u \in \mathbb{R}\}$, see, e.g., Baraud (2016). We denote

$$ \mathcal{V} = \text{VC dimension of } \{\ell_q(., \theta), \theta \in \Theta, 1 \le q \le Q\}; \quad (16) $$

$$ \mathcal{V}_{-} = \text{VC dimension of } \{1\{\ell_q(., \theta) - \ell_{q'}(., \theta) \ge 0\}, \theta \in \Theta, 1 \le q, q' \le Q\}. \quad (17) $$

**Assumption 4.** The VC dimensions $\mathcal{V}$ and $\mathcal{V}_{-}$ are finite.

**Example 2.** In the two-component case, $Q = 2$ where $P_{\theta}$ belongs to an exponential family, we have that $\mathcal{V}, \mathcal{V}_{-} = O(k^2 \log(k))$ (see Lemma 13) with $k$ the dimension of the sufficient statistic vector. For instance, $k = d + d^2$ for the Gaussian family, hence $\mathcal{V}, \mathcal{V}_{-} = O(d^4 \log(d))$ in that case. (For the specific case of the homoscedastic Gaussian family, we have that $\mathcal{V}, \mathcal{V}_{-} = O(d)$, see Lemma 18).

Let us now discuss conditions on the estimator $\hat{\theta}$ on which the plug-in procedure is based. We start by introducing the following assumption (used in the concentration part of the proof, see Lemma 11).

**Assumption 5.** The estimator $\hat{\theta}$ is assumed to take its values in a countable subset $\mathcal{D}$ of $\Theta$.

This assumption is a minor restriction, because we can always choose $\mathcal{D} \subset \mathbb{Q}^K$ (recall $\Theta \subset \mathbb{R}^K$). Next, we additionally define a quantity measuring the quality of the estimator: for all $\epsilon > 0$,

$$ \eta(\epsilon, \theta^*) = \mathbb{P}_{\theta^*} \left( \min_{\sigma \in [Q]} \|\hat{\theta}^{\sigma} - \theta^*\|_2 \ge \epsilon \right). \quad (18) $$
<figure><img src="image_4.png" /><figcaption>Figure 3: Same as Figure 2 in the diagonal covariance setting.</figcaption></figure>

<figure><img src="image_5.png" /><figcaption>Figure 4: FCR (left panel) and selection frequency (right panel) as a function of the sample size $n$. Diagonal covariances setting with $Q = 2$, $d = 2$, $\epsilon = \sqrt{2}$, $\alpha = 10\%$.</figcaption></figure>

**Diagonal covariances** In this setting, the true parameters are the same as in previous paragraph, but the true mixture proportions and covariance matrices are unknown and not used in the EM algorithm. However, to help the estimation, the diagonal structure of $\Sigma_1$ and $\Sigma_2$ is used in this algorithm.

Figure 3 displays the FCR and the selection frequency in function of the separation $\epsilon$ for this case. The conclusion is qualitatively the same as in the previous case, but with larger FCR values for a weak separation. Overall, it shows that the plug-in procedure is anti-conservative and that the bootstrap corrections are able to recover an FCR/selection frequency close to the one of the oracle. However, for a weak separation $\epsilon = 1$, note that the parametric bootstrap correction is not enough and that the latter procedure still overshoots the level $\alpha$. Indeed, in our simulations, it appears that $P_{\hat{\theta}}$ is typically a distribution that is more favorable than $P_{\theta^*}$ from a statistical point of view (for instance, with more separated clusters). Note that these conclusions also hold when making the FCR/selection frequency vary with $n$, see Figure 4.

Figure 5 displays the FCR and the selection frequency for varying nominal level $\alpha$, with $\epsilon = \sqrt{2}$ and $n \in \{200, 1000\}$. The plug-in is still anti-conservative while the bootstrap procedures have an FCR that is close to $\alpha$ for both sample sizes and uniformly on the considered $\alpha$ range. Moreover, we note that for all our procedures (including the plug-in), the gap between the FCR and the nominal level is roughly constant with $\alpha$: this illustrates the adaptive aspect of our procedures. This is in contrast with the baseline procedure, for which this gap highly depends on $\alpha$, and which may be either anti-conservative or sub-optimal depending on the $\alpha$ value.
Above, $\hat{\theta}^{\sigma} = (\hat{\pi}^{\sigma}, \hat{\phi}^{\sigma})$ denotes the estimator $\theta^{\sigma} = (\hat{\pi}, \hat{\phi})$ after having permuted the labels according to the permutation $\sigma$, that is, with $\hat{\pi}_{q}^{\sigma} = \hat{\pi}_{\sigma(q)}$ and $\hat{\phi}_{q}^{\sigma} = \hat{\phi}_{\sigma(q)}$ for $1 \le q \le Q$.

**Example 3.** The literature provides several results regarding the estimation of Gaussian mixtures, see e.g. Regev and Vijayaraghavan (2017) for a review. Proposition 1 revisits some of these results, for the estimator derived from the EM algorithm (Balakrishnan et al., 2017; Dempster et al., 1977) and the constrained MLE (Ho & Nguyen, 2016).

## 4.2 | Results

We now state our main results, starting with the consistency of the plug-in procedure.

**Theorem 1** (Asymptotic optimality of the plug-in procedure). Let Assumptions 1, 2, and 4 be true. Consider an estimator $\hat{\theta}$ satisfying Assumption 5 that is assumed to be consistent in the sense that for all $\epsilon > 0$, the probability $\eta(\epsilon, \theta^*)$ given by (18) tends to 0 as $n$ tends to infinity. Then the corresponding plug-in procedure $\hat{C}_{\alpha}^{\text{PI}}$ (Algorithm 2) satisfies the following: for any $\alpha \in (\alpha_c, \bar{\alpha})$, we have

$$ \limsup_n \text{FSR}(\hat{C}_{\alpha}^{\text{PI}}) \le \alpha, \quad \limsup_n \text{mFSR}(\hat{C}_{\alpha}^{\text{PI}}) \le \alpha, \qquad (19) $$

and for any procedure $C = (\hat{Z}, S)$ that controls the mFSR at level $\alpha$, we have

$$ \liminf_n \{n^{-1}\mathbb{E}_{\theta^*}(|\hat{S}_{\alpha}^{\text{PI}}|) - n^{-1}\mathbb{E}_{\theta^*}(|S|)\} \ge 0. \qquad (20) $$

Next, we derive convergence rates under the additional regularity conditions given by Assumption 3.

**Theorem 2** (Optimality of the plug-in procedure with rates). Consider the setting of Theorem 1, where in addition Assumption 3 holds. Take $\alpha \in (\alpha_c, \bar{\alpha})$ and recall $\eta(\epsilon, \theta^*)$ defined by (18). With constants $A > 0$ and $B > 0$ only depending on $Q, C_1, C_2, C_3, \mathcal{V}$ (16), $\mathcal{V}_-(17)\alpha$ and $\theta^*$, we have for any sequence $\epsilon_n > 0$ tending to zero, for $n$ larger than a constant only depending on $\alpha$ and $\theta^*$,

$$ \text{FSR}(\hat{C}_{\alpha}^{\text{PI}}) \le \alpha + A\sqrt{\epsilon_n} + B\sqrt{\log n/n} + 5/n^2 + \eta(\epsilon_n, \theta^*) \qquad (21) $$

$$ n^{-1}\mathbb{E}_{\theta^*}(|\hat{S}_{\alpha}^{\text{PI}}|) - n^{-1}\mathbb{E}_{\theta^*}(|S|) \ge -A\sqrt{\epsilon_n} - B\sqrt{\log n/n} - \eta(\epsilon_n, \theta^*), \qquad (22) $$

for any procedure $C = (\hat{Z}, S)$ that controls the mFSR at level $\alpha$.

The proofs of Theorems 1 and 2 are based on a more general nonasymptotical result, for which the remainder terms are more explicit, see Theorem 3 in Appendix A. A sketch of the proof is as follows: first, we identify the optimal procedure, which corresponds to invert the mFSR$_{\theta^*}$ functional with respect to the threshold of the procedure, which is possible thanks to Assumption 1. It is optimal in the sense that it maximizes the expected size of the selection among thresholding procedures that control the mFSR at level $\alpha$ (see Lemmas 6 and 8). Second, we quantify how the plug-in procedure (which can be written as a thresholding-based procedure, see Lemma 9) gets to
<figure><img src="image_6.png" /><figcaption>Figure 5: FCR (left panel) and selection frequency (right panel) as a function of the nominal level $\alpha$. Diagonal covariances setting with $Q = 2$, $d = 2$, $\epsilon = \sqrt{2}$, $n \in \{200, 1000\}$.</figcaption></figure>
the optimal procedure when $n$ grows. For this, we investigate how the expected quantities vary in $\theta$ (Lemma 10), by using a regularity assumption (specifically, Assumption 2 or Assumption 3) and how the empirical quantities approximate the expected ones (Lemma 11), by using Assumption 4. This allows to approximate the FSR (resp. expected selection size) of the plug-in procedure by the one of the optimal procedure, which leads to (19), (21) (resp. 20, 22). Overall, this proof uses techniques that share similarities with the work of Rebafka et al. (2022) developed in a different context. Here, an additional difficulty is to handle the new statistic $T(X_i, \hat{\theta})$ which is defined as an extremum, see (5).

Theorem 2 establishes that, given a model which is regular enough and a consistent estimator, the plug-in procedure controls the FSR and is asymptotically optimal up to remainder terms which are of the order of $\sqrt{\epsilon_n} + \sqrt{\log n/n} + \eta(\epsilon_n, \theta^*)$. Here, $\epsilon_n$ dominates the convergence rate of the parameter estimate, and is taken large enough to ensure that $\eta(\epsilon_n, \theta^*)$ vanishes.

For instance, in the multivariate Gaussian mixture model (with further assumptions) and by considering either the EM estimator or the constrained MLE, we have $\eta(\epsilon_n, \theta^*) \le 1/n$ for $\epsilon_n = C\sqrt{\log(n)/n}$, see Proposition 1. This implies that the remainder terms in (21) and (22) are at most of order $((\log n)/n)^{1/4}$.

# 5 | EXPERIMENTS

In this section, we evaluate the performance of the new procedures: plug-in (Algorithm 2), parametric bootstrap and nonparametric bootstrap (Algorithm 3). For this, we use both synthetic and real data. We publicly release the code of the numerical experiments at https://github.com/arianemarandon/fsrcontrol.

## 5.1 | Synthetic data

The performance of our procedures is studied via simulations in different settings with various difficulties. In this section, all mixtures are Gaussian, whereas Appendix E2 provides additional results for Student's t-mixture models. For parameter estimation, the classical EM algorithm is applied with 100 iterations and 10 starting points chosen with Kmeans++ (Arthur & Vassilvitskii, 2007). To help avoiding degenerate likelihood maxima when estimating covariance matrices (Chen, 2017), the covariances are regularized by adding 1e-6 to the diagonal after each M-step. In the bootstrap procedures, $B = 1000$ bootstrap samples are generated. The performance of all procedures is assessed via the sample FSRand the proportion of classified data points, which is referred to as the selection frequency. For every setting and every set of parameters, depicted results display the mean over 100 simulated datasets. As a baseline, we consider the fixed threshold procedure in which one selects data points that have a maximum posterior group membership probability that exceeds $1 - \alpha$. The oracle procedure (Algorithm 1) is also considered in our experiments for comparison.

### 5.1.1 | Degree of cluster separation

In the first setting, the true mixture proportions and covariance matrices are known and used in the EM algorithm. We consider the case $Q = 2$, $\pi_1 = \pi_2 = 1/2$ and $\Sigma_1 = \Sigma_2 = I_d$ with $I_d$ the
**Other settings** We also consider the following settings:

* Larger dimension (Figure 9). Same as in the previous paragraph but for $d = 20$.

* Three components (Figure 10): Same as in the previous paragraph but for $Q = 3$, $\pi_1 = \pi_2 = \pi_3 = 1/3$, $\Sigma_1 = \Sigma_2 = \Sigma_3 = I_2$, $\mu_1 = 0$, $\mu_2 = (\epsilon/\sqrt{2}, \epsilon/\sqrt{2})$, $\mu_3 = (0, \sqrt{\epsilon})$ (so that $\mu_3$ is at equal L2-distances from $\mu_1$ and $\mu_2$).

* Unconstrained model (Figure 11): In this case, the EM algorithm is agnostic in the sense that does not use any knowledge, whatsoever, of the structure or values of the parameters. The true parameters are the same as the previous paragraph.

* An additional “typical” setting (Figure 6): an unconstrained model with three components $Q = 3$, dimension $d = 4$ and a separation $\epsilon = 2$.

In the three first settings, the difficulty of the estimation problem is increased with respect to the previous section. Hence, the behaviors the aforementioned procedures is deteriorated. In particular, it puts forward that the non parametric bootstrap can be anti-conservative (resp. over-conservative) in some cases when $n$ is too large (resp. small). Finally, the last setting is intended to be a faithful summary of the procedure behaviors in a “typical” situation. For a fairly large sample size $n = 1000$, the new bootstrap procedures are correctly mimicking the oracle. The more challenging case of a small sample size $n = 100$ is handled differently by these two corrections: the parametric one is anti-conservative while the non-parametric is very conservative.

In conclusion, in our experiments, the parametric bootstrap procedure shows overall the more “stable” behavior: it uniformly improves the plug-in procedure across all the explored parameter ranges. In addition, it gets an FCR and a selection frequency close to the one of the oracle when the sample size $n$ is fairly large. For more challenging cases where the sample size is small and we seek for a stricter FCR control, the non-parametric bootstrap is an interesting alternative.

## 5.2 Real data set

We consider the Wisconsin Breast Cancer Diagnosis (WCBD) dataset from the UCI ML repository. The data consists of features computed from a digitalized image of a fine needle aspirate (FNA) of a breast mass, on a total of 569 patients (each corresponds to one FNA sample) of which 212 are diagnosed as Benign and 357 as Malignant. Ten real-valued measures were computed for each cell nucleus present in any one image (e.g. radius, perimeter, texture, etc.). Then, the mean, standard error and mean of the three largest values of these measures were computed for each image, resulting in a total of 30 features. Here, we restrict the analysis to the variables that correspond to the mean of these measures.

We choose to model the data as a mixture of Students as proposed in Peel and McLachlan (2000). Student mixtures are appropriate for data containing observations with longer than normal tails or atypical observations leading to overlapping clusters. Compared to Gaussian mixtures, Students are less concentrated and thus produce estimates of the posterior probabilities of class memberships that are less extreme, which is favorable for our selection procedures. In our study, the degree of freedom of each component is set to 4, and no constraints are put on the rest of the parameters. The t-mixture is fit via the EM algorithm provided by the Python package `studenttmixture` (Peel and McLachlan, 2000).

For illustration, Figure 7 (panel (a)) displays the first two variables of the dataset, the mean radius and the mean texture of the images. Different colors indicate the ground truth labels
<figure><img src="image_4.png" /><figcaption>FIGURE 2 Cluster separation: (a) FSR and (b) selection frequency as a function of the mean separation $\epsilon$. Known mixture proportions and covariances setting with $Q = 2$, $n = 100$, $d = 2$, $\alpha = 0.1$.</figcaption></figure>

$(d \times d)$-identity matrix. For the mean vectors, we set $\mu_1 = 0$ and $\mu_2 = (\epsilon/\sqrt{d}, \dots, \epsilon/\sqrt{d})$. The quantity $\epsilon$ corresponds to the mean separation, that is, $\|\mu_1 - \mu_2\|_2 = \epsilon$ and accounts for the difficulty of the clustering problem.

Figure 2a displays the FSR for nominal level $\alpha = 0.1$, sample size $n = 100$, dimension $d = 2$ and varying mean separation $\epsilon \in \{1, \sqrt{2}, 2, 4\}$. Globally, our procedures all have an FSR close to the target level $\alpha$ (excepted for the very well separated case $\epsilon = 4$ for which the FSR is much smaller because a large part of the items can be trivially classified). In addition, the selection rate (b) is always close to the one of the oracle procedure. On the other hand, the baseline procedure is rather conservative: its FSR can be well below the nominal level and it selects up to 50% less than the other procedures. This is well expected, because unlike our procedures, the baseline has a fixed threshold and thus does not adapt to the difficulty of the problem.

We also note that the FSR of the plug-in approach is slightly inflated for a weak separation ($\epsilon = 1$). This comes from the parameter estimation, which is difficult in that case. This also illustrates the interest of the bootstrap methods, that allow to recover the correct level in that case, by appropriately correcting the plug-in approach.

### 5.1.2 | Unknown diagonal covariance matrices

In this setting, the true parameters are the same as in the previous paragraph, but the true mixture proportions and covariance matrices are unknown when fitting the mixture model. However, to help the estimation, we suppose a diagonal structure for $\Sigma_1$ and $\Sigma_2$, which is used in the EM algorithm.

Figure 3a displays the FSR and the selection frequency as a function of the separation $\epsilon$. The conclusion is qualitatively the same as in the previous case, but with larger FSR values for a weak separation. Overall, it shows that the plug-in procedure is anticonservative and that the bootstrap corrections are able to recover an FSR and a selection frequency close to the one of the oracle. However, for a weak separation, namely $\epsilon = 1$, the parametric bootstrap correction is not enough and the latter procedure still overshoots the nominal level $\alpha$. Indeed, in our simulations, it appears that $P_{\theta}$ is often a distribution that is more favorable than $P_{\theta^*}$ from a statistical point of view (for instance, with more separated clusters). These conclusions also hold for varying sample size $n$, see Figure 3b.
<figure><img src="image_7.png" /><figcaption>Figure 6: Same as Figure 5 with dimension $Q = 3$, $d = 4$, $\epsilon = 2$, for an unconstrained model.</figcaption></figure>
<figure><img src="image_5.png" /><figcaption>FIGURE 3 Unknown covariance matrices: FSR (left panel) and selection frequency (right panel) as a function of: (a) the mean separation; (b) the sample size *n*; and (c) the nominal level *α*. Diagonal covariance matrix setting with *Q* = 2, *d* = 2. Default settings are: *n* = 200, *α* = 0.1, *ε* = √2.</figcaption></figure>

Figure 3c displays the FSR and the selection frequency for varying nominal level α, with ε = √2 and n = 200. The plug-in is still anticonservative, while the bootstrap procedures have an FSR that is close to α uniformly on the considered α range. Moreover, we note that for all our procedures (including the plug-in), the gap between the FSR and the nominal level is roughly constant with α: this illustrates the adaptive aspect of our procedures. This is in contrast with the baseline procedure, for which this gap highly depends on α, and which may be either anticonservative or suboptimal depending on the α value.

5.1.3 | Larger dimension

We now increase the dimension to *d* = 10. As in the previous setting, we consider a two-component mixture with balanced proportions and covariance matrices Σ<sub>q</sub> = I<sub>d</sub>, for any *q*, the diagonal structure being used as a constraint in the estimation. The mean separation ||μ<sub>1</sub> - μ<sub>2</sub>|| = ε is set to √2. Figure 4 displays the FSR and the selection frequency for varying α. In that case, the parameter estimation is more challenging, because the maximum posterior probability for any point tends to be over-estimated. Hence, taking *n* = 200 is not enough to obtain
<figure><img src="image_8.png" /><figcaption>Figure 7: Comparison of the clustering result with t-mixture modelling with ground truth labels on the WCBD dataset, restricted to the variables *radius* and *texture*, with and without selection. With the parametric bootstrap procedure applied at $\alpha = 5\%$, the FCR w.r.t. the ground truth labels is of 3% versus 14% without selection.</figcaption></figure>

(this information is not used in the clustering). One can see that the Student approximation is fairly good for each of the groups, and there is some overlap between them. Figure 7 (panel (b)) displays the MAP clustering result for the t-mixture model without any selection. The FCR is computed with respect to the ground truth labels and amounts to 14%. Finally, Figure 7 (panel (c)) provides the result of our parametric bootstrap procedure with nominal level $\alpha = 5\%$. The procedure does not classify points that are at the intersection of the clusters and the FCR equals 3%, which is below the targeted level.

Figure 8 (panel (a)) provides a comparison of all our procedures and the baseline for varying nominal level $\alpha$. The new procedures yield quite similar results. In particular, we see that the false clustering rate is always close to the nominal level, which illustrates once again the adaptive behavior of our procedures. By contrast, the false clustering rate of the baseline procedure is significantly below the nominal level and selects far less data points than our procedures, especially for large values of $\alpha$. In addition, Figure 8 (panel (b)) displays results for the dataset where seven variables ($d = 7$) are analyzed and the conclusions are qualitatively the same.

# 6 Conclusion and discussion

We have presented new data-driven methods providing both clustering and selection that ensure an FCR control guarantee in a mixture model. The plug-in approach was shown to be theoretically valid both when the parameter estimation is accurate and the sample size is large enough. When this is not necessarily the case, we proposed two second-order bootstrap corrections that have been shown to increase the FCR control ability on numerical experiments. Finally, applying our unsupervised methods to a supervised data set, our approach has been qualitatively validated by considering the attached labels as revealing the true clusters: our approach produces a selection/clustering close to what one would have done by observing the labels with some misclassification error range.

We would like to underline that the cluster number $Q$ is assumed to be fixed and known throughout the study. In practice, it can be fitted from the data by using the standard AIC or BIC criteria. In addition, if several values of $Q$ make sense from a practical viewpoint, we recommend to provide to the practitioner the collection of the corresponding outputs.

Concerning the pure task of controlling the FCR in the mixture model, our methods provide a correct FCR control in some area of the parameter space, leaving other less favorable
<figure><img src="image_6.png" /><figcaption>FIGURE 4 Large dimension: FSR (left panel) and selection frequency (right panel) as a function of the nominal level $\alpha$. Diagonal covariances setting with $Q = 2$, $d = 10$, $\epsilon = \sqrt{2}$.</figcaption></figure>

<figure><img src="image_7.png" /><figcaption>FIGURE 5 5-component mixture model: FSR (left panel) and selection frequency (right panel) as a function of the nominal level $\alpha$. Diagonal covariances setting with $Q = 5$, $d = 2$, $\epsilon = 2$.</figcaption></figure>

satisfactory performances. For $n = 1000$, however, the conclusions are qualitatively similar to dimension $d = 2$, which becomes clear by comparing Figures 4b and 3c.

### 5.1.4 | Five-component mixture

We next increase the number of classes to $Q = 5$, and set $\mu_1 = (0, 0)$, $\mu_2 = (\epsilon, 0)$, $\mu_3 = (0, \epsilon)$, $\mu_4 = (-\epsilon, 0)$, $\mu_5 = (0, -\epsilon)$ with $\epsilon = 2$. The mean separation is chosen so that the selection frequency of the oracle rule is approximately the same as in the previous cases. The mixture proportions are $\pi_q = 1/Q$, for any $q$ and the covariance matrices are $\Sigma_q = I_d$, for any $q$, the diagonal structure being used as a constraint in the estimation. Figure 5 displays the FSR and the selection frequency for varying $\alpha$. While the problem is more challenging in itself, the conclusions are qualitatively the same as in the previous settings.

## 5.2 | Real dataset

We consider the Vertebral Column dataset from the UCI ML repository. The data consists of $d = 6$ biomechanical features measured for $n = 310$ orthopaedic patients. Each attribute is derived from the shape and orientation of the pelvis and lumbar spine. Patients are classified
<figure><img src="image_9.png" /><figcaption>(a) Variables: radius, texture</figcaption></figure>

<figure><img src="image_10.png" /></figure>

<figure><img src="image_11.png" /></figure>

<figure><img src="image_12.png" /><figcaption>(b) Variables: radius, texture, smoothness, compactness, concavity, symmetry, fractal dimension</figcaption></figure>

Figure 8: FCR and selection frequency as a function of the nominal level $\alpha$.
<figure><img src="image_8.png" /><figcaption>**FIGURE 6** Relationship between each pair of features in the Vertebral Column dataset: plot of $X_i$ as a function of $X_j$ for each $i \neq j = 1, \dots, 6$. Data points are colored according to which class they belong to: “normal” (green), “disk hernia” (purple), and “spondilolysis” (yellow).</figcaption></figure>

into three classes: “normal” (100 patients), “disk hernia” (60 patients), or “spondilolysis” (150 patients). Figure 6 displays all pairs of variables ($X_i, X_j$) with $i \neq j$ for all six variables.

We choose to model the data as a mixture of Student's t-distributions as proposed in Peel and McLachlan (2000). Student's t-mixtures are appropriate for data containing observations with longer than normal tails or atypical observations leading to overlapping clusters. No constraints are imposed on the parameters, including the covariance matrices, for which no structural assumptions are made. The t-mixture is fit via the EM algorithm for t-mixtures (Peel & McLachlan, 2000) provided by the Python package `studenttmixture` (Parkinson, 2018). While the number of ground truth classes is $Q = 3$, the EM solution for 3 components leads to a cluster that is almost empty. Furthermore, one cluster essentially corresponds to the “spondilolysis” class, whereas members of the two classes “normal” and “disk hernia” are grouped into a unique cluster (which is clear from Figure 6, where the “normal” and “disk hernia” classes are hardly distinguishable). This leads us to run again the procedure for a mixture with only $Q = 2$ components.

For $Q = 2$ Table 1 reports the sample FSR and the selection frequency for varying $\alpha$, where the FSR is computed with respect to the ground truth classes, by merging the classes “normal” and “disk hernia” into a single class. However, given that the assumption of a Student's t-mixture does not necessarily hold for these ground truth groups, the sample FSR value computed here should be interpreted loosely, and a precise control of this value is not expected. Here the ground truth groups are well separated and all procedures display a FSR value below the nominal level $\alpha$. In terms of the selection frequency, the parametric bootstrap is close to the plug-in, as in the simulations, whereas the baseline procedure (which is necessarily more conservative than the plug-in) selects less data samples. In contrast to the other procedures, the nonparametric bootstrap displays adaptation to the level, selecting few examples when $\alpha$ is small, and selecting all examples when $\alpha$ is large.
parameter configurations with a slight inflation in the FCR level. This phenomenon is well known for FDR control in the two group mixture multiple testing model (Sun and Cai, 2007; Roquain and Verzelen, 2019), and facing a similar problem in our framework is well expected. On the one hand, in some cases, this problem can certainly be solved by improving on parameter estimation: here the EM algorithm seems to over-estimate the extreme posterior probabilities, which makes the plug-in procedure too anti-conservative. On the other hand, it could be hopeless to expect a robust FCR control uniformly valid over all configurations, while being optimal in the favorable cases. To illustrate that point, we refer to the work Roquain and Verzelen (2019) that shows that such a procedure does not exist in the FDR controlling case, when the null distribution is Gaussian with an unknown scaling parameter (which is a framework sharing similarities with the one considered here). Investigating such a “lower bound” result in the current setting would provide better guidelines for the practitioner and is therefore an interesting direction for future research. In addition, in these unfavorable cases, adding labeled samples and considering a semi-supervised framework can be an appropriate alternative for practical use. This new sample is likely to considerably improve the inference. Studying the FCR control in that setting is another promising avenue.

## Acknowledgments

This work has been supported by ANR-16-CE40-0019 (SansSouci), ANR-17-CE40-0001 (BA-SICS) and by the GDR ISIS through the “projets exploratoires” program (project TASTY). A. Marandon has been supported by a grant from Région Île-de-France (“DIM Math Innov”). We would like to thank Gilles Blanchard and Stéphane Robin for interesting discussions. We are also grateful to Eddie Aamari and Nhat Ho for their help for proving Lemma 17.

## References

* Kweku Abraham, Ismael Castillo, and Etienne Roquain. Empirical bayes cumulative $\ell$-value multiple testing procedure for sparse sequences. arXiv preprint arXiv:2102.00929, 2021.

* Anastasios Nikolas Angelopoulos, Stephen Bates, Michael I. Jordan, and Jitendra Malik. Uncertainty sets for image classifiers using conformal prediction. In 9th International Conference on Learning Representations, ICLR 2021, Virtual Event, Austria, May 3-7, 2021. OpenReview.net, 2021. URL https://openreview.net/forum?id=eNdiU_DbM9.

* Sanjeev Arora and Ravi Kannan. Learning mixtures of separated nonspherical Gaussians. The Annals of Applied Probability, 15(1A):69–92, 2005. doi: 10.1214/105051604000000512. URL https://doi.org/10.1214/105051604000000512.

* David Arthur and Sergei Vassilvitskii. k-means++: The advantages of careful seeding. Technical report, Stanford, 2006.

* Martin Azizyan, Aarti Singh, and Larry Wasserman. Minimax theory for high-dimensional gaussian mixtures with sparse mean separation. In Proceedings of the 26th International Conference on Neural Information Processing Systems - Volume 2, NIPS'13, page 2139–2147, Red Hook, NY, USA, 2013. Curran Associates Inc.

* Sivaraman Balakrishnan, Martin J. Wainwright, and Bin Yu. Statistical guarantees for the EM algorithm: From population to sample-based analysis. The Annals of Statistics, 45(1): 77–120, 2017. doi: 10.1214/16-AOS1435. URL https://doi.org/10.1214/16-AOS1435.
<table><caption>TABLE 1 Vertebral Column dataset: sample FSR and selection frequency of each method for varying $\alpha$.</caption><thead><tr><th>α</th><th>Baseline</th><th>Plug-in</th><th>Parametric bootstrap</th><th>Nonparametric bootstrap</th></tr></thead><tbody><tr><td>Sample FSR</td><td></td><td></td><td></td><td></td></tr><tr><td>1%</td><td>0.009</td><td>0.017</td><td>0.015</td><td>0</td></tr><tr><td>5%</td><td>0.015</td><td>0.032</td><td>0.032</td><td>0.009</td></tr><tr><td>10%</td><td>0.018</td><td>0.032</td><td>0.032</td><td>0.032</td></tr><tr><td>Selection frequency</td><td></td><td></td><td></td><td></td></tr><tr><td>1%</td><td>0.739</td><td>0.926</td><td>0.861</td><td>0</td></tr><tr><td>5%</td><td>0.848</td><td>1.000</td><td>1.000</td><td>0.745</td></tr><tr><td>10%</td><td>0.919</td><td>1.000</td><td>1.000</td><td>1.000</td></tr></tbody></table>

Note: The number of components *Q* is set to 2 and the sample FSR is computed with respect to the ground truth classes "normal/disk hernia" (merged into a single class), and "spondilolysthesis".

## 6 | CONCLUSION AND DISCUSSION

We have presented new data-driven methods providing both clustering and selection that ensure an FSR control guarantee in a mixture model. The plug-in approach was shown to be theoretically valid both when the parameter estimation is accurate and the sample size is large enough. When this is not necessarily the case, we proposed two second-order bootstrap corrections that have been shown to increase the FSR control ability on numerical experiments. Finally, applying our unsupervised methods to a supervised dataset, our approach has been qualitatively validated by considering the attached labels as revealing the true clusters.

We underline that the cluster number *Q* is assumed to be fixed and known throughout the study. In practice, the cluster number can be unknown, in which case it is estimated from the data. In that case, the interpretation of the FSR control is still meaningful, in the sense that it provides error control within a *specified* model. To obtain a more general guarantee, a solution could be to consider a definition of the clustering error that is able to compare partitions with different levels of granularity (and hence cluster numbers). There are several such measures in use in clustering (not restricted to a selection), such as the Adjusted Rand Index (Hubert & Arabie, 1985), among others (Rosenberg & Hirschberg, 2007; Nguyen et al., 2009). Building a selection rule that controls the clustering error in restriction to a part of the sample is expected to be much more challenging for such a criterion, but represents an interesting research direction.

Concerning the pure task of controlling the FSR in the mixture model, our methods provide a correct FSR control in some region of the parameter space, leaving other less favorable parameter configurations with a slight inflation in the FSR level. This phenomenon is well known for FDR control in the two-component mixture multiple testing model (Roquain & Verzelen, 2022; Sun & Cai, 2007), and facing a similar problem in our framework is well expected. On the one hand, in some cases, this problem can certainly be solved by improving on parameter estimation: here the EM algorithm seems to over-estimate the extreme posterior probabilities, which makes the plug-in procedure too anticonservative. On the other hand, it could be hopeless to expect a robust FSR control uniformly valid over all configurations, while being optimal in the favorable cases. To illustrate that point, we refer to Roquain and Verzelen (2022) who show that such a procedure does not exist in the FDR controlling case, when the null distribution is Gaussian with
Peter L. Bartlett and Marten H. Wegkamp. Classification with a reject option using a hinge loss. Journal of Machine Learning Research, 9(59):1823–1840, 2008. URL http://jmlr.org/papers/v9/bartlett08a.html.

Yoav Benjamini and Yosef Hochberg. Controlling the false discovery rate: a practical and powerful approach to multiple testing. J. Roy. Statist. Soc. Ser. B, 57(1):289–300, 1995. ISSN 0035-9246.

Peter J. Bickel and Aiyou Chen. A nonparametric view of network models and newman-girvan and other modularities. Proceedings of the National Academy of Sciences, 106(50):21068–21073, 2009. ISSN 0027-8424. doi: 10.1073/pnas.0907096106. URL https://www.pnas.org/content/106/50/21068.

T. Cai, Wenguang Sun, and Weinan Wang. Covariate-assisted ranking and screening for large-scale two-sample inference. Journal of the Royal Statistical Society: Series B (Statistical Methodology), 81(2):187–234, 2019.

Jiahua Chen. Consistency of the MLE under Mixture Models. Statistical Science, 32(1):47–63, 2017. doi: 10.1214/16-STS578. URL https://doi.org/10.1214/16-STS578.

C. Chow. On optimum recognition error and reject tradeoff. IEEE Transactions on Information Theory, 16(1):41–46, 1970. doi: 10.1109/TIT.1970.1054406.

A. P. Dempster, N. M. Laird, and D. B. Rubin. Maximum likelihood from incomplete data via the em algorithm. Journal of the Royal Statistical Society. Series B (Methodological), 39(1):1–38, 1977.

Christophe Denis and Mohamed Hebiri. Consistency of plug-in confidence sets for classification in semi-supervised learning. Journal of Nonparametric Statistics, 32(1):42–72, 2020. doi: 10.1080/10485252.2019.1689241. URL https://doi.org/10.1080/10485252.2019.1689241.

Bradley Efron, Robert Tibshirani, John D Storey, and Virginia Tusher. Empirical bayes analysis of a microarray experiment. Journal of the American Statistical Association, 96 (456):1151–1160, 2001. doi: 10.1198/016214501753382129.

Radu Herbei and Marten H. Wegkamp. Classification with reject option. The Canadian Journal of Statistics / La Revue Canadienne de Statistique, 34(4):709–721, 2006. ISSN 03195724. URL http://www.jstor.org/stable/20445230.

Nhat Ho and XuanLong Nguyen. On strong identifiability and convergence rates of parameter estimation in finite mixtures. Electronic Journal of Statistics, 10(1):271–307, 2016. doi: 10.1214/16-EJS1105. URL https://doi.org/10.1214/16-EJS1105.

Jing Lei and Alessandro Rinaldo. Consistency of spectral clustering in stochastic block models. The Annals of Statistics, 43(1):215–237, 2015. doi: 10.1214/14-AOS1274. URL https://doi.org/10.1214/14-AOS1274.

Yu Lu and Harrison H Zhou. Statistical and computational guarantees of lloyd's algorithm and its variants. arXiv preprint arXiv:1612.02099, 2016.

Tristan Mary-Huard, Vittorio Perduca, Gilles Blanchard, and Martin-Magniette Marie-Laure. Error rate control for classification rules in multiclass mixture models, 2021.
an unknown scaling parameter (which is a framework sharing similarities with the one consid-
ered here). Investigating such a “lower bound” result in the current setting would provide better
guidelines for the practitioner and is therefore an interesting direction for future research.

Another avenue would be to assess the theoretical performances of the bootstrap procedures.
While a consistency result for FSR/Power similar to Theorem 1 might be possible to obtain from
classical theory (see, e.g., Politis et al., 1999)—but would require a fully devoted work—it seems
even more challenging to theoretically justify why the corresponding remainder terms are smaller
than those obtained for the plug-in procedure, hence to fully support our experimental findings.

Next, the continuity of the map $t \mapsto \mathbb{P}_{\theta^{-}}(T(X, \theta) < t)$ (see Assumption 1) is a key property to apply our methodology (see in particular Lemma 6). For instance, Assumption 1 is not true for discrete data (e.g., when considering a Poisson model). However, one idea is to circumvent this difficulty by using a randomization technique, as for the original Neyman-Pearson test. Investigating this point further is left for future research.

Finally, to attempt to overcome problems of model misspecification related to strong dis-
tributional assumptions of the mixture model, one may explore the extension of our approach
to semiparametric mixture models. In recent years, several new semiparametric mixtures have
been proposed in the literature, see Xiang et al. (2019) for an overview. For univariate obser-
vations, there are mixtures over a location parameter (Bordes et al., 2006), and others that use
shape constraints, namely log concavity, on the densities of the mixture components (Chang &
Walther, 2007). For multivariate data, there also exist some approaches, namely mixtures where
conditional independence of the covariates is assumed (Benaglia et al., 2009; Hall & Zhou, 2003).
In most of these models, identifiability is an issue and the construction of consistent and com-
putationally efficient estimators is challenging. Most often, inference is performed using some
EM-type algorithm and kernel-based estimators for the component densities. Adapting our clus-
tering approach to such approaches is challenging both from a practical point of view (in regard
of the dimension of the observations and the constraints on the component densities) and from a
theoretical point of view (given the limited existing literature proving consistency).

ACKNOWLEDGMENTS

This work has been supported by ANR-16-CE40-0019 (SansSouci), ANR-17-CE40-0001 (BASICS) and by the GDR ISIS through the “projets exploratoires” program (project TASTY). A. Marandon has been supported by a grant from Région Île-de-France (“DIM Math Innov”). We would like to thank Gilles Blanchard and Stéphane Robin for interesting discussions. We are also grateful to Eddie Aamari and Nhat Ho for their help for proving Lemma 16. Finally, we would like to thank anonymous referees and an associate editor for helpful comments.

ENDNOTE

ORCID

Ariane Marandon ID https://orcid.org/0000-0002-9114-2527

Etienne Roquain ID https://orcid.org/0000-0003-0214-5185

REFERENCES

Abbe, E. (2018). Community detection and stochastic block models: Recent developments. *Journal of Machine Learning Research*, *18*(177), 1–86.

Abraham, K., Castillo, I., & Roquain, É. (2022). Empirical bayes cumulative *ℓ*-value multiple testing procedure for sparse sequences. *Electronic Journal of Statistics*, *16*(1), 2033–2081.

<footnote n="1">¹Here, $\lambda_1(\Sigma)$ (resp. $\lambda_d(\Sigma)$) denotes the smallest (resp. largest) eigenvalue of $\Sigma$.</footnote>
Pascal. Massart. *Concentration Inequalities and Model Selection Ecole d'Eté de Probabilités de Saint-Flour XXXIII - 2003*. École d'Eté de Probabilités de Saint-Flour, 1896. Springer Berlin Heidelberg, Berlin, Heidelberg, 1st ed. 2007. edition, 2007. ISBN 1-280-85333-6.

Volodymyr Melnykov. On the distribution of posterior probabilities in finite mixture models with application in clustering. *Journal of Multivariate Analysis*, 122:175–189, 2013.

Amir Najafi, Seyed Abolfazl Motahari, and Hamid R. Rabiee. Reliable clustering of Bernoulli mixture models. *Bernoulli*, 26(2):1535 – 1559, 2020. doi: 10.3150/19-BEJ1173. URL https://doi.org/10.3150/19-BEJ1173.

D. Peel and G. J. McLachlan. Robust mixture modelling using the t distribution. *Statistics and Computing*, 10(4):339–348, 2000. doi: 10.1023/A:1008981510081. URL https://doi.org/10.1023/A:1008981510081.

Tabea Rebafka, Etienne Roquain, and Fanny Villers. Graph inference with clustering and false discovery rate control, 2019.

Etienne Roquain and Nicolas Verzelen. False discovery rate control with unknown null distribution: is it possible to mimic the oracle? *arXiv preprint arXiv:1912.03109*, 2019.

John D Storey. The positive false discovery rate: a bayesian interpretation and the q-value. *The Annals of Statistics*, 31(6):2013–2035, 2003.

Wenguang Sun and T. Tony Cai. Oracle and adaptive compound decision rules for false discovery rate control. *Journal of the American Statistical Association*, 102(479):901–912, 2007.

Vladimir Vovk, Alexander Gammerman, and Glenn Shafer. *Algorithmic learning in a random world*. Springer Science & Business Media, 2005.

Volodya Vovk, Alexander Gammerman, and Craig Saunders. Machine-learning applications of algorithmic randomness. In *Proceedings of the Sixteenth International Conference on Machine Learning*, ICML '99, page 444–453, 1999.

Marten Wegkamp and Ming Yuan. Support vector machines with a reject option. *Bernoulli*, 17(4):1368–1385, 2011.

# A Proof of our main result

## A.1 An optimal procedure

We consider in this section the procedure that serves as an optimal procedure in our theory. For $t \in [0, 1]$, let $C_t^* = (\hat{Z}^*, S_t^*)$ be the procedure using the Bayes clustering $\hat{Z}^*$ (6) and the selection rule $S_t^* = \{i \in \{1, \dots, n\} : T_i^* < t\}$. Let us consider the map $t \in [0, 1] \mapsto \text{mFCR}(C_t^*)$ and note that $\text{mFCR}(C_t^*) = \text{mFCR}_t^*$ as defined by (12). Lemma 11 below provides the key properties for this function.

**Definition 3.** The optimal procedure at level $\alpha$ is defined by $C_{t^*(\alpha)}^*$ where $t^*(\alpha)$ is defined by (13).
Angelopoulos, A. N., & Bates, S. (2021). A gentle introduction to conformal prediction and distribution-free uncertainty quantification. arXiv preprint arXiv:2107.07511.

Angelopoulos, A. N., Bates, S., Candès, E. J., Jordan, M. I., & Lei, L. (2025). Learn then test: Calibrating predictive algorithms to achieve risk control. *The Annals of Applied Statistics*, 19(2), 1641–1662.

Arora, S., & Kannan, R. (2005). Learning mixtures of separated nonspherical Gaussians. *Annals of Applied Probability*, 15(1A), 69–92.

Arthur, D., & Vassilvitskii, S. (2007). *K-means++: The advantages of careful seeding*. In *Proceedings of the Eighteenth Annual ACM-SIAM Symposium on Discrete Algorithms* (pp. 1027–1035). Society for Industrial and Applied Mathematics.

Azizyan, M., Singh, A., & Wasserman, L. (2013). *Minimax theory for high-dimensional gaussian mixtures with sparse mean separation*. In *Proceedings of the 26th international conference on neural information processing systems*—*volume 2*, NIPS'13 (pp. 2139–2147). Curran Associates Inc., USA.

Balakrishnan, S., Wainwright, M. J., & Yu, B. (2017). Statistical guarantees for the EM algorithm: From population to sample-based analysis. *Annals of Statistics*, 45(1), 77–120.

Bao, H., Shimada, T., Xu, L., Sato, I., & Sugiyama, M. (2020). Pairwise supervision can provably elicit a decision boundary. arXiv preprint arXiv:2006.06207.

Baraud, Y. (2016). Bounding the expectation of the supremum of an empirical process over a (weak) vc-major class. *Electronic Journal of Statistics*, 10(2), 1709–1728.

Bartlett, P. L., & Wegkamp, M. H. (2008). Classification with a reject option using a hinge loss. *Journal of Machine Learning Research*, 9(59), 1823–1840.

Benaglia, T., Chauveau, D., & Hunter, D. R. (2009). An EM-like algorithm for semi- and non-parametric estimation in multivariate mixtures. *Journal of Computational and Graphical Statistics*, 18(2), 505–526.

Benjamini, Y., & Hochberg, Y. (1995). Controlling the false discovery rate: A practical and powerful approach to multiple testing. *Journal of the Royal Statistical Society. Series B, Statistical Methodology*, 57(1), 289–300.

Bordes, L., Mottelet, S., & Vandekerkhove, P. (2006). Semiparametric estimation of a two-component mixture model. *Annals of Statistics*, 34(3), 1204–1232.

Cai, T., Sun, W., & Wang, W. (2019). Covariate-assisted ranking and screening for large-scale two-sample inference. *Journal of the Royal Statistical Society: Series B*, 81(2), 187–234.

Cao, Y., Cai, T., Feng, L., Gu, L., Gu, J., An, B., Niu, G., & Sugiyama, M. (2022). Generalizing consistent multi-class classification with rejection to be compatible with arbitrary losses. *Advances in Neural Information Processing Systems*, 35, 521–534.

Chang, G. T., & Walther, G. (2007). Clustering with mixtures of log-concave distributions. *Computational Statistics & Data Analysis*, 51(12), 6242–6251.

Chen, J. (2017). Consistency of the MLE under mixture models. *Statistical Science*, 32(1), 47–63.

Chow, C. (1970). On optimum recognition error and reject tradeoff. *IEEE Transactions on Information Theory*, 16(1), 41–46.

Chretien, S., Dombry, C., & Faivre, A. (2019). The Guedon-Vershynin semi-definite programming approach to low dimensional embedding for unsupervised clustering. *Frontiers in Applied Mathematics and Statistics*, 5, 41.

Cortes, C., DeSalvo, G., & Mohri, M. (2016). *Learning with rejection*. In *Algorithmic learning theory: 27th international conference* (pp. 67–82). Springer.

Dempster, A. P., Laird, N. M., & Rubin, D. B. (1977). Maximum likelihood from incomplete data via the em algorithm. *Journal of the Royal Statistical Society. Series B, Statistical Methodology*, 39(1), 1–38.

Denis, C., & Hebiri, M. (2020). Consistency of plug-in confidence sets for classification in semi-supervised learning. *Journal of Nonparametric Statistics*, 32(1), 42–72.

Eckardt, J.-N., Röllig, C., Metzeler, K., Heisig, P., Stasik, S., Georgi, J.-A., Kroschinsky, F., Stölzel, F., Platzbecker, U., Spiekermann, K., Krug, U., Braess, J., Görlich, D., Sauerland, C., Woermann, B., Herold, T., Hiddemann, W., Müller-Tidow, C., Serve, H., ... Middeke, J. M. (2023). Unsupervised meta-clustering identifies risk clusters in acute myeloid leukemia based on clinical and genetic profiles. *Communications Medicine*, 3(1), 68.

Efron, B., Tibshirani, R., Storey, J. D., & Tusher, V. (2001). Empirical bayes analysis of a microarray experiment. *Journal of the American Statistical Association*, 96(456), 1151–1160.

Fahad, A., Alshatri, N., Tari, Z., Alamri, A., Khalil, I., Zomaya, A. Y., Foufou, S., & Bouras, A. (2014). A survey of clustering algorithms for big data: Taxonomy and empirical analysis. *IEEE Transactions on Emerging Topics in Computing*, 2(3), 267–279.
Note that the optimal procedure is not the same as the oracle procedure defined in Section 3.1, although these two procedures are expected to behave roughly in the same way (at least for a large $n$).

Under Assumption 1, Lemma 11 entails that, for $\alpha > \alpha_c$, mFCR($C_{t^*(\alpha)}^*$) $\le \alpha$. Hence, $C_{t^*(\alpha)}^*$ controls the mFCR at level $\alpha$. In addition, it is optimal in the following sense: any other mFCR controlling procedure should select less items than $C_{t^*(\alpha)}^*$.

**Lemma 9 (Optimality of $C_{t^*(\alpha)}^*$).** Let Assumption 1 be true and choose $\alpha \in (\alpha_c, \bar{\alpha}]$. Then the oracle procedure $C_{t^*(\alpha)}^* = (\hat{\mathbf{Z}}^*, S_{t^*(\alpha)}^*)$ satisfies the following:

$$ (i) \text{ mFCR}(C_{t^*(\alpha)}^*) = \alpha; $$

$$ (ii) \text{ for any procedure } \mathcal{C} = (\hat{\mathbf{Z}}, S) \text{ such that mFCR}(\mathcal{C}) \le \alpha, \text{ we have } \mathbb{E}_{\theta^*}(|S|) \le \mathbb{E}_{\theta^*}(|S_{t^*(\alpha)}^*|). $$

## A.2 Preliminary steps for proving Theorem 6

To keep the main proof concise, we need to define several additional notation. Let for $t \in [0, 1]$ and $\theta \in \Theta$ (recall (5))

$$ \hat{\mathbf{L}}_0(\theta, t) = \frac{1}{n} \sum_{i=1}^n T(X_i, \theta) \mathbb{1}_{T(X_i, \theta) < t}; $$

$$ \hat{\mathbf{L}}_1(\theta, t) = \frac{1}{n} \sum_{i=1}^n \mathbb{1}_{T(X_i, \theta) < t}. $$

Denote $\hat{\mathbf{L}} = \hat{\mathbf{L}}_0/\hat{\mathbf{L}}_1$, $\mathbf{L}_0 = \mathbb{E}_{\theta^*} \hat{\mathbf{L}}_0$, $\mathbf{L}_1 = \mathbb{E}_{\theta^*} \hat{\mathbf{L}}_1$, $\mathbf{L} = \mathbf{L}_0/\mathbf{L}_1$ (with the convention $0/0 = 0$). Note that for any $\alpha > \alpha_c$, the mFCR of the optimal procedure $C_{t^*(\alpha)}^*$ defined in Section A.1 is given by mFCR($C_{t^*(\alpha)}^*$) = $\mathbf{L}(\theta^*, t^*(\alpha)) = \alpha$.

Also, we denote from now on $\ell_{i,q}^* = \mathbb{P}_{\theta^*}(Z_i = q|X_i)$ for short and introduce for any parameter $\theta \in \Theta$ (recall (4) and (5))

$$ \bar{q}(X_i, \theta) \in \underset{q \in \{1, \dots, Q\}}{\operatorname{argmax}} \ell_q(X_i, \theta), \quad 1 \le i \le n; \qquad (23) $$

$$ U(X_i, \theta) = 1 - \ell_{i, \bar{q}}^*(X_i, \theta), \quad 1 \le i \le n; \qquad (24) $$

$$ \hat{\mathbf{M}}_0(\theta, t) = \frac{1}{n} \sum_{i=1}^n U(X_i, \theta) \mathbb{1}_{T(X_i, \theta) < t}, \quad t \in [0, 1], \qquad (25) $$

Note that $\hat{\mathbf{M}}_0(\theta^*, t) = \hat{\mathbf{L}}_0(\theta^*, t)$ but in general $\hat{\mathbf{M}}_0(\theta, t)$ is different from $\hat{\mathbf{L}}_0(\theta, t)$. We denote $\hat{\mathbf{M}} = \hat{\mathbf{M}}_0/\hat{\mathbf{L}}_1$, $\mathbf{M}_0 = \mathbb{E}_{\theta^*} \hat{\mathbf{M}}_0$ and $\mathbf{M} = \mathbf{M}_0/\mathbf{L}_1$ (with the convention $0/0 = 0$).

When $\alpha \in (\alpha_c, \bar{\alpha}]$ (recall (14) and (15)), we also let

$$ s^* = s^*(\alpha, \theta^*) = n^{-1} \mathbb{E}_{\theta^*} \left( |S_{t^*}^{*\left(\frac{\alpha+\alpha_c}{2}\right)}| \right) = \mathbf{L}_1(\theta^*, t^*\left(\frac{\alpha+\alpha_c}{2}\right)) > 0. \qquad (26) $$

We easily see that the latter is positive: if it was zero then $S_{t^*\left(\frac{\alpha+\alpha_c}{2}\right)}^*$ would be empty which would entails that mFCR($C_{t^*\left(\frac{\alpha+\alpha_c}{2}\right)}^*$) is zero. This is excluded by definition (14) of $\alpha_c$ because $(\alpha + \alpha_c)/2 > \alpha_c$.

Also, we are going to extensively use the event

$$ \Omega_{\epsilon} = \left\{ \min_{\sigma \in [Q]} \|\hat{\theta}^{\sigma} - \theta^{*}\|_{2} < \epsilon \right\}. $$
Geifman, Y., & El-Yaniv, R. (2017). Selective classification for deep neural networks. In Proceedings of the 31st International Conference on Neural Information Processing Systems, NIPS'17. (pp. 4885-4894). Curran Associates Inc. USA.

Giraud, C., & Verzelen, N. (2018). Partial recovery bounds for clustering with the relaxed K-means. *Mathematical Statistics and Learning, **1**(3), 317-374.

Grün, B. (2019). *Model-based clustering*. In *Handbook of mixture analysis* (pp. 157-192). Chapman and Hall/CRC.

Hall, P., & Zhou, X.-H. (2003). Nonparametric estimation of component distributions in a multivariate mixture. *Annals of Statistics, **31*** (1), 201-224.

Herbei, R., & Wegkamp, M. H. (2006). Classification with reject option. *Canadian Journal of Statistics, **34*** (4), 709-721.

Ho, N., & Nguyen, X. (2016). On strong identifiability and convergence rates of parameter estimation in finite mixtures. *Electronic Journal of Statistics, **10*** (1), 271-307.

Hubert, L., & Arabie, P. (1985). Comparing partitions. *Journal of Classification, **2*** (1), 193-218.

Lei, J. (2014). Classification with confidence. *Biometrika, **101*** (4), 755-769.

Lei, J., & Rinaldo, A. (2015). Consistency of spectral clustering in stochastic block models. *Annals of Statistics, **43*** (1), 215-237.

Lu, Y., & Zhou, H. H. (2016). Statistical and computational guarantees of lloyd's algorithm and its variants. *arXiv preprint arXiv:1612.02099*.

Mary-Huard, T., Perduca, V., Blanchard, G., & Marie-Laure, M.-M. (2021). Error rate control for classification rules in multiclass mixture models.

Massart, P. (2007). *Concentration inequalities and model selection Ecole d'Eté de Probabilités de saint-flour XXXIII - 2003 Tex2Doc conversion from DLF (18-08-2025)-Batch-III (1st ed.). Springer*.

Melnykov, V. (2013). On the distribution of posterior probabilities in finite mixture models with application in clustering. *Journal of Multivariate Analysis, **122***, 175-189.

Mohri, M., Rostamizadeh, A., & Talwalkar, A. (2012). *Foundations of machine learning*. The MIT Press.

Najafi, A., Motahari, S. A., & Rabiee, H. R. (2020). Reliable clustering of Bernoulli mixture models. *Bernoulli, **26*** (2), 1535-1559.

Ng, S.-K. (2019). *Mixture modelling for medical and health sciences*. Chapman and Hall/CRC.

Nguyen, X. V., Epps, J., & Bailey, J. (2009). Information theoretic measures for clusterings comparison: Is a correction for chance necessary? In A. P. Danyluk, L. Bottou, & M. L. Littman (Eds.), *Proceedings of the 26th Annual international conference on machine learning, ICML 2009 International Conference Proceeding Series* (Vol. 382, pp. 1073-1080). ACM. https://doi.org/10.1145/1553374.1553511

Ni, C., Charoenphakdee, N., Honda, J., & Sugiyama, M. (2019). On the calibration of multiclass classification with rejection. In H. M. Wallach, H. Larochelle, A. Beygelzimer, F. d'Alché-Buc, E. B. Fox, & R. Garnett (Eds.), *Advances in neural information processing systems 32: Annual conference on neural information processing systems 2019* (pp. 2582-2592). NeurIPS. https://proceedings.neurips.cc/paper/2019/hash/571d3a9420bfd9219f65b643d0003bf4-Abstract.html

O'Hagan, A., Murphy, T. B., Scrucca, L., & Gormley, I. C. (2019). Investigation of parameter uncertainty in clustering using a gaussian mixture model via jackknife, bootstrap and weighted likelihood bootstrap. *Computational Statistics, **34*** (4), 1779-1813.

Parkinson, J. (2018). Python package studenttmixture. https://github.com/jlparki/mix_T

Peel, D., & McLachlan, G. J. (2000). Robust mixture modelling using the t distribution. *Statistics and Computing, **10*** (4), 339-348.

Politis, D. N., Romano, J. P., Wolf, M., Politis, D. N., Romano, J. P., & Wolf, M. (1999). *Subsampling in the IID case*. Springer.

Rebafka, T., Roquain, É., & Villers, F. (2022). Powerful multiple testing of paired null hypotheses using a latent graph model. *Electronic Journal of Statistics, **16*** (1), 2796-2858.

Regev, O., & Vijayaraghavan, A. (2017). On learning mixtures of well-separated Gaussians. In C. Umans (Ed.), *58th IEEE annual symposium on foundations of computer science, FOCS 2017* (pp. 85-96). IEEE Computer Society. https://doi.org/10.1109/FOCS.2017.17

Romano, Y., Sesia, M., & Candès, E. J. (2020). Classification with valid and adaptive coverage. In H. Larochelle, M.' A. Ranzato, R. Hadsell, M.{-.}, F. Balcan, & H.{-.}. T. Lin (Eds.), *Advances in neural information
On this event, we fix any permutation $\sigma \in [Q]$ (possibly depending on $X$) such that $\|\hat{\theta}^\sigma - \theta^*\|_2 < \epsilon$. Now using Lemma 12, the plug-in selection rule can be rewritten as $\hat{S}_\alpha^{\text{PI}} = \{i \in \{1, \dots, n\} : \hat{T}_i < \hat{t}(\alpha)\}$ (denoted by $\hat{S}$ in the sequel for short), where

$$ \hat{t}(\alpha) = \sup\{t \in [0, 1] : \hat{\mathbf{L}}(\hat{\theta}, t) \le \alpha\}. \qquad (27) $$

With the above notation, we can upper bound what is inside the brackets of FCR($\hat{\mathcal{C}}^{\text{PI}}$) and mFCR($\hat{\mathcal{C}}^{\text{PI}}$) as follows.

**Lemma 10.** For the permutation $\sigma$ in $\Omega_\epsilon$ realizing $\|\hat{\theta}^\sigma - \theta^*\|_2 < \epsilon$, we have on the event $\Omega_\epsilon$ the following relations:

$$ |\hat{S}| = \hat{\mathbf{L}}_1(\hat{\theta}^\sigma, \hat{t}(\alpha)); $$

$$ \min_{\sigma' \in [Q]} \mathbb{E}_{\theta^*} \left( \varepsilon_{\hat{S}}(\sigma'(\hat{\mathbf{Z}}), \mathbf{Z}) \middle| \mathbf{X} \right) \le \hat{\mathbf{M}}_0(\hat{\theta}^\sigma, \hat{t}(\alpha)); $$

$$ \min_{\sigma' \in [Q]} \mathbb{E}_{\theta^*} \left( \frac{\varepsilon_{\hat{S}}(\sigma'(\hat{\mathbf{Z}}), \mathbf{Z})}{\max(|\hat{S}|, 1)} \middle| \mathbf{X} \right) \le \hat{\mathbf{M}}(\hat{\theta}^\sigma, \hat{t}(\alpha)). $$

Finally, we make use of the concentration of the empirical processes $\hat{\mathbf{L}}_0(\theta, t)$, $\hat{\mathbf{L}}_1(\theta, t)$, and $\hat{\mathbf{M}}_0(\theta, t)$, uniformly with respect to $\theta \in \mathcal{D}$ (where $\mathcal{D}$ is defined in Assumption 4). Thus, we define the following events, for $\delta, \epsilon > 0$ (recall $s^*$ defined by (26)):

$$ \Gamma_{0,\delta,t} = \left\{ \sup_{\theta \in \mathcal{D}} \left| \hat{\mathbf{L}}_0(\theta, t) - \mathbf{L}_0(\theta, t) \right| \le \delta \right\}; $$

$$ \Gamma_{1,\delta,t} = \left\{ \sup_{\theta \in \mathcal{D}} \left| \hat{\mathbf{L}}_1(\theta, t) - \mathbf{L}_1(\theta, t) \right| \le \delta \right\}; $$

$$ \Gamma_{\delta,t} = \left\{ \sup_{\substack{\theta \in \mathcal{D}, \\ \mathbf{L}_1(\theta,t) \ge s^*}} \left| \hat{\mathbf{L}}(\theta,t) - \mathbf{L}(\theta,t) \right| \le \delta \right\}; $$

$$ \Upsilon_{0,\delta,t} = \left\{ \sup_{\theta \in \mathcal{D}} \left| \hat{\mathbf{M}}_0(\theta, t) - \mathbf{M}_0(\theta, t) \right| \le \delta \right\}. $$

Note that the following holds:

$$ \Gamma_{0,\delta s^*/2,t} \cap \Gamma_{1,\delta s^*/2,t} \subset \Gamma_{\delta,t}. \qquad (28) $$

Indeed, on the event $\Gamma_{0,\delta s^*/2,t} \cap \Gamma_{1,\delta s^*/2,t}$, provided that $\mathbf{L}_1(\theta, t) \ge s^*$, we have

$$ \begin{aligned} & \left| \frac{\hat{\mathbf{L}}_0(\theta, t)}{\hat{\mathbf{L}}_1(\theta, t)} - \frac{\mathbf{L}_0(\theta, t)}{\mathbf{L}_1(\theta, t)} \right| \\ & \le \left| \frac{\mathbf{L}_0(\theta, t) - \hat{\mathbf{L}}_0(\theta, t)}{\mathbf{L}_1(\theta, t)} \right| + \hat{\mathbf{L}}_0(\theta, t) \left| \frac{1}{\hat{\mathbf{L}}_1(\theta, t)} - \frac{1}{\mathbf{L}_1(\theta, t)} \right| \\ & \le (\delta s^*/2)/s^* + (\delta s^*/2)/s^* = \delta, \end{aligned} $$

because $\hat{\mathbf{L}}_0(\theta, t) \le \hat{\mathbf{L}}_1(\theta, t)$. This proves the desired inclusion.

## A.3 Proof of Theorem 6

Let us now provide a proof for Theorem 6.
processing systems 33: Annual conference on neural information processing systems 2020. NeurIPS. https://proceedings.neurips.cc/paper/2020/hash/244edd7e85dc81602b7615cd705545f5-Abstract.html

Roquain, E., & Verzelen, N. (2022). False discovery rate control with unknown null distribution: Is it possible to mimic the oracle? *Annals of Statistics*, 50(2), 1095–1123.

Rosenberg, A., & Hirschberg, J. (2007). V-measure: A conditional entropy-based external cluster evaluation measure. In J. Eisner (Ed.), *EMNLP-CoNLL 2007, Proceedings of the 2007 joint conference on empirical methods in natural language processing and computational natural language learning* (pp. 410–420). ACL. https://aclanthology.org/D07-1043/

Sadinle, M., Lei, J., & Wasserman, L. (2019). Least ambiguous set-valued classifiers with bounded error levels. *Journal of the American Statistical Association*, 114(525), 223–234.

Shalev-Shwartz, S., & Ben-David, S. (2014). *Understanding machine learning: From theory to algorithms*. Cambridge University Press.

Stoitsas, K., Bahulikar, S., de Munter, L., de Jongh, M. A., Jansen, M. A., Jung, M. M., van Wingerden, M., & Van Deun, K. (2022). Clustering of trauma patients based on longitudinal data and the application of machine learning to predict recovery. *Scientific Reports*, 12(1), 16990.

Storey, J. D. (2003). The positive false discovery rate: A bayesian interpretation and the q-value. *Annals of Statistics*, 31(6), 2013–2035.

Sun, W., & Cai, T. T. (2007). Oracle and adaptive compound decision rules for false discovery rate control. *Journal of the American Statistical Association*, 102(479), 901–912.

Vempala, S., & Wang, G. (2004). A spectral algorithm for learning mixture models. *Journal of Computer and System Sciences*, 68(4), 841–860. Special Issue on FOCS 2002.

Vovk, V., Gammerman, A., & Shafer, G. (2005). *Algorithmic learning in a random world*. Springer Science & Business Media.

Wegkamp, M., & Yuan, M. (2011). Support vector machines with a reject option. *Bernoulli*, 17(4), 1368–1385.

Xiang, S., Yao, W., & Yang, G. (2019). An overview of semiparametric extensions of finite mixture models. *Statistical Science*, 34(3), 391–404.

**How to cite this article:** Marandon, A., Rebafka, T., Roquain, E., & Sokolovska, N. (2025). False selection rate control in mixture models. *Scandinavian Journal of Statistics*, 1–47. https://doi.org/10.1111/sjos.70017

## APPENDIX A. PROOF OF THEOREMS 1 AND 2

### A1 A general result

In this section, we establish a general result, from which Theorems 1 and 2 can be deduced. It provides nonasymptotic bounds on the mFSR and the FSR of the plug-in procedure and on its average selection number, by relying only on Assumption 1. To state the result, we introduce some additional quantities measuring the regularity of the model which will appear in our remainder terms. Recall definitions (4), (5), and (13) of $\ell_q(X, \theta)$, $T(X, \theta)$, and $t^*(\alpha)$ respectively, and let for $\epsilon, \delta, \nu > 0$,

$$W_{\ell}(\epsilon) = \sup \left\{ \sum_{q=1}^{Q} \mathbb{E}_{\theta^*} \left[ |\ell_q(X, \theta^*) - \ell_q(X, \theta)| \right], \|\theta - \theta^*\|_2 \le \epsilon, \theta \in \Theta \right\}; \quad (A1)$$

$$W_T(\delta) = \sup \{ |\mathbb{P}_{\theta^*}(T(X, \theta^*) < t') - \mathbb{P}_{\theta^*}(T(X, \theta^*) < t)|, \quad (A2)$$

$$t, t' \in [0, 1], |t' - t| \le \delta; \quad (A3)$$

$$\Psi(\epsilon) = W_T(W_{\ell}(\epsilon)^{1/2}) + W_{\ell}(\epsilon)^{1/2}; \quad (A4)$$
**Step 1: bounding $\hat{t}(\alpha)$ w.r.t. $t^*(\alpha)$** Recall (13), (27) and (26). In this part, we only consider realizations on the event $\Omega_\epsilon$. Let $\beta \in [\frac{2\alpha+\alpha_c}{3}, \frac{\alpha+\bar{\alpha}}{2}]$. By Lemma 13, we have

$$ \mathbf{L}_1(\hat{\theta}^{\sigma}, t^*(\beta)) \geq \mathbf{L}_1(\theta^*, t^*(\beta)) - \Psi(\|\hat{\theta}^{\sigma} - \theta^*\|_2) \geq \mathbf{L}_1(\theta^*, t^*(2\alpha + \alpha_c)/3) - \Psi(\epsilon), $$

because $t^*(\beta) \geq t^*(\frac{2\alpha+\alpha_c}{3})$ since $t^*(\cdot)$ is non decreasing by Lemma 11. Hence $\mathbf{L}_1(\hat{\theta}^{\sigma}, t^*(\beta)) \geq s^*$ for $\epsilon$ smaller than a threshold only depending on $\theta^*$ and $\alpha$. Hence, we have on $\Gamma_{\delta,t^*(\beta)}$ that

$$ \mathbf{L}(\hat{\theta}^{\sigma}, t^*(\beta)) - \delta \leq \hat{\mathbf{L}}(\hat{\theta}^{\sigma}, t^*(\beta)) \leq \delta + \mathbf{L}(\hat{\theta}^{\sigma}, t^*(\beta)). $$

By using again Lemma 13, we have

$$ \mathbf{L}(\theta^*, t^*(\beta)) - 3\Psi(\epsilon)/s^* \leq \mathbf{L}(\hat{\theta}^\sigma, t^*(\beta)) \leq \mathbf{L}(\theta^*, t^*(\beta)) + 3\Psi(\epsilon)/s^*. $$

Given that $\mathbf{L}(\theta^*, t^*(\beta)) = \text{mFCR}(C_{t^*(\beta)}^*) = \beta$ (see Lemma 9 (i)), it follows that for $\gamma = \gamma(\epsilon, \delta) = \delta + 4\Psi(\epsilon)/s^*$, on the event $\Gamma_{\delta,t^*(\alpha-\gamma)} \cap \Gamma_{\delta,t^*(\alpha+\gamma)}$,

$$ \hat{\mathbf{L}}(\hat{\theta}^{\sigma}, t^*(\alpha - \gamma)) \leq \alpha, \quad \hat{\mathbf{L}}(\hat{\theta}^{\sigma}, t^*(\alpha + \gamma)) > \alpha, $$

where we indeed check that $\alpha - \gamma \geq \frac{2\alpha+\alpha_c}{3}$ and $\alpha + \gamma \leq \frac{\alpha+\bar{\alpha}}{2}$ for $\delta$ and $\epsilon$ smaller than some threshold only depending on $\theta^*$ and $\alpha$. In a nutshell, we have established

$$ \Gamma_{\delta,t^*(\alpha-\gamma)} \cap \Gamma_{\delta,t^*(\alpha+\gamma)} \cap \Omega_{\epsilon} \subset \{t^*(\alpha-\gamma) \leq \hat{t}(\alpha) \leq t^*(\alpha+\gamma)\}. \quad (29) $$

**Step 2: upper-bounding the FCR** Let us consider the event

$$ \Lambda_{\alpha, \delta, \epsilon} := \Gamma_{0, \delta s^*/2, t^*(\alpha-\gamma)} \cap \Gamma_{1, \delta s^*/2, t^*(\alpha-\gamma)} \cap \Gamma_{0, \delta s^*/2, t^*(\alpha+\gamma)} \cap \Gamma_{1, \delta s^*/2, t^*(\alpha+\gamma)} \cap \Upsilon_{0, \delta, t^*(\alpha+\gamma)} \cap \Omega_{\epsilon}, $$

where the different events have been defined in the previous section.

Let us prove (21). By using Lemma 10 and (29),

$$ \begin{aligned} \text{FCR}(\hat{C}) &\leq \mathbb{E}_{\theta^*} [\widehat{\mathbf{M}}(\hat{\theta}^{\sigma}, \hat{t}(\alpha)) \mathbb{1}_{\Lambda_{\alpha, \delta, \epsilon}}] + \mathbb{P}((\Lambda_{\alpha, \delta, \epsilon})^c) \\ &\leq \mathbb{E}_{\theta^*} \left[ \frac{\widehat{\mathbf{M}}_0(\hat{\theta}^{\sigma}, t^*(\alpha + \gamma))}{\widehat{\mathbf{L}}_1(\hat{\theta}^{\sigma}, t^*(\alpha - \gamma))} \mathbb{1}_{\Lambda_{\alpha, \delta, \epsilon}} \right] + \mathbb{P}((\Lambda_{\alpha, \delta, \epsilon})^c). \end{aligned} $$

Now using a concentration argument on the event $\Lambda_{\alpha, \delta, \epsilon} \subset \Gamma_{1, \delta, t^*(\alpha-\gamma)} \cap \Upsilon_{0, \delta, t^*(\alpha+\gamma)}$, we have

$$ \begin{aligned} \text{FCR}(\hat{C}) &\leq \mathbb{E}_{\theta^*} \left[ \frac{\mathbf{M}_0(\hat{\theta}^{\sigma}, t^*(\alpha + \gamma)) + \delta}{\mathbf{L}_1(\hat{\theta}^{\sigma}, t^*(\alpha - \gamma)) - \delta} \mathbb{1}_{\Lambda_{\alpha, \delta, \epsilon}} \right] + \mathbb{P}((\Lambda_{\alpha, \delta, \epsilon})^c) \\ &\leq \frac{\mathbf{M}_0(\theta^*, t^*(\alpha + \gamma)) + 3\Psi(\epsilon) + \delta}{\mathbf{L}_1(\theta^*, t^*(\alpha - \gamma)) - \Psi(\epsilon) - \delta} + \mathbb{P}((\Lambda_{\alpha, \delta, \epsilon})^c) \\ &= \frac{\mathbf{L}_0(\theta^*, t^*(\alpha + \gamma)) + 3\Psi(\epsilon) + \delta}{\mathbf{L}_1(\theta^*, t^*(\alpha - \gamma)) - \Psi(\epsilon) - \delta} + \mathbb{P}((\Lambda_{\alpha, \delta, \epsilon})^c), \end{aligned} \quad (30) $$

by using Lemma 13 and that $\mathbf{M}_0(\theta^*, t) = \mathbf{L}_0(\theta^*, t)$ for all $t$ by definition. Now, using again Lemma 13, we have

$$ \begin{aligned} \mathbf{L}_0(\theta^*, t^*(\alpha + \gamma)) &\leq \mathbf{L}_0(\theta^*, t^*(\alpha - \gamma)) + \mathcal{W}_T(t^*(\alpha + \gamma) - t^*(\alpha - \gamma)) \\ &\leq \mathbf{L}_0(\theta^*, t^*(\alpha - \gamma)) + \mathcal{W}_T(\mathcal{W}_{t^*(\alpha, 2\gamma)}) \end{aligned} $$
$$W_{t^*, \alpha}(v) = \sup \{|t^*(\alpha + \beta) - t^*(\alpha)|, |\beta| \le v\}. \quad (A5)$$

**Theorem 3.** Let Assumption 1 be true. For any $\alpha \in (\alpha_c, \bar{\alpha})$ and constants $s^* = s^*(\alpha, \theta^*) \in (0, 1)$ and $e^* = e(\alpha, \theta^*) > 0$ depending only on $\alpha$ and $\theta^*$, the following holds. Consider the plug-in procedure $\hat{C}_\alpha^{\text{PI}} = (\hat{\mathbf{Z}}_\alpha^{\text{PI}}, \hat{S}_\alpha^{\text{PI}})$ introduced in Algorithm 2 and based on an estimator $\hat{\theta}$ satisfying Assumption 5, with $\eta(\epsilon, \theta^*)$ defined by (18). Then for $\epsilon \le e^*$ and $n \ge (2e)^3$, letting

$$\Delta_n(\epsilon) = 2(W_T(W_{t^*,\alpha}(2\delta_n + 8\Psi(\epsilon)/s^*)) + 4\Psi(\epsilon) + 2\delta_n),$$

for $\delta_n = C\sqrt{(\log n)/n}/s^*$ where $C = 2 + 56Q\sqrt{\mathcal{V}} + 28Q^2\sqrt{\mathcal{V}}$ and with the quantities $W_T$, $W_\ell$, $\Psi$, $W_{t^*,\alpha}$ defined by (A3), (A1), (A4), (A5), respectively, it holds:

* The procedure $\hat{C}_\alpha^{\text{PI}}$ controls both the FSR and the mFSR at level close to $\alpha$ in the following sense:

$$\begin{aligned} \text{FSR}(\hat{C}_\alpha^{\text{PI}}) &\le \alpha + \Delta_n(\epsilon)/s^* + 5/n^2 + \eta(\epsilon, \theta^*); \\ \text{mFSR}(\hat{C}_\alpha^{\text{PI}}) &\le \alpha + \Delta_n(\epsilon)/s^* + s^{-1}[50/n^2 + 10\eta(\epsilon, \theta^*)]. \end{aligned}$$

* The procedure $\hat{C}_\alpha^{\text{PI}}$ is nearly optimal in the following sense: for any other procedure $C = (\hat{\mathbf{Z}}, S)$ that controls the mFSR at level $\alpha$,

$$n^{-1}\mathbb{E}_{\theta^*}(|\hat{S}_\alpha^{\text{PI}}|) \ge n^{-1}\mathbb{E}_{\theta^*}(|S|) - \Delta_n(\epsilon).$$

Before proving this result (which will be done in the next subsections, see Sections A.5 and A.6), let us first show that Theorem 3 implies Theorems 1 and 2.

## A2 Proof of Theorem 1

By Lemma 4 below, $\Delta_n(\epsilon)$ tends to 0 when $n$ tends to infinity and $\epsilon$ tends to 0. Moreover, by consistency of $\hat{\theta}$, $\eta(\epsilon, \theta^*)$ tends to 0 for all $\epsilon > 0$. This implies the result.

**Lemma 4.** Under Assumption 1, we have $\lim_{\hat{\delta} \to 0} W_T(\delta) = 0$, $\lim_{v \to 0} W_{t^*, \alpha}(v) = 0$. Under Assumption 2, we have $\lim_{\epsilon \to 0} W_\ell(\epsilon) = 0$. Under both assumptions, we have $\lim_{\epsilon \to 0} \Psi(\epsilon) = 0$.

Lemma 4 is proved in Section B.4.

## A3 Proof of Theorem 2

Using Assumption 3 (with the notation therein) and Lemma 5, we have

$$\begin{aligned} \Delta_n(\epsilon) &= 2(W_T(W_{t^*,\alpha}(2\delta_n + 8\Psi(\epsilon)/s^*)) + 4\Psi(\epsilon) + 2\delta_n) \\ &\le 2C_2C_3(2\delta_n + (8/s^*)\sqrt{C_1}(C_2+1)\sqrt{\epsilon}) + 8\sqrt{C_1}(C_2+1)\sqrt{\epsilon} + 4\delta_n \\ &= 8\sqrt{C_1}(C_2+1)(1+2C_2C_3)\sqrt{\epsilon}/s^* + 4(C_2C_3+1)C\sqrt{\log n/n}/s^*, \end{aligned}$$

because $s^* \le 1$ and by definition of $\delta_n$. This gives (21) and (22) with $A = 8\sqrt{C_1}(C_2+1)(1+2C_2C_3)/s^{*2}$ and $B = 4(C_2C_3+1)C/s^{*2}$.
This entails

$$
\begin{align*}
\mathrm{FCR}(\hat{C}) &\le \frac{\mathbf{L}_0(\theta^*, t^*(\alpha - \gamma)) + \mathcal{W}_T(\mathcal{W}_{t^*, \alpha}(2\gamma)) + 3\Psi(\epsilon) + \delta}{\mathbf{L}_1(\theta^*, t^*(\alpha - \gamma)) - \Psi(\epsilon) - \delta} + \mathbb{P}((\Lambda_{\alpha, \delta, \epsilon})^c) \\
&\le \frac{\mathbf{L}_0(\theta^*, t^*(\alpha - \gamma))}{\mathbf{L}_1(\theta^*, t^*(\alpha - \gamma)) - \Psi(\epsilon) - \delta} + (s^*/2)^{-1} (\mathcal{W}_T(\mathcal{W}_{t^*, \alpha}(2\gamma)) + 3\Psi(\epsilon) + \delta) + \mathbb{P}((\Lambda_{\alpha, \delta, \epsilon})^c),
\end{align*}
$$

by choosing $\epsilon, \delta$ smaller than a threshold (only depending on $\theta^*$ and $\alpha$) so that $\mathbf{L}_1(\theta^*, t^*(\alpha - \gamma)) - \Psi(\epsilon) - \delta \ge s^*/2$. Now using $\mathbf{L}_0(\theta^*, t^*(\alpha - \gamma)) = (\alpha - \gamma)\mathbf{L}_1(\theta^*, t^*(\alpha - \gamma))$, we have

$$
\begin{align*}
\frac{\mathbf{L}_0(\theta^*, t^*(\alpha - \gamma))}{\mathbf{L}_1(\theta^*, t^*(\alpha - \gamma)) - \Psi(\epsilon) - \delta} &= (\alpha - \gamma) \left( 1 + \frac{\Psi(\epsilon) + \delta}{\mathbf{L}_1(\theta^*, t^*(\alpha - \gamma)) - \Psi(\epsilon) - \delta} \right) \\
&\le \alpha (1 + (s^*/2)^{-1} (\Psi(\epsilon) + \delta)).
\end{align*}
$$

This leads to

$$
\mathrm{FCR}(\hat{C}) \le \alpha + (2/s^*) (\mathcal{W}_T(\mathcal{W}_{t^*, \alpha}(2\delta + 8\Psi(\epsilon)/s^*)) + 4\Psi(\epsilon) + 2\delta) + \mathbb{P}((\Lambda_{\alpha, \delta, \epsilon})^c),
$$

which holds true for $\delta, \epsilon$ smaller than a threshold only depending on $\theta^*$ and $\alpha$.

**Step 3: upper-bounding the mFCR** We apply a similar technique as for step 2. By using Lemma 10 and (29),

$$
\begin{align*}
\mathrm{mFCR}(\hat{C}) &\le \frac{\mathbb{E}_{\theta^*}[\widehat{\mathbf{M}}_0(\hat{\theta}^\sigma, \hat{t}(\alpha))1_{\Lambda_{\alpha, \delta, \epsilon}}] + \mathbb{P}((\Lambda_{\alpha, \delta, \epsilon})^c)}{\mathbb{E}_{\theta^*}[\widehat{\mathbf{L}}_1(\hat{\theta}^\sigma, \hat{t}(\alpha))1_{\Lambda_{\alpha, \delta, \epsilon}}]} \\
&\le \frac{\mathbb{E}_{\theta^*}[\widehat{\mathbf{M}}_0(\hat{\theta}^\sigma, t^*(\alpha + \gamma))1_{\Lambda_{\alpha, \delta, \epsilon}}] + \mathbb{P}((\Lambda_{\alpha, \delta, \epsilon})^c)}{\mathbb{E}_{\theta^*}[\widehat{\mathbf{L}}_1(\hat{\theta}^\sigma, t^*(\alpha - \gamma))1_{\Lambda_{\alpha, \delta, \epsilon}}]}.
\end{align*}
$$

Now using a concentration argument on $\Lambda_{\alpha, \delta, \epsilon} \subset \Gamma_{1, \delta, t^*(\alpha-\gamma)} \cap \Upsilon_{0, \delta, t^*(\alpha+\gamma)}$, we have

$$
\begin{align*}
\mathrm{mFCR}(\hat{C}) &\le \frac{\mathbb{E}_{\theta^*}[(\mathbf{M}_0(\hat{\theta}^\sigma, t^*(\alpha + \gamma)) + \delta)1_{\Lambda_{\alpha, \delta, \epsilon}}] + \mathbb{P}((\Lambda_{\alpha, \delta, \epsilon})^c)}{\mathbb{E}_{\theta^*}[(\mathbf{L}_1(\hat{\theta}^\sigma, t^*(\alpha - \gamma)) - \delta)1_{\Lambda_{\alpha, \delta, \epsilon}}]} \\
&\le \frac{\mathbf{M}_0(\theta^*, t^*(\alpha + \gamma)) + 3\Psi(\epsilon) + \delta + \mathbb{P}((\Lambda_{\alpha, \delta, \epsilon})^c)}{\mathbf{L}_1(\theta^*, t^*(\alpha - \gamma)) - \Psi(\epsilon) - \delta - \mathbb{P}((\Lambda_{\alpha, \delta, \epsilon})^c)} \\
&= \frac{\mathbf{L}_0(\theta^*, t^*(\alpha + \gamma)) + 3\Psi(\epsilon) + \delta + \mathbb{P}((\Lambda_{\alpha, \delta, \epsilon})^c)}{\mathbf{L}_1(\theta^*, t^*(\alpha - \gamma)) - \Psi(\epsilon) - \delta - \mathbb{P}((\Lambda_{\alpha, \delta, \epsilon})^c)},
\end{align*}
$$

by using Lemma 13 and that $\mathbf{M}_0(\theta^*, t) = \mathbf{L}_0(\theta^*, t)$ by definition. Letting $x = \mathbf{L}_0(\theta^*, t^*(\alpha + \gamma)) + 3\Psi(\epsilon) + \delta$, $y = \mathbf{L}_1(\theta^*, t^*(\alpha - \gamma)) - \Psi(\epsilon) - \delta$ and $u = \mathbb{P}((\Lambda_{\alpha, \delta, \epsilon})^c)$, we have obtained the bound $(x + u)/(y - u)$, which has to be compared with the FCR bound (30), which reads $x/y + u$. Now, when $y \in [0, 1]$, $x \ge 0$, $x/y \le 2$, $u/y \le 1/2$, $y - u \ge s^*/2$, we have

$$
(x+u)/(y-u) \le \frac{x/y}{1-u/y} + (2/s^*)u \le x/y(1+2u/y) + (2/s^*)u \le x/y + (10/s^*)u.
$$

As a result, for $\epsilon, \delta$ small enough, and $\mathbb{P}((\Lambda_{\alpha, \delta, \epsilon})^c) \le s^*/4$, we obtain the same bound as for the FCR, with $\mathbb{P}((\Lambda_{\alpha, \delta, \epsilon})^c)$ replaced by $(10/s^*)\mathbb{P}((\Lambda_{\alpha, \delta, \epsilon})^c)$.
**Lemma 5.** Under Assumption 3, we have $W_c(\epsilon) \le C_1\epsilon$, $W_T(\delta) \le C_2\delta$, $W_{t^*,\alpha}(\nu) \le C_3\nu$ and $\Psi(\epsilon) \le \sqrt{C_1(C_2+1)}\sqrt{\epsilon}$ for $\epsilon, \delta, \nu$ small enough.

Lemma 5 is proved in Section B.5

## A4 An optimal procedure

We consider in this section the procedure that serves as an optimal procedure in our theory. For $t \in [0, 1]$, let $C_t^* = (\hat{Z}^*, S_t^*)$ be the procedure using the Bayes clustering $\hat{Z}^*$ (6) and the selection rule $S_t^* = \{i \in \{1, \dots, n\} : T_i^* < t\}$. Let us consider the map $t \in [0, 1] \mapsto \text{mFSR}(C_t^*)$ and note that $\text{mFSR}(C_t^*) = \text{mFSR}_t^*$ as defined by (12). Lemma 8 below provides the key properties for this function.

**Definition 3.** The optimal procedure at level $\alpha$ is defined by $C_{t^*(\alpha)}^*$ where $t^*(\alpha)$ is defined by (13).

Under Assumption 1, Lemma 8 entails that, for $\alpha > \alpha_c$, mFSR($C_{t^*(\alpha)}^*$) $\le \alpha$. Hence, $C_{t^*(\alpha)}^*$ controls the mFSR at level $\alpha$. In addition, it is optimal in the following sense: any other mFSR controlling procedure should select less items than $C_{t^*(\alpha)}^*$.

**Lemma 6 (Optimality of $C_{t^*(\alpha)}^*$).** Let Assumption 1 be true and choose $\alpha \in (\alpha_c, \bar{\alpha}]$. Then the oracle procedure $C_{t^*(\alpha)}^* = (\hat{Z}^*, S_{t^*(\alpha)}^*)$ satisfies the following:

(i) mFSR($C_{t^*(\alpha)}^*$) = $\alpha$;

(ii) for any procedure $C = (\hat{Z}, S)$ such that mFSR($C$) $\le \alpha$, we have $\mathbb{E}_{\theta^*}(|S|) \le \mathbb{E}_{\theta^*}(|S_{t^*(\alpha)}^*|)$.

Lemma 6 is proved in Section B.6.

## A5 Preliminary steps for proving Theorem 3

To keep the main proof concise, we need to define several additional notation. Let for $t \in [0, 1]$ and $\theta \in \Theta$ (recall 5)

$$ \hat{\mathbf{L}}_0(\theta, t) = \frac{1}{n} \sum_{i=1}^{n} T(X_i, \theta) \mathbf{1}_{T(X_i, \theta) < t}; \quad (A6) $$

$$ \hat{\mathbf{L}}_1(\theta, t) = \frac{1}{n} \sum_{i=1}^{n} \mathbf{1}_{T(X_i, \theta) < t}. \quad (A7) $$

Denote $\hat{\mathbf{L}} = \hat{\mathbf{L}}_0/\hat{\mathbf{L}}_1$, $\mathbf{L}_0 = \mathbb{E}_{\theta^*} \hat{\mathbf{L}}_0$, $\mathbf{L}_1 = \mathbb{E}_{\theta^*} \hat{\mathbf{L}}_1$, $\mathbf{L} = \mathbf{L}_0/\mathbf{L}_1$ (with the convention $0/0 = 0$). Note that for any $\alpha > \alpha_c$, the mFSR of the optimal procedure $C_{t^*(\alpha)}^*$ defined in Section A.4 is given by $\text{mFSR}(C_{t^*(\alpha)}^*) = \mathbf{L}(\theta^*, t^*(\alpha)) = \alpha$.

Also, we denote from now on $\ell_{i,q}^* = \mathbb{P}_{\theta^*}(Z_i = q|X_i)$ for short and introduce for any parameter $\theta \in \Theta$ (recall 4 and 5)

$$ \bar{q}(X_i, \theta) \in \underset{q \in \{1, \dots, Q\}}{\operatorname{argmax}} \ell_q(X_i, \theta), \quad 1 \le i \le n; \quad (A8) $$

$$ U(X_i, \theta) = 1 - \ell_{i, \bar{q}(X_i, \theta)}^*, \quad 1 \le i \le n; \quad (A9) $$

$$ \hat{\mathbf{M}}_0(\theta, t) = \frac{1}{n} \sum_{i=1}^{n} U(X_i, \theta) \mathbf{1}_{T(X_i, \theta) < t}, \quad t \in [0, 1], \quad (A10) $$
**Step 4: lower-bounding the selection rate** In Step 3, when bounding the mFCR, we derived a lower bound for the denominator of the mFCR, that is, $\mathbb{E}_{\theta^*}(|\hat{S}|)$. It reads

$$ 
\begin{aligned}
n^{-1} \mathbb{E}_{\theta^*}(|\hat{S}|) &\geq \mathbf{L}_1(\theta^*, t^*(\alpha - \gamma)) - \Psi(\epsilon) - \delta - \mathbb{P}((\Lambda_{\alpha, \delta, \epsilon})^c) \\
&\geq \mathbf{L}_1(\theta^*, t^*(\alpha)) - \mathcal{W}_T(t^*(\alpha) - t^*(\alpha - \gamma)) - \Psi(\epsilon) - \delta - \mathbb{P}((\Lambda_{\alpha, \delta, \epsilon})^c) \\
&\geq n^{-1} \mathbb{E}_{\theta^*}(|S_{t^*(\alpha)}^*|) - \mathcal{W}_T(\mathcal{W}_{t^*, \alpha}(\gamma)) - \Psi(\epsilon) - \delta - \mathbb{P}((\Lambda_{\alpha, \delta, \epsilon})^c),
\end{aligned}
$$

by using (17) and (19). Now consider another procedure $\mathcal{C} = (\hat{\mathbf{Z}}, S)$ that controls the mFCR at level $\alpha$, that is, mFCR($\mathcal{C}$) $\leq \alpha$. By Lemma 9, we then have $\mathbb{E}_{\theta^*}(|S_{t^*(\alpha)}^*|) \geq \mathbb{E}_{\theta^*}(|S|)$.

**Step 5: concentration** Finally, we bound $\mathbb{P}((\Lambda_{\alpha, \delta, \epsilon})^c)$ by using Lemma 14 with $u = \sqrt{(\log n)/n}$ and $x = 2\sqrt{(\log n)/n}$. This gives for $\delta = 2x/s^*$, and $n$ such that $\log(n) \geq 1$,

$$ \mathbb{P}((\Lambda_{\alpha, \delta, \epsilon})^c) \leq 5/n + \mathbb{P}(\Omega_{\epsilon}^c). $$

## B Proofs of lemmas

**Proof of Lemma 1** The clustering risk of $\hat{\mathbf{Z}}$ is given by

$$ 
\begin{aligned}
R(\hat{\mathbf{Z}}) &= \mathbb{E}_{\theta^*} \left( \min_{\sigma \in [Q]} \mathbb{E}_{\theta^*} \left( n^{-1} \sum_{i=1}^n \mathbb{1}\{Z_i \neq \sigma(\hat{Z}_i)\} \middle| \mathbf{X} \right) \right) \\
&= \mathbb{E}_{\theta^*} \left( \min_{\sigma \in [Q]} n^{-1} \sum_{i=1}^n \mathbb{P}_{\theta^*}(Z_i \neq \sigma(\hat{Z}_i) \middle| \mathbf{X}) \right) \\
&\geq \mathbb{E}_{\theta^*} \left( \min_{\hat{\mathbf{Z}}} n^{-1} \sum_{i=1}^n \mathbb{P}_{\theta^*}(Z_i \neq \hat{Z}_i \middle| \mathbf{X}) \right),
\end{aligned}
$$

where, by independence, the minimum in the lower bound is achieved for the Bayes clustering. Thus, $R(\hat{\mathbf{Z}}) \geq n^{-1} \sum_{i=1}^n \mathbb{E}_{\theta^*}(T_i^*)$. Moreover, $n^{-1} \sum_{i=1}^n \mathbb{E}_{\theta^*}(T_i^*) \geq R(\hat{\mathbf{Z}}^*)$, since

$$ R(\hat{\mathbf{Z}}^*) = \mathbb{E}_{\theta^*} \left( \min_{\sigma \in [Q]} n^{-1} \sum_{i=1}^n \mathbb{P}_{\theta^*}(Z_i \neq \sigma(\hat{Z}_i^*) | \mathbf{X}) \right) \leq \mathbb{E}_{\theta^*} \left( n^{-1} \sum_{i=1}^n \mathbb{P}_{\theta^*}(Z_i \neq \hat{Z}_i^* | \mathbf{X}) \right). $$

Thus, $\min_{\hat{\mathbf{Z}}} R(\hat{\mathbf{Z}}) = R(\hat{\mathbf{Z}}^*)$ and the proof is completed.

**Proof of Lemma 2** Following the reasoning of the proof of Lemma 1 that, we have

$$ FCR_{\theta^*}(\mathcal{C}) = \mathbb{E}_{\theta^*} \left( \min_{\sigma \in [Q]} \mathbb{E}_{\theta^*} \left( \frac{\sum_{i \in S} \mathbb{1}\{Z_i \neq \sigma(\hat{Z}_i^*)\}}{\max(|S|, 1)} \middle| \mathbf{X} \right) \right) = \mathbb{E}_{\theta^*} \left( \frac{\sum_{i \in S} T_i^*}{\max(|S|, 1)} \right). $$

**Proof of Lemma 3** By definition, we have

$$ FCR(\hat{\mathcal{C}}_{\alpha}^{\mathrm{PI}}) = \mathbb{E}_{\theta^*} \left( \min_{\sigma \in [Q]} \mathbb{E}_{\theta^*} \left( \frac{\sum_{i=1}^n \mathbb{1}_{Z_i \neq \sigma(\hat{Z}_i^{\mathrm{PI}}(\mathbf{X}))} \mathbb{1}\{i \in \hat{S}^{\mathrm{PI}}(\mathbf{X})\}}{\max(|\hat{S}^{\mathrm{PI}}(\mathbf{X})|, 1)} \middle| \mathbf{X} \right) \right) $$

so that (9) follows by a direct integration w.r.t. the latent variable Z.
Note that $\hat{\mathbf{M}}_0(\theta^*, t) = \hat{\mathbf{L}}_0(\theta^*, t)$ but in general $\hat{\mathbf{M}}_0(\theta, t)$ is different from $\hat{\mathbf{L}}_0(\theta, t)$. We denote $\hat{\mathbf{M}} = \hat{\mathbf{M}}_0/\hat{\mathbf{L}}_1$, $\mathbf{M}_0 = \mathbb{E}_{\theta^*}\hat{\mathbf{M}}_0$ and $\mathbf{M} = \mathbf{M}_0/\mathbf{L}_1$ (with the convention $0/0 = 0$).

When $\alpha \in (\alpha_c, \bar{\alpha}]$ (recall 14 and 15), we also let

$$s^* = s^*(\alpha, \theta^*) = n^{-1}\mathbb{E}_{\theta^*}\left(\left|S_{t^*\left(\frac{\alpha+\alpha_c}{2}\right)}^*\right|\right) = \mathbf{L}_1(\theta^*, t^*\left(\frac{\alpha+\alpha_c}{2}\right)) > 0. \quad (A11)$$

We easily see that the latter is positive: if it was zero then $S_{t^*\left(\frac{\alpha+\alpha_c}{2}\right)}^*$ would be empty which would entails that mFSR($C_{t^*\left(\frac{\alpha+\alpha_c}{2}\right)}^*$) is zero. This is excluded by definition (14) of $\alpha_c$ because $(\alpha + \alpha_c)/2 > \alpha_c$.

Also, we are going to extensively use the event

$$\Omega_{\epsilon} = \left\{ \min_{\sigma \in [Q]} \|\hat{\theta}^{\sigma} - \theta^{*}\|_{2} < \epsilon \right\}.$$

On this event, we fix any permutation $\sigma \in [Q]$ (possibly depending on $X$) such that $\|\hat{\theta}^{\sigma} - \theta^{*}\|_{2} < \epsilon$. Now using Lemma 9, the plug-in selection rule can be rewritten as $\hat{S}_{\alpha}^{\text{PI}} = \{i \in \{1, \dots, n\} : \hat{T}_{i} < \hat{t}(\alpha)\}$ (denoted by $\hat{S}$ in the sequel for short), where

$$\hat{t}(\alpha) = \sup\{t \in [0, 1] : \hat{\mathbf{L}}(\hat{\theta}, t) \le \alpha\}. \quad (A12)$$

With the above notation, we can upper bound what is inside the brackets of FSR($\hat{C}^{\text{PI}}$) and mFSR($\hat{C}^{\text{PI}}$) as follows.

**Lemma 7.** For the permutation $\sigma$ in $\Omega_{\epsilon}$ realizing $\|\hat{\theta}^{\sigma} - \theta^{*}\|_{2} < \epsilon$, we have on the event $\Omega_{\epsilon}$ the following relations:

$$\begin{aligned} |\hat{S}| &= \hat{\mathbf{L}}_1(\hat{\theta}^{\sigma}, \hat{t}(\alpha)); \\ \min_{\sigma' \in [Q]} \mathbb{E}_{\theta^{*}} \left( \varepsilon_{\hat{S}}(\sigma'(\hat{\mathbf{Z}}), \mathbf{Z}) \middle| \mathbf{X} \right) &\le \hat{\mathbf{M}}_0(\hat{\theta}^{\sigma}, \hat{t}(\alpha)); \\ \min_{\sigma' \in [Q]} \mathbb{E}_{\theta^{*}} \left( \frac{\varepsilon_{\hat{S}}(\sigma'(\hat{\mathbf{Z}}), \mathbf{Z})}{\max(|\hat{S}|, 1)} \middle| \mathbf{X} \right) &\le \hat{\mathbf{M}}(\hat{\theta}^{\sigma}, \hat{t}(\alpha)). \end{aligned}$$

Lemma 7 is proved in Section B.7.

Finally, we make use of the concentration of the empirical processes $\hat{\mathbf{L}}_0(\theta, t)$, $\hat{\mathbf{L}}_1(\theta, t)$, and $\hat{\mathbf{M}}_0(\theta, t)$, uniformly with respect to $\theta \in D$ (where $D$ is defined in Assumption 5). Thus, we define the following events, for $\delta > 0$ (recall $s^*$ defined by A11):

$$\begin{aligned} \Gamma_{0,\delta,t} &= \left\{ \sup_{\theta \in D} \left| \hat{\mathbf{L}}_0(\theta, t) - \mathbf{L}_0(\theta, t) \right| \le \delta \right\}; \\ \Gamma_{1,\delta,t} &= \left\{ \sup_{\theta \in D} \left| \hat{\mathbf{L}}_1(\theta, t) - \mathbf{L}_1(\theta, t) \right| \le \delta \right\}; \\ \Gamma_{\delta,t} &= \left\{ \sup_{\theta \in D, \mathbf{L}_1(\theta,t) \ge s^*} \left| \hat{\mathbf{L}}(\theta, t) - \mathbf{L}(\theta, t) \right| \le \delta \right\}; \\ \Upsilon_{0,\delta,t} &= \left\{ \sup_{\theta \in D} \left| \hat{\mathbf{M}}_0(\theta, t) - \mathbf{M}_0(\theta, t) \right| \le \delta \right\}. \end{aligned}$$
**Proof of Lemma 4** The only non-trivial fact is for $W_{t^*,\alpha}(v)$. Assumption 1 and Lemma 11 provide that $t \mapsto \text{mFCR}_t^*$ is a one-to-one continuous increasing map from $(t^*(\alpha_c), t^*(\bar{\alpha}))$ to $(\alpha_c, \bar{\alpha})$. Hence, for $\alpha \in (\alpha_c, \bar{\alpha})$, $\beta \mapsto t^*(\alpha + \beta)$ is continuous in 0 and $\lim_{v \to 0} W_{t^*,\alpha}(v) = 0$.

**Proof of Lemma 9** By Lemma 11, we have that mFCR($C_t^*$) is monotonous in $t$ and continuous w.r.t. $t$ on $(t^*(\alpha_c), 1]$, thus for $\alpha \in (\alpha_c, \bar{\alpha}]$, mFCR($C_{t^*(\alpha)}^*$) = $\alpha$ which gives (i). For (ii), let $C = (\hat{\mathbf{Z}}, S)$ be a procedure such that mFCR($C$) $\le \alpha$. Let us consider the procedure $C'$ with the Bayes clustering $\hat{\mathbf{Z}}^*$ and the same selection rule $S$. Since $C'$ is based on a Bayes clustering, by the same reasoning leading to $R(\hat{\mathbf{Z}}^*) \le R(\hat{\mathbf{Z}})$ in Section 3.1, we have that mFCR($C'$) $\le mFCR(C) \le \alpha$ with

$$ \text{mFCR}(C') = \frac{\mathbb{E}_{\theta^*} \left( \sum_{i \in S} T_i^* \right)}{\mathbb{E}_{\theta^*}(|S|)}. $$

Hence,

$$ \mathbb{E}_{\theta^*} \left( \sum_{i \in S} T_i^* \right) \le \alpha \mathbb{E}_{\theta^*}(|S|). \qquad (31) $$

Now we use an argument similar to the proof of Theorem 1 in Cai et al. (2019). By definition of $S_{t^*(\alpha)}^*$, we have that

$$ \sum_{i=1}^{n} \left( \mathbb{1}_{i \in S_{t^*(\alpha)}^*}(\mathbf{x}) - \mathbb{1}_{i \in S}(\mathbf{x}) \right) (T_i^* - t^*(\alpha)) \le 0 $$

which we can rewrite as

$$ \sum_{i=1}^{n} \left( \mathbb{1}_{i \in S_{t^*(\alpha)}^*}(\mathbf{x}) - \mathbb{1}_{i \in S}(\mathbf{x}) \right) (T_i^* - t^*(\alpha) + \alpha - \alpha) \le 0 $$

and so

$$ \mathbb{E}_{\theta^*} \left( \sum_{i=1}^{n} \left( \mathbb{1}_{i \in S_{t^*(\alpha)}^*}(\mathbf{x}) - \mathbb{1}_{i \in S}(\mathbf{x}) \right) (T_i^* - \alpha) \right) \le (t^*(\alpha) - \alpha) \mathbb{E}_{\theta^*} \left( \sum_{i=1}^{n} \left( \mathbb{1}_{i \in S_{t^*(\alpha)}^*}(\mathbf{x}) - \mathbb{1}_{i \in S}(\mathbf{x}) \right) \right) \\ = (t^*(\alpha) - \alpha) (\mathbb{E}_{\theta^*}(|S_{t^*(\alpha)}^*|) - \mathbb{E}_{\theta^*}(|S|)). $$

On the other hand, mFCR($C_{t^*(\alpha)}^*$) = $\alpha$ together with (31) implies that

$$ \mathbb{E}_{\theta^*} \left( \sum_{i=1}^{n} \left( \mathbb{1}_{i \in S_{t^*(\alpha)}^*}(\mathbf{x}) - \mathbb{1}_{i \in S}(\mathbf{x}) \right) (T_i^* - \alpha) \right) \\ = \mathbb{E}_{\theta^*} \left( \sum_{i \in S_{t^*(\alpha)}^*} T_i^* - \alpha |S_{t^*(\alpha)}^*| - \sum_{i \in S} T_i^* + \alpha |S| \right) \ge 0. $$

Combining, the relations above provides

$$ (t^*(\alpha) - \alpha)(\mathbb{E}_{\theta^*}(|S_{t^*(\alpha)}^*|) - \mathbb{E}_{\theta^*}(|S|)) \ge 0. $$

Finally, noting that $t^*(\alpha) - \alpha > 0$ since $\alpha = \text{mFCR}(C_{t^*(\alpha)}^*) < t^*(\alpha)$ by (ii) Lemma 11, this gives $\mathbb{E}_{\theta^*}(|S_{t^*(\alpha)}^*|) - \mathbb{E}_{\theta^*}(|S|) \ge 0$ and concludes the proof.
Note that the following holds:

$$ \Gamma_{0,\delta s^*/2,t} \cap \Gamma_{1,\delta s^*/2,t} \subset \Gamma_{\delta,t}. \quad (A13) $$

Indeed, on the event $\Gamma_{0,\delta s^*/2,t} \cap \Gamma_{1,\delta s^*/2,t}$, provided that $\mathbf{L}_1(\theta, t) \ge s^*$, we have

$$ \begin{aligned} & \left| \frac{\hat{\mathbf{L}}_0(\theta, t)}{\hat{\mathbf{L}}_1(\theta, t)} - \frac{\mathbf{L}_0(\theta, t)}{\mathbf{L}_1(\theta, t)} \right| \\ & \le \left| \frac{\mathbf{L}_0(\theta, t) - \hat{\mathbf{L}}_0(\theta, t)}{\mathbf{L}_1(\theta, t)} \right| + \hat{\mathbf{L}}_0(\theta, t) \left| \frac{1}{\hat{\mathbf{L}}_1(\theta, t)} - \frac{1}{\mathbf{L}_1(\theta, t)} \right| \\ & \le (\delta s^*/2)/s^* + (\delta s^*/2)/s^* = \delta, \end{aligned} $$

because $\hat{\mathbf{L}}_0(\theta, t) \le \hat{\mathbf{L}}_1(\theta, t)$. This proves the desired inclusion.

## A6 Proof of Theorem 3

Let us now provide a proof for Theorem 3.

### A6.1 Step 1: Bounding $\hat{t}(\alpha)$ w.r.t. $t^*(\alpha)$

Recall (13), (A12) and (A11). In this part, we only consider realizations on the event $\Omega_c$. Let $\beta \in [\frac{2\alpha+\alpha_c}{3}, \frac{\alpha+\bar{\alpha}}{2}]$. By Lemma 10, we have

$$ \mathbf{L}_1(\hat{\theta}^{\sigma}, t^*(\beta)) \ge \mathbf{L}_1(\theta^*, t^*(\beta)) - \Psi(\|\hat{\theta}^{\sigma} - \theta^*\|_2) \ge \mathbf{L}_1(\theta^*, t^*(2\alpha + \alpha_c)/3) - \Psi(\epsilon), $$

because $t^*(\beta) \ge t^*(\frac{2\alpha+\alpha_c}{3})$ since $t^*(\cdot)$ is nondecreasing by Lemma 8. Hence $\mathbf{L}_1(\hat{\theta}^{\sigma}, t^*(\beta)) \ge s^*$ for $\epsilon$ smaller than a threshold only depending on $\theta^*$ and $\alpha$. Hence, we have on $\Gamma_{\delta,t^*(\beta)}$ that

$$ \mathbf{L}(\hat{\theta}^{\sigma}, t^*(\beta)) - \delta \le \hat{\mathbf{L}}(\hat{\theta}^{\sigma}, t^*(\beta)) \le \delta + \mathbf{L}(\hat{\theta}^{\sigma}, t^*(\beta)). $$

By using again Lemma 10, we have

$$ \mathbf{L}(\theta^*, t^*(\beta)) - 3\Psi(\epsilon)/s^* \le \mathbf{L}(\hat{\theta}^{\sigma}, t^*(\beta)) \le \mathbf{L}(\theta^*, t^*(\beta)) + 3\Psi(\epsilon)/s^*. $$

Given that $\mathbf{L}(\theta^*, t^*(\beta)) = \text{mFSR}(C_{t^*(\beta)}^*) = \beta$ (see Lemma 6 (i)), it follows that for $\gamma = \gamma(\epsilon, \delta) = \delta + 4\Psi(\epsilon)/s^*$, on the event $\Gamma_{\delta,t^*(\alpha-\gamma)} \cap \Gamma_{\delta,t^*(\alpha+\gamma)}$,

$$ \hat{\mathbf{L}}(\hat{\theta}^{\sigma}, t^*(\alpha - \gamma)) \le \alpha, \quad \hat{\mathbf{L}}(\hat{\theta}^{\sigma}, t^*(\alpha + \gamma)) > \alpha, $$

where we indeed check that $\alpha - \gamma \ge \frac{2\alpha+\alpha_c}{3}$ and $\alpha + \gamma \le \frac{\alpha+\bar{\alpha}}{2}$ for $\delta$ and $\epsilon$ smaller than some threshold only depending on $\theta^*$ and $\alpha$. In a nutshell, we have established

$$ \Gamma_{\delta,t^*(\alpha-\gamma)} \cap \Gamma_{\delta,t^*(\alpha+\gamma)} \cap \Omega_c \subset \{t^*(\alpha-\gamma) \le \hat{t}(\alpha) \le t^*(\alpha+\gamma)\}. \quad (A14) $$

### A6.2 Step 2: Upper-bounding the FSR

Let us consider the event
**Proof of Lemma 10** First, we have by definition $\ell_q(X_i, \theta^\sigma) = \ell_{\sigma(q)}(X_i, \theta)$ and thus $T(X_i, \hat{\theta}) = T(X_i, \hat{\theta}^\sigma)$ by taking the maximum over $q$. This gives $\hat{S}^\sigma = \hat{S}$ and yields the first equality. Next, we have on $\Omega_\epsilon$,

$$ \min_{\sigma' \in [Q]} \mathbb{E}_{\theta^*} \left( \varepsilon_{\hat{S}}(\sigma'(\hat{\mathbf{Z}}), \mathbf{Z}) \middle| \mathbf{X} \right) \le \mathbb{E}_{\theta^*} \left( \varepsilon_{\hat{S}}(\sigma(\hat{\mathbf{Z}}), \mathbf{Z}) \middle| \mathbf{X} \right) \\ \le \mathbb{E}_{\theta^*} \left( \varepsilon_{\hat{S}^\sigma}(\sigma(\hat{\mathbf{Z}}), \mathbf{Z}) \middle| \mathbf{X} \right), $$

still because $\hat{S}^\sigma = \hat{S}$. Now observe that,

$$ \begin{aligned} \mathbb{E}_{\theta^*} \left( \varepsilon_{\hat{S}^\sigma}(\sigma(\hat{\mathbf{Z}}), \mathbf{Z}) \middle| \mathbf{X} \right) &= \frac{1}{n} \sum_{i=1}^{n} \mathbb{P}_{\theta^*} (Z_i \neq \sigma(\bar{q}(X_i, \hat{\theta})) \mid \mathbf{X}) \mathbb{1}_{T(X_i, \hat{\theta}^\sigma) < \hat{t}(\alpha)} \\ &= \frac{1}{n} \sum_{i=1}^{n} (1 - \ell_{i, \sigma(\bar{q}(X_i, \hat{\theta}))}^*) \mathbb{1}_{T(X_i, \hat{\theta}^\sigma) < \hat{t}(\alpha)} \\ &= \hat{M}_0(\hat{\theta}^\sigma, \hat{t}(\alpha)), \end{aligned} $$

because $\sigma(\bar{q}(X_i, \hat{\theta})) = \bar{q}(X_i, \hat{\theta}^\sigma)$. This proves the result.

## C Auxiliary results

**Lemma 11.** Let us consider the procedure $C_t^*$ defined in Section A.1 and the functional mFCR$_t^*$ defined by (12). Then we have

$$ \text{mFCR}(C_t^*) = \frac{\mathbb{E}_{\theta^*} \left( \sum_{i=1}^n T_i^* \mathbb{1}_{T_i^* < t} \right)}{\mathbb{E}_{\theta^*} \left( \sum_{i=1}^n \mathbb{1}_{T_i^* < t} \right)} = \text{mFCR}_t^*, \quad t \in [0, 1]. \quad (32) $$

Moreover, the following properties for the function $t \in [0, 1] \mapsto \text{mFCR}(C_t^*)$:

(i) mFCR($C_t^*$) is non-decreasing in $t \in [0, 1]$ and, under Assumption 1, it is increasing in $t \in (t^*(\alpha_c), t^*(\bar{\alpha}))$;

(ii) mFCR($C_t^*$) $< t$ for $t \in (0, 1]$;

(iii) Under Assumption 1, mFCR($C_t^*$) is continuous w.r.t. $t$ on $(t^*(\alpha_c), 1]$, where $t^*(\alpha_c)$ is given by (14).

**Proof.** First, (32) is obtained similarly than (8). For proving (i), let $t_1, t_2 \in [0, 1]$ such that $t_1 < t_2$. We show that mFCR($C_{t_1}^*$) $\le$ mFCR($C_{t_2}^*$). Remember here the convention $0/0 = 0$ and that mFCR($C_t^*$) = $\mathbb{E}_{\theta^*} (T(X, \theta^*) \mid T(X, \theta^*) < t)$. First, if $\mathbb{P}_{\theta^*} (T(X, \theta^*) < t_1) = 0$ then the result is immediate. Otherwise, we have that

$$ \text{mFCR}(C_{t_1}^*) - \text{mFCR}(C_{t_2}^*) = (\mathbb{P}_{\theta^*} (T(X, \theta^*) < t_1))^{-1} \mathbb{E}_{\theta^*} \left( T(X, \theta^*) \left\{ \mathbb{1}_{T(X, \theta^*) < t_1} - \frac{\mathbb{P}_{\theta^*} (T(X, \theta^*) < t_1)}{\mathbb{P}_{\theta^*} (T(X, \theta^*) < t_2)} \mathbb{1}_{T(X, \theta^*) < t_2} \right\} \right) $$

where, given that $t_1 < t_2$, the quantity in the brackets is positive when $T(X, \theta^*) < t_1$ and is negative or zero otherwise. Hence,

$$ T(X, \theta^*) \left\{ \mathbb{1}_{T(X, \theta^*) < t_1} - \frac{\mathbb{P}_{\theta^*} (T(X, \theta^*) < t_1)}{\mathbb{P}_{\theta^*} (T(X, \theta^*) < t_2)} \mathbb{1}_{T(X, \theta^*) < t_2} \right\} \\ \le t_1 \left\{ \mathbb{1}_{T(X, \theta^*) < t_1} - \frac{\mathbb{P}_{\theta^*} (T(X, \theta^*) < t_1)}{\mathbb{P}_{\theta^*} (T(X, \theta^*) < t_2)} \mathbb{1}_{T(X, \theta^*) < t_2} \right\}. $$
$$ \Lambda_{\alpha, \delta, \epsilon} := \Gamma_{0, \delta s^*/2, t^*(\alpha-\gamma)} \cap \Gamma_{1, \delta s^*/2, t^*(\alpha-\gamma)} \cap \Gamma_{0, \delta s^*/2, t^*(\alpha+\gamma)} \\ \cap \Gamma_{1, \delta s^*/2, t^*(\alpha+\gamma)} \cap \Upsilon_{0, \delta, t^*(\alpha+\gamma)} \cap \Omega_{\epsilon}, $$

where the different events have been defined in the previous section.

Let us prove (21). Using Lemma 7 and (A14),

$$ \begin{aligned} \mathrm{FSR}(\hat{C}) &\le \mathbb{E}_{\theta^*} [\hat{\mathbf{M}}(\hat{\theta}^{\sigma}, \hat{t}(\alpha)) 1_{\Lambda_{\alpha, \delta, \epsilon}}] + \mathbb{P}((\Lambda_{\alpha, \delta, \epsilon})^c) \\ &\le \mathbb{E}_{\theta^*} \left[ \frac{\hat{\mathbf{M}}_0(\hat{\theta}^{\sigma}, t^*(\alpha+\gamma))}{\hat{\mathbf{L}}_1(\hat{\theta}^{\sigma}, t^*(\alpha-\gamma))} 1_{\Lambda_{\alpha, \delta, \epsilon}} \right] + \mathbb{P}((\Lambda_{\alpha, \delta, \epsilon})^c). \end{aligned} $$

Now using a concentration argument on the event $\Lambda_{\alpha, \delta, \epsilon} \subset \Gamma_{1, \delta, t^*(\alpha-\gamma)} \cap \Upsilon_{0, \delta, t^*(\alpha+\gamma)}$, we have

$$ \begin{aligned} \mathrm{FSR}(\hat{C}) &\le \mathbb{E}_{\theta^*} \left[ \frac{\mathbf{M}_0(\hat{\theta}^{\sigma}, t^*(\alpha+\gamma)) + \delta}{\mathbf{L}_1(\hat{\theta}^{\sigma}, t^*(\alpha-\gamma)) - \delta} 1_{\Lambda_{\alpha, \delta, \epsilon}} \right] + \mathbb{P}((\Lambda_{\alpha, \delta, \epsilon})^c) \\ &\le \frac{\mathbf{M}_0(\theta^*, t^*(\alpha+\gamma)) + 3\Psi(\epsilon) + \delta}{\mathbf{L}_1(\theta^*, t^*(\alpha-\gamma)) - \Psi(\epsilon) - \delta} + \mathbb{P}((\Lambda_{\alpha, \delta, \epsilon})^c) \\ &= \frac{\mathbf{L}_0(\theta^*, t^*(\alpha+\gamma)) + 3\Psi(\epsilon) + \delta}{\mathbf{L}_1(\theta^*, t^*(\alpha-\gamma)) - \Psi(\epsilon) - \delta} + \mathbb{P}((\Lambda_{\alpha, \delta, \epsilon})^c), \end{aligned} \quad (A15) $$

by using Lemma 10 and that $\mathbf{M}_0(\theta^*, t) = \mathbf{L}_0(\theta^*, t)$ for all $t$ by definition. Now, using again Lemma 10, we have

$$ \begin{aligned} \mathbf{L}_0(\theta^*, t^*(\alpha + \gamma)) &\le \mathbf{L}_0(\theta^*, t^*(\alpha - \gamma)) + \mathcal{W}_T(t^*(\alpha + \gamma) - t^*(\alpha - \gamma)) \\ &\le \mathbf{L}_0(\theta^*, t^*(\alpha - \gamma)) + \mathcal{W}_T(\mathcal{W}_{t^*, \alpha}(2\gamma)) \end{aligned} $$

This entails

$$ \begin{aligned} \mathrm{FSR}(\hat{C}) &\le \frac{\mathbf{L}_0(\theta^*, t^*(\alpha - \gamma)) + \mathcal{W}_T(\mathcal{W}_{t^*, \alpha}(2\gamma)) + 3\Psi(\epsilon) + \delta}{\mathbf{L}_1(\theta^*, t^*(\alpha - \gamma)) - \Psi(\epsilon) - \delta} + \mathbb{P}((\Lambda_{\alpha, \delta, \epsilon})^c) \\ &\le \frac{\mathbf{L}_0(\theta^*, t^*(\alpha - \gamma))}{\mathbf{L}_1(\theta^*, t^*(\alpha - \gamma)) - \Psi(\epsilon) - \delta} \\ &\quad + (s^*/2)^{-1}(\mathcal{W}_T(\mathcal{W}_{t^*, \alpha}(2\gamma)) + 3\Psi(\epsilon) + \delta) + \mathbb{P}((\Lambda_{\alpha, \delta, \epsilon})^c), \end{aligned} $$

by choosing $\epsilon, \delta$ smaller than a threshold (only depending on $\theta^*$ and $\alpha$) so that $\mathbf{L}_1(\theta^*, t^*(\alpha - \gamma)) - \Psi(\epsilon) - \delta \ge s^*/2$. Now using $\mathbf{L}_0(\theta^*, t^*(\alpha - \gamma)) = (\alpha - \gamma)\mathbf{L}_1(\theta^*, t^*(\alpha - \gamma))$, we have

$$ \begin{aligned} \frac{\mathbf{L}_0(\theta^*, t^*(\alpha - \gamma))}{\mathbf{L}_1(\theta^*, t^*(\alpha - \gamma)) - \Psi(\epsilon) - \delta} &= (\alpha - \gamma) \left( 1 + \frac{\Psi(\epsilon) + \delta}{\mathbf{L}_1(\theta^*, t^*(\alpha - \gamma)) - \Psi(\epsilon) - \delta} \right) \\ &\le \alpha (1 + (s^*/2)^{-1} (\Psi(\epsilon) + \delta)). \end{aligned} $$
Taking the expectation makes the right-hand-side equal to zero, from which the result follows. Now, to show the increasingness, if mFCR($C_{t_1}^*$) = mFCR($C_{t_2}^*$) for $t^*(\alpha_c) < t_1 < t_2 < t^*(\bar{\alpha})$, then the above reasoning shows that

$$ (T(X, \theta^*) - t_1) \left\{ 1_{T(X, \theta^*) < t_1} - \frac{\mathbb{P}_{\theta^*} (T(X, \theta^*) < t_1)}{\mathbb{P}_{\theta^*} (T(X, \theta^*) < t_2)} 1_{T(X, \theta^*) < t_2} \right\} \le 0 $$

and has an expectation equal to 0. Hence, given that $T(X, \theta^*)$ is continuous, we derive that almost surely

$$ \mathbb{P}_{\theta^*} (T(X, \theta^*) < t_2) 1_{T(X, \theta^*) < t_1} = \mathbb{P}_{\theta^*} (T(X, \theta^*) < t_1) 1_{T(X, \theta^*) < t_2}, $$

that is, $\mathbb{P}_{\theta^*} (t_1 \le T_i^* < t_2) = 0$, which is excluded by Assumption 1. This entails mFCR($C_{t_1}^*$) < mFCR($C_{t_2}^*$).

For proving (ii), let $t > 0$. If $\mathbb{P}_{\theta^*} (T(X, \theta^*) < t) = 0$ then the result is immediate. Otherwise, we have that mFCR($C_t^*$) - $t = (\mathbb{P}_{\theta^*} (T(X, \theta^*) < t))^{-1} \mathbb{E}_{\theta^*} ((T(X, \theta^*) - t)1\{T(X, \theta^*) < t\})$. The latter is clearly not positive, and is moreover negative because $(T(X, \theta^*) - t)1\{T(X, \theta^*) < t\} \le 0$ and $\mathbb{P}_{\theta^*} (T(X, \theta^*) = t) = 0$ by Assumption 1.

For proving (iii), let $\psi_0(t) = \mathbb{E}_{\theta^*} (T(X, \theta^*) 1\{T(X, \theta^*) < t\})$ and $\psi_1(t) = \mathbb{P}_{\theta^*} (T(X, \theta^*) < t)$, the numerator and denominator of mFCR($C_t^*$) = mFCR$_t^*$, respectively. $\psi_1(t)$ is non-decreasing in $t$, with $\psi_1(0) = 0$ and $\psi_1(1) > 0$. Moreover, $\psi_0$ and $\psi_1$ are both continuous under Assumption 1. Then denote by $t_c$ the largest $t$ s.t. $\psi_1(t) = 0$. $\psi_1$ is zero on $[0, t_c]$ then strictly positive and non-decreasing on $(t_c, 1]$, and we have that $t_c = t^*(\alpha_c)$. Hence, mFCR($C_t^*$) is zero on $[0, t_c]$ then strictly positive and continuous on $(t_c, 1]$. $\square$

**Remark 3.** With the notation of the above proof, $t \mapsto \text{mFCR}(C_t^*)$ may have a discontinuity point at $t_c$ since for $t_n \xrightarrow{t_n \to t_c} t_c$, as $\psi_1(t_n) \to 0$, one does not necessarily have that $\text{mFCR}(C_t^*) \to 0$.

**Lemma 12** (Expression of plug-in procedure as a thresholding rule). For any $\alpha \in (0, 1)$, let us consider the plug-in procedure $\hat{C}_\alpha^{PI} = (\hat{\mathbf{Z}}_\alpha^{PI}, \hat{S}_\alpha^{PI})$ defined by Algorithm 2 and denote $K = |\hat{S}_\alpha^{PI}|$ the maximum of the $k \in \{0, \dots, n\}$ such that $\max(k, 1)^{-1} \sum_{j=1}^k \hat{T}_j \le \alpha$ for $\hat{T}_i = 1 - \max_q \ell_q(X_i, \hat{\theta})$, $1 \le i \le n$. Consider also $\hat{t}(\alpha)$ defined by (27). Let Assumption 1 be true and consider an estimator $\hat{\theta}$ satisfying Assumption 4. Then it holds that $\hat{t}(\alpha) = \hat{T}_{(K+1)}$ and

$$ \hat{S}_\alpha^{PI} = \{i \in \{1, \dots, n\} : \hat{T}_i < \hat{t}(\alpha)\}. $$

*Proof.* If $\hat{T}_{(K)} < \hat{T}_{(K+1)}$ then the result is immediate. Thus it suffices to show that $\hat{T}_{(K)} = \hat{T}_{(K+1)}$ occurs with probability 0. From Assumption 4 (with the countable set $\mathcal{D}$ defined therein), we have

$$ \mathbb{P}_{\theta^*} (\hat{T}_{(K)} = \hat{T}_{(K+1)}) \le \mathbb{P}_{\theta^*} \left( \bigcup_{i \ne j} \{\hat{T}_i = \hat{T}_j\} \right) \le \mathbb{P}_{\theta^*} \left( \bigcup_{\theta \in \mathcal{D}} \bigcup_{i \ne j} \{T(X_i, \theta) = T(X_j, \theta)\} \right). $$

Now, the right term is a countable union of events which are all of null probability under Assumption 1. The result follows.
This leads to

$$\mathrm{FSR}(\hat{C}) \leq \alpha + (2/s^*)[W_T(W_{t^*,\alpha}(2\delta + 8\Psi(\epsilon)/s^*)) + 4\Psi(\epsilon) + 2\delta] + \mathbb{P}((\Lambda_{\alpha,\delta,\epsilon})^c),$$

which holds true for $\delta, \epsilon$ smaller than a threshold only depending on $\theta^*$ and $\alpha$.

### A6.3 Step 3: Upper-bounding the mFSR

We apply a similar technique as for step 2. Using Lemma 7 and (A14),

$$\begin{aligned} \mathrm{mFSR}(\hat{C}) &\leq \frac{\mathbb{E}_{\theta^*}[\hat{\mathbf{M}}_0(\hat{\theta}^\sigma, \hat{t}(\alpha))1_{\Lambda_{\alpha,\delta,\epsilon}}] + \mathbb{P}((\Lambda_{\alpha,\delta,\epsilon})^c)}{\mathbb{E}_{\theta^*}[\hat{\mathbf{L}}_1(\hat{\theta}^\sigma, \hat{t}(\alpha))1_{\Lambda_{\alpha,\delta,\epsilon}}]} \\ &\leq \frac{\mathbb{E}_{\theta^*}[\hat{\mathbf{M}}_0(\hat{\theta}^\sigma, t^*(\alpha + \gamma))1_{\Lambda_{\alpha,\delta,\epsilon}}] + \mathbb{P}((\Lambda_{\alpha,\delta,\epsilon})^c)}{\mathbb{E}_{\theta^*}[\hat{\mathbf{L}}_1(\hat{\theta}^\sigma, t^*(\alpha - \gamma))1_{\Lambda_{\alpha,\delta,\epsilon}}]} \end{aligned}$$

Now using a concentration argument on $\Lambda_{\alpha,\delta,\epsilon} \subset \Gamma_{1,\delta,t^*(\alpha-\gamma)} \cap \Upsilon_{0,\delta,t^*(\alpha+\gamma)}$, we have

$$\begin{aligned} \mathrm{mFSR}(\hat{C}) &\leq \frac{\mathbb{E}_{\theta^*}[(\mathbf{M}_0(\hat{\theta}^\sigma, t^*(\alpha + \gamma)) + \delta)1_{\Lambda_{\alpha,\delta,\epsilon}}] + \mathbb{P}((\Lambda_{\alpha,\delta,\epsilon})^c)}{\mathbb{E}_{\theta^*}[(\mathbf{L}_1(\hat{\theta}^\sigma, t^*(\alpha - \gamma)) - \delta)1_{\Lambda_{\alpha,\delta,\epsilon}}]} \\ &\leq \frac{\mathbf{M}_0(\theta^*, t^*(\alpha + \gamma)) + 3\Psi(\epsilon) + \delta + \mathbb{P}((\Lambda_{\alpha,\delta,\epsilon})^c)}{\mathbf{L}_1(\theta^*, t^*(\alpha - \gamma)) - \Psi(\epsilon) - \delta - \mathbb{P}((\Lambda_{\alpha,\delta,\epsilon})^c)} \\ &= \frac{\mathbf{L}_0(\theta^*, t^*(\alpha + \gamma)) + 3\Psi(\epsilon) + \delta + \mathbb{P}((\Lambda_{\alpha,\delta,\epsilon})^c)}{\mathbf{L}_1(\theta^*, t^*(\alpha - \gamma)) - \Psi(\epsilon) - \delta - \mathbb{P}((\Lambda_{\alpha,\delta,\epsilon})^c)}, \end{aligned}$$

by using Lemma 10 and that $\mathbf{M}_0(\theta^*, t) = \mathbf{L}_0(\theta^*, t)$ by definition. Letting $x = \mathbf{L}_0(\theta^*, t^*(\alpha + \gamma)) + 3\Psi(\epsilon) + \delta, y = \mathbf{L}_1(\theta^*, t^*(\alpha - \gamma)) - \Psi(\epsilon) - \delta$ and $u = \mathbb{P}((\Lambda_{\alpha,\delta,\epsilon})^c)$, we have obtained the bound $(x + u)/(y - u)$, which has to be compared with the FSR bound (A15), which reads $x/y + u$. Now, when $y \in [0, 1], x \geq 0, x/y \leq 2, u/y \leq 1/2, y - u \geq s^*/2$, we have

$$(x+u)/(y-u) \leq \frac{x/y}{1-u/y} + (2/s^*)u \leq x/y(1+2u/y) + (2/s^*)u \leq x/y + (10/s^*)u.$$

As a result, for $\epsilon, \delta$ small enough, and $\mathbb{P}((\Lambda_{\alpha,\delta,\epsilon})^c) \leq s^*/4$, we obtain the same bound as for the FSR, with $\mathbb{P}((\Lambda_{\alpha,\delta,\epsilon})^c)$ replaced by $(10/s^*)\mathbb{P}((\Lambda_{\alpha,\delta,\epsilon})^c)$.

### A6.4 Step 4: Lower-bounding the selection rate

In Step 3, when bounding the mFSR, we derived a lower bound for the denominator of the mFSR, that is, $\mathbb{E}_{\theta^*}(|\hat{S}|)$. It reads

$$\begin{aligned} n^{-1}\mathbb{E}_{\theta^*}(|\hat{S}|) &\geq \mathbf{L}_1(\theta^*, t^*(\alpha - \gamma)) - \Psi(\epsilon) - \delta - \mathbb{P}((\Lambda_{\alpha,\delta,\epsilon})^c) \\ &\geq \mathbf{L}_1(\theta^*, t^*(\alpha)) - W_T(t^*(\alpha - t^*(\alpha - \gamma)) - \Psi(\epsilon) - \delta - \mathbb{P}((\Lambda_{\alpha,\delta,\epsilon})^c) \\ &\geq n^{-1}\mathbb{E}_{\theta^*}(|S^*_{t^*(\alpha)}|) - W_T(W_{t^*,\alpha}(\gamma)) - \Psi(\epsilon) - \delta - \mathbb{P}((\Lambda_{\alpha,\delta,\epsilon})^c), \end{aligned}$$

by using (A3) and (A5). Now consider another procedure $C = (\hat{Z}, S)$ that controls the mFSR at level $\alpha$, that is, mFSR($C$) $\leq \alpha$. By Lemma 6, we then have $\mathbb{E}_{\theta^*}(|S^*_{t^*(\alpha)}|) \geq \mathbb{E}_{\theta^*}(|S|)$.
**Lemma 13.** We have for all $\theta \in \Theta$,

$$ \sup_{t \in [0,1]} |\mathbf{L}_1(\theta, t) - \mathbf{L}_1(\theta^*, t)| \le \Psi(\|\theta^* - \theta\|); \quad (33) $$

$$ \sup_{t \in [0,1]} |\mathbf{L}_0(\theta, t) - \mathbf{L}_0(\theta^*, t)| \le 2\Psi(\|\theta^* - \theta\|); \quad (34) $$

$$ \sup_{t \in [t^*((\alpha+\alpha_c)/2),1]} |\mathbf{L}(\theta, t) - \mathbf{L}(\theta^*, t)| \le 3\Psi(\|\theta^* - \theta\|)/s^*; \quad (35) $$

$$ \sup_{t \in [0,1]} |\mathbf{M}_0(\theta, t) - \mathbf{M}_0(\theta^*, t)| \le 3\Psi(\|\theta^* - \theta\|); \quad (36) $$

$$ \sup_{t \in [t^*((\alpha+\alpha_c)/2),1]} |\mathbf{M}(\theta, t) - \mathbf{M}(\theta^*, t)| \le 4\Psi(\|\theta^* - \theta\|)/s^*; \quad (37) $$

where $\alpha \in (\alpha_c, \bar{\alpha}]$ and $s^* > 0$ is given by (26). In addition, for all $\theta \in \Theta$ and $t, t' \in [0, 1]$,

$$ |\mathbf{L}_0(\theta, t) - \mathbf{L}_0(\theta, t')| \le 4\Psi(\|\theta^* - \theta\|) + W_T(|t - t'|). \quad (38) $$

*Proof.* Fix $\theta \in \Theta$ and $t \in [0, 1]$. We have for any $\delta > 0$,

$$ \begin{aligned} & |\mathbb{P}_{\theta^*}(T(X, \theta) < t) - \mathbb{P}_{\theta^*}(T(X, \theta^*) < t)| \\ & \le (\mathbb{P}_{\theta^*}(T(X, \theta^*) < t + \delta) - \mathbb{P}_{\theta^*}(T(X, \theta^*) < t)) \vee (\mathbb{P}_{\theta^*}(T(X, \theta^*) < t) - \mathbb{P}_{\theta^*}(T(X, \theta^*) < t - \delta)) \\ & \quad + \mathbb{P}_{\theta^*}(|T(X, \theta^*) - T(X, \theta)| > \delta) \\ & \le W_T(\delta) + \mathbb{E}_{\theta^*}(|T(X, \theta^*) - T(X, \theta)|)/\delta. \end{aligned} $$

In addition, by definition (5),

$$ \begin{aligned} |T(X, \theta^*) - T(X, \theta)| &\le \max_{1 \le q \le Q} \ell_q(X, \theta^*) - \max_{1 \le q \le Q} \ell_q(X, \theta)| \\ &\le \max_{1 \le q \le Q} |\ell_q(X, \theta^*) - \ell_q(X, \theta)|. \end{aligned} $$

Hence,

$$ |\mathbb{P}_{\theta^*}(T(X, \theta) < t) - \mathbb{P}_{\theta^*}(T(X, \theta^*) < t)| \le \inf_{\delta \in (0,1)} \{W_T(\delta) + W_\ell(\|\theta^* - \theta\|)/\delta\} = \Psi(\|\theta^* - \theta\|), $$

which establishes (33).

Next, we have

$$ \begin{aligned} & \mathbf{L}_0(\theta, t) - \mathbf{L}_0(\theta^*, t) \\ &= \mathbb{E}_{\theta^*}[T(X, \theta)(1_{T(X,\theta)<t} - 1_{T(X,\theta^*)<t}) + 1_{T(X,\theta^*)<t}(T(X, \theta) - T(X, \theta^*)]| \\ &\le t|\mathbb{P}_{\theta^*}(T(X, \theta) < t) - \mathbb{P}_{\theta^*}(T(X, \theta^*) < t)| + |\mathbb{E}_{\theta^*}[1_{T(X,\theta^*)<t}(T(X, \theta) - T(X, \theta^*))]| \\ &\le |\mathbb{P}_{\theta^*}(T(X, \theta) < t) - \mathbb{P}_{\theta^*}(T(X, \theta^*) < t)| + \mathbb{E}_{\theta^*}|T(X, \theta) - T(X, \theta^*)| \le 2\Psi(\|\theta^* - \theta\|) \end{aligned} $$

By exchanging the role of $\theta$ and $\theta^*$ in the above reasoning, the same bound holds for $\mathbf{L}_0(\theta^*, t) - \mathbf{L}_0(\theta, t)$, which gives (34). To prove (35), we use for any $t \in [t^*((\alpha+\alpha_c)/2), 1]$,

$$ \begin{aligned} & \left|\frac{\mathbf{L}_0(\theta, t)}{\mathbf{L}_1(\theta, t)} - \frac{\mathbf{L}_0(\theta^*, t)}{\mathbf{L}_1(\theta^*, t)}\right| \\ &\le \left|\frac{\mathbf{L}_0(\theta, t) - \mathbf{L}_0(\theta^*, t)}{\mathbf{L}_1(\theta^*, t)}\right| + \mathbf{L}_0(\theta, t) \left|\frac{1}{\mathbf{L}_1(\theta^*, t)} - \frac{1}{\mathbf{L}_1(\theta, t)}\right| \\ &\le 2\Psi(\|\theta^* - \theta\|)/s^* + \frac{1}{\mathbf{L}_1(\theta^*, t)} \frac{\mathbf{L}_0(\theta, t)}{\mathbf{L}_1(\theta, t)} |\mathbb{P}_{\theta^*}(T(X, \theta^*) < t) - \mathbb{P}_{\theta^*}(T(X, \theta) < t)| \\ &\le 3\Psi(\|\theta^* - \theta\|)/s^*, \end{aligned} $$
## A6.5 Step 5: Concentration

Finally, we bound $P((\Lambda_{\alpha, \delta, \epsilon})^c)$ by using Lemma 11 with $x = (1 + 2c) \sqrt{\frac{\log n}{n}}$ (with $c$ defined in Lemma 11). This gives for $\delta = 2x/s^*$, and $n \ge (2e)^3$

$$P((\Lambda_{\alpha, \delta, \epsilon})^c) \le 5/n^2 + P(\Omega_{\epsilon}^c).$$

# APPENDIX B. PROOFS OF LEMMAS

## B1 Proof of Lemma 1

The clustering risk of $\hat{Z}$ is given by

$$ 
\begin{aligned}
R(\hat{Z}) &= E_{\theta^*} \left( \min_{\sigma \in [Q]} E_{\theta^*} \left( n^{-1} \sum_{i=1}^n 1\{Z_i \neq \sigma(\hat{Z}_i)\} \middle| \mathbf{X} \right) \right) \\
&= E_{\theta^*} \left( \min_{\sigma \in [Q]} n^{-1} \sum_{i=1}^n P_{\theta^*}(Z_i \neq \sigma(\hat{Z}_i) \middle| \mathbf{X}) \right) \\
&\ge E_{\theta^*} \left( \min_{\hat{Z}} n^{-1} \sum_{i=1}^n P_{\theta^*}(Z_i \neq \hat{Z}_i \middle| \mathbf{X}) \right),
\end{aligned}
$$

where, by independence, the minimum in the lower bound is achieved for the Bayes clustering. Thus, $R(\hat{Z}) \ge n^{-1} \sum_{i=1}^n E_{\theta^*}(T_i^*)$. Moreover, $n^{-1} \sum_{i=1}^n E_{\theta^*}(T_i^*) \ge R(\hat{Z}^*)$, since

$$ 
\begin{aligned}
R(\hat{Z}^*) &= E_{\theta^*} \left( \min_{\sigma \in [Q]} n^{-1} \sum_{i=1}^n P_{\theta^*}(Z_i \neq \sigma(\hat{Z}_i^*) \middle| \mathbf{X}) \right) \\
&\le E_{\theta^*} \left( n^{-1} \sum_{i=1}^n P_{\theta^*}(Z_i \neq \hat{Z}_i^* \middle| \mathbf{X}) \right).
\end{aligned}
$$

Thus, $\min_{\hat{Z}} R(\hat{Z}) = R(\hat{Z}^*)$ and the proof is completed.

## B2 Proof of Lemma 2

Following the reasoning of the proof of Lemma 1, we have

$$ 
\begin{aligned}
\text{FSR}_{\theta^*}(C) &= E_{\theta^*} \left( \min_{\sigma \in [Q]} E_{\theta^*} \left( \frac{\sum_{i \in S} 1\{Z_i \neq \sigma(\hat{Z}_i^*)\}}{\max(|S|, 1)} \middle| \mathbf{X} \right) \right) \\
&= E_{\theta^*} \left( \frac{\sum_{i \in S} T_i^*}{\max(|S|, 1)} \right).
\end{aligned}
$$

## B3 Proof of Lemma 3

By definition, we have

$$ \text{FSR}(\hat{C}_{\alpha}^{\text{PI}}) = E_{\theta^*} \left( \min_{\sigma \in [Q]} E_{\theta^*} \left( \frac{\sum_{i=1}^n 1_{Z_i \neq \sigma(\hat{Z}_i^{\text{PI}}(\mathbf{X}))} 1\{i \in \hat{S}^{\text{PI}}(\mathbf{X})\}}{\max(|\hat{S}^{\text{PI}}(\mathbf{X})|, 1)} \middle| \mathbf{X} \right) \right). $$