inst/doc/manual.R

In neldermead: R Port of the 'Scilab' Neldermead Module

## ----setup, include = FALSE---------------------------------------------------
knitr::opts_chunk$set(echo = TRUE)

## ---- echo = FALSE, message = FALSE--------------------------------------
options(width = 75)
require(neldermead)

## ----eval = FALSE--------------------------------------------------------
#    costf <- function(x, index, fmsfundata){
#      # Define f and c here #
#      return(list(f, g = NULL, c, gc = NULL, index = index,
#                  this = list(costfargument = fmsfundata)))
#    }

## ----eval = FALSE--------------------------------------------------------
#    myalgorithm <- function( this ){
#      ...
#      return(this)
#    }

## ----eval = FALSE--------------------------------------------------------
#    mystoppingrule <- function( this , simplex ){
#      ...
#      return(list(this = this, terminate = terminate, status = status))
#    }

## ------------------------------------------------------------------------
  quadratic <- function(x = NULL,index = NULL,fmsfundata = NULL){
    return(list(f = x[1]^2 + x[2]^2,
                g = c(),
                c = c(), 
                gc = c(),
                index = index,
                this = list(costfargument = fmsfundata)))
  }

  x0 <- transpose( c(1.0,1.0) )
  nm <- neldermead()
  nm <- neldermead.set(nm,'numberofvariables',2)
  nm <- neldermead.set(nm,'function',quadratic)
  nm <- neldermead.set(nm,'x0',x0)
  nm <- neldermead.search(nm)
  summary(nm)

## ------------------------------------------------------------------------
  rosenbrock <- function(x = NULL,index = NULL,fmsfundata = NULL){
    return(list(f = 100*(x[2]-x[1]^2)^2+(1-x[1])^2,
                g = c(),
                c = c(),
                gc = c(),
                index = index,
                this = list(costfargument = fmsfundata)))
  }
  x0 <- transpose(c(-1.2,1.0))
  nm <- neldermead()
  nm <- neldermead.set(nm,'numberofvariables',2)
  nm <- neldermead.set(nm,'function',rosenbrock)
  nm <- neldermead.set(nm,'x0',x0)
  nm <- neldermead.set(nm,'maxiter',200)
  nm <- neldermead.set(nm,'maxfunevals',300)
  nm <- neldermead.set(nm,'tolfunrelative',10*.Machine$double.eps)
  nm <- neldermead.set(nm,'tolxrelative',10*.Machine$double.eps)
  nm <- neldermead.set(nm,'simplex0method','axes')
  nm <- neldermead.set(nm,'simplex0length',1.0)
  nm <- neldermead.set(nm,'method','variable')
  nm <- neldermead.set(nm,'verbose',FALSE)
  nm <- neldermead.set(nm,'storehistory',TRUE)
  nm <- neldermead.set(nm,'verbosetermination',FALSE)
  nm <- neldermead.search(nm)
  
  xmin <- ymin <- -2.0 
  xmax <- ymax <- 2.0 
  nx <- ny <- 100
  stepy <- stepx <- (xmax - xmin)/nx
  ydata <- xdata <- seq(xmin,xmax,stepx)
  zdata <- apply(expand.grid(xdata,ydata),1,
                 function(x) neldermead.function(nm,transpose(x)))
  zdata <- matrix(zdata,ncol = length(ydata))
  optimpath <- matrix(unlist((neldermead.get(nm,'historyxopt'))),
                      nrow = 2)
  optimpath <- data.frame(x = optimpath[1,],y = optimpath[2,])

  contour(xdata,ydata,zdata,levels = c(1,10,100,500,1000,2000))
  par(new = TRUE,ann = TRUE)
  plot(c(x0[1],optimpath$x[158]), c(x0[2],optimpath$y[158]),
       col = c('red','green'),pch = 16,xlab = 'x[1]',ylab = 'x[2]',
       xlim = c(xmin,xmax),ylim = c(ymin,ymax))
  par(new = TRUE,ann = FALSE)  
  plot(optimpath$x,optimpath$y,col = 'blue',type = 'l',
       xlim = c(xmin,xmax),ylim = c(ymin,ymax))

## ------------------------------------------------------------------------
  quadratic <- function(x = NULL,index = NULL,fmsfundata = NULL){
    return(list(f = x[1]^2 + x[2]^2,
                g = c(),
                c = c(), 
                gc = c(),
                index = index,
                this = list(costfargument = fmsfundata)))
  }
  set.seed(0)
  x0 <- transpose(c(1.2,1.9))
  nm <- neldermead()
  nm <- neldermead.set(nm,'numberofvariables',2)
  nm <- neldermead.set(nm,'function',quadratic)
  nm <- neldermead.set(nm,'x0',x0)
  nm <- neldermead.set(nm,'verbose',FALSE)
  nm <- neldermead.set(nm,'storehistory',TRUE)
  nm <- neldermead.set(nm,'verbosetermination',FALSE)
  nm <- neldermead.set(nm,'method','box')
  nm <- neldermead.set(nm,'boundsmin',c(1,1))
  nm <- neldermead.set(nm,'boundsmax',c(2,2))
  nm <- neldermead.search(nm)
  summary(nm)

## ------------------------------------------------------------------------
  michalewicz <- function(x = NULL,index = NULL,fmsfundata = NULL){
    f <- c()
    c <- c()
    if (index == 2 | index == 6) 
      f <- (x[1]-10)^3+(x[2]-20)^3
    
    if (index == 5 | index == 6)
      c <- c((x[1]-5)^2+(x[2]-5)^2 -100, 
          82.81-((x[1]-6)^2+(x[2]-5)^2))
    varargout <- list(f = f,
        g = c(),
        c = c, 
        gc = c(),
        index = index,
        this = list(costfargument = fmsfundata))
    return(varargout)
  }
  set.seed(0)
  x0 <- transpose(c(15,4.99))
  nm <- neldermead()
  nm <- neldermead.set(nm,'numberofvariables',2)
  nm <- neldermead.set(nm,'nbineqconst',2)
  nm <- neldermead.set(nm,'function',michalewicz)
  nm <- neldermead.set(nm,'x0',x0)
  nm <- neldermead.set(nm,'maxiter',300)
  nm <- neldermead.set(nm,'maxfunevals',1000)
  nm <- neldermead.set(nm,'simplex0method','randbounds')
  nm <- neldermead.set(nm,'boxnbpoints',3)
  nm <- neldermead.set(nm,'storehistory',TRUE)
  nm <- neldermead.set(nm,'method','box')
  nm <- neldermead.set(nm,'boundsmin',c(13,0))
  nm <- neldermead.set(nm,'boundsmax',c(20,10))
  nm <- neldermead.search(nm)
  summary(nm)

## ------------------------------------------------------------------------
  negLL <- function(x = NULL, index = NULL, fmsfundata = NULL){
    mn <- x[1]
    sdv <- x[2]
    out <- -sum(dnorm(fmsfundata$data, mean = mn, sd = sdv, log = TRUE))
  
    return(list(f = out, 
           index = index,
           this = list(costfargument = fmsfundata)))
  }

  set.seed(12345)
  fmsfundata <- structure(
    list(data = rnorm(500,mean = 50,sd = 2)),
    class = 'optimbase.functionargs')

  x0 <- transpose(c(45,3))
  nm <- neldermead()
  nm <- neldermead.set(nm,'numberofvariables',2)
  nm <- neldermead.set(nm,'function',negLL)
  nm <- neldermead.set(nm,'x0',x0)
  nm <- neldermead.set(nm,'costfargument',fmsfundata)
  nm <- neldermead.set(nm,'maxiter',500)
  nm <- neldermead.set(nm,'maxfunevals',1500)
  nm <- neldermead.set(nm,'method','box')
  nm <- neldermead.set(nm,'storehistory',TRUE)
  nm <- neldermead.set(nm,'boundsmin',c(-100, 0))
  nm <- neldermead.set(nm,'boundsmax',c(100, 100))
  nm <- neldermead.search(this = nm)
  summary(nm)

## ------------------------------------------------------------------------
  rosenbrock <- function(x = NULL){
    f <- 100*(x[2]-x[1]^2)^2+(1-x[1])^2
  }
  x0 <- c(-1.2,1.0)
  npts <- 6
  xmin <- c(-2,-2)
  xmax <- c(2,2)
  grid <- fmin.gridsearch(fun = rosenbrock,x0 = x0,xmin = xmin,xmax = xmax,npts = npts,alpha = alpha)
  grid