4
$\begingroup$

I am trying to write a simple piece of code to produce the sample space for the Kolmogorov-Smirnov statistic $D_{m,n} = \underset{x}{\max}\left|S_{m}(x)-S_{n}(x)\right|$, where $m$ and $n$ are the sample sizes, and $S_{m}(x)$ and $S_{n}(x)$ the two empirical distribution functions of the two samples. In plain language, $D_{m,n}$ is the absolute value of the largest difference between the two empirical distribution functions at points $x$ (Gibbons & Chakraborti, 2010).

My thinking is that the possible values of $D_{m,n}$ range from $\min\left(\frac{1}{m},\frac{1}{n}\right)$ at the lower bound (because at the very least the difference between $S_{m}(x)$ and $S_{n}(x)$ at the smallest value of $x$ will be either $\frac{1}{m}$ or $\frac{1}{n}$), and $1$ at the upper bound (because it is possible one empirical distribution function reaches its maximum while the other is still at its minimum).

Intermediate values are all the differences in the set comprising $\min\left(\left|\frac{i}{m}-\frac{j}{n}\right|,\min\left(\frac{1}{m},\frac{1}{n}\right)\right)$ for $i=1,2,\dots,m$ and $j=1,2,\dots,n$.

Here is my R code:

KSSS <- function(m, n) {
  # start with steps based on the inverse of each sample size
  values <- sort(unique(c((1:m)/m, (1:n)/n)))
  nv <- length(values)
  additional.values <- c()  # list the magnitude of differences in values
  for (i in 1:(nv-1) ) {
    for (j in (i+1):nv) {
      additional.values <- c(additional.values, abs(values[i]-values[j]))
      }
    }
  if (length(additional.values[additional.values >= min(1/m,1/n)]) == 0) {
    # If no additional.values are greater than the smallest possible 
    # value of D, return values
    return(values)
    }
   else{
    # Include any additional.values greater than the smallest 
    # possible value of D
    values <- sort(c(additional.values, values))
    values <- values[values >= min(1/nx,1/ny)]
    return(values[!duplicated(round(values,12))])
    }
  }

My questions:

  1. Is my reasoning about the set of possible values of $D_{m,n}$ correct?

  2. Does my code accomplish this? (General improvements to the algorithm welcome :).

References:
Gibbons, J. D. and Chakraborti, S. (2010). Nonparametric Statistical Inference. Chapman & Hall/CRC, 5th edition edition.

$\endgroup$
4
  • 1
    $\begingroup$ It is possible for $D_{m,n}$ to be zero: this will happen for two identical datasets (up to permutation). Since $S_m \in \{0, 1/m, 2/m, \ldots, 1\}$ and $S_n\in\{0,1/n,2/n, \ldots, 1\},$ wouldn't it be enough to say that the sample space is a subset of the set of their absolute differences? That could be generated (in ascending order) in R with a single line, function(m,n) sort(unique(c(outer(0:m, 0:n, function(i,j) abs(i*n-j*m))))) / (m*n) Indeed, by leaving off the final /(m*n) you will generate the numerators of all the fractions, producing exact results. $\endgroup$
    – whuber
    Jun 29 '18 at 18:08
  • $\begingroup$ Thank you, whuber (for both insights)! I would be delighted to accept that as an answer. :) $\endgroup$
    – Alexis
    Jun 29 '18 at 18:18
  • 1
    $\begingroup$ Hmm... I was really just feeling you out for clarification of your question :-). $\endgroup$
    – whuber
    Jun 29 '18 at 18:25
  • $\begingroup$ @whuber too bad you are such a natural with the bullseye answers. ;) $\endgroup$
    – Alexis
    Jun 29 '18 at 18:27
3
$\begingroup$

It is possible for $D_{m,n}$ to be zero: this will happen when the K-S test is applied to two identical datasets (up to permutation).

Since $S_m \in \{0, 1/m, 2/m, \ldots, 1\}$ and $S_n\in\{0,1/n,2/n, \ldots, 1\},$ wouldn't it be enough to say that the sample space is a subset of the set of their absolute differences? A function to do this could be written in R with a single line,

function(m,n) unique(c(outer(0:m, 0:n, function(i,j) abs(i*n-j*m)))) / (m*n)

outer generates all ordered pairs $(i,j)$ in $\{0,1,\ldots, m\}\times \{0,1,\ldots,n\}$ and for each one computes $|in-jm|.$ c flattens this from a matrix to an array and unique converts that to the set of its elements (in some arbitrary order). The final division by /(m*n) yields all the the distinct values

$$\frac{|in-jm|}{mn}= \left| \frac{i}{m} - \frac{j}{n} \right|.$$

You could sort these (using sort) if you like. Remove the $0$ when $m\ne n,$ too, if that helps. By leaving off the final /(m*n) you will generate the numerators of all the fractions, producing exact results (at least when $mn$ doesn't overflow the long integer storage mode in R) if that's desired.

$\endgroup$
1
  • 1
    $\begingroup$ Just noticed that if $m=n$, then the result should just be c(0:m)/m. $\endgroup$
    – Alexis
    Nov 26 '19 at 17:26

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.