Rで遺伝的アルゴリズムの勉強

遺伝的アルゴリズム(GA)をPackageを使わずにRでスクリプトを作成する例です。
基本的なGAの流れに沿って、各ステップを関数として実装し、最後にそれらを組み合わせてGAを実行するスクリプトを作成しました。

1. 目的関数の設定

まずは、GAで最適化したい目的関数を定義します。
ここでは、簡単な例として、与えられたバイナリベクトル（0と1の配列）中の1の数を最大化する という問題を考えます。

# 目的関数：バイナリベクトル中の1の数を数える
objective_function <- function(chromosome) {
  sum(chromosome)
}

2. パラメータ設定

GAの動作を制御するためのパラメータを設定します。

# パラメータ設定
population_size <- 50       # 個体数
chromosome_length <- 20     # 染色体の長さ（バイナリベクトルの長さ）
generations <- 100          # 世代数
mutation_rate <- 0.01       # 突然変異率
crossover_rate <- 0.7       # 交叉率

3. 個体群の初期化

ランダムなバイナリベクトルを生成して、初期個体群を作成します。

# 初期個体群の生成
initialize_population <- function(population_size, chromosome_length) {
  population <- list()
  for (i in 1:population_size) {
    # ランダムな0と1のベクトルを生成
    chromosome <- sample(c(0, 1), chromosome_length, replace = TRUE)
    population[[i]] <- chromosome
  }
  population
}

population <- initialize_population(population_size, chromosome_length)

4. 評価関数

各個体の適応度を評価する関数です。ここでは、objective_function をそのまま使います。

# 適応度評価
evaluate_fitness <- function(population, objective_function) {
  fitness_scores <- numeric(length(population))
  for (i in 1:length(population)) {
    fitness_scores[i] <- objective_function(population[[i]])
  }
  fitness_scores
}

5. 選択

次世代に残す個体を選択します。ここでは、ルーレット選択 を実装します。
ルーレット選択は、適応度が高い個体ほど選択されやすいように確率を設定する方式です。

# ルーレット選択
selection <- function(population, fitness_scores) {
  # 適応度に基づいて選択確率を計算
  selection_probabilities <- fitness_scores / sum(fitness_scores)
  # 次世代の個体群を選択
  parent_indices <- sample(1:length(population), 
                            size = length(population), 
                            replace = TRUE, 
                            prob = selection_probabilities)
  
  # 選択された親個体から次世代の個体群を生成
  next_generation <- list()
  for(i in parent_indices){
    next_generation[[length(next_generation) + 1]] <- population[[i]]
  }
  next_generation
}

6. 交叉

選択された個体間で交叉を行い、新たな個体を生成します。ここでは、一点交叉 を実装します。

# 一点交叉
crossover <- function(parent1, parent2, crossover_rate) {
  if (runif(1) < crossover_rate) {
    # 交叉点をランダムに決定
    crossover_point <- sample(1:(length(parent1) - 1), 1)
    # 交叉を実行
    child1 <- c(parent1[1:crossover_point], parent2[(crossover_point + 1):length(parent2)])
    child2 <- c(parent2[1:crossover_point], parent1[(crossover_point + 1):length(parent1)])
    return(list(child1, child2))
  } else {
    # 交叉を行わない場合は親をそのまま返す
    return(list(parent1, parent2))
  }
}

7. 突然変異

生成された個体に突然変異を導入します。ここでは、ビット反転突然変異 を実装します。

# ビット反転突然変異
mutation <- function(chromosome, mutation_rate) {
  for (i in 1:length(chromosome)) {
    if (runif(1) < mutation_rate) {
      # ビットを反転
      chromosome[i] <- 1 - chromosome[i]
    }
  }
  chromosome
}

8. GAメインループ

上記の関数を組み合わせて、GAのメインループを実装します。

# 遺伝的アルゴリズムの実行
genetic_algorithm <- function(population_size, chromosome_length, generations, 
                              mutation_rate, crossover_rate, objective_function) {
  
  population <- initialize_population(population_size, chromosome_length)
  best_fitness_history <- numeric(generations) # 各世代の最良のfitnessを記録
  best_chromosome <- NULL # 最良の染色体
  max_fitness <- -Inf # 最大fitnessの初期値

  for (generation in 1:generations) {
    fitness_scores <- evaluate_fitness(population, objective_function)
    
    # 各世代の最良のfitnessを記録
    current_best_fitness_index <- which.max(fitness_scores)
    current_max_fitness <- fitness_scores[current_best_fitness_index]
    best_fitness_history[generation] <- current_max_fitness
    
    # 最良個体の更新
    if (current_max_fitness > max_fitness) {
      max_fitness <- current_max_fitness
      best_chromosome <- population[[current_best_fitness_index]]
    }

    # 選択
    parents <- selection(population, fitness_scores)
    
    # 交叉と突然変異
    next_generation <- list()
    for (i in seq(1, length(parents), by = 2)) {
      # 偶数個体群を処理するため、範囲チェック
      if (i+1 <= length(parents)) {
        parent1 <- parents[[i]]
        parent2 <- parents[[i+1]]
        
        # 交叉
        children <- crossover(parent1, parent2, crossover_rate)
        child1 <- children[[1]]
        child2 <- children[[2]]
        
        # 突然変異
        child1 <- mutation(child1, mutation_rate)
        child2 <- mutation(child2, mutation_rate)
        
        next_generation[[length(next_generation) + 1]] <- child1
        next_generation[[length(next_generation) + 1]] <- child2
      } else {
        # 親の数が奇数の場合、最後の親はそのまま次世代へ（または何らかの処理を追加）
        next_generation[[length(next_generation) + 1]] <- mutation(parents[[i]], mutation_rate) # 突然変異のみ適用
      }
    }
    
    population <- next_generation
  }
  
  # 結果をリストで返す：最良の染色体とfitnessの履歴
  list(best_chromosome = best_chromosome, 
       best_fitness = max_fitness,
       fitness_history = best_fitness_history)
}

9. GAの実行と結果の確認

作成した genetic_algorithm 関数を実行し、結果を確認します。

# GA実行
result <- genetic_algorithm(population_size, chromosome_length, generations, 
                            mutation_rate, crossover_rate, objective_function)

# 結果の表示
cat("最良の染色体:", result$best_chromosome, "\n")
cat("最良の適応度:", result$best_fitness, "\n")

# Fitnessの履歴をプロット
plot(1:generations, result$fitness_history, type="l", 
     xlab="世代", ylab="最良適応度", 
     main="遺伝的アルゴリズムの最適化過程")

スクリプト全体

上記のコードを全てまとめたスクリプトは以下のようになります。

# 1. 目的関数の設定
objective_function <- function(chromosome) {
  sum(chromosome)
}

# 2. パラメータ設定
population_size <- 50
chromosome_length <- 20
generations <- 100
mutation_rate <- 0.01
crossover_rate <- 0.7

# 3. 個体群の初期化
initialize_population <- function(population_size, chromosome_length) {
  population <- list()
  for (i in 1:population_size) {
    chromosome <- sample(c(0, 1), chromosome_length, replace = TRUE)
    population[[i]] <- chromosome
  }
  population
}

# 4. 評価関数
evaluate_fitness <- function(population, objective_function) {
  fitness_scores <- numeric(length(population))
  for (i in 1:length(population)) {
    fitness_scores[i] <- objective_function(population[[i]])
  }
  fitness_scores
}

# 5. 選択 (ルーレット選択)
selection <- function(population, fitness_scores) {
  selection_probabilities <- fitness_scores / sum(fitness_scores)
  parent_indices <- sample(1:length(population), 
                            size = length(population), 
                            replace = TRUE, 
                            prob = selection_probabilities)
  next_generation <- list()
  for(i in parent_indices){
    next_generation[[length(next_generation) + 1]] <- population[[i]]
  }
  next_generation
}

# 6. 交叉 (一点交叉)
crossover <- function(parent1, parent2, crossover_rate) {
  if (runif(1) < crossover_rate) {
    crossover_point <- sample(1:(length(parent1) - 1), 1)
    child1 <- c(parent1[1:crossover_point], parent2[(crossover_point + 1):length(parent2)])
    child2 <- c(parent2[1:crossover_point], parent1[(crossover_point + 1):length(parent1)])
    return(list(child1, child2))
  } else {
    return(list(parent1, parent2))
  }
}

# 7. 突然変異 (ビット反転突然変異)
mutation <- function(chromosome, mutation_rate) {
  for (i in 1:length(chromosome)) {
    if (runif(1) < mutation_rate) {
      chromosome[i] <- 1 - chromosome[i]
    }
  }
  chromosome
}

# 8. GAメインループ
genetic_algorithm <- function(population_size, chromosome_length, generations, 
                              mutation_rate, crossover_rate, objective_function) {
  population <- initialize_population(population_size, chromosome_length)
  best_fitness_history <- numeric(generations)
  best_chromosome <- NULL
  max_fitness <- -Inf

  for (generation in 1:generations) {
    fitness_scores <- evaluate_fitness(population, objective_function)
    current_best_fitness_index <- which.max(fitness_scores)
    current_max_fitness <- fitness_scores[current_best_fitness_index]
    best_fitness_history[generation] <- current_max_fitness
    
    if (current_max_fitness > max_fitness) {
      max_fitness <- current_max_fitness
      best_chromosome <- population[[current_best_fitness_index]]
    }

    parents <- selection(population, fitness_scores)
    next_generation <- list()
    for (i in seq(1, length(parents), by = 2)) {
      if (i+1 <= length(parents)) {
        parent1 <- parents[[i]]
        parent2 <- parents[[i+1]]
        children <- crossover(parent1, parent2, crossover_rate)
        child1 <- children[[1]]
        child2 <- children[[2]]
        child1 <- mutation(child1, mutation_rate)
        child2 <- mutation(child2, mutation_rate)
        next_generation[[length(next_generation) + 1]] <- child1
        next_generation[[length(next_generation) + 1]] <- child2
      } else {
        next_generation[[length(next_generation) + 1]] <- mutation(parents[[i]], mutation_rate)
      }
    }
    population <- next_generation
  }
  list(best_chromosome = best_chromosome, 
       best_fitness = max_fitness,
       fitness_history = best_fitness_history)
}

# 9. GAの実行と結果の確認
result <- genetic_algorithm(population_size, chromosome_length, generations, 
                            mutation_rate, crossover_rate, objective_function)

cat("最良の染色体:", result$best_chromosome, "\n")
cat("最良の適応度:", result$best_fitness, "\n")

plot(1:generations, result$fitness_history, type="l", 
     xlab="世代", ylab="最良適応度", 
     main="遺伝的アルゴリズムの最適化過程")

実行方法

Rスクリプトとして保存: 上記のコードを ga_script.R などの名前で保存します。
Rで実行: RコンソールまたはRStudioで、保存したスクリプトのディレクトリに移動し、以下を実行します。
```
source("ga_script.R")
```

解説

目的関数 (objective_function): 最適化したい問題を関数として定義します。ここではバイナリベクトルの1の数を最大化する簡単な例です。
パラメータ設定: population_size, chromosome_length, generations, mutation_rate, crossover_rate などのGAの動作を調整するパラメータを設定します。
初期化 (initialize_population): 最初にランダムな個体群を生成します。
評価 (evaluate_fitness): 各個体の適応度を計算します。
選択 (selection): ルーレット選択で次世代の親個体を選択します。適応度が高いほど選択されやすくなります。
交叉 (crossover): 一点交叉により、親個体から新しい個体を生成します。
突然変異 (mutation): ビット反転突然変異により、個体の遺伝子をランダムに変化させます。
GAメインループ (genetic_algorithm): これらのステップを世代数だけ繰り返します。各世代で、評価、選択、交叉、突然変異を行い、個体群を進化させていきます。
結果表示: 最終的に得られた最良の個体（染色体）とその適応度、および世代ごとの最良適応度の変化をプロットで表示します。

補足

Packageを使わない理由: Packageを使わずにスクリプトを作成することで、遺伝的アルゴリズムの基本的な仕組みをより深く理解することができます。
カスタマイズ: このスクリプトは基本的なGAの枠組みを提供しています。目的関数、選択方法、交叉方法、突然変異方法などを変更することで、様々な問題に対応できます。
発展:
- 選択方法の変更: ルーレット選択以外にも、トーナメント選択、ランキング選択など様々な選択方法があります。
- 交叉方法の変更: 一点交叉以外にも、二点交叉、一様交叉などがあります。
- 突然変異方法の変更: ビット反転突然変異以外にも、挿入、削除、逆位などの突然変異 оператор があります。
- エリート戦略: 各世代で最も適応度の高い個体を次世代に必ず残すエリート戦略を導入することで、解の探索性能が向上する場合があります。
- より複雑な問題への適用: このスクリプトを基に、より複雑な最適化問題（ナップサック問題、巡回セールスマン問題など）に挑戦してみてください。

このスクリプトが、Rで遺伝的アルゴリズムを学ぶための一助となれば幸いです。

このブログを検索/Search this blog

Translate

石川昌和 | Masakazu Ishikawa | いしかわまさかず