code/draft_R_pkg/R/main.R

   1 #' @include defaults.R
   2
   3 #' @title Cluster power curves with PAM in parallel
   4 #'
   5 #' @description Groups electricity power curves (or any series of similar nature) by applying PAM
   6 #' algorithm in parallel to chunks of size \code{nbSeriesPerChunk}
   7 #'
   8 #' @param data Access to the data, which can be of one of the three following types:
   9 #' \itemize{
  10 #'   \item data.frame: each line contains its ID in the first cell, and all values after
  11 #'   \item connection: any R connection object (e.g. a file) providing lines as described above
  12 #'   \item function: a custom way to retrieve the curves; it has two arguments: the start index
  13 #'     (start) and number of curves (n); see example in package vignette.
  14 #' }
  15 #' @param K Number of clusters
  16 #' @param nbSeriesPerChunk Number of series in each group
  17 #' @param writeTmp Function to write temporary wavelets coefficients (+ identifiers);
  18 #'   see defaults in defaults.R
  19 #' @param readTmp Function to read temporary wavelets coefficients (see defaults.R)
  20 #' @param WER "end" to apply stage 2 after stage 1 has iterated and finished, or "mix"
  21 #'   to apply it after every stage 1
  22 #' @param ncores number of parallel processes; if NULL, use parallel::detectCores()
  23 #'
  24 #' @return A data.frame of the final medoids curves (identifiers + values)
  25 epclust = function(data, K, nbSeriesPerChunk, writeTmp=ref_writeTmp, readTmp=ref_readTmp,
  26     WER="end", ncores=NULL)
  27 {
  28     #TODO: setRefClass(...) to avoid copy data:
  29     #http://stackoverflow.com/questions/2603184/r-pass-by-reference
  30
  31     #0) check arguments
  32     if (!is.data.frame(data) && !is.function(data))
  33         tryCatch({dataCon = open(data)},
  34             error="data should be a data.frame, a function or a valid connection")
  35     if (!is.integer(K) || K < 2)
  36         stop("K should be an integer greater or equal to 2")
  37     if (!is.integer(nbSeriesPerChunk) || nbSeriesPerChunk < K)
  38         stop("nbSeriesPerChunk should be an integer greater or equal to K")
  39     if (!is.function(writeTmp) || !is.function(readTmp))
  40         stop("read/writeTmp should be functional (see defaults.R)")
  41     if (WER!="end" && WER!="mix")
  42         stop("WER takes values in {'end','mix'}")
  43     #concerning ncores, any non-integer type will be treated as "use parallel:detectCores()"
  44
  45     #1) acquire data (process curves, get as coeffs)
  46     index = 1
  47     nbCurves = nrow(data)
  48     while (index < nbCurves)
  49     {
  50         if (is.data.frame(data))
  51         {
  52             #full data matrix
  53             writeTmp( getCoeffs( data[index:(min(index+nbSeriesPerChunk-1,nbCurves)),] ) )
  54         } else if (is.function(data))
  55         {
  56             #custom user function to retrieve next n curves, probably to read from DB
  57             writeTmp( getCoeffs( data(index, nbSeriesPerChunk) ) )
  58         } else
  59         {
  60             #incremental connection
  61             #TODO: find a better way to parse than using a temp file
  62             ascii_lines = readLines(dataCon, nbSeriesPerChunk)
  63             seriesChunkFile = ".tmp/seriesChunk"
  64             writeLines(ascii_lines, seriesChunkFile)
  65             writeTmp( getCoeffs( read.csv(seriesChunkFile) ) )
  66         }
  67         index = index + nbSeriesPerChunk
  68     }
  69     if (exists(dataCon))
  70         close(dataCon)
  71
  72     library(parallel)
  73     ncores = ifelse(is.integer(ncores), ncores, parallel::detectCores())
  74     cl = parallel::makeCluster(ncores)
  75     parallel::clusterExport(cl=cl, varlist=c("X", "Y", "K", "p"), envir=environment())
  76     li = parallel::parLapply(cl, 1:B, getParamsAtIndex)
  77
  78     #2) process coeffs (by nbSeriesPerChunk) and cluster in parallel (just launch async task, wait for them to complete, and re-do if necessary)
  79     repeat
  80     {
  81         completed = rep(FALSE, ............)
  82         #while there is jobs to do (i.e. size of tmp "file" is greater than nbSeriesPerChunk),
  83         #A) determine which tasks which processor will do (OK)
  84         #B) send each (sets of) tasks in parallel
  85         #C) flush tmp file (current parallel processes will write in it)
  86         #always check "complete" flag (array, as I did in MPI) to know if "slaves" finished
  87     }
  88
  89     parallel::stopCluster(cl)
  90
  91     #3) readTmp last results, apply PAM on it, and return medoids + identifiers
  92
  93     #4) apply stage 2 (in parallel ? inside task 2) ?)
  94     if (WER == "end")
  95     {
  96         #from center curves, apply stage 2...
  97     }
  98 }
  99
 100 getCoeffs = function(series)
 101 {
 102     #... return wavelets coeffs : compute in parallel !
 103 }