chi_sq_test.R

library(tibble)
library(dplyr)
library(purrr)

chisq_test_distr <- function(sample, fit_p, nbins = floor(length(sample) / 8)) {
  cut_quantiles <- seq(0, 1, length.out = nbins + 1)
  cut_values <- quantile(sample, cut_quantiles)
  cut_values[1] <- -Inf

  cuts <- tibble(
    bin = 1:(nbins+1),
    min_quantile = cut_quantiles,
    max_quantile = lead(min_quantile, default = Inf),
    min_value = cut_values,
    max_value = lead(cut_values, default = Inf),
    n_sample = map2_dbl(min_value, max_value, ~sum(sample > .x & sample <= .y)),
    n_fit = map2_dbl(
      min_value, max_value, 
      ~fit_p(.y) - fit_p(.x)
    ) * length(sample)
  )[1:nbins, ]
  
  chisq_stat <- sum((cuts$n_sample - cuts$n_fit)^2 / cuts$n_fit)
  degrees_of_freedom <- nbins - 1 - 1
  p_value <- 1 - pchisq(
    chisq_stat, 
    df = degrees_of_freedom
  )
  
  list(
    cuts = cuts, 
    chisq_stat = chisq_stat, 
    df = degrees_of_freedom,
    p_value = p_value
  )
}

local({
  # tests
  norm_sample <- rnorm(1000)
  fit <- chisq_test_distr(norm_sample, function(x) pnorm(x, mean = 0, sd = 1))
  stopifnot(fit$p_value > 0.001)
  
  norm_sample_mean_2 <- rnorm(100, mean = 2)
  fit <- chisq_test_distr(norm_sample_mean_2, function(x) pnorm(x, mean = 0, sd = 1))
  stopifnot(fit$p_value < 0.05)
})