LSopt: Stochastic Local Search

Description Usage Arguments Details Value Author(s) References See Also Examples

Description

Performs a simple stochastic Local Search.

Usage

1
LSopt(OF, algo = list(), ...)

Arguments

OF

The objective function, to be minimised. Its first argument needs to be a solution; ... arguments are also passed.

algo

List of settings. See Details.

...

Other variables to be passed to the objective function and to the neighbourhood function. See Details.

Details

Local Search (LS) changes an initial solution for a number of times, accepting only such changes that lead to an improvement in solution quality (as measured by the objective function OF). More specifically, in each iteration, a current solution xc is changed through a function algo$neighbour. This function takes xc as an argument and returns a new solution xn. If xn is not worse than xc, ie, if OF(xn,...)<=OF(xc,...), then xn replaces xc.

The list algo contains the following items:

nS

The number of steps. The default is 1000; but this setting depends very much on the problem.

x0

The initial solution. This can be a function; it will then be called once without arguments to compute an initial solution, ie, x0 <- algo$x0(). This can be useful when LSopt is called in a loop of restarts and each restart is to have its own starting value.

neighbour

The neighbourhood function, called as neighbour(x, ...). Its first argument must be a solution x; it must return a changed solution.

printDetail

If TRUE (the default), information is printed. If an integer i greater then one, information is printed at very ith step.

printBar

If TRUE (the default), a txtProgressBar (from package utils) is printed).

storeF

if TRUE (the default), the objective function values for every solution in every generation are stored and returned as matrix Fmat.

storeSolutions

default is FALSE. If TRUE, the solutions (ie, decision variables) in every generation are stored and returned in list xlist (see Value section below). To check, for instance, the current solution at the end of the ith generation, retrieve xlist[[c(2L, i)]].

At the minimum, algo needs to contain an initial solution x0 and a neighbour function.

LS works on solutions through the functions neighbour and OF, which are specified by the user. Thus, a solution need not be a numeric vector, but can be any other data structure as well (eg, a list or a matrix).

To run silently (except for warnings and errors), algo$printDetail and algo$printBar must be FALSE.

Value

A list:

xbest

best solution found.

OFvalue

objective function value associated with best solution.

Fmat

a matrix with two columns. Fmat[ ,1L] contains the proposed solution over all iterations; Fmat[ ,2L] contains the accepted solutions.

xlist

if algo$storeSolutions is TRUE, a list; else NA. Contains the neighbour solutions at a given iteration (xn) and the current solutions (xc). Example: Fmat[i, 2L] is the objective function value associated with xlist[[c(2L, i)]].

initial.state

the value of .Random.seed when the function was called.

Author(s)

Enrico Schumann

References

Gilli, M., Maringer, D. and Schumann, E. (2011) Numerical Methods and Optimization in Finance. Elsevier. http://www.elsevierdirect.com/product.jsp?isbn=9780123756626

Schumann, E. (2013) The NMOF Manual. http://enricoschumann.net/NMOF.htm

See Also

TAopt, restartOpt

Examples

  1
  2
  3
  4
  5
  6
  7
  8
  9
 10
 11
 12
 13
 14
 15
 16
 17
 18
 19
 20
 21
 22
 23
 24
 25
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
## Aim: find the columns of X that, when summed, give y

## random data set
nc <- 25L          ## number of columns in data set
nr <- 5L           ## number of rows in data set
howManyCols <- 5L  ## length of true solution
X <- array(runif(nr*nc), dim = c(nr, nc))
xTRUE <- sample(1L:nc, howManyCols)
Xt <- X[ , xTRUE, drop = FALSE]
y <- rowSums(Xt)

## a random solution x0 ...
makeRandomSol <- function(nc) {
    ii <- sample.int(nc, sample.int(nc, 1L))
    x0 <- logical(nc); x0[ii] <- TRUE
    x0
}
x0 <- makeRandomSol(nc)

## ... but probably not a good one
sum(y - rowSums(X[ , xTRUE, drop = FALSE])) ## should be 0
sum(y - rowSums(X[ , x0, drop = FALSE])) 

## a neighbourhood function: switch n elements in solution
neighbour <- function(xc, Data) {
    xn <- xc
    p <- sample.int(Data$nc, Data$n)
    xn[p] <- !xn[p]
    if (sum(xn) < 1L)
        xn <- xc
    xn
}

## a greedy neighbourhood function
neighbourG <- function(xc, Data) {
    of <- function(x)
        abs(sum(Data$y - rowSums(Data$X[ ,x, drop = FALSE])))
    xbest <- xc
    Fxbest <- of(xbest)
    for (i in 1L:Data$nc) {
        xn <- xc; p <- i
        xn[p] <- !xn[p]
        if (sum(xn) >= 1L) {
            Fxn <- of(xn)
            if (Fxn < Fxbest) {
                xbest <- xn
                Fxbest <- Fxn
            }
        }
    }
    xbest
}

## an objective function
OF <- function(xn, Data)
    abs(sum(Data$y - rowSums(Data$X[ ,xn, drop = FALSE])))


## (1) GREEDY SEARCH
## note: this could be done in a simpler fashion, but the
##       redundancies/overhead here are small, and the example is to
##       show how LSopt can be used for such a search
Data <- list(X = X, y = y, nc = nc, nr = nr, n = 1L)
algo <- list(nS = 500L, neighbour = neighbourG, x0 = x0,
             printBar = FALSE, printDetail = FALSE)
solG <- LSopt(OF, algo = algo, Data = Data)

## after how many iterations did we stop?
iterG <- min(which(solG$Fmat[ ,2L] == solG$OFvalue))
solG$OFvalue  ## the true solution has OF-value 0

## (2) LOCAL SEARCH
algo$neighbour <- neighbour
solLS <- LSopt(OF, algo = algo, Data = Data)
iterLS <- min(which(solLS$Fmat[ ,2L] == solLS$OFvalue))
solLS$OFvalue  ## the true solution has OF-value 0

## (3) *Threshold Accepting*
algo$nT <- 10L
algo$nS <- ceiling(algo$nS/algo$nT)
algo$q <- 0.99
solTA <- TAopt(OF, algo = algo, Data = Data)
iterTA <- min(which(solTA$Fmat[ ,2L] == solTA$OFvalue))
solTA$OFvalue  ## the true solution has OF-value 0

## look at the solution
all <- sort(unique(c(which(solTA$xbest),
                     which(solLS$xbest),
                     which(solG$xbest),
                     xTRUE)))
ta <- ls <- greedy <- true <- character(length(all))
true[  match(xTRUE, all)] <- "o"
greedy[match(which(solG$xbest),  all)] <- "o"
ls[    match(which(solLS$xbest), all)] <- "o"
ta[    match(which(solTA$xbest), all)] <- "o"
data.frame(true = true, greedy = greedy, LS = ls , TA = ta,
           row.names=all)

## plot results
par(ylog = TRUE, mar = c(5,5,1,6), las = 1)
plot(solTA$Fmat[seq_len(iterTA) ,2L],type = "l", log = "y",
     ylim = c(1e-4,
              max(pretty(c(solG$Fmat,solLS$Fmat,solTA$Fmat)))),
     xlab = "iterations", ylab = "OF value", col = grey(0.5))
lines(cummin(solTA$Fmat[seq_len(iterTA), 2L]), type = "l")
lines(solG$Fmat[ seq_len(iterG),  2L], type = "p", col = "blue")
lines(solLS$Fmat[seq_len(iterLS), 2L], type = "l", col = "goldenrod3")
legend(x = "bottomleft",
       legend = c("TA best solution", "TA current solution",
                  "Greedy", "LS current/best solution"),
       lty = c(1,1,0,1),
       col = c("black",grey(0.5),"blue","goldenrod2"),
       pch = c(NA,NA,21,NA))
axis(4, at = c(solG$OFvalue, solLS$OFvalue, solTA$OFvalue),
        labels = NULL, las = 1)
lines(x = c(iterG, par()$usr[2L]), y = rep(solG$OFvalue,2),
      col = "blue", lty = 3)
lines(x = c(iterTA, par()$usr[2L]), y = rep(solTA$OFvalue,2),
      col = "black", lty = 3)
lines(x = c(iterLS, par()$usr[2L]), y = rep(solLS$OFvalue,2),
      col = "goldenrod3", lty = 3)

Example output

[1] 0
[1] -33.59116
[1] 0.309139
[1] 0.2291871
[1] 0.002299774
   true greedy LS TA
1     o      o     o
2                  o
3            o  o  o
4            o      
6               o   
7            o      
8                  o
10    o         o   
11              o   
12    o             
14           o  o   
15                 o
17           o     o
22           o  o  o
23    o      o      
25    o             

NMOF documentation built on May 2, 2019, 6:39 p.m.