started to look at EMGLLF.c
[valse.git] / reports / simulData_17mars.R
... / ...
CommitLineData
1simulData_17mars = function(ite){
2 set.seed = 22021989+ite
3
4 ###########
5 ## Modele
6 ###########
7 K = 2
8 p = 48
9 T = seq(0,1.5,length.out = p)
10 T2 = seq(0,3, length.out = 2*p)
11 n = 100
12 x1 = cos(2*base::pi*T) + 0.2*cos(4*2*base::pi*T) + 0.3*c(rep(0,round(length(T)/7)),rep(1,round(length(T)*(1-1/7))))+1
13 sigmaX = 0.12
14 sigmaY = 0.12
15 beta = list()
16 p1= 0.5
17 beta[[1]] =diag(c(rep(p1,5),rep(1,5), rep(p1,5), rep(1, p-15)))
18 p2 = 2
19 beta[[2]] = diag(c(rep(p2,5),rep(1,5), rep(p2,5), rep(1, p-15)))
20 ARI1 = ARI2 = ARI3 = 0
21
22 ###########
23 ## Data + Projection
24 ###########
25 require(wavelets)
26 XY = array(0, dim = c(2*p,n))
27 XYproj = array(0, dim=c(96,n))
28 x = x1 + matrix(rnorm(n*p, 0, sigmaX), ncol = n)
29 affec = sample(c(1,2), n, replace = TRUE)
30 y = x
31 xy = matrix(0,ncol=n, nrow= 2*p)
32 for (i in c(1:n)){
33 y[,i] = x[,i] %*% beta[[affec[i]]] + rnorm(p, 0, sigmaY)
34 xy[,i] = c(x[,i],y[,i])
35 XY[,i] = xy[,i] - mean(xy[,i])
36 Dx = dwt(x[,i], filter='haar')@W
37 Dx = rev(unlist(Dx))
38 Dx = Dx[2:(1+3+6+12+24)]
39 Ax = dwt(x[,i], filter='haar')@V
40 Ax = rev(unlist(Ax))
41 Ax = Ax[2:(1+3)]
42 Dy = dwt(y[,i], filter='haar')@W
43 Dy = rev(unlist(Dy))
44 Dy = Dy[2:(1+3+6+12+24)]
45 Ay = dwt(y[,i], filter='haar')@V
46 Ay = rev(unlist(Ay))
47 Ay = Ay[2:(1+3)]
48 XYproj[,i] = c(Ax,Dx,Ay,Dy)
49 }
50
51 res_valse = valse(t(x),t(y), kmax=2, verbose=TRUE, plot=FALSE, size_coll_mod = 1000)
52 res_valse_proj = valse(t(XYproj[1:p,]),t(XYproj[(p+1):(2*p),]), kmax=2, verbose=TRUE, plot=FALSE, size_coll_mod = 1000)
53
54 save(res_valse,file=paste("Res_",ite, ".RData",sep=""))
55 save(res_valse_proj,file=paste("ResProj_",ite, ".RData",sep=""))
56}