#  READING CREDSCO_BIN

load("E:/Curs MD english/R/credsco_bin")
objects()

#
# VISUALISATION OF DATA
#
# PRINCIPAL COMPONENT ANALYSIS OF CONTINUOUS VARIABLES, WITH Dictamen PROJECTED AS ILLUSTRATIVE
#

# CREATION OF THE DATA FRAME OF CONTINUOUS VARIABLES

dcon <- data.frame(Antig_feina,Plaç,Edat,Despeses,Ingressos,Patrimoni,Carrecs_pat,Import_sol,Preu_finan,Rati_fin,Estalvi)



# PRINCIPAL COMPONENT ANALYSIS OF dcon

pc1 = prcomp(dcon, scale=T)

print(pc1)

screeplot(pc1,type="lines")


# WHICH PERCENTAGE OF THE TOTAL INERTIA IS REPRESENTED IN SUBSPACES?

100*sum(pc1$sdev[1:nd]^2)/dim(dcon)[2]
barplot(100*cumsum(pc1$sdev[1:dim(dcon)[2]]^2)/dim(dcon)[2])

# SELECTION OF THE SIGIFICATIVE DIMENSIONS (LAST ELBOW TECHNIQUE) (DIFFICULT THIS TIME EH!)

nd = 7

# STORAGE OF THE EIGENVALUES, EIGENVECTORS AND PROJECTIONS IN THE nd DIMENSIONS

lbd = pc1$sdev[1:nd]^2
U = pc1$rotation[,1:nd]
Psi = pc1$x[,1:nd]

# STORAGE OF LABELS FOR INDIVIDUALS AND VARIABLES

iden = row.names(dcon)
etiq = names(dcon)
ze = rep(0,length(etiq)) # WE WILL NEED THIS VECTOR AFTERWARDS FOR THE GRAPHICS

# PLOT OF INDIVIDUALS

plot(Psi[,1],Psi[,2],type="n")
text(Psi[,1],Psi[,2],labels=iden)
axis(side=1, pos= 0, labels = F, col="cyan")
axis(side=3, pos= 0, labels = F, col="cyan")
axis(side=2, pos= 0, labels = F, col="cyan")
axis(side=4, pos= 0, labels = F, col="cyan")

# PLOT OF VARIABLES

Phi = cor(dcon,Psi)
plot(Phi[,1],Phi[,2],type="none",xlim=c(min(Phi[,1],0),max(Phi[,1],0)))
axis(side=1, pos= 0, labels = F)
axis(side=3, pos= 0, labels = F)
axis(side=2, pos= 0, labels = F)
axis(side=4, pos= 0, labels = F)
arrows(ze, ze, Phi[,1], Phi[,2], length = 0.07,col="blue")
text(Phi[,1],Phi[,2],labels=etiq,col="darkblue")


# BIPLOT DE RP

#biplot(Psi,U)

#fm = round(max(abs(Psi[,1]))) # ES UN FACTOR D'ESCALA PER DIBUIXAR LES FLETXES MES VISIBLES EN EL GRAFIC

#plot(Psi[,1],Psi[,2],col="gray")
#arrows(ze, ze, fm*U[,1], fm*U[,2], length = 0.07,col="red")
#text(fm*U[,1],fm*U[,2],labels=etiq,col="red")


# 
# INTERPRET THE TWO FIRST AXES?
#
cor(dcon,Psi)

# CALCULATE THE PROJECTIONS OF THE VARIABLES PHI

U %*% diag(sqrt(lbd))

# WRITE THE FIRST PRINCIPAL COMPONENT FUNCTION OF THE OBSERVED VARIABLES
#

# WRITE THE Ingressos AS FUNCTION OF THE SELECTED FACTORS
# EVERY PRINCIPAL COMPONENT IS A COMBINATION OF THE ORIGINAL VARIABLESS (STANDARDIZED SICE SCALE=T)
# 

Z    = scale(dcon)
reg1 = lm(Psi[,1]~ Z)
print(reg1)
U[,1]

# LIKEWISE EVERY VARIABLE CAN BE EXPLAINED BY THE PRINCIPAL COMPONENTS (LATENT FACTORS) 

#reg2 = lm(Z[,1]~ pc1$x[,1:nd])
#print(reg2)
#U[1,]

# HOWEVER TO IDENTIFY THE LATENT FACTORS IS BETTER TO PERFORM A ROTATION OF THE DIRECTIONS TO FACILITATE THIR INTERPRETABILITY
# RUN VARIMAX ROTATION AND LOOK FOR LOADING (=CORRELATION) > 0.7

#pcrot = varimax(Phi)
#print(pcrot)

# PROJECTION OF dICATAMEN AS ILLUSTRATIVE

# PROJECCIÓ OF INDIVIDUALS DIFFERENTIATING THE Dictamen
# (we need a numeric Dictamen to color)

idict=dd[,1]
plot(Psi[,1],Psi[,2],col=idict)
axis(side=1, pos= 0, labels = F, col="darkgray")
axis(side=3, pos= 0, labels = F, col="darkgray")
axis(side=2, pos= 0, labels = F, col="darkgray")
axis(side=4, pos= 0, labels = F, col="darkgray")
legend("topright",c("pos","neg"),pch=1,col=c(1,2))

# NOW WE PROJECT THE CDG OF BOTH LEVELS OF Dictament

fdic1 = tapply(Psi[,1],Dictamen,mean)
fdic2 = tapply(Psi[,2],Dictamen,mean) 

points(fdic1,fdic2,pch=16,col="yellow")
text(fdic1,fdic2,labels=levels(Dictamen),col="yellow")


# 
# LETS TAKE ALL THE INFORMATION AVAILABLE: ALL CATEGORICAL VARIABLES 
# WE WILL REPEAT THE SAME PCA ANALYSIS BUT WITH MULTIPLE CORRESPONDENCE ANALYSIS

# CREATION OF THE DATA FRAME OF CATEGORICAL VARIABLES

dcat <- data.frame(edatR,antigR,plaçR,despesesR,ingressosR,patrimoniR,carrecsR,importR,preuR,ratfinR,estalviR,Vivenda,Estat_civil,Registres,Tipus_feina)

# CREATION OF THE DATA FRAME WITH THE ILLUSTRATIVE VARIABLE

dict <- data.frame(Dictamen)


# LOADING ACM FUNTION

source("C:/Users/aluja/Documents/txt/docent/MD/R/acm.r") # LA TROBAREU A LA PAGINA WEB DE L'ASSIGNATURA

ac1 <- acm(dcat,dict)

attributes(ac1)

# HISTOGRAM OF EIGENVALUES OF THE MCA (WHY ARE SO DIFFERENT?)

barplot(ac1$vaps)

# WHAT IS THE SIGNIFICANT DIMENSION IN MCA?

i <- 1

while (ac1$vaps[i] > 1/ncol(dcat)) i = i+1

nd = i-1

# PLOT OF INDIVIDUALS (USELESS)

plot(ac1$rs[,1],ac1$rs[,2],type="n")
text(ac1$rs[,1],ac1$rs[,2],labels=row.names(ac1$rs))
axis(side=1, pos= 0, labels = F, col="cyan")
axis(side=3, pos= 0, labels = F, col="cyan")
axis(side=2, pos= 0, labels = F, col="cyan")
axis(side=4, pos= 0, labels = F, col="cyan")

# PLOT OF MODALITIES (LEVELS OF CATEGORICAL VARIABLES) 

ivar = gl(11,5)
ivar = c(ivar,rep(12,6),rep(13,5),rep(14,2),rep(15,4)) # ivar IDENTIFIES TO WHICH CATEGORICAL VARIABLE BELONGS EACH MODALITY, TO BE USED FOR COLOR PURPOSES)

plot(ac1$cs[,1],ac1$cs[,2],type="none")
axis(side=1, pos= 0, labels = F, col="gray")
axis(side=3, pos= 0, labels = F, col="gray")
axis(side=2, pos= 0, labels = F, col="gray")
axis(side=4, pos= 0, labels = F, col="gray")
text(ac1$cs[,1],ac1$cs[,2],labels=row.names(ac1$cs),col=ivar,cex=0.6)


# PLOT OF THE LEVELS OF Dictamen

text(ac1$csup[,1],ac1$csup[,2],labels=row.names(ac1$csup))

# PLOT of INDIVIDUALS ACCODING THE LEVEL OF dICTAMEN

plot(ac1$rs[,1],ac1$rs[,2],col=idict)