GenSA和SA对Knapsack问题给出了无意义的输出。

Question

我有以下CSV。

Knapsack.CSV

,gewichten(gr),waarde
Voorwerp 1,70,135
Voorwerp 2,73,139
Voorwerp 3,77,149
Voorwerp 4,80,150
Voorwerp 5,82,156
Voorwerp 6,87,163
Voorwerp 7,90,173
Voorwerp 8,94,184
Voorwerp 9,98,192
Voorwerp 10,106,201
Voorwerp 11,110,210
Voorwerp 12,113,214
Voorwerp 13,115,221
Voorwerp 14,118,229
Voorwerp 15,120,240

我正试图用GenSA和GA来解决Knapsack问题。这组数据的解法应该是1458左右。

然而，用这个代码。

install.packages("GenSA")
install.packages("GA")
require(GenSA)
library(GenSA)
require(GA)
library(GA)

#Loading data
df <- read.csv("knapsack.csv", header=TRUE, sep=",")


#Define function
knapsack <- function(x) {
  f <- sum(x * df[3])
  penalty <- sum(df[2]) * abs(sum(x*df[2]) - 750)
  f - penalty
}

init <- runif(1, -5000, 5000)

onder <- rep(-5000, length(init))
boven <- rep(5000, length(init)) 



controlelijst <- list(max.time=25, nb.stop.improvement = 100)

resultaatSA <- GenSA(par=init, lower = onder, upper = boven, fn=knapsack, control=controlelijst)

resultaatSA$par




# Solution num 2
SGA <- ga(type="binary", fitness=knapsack, nBits=length(df[1]), maxiter=150, run=250, popSize=100, seed=101)

SGA
SGA@solution

我得到了很多无意义的输出。比如GenSA说解是-5000，或者有时是5000。这些都是我设置的边界约束（boundariesconstraints）。

SA给出的解是1。

我到底做错了什么，我需要如何正确使用这两个函数？

Answer 1

我认为John Coleman是正确的，你必须简单地改变成本函数的符号，使其最小化。下面是一个（略显夸张）的例子，说明在你的函数上加一个小减号将如何导致一个非常不同的解决方案（希望是正确的）。我以类似的方式绘制了结果，如下图所示 GA 的，因为我从你的代码中看到你也在研究这个选项（事实上是一个很好的选项）。对于 GA，成本函数为 最大化所以你需要去掉减号。

library(GenSA)
library(GA)

df <- read.table(text = ",gewichten(gr),waarde
Voorwerp 1,70,135
Voorwerp 2,73,139
Voorwerp 3,77,149
Voorwerp 4,80,150
Voorwerp 5,82,156
Voorwerp 6,87,163
Voorwerp 7,90,173
Voorwerp 8,94,184
Voorwerp 9,98,192
Voorwerp 10,106,201
Voorwerp 11,110,210
Voorwerp 12,113,214
Voorwerp 13,115,221
Voorwerp 14,118,229
Voorwerp 15,120,240", sep = ",", header = T
)


#Define function
knapsack <- function(x) {
  f <- sum(x * df[3])
  penalty <- sum(df[2]) * abs(sum(x*df[2]) - 750)
  -(f - penalty) # SIMPLY ADDED A MINUS SIGN
}

init <- runif(1, -5000, 5000)

onder <- rep(-5000, length(init))
boven <- rep(5000, length(init)) 


controlelijst <- list(max.time=25, nb.stop.improvement = 100)

resultaatSA <- GenSA(par=init, lower = onder, upper = boven, fn=knapsack, control=controlelijst)

resultaatSA$par # 0.5233775
head(resultaatSA$trace.mat)


# summarize results
tmp <- as.data.frame(resultaatSA$trace.mat)
meani <- aggregate(tmp$function.value, list(step = tmp$nb.steps),mean, na.rm = TRUE)
exe <- aggregate(tmp$current.minimum, list(step = tmp$nb.steps),mean, na.rm = TRUE)
medi <- aggregate(tmp$function.value, list(step = tmp$nb.steps),median, na.rm = TRUE)
ylim <- c(min(range(exe$x,na.rm = TRUE, finite = TRUE)),
          max(range(meani$x, na.rm = TRUE, finite = TRUE)))

# plot
op <- par(mar=c(5.1, 4.1, 1, 4.1))
plot(tmp$nb.steps, tmp$function.value, type = "n", ylim = ylim, xlab = "Iteration",
     ylab = "Cost value")
graphics::grid(equilogs = FALSE)
points(tmp$nb.steps, tmp$current.minimum, type = "o", pch = 16, lty = 1,
       col = "green3", cex = 0.7)
points(meani$step, meani$x, type = "o", pch = 1, lty = 2,
       col = "dodgerblue3", cex = 0.7)
polygon(c(meani$step, rev(meani$step)),
        c(exe$x, rev(medi$x)),
        border = FALSE, col = adjustcolor("green3", alpha.f = 0.1))
par(new=TRUE)
plot(tmp$nb.steps, tmp$temperature, t="l", col=2, lty=2, log="y", axes = FALSE, xlab = "", ylab = "")
axis(4, col=2, col.axis=2); mtext(text = "Temperature", side = 4, line = par()$mgp[1], col=2)
legend("topright", legend = c("Best", "Mean", "Median", "Temperature"),
       col = c("green3", "dodgerblue3", adjustcolor("green3", alpha.f = 0.1), 2),
       pch = c(16, 1, NA, NA), lty = c(1,2,1,2),
       lwd = c(1, 1, 10, 1), pt.cex = c(rep(0.7,2), 2, NA),
       inset = 0.02)
par(op)

Answer 2

在GA调用中，你通过了错误的 nBits 如上所述 用户2957945

基于。用简单的遗传算法解决背包问题；我得到了两个解决方案，利润为1449，重量为750。

编辑我得到了一个1456个利润和750个重量的解决方案。

library(GA)
# --------------------------------------------------------------------
# Read Data
# --------------------------------------------------------------------
my_df <- read.table(text=',gewichten(gr),waarde
Voorwerp 1,70,135
Voorwerp 2,73,139
Voorwerp 3,77,149
Voorwerp 4,80,150
Voorwerp 5,82,156
Voorwerp 6,87,163
Voorwerp 7,90,173
Voorwerp 8,94,184
Voorwerp 9,98,192
Voorwerp 10,106,201
Voorwerp 11,110,210
Voorwerp 12,113,214
Voorwerp 13,115,221
Voorwerp 14,118,229
Voorwerp 15,120,240', sep=',', header=T)

# --------------------------------------------------------------------
# Define  profit, weights, Knapsack limit, and fitness function
# --------------------------------------------------------------------
p <- my_df$waarde
w <- my_df$gewichten.gr.
W <- 750
n <- length(p)

# Define fitness function 
knapsack <- function(x) { 
  f <- sum(x * p) 
  penalty <- sum(w) * abs(sum(x * w) - W) 
  f - penalty 
}

# --------------------------------------------------------------------
# Run SGA
# --------------------------------------------------------------------
SGA <- ga(type="binary", 
          fitness=knapsack , 
          nBits=n, 
          maxiter=500, # Maximum number of generations 
          run=200,     # Stop if the best-so-far fitness
          # hasn't improved for 'run' generations 
          popSize=200)


# --------------------------------------------------------------------
# see results
# --------------------------------------------------------------------
x.star <- SGA@solution 
# solutions
x.star
# number of elements in each solution
rowSums(x.star) 
# profit in each solution
rowSums(sweep(x.star, MARGIN=2, p, `*`))
# weight of each solution
rowSums(sweep(x.star, MARGIN=2, w, `*`))

关于 GenSA 函数，我不是专家，但我认为它是为了更复杂的优化，而且它能使函数最小化（不同于 ga 的函数最大化），所以你不能对两种方法使用相同的函数。另一个需要注意的是，knapsack问题通常被认为是一个二进制问题（选择0，1为nBits），但我不确定你是否可以强求 GenSA 来实现，所以你需要找到一个变通的方法。

下面是我的尝试 GenSA，返回1456利润和750权重的解，我用了 round 作为获取二进制输出的变通方法。

library(GenSA)

knapsack_gensa <- function(x) { 
  f <- sum(round(x) * p)
  penalty <- sum(w) * abs(sum(round(x) * w) - W) 
  penalty - f 
}
gensa <- GenSA(lower = rep(0, n), 
               upper = rep(1, n), 
               fn=knapsack_gensa)

solution <- round(gensa$par)
sum(solution)
sum(solution*p)
sum(solution*w)