R/selectVariables.R

   1 #' selectVaribles
   2 #' It is a function which construct, for a given lambda, the sets of
   3 #' relevant variables and irrelevant 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 thres  threshold to consider a coefficient to be equal to 0
  16 #' @param tau        threshold to say that EM algorithm has converged
  17 #'
  18 #' @return
  19 #' @export
  20 #'
  21 #' @examples
  22 selectVariables <- function(phiInit,rhoInit,piInit,gamInit,
  23     mini,maxi,gamma,glambda,X,Y,thres,tau)
  24 {
  25     dimphi <- dim(phiInit)
  26     p <- dimPhi[1]
  27     m <- dimPhi[2]
  28     k <- dimPhi[3]
  29     L <- length(glambda);
  30     A1 <- array(0, dim <- c(p,m+1,L))
  31     A2 <- array(0, dim <- c(p,m+1,L))
  32     Rho <- array(0, dim <- c(m,m,k,L))
  33     Pi <- array(0, dim <- c(k,L));
  34
  35     # For every lambda in gridLambda, comutation of the coefficients
  36     for (lambdaIndex in c(1:L))
  37     {
  38         Res <- EMGLLF(phiInit,rhoInit,piInit,gamInit,mini,maxi,
  39             gamma,glambda[lambdaIndex],X,Y,tau);
  40         phi <- Res$phi
  41         rho <- Res$rho
  42         pi <- Res$pi
  43
  44         # If a coefficient is larger than the threshold, we keep it
  45         selectedVariables <- array(0, dim = c(p,m))
  46         discardedVariables <- array(0, dim = c(p,m))
  47         atLeastOneSelectedVariable <- false
  48         for (j in c(1:p))
  49         {
  50             cpt <- 1
  51             cpt2 <-1
  52             for (mm in c(1:m))
  53             {
  54                 if (max(abs(phi[j,mm,])) > thres)
  55                 {
  56                     selectedVariables[j,cpt] <- mm
  57                     cpt <- cpt+1
  58                     atLeastOneSelectedVariable <- true
  59                 } else
  60                 {
  61                     discardedVariables[j,cpt2] <- mm
  62                     cpt2 <- cpt2+1
  63                 }
  64             }
  65         }
  66
  67         # If no coefficients have been selected, we provide the zero matrix
  68         # We delete zero coefficients: vec = indices of zero values
  69         if (atLeastOneSelectedVariable)
  70         {
  71             vec <- c()
  72             for (j in c(1:p))
  73             {
  74                 if (selectedVariables(j,1) != 0)
  75                     vec <- c(vec,j)
  76                 # Else ( NOTE: [auder] else ?! TODO: explain? )
  77                 # we provide the indices of relevant coefficients
  78                 A1[,1,lambdaIndex] <- c(vec,rep(0,p-length(vec)))
  79                 A1[1:length(vec),2:(m+1),lambdaIndex] <- selectedVariables[vec,]
  80                 A2[,1,lambdaIndex] <- 1:p
  81                 A2[,2:(m+1),lambdaIndex] <- discardedVariables
  82                 Rho[,,,lambdaIndex] <- rho
  83                 Pi[,lambdaIndex] <- pi
  84             }
  85         }
  86     }
  87
  88     return(res = list(A1 = A1, A2 = A2 , Rho = Rho, Pi = Pi))
  89 }