R 语言手写差分进化算法(待完善)

R
启发式优化算法
Author

Rui

Published

September 9, 2022

文字部分回头补上

参数设置

NP = 50     # 种群的数量(染色体总数)
D = 2      # 变量的维度(基因数)
G = 500     # 最大进化代数(迭代次数)
F0 = 0.4    # 初始变异因子(缩放因子)
CR = 0.25    # 交叉概率
xmax = 10     # 变量上限
xmin = -10    # 变量下限
theta = 1e-7   # 迭代停止阈值

随机初始化种群

data = matrix(runif(NP*D, min = xmin, max = xmax), nrow = NP, ncol = D) # 初始种群

v = matrix(vector(length = NP*D), nrow = NP, ncol = D) # 变异种群(空矩阵)

u = matrix(vector(length = NP*D), nrow = NP, ncol = D) # 选择种群(空矩阵)

目标函数

Ob <- function (x) {
  # 将这个复杂的函数拆解成
  term1 <- (sin(sqrt(sum(x^2))) + 0.001) / (sqrt(sum(x^2)) + 0.001) 
  term2 <- exp(0.5*sum(cos((2*pi*x))))
  
  # 将各部分合并起来得到最终的目标函数
  term1 + term2 + 10
}
# 这是向量化的目标函数(x是向量)

适应度与适应度增量

fitness_trace <- c()
FITNESS <- apply(data, MARGIN = 1, FUN = Ob)
fitness_trace[1] <- min(FITNESS)

delta_fitness <- c()

差分进化循环

for (i in 1:G) {
  
  iteration <- i
  
  
  # 变异操作
  # 自适应变异因子
  lamda = exp(1 - G/(G+1-i))
  FN = F0*(2^lamda)
  # 提取3个不同行的染色体(行索引也不能等于m)进行变异
  temp <- c(1:NP)
  for (m in 1:NP) {
    indx <- sample(temp[-m], size = 3, replace = FALSE)
    v[m, ] <- data[indx[1], ] + FN*(data[indx[2], ] - data[indx[3], ]) # 差分变异
  }
  
  # 交叉操作
  r <- sample(1:NP, size = 1) # 整数随机数
  for (n in 1:D) {
    cr <- runif(1) # 1个介于0-1之间的随机数
    
    if (cr <= CR | n == r) {
      u[ , n] = v[ , n] # 所有个体的第n维变量
    } else {
      u[ , n] = data[ , n]
    }
  }
  
  # 边界值判定
  for (n in 1:D) {
    for (m in 1:NP) {
      if (u[m, n] < xmin | u[m, n] > xmax) {
        u[m, n] = runif(1)*(xmax - xmin) + xmin
      }
    }
  }
  
  # 选择操作
  fitness <- apply(u, MARGIN = 1, FUN = Ob)
  for (m in 1:NP) {
    if (fitness[m] > FITNESS[m]) {
      data[m, ] = u[m, ]
    }
  }
  
  
  # 监视适应度与适应度增量
  fitness_trace[i+1] <- min(FITNESS)
  delta_fitness[i] <- fitness_trace[i+1] - fitness_trace[i]
  
  FITNESS <- apply(data, MARGIN = 1, FUN = Ob)
  
  # 迭代停止判定
  #if (!is.null(delta_fitness[i])) {
    #if (delta_fitness[i] < theta) break
  #}
  
}

我很想在迭代的最后加入代码:

if (!is.null(delta_fitness[i])) {
    if (delta_fitness[i] < theta) break
  }

用以判定停止条件。即当前后两次迭代的适应度的增量(差值)小于给定的阈值 theta 时终止循环。但在几次试验后发现效果并不好:迭代不超过 10 次就停止了,而输出值并不理想。

可视化

适应度

plot(x = 0:iteration, y = fitness_trace, type = "l")

适应度增量

plot(x = 1:iteration, y = delta_fitness, type = "l")