1
$\begingroup$

We've got a sample $X_1,\dots,X_n$ from an unknown distribution $F$, and we want to estimate it. The problem is that it is discrete and with support $S=\{0,1/100,1/99,\dots,1/2,1\}$, and quite many of the possible values are not observed in the sample.

In the continuous case I'd try kernel density estimation or a mixture model, but for discrete data with that support I'm not sure of what to do. I could try KDE anyway and "discretize" it, but we are in trouble with $\hat{F}(0)$:

X <- c(0.08333, 0.33333, 1, 0.14286, 1, 0.03571, 0.5, 0.11111, 1, 
       0, 1, 1, 0, 0.16667, 1, 0.16667, 1, 0.2, 0.5, 1, 0.5, 0, 0.05556, 
       0.5, 0.125, 0.07143, 0.5, 0.2, 0.5, 1, 1, 0, 1, 1, 0.5, 0.05263, 
       1, 1, 0.25, 1, 0, 1, 0.33333, 0.5, 0.2, 0.09091, 1, 1, 0.5, 0.03571)
S <- round(c(0, 1/100:1), 5) # rounded so it's easier to use as an example

library(ks)
k <- kde(X)
plot(ecdf(X))
xx <- seq(-2,1.2,.01)
lines(xx, pkde(xx, k))

Fhat <- function(x) {
  F01 <- pkde(0:1, k)
  p <- (pkde(x, k)-F01[1]) / (F01[2]-F01[1])
  return( pmax(0, pmin(1, p)) )
}
lines(xx, Fhat(xx), col = "red")

enter image description here

If $F$ were continuous in $[0,1]$ I could bound $\hat{F}$ so that $\hat{F}(0)=0$ and $\hat{F}(1)=1$, but in our discrete case that would make $P(X=0)=0$ no matter what. In fact, the support seems tricky because most of the possible values are close to 0.

So what is the best way to approach this?

Where each value comes from: imagine a list of $n>>100$ elements, some of which meet some criteria. A function returns the inverse of the rank of the first element meeting the criteria if it is in the top 100, or 0 if it is below the top 100 (this is why there are many zeros). So if the first such element is the third, it returns $1/3$; if it's the 106th, it returns $0$. There are always elements meeting the criteria, but that the first of them might be ranked below 100. In this sense, maybe it makes sense to have $F(0)=0$?

However, there can be other functions that don't necessarily return something like $1/i$, but rather have some arbitrary support within $[0,1]$, such as $\{0, 0.2, 0.4, ..., 1\}$. It's also important to note that there is a notion of distance between two values; it's just that not all real numbers are observable.

I didn't want to use any parametric assumption because I don't really see one and because in the general case we can have an arbitrary support, so I thought that I had to use some kind of smoothing.

$\endgroup$
  • 2
    $\begingroup$ Could you tell us more about your data? You say it's discrete, but you present us a bunch of real numbers -- it's confusing... (I'm not saying that discrete data cannot be like this, but it is not perfectly clear what is your data.) $\endgroup$ – Tim Dec 19 '16 at 16:19
  • $\begingroup$ These data are the output of a function that computes a bunch of things for a list of elements, and returns the inverse of the rank of the first element which meets some criteria. So if the first such element is the third, it returns $1/3$. However, there can be other functions that don't necessarily return something like $1/i$, but rather have some arbitrary support within $[0,1]$. $\endgroup$ – Julián Urbano Dec 19 '16 at 16:37
  • $\begingroup$ In your example you show ranks only going up to 100. To estimate its cdf what is your sample size? Why do you want to fill in values outside the support of the distribution. $\endgroup$ – Michael R. Chernick Dec 19 '16 at 19:03
  • $\begingroup$ @MichaelChernick the sample size is usually around 50, 100 at the most. I don't understand what you mean by "values outside the distribution". Where did I use such values? All values in X are in the support S, up to rounding. $\endgroup$ – Julián Urbano Dec 19 '16 at 21:18
  • 1
    $\begingroup$ We can give better answers if you provide context. $\endgroup$ – Kodiologist Dec 22 '16 at 19:47
1
+50
$\begingroup$

You'll have many more tools at your disposal if you choose to model the rank itself rather than its reciprocal. There are a few reasons for this: Almost no common parametric families exist with support over just the rational numbers (or a subset thereof). Additionally, KDE is complicated by the discontinuous support and the fact that it's bounded. Lastly, your data should be viewed as truncated data. It's easier to model a distribution whose tail is truncated rather than a bounded truncated region within its support.

Fortunately, there's an easy workaround for these issues. If you choose to model the rank and you get distribution estimates via some convenient method, you can easily transform the results at the end to model the reciprocal. Let's say you have the true (but unobserved) rank data given by $R_1,\dots,R_n$. Define the observed rank data as $Z_i = R_i$ if $R_i \le 100$ and $Z_i := ``Tail"$ otherwise. Parametric modeling options include truncated versions of any semi-infinite discrete distribution.

For instance, if you went with the appropriately heavy-tailed Zipf distribution, the original pdf would be $\hat{P}(R=k|s) = k^{-s}/\zeta(s)$, for $k=1,2,...$. But to get the truncated version, you'd use $\hat{P}(Z=k|s,\theta)=(1-\theta)k^{-s}/\zeta(s)$, for $k=1,2,...,100$, and $\hat{P}(Z=``Tail"|s,\theta)=\theta$.

If you'd prefer a nonparametric method, you could still use a KDE. However, you'll have to use a discrete kernel with an adaptive bandwidth. Something like this: $$\hat{P}(Z=x) = \dfrac 1 n \sum_{i=1}^n K(x,X_i|h(x))$$ where $K$ is a pmf over the integers which is parameterized by a nonconstant $h(x)$. For instance, you could set $K(x,X_i|h(x)) = \dfrac {[h(x)]^{|X_i - x|} e^{-h(x)}} {|X_i - x| !} \dfrac {(1+I(X_i=x))} 2$, where perhaps $h(x)=\beta x^2$ for some $\beta>0$. Defining the bandwidth function this way would then treat two equally spaced large ranks as closer together than two equally spaced small ranks. This might be useful since you mentioned that the distance between rank reciprocals was meaningful.

$\endgroup$
  • $\begingroup$ Thanks for the idea. In the first option, I guess I should estimate $\theta$ as the frequency of zeros in my sample, right? In the second option, do you have any indication as to how to implement this custom kernel in R? $\endgroup$ – Julián Urbano Dec 29 '16 at 2:47

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.