Commit | Line | Data |
---|---|---|
3453829e BA |
1 | ### Regression matrices |
2 | model = Res | |
3 | K = dim(model$phi)[3] | |
4 | valMax = max(abs(model$phi)) | |
5 | ||
6 | require(fields) | |
7 | ||
8 | if (K<4){ | |
9 | par(mfrow = c(1,K)) | |
10 | } else op = par(mfrow = c(2, (K+1)/2)) | |
11 | ||
12 | ## Phi | |
13 | ||
14 | for (r in 1:K){ | |
15 | image.plot(t(abs(model$phi[,,r])), | |
16 | col=gray(rev(seq(0,64,length.out=65))/65),breaks=seq(0,valMax,length.out=66)) | |
17 | } | |
18 | par(mfrow = c(1,K),oma = c(0,0,3,0)) | |
19 | mtext("Regression matrices in each cluster", side=3, line=4, font=2, cex=2, col='red') | |
20 | ||
21 | par(mfrow = c(1,2), oma=c(0,0,3,0)) | |
22 | for (i in 1:4) | |
23 | plot(runif(20), runif(20), | |
24 | main=paste("random plot (",i,")",sep='')) | |
25 | par(op) | |
26 | mtext("Four plots", | |
27 | side=3, line=4, font=2, cex=2, col='red') | |
28 | ||
29 | ### Zoom onto two classes we want to compare | |
30 | kSel = c(1,2) | |
31 | par(mfrow = c(1,3)) | |
32 | ||
33 | for (r in kSel){ | |
34 | image.plot(t(abs(model$phi[,,r])),xaxt="n",yaxt="n", | |
35 | col=gray(rev(seq(0,64,length.out=65))/65),breaks=seq(0,valMax,length.out=66)) | |
36 | } | |
37 | image.plot(t(abs(model$phi[,,kSel[1]]-model$phi[,,kSel[2]])), | |
38 | col=gray(rev(seq(0,64,length.out=65))/65),breaks=seq(0,valMax,length.out=66)) | |
39 | ||
40 | ### Covariance matrices | |
41 | par(mfrow = c(K, 1)) | |
42 | for (r in 1:K){ | |
43 | image.plot(matrix(diag(model$rho[,,r]), ncol= 1), | |
44 | col=gray(rev(seq(0,64,length.out=65))/65),breaks=seq(0,valMax,length.out=66)) | |
45 | } | |
46 | ||
47 | ### proportions | |
48 | Gam = matrix(0, ncol = K, nrow = n) | |
49 | gam = Gam | |
50 | for (i in 1:n){ | |
51 | for (r in 1:K){ | |
52 | sqNorm2 = sum( (Y[i,]%*%model$rho[,,r]-X[i,]%*%model$phi[,,r])^2 ) | |
53 | Gam[i,r] = model$pi[r] * exp(-0.5*sqNorm2)* det(model$rho[,,r]) | |
54 | } | |
55 | gam[i,] = Gam[i,] / sum(Gam[i,]) | |
56 | } | |
57 | affec = apply(gam, 1,which.max) | |
58 | gam2 = matrix(NA, ncol = K, nrow = n) | |
59 | for (i in 1:n){ | |
60 | gam2[i, affec[i]] = gam[i, affec[i]] | |
61 | } | |
62 | boxplot(gam2) | |
63 | ||
64 | ### Mean in each cluster | |
65 | XY = cbind(X,Y) | |
66 | XY_class= list() | |
67 | meanPerClass= matrix(0, ncol = K, nrow = dim(XY)[2]) | |
68 | for (r in 1:K){ | |
69 | XY_class[[r]] = XY[affec == r, ] | |
70 | meanPerClass[,r] = apply(XY_class[[r]], 2, mean) | |
71 | } | |
72 | ||
73 | matplot(meanPerClass, type='l') |