Create a gist now

Instantly share code, notes, and snippets.

What would you like to do?
Nelder and Mead 1965 for R
##### definition #####
Nelder_mead <- function(x1, lambda=1, obj_func, alpha=1, gamma=2, rho=1/2, sigma=1/2){
# dimmension
dim_x <- length(x1)
# generate n point using given x
# and summarize n+1 point (x1, ... , xn+1)
simplex_point <- cbind(x1, x1 + lambda*diag(dim_x))
# initial criterion value
criterion <- 1
# iteration number
iter <- 1
while(!criterion<10^{-8}){
# calculate function values
if(iter==1){
f_value <- apply(simplex_point, 2, obj_func)
asc_num <- order(f_value)
asc_point <- simplex_point[,asc_num]
}else{
f_value <- apply(asc_point, 2, obj_func)
asc_num <- order(f_value)
asc_point <- asc_point[,asc_num]
}
# calculate centroid without max(f(xi))
x_cent <- {1/dim_x}*{apply(asc_point[,-{dim_x+1}], 1, sum)}
# calculate refrection
x_ref <- x_cent + alpha*(x_cent - asc_point[,{dim_x+1}])
f_1 <- obj_func(asc_point[,1])
f_ref <- obj_func(x_ref)
f_n <- obj_func(asc_point[,{dim_x}])
f_np1 <- obj_func(asc_point[,{dim_x+1}])
if({f_1<=f_ref}&&{f_ref<f_n}){
## replace by reflection
asc_point[,{dim_x+1}] <- x_ref
}else if(f_ref < f_1){
## replace by expansion
x_exp <- x_cent + gamma*(x_cent - asc_point[,{dim_x+1}])
f_exp <- obj_func(x_exp)
if(f_exp < f_ref){
asc_point[,{dim_x+1}] <- x_exp
}else{
asc_point[,{dim_x+1}] <- x_ref
}
}else{
## replace by contraction
x_cont <- x_cent + rho*(x_cent - asc_point[,{dim_x+1}])
f_cont <- obj_func(x_cont)
if(f_cont < f_np1){
asc_point[,{dim_x+1}] <- x_cont
}else{
tmp1 <- asc_point[,1]
tmp2 <- asc_point[,1] + sigma*(asc_point[,-1] - asc_point[,1])
asc_point <- cbind(tmp1, tmp2)
}
}
f_value <- apply(asc_point, 2, obj_func)
criterion <- sqrt(sum({f_value-mean(f_value)}^2))
min_num <- which.min(apply(asc_point, 2, obj_func))
x_min <- asc_point[,min_num]
iter <- iter + 1
# print(c(x_min, iter))
}
return(list(x=x_min, iteration=round(iter)))
}
##### definition #####
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment