# R/plot.fem.R In FisherEM: The FisherEM Algorithm to Simultaneously Cluster and Visualize High-Dimensional Data

#### Documented in plot.fem

```plot.fem <- function(x,frame=0,crit=c(),...){
# Color palette
palette(c("#A6CEE3", "#1F78B4", "#B2DF8A", "#33A02C", "#FB9A99", "#E31A1C",
"#FDBF6F", "#FF7F00", "#CAB2D6", "#6A3D9A", "#FFFF99", "#B15928"))
Y = eval.parent(x\$call[[2]])

# Model selection
if ((frame==0 || frame==1) & (length(x\$allCriteria\$K)>1 || length(x\$allCriteria\$models)>1)){
K = x\$allCriteria\$K; models = x\$allCriteria\$model
if (length(crit)==0) crit = x\$crit
if (crit=='bic') val = x\$allCriteria\$bic; if (crit=='aic') val = x\$allCriteria\$aic; if (crit=='icl') val = x\$allCriteria\$icl
val[val==-Inf] = NA
if (length(models)>1){
plot(K,val,type='n',xlab='K',ylab=crit,main='Model selection')
for (m in 1:length(levels(models))){
modl = levels(models)[m]
lines(K[models==modl],val[models==modl],col=m,lty=m,lwd=2)
}
legend('bottomleft',levels(models),col=1:length(models),lty=1:length(models),
ncol=round(length(models)/3),cex=1)
}
else plot(K,val,type='b',xlab='K',ylab=crit,col=1:length(models),lty=1:length(models),main='Selection of the number of groups',...)
if (x\$call[[1]]=='sfem'){
if (crit=='bic') val = x\$allCriteria\$l1\$bic; if (crit=='aic') val = x\$allCriteria\$l1\$aic; if (crit=='icl') val = x\$allCriteria\$l1\$icl
plot(x\$allCriteria\$l1\$l1,val,type='b',xlab='l1 value',ylab=crit,main='Selection of the sparsity penalty',...)
}
}

# Log-likelihood
if (frame==0) Sys.sleep(0.5)
if (frame==0 || frame==2) if (x\$call[[1]]!='sfem') plot(x\$loglik.all,type='b',xlab='Iterations',
ylab='Log-likelihood',main='Log-likelihood',col=2,pch=20,cex=0.5,...)

# Discriminative subspace
if (frame==0) Sys.sleep(0.5)
if (frame==0 || frame==3){
Y = as.matrix(Y)
p = ncol(Y)
n = nrow(Y)
if (ncol(x\$U)>1){
xx = as.matrix(Y)%*%x\$U[,1:2]
cls = x\$cls
min1= round(min(xx[,1]),1)-1
max1= round(max(xx[,1],1))+1
min2= round(min(xx[,2]),1)-1
max2= round(max(xx[,2]),1)+1
topX = topY = 0
xhist = yhist = list()
for (k in 1:max(cls)){
xhist[[k]] = hist(xx[cls==k,1],breaks=seq(min1,max1,0.1), plot=FALSE)
yhist[[k]] = hist(xx[cls==k,2],breaks=seq(min2,max2,0.1), plot=FALSE)
topX = max(topX,xhist[[k]]\$counts); topY = max(topY,yhist[[k]]\$counts)
}
nf <- layout(matrix(c(2,0,1,3),2,2,byrow=TRUE), c(3,1), c(1,3), TRUE)
xrange <- c(min1,max1); yrange <- c(min2,max2)

par(mar=c(3,3,1,1))
plot(xx,col=cls,pch=cls,xlim=xrange, ylim=yrange,cex=1.5)
par(mar=c(0,3,1,1))
par(mar=c(3,0,1,1))
par(old.par)
}
else {
cat('Since K=2, the data and the discriminative subspace have been projected on the 2 first PCs','\n')
if (ncol(Y)<nrow(Y))  Z = -eigen(cov(Y),symmetric=T)\$vectors[,1:2]
else {
z = eigen(Y%*%t(Y),symmetric=T)
Z = matrix(NA,p,2)
for (i in 1:2) Z[,i] = t(Y)%*%z\$vectors[,i]/sqrt(n*z\$values[i])
}
MU	= colMeans(Y)
proj= matrix(NA,2,p)
axU = matrix(x\$U,p,1)%*%matrix(x\$U,1,p)%*%matrix(10, p, 1)
proj= axU+matrix(MU,p,1)
Yproj = Y%*%Z
u 	= matrix(proj,1,p)%*%Z # projection sur les 2 pc
ybar 	= matrix(MU,1,p) %*% Z # proj des moyennes des cls sur les 2 pc
plot(Yproj,col=x\$cls+1,xlab='comp.1',ylab='comp.2',pch=x\$cls+1)
pente=(u[1,2]-ybar[1,2])/(u[1,1]-ybar[1,1])
oo=u[1,2]-pente*u[1,1]
xb=(2*ybar[1,1]-sqrt(50^2/(pente^2+1)))/2
xa=(2*ybar[1,1]+sqrt(50^2/(pente^2+1)))/2
#cat(xa,xb,oo+pente*c(xa,xb),'\n')
lines(c(xa,xb),oo+pente*c(xa,xb),col=1,type='l',lwd=2)
title(main='Estimated discriminative subspace (projected onto PCA axes)')
}
par(old.par)
}

# Plot of the group means
if (frame==0) Sys.sleep(0.5)
if (frame==0 || frame==4){
matplot(t(x\$mean),type='l',lwd=1.5,xaxt='n',ylab='')
axis(1,at=1:ncol(Y),labels=colnames(Y),las=2)
title(main='Group means')
legend('bottomleft',paste('Group ',1:max(x\$cls),sep=''),col=1:max(x\$cls),lty=1:max(x\$cls),
ncol=round(max(cls)/3),cex=0.8,...)
}