Uniformly sample from Gaussian distributed data

Question

I have data that is roughly Gaussian distributed, bounded on a range of [x0, x1], w/ mean m and standard deviation s. I want to sample from the data so that I am equally likely to sample from the values far from the mean, m as I am from values close to the mean m. This in practice would transform the distribution of the data from N(m, s) to a uniform distribution U[x0, x1]. I can imagine several heuristic ways to doing this, but I am wondering what the best way to do this is?

Initial attempt: discretize the data into N bins bins between x0 and x1. Assign each point to a bin. Uniformly sample each bin.

This seems like a very straighforward way of doing this. Is there a better or more optimal way?

Answer 1

It doesn't matter what the distribution of the data might be. You ask to create samples in which only the data values appear, yet are as uniformly distributed as possible.

Remove any data not in the interval $[x_l, x_u].$ Let the ordered remaining values be $X_1 \le X_2 \le X_n.$ Set $X_0 = x_l$ and $X_{n+1} = x_u.$

Sample uniformly within $[x_l, x_u],$ producing a batch of values $U_1, U_2, \ldots, U_N$. Randomly replace each $U_i$ with one of the two neighboring values $X_{j(i)} \le U_i \le X_{j(i)+1}$ by selecting the larger with probability proportional to $U_i - X_{j(i)}.$

This is equivalent to putting each $X_i$ into its own bin with boundaries erected halfway between it and its nearest neighbors and sampling the $X_i$ with probabilities proportional to these bin widths. It minimizes the Kolmogorov-Smirnov distance between this distribution and the uniform distribution on $[x_l, x_u].$

Here is an example where there are $n=33$ data values in the interval $[1,3]$ and a random sample of $N=1000$ values has been created.

The left panel shows "rug plots" of the data and their sample, with vertical jittering to help resolve overlaps. You can see the $U_i$ are taken from the set of values of the $X_j.$

The middle and right panels graph the empirical cumulative distribution functions (ecdfs) of the data and the sample, respectively. For reference, as black curves, are the truncated Normal distribution from which the data were sampled (middle) and the uniform distribution on the interval $[1,3]$ (right).

An R implementation of this method is given by the function ru. Its arguments are the sample size n, the data x, and the interval $[x_l,x_u]$ as an array bounds. (It expands this interval if needed to accommodate all the values in x.) Its output is the array of $U_i.$

It is reasonably efficient, requiring $O(N+n)$ storage and $O((N+n)\log(N+n))$ time.

ru <- function(n, x, bounds = range(x)) {
  if(!missing(bounds)) bounds <- range(c(bounds, range(x)))
  r <- diff(bounds) + 1
  x <- c(bounds[1]-r, sort(x), bounds[2]+r) # Place sentinels around `x`
  u <- runif(n, bounds[1], bounds[2])       # Sample range of `bounds`
  #
  # Algorithm to find nearest `x` values to each `u` value.
  #
  a <- c(x, u)
  i <- c(rep(TRUE, length(x)), rep(FALSE, length(u)))[order(a)]
  j <- cumsum(i)[!i]
  k <- length(x) + 1 - rev(cumsum(rev(i)))[!i]
  ifelse(runif(n, a[j], a[k]) > u, a[j], a[k])
}