# Invidious Comparisons:
Ranking and Selection as Compound Decisions

Jiaying Gu  and  Roger Koenker

###### Abstract.

There is an innate human tendency, one might call it the “league table
mentality,” to construct rankings. Schools, hospitals, sports teams, movies,
and myriad other objects are ranked even though their inherent multi-
dimensionality would suggest that – at best – only partial orderings were
possible. We consider a large class of elementary ranking problems in which we
observe noisy, scalar measurements of merit for $n$ objects of potentially
heterogeneous precision and are asked to select a group of the objects that
are “most meritorious.” The problem is naturally formulated in the compound
decision framework of Robbins’s (1956) empirical Bayes theory, but it also
exhibits close connections to the recent literature on multiple testing. The
nonparametric maximum likelihood estimator for mixture models (Kiefer and
Wolfowitz (1956)) is employed to construct optimal ranking and selection
rules. Performance of the rules is evaluated in simulations and an application
to ranking U.S kidney dialysis centers.

Version: . This paper was presented as the Walras-Bowley Lecture at the 2020
World Congress of the Econometric Society, and is dedicated to the memory of
Larry Brown who introduced us to empirical Bayes methods. We thank Michael
Gilraine, Keisuke Hirano, Robert McMillan, Stanislav Volgushev and Sihai Dave
Zhao for useful discussions. Jiaying Gu acknowledges financial support from
Social Sciences and Humanities Research Council of Canada.

##  1\. Introduction

In the wake of Wald’s seminal monograph on statistical decision theory there
was a growing awareness that the Neyman-Pearson testing apparatus was
inadequate for many important statistical tasks. Ranking and selection
problems featured prominently in this perception. Motivated by a suggestion of
Harold Hotelling, Bahadur50 studied selection of the best of several Gaussian
populations. Assuming that sample means were observed for each of $K$
populations with means, $\theta_{k}$ and common variance, the problem of
selecting the best population, $\theta^{*}$, was formulated as choosing
weights $z_{1},\cdots,z_{K}$ to minimize,

$L(\theta,z)=\theta^{*}-\sum_{k=1}^{K}z_{k}\theta_{k}/\sum_{k=1}^{K}z_{k}.$

Bahadur showed that it was minimax to select only the population with the
largest sample mean, that is to choose $z_{i}^{*}=1$ if
$\bar{X}_{i}=\max\\{\bar{X}_{1},\cdots,\bar{X}_{K}\\}$ and $z_{i}^{*}=0$
otherwise, thereby clearly demonstrating that procedures that did preliminary
tests of equality of means and then chose $z_{i}>0$ for several or even all of
the populations when tests failed to reject were inadmissible. This finding
was reinforced in BR50 who focused on the two-sample setting but relaxed the
common variance assumption. In related work, Bechhofer and Gupta56 sought to
optimize the number of selected populations as well as their identities, see
GuptaPan and BKS for extensive reviews of subsequent developments.

GoelRubin pioneered the hierarchical Bayesian approach to selection that has
been adopted by numerous authors in the ensuing decades, early on by
BergerDeeley and LairdLouis. Portnoy showed that rankings based on best linear
predictors were optimal in Gaussian multivariate variance components models,
but cautioned that departures from normality could easily disrupt this
optimality. A notable feature of the hierarchical model paradigm is the
recognition that sample observations may exhibit heterogeneous precision; this
is typically accounted for by assuming known variances for observed sample
means. As ranking and selection methods became increasingly relevant in
genomic applications there has been renewed interest in loss functions and
linkages to the burgeoning literature on multiple testing. Our perspective is
informed by recent developments in the nonparametric estimation of mixture
models and its relevance for a variety of compound decision problems. This
approach seeks to reduce the reliance on Gaussian distributional assumptions
that pervades the earlier literature. As we have argued elsewhere, ProbRob,
and Minimalist nonparametric empirical Bayes methods offer powerful
complementary methods to more conventional parametric hierarchical Bayes for
multiple testing and compound decision problems. Our primary objective in this
paper is to elaborate this assertion for ranking and selection applications.
Throughout we try to draw parallels and contrasts with the literature on
multiple testing. We will restrict our attention to settings where we observe
a scalar estimate of an unobserved latent quality measure accompanied by some
measure of its precision, thereby evading more complex multivariate settings,
as in accuracy who employ quantile regression methods.

An important motivation for revived interest in ranking and selection problems
in econometrics has been the influential work of Chetty and his collaborators
on teacher evaluation and geographic mobility in the U.S. This has stimulated
the important recent work of MRSW proposing new resampling methods for
constructing confidence sets for ranking and selectionfor a finite population.
In contrast to this inferential approach we focus instead on the complementary
perspective of compound decision making, constructing decision rules for
selecting the best, or worst, populations subject to control of the expected
number of elements selected and among those selected, the expected proportion
of false discovery.

Before proceeding it is important to acknowledge that despite its universal
appeal and application there is something inherently futile about many ranking
and selection problems as intimated by our title. If our latent measure of
true quality is Gaussian, as assumed in virtually all of the econometric
applications of the selection problem, and we wish to select the top ten
percent of individuals given that their true quality is contaminated by
Gaussian noise. We will see that conventional linear shrinkage as embodied in
the classical James-Stein formula can improve performance considerably over
naive maximum likelihood (fixed effects) procedures, and some further
improvement is possible by carefully tailoring the decision rules for tail
probability loss, however we find that even oracle decision rules that
incorporate complete knowledge of the precise distributional features of the
problem cannot achieve better than about even odds that selected individuals
have latent ability above the selection thresholds. When the latent
distribution of ability is heavier tailed then selection becomes somewhat
easier, and more refined selection rules are more advantageous, but as we will
show the selection problem still remains quite challenging.

Thus, a secondary objective of the paper is to add another cautionary voice to
those who have already questioned the reliability of existing ranking and
selection methods. A critical overview of the role of ranking and selection in
public policy applications is provided by LeagueTables. It is widely
acknowledged that league tables as currently employed can be a pernicious
influence on policy, a viewpoint underscored in GelmanPrice. While much of
this criticism can be attributed to inadequate data collection, we believe
that there is also room for methodological improvements.

Section 2 provides a brief overview of compound decision theory and describes
nonparametric methods for estimation of Gaussian mixture models. Section 3
introduces a basic framework for our approach to ranking and selection in a
setting with homogeneous precision of the observed measurements. In Section 4
we introduce heterogeneous precision of known form, and Section 5 considers
settings in which the joint distribution of the observed measurements and
their precision determines the form of the ranking and selection rules.
Optimal ranking and selection rules are derived in each of these sections
under the assumption that the form of the mixing distribution of the
unobserved, latent quality of the observations is known. Section 6 introduces
feasible ranking and selection rules and conditions under which they attain
the same asymptotic performance as the optimal rules. Section 7 then compares
several feasible ranking and selection methods, some that ignore the compound
decision structure of the problem, some that employ parametric empirical Bayes
methods and some that rely on nonparametric empirical Bayes methods. Finally,
Section 8 describes an empirical application on evaluating the performance of
medical dialysis centers in the United States. Proofs of all formal results
are collected in Appendix A.

##  2\. The Compound Decision Framework

Robbins.51 posed a challenge to the nascent minimax decision theory of Wald50:
Suppose we observe independent Gaussian realizations,
$Y_{i}\sim\mathcal{N}(\theta_{i},1),\;i=1,\cdots,n$ with means $\theta_{i}$
taking either the value $+1$ or $-1$. We are asked to estimate the $n$-vector
$\theta=(\theta_{1},\cdots,\theta_{n})$ subject to mean absolute error loss,

$L(\hat{\theta},\theta)=n^{-1}\sum_{i=1}^{n}|\hat{\theta}_{i}-\theta_{i}|.$

When $n=1$ Robbins shows that the minimax decision rule is
$\delta(y)=\mbox{sgn}(y)$; in the least favorable variant of the problem
malevolent nature chooses $\pm 1$ with equal probability, and the optimal
response is to estimate $\theta_{i}=+1$ when $Y_{i}$ is positive, and
$\theta_{i}=-1$ otherwise. Robbins goes on to show that when $n>1$, this rule
remains minimax, each coordinate is treated independently as if viewed in
complete isolation. This is also the maximum likelihood estimator, and may be
viewed in econometrics terms as a classical fixed-effects estimator. But is it
at all reasonable?

Doesn’t our sample convey information about the relative frequency of $\pm 1$
that might potentially contradict the pessimistic presumption of the minimax
rule? If we happened to know the unconditional probability,
$p=\mathbb{P}(\theta_{i}=1)$, then the conditional probability that $\theta=1$
given $Y_{i}=y$, is given by,

$P(\theta=1|y)=\frac{p\varphi(y-1)}{p\varphi(y-1)+(1-p)\varphi(y+1)}.$

We should guess $\hat{\theta}_{i}=1$ if this probability exceeds 1/2, giving
us the revised decision rule,

$\delta_{p}(y)=\mbox{sgn}(y-\tfrac{1}{2}\log((1-p)/p)).$

Each observed, $Y_{i}$, is modified by a simple logistic perturbation before
computing the sign. Our observed random sample, $y=(y_{1},\cdots,y_{n})$, is
informative about $p$, we have the log likelihood,

$\ell_{n}(p|y)=\sum_{i=1}^{n}\log(p\varphi(y_{i}-1)+(1-p)\varphi(y_{i}+1)),$

which could be augmented by a prior of some form, if desired, to obtain a
posterior mean for $p$ and a plug-in Bayes rule for estimating each of the
$\theta_{i}$’s. The Bayes risk of this procedure is substantially less than
the minimax risk when $p\neq 1/2$ and is asymptotically equivalent to the
minimax risk when $p=1/2$. This is the first principle of compound decision
theory: borrowing strength across an entire ensemble of related decision
problems yields improved collective performance.

What happens when we relax the restriction on the support of the $\theta$’s
and allow support on the whole real line? We now have a general Gaussian
mixture setting where the observed $Y_{i}$’s have marginal density given by
the convolution, $f=\varphi*G$, that is,

$f(y)=\int\varphi(y-\theta)dG(\theta),$

and instead of merely needing to estimate one probability we need an estimate
of an entire distribution function, $G$. KW, anticipated by an abstract of
R50, established that the nonparametric maximum likelihood estimator (NPMLE),

$\hat{G}=\mbox{argmin}_{G\in\mathcal{G}}\\{-\sum_{i=1}^{n}\log
f(y_{i})\;|\;f(y_{i})=\int\varphi(y_{i}-\theta)dG(\theta)\\}$

where $\mathcal{G}$ is the space of probability measures on $\mathbb{R}$, is a
consistent estimator of $G$. This is an infinite dimensional convex
optimization problem with a strictly convex objective subject to linear
constraints. See Lindsay.95 and KM for further details on the geometry and
computational aspects of the NPMLE problem. Solutions are atomic as a
consequence of the Carathéodory theorem, but until quite recently little was
known about the growth rate of the number of atoms characterizing the
solutions. PW20 have established that for $G$ with sub-Gaussian tails the
cardinality of its support, i.e. the number of atoms, of $\hat{G}$, grows like
$\mathcal{O}(\log n)$. In this respect the NPMLE shares a property of several
shape-constrained density estimators. Polyanskiy and Wu call such procedures
self-regularizing: without any further penalization, the maximum likelihood
estimator in such cases automatically selects a highly parsimonious estimator.
Thus, in sharp contrast to the well-known difficulties with maximum likelihood
for finite dimensional mixture models, or with Gaussian deconvolution
employing Fourier methods, there are no awkward problems of tuning parameter
selection for the NPMLE.

Having seen that the upper bound on the complexity of the NPMLE $\hat{G}$ was
only $\mathcal{O}(\log n)$, one might begin to wonder whether
$\mathcal{O}(\log n)$ mixtures are “complex enough” to adequately represent
the process that generated our observed data. PW20 also address this concern:
they note that for any sub-Gaussian $G$, there exists a discrete distribution,
$G_{k}$, with $k=\mathcal{O}(\log n)$ atoms, such that for
$f_{k}=\varphi*G_{k}$, the total variation distance, $TV(f,f_{k})=o(1/n)$, and
consequently there is no statistical justification for considering estimators
of $G$ whose complexity grows more rapidly than $\mathcal{O}(\log n)$. This
observation is related to recent literature on generative adversarial
networks, e.g. AIMM20, that target models and estimators that, when simulated,
successfully mimic observed data. This viewpoint also underlies the innovative
monograph of Davies.

Other nonparametric maximum likelihood estimators for $G$ are potentially also
of interest. E16 and E19 has proposed an elegant log-spline sieve approach
that yields a smooth estimate of $G$; this has advantages especially from an
inferential perspective, at the cost of reintroducing the task of selecting
tuning parameters. An early proposal of LL91 merged parametric empirical Bayes
estimation of $G$ with an EM step that pulled the parametric estimate back
toward the NPMLE.

Given an estimate, $\hat{G}$, it is straightforward to compute posterior
distributions for each of sample observations, or for that matter, for out-of-
sample observations. In effect, we have estimated the prior, as in Robbins.51
binary means problem, but we have ignored the variability of $\hat{G}$ when we
adopt plug-in procedures that use it. This may account for the improved
performance of smoothed estimates of $G$ in certain inferential problems, as
conjectured in K20. In the sequel we will compare ranking and selection
procedures based on various functionals of these posterior distributions. A
leading example is the posterior mean, but ranking and selection problems
suggest other functionals of potential interest.

If we are asked to estimate the $\theta_{i}$’s subject to quadratic loss, and
assuming standard Gaussian noise, the Bayes rule is given by the posterior
mean,

(2.1) $\delta(y)=\mathbb{E}(\theta|y)=y+f^{\prime}(y)/f(y).$

efron.11 refers to this as Tweedie’s formula, it appears in Robbins.56
credited to M.C.K. Tweedie. Appendix A of Haydn provides an elementary
derivation. The nonlinear shrinkage term takes a particularly simple affine
form when $G$ happens to be Gaussian, since in this case $f$ is itself also
Gaussian and the formula reduces to well-known linear shrinkage variants of
classical Stein rules.

A striking feature of the Tweedie formula is that the dependence on $G$ of the
nonlinear shrinkage term is hidden in the log derivative which depends only on
the marginal density of the $Y_{i}$’s. This creates a temptation to believe
that estimation of $G$ is really superfluous, that it would suffice to simply
estimate $f$ instead. This $f$-modeling temptation should be resisted for a
variety of reasons elucidated in E19 not the least of which is that it is
difficult to account for the full probabilistic structure of the problem when
estimating $f$ directly. For example, it is known that the Bayes rule for the
posterior mean in the Gaussian case is monotone increasing in $y$, a fact not
easily incorporated into the conventional kernel density estimation approach.
KM consider this point in more detail. More importantly, estimation of $G$
opens the way to improved methods for attacking many practical decision
problems including ranking and selection.

We have focused in this brief overview on compound decision problems for
Gaussian location mixtures, however the NPMLE is adaptable to a wide variety
of other mixture problems. Some of these other options are described in
REBayesVig and are implemented in the R package REBayes.

##  3\. Homogeneous Variances

Suppose that you are given real-valued measurements, $y_{i}:i=1,2,\cdots,n$ of
some attribute like test score performance for students or their teachers,
survival rates for hospital surgical procedures, etc., and are told that the
measurements are exchangeable and approximately Gaussian with unknown means
$\theta_{i}$ and known variances $\sigma_{i}^{2}$ assumed provisionally to
take the same value $\sigma^{2}$. Your task, should you decide to accept it,
is to choose a group of size not to exceed $\alpha n$ of the elements with the
largest $\theta_{i}$’s. One’s first inclination might be to view each $y_{i}$
as the maximum likelihood estimate for the corresponding $\theta_{i}$, and
select the $\alpha n$ largest observed values, but the compound decision
framework suggests that it would be better to treat the problems as an
ensemble. A second natural inclination might be to compute posterior means of
the $\theta$’s with some linear or nonlinear shrinkage rule, rank them and
select the $\alpha$ best, but we will see that this too may be questionable.

###  3.1. Posterior Tail Probability

A natural alternative to ranking by the posterior means is to rank by
posterior tail probabilities. Let $\theta_{\alpha}=G^{-1}(1-\alpha)$, and
define, $v_{\alpha}(y):=\mathbb{P}(\theta\geq\theta_{\alpha}|Y=y)$, then
ranking by posterior tail probability gives the decision rule,

$\delta(y)=\mymathbb 1\\{v_{\alpha}(y)\geq\lambda_{\alpha}\\}$

where $\lambda_{\alpha}$ is chosen so that
$\mathbb{P}(v_{\alpha}(Y)\geq\lambda_{\alpha})=\alpha$. This ranking criterion
has been proposed by HN motivated as a ranking device for a fixed quantile
level $\alpha$. It can be interpreted in multiple testing terms:
$1-v_{\alpha}(y)$ is the local false discovery rate of ETST and storey02, for
testing the hypothesis $H_{0}:\theta<\theta_{\alpha}$ vs.
$H_{A}:\theta\geq\theta_{\alpha}$. To see this, let $h_{i}$ be a binary random
variable $h_{i}=\mymathbb 1\\{\theta_{i}\geq\theta_{\alpha}\\}$, the loss
function for observation $i$ is

$L(\delta_{i},\theta_{i})=\lambda\mymathbb 1\\{h_{i}=0,\delta_{i}=1)+\mymathbb
1\\{h_{i}=1,\delta_{i}=0\\}.$

The compound Bayes risk is,

$\mathbb{E}[\sum_{i=1}^{n}L(\delta_{i},\theta_{i})]=n[\alpha+\int\delta(y)[(1-\alpha)\lambda
f_{0}(y)-\alpha f_{1}(y)]dy]$

where
$f_{0}(y)=(1-\alpha)^{-1}\int_{-\infty}^{\theta_{\alpha}}\varphi(y|\theta,\sigma^{2})dG(\theta)$
and
$f_{1}(y)=\alpha^{-1}\int_{\theta_{\alpha}}^{+\infty}\varphi(y|\theta,\sigma^{2})dG(\theta)$,
$\varphi(y|\theta,\sigma^{2})=\phi((y-\theta)/\sigma)/\sigma$ and
$\phi(\cdot)$ is the standard normal density. The Bayes rule for a fixed
$\lambda$ is

$\delta(y_{i})=\mymathbb
1\Big{\\{}v_{\alpha}(y_{i})\geq\frac{\lambda}{1+\lambda}\Big{\\}}$

where $v_{\alpha}(y)=\alpha
f_{1}(y)/f(y)=\mathbb{P}(\theta\geq\theta_{\alpha}|Y=y)$, and
$f(y)=(1-\alpha)f_{0}(y)+\alpha f_{1}(y).$ Provided that $v_{\alpha}(y)$ is
monotone in $y$ a unique $\lambda_{\alpha}$ can be found such that
$\mathbb{P}(\delta(Y)=1)=\mathbb{P}(v_{\alpha}(Y)\geq\lambda_{\alpha}/(1+\lambda_{\alpha}))=\alpha$.

######  Lemma 3.1.

For fixed $\alpha$, assuming
$\mathbb{E}_{\theta|Y}[\nabla_{y}\log\varphi(y|\theta,\sigma^{2})|Y]<\infty$,
$v_{\alpha}(y)$ is monotone in $y$ and the sets
$\Omega_{\alpha}:=\\{Y:v_{\alpha}(Y)\geq\lambda_{\alpha}/(1+\lambda_{\alpha})\\}$
have a nested structure, that is if $\alpha_{1}>\alpha_{2}$, then
$\Omega_{\alpha_{2}}\subseteq\Omega_{\alpha_{1}}$.

Any implementation of such a Bayes rule requires an estimate of the mixing
distribution, $G$, or something essentially equivalent that would enable us to
compute the local false discovery rates $v_{\alpha}(y)$ and the cut-off
$\theta_{\alpha}$. The NPMLE, or perhaps a smoothed version of it, will
provide a natural $\hat{G}$ for this task.

###  3.2. Posterior Tail Expectation and Other Losses

Rather than assessing loss by simply counting misclassifications we might
consider weighting such misclassifications by the magnitude of $\theta$, for
example,

$L(\delta_{i},\theta_{i})=\sum_{i=1}^{n}(1-\delta_{i})\mymathbb
1(\theta_{i}\geq\theta_{\alpha})\theta_{i}.$

This presumes, of course, that we have centered the distribution $G$ in some
reasonable way, perhaps by forcing the mean or median to be zero. Minimizing
with respect to $\delta$ subject to the constraint that
$\mathbb{P}(\delta(Y)=1)=\alpha$ leads to the Lagrangian,

$\underset{\delta}{\min}\int\int(1-\delta(y))\mymathbb
1\\{\theta\geq\theta_{\alpha}\\}\theta\varphi(y|\theta,\sigma^{2})dG(\theta)dy+\lambda\Big{[}\int\int\delta(y)\varphi(y|\theta,\sigma^{2})dG(\theta)dy-\alpha\Big{]}$

which is equivalent to

$\displaystyle\underset{\delta}{\min}$ $\displaystyle\int\int\mymathbb
1\\{\theta\geq\theta_{\alpha}\\}(\theta-\lambda)\varphi(y|\theta,\sigma^{2})dG(\theta)dy$
$\displaystyle-\int\delta(y)\Big{[}\int\mymathbb
1\\{\theta\geq\theta_{\alpha}\\}(\theta-\lambda)\varphi(y|\theta,\sigma^{2})dG(\theta)-\int\lambda\mymathbb
1\\{\theta<\theta_{\alpha}\\}\varphi(y|\theta,\sigma^{2})dG(\theta)\Big{]}dy.$

Ignoring the first term since it doesn’t depend upon $\delta$, the Oracle
Bayes rule, becomes, choose $\delta(y)=1$ if,

$\frac{\int\mymathbb
1\\{\theta\geq\theta_{\alpha}\\}\theta\varphi(y|\theta,\sigma^{2})dG(\theta)}{\int\varphi(y|\theta,\sigma^{2})dG(\theta)}\geq\lambda,$

with $\lambda$ chosen so that $\mathbb{P}(\delta(Y)=1)=\alpha$. Such criteria
are closely related to expected shortfall criteria appearing in the literature
on risk assessment. Again, the NPMLE can be employed to construct feasible
posterior ranking criteria.

Several other loss functions are considered by LLPR06 including some based on
global alignment of the ranks. However, if the objective of the ranking
exercise is eventually to select the best, or worst, elements it seems
difficult to rationalize global loss functions of this type.

###  3.3. False Discovery and the $\alpha$-Level

Although our loss functions yield distinct criteria for ranking, their
decision rules lead to the same selections when the precision of the
measurements is homogeneous. When variances are homogeneous there is a cut-
off, $\eta_{\alpha}$, and a decision rule,
$\delta_{\alpha}(Y)=1(Y\geq\eta_{\alpha})$, determining a common selection for
all decision rules.

######  Lemma 3.2.

For fixed $\alpha$ and homogeneous variance, posterior mean, posterior tail
probability and posterior tail expectation all yield the same ranking and
therefore the same selection.

The marginal false discovery rate for this selection in our Gaussian mixture
setting is,

$\displaystyle mFDR$
$\displaystyle=\mathbb{P}(\theta<\theta_{\alpha}|\delta_{\alpha}(Y)=1)$
$\displaystyle=\alpha^{-1}\mathbb{P}(\theta<\theta_{\alpha},\theta+\sigma
Z\geq\eta_{\alpha})$
$\displaystyle=\alpha^{-1}\int_{-\infty}^{\theta_{\alpha}}\Phi((\theta-\eta_{\alpha})/\sigma)dG(\theta),$

where $Z$ is standard Gaussian. The marginal false non-discovery rate is,

$\displaystyle mFNR$
$\displaystyle=\mathbb{P}(\theta\geq\theta_{\alpha}|\delta_{\alpha}(Y)=0)$
$\displaystyle=(1-\alpha)^{-1}\int_{\theta_{\alpha}}^{\infty}\Phi((\eta_{\alpha}-\theta)/\sigma)dG(\theta),$

Figure 3.1 shows the false discovery rate, and false non-discovery rate for a
range of capacity constraints, $\alpha$, when the mixing distribution, $G$, is
standard Gaussian and $\sigma^{2}=1$. In this case, the cut-off value
$\eta_{\alpha}$ is the $(1-\alpha)$ quantile of $\mathcal{N}(0,2)$. For
selecting individuals at the top $\alpha$ quantile, the false discovery rate
remains alarmingly high especially for smaller $\alpha$, implying that the
selected set may consist of a very high proportion of false discoveries. For
example, for $\alpha=0.10$ the proportion of selected observations with
$\theta$ below the threshold $\theta_{\alpha}$ is slightly greater than 50
percent.

Figure 3.1. False Discovery Rates and False Non-Discovery Rates for a Standard
Gaussian Mixing Distribution.

Thus far we have implicitly assumed that the size of the selected set is
predetermined by the parameter $\alpha$. Having established a ranking based on
a particular loss function, we simply select a subset of size $\lceil{\alpha
n}\rceil$ consisting of the highest ranked observations. In the next
subsection we begin to consider modifying this strategy by constraining the
probability of false discoveries. This will allow the size of the selected set
to adapt to the difficulty of the selection task.

###  3.4. Guarding against false discovery

Recognizing the risk of false “discoveries” among those selected, we will
consider an expanded loss function,

(3.1)
$L(\bm{\delta},\bm{\theta})=\sum_{i=1}^{n}h_{i}(1-\delta_{i})+\tau_{1}\Big{(}\sum_{i=1}^{n}\Big{\\{}(1-h_{i})\delta_{i}-\gamma\delta_{i}\Big{\\}}\Big{)}+\tau_{2}\Big{(}\sum_{i=1}^{n}\delta_{i}-\alpha
n\Big{)}$

where $h_{i}=\mymathbb 1\\{\theta_{i}\geq\theta_{\alpha}\\}$. If we set
$\tau_{1}$ to zero, then minimizing the expected loss leads to the Bayes rule
discussed in Section 3.1. On the other hand, if we set $\tau_{2}$ to zero,
then minimizing expected loss leads to a decision rule that is equivalent to a
multiple testing problem with null hypothesis
$H_{0i}:\theta_{i}\leq\theta_{\alpha}$; the goal is to minimize the expected
number of over-looked discoveries subject to the constraint that the marginal
FDR rate is controlled at level $\gamma$, that is,
$\mathbb{E}[\sum_{i=1}^{n}(1-h_{i})\delta_{i}]/\mathbb{E}[\sum_{i=1}^{n}\delta_{i}]\leq\gamma.$
When $\tau_{1}=0$, the risk can be expressed as,

$\mathbb{E}_{\bm{\theta}|Y}\Big{[}L(\bm{\delta},\bm{\theta})\Big{]}=\sum_{i=1}^{n}(1-\delta_{i})v_{\alpha}(Y_{i})+\tau_{2}\Big{(}\sum_{i=1}^{n}\delta_{i}-\alpha
n\Big{)}$

where
$v_{\alpha}(y_{i})=\mathbb{P}(\theta_{i}\geq\theta_{\alpha}|Y_{i}=y_{i})$.
Taking another expectation over $Y$, and minimizing over both $\bm{\delta}$
and $\tau_{2}$, leads to the decision rule,

$\delta_{i}^{*}=\begin{cases}1,&\text{if }v_{\alpha}(y_{i})\geq\tau_{2}^{*}\\\
0,&\text{if }v_{\alpha}(y_{i})<\tau_{2}^{*}.\end{cases}$

The Lagrange multiplier is chosen so that the constraint
$\mathbb{P}(\delta_{i}=1)\leq\alpha$ holds with equality:

$\tau_{2}^{*}=\text{min}\\{\tau_{2}:\mathbb{P}(v_{\alpha}(y_{i})\geq\tau_{2})\leq\alpha\\}$

Each selection improves the objective function by $v_{\alpha}(y_{i})$, but
incurs a cost of $\tau_{2}$. Since all selections incur the same cost, we may
rank according to $v_{\alpha}(y_{i})$, selecting units until the capacity
constraint $\alpha n$ is achieved. Selection of the last unit may need to be
randomized to exactly satisfy the constraint, as we note below.

When $\tau_{2}=0$ the focus shifts to the marginal FDR, the ratio of the
expected number of false discoveries, to the expected number of selections.
This is slightly different from the original FDR as defined in BH. However,
when $n$ is large the two concepts are asymptotically equivalent as shown by
GenoWass. Our objective becomes,

$\mathbb{E}_{\bm{\theta}|\bm{Y}}\Big{[}L(\bm{\delta},\bm{\theta})\Big{]}=\sum_{i=1}^{n}(1-\delta_{i})v_{\alpha}(Y_{i})+\tau_{1}\Big{(}\sum_{i=1}^{n}\\{\delta_{i}(1-v_{\alpha}(Y_{i}))-\gamma\delta_{i}\\}\Big{)}.$

Taking expectations again over $Y$ and minimizing over both $\bm{\delta}$ and
$\tau_{1}$ yields,

$\delta_{i}^{*}=\begin{cases}1,&\text{if
}v_{\alpha}(y_{i})>\tau_{1}^{*}(1-v_{\alpha}(y_{i})-\gamma)\\\ 0,&\text{if
}v_{\alpha}(y_{i})\leq\tau_{1}^{*}(1-v_{\alpha}(y_{i})-\gamma)\end{cases}$

and the Lagrange multiplier takes a value $\tau_{1}^{*}$ to make the marginal
FDR constraint hold with equality.

When both constraints are incorporated we must balance the power gain from
more selections and the cost that occurs from both the capacity constraint and
FDR control. The Bayes rule solves,

$\min_{\bm{\delta}}\mathbb{E}\Big{[}\sum_{i=1}^{n}(1-\delta_{i})v_{\alpha}(y_{i})\Big{]}+\tau_{1}\Big{(}\mathbb{E}\Big{[}\sum_{i=1}^{n}\Big{\\{}(1-v_{\alpha}(y_{i}))\delta_{i}-\gamma\delta_{i}\Big{\\}}\Big{]}\Big{)}+\tau_{2}\Big{(}\mathbb{E}\Big{[}\sum_{i=1}^{n}\delta_{i}\Big{]}-\alpha
n\Big{)}.$

Given the discrete nature of the decision function, this problem appears to
take the form of a classical napsack problem, however following the approach
of Basu we will consider a relaxed version of the problem in which units are
selected sequentially until one or the other constraint would be violated,
with the final selection randomized to satisfy the constraint exactly.

######  Remark 3.3.

Given the Lagrangian form of our loss function it is natural to consider an
optimization perspective for the selection problem. Minimizing the expectation
of the loss defined in (3.1) is equivalent to minimizing
$\mathbb{P}[\delta_{i}=0,\theta_{i}\geq\theta_{\alpha}]$ subject to the
constraint that
$\mathbb{P}[\delta_{i}=1,\theta_{i}<\theta_{\alpha}]/\mathbb{P}[\delta_{i}=1]\leq\gamma$
and $\mathbb{P}[\delta_{i}=1]\leq\alpha$. So we are looking for a thresholding
rule that minimizes the expected number of missed discoveries subject to the
capacity constraint and the constraint that the marginal FDR rate of the
decision rule is below level $\gamma$. This minimization problem is also
easily seen – from a testing perspective – to be equivalent to maximizing
power of the decision rule $\delta$,
$\mathbb{P}[\delta_{i}=1|\theta_{i}\geq\theta_{\alpha}]$, subject to the same
two constraints.

######  Proposition 3.4.

For any pair, $(\alpha,\gamma)$ such that $\gamma<1-\alpha$, the optimal Bayes
rule takes the form of

$\delta_{i}^{*}=\mymathbb
1\\{v_{\alpha}(y_{i})\geq\lambda^{*}(\alpha,\gamma)\\}$

where $\lambda^{*}(\alpha,\gamma)=v_{\alpha}(t^{*})$ with
$t^{*}=\max\\{t_{1}^{*},t_{2}^{*}\\}$ and

$\displaystyle t_{1}^{*}$
$\displaystyle=\min\Big{\\{}t:\frac{\int_{-\infty}^{\theta_{\alpha}}\tilde{\Phi}((t-\theta)/\sigma)dG(\theta)}{\int_{-\infty}^{+\infty}\tilde{\Phi}((t-\theta)/\sigma)dG(\theta)}-\gamma\leq
0\Big{\\}}$ $\displaystyle t_{2}^{*}$
$\displaystyle=\min\Big{\\{}t:\int_{-\infty}^{+\infty}\tilde{\Phi}((t-\theta)/\sigma)dG(\theta)-\alpha\leq
0\Big{\\}}$

with $\tilde{\Phi}$ being the survival function of a standard normal random
variable.

######  Remark 3.5.

The optimal cutoff $t^{*}$ depends on the data generating process and also the
choice of $\alpha$ and $\gamma$. When data is noisy, the FDR control
constraint may be binding before the capacity constraint is reached, and
consequently the selected set may be strictly smaller than the pre-specified
$\alpha$ proportion. On the other hand, when the signal is strong, the FDR
control constraint is unlikely to be binding before the capacity constraint is
attained.

We have seen that when variances are homogeneous, the optimal selection rule
thresholds on $Y$, so it is clear then that any ranking that is based on a
monotone transformation of $Y$ will lead to an equivalent selected set. We
should also highlight that we have focused on a null hypothesis that depends
on $\alpha$, while the multiple testing literature, for example ETST, SunCai
and Basu, typically focuses on the null hypothesis of $H_{0i}:\theta_{i}=0$.
When variances are homogeneous, it doesn’t matter whether we use an $\alpha$
dependent null or the conventional zero null, because the transformation based
on the conventional null, $\mathbb{P}(\theta>0|Y=y)$, is also a monotone
function of $Y$, and therefore yields an equivalent decision rule. However,
when variances are heterogeneous, this invariance no longer holds; different
transformations of the pair $(y,\sigma)$ lead to distinct decision rules that
lead to distinct performance, and using the conventional null hypothesis is no
longer advisable for the ranking and selection problem as we will show in the
next section.

##  4\. heterogeneous known variances

The homogeneous variance assumption of the preceding section is unsustainable
in most applications. Batting averages are accompanied by a number of “at
bats” and mean test score performances are accompanied by student sample
sizes. In this section we will consider the expanded model,

$Y_{i}\sim\mathcal{N}(\theta_{i},\sigma_{i}^{2}),\quad\text{and
}\quad\theta_{i}\sim G,\quad\sigma_{i}\sim
H,\quad\sigma_{i}\perp\\!\\!\\!\perp\theta_{i}$

We will assume that we observe $\sigma_{i}$, an assumption that will be
relaxed in the next section.

###  4.1. Posterior Tail Probability

With the same alternative hypothesis as $H_{A}:\theta\geq\theta_{\alpha}$, it
is natural to consider the posterior tail probability again, now as a function
of the pair, $(y_{i},\sigma_{i})$,

$v_{\alpha}(y_{i},\sigma_{i})=\mathbb{P}(\theta_{i}\geq\theta_{\alpha}|y_{i},\sigma_{i})=\frac{\int_{\theta_{\alpha}}^{+\infty}\frac{1}{\sigma_{i}}\phi((y_{i}-\theta)/\sigma_{i})dG(\theta)}{\int_{-\infty}^{+\infty}\frac{1}{\sigma_{i}}\phi((y_{i}-\theta)/\sigma_{i})dG(\theta)}.$

Solving the same decision problem with the loss function specified in (3.1),
we have the conditional risk,

$\mathbb{E}_{\bm{\theta}|\bm{Y},\bm{\sigma}}\Big{[}L(\bm{\delta},\bm{\theta})\Big{]}=\sum_{i=1}^{n}(1-\delta_{i})v_{\alpha}(Y_{i},\sigma_{i})+\tau_{1}\Big{(}\sum_{i=1}^{n}\\{\delta_{i}(1-v_{\alpha}(Y_{i},\sigma_{i}))-\gamma\delta_{i}\\}\Big{)}+\tau_{2}\Big{(}\sum_{i=1}^{n}\delta_{i}-\alpha
n\Big{)}.$

Taking another expectation with respect to the joint distribution of the
$(Y_{i},\sigma_{i})$’s, the Bayes rule solves

$\min_{\bm{\delta}}\mathbb{E}\Big{[}\sum_{i=1}^{n}(1-\delta_{i})v_{\alpha}(y_{i},\sigma_{i})\Big{]}+\tau_{1}\Big{(}\mathbb{E}\Big{[}\sum_{i=1}^{n}\Big{\\{}(1-v_{\alpha}(y_{i},\sigma_{i}))\delta_{i}-\gamma\delta_{i}\Big{\\}}\Big{]}\Big{)}+\tau_{2}\Big{(}\mathbb{E}\Big{[}\sum_{i=1}^{n}\delta_{i}\Big{]}-\alpha
n\Big{)}$

The optimal selection rule can again be characterized as a thresholding rule
on $v_{\alpha}(y_{i},\sigma_{i})$.

######  Proposition 4.1.

For a pre-specified pair $(\alpha,\gamma)$ such that $\gamma<1-\alpha$, the
optimal Bayes rule takes the form

$\delta^{*}(y,\sigma)=\mymathbb
1\\{v_{\alpha}(y,\sigma)\geq\lambda^{*}(\alpha,\gamma)\\}$

where
$\lambda^{*}(\alpha,\gamma)=\max\\{\lambda_{1}^{*}(\alpha,\gamma),\lambda_{2}^{*}(\alpha)\\}$
and

$\displaystyle\lambda_{1}^{*}(\alpha,\gamma)$
$\displaystyle=\min\Big{\\{}\lambda:\frac{\int\int_{-\infty}^{\theta_{\alpha}}\tilde{\Phi}((t_{\alpha}(\lambda,\sigma)-\theta)/\sigma)dG(\theta)dH(\sigma)}{\int\int_{-\infty}^{+\infty}\tilde{\Phi}((t_{\alpha}(\lambda,\sigma)-\theta)/\sigma)dG(\theta)dH(\sigma)}-\gamma\leq
0\Big{\\}}$ $\displaystyle\lambda_{2}^{*}(\alpha)$
$\displaystyle=\min\Big{\\{}\lambda:\int\int_{-\infty}^{+\infty}\tilde{\Phi}((t_{\alpha}(\lambda,\sigma)-\theta)/\sigma)dG(\theta)dH(\sigma)-\alpha\leq
0\Big{\\}}$

denoting $\tilde{\Phi}=1-\Phi$ and with $t_{\alpha}(\lambda,\sigma)$ defined
as $v_{\alpha}(t_{\alpha}(\lambda,\sigma),\sigma)=\lambda$ for all
$\lambda\in[0,1]$.

######  Remark 4.2.

Note that although the thresholding value $\lambda^{*}$ does not depend on the
value of $\sigma$, the ranking does depend on $\sigma$. One way to see this is
that since $v_{\alpha}(y,\sigma)$ is monotone in $y$ for all $\sigma>0$, the
optimal rule is equivalent to $\mymathbb
1\\{y_{i}>t_{\alpha}(\lambda^{*},\sigma)\\}$, where
$t_{\alpha}(\lambda,\sigma)$ is a function of $\sigma$. For a fixed value of
$\lambda^{*}$, the selection region for $Y$ depends on $\sigma$ in a nonlinear
way. Comparing individuals $i$ and $j$, it may be the case that $y_{i}>y_{j}$,
but $y_{j}$ belongs to the selection region while $y_{i}$ does not. An example
to illustrate this appears below. It should also be emphasized that when
variances are heterogeneous, different loss functions, need not lead to
equivalent selections.

###  4.2. The Conventional Null Hypothesis

The posterior tail probability criterion is motivated by viewing the ranking
and selection problems as hypothesis testing while allowing the null
hypothesis to be $\alpha$ dependent. The particular construction of the null
hypothesis turns out to be critical for the ranking exercise. In this
subsection we present a simple example to illustrate that tail probability
based on the conventional null hypothesis of zero effect does not lead to a
powerful ranking device. Consider data generated from a three component normal
mixture model,

(4.1) $Y_{i}|\sigma_{i}\sim
0.85\mathcal{N}(-1,\sigma_{i}^{2})+0.1\mathcal{N}(0.5,\sigma_{i}^{2})+0.05\mathcal{N}(5,\sigma_{i}^{2}),\quad\sigma_{i}\sim
U[0.5,4]$

Instead of transforming the data by $v_{\alpha}$, we consider the
transformation,

$\texttt{T}(y_{i},\sigma_{i}^{2})=\mathbb{P}(\theta_{i}>0|y_{i},\sigma_{i}^{2})=\frac{\int_{0}^{+\infty}\frac{1}{\sigma_{i}}\phi((y_{i}-\theta)/\sigma_{i})dG(\theta)}{\int_{-\infty}^{+\infty}\frac{1}{\sigma_{i}}\phi((y_{i}-\theta)/\sigma_{i})dG(\theta)}$

and rank individuals accordingly. This transformation corresponds to the
procedure proposed in SM, and is motivated for multiple testing problems under
the conventional null hypothesis $H_{0}:\theta\leq 0$. The decision rule
$\delta_{i}^{\texttt{T}}=\mymathbb
1\\{\texttt{T}(y_{i},\sigma_{i})\geq\lambda\\}$ then chooses the cutoff value
$\lambda$ that respects both the capacity constraint and the FDR control
constraint for selecting the top $\alpha$ proportion.

Figure 4.1. Selection boundaries based on the model (4.1) with $\alpha=0.05$
and $\gamma=0.1$. The solid black curve corresponds to the boundary of the
selection region based on transformation $v_{\alpha}$. The dash red curve
corresponds to the boundary of the selection region based on transformation
$T$. Density of $\sigma$ is assumed to be uniform on the interval $[0.5,4]$.

Figure 4.1 compares the selection region for the two ranking procedures with
$\alpha=5\%$ and marginal FDR control at level $10\%$. The solid black line
corresponds to the selection boundary using ranking based on transformation
$v_{\alpha}$ and the dashed red line corresponds to the selection boundary
using ranking based on the transformation T. The black highlighted area below
the black selection boundary corresponds to a region where the ranking method
based on T will select but the ranking method based on $v_{\alpha}$ does not.
On the other hand, the blue highlighted area corresponds to a region selected
by $v_{\alpha}$, but not for T. The transformation T ranks those in the black
region higher than those in the blue region because although they have a
relatively smaller mean effect $y$, their associated variances are also
smaller indicating stronger evidence that such individuals have a positive
$\theta$ than those located in the blue area. However, our task is to find
individuals with true effects, $\theta_{i}$, in the upper tail. For
$\alpha=5\%$, we aim to select all individuals with $\theta=5$, individuals in
the black region present strong evidence that their true effect can not be too
large because their observed effect $y$ is small and their associated variance
is also small, while those in the blue region, although their observed mean
effects are associated with larger variances offer reasonable evidence that
their associated true effect $\theta$ may be large. This evidence is not
apparent in transformation T, but is captured in the transformation
$v_{\alpha}$.

Indeed, the average power of ranking based on the two different transformation
$v_{\alpha}$ and T differ significantly. Defining the power of the selection
rule as
$\beta(\delta):=\mathbb{P}(\theta_{i}\geq\theta_{\alpha},\delta_{i}=1)/\mathbb{P}(\theta_{i}\geq\theta_{\alpha})$,
the proportion of true top $\alpha$ cases selected based on decision rule
$\bm{\delta}$, then $\beta(\bm{\delta^{\texttt{T}}})=39\%$ and
$\beta(\bm{\delta^{*}})=69\%$. Thus, although much of the literature relies on
ranking and selection rules based on some form of posterior means and
conventional hypothesis testing apparatus we would caution that such methods
can be quite misleading and inefficient.

###  4.3. Nestedness of Selection Sets

If we were to relax the capacity constraint to allow a larger proportion,
$\alpha_{1}>\alpha_{0}$ be be selected, while maintaining our initial false
discovery control, we would expect that members selected under the more
stringent capacity constraint should remain selected under the relaxed
constraint. We now discuss sufficient conditions under which we obtain this
nestedness of the selection sets when using the posterior tail probability
rule.

The optimal Bayes rule defines the selection set for each pair of
$(\alpha,\gamma)$ as

$\Omega_{\alpha,\gamma}:=\\{(y,\sigma):v_{\alpha}(y,\sigma)\geq\lambda^{*}(\alpha,\gamma)\\}$

and when $\sigma$ is known, $v_{\alpha}(y,\sigma)$ is monotone in $y$ as shown
in Lemma 3.2 for each fixed $\sigma$, hence the selection set can also be
represented as

$\Omega_{\alpha,\gamma}=\\{(y,\sigma):y\geq
t_{\alpha}(\lambda^{*}(\alpha,\gamma),\sigma)\\}$

It is also convenient for later discussion to define

$\displaystyle\Omega_{\alpha,\gamma}^{FDR}:$
$\displaystyle=\\{(y,\sigma):v_{\alpha}(y,\sigma)\geq\lambda_{1}^{*}(\alpha,\gamma)\\}=\\{(y,\sigma):y\geq
t_{\alpha}(\lambda_{1}^{*}(\alpha,\gamma),\sigma)\\}$
$\displaystyle\Omega_{\alpha}^{C}$
$\displaystyle:=\\{(y,\sigma):v_{\alpha}(y,\sigma)\geq\lambda_{2}^{*}(\alpha)\\}=\\{(y,\sigma):y\geq
t_{\alpha}(\lambda_{2}^{*}(\alpha),\sigma)\\}$

which are respectively the selection sets when the false discovery rate
constraint or the capacity constraint is binding. It is easy to see that
$\Omega_{\alpha,\gamma}=\Omega_{\alpha,\gamma}^{FDR}\cap\Omega_{\alpha}^{C}$.

######  Lemma 4.3.

Let the density function of $v_{\alpha}(y_{i},\sigma_{i})$ be denoted as
$f_{v}(v;\alpha)$ and let

$\lambda_{1}^{*}(\alpha,\gamma)=\min\Big{\\{}\lambda:\frac{\int\int_{-\infty}^{\theta_{\alpha}}\tilde{\Phi}((t_{\alpha}(\lambda,\sigma)-\theta)/\sigma)dG(\theta)dH(\sigma)}{\int\int_{-\infty}^{+\infty}\tilde{\Phi}((t_{\alpha}(\lambda,\sigma)-\theta)/\sigma)dG(\theta)dH(\sigma)}-\gamma\leq
0\Big{\\}}$

with $t_{\alpha}(\lambda,\sigma)$ defined as
$v_{\alpha}(t_{\alpha}(\lambda,\sigma),\sigma)=\lambda$ and $\tilde{\Phi}$ be
the survival function of the standard normal random variable. If
$\nabla_{\alpha}\log f_{v}(v;\alpha)$ is non-decreasing in $v$, then for fixed
$\gamma$, if $\alpha_{1}>\alpha_{2}$, we have
$\lambda_{1}^{*}(\alpha_{1},\gamma)\leq\lambda_{1}^{*}(\alpha_{2},\gamma)$.

######  Remark 4.4.

The density function $f_{v}(v;\alpha)$ can be viewed as a function of $v$
indexed by the parameter $\alpha$. An explicit form for $f_{v}(v;\alpha)$
appears in Section 4.4 for the normal-normal model). The condition imposed in
Lemma 4.3 is equivalent to a monotone likelihood ratio condition, that is that
the likelihood ratio $f_{v}(v;\alpha_{1})/f_{v}(v;\alpha_{2})$ is non-
decreasing in $v$ if $\alpha_{1}>\alpha_{2}$. The monotone likelihood ratio
condition is stronger in that it implies stochastic dominance between the
distributions of the random variables $v_{\alpha_{1}}$ and $v_{\alpha_{2}}$,
which holds by construction.

######  Corollary 4.5.

If the condition in Lemma 4.3 holds, then
$\Omega_{\alpha_{2},\gamma}^{FDR}\subseteq\Omega_{\alpha_{1},\gamma}^{FDR}$
for any $\alpha_{1}>\alpha_{2}$.

######  Remark 4.6.

The condition in Lemma 4.3 is sufficient but not necessary for nestedness of
$\Omega_{\alpha,\gamma}^{FDR}$ because even when
$\lambda_{1}^{*}(\alpha_{1},\gamma)>\lambda_{1}^{*}(\alpha_{2},\gamma)$, we
can still have
$t_{\alpha_{1}}(\lambda_{1}^{*}(\alpha_{1},\gamma),\sigma)<t_{\alpha_{2}}(\lambda^{*}_{1}(\alpha_{2},\gamma),\sigma)$
because the function $v_{\alpha}(Y,\sigma)$ depends on $\alpha$, so does its
inverse function $t_{\alpha}$.

######  Lemma 4.7.

Let

$\lambda_{2}^{*}(\alpha)=\min\\{\lambda:\int\int_{-\infty}^{+\infty}\tilde{\Phi}((t_{\alpha}(\lambda,\sigma)-\theta)/\sigma)dG(\theta)dH(\sigma)\leq\alpha\\}$

If for any $\alpha_{1}>\alpha_{2}$,
$t_{\alpha_{1}}(\lambda_{2}^{*}(\alpha_{1}),\sigma)\leq
t_{\alpha_{2}}(\lambda_{2}^{*}(\alpha_{2}),\sigma)$ for each $\sigma$, then
$\Omega_{\alpha_{2}}^{C}\subseteq\Omega_{\alpha_{1}}^{C}$.

######  Remark 4.8.

Monotonicity coincides with the condition in Theorem 3 in HN. They demonstrate
that it holds when $G$ follows a normal distribution. However, it need not
hold as shown in the counter-example in Section 4.5.

######  Lemma 4.9.

If $\nabla_{\alpha}\log f_{v}(v;\alpha)$ is non-decreasing in $v$ and the
condition in Lemma 4.7 holds, then for a fixed $\gamma$, the selection region
has a nested structure: if $\alpha_{1}>\alpha_{2}$ then
$\Omega_{\alpha_{2},\gamma}\subseteq\Omega_{\alpha_{1},\gamma}$.

###  4.4. Examples

In this section we consider several examples beginning with the simplest
classical case in which the $\theta_{i}$ constitute a random sample from the
standard Gaussian distribution. This Gaussian assumption on the form of the
mixing distribution $G$ underlies almost all of the empirical Bayes literature
in applied economics; it is precisely what justifies the linear shrinkage
rules that are typically employed.

######  Example 4.10 (Gaussian $G$).

Consider the normal-normal model, where
$y|\theta,\sigma^{2}\sim\mathcal{N}(\theta,\sigma^{2})$ and
$\theta\sim\mathcal{N}(0,\sigma_{\theta}^{2})$ and $\sigma\sim H$ with density
function $h(\sigma)$. The marginal distribution of $y$ given $\sigma^{2}$ is
$\mathcal{N}(0,\sigma^{2}+\sigma_{\theta}^{2})$ and the joint density of
$(y,\sigma)$ takes the form

$f(y,\sigma)=\frac{1}{\sqrt{2\pi(\sigma^{2}+\sigma_{\theta}^{2})}}\exp\Big{\\{}-\frac{y^{2}}{2(\sigma^{2}+\sigma_{\theta}^{2})}\Big{\\}}h(\sigma)$

Given the normal conjugacy, the posterior distribution of
$\theta|y,\sigma^{2}$ follows $\mathcal{N}(\rho y,\rho\sigma^{2})$ where
$\rho=\sigma_{\theta}^{2}/(\sigma_{\theta}^{2}+\sigma^{2})$. The random
variable $v$ is a transformation of the pair $(Y,\sigma^{2})$, defined as

$v=\psi(y,\sigma^{2}):=\mathbb{P}(\theta\geq\theta_{\alpha}|y,\sigma^{2})=\Phi((\rho
y-\theta_{\alpha})/\sqrt{\rho\sigma^{2}})$

For fixed $\sigma^{2}$, $\psi$ is monotone increasing in $y$ and
$\psi^{-1}(v)=\theta_{\alpha}/\rho+\sqrt{\sigma^{2}/\rho}\Phi^{-1}(v)$ with
$\nabla_{v}\psi^{-1}(v)=\sqrt{\sigma^{2}/\rho}/\phi(\Phi^{-1}(v))$. The joint
density of $v$ and $\sigma$ is thus,

$\displaystyle g(v,\sigma)$
$\displaystyle=f(\psi^{-1}(v),\sigma)|\nabla_{v}\psi^{-1}(v)|$
$\displaystyle=\frac{1}{\sqrt{2\pi(\sigma^{2}+\sigma_{\theta}^{2})}}\exp\Big{\\{}-\frac{(\theta_{\alpha}/\rho+\sqrt{\sigma^{2}/\rho}\Phi^{-1}(v))^{2}}{2(\sigma^{2}+\sigma_{\theta}^{2})}\Big{\\}}\frac{\sqrt{\sigma^{2}/\rho}}{\phi(\Phi^{-1}(v))}h(\sigma).$

Integrating out $\sigma$ we have the marginal density of $v$,

$f_{v}(v;\alpha)=\int\frac{1}{\sqrt{2\pi(\sigma^{2}+\sigma_{\theta}^{2})}}\exp\Big{\\{}-\frac{(\theta_{\alpha}/\rho+\sqrt{\sigma^{2}/\rho}\Phi^{-1}(v))^{2}}{2(\sigma^{2}+\sigma_{\theta}^{2})}\Big{\\}}\frac{\sqrt{\sigma^{2}/\rho}}{\phi(\Phi^{-1}(\tilde{v}))}dH(\sigma)$

The capacity constraint is $\mathbb{P}(v\geq\lambda^{*}_{2})=\alpha$, with
cut-off value, $\lambda^{*}_{2}$, satisfying,

$\alpha=\mathbb{P}(v\geq\lambda_{2}^{*})=1-\int\Phi\Big{(}\theta_{\alpha}\frac{\sqrt{\sigma^{2}+\sigma_{\theta}^{2}}}{\sigma_{\theta}^{2}}-\Phi^{-1}(1-\lambda_{2}^{*})\sqrt{\sigma^{2}/\sigma_{\theta}^{2}}\Big{)}dH(\sigma).$

To find $\lambda_{1}^{*}$, we can use the formula provided in Proposition 4.1.
A more direct approach is to recognize, see Section 6, that the FDR control
constraint is defined as $\gamma=\mathbb{E}[(1-v)\mymathbb
1\\{v\geq\lambda_{1}^{*}\\}]/\mathbb{P}(v\geq\lambda_{1}^{*})$, where the cut-
off value $\lambda_{1}^{*}$ is defined through

$\gamma=\int_{\lambda_{1}^{*}}^{1}(1-v)f_{v}(v;\alpha)dv/\int_{\lambda_{1}^{*}}^{1}f_{v}(v;\alpha)dv.$

Let $\lambda^{*}=\max\\{\lambda_{1}^{*},\lambda_{2}^{*}\\}$, the selection
region is then $\\{(y,\sigma):y\geq t_{\alpha}(\lambda^{*},\sigma)\\}$ with

$t_{\alpha}(\lambda^{*},\sigma)=\theta_{\alpha}/\rho-\Phi^{-1}(1-\lambda^{*})\sqrt{\sigma^{2}/\rho}.$

Suppose we use the posterior mean of $\theta$ as a ranking device, so
$\delta_{i}^{\texttt{PM}}=\mymathbb 1\\{y\rho\geq\omega^{*}\\}$ for some
suitably chosen $\omega^{*}$ that guarantees both capacity and FDR control.
For the capacity constraint, the thresholding value solves,

$\displaystyle 1-\alpha$
$\displaystyle=\int\mathbb{P}(y\rho<\omega_{2}^{*})dH(\sigma)$
$\displaystyle=\int\Phi\Big{(}\omega_{2}^{*}/(\sigma_{\theta}^{2}\sqrt{\sigma_{\theta}^{2}+\sigma^{2}})\Big{)}dH(\sigma),$

while FDR control requires a thresholding value that solves,

$\displaystyle\gamma$
$\displaystyle=\int\mathbb{P}(y\geq\omega_{1}^{*}/\rho,\theta<\theta_{\alpha})dH(\sigma)/\int\mathbb{P}(y\geq\omega_{1}^{*}/\rho)dH(\sigma)$
$\displaystyle=\int\int_{[\omega_{1}^{*}/\rho,+\infty)}(1-\alpha)f_{0}(y)dydH(\sigma)/\int
1-\Phi\Big{(}\omega_{1}^{*}/(\sigma_{\theta}^{2}\sqrt{\sigma_{\theta}^{2}+\sigma^{2}}\Big{)}dH(\sigma),$

with

$f_{0}(y)=\frac{1}{1-\alpha}\frac{1}{\sqrt{2\pi(\sigma_{\theta}^{2}+\sigma^{2})}}\exp\Big{\\{}-\frac{y^{2}}{2(\sigma_{\theta}^{2}+\sigma^{2})}\Big{\\}}\Phi\Big{(}\frac{(\theta_{\alpha}-y\rho)}{\sqrt{\rho\sigma^{2}}}\Big{)},$

denoting the density of $y$ under the null $\theta<\theta_{\alpha}$. Setting
$\omega^{*}=\max\\{\omega_{1}^{*},\omega_{2}^{*}\\}$, the selection region is
then $\\{(y,\sigma):y\geq\omega^{*}/\rho\\}$.

Figure 4.2. The left panel plots the selection boundaries for the normal-
normal model with $\sigma_{\theta}^{2}=1$ and $\alpha=0.05$ and $\gamma=0.2$.
The density of $\sigma$ is assumed to be uniform on the range $[0.5,1]$.
Selected units must have $(y_{i},\sigma_{i})$ above the curves. The red curves
correspond to the selection region boundaries with FDR controlled at level
$0.2$; solid lines for posterior mean ranking and dash line for posterior tail
probability ranking. The black curves correspond to the selection boundaries
with capacity control at level $0.05$. The middle and right panels illustrate
the selected set from a realized sample of size 10,000. The grey circles
correspond to individuals selected by both the posterior tail probability rule
and the posterior mean rule. The green crosses depict individuals selected by
the posterior mean rule but not the tail probability and the red crosses
indicate individuals selected by the tail probability rule, but not by the
posterior mean rule.

Figure 4.2 plots the selection boundaries for both constraints with
$\theta\sim\mathcal{N}(0,1)$ and $\sigma\sim U[0.5,1]$. With $\alpha=0.05$ and
$\gamma=0.2$, the FDR constraint is binding, but not the capacity constraint.
In this example, if we only impose the capacity constraint to be 5 percent,
even an Oracle totally aware of $G$, will face a false discovery rate of
nearly 52 percent. In other words more than half of those selected to be in
the right tail will be individuals with $\theta<\theta_{\alpha}$ rather than
from the intended $\theta\geq\theta_{\alpha}$ group. This fact motivates our
more explicit incorporation of FDR into the selection constraints. We may
recall that in the homogeneous variance Gaussian setting we saw in Figure 3.1
that FDR was also very high when $\alpha$ is set at 0.05. Figure 4.2 also
depicts the selected set with a realized sample of 10,000 from the normal-
normal model. With capacity constraint alone, the posterior mean criteria
favours individuals with smaller variances. When the FDR constraint is
implemented, with $\gamma=0.2$, it becomes the binding constraint in this
setting, both criteria become more stringent and only a much smaller set of
individuals are selected, and there is less conflict in the selections. The
corresponding selected sets are plotted in the right panel of Figure 4.2. When
the variance parameter $\sigma_{\theta}^{2}$ in $G$ is not observed, we can
estimate it via the MLE based on the marginal likelihood of $Y$. This leads to
the generalized James-Stein estimator proposed in EfronMorris.

######  Example 4.11 (A Discrete $G$).

Suppose $\theta\sim 0.85\delta_{-1}+0.1\delta_{2}+0.05\delta_{5}$. Then the
marginal density of $y$ given $\sigma^{2}$ takes the form,

$\displaystyle f(y|\sigma^{2})$
$\displaystyle=\int\frac{1}{\sigma}\phi((y-\theta)/\sigma)dG(\theta)$
$\displaystyle=\frac{1}{\sigma}\Big{(}0.85\phi((y+1)/\sigma)+0.1\phi((y-2)/\sigma)+0.05\phi((y-5)/\sigma))$

And the random variable $v$ is a transformation of the pair $(y,\sigma^{2})$,
defined as,

$\displaystyle v$
$\displaystyle=\psi(y,\sigma^{2}):=\mathbb{P}(\theta\geq\theta_{\alpha}|y,\sigma^{2})=\frac{\int_{\theta_{\alpha}}^{+\infty}\frac{1}{\sigma}\phi((y-\theta)/\sigma)dG(\theta)}{\int_{-\infty}^{\infty}\frac{1}{\sigma}\phi((y-\theta)/\sigma)dG(\theta)}.$

The capacity constraint leads to a thresholding rule on $v$ such that
$\mathbb{P}(v\geq\lambda_{2}^{*})=\alpha$, while the FDR controls leads to a
cutoff value $\lambda_{1}^{*}$ defined through
$\gamma=\mathbb{E}[(1-v)\mymathbb
1\\{v\geq\lambda_{1}^{*}\\}]/\mathbb{P}(v\geq\lambda_{1}^{*})$. Let
$\lambda^{*}=\max\\{\lambda_{1}^{*},\lambda_{2}^{*}\\}$, the selection region
is then defined by $\\{(y,\sigma):y\geq t_{\alpha}(\lambda^{*},\sigma)\\},$
and can be found easily numerically.

Figure 4.3. The left panel plots the selection boundaries for the normal-
discrete model with $\theta\sim
G=0.85\delta_{-1}+0.1\delta_{2}+0.05\delta_{5}$ and $\alpha=0.05$ and
$\gamma=0.1$. The density of $\sigma$ is uniform on the range $[0.5,4]$. The
red curves correspond to the selection region with FDR controlled at level
$0.2$, solid lines for posterior mean ranking and dashed lines for posterior
tail probability ranking. The black curves correspond to the selection region
with capacity control at level $0.05$. The other panels are structured as in
the previous figure.

Figure 4.3 plots the selection boundaries for both constraints when $\theta$
follows this discrete distribution. We again set $\alpha=0.05$ and
$\gamma=0.2$, so we would like to select all the individuals associated with
the largest effect size, $\\{\theta=5\\}$, while controlling the FDR rate
below $20\%$. The red curves again correspond to FDR control with the two
ranking procedure, while the black curves correspond to capacity control. For
the two regions to overlap with $\alpha$ fixed at $0.05$, we must be willing
to tolerate $\gamma\approx 37\%$. In this case, we see that the posterior
probability ranking procedure prefers individuals with larger variances, while
the posterior mean ranking procedure prefers smaller variances. Based on a
realized sample of 10,000, Figure 4.3 again shows the selected observations
and once more we see that the posterior mean criteria favours individuals with
smaller variances, both under the capacity constraint and the FDR constraint.
In contrast to the normal-normal setting, now the FDR constraint is much less
severe and allows us to select considerably more individuals.

###  4.5. A Counterexample: Non-nested Selection Regions

Thus far we have stressed conditions under which selection regions are nested
with respect to $\alpha$, that is, for
$\alpha_{1}<\alpha_{2}<\cdots<\alpha_{m}$, we have
$\Omega_{\alpha_{1}}\subseteq\Omega_{\alpha_{2}}\subseteq\cdots\subseteq\Omega_{\alpha_{m}}$,
for the selection regions. However, this need not hold when there is variance
heterogeneity and when nesting fails we can have seemingly anomalous
situations in which units are selected by the tail probability rule at some
stringent, low $\alpha$, but are then rejected for some less stringent, larger
$\alpha$’s. To illustrate this phenomenon we will neglect the FDR constraint
and focus on our discrete mixing distribution,
$G=0.85\delta_{-1}+0.10\delta_{2}+0.05\delta_{5}$, with $\sigma\sim U[1/2,4]$.
The selection regions are depicted in Figure 4.4 for
$\alpha\in\\{0.04,0.05,0.06,0.08\\}$. Units are selected when their observed
pair, $(y_{i},\sigma_{i})$ lies above these curves for various $\alpha$’s.
When $\sigma$ is small we see, as expected, that selection is nested: if a
unit is selected at low $\alpha$ it stays selected at larger $\alpha$’s.
However, when $\sigma=3$, we see that there are units selected at
$\alpha=0.05$ and even $\alpha=0.04$ and yet they are rejected for
$\alpha=0.06$. How can this be?

Figure 4.4. Oracle selection boundaries (with just capacity constraint) for
several $\alpha$ levels for the tail probability criterion (left panel) and
posterior mean criterion (right panel) with a discrete example with
$G=0.85\delta_{-1}+0.10\delta_{2}+0.05\delta_{5}$, and $\sigma\sim U[1/2,4]$.
Crossing of the boundaries implies that the selection regions are non-nested
as explained in the text.

Imagine you are the oracle, so you know $G$, and when you decide to select
with $\alpha=0.06$ you know that you will have to select a few $\theta=2$
types, since there are only 5 percent of the $\theta=5$ types. Your main worry
at that point is to try to avoid selecting any $\theta=-1$ types; this can be
accomplished but only by avoiding the high $\sigma$ types. In contrast when
$\alpha=0.05$ so we are trying to vacuum up all of the $\theta=5$ types it is
worth taking more of a risk with high $\sigma$ types as long as their $y_{i}$
is reasonably large.

The crossing of the selection boundaries and non-nestedness of the selection
regions is closely tied up with the tail probability criterion and the
$\alpha$ dependent feature of the hypothesis. If we repeat our exercise with
the same $G$, and $\sigma$ distribution, but select according to posterior
means, we get the nested selection boundaries illustrated in Figure 4.4.
Proposition 4.12 establishes this to be a general phenomenon for any
distributution $G$.

It should be noted that the crossing of selection boundaries we have
illustrated seems to have been anticipated in HN, who consider similar tail
criteria. They propose a ranking scheme that assigns rank equal to the
smallest $\alpha$ for which a unit would be selected as a way to resolve the
ambiguities generated by crossing. We don’t see a compelling decision
theoretic rationale for this revised ranking rule, instead we prefer to
maintain some separation between the ranking and selection problems and focus
on risk assessment as a way to reconcile them.

Figure 4.5. Oracle risk evaluation for several $\alpha$ levels for the tail
probability criterion and a discrete example with
$G=0.85\delta_{-1}+0.10\delta_{2}+0.05\delta_{5}$, and $\sigma\sim U[1/2,4]$.
The solid lines correspond to the evaluation of the loss function specified in
(3.1) with both capacity and FDR control constraints. The dotted line
corresponds to oracle risk when only the capacity constraint is imposed.

The risk based on the loss function defined in (3.1) clearly depends on
$\alpha$ and $\gamma$. More specifically it consists of three pieces, the
leading term has the interpretation of “missed discovery” probability, which
we try to minimize and the second and third pieces correspond to the FDR and
capacity constraints respectively each weighted by a Lagrangian multiplier.
Focusing on the first term, we have

$\displaystyle\mathbb{E}[H_{i}(1-\delta_{i})]$
$\displaystyle=\mathbb{E}[\mymathbb
1\\{\theta_{i}\geq\theta_{\alpha}\\}(1-\delta_{i})]$
$\displaystyle=\mathbb{P}[\theta_{i}\geq\theta_{\alpha}]-\mathbb{E}[\delta_{i}\mymathbb
1\\{\theta_{i}\geq\theta_{\alpha}\\}]$
$\displaystyle=\mathbb{P}[\theta_{i}\geq\theta_{\alpha}]-\int\int_{\theta_{\alpha}}^{+\infty}\Big{[}1-\Phi((t_{\alpha}(\lambda,\sigma)-\theta)/\sigma)\Big{]}dG(\theta)dH(\sigma)$

where $\lambda$, depends on $\alpha$ and $\gamma$, is determined by either the
false discovery rate control or the capacity constraint, whichever binds.

A feature of discrete mixing distributions, $G$, is that the first term in the
loss, $\mathbb{P}(\theta_{i}\geq\theta_{\alpha})$, is piece-wise constant,
with jumps occurring only at discontinuity points of $G$, while the second
term depends on both $\theta_{\alpha}$ and the cut-off values
$t_{\alpha}(\lambda,\sigma)$. When the capacity constraint binds there exist
ranges of $\alpha$ such that $\theta_{\alpha}$ remains constant, while
$t_{\alpha}(\lambda,\sigma)$ decreases for each $\sigma$, hence the risk with
just capacity constraint binding is a decreasing function for $\alpha$ in the
interval $(0.05,0.15)$. On the other hand, when the FDR constraint binds, it
can be shown that the cut-off value $t_{\alpha}(\lambda,\sigma)$ is constant.
To see this recall that the cutoff $\lambda$ determined by the FDR constraint
is defined as $\mathbb{E}[(1-v_{\alpha}(Y,\sigma))\mymathbb
1\\{v_{\alpha}(Y,\sigma)\geq\lambda\\}]=\gamma\mathbb{P}(v_{\alpha}(Y,\sigma)\geq\lambda)$,
so when $\theta_{\alpha}$ is constant over a range of $\alpha$ the
distribution of $v_{\alpha}(Y,\sigma)$ does not change, and consequently the
value $\lambda$ is constant over that range of $\alpha$.

Figure 4.5 evaluates risk based on the optimal selection rule for various
$\alpha$ and FDR levels, $\gamma\in\\{0.01,0.05,0.1,0.15,0.3\\}$. The solid
curves correspond to risk evaluated at the optimal Bayes rule defined in
Proposition 4.1. The dotted line corresponds to the risk evaluated at the
Bayes rule when only the capacity constraint is imposed. As $\gamma$
increases, the risk decreases as expected. For FDR levels as stringent as
$\gamma=0.01$, the FDR constraint binds and the risk is piece-wise constant.
As the FDR level is relaxed, there are range of $\alpha$ such that capacity
constraint becomes binding, and the risk decreases after the initial jump at
$\alpha=0.05$.

The evaluation of the risk for this particular example indicates that it is
easier to select the top 5% individuals, those with $\theta=5$. As we intend
to select more in the right tail, we are facing more uncertainty. This also
motivates a more systematic choice of $(\alpha,\gamma)$. Although selection
based on the tail probability criterion can lead to non-nested selection
regions, we conclude this sub-section by demonstrating that posterior mean
selection is necessarily nested.

######  Proposition 4.12.

Let the density function of $y$ conditional on $\theta$ and $\sigma$ be
denoted as $f(y|\theta,\sigma)$. If selection is based on the posterior mean,
$\delta_{i}=\\{M(y,\sigma)\geq c(\alpha,\gamma)\\}$ with

$M(y,\sigma):=\frac{\int\theta f(y|\theta,\sigma)dG(\theta)}{\int
f(y|\theta,\sigma)dG(\theta)}$

and $c(\alpha,\gamma)$ is chosen to satisfy both the capacity constraint at
level $\alpha$ and the FDR constraint at level $\gamma$, then the selection
regions, defined as $\Lambda_{\alpha,\gamma}=\\{(y,\sigma):M(y,\sigma)\geq
c(\alpha,\gamma)\\}$, are nested, that is, for any $\alpha_{1}>\alpha_{2}$,
$\Lambda_{\alpha_{2},\gamma}\subseteq\Lambda_{\alpha_{1},\gamma}$.

##  5\. Heterogeneous unknown variances

Assuming that the $\sigma^{2}$’s are known, up to a common scale parameter,
may be plausible in some applications such as baseball batting averages, but
it is frequently more plausible to adopt the view that we are simply
confronted with an estimate available perhaps from longitudinal data. In such
cases we need to consider the pairs, $(y_{i},S_{i})$ as potentially jointly
dependent random variables arising from the longitudinal model,

$Y_{it}=\theta_{i}+\sigma_{i}\epsilon_{it},\quad\epsilon_{it}\sim_{iid}\mathcal{N}(0,1),\quad(\theta_{i},\sigma_{i}^{2})\sim
G,$

with sufficient statistics, $Y_{i}=T_{i}^{-1}\sum_{t=1}^{T_{i}}Y_{it}$ and
$S_{i}=(T_{i}-1)^{-1}\sum_{t=1}^{T_{i}}(Y_{it}-Y_{i})^{2}$, for
$(\theta_{i},\sigma_{i}^{2})$. Conditional on $(\theta_{i},\sigma_{i}^{2})$,
we have
$Y_{i}|\theta_{i},\sigma_{i}^{2}\sim\mathcal{N}(\theta_{i},\sigma_{i}^{2}/T_{i})$
and $S_{i}|\sigma_{i}^{2}$ is distributed as Gamma with shape parameter
$r_{i}=(T_{i}-1)/2$, scale parameter, $\sigma_{i}^{2}/r_{i}$, and density
function denoted as $\Gamma(S_{i}|r_{i},\sigma_{i}^{2}/r_{i})$.

Given the loss function (3.1) and defining $\theta_{\alpha}$ as
$\alpha=\mathbb{P}(\theta_{i}\geq\theta_{\alpha})=\int\int_{\theta_{\alpha}}^{+\infty}dG(\theta,\sigma^{2})$,
the conditional risk is,

$\displaystyle\mathbb{E}_{\bm{\theta}|\bm{Y},\bm{S}}\Big{[}L(\bm{\delta},\bm{\theta})\Big{]}$
$\displaystyle=\sum_{i=1}^{n}(1-\delta_{i})v_{\alpha}(Y_{i},S_{i})$
$\displaystyle+\tau_{1}\Big{(}\sum_{i=1}^{n}\\{\delta_{i}(1-v_{\alpha}(Y_{i},S_{i}))-\gamma\delta_{i}\\}\Big{)}+\tau_{2}\Big{(}\sum_{i=1}^{n}\delta_{i}-\alpha
n\Big{)}$

with

$\displaystyle v_{\alpha}(y_{i},s_{i})$
$\displaystyle=\mathbb{P}(\theta_{i}\geq\theta_{\alpha}|Y_{i}=y_{i},S_{i}=s_{i})$
$\displaystyle=\frac{\int\int_{\theta_{\alpha}}^{+\infty}\Gamma(s_{i}|r_{i},\sigma_{i}^{2}/r_{i})\phi((y_{i}-\theta)/\sqrt{\sigma^{2}/T_{i}})/\sqrt{\sigma^{2}/T_{i}}dG(\theta,\sigma^{2})}{\int\int\Gamma(s_{i}|r_{i},\sigma_{i}^{2}/r_{i})\phi((y_{i}-\theta)/\sqrt{\sigma^{2}/T_{i}})/\sqrt{\sigma^{2}/T_{i}}dG(\theta,\sigma^{2})}.$

Taking expectations with respect to the joint distribution of $(Y,S)$, the
optimal Bayes rule solves,

$\underset{\bm{\delta}}{\min}\mathbb{E}\Big{[}\sum_{i=1}^{n}(1-\delta_{i})v_{\alpha}(y_{i},s_{i})\Big{]}+\tau_{1}\Big{(}\mathbb{E}\Big{[}\sum_{i=1}^{n}\Big{\\{}(1-v_{\alpha}(y_{i},s_{i}))\delta_{i}-\gamma\delta_{i}\Big{\\}}\Big{]}\Big{)}+\tau_{2}\Big{(}\mathbb{E}\Big{[}\sum_{i=1}^{n}\delta_{i}\Big{]}-\alpha
n\Big{)}.$

Before characterizing the Bayes rule, we should observe that when variances
$\sigma^{2}$ are not directly observed, the tail probability $v_{\alpha}(Y,S)$
may no longer have the monotonicity property we have described above.

######  Lemma 5.1.

Consider the transformation
$v_{\alpha}(Y,S)=\mathbb{P}(\theta\geq\theta_{\alpha}|Y,S]$, then for fixed
$S=s$, the function $v_{\alpha}(Y,s)$ may not be monotone in $Y$; and for
fixed $Y=y$, the function $v_{\alpha}(y,S)$ may not be monotone in $S$.

######  Proposition 5.2.

For a pair of pre-specified $(\alpha,\gamma)$ such that $\gamma<1-\alpha$, the
optimal Bayes selection rule takes the form

$\delta_{i}^{*}=\mymathbb
1\\{v_{\alpha}(Y,S)\geq\lambda^{*}(\alpha,\gamma)\\}$

where
$\lambda^{*}(\alpha,\gamma)=\max\\{\lambda_{1}^{*}(\alpha,\gamma),\lambda_{2}^{*}(\alpha)\\}$
with

$\lambda_{1}^{*}(\alpha,\gamma)=\min\Big{\\{}\lambda:\mathbb{E}\Big{[}(1-v_{\alpha}(Y,S)-\gamma)\mymathbb
1\\{v_{\alpha}(Y,S)\geq\lambda\\}\Big{]}\leq 0\Big{\\}}$

and

$\lambda_{2}^{*}(\alpha)=\min\Big{\\{}\lambda:\mathbb{P}(v_{\alpha}(Y,S)\geq\lambda)-\alpha\leq
0\Big{\\}}$

Based on the optimal Bayes rule, the selected set is defined as

$\Omega_{\alpha,\gamma}=\\{(Y,S):v_{\alpha}(Y,S)\geq\lambda^{*}(\alpha,\gamma)\\}.$

Note that for each prespecified pair $(\alpha,\gamma)$,
$\Omega_{\alpha,\gamma}$ is just the $\lambda^{*}(\alpha,\gamma)$-superlevel
set of the function $v_{\alpha}(Y,S)$. For any $\alpha_{1}>\alpha_{2}$,
nestedness of the selected sets would mean that the
$\lambda^{*}(\alpha_{2},\gamma)$-superlevel set of the function
$v_{\alpha_{2}}$ must be a subset of the
$\lambda^{*}(\alpha_{1},\gamma)$-superlevel set of the function
$v_{\alpha_{1}}$.

######  Remark 5.3.

The construction and the form of the optimal selection rule may appear to be
very similar to the case where $\sigma_{i}^{2}$ is observed. However, the
crucial difference is that we no longer require the independence between
$\theta$ and $\sigma^{2}$ in this section. In contrast, when $\sigma_{i}^{2}$
is assumed to be directly observed, the independence assumption is critical
for all the derivations. For instance, the non-null proportion, defined as
$\mathbb{P}(\theta_{i}\geq\theta_{\alpha})$, must change for different values
of $\sigma_{i}$ if we allow the distribution of $\theta$ to depend on
$\sigma$.

###  5.1. Two Examples

In this subsection we discuss two explicit examples. In the first example,
$G$, is the classical conjugate prior for the Gaussian model with unknown
means and variances.

######  Example 5.4.

Suppose we have balanced panel data
$y_{i1},\dots,y_{iT}\sim\mathcal{N}(\theta,\sigma^{2})$ with sample means
$Y_{i}=\frac{1}{T}\sum_{t}y_{it}$ and sample variances
$S_{i}=\frac{1}{T-1}\sum_{t}(y_{it}-Y_{i})^{2}$. Further, suppose that
$G(\theta,\sigma^{2})$ takes the normal-inverse-chi-squared form,
$\mbox{NIX}(\theta_{0},\kappa_{0},\nu_{0},\sigma_{0}^{2})=\mathcal{N}(\theta|\theta_{0},\sigma^{2}/\kappa_{0})\chi^{-2}(\sigma^{2}|\nu_{0},\sigma_{0}^{2})$.
Integrating out $\sigma^{2}$, the marginal distribution of $\theta$ is a
generalized $t$ distribution

$\frac{\theta-\theta_{0}}{\sigma_{0}/\sqrt{\kappa_{0}}}\sim t_{\nu_{0}}$

where $t_{\nu_{0}}$ is the $t$-distribution with degree of freedom $\nu_{0}$.
Therefore, the $1-\alpha$ quantile of $\theta$, denoted $\theta_{\alpha}$ is
simply,

$\theta_{\alpha}=\theta_{0}+\frac{\sigma_{0}}{\sqrt{\kappa_{0}}}F^{-1}_{t_{\nu_{0}}}(1-\alpha)$

where $F^{-1}_{t_{\nu_{0}}}$ denotes the quantile function of $t_{\nu_{0}}$.

Conjugacy of the distribution $G$ implies that the posterior distribution of
$(\theta,\sigma^{2}|Y,S)$ follows
$NIX(\theta_{T},\kappa_{T},\nu_{T},\sigma_{T}^{2})=\mathcal{N}(\theta|\theta_{T},\sigma_{T}^{2}/\kappa_{T})\chi^{-2}(\sigma^{2}|\nu_{T},\sigma_{T}^{2})$
with

$\displaystyle\nu_{T}$ $\displaystyle=\nu_{0}+T$ $\displaystyle\kappa_{T}$
$\displaystyle=\kappa_{0}+T$ $\displaystyle\theta_{T}$
$\displaystyle=\frac{\kappa_{0}\theta_{0}+TY}{\kappa_{T}}$
$\displaystyle\sigma_{T}^{2}$
$\displaystyle=\frac{1}{\nu_{T}}\Big{(}\nu_{0}\sigma_{0}^{2}+(T-1)S+\frac{T\kappa_{0}}{\kappa_{0}+T}(\theta_{0}-Y)^{2}\Big{)}.$

Integrating out $\sigma^{2}$, the marginal posterior of $\theta$ again follows
a generalized $t$-distribution,

$\frac{\theta-\theta_{T}}{\sigma_{T}/\sqrt{\kappa_{T}}}\sim t_{\nu_{T}},$

where $t_{\nu_{T}}$ is the $t$-distribution with degree of freedom $\nu_{T}$.
It is then clear that the posterior mean of $\theta$, is simply a linear
function of $Y$ and independent of $S$,

$\mathbb{E}[\theta|Y,S]=\theta_{T}=\frac{\kappa_{0}\theta_{0}+TY}{\kappa_{T}},$

and the posterior tail probability is given by,

$v_{\alpha}(Y,S)=\mathbb{P}(\theta\geq\theta_{\alpha}|Y,S)=\mathbb{P}\Big{(}\frac{\theta-\theta_{T}}{\sigma_{T}/\sqrt{\kappa_{T}}}\geq\frac{\theta_{\alpha}-\theta_{T}}{\sigma_{T}/\sqrt{\kappa_{T}}}|Y,S\Big{)}=1-F_{t_{\nu_{T}}}\Big{(}\frac{\theta_{\alpha}-\theta_{T}}{\sigma_{T}/\sqrt{\kappa_{T}}}\Big{)}.$

Figure 5.1. The left panel shows level curves of the posterior mean (marked as
red dashed lines) and posterior tail probability (marked as black solid lines)
for the normal model with $(\theta,\sigma^{2})\sim NIX(0,1,6,1)$ and panel
time dimension $T=9$. The right plot shows the boundary of the selection
region based on posterior mean ranking (marked as the red dashed line) and the
posterior tail probability ranking (marked as the solid black line) with
$\alpha=5\%$ and $\gamma=10\%$.

Figure 5.2. Selection set comparison for one sample realization from the
normal model with $(\theta,\sigma^{2})\sim NIX(0,1,6,1)$ and panel time
dimension $T=9$: The left panel shows in grey circles the agreed selected
elements by both the posterior mean and posterior tail probability criteria
under the capacity constraint, extra elements selected by the posterior mean
are marked in green and extra elements selected by the posterior tail
probability rule are marked in red. The right panel shows the comparison of
the selected sets under both the capacity and FDR constraint with $\alpha=5\%$
and $\gamma=10\%$.

To illustrate this case, suppose $\theta_{0}=0$, $\kappa_{0}=1$ ,
$\sigma_{0}^{2}=1$ and $\nu_{0}=6$ and $T=9$, it can be verified that
$v_{\alpha}(Y,S)$ is in fact a monotone function of $Y$ for each fixed $S$ and
any $\alpha>0$, hence in this example we can invert the function
$v_{\alpha}(y,s)$ to obtain the level curves. The left panel of Figure 5.1
shows the level curves for $v_{\alpha}(Y,S)$ and $\mathbb{E}(\theta|Y,S)$ for
$\alpha=5\%$. It is clear that the posterior mean is a constant function of
$S$, while the posterior tail probability exhibits more exotic behaviour with
respect to $S$, especially for more extreme values of $Y$. If we fix
$S=s_{0}$, then $v_{\alpha}(Y,s_{0})$ is an increasing function of $Y$. On the
other hand, fixing $Y=y_{0}$, for small $y_{0}$ implies that
$v_{\alpha}(y_{0},S)$ is a increasing function of $S$, while for $y_{0}$
large, $v_{\alpha}(y_{0},S)$ becomes a decreasing function of $S$.

A capacity constraint of size $\alpha$ implies the thresholding rule,

$\mathbb{P}(v_{\alpha}(Y,S)\geq\lambda_{2}^{*})=\alpha,$

while FDR control at level $\gamma$ leads to a cutoff value $\lambda_{1}^{*}$
defined as

$\gamma=\mathbb{E}[(1-v_{\alpha}(Y,S))\mymathbb
1\\{v_{\alpha}(Y,S)\geq\lambda_{1}^{*}\\}]/\mathbb{P}(v_{\alpha}(Y,S)\geq\lambda_{1}^{*}).$

The larger of the two thresholds, denoted
$\lambda^{*}=\max\\{\lambda_{1}^{*},\lambda_{2}^{*}\\}$ defines the selection
region based on posterior tail probability ranking
$\Omega_{\alpha,\gamma}=\\{(Y,S):v_{\alpha}(Y,S)\geq\lambda^{*}\\}$. For
$\alpha=5\%$ and $\gamma=10\%$, the selection region based on the tail
probability rule is $\\{(Y,S):v_{\alpha}(Y,S)\geq 0.72\\}$. The posterior mean
ranking is defined as $\\{(Y,S):\mathbb{E}[\theta|Y,S]\geq 2.2\\}$. These
selection boundaries are depicted as the red dashed line and black solid line
respectively in the right panel of Figure 5.1. In this case, the FDR
constraint binds. If only the capacity constraint were in place, we would have
a cut-off for tail probability at $0.40$ and the cut-off for posterior mean at
$1.84$. Figure 5.2 further shows the comparison of the selected set based on a
sample realization from the model.

In Appendix B we consider a more complex bivariate discrete example that
illustrates somewhat more exotic behavior of the decision boundaries and
compares performance of several different ranking and selection rules.

###  5.2. Variants of the unknown variance model

We’ve assumed that the only scale heterogeneity is driven by $\sigma_{i}$ in
the above model, but often there may be more heteroskedasticity that should be
allowed in $\epsilon_{it}$. Here we consider the variant where

$Y_{it}=\theta_{i}+\sigma_{i}\epsilon_{it},\quad\epsilon_{it}\sim\mathcal{N}(0,1/w_{it}),\quad(\theta_{i},\sigma_{i}^{2})\sim
G$

we also assume that $w_{it}\sim H$ are known quantity and are independent from
$(\theta_{i},\sigma_{i}^{2})$. Denoting $w_{i}=\sum_{t=1}^{T_{i}}w_{it}$, the
sufficient statistics now takes the form
$Y_{i}=\sum_{t=1}^{T_{i}}w_{it}Y_{it}/w_{i}$ and
$S_{i}=(T_{i}-1)^{-1}\sum_{t=1}^{T_{i}}(Y_{it}-Y_{i})^{2}$. It is easy to show
that
$Y_{i}|\theta_{i},\sigma_{i}^{2}\sim\mathcal{N}(\theta_{i},\sigma_{i}^{2}/w_{i})$
and $S_{i}|\sigma_{i}^{2}$ follows Gamma distribution with shape parameter
$r_{i}=(T_{i}-1)/2$ and scale parameter $\sigma_{i}^{2}/r_{i}$. The decision
rules now become a function of the tuple $(Y_{i},S_{i},w_{i})$, for instance
the tail probability can be specified as

$v_{\alpha}(y_{i},s_{i},w_{i})=\mathbb{P}(\theta_{i}\geq\theta_{\alpha}|y_{i},s_{i},w_{i})=\frac{\int_{\theta_{\alpha}}^{+\infty}f(y_{i}|\theta,\sigma^{2}/w_{i})\Gamma(s_{i}|r_{i},\sigma^{2}/r_{i})dG(\theta,\sigma^{2})}{\int_{-\infty}^{+\infty}f(y_{i}|\theta,\sigma^{2}/w_{i})\Gamma(s_{i}|r_{i},\sigma^{2}/r_{i})dG(\theta,\sigma^{2})}.$

and posterior mean takes the form

$\mathbb{E}[\theta_{i}|y_{i},s_{i},w_{i}]=\int\theta
f(y_{i},s_{i}|\theta,\sigma^{2},w_{i})dG(\theta,\sigma^{2}).$

The threshold values under either the capacity or the FDR constraint can be
worked out in a similar fashion. For any ranking statistics
$\delta(Y_{i},S_{i},w_{i})$ together with a decision rule $\mymathbb
1\\{\delta(Y_{i},S_{i},w_{i})\geq\lambda\\}$, the capacity constraint requires
choosing a thresholding value $\lambda_{2}^{*}(\alpha)$ such that

$\alpha=\int\int\mymathbb
1\\{\delta(y,s,w)\geq\lambda_{2}^{*}(\alpha)\\}f(y,s|\theta,\sigma^{2},w)dG(\theta,\sigma^{2})dH(w),$

while the thresholding value to control in addition the FDR rate under size
$\gamma$ requires solving for $\lambda_{1}^{*}(\alpha,\gamma)$ such that

$\gamma=\frac{\mathbb{P}(\delta(y,s,w)\geq\lambda_{1}^{*}(\alpha,\gamma);\theta<\theta_{\alpha})}{\mathbb{P}(\delta(y,s,w)\geq\lambda_{1}^{*}(\alpha,\gamma))}$

which can be further represented as

$\gamma=\frac{\int\int\mymathbb
1\\{\delta(y,s,w)\geq\lambda_{1}^{*}(\alpha,\gamma)\\}(1-\alpha)f_{0}(y,s|\theta,\sigma^{2},w)dG(\theta,\sigma^{2})dH(w)}{\int\int\mymathbb
1\\{\delta(t,s,w)\geq\lambda_{1}^{*}(\alpha,\gamma)\\}f(y,s|\theta,\sigma^{2},w)dG(\theta,\sigma^{2})dH(w)}$

where $f_{0}(y,s|\theta,\sigma^{2},w)$ is the density of $(y,s)$ under the
null hypothesis $\theta<\theta_{\alpha}$.

We can again consider selection region as those plotted in Figure 5.1 and
Figure B.1 to appreciate how different decision criteria determines the
selection. As soon as the ranking statistics depends on $w$, the selection
region of the thresholding rule $\mymathbb
1\\{\delta(y,s,w)\geq\lambda^{*}\\}$ will also depend on the magnitude of $w$.

##  6\. Asymptotic Adaptivity

The previous sections propose the optimal Bayes rule for minimizing the
expected number of missed discoveries subject to both capacity and FDR
constraint under various modeling environments. In each of these environments,
the optimal Bayes rule takes the form
$\delta^{*}=1\\{v_{\alpha}\geq\lambda^{*}\\}$ where $v_{\alpha}$ is defined as
the posterior probability of $\theta\geq\theta_{\alpha}$ conditional on the
data. The thresholding value $\lambda^{*}$ is defined to satisfy both the
capacity and FDR constraints. The optimal Bayes rule involves several unknown
quantities, in particular the $v_{\alpha}$’s and the optimal thresholding
value, $\lambda^{*}$, that require knowledge on the distribution of
$\theta_{i}$ or the joint distribution of $(\theta_{i},\sigma_{i}^{2})$ when
the variances are latent variables. For estimating this distribution of the
latent variables, we propose a plug-in procedure that is very much in the
spirit of empirical Bayes methods pioneered by Robbins.56. In this Section we
also establish that the resulting feasible rules achieve asymptotic validity
and asymptotically attain the same performance as the infeasible Bayes rule.

We begin by discussing properties of the Oracle procedure assuming that
$v_{\alpha}$ is known and we only need to estimate the optimal thresholding
value. We establish asymptotic validity of this Oracle procedure and then
propose a plug-in method for both $v_{\alpha}$ and the thresholding value
thereby establishing the asymptotic validity of the empirical rule. Before
presenting the formal results, we introduce regularity conditions that will be
required. We distinguish two cases depending on whether the $\sigma_{i}^{2}$’s
are observed.

######  Assumption 1.

  1. (1)

(Variances observed) $\\{Y_{i},\sigma_{i}^{2},\theta_{i}\\}$ are independent
and identically distributed with a joint distribution with $\sigma_{i}^{2}$
and $\theta_{i}$ independent. The random variables $\theta_{i}$ and
$\sigma_{i}^{2}$ have positive densities with respect to Lebesgue measure on a
compact set $\Theta\subset\mathbb{R}$ and
$[\underline{\sigma}^{2},\bar{\sigma}^{2}]$ respectively for some
$\underline{\sigma}^{2}>0$ and $\bar{\sigma}^{2}<+\infty$.

  2. (2)

(Variance unobserved) Let $S_{i}$ be an individual sample variance based on
$T$ repeated measurements and $Y_{i}$ be the sample means with $T\geq 4$.
Suppose further that $\\{Y_{i},S_{i},\theta_{i},\sigma_{i}\\}$ are independent
and identically distributed and that the random variables
$\\{\theta_{i},\sigma_{i}^{2}\\}$ have a joint distribution $G$ with a joint
density positive everywhere on its support.

###  6.1. Optimal thresholding

Whether $\sigma_{i}^{2}$ is observed or estimated, the optimal thresholding
value can be defined in a unified manner by
$\lambda^{*}=\max\\{\lambda_{1}^{*},\lambda_{2}^{*}\\}$ with

$\displaystyle\lambda_{1}^{*}$ $\displaystyle=\inf\\{t\in(0,1),H_{v}(t)\geq
1-\alpha\\}$ $\displaystyle\lambda_{2}^{*}$
$\displaystyle=\inf\\{t\in(0,1),Q(t)\leq\gamma\\}$

where $H_{v}$ denotes the cumulative distribution of either
$v_{\alpha}(y_{i},\sigma_{i})$ or $v_{\alpha}(y_{i},s_{i})$, induced by the
marginal distribution of the data, either as the pair $\\{y_{i},\sigma_{i}\\}$
when variances are observed or the pair $\\{y_{i},s_{i}\\}$ otherwise. Hence
$\lambda_{1}^{*}$ is the $1-\alpha$ quantile of $H_{v}$.

The function $Q(t)$ is defined as
$Q(t)=\mathbb{E}[(1-v_{\alpha})1\\{v_{\alpha}\geq
t\\}]/\mathbb{E}[1\\{v_{\alpha}\geq t\\}]$. Its formulation recalls
Proposition 5.2 and the existence of $\lambda_{2}^{*}$ is guaranteed as long
as $\alpha<1-\gamma$. The thresholding value is also equivalent to those
defined in Proposition 3.4 and Proposition 4.1. In particular, the
thresholding values $t_{1}^{*}$ and $t_{2}^{*}$ in Proposition 3.4 are cast in
terms of $Y$ directly and it is easy to see
$\lambda_{j}^{*}=v_{\alpha}(t_{j}^{*})$ for $j=1,2$ when variances are
homogeneous. In a similar spirit, the explicit formulae for $\lambda_{1}^{*}$
and $\lambda_{2}^{*}$ in Proposition 4.1 are a result of invoking the
monotonicity of $v_{\alpha}(y,\sigma)$ with respect to $y$ for each fixed
value of $\sigma$. The function $Q(t)$ is the mFDR of the procedure
$\delta=1\\{v_{\alpha}\geq t\\}$ for any $\alpha\in(0,1)$, and is
monotonically decreasing in $t$. Monotonicity of $Q(t)$ is crucial to justify
this thresholding procedure insuring that either the capacity constraint or
the mFDR constraint must be binding. cao2013optimal have observed that a
sufficient condition for monotonicity for a broad class of multiple testing
procedures is that the ratio of densities under the null and alternative of
the test statistics employed for ranking be monotone and discuss the
consequences of the violation of this condition. For the posterior tail
probability criterion this monotone likelihood ratio condition, as we will
see, can be verified directly.

To see this, recall that mFDR is defined as
$\sum_{i=1}^{n}\mathbb{P}[\delta_{i}=1,\theta_{i}<\theta_{\alpha}]/\sum_{i=1}^{n}\mathbb{P}(\delta_{i}=1)$.
It suffices to show that
$\mathbb{P}[\delta_{i}=1,\theta_{i}<\theta_{\alpha}]=\mathbb{E}[(1-v_{\alpha,i})\delta_{i}]$.
Since $v_{\alpha,i}=P[\theta_{i}\geq\theta_{\alpha}|D_{i}]=\alpha
f_{1}(D_{i})/f(D_{i})$ where $D_{i}$ is the individual data being either
$\\{y_{i},\sigma_{i}\\}$ or $\\{y_{i},s_{i}\\}$ depending on the model and
$f_{1}$ is the marginal density of the data when
$\theta_{i}\geq\theta_{\alpha}$ and $f$ is the marginal density of $D_{i}$.
Then it is easy to see that
$\mathbb{P}[\delta_{i}=1,\theta_{i}<\theta_{\alpha}]=(1-\alpha)\int
1\\{v_{\alpha,i}\geq t\\}f_{0}(D_{i})dD_{i}=\int 1\\{v_{\alpha,i}\geq
t\\}(1-v_{\alpha,i})f(D_{i})dD_{i}=\mathbb{E}[(1-v_{\alpha,i})1\\{v_{\alpha,i}\geq
t\\}]$. Then $Q(t)=\int_{t}^{1}(1-v)h_{v}dv/\int_{t}^{1}h_{v}dv$ where $h_{v}$
is the density function of $v_{\alpha}$. Monotonicity of $Q(t)$ can be easily
verified by showing that the derivative with respect to $t$ of the right hand
side quantity is nonpositive.

###  6.2. Oracle Procedure

The only unknown quantity in the Oracle procedure is the thresholding value
and we now discuss how to estimate it to achieve asymptotic validity. $H_{v}$
and $Q$ can be estimated by the following quantities:

$\displaystyle H_{n}(t)$
$\displaystyle=\frac{1}{n}\sum_{i=1}^{n}1\\{v_{\alpha,i}\leq t\\}$
$\displaystyle Q_{n}(t)$
$\displaystyle=\frac{\sum_{i=1}^{n}(1-v_{\alpha,i})1\\{v_{\alpha,i}\geq
t\\}}{\sum_{i=1}^{n}1\\{v_{\alpha,i}\geq t\\}}$

and the associated thresholding values are then defined as

$\displaystyle\lambda_{1n}$ $\displaystyle=\inf\\{t\in[0,1],H_{n}(t)\geq
1-\alpha\\}$ $\displaystyle\lambda_{2n}$
$\displaystyle=\inf\\{t\in[0,1],Q_{n}(t)\leq\gamma\\}$
$\displaystyle\lambda_{n}$ $\displaystyle=\max\\{\lambda_{1n},\lambda_{2n}\\}$

######  Theorem 6.1.

(Asymptotic validity of the Oracle procedure) Under Assumption 1, the
procedure $\delta_{i}=1\\{v_{\alpha,i}\geq\lambda_{n}\\}$ asymptotically
controls the false discovery rate below $\gamma$ and the expected proportion
of rejections below $\alpha$ for any $(\alpha,\gamma)\in[0,1]^{2}$ and
$\gamma<1-\alpha$ when $n\to\infty$, more specifically

$\displaystyle\underset{n\to\infty}{\limsup}\mathbb{E}\Big{[}\frac{\sum_{i=1}^{n}1\\{\theta_{i}<\theta_{\alpha},v_{\alpha,i}\geq\lambda_{n}\\}}{\sum_{i=1}^{n}1\\{v_{\alpha,i}\geq\lambda_{n}\\}\bigvee
1}\Big{]}\leq\gamma$
$\displaystyle\underset{n\to\infty}{\limsup}\mathbb{E}\Big{[}\frac{1}{n}\sum_{i=1}^{n}1\\{v_{\alpha,i}\geq\lambda_{n}\\}\Big{]}\leq\alpha$

###  6.3. Adaptive procedure

In practise the posterior tail probability also involves the unknown quantity
$\theta_{\alpha}=G^{-1}(1-\alpha)$ that needs to be estimated. We propose a
plug-in estimator in the spirit of the empirical Bayes method: estimating $G$
by its nonparametric maximum likelihood estimator $\hat{G}_{n}$ and estimating
$\theta_{\alpha}$ as its $1-\alpha$ quantile.

Consistency of the nonparametric maximum likelihood estimator, $\hat{G}_{n}$,
was first proven by KW using Wald type arguments. A Hellinger risk bound for
the associated marginal density estimate and adaptivity of $\hat{G}_{n}$ and a
self-regularization property have been recently established in
saha2020nonparametric and PW20. In particular, the following established
result, stated here as an assumption, is crucial for establishing the
asymptotic validity of the adaptive procedure.

######  Assumption 2.

The nonparametric maximum likelihood estimator $\hat{G}_{n}$ is strongly
consistent for $G$. That is, for all continuity point $k$ of $G$,
$\hat{G}_{n}(k)\to G(k)$ almost surely when $n\to\infty$. Furthermore, the
estimated marginal (mixture) density converges almost surely in Hellinger
distance.

When variances are homogeneous or when variances are unknown but we have
longitudinal data so that we have a mixture model for the pair
$\\{Y_{i},S_{i}\\}$, the Hellinger convergence is established in
van1993hellinger. When variances are heterogeneous but known, the Hellinger
bound for marginal density is established recently in jiang2020general.

The plug-in estimators for the posterior tail probability,
$v_{\alpha}(y_{i},\sigma_{i})$ when variances are known or
$v_{\alpha}(y_{i},s_{i})$ when variances are unknown is then defined
respectively as

$\displaystyle\hat{v}_{\alpha}(y_{i},\sigma_{i})$
$\displaystyle=\int_{\hat{\theta}_{\alpha}}^{+\infty}\frac{1}{\sigma_{i}}\phi((y_{i}-\theta)/\sigma_{i})d\hat{G}_{n}(\theta)/\int_{-\infty}^{+\infty}\frac{1}{\sigma_{i}}\phi((y_{i}-\theta)/\sigma_{i})d\hat{G}_{n}(\theta)$
$\displaystyle\hat{v}_{\alpha}(y_{i},s_{i})$
$\displaystyle=\int_{\hat{\theta}_{\alpha}}^{+\infty}f(y_{i},s_{i}|\theta,\sigma)d\hat{G}_{n}(\theta,\sigma^{2})/\int_{-\infty}^{+\infty}f(y_{i},s_{i}|\theta,\sigma)d\hat{G}_{n}(\theta,\sigma^{2})$

where $\phi$ is the standard normal density and $f$ is the density function
for $(y_{i},s_{i})$ which is a product of Gaussian and gamma density.
Abbreviating the estimated posterior tail probability by $\hat{v}_{\alpha,i}$,
we mimic the Oracle procedure and estimate the thresholding value by
$\hat{\lambda}_{n}=\max\\{\hat{\lambda}_{1n},\hat{\lambda}_{2n}\\}$, where,

$\displaystyle\hat{\lambda}_{1n}$
$\displaystyle=\inf\\{t\in[0,1]:\frac{1}{n}\sum_{i=1}^{n}1\\{\hat{v}_{\alpha,i}\leq
t\\}\geq 1-\alpha\\}$ $\displaystyle\hat{\lambda}_{2n}$
$\displaystyle=\inf\\{t\in[0,1]:\frac{\sum_{i=1}^{n}(1-\hat{v}_{\alpha,i})1\\{\hat{v}_{\alpha,i}\geq
t\\}}{\sum_{i=1}^{n}1\\{\hat{v}_{\alpha,i}\geq t\\}}\geq\gamma\\}$

######  Theorem 6.2.

(Asymptotic validity of adaptive procedure) Under Assumptions 1 and 2, the
adaptive procedure $\delta_{i}=1\\{\hat{v}_{\alpha,i}\geq\hat{\lambda}_{n}\\}$
asymptotically controls the false discovery rate below $\gamma$ and the
expected proportion of rejections below $\alpha$ for any
$(\alpha,\gamma)\in[0,1]^{2}$ with $\alpha<1-\gamma$ when $n\to\infty$, more
specifically

$\displaystyle\underset{n\to\infty}{\limsup}\mathbb{E}\Big{[}\frac{\sum_{i=1}^{n}1\\{\theta_{i}<\theta_{\alpha},\hat{v}_{\alpha,i}\geq\hat{\lambda}_{n}\\}}{\sum_{i=1}^{n}1\\{\hat{v}_{\alpha,i}\geq\hat{\lambda}_{n}\\}\bigvee
1}\Big{]}\leq\gamma$
$\displaystyle\underset{n\to\infty}{\limsup}\mathbb{E}\Big{[}\sum_{i=1}^{n}1\Big{\\{}\frac{1}{n}\hat{v}_{\alpha,i}\geq\hat{\lambda}_{n}\Big{\\}}\Big{]}\leq\alpha$

It is clear that given the Lagrangian formulation of the compound decision
problem, it can be viewed equivalently as a constrained optimization problem.
See also the discussion in Remark 3.3. We seek to maximize power defined as
$\beta(t):=\mathbb{P}(\theta_{i}\geq\theta_{\alpha},\delta_{i}=1)/\alpha$
subject to two constraints: the first marginal FDR rate and the other on the
selected proportion. For each fixed pair of $\\{\alpha,\gamma\\}$, the optimal
Bayes rule achieves the best power among all thresholding procedures that
respect the two constraints. The next theorem establishes that our feasible,
adaptive procedure achieves the same power as the oracle rule asymptotically.

######  Theorem 6.3.

Under Assumption 1 and 2, the adaptive procedure
$\delta_{i}=1\\{\hat{v}_{\alpha,i}\geq\hat{\lambda}_{n}\\}$ attains the same
power as the optimal Bayes rule asymptotically. In particular, as
$n\to\infty$,

$\displaystyle\frac{\sum_{i=1}^{n}1\\{\theta_{i}\geq\theta_{\alpha},\hat{v}_{\alpha,i}\geq\hat{\lambda}_{n}\\}}{\sum_{i=1}^{n}1\\{\theta_{i}\geq\theta_{\alpha}\\}}\overset{p}{\to}\beta(\lambda^{*})$

This result is supported by the simulation evidence presented in the next
section.

##  7\. Simulation Evidence

In this Section we describe two small simulation exercises designed to
illustrate performance of several competing methods for ranking and selection.
As a benchmark for evaluating performance we consider several Oracle methods
that presume knowledge of the true distribution, $G$, generating the
$\theta$’s as well as several feasible methods that rely on estimation of $G$.
These are contrasted with more traditional methods that are based on linear
shrinkage rules of the Stein type. The linear shrinkage rule is the posterior
mean of $\theta$ under the assumption that $G$ follows a Gaussian distribution
with unknown mean and variance parameters. This is the commonly used estimator
for ranking and selection in applied work, notably Chetty, Friedman and
Rockoff (2014a, 2014b) for teacher evaluation and chetty2018impacts for
studying intergenerational mobility. Typically the linear shrinkage estimator
is used in the context of heterogeneous known variances, this will be the
model we focus on in our simulation experiments. The linear shrinkage formula
defined in (2.1) easily adapts to the heterogeneous variances case and leads
to the Jame-Stein shrinkage rule with heterogeneous known variances.
EfronMorris introduced some further modifications. As we have already
demonstrated, when variances are heterogeneous, the linear shrinkage estimator
provides a different ranking compared to posterior tail probability rules. A
further complication arises when we seek procedures that also control for
false discovery. To estimate the false discovery rate for different
thresholding values we requires knowledge of $G$. If the Gaussian assumption
on $G$ underlying the linear shrinkage rules is misplaced, it may lead to an
inaccurate estimator for FDR, and consequently to procedures that fail to
control for false discovery.

Performance will be evaluated primarily on the basis of power, which we define
as the proportion of the individuals whose true $\theta_{i}$ exceeds the
cutoff $\theta_{\alpha}=G^{-1}(1-\alpha)$, who are actually selected. This is
the sample counterpart of
$\mathbb{P}(\delta_{i}=1,\theta_{i}\geq\theta_{\alpha})/\mathbb{P}(\theta_{i}\geq\theta_{\alpha})$.
FDR is calculated as the sample counterpart of
$\mathbb{P}(\delta_{i}=1,\theta_{i}<\theta_{\alpha})/\mathbb{P}(\delta_{i}=1)$,
that is the proportion of selected individuals whose true $\theta_{i}$ falls
below the threshold. While our selection rules are intended to constrain FDR
below the $\gamma$ threshold, as in other testing problems they are not always
successful in this objective in finite samples so empirical power comparisons
must be interpreted cautiously in view of this. Nonetheless, asymptotic
validity is assured by the results in Section 6. We compare performance for
three distinct $\alpha$ levels, $\\{0.05,0.10,0.15\\}$ and three $\gamma$
levels $\\{0.05,0.10,0.15\\}$.

###  7.1. The Student $t$ Setting

Figure 7.1. Power Performance for Several Selection Rules with Student $t$
Signal. Capacity and FDR constraints are indicated at the top of each panel in
the Figure.

Our first simulation setting focuses on the effect of tail behavior of the
distribution on performance of competing rules. For these simulations we take
$G$ to be a discrete approximation to Student $t$ distributions with degrees
of freedom in the set $\\{1,2,3,5,10\\}$, and supported on the interval
$[-20,20]$. The scale parameters of the Gaussian noise contribution are
independent and uniformly distributed on the interval $[0.5,1.5]$. We report
power performance for several alternative ranking and selection rules:

OTP:

Oracle Tail Probability Rule

OPM:

Oracle Posterior Mean Rule

Efron:

Efron Tail Probability Rule

KWs:

Kiefer-Wolfowitz Smoothed Tail Probability Rule

EM:

Efron and Morris (1973) Linear Shrinkage Rule

We illustrate the results in Figure 7.1, where we plot empirical power against
degrees of freedom of the $t$ distribution for a selected set of values for
the capacity constraint, $\alpha\in\\{0.05,0.10,0.15\\}$ and FDR constraint,
$\gamma\in\\{0.05,0.10,0.15\\}$ as indicated at the top of each panel of the
figure. The most striking conclusion from this exercise is the dramatic
decrease in power as we move toward the Gaussian distribution. At the Cauchy,
$t_{1}$, power is quite respectable for all choices of $\alpha$ and $\gamma$,
but power declines rapidly as the degrees of freedom increases, reenforcing
our earlier conclusion that the Gaussian case is extremely difficult. We would
stress, in view of this finding, that classical linear shrinkage procedures
designed for the Gaussian setting are poorly adapted to heavy tailed settings
in which the reliability of selection procedures is potentially greatest.

Careful examination of this figure reveals that there is a slight advantage to
the posterior tail probability rules over the posterior mean procedures, both
for the Oracle rules and for our feasible procedures. There is surprisingly
little sacrifice in power in moving from the Oracle methods to the Efron or
Kiefer-Wolfowitz rules. The Efron and Morris selection rule is very
competitive in the almost Gaussian, $t_{10}$ setting but sacrifices
considerable power in the lower degrees of freedom settings due to the
misspecification of the distribution $G$ and consequent inaccurate estimation
of the false discovery rate.

###  7.2. A Teacher Value-Added Setting

Our second simulation setting is based on a discrete approximation of the data
structure employed in GGM to study teacher value-addded methods. Several
longitudinal waves of student test scores from the Los Angeles Unified School
District were combined in this study. Here we abstract from many features of
the full longitudinal structure of this data, and focus instead on simulating
performance of several selection methods. We maintain our standard known
variance model in which we observe
$Y_{i}\sim\mathcal{N}(\theta_{i},\sigma_{i}^{2})$ with $\theta_{i}$’s drawn
iidly from a distribution $\tilde{G}$ estimated by GGM. This distribution was
estimated from the full longitudinal LA sample using the nonparametric maximum
likelihood estimator of Kiefer and Wolfowitz and then smoothed slightly by
convolution with a biweight kernel and illustrated in the left panel of Figure
7.2. Variances, in keeping with our hypothesis in Section 4, are drawn from a
distribution with density illustrated in the right panel of Figure 7.2. We
focus on selection from the left tail of the resulting distribution since it
is those teachers whose jobs are endangered by recent policy recommendations
in the literature. (see for instance hanushek2011economic).

Figure 7.2. Densities of “True” (Mean) Ability and Standard Deviation for the
Teacher Value Added Simulations

We draw samples of size 10,000 from the foregoing distribution and compute
performance measures based on 100 replications. The fitted densities for this
simulation exercise are based on a sample of roughly 11,000 teachers, so the
simulation sample size is chosen to be commensurate with this. In Table 7.1 we
report power, FDR and the proportion selected by ten selection rules. The
Oracle rules, OTP and OPM, based ranking by the tail probability and posterior
mean criteria can be considered benchmarks for the remaining feasible
procedures. Only the Oracle procedures can be considered reliable from the
perspective of adhering to the capacity and FDR constraints. Consequently,
some caution is required in the interpretation of the power comparisons since
feasible procedures can exhibit good power at the expense of violating these
constraints. This is analogous to the common difficulty in interpreting power
in testing problems when different procedures have differing size. When FDR is
constrained to 5%, even the Oracle is only able to select about half of the
deserving individuals; OTP is consistently preferable to OPM as expected and
power performance improves somewhat as the capacity constraint is relaxed.
Among the feasible $G$-modeling selection procedures, the Efron rules have
good power performance, but fail to meet the FDR constraints. We conjecture
that somewhat less aggressive smoothing than the default, $df=5,c_{0}=0.1$
might help to rectify this. In contrast the smoothed Kiefer-Wolfowitz rules
are somewhat overly conservative in meeting the FDR constraints and might
benefit from somewhat less aggressive smoothing.

Among the other procedures the linear posterior mean rule, LPM, as employed by
Chetty, Friedman and Rockoff (2014a, 2014b), and the linear posterior mean
rule, EM, of EfronMorris behaves identically and exhibit somewhat erratic FDR
control due to the misspecified Gaussian assumption on $G$, and thus leads to
weaker power performance. As a further comparison, when the linear shrinkage
rules are implemented without any FDR constraint, denoted LPM* and EM* in the
table, as they typically would be used in practice, the false discovery
proportion is considerably higher than the targeted $\gamma$. We also report
the performance of MLE and P-value rules, implemented without FDR control;
again both yield a higher FDR rate, making it futile to evaluate their power
performance.

| $\gamma=5\%$ |  | $\gamma=10\%$
---|---|---|---
| $\alpha=1\%$ | $\alpha=3\%$ | $\alpha=5\%$ | $\alpha=10\%$ |  | $\alpha=1\%$ | $\alpha=3\%$ | $\alpha=5\%$ | $\alpha=10\%$
Power |  |  |  |  |  |  |  |  |
OTP | $0.394$ | $0.520$ | $0.554$ | $0.626$ |  | $0.494$ | $0.625$ | $0.661$ | $0.733$
OPM | $0.365$ | $0.492$ | $0.521$ | $0.599$ |  | $0.484$ | $0.620$ | $0.654$ | $0.731$
ETP | $0.435$ | $0.540$ | $0.580$ | $0.657$ |  | $0.540$ | $0.647$ | $0.688$ | $0.759$
KWTP | $0.355$ | $0.477$ | $0.521$ | $0.614$ |  | $0.452$ | $0.583$ | $0.631$ | $0.723$
EPM | $0.398$ | $0.511$ | $0.552$ | $0.632$ |  | $0.528$ | $0.642$ | $0.683$ | $0.758$
KWPM | $0.325$ | $0.447$ | $0.492$ | $0.588$ |  | $0.440$ | $0.576$ | $0.624$ | $0.719$
LPM | $0.162$ | $0.341$ | $0.418$ | $0.689$ |  | $0.246$ | $0.460$ | $0.542$ | $0.805$
LPM* | $0.726$ | $0.781$ | $0.796$ | $0.829$ |  | $0.726$ | $0.781$ | $0.796$ | $0.829$
EM | $0.162$ | $0.341$ | $0.418$ | $0.689$ |  | $0.246$ | $0.460$ | $0.542$ | $0.805$
EM* | $0.726$ | $0.781$ | $0.796$ | $0.829$ |  | $0.726$ | $0.781$ | $0.796$ | $0.829$
MLE | $0.699$ | $0.768$ | $0.787$ | $0.824$ |  | $0.699$ | $0.768$ | $0.787$ | $0.824$
P-val | $0.374$ | $0.478$ | $0.535$ | $0.635$ |  | $0.374$ | $0.478$ | $0.535$ | $0.635$
FDR |  |  |  |  |  |  |  |  |
OTP | $0.050$ | $0.050$ | $0.051$ | $0.051$ |  | $0.103$ | $0.103$ | $0.100$ | $0.101$
OPM | $0.047$ | $0.050$ | $0.051$ | $0.053$ |  | $0.103$ | $0.101$ | $0.101$ | $0.102$
ETP | $0.070$ | $0.059$ | $0.061$ | $0.064$ |  | $0.128$ | $0.115$ | $0.117$ | $0.119$
KWTP | $0.035$ | $0.037$ | $0.041$ | $0.048$ |  | $0.082$ | $0.081$ | $0.085$ | $0.096$
EPM | $0.063$ | $0.057$ | $0.059$ | $0.062$ |  | $0.129$ | $0.115$ | $0.116$ | $0.118$
KWPM | $0.038$ | $0.040$ | $0.045$ | $0.051$ |  | $0.081$ | $0.084$ | $0.087$ | $0.097$
LPM | $0.016$ | $0.025$ | $0.033$ | $0.083$ |  | $0.031$ | $0.048$ | $0.061$ | $0.151$
LPM* | $0.276$ | $0.226$ | $0.207$ | $0.172$ |  | $0.276$ | $0.226$ | $0.207$ | $0.172$
EM | $0.016$ | $0.025$ | $0.033$ | $0.083$ |  | $0.031$ | $0.048$ | $0.061$ | $0.151$
EM* | $0.276$ | $0.225$ | $0.207$ | $0.172$ |  | $0.276$ | $0.225$ | $0.207$ | $0.172$
MLE | $0.304$ | $0.238$ | $0.216$ | $0.177$ |  | $0.304$ | $0.238$ | $0.216$ | $0.177$
P-val | $0.627$ | $0.526$ | $0.467$ | $0.365$ |  | $0.627$ | $0.526$ | $0.467$ | $0.365$
Selected |  |  |  |  |  |  |  |  |
OTP | $0.004$ | $0.016$ | $0.029$ | $0.066$ |  | $0.006$ | $0.021$ | $0.037$ | $0.082$
OPM | $0.004$ | $0.015$ | $0.027$ | $0.063$ |  | $0.005$ | $0.021$ | $0.036$ | $0.081$
ETP | $0.005$ | $0.017$ | $0.031$ | $0.070$ |  | $0.006$ | $0.022$ | $0.039$ | $0.086$
KWTP | $0.004$ | $0.015$ | $0.027$ | $0.064$ |  | $0.005$ | $0.019$ | $0.034$ | $0.080$
EPM | $0.004$ | $0.016$ | $0.029$ | $0.067$ |  | $0.006$ | $0.022$ | $0.038$ | $0.086$
KWPM | $0.003$ | $0.014$ | $0.026$ | $0.062$ |  | $0.005$ | $0.019$ | $0.034$ | $0.080$
LPM | $0.002$ | $0.010$ | $0.022$ | $0.075$ |  | $0.003$ | $0.014$ | $0.029$ | $0.095$
LPM* | $0.010$ | $0.030$ | $0.050$ | $0.100$ |  | $0.010$ | $0.030$ | $0.050$ | $0.100$
EM | $0.002$ | $0.010$ | $0.022$ | $0.075$ |  | $0.003$ | $0.014$ | $0.029$ | $0.095$
EM* | $0.010$ | $0.030$ | $0.050$ | $0.100$ |  | $0.010$ | $0.030$ | $0.050$ | $0.100$
MLE | $0.010$ | $0.030$ | $0.050$ | $0.100$ |  | $0.010$ | $0.030$ | $0.050$ | $0.100$
P-val | $0.010$ | $0.030$ | $0.050$ | $0.100$ |  | $0.010$ | $0.030$ | $0.050$ | $0.100$
Table 7.1. Comparison of Performance of Several Selection Rules for the
Teacher Value Added Simulation.

##  8\. Ranking and Selection of U.S. Dialysis Centers

Motivated by important prior work on ranking and selection by Lin, Louis,
Paddock and Ridgeway (2006, 2009) illustrated by applications to ranking U.S.
dialysis centers, we have chosen to maintain this focus to illustrate our own
approach. Kidney disease is a growing medical problem in the U.S and
considerable effort has been devoted to data collection and evaluation of the
relative performance of the more than 6000 dialysis centers serving the
afflicted population. Centers are evaluated on multiple criteria, but the
primary focus of center ranking is their standardized mortality rate, or SMR,
the ratio of observed deaths to expected deaths for center patients.
Allocating patients to centers is itself a complex task since patients may
move from one center to another in the course of a year. Centers also vary
considerably in the mix of patients they serve. Predictions from an estimated
Cox proportional hazard model that attempts to account for this heterogeneity
are employed to estimate expected deaths for each center.

Our analysis focuses exclusively on the SMR evaluation of centers using
longitudinal data from 2004-18 as reported in DataSource. We restrict
attention to 3230 centers that have consistently reported SMR data over this
sample period. Observed deaths, denoted $y_{it}$ for center $i$ in year $t$
are conventionally modeled as Poisson,

$y_{it}\sim\mbox{Pois}(\rho_{i}\mu_{it})$

where $\mu_{it}$ is center $i$’s expected deaths as predicted by the Cox model
in year $t$ and $\rho_{i}$ is the center’s unobserved mortality rate. We view
$\mu_{it}$ as the effective sample size for the center, after adjustment for
patient characteristics of the center. Center characteristics are expliciitly
excluded from the Cox model. The classical variance stabilizing transformation
for the Poisson brings us back to the Gaussian model,

$z_{it}=\sqrt{y_{it}/\mu_{it}}\sim\mathcal{N}(\theta_{i},1/w_{it}),$

where $\theta_{i}=\sqrt{\rho_{i}}$ and $w_{it}=4\mu_{it}$. Exchangeability of
the centers yields a mixture model in which the parameter $\theta_{i}$, is
effectively assumed to be drawn iidly from a distribution, $G$. The
predictions of expected mortality, $\mu_{it}$, are assumed to be sufficiently
accurate that we treat $w_{it}$ as known, and independent of $\theta_{i}\sim
G$.

Over short time horizons like 3 years we assume that centers have a fixed draw
of $\theta_{i}$ from $G$, and thus we have sufficient statistics for
$\theta_{i}$, as,

$T_{i}=\sum_{t\in\mathcal{T}}w_{it}z_{it}/w_{i}\sim\mathcal{N}(\theta_{i},1/w_{i}),$

where the set $\mathcal{T}$ is the corresponding three year window and
$w_{i}=\sum_{t}w_{it}$. Given these ingredients it is straightforward to
construct a likelihood for the mixing distribution, $G$, and proceed with
estimation of it.

Our objective is then to select centers based on the posterior distributions
of their $\theta_{i}$’s. For example, the posterior tail probability of center
$i$ is given by,

$v_{\alpha}(t_{i},w_{i})=\mathbb{P}(\theta_{i}\geq\theta_{\alpha}|t_{i},w_{i})=\frac{\int_{\theta_{\alpha}}^{+\infty}f(t_{i}|\theta,w_{i})dG(\theta)}{\int_{-\infty}^{+\infty}f(t_{i}|\theta,w_{i})dG(\theta)},$

where $f$ is the density function of $T_{i}$ conditional on $\theta_{i}$ and
$w_{i}$. The capacity constraint requires choosing a thresholding value
$\lambda_{2}^{*}(\alpha)$ such that

$\alpha=\int\int\mymathbb
1\\{v_{\alpha}(t,w)\geq\lambda_{2}^{*}(\alpha)\\}\varphi(t|\theta,w)dG(\theta)dH(w),$

which can be approximated by $\frac{1}{n}\sum_{i}\mymathbb
1\\{v_{\alpha}(t_{i},w_{i})\geq\lambda_{2}^{*}(\alpha)\\}$, and inverted to
obtain the threshold. Based on the discussion in Section 6, for the FDR
constraint we chose a thresholding value $\lambda_{1}^{*}(\alpha,\gamma)$ such
that

(8.1) $\gamma=\frac{\int\int\mymathbb
1\\{v_{\alpha}(t,w)\geq\lambda_{1}^{*}(\alpha,\gamma)\\}(1-v_{\alpha}(t,w))f(t|\theta,w)dG(\theta)dH(w)}{\int\int\mymathbb
1\\{v_{\alpha}(t,w)\geq\lambda_{1}^{*}(\alpha,\gamma)\\}f(t|\theta,w)dG(\theta)dH(w)}$

where $H$ is the marginal distribution of the observed portion of the variance
effect. The numerator can be approximated by
$\frac{1}{n}\sum_{i}(1-v_{\alpha}(t_{i},w_{i}))\mymathbb
1\\{v_{\alpha}(t_{i},w_{i})\geq\lambda_{1}^{*}(\alpha,\gamma)\\}$ and the
denominator can be approximated by, $\frac{1}{n}\sum_{i}\mymathbb
1\\{v_{\alpha}(t_{i},w_{i})\geq\lambda_{1}^{*}(\alpha,\gamma)\\}$.

Posterior mean ranking, in contrast, is based on,

$\delta(t_{i},w_{i})=\mathbb{E}[\theta_{i}|t_{i},w_{i}]=\int\theta
f(t_{i}|\theta,w_{i})dG(\theta).$

For the capacity constraint we choose a thresholding value $C_{2}^{*}(\alpha)$
such that

$\alpha=\int\int\mymathbb 1\\{\delta(t,w)\geq
C_{2}^{*}(\alpha)\\}f(t|\theta,w)dG(\theta)dH(w).$

For FDR constraint, we pick a thresholding value $C_{1}^{*}(\alpha,\gamma)$
such that

$\gamma=\frac{\mathbb{P}(\delta(t,w)\geq
C_{1}^{*}(\alpha,\gamma);\theta<\theta_{\alpha})}{\mathbb{P}(\delta(t,w)\geq
C_{1}^{*}(\alpha,\gamma))}.$

The right hand side of the FDR constraint can be approximated by

$\frac{1}{n}\sum_{i}\mymathbb 1\\{\delta(t_{i},s_{i},w_{i})\geq
C_{1}^{*}(\alpha,\gamma)\\}(1-v_{\alpha}(t_{i},w_{i}))/\frac{1}{n}\sum_{i}\mymathbb
1\\{\delta(t_{i},w_{i})\geq C_{1}^{*}(\alpha,\gamma)\\},$

while the right hand side of the capacity constraint can be approximated by

$\frac{1}{n}\sum_{i}\mymathbb 1\\{\delta(t_{i},w_{i})\geq
C_{2}^{*}(\alpha)\\},$

so $C_{2}^{*}(\alpha)$ is simply the empirical quantile of the
$\delta(t_{i},w_{i})$.

We will compare the foregoing ranking and selection rules with more naive
rules based upon the Poisson and Gaussian MLEs,
$\sum_{t\in\mathcal{T}}y_{it}/\sum_{t\in\mathcal{T}}\mu_{it}$, and $T_{i}$,
respectively, a variant of the much maligned P-value, as well as a linear
shrinkage procedure. For these rules we do not attempt to control for FDR
since this is how they are typically implemented in practice.

To help appreciate the difficulty of the selection task, Table 8.1 reports
estimated FDR rates for several selection rules under a range of capacity
constraints $\alpha$ for both right and left tail selection based on the data
from 2004 - 2006. Right tail selection corresponds to identifying centers
whose mortality rate is higher than expected; left tail selection to centers
with mortality lower than expected. To estimate FDR we require an estimate of
the distribution of distribution, $G$. For this purpose we use the smoothed
version of the Kiefer-Wolfowitz NPMLE introduced in Section 2. The biweight
bandwidth for the smoothing was chosen as the mean absolute deviation from the
median of the discrete NPMLE, $\hat{G}$. The assessment of FDR reported in
Table 8.1 reflects the considerable uncertainty associated with the selected
set of centers deemed by the capacity constraint to be in the upper (or lower)
$\alpha$ quantile based upon our estimate of the distribution, $G$, of
unobserved quality.

The MLE rule ranks centers based on their Gaussian MLE, $T_{i}$, while the
Poisson-MLE rule ranks on $\sum_{t}y_{it}/\sum_{t}\mu_{it}$, which is the MLE
of $\rho_{i}$ from the Poisson model. P-value ranks centers based on the
variance stabilizing transformation from the Poisson model under the null
hypothesis $\rho_{i}=1$ and $\rho_{i}>1$ as the alternative hypothesis for
right selection and $\rho_{i}<1$ for the left selection. All these rules
ignore the compound decision perspective of the problem entirely.

Among the compound decision rules, we consider the linear (James-Stein)
shrinkage rule,
$\hat{\mu}_{\theta}+(T_{i}-\hat{\mu}_{\theta})\hat{\sigma}_{\theta}^{2}/(\hat{\sigma}_{\theta}^{2}+1/w_{i})$
which is the posterior mean of $\theta_{i}$ based on the model
$T_{i}\sim\mathcal{N}(\theta_{i},1/w_{i})$ assuming that the latent variable
$\theta_{i}$ follows a Gaussian distribution with mean $\mu_{\theta}$ and
variance $\sigma_{\theta}^{2}$. We also consider the EfronMorris estimator
which is a slight modification of the James-Stein estimator.

Finally, PM and TP are the posterior mean of $\theta$ and posterior tail
probability of $\theta\geq\theta_{\alpha}$, for right selection, and
$\theta\leq\theta_{\alpha}$ for left selection based on our estimated
$\hat{G}$. For both left and right tail selection, as $\alpha$ increases, the
FDR rate decreases, indicating the selection task becomes easier. All rules
that account for the compound decision perspective of the problem have
slightly lower FDRs than those that consider each center individually.

| $\alpha=4\%$ | $\alpha=10\%$ | $\alpha=15\%$ | $\alpha=20\%$ | $\alpha=25\%$
---|---|---|---|---|---
Right Selection |  |  |  |  |
MLE | $0.564$ | $0.481$ | $0.432$ | $0.399$ | $0.353$
Poisson-MLE | $0.564$ | $0.485$ | $0.436$ | $0.402$ | $0.356$
Pvalue | $0.559$ | $0.479$ | $0.430$ | $0.396$ | $0.351$
Efron-Morris | $0.547$ | $0.476$ | $0.427$ | $0.395$ | $0.350$
James-Stein | $0.547$ | $0.476$ | $0.427$ | $0.395$ | $0.350$
PM | $0.543$ | $0.474$ | $0.426$ | $0.395$ | $0.350$
TP | $0.542$ | $0.473$ | $0.425$ | $0.394$ | $0.350$
Left Selection |  |  |  |  |
MLE | $0.643$ | $0.559$ | $0.479$ | $0.437$ | $0.385$
Poisson-MLE | $0.634$ | $0.556$ | $0.478$ | $0.436$ | $0.384$
Pvalue | $0.660$ | $0.566$ | $0.482$ | $0.440$ | $0.386$
Efron-Morris | $0.635$ | $0.550$ | $0.475$ | $0.434$ | $0.384$
James-Stein | $0.635$ | $0.550$ | $0.475$ | $0.434$ | $0.384$
PM | $0.632$ | $0.549$ | $0.474$ | $0.434$ | $0.383$
TP | $0.627$ | $0.547$ | $0.472$ | $0.432$ | $0.382$
Table 8.1. FDR rate: 2004-2006

The Stars assigns ratings of five stars down to one star to centers in the
proportions $\\{0.22,0.30,0.35,0.09,0.04\\}$ respectively. We will abbreviate
these ratings to the conventional academic scale of A-F. To illustrate the
conflict between the selection criteria we plot in Figure 8.1 the centers
selected for the grade A (Five star, which consists 22% of the centers that
suppose to have their true mortality rate being the lowest) category with and
without FDR control. Centers are characterized by pairs, $(T_{i},w_{i})$,
consisting of their weighted mean standarized mortality, $T_{i}$, and their
estimate of the precision, $w_{i}$, of these mortality estimates. In each plot
the solid curves represent the decision boundaries of the selection rule under
comparison. Centers with low mortality and relatively high precision appear
toward the northwest in each figure.

Panel (a) of the figure compares the posterior tail probability selection with
the MLE, or fixed effect, selection. The selection boundary for the MLE is the
(red) vertical line, since the MLE ignores the precision of the estimates
entirely. The selection boundary for the tail probability rule is indicated by
the (blue) curve. A few centers with high precision excluded by the MLE rule
are selected by the TP rule, and on the contrary a few centers with low
precision are selected by the MLE rule, but excluded by the TP criterion.
Panel (b) imposes FDR control with $\gamma=0.20$ on the TP selection with an
estimated thresholding value implied by the FDR constraint using the smoothed
NPMLE. The MLE selection is the same as in Panel (a) without the FDR control.
We see that under TP rule with FDR control, the number of selected centers is
reduced considerably. Instead of selecting 711 centers allowed by the
capacity, it selects only 230 centers. In comparison, the MLE rule under
capacity constraint has an estimated FDR rate at 0.431. Panel (c) compares
selected centers by the TP rule with those selected by a James-Stein linear
shrinkage rule. Now the TP rule tolerates a few more low precision centers,
while it is the James-Stein rule that demands higher precision to be selected.
Finally, in Panel (d) we again subject the TP rule to FDR control of 20
percent, while the James-Stein rule continues to adhere only to the capacity
constraint. The TP boundary scales back substantially, suggesting that a large
proportion of the extra selections made by James-Stein linear shrinkage rules
are likely to be false discoveries. In fact, the estimated FDR rate of the
James-Stein rule under just capacity constraint is also 0.431, same as that of
the MLE rule.

Figure 8.1. Contrasting Selections for A-rated centers: The two upper panels
compare posterior tail probability selection with MLE (fixed effects)
selection, while the lower panels compare TP selection with James-Stein
(linear shrinkage) selection. Left panels impose capacity control only, while
the right panels impose 20 percent FDR control for the TP rule. The estimated
FDR rate for both the MLE and James-Stein selection under capacity constraint,
using the smoothed NPMLE estimator for $G$, is 0.431. Comparisons are based on
the 2004-2006 data.

Given the longtitudinal structure of the Dialysis data, it would be possible
to consider the models in Section 5 that allow for unobserved variance
heterogeneity. We refrain from doing so partly due to space considerations and
because we are reluctant to assume stationarity of random effects over longer
time horizons.

###  8.1. Temporal Stability, Ranking and Selection

Given the longitudinal nature of the data, it is natural to ask, “How stable
are rankings over time, and isn’t there some temporal dependence in the
observed data that should be accounted for?” Surprisingly, the year-to-year
dependence in the observed mortality is quite weak. In Figure 8.2 we plot a
histogram of estimated AR(1) coefficients for the 3230 centers; it is roughly
centered at zero and slightly skewed to the left. We do not draw the
conclusion from this that there is no temporal dependence in the observed
$y_{it}$, but only that there is considerable heterogeneity in the nature of
this dependence with roughly as many centers exhibiting negative serial
dependence as those with positive dependence. Our approach of considering
brief, 3-5 year, windows of presumed stability in center performance is
consistent with the procedures of the official ranking agency. In each of
these windows we can compute a ranking according to one of the criteria
introduced above, and it is of interest to know how much stability do we see
in these rankings?

Figure 8.2. Histogram of estimated AR(1) coefficients for 3230 Dialysis
centers based on annual data 2004-2017

To address this question we consider rankings based on the posterior tail
probability criterion for three year windows. In each of the 5 3-year windows
we assign centers letter grades, A-F, with proportions
$\\{0.22,0.30,0.35,0.09,0.04\\}$ respectively. Table 8.2 reports the estimated
transition matrix between these categories, so entry $i,j$ in the matrix
represents the estimated probability of a center in state $i$ moving to state
$j$ in the next period.

| A | B | C | D | F
---|---|---|---|---|---
A | $0.440$ | $0.330$ | $0.200$ | $0.024$ | $0.006$
B | $0.248$ | $0.357$ | $0.328$ | $0.059$ | $0.007$
C | $0.122$ | $0.286$ | $0.440$ | $0.113$ | $0.039$
D | $0.060$ | $0.188$ | $0.436$ | $0.208$ | $0.108$
F | $0.021$ | $0.081$ | $0.352$ | $0.217$ | $0.329$
Table 8.2. Estimated First Order Markov Transition Matrix: Entry $i,j$ of the
matrix estimates the probability of a transition from state $i$ to state $j$
based on posterior tail probability rankings for 3-year longitudinal grouping
of the center data.

It is obviously difficult to maintain an “A” rating for more than a couple of
periods, but centers with poor performance are also likely to move into the
middle of the rankings. Although, as we have seen, there is no guarantee that
the posterior tail probability criterion yields a nested ranking, nestedness
does hold in this particular application. Posterior mean ranking yields
similar transition behavior. The high degree of mobility between rating
categories reenforces our conclusion that ranking and selection into rating
categories is subject to considerable uncertainty.

##  9\. Conclusions

Robbins’s compound decision framework is well suited to ranking and selection
problems, and nonparametric maximum likelihood estimation of mixture models
offers a powerful tool for implementing empirical Bayes rules for such
problems. Posterior tail probability selection rules perform better than
posterior mean rules when precision is heterogeneous. Ranking and selection is
especially difficult in Gaussian settings where classical linear shrinkage
methods are most appropriate. Nonparametric empirical Bayes methods can
substantially improve upon selection methods based on linear shrinkage and
traditional $p$-values when the latent mixing distribution is not Gaussian
both in terms of power and false discovery rate.

##  Appendix A Proofs

Proof. [Lemma 3.1] We can write:

$\displaystyle\nabla_{y}\frac{\alpha f_{1}(y)}{f(y)}$
$\displaystyle=\frac{\int_{[\theta_{\alpha},+\infty)}\nabla_{y}\log\varphi(y|\theta,\sigma^{2})\varphi(y|\theta,\sigma^{2})dG(\theta)}{\int_{(-\infty,+\infty)}\varphi(y|\theta,\sigma^{2})dG(\theta)}$
$\displaystyle-\frac{\int_{[\theta_{\alpha},+\infty)}\varphi(y|\theta,\sigma^{2})dG(\theta)}{\int_{(-\infty,+\infty))}\varphi(y|\theta,\sigma^{2})dG(\theta)}\frac{\int_{(-\infty,+\infty)}\nabla_{y}\log\varphi(y|\theta,\sigma^{2})\varphi(y|\theta,\sigma^{2})dG(\theta)}{\int_{(-\infty,+\infty)}\varphi(y|\theta,\sigma^{2})dG(\theta)}$
$\displaystyle=\mathbb{E}\Big{[}\mymathbb
1\\{\theta\geq\theta_{\alpha}\\}\nabla_{y}\log\varphi(y|\theta,\sigma^{2})|Y\Big{]}-\mathbb{E}\Big{[}\mymathbb
1\\{\theta\geq\theta_{\alpha}\\}|Y]\mathbb{E}[\nabla_{y}\log\varphi(y|\theta,\sigma^{2})|Y\Big{]}$
$\displaystyle=\text{Cov}\Big{[}\mymathbb
1\\{\theta\geq\theta_{\alpha}\\},\nabla_{y}\log\varphi(y|\theta,\sigma^{2})|Y\Big{]}\geq
0$

The last inequality holds since $\nabla_{y}\log\varphi(y|\theta,\sigma^{2})$
is increasing in $\theta$ for each fixed $\sigma^{2}$, by the Gaussian
assumption and the fact that covariance of monotone functions of $\theta$ is
non-negative (see schmidt2014inequalities) assuming the existence of
$\mathbb{E}_{\theta|Y}[\nabla_{y}\log\varphi(Y|\theta,\sigma^{2})|Y]$ which we
assume. Nesting follows from the monotonicity of the $v_{\alpha}(y)$
criterion: monotonicity of $v_{\alpha}(y)$ implies that there exists
$t_{\alpha}$ such that $\mymathbb
1\\{v_{\alpha}(y)\geq\lambda_{\alpha}/(1+\lambda_{\alpha})\\}=\mymathbb
1\\{y\geq t_{\alpha}\\}$, hence if $\alpha_{1}>\alpha_{2}$ and
$\mathbb{P}(y\geq t_{\alpha_{1}})=\alpha_{1}$ and $\mathbb{P}(y\geq
t_{\alpha_{2}})=\alpha_{2}$, then it must be that $t_{\alpha_{1}}\leq
t_{\alpha_{2}}$, implying nestedness.

Proof. [Lemma 3.2] Denote the three decision criteria,
$v_{1}(y)=\mathbb{E}(\theta|Y=y)$, $v_{2}(y)=\mathbb{P}(\theta\geq
G^{-1}(1-\alpha)|Y=y)$ and $v_{3}(y)=\mathbb{E}(\theta 1(\theta\geq
G^{-1}(1-\alpha))|Y=y)$. Assuming that $\mathbb{E}[\theta|Y]<\infty$,
$\mathbb{E}_{\theta|Y}[\nabla_{y}\log\varphi(y|\theta,\sigma^{2})|Y]<\infty$
and
$\mathbb{E}_{\theta|Y}[\theta\nabla_{y}\log\varphi(y|\theta,\sigma^{2})|Y]<\infty$,
the calculation leading to the proof of Lemma 1 shows:

$\displaystyle\nabla_{y}v_{1}(y)$
$\displaystyle=\text{Cov}(\theta,\nabla_{y}\log\varphi(y|\theta)|Y=y)$
$\displaystyle\nabla_{y}v_{2}(y)$ $\displaystyle=\text{Cov}(1(\theta\geq
G^{-1}(1-\alpha)),\nabla_{y}\log\varphi(y|\theta)|Y=y)$
$\displaystyle\nabla_{y}v_{3}(y)$ $\displaystyle=\text{Cov}(\theta
1(\theta\geq G^{-1}(1-\alpha)),\nabla_{y}\log\varphi(y|\theta)|Y=y).$

Thus, the monotonicity of $\nabla_{y}\log\varphi(y|\theta,\sigma^{2})$ implies
they all yield identical rankings.

Proof. [Proposition 3.4] The Bayes rule for the non-randomized selections can
be characterized as,

$\delta_{i}^{*}=\begin{cases}1,&\text{if
}v_{\alpha}(y_{i})\geq\tau_{1}^{*}(1-v_{\alpha}(y_{i})-\gamma)+\tau_{2}^{*}\\\
0,&\text{if
}v_{\alpha}(y_{i})<\tau_{1}^{*}(1-v_{\alpha}(y_{i})-\gamma)+\tau_{2}^{*}\end{cases}$

with Karush-Kuhn-Tucker conditions,

(A.1)
$\displaystyle\tau_{1}^{*}\Big{(}\mathbb{E}\Big{[}\sum_{i=1}^{n}\Big{\\{}(1-v_{\alpha}(y_{i}))\delta_{i}^{*}-\gamma\delta_{i}^{*}\Big{\\}}\Big{]}\Big{)}$
$\displaystyle=0$ (A.2)
$\displaystyle\tau_{2}^{*}\Big{(}\mathbb{E}\Big{[}\sum_{i=1}^{n}\delta_{i}^{*}\Big{]}-\alpha
n\Big{)}$ $\displaystyle=0$ (A.3)
$\displaystyle\mathbb{E}\Big{[}\sum_{i=1}^{n}\Big{\\{}(1-v_{\alpha}(y_{i}))\delta_{i}^{*}-\gamma\delta_{i}^{*}\Big{\\}}\Big{]}$
$\displaystyle\leq 0$ (A.4) $\displaystyle\frac{\mymathbb
1}{n}\mathbb{E}\Big{[}\sum_{i=1}^{n}\delta_{i}^{*}\Big{]}-\alpha$
$\displaystyle\leq 0$ (A.5) $\displaystyle\tau_{1}^{*}$ $\displaystyle\geq 0$
(A.6) $\displaystyle\tau_{2}^{*}$ $\displaystyle\geq 0$

The Bayes rule takes the form of thresholding on the posterior tail
probability $v_{\alpha}(y)$ and since $v_{\alpha}(y)$ is monotone in $y$ as
shown in Lemma 3.2, it is therefore a thresholding rule on $Y$,
$\delta_{i}^{*}=\mymathbb 1\\{y_{i}\geq t^{*}\\}$ with cutoff $t^{*}$
depending on the values of $(\tau_{1}^{*},\tau_{2}^{*},\alpha,\gamma)$.
Condition (A.3) is equivalent to the condition that the marginal false
discovery rate, mFDR, since it requires

$\mathbb{E}\Big{[}\sum_{i=1}^{n}\\{(1-v_{\alpha}(y_{i})\delta_{i}^{*}\\}\Big{]}/\mathbb{E}\Big{[}\sum_{i=1}^{n}\delta_{i}^{*}\Big{]}\leq\gamma$

and we can show that the left hand side quantity is precisely the mFDR since

$\displaystyle\text{mFDR}(t^{*})=$
$\displaystyle\mathbb{P}(\delta_{i}^{*}=1,\theta_{i}\leq\theta_{\alpha})/\mathbb{P}(\delta_{i}^{*}=1)$
$\displaystyle=$ $\displaystyle(1-\alpha)\int\mymathbb 1\\{y\geq
t^{*}\\}f_{0}(y)dy/\int\mymathbb 1\\{y\geq t^{*}\\}f(y)dy$ $\displaystyle=$
$\displaystyle\int\mymathbb 1\\{y\geq
t^{*}\\}(1-v_{\alpha}(y))f(y)dy/\int\mymathbb 1\\{y\geq t^{*}\\}f(y)dy$
$\displaystyle=$
$\displaystyle\mathbb{E}\Big{[}\sum_{i=1}^{n}\\{(1-v_{\alpha}(y_{i}))\delta_{i}^{*}\\}\Big{]}/\mathbb{E}\Big{[}\sum_{i=1}^{n}\delta_{i}^{*}\Big{]}$
$\displaystyle=$
$\displaystyle\frac{\int_{-\infty}^{\theta_{\alpha}}\tilde{\Phi}((t^{*}-\theta)/\sigma)dG(\theta)}{\int_{-\infty}^{+\infty}\tilde{\Phi}((t^{*}-\theta)/\sigma)dG(\theta)}\leq\gamma.$

For any mixing distribution $G$, as $t^{*}$ increases, it becomes less likely
for condition (A.3) to bind. And as $t^{*}$ approaches $-\infty$, left side of
(A.3) approaches $1-\alpha-\gamma$, and hence we’ve restricted
$\gamma<1-\alpha$ to avoid cases where the condition (A.3) never binds. On the
other hand, condition (A.4) is equivalent to

$\mathbb{P}(\delta_{i}^{*}=1)-\alpha=\int_{-\infty}^{+\infty}\tilde{\Phi}((t^{*}-\theta)/\sigma)dG(\theta)-\alpha\leq
0.$

As $t^{*}$ increases, it also becomes less likely that condition (A.4) binds.
Therefore, we can define,

$\displaystyle t_{1}^{*}$ $\displaystyle=\min\\{t:\text{mFDR}(t)-\gamma\leq
0\\}$ $\displaystyle t_{2}^{*}$
$\displaystyle=\min\\{t:\int_{-\infty}^{+\infty}\tilde{\Phi}((t-\theta)/\sigma)dG(\theta)-\alpha\leq
0\\}$

When $t_{1}^{*}<t_{2}^{*}$, the feasible region for $Y$ defined by inequality
(A.4) is a strict subset of that defined by inequality (A.3). When
$t_{1}^{*}>t_{2}^{*}$, then the feasible region defined by inequality (A.3) is
a strict subset of that defined by inequality (A.4). When
$t_{1}^{*}=t_{2}^{*}$, the feasible regions coincide. This case occurs when
$\text{mFDR}(t_{1}^{*})=\gamma$ and $\mathbb{P}(y\geq t_{2}^{*})=\alpha$, so
$\mathbb{E}[v_{\alpha}(Y)\mymathbb
1\\{v_{\alpha}(Y)\geq\lambda^{*}\\}]=\alpha-\alpha\gamma$ with
$v_{\alpha}(t_{1}^{*})=\lambda^{*}(\alpha,\gamma)$. Again, the strict
thresholding enforced by the statement of the proposition can be relaxed
slightly by randomizing the selection probability of the last unit so that the
active constraint is satisfied exactly.

Proof. [Proposition 4.1] In the proof we will suppress the dependence of
$\lambda^{*}$ on the $(\alpha,\gamma)$. The argument is very similar to the
proof for Proposition 3.4, except that now the feasible region defined by
constraint (A.3) and (A.4) is a two-dimensional region for
$(y_{i},\sigma_{i})$. Since the posterior tail probability
$v_{\alpha}(y,\sigma)$ is monotone in $y$ for any fixed $\sigma$ as a result
of Lemma 3.2, the optimal rule can be reformulated as a thresholding rule on
$Y$ again, $\delta_{i}^{*}=\mymathbb
1\\{y_{i}>t_{\alpha}(\lambda^{*},\sigma_{i})\\}$ except now the threshold
value also depends on $\sigma_{i}$.

Now consider the constraint (A.3) and (A.4). Condition (A.3) is equivalent to
the condition that,

$\frac{\mathbb{E}\Big{[}\sum_{i=1}^{n}\\{(1-v_{\alpha}(y_{i},\sigma_{i}))\delta_{i}^{*}\\}\Big{]}}{\mathbb{E}\Big{[}\sum_{i=1}^{n}\delta_{i}^{*}\Big{]}}-\gamma=\frac{\int\int_{-\infty}^{\theta_{\alpha}}\tilde{\Phi}((t_{\alpha}(\lambda^{*},\sigma)-\theta)/\sigma)dG(\theta)dH(\sigma)}{\int\int_{-\infty}^{+\infty}\tilde{\Phi}((t_{\alpha}(\lambda^{*},\sigma)-\theta)/\sigma)dG(\theta)dH(\sigma)}-\gamma\leq
0.$