pirika logo

ホームページ Pirikaで化学 ブログ 業務リスト お問い合わせ
Pirikaで化学トップ 情報化学+教育 HSP 化学全般
情報化学+教育トップ 情報化学 MAGICIAN MOOC プログラミング
MAGICIANトップ MAGICIAN-Jr. MAGICIAN-講義 過去の資料

MAGICIANとは、材料ゲノム(Materials Genome)、材料情報学(Materials Informatics)、情報化学(Chemo-Informatics)とネットワーク(Networks)を結びつけて(Associate)いかれる人材です。

MAGICIANトップ > MAGICIAN講義資料 > 第7回 遺伝的アルゴリズムを用いた最適解探索

2018.11.8

MAGICIAN(MAterials Genome/Informatics and Chemo-Informatics Associate Network)養成講座

Materials Genomeという言い方をした時、Genome:遺伝子はどのように扱うのでしょうか?

医薬品でゲノム解析と言った場合には、ある遺伝子特有の病気に対して、オーダーメード医薬品を作るようなイメージで使われます。

マテリアルでゲノム解析と使う場合、例えば触媒を考える時に、触媒の性能を司る“遺伝子”を特定して、性能を向上させるためには“どう進化”させたら良いのかを考えます。

例えば、馬からサラブレッドを作る、美味しいお米を作る、病害虫に強い作物を作るのと同じようなものです。

サラブレッドを育てる、適者生存の法則を使います。

進化した生物は雌雄を持つものが多いです。
これは親の遺伝子を組み合わせて、性質の良いとこ取りするのに適しています。

触媒の良い遺伝子と言った場合に、例えば、転化率が高い、選択率が高い、寿命が長い、コストが安いなどを目指して、処方、プロセスを最適化していきます。

このように、目的変数が複数ある場合の最適化を考える時、昔から使われているのが、遺伝的アルゴリズム法というやり方です。
これは、簡単にいえば“適者生存の法則”をコンピュータ上に載せたものです。
レースに勝つ馬の遺伝子がどうなっているかは置いておいて、レースに勝つ馬同士を掛け合わせて子供を産ませればレースに強いサラブレッドができます。

その際には、レースの定義が重要になります。
触媒の転化率がいくら高くても、選択率が低ければ触媒としてのレースの適者ではありません。
そこで、転化率や選択率を予測する式をニューラルネットワーク法で作成し、遺伝的アルゴリズム法で最適化するという方法が成立します。

ニューロンや遺伝などという言葉が入っているので、まるでバイオの世界(実際、昔その手の書籍は生物のコーナーに置いてあった)ですが、医薬品の設計から、マテリアルの設計まで広く使われています。

今回はこの遺伝的アルゴリズム法を学びましょう。

遺伝的アルゴリズム法

遺伝的アルゴリズム法(GA: Genetic Algorithm)は1975年に提案された近似解を探索するアルゴリズムです。

アルゴリズムというのは、“考え方”のことです。
そこで、汎用のニューラルネットワークのシミュレータは存在しますが、汎用のGAのソフトは、私は知りません。

“考え方”に従って、自分でプログラムを実装していかなければならないので、広く一般的には使われていないのでしょう。

GAの考え方の例題として、一つの山の頂点を探すことを考えて見ましょう。

自分のいる地点の傾斜を見て、登りの方向に常に進めば、何時かは頂点(どの方向へ進んでも下りになる)にたどり着きます。

しかし、そうした微分係数を用いた探索方法では、多峰性の山登りには使えません。
初期値をどこに取ったかに依存してしまいます。

また、評価関数が微分できないと傾きがわからないので進む方向がわかりません。

t

t

GAでは、多次元に散りばめられた多くの評価点から傾きを使わずに探索を行います。

完全にランダムな方向に動かすランダム探索では、特に次元数が増えてきたときに近似解の探索が非常に困難になります。

例えば、宅配便のトラックが50個の荷物を配達する事を考えてみましょう。
各家を回るルートの数は50!(50*49* *2*1)存在します。

そのルートを最短に(ガソリン・コストを最小に)するために、50!の距離を全部計算して最短距離を探すのはドライバーには不可能ですし、スパコン使っても難しいのです。

ましてや、その1軒が再配達となったら、その1軒に戻るルートも含めて最短距離を再計算しなくてはならなくなります。

囲碁や将棋の次の数手先を読むなども同じような問題を抱えています。
そこでランダムに動かすのでは無く、”適者生存の法則”で動かすのがGAです。

例えば、触媒の組成を遺伝子に見立てます。(本来はスキーマ理論というのがあって遺伝子は0か1のビットで表します。しかし、化学の現象はビットでは表しにくいので、私は%の実数で計算を行っています。)

そうした遺伝子を多数用意し遺伝子のプールに入れます。

t

そして、各遺伝子の評価値を計算します。

例えば、次のような評価式を自分のニーズに合わせて作ります。 E= 30*転化率 + 50*選択率 + 100*寿命 これらの各係数は自分の考えた“適者”によって様々な値になり得ます。
これがレースを定義するという意味で、GAの設計の一番難しいところで、一番楽しいところです。

転化率や選択率は少し余裕があるので多少落ちても、寿命が長くなってほしいなどによって評価式を作ります。
そして各遺伝子の評価値を計算したら、ルーレットを作成します。

このルーレットのポケットのサイズは、評価値が高いものほど広くとります。
つまり適者ほど親になる権利が高くなります。

そして、ルーレットを2回まわして、遺伝子を二つ選択します。

そして、ある確率で突然変異と交叉を行って、次世代の遺伝子プールにためます。
次世代の遺伝子プールに初期と同じ数の遺伝子が入ったら、世代交代をして、同じ事を繰り返します。

突然変異とは、ある元素の組成が変わることに相当します。
メキシカン・ハットの山登りの例では、Yの値は良いが、Xの値が悪いので評価が低いが、Xの値が突然変異すると高い評価になったりします。

t

交叉は、遺伝子のある部分を切ってもう一つの遺伝子とつなぎ変えてあげます。
X◯Y✖️とX✖️Y◯を繋ぎかえるとX◯Y◯とX✖️ Y✖️の2つの遺伝子になります。

どのくらいの頻度で突然変異と交叉を行うか?は難しい問題です。

基本的には突然変異と交叉の操作を行うと遺伝子的には悪い遺伝子になります。
生物の進化が何百万年もかかるゆっくりとした進化なことからも明らかでしょう。

詳しい事は自分もわかりませんが、癌細胞のような悪い遺伝子を持った細胞が優性遺伝子なら人類はもっと早く滅亡していたでしょう。

そこで、一番性能の高い遺伝子は、突然変異と交叉を行わず必ず次世代に残すという操作を取り入れることもあります。
これをエリート保存といいますが、必ずしも良い結果になるという保障はありません。
そのエリートに引きずられてローカルミニマムから抜け出せない事もあります。

本来、ビットで表せば
突然変異  11100101 → 11110101
交叉    10100101 → 11110101
       11110100 → 10100100

となり余計な操作は必要ないのですが、組成を実数の%で表しているので、突然変異や交叉の後、%に直します。

その際には他の元素の組成も変わってしまうので、収束は遅くなります。世代を繰り返すごとに動かす絶対量を小さくするなどの工夫が必要ですが、余り気にしなくても、計算負荷は大きくないので、取り敢えずこれで良しとします。

この部分にこだわるより、この元素は高価なので、使用量は30%まで(30%を超えたら評価値を下げる)とか、レアーアースは使わないでとか、化学者ならではの評価関数の作成に労力を集中した方が良いと思います。

言葉で表すと以上のように難しいので、実例を示しましょう。
本来“遺伝的アルゴリズム”と言ったときに、アルゴリズムとは考え方を示すものなので、汎用の遺伝的アルゴリズム(GA)・ソフトウエアーはありません。
現象ごとにGAソフトを自分で組み立てることになります。

これがGAを勉強する上でネックになっているようです。

GA検討の実例

今回は、まず、酸化銅のナノ粒子の合成条件の最適化を、GAを用いて行ってみましょう。これは、次の反応で酸化銅のナノ粒子を作成します。

t

日本金属学会誌, 62, no.6, 1998, 519-526

様々な実験条件の時に、ナノ粒子の収量と粒子径が記載されています。
粒子径は範囲で与えられているので、小さい方をL、大きい方をHとして、沈殿量、粒径L、粒径Hの3つを予測する式を作成します。

式の作成自体は授業でも取り上げ、じっくりやったのでここでは省略します。
そのような式ができると、ある実験条件の時に目的変数である沈殿量、粒径L、粒径Hがどうなるかを予測できるようになります。

t

授業の際には上記のようなシミュレータを構築しました。

実験条件をエクセルのスライドバーで動かし、“沈殿量が80以上で一番小さなH/Lを求めなさい”が課題だした。

H/Lが一番小さいというのは、粒径が揃っている事の指標です。
こうした酸化銅のナノ粒子は、インクジェットプリンターなどでナノ回路を印刷し水素で還元して金属銅にして使われます。

粒径が揃っているというのは、ノズルが詰まらないために重要な指標となります。

粒径が揃うと沈殿量が少なくなったりしますので、両方の条件を満足する条件を手作業でやるのは大変だろうというのを実感してもらうために学生に競ってもらったのですがが、自分が事前に見つけておいた実験条件と結果はあっさりクリアーされてしまいました。

そうすると、GAをわざわざ組む必要も無いので、今回はさらに、粒径Lが2.5μm付近という条件をつけます。

実験条件5つを動かし、粒径Lが2.5μmで、沈殿量80以上、L/Hがなるべく1に近い実験条件を求めなさいが課題です。
(実際の製造では、製造コストが最小にとか、製造サイクルとか、さらに色々な制約条件がかかってきます。) 目的変数が3つ以上になると手作業で行うのは段々と難しくなります。
(それにしても、真面目で優秀な横浜国大の学生はできてしまうのが実に問題です。GAなどは、無精者の科学で、めんどくさいという感覚を持たない真面目な学生(やAI)には縁遠いのかもしれません。)

それでは、実際にGAソフトをエクセル上に組んでみましょう。

まず、最初の遺伝子のプールを作成します。
遺伝子を何個用意したらいいでしょうか? 

さしあたって、論文に記載されている実験条件48種類を初期遺伝子にします。
実験条件範囲内で、乱数で初期遺伝子を作成しても良いです。

探索する実験条件の上限値、下限値は定めておいたほうが良いでしょう。
それは、推算式の適用範囲の問題があるからです。

基本的には推算式を構築した際の実験条件を上下20%以上逸脱する場合には余り高い推算精度は期待しないほうがいいです。

実際には個々の説明変数の上限、下限だけでなく、組み合わせの上限、下限の問題もあるので複雑ではありますが、ここでは論文にある実験条件の範囲内で動かすことにします。

t

黄色でマークした部分が実験条件になります。
その時の粒径L、沈殿量、粒径H(緑でマーク)は予測値になります。
差し当たってこの予測値は十分な精度を持つとします。

すると、今回評価したい物性値(オレンジ色でマーク)L/H(粒径L/粒径H), Pre/80(沈殿量/80), (DL-2.5)2 (粒径Lから2.5を引いた値の2乗)が計算できます。

L/Hは1に近いほど優秀です。
Pre/80は1以上であるほど良いです。
(DL-2.5)2 は粒径Lが2.5μmから外れるほど急に大きな値になるので、なるべく小さな値になれば良いです。
この3つの評価値を組み合わせてEvaluation(青でマーク)を定義します。

このEvaluationの定義が一番難しく、一番面白いところです。
3つの評価値をどう捉えるかは実験を行う者の主観に基づきます。
Pre/80は1以上であればそれ以上大きな値を探索する必要がなく、L/Hを重視したいとか、粒径Lが2.5μmに近いものを優先したいとか主観に基づいて評価関数を定義します。

ここでは、 3*L/H+1.15*Pre/80-0.4*(DL-2.5)2 で評価します。
自分で高く評価したい評価値の係数を調整する。(DL-2.5)2 は小さいほど良いので係数はマイナスになります。

ここまで出来たら、Evaluationの列をキーにソート(並び替え)をかけます。
大きい順にソートをかけると、自分の目的にあった実験操作条件が上位にくることになります。

本来のGAでは、このEvaluationの値をルーレットのポケットの面積に変換して、ルーレットを回して次世代の遺伝子のプールに入れていく遺伝子を決めます。
エクセルでそれを行うのは大変なので、そこは簡略化します。

t

Evaluationのカラム(N列)の後ろには空のカラム(O列)があるのでソートをかけても、P列以降は変更されません。

ソートがすむと、二行目には評価関数が一番大きいものが入ってきます。
その値は2.796であることがわかります。
L/Hは0.5, Pre/80は1.128, (DL-2.5)2 は0.005です。

この実験条件はエリート(一番優れている)なので、これは保存することにします。

具体的には、Q2からU2には、C2からG2の値を入れます。その際には、=$C$2のようにセルの絶対位置を指定します。
3行目以降は、交叉(Cross Over), 突然変異(Mutation)の操作を加えます。

本来はルーレットを回して、交叉、突然変異をおこなう行を決めるのですが、ここでは固定行を使います。
Q3とQ4のセルを見てみましょう。

Q3=$C$2, Q4=$C$29となっている。

つまり、三行目には一番評価関数の高い二行目が代入され、四行目には29行目が代入されています。
29というのは適当に選んだ数字です。

上位にあるもの同士を交叉すると、多様性が確保できなくなり、ローカルミニマムに収束してしまう事が多いので、適当に悪い遺伝子(行の下の方)と交叉させます。

ここで、R列までは2行目、29行目を絶対位置指定で持ってきますが、S-U列では2行目、29行目を入れ替えます。
入れ替える位置によってオレンジ色でマークしてあります。
5列ある所のどこで入れかえるかは、本来は乱数で決めるのですが、固定位置にしてあります。

Q-U列では2行目、29行目の遺伝子が入り混じった遺伝子が作成されていることになります。
上の画面には表示されていませんが、27行目以降には突然変異の操作が示されています。

Q27=$C$2+($C$53-$C$2)*RAND()*0.3
Q28=$C$2+($C$52-$C$2)*RAND()*0.3

$C$2はエリートの値になっています。$C$53はC列の最大値、$C$52は最小値が入っています。

($C$53-$C$2)*RAND()*0.3はC列の最大値からエリート値を引いたものに0.000-1.000の乱数をかけ、さらに0.3を掛けた値分だけエリート値を大きい方へ動かすことに相当します。

逆にQ28では、小さい方に動かすことに相当します。

どの列を突然変異させるかは、本来は乱数で決めるが、交叉と同様、固定位置で突然変異を行います。

ソートをかけた段階で、最初のプールの値(C-G列)は確定しますので、次世代プールの値も確定されます。

そこで、次世代プール(Q-U列)にある48組の遺伝子をコピーし、C-G列に“値のみ”のペーストを行います。
すると、H―Jの物性予測値はすぐに計算されます。
そこで、再び、Evaluationの値でソートをかけます。

評価関数の値は、2.796から2.867に改善されます。
(乱数を用いているので、ここから先は計算毎に異なるので同じ値にならなくても気にすることはありません。)

L/Hは0.5から0.517, Pre/80は1.128から1.145, (DL-2.5)2 は0.005から0.002になります。実験条件的には、アンモニアの量と1/ksの値が変わっています。
これが第1世代になります。
そしてその時の「Q-U列をコピーしてC-G列に値のみのペーストをし、ソート」を繰り返します。

t

4世代ぐらい繰り返すと大体収束してしまいます。

10世代目には、L/Hは0.5から0.571, Pre/80は1.128から1.035, (DL-2.5)2 は0.005から0.037になります。

世代交代を繰り返すと、上位の実験条件は非常に似たものになってきます。
(10世代目ではトップの5つが同じになる。)

そこで、解の多様性を確保するためには、交叉をするペアの一つは行番号の大きなものを選ぶ方が良い結果を与えます。

評価関数の係数を変化させたり(この例では、L/Hが大きいものが優先されています)、交叉、突然変異の配分などを変えて、GAをどう設計すれば自分の欲しい答えが見つかるか、試してみましょう。 粒径の評価値を(DL-2.5)2 からABS((DL-2.5)に変えて評価してみましょう。

エクセルのマクロを設定できるなら、ボタンを作って、ボタンが押されたら、コピー、値をペースト、ソートまでの自動化をしてみましょう。ボタンが押されるたびに世代を交代できるようになるので、計算は非常に楽になります。

GAの基本形はこのエクセルでやったことが全てなので、目的変数、評価関数を書き換えるだけで他の事象に応用できるので、色々試してみましょう。

言葉で説明すると、とても複雑で理解できるとは思えません。
最初は実際にエクセルで走らせて、何をやっているかをじっくり理解するのが大事です。

GAの例その2

次はゼオライトによるナフサのリフォーミングを検討してみましょう。
実験条件(Temp., LHSV, Feed Rate, P, Time)を変えた時の転嫁率(Conversion), オクタン価が論文に記載されています。
Pak. J. Chem. 4(4): 167-184, 2014

ナフサはオクタン価が低いので高機能ガソリンを作る上で重要なプロセスです。 授業では、転嫁率とオクタン価の予測式を作成しました。

t

目的変数(緑でマーク)は2つあるので入れ替えながら、重回帰計算を行い、計算値(青でマーク)を得ました。

t

t

グラフ表示してみると、実験条件と実験結果は重回帰式で良好に記述できている事がわかります。

t

転嫁率とオクタン価には上図のような関係があります。

そこで、GAで求めたい答えは、重回帰式は正しいとして、転嫁率とオクタン価の両方が高くなるような実験条件を求めることになります。

説明変数は5つなので、例題1で作ったGA計算のシートがそのまま流用できます。
実験件数は14件なので、実験条件の範囲内で48件まで水増ししてあげます。
目的変数の数は2つなので、1列削除します。

さて、問題は評価関数です。
” 転嫁率とオクタン価の両方が高くなる”を関数で表すにはどうしたら良いでしょうか?

自分で考えてみましょう。

GAの例その3

この講座のNo.5 データのクレンジングでLiイオンの拡散係数を取り扱いました。

そこではニューラルネットワーク法で予測した拡散係数をGAでさらに拡散係数の高いものを探索しました。
拡散係数自体は組成に対して非線形であるので、ニューラルネットワーク法を用いないと正しく記述できません。

ニューラルネットワーク法を用いた推算式は、GAで用いようとすると、プログラムを書けるか、マクロを書けるかしないと扱えないので例題としては使いにくいです。

そこで、原子種を絞って(O, Li, Al, Si, Nb)YSBを用いてクロスタームを導入した予測式を作成します。

t

t

すると、上図のように拡散係数が高い領域では、そこそこ精度の高い相関式が得られます。

説明変数は5つであるので、前のGAのエクセルがそのまま使えます。
評価関数は目的変数がLog DLiの1つで、大きいほど良いので設定は簡単です。

この場合、問題になるのは、5つの原子種の%を全部たすと100%になるように調整しなければならない事です。

交叉や突然変異させた後にトータルが100になるように変換します。
(この場合は、酸素欠陥の問題は端折ってあります)

具体的には、交叉や突然変異を行った結果の入っている列のさらに後ろの列に、100%で規格化した列を作り、そちらをコピーして、C-G列に値のみのペーストを行います。

もう1つの問題は、5種類全部を使ってはいけない事です。
さらに、Nbを使った場合には、Al, Siは必ず0%にしなければいけません。
Li, O以外に1つ以上は選択しなければいけません。

この条件を満たすように改造して、シートを完成させましょう。

雑記

遺伝的アルゴリズム(GA)の汎用シミュレータはありません。
  それは評価関数や制限事項がケース・バイ・ケースになってしまうからです。
しかし、それを自分でなんとかできるのであれば、GAの動作原理を簡単に学ぶ事ができます。

実際にはプログラムを書ければさらに高度な事が可能です。
今回の例のように一つシートを作れば、使い回しができるのと同じで、プログラムに関しても一つ雛形を作れば後の改造は容易です。ぜひチャレンジしてください。

自分がGAのプログラムを初めて書いたのは2000年頃です。化学用のプログラムではありません。(出典が何だったか、色々探しているが見つかりません。きちんと引用文献をつけたいのですが。。。 誰かご存知の方がおられたらお知らせください。)

GAで解いた最初の問題は次の画像認識の問題だ。320*240の領域に色々な図形が書かれています。コンピュータの認識できるのは、各ドットが0(白)か1(黒)かだけです。その中に85*73の探索図形があるかどうか、その位置はどこかを回答させる問題です。人間であれば簡単にどこにその図形があるか認識できるでしょう。

t

最近のGoogleのディープ・ラーニングによる画像認識の走りのような技術開発でした。

簡単のため、図形の拡大、縮小は無視できるとする。探索図形の点をθ回転し,(Xm、Ym)に平行移動する計算は、

t

で表す事ができます。
探索空間は、320(横)240(縦)360(回転)で27,648,000になります。

これを全て計算して一致度が一番高いものを答えとするのは大変な事なのはわかるでしょう。

遺伝子型を[θ, Xm, Ym]として、乱数で作った遺伝子を1000個ばらまいて、一致度が高いものほどルーレットのポケットを大きくして遺伝的アルゴリズムで探索させると直ぐに答えは見つかりました。

図として保存された構造式やグラフをコンピュータに認識させる事ができるので必死になってGAを勉強しました。

勉強が進むにつれて、画像認識よりは利用価値の高い化学の分野が多いことに気がついて、画像認識自体は余りやらなくなりましたが、先に述べたようにプログラム自体は使い回しが効くので、未だに変数名などはこのころの物を使い回しています。

遺伝的アルゴリズムを使うか、ディープ・ラーニングを使うかは別として、「コンピュータは画像を0,1でしか認識できない」というのは非常な驚きと言えます。

子供にも簡単にできるような図形認識をさせるのに、膨大な計算をするか、膨大なタグ付けされた学習データを用意して学習させるしかありません。しかも、ちょっとでも変化させると(デフォルメさせると)人間と違って応用は全く効きません。

拡大、縮小などは行列計算に簡単に組み込めます。結局、コンピュータのできる認識は、0,1と行列計算のような論理計算だけということになります。

最近、人工知能(AI)が非常に流行りだが、脳科学者や数学者は、機械には知能は無いと言い切っています。
あるのは、統計と確率、論理計算と検索だけで、これは知能では無いと。

このAIブームには限界が見えているし、直ぐに廃れるでしょう。
でも廃れるというのも間違いで、熟成期に入ると考えた方が良いかもしれません。

自分がGAをゴリゴリやっていた2000年頃は、第2次AIブームの終わりの頃で、そのあとは下火になりました。

その頃から画像認識やディープ・ラーニングの基礎が熟成し始め、昨今のインターネット上のビッグデータが利用できるようになって一気に花開いたのが、現在の第3次AIブームでしょう。

今ここで撒かれたタネが熟成期にどう育って行くのか楽しみなところです。
少なくとも人間の子供程度の認識力は持って欲しいものです。

人間の脳が図形を見たときには、図形の認識をする脳部位、記憶を司る脳部位、推論をする脳部位など様々なNeural Networkが同時に活性化される。Activate Networksがキーになると自分は思っています。

では、その間、化学系の研究者は何もしないでただ待っているのでしょうか?

今回の例のように、論理的な評価関数を主観的に作る。化学者にしか理解できない制限事項の論理式を作る。それができれば、MAGICIAN、手品師のように問題解決できるようになるでしょう。

MAGICIANトップ > MAGICIAN講義資料 > 第7回 遺伝的アルゴリズムを用いた最適解探索(2018.11.8)


Copyright pirika.com since 1999-
Mail: yamahiroXpirika.com (Xを@に置き換えてください) メールの件名は[pirika]で始めてください。