pkg/R/selectVariables.R

   1 #' selectVariables
   2 #'
   3 #' It is a function which construct, for a given lambda, the sets of relevant variables.
   4 #'
   5 #' @param phiInit an initial estimator for phi (size: p*m*k)
   6 #' @param rhoInit an initial estimator for rho (size: m*m*k)
   7 #' @param piInit    an initial estimator for pi (size : k)
   8 #' @param gamInit an initial estimator for gamma
   9 #' @param mini      minimum number of iterations in EM algorithm
  10 #' @param maxi      maximum number of iterations in EM algorithm
  11 #' @param gamma  power in the penalty
  12 #' @param glambda grid of regularization parameters
  13 #' @param X          matrix of regressors
  14 #' @param Y          matrix of responses
  15 #' @param thresh real, threshold to say a variable is relevant, by default = 1e-8
  16 #' @param eps        threshold to say that EM algorithm has converged
  17 #' @param ncores Number or cores for parallel execution (1 to disable)
  18 #'
  19 #' @return a list of outputs, for each lambda in grid: selected,Rho,Pi
  20 #'
  21 #' @examples TODO
  22 #'
  23 #' @export
  24 #'
  25 selectVariables = function(phiInit,rhoInit,piInit,gamInit,mini,maxi,gamma,glambda,
  26                            X,Y,thresh=1e-8,eps, ncores=3, fast=TRUE)
  27 {
  28   if (ncores > 1)
  29   {
  30     cl = parallel::makeCluster(ncores, outfile='')
  31     parallel::clusterExport(cl=cl,
  32                             varlist=c("phiInit","rhoInit","gamInit","mini","maxi","glambda","X","Y","thresh","eps"),
  33                             envir=environment())
  34   }
  35
  36   # Computation for a fixed lambda
  37   computeCoefs <- function(lambda)
  38   {
  39     params = EMGLLF(phiInit,rhoInit,piInit,gamInit,mini,maxi,gamma,lambda,X,Y,eps,fast)
  40
  41     p = dim(phiInit)[1]
  42     m = dim(phiInit)[2]
  43
  44     #selectedVariables: list where element j contains vector of selected variables in [1,m]
  45     selectedVariables = lapply(1:p, function(j) {
  46       #from boolean matrix mxk of selected variables obtain the corresponding boolean m-vector,
  47       #and finally return the corresponding indices
  48       seq_len(m)[ apply( abs(params$phi[j,,]) > thresh, 1, any ) ]
  49     })
  50
  51     list("selected"=selectedVariables,"Rho"=params$rho,"Pi"=params$pi)
  52   }
  53
  54   # For each lambda in the grid, we compute the coefficients
  55   out <-
  56     if (ncores > 1)
  57       parLapply(cl, glambda, computeCoefs)
  58   else
  59     lapply(glambda, computeCoefs)
  60   if (ncores > 1)
  61     parallel::stopCluster(cl)
  62   # Suppress models which are computed twice
  63   #En fait, ca ca fait la comparaison de tous les parametres
  64   #On veut juste supprimer ceux qui ont les memes variables sélectionnées
  65   #sha1_array <- lapply(out, digest::sha1)
  66   #out[ duplicated(sha1_array) ]
  67   selec = lapply(out, function(model) model$selected)
  68   ind_dup = duplicated(selec)
  69   ind_uniq = which(!ind_dup)
  70   out2 = list()
  71   for (l in 1:length(ind_uniq)){
  72     out2[[l]] = out[[ind_uniq[l]]]
  73   }
  74   out2
  75 }