##
##  Copyright (C) 2000-2005 CEA
##
##  This software and supporting documentation were developed by
##      CEA/DSV/SHFJ
##      4 place du General Leclerc
##      91401 Orsay cedex
##      France
##

.PRINT.ERROR<<-TRUE
.DEBUG2<<-FALSE
##- -------------------------------------------------------------------------
##- Macro functions for wrappers : Do
##  - load precomputed data (X)
##  - run wrapper if nedded
##  - apply dimension selection
## The only difference is that one should not provide an X but rather an
## number of best features in Y to take
.noWrapperCACHE<<-list()
noWrapper<-function(
                    nX=NULL,Y,
                    T.in=Inf,
                    callBack=NULL,
                    DEBUG=FALSE,
                    var='',file=''
                    )
  {
    if(is.null(.noWrapperCACHE[[ var ]] ))
       {
         YT=fixedThreshold(Y,T.in, bigValsAreGood=TRUE)
         Y=featNamedVec(YT)
         ## 1. Build X list of preselected features
         if(!is.null(nX))
           {
             if(is.numeric(nX)) preSelectY=Y[1:nX] else preSelectY=nX
             default=preSelectY#Y[1:nX]
             
             x= list(
               selectedFeature=preSelectY[[length(preSelectY)]],#Y[[nX]],
               score          =1,
               features       =preSelectY#Y[1:nX]
               )
             X=list()
             X[[1]]=x
             Y=setdiff(Y,preSelectY)
             ##        Y = Y[(nX+1):length(Y)]
           }
         else
           {
             X=list()
             default=c()
           }
         
         for(i in 1:length(Y))
           {
             x= list(
               selectedFeature=Y[[i]],
               score          =1,
               features       =c(default,Y[1:i])
               )
             X[[length(X)+1]]=x
           }
         attr(X,'nested')=TRUE
         selection=list(featSets=X, ite.counter=0)
         .noWrapperCACHE[[ var ]]<<-selection
       }
       else
       {
         selection=.noWrapperCACHE[[ var ]]
       }
    return(selection)
  }

LFS<-function(
              iterations=NULL, nX=NULL, Y, J,
              T.in=Inf,improve=0,
              callBack=NULL,
              DEBUG=FALSE,
              var='',file=''
              )
  {
    YT=fixedThreshold(Y,T.in, bigValsAreGood=TRUE)
    Y=featNamedVec(YT)
    ## 1. Build X list of preselected features
    if(!is.null(nX))
      {
        if(is.numeric(nX)) preSelectY=Y[1:nX] else preSelectY=nX
        score=evalOnFeatures(self=J, features=preSelectY)
##        if(inherits(score,"try-error"))stop('Cannot evaluate J on features')
        x= list(
          selectedFeature=preSelectY[[length(preSelectY)]],
          score          =score,
          features       =preSelectY
          )
        X=list()
        X[[1]]=x
        Y=setdiff(Y,preSelectY)
##        Y = Y[(nX+1):length(Y)]
      }
    else
        X=list()

    ## 2. Call the wrapper if needed
    func=list(
      'wrapper.LFS',
      list(
           iterations=iterations,
           X=X, Y=Y, J=J, improve=improve,
           callBack=callBack,
           DEBUG=DEBUG
           )
      )
    selection=ifLoadElseComputeAndStore(var,file,func,DEBUG=DEBUG)
    if(DEBUG&&.DEBUG2)
      {
        cat('#-- LFS, before dim selection X:')
        print(selection$featSets)
      }
    return(selection)
  }
SFS<-function(
              iterations=10, nX=NULL, Y, J,improve=NULL,
              T.in=Inf,
              callBack=NULL,
              DEBUG=FALSE,
              var='',file=''
              )
  {
    YT=fixedThreshold(Y,T.in, bigValsAreGood=TRUE)
    Y=featNamedVec(YT)
    ## 1. Build X list of preselected features
    if(!is.null(nX))
      {
        if(is.numeric(nX)) preSelectY=Y[1:nX] else preSelectY=nX
        score=evalOnFeatures(self=J, features=preSelectY)
##        if(inherits(score,"try-error"))stop('Cannot evaluate J on features')
        x= list(
          selectedFeature=preSelectY[[length(preSelectY)]],
          score          =score,
          features       =preSelectY
          )
        X=list()
        X[[1]]=x
        Y=setdiff(Y,preSelectY)
      }
    else
        X=list()

    ## 2. Call the wrapper if needed
    func=list(
      'wrapper.SFS',
      list(
           iterations=iterations,
           X=X,Y=Y,J=J,improve=improve,callBack=callBack,
           DEBUG=DEBUG
           )
      )
    selection=ifLoadElseComputeAndStore(var,file,func,DEBUG=DEBUG)
    if(DEBUG&&.DEBUG2)
      {
        cat('#-- SFS, before dim selection X:')
        print(selection$featSets)
      }
    return(selection)
  }
SBS<-function(
                                     iterations=NULL, Y, J,
                                     T.in=Inf,
                                     callBack=NULL,
                                     DEBUG=FALSE,
                                     var='',file=''
                                     )
  {
    YT=fixedThreshold(Y,T.in, bigValsAreGood=TRUE)
    Y=featNamedVec(YT)
    if(is.null(iterations))iterations=length(Y)-1
    ## 2. Call the wrapper if needed
    func=list(
      'wrapper.SBS',
      list(
           iterations=iterations,
           Y=Y, J=J,
           callBack=callBack,
           DEBUG=DEBUG
           )
      )
    selection=ifLoadElseComputeAndStore(var,file,func,DEBUG=DEBUG)
    if(DEBUG&&.DEBUG2)
      {
        cat('#-- SBS, before dim selection X:')
        print(selection$featSets)
      }
    return(selection)
  }

SFFS<-function(
               iterations=10, nX=NULL, Y, J, improve=NULL,
               T.in=Inf,
               callBack=NULL,
               DEBUG=FALSE,
               var='',file=''
               )
  {
    YT=fixedThreshold(Y,T.in, bigValsAreGood=TRUE)
    Y=featNamedVec(YT)
    ## 1. Build X list of preselected features
    if(!is.null(nX))
      {
        if(is.numeric(nX)) preSelectY=Y[1:nX] else preSelectY=nX
        score=evalOnFeatures(self=J, features=preSelectY)
##        if(inherits(score,"try-error"))stop('Cannot evaluate J on features')
        x= list(
          selectedFeature=preSelectY[[length(preSelectY)]],
          score          = score,
          features       = preSelectY
          )
        X=list()
        X[[1]]=x
        Y=setdiff(Y,preSelectY)
      }
    else
        X=list()

    ## 2. Call the wrapper if needed
    func=list(
      'wrapper.SFFS',
      list(
           iterations=iterations,
           X=X, Y=Y, J=J, improve=improve,
           callBack=callBack,
           DEBUG=DEBUG
           )
      )
    selection=ifLoadElseComputeAndStore(var,file,func,DEBUG=DEBUG)
    if(DEBUG&&.DEBUG2)
      {
        cat('#-- SFFS, before dim selection X:')
        print(selection$featSets)
      }
    return(selection)
  }

linearSequentialForwardFloatingSelection<-
  function(
           iterations=10,nX=NULL,Y,J,improve=0,
           linSelectDimMax=Inf,T.in=Inf,
           callBack=NULL,
           DEBUG=FALSE,
           var='',file=''
           )
  {
    Ynamed=featNamedVec(Y)
 
    ## 1. Build X list of preselected features
    if(!is.null(nX))
      {
        if(is.numeric(nX)) preSelectY=Ynamed[1:nX] else preSelectY=nX
        score=evalOnFeatures(self=J, features=preSelectY)
##        if(inherits(score,"try-error"))stop('Cannot evaluate J on features')
        x= list(
          selectedFeature= preSelectY[[length(preSelectY)]],
          score          = score,
          features       = preSelectY
          )
        X=list()
        X[[1]]=x
        Ynamed=setdiff(Ynamed,preSelectY)
##        Ynamed = Ynamed[(nX+1):length(Y)]
      }
    else
        X=list()

    ## 2. Linear wrapper on all features
    func=list(
      'wrapper.LFS',
      list(
           iterations=Inf,
           X=X, Y=Ynamed, J=J, dimMax=linSelectDimMax, improve=improve,
           callBack=callBack,
           DEBUG=DEBUG
           )
      )
    ## Remove iterations=blabla and T.in=blabla in the var name
    var1=sub('iterations=[^,]+,','',var)
    var1=sub('T.in=[^,]+,','',var1)

    lin.selection=ifLoadElseComputeAndStore(paste(var1,'_lin',sep=''),
      file,func,DEBUG=DEBUG)
    Xlin=lin.selection$featSets
    
    ## 3. SFFS on thresholded features with feat sets linearly selected
    ## Feature = thresholded features + linearly selected features
    YT=fixedThreshold(Y,T.in, bigValsAreGood=TRUE)
    Ynamed=featNamedVec(YT)
    Xfeats=Xlin[[length(Xlin)]]$features
    Ynamed=setdiff(Ynamed,Xfeats)
    
    ## 4. Call the SFFS wrapper if needed
    func=list(
      'wrapper.SFFS',
      list(
           iterations=iterations,
           X=Xlin, Y=Ynamed, J=J, improve=improve,
           callBack=callBack,
           DEBUG=DEBUG
           )
      )
    selection=ifLoadElseComputeAndStore(paste(var,'_sffs',sep=''),
      file,func,DEBUG=DEBUG)
    if(DEBUG&&.DEBUG2)
      {
        cat('#-- linearSequentialForwardFloatingSelection, before dim selection X:')
        print(selection$featSets)
      }
    return(selection)
  }

combinatorialSelection<-function(
                                 iterations=c(0.6, 0.7), nX=0, Y, J,
                                 callBack=NULL,
                                 DEBUG=FALSE,
                                 var='',file=''
                                 )
  {  
    X=list()

    ## 2. Call the wrapper if needed
    func=list(
      'wrapper.combinatorialSelection',
      list(
           iterations=iterations,#20
           X=X, Y=Y, J=J,
           callBack=callBack,
           DEBUG=DEBUG
           )
      )
    selection=ifLoadElseComputeAndStore(var,file,func,DEBUG=DEBUG)
    if(DEBUG&&.DEBUG2)
      {
        cat('#-- combinatorialSelection, before dim selection X:')
        print(selection$featSets)
      }
    return(selection)
  }
## - -------------------------------------------------------------------------
## iterations    : Number of forward feature selection
## X             : Already selected features and associated scores list
##                 of names vector
##               : X:= list(of X[[k]] ), where each X[[k]]$features
##                 X[[k]]$score
## k             : The current step
## Y             : Available features
## J             : Function to be optimized, it shuld return a list as
##                 list(features=c(), score=value, ....)
wrapper.LFS <- function
(
 iterations=Inf, X=list(), Y=c(), J, dimMax=Inf,improve=0,
 callBack=NULL,
 DEBUG=FALSE
 )
{
  Y=featNamedVec(Y)
  if(is.null(iterations)||(iterations==Inf))iterations=length(Y)
  if(length(X)!=0)
    {
      k=length(X)
      scoremax = X[[k]]$score
      dim=length(X[[k]]$features)
    }
  else
    {
      k=0
      dim=0
      scoremax = -1
    }

  ite = 1
  while (ite <= iterations && dim<=dimMax){
    if(DEBUG)  cat('# - LFS, ite:',ite,'\n')
    xi = Y[[1]]
    
    if (k!=0) Xk = X[[k]]$features
    else Xk = c()
    
    features = c(Xk, xi)
    score = try( evalOnFeatures(self=J, features=features) )
    if(inherits(score,"try-error")){
      if(.PRINT.ERROR)
        {
          cat('# - ERROR in LFS\n')
          cat('------------- Error:\n')
          print(score)
          cat('------------- Features:\n')
          print(features)
          cat('-------------\n')
        }
      if(length(grep('break forward loop',as.character(score))))break
    }
    else{
      if(DEBUG) cat('\t',xi,score,'max',scoremax,scoremax+improve,'\n')
      Xcurrent=list(selectedFeature=xi, score=score, features=features)
      if(score>=(scoremax+improve)){
        cat('\t',score,'-',scoremax,'=',score-scoremax,'%\t=>+',(score-scoremax)*.NBSUBJECTS,'\n',sep='')
        scoremax    = score
        Xmost = Xcurrent

        if(DEBUG) cat('\tadd feature:',
                      Xmost$selectedFeature,' ',Xmost$score,'\n')
                                        # - Add the feature in X[[k+1]]
        X[[k+1]] = Xmost
        k        = k+1
        dim      = dim+1
      }
                                        # - Call back with best Xk
      if(!is.null(callBack))
        do.call(callBack[[1]],c(callBack[[2]],list(Xcurrent)))
                                        # - Pop feature from Y
      Y=Y[2:length(Y)]
      ite = ite+1
    }
  }
  selection=list(featSets=X, ite.counter=(ite-1))
  if(DEBUG&&.DEBUG2){
    cat('# - LFS END, X:\n')
    print(selection)
  }
  return(selection)
}

## - -------------------------------------------------------------------------
## iterations    : Number of forward feature selection
## X             : Already selected features and associated scores list
##                 of names vector
##               : X:= list(of X[[k]] ), where each X[[k]]$features
##                 X[[k]]$score
## k             : The current step
## Y             : Available features
## J             : Function to be optimized, it shuld return a list as
##                 list(features=c(), score=value, ....)
wrapper.SFS <- function
(
 iterations=Inf,X=list(),Y=c(),J,improve=NULL,
 callBack=NULL,
 DEBUG=FALSE
 )
{
  Y=featNamedVec(Y)
  iterations=min(iterations,length(Y))
  if(length(X)!=0)
      k=length(X)
  else
      k=0
  
  ite.counter=0
  ite = 1
  while (ite <= iterations){
    if(DEBUG)  cat('# - SFS, ite:',ite,'\n')

    if (k!=0){ Xk = X[[k]]$features; currentScore=X[[k]]$score}
    else { Xk = c();currentScore=0}

### FIX A BUG HERE-----------------------------
###    ## Send the restart iteration signal
###    ## (may be ignored or use to relearn model)
    if(k!=0) evalOnFeatures(self=J, Xk,reStart=TRUE)
    
    Xmost = seqFeatSelect.most(Xk=Xk, Y=Y, J=J)

    if(!is.null(improve) && Xmost$score<(currentScore+improve))
      {
        if(DEBUG) cat('# - SFS Stop:',
                      Xmost$score,'<',currentScore,'+',improve,'\n')
        break
      }
    currentScore=Xmost$score
    ite.counter=ite.counter+Xmost$ite.counter
    
    if(DEBUG) cat('\t add feature:',
                  Xmost$selectedFeature,' ',Xmost$score,'\n')
                                        # - Call back with best Xk
    if(!is.null(callBack))do.call(callBack[[1]], c(callBack[[2]],list(Xmost)))
                                        # - Add the mas feature in X[[k+1]]
    X[[k+1]] = Xmost
    k        = k+1
                                        # - Remove best new feature from Y
    Y = setdiff(Y, Xmost$selectedFeature)
    ite = ite+1
  }
  if(DEBUG&&.DEBUG2){
    cat('# - SFS END, X:\n')
    print(X)
  }
  return(list(featSets=X,ite.counter=ite.counter))
}

wrapper.SBS <- function
(
 iterations=NULL, Y=c(), J,
 callBack=NULL,
 DEBUG=FALSE
 )
{
  Y=featNamedVec(Y)
  iterations=min(iterations,(length(Y)-1))
  k=length(Y)
  X=list()
  ## Compute X[[k]]
  X[[k]]=list()
  X[[k]]$features = Y
  score=evalOnFeatures(self=J, X[[k]]$features)
##  if(inherits(score,"try-error"))stop('Cannot evaluate J on features')
  X[[k]]$score=score
  
  ite.counter=0
  ite = 1
  while (ite <= iterations){
    ## Send the restart iteration signal
    ## (may be ignored or use to relearn model)
    evalOnFeatures(self=J, X[[k]]$features,reStart=TRUE)
    if(DEBUG)  cat('# - SBS, ite:',ite,'\n')
    
    Xleast = seqFeatSelect.least(Xk=X[[k]]$features, J=J)
    ite.counter=ite.counter+Xleast$ite.counter
    
    if(DEBUG) cat('# - SBS remove feature:',
                  Xleast$selectedFeature,' ',Xleast$score,'\n')
                                        # - Call back with best Xk
    if(!is.null(callBack))do.call(callBack[[1]], c(callBack[[2]],list(Xleast)))
                                        #~- Add the mas feature in X[[k+1]]
    X[[k-1]] = Xleast
    k        = k-1
    ite = ite+1
  }
  if(DEBUG&&.DEBUG2){
    cat('# - SBS END, X:\n')
    print(X)
  }
  return(list(featSets=X, ite.counter=ite.counter))
}

## - -------------------------------------------------------------------------
## Stop :
## <iterations> : to fix the number of max iteration
## <l>          : to fix the number of max feature ie. : k
## Nothing     :  will run until there are no more features in Y
wrapper.SFFS<-function
(
 iterations=Inf, X=list(), Y=c(), J,improve=NULL,
 callBack=NULL,
 DEBUG=FALSE
 )
{
  Y=featNamedVec(Y)
  ite.counter=0
  
  if(length(X)!=0)
    k=length(X)
  else
    {
                                        # - if X is empty run twice forward
      most = wrapper.SFS(iterations=2, X=X, Y=Y, J=J,
        DEBUG=DEBUG)
      X = most$featSets
      Y = most$Y
      k = 2
      ite.counter=ite.counter+most$ite.counter
    }
  
  ite=1
  while( length(Y)>0 && ite<=iterations)
    {
      if(DEBUG){
        cat('# - SFFS --------------------------------------------- X[',k,']=',
            X[[k]]$score,', Features :\n',sep='')
        #print(X[[k]]$features)
      }
                                        # - ----------------------------------
                                        # - Step I  : Inclusion
      currentScore=X[[k]]$score
      Xmost = seqFeatSelect.most(Xk=X[[k]]$features, Y=Y, J=J)
      if(!is.null(improve)&&Xmost$score<(currentScore+improve))
        {
          if(DEBUG) cat('# - SFFS Stop:',
                        Xmost$score,'<',currentScore,'+',improve,'\n')
          break
        }
      ite.counter=ite.counter+Xmost$ite.counter
                                        # - Add the most feature in X[[k+1]]
      X[[k+1]] = Xmost
      k        = k+1
                                        # - Remove most new feature from Y
      Y = setdiff(Y, Xmost$selectedFeature)
      
      if(DEBUG) cat('\tSFFS most   :',
                    X[[k]]$selectedFeature,' ** ',X[[k]]$score,'\n')
      
                                        # - ----------------------------------
                                        # Step II : Conditionnal Exclusion
                                        # (Sequential Backward Selection)
                                        # Remove features while it improves
                                        # the score
      continueBackward = TRUE
      firstBackward    = TRUE
      while(continueBackward){
        least = seqFeatSelect.least(Xk=X[[k]]$features, J=J)
        ite.counter=ite.counter+least$ite.counter
        if(DEBUG) cat('\tSFFS least  :',least$selectedFeature,
                      ' ** ',least$score,'\n')
        
        if(firstBackward)firstBackward = FALSE
                                        # Removal the least if it does
                                        # improve the score
                                        # AND it is not the feature
                                        # that we'eve juste add
                                        # (in order to avoid infinite loop)
        if( (k>=2) && (least$score > X[[k-1]]$score) &&
           !(firstBackward && (least$selectedFeature==X[[k]]$selectedFeature)))
          {
            oldKminusOne = X[[k-1]]$score
            X[[k-1]]$score    = least$score
            X[[k-1]]$features = setdiff(X[[k]]$features, least$selectedFeature)
            Y = c(Y, least$selectedFeature)
            k = k-1
            if(DEBUG)
              cat('\tSFFS remove :',least$selectedFeature,
                  ' ** ',least$score,' -> X[',k,']=',X[[k]]$score,
                  ' ** Delta=',(least$score-oldKminusOne),'\n')
          }
                                        # Removal the least does not improve
                                        # the score
                                        # no no further backward search is
                                        # performed
        else continueBackward = FALSE
      }
      if(!is.null(callBack))
        do.call(callBack[[1]], c(callBack[[2]],list(X[[k]])))
      ite=ite+1
    }
  return(list(featSets=X, ite.counter=ite.counter))
}


## - -------------------------------------------------------------------------
## Least : given a set of <Xk>, try to remove each of them and return least
## significant feature,
## ie. : the one when removed that affect the least the Criterion J
## feature = arg max J(Xk \ xi)
seqFeatSelect.least <- function(Xk,J,DEBUG=FALSE){
  scoremax=-1                           # Best set of features
                                        # Iterate over already selected
                                        # features
  ite.counter=0
  for(xi in Xk){
##    cat('----------------------------------\n')
    subSetFeatures = setdiff(Xk, xi)
    score = try(evalOnFeatures(self=J, features=subSetFeatures))
##    if(DEBUG)
##    cat(xi,score,'\n')

    if(inherits(score,"try-error")){
      if(.PRINT.ERROR)
        {
          cat('# - ERROR in seqFeatSelect.least\n')
          cat('------------- Error:\n')
          print(score)
          cat('------------- Features:\n')
          print(subSetFeatures)
          cat('-------------\n')
        }
    }
    else{
      ite.counter=ite.counter+1
      if(DEBUG) cat('# - seqFeatSelect.least try remove feature:',xi,
                    ' ',score,'\n')
      if(score>scoremax){
        xmax        = xi                         # - Best new feature
        scoremax    = score
        featuresmax = subSetFeatures
      }
    }
  }
  return( list(selectedFeature = xmax,
               score           = scoremax,
               features        = featuresmax,
               ite.counter     = ite.counter)
         )
}

## - -------------------------------------------------------------------------
## Most : given a set of features <Xk>, try to add each feature in <Y> and
## return most significant feature,
## ie. : the one when added that improve the most the Criterion J
## feature = arg max J(Xk U {xi in Y})
seqFeatSelect.most <- function(Xk,Y,J,DEBUG=FALSE){
  scoremax=-1                           # - Best set of features
                                        # - Iterate over availables features
  ite.counter=0
  for(xi in Y){
    features = c(Xk, xi)
    score = try( evalOnFeatures(self=J, features=features) )
    if(inherits(score,"try-error")){
      if(.PRINT.ERROR)
        {
          cat('# - ERROR in seqFeatSelect.most\n')
          cat('------------- Error:\n')
          print(score)
          cat('------------- Features:\n')
          print(features)
          cat('-------------\n')
        }
    }
    else{
      ite.counter=ite.counter+1
      if(DEBUG) cat('# - seqFeatSelect.most:',xi,' ',score,'\n')
      if(score>scoremax){
        xmax        = xi                         # - Best new feature
        scoremax    = score
        featuresmax = features                  
      }
    }
  }
  return( list(selectedFeature= xmax,
               score          = scoremax,
               features       = featuresmax,
               ite.counter    = ite.counter)
         )
}
##- -------------------------------------------------------------------------
## For each item in X (witk k feature) ad one feature from Y, and repeat it for all Y.
## This will form a new X list with k+1 feature
##
## iterations is a set of thresholds to apply in each dimension before computing the next dimension
## Exemples : c(0.6, 0.7, 0.72)
## Threshold 1D with with 0.6  to get -> 2D
## Threshold 2D with with 0.7  to get -> 3D
## Threshold 3D with with 0.72 to get -> 4D
wrapper.combinatorialSelection <- function
(
 iterations=c(0.6, 0.7), X=list(), Y=c(), J,
 callBack=NULL,
 DEBUG=FALSE
 )
{
  DIM=list()

  if(length(X)==0)##Build X from Y
      DIM[[1]]=featWeightedVec2featScoreList(Y)
  else
    DIM[[1]]=X

  Y=featNamedVec(Y)
##  print(DIM[[1]])
##  print(Y)
  for(dim in 2:(length(iterations)+1))
    {
      cat('--------------- ite=',iterations[[dim-1]],'\n')
      Xnew=list()
      item=1
      X=wrapper.fixedThreshold(DIM[[dim-1]],iterations[[dim-1]])
      print(X)
      ## Do some cleaning : remove redondant set of features
##      cat('CLEANING START-------------------------------------------- X:\n')
      if(length(X)>1)
        {
          index=rep(TRUE,length(X))
          for(i in 1:(length(X)-1))
            {
              xFeat=sort(X[[i]]$features)
              for(j in (i+1):length(X))
                {
                  if(index[[j]] && (all(xFeat==sort(X[[j]]$features))))index[[j]]=FALSE
                }
            }
##          cat('CLEANING END--------------------------------------------\n')
          X=X[index]
          ## Store the cleaned list in the previous dimension (for further analysis)
           DIM[[dim-1]]=X
##          storeObject('X','toto.R')
          if(DEBUG){cat('#- wrapper.combinatorialSelection found',length(X),'feat. set of dim',(dim-1),'>',iterations[[dim-1]],'\n')}
        }
      
      for(x in X)
        {
          xFeat=x$features
          cat('--- x',xFeat,'\n')
          for(y in Y)
            {
##              cat('- y',y,'\n')
              if(!(y %in% xFeat))
                {
                  features=c(xFeat,y)
                  score=try(evalOnFeatures(self=J, features=features))
                  if(inherits(score,"try-error")){
                    if(.PRINT.ERROR)
                      {
                        cat('# - ERROR in combinatorialSelection\n')
                        cat('------------- Error:\n')
                        print(score)
                        cat('------------- Features:\n')
                        print(features)
                        cat('-------------\n')
                      }
                  }
                  else
                    {
                      if(score>x$score)
                        {
                          Xnew[[item]]=list(features=features,score=score)
                          item=item+1
                        }
                    }
                }
            }
        }
      if(length(Xnew)==0) break
      DIM[[dim]]=Xnew
    }
  ## Sort & Merge
  DIM2=list()
    for ( dim in 1:length(DIM))
      {
        X=DIM[[dim]]
        if(length(DIM)==0)next
        X=wrapper.sort(X)
        DIM2=c(DIM2,X)
      }
    XY=list()
    XY$X=DIM2
##  cat('*** wrapper.combinatorialSelection S *******\n')
##  print(l)
##  print(XY)
##  cat('*** wrapper.combinatorialSelection E *******\n')

  return(XY)
}
wrapper.combinatorialSelection.exemples <- function()
  {
    source('~duchesnay/p4/brainvisa-main/lib/dataMind/modules/PatternRecognition/mlInit.R')
    mlInit()

    ## Exemples
    D=rt('arti_correlated_100.dat')
    class='class'
    
    D=rt('df_50Neg_r0.csv')
    class='SEX'
 
    
    data.y = D[[class]]
    data=D
    data.x=rmcol(data,c(class,'SUBJECT'))
    data.x=rmcol(data.x,c(class,'SUBJECT',  'x0.5', 'x0.5c', 'x0.5ci', 'x0.4', 'x0.4c', 'x0.4ci', 'x0.3', 'x0.3c', 'x0.3ci' ))

    data.x=as.matrix(data.x)
    X.norm=filter.nonNormality(data.x,data.y)
    XT.norm=fixedThreshold(X.norm,0.95)
    X.looLda   =filter.looRanking(data.x,data.y,Y=XT.norm,classifier='lda',DEBUG=TRUE)
    X.looLdaScore=pvalOrErr2Score(X.looLda)
    J=ClassSepLdaLoo(data.x=data.x, data.y=data.y)
    res=wrapper.combinatorialSelection(iterations=c(0.6, 0.7, 0.75, 0.8, 0.82, 0.84, 0.86, 0.88, 0.9, 0.92), Y=X.looLdaScore, J=J,DEBUG=TRUE)
    
    t1=sort.dataframe(featScoreList2Tab(res[[1]]),'score',decreasing=TRUE);wt(t1,'t1.csv')
    t2=sort.dataframe(featScoreList2Tab(res[[2]]),'score',decreasing=TRUE);wt(t2,'t2.csv')
    t3=sort.dataframe(featScoreList2Tab(res[[3]]),'score',decreasing=TRUE);wt(t3,'t3.csv')
    t4=sort.dataframe(featScoreList2Tab(res[[4]]),'score',decreasing=TRUE);wt(t4,'t4.csv')
    t5=sort.dataframe(featScoreList2Tab(res[[5]]),'score',decreasing=TRUE);wt(t5,'t5.csv')
    t6=sort.dataframe(featScoreList2Tab(res[[6]]),'score',decreasing=TRUE);wt(t6,'t6.csv')
    t7=sort.dataframe(featScoreList2Tab(res[[7]]),'score',decreasing=TRUE);wt(t7,'t7.csv')
    t8=sort.dataframe(featScoreList2Tab(res[[8]]),'score',decreasing=TRUE);wt(t8,'t8.csv')
    t9=sort.dataframe(featScoreList2Tab(res[[9]]),'score',decreasing=TRUE);wt(t9,'t9.csv')
    t10=sort.dataframe(featScoreList2Tab(res[[10]]),'score',decreasing=TRUE);wt(t10,'t10.csv')

  }

##- -------------------------------------------------------------------------
## Utils
featScoreList2Tab<-function(list,featAsStr=FALSE)
  {
    RES=NULL
    for(x in list)
      {
        if(featAsStr)
          {
            features=vec2str(x$features)
            names(features)='features'
          }
        else
          {
            features=as.list(x$features)
            names(features)=paste('x',1:length(x$features),sep='')
          }
        RES=append.dataframe(RES,c(features,list(score=x$score)))
      }
    rownames(RES)=1:nrow(RES)
    return(RES)
  }

##- -------------------------------------------------------------------------
## Utils convert from filters
featWeightedVec2featScoreList<-function(vec)
  {
    l=list()
    i=1
    for(n in names(vec))
        {
          t=list(features=n, score=vec[[n]])
          l[[i]]=t
          i=i+1
        }
    return(l)
  }

pvalOrErr2Score<-function(v)
  {
    return(1-v)
  }

## = ========================================================================
## =
## =  Feature subsets manipulation
## =
## = ========================================================================

##- -------------------------------------------------------------------------
## 2. Dimension selection
## dimension selection take as input : X a list of several feature subset
## return a list of feature set

selectDim<-function(X,dim)
  {
    if(dim<length(X[[1]]$features))
      {
        X=X[[1]]
        X$features=X$features[1:dim]
        return(list(X))
      }
    if(!is.null(attr(X,'nested')) )
      {
        firstDim=length(X[[1]]$features)
        k=dim-firstDim+1
        cat(dim, length(X[[k]]$features), '-------\n')
        return(list(X[[k]]))
      }
    for(k in 1:length(X)){
      d=length(X[[k]]$features)
##      cat(d,dim,'\n')
      if(d==dim)return(list(X[[k]]))
      }
      stop('Unavailable dimension')
    }

selectDimBests<-function(X,nb=1)
  {
    X=wrapper.sort(X)
    return( X[1:nb] )
  }

selectDimBest<-function(X)
  {
    bestScore = -Inf
    for(k in 1:length(X)){
      score     = X[[k]]$score
      if(score>bestScore){
        bestK = k
        bestScore = score
      }
    }
    return( list(X[[bestK]]) )
  }


##- -------------------------------------------------------------------------
## Thresholding & Sorting
wrapper.fixedThreshold<-function(X,T)
  {
    XT=list()
    i=1
    for(x in X)
      {
        if(x$score>T)
          {
            XT[[i]]=x
            i=i+1
          }
      }
    return(XT)
  }

wrapper.sort<-function(X)
  {
    XS=NULL
    for(x in X)
      {
        if(is.null(XS))XS=list(x)
        else
          {
            len=length(XS)
            score=x$score
            i=1
###            cat((i<=len), score, XS[[i]]$score,'\n')
            while((i<=len) && (XS[[i]]$score>=score))i=i+1
            if(i>len)XS[[i]]=x
            else if(i==1)XS=c(list(x),XS)
            else
              {
                XS=c(XS[1:(i-1)], list(x), XS[(i):len])
              }
          }
      }
    return(XS)
  }