Nothing
plot.MIPCA <- function(x,choice="all",axes=c(1,2),new.plot=TRUE,main=NULL,level.conf=0.95, ...){
####
procrustes <- function(amat, target, orthogonal = FALSE, translate = FALSE,
magnify = FALSE) {
for (i in nrow(amat):1) {
if (any(is.na(amat)[i, ]) | any(is.na(target)[i,
])) {
amat <- amat[-i, ]
target <- target[-i, ]
}
}
dA <- dim(amat)
dX <- dim(target)
if (length(dA) != 2 || length(dX) != 2)
stop("arguments amat and target must be matrices")
if (any(dA != dX))
stop("dimensions of amat and target must match")
if (length(attr(amat, "tmat")))
stop("oblique loadings matrix not allowed for amat")
if (orthogonal) {
if (translate) {
p <- dX[1]
target.m <- (rep(1/p, p) %*% target)[, ]
amat.m <- (rep(1/p, p) %*% amat)[, ]
target.c <- scale(target, center = target.m,
scale = FALSE)
amat.c <- scale(amat, center = amat.m, scale = FALSE)
j <- svd(crossprod(target.c, amat.c))
}
else {
amat.c <- amat
j <- svd(crossprod(target, amat))
}
rot <- j$v %*% t(j$u)
if (magnify)
beta <- sum(j$d)/sum(amat.c^2)
else beta <- 1
B <- beta * amat.c %*% rot
if (translate)
B <- B + rep(as.vector(target.m), rep.int(p,
dX[2]))
value <- list(rmat = B, tmat = rot, magnify = beta)
if (translate)
value$translate <- target.m - (rot %*% amat.m)[,
]
}
else {
b <- solve(amat, target)
gamma <- sqrt(diag(solve(crossprod(b))))
rot <- b * rep(gamma, rep.int(dim(b)[1], length(gamma)))
B <- amat %*% rot
fcor <- solve(crossprod(rot))
value <- list(rmat = B, tmat = rot, correlation = fcor)
}
return(value)
}
####
res <- x
if (!inherits(res, "MIPCA")) stop("non convenient data")
ncp <- res$call$ncp
reference <- PCA(res$res.imputePCA,scale=res$call$scale,graph=FALSE,ncp=ncp)
rec.pca <- res$res.imputePCA
# rec <- reconst(reference,ncp)
# rec.pca <- as.matrix(res$call$X)
# rec.pca[res$call$missing] <- rec[res$call$missing]
res.var <- res.supp <- rec.pca
res.procrustes <- reference$ind$coord[,1:ncp]
res.dim <- as.matrix(res$res.imputePCA)
for (i in 1:dim(res$res.MI)[3]){
rec.pca <- res$res.MI[,,i]
acpfin <- PCA(rec.pca, scale=res$call$scale,graph=FALSE,ncp=ncp)
tourne <- procrustes(acpfin$ind$coord[,1:ncp], reference$ind$coord[,1:ncp],orthogonal = TRUE, translate = TRUE, magnify = TRUE)$rmat
colnames(tourne) <- colnames(res.procrustes)
colnames(rec.pca) <- colnames(res.supp)
res.procrustes <- rbind.data.frame(res.procrustes,tourne)
res.supp <- rbind.data.frame(res.supp,rec.pca)
res.var <- cbind.data.frame(res.var,rec.pca)
res.dim <- cbind.data.frame(res.dim,acpfin$ind$coord[,1:ncp])
}
####
if (!is.null(main)) title <- main
if ((choice=="all")|(choice=="ind.proc")){
if (new.plot) dev.new()
oo=PCA(res.procrustes,ind.sup=c((nrow(res$call$X)+1):nrow(res.procrustes)),scale=FALSE,graph=FALSE)
oo$eig=reference$eig
el=coord.ellipse(cbind.data.frame(as.factor(rep(rownames(res$call$X),res$call$nboot)),oo$ind.sup$coord[,axes]),level.conf=level.conf)
if (is.null(main)) title="Multiple imputation using Procrustes"
plot(oo,axes=axes,col.ind.sup=rep(1:nrow(res$call$X),res$call$nboot),label="ind",ellipse=el,col.quali="black", title=title,invisible="ind.sup",new.plot=FALSE)
# if (!is.null(add.tab)){
# vrai = PCA(add.tab,graph=FALSE,scale=res$call$scale)
# tourne <- procrustes(vrai$ind$coord[,axes], reference$ind$coord[,axes],orthogonal = TRUE, translate = TRUE, magnify = TRUE)$rmat
# points(tourne[,axes],cex=0.9,col=2)
# }
}
if ((choice=="all")|(choice=="dim")){
if (new.plot) dev.new()
colnames(res.dim)=paste("V",1:ncol(res.dim))
ooo=PCA(res.dim,quanti.sup=(ncol(res$call$X)+1):ncol(res.dim),scale=res$call$scale,graph=FALSE)
ooo$eig=reference$eig
if (is.null(main)) title="Projection of the Principal Components"
plot(ooo,choi="var",axes=axes,title=title,label="none",new.plot=FALSE,invisible="var")
}
if ((choice=="all")|(choice=="ind.supp")){
if (new.plot) dev.new()
oo=PCA(res.supp,ind.sup=c((nrow(res$call$X)+1):nrow(res.supp)),scale=res$call$scale,graph=FALSE)
el=coord.ellipse(cbind.data.frame(as.factor(rep(rownames(res$call$X),res$call$nboot)),oo$ind.sup$coord[,1:2]),level.conf = level.conf)
if (is.null(main)) title="Supplementary projection"
plot(oo,axes=axes,col.ind.sup=rep(1:nrow(res$call$X),res$call$nboot),label="ind",ellipse=el,col.quali="black",
title=title,invisible="ind.sup",new.plot=FALSE)
# if (!is.null(add.tab)){
# dele = PCA(rbind.data.frame(rec.pca,add.tab),ind.sup=c((nrow(res$call$X)+1):(2*nrow(res$call$X))),scale=scale,graph=FALSE)
# points(dele$ind.sup$coord[,axes],col=2)
# }
}
if ((choice=="all")|(choice=="var")){
if (new.plot) dev.new()
color = c("black", "red", "green3", "blue", "cyan", "magenta",
"darkgray", "darkgoldenrod", "darkgreen", "violet",
"turquoise", "orange", "lightpink", "lavender", "yellow",
"lightgreen", "lightgrey", "lightblue", "darkkhaki",
"darkmagenta", "darkolivegreen", "lightcyan", "darkorange",
"darkorchid", "darkred", "darksalmon", "darkseagreen",
"darkslateblue", "darkslategray", "darkslategrey",
"darkturquoise", "darkviolet", "lightgray", "lightsalmon",
"lightyellow", "maroon")
colnames(res.var)=paste("V",1:ncol(res.var))
colnames(res.var)[1:ncol(res$call$X)]=colnames(res$call$X)
oo=PCA(res.var,quanti.sup=c((ncol(res$call$X)+1):ncol(res.var)),scale=res$call$scale,graph=FALSE)
if (is.null(main)) title="Variable representation"
plot(oo, axes=axes, choix = "var", title=title,invisible = "quanti.sup", fill = color[1:ncol(res$call$X)],new.plot=FALSE)
for (k in 1:res$call$nboot) points(oo$quanti.sup$coord[((k-1)*ncol(res$call$X)+1):(k*ncol(res$call$X)),axes[1]], oo$quanti.sup$coord[((k-1)*ncol(res$call$X)+1):(k*ncol(res$call$X)),axes[2]], col = color[1:ncol(res$call$X)], pch = 15, cex = 0.3)
}
}
Any scripts or data that you put into this service are public.
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.