fix initialization and made some update

diff --git a/pkg/R/computeGridLambda.R b/pkg/R/computeGridLambda.R
index 8ec4d66..9d06aed 100644 (file)
--- a/pkg/R/computeGridLambda.R
+++ b/pkg/R/computeGridLambda.R
@@ -33,5 +33,5 @@ computeGridLambda = function(phiInit, rhoInit, piInit, gamInit, X, Y,
                        grid[i,j,] = abs(list_EMG$S[i,j,]) / (n*list_EMG$pi^gamma)
        }
        grid = unique(grid)
                        grid[i,j,] = abs(list_EMG$S[i,j,]) / (n*list_EMG$pi^gamma)
        }
        grid = unique(grid)

-       grid
+       sort(grid)

 }
 }

diff --git a/pkg/R/initSmallEM.R b/pkg/R/initSmallEM.R
index 5dcafb8..bd5ce17 100644 (file)
--- a/pkg/R/initSmallEM.R
+++ b/pkg/R/initSmallEM.R
@@ -13,21 +13,21 @@ initSmallEM = function(k,X,Y, fast=TRUE)
        n = nrow(Y)
        m = ncol(Y)
        p = ncol(X)
        n = nrow(Y)
        m = ncol(Y)
        p = ncol(X)

-  
-       Zinit1 = array(0, dim=c(n,20))
-       betaInit1 = array(0, dim=c(p,m,k,20))
-       sigmaInit1 = array(0, dim = c(m,m,k,20))
-       phiInit1 = array(0, dim = c(p,m,k,20))
-       rhoInit1 = array(0, dim = c(m,m,k,20))
+  nIte = 20
+       Zinit1 = array(0, dim=c(n,nIte))
+       betaInit1 = array(0, dim=c(p,m,k,nIte))
+       sigmaInit1 = array(0, dim = c(m,m,k,nIte))
+       phiInit1 = array(0, dim = c(p,m,k,nIte))
+       rhoInit1 = array(0, dim = c(m,m,k,nIte))

        Gam = matrix(0, n, k)
        Gam = matrix(0, n, k)

-       piInit1 = matrix(0,20,k)
-       gamInit1 = array(0, dim=c(n,k,20))
+       piInit1 = matrix(0,nIte,k)
+       gamInit1 = array(0, dim=c(n,k,nIte))

        LLFinit1 = list()
 
        #require(MASS) #Moore-Penrose generalized inverse of matrix
        LLFinit1 = list()
 
        #require(MASS) #Moore-Penrose generalized inverse of matrix

-       for(repet in 1:20)
+       for(repet in 1:nIte)

        {
        {

-         distance_clus = dist(X)
+         distance_clus = dist(cbind(X,Y))

          tree_hier = hclust(distance_clus)
          Zinit1[,repet] = cutree(tree_hier, k)
 
          tree_hier = hclust(distance_clus)
          Zinit1[,repet] = cutree(tree_hier, k)
 


@@ -62,13 +62,12 @@ initSmallEM = function(k,X,Y, fast=TRUE)
                miniInit = 10
                maxiInit = 11
                
                miniInit = 10
                maxiInit = 11
                

-               new_EMG = EMGLLF(phiInit1[,,,repet], rhoInit1[,,,repet], piInit1[repet,],
+               init_EMG = EMGLLF(phiInit1[,,,repet], rhoInit1[,,,repet], piInit1[repet,],

                        gamInit1[,,repet], miniInit, maxiInit, gamma=1, lambda=0, X, Y, eps=1e-4, fast)
                        gamInit1[,,repet], miniInit, maxiInit, gamma=1, lambda=0, X, Y, eps=1e-4, fast)

-               LLFEessai = new_EMG$LLF
+               LLFEessai = init_EMG$LLF

                LLFinit1[repet] = LLFEessai[length(LLFEessai)]
        }
                LLFinit1[repet] = LLFEessai[length(LLFEessai)]
        }

-
-       b = which.max(LLFinit1)
+       b = which.min(LLFinit1)

        phiInit = phiInit1[,,,b]
        rhoInit = rhoInit1[,,,b]
        piInit = piInit1[b,]
        phiInit = phiInit1[,,,b]
        rhoInit = rhoInit1[,,,b]
        piInit = piInit1[b,]

diff --git a/pkg/R/main.R b/pkg/R/main.R
index 634c273..238160c 100644 (file)
--- a/pkg/R/main.R
+++ b/pkg/R/main.R
@@ -123,36 +123,39 @@ valse = function(X, Y, procedure='LassoMLE', selecMod='DDSE', gamma=1, mini=10,
   
   print(tableauRecap)
   tableauRecap = tableauRecap[which(tableauRecap[,4]!= Inf),]
   
   print(tableauRecap)
   tableauRecap = tableauRecap[which(tableauRecap[,4]!= Inf),]

-  modSel = capushe::capushe(tableauRecap, n)
-  indModSel <-
-    if (selecMod == 'DDSE')
-      as.numeric(modSel@DDSE@model)
-  else if (selecMod == 'Djump')
-    as.numeric(modSel@Djump@model)
-  else if (selecMod == 'BIC')
-    modSel@BIC_capushe$model
-  else if (selecMod == 'AIC')
-    modSel@AIC_capushe$model

   
   

-  mod = as.character(tableauRecap[indModSel,1])
-  listMod = as.integer(unlist(strsplit(mod, "[.]")))
-  modelSel = models_list[[listMod[1]]][[listMod[2]]]
+  return(tableauRecap)

   
   

-  ##Affectations
-  Gam = matrix(0, ncol = length(modelSel$pi), nrow = n)
-  for (i in 1:n){
-    for (r in 1:length(modelSel$pi)){
-      sqNorm2 = sum( (Y[i,]%*%modelSel$rho[,,r]-X[i,]%*%modelSel$phi[,,r])^2 )
-      Gam[i,r] = modelSel$pi[r] * exp(-0.5*sqNorm2)* det(modelSel$rho[,,r])
-    }
-  }
-  Gam = Gam/rowSums(Gam)
-  modelSel$affec = apply(Gam, 1,which.max)
-  modelSel$proba = Gam
-  
-  if (plot){
-    print(plot_valse(X,Y,modelSel,n))
-  }
-  
-  return(modelSel)
+  # modSel = capushe::capushe(tableauRecap, n)
+  # indModSel <-
+  #   if (selecMod == 'DDSE')
+  #     as.numeric(modSel@DDSE@model)
+  # else if (selecMod == 'Djump')
+  #   as.numeric(modSel@Djump@model)
+  # else if (selecMod == 'BIC')
+  #   modSel@BIC_capushe$model
+  # else if (selecMod == 'AIC')
+  #   modSel@AIC_capushe$model
+  # 
+  # mod = as.character(tableauRecap[indModSel,1])
+  # listMod = as.integer(unlist(strsplit(mod, "[.]")))
+  # modelSel = models_list[[listMod[1]]][[listMod[2]]]
+  # 
+  # ##Affectations
+  # Gam = matrix(0, ncol = length(modelSel$pi), nrow = n)
+  # for (i in 1:n){
+  #   for (r in 1:length(modelSel$pi)){
+  #     sqNorm2 = sum( (Y[i,]%*%modelSel$rho[,,r]-X[i,]%*%modelSel$phi[,,r])^2 )
+  #     Gam[i,r] = modelSel$pi[r] * exp(-0.5*sqNorm2)* det(modelSel$rho[,,r])
+  #   }
+  # }
+  # Gam = Gam/rowSums(Gam)
+  # modelSel$affec = apply(Gam, 1,which.max)
+  # modelSel$proba = Gam
+  # 
+  # if (plot){
+  #   print(plot_valse(X,Y,modelSel,n))
+  # }
+  # 
+  # return(modelSel)

 }
 }