{"id":931,"date":"2024-05-25T20:21:17","date_gmt":"2024-05-26T00:21:17","guid":{"rendered":"https:\/\/www.clayford.net\/statistics\/?p=931"},"modified":"2024-05-29T21:48:02","modified_gmt":"2024-05-30T01:48:02","slug":"parametric-bootstrap-of-kolmogorov-smirnov-test","status":"publish","type":"post","link":"https:\/\/www.clayford.net\/statistics\/parametric-bootstrap-of-kolmogorov-smirnov-test\/","title":{"rendered":"Parametric Bootstrap of Kolmogorov\u2013Smirnov Test"},"content":{"rendered":"

Zeimbekakis, et al.\u00a0recently published an article in The American Statistician titled On Misuses of the Kolmogorov\u2013Smirnov Test for One-Sample Goodness-of-Fit<\/a>. One of the misues they discuss is using the KS test with parameters estimated from the sample. For example, let\u2019s sample some data from a normal distribution.<\/p>\n

x <- rnorm(200, mean = 8, sd = 8)\r\nc(xbar = mean(x), s = sd(x))<\/code><\/pre>\n
##     xbar        s \r\n## 8.333385 7.979586<\/code><\/pre>\n
If we wanted to assess the goodness-of-fit of this sample to a normal distribution, the following is a bad way to use the KS test:<\/p>\n
ks.test(x, "pnorm", mean(x), sd(x))<\/code><\/pre>\n
## \r\n##  Asymptotic one-sample Kolmogorov-Smirnov test\r\n## \r\n## data:  x\r\n## D = 0.040561, p-value = 0.8972\r\n## alternative hypothesis: two-sided<\/code><\/pre>\n
The appropriate way to use the KS test is to actually supply hypothesized parameters. For example:<\/p>\n
ks.test(x, "pnorm", 8, 8)<\/code><\/pre>\n
## \r\n##  Asymptotic one-sample Kolmogorov-Smirnov test\r\n## \r\n## data:  x\r\n## D = 0.034639, p-value = 0.9701\r\n## alternative hypothesis: two-sided<\/code><\/pre>\n
The results of both tests are the same. We fail to reject the null hypothesis that the sample is from a Normal distribution with the stated mean and standard deviation. However, the former test is very conservative. Zeimbekakis, et al.\u00a0show this via simulation. I show a simplified version of this simulation. The basic idea is that if the test were valid, the p-values would be uniformly distributed and the points in the uniform distribution QQ-plot would fall along a diagonal line. Clearly that\u2019s not the case.<\/p>\n
n <- 200\r\nrout <- replicate(n = 1000, expr = {\r\n  x <- rnorm(n, 8 , 8)\r\n  xbar <- mean(x)\r\n  s <- sd(x)\r\n  ks.test(x, "pnorm", xbar, s)$p.value\r\n})\r\nhist(rout, main = "Histogram of p-values")<\/code><\/pre>\n
<\/p>\n
qqplot(x = ppoints(n), y = rout, main = "Uniform QQ-plot")\r\nqqline(rout, distribution = qunif)<\/code><\/pre>\n
<\/p>\n
Conclusion: using fitted parameters in place of the true parameters in the KS test yields conservative results. The authors state in the abstract that this \u201chas been \u2018discovered\u2019 multiple times.\u201d<\/p>\n
When done the right way, the KS test yields uniformly distributed p-values.<\/p>\n
rout2 <- replicate(n = 1000, expr = {\r\n  x <- rnorm(n, 8 , 8)\r\n  ks.test(x, "pnorm", 8, 8)$p.value\r\n})\r\nhist(rout2)<\/code><\/pre>\n
<\/p>\n
qqplot(x = ppoints(n), y = rout2, main = "Uniform QQ-plot")\r\nqqline(rout2, distribution = qunif)<\/code><\/pre>\n
<\/p>\n
Obviously it\u2019s difficult to know which parameters to supply to the KS test. Above we knew to supply 8 as the mean and standard deviation because that\u2019s what we used to generate the data. But what to do in real life? Zeimbekakis, et al.\u00a0propose a parametric bootstrap to approximate the null distribution of the KS test statistic. The steps to implement the bootstrap are as follows:<\/p>\n
    \n
  1. draw a random sample from the fitted distribution<\/li>\n
  2. get estimates of parameters of random sample<\/li>\n
  3. obtain the empirical distribution function<\/li>\n
  4. calculate the bootstrapped KS statistic<\/li>\n
  5. repeat steps 1 – 4 many times<\/li>\n<\/ol>\n
    Let\u2019s do it. The following code is a simplified version of what the authors provide with the paper. Notice they use MASS::fitdistr()<\/code> to obtain MLE parameter estimates. This returns the same mean for the normal distribution but a slightly smaller (i.e.\u00a0biased) estimated standard deviation.<\/p>\n
    param  <- MASS::fitdistr(x, "normal")$estimate\r\nks <- ks.test(x, function(x)pnorm(x, param[1], param[2]))\r\nstat <- ks$statistic\r\nB <- 1000\r\nstat.b <- double(B)\r\nn <- length(x)\r\n\r\n## bootstrapping\r\nfor (i in 1:B) {\r\n  # (1) draw a random sample from a fitted dist\r\n  x.b <- rnorm(n, param[1], param[2])\r\n  # (2) get estimates of parameters of random sample\r\n  fitted.b <- MASS::fitdistr(x.b, "normal")$estimate\r\n  # (3) get empirical distribution function\r\n  Fn <- function(x)pnorm(x, fitted.b[1], fitted.b[2])\r\n  # (4) calculate bootstrap KS statistic\r\n  stat.b[i] <- ks.test(x.b, Fn)$statistic\r\n}\r\nmean(stat.b >= stat)<\/code><\/pre>\n
    ## [1] 0.61<\/code><\/pre>\n
    The p-value is the proportion of statistics greater than or equal to the observed statistic calculated with estimated parameters.<\/p>\n
    Let\u2019s turn this into a function and show that it returns uniformly distributed p-values when used with multiple samples. Again this is a simplified version of the R code the authors generously shared with their paper.<\/p>\n
    ks.boot <- function(x, B = 1000){\r\n  param  <- MASS::fitdistr(x, "normal")$estimate\r\n  ks <- ks.test(x, function(k)pnorm(k, param[1], param[2]))\r\n  stat <- ks$statistic\r\n  stat.b <- double(B)\r\n  n <- length(x)\r\n  for (i in 1:B) {\r\n    x.b <- rnorm(n, param[1], param[2])\r\n    fitted.b <- MASS::fitdistr(x.b, "normal")$estimate\r\n    Fn <- function(x)pnorm(x, fitted.b[1], fitted.b[2])\r\n    stat.b[i] <- ks.test(x.b, Fn)$statistic\r\n  }\r\n  mean(stat.b >= stat)\r\n}<\/code><\/pre>\n
    Now replicate the function with many samples. This takes a moment to run. It took my Windows 11 PC with an Intel i7 chip about 100 seconds to run.<\/p>\n
    rout_boot <- replicate(n = 1000, expr = {\r\n  x <- rnorm(n, 8 , 8)\r\n  ks.boot(x)\r\n})<\/code><\/pre>\n
    hist(rout_boot)<\/code><\/pre>\n
    <\/p>\n
    qqplot(x = ppoints(n), y = rout_boot, main = "Uniform QQ-plot")\r\nqqline(rout_boot, distribution = qunif)<\/code><\/pre>\n
    <\/p>\n","protected":false},"excerpt":{"rendered":"
    Zeimbekakis, et al.\u00a0recently published an article in The American Statistician titled On Misuses of the Kolmogorov\u2013Smirnov Test for One-Sample Goodness-of-Fit…. Read more<\/a><\/p>\n","protected":false},"author":1,"featured_media":0,"comment_status":"open","ping_status":"open","sticky":false,"template":"","format":"standard","meta":{"footnotes":""},"categories":[22,8],"tags":[31,83,84],"class_list":["post-931","post","type-post","status-publish","format-standard","hentry","category-hypothesis-testing","category-simulation","tag-bootstrap","tag-kolmogorov-smirnov-test","tag-regression-and-other-stories"],"_links":{"self":[{"href":"https:\/\/www.clayford.net\/statistics\/wp-json\/wp\/v2\/posts\/931","targetHints":{"allow":["GET"]}}],"collection":[{"href":"https:\/\/www.clayford.net\/statistics\/wp-json\/wp\/v2\/posts"}],"about":[{"href":"https:\/\/www.clayford.net\/statistics\/wp-json\/wp\/v2\/types\/post"}],"author":[{"embeddable":true,"href":"https:\/\/www.clayford.net\/statistics\/wp-json\/wp\/v2\/users\/1"}],"replies":[{"embeddable":true,"href":"https:\/\/www.clayford.net\/statistics\/wp-json\/wp\/v2\/comments?post=931"}],"version-history":[{"count":5,"href":"https:\/\/www.clayford.net\/statistics\/wp-json\/wp\/v2\/posts\/931\/revisions"}],"predecessor-version":[{"id":937,"href":"https:\/\/www.clayford.net\/statistics\/wp-json\/wp\/v2\/posts\/931\/revisions\/937"}],"wp:attachment":[{"href":"https:\/\/www.clayford.net\/statistics\/wp-json\/wp\/v2\/media?parent=931"}],"wp:term":[{"taxonomy":"category","embeddable":true,"href":"https:\/\/www.clayford.net\/statistics\/wp-json\/wp\/v2\/categories?post=931"},{"taxonomy":"post_tag","embeddable":true,"href":"https:\/\/www.clayford.net\/statistics\/wp-json\/wp\/v2\/tags?post=931"}],"curies":[{"name":"wp","href":"https:\/\/api.w.org\/{rel}","templated":true}]}}