Consistent estimation of the missing mass for feature models

TL;DR

This paper investigates the challenge of estimating the unseen features in feature models, proving the non-existence of a universal consistent estimator and demonstrating the consistency of a specific estimator under heavy-tailed probabilities.

Contribution

It establishes the impossibility of universal consistent estimation of the missing mass and shows the consistency of a particular estimator for heavy-tailed feature probabilities.

Findings

No universally consistent estimator exists for the missing mass.

The estimator by Ayed et al. (2017) is strongly consistent for heavy-tailed probabilities.

Derived concentration inequalities for missing mass and feature frequency counts.

Abstract

Feature models are popular in machine learning and they have been recently used to solve many unsupervised learning problems. In these models every observation is endowed with a finite set of features, usually selected from an infinite collection $(F_{j})_{j\geq 1}$. Every observation can display feature $F_{j}$ with an unknown probability $p_{j}$. A statistical problem inherent to these models is how to estimate, given an initial sample, the conditional expected number of hitherto unseen features that will be displayed in a future observation. This problem is usually referred to as the missing mass problem. In this work we prove that, using a suitable multiplicative loss function and without imposing any assumptions on the parameters $p_{j}$, there does not exist any universally consistent estimator for the missing mass. In the second part of the paper, we focus on a special class of heavy-tailed probabilities $(p_{j})_{j\geq 1}$, which are common in many real applications, and we show that, within this restricted class of probabilities, the nonparametric estimator of the missing mass suggested by Ayed et al. (2017) is strongly consistent. As a byproduct result, we will derive concentration inequalities for the missing mass and the number of features observed with a specified frequency in a sample of size $n$.

Equations

$$M_n(Y_1,\dots ,Y_n;(p_j{)}_{j\ge 1})=\mathbb{E}\left[\sum _{j\ge 1}11_{\{X_{n,j}=0,Y_{n+1,j}=1\}}|Y_1,\dots ,Y_n\right]=\sum _{j\ge 1}{p}_{j}11_{\{X_{n,j}=0\}},$$

$$L({\hat{T}}_n,M_n):=\left|\frac{{\hat{T}}_n}{M_n}-1\right|.$$

$$\mathcal{P}:=\left\{\sum _{j\ge 1}{p}_j{\delta }_{F_j}:F_j,p_j\in [0,1],\forall j\ge 1,\sum _{j\ge 1}{p}_j<+\infty \right\}.$$

$$\underset{n\to +\infty }{\lim }{\mathbb{P}}_{{\mathbf{Y}}_n|p}\left(\left|\frac{{\hat{T}}_n}{M_n}-1\right|\ge \epsilon \right)=0,$$

$$\underset{n\to +\infty }{\mathrm{limsup}}{\mathbb{P}}_{{\mathbf{Y}}_n|p}\left(\left|\frac{{\hat{T}}_n}{M_n}-1\right|\ge \epsilon \right)>C.$$

$$\underset{p\in \mathcal{P}}{\sup }\underset{n\to +\infty }{\mathrm{limsup}}{\mathbb{P}}_{{\mathbf{Y}}_n|p}\left(\left|\frac{{\hat{T}}_n}{M_n}-1\right|\ge ϵ\right)>C,$$

$$\underset{p\in \mathcal{P}}{\sup }\underset{n\to +\infty }{\mathrm{limsup}}{\mathbb{P}}_{{\mathbf{Y}}_n|p}\left(\left|\frac{{\hat{T}}_n}{M_n}-1\right|\ge ϵ\right)\ge \underset{p\in \mathcal{P}}{\sup }\underset{n\to +\infty }{\mathrm{limsup}}{\mathbb{P}}_{{\mathbf{Y}}_n|p}\left(\left|\frac{M_n}{{\hat{T}}_n}-1\right|\ge 2ϵ\right).$$

$$-M_nϵ<{\hat{T}}_n-M_n<M_nϵ,$$

$$\left|\frac{M_n}{{\hat{T}}_n}-1\right|\ge 2ϵ\Rightarrow \left|\frac{{\hat{T}}_n}{M_n}-1\right|\ge ϵ,$$

$$\underset{p\in \mathcal{P}}{\sup }\underset{n\to +\infty }{\mathrm{limsup}}{\mathbb{P}}_{{\mathbf{Y}}_n|p}\left(\left|\frac{M_n}{{\hat{T}}_n}-1\right|\ge 2ϵ\right)>C,$$

$$\underset{p\in \mathcal{P}}{\sup }\underset{n\to +\infty }{\mathrm{limsup}}{\mathbb{P}}_{{\mathbf{Y}}_n|p}\left(\left|\frac{M_n}{{\hat{T}}_n}-1\right|\ge \epsilon \right)\ge {\mathbb{E}}_p\left[\underset{n\to +\infty }{\mathrm{limsup}}{\mathbb{P}}_{{\mathbf{Y}}_n|p}\left(\left|\frac{M_n}{{\hat{T}}_n}-1\right|\ge \epsilon \right)\right]$$

$$\ge \underset{n\to +\infty }{\mathrm{limsup}}{\mathbb{E}}_{{\mathbf{Y}}_n}\left[{\mathbb{P}}_{p|{\mathbf{Y}}_n}\left(\left|\frac{M_n}{{\hat{T}}_n}-1\right|\ge \epsilon \right)\right].$$

$$\mathbb{E}\left[e^{-{\int }_{[0,1]}f(F)\tilde{N}(\mathrm{d}F)}\right]=\exp \left\{-{\int }_{\text{mathds}{R}^{+}\times [0,1]}(1-e^{-sf(F)})\nu (\mathrm{d}s,\mathrm{d}F)\right\},$$

$$\tilde{N}|{\mathbf{Y}}_n\stackrel{d}{=}{\tilde{N}}_n+\sum _{\ell =1}^{k_n}{J}_{\ell }{\delta }_{F_{\ell }^{\ast }}$$

$$\begin{array}{rl}{\mathbb{P}}_{p|{\mathbf{Y}}_n}(M_n\in B)& ={\mathbb{P}}_{p|{\mathbf{Y}}_n}\left(\sum _{j\ge 1}{p}_j{\delta }_{F_j}(A_n^c)\in B\right)\\ & ={\mathbb{P}}_{\tilde{N}|{\mathbf{Y}}_n}\left(\sum _{j\ge 1}(1-e^{-s_j}){\delta }_{F_j}(A_n^c)\in B\right)\\ & ={\mathbb{P}}_{\tilde{N}|{\mathbf{Y}}_n}\left({\int }_{A_n^c}(1-e^{-s})\tilde{N}(\mathrm{d}F)\in B\right)\\ & \stackrel{\text{eqref}eq:posterior}{=}{\mathbb{P}}_{{\tilde{N}}_n}\left({\int }_{[0,1]}(1-e^{-s}){\tilde{N}}_n(\mathrm{d}F)\in B\right)\end{array}$$

$$S_n:={\tilde{N}}_n([0,1])={\int }_{[0,1]}s{\tilde{N}}_n(\mathrm{d}F)$$

$$\mathbb{E}[e^{x{S}_n}]$$

$$=\exp \left\{-{\int }_0^1{\int }_0^{+\infty }(1-e^{xs})e^{-sn}\nu (\mathrm{d}s,\mathrm{d}\theta )\right\}$$

$$=\exp \left\{-{\int }_0^{+\infty }(1-e^{xs})\frac{e^{-s(n+1)}}{s}\mathrm{d}s\right\}={\left(1-\frac{x}{n+1}\right)}^{-1},$$

$$\begin{array}{rl}\left|\frac{M_n}{{\hat{T}}_n}-1\right|& =\left|\frac{M_n}{S_n}\left(\frac{S_n}{{\hat{T}}_n}-1+1\right)-1\right|\ge \left|\frac{M_n}{S_n}\left|\frac{S_n}{{\hat{T}}_n}-1\right|-\left|\frac{M_n}{S_n}-1\right|\right|\\ & \ge \frac{M_n}{S_n}\left|\frac{S_n}{{\hat{T}}_n}-1\right|-\left|\frac{M_n}{S_n}-1\right|\end{array}$$

$$\mathbb{P}_{p|\mathbf{Y}_{n}}\left(1-\frac{\varepsilon}{2}\leq\frac{M_{n}}{S_{n}}\leq 1\,,\,\Big{|}\frac{S_{n}}{\hat{T}_{n}}-1\Big{|}>3\varepsilon\right)\leq\mathbb{P}_{p|\mathbf{Y}_{n}}\left(\Big{|}\frac{M_{n}}{\hat{T}_{n}}-1\Big{|}>\varepsilon\right)$$

$$1-\frac{\varepsilon}{2}\leq\frac{M_{n}}{S_{n}}\leq 1\,,\,\Big{|}\frac{S_{n}}{\hat{T}_{n}}-1\Big{|}>3\varepsilon$$

$$\left|\frac{M_n}{{\hat{T}}_n}-1\right|\ge \frac{M_n}{S_n}\left|\frac{S_n}{{\hat{T}}_n}-1\right|-\left|\frac{M_n}{S_n}-1\right|\ge \left(1-\frac{\epsilon }{2}\right)3\epsilon -\frac{\epsilon }{2}=\epsilon (5-3\epsilon )/2>\epsilon$$

$${\mathbb{P}}_{p|{\mathbf{Y}}_n}\left(\left|\frac{M_n}{{\hat{T}}_n}-1\right|>\epsilon \right)\ge {\mathbb{P}}_{p|{\mathbf{Y}}_n}\left(1-\frac{\epsilon }{2}\le \frac{M_n}{S_n}\le 1\right)-1+{\mathbb{P}}_{p|{\mathbf{Y}}_n}\left(\left|\frac{S_n}{{\hat{T}}_n}-1\right|>3\epsilon \right)$$

$$\begin{array}{rl}& \underset{p\in \mathcal{P}}{\sup }\underset{n\to +\infty }{\mathrm{limsup}}{\mathbb{P}}_{{\mathbf{Y}}_n|p}\left(\left|\frac{{\hat{T}}_n}{M_n}-1\right|\ge \epsilon \right)\\ & \ge \underset{n\to +\infty }{\mathrm{limsup}}{\mathbb{E}}_{{\mathbf{Y}}_n}\left[{\mathbb{P}}_{p|{\mathbf{Y}}_n}\left(1-\frac{\epsilon }{2}\le \frac{M_n}{S_n}\le 1\right)-1\right]\\ & +\underset{x>0}{\inf }\underset{{\mathbf{Y}}_n}{\inf }{\mathbb{P}}_{p|{\mathbf{Y}}_n}\left(\left|\frac{S_n}{x}-1\right|>3\epsilon \right).\end{array}$$

$$s_j-\frac{1}{2}{s}_j^2\le 1-e^{-s_j}\le s_j.$$

$$S_n-\frac{1}{2}{S}_n^2\le S_n-\frac{1}{2}\sum _{j\ge 1}{s}_j^2{\delta }_{F_j}(A_n^c)\le M_n\le S_n,$$

$$1-\frac{1}{2}{S}_n\le \frac{M_n}{S_n}\le 1.$$

$$\begin{array}{rl}& {\mathbb{P}}_{p|{\mathbf{Y}}_n}\left(1-\frac{\epsilon }{2}\le \frac{M_n}{S_n}\le 1\right)-1\ge {\mathbb{P}}_{p|{\mathbf{Y}}_n}\left(S_n\le \epsilon \right)-1\\ & =(n+1){\int }_0^{\epsilon }{e}^{-x(n+1)}\mathrm{d}x-1=-e^{-\epsilon (n+1)},\end{array}$$

$${\mathbb{P}}_{p|{\mathbf{Y}}_n}\left(\left|\frac{S_n}{x}-1\right|>3\epsilon \right)$$

$$=1+e^{-(1+3\epsilon )x(n+1)}-e^{-(1-3\epsilon )x(n+1)}$$

$$\ge \underset{y>0}{\inf }\left[1+e^{-(1+3\epsilon )y}-e^{-(1-3\epsilon )y}\right]$$

$${\mathbb{P}}_{p|{\mathbf{Y}}_n}\left(\left|\frac{S_n}{x}-1\right|>3\epsilon \right)\ge f(y^{\ast })=:C>0.$$

Full text

Keywords: Feature models; missing mass; multiplicative consistency; regular variation; nonparametric estimator.

1 Introduction

Feature models generalize species sampling models by allowing every observation to belong to more than one species, now called features. In particular, every observation is endowed with a finite set of features selected from a (possibly infinite) collection of features $(F_{j})_{j\geq 1}$. Every feature $F_{j}$ is associated with an unknown probability $p_{j}$, and each observation displays feature $F_{j}$ with probability $p_{j}$. We may conveniently represent each observation with a binary sequence, whose entries indicate the presence (1) or absence (0) of each feature. Feature models have been first applied in ecology for modeling incidence vectors collecting the presence or absence of species traps (Colwell et al. (2012) and Chao et al. (2014)), and more recently in several fields of biosciences, such as the study of genetic variation and protein interactions (Chu et al. (2006), Ionita-Laza et al. (2009), Ionita-Laza et al. (2010) and Zou et al. (2016)). They also found applications in the analysis of choice behaviour arising from psychology, marketing and computer science (Görür et al. (2006)); in the context of binary matrix factorization for modeling dyadic data to design recommender system (Meeds et al. (2007)); in graphical models (Wood et al. (2006) and Wood & Griffiths (2007)); in cognitive psychology for the analysis of similarity judgement matrices (Navarro & Griffiths (2007)); in the context of independent component analysis and sparse factor analysis (Knowles & Ghahramani (2007)); in link prediction using network data (Miller et al. (2010)).

The Bernoulli product model is arguably the most popular feature model. It assumes that the $i$–th observation is a sequence $Y_{i}=(Y_{i,j})_{j\geq 1}$ of independent Bernoulli random variables with unknown success probabilities $(p_{j})_{j\geq 1}$, and that $Y_{r}$ is independent of $Y_{s}$ for any $r\neq s$. Therefore $X_{n,j}:=\sum_{1\leq i\leq n}Y_{i,j}$, namely the number of times that feature $F_{j}$ has been observed in a sample $(Y_{1},\ldots,Y_{n})$, is a Binomial random variable with parameter $(n,p_{j})$ for any $j\geq 1$. Recently, the Bernoulli product model has been extensively applied to the fundamental problem of discovering genetic variation in human populations. See, e.g., Ionita-Laza et al. (2009), Zou et al. (2016)) and references therein. In such a context, interest is in estimating the conditional expected number, given a sample $(Y_{1},\ldots,Y_{n})$, of hitherto unseen features that would be observed if an additional sample $Y_{n+1}$ was collected, namely

[TABLE]

where $1\!\!1$ is the indicator function. The statistic $M_{n}(Y_{1},\ldots,Y_{n};(p_{j})_{j\geq 1})$ is referred to as the missing mass, i.e. the sum of the probability masses of unobserved features in a sample of size $n$. In genetics, interest in estimating (1) is motivated by the ambitious prospect of growing databases to encompass hundreds of thousands of genomes, which makes important to quantify the power of large sequencing projects to discover new genetic variants (Auton et al. (2015)). An accurate estimate of the missing mass provides a quantitative evaluation of the potential and limitations of these datasets, providing a roadmap for large-scale sequencing projects.

Let $\hat{T}_{n}(Y_{1},\ldots,Y_{n})$ denote an arbitrary estimator of $M_{n}(Y_{1},\ldots,Y_{n};(p_{j})_{j\geq 1})$. For easiness of notation, in the rest of the paper we will not highlight the dependence on $(Y_{1},\ldots,Y_{n})$ and $(p_{j})_{j\geq 1}$, and we simply write $M_{n}$ and $\hat{T}_{n}$. Motivated by the recent works of Ohannessian & Dahleh (2012), Mossel & Ohannessian (2015), Ben-Hamou et al. (2017) and Ayed et al. (2018) on the estimation of the missing mass in species sampling models, in this paper we consider the problem of consistent estimation of $M_{n}$ under the Bernoulli product model. The classical notion of additive consistency, involving the large $n$ limiting behaviour of $\hat{T}_{n}-M_{n}$, is not suitable in the context of the estimation of $M_{n}$. This is because $M_{n}\to 0$, as $n\to+\infty$, which implies that [math] is a consistent estimator of the missing mass for any sequence $(p_{j})_{j\geq 1}$. Hence, in such a framework, one should invoke a more adequate notion of consistency, which allows to achieve more informative results. This notion of consistency is based on the limiting behaviour of the multiplicative loss function

[TABLE]

More precisely we say that the estimator $\hat{T}_{n}$ is multiplicative consistent for $M_{n}$ if $\hat{T}_{n}/M_{n}\to 1$ as $n\rightarrow+\infty$, either almost surely or in probability. The multiplicative loss function has been already used in statistics, e.g. for the estimation of small value probabilities using importance sampling (Chatterjee & Diaconis (2018)) and for the estimation of tail probabilities in extreme value theory (Beirlant & Devroye (1999)). We show that there do not exist universally consistent estimators, in the multiplicative sense, of the missing mass $M_{n}$. That is, under the Bernoulli product model and the loss function (2), we prove that for any estimator $\hat{T}_{n}$ of $M_{n}$ there exists at least a choice of $(p_{j})_{j\geq 1}$ for which $\hat{T}_{n}/M_{n}$ does not converge to $1$ in probability, as $n\rightarrow+\infty$. The proof relies on non-trivial extensions of Bayesian nonparametric ideas and techniques developed by Ayed et al. (2018) for the estimation of the missing mass in species sampling models. In particular, the key argument makes use of a generalized Indian Buffet construction (James (2017)), which allows to prove inconsistency by exploiting properties of the posterior distribution of $M_{n}$. Our inconsistency result is the natural counterpart for feature models of the work of Mossel & Ohannessian (2015), showing the impossibility of estimating the missing mass without imposing any structural (distributional) assumption on the $p_{j}$’s. We complete our study by investigating the consistency of an estimator of $M_{n}$ recently proposed by Ayed et al. (2017). To the best of our knowledge this is the first nonparametric estimator of $M_{n}$, in the sense that its derivation does not rely on any distributional assumption on the $p_{j}$’s. We show that the estimator of Ayed et al. (2017) is strongly consistent, in the multiplicative sense, under the assumption that the tail of $(p_{j})_{j\geq 1}$ decays to zero as a regularly varying function (Bingham et al. (1987)). The proof relies on novel concentration inequalities for $M_{n}$, as well as for related statistics, which are of independent interest.

The paper is structured as follows. In Section 2 we prove that for the Bernoulli product model there do not exist universally consistent estimators, in the multiplicative sense, of the missing mass $M_{n}$. Section 3 introduces some exponential tail bounds for $M_{n}$, as well as for related statistics, which are then applied in Section 4 to show that the estimator of $M_{n}$ in Ayed et al. (2018) is consistent under the assumption of regularly varying probabilities $p_{j}$’s.

2 Non existence of universally consistent estimators of the missing mass

Consider the Bernoulli product model described in the Introduction. Without loss of generality, we assume that each feature $F_{j}$ is labeled by a value in $[0,1]$ and therefore $(F_{j})_{j\geq 1}$ is a sequence of distinct points in $[0,1]$. Furthermore, the probabilities $(p_{j})_{j\geq 1}$ are assumed to be summable, i.e. $\sum_{j\geq 1}p_{j}<+\infty$; this condition is needed in order to guarantee that every observation $Y_{i}$ will display only a finite number of features almost surely. Indeed, $\sum_{j\geq 1}p_{j}<+\infty$ is equivalent to $\sum_{j\geq 1}\mathbb{P}(F_{j}\in Y_{i})=\sum_{j\geq 1}\mathbb{E}[1\!\!1_{\{F_{j}\in Y_{i}\}}]<+\infty$, which in turns implies $\sum_{j\geq 1}1\!\!1_{\{F_{j}\in Y_{i}\}}<+\infty$ almost surely, by Tonelli-Fubini Theorem. The two unknown sequences $(F_{j})_{j\geq 1}$ and $(p_{j})_{j\geq 1}$ can be uniquely encoded in a finite measure on $[0,1]$, $\sum_{j\geq 1}p_{j}\delta_{F_{j}}(\cdot)$, with all masses smaller than one. We can therefore consider as parameter space the set

[TABLE]

Recall that $X_{n,j}$ denotes the number of times that feature $F_{j}$ has been observed in the sample $(Y_{1},\ldots,Y_{n})$, that is $X_{n,j}=\sum_{1\leq i\leq n}Y_{i,j}=\sum_{1\leq i\leq n}1\!\!1_{\{F_{j}\in Y_{i}\}}$ is a Binomial random variable with parameter $(n,p_{j})$. For a fixed $n\geq 1$, an estimator $\hat{T}_{n}:[0,1]^{n}\rightarrow\mathds{R}_{+}$ of the missing mass $M_{n}$ is a measurable map which argument is the observed sample $\mathbf{Y}_{n}=(Y_{1},\ldots,Y_{n})$. We say that the estimator $\hat{T}_{n}$ is multiplicative consistent under the parameter space $\mathscr{P}$ if for every $\epsilon>0$ and every $p\in\mathscr{P}$,

[TABLE]

where $\mathbb{P}_{\mathbf{Y}_{n}|p}$ denotes the law of the observations $\mathbf{Y}_{n}$ under a feature allocation model of parameter $p$. Theorem 2.1 shows that there are no universally multiplicative consistent estimators of $M_{n}$ for the class $\mathscr{P}$. This means that for any estimator $\hat{T}_{n}$ of the missing mass, there exists at least one element $p\in\mathscr{P}$ for which $\hat{T}_{n}/M_{n}$ does not converge to $1$ in probability, as $n\rightarrow+\infty$.

Theorem 2.1

Under the feature allocation model, there are no universally consistent estimators, i.e. there are no estimators satisfying (4). In particular, for every estimator $\hat{T}_{n}$, it is possible to find an element $p\in\mathscr{P}$ such that for any $\varepsilon\in(0,1/6)$

[TABLE]

for some strictly positive constant $C$.

2.1 Proof of Theorem 2.1

In order to prove Theorem 2.1, it is enough to show that for every estimator $\hat{T}_{n}$ and every $\epsilon\in(0,1/6)$,

[TABLE]

and therefore there exists a $p\in\mathscr{P}$ for which $\hat{T}_{n}$ is not consistent.

First, let us notice that, for every $\epsilon\in(0,1/6)$,

[TABLE]

Indeed, if $\left|\frac{\hat{T}_{n}}{M_{n}}-1\right|<\epsilon$, then

[TABLE]

and, from the lower bound of (8), $\hat{T}_{n}<(1-\epsilon)M_{n}$. Because $\epsilon<1/2$, it follows that $\frac{1}{\hat{T}_{n}}<\frac{2}{M_{n}}$. This last inequality together with (8) leads to $\left|\frac{M_{n}}{\hat{T}_{n}}-1\right|<2\epsilon$. Considering the complements of the two events, it follows that

[TABLE]

and, as a consequence, $\mathbb{P}(|\frac{\hat{T}_{n}}{M_{n}}-1|\geq\epsilon)\geq\mathbb{P}(|\frac{M_{n}}{\hat{T}_{n}}-1|\geq 2\epsilon)$, proving (7). From now on, we will denote $\varepsilon=2\epsilon\in(0,1/3)$ and prove that

[TABLE]

for some strictly positive constant $C$.

The main idea of the proof is in the following formula and works as follows: we lower bound the supremum over $\mathscr{P}$ in (7) by an average with respect to a (carefully chosen) prior for $p$; we swap the conditional distribution of $\mathbf{Y}_{n}|p$ and the marginal of $p$ with the conditional of $p|\mathbf{Y}_{n}$ and the marginal of $\mathbf{Y}_{n}$; we lower bound the event probability with respect to the posterior of $p$ given $\mathbf{Y}_{n}$. Formally,

[TABLE]

where we have applied reverse Fatou’s lemma to take the $\limsup$ outside the expectation. In (10), $\mathbb{E}_{p}$ denotes the expectation with respect to the prior for $p$, $\mathbb{E}_{\mathbf{Y}_{n}}$ the expectation with respect to the marginal distribution of $\mathbf{Y}_{n}$ and $\mathbb{P}_{p|\mathbf{Y}_{n}}$ the probability under the posterior of $p$ given $\mathbf{Y}_{n}$.

Our choice of the nonparametric prior for $p$ is based on completely random measures (see Daley and Vere-Jones (2008)) and the generalized Indian Buffet process prior of James (2017). In particular, a prior for $p\in\mathscr{P}$ can defined through a completely random measure $\tilde{N}(\cdot)=\sum_{j}s_{j}\delta_{F_{j}}(\cdot)$ on $[0,1]$, where $(\{s_{j},F_{j}\})_{j\geq 1}$ is a Poisson Point Process on $\mathds{R}^{+}\times[0,1]$, by setting $p(\cdot)=\sum_{j}(1-e^{-s_{j}})\delta_{F_{j}}(\cdot)\in\mathscr{P}$. We select $\tilde{N}$ to be a completely random measure with Lévy intensity $\nu({\rm d}s,{\rm d}F)=e^{-s}/s\,{\rm d}s1\!\!1_{(0,1)}(F){\rm d}F$. The distribution of $\tilde{N}$ is completely characterized by its Laplace functional defined as follows,

[TABLE]

for any measurable function $f:[0,1]\to\mathds{R}^{+}$. See also Kingman (1993).

Theorem 3.1 of James (2017) provides with a distributional equality for the posterior of $\tilde{N}$ given $\mathbf{Y}_{n}$. Denoting by $F_{1}^{*},\ldots,F_{k_{n}}^{*}$ the $k_{n}$ distinct features observed in $\mathbf{Y}_{n}$, we have the following distributional equality

[TABLE]

where the $J_{\ell}$’s are non-negative random jumps and $\tilde{N}_{n}$ is an independent completely random measure with updated Lévy intensity $\nu_{n}({\rm d}s,{\rm d}F)=e^{-sn}\nu({\rm d}s,{\rm d}F)$.

Defining $A_{n}:=\{F_{1}^{*},\ldots,F_{k_{n}}^{*}\}$, from (11) we have that, for any Borel set $B$ in $\mathds{R}^{+}$, the missing mass $M_{n}$ satisfies

[TABLE]

showing that the posterior distribution of the missing mass $M_{n}$ is equal in distribution to the random variable $\int_{[0,1]}(1-e^{-s})\tilde{N}_{n}({\rm d}F)$. Besides, it is worth to introduce the random variable

[TABLE]

whose distribution can be computed exactly and turns out to be a Gamma random variable of parameters $(1,n+1)$. Indeed, from the Laplace functional, for every $x\in\mathds{R}$ we have

[TABLE]

which is the characteristic function of a Gamma$(1,n+1)$ random variable.

We now have all the necessary ingredients to prove the lower bound (10). Fix $\varepsilon\in(0,1/3)$. First note that, the inverse triangular inequality entails

[TABLE]

which implies

[TABLE]

indeed, thanks to (13), the two events together

[TABLE]

imply that

[TABLE]

where the last inequality follows from the fact that $\varepsilon<1$. Hence, from (14), we have that

[TABLE]

which may be plugged into (10) to obtain

[TABLE]

We are going to lower bound separately the two terms on the r.h.s. of (15). With regard to the first term, let us observe that the elementary inequality $x-x^{2}/2\leq 1-e^{-x}\leq x$, for $x>0$, implies that for all $j\geq 1$

[TABLE]

Summing over $j$,

[TABLE]

and therefore,

[TABLE]

As a simple consequence of the last inequality, for any $\varepsilon>0$, the event $\{S_{n}\leq\varepsilon\}$ implies the validity of $\{1-\varepsilon/2\leq M_{n}/S_{n}\leq 1\}$ and therefore we can upper bound the first term in (15) as follows

[TABLE]

where we have used the fact that the posterior distribution of $S_{n}$ is ${\rm Gamma}(1,n+1)$.

Let us now consider the second term on the r.h.s. of (15). Using again the fact that $S_{n}$ is Gamma distributed and $\varepsilon<1/3$, we have

[TABLE]

it is now easy to see that the function $f(y):=1+e^{-(1+3\varepsilon)y}-e^{-(1-3\varepsilon)y}$ is strictly positive, continuous and admits a global minimum on $\mathds{R}^{+}$ at the point $y^{*}=\log((1+3\varepsilon)/(1-3\varepsilon))/(6\varepsilon)$, therefore

[TABLE]

Using the two bounds (16) and (17) in (15), for any $\varepsilon\in(0,1/3)$ we get

[TABLE]

which completes the proof.

3 Concentration inequalities for feature models

In this section we will establish exponential tail bounds for the missing mass $M_{n}$ and the statistic $K_{n,r}$ defined by

[TABLE]

which counts the number of features observed with frequency $r$ in the sample $\mathbf{Y}_{n}$. The statistic $K_{n,r}$ is of interest in different applications of feature allocation models and its analysis will be important for the study of the estimator of missing mass considered in Section 4, which involves $K_{n,1}$. The tail bounds we present in this Section are valid in full generality, i.e. without any assumptions on the probability masses $(p_{j})_{j\geq 1}$. In Section 4, we will use these results to prove consistency results under the assumption of regularly varying heavy tails $(p_{j})_{j\geq 1}$.

In order to derive the concentration inequalities for $K_{n,r}$ we will use Chernoff bounds, which require suitable bounds on the log-Laplace transform. First, let us recall some definitions from Boucheron et al. (2013) and Ben-Hamou et al. (2017).

Definition 3.1

Let $X$ be a real–valued random variable defined on some probability space, then:

• i.

$X$* is sub-Gaussian on the right tail (resp. on the left tail) with variance factor $v$ if for any $\lambda\geq 0$ (resp. $\lambda\leq 0$)*

[TABLE]

• ii.

$X$* is sub-Gamma on the right tail with variance factor $v$ and scale parameter $c$ if*

[TABLE]

• iii.

$X$* is sub-Gamma on the left tail with variance factor $v$ and scale parameter $c$ if $-X$ is sub-gamma on the right tail with variance factor $v$ and scale parameter $c$;*

• iv.

$X$* is sub-Poisson with variance factor $v$ if for all $\lambda\in\mathds{R}$*

[TABLE]

being $\phi(\lambda)=e^{\lambda}-1-\lambda$.

Note that a sub-Gaussian random variable is also sub-Gamma for any choice of the scale parameter $c$, but in general the inverse is not true. As we will see in the sequel, the bounds on the log-Laplace (18)–(19) imply exponential tails bounds by means of the Chernoff inequality. See Boucheron et al. (2013) for the details.

The following proposition shows that the missing mass $M_{n}$ is sub-Gaussian on the left tail and sub-Gamma on the right one.

Proposition 3.1

Let $n>2$. On the left tail, the random variable $M_{n}$ is sub-Gaussian with variance factor $v_{n}^{-}:=2\mathbb{E}[K_{n+2,2}]/((n+2)\cdot(n+1))$, i.e. for any $\lambda\leq 0$ it holds

[TABLE]

On the right tail, the random variable $M_{n}$ is sub-Gamma with variance factor $v_{n}^{+}:=2\mathbb{E}[K_{n}]/(n^{2}-2n)$ and scale parameter $1/n$, i.e. for any $0\leq\lambda<1/n$ one has

[TABLE]

Proof.

We first focus on the proof of (21). Let $\lambda\leq 0$, exploiting the independence of the random variables $X_{n,j}$’s and the elementary inequality $\log(z)\leq z-1$, valid for any $z>0$, we obtain

[TABLE]

We observe that, being $\lambda\leq 0$, one has:

[TABLE]

hence (21) has been proven.

We now concentrate on the proof of (22), arguing exactly as before we obtain that

[TABLE]

where we have used the infinite series representation for the exponential function. Fixing the useful notation

[TABLE]

and observing that $\Phi_{n,k}\leq\Phi_{n}$, for any $k\geq 1$, we get

[TABLE]

for any $0<\lambda<1/n$. Proceeding along similar lines as in (Gnedin et al., 2007, Lemma 1), it is not difficult to see that

[TABLE]

which entails $\Phi_{n}\leq\mathbb{E}[K_{n}]/(1-2/n)$, for any $n>2$. The last inequality can be used to provide an upper bound for the r.h.s. of (23) as follows

[TABLE]

and (22) has been now proved. ∎

As already mentioned at the beginning of this section, the sub-Gaussian and sub-Gamma bounds obtained in Proposition 3.1 imply useful exponential tail bounds for $M_{n}$ (see Boucheron et al. (2013)). More specifically we have that:

Corollary 3.1

For any $n>2$ and $x\geq 0$, the following hold

[TABLE]

Proof.

The two inequalities follow by the Chernoff bound and the log-Laplace bound proved in Proposition 3.1. This is a standard argument, see Boucheron et al. (2013) for details. ∎

Proceeding along similar lines as before we show that $K_{n,r}$ is a sub-Poisson random variable, this result is implicitly proved in the Supplementary material by Ayed et al. (2017), but for the sake of completeness we report it also here.

Proposition 3.2

For any $r\geq 1$ and $n\geq 1$, the random variable $K_{n,r}$ is sub-Poisson with variance factor $\mathbb{E}[K_{n,r}]$. Indeed, for any $\lambda\in\mathds{R}$ the following bound holds true

[TABLE]

where $\phi(\lambda):=e^{\lambda}-1-\lambda$.

Proof.

Exploiting the independence of the random variables $X_{n,j}$’s, for any $\lambda\in\mathds{R}$ we can write:

[TABLE]

where we have used the inequality $\log(z)\leq z-1$, for any $z>0$. ∎

The previous proposition and the Chernoff bounds imply an exponential tail bound for $K_{n,r}$, indeed one can prove that

Corollary 3.2

For any $n\geq 1$, $r\geq 1$ and $x\geq 0$ the following holds true

[TABLE]

Corollary 3.1 and 3.2 provide us with concentration inequalities of the missing mass and the statistic $K_{n,r}$, respectively, around their mean. These results have been derived without any assumption on the probabilities $(p_{j})_{j\geq 1}$ and hold for all elements of $\mathscr{P}$. In the next Section, we will focus on the class of regularly varying probabilities and, after recalling the nonparametric estimator proposed by Ayed et al. (2018) we will prove that this estimator is consistent within such a subset of $\mathscr{P}$.

4 A consistent estimator for regularly varying feature probabilities

Ayed et al. (2018) have introduced a nonparametric estimator of the missing mass, defined as follows

[TABLE]

Namely, $\hat{M}_{n}$ is the number of features having frequency one divided by the sample size $n$. Such an estimator is attractive both from a theoretical and a computational standpoint. Indeed, on the one side, it admits two different interpretations as a Jackknife estimator in the sense of Quenouille (1956) and as a non-parametric empirical Bayes estimator in the same spirit as Efron and Morris (1973); on the other side, it is feasible and easy to implement. See Ayed et al. (2018) for details. Here we want to study the consistency of (26). In order to do this we have seen that, without assumptions on the features’ proportions, any estimator of the missing mass is always inconsistent (Theorem 2.1), hence we study the consistency of (26) under the ubiquitous assumption of heavy tailed probabilities $(p_{j})_{j\geq 1}$. We rely on the theory of regular variation by Karamata, J. (1930, 1933) (see also Karlin (1967)) to define a suitable class of heavy-tailed $(p_{j})_{j\geq 1}$, showing that, under this class, $\hat{M}_{n}$ turns out to be multiplicative consistent.

We use the limiting notation $f\simeq g$ to mean $f/g\rightarrow 1$; we further write $f\lesssim g$ if there exists a fixed constant $C>0$ such that $f\leq Cg$. Then, similarly as done by Karlin (1967) we give the following

Definition 4.1

Let $\nu(\mathrm{d}x):=\sum_{i\geq 1}\delta_{p_{i}}(\mathrm{d}x)$ and define the measure $\overline{\nu}(x):=\nu[x,1]$, which is the cumulative count of all features having no less than a certain probability mass. We say that $(p_{j})_{j\geq 1}$ is regularly varying with regular variation index $\alpha\in(0,1)$ if $\overline{\nu}(x)\simeq x^{-\alpha}\ell(1/x)$ as $x\downarrow 0$, where $\ell(t)$ is a slowly varying function, that is $\ell(ct)/\ell(t)\to 1$ as $t\to+\infty$ for all $c>0$.

Let us remark that if we denote $(p_{[j]})_{j\geq 1}$ the sorted probabilities in decreasing order, definition 4.1 is equivalent to

[TABLE]

as $j\rightarrow\infty$, where $\ell_{*}$ is another slowly varying function. For simplicity, the relation between $\ell$, $\ell_{*}$ and $\alpha$ is skipped here, interested readers can refer to Lemma 22 and Proposition 23 of Gnedin et al. (2007). Definition 4.1 is in the same spirit as Karlin (1967), but for our purposes here we consider the case $\sum_{j\geq 1}p_{j}<+\infty$, while in Karlin (1967) the $p_{j}$’s satisfy the more restrictive condition $\sum_{j\geq 1}p_{j}=1$. The next theorem is similar to a result proved by Karlin (1967) and provides the first order asymptotic of $\mathbb{E}K_{n,r}$.

Theorem 4.1

Let $(p_{j})_{j\geq 1}$ be regularly varying with $\alpha\in(0,1)$. If $\Gamma(\cdot)$ denotes the Gamma function, then as $n\to+\infty$, $\mathbb{E}[K_{n,r}]\simeq\frac{\alpha\Gamma(r-\alpha)}{r!}n^{\alpha}\ell(n)$.

Proof.

It is worth to recall the notation already used in Section 3

[TABLE]

roughly speaking $\Phi_{n,r}$ can be considered an asymptotic approximation of $\mathbb{E}K_{n,r}$. Indeed, in order to prove the theorem, we first show that $\Phi_{n,r}\simeq\frac{\alpha\Gamma(r-\alpha)}{r!}n^{\alpha}\ell(n)$ as $n\to+\infty$, and then we prove that $\Phi_{n,r}\simeq\mathbb{E}[K_{n,r}]$. In order to prove the former asymptotic equivalence it is worth noticing that (Gnedin et al., 2007, Proposition 13) applies also for the feature setting under regularly varying heavy tails, indeed the measure defined by $\nu_{r}({\rm d}p):=p^{r}\nu({\rm d}p)$ is such that

[TABLE]

Since $\Phi_{n,r}=n^{r}/r!\int_{0}^{1}e^{-np}\nu_{r}({\rm d}p)$ is the Laplace transform of $\Phi_{n,r}$ multiplied by a suitable quantity, we can apply Tauberian theorems to connect the asymptotic behaviour of the cumulative distribution function of $\nu_{r}$ given in (27) to that of $\Phi_{n,r}$. In particular, from Tauberian theorems (see Feller (1971)), we obtain

[TABLE]

as $n\to+\infty$. As a byproduct of (28), we get $\Phi_{n,r}\to+\infty$. Finally to show $\Phi_{n,r}\simeq\mathbb{E}[K_{n,r}]$, we can easily observe that (Gnedin et al., 2007, Lemma 1) applies in this setting as well, hence there exists a constant $c$ such that

[TABLE]

as $n\to+\infty$. From (29), along with $\Phi_{n,r}\to+\infty$, we obtain

[TABLE]

in other words we have shown that $\Phi_{n,r}\simeq\mathbb{E}[K_{n,r}]$ as $n\to+\infty$.

∎

We are now ready to prove that $\hat{M}_{n}$ is multiplicative consistent, when the feature probabilities $(p_{j})_{j\geq 1}$ are regularly varying. In the proof we will use the concentration inequalities of Section 3 along with Theorem 4.1 to tune the concentration inequalities under the assumption of regular variation.

Proposition 4.1

Let $(p_{j})_{j\geq 1}$ be regularly varying with index $\alpha\in(0,1)$. Let $\hat{M}_{n}:=K_{n,1}/n$ be the nonparametric estimator of the missing mass in a sample of size $n$, then $\hat{M}_{n}$ is strongly multiplicative consistent, i.e. $M_{n}/\hat{M}_{n}\stackrel{{\scriptstyle a.s.}}{{\longrightarrow}}1$.

Proof.

In order to prove the multiplicative consistency we first show that $K_{n,1}/\mathbb{E}[K_{n,1}]\stackrel{{\scriptstyle a.s.}}{{\longrightarrow}}1$ and that $M_{n}/\mathbb{E}[M_{n}]\stackrel{{\scriptstyle a.s.}}{{\longrightarrow}}1$. As for the former convergence, we can use the concentration inequality (25) given in Corollary 3.2 when $r=1$, which, for any $\varepsilon>0$, gives

[TABLE]

When $\varepsilon>0$ is fixed, we can use the asymptotic $\mathbb{E}[K_{n,1}]\simeq\alpha\Gamma(1-\alpha)n^{\alpha}\ell(n)$ in Theorem 4.1 to say that

[TABLE]

which implies that for any $\varepsilon>0$, $\mathbb{P}(\limsup_{n}(|K_{n,1}/\mathbb{E}[K_{n,1}]-1|\geq\varepsilon))=0$ by the first Borel-Cantelli lemma, hence $K_{n,1}/\mathbb{E}[K_{n,1}]\stackrel{{\scriptstyle a.s.}}{{\longrightarrow}}1$.

Analogously we may use Corollary 3.1 to prove the almost sure convergence to $1$ of the ratio $M_{n}/\mathbb{E}[M_{n}]$. Indeed, for any $\varepsilon>0$, we have

[TABLE]

By observing that $\mathbb{E}[M_{n}]=\mathbb{E}[K_{n+1,1}]/(n+1)$, the previous upper bound boils down to

[TABLE]

Now, using again Theorem 4.1, it is not difficult to see that for any fixed $\varepsilon>0$

[TABLE]

then, by the first Borel-Cantelli lemma, we get $M_{n}/\mathbb{E}[M_{n}]\stackrel{{\scriptstyle a.s.}}{{\longrightarrow}}1$, as well.

By the previous results the consistency of $\hat{M}_{n}$ easily follows, indeed

[TABLE]

since all the ratios on the r.h.s. converge to $1$ almost surely. ∎