# R/algorithms.R In zoomgrid: Grid Search Algorithm with a Zoom

#################################################################################
## functions implementing the algorithms
#################################################################################

#' Build the grid for the grid search algorithm with a zoom.
#'
#' This function builds the grid for the grid search algorithm with a zoom.
#'
#' The argument \code{...} is a sequence of vectors or lists contatining the information about the grid to be built.
#' Each element in the sequence is either a vector or a list taking one of the following forms
#'
#' - \code{x}, if \code{x} is already an sequence of the grid points for the corresponding argument.
#'
#' - \code{c(from=, to=, by=)}
#'
#' - \code{c(from=, to=, length=)}
#'
#' - \code{list(from=, to=, by=)}
#'
#' - \code{list(from=, to=, length=)}
#'
#' where
#'
#' - \code{from}: the min of the argument of the target function
#'
#' - \code{to}: the max of the argument of the target function
#'
#' - \code{by}: the increment of the sequence
#'
#' - \code{length}: desired length.
#'
#'
#' There are many different ways to organize the points on the grid for certain argument of the target function,
#' the user can make them freely and input directly by \code{build_grid(x, ...)}.
#' Notice that \code{x} does not need to be increasing, as the function will sort it.
#' The design that \code{x} does not need to be increasing makes it convenient for the user
#' to interpolate more points at some region without considering to sort it all the time.
#'
#' When \code{by} is provided, the \code{length} will be ignored.
#' So if the user wanna specify the \code{length}, please do not use \code{by}.
#'
#' The order of the sequence \code{...} matters as it represents the order of the corresponding arguments of the target function to be optimized.
#'
#' @param ... a sequence of vectors or lists containing the information about the grid to be built, see Usage and Details.
#'
#' @return a new object of the class GRID with the grid ready for the grid search with a zoom.
#'
#' The object contains the following components:
#' \item{grid}{the grid}
#' \item{size}{number of points in the grid}
#' \item{npar}{number of arguments or parameters}
#'
#' @author Yukai Yang, \email{yukai.yang@@statistik.uu.se}
#' @keywords algorithms
#'
#' @examples
#' vx = 1:5
#' build_grid(vx, c(from=1, to=2, by=.2), list(from=3, to=4, length=5))
#'
#' @export
build_grid <- function(...){
args = list(...)
ret = list()

grid_base = list(); iter = 1
for(arg in args){
if(is.null(arg)) next

if(is.list(arg)){
nfrom = arg$from; if(is.null(nfrom)) next nto = arg$to; if(is.null(nto)) next
if(nfrom > nto)
stop(simpleError("'from' is bigger than 'to'! Please check what you have input..."))

nby = arg$by if(!is.null(nby)){ if(nby <= 0) stop(simpleError("'by' is not positive...")) # with from, to, by tmp = seq(from=nfrom, to=nto, by=nby) }else{ nlength = arg$length; if(is.null(nlength)) next
if(nlength <= 0) stop(simpleError("'length' is not positive..."))

# with from, to, length
tmp = seq(from=nfrom, to=nto, length.out=nlength)
}
}else if(is.vector(arg)){
if(!is.numeric(arg)) next

nfrom = arg["from"]
if(is.na(nfrom)){
# x case
tmp = sort(arg)
}else{
nto = arg["to"]; if(is.na(nto)) next
if(nfrom > nto)
stop(simpleError("'from' is bigger than 'to'! Please check what you have input..."))

nby = arg["by"]
if(!is.na(nby)){
if(nby <= 0) stop(simpleError("'by' is not positive..."))

# with from, to, by
tmp = seq(from=nfrom, to=nto, by=nby)
}else{
nlength = arg["length"]; if(is.na(nlength)) next
if(nlength <= 0) stop(simpleError("'length' is not positive..."))

# with from, to, length
tmp = seq(from=nfrom, to=nto, length.out=nlength)
}
}
}else next

grid_base[[iter]] = tmp; iter = iter+1
}

ret$grid = expand.grid(grid_base, KEEP.OUT.ATTRS=FALSE) ret$grid_base = grid_base
ret$size = nrow(ret$grid)
ret$npar = length(grid_base) class(ret) = "GRID" return(ret) } #' Carry out the grid search algorithm with a zoom. #' #' This function carries out the grid search algorithm with a zoom. #' #' The target function \code{FUN} to be minimized is a scalar real valued function with multiple arguments. #' The maximization can be achieved by multiplying -1 to the original function, #' and then input the new function to \code{FUN}. #' #' The \code{grid} must be created by using the function \code{\link{build_grid}} function in the package. #' #' Any other invariant arugments to the function \code{FUN} can be specified in \code{MoreArgs} by using a list with variable names. #' #' The common grid search first build a grid within some bounded area. #' And then the target function will be evaluated at each point in the grid. #' The points that produce the smallest \code{num} target function values shall be returned. #' #' \code{zoom = 0} implies that no zoom-in will be applied and the corresponding grid search is the common one introduced above, #' while any integer \code{zoom > 0} implies that the zoom-in will be applied in the grid search. #' When a grid search with a zoom is applied, \code{zoom > 0} is actually the number of rounds or layers of the algorithm, #' and therefore the grid search algorithm with a zoom consists of #' \deqn{ n^{0} + n^{1} + n^{2} + ... + n^{z} } #' grid searches, where \eqn{n} = \code{num} and \eqn{z} = \code{zoom}. #' #' As mentioned above, in each grid search, \code{num} is the number of points that will be returned. #' And therefore, in the end, there will be #' \deqn{ n^{1} + n^{2} + n^{3} + ... + n^{z+1} } #' points returned, where \eqn{n} = \code{num} and \eqn{z} = \code{zoom}. #' #' Each time when the algorithm zooms in, it will automatically build subgrids #' based on the points that have been found in the super gird search. #' Due to the exhaustive property of the grid search algorithm, #' it is desirable to make fewer points in the subgrid. #' The decay rate \code{decay} provides the opportunity to control the number of points in the subgrids. #' The number of points for each argument of the target function in the subgrid will be #' \code{max} [ \code{Int} (\code{decay} \eqn{*} \eqn{N}), 3 ]. #' #' Parallel computation is implemented in the function, which can be activated by setting \code{parallel = TRUE}. #' #' \code{cores}, which represents the number of cores, works only when \code{parallel = TRUE}. #' By default \code{cores=NULL} implies that the function will detect the number of cores and use it. #' #' The boolean \code{silent} controls if there will be output in the console. #' #' #' @param FUN the target function to be minimized. #' @param grid an object of the class GRID from \code{\link{build_grid}}. #' @param MoreArgs a list of other arguments to \code{FUN}, see \code{\link{mapply}}. #' @param zoom number of (additional) rounds or layers of the zoom-in, 0 by default. #' @param decay a number in between 0 and 1 representing the decay rate of the grid sizes of the zoom. #' @param num number of points to return, i.e. the smallest \code{num} points, 1 by default the minimum. #' @param parallel a boolean indicating if the parallel computation is carried out, by default \code{FALSE}. #' @param cores The number of cores to use, i.e. at most how many child processes will be run simultaneously. For details, see \code{mcmapply} in \code{parallel} package. #' @param silent a boolean indicating if the information regarding the omputation is printed. #' #' @return a list containing the results from the grid search with a zoom. #' #' The list contains the following components: #' \item{par}{the approximate global minimizer} #' \item{points}{all the local minimizer points found by the grid search with a zoom} #' #' @author Yukai Yang, \email{yukai.yang@@statistik.uu.se} #' @seealso \code{\link{build_grid}}, \code{\link{grid_search_check}} #' @keywords algorithms #' #' @examples #' # Rastrigin function #' ndim = 2 # number of dimension #' nA = 10 # parameter A #' # vx in [-5.12, 5.12] #' #' # minimizer = rep(0, ndim) #' # minimum = 0 #' Rastrigin <- function(vx) return(nA * ndim + sum(vx*vx - nA * cos(2*pi*vx))) #' #' \donttest{ #' # set seed and initialize the initial or starting value #' set.seed(1) #' par = runif(ndim, -5.12, 5.12) #' cat("start from", par) #' #' # results from different optimization algorithms #' optim(par = par, Rastrigin, method='Nelder-Mead') #' optim(par = par, Rastrigin, method='BFGS') #' optim(par = par, Rastrigin, method='L-BFGS-B') #' optim(par = par, Rastrigin, method='SANN') #' } #' #' # a toy example #' # build the grid first #' bin = c(from=-5.12, to=5.12, by=.5) #' grid = build_grid(bin,bin) #' # so this is a relatively sparse grid #' #' # serial computation #' ret0 = grid_search(Rastrigin, grid, silent=FALSE) #' ret0$par
#'
#' \donttest{
#'
#' # We can build a finer grid
#' bin = c(from=-5.12, to=5.12, by=.1)
#' grid = build_grid(bin,bin)
#'
#' # serial computation
#' ret1 = grid_search(Rastrigin, grid, silent=FALSE)
#' ret1$par #' #' # parallel computation #' ret2 = grid_search(Rastrigin, grid, num=2, parallel=TRUE, silent=FALSE) #' ret2$par
#'
#' # grid search with a zoom!
#' ret3 = grid_search(Rastrigin, grid, zoom=2, num=2, parallel=TRUE, silent=FALSE)
#' ret3$par #' } #' #' @export grid_search <- function(FUN, grid, MoreArgs=NULL, zoom=0, decay=0.5, num=1, parallel=FALSE, cores=NULL, silent=TRUE){ if(class(grid)!="GRID") stop(simpleError("The argument 'grid' is not an object of class 'GRID'.")) if(!silent){ cat0(paste0(rep("#",getOption("width")),collapse='')) cat0("zoomgrid version ",vnum," ",packname) cat0(paste0(rep("-",getOption("width")),collapse='')) } # tmp = split(grid$grid, seq_len(grid$size)) if(is.null(cores)){ cores = parallel::detectCores() if(is.na(cores)){ if(!silent) cat0("No cores are detected! Anyway, let's use one core...") parallel=FALSE } } if(parallel){ ftmp = grid_peval if(!silent) cat0("Parallel computation runs with ",cores," cores.") }else ftmp = grid_seval if(!silent) ptm = proc.time() ret = recursive_search(ftmp=ftmp,FUN=FUN,grid=grid,MoreArgs=MoreArgs,zoom=zoom,decay=decay,num=num,cores=cores) if(!silent) cat0("The Grid Search of ",zoom," zoom-in layers with ",num," points each gives ",length(ret), " results.") par = ftmp(FUN=FUN,grid=ret,MoreArgs=MoreArgs,num=1,cores=cores)[[1]] if(!silent) cat0("The minimizer is believed to be in the neighbourhood of ",par,".") if(!silent){ cat0(paste0(rep("-",getOption("width")),collapse='')) print(proc.time() - ptm) cat0(paste0(rep("#",getOption("width")),collapse='')) } return(list(par=par, points=ret)) } #' Check the time consumed by running the grid search algorithm with a zoom. #' #' This function checks the time consumed by running the grid search algorithm with a zoom as well as some other conditions. #' #' The running of this function takes only several seconds. #' So it is recommended to run this function before \code{\link{grid_search}} to check the approximate time consumed #' by \code{\link{grid_search}} by using exactly the same arguments. #' #' This function is extremely useful when the user is going to run \code{\link{grid_search}} on some super-computing server #' and need to know approximately how long time it will take in order to specify the corresponding settings #' according to some batch system like SLURM for example. #' #' The boolean \code{silent} controls if there will be output in the console. #' #' For details, see \code{\link{grid_search}}. #' #' @param FUN the target function to be minimized. #' @param grid an object of the class GRID from \code{\link{build_grid}}. #' @param MoreArgs a list of other arguments to \code{FUN}, see \code{\link{mapply}}. #' @param zoom number of (additional) rounds or layers of the zoom-in, 0 by default. #' @param decay a number in between 0 and 1 representing the decay rate of the grid sizes of the zoom. #' @param num number of points to return, i.e. the smallest \code{num} points, 1 by default the minimum. #' @param parallel a boolean indicating if the parallel computation is carried out, by default \code{FALSE}. #' @param cores The number of cores to use, i.e. at most how many child processes will be run simultaneously. For details, see \code{mcmapply} in \code{parallel} package. #' @param silent a boolean indicating if the information regarding the computation is printed. #' #' @return a number of the time in seconds. #' #' @author Yukai Yang, \email{yukai.yang@@statistik.uu.se} #' @seealso \code{\link{build_grid}}, \code{\link{grid_search}} #' @keywords algorithms #' #' @examples #' # Rastrigin function #' ndim = 2 # number of dimension #' nA = 10 # parameter A #' # vx in [-5.12, 5.12] #' #' # minimizer = rep(0, ndim) #' # minimum = 0 #' Rastrigin <- function(vx) return(nA * ndim + sum(vx*vx - nA * cos(2*pi*vx))) #' #' # a toy example #' # build the grid first #' bin = c(from=-5.12, to=5.12, by=.5) #' grid = build_grid(bin,bin) #' # so this is a relatively sparse grid #' #' # serial computation #' ret0 = grid_search(Rastrigin, grid, silent=FALSE) #' ret0$par
#'
#' \donttest{
#' # If we expand the grid to allow for more points
#' bin = c(from=-5.12, to=5.12, by=.1)
#' grid = build_grid(bin,bin)
#'
#' # run the check before the grid search
#' ret1 = grid_search_check(Rastrigin, grid, silent=FALSE)
#' ret1 = grid_search(Rastrigin, grid, silent=FALSE)
#' }
#'
#' @export
grid_search_check <- function(FUN, grid, MoreArgs=NULL, zoom=0, decay=0.5, num=1, parallel=FALSE, cores=NULL, silent=TRUE){
if(class(grid)!="GRID")
stop(simpleError("The argument 'grid' is not an object of class 'GRID'."))

if(!silent){
cat0(paste0(rep("#",getOption("width")),collapse=''))
cat0("zoomgrid version ",vnum," ",packname)
cat0(paste0(rep("-",getOption("width")),collapse=''))
}

if(is.null(cores)){
cores = parallel::detectCores()
if(is.na(cores)){
if(!silent) cat0("No cores are detected! Anyway, let's use one core...")
parallel=FALSE
}
}

if(parallel){
ftmp = grid_peval
if(!silent) cat0("Parallel computation runs with ",cores," cores.")
}else ftmp = grid_seval

set.seed(1)

# test 1
np = 2000**(1/grid$npar) tmp = NULL for(base in grid$grid_base){
tmp = paste(tmp, paste0("c(",paste(sort(sample(base,min(np,length(base)))),collapse=','),")"), sep=", ")
}
tmp = substr(tmp, 3, nchar(tmp))
new_grid = eval(parse(text=paste0("build_grid(",tmp,")")))

ptm = proc.time()
ret = recursive_search(ftmp=ftmp,FUN=FUN,grid=new_grid,MoreArgs=MoreArgs,zoom=0,num=num,cores=cores)
ptm = proc.time() - ptm

tx1 = ptm[3]; nn1 = new_grid$size/cores # test 2 np = 3000**(1/grid$npar)
tmp = NULL
for(base in grid$grid_base){ tmp = paste(tmp, paste0("c(",paste(sort(sample(base,min(np,length(base)))),collapse=','),")"), sep=", ") } tmp = substr(tmp, 3, nchar(tmp)) new_grid = eval(parse(text=paste0("build_grid(",tmp,")"))) ptm = proc.time() ret = recursive_search(ftmp=ftmp,FUN=FUN,grid=new_grid,MoreArgs=MoreArgs,zoom=0,num=num,cores=cores) ptm = proc.time() - ptm tx2 = ptm[3]; nn2 = new_grid$size/cores

# forecast the time

tmp = solve(matrix(c(nn1,nn2,1,1),2,2), c(tx1,tx2))
rr = tmp[1]; aa = tmp[2]
nn = grid$size/cores tmp = 0:zoom ret = sum( (nn*(decay**(tmp*grid$npar)) * rr + aa) * (num**tmp) )
if(!silent){
cat0("The expected time consumed by running the grid search is around ",ret," seconds.")
cat0(paste0(rep("#",getOption("width")),collapse=''))
}

return(ret)
}


## Try the zoomgrid package in your browser

Any scripts or data that you put into this service are public.

zoomgrid documentation built on May 2, 2019, 2:41 a.m.