Excelで進化計算

※こちらは，会社の技術ブログとのクロスポスト記事です．元の記事 ➡ [(1)](https://www.timedia.co.jp/tech/20220606-tech)，[(2)](https://www.timedia.co.jp/tech/20220608-tech)，[(3)](https://www.timedia.co.jp/tech/20220623-tech)，[(4)](https://www.timedia.co.jp/tech/20230112-tech)

この記事では，Microsoft社の製品であるExcelの進化計算を紹介します．
* Excelで進化計算を実装しようという話ではありません
* VBAは利用せず，プログラムも書きません
* 初めから，Excelには進化計算が実装されています
* 進化計算は，Excelのソルバー機能の1つです

Excelの隠れた機能である進化計算を呼び出し，実行する方法を紹介します．
Windows 10，Excel 2016 環境で説明しますが，他の環境でも実行可能なはずです．

# (1) Excelで進化計算を動かす
まずは，Excelで進化計算を動かしてみましょう．

## ソルバーを有効にする
最初に，以下の手順でソルバーを有効にします．
1. Excelを起動し，空白のブックを開きます．
2. 「ファイル」，「オプション」，「アドイン」を順番にクリックします．
3. 「ソルバーアドイン」を選択し，「設定(G)…」をクリックします．
![](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/308259/dbd05a17-d240-8d16-e4fd-0c4ccdff7c99.png)
4. 「ソルバーアドイン」にチェック☑を入れて，「OK」をクリックします．
![](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/308259/6e44877b-b53f-3d40-8210-98937c893b8f.png)
5. 「データ」をクリックしたとき，一番右に「ソルバー」「分析」の項目があることを確認します．
![](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/308259/524aa2d2-220c-6a10-aefb-1f8896652637.png)

## 問題を準備する
ソルバーで解くための最適化問題を準備します．
A2のセルに変数の下限値「-5.12」，B2のセルに変数の上限値「5.12」を書き込みます．
E2のセルに目的関数「=20+C2^2+D2^2-10\*(COS(2\*PI()\*C2)+COS(2\*PI()\*D2))」を書きます．
これは，2変数の[Rastrigin function](https://en.wikipedia.org/wiki/Rastrigin_function)と呼ばれる最適化のベンチマーク問題です．
![](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/308259/6d64ca86-3a44-7437-475d-5f5db2b01433.png)
Rastrigin functionは，次式で定義されます
　　$\text{Minimize  }\qquad f(\boldsymbol{x}) = 10n + \sum_{i=1}^n [x_i^2 - 10\cos(2\pi x_i)]$  
　　$\text{Subject to}\quad-5.12\leq x_i \leq 5.12,\quad$ for $i = 1,2,\dots,n.$

この問題は，目的関数値$f(\boldsymbol{x})$が最少となる変数$\boldsymbol{x}$を見出すことが目標となります．
$\boldsymbol{x}=(0,0,\dots,0)$のとき，最小値（最適値）$f(\boldsymbol{x})=0$をとります．

今回は，変数の次元数$n$＝2，$x_1$＝C2セルの値，$x_2$＝D2セルの値，$f(\boldsymbol{x})$＝E2セルの値としています．

変数の次元数$n$＝1のとき，変数$\boldsymbol{x}$と目的関数値$f(\boldsymbol{x})$は，以下の関係性になります．
![](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/308259/c3b64ccc-d57d-5fc7-6da0-e43531bde650.png)
変数の次元数$n$＝2のとき，変数$\boldsymbol{x}$と目的関数値$f(\boldsymbol{x})$は，以下の関係性になります．
![](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/308259/4ea3af55-dcc8-9e26-f959-82ad56df26eb.png)

## 進化計算を実行する
実際に進化計算を実行します．
1. 「データ」，「ソルバー」を順番にクリックします．
2. 以下の画像を参考に，「目的セル」，「目標値」，「変数セル」，「制約条件」を設定します．
![](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/308259/389d4257-50e2-13ab-9ad7-08bcedad6ed4.png)
3. 解決方法をエボリューショナリー (= evolutionary)に変更したうえで，「解決(S)」をクリックします．

以上で進化計算が実行されます．90秒くらい待つと，実行が終了します．

手元の実行では，極めて最適解に近い値を獲得できました.
上の図，C2，D2，E2セルの値が最適化の実行結果です．
変数は小数第6位，目的関数値は小数第10位まで，最適解と同じ値を得ることができました．

# (2) 3つのエンジンと3つの最適化問題
ここまで，Excelで進化計算を実行する方法を紹介しました．
ここからは，ソルバーのパラメーターについて紹介します．
変数$\boldsymbol{x}$と目的関数値$f(\boldsymbol{x})$の関係性を表す画像は，MATLABで作成しています．

## 3つのエンジン
ソルバーのパラメーターには，「解決方法の選択」というタブがあり，3つのエンジンの選択肢があります．
「GRG 非線形」「シンプレックス LP」「エボリューショナリー」の3つです．
エンジンの日本語，英語対応は以下になります．
* GRG 非線形 = Generalized reduced gradient method, nonlinear = 一般化簡約勾配法，非線形
* シンプレックス LP = Simplex method, linear programming = 単体法，線型計画法
* エボリューショナリー = Evolutionary computation, evolutionary algorithm = 進化計算，進化アルゴリズム

問題側から見た場合，説明の記載通りにエンジンを選択することが望ましいです．すなわち，
* 滑らかな非線形性を示す問題　　　 → GRG 非線形 エンジン
* 線形を示す問題　　　　　　　　　 → シンプレックス LP エンジン
* 滑らかではない非線形性を示す問題 → エボリューショナリー エンジン
という選択になります．

では，問題と最適化エンジンの組み合わせを変えるとどうなるか見ていきましょう．

## 滑らかではない非線形性を示す問題
先ほどの[Rastrigin function](https://en.wikipedia.org/wiki/Rastrigin_function)は，滑らかではない非線形性を示す問題になります．
![](https://github.com/tomtkg/Utility-Functions/blob/main/Data/function/rastrigin.png?raw=true)
この問題に対して「GRG 非線形 エンジン」を適用した場合，下図のような，役に立たない解を獲得します．
C2セル=「3」，D2セル=「4」の状態でエンジンを実行したため，その近傍の局所最適解が結果になっています．
（一応，C2セル=D2セル=「0付近」でエンジンを実行すれば良い解を獲得できます）
![](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/308259/68cdbdad-295e-958c-0aaf-8d8935e0f4e0.png)
この問題に対して「シンプレックス LP エンジン」を適用した場合，下図のエラーが出力されます．
![](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/308259/b844d903-3146-75a5-2950-5b3951025d82.png)
すなわち，滑らかではない非線形性を示す問題は，**「エボリューショナリー エンジン」しか対応できません．**

## 滑らかな非線形性を示す問題
問題として，[Sphere function](https://en.wikipedia.org/wiki/Test_functions_for_optimization)を取り上げます．E2のセルを「=C2^2+D2^2」に書き換えます．
![](https://github.com/tomtkg/Utility-Functions/blob/main/Data/function/sphere.png?raw=true)
この問題に対して各エンジンを適用した場合，
* GRG 非線形 エンジン　　　　　→ 最適解 (変数1=変数2=目的関数値=0) を瞬時に獲得します
* シンプレックス LP エンジン　　→ エラーが出力されます
* エボリューショナリー エンジン → 下図の結果が得られます
![excel23.PNG](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/308259/a3272232-4186-fe4a-4523-cddddba42ab7.png)
実行時間および精度の面で，「GRG 非線形 エンジン」が適切な選択になります．
一方で，90秒程度の実行時間，10万分の1以下の誤差を無視できる場合，**「エボリューショナリー エンジン」でも対応可能です．**

## 線形を示す問題
問題として，単純な足し算を取り上げます．E2のセルを「=C2+D2」に書き換えます．
![](https://github.com/tomtkg/Utility-Functions/blob/main/Data/function/linear.png?raw=true)
この問題に対して各エンジンを適用した場合，
* GRG 非線形 エンジン　　　　　→ 最適解 (変数1=変数2=-5.12, 目的関数値=-10.24) を瞬時に獲得します
* シンプレックス LP エンジン　　→ 最適解 (変数1=変数2=-5.12, 目的関数値=-10.24) を瞬時に獲得します
* エボリューショナリー エンジン → 最適解 (変数1=変数2=-5.12, 目的関数値=-10.24) を獲得します
精度の面で，3つのエンジンに差は見られませんでした．
この問題においても，90秒程度の実行時間と僅かな誤差を無視できる場合，**「エボリューショナリー エンジン」で対応可能です．**

# ここまでのまとめ
Excelで進化計算を実行する方法，Excelソルバーの3つのエンジンと3つの最適化問題を紹介しました．

最初に取り上げた最適化問題[Rastrigin function](https://en.wikipedia.org/wiki/Rastrigin_function)は，エボリューショナリー エンジンでしか解くことができませんでした．
実世界の最適化問題も滑らかでも線形でもない複雑な問題であることが多く，エボリューショナリー エンジンしか解けないかもしれません．
エボリューショナリー エンジンは，今回紹介した3つの最適化問題を全て解くことができました．
最適化問題に出会ったら，「最初に進化計算を実行して，結果を確認する」というのが良い戦略だと思います．

# (3) オプションの利用
ここまで，オプション設定を利用してきませんでした．そのため，
* エボリューショナリー エンジンは，90秒の実行時間が必要だが，紹介した最適化問題を全て解ける
* 紹介した最適化問題を全て解けるのは，エボリューショナリー エンジンのみ

ということが起きていました．これは，正確には基本(デフォルト)設定を採用した場合の話です．

ここからは，オプション設定を変更することで，以下の2つの課題を解決する方法を紹介します．
* GRG非線形 エンジンは，滑らかでも線形でもない問題を解けない
* エボリューショナリー エンジンは，実行時間が90秒必要

## オプションについて
まず，ソルバーのパラメーターを開きます．
最適化エンジンを選択するタブの右隣にある「オプション(P)」をクリックします．
ここから，設定を変更することができます．
![](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/308259/43cb267d-6159-40f6-bcdd-07647c377b68.png)

以下が，オプションで変更可能な設定の全てです．
![excel32.PNG](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/308259/d0008cda-2f44-1679-2a06-7481b9927e6b.png)
![excel33.PNG](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/308259/0058793f-6677-b07b-d52c-ab29abedf1e6.png)
![excel34.PNG](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/308259/1a6b2262-76b5-267e-01c4-00d0f20f9fca.png)

## GRG非線形 エンジンで，滑らかでも線形でもない問題を解く
GRG非線形のオプションである「マルチスタートを使用する」にチェックを入れてエンジンを動かすと，滑らかでも線形でもない問題を解くことができます．
この方法は，進化計算と同じく，乱数を利用する確率的な手法です．
母集団サイズを大きくする or ランダムシードを変えて複数回試行することで，より良い解を獲得できる可能性が高まります．

左の図が母集団サイズ100，右の図が母集団サイズ1000で問題（[Rastrigin function](https://en.wikipedia.org/wiki/Rastrigin_function)）を解いた結果です．
|![excel35.PNG](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/308259/24fedf20-ebff-20bd-fb28-cd9d6d9cc0de.png)|![excel36.PNG](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/308259/35a9530e-9110-03a3-f7e0-8eb1623eb135.png)|
|:---:|:---:|

マルチスタート無しの通常のGRG非線形 エンジンは，母集団サイズ1で入力値が初期解です．
母集団サイズ100の場合，通常のエンジンより良い結果ですが，変数値も目的関数値も最適解から離れています．
母集団サイズ1000の場合，変数値が最適解に近く，目的関数値0の解を獲得できています．
マルチスタートを使用した場合，通常のエンジンと比較して計算時間は100倍，1000倍に増えますが，それでも5秒もかからずに問題を解くことができます．
また，母集団サイズ1000の結果から，Excelの数値計算の精度では，これ以上に最適な解の獲得が難しいことが分かります．（目的関数値で差がつかないため）

**「マルチスタート」を使用すれば，GRG非線形 エンジンで滑らかでも線形でもない問題を解けます．**

## マルチスタート（戦略）とは
マルチスタート（戦略）は，文字通り，複数の解から最適化を開始する戦略です．
最適化の最初に，一様乱数や[LHS](https://qiita.com/tomtkg/items/759126576184750a9971)を利用して解を用意します．
進化計算やベイズ最適化は，マルチスタートを利用するのが基本です．
マルチスタートのGRG非線形は，「ランダムに解を1つ生成して，それを最適化する」操作を繰り返します．
すべての計算が独立しているため，他の解の影響を受けません．
一方で，進化計算やベイズ最適化は，解同士が相互に影響し合いながら最適化します．
交叉や環境選択，確率分布の生成に複数の解の情報を利用します．

以下の画像は，マルチスタートを説明するイメージ図です．（MATLABで作成）
|![](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/308259/17b52440-51e9-935f-5a14-399ea5593021.png)|![](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/308259/ed546974-f61a-101c-b290-f66b593d7fed.png)|![](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/308259/9e63208b-a16a-11a9-e8d1-49e972c86e92.png)|
|:---:|:---:|:---:|
|図1．1個の解|図2．ランダムな100個の解|図3ランダムな1000個の解|

図は，最適化開始時点の様子を表しています．
$n=2$次元のRastrigin functionの等高線図を共通してプロットしています．
横軸が$x_1$，縦軸が$x_2$，色がその地点の目的関数値$f(\boldsymbol{x})$の値を表します．
赤点が解を表します．
マルチスタート無しの通常のGRG非線形 エンジンは，図1のように与えられた1個の初期解から最適化を行い，近傍の局所最適解を獲得します．
一方で，マルチスタート有りの場合，図2, 3のように，ランダムに生成した複数の初期解が与えられます．
図2は，$(x_1,x_2)=(0,0)$付近に初期解が無いため，$(x_1,x_2)=(-1,0)$付近の局所最適解を最良解として獲得すると考えられます．
図3は，$(x_1,x_2)=(0,0)$付近に初期解があるため，最適解に極めて近い解を獲得できると考えられます．

## エボリューショナリー エンジンで，短時間で問題を解く
エボリューショナリー エンジンの実行時間が長い一番の原因は，「改善が見られない最大時間」が30(秒)と設定されていることにあります．
この値を3などに小さくすることで，劇的にエンジンの実行時間を短くできます．
他にも，「最大時間（秒）」や「最大実行可能解数」などの「解決の制限」をすることで，エンジンの実行時間を短くできます．

プログラムの実行時間を自由に設定できるのは，進化計算の魅力の1つです．
他の最適化手法では，実行時間が決まっており，実行時間が短すぎると何も出力できず，実行時間が長くても出来ることが無い場合があります．
また，短時間の実行で良い解を獲得できますが，長時間実行できても，より良い解の獲得が困難な場合があります．
進化計算は，短時間の実行でも良い解を獲得できますし，長時間実行できるなら，その分より良い解を獲得できる可能性が高い利点があります．

最適化問題を解くにあたって課題が見つかったときには，別の方法を試す前にオプションの変更を検討しましょう．
特に進化計算は，パラメータやオプションが多く，様々な最適化問題に対して通用する十分なポテンシャルを秘めています．

# (4) 制約条件の利用
Excelソルバーのダイアログにおいて，最も大きな面積を占めているのが制約条件に関する表示です．
逆説的に，最適化において最も大切なことは，制約条件の設定と考えることができます．

最初の進化計算実行時から，変数に上限値，下限値を設定する4つの制約条件を追加していました．
次は，この制約条件を変更してみます．

## 制約条件について
Excelで設定できる制約条件は，全部で6つです．
上が設定画面，下が設定例です．
![](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/308259/d2c50c02-800e-70b2-90e2-26366ff23427.png)
![](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/308259/fecb84b4-6754-8d38-b823-07789aee75bd.png)

上の3つ`<=`，`=`，`>=`は，左辺と右辺の大小関係を表す比較演算子です．
下の3つ`int`．`bin`，`dif`は，設定すると`=`を含む条件に自動で変更されます．
設定すると，`int`は整数，`bin`は0か1，`dif`は１以上の互いに異なる整数値しか取れなくなります．
`dif`は，セルを範囲指定して使います．
右上の設定例では，`$A$6:$A$9`セルには`1,2,3,4`や`3,1,4,2`が入ります．

## `=`の使用例
総和が100という制限のもと，3変数の総乗を最大化する問題は，以下になります．
　　　　$\text{Maximize  }\quad f(\boldsymbol{x}) = \prod_{i=1}^3 x_i = x_1 \times x_2 \times x_3$　　E2セル=A2\*B2\*C2
　　　　$\text{Subject to}\quad \sum_{i=1}^3 x_i = x_1 + x_2 + x_3 = 100$　　 D2セル=A2+B2+C2
![](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/308259/c6093ee5-54ab-be46-971b-80af5369a85b.png)
答えは，$x_1=x_2=x_3=100/3$です．

## `int`の使用例
[Rastrigin function](https://en.wikipedia.org/wiki/Rastrigin_function)の変数を整数に限定します．
　　　　$\text{Minimize  }\quad f(\boldsymbol{x}) = 10n + \sum_{i=1}^n [x_i^2 - 10\cos(2\pi x_i)]$!
　　　　$\text{Subject to}\quad-5.12\leq x_i \leq 5.12,\quad x_i\in\mathbb{Z},\quad$ for $i = 1,2,\dots,n.$
![](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/308259/d5ab712e-9c86-b696-3e29-912195e281a3.png)
実数値を離散化して整数にすると，最適化問題が解きやすくなる場合があります

## `bin`の使用例
0-1ナップサック問題で考えます．
重み$w_i$，価値$v_i$を持つ$n$個のアイテムを取捨選択する問題です．
重量の合計が$W$以下で，価値の合計が最大となる組み合わせを求めます．
　　　　$\text{Maximize  }\quad f(\boldsymbol{x}) = \sum_{i=1}^n v_ix_i$
　　　　$\text{Subject to}\quad \sum_{i=1}^n w_ix_i \leq W,\quad x_i\in\{0,1\},\quad$ for $i = 1,2,\dots,n.$
![](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/308259/e08f6d9f-f5f6-b569-218d-16f67116a5e2.png)
$W=67$の場合，アイテム1, 4, 8を選択する，重量の合計が67で価値の合計が1270の結果が得られました．

## `dif`の使用例
$1,2,\dots,n$までの数字をランダムに並べ替えることを考えます．
第$i$要素を$10^i$と掛け合わせた合計の最大化は，以下になります．
　　　　$\text{Maximize  }\quad f(\boldsymbol{x}) = \sum_{i=1}^n x_i\times10^i$
　　　　$\text{Subject to}\quad\boldsymbol{x}=\{x\in\mathbb{N^+}|x\leq n\},\quad n(\boldsymbol{x}) = n.$
![](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/308259/8c3609c3-1a06-0061-98fa-e0d6c01951bd.png)
[巡回セールスマン問題](https://ja.wikipedia.org/wiki/%E5%B7%A1%E5%9B%9E%E3%82%BB%E3%83%BC%E3%83%AB%E3%82%B9%E3%83%9E%E3%83%B3%E5%95%8F%E9%A1%8C)は，このような順列最適化の代表的な問題になります．

# 終わりに
Excelで進化計算を利用する方法，Excelのソルバー機能を紹介しました．
記事では，進化計算によって色々な最適化問題を解きました．
皆様も日常のちょっとした問題をExcelで解決してみてはいかかでしょうか．

※この記事の「問題を解ける」は，「最適解に近い近似解を獲得できる」ことを意味します．
※厳密解，精度保証付き近似解が要求されない環境を想定しています．