2c74554899dfbeca72187047b7f6f1f393ae3fda
[valse.git] / pkg / R / plot_valse.R
1 #' Plot
2 #'
3 #' It is a function which plots relevant parameters
4 #'
5 #' @param model the model constructed by valse procedure
6 #' @param n sample size
7 #' @return several plots
8 #'
9 #' @examples TODO
10 #'
11 #' @export
12 #'
13 plot_valse = function(model,n){
14 require("gridExtra")
15 require("ggplot2")
16 require("reshape2")
17 require("cowplot")
18
19 K = length(model$pi)
20 ## regression matrices
21 gReg = list()
22 for (r in 1:K){
23 Melt = melt(t((model$phi[,,r])))
24 gReg[[r]] = ggplot(data = Melt, aes(x=Var1, y=Var2, fill=value)) + geom_tile() +
25 scale_fill_gradient2(low = "blue", high = "red", mid = "white", midpoint = 0, space = "Lab") +
26 ggtitle(paste("Regression matrices in cluster",r))
27 }
28 print(gReg)
29
30 ## Differences between two clusters
31 k1 = 1
32 k2 = 2
33 Melt = melt(t(model$phi[,,k1]-model$phi[,,k2]))
34 gDiff = ggplot(data = Melt, aes(x=Var1, y=Var2, fill=value)) + geom_tile() +
35 scale_fill_gradient2(low = "blue", high = "red", mid = "white", midpoint = 0, space = "Lab") +
36 ggtitle(paste("Difference between regression matrices in cluster",k1, "and", k2))
37 print(gDiff)
38
39 ### Covariance matrices
40 matCov = matrix(NA, nrow = dim(model$rho[,,1])[1], ncol = K)
41 for (r in 1:K){
42 matCov[,r] = diag(model$rho[,,r])
43 }
44 MeltCov = melt(matCov)
45 gCov = ggplot(data =MeltCov, aes(x=Var1, y=Var2, fill=value)) + geom_tile() +
46 scale_fill_gradient2(low = "blue", high = "red", mid = "white", midpoint = 0, space = "Lab") +
47 ggtitle("Covariance matrices")
48 print(gCov )
49
50 ### proportions
51 Gam = matrix(0, ncol = K, nrow = n)
52 gam = Gam
53 for (i in 1:n){
54 for (r in 1:K){
55 sqNorm2 = sum( (Y[i,]%*%model$rho[,,r]-X[i,]%*%model$phi[,,r])^2 )
56 Gam[i,r] = model$pi[r] * exp(-0.5*sqNorm2)* det(model$rho[,,r])
57 }
58 gam[i,] = Gam[i,] / sum(Gam[i,])
59 }
60 affec = apply(gam, 1,which.max)
61 gam2 = matrix(NA, ncol = K, nrow = n)
62 for (i in 1:n){
63 gam2[i, ] = c(gam[i, affec[i]], affec[i])
64 }
65 bp <- ggplot(data.frame(gam2), aes(x=X2, y=X1, color=X2, group = X2)) +
66 geom_boxplot() + theme(legend.position = "none")+ background_grid(major = "xy", minor = "none")
67 print(bp )
68
69 ### Mean in each cluster
70 XY = cbind(X,Y)
71 XY_class= list()
72 meanPerClass= matrix(0, ncol = K, nrow = dim(XY)[2])
73 for (r in 1:K){
74 XY_class[[r]] = XY[affec == r, ]
75 meanPerClass[,r] = apply(XY_class[[r]], 2, mean)
76 }
77 data = data.frame(mean = as.vector(meanPerClass), cluster = as.character(rep(1:K, each = dim(XY)[2])), time = rep(1:dim(XY)[2],K))
78 g = ggplot(data, aes(x=time, y = mean, group = cluster, color = cluster))
79 print(g + geom_line(aes(linetype=cluster, color=cluster))+ geom_point(aes(color=cluster)) + ggtitle('Mean per cluster'))
80
81 }