RでGA(遺伝的アルゴリズム)

(工事中)

パッケージ

今のところ Ca' Foscali University(Venezia)のga package、Sydney工科大学のgafit packageおよび rgenoud そして、seaoが存在する。

汎用

特定用途

使い方(ga : 汎用)

パッケージのロード

library(ga)

DecodeとEncodeについて

Decode (binary文字列(gene)を数値にする)

binary2integer
binary2real

Encode (数値を binary文字列(gene)にする)

integer2binary
real2binary

扱う対象の数値範囲などを考えながら binary文字列長あるいは Decode, Encode方式を考える。

例えば、関数に与える数値が整数でもよいなら integer系の encode/decodeでよいが、 実数にしたいなら 変数の区間(from, to)や binary長(size)を定義した上で real系の encode/decodeを 使う必要がある。

メインオブジェクトの作成

関数の最小/大値を探査するためには genetic関数を使って geneticオブジェクトを 生成しないといけない。以降、このオブジェクトを使って、 新世代の生成などを実施することになる。

例えば binary文字列長が5でよいのであれば

exa <- genetic(encode=integer2binary, par.encode=list(size=5), 
decode=binary2integer, length.string=5, fitness=fit)

とする。

対象とする関数は

fit <- function(x) x^2

などとして*1定義する。

> exa
Alphabet           0 1 
Wildcard           * 
Length String      5 
Size Generation    100 
Size Population    20 
Encode             integer2binary 
Parameter Encode  
$size
[1] 5


Decode             binary2integer 
Fitness            fit 
Type of Selection  roulette 
Prob. Crossover    0.8 
Prob. Mutation     0.01 

交配確率(Prob. Crossover = 0.8)や突然変異確率(Prob. Mutation = 0.01)、 一世代中の個体数(Size Population = 20)や世代数(Size Generation = 100) などは指定しないとDefault値を充てられるので、気に入らなければ適宜 指定すること。

問題を解く

ga packageは学習用だけあって、内部でどのようなことをなしているのかを こちらに見えるようにすることも可能である。

が、今はとりあえず手っ取り早く解いてみて、このパッケージの適用性を 示してみる*2

garesult <- ga.simple(exa)
plot(garesult)

結果は plot.ga関数を使用して (上では plot(garesult)となっているが やっていることは同じこと)プロットできる。

内部の関数たち

このパッケージに附属している関数を使えば 遺伝的アルゴリズム内部で何をやっているのかを 見ることができる。いくつかを紹介しておく。

交差 (crossover)

以下は 値21の示すbinary文字列と 値13の示すそれを交差させたものである。 交差する位置ははじめは1-2の位置となっているが 2回目の試行では3-4となっている。この位置はランダムに決定される。

> x <- encode(21,exa)
> x
 1~0~1~0~1 
> y <- encode(13,exa)
> y
  1~0~1~1~0 
> crossover(x,y,exa)
  Parents:
  1|0~1~0~1 
  1|0~1~1~0 
  Offspring:
  1|0~1~1~0 
  1|0~1~0~1 
> crossover(x,y,exa)
  Parents:
  1~0~1|0~1 
  1~0~1|1~0 
  Offspring:
  1~0~1|1~0 
  1~0~1|0~1 

ちなみに、複数の位置で交差させるのは現時点では無理かと思われる。

突然変異 (mutation)

例えば一世代 5個体とする。世代中の個体の組 x は mutation関数により変異させられる。

> x
1    1    0~0~0~0~1 
2    2    1~1~1~1~0 
3    3    0~0~1~0~1 
4    4    0~1~1~1~1
5    5    0~1~1~1~0
> mutation(x,exa)
String before mutation:
  0~1~0~0~0~0~1~0~1~1~0~1~1~1~1~0~1~0~1~1~1~0~1~1~0 
String after mutation:
  0~0~0~0~1~0~1~0~0~1~0~0~1~1~1~0~1~0~1~1~1~0~1~1~0 

親の選択 (selection)

xをある世代の個体集合、fitを対象とする関数、exaをgeneticオブジェクトとして 与えると、世代の中から比較的優位な個体二つをランダムに選択して親が決定される。

selection(x,fit,exa)
$x
  0~1~1~1~0 
$y
  1~0~0~0~0 
$parents.position
[1] 17  8

彼らの子どもは前述の crossoverを使って作られる。

世代の作成 (generation)

まずは AdamとEve(たち)がいないとはじまらないので generate.population関数を使って 第一世代を作る。

gen00 <- generate.population(exa)

以降はこの第一世代を用いて generation関数を使って 新たな世代を作っていくことが可能となる。

gen01 <- generation(gen00, exa)

世代の評価 (fitness)

世代に含まれる個体を評価するのはfitness.population関数を使う

fitness.population(gen01, exa)

世代の平均値を知りたい場合は以下を用いる。

fitness.mean(gen01, exa)

なお、先に示した ga.simpleの結果のプロットのうち左図は世代の fitnessのうち最大値を実線、平均を点線で示している。

系譜 (genealogy)

系譜を取得することが可能である。 優良な個体の血統を知ることができる。 例えば、先ほど実施した simple.gaの結果のうち、3世代から(to) 5世代まで(from) の系譜を調べたいのなら以下のようにする。

この時 fromの方が toより大きくないといけないことに注意(遡る感じで)

> y <- genealogy(,offspring=1,from=5,to=3,exa)
> plot(y)

この場合の plotには、plot.genealogyが使われている(が意識する必要はない)

使い方(genalg)

用途

パッケージのロード

使い方の例

使い方(gafit : 曲線当てはめ用)

用途

探査対象関数の最小値を到達可能な範囲で捜査する。 探査対象の関数(target)および、 探索開始点(start)、thermal(突然変異発生確率)、maxiter(世代数上限)、 samples(一世代辺りの個体数)、step(step size?)を与える。

こちらは ga packageとは違って収束過程、個体の性質などは 内部で処理しており R上で見えるわけではない。

パッケージのロード

library (gafit)

使い方の例

まずは関数の定義から。 まず、リファレンスにあるように cos(θ) + sin(θ)から。

e <- expression( cos( theta ) + sin( theta ))
gafit( e, list( theta=3 ), step=0.1 )
$theta
[1] 3.922567

attr(,"score")
[1] -1.414200

この場合 thetaが 3.926991のとき -1.414214で最小値を取る。

x <- seq(3,5,by=0.001)
f <- function(x) (cos(x) + sin(x))
plot(x, f(x), type="l")
s <- sort(f(x),index=T)
x[s$ix[1]]

などで検証可能である。

なお、パラメータが二つになっても適用可能 (ga packageは基本的に一つ)

ここでは有名なRosenbrockのbanana関数に適用することにする。 この関数は (x, y) = (1, 1)で最小値 0となる。 最小値近傍で傾斜が複雑になっているため、最適化ルーチンのベンチマークに 使われることがよくある。

e_banana <- expression( 100 * ( Y - X^2 )^2 + ( 1 - X )^2 )
gafit(e_banana, list( X=0, Y=0 ), step=0.1, maxiter=100)
$X
[1] 1.000984

$Y
[1] 1.002066

attr(,"score")
[1] 1.897574e-06

この例*3では良い結果を出すことができた。

他に 複素数も使うことができる。これは ?gafitしてExampleを参照下さい。

使い方(rgenoud)

用途

最適化

パッケージのロード

library(rgenoud)

使い方の例

使い方(seao)

用途

最適化

パッケージのロード

library(seao)

使い方の例

まずは遺伝子の構成を決める。seaoの特徴的なところは 遺伝子の構成が上記のライブラリよりも自由に決定できることだと思う。

まずはgenomestruc (Genome Structureの略?)を使って 対話的に 遺伝子の構成を決めていく。

今回は 遺伝子の持つ意味を3種類とした。

> gs <- genomestruc(n.genes=3)
---- GENE  1 
Name of the gene
x
Give the minimum value
-10
Give the maximum value
10
Give the step (precision)
0.1
Maximum set to: 10 
There are 201 alleles
Optional: Give a list of names for the different alleles (seperated by a ';')

WARNING: Number of allelenames not equal to number of alleles; namevector set to NULL 
---- GENE  2 
Name of the gene
y
Give the minimum value
-10
Give the maximum value
10
Give the step (precision)
0.05
Maximum set to: 10 
There are 401 alleles
Optional: Give a list of names for the different alleles (seperated by a ';') 

WARNING: Number of allelenames not equal to number of alleles; namevector set to NULL 


---- GENE  3 
Name of the gene
flag
Give the minimum value
1
Give the maximum value
2
Give the step (precision)
1
Maximum set to: 2 
There are 2 alleles
Optional: Give a list of names for the different alleles (seperated by a ';')
f1;f2

上記で定義した遺伝子構成を表示するには以下。

> str(gs)
List of 1
 $ genes:List of 3
  ..$ :List of 5
  .. ..$ name     : chr "x"
  .. ..$ min      : num -10
  .. ..$ max      : num 10
  .. ..$ n.alleles: num 201
  .. ..$ names    : NULL
  ..$ :List of 5
  .. ..$ name     : chr "y"
  .. ..$ min      : num -10
  .. ..$ max      : num 10
  .. ..$ n.alleles: num 401
  .. ..$ names    : NULL
  ..$ :List of 5
  .. ..$ name     : chr "flag"
  .. ..$ min      : num 1
  .. ..$ max      : num 2
  .. ..$ n.alleles: num 2
  .. ..$ names    : chr [1:2] "f1" "f2" 

Gene1 (x) は -10から10まで。分割は0.1毎なので 201種。
Gene2 (y) は -10から10まで。分割は0.05毎なので401種。
Gene3 (flag) は f1 (=1)もしくは f2 (=2)である。

この構成をもとに 10個体を生成する

> g00 <- newgen(gs, n.ind=10,method="random")
> g00
# 中略
$generations[[1]]$allele
      [,1]  [,2] [,3]
 [1,]  8.5 -9.25    2
 [2,] -8.4  9.00    1
 [3,]  9.2  3.55    1
 [4,] -4.9 -6.30    2
 [5,]  3.3 -6.25    1
 [6,] -0.5 -5.45    1
 [7,] -6.7  1.05    2
 [8,] -2.1  5.10    2
 [9,]  1.8 -4.75    1
[10,] -8.8  7.25    2 

$generations[[1]]$parents
[1] "-1"     "random"

[,1]が x, [,2]がyそして、[,3]がflagになっていることがわかる。

なお、このseaoパッケージだが最適化の目的関数は自前で定義して generation (g00)に与えないといけない。

例えば

> g00$generations[[1]]$fit <- c(1,2,3,4,5,6,7,8,9,10)

以下、 g00(fit定義済)から発生させた g01と、その後 g01のfitを定義せずに g02を発生させた場合に示す挙動を書いておく。 fitがないと収束もへったくれもないので"怒る"ようである。

> g01 <- nextgen(g00,n.ind=10,gen.parent=NULL,selection=list(base="fit", rescale=0),
                 corate=90,mutation=list(base="unif", spread=1, rate=15))
> g02 <- nextgen(g01,n.ind=10,gen.parent=NULL,selection=list(base="fit", rescale=0),
                 corate=90,mutation=list(base="unif", spread=1, rate=15))
Error in if (d.maxfit - d.minfit < 1e-12) { : 
        missing value where TRUE/FALSE needed
> 

例えば以下のような目的関数を与えて100世代の計算をさせてみた。

ちなみに、g00は Generation, fitが目的関数、

g00 <- newgen(gs, n.ind=10,method="random");
fit <- function(x) ( 100 * (x[2]-x[1]^2)^2 + (1-x[1])^2 ) * x[3] * -1  + 1000;
tmpfit = c(1:10); # 初期化のつもり 
for(j in 1:100) {  # 100世代繰り返す
      for(i in 1:10){ # 10は個体数
	    tmpfit[i]=fit(g00$generations[[j]]$allele[i,]);
      }
      g00$generations[[j]]$fit <- tmpfit ;
      g00 <- nextgen(g00,n.ind=10,gen.parent=NULL,selection=list(base="fit", rescale=0),
                     corate=100,mutation=list(base="unif", spread=1, rate=15));
}

こうして 100世代まで繰り返すと

> g00$generations[[2]]
$fit
 [1] 984.724659 906.045219 984.724659 984.724659 984.724659 984.724659
 [7] 812.090438 984.724659 842.146659   2.241219 

$allele
       [,1]  [,2] [,3]
 [1,] -0.89  0.45    1
 [2,] -0.71 -0.45    1
 [3,] -0.89  0.45    1
 [4,] -0.89  0.45    1
 [5,] -0.89  0.45    1
 [6,] -0.89  0.45    1
 [7,] -0.71 -0.45    2
 [8,] -0.89  0.45    1
 [9,] -0.89 -0.45    1
[10,] -0.71 -2.65    1

... 

$generations[[100]]$fit
 [1]  -10395.395     999.830     999.830     999.830     661.280     999.830
 [7]   -5146.720   -9894.494   -3108.970 -437799.376

$generations[[100]]$allele
       [,1]  [,2] [,3]
 [1,]  3.32  0.35    1
 [2,]  0.60  0.35    1
 [3,]  0.60  0.35    1
 [4,]  0.60  0.35    1
 [5,]  0.60  2.20    1
 [6,]  0.60  0.35    1
 [7,]  0.60  8.20    1
 [8,]  2.78  0.35    2
 [9,]  0.60 -6.05    1
[10,] -8.16  0.35    1

なんとなく最大値 1000に近い値を算出できていることがわかる。

使い方(subselect)

用途

パッケージのロード

library(subselect)

使い方の例

応用例

maptools パッケージの pointLabel 関数

 地図のラベル配置にGAを利用




*1 この場合とてつもなくつまらない関数
*2 本来結果も示すべきですが 例示の関数が面白くないので今日は示しません。なんか面白そうなのを適用して改訂させて頂きますデス
*3 当然合わない場合もあります example x=1.275, y=1.623で 0.077

トップ   編集 凍結 差分 バックアップ 添付 複製 名前変更 リロード   新規 一覧 検索 最終更新   ヘルプ   最終更新のRSS
Last-modified: 2023-03-25 (土) 11:19:17