@yuichirominato 2018.06.19更新 73views

D-Waveで深層学習の基礎となるRBMのボルツマン学習を実行してみた

D-Wave QUBO RBM イジング ディープラーニング 深層学習 量子アニーリング

はじめに

普段量子コンピュータの勉強会をしているのですが、D-Waveを使用した機械学習に関して興味がある方が多いのと、質問が多いのでまず基本的な学習の過程のおさらいとD-Waveを活用してまずは簡単な例題を解いて見たいと思います。

参考資料や記事

これまで様々な問題を検証して記事にしてきました。基本的な用語などは過去の記事から参照してしまっているので、そちらもご覧ください。

イジングモデルや量子アニーリングについては過去の記事
「量子アニーリング、イジングモデルとフレームワーク」
http://blog.mdrft.com/post/6

一般的なイジングモデルのNP問題の解法はこちら
「NP問題のイジング」
http://blog.mdrft.com/post/42

また、Yoshua Bengio先生の論文にもある通り、ソフトウェアシミュレーションが活用できます。自社でもGPUやスパコンを活用したシミュレーションを活用していましたので、実機が活用できない環境でもシミュレーションでの実装をお勧めします。

RBMとは

RBMは制限ボルツマンマシン、制限付きボルツマンマシンなどと呼ばれ、ボルツマンマシンと呼ばれるタイプのNNに層を制限して2層の可視層と隠れ層に制限したものです。

参考:
Restricted Boltzmann machine
https://en.wikipedia.org/wiki/Restricted_Boltzmann_machine

D-WaveマシンとRBMとの関係とはじまり

モントリオール大学のこの論文が始まりです。ディープラーニングでの有名な先生がD-Waveでの深層学習の実装可能性として下記の論文を発表し、D-Waveマシンへの実装が検討されました。

On the Challenges of Physical Implementations of RBMs
Vincent Dumoulin and Ian J. Goodfellow and Aaron Courville and Yoshua Bengio
https://arxiv.org/pdf/1312.5258.pdf

この際にはD-Waveを模したソフトウェアシミュレータでNLLとよばれる学習の際に必要なRBMの勾配の推定にサンプリングが活用できるという方針は模索されましたが、実機は使われませんでした。

かなりいい論文で、ボルツマンマシンの基本が書かれています。この論文とGeoffrey Hinton先生の下記の論文はRBMを学ぶ上で大変役立ちました。

A Practical Guide to Training Restricted Boltzmann Machines
Geoffrey Hinton
https://www.cs.toronto.edu/~hinton/absps/guideTR.pdf

量子ボルツマンサンプリングの提唱

次に、実際に上記の問題をD-Wave実機で行なった結果として、

Application of Quantum Annealing to Training of Deep Neural Networks
Steven H. Adachi, Maxwell P. Henderson
https://arxiv.org/abs/1510.06356

があります。これは量子アニーラをボルツマンサンプリングマシンとして、上記のNLLの推定に役立てようという具体的な方針が示され、実機で学習が実際にされています。

引用:https://arxiv.org/pdf/1510.06356.pdf

D-Waveクラウドシステムへの実装

上記のボルツマンサンプリングを活用した機能は実はD-Waveマシンに搭載がされています。論文中に出て来る逆温度パラメータ$\beta$とともに、機械学習への応用ということで近年のD-Waveの主要機能としてD-Wave社は提供を開始し、活用を期待しているようです。

ボルツマンサンプリング、RBM、DBMを活用した強化学習

さらに、上記のサンプリング機能を使って、RBMを複数接続したDBMを学習し、簡単な強化学習を解くという論文も出ています。

Reinforcement Learning Using Quantum Boltzmann Machines
Daniel Crawford, Anna Levit, Navid Ghadermarzy, Jaspreet S. Oberoi, Pooya Ronagh
https://arxiv.org/abs/1612.05695

詳細なアルゴリズムも掲載されていますので、実装に関してはとても有用かと思います。

引用:https://arxiv.org/pdf/1612.05695.pdf

今回やること・学ぶこと

今回は上記2つの論文を確認しながら、自社で借りている実機のD-Waveを使いながら、RBMの基礎とサンプリング学習の実装をみていきたいと思います。また、D-Waveのマニュアルにも最近は正式にボルツマンマシンのドキュメントが搭載されていますので、それを参照しながら機械学習の例題を行なってみたいと思います。

Application of Quantum Annealing to Training of Deep Neural Networks
Steven H. Adachi, Maxwell P. Henderson
https://arxiv.org/abs/1510.06356

Reinforcement Learning Using Quantum Boltzmann Machines
Daniel Crawford, Anna Levit, Navid Ghadermarzy, Jaspreet S. Oberoi, Pooya Ronagh
https://arxiv.org/abs/1612.05695

使用するツールとマシンについて

ツール類は今回はpythonを使います。マシンはD-Waveを使用します。

引用元:https://www.dwavesys.com/resources/media-resources

引用元:https://www.dwavesys.com/resources/media-resources

詳細はカナダのホームページを参照してくださいませ。
https://www.dwavesys.com/

日本語サイトもあるようです。
http://dwavejapan.com/

モデルと結合について

実は最初に前提としてモデルを構成する結合数に関しての考察がされているので、それは遵守したいと思います。

On the Challenges of Physical Implementations of RBMs
Vincent Dumoulin and Ian J. Goodfellow and Aaron Courville and Yoshua Bengio
https://arxiv.org/pdf/1312.5258.pdf

この論文中で、D-Waveの結合数とRBMの結合に関して、RBMの結合を簡略化して学習ができないかどうかという試みがされました。実際にはD-Waveはキメラグラフと呼ばれる結合を持っており、1量子ビットからの接続数が少ないので、RBMの結合を簡略化してD-Waveに合わせて学習ができないかということが検討されました。

しかし、上記の論文での結論は結合数を減らすこと(モデルを簡略化すること)では学習がうまくいかないという結論でしたので、以降はRBMの結合数をきちんと守った形で結合を遵守するという方針があります。

下記にRBMとchimera restricted RBMの模式図を描いてみました。

8量子ビットある場合、上記の完全二部グラフは結合数が、N*Nありますが、下のキメラグラフは元々の結合が少ないので、二部グラフがつくりにくいという事情があります。

ただ、下のChimera Restricted RBMはうまくはいかないので、上記のRBMモデルを忠実に再現するということを、D-Waveのキメラグラフを使いながら実現していく必要があります。

早速仕組みを少しずつみていきたいと思います。

モデル概略

まず、モデルは可視層と呼ばれる層と隠れ層と呼ばれる2層からなるモデルです。結合に方向性がなく無向グラフです。

引用:https://arxiv.org/pdf/1510.06356.pdf

結合分布がギブス分布、ボルツマン分布に従うようです。

$p(v,h) = \frac{1}{Z}exp(-E(v,h))$

結合分布
確率論において、同時分布(どうじぶんぷ)または結合分布(けつごうぶんぷ, joint distribution)とは、確率変数が複数個ある場合に、複数の確率変数がとる値の組に対して、その発生の度合いを確率を用いて記述するもので、確率分布の一種である。
引用:https://ja.wikipedia.org/wiki/%E5%90%8C%E6%99%82%E5%88%86%E5%B8%83

ボルツマン分布
ボルツマン分布(ボルツマンぶんぷ、英語: Boltzmann distribution)は、一つのエネルギー準位にある粒子の数(占有数)の分布を与える理論式の一つである。ギブス分布とも呼ばれる。気体分子の速度の分布を与えるマクスウェル分布をより一般化したものに相当する。

量子統計力学においては、占有数の分布がフェルミ分布に従うフェルミ粒子と、ボース分布に従うボース粒子の二種類の粒子に大別できる。ボルツマン分布はこの二種類の粒子の違いが現れないような条件におけるフェルミ分布とボーズ分布の近似形(古典近似)である。ボルツマン分布に従う粒子は古典的粒子とも呼ばれる。
引用:https://ja.wikipedia.org/wiki/%E5%90%8C%E6%99%82%E5%88%86%E5%B8%83

ここで、この結合分布はエネルギー関数で規定され、可視層のノード数$n$と隠れ層のノード数$m$とで下記のように示されます。

$
E(v,h) = -\sum_{i=1}^nb_iv_i-\sum_{j=1}^mc_jh_j-\sum_{i=1}^n\sum_{j=1}^mW_{ij}v_ih_j\\
ここで、v_i,h_j\in\{0,1\}
$

また、正規化のための規格化定数、分配関数は下記の通りになります。

$
Z = \sum_{v_k}\sum_{h_l}exp\left(\sum_kb_kv_k + \sum_lc_lh_l + \sum_{k,l}W_{kl}v_kh_l \right)
$

また、完全二部グラフより、条件付き確率分布はそれぞれ、vとhについてシグモイド関数を用いて下記の通りになります。

$
P(h_j=1|v) = sigm(c_j+\sum_iW_{ij}v_i)\\
P(v_i=1|h) = sigm(b_i+\sum_jW_{ij}h_j)
$

学習について

次に上記の確率分布からなるNNのモデルの学習方法を確認したいと思います。複層のDBN(DBM)でも、RBMの形にして学習をします。これらの学習は$logP$を最大にするように、トレーニングデータとモデルの誤差計算をします。結合係数やバイアスの勾配計算は$logP$を用いて下記のようになります。

$\frac{\partial logP}{\partial W_{ij}} = < v_ih_j > _{data}- < v_ih_j > _{model}\\
\frac{\partial logP}{\partial b_i} = < v_i > _{data}- < v_i > _{model}\\ \frac{\partial logP}{\partial c_j} = < h_j > _{data}- < h_j > _{model}$

ここで、上記の結合係数の勾配計算で、$ _{model}$の期待値の計算はあまり効率的な計算がありません。

$ < v_ih_j > _{model} = \frac{1}{Z}\sum_{\{v_k\}}\sum_{\{h_l\}}v_ih_j exp\left(\sum_kb_kv_k + \sum_lc_lh_l+\sum_{kl}W_{kl}v_kh_l\right)$

実際には、この値を直接求めるのが大変なので、近似計算でCD法など、Gibbsサンプリングを応用した手法で、順番に隠れ層と可視層の計算を行なって値を取るという方法が行われます。

参考記事:https://qiita.com/t_Signull/items/f776aecb4909b7c5c116

なので、CD法などでは計算量を節約する意味で、かなり近似的な計算が行われ、これを精度改善しようとすると計算量と時間がかかります。

$\Delta W_{ij} = \epsilon[ < v_ih_j > _{data}- < v_ih_j > _{recon}]$

上記は結合係数の更新ですが、バイアスの更新でも同様です。

量子ボルツマンサンプリングを活用したパラメータ更新

ここで、このCD法の代わりに量子アニーリングマシンを使って元の勾配計算どおりしまおうというのが実機を使ったボルツマンサンプリング学習の基本です。

引用:https://arxiv.org/pdf/1510.06356.pdf

ボルツマンサンプリングを使った際の結合係数やバイアスの更新式は、下記の通りです。$\alpha$はモーメンタムで$\epsilon$は学習率です。

$W_{ij}^{(t+1)} = \alpha W_{ij}^{(t)} + \epsilon[ < v_i h_j > _{data} – < v_ih_j > _{model}]\\
b_{i}^{(t+1)} = \alpha b_i^{(t)} + \epsilon[ < v_i > _{data} – < v_i > _{model}]\\
c_{j}^{(t+1)} = \alpha c_j^{(t)} + \epsilon[ < h_j > _{data} – < h_j > _{model}]$

RBMを学習した後は古典計算機でバックプロパゲーションで仕上げの学習をするようです。

量子アニーラを活用したサンプリング手法

D-Waveマシンは量子アニーリングという理論をベースに作られており、基本的には最適化問題を解くために最小値問題を求めるような設計になっています。しかし実際には、外部環境の影響やその他の理由で、最適解に落ち着かない場合が多いです。そのため、この局所解へ落ちる性質を利用して、分布を求めるサンプリングマシンとして活用するというのが、サンプリング学習の基本です。

この励起状態のばらつきがボルツマン分布として仮定できると、下記式に近似できこれと初期に出てきた分布の式を対応させることで、サンプリング手法を更新式に導入できそうです。

$P(x) = \frac{1}{Z}exp\left(-\beta_{eff}H_f(x)\right)$

$H_f$は最終的に求めるハミルトニアン(エネルギーコスト関数)で、$\beta_{eff}$はサンプリングの分布を調整する変数です。

これを使用することで、一番計算量のかかる部分を下記のように近似します。

$\overline{v_ih_j} = \frac{1}{N}\sum_{n=1}^N v_i^{(n)}h_j^{(n)}$

これを$ < v_ih_j > _{model}$に適用します。他の可視層や隠れ層のモデル期待値も同様です。

D-Waveマシンを想定した結合方法

D-Waveはキメラグラフと呼ばれる結合方式を採用しているので、RBMの完全二部グラフを再現するには少し工夫が必要です。

キメラグラフ参考資料:https://qiita.com/YuichiroMinato/items/e6952fec1a9965156873

具体的には、量子ビット同士を繋げてコピー操作ができるので、それを活用してRBMを実現します。繋ぎ方は下記のようです。実際には完全二部グラフはとてもわかりやすい、縦横の交差になります。この交差部分に結合係数の値を設定していきます。

上の図では左側はハードウェアを想定した結合、実際にハードは量子ビットの形状は細長くなっています。右側は論理的に書いた図です。

QUBOの作成と、逆温度パラメータ

D-Waveで特には、QUBOmatrixというものを作ります。QUBOmatrixは問題をバイナリ値でかんがえた時に作る行列で、これを少し変形してD-Waveにかけます。

引用:https://arxiv.org/pdf/1510.06356.pdf

結果として、対角項にバイアスbとcがはいり、非対角項に結合係数wが入ります。バイアスは局所磁場とか磁場項とよび、結合係数は相互作用項と呼んだりします。

QUBOからイジングへ

次に、逆温度のパラメータはここで見ると、行列の外にだせるので全体のスケールを調整するパラメータになってます。QUBOからイジングに直す時にはまだかんがえなくていいので、イジングmatrixに直します。

あとはスケーリングパラメータを決めてD-Waveでサンプリングを実行

あとは$\beta_{eff}$を決めて実行です。どうやら論文によると、$\beta_{eff}$のあたいは実行するサンプルによって依存するようです。MNISTの画像認識の学習では、画像サイズが大きくなるにつれて、$\beta_{eff}$は小さくするのが適切のようでした。

また、D-Waveの量子ビット数は多くないので、完全二部グラフを作って学習をさせる際には、MNISTをダウンサンプリングして実装しています。

MNISTの28*28pixelの画像を6*6まで下げています。

引用:https://arxiv.org/pdf/1510.06356.pdf

実行結果などの学習具合は詳しくは論文を参考にしてもらって、だいたい下記のように収束が良いようです。

引用:https://arxiv.org/pdf/1510.06356.pdf

続いて実機を使ってこれを実装する

D-Waveには前処理を変更して、逆温度パラメータを設定する項目があるので、これを使用してみます。

簡単な例題で実行する

簡単な例題で上記のRBMのサンプリングに関して学んでみたいと思います。
まずはD-Waveの例題でも載っている、簡単な問題で。

$Q = \left[
\begin{array}{rr}
-1 & 2 \\
& -1 \\
\end{array}
\right]$

このQUBOmatrixの場合、エネルギーコスト関数は、

$E(x) = -x_1-x_2+2x_1x_2$

となります。$x$の場合の数でエネルギーを求めて見ると、

解析的にかんがえて見ると、ここではすべての場合の数が出せているので、まず規格化定数$Z$を求めます。

$Z = \sum exp(E(x)) = exp(0) + exp(1) + exp(1) + exp(0) = 1+2.718+2.718+1 = 7.44$

そして、

$P = \frac{1}{Z}exp(E)$

より、確率は0.13と0.37がでます。

これを手作業で出すのは大変なので、サンプリングでかんがえてみます。

D-Waveでサンプリングする

まずこれをQUBOからイジングに直します。

$
Q = \left[
\begin{array}{rr}
0 & 0.5 \\
& 0 \\
\end{array}
\right]
$

こうなりましたので、これをいれます。
係数のJijに-0.5いれて、1000回のサンプリングでかけると、、、

結構衝撃的です。。。大体の確率でサンプリングできてる。。。まじですか、、、

解を確認したところ大丈夫そうです。


速度も速い。

ということで、サンプリングを行なった結果、モデルの分布が取れています。先ほどは解析的にときましたが、サンプリングから逆算してモデルの期待値が取れますので、高速に誤差計算できそうです。サンプリングした結果は古典計算機で値の更新をします。

大成功でした。

あわせて読みたい

SERVICE

汎用量子コンピュータ向けアプリケーション開発SDK

詳しく見る Githubで入手する(無料)

汎用量子コンピュータ向け高速シミュレータ

詳しく見る

量子コンピュータ向けクラウドサービス(準備中)

詳しく見る

イジングマシン向けアプリケーション開発SDK

詳しく見る Githubで入手する(無料)

COMMUNITY

量子コンピュータのことを学びたい、仕事としたいなどの情報交換の場を設け、コミュニティの形成を進めています。オフラインの勉強会と、オンラインのチャットコミュニティの2種類あります。オフラインのConnpassは1400名、オンラインのSlackは880名を超える参加があります。どちらも無料ですのでお気軽にご参加ください。

CONNPASS SLACK

CONTACT

弊社あての連絡はinfo@mdrft.comより用件を明記の上、メールで連絡を頂けますと幸いです。

ブログトップに戻る

Recent Posts