plot.MANOVA <- function(x, A1 = 1, A2 = 2, ScaleGraph = TRUE, PlotGroups = TRUE, PlotVars = TRUE, PlotInd = TRUE, WhatInds=NULL,
WhatVars=NULL, WhatGroups=NULL, IndLabels=NULL, VarLabels=NULL, GroupLabels=NULL, AbbreviateLabels=FALSE,
LabelInd=TRUE, LabelVars = TRUE, CexGroup=1, PchGroup=16, margin=0.1, AddLegend=FALSE, ShowAxes=FALSE, LabelAxes=FALSE, LabelGroups=TRUE,
PlotCircle = TRUE, ConvexHulls = FALSE, TypeCircle = "M", ColorGroups = NULL, ColorVars = NULL, LegendPos="topright",
ColorInd = NULL, voronoi = TRUE, mode="a", TypeScale = "Complete", ValuesScale = "Original",
MinQualityVars = 0, dpg = 0, dpi=0, dp=0, PredPoints=0, PlotAxis = FALSE, CexInd = NULL, CexVar = NULL,
PchInd = NULL, PchVar = NULL, ColorVar=NULL, ShowAxis=FALSE, VoronoiColor="black", ShowBox=FALSE, ShowTitle=TRUE,
PlotClus = FALSE, TypeClus = "ch", ClustConf = 1, ClustCenters = FALSE, UseClusterColors = TRUE, CexClustCenters=1, ...) {
modes=c("p", "a", "b", "h", "ah", "s")
if (is.numeric(mode))
mode = modes[mode]
TypeScales=c("Complete", "StdDev", "BoxPlot")
if (is.numeric(TypeScale))
TypeScale = TypeScales[TypeScale]
ValuesScales=c("Original", "Transformed")
if (is.numeric(ValuesScale))
ValuesScale = ValuesScales[ValuesScale]
if (is.null(ColorGroups))
ColorGroups = x$ClusterColors
if (is.null(ColorInd))
ColorInd = x$ClusterColors[as.integer(x$groups)]
if (LabelAxes){
xlabel = paste("CV", A1, "(", round(x$Inertia[A1], digits=2),"%)")
ylabel = paste("CV", A2, "(", round(x$Inertia[A2], digits=2),"%)")
Title="Canonical/MANOVA Biplot"}
else{
xlabel = ""
ylabel = ""
Title=paste("Canonical Biplot", " / ", A1, "-", A2, " (", round(x$Inertia[A1]+x$Inertia[A2], digits=2),"%)")
}
scf = 1
if (ScaleGraph) {
sca = sum(x$GroupCoordinates^2)
scb = sum(x$ColCoordinates^2)
sca = sca/x$g
scb = scb/x$p
scf = sqrt(sqrt(scb/sca))
}
J = x$GroupCoordinates[, c(A1, A2)] * scf
H = x$ColCoordinates[, c(A1, A2)]/scf
V = x$RowCoordinates[, c(A1, A2)] * scf
g = dim(J)[1]
p = dim(H)[1]
n = dim(V)[1]
if (is.null(GroupLabels)) GroupLabels=rownames(J)
if (is.null(IndLabels)) IndLabels=rownames(V)
if (is.null(VarLabels)) VarLabels=rownames(H)
if (AbbreviateLabels){
GroupLabels=abbreviate(GroupLabels, minlength = 5L)
IndLabels=abbreviate(IndLabels, minlength = 5L)
VarLabels=abbreviate(VarLabels, minlength = 5L)
}
# Determining what rows to plot
if (is.null(WhatInds))
WhatInds = matrix(1, n, 1)
else
if (!CheckBinaryVector(WhatInds)) {
AllRows = matrix(0, n, 1)
AllRows[WhatInds]=1
WhatInds=AllRows
}
WhatInds=as.logical(WhatInds)
# Determining what columns to plot
if (is.null(WhatVars))
WhatVars = matrix(1, p, 1)
else
if (!CheckBinaryVector(WhatVars)){
AllCols = matrix(0, p, 1)
AllCols[WhatVars]=1
WhatVars=AllCols
}
WhatVars=as.logical(WhatVars)
# Determining what groups to plot
if (is.null(WhatGroups))
WhatGroups = matrix(1, p, 1)
else
if (!CheckBinaryVector(WhatGroups)){
AllCols = matrix(0, p, 1)
AllCols[WhatGroups]=1
WhatGroups=AllCols
}
WhatVars=as.logical(WhatVars)
if (is.null(CexInd))
CexInd = 0.8
if (is.null(CexVar))
CexVar = rep(0.8, x$p)
else if (length(CexVar == 1))
CexVar = rep(CexVar, p)
if (is.null(PchInd))
PchInd = 1
if (is.null(PchVar))
PchVar = rep(16, p)
if (is.null(ColorVar))
ColorVar = rep("black", p)
if (!LabelVars) VarLabels=rep(" ", p)
if (!is.null(x$Sup_Individual_Coord)){VS=x$Sup_Individual_Coord*scf}
if (ShowAxis) {
xaxt = "s"
yaxt = "s"
} else {
xaxt = "n"
yaxt = "n"
}
if (PlotInd | PlotClus)
PP = rbind(J, H, V)
else
PP = rbind(J, H)
xmin = min(PP[, 1])
xmax = max(PP[, 1])
ymin = min(PP[, 2])
ymax = max(PP[, 2])
xrang=abs(xmax-xmin)
yrang=abs(ymax-ymin)
fact=0.03
xmin=xmin-fact*xrang
xmax=xmax+fact*xrang
ymin=ymin-fact*yrang
ymax=ymax+fact*yrang
if (xmax <0 ) xmax=xmax*(-1)
PP = rbind(PP, c(xmin - (xmax - xmin) * margin, ymin - (ymax - ymin) * margin))
PP = rbind(PP, c(xmax + (xmax - xmin) * margin, ymax + (ymax - ymin) * margin))
plot(PP[, 1], PP[, 2], cex = 0, asp = 1, main = Title, xlab = xlabel, ylab = ylabel, xaxt = xaxt, yaxt = yaxt, axes=ShowAxes, ...)
#op=par(mai=c(0,0,0.5,0))
#op=par(mar=c(1, 1, 1, 1) + 0.1)
if (ShowBox) rect(xmin, ymin, xmax, ymax)
if (PlotClus) {
ColorInd=PlotBiplotClusters(V, x$Clusters, TypeClus = TypeClus, ClusterColors = x$ClusterColors, ClusterNames=x$ClusterNames, centers = ClustCenters, ClustConf=ClustConf, CexClustCenters=CexClustCenters, ...)
if (x$ClusterType=="gm"){
ColorInd2=rgb((x$P %*% t(col2rgb(x$ClusterColors)))/255)
if (UseClusterColors) ColorInd = ColorInd2
PchInd=rep(16,n)
}
}
if (PlotInd) {
points(V[, 1], V[, 2], col = ColorInd, pch=PchInd, cex=CexInd)
if (LabelInd)
text(V[, 1], V[, 2], IndLabels, adj = 0, col = ColorInd, cex = 0.8)
if (!is.null(x$Sup_Individual_Coord)){
points(VS[, A1], VS[, A2], cex = 0.5, col = "Black")
text(VS[, A1], VS[, A2], paste(rownames(VS),"*",sep=""), adj = 0, col = "Black", cex = 0.8)
}
}
if (PlotGroups) {
points(J[, 1], J[, 2], cex = CexGroup, pch = PchGroup, col = ColorGroups)
if (LabelGroups)
text(J[, 1], J[, 2], GroupLabels, adj = 0, col = ColorGroups, cex = 1.5)
}
if (AddLegend) legend(LegendPos, rownames(J), fill = ColorGroups)
if (PlotVars) {
qlrcols = x$ColContributions[, A1] + x$ColContributions[, A2]
if (mode=="s")
Scales = GetBiplotScales(x, TypeScale = TypeScale, ValuesScale = ValuesScale)
for (j in 1:p) if (qlrcols[j] > MinQualityVars)
VarBiplot(H[j, 1], H[j, 2], xmin = xmin, xmax = xmax, ymin = ymin, ymax = ymax, label = VarLabels[j], mode = mode, CexPoint = CexVar[j], Color = ColorVar[j],
ticks = Scales$Ticks[[j]], ticklabels = Scales$Labels[[j]], ts = TypeScale, PchPoint = PchVar[j], ...)
for (idp in dpg)
if ((idp > 0) & (idp < (p + 1))) {
g = H[idp, ]
nn = (t(g) %*% g)
scal <- (J %*% g)/nn[1, 1]
Dscal <- diag(as.vector(scal))
Fpr <- Dscal %*% matrix(rep(1, nrow(J)), ncol = 1) %*% t(g)
nrFpr <- nrow(Fpr)
dlines(J, Fpr, color=ColorGroups)
}
if (PlotInd){
for (idp in dp)
if ((idp > 0) & (idp < (p + 1))) {
g = H[idp, ]
nn = (t(g) %*% g)
scal <- (V %*% g)/nn[1, 1]
Dscal <- diag(as.vector(scal))
Fpr <- Dscal %*% matrix(rep(1, nrow(V)), ncol = 1) %*% t(g)
nrFpr <- nrow(Fpr)
dlines(V, Fpr, color=ColorInd)
}
}
for (idp in PredPoints)
if ((idp > 0) & (idp < (n + 1)))
for (j in 1:p){
g = H[j, ]
nn = (t(g) %*% g)
scal <- (J[idp,] %*% g)/nn[1, 1]
Fpr <- scal %*% t(g)
nrFpr <- nrow(Fpr)
dlines(matrix(J[idp,],1,2) , Fpr, color=ColorVar[j])
}
}
if (PlotCircle) {
switch(TypeCircle, U = {
radius = x$UnivRad
}, B = {
radius = x$BonfRad
}, M = {
radius = x$MultRad
}, C = {
radius = x$ChisRad
})
radius=radius*scf;
for (i in 1:x$g) {
Circle(radius[i], c(J[i, 1], J[i, 2]), col = ColorGroups[i])}
}
if (ConvexHulls) {
lev = levels(x$groups)
for (i in 1:nlevels(x$groups)) {
XP = V[which(x$groups == lev[i]), ]
XP = cbind(XP[, 1], XP[, 2])
hpts <- chull(XP)
hpts <- c(hpts, hpts[1])
lines(XP[hpts, ], col = ColorGroups[i])
}
}
if (voronoi) {
tv = deldir(J[, 1], J[, 2], rw = c(xmin, xmax, ymin, ymax))
plot(tv, add = TRUE, wlines = "tess", xlim = c(xmin, xmax), ylim = c(ymin, ymax),col=VoronoiColor)
}
if (!is.null(x$ContSupVarsBiplot))
plot(x$ContSupVarsBiplot, F1=A1, F2=A2, xmin = xmin, xmax = xmax, ymin = ymin, ymax = ymax, mode=mode, TypeScale=TypeScale)
#par(op)
}
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.