[This article was first published on Statistic on aiR, and kindly contributed to R-bloggers]. (You can report issue about the content on this page here)
Want to share your content on R-bloggers? click here if you have a blog, or here if you don't.
This is the code to perform the Bhapkar V test. I’ve rapidly wrote it, in 2 hours. The code is then quite brutal and it could be done better. As soon as possible, I will correct it.Want to share your content on R-bloggers? click here if you have a blog, or here if you don't.
WARNING: it works *ONLY* with 3 groups, for now!
bhapkar.test.3g <- function(data1=list){ sample <- c() for(i in 1:length(data1)){ sample <- c(sample, rep(i, length(data1[[i]]))) } obs <- c() for(i in 1:length(data1)){ obs <- c(obs, data1[[i]]) } rank <- rank(obs) cplets <- list() vec <- c() for(i in 1:length(data1[[1]])){ vec < - c(vec, (length(data1[[2]][data1[[2]]>data1[[1]][i]]) * length(data1[[3]][data1[[3]]>data1[[1]][i]]))) } cplets[[1]] <- vec vec <- c() for(i in 1:length(data1[[2]])){ vec < - c(vec, (length(data1[[1]][data1[[1]]>data1[[2]][i]]) * length(data1[[3]][data1[[3]]>data1[[2]][i]]))) } cplets[[2]] <- vec vec <- c() for(i in 1:length(data1[[3]])){ vec < - c(vec, (length(data1[[2]][data1[[2]]>data1[[3]][i]]) * length(data1[[1]][data1[[1]]>data1[[3]][i]]))) } cplets[[3]] <- vec cplets1 <- c(cplets[[1]], cplets[[2]], cplets[[3]]) mydata <- data.frame(obs=obs, sample=sample, rank=rank, cplets=cplets1) v1 <- sum(cplets[[1]]) v2 <- sum(cplets[[2]]) v3 <- sum(cplets[[3]]) vtot <- v1+v2+v3 u1 <- v1/vtot u2 <- v2/vtot u3 <- v3/vtot u <- c(u1,u2,u3) lengths <- c(length(data1[[1]]), length(data1[[2]]), length(data1[[3]])) N <- sum(lengths) P <- c(lengths / N) ngroup <- length(data1) V <- N * (2*length(data1)-1)* (sum(P*((u-1/ngroup)^2)) - (sum(P*((u-1/ngroup))))^2) prop <- pchisq(V, df=length(data1)-1) names(V) = "V = " method = "Bhapkar V-test" rval <- list(method = method, statistic = V, p.value = prop) class(rval) = "htest" return(rval) }
An example:
a <- c(42, 46, 48.5, 49, 68, 51) b <- c(70.5, 54, 60,72) c <- c(66, 54, 43, 105, 94) mydata <- list(a,b,c) bhapkar.test.3g(mydata) Bhapkar V-test data: V = 6.7713, p-value = 0.9661
REFERENCES:
Statistical analysis of nonnormal data
By J. V. Deshpande, A. P. Gore, A. Shanubhogue
pag. 61
To leave a comment for the author, please follow the link and comment on their blog: Statistic on aiR.
R-bloggers.com offers daily e-mail updates about R news and tutorials about learning R and many other topics. Click here if you're looking to post or find an R/data-science job.
Want to share your content on R-bloggers? click here if you have a blog, or here if you don't.