R Script for the significance of CA's Dimensions
In order to perform the test in R you can:
-(if you are familiar with R) copy/paste the following lines in the R script Editor and save it as a Script
-copy/paste into the R console and press <return> in order to run the script
Note: the required Package (i.e., FactoMineR) has to be installed first. In any event, my script for CA in R already implements this test. Also note that my R package 'CAinterprTools' (described in this site) implements the Malinvaud test as well as a permutation test for the significance of the CA dimensions.
-(if you are familiar with R) copy/paste the following lines in the R script Editor and save it as a Script
-copy/paste into the R console and press <return> in order to run the script
Note: the required Package (i.e., FactoMineR) has to be installed first. In any event, my script for CA in R already implements this test. Also note that my R package 'CAinterprTools' (described in this site) implements the Malinvaud test as well as a permutation test for the significance of the CA dimensions.
Code:
library(FactoMineR)
data <- read.table(file=file.choose(), header=TRUE)
grandtotal <- sum(data)
nrows <- nrow(data)
ncols <- ncol(data)
numb.dim.cols<-ncol(data)-1
numb.dim.rows<-nrow(data)-1
a <- min(numb.dim.cols, numb.dim.rows) #dimensionality of the table
labs<-c(1:a) #set the number that will be used as x-axis' labels on the scatterplots
res.ca<-CA(data, ncp=a, graph=FALSE)
malinv.test.rows <- a
malinv.test.cols <- 6
malinvt.output <-matrix(ncol= malinv.test.cols, nrow=malinv.test.rows)
colnames(malinvt.output) <- c("K", "Dimension", "Eigen value", "Chi-square", "df", "p value")
malinvt.output[,1] <- c(0:(a-1))
malinvt.output[,2] <- c(1:a)
for(i in 1:malinv.test.rows){
k <- -1+i
malinvt.output[i,3] <- res.ca$eig[i,1]
malinvt.output[i,5] <- (nrows-k-1)*(ncols-k-1)
}
malinvt.output[,4] <- rev(cumsum(rev(malinvt.output[,3])))*grandtotal
malinvt.output[,6] <- round(pchisq(malinvt.output[,4], malinvt.output[,5], lower.tail=FALSE), digits=6)
print(malinvt.output)
dev.new()
plot(malinvt.output[,6], , xaxt="n", xlim=c(1, a), xlab="Dimensions", ylab="p value")
axis(1, at=labs, labels=sprintf("%.0f",labs))
title(main="Malinvaud's test Plot", sub="dashed line: alpha 0.05 threshold", col.sub="RED", cex.sub=0.80)
abline(h=0.05, lty=2, col="RED")
Have you found this website helpful? Consider to leave a comment in this page.