R/CompClassicDist.R

In musto101/wilcox_R: What the package does (short line)

CompClassicDist <-
function(svd, out) {
	out$classic$P <- as.matrix(svd$loadings[,1:out$k])
	out$classic$L <- as.vector(svd$eigenvalues[1:out$k])
	out$classic$M <- as.vector(svd$centerofX)
	out$classic$T <- as.matrix(svd$scores[,1:out$k])

	out$classic$k <- out$k
	out$classic$Xc <- as.matrix(svd$Xcentered)
	Tclas <- out$classic$Xc %*% out$classic$P
	out$classic$sd <- sqrt(mahalanobis(Tclas, center=rep(0, length=ncol(Tclas)), cov=diag(x=out$classic$L, nrow=out$classic$k)))
	out$classic$cutoff$sd <- sqrt(qchisq(0.975, out$classic$k))
	Xtilde <- Tclas %*% t(out$classic$P)
	Cdiff <- out$classic$Xc - Xtilde
	out$classic$od <- vector(mode="numeric", length=nrow(out$classic$Xc))
	if(is.list(dimnames(out$classic$T))) {
		names(out$classic$od) <- dimnames(out$classic$T)[[1]]
	}
	out$classic$od<-apply(Cdiff,1,vecnorm)
	if(out$k != svd$rank) {
		m <- mean(out$classic$od^(2/3))
		s <- sqrt(var(out$classic$od^(2/3)))
		out$classic$cutoff$od <- sqrt((qnorm(0.975, m, s))^3)
	}
	else {
		out$classic$cutoff$od <- 0
	}
	out$classic$flag <- (out$classic$od <= rep(x=out$classic$cutoff$od, times=length(out$classic$od))) &
								(out$classic$sd <= rep(x=out$classic$cutoff$sd, times=length(out$classic$sd)))
	out$classic$class <- "CPCA"
	return(out)
}