#-------------------------------------------------------------
#
# Licensed to the Apache Software Foundation (ASF) under one
# or more contributor license agreements.  See the NOTICE file
# distributed with this work for additional information
# regarding copyright ownership.  The ASF licenses this file
# to you under the Apache License, Version 2.0 (the
# "License"); you may not use this file except in compliance
# with the License.  You may obtain a copy of the License at
#
#   http://www.apache.org/licenses/LICENSE-2.0
#
# Unless required by applicable law or agreed to in writing,
# software distributed under the License is distributed on an
# "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY
# KIND, either express or implied.  See the License for the
# specific language governing permissions and limitations
# under the License.
#
#-------------------------------------------------------------

m_bandit = function(Matrix[Double] X_train, Matrix[Double] Y_train, List[Unknown] metaList, List[Unknown] targetList,
  Frame[Unknown] lp, Frame[Unknown] primitives, Frame[Unknown] param,  Integer k = 3, Integer R=50,
  Boolean verbose = TRUE)
  return (Frame[Unknown] bestPipeline, Matrix[Double] bestHyperparams,  Matrix[Double] bestAccuracy, Frame[String] feaFrameOuter) 
{
  print("Starting optimizer")
  NUM_FEATURES = 14
  HYPERPARAM_LENGTH = 110
  print("null in data "+sum(is.na(X_train)))
  bestPipeline = frame("", rows=1, cols=1)
  bestHyperparams = as.matrix(0)
  bestAccuracy = as.matrix(0)
  # initialize bandit variables
  # variable names follow publication where algorithm is introduced
  eta = 2  # the halving ratio is fixed to 2
  s_max = floor(log(R,eta));
  B = (s_max + 1) * R;
  
  # initialize output variables
  hparam = matrix(0, rows=k*(s_max+1), cols=HYPERPARAM_LENGTH)
  pipeline = frame(0, rows=k*(s_max+1), cols=ncol(lp)+1)
  startOut=0; endOut=0;
  feaFrameOuter = frame(data=["#MissingValues", "MinVla", "MaxVal", "AverageMin", "AverageMax", 
  "#CategoricalFeatures", "#NumericFeatures", "Mean", "#Outliers", "#OHEfeatures", "#Classes",
  "Imbalance", "#rows", "#cols", "pipelines", "accuracy", "execution time in ms", "CV time in ms"],
  rows = 1, cols = NUM_FEATURES + 4 )

  for(s in s_max:0) {
    
   # result variables
    bracket_hp = matrix(0, rows=k*(s+1)+k, cols=HYPERPARAM_LENGTH)
    bracket_pipel = matrix(0, rows=k*(s+1)+k, cols=3)
    start=1; end=0;
    
    # # compute the number of initial pipelines n
    n = ceil(floor(B/R/(s+1)) * eta^s);
    r = R * eta^(-s);
    # get the physical pipelines, the pipelines, pipelines are recoded
    [configurations, n] = get_physical_configurations(lp, n, primitives)
    
    # append configuration keys for extracting the pipeline later on
    id = seq(1, nrow(configurations))
    configurations = cbind(as.frame(id), configurations)
    # save the original configuration as a lookup table
    lookup = configurations
    
    if(verbose)
      print("n "+ n +"\nR "+ R +"\ns_max "+ s_max +"\nB "+ B +"\nn "+ n +"\nr "+ r)
    
    for( i in 0:s) {
      # successive halving
      n_i = min(max(as.integer(floor(n * eta^(-i))), 1), nrow(configurations));
      r_i = as.integer(floor(r * eta^i));
      
      if(verbose) {
        print("no of configurations ---------"+n_i)
        print("no of resources --------------"+r_i)
        print("iteration  ---------------------"+i+" out of "+s)
      }
      
      configurations = configurations[1:n_i, ]
      [outPip,outHp, feaFrameOuter] = run_with_hyperparam(configurations, r_i, X_train, Y_train, metaList,
        targetList, param, feaFrameOuter, verbose)
      # sort the pipelines by order of accuracy decreasing
      a = order(target = outPip, by = 1, decreasing=TRUE, index.return=FALSE)
      b = order(target = outHp, by = 1, decreasing=TRUE, index.return=FALSE)
      rowIndex = min(k, nrow(a))

      # maintain the brackets results
      end = end + rowIndex
      bracket_pipel[start:end, ] =  a[1:rowIndex,]
      bracket_hp[start:end, 1:ncol(b)] =  b[1:rowIndex,]
      start = end + 1

      # sort the configurations for successive halving
      avergae_perf =  getMaxPerConf(outPip, nrow(configurations)) 
      configurations = frameSort(cbind(avergae_perf, configurations))
      configurations = configurations[, 2:ncol(configurations)]
    }
    bracket_pipel = removeEmpty(target=bracket_pipel, margin="rows")
    bracket_hp = removeEmpty(target=bracket_hp, margin="rows")
    # keep the best k results for each bracket
    [bracket_bestPipeline, bracket_bestHyperparams] = extractBracketWinners(bracket_pipel, bracket_hp, k, lookup)
    # optimize by the features

    startOut = endOut + 1
    endOut = endOut + nrow(bracket_bestPipeline)
    pipeline[startOut: endOut, ] = bracket_bestPipeline
    hparam[startOut:endOut, 1:ncol(bracket_bestHyperparams)] = bracket_bestHyperparams
  }

  [bestPipeline, bestHyperparams] = extractTopK(pipeline, hparam, as.scalar(targetList['dirAcc']), k)

  bestAccuracy = as.matrix(bestPipeline[,1])
  bestPipeline = bestPipeline[,2:ncol(bestPipeline)]
  bestHyperparams = bestHyperparams[,2:ncol(bestHyperparams)]
  
  if(verbose) {
    print("best pipeline"+ toString(bestPipeline))
    print("best hyper-parameters \n"+ toString(bestHyperparams))
    print("best accuracy \n"+ toString(bestAccuracy))
    if(as.scalar(targetList['target']) != "compare")
      print("dirty accuracy "+as.scalar(targetList['dirAcc']))
  }
}

# this method will extract the physical pipelines for a given logical pipelines
get_physical_configurations = function(Frame[String] logical, Scalar[int] numConfigs, 
  Frame[Unknown] primitives)
  return(Frame[String] physical, Double min)
{
  # load the primitives
  physical = as.frame("NaN")
  outliers = primitives[,1]
  mvi = primitives[,2]
  noise = primitives[,3]
  ci = primitives[,4]
  dim = primitives[,5]
  dummy = primitives[,6]
  scale = primitives[,7]
 
  operator = frame(0, rows=nrow(primitives), cols=ncol(logical))  # combine all logical primitives
  for(j in 1:ncol(logical))
  {
    # extract the physical primitives
    if(as.scalar(logical[1,j]) == "OTLR")
      operator[, j] = outliers;
    else if(as.scalar(logical[1,j]) == "MVI")
      operator[, j] = mvi;
    else if(as.scalar(logical[1,j]) == "NR")
      operator[, j] = noise;
    else if(as.scalar(logical[1,j]) == "CI")
      operator[, j] = ci;
    else if(as.scalar(logical[1,j]) == "DIM")
      operator[, j] =  dim;
    else if(as.scalar(logical[1,j]) == "DUMMY")
      operator[, j] =  dummy;
    else if(as.scalar(logical[1,j]) == "SCALE")
      operator[, j] = scale;
    else stop("invalid operation "+as.scalar(logical[1,j]))
  }

  idx = matrix(1, rows=1, cols=ncol(logical))
  # get the indexes of columns for recode transformation
  index = vectorToCsv(idx)
  # recode logical pipelines for easy handling
  jspecR = "{ids:true, recode:["+index+"]}";
  [X, M] = transformencode(target=operator, spec=jspecR);
  X = replace(target= X, pattern = NaN, replacement = 0)
  
  paramLens = matrix(0, ncol(logical), 1);
  for( j in 1:ncol(logical)) {
    vect = removeEmpty(target = X[,j], margin = "rows");
    paramLens[j,1] = nrow(vect);
  }
  min = prod(paramLens)
  sample = ifelse(min > numConfigs, TRUE, FALSE)
  paramVals = matrix(0, ncol(logical), max(paramLens));
  for( j in 1:ncol(logical) ) {
    vect = removeEmpty(target = X[,j], margin = "rows");
    paramVals[j,1:nrow(vect)] = t(vect);
  }
  cumLens = rev(cumprod(rev(paramLens))/rev(paramLens));
  # materialize hyper-parameter combinations 
  HP = matrix(0, min(numConfigs, min), ncol(logical));
  if(sample) 
    pip = sample(numConfigs,numConfigs)
  else pip = seq(1,nrow(HP))
  for( i in 1:nrow(HP) ) {
    for( j in 1:ncol(logical) ) {
      HP[i,j] = paramVals[j,as.scalar((as.scalar(pip[i,1])/cumLens[j,1])%%paramLens[j,1]+1)];
    }
  }
  
  physical = transformdecode(target=HP, spec=jspecR, meta=M);
  print("physical pipeline "+toString(physical))
}

# this method will call the execute pipelines with their hyper-parameters
run_with_hyperparam = function(Frame[Unknown] ph_pip, Integer r_i, Matrix[Double] X, Matrix[Double] Y,
  List[Unknown] metaList, List[Unknown] targetList, Frame[Unknown] param, Frame[Unknown] featureFrameOuter,
  Boolean verbose)
  return (Matrix[Double] output_operator, Matrix[Double] output_hyperparam, Frame[Unknown] featureFrameOuter) {

  output_hp = matrix(0, nrow(ph_pip)*r_i, 100)
  output_accuracy = matrix(0, nrow(ph_pip)*r_i, 1)
  output_pipelines = matrix(0, nrow(ph_pip)*r_i, 2)
  
  # rows in validation set
  clone_X = X
  clone_Y = Y
  index = 1
  id = as.matrix(ph_pip[, 1])
  ph_pip = ph_pip[, 2:ncol(ph_pip)]
  
  feaVec = gatherStats(X, Y, as.matrix(metaList['mask']), as.scalar(targetList['target']))

  for(i in 1:nrow(ph_pip))
  {
    # execute configurations with r resources
    [hp, no_of_res, no_of_flag_vars] = getHyperparam(ph_pip[i], param, r_i)
    if(ncol(featureFrameOuter) > 1)
      feaFrame = frame("", rows = no_of_res, cols = ncol(featureFrameOuter))
    pip_toString = pipToString(ph_pip[i])
    for(r in 1:no_of_res)
    {
      # as the matrix first block of r rows belongs to first operator and r+1 block of rows to second operator 
      # we need to extract a row from each block
      indexes = matrix(no_of_res, rows=ncol(ph_pip), cols=1)
      indexes[1, 1] = r
      indexes = cumsum(indexes)
      indexes = table(indexes, 1, 1, nrow(hp), 1)
      hp_matrix = removeEmpty(target = hp, margin="rows", select = indexes)
      [X, Y, T] = executePipeline(ph_pip[i], X, Y, as.matrix(metaList['mask']), as.matrix(metaList['fd']), hp_matrix, no_of_flag_vars, FALSE)
      if(as.scalar(targetList['target']) == "compare")
        [accuracy, accT] = compareValue(clone_X, X, as.matrix(targetList['cleanData']), as.matrix(metaList['mask']))
      else
        [accuracy, accT] = fclassify(X, Y, as.matrix(metaList['mask']), as.matrix(targetList['mlHp']), as.scalar(targetList['dirAcc']), 
        as.scalar(targetList['wAccuracy']), as.scalar(targetList['cv']))
      matrix_width = as.matrix(nrow(hp_matrix) * ncol(hp_matrix))
      hp_vec = cbind(matrix_width, matrix(hp_matrix, rows=1, cols=nrow(hp_matrix)*ncol(hp_matrix), byrow=TRUE))
      output_accuracy[index, 1] = accuracy
      output_hp[index, 1:ncol(hp_vec)] = hp_vec
      output_pipelines[index, ] = cbind(as.matrix(index), id[i,1])
      X = clone_X
      Y = clone_Y
      index = index + 1
      
      if(ncol(featureFrameOuter) > 1) {
        feaFrame[r, 1:ncol(feaVec)] = as.frame(feaVec)
        feaFrame[r, (ncol(feaVec)+1)] = pip_toString
        feaFrame[r, (ncol(feaVec)+2)] = accuracy
        feaFrame[r, (ncol(feaVec)+3)] = T
        feaFrame[r, (ncol(feaVec)+4)] = accT
      }
    }
    
    X = clone_X
    Y = clone_Y
    if(ncol(featureFrameOuter) > 1)
      featureFrameOuter = rbind(featureFrameOuter, feaFrame)
  }
  output_hyperparam = removeEmpty(target=cbind(output_accuracy, output_hp), margin="rows")
  output_operator = removeEmpty(target=cbind(output_accuracy, output_pipelines) ,margin="rows")
}

# extract the hyper-parameters for pipelines
getHyperparam = function(Frame[Unknown] pipeline, Frame[Unknown]  hpList, Integer no_of_res)
  return (Matrix[Double] paramMatrix, Integer no_of_res, Integer NUM_META_FLAGS)
{

  allParam = 0;
  START_INDEX = 8 # value from where the hyper-params starts after skipping meta flags
  NUM_META_FLAGS = 5
  # load the hyper-parameters values
  paramList = list()
  # store the row indexes of the operator matches
  indexes = matrix(0, rows= ncol(pipeline), cols=1)
  paramCount = matrix(0, rows= ncol(pipeline), cols=1)
  for(k in 1:ncol(pipeline))
  {
    op = as.scalar(pipeline[1,k])
    hasParam = map(hpList[,1], "x->x.split(\",\")[0].equals(\""+op+"\")")
    # convert the boolean vector to 0/1 matrix representation
    m_hasParam = hasParam == frame("true", rows=nrow(hasParam), cols=1)
    m_hasParam = as.matrix(m_hasParam)
    # compute the relevant index
    index = m_hasParam * seq(1, nrow(m_hasParam))
    index = as.scalar(removeEmpty(target = index, margin = "rows"))
    indexes[k] = index
    paramCount[k] = as.integer(as.scalar(hpList[index, 2]))
  }
  # if there are no hyper-parameters than change the values of resources
  # so that the pipeline is only executed once and no resource are wasted, saving looping
  no_of_res = ifelse(sum(paramCount) > 0, no_of_res, 1)
  # the below matrix stores the different combinations of hyper-parameter value for each pipeline
  # if the resource value is greater than zero this means for 1 pipeline it will store r rows where each row store set
  # of hyperparameter values for ith pipeline. If resource value rv = 10 and ncol(pip) = 3 then the output matrix will have
  # 10*3= 30 rows and 1:10 hyper-paramters for i-the pipeline 11:20 for (i+1)-th pipeline and so on
  # this matrix stores no. of hps, values of hps, and flags
  paramMatrix = matrix(0, rows=ncol(pipeline)*no_of_res, cols=max(paramCount)+NUM_META_FLAGS+1)

  for(i in 1:ncol(pipeline)) {
    index = as.scalar(indexes[i])
    no_of_param = as.integer(as.scalar(paramCount[i]))
    # extract hasY and verbose flags
    attachMask = matrix(as.scalar(hpList[index, 3]), rows=no_of_res, cols=1)
    attachFD = matrix(as.scalar(hpList[index, 4]), rows=no_of_res, cols=1)
    attachY = matrix(as.scalar(hpList[index, 5]), rows=no_of_res, cols=1)
    isVerbose = matrix(as.scalar(hpList[index, 6]), rows=no_of_res, cols=1)
    dataFlag = matrix(as.scalar(hpList[index, 7]), rows=no_of_res, cols=1)
    
    if(no_of_param > 0) {
      paramIdx = START_INDEX
      typeIdx = START_INDEX
      OpParam = matrix(0, rows=no_of_res, cols=max(paramCount))
      
      for(j in 1:no_of_param) {
        type = as.scalar(hpList[index, typeIdx])
        paramValIndex = (no_of_param) + paramIdx
        minVal =  as.scalar(hpList[index, paramValIndex])
        maxVal = as.scalar(hpList[index, paramValIndex + 1])
        if(type == "FP") {
          val = rand(rows=no_of_res, cols=1, min=minVal,max=maxVal, pdf="uniform");
          OpParam[, j] = val;
        }
        else if(type == "INT") {
          val = sample(maxVal, no_of_res, TRUE);
          less_than_min = val < minVal;
          val = (less_than_min * minVal) + val;
          OpParam[, j] = val;
        }
        else if(type == "BOOL") {
          if(maxVal == 1) {
            s = sample(2, no_of_res, TRUE);
            b = s - 1;
            OpParam[, j] = b;
          } 
          else
            OpParam[, j] = matrix(0, rows=no_of_res, cols=1)
        }
        else
          print("invalid data type")  # TODO handle string set something like {,,}
          
        paramIdx = paramIdx + 2
        typeIdx = typeIdx + 1
      }
      # hyper-parameter vector contains no. of hp, values of hp, and flag values
      OpParam = cbind(matrix(no_of_param, rows=nrow(OpParam), cols=1),OpParam, attachMask,
        attachFD, attachY, isVerbose, dataFlag)
    }
    else {
      # no hyper-parameters, so create a dummy matrix of zeros so flags are always aligned
      dummy = matrix(0, rows=no_of_res, cols=max(paramCount)+1)
      OpParam = cbind(dummy, attachMask, attachFD, attachY)
      OpParam = cbind(OpParam, isVerbose, dataFlag)
    }
    paramMatrix[((i-1)*no_of_res)+1:i*no_of_res, 1:ncol(OpParam)] = OpParam
  }
}


# extract the top k pipelines as a final result after deduplication and sorting
extractTopK = function(Frame[Unknown] pipeline, Matrix[Double] hyperparam, 
  Double testAccuracy, Integer k)
  return (Frame[Unknown] bestPipeline, Matrix[Double] bestHyperparams)
{

  idx = vectorToCsv(seq(1, ncol(pipeline)))
  jspecDC = "{ids:true, recode:["+idx+"]}";
  # OHE of categorical features
  [dpipeline, dM] = transformencode(target=pipeline, spec=jspecDC);
  # bind the pipelines and hyper-parameters into one matrix
  forDedup = cbind(dpipeline, hyperparam) 
  # perform the similarity based deduplication
  dup = mdedup(cbind(pipeline, as.frame(hyperparam)), matrix(seq(2, ncol(forDedup)), 1,
    ncol(forDedup)-1), matrix(1,1,ncol(forDedup)-1), as.matrix(1), as.matrix(1), FALSE)

  if(sum(dup) > 0)
  {
    # take out the unique tuples
    uniqueTuples =  removeEmpty(target= forDedup, margin="rows", select = (dup ==0))
    # remove the zero rows, identifiers of unique records
    dup =  removeEmpty(target = dup, margin="rows")
    # get the counts of duplicate tuples with their tuple id
    countDist = table(dup, 1) > 0
    countDist = countDist * seq(1, nrow(countDist))
    countsVal = removeEmpty(target= countDist, margin="rows")
    indexes = table(seq(1, nrow(countsVal)),countsVal,1,nrow(countsVal), cols=nrow(forDedup))

    # for each duplicate record just take the one reocrd and strip the others
    deduplicates = indexes %*% forDedup
  
    # combine the deduplicated tuples and unique tuples again
    forDedup = rbind(uniqueTuples, deduplicates)
  }
  
  # decode the pipelines
  decoded = transformdecode(target=forDedup[, 1:ncol(pipeline)], meta=dM, spec=jspecDC)
  
  # separate the pipelines and hyper-parameters
  pipeline = decoded[, 1:ncol(pipeline)]
  hyperparam = forDedup[, ncol(pipeline)+1:ncol(forDedup)]

  # sort results
  hyperparam = order(target = hyperparam, by = 1, decreasing=TRUE, index.return=FALSE)
  pipeline = frameSort(pipeline)


  # remove the row with accuracy less than test accuracy 
  mask = (hyperparam[, 1] < testAccuracy) == 0
  hyperparam = removeEmpty(target = hyperparam, margin = "rows", select = mask)
  rowIndex = ifelse(nrow(hyperparam) > k, k, nrow(hyperparam))
  # select the top k
  bestPipeline = pipeline[1:rowIndex,]
  bestHyperparams = hyperparam[1:rowIndex,]
  
}

# extract the top k pipelines for each bracket, the intermediate results
extractBracketWinners = function(Matrix[Double] pipeline, Matrix[Double] hyperparam, 
  Integer k, Frame[String] conf)
  return (Frame[Unknown] bestPipeline, Matrix[Double] bestHyperparams)
{
  # bestPipeline = frameSort(bestPipeline)
  hyperparam = order(target = hyperparam, by = 1, decreasing=TRUE, index.return=FALSE)
  pipeline = order(target = pipeline, by = 1, decreasing=TRUE, index.return=FALSE)
  rowIndex = min(k, nrow(pipeline))

  pipeline = pipeline[1:rowIndex,]
  bestHyperparams = hyperparam[1:rowIndex,]
  bestPipeline = frame(data="|", rows=nrow(pipeline), cols=ncol(conf))
  for(i in 1: nrow(pipeline)) {
    index = as.scalar(pipeline[i, 3])
    out = conf[index, 2:ncol(conf)]
    bestPipeline[i, 1] = as.frame(pipeline[i, 1])
    bestPipeline[i, 2:ncol(bestPipeline)] = out
  }  
}

###########################################################################
# The function will return the max performance by each individual pipeline
############################################################################
getMaxPerConf = function(Matrix[Double] pipelines, Double size)
return (Frame[Unknown] maxperconf)
{
  tab = removeEmpty(target=table(pipelines[, 2], pipelines[, 3], pipelines[, 1]), margin="cols")
  maxperconf = frame(0, rows=size, cols=1)
  maxperconf[1:ncol(tab),] = as.frame(t(colMaxs(tab)))
}

###########################################################################
## function to classify the data using cross validation
############################################################################
fclassify = function(Matrix[Double] X, Matrix[Double] Y, Matrix[Double] mask, Matrix[Double] MLhp,
  Double testAccuracy, Boolean isWeighted, Integer cv=3)
  return (Double accuracy, Double T)
{
 
  if(max(Y) == min(Y)) {
    print("Y contains only one class")
    accuracy = as.double(0)
  }
  else {
    print("STARTING "+cv+" CROSS VALIDATIONS")
    # do the k = 3 cross validations
    [accuracyMatrix, T] = crossV(X, Y, cv, mask, MLhp, isWeighted)
    accuracyMatrix = removeEmpty(target=accuracyMatrix, margin="rows")
    acc = colMeans(accuracyMatrix)
    accuracy = as.scalar(acc[1,1])
    print(cv +" validation accuracy "+accuracy+" in "+T+" ms\n\n")
  }
}

# # ######################################################################
# # # # Function for cross validation using hold out method
# # # # Inputs: The input dataset X, Y and the value of k validation, mask of the 
# # # # dataset for OHE of categorical columns, vector of ML hyper-parameters identified 
# # # # via gridsearch and a boolean value of (un)weighted accuracy.
# # # # Output: It return a matrix having the accuracy of each fold.
# # ######################################################################

crossV = function(Matrix[double] X, Matrix[double] y, Integer k, Matrix[Double] mask,
  Matrix[Double] MLhp, Boolean isWeighted) 
return (Matrix[Double] accuracyMatrix, Double T)
{
  t1 = time()
  accuracyMatrix = matrix(0, k, 1)
  dataList = list()
  testL = list()
  data = order(target = cbind(y, X),  by = 1, decreasing=FALSE, index.return=FALSE)
  classes = table(data[, 1], 1)
  ins_per_fold = classes/k
  start_fold = matrix(1, rows=nrow(ins_per_fold), cols=1)
  fold_idxes = cbind(start_fold, ins_per_fold)

  start_i = 0; end_i = 0; idx_fold = 1;;
  for(i in 1:k)
  {
    fold_i = matrix(0, 0, ncol(data))
    start=0; end=0; 
    for(j in 1:nrow(classes))
    {
      idx = as.scalar(classes[j, 1])
      start = end + 1;
      end = end + idx
      class_j =  data[start:end, ]
      start_i = as.scalar(fold_idxes[j, 1]);
      end_i = as.scalar(fold_idxes[j, 2])
      fold_i = rbind(fold_i, class_j[start_i:end_i, ])
    }
    dataList = append(dataList, fold_i)
    fold_idxes[, 1] = fold_idxes[, 2] + 1
    fold_idxes[, 2] += ins_per_fold
  }

  parfor(i in seq(1,k))
  {
    [trainList, hold_out] = remove(dataList, i)
    trainset = rbind(trainList)
    testset = as.matrix(hold_out)
    trainX = trainset[, 2:ncol(trainset)]
    trainy = trainset[, 1]
    testX = testset[, 2:ncol(testset)]
    testy = testset[, 1]
    beta = multiLogReg(X=trainX, Y=trainy, icpt=as.scalar(MLhp[1,1]), reg=as.scalar(MLhp[1,2]), tol= 1e-9, 
    maxi=as.scalar(MLhp[1,3]), maxii= 50, verbose=FALSE);
    [prob, yhat, a] = multiLogRegPredict(testX, beta, testy, FALSE)
    accuracy = getAccuracy(testy, yhat, isWeighted)
    accuracyMatrix[i] = accuracy
  }
  T = floor((time() - t1) / 1e+6)
}


###############################################################################################
# The function will collect the features like statistics and pipelines and accuracy 
# so that they could be used for training a model and predicting pipelines without enumeration
###############################################################################################
gatherStats = function(Matrix[Double] X, Matrix[Double] Y, Matrix[Double] mask, String target)
return (Matrix[Double] features)
{

  features = matrix(0, rows = 1, cols= 14)
  features[1, 1]=  sum(is.na(X)) # number of missing values
  X = replace(target= X, pattern = NaN, replacement = 0)
  num = removeEmpty(target=X, margin="cols", select=(mask == 0))
  # get the stats
  features[1, 2] =  min(num) # minimum value
  features[1, 3] = max(num)
  features[1, 4] = mean(colMins(num)) # average minimum value
  features[1, 5] = mean(colMaxs(num)) # average maximum value
  features[1, 6] = sum(mask) # number of categorical features
  features[1, 7] = sum(mask == 0) # number of numerical features
  features[1, 8] = mean(num) # mean value
  colSd = colSds(num)
  count3sdplus = sum(num > (colMeans(num) + 3*colSd ))
  count3sdminus = sum(num < (colMeans(num) - 3*colSd ))
  outliers = count3sdplus + count3sdminus
  features[1, 9] = outliers
  # OHE features 
  OHE = sum(colMaxs(X) * mask)
  features[1, 10] = OHE
  if(target != "compare")
  {
    ctab = table(Y, 1)
    features[1, 11] = nrow(ctab) # number of classes
    minCat = min(ctab) / nrow(ctab)
    maxCat = max(ctab) / nrow(ctab)
    # class imabalance 1=YES, 0=NO
    features[1, 12]= ifelse((maxCat - minCat) > 0.3, 1, 0)
  }
  else 
  {
    features[1, 11] = 0
    features[1, 12] = 0
  }
  features[1, 13] = nrow(X)
  features[1, 14] = ncol(X)
  
}


######################################################################
# # Function for cross validation using hold out method
# # Inputs: The input dataset X, Y and the value of k validation, mask of the
# # dataset for OHE of categorical columns, vector of ML hyper-parameters identified 
# # via grid-search and a boolean value of (un)weighted accuracy.
# # Output: It return a matrix having the accuracy of each fold.
######################################################################

compareValue = function(Matrix[double] dirtyX, Matrix[double] fixedX,  Matrix[Double] cleanX, Matrix[Double] mask) 
return (Double precision, Double T)
{
  t1 = time()
  DEFAULT = 404
  mv = is.na(dirtyX)
  correctionsRequired = 0
  mv = is.na(fixedX)
  dirtyX = replace(target= dirtyX, pattern=NaN, replacement=DEFAULT)
  cleanX = replace(target= cleanX, pattern=NaN, replacement=DEFAULT)
  fixedX = replace(target= fixedX, pattern=NaN, replacement=DEFAULT)
  diffCleanDirty =  sum((abs(cleanX - dirtyX) < 0.001) < 1) #sum(cleanX == dirtyX) #
  print("dirty != clean: "+diffCleanDirty)
  correctionsRequired =  (abs(cleanX - dirtyX) < 0.001) < 1#dirtyX != cleanX
  print("corrections required: "+sum(correctionsRequired))
  correctionsMade =  sum(dirtyX != fixedX)
  print("corrections made: "+correctionsMade)
  dim = nrow(dirtyX) * ncol(dirtyX) 
  match = (abs(cleanX - fixedX) < 0.001) * correctionsRequired
  print("total matches "+sum(match))
  # print("total matches \n"+toString(match))
  precision = max(0.001, sum(match) / max(1, correctionsMade))
  T = floor((time() - t1) / 1e+6)
  print("Precision: "+toString(precision) + " in "+T+" ms")


}

pipToString = function(Frame[String] F)
return (String s)
{
  s = ""
  for(i in 1:ncol(F))
    s = s + as.scalar(F[,i])+";"

}