R/TestUNey.R

In MissMech: Testing Homoscedasticity, Multivariate Normality, and Missing Completely at Random

#' Test of Goodness of Fit (Uniformity)
#' @rdname TestUNey
#' @aliases TestUNey SimNey
#' 
#' @description
#' This routine tests whether the values in a vector x is distributed as uniform (0,1).  
#' The Neyman's smooth test of fit, as described by Ladwina (1994) is used. 
#' The p-values are obtained based on a resampling method from uniform (0,1).
#' 
#' @param x A vector of values, each in the interval [0,1].
#' @param nrep The number of replications used to simulate the Neyman distribution.
#' @param sim A vector of simulated values from the Neyman distribution. If sim = NA this vector is generated by the function SimNev, otherwise the vector inputted is used.
#' @param n.min The minimum sample size that triggers the use of asymptotic Chi distribution in place of the emprical distribution in the Neyman test of uniformity.
#' 
#' @returns
#' \item{pn }{The p-value for the test.}
#' \item{n4 }{The value of the test statistics.}
#' 
#' @references 
#' Ledwina, T. (1994). ``Data-driven version of neyman's smooth test of fit,'' \emph{Journal of the American Statistical Association,} 89, 1000-1005, \doi{10.2307/2290926}.
#' @author 
#' Mortaza Jamshidian, Siavash Jalal, and Camden Jansen
#' 
#' @examples
#' # Example 1
#' x <- runif(100)
#' TestUNey(x, nrep = 10000, sim = NA)
#' 
#' # Example 2
#' x <- runif(30,2,5)
#' x <- (x-min(x))/(max(x)-min(x))
#' TestUNey(x, nrep = 10000, sim = NA)
#' 
#' # Example 3
#' x <- c(0.6,0.6,0.5,0.7,0.3,0.4,0.5,0.4,0.2,0.4,0.2,0.5,0.7,0.1,0.7,0.1,0.5,0.5,0.4,0.6,0.3)
#' TestUNey(x, nrep = 10000, sim = NA)
#' 
#' @export
TestUNey <- function(x, nrep = 10000, sim = NA, n.min = 30)
{
  # This routine tests whether the values in each row of x are unif(0,1). It
  # uses the Neyman's smooth test (see e.g., Ledwina 1994, TAS)
  # x is a vector
  # P-values are computed based on a
  # resampling method from unif(0,1).
  # All values of $x$ are between 0 and 1
  
  n <- length(x)
  pi <- LegNorm(x)
  n4 <- (apply(pi$p1, 2, sum) ^ 2 + apply(pi$p2, 2, sum) ^ 2 + 
           apply(pi$p3, 2, sum) ^ 2 + apply(pi$p4, 2, sum) ^ 2) / n 
  if (n < n.min){
    if(is.na(sim[1])) {
      sim <- SimNey(n, nrep)
    }
    pn <- length(which(sim > n4)) / nrep
  } else {
    pn <- stats::pchisq(n4, 4, lower.tail = FALSE)
  }
  list(pn = pn, n4 = n4)
}

#' @rdname TestUNey
#' @aliases TestUNey SimNey
#' @usage NULL
#' @export
SimNey <- function(n, nrep)
{
  x <- matrix(stats::runif(nrep * n), ncol = nrep)
  pi <- LegNorm(x)
  n4sim <- (apply(pi$p1, 2, sum) ^ 2 + apply(pi$p2, 2, sum) ^ 2 + 
              apply(pi$p3, 2, sum) ^ 2 + apply(pi$p4, 2, sum) ^ 2) / n 
  n4sim <- sort(n4sim)
}