ホップフィールドネットワークの想起過程を再現してみた

はじめに

甘利先生らは、ホップフィールドネットワーク（以降 Hopfieldモデル）について、想起過程のダイナミクスをかつて研究されていた [1]。甘利先生の著書 [2] にもまた、想起過程の計算機シミュレーションの図が示されている。

先日、HopfieldモデルをC言語で実装した（Link）ので、この想起過程の図も再現できるだろうと考えた。実際にやってみたらうまくいったので、その結果を報告するのが本記事の主旨である。

結果だけを知りたいひとは目次から【実験結果】へ
プログラムを動かしてみたいひとは Google Colab の Notebook へ → Link

はじめに
自己相関型の連想記憶モデル
実装
実験
- 実験結果その1：ニューロンの自己結合を許容しない場合
- 実験結果その2：ニューロンの自己結合を許容する場合
おわりに
参考文献

自己相関型の連想記憶モデル

モデルの学習

まずHopfieldモデル（自己相関型の連想記憶モデル）にかかる記号の準備を兼ねて、学習について簡単に説明しておく。本記事を書くにあたり、文献 [3] と [4] を参考にした。

いまHopfieldモデルが $n$ 個のニューロン（素子）からなるとする。各ニューロンは $\pm 1$ の値をランダムにとる。それらをまとめてベクトルで表示すれば、 $\mathbf{x} = [ x_1, x_2, \ldots, x_n ]^\top \in \{-1, 1 \}^{n}$ となり、これがモデルの特定の状態を1つ指定する。

モデルに記憶させるパタンを $\{ \mathbf{s}^{(1)}, \mathbf{s}^{(2)}, \ldots, \mathbf{s}^{(m)} \}$ により与える。 $\mathbf{s}^{(p)} = [ s^{(p)}_1, s^{(p)}_2, \ldots, s^{(p)}_n ]^\top$ であり、各要素は $\pm 1$ の値をとる。記憶対象パタンをモデルが取りうる特定の状態と見立てているわけである。

記憶行列 $\mathbf{W} = [ w_{ij}]$ を次式で定義する。

$\begin{align*} w_{ij} := \frac{1}{n} \sum_{p=1}^{m}s^{(p)}_{i}s^{(p)}_{j} \end{align*}$

$\mathbf{W}$ は $n \times n$ 行列である。 $w_{ij}$ は第 $i$ ニューロンと第 $j$ ニューロンの結合の強さを表す。ネットワークの「学習」はこの $\mathbf{W}$ の計算で尽きている（Hebb則による学習）。

上式をベクトルと行列で表示してみると、

$\begin{align*} \mathbf{W} := \frac{1}{n} \sum_{p=1}^{m}\mathbf{s}^{(p)}{\mathbf{s}^{(p)}}^{\top} = \frac{1}{n} \mathbf{S} \mathbf{S}^\top \end{align*}$

となる。ここで $\mathbf{S} = [ \mathbf{s}^{(1)}, \mathbf{s}^{(2)}, \ldots, \mathbf{s}^{(m)} ]$ と定義した。 $\mathbf{S}$ は $n \times m$ 行列である。

ニューロンの自己結合はないと仮定するならば、その対角成分は0となる。すなわち $i = j$ のとき、 $w_{ii} = 0$ である。この仮定のもとで改めて $\mathbf{W}$ を定義すれば、

$\begin{align*} \mathbf{W} := \frac{1}{n} \mathbf{S} \mathbf{S}^\top - r \mathbf{I}_n \end{align*}$

となる。ただし $\mathbf{I}_n$ は $n$ 次の正方行列である。ニューロンの自己結合を許容しない制約は、国内外の多くの文献で見られるものであり、記憶行列も通常はこちらの定義が採用される。 $r$ は記憶率と呼ばれ、次式で定義される。

$\begin{align*} r := \frac{m}{n} \end{align*}$

記憶率とはすなわちニューロン数に対する記憶パタン数の割合である。よく知られた結果により、Hopfieldモデルが記憶可能なパタン数は約 $0.15 \:n$ が一つの目安である。ただしこの $0.15$ という数値は、想起されたパタンにある程度のビット誤りを許容する、緩やかな条件に基づいて算出されている。より厳しく、1ビットの誤りなく正確に想起できる条件を課すならば、記憶率の限界値として $r_c = 1/(2 \ln (n))$ を導出できる [2]。例えば $n = 100$ のとき $r_c \approx 0.109$ 、 $n = 1000$ のとき $r_c \approx 0.0724$ となり、 $0.15$ と比べてずっと小さくなる。

想起のダイナミクスと連想記憶

時刻 $t$ における第 $i$ ニューロンの出力を $x_i(t)$ とする。時刻 $t+1$ におけるニューロンの値（状態）は次式で与えられる。

$\begin{align*} x_i(t+1)= \text{sgn}\left(\sum_{j=1}^{n} w_{ij} x_j(t) -\theta \right) \end{align*}$

ただし $\text{sgn} (u)$ は $u \gt 0$ のとき $1$ 、 $u \leq 0$ のとき $-1$ を取る関数である。また $\theta$ はしきい値であり、しばしば 0 に設定される。

あるいは時刻 $t$ におけるモデルの状態を $\mathbf{x}(t) = [ x_1(t), x_2(t), \ldots, x_n(t)]^{\top}$ とすれば、ベクトルと行列を用いた表記が得られる。

$\begin{align*} \mathbf{x} (t+1)= \text{sgn}\left(\mathbf{W} \mathbf{x}(t) -\theta \mathbf{1}_n \right) \end{align*}$

ただし $\text{sgn} (\cdot)$ はベクトルの要素ごとに適用するものとする。また $\mathbf{1}_n$ は全ての要素が1からなる $n$ 次元ベクトルである。

想起の際は入力パタン $\mathbf{x}(0) = \mathbf{x}_{0}$ を与え、上式に従って状態を更新する。全てのニューロンの状態が同じタイミングで一斉に更新されるので、同期型の状態更新とも呼ばれる。この状態更新則から生じる想起のダイナミクスこそが考察の対象となる。

状態更新を十分に繰り返すと、状態 $\mathbf{x}(t)$ はある平衡状態 $\mathbf{x}_e$ へと到達する。平衡状態においては更新則の適用前後で状態の変化が見られない。

$\begin{align*} \mathbf{x}_e = \text{sgn}\left(\mathbf{W} \mathbf{x}_e - \theta \mathbf{1}_n\right) \end{align*}$

$\mathbf{x}_e$ を初期パタン $\mathbf{x}_{0}$ に対する想起パタンと呼ぶ。もし初期値 $\mathbf{x}_{0}$ がある程度、記憶パタン $\mathbf{s} \in \{ \mathbf{s}^{(1)}, \ldots, \mathbf{s}^{(m)} \}$ に類似していれば、 $\mathbf{x}_e = \mathbf{s}$ への収束が期待できる。まさにこれが連想記憶である。特にネットワークの重みを表す記憶行列を状態ベクトルたちの（標本）自己相関行列により与えているため、このHopfieldモデルはまた自己相関型の連想記憶モデルと呼ばれる。

エネルギー関数と極小点

本記事ではこれまでエネルギー関数を定義してこなかったので、定義を書いておこう。

$\begin{align*} E (\mathbf{x})= -\frac{1}{2} \sum_{i=1}^{n} \sum_{j=1}^{n} w_{ij} x_i x_j + \theta \sum_{i=1}^{n} x_i \end{align*}$

同期型の状態更新則を上で述べたが、一度に1つのニューロンの状態のみが更新されるという、非同期型の状態更新則も存在する。本記事のように記憶行列が対称行列ならば、非同期型の更新則の場合はエネルギー関数の値が必ず減少する（ $\Delta E \lt 0$ ）ことが保証されている。しかし同期型の更新則の場合、その保証はない。

平衡状態はエネルギー関数の極小点（ $\nabla E = 0$ ）にほぼ対応している。極小点は唯一存在するとは限らず、一般には複数存在する。エネルギー関数の「地形」は一様な形状ではなく、極小点の周辺は山と谷が入り組んだ形状をしていることが知られている。それゆえ、収束に至るまでの状態の挙動もまた一様ではなく、初期パタンに最も近い極小点に収束していくとも限らない。本記事ではこれ以上の詳細な性質には立ち入らない。

状態間の類似度を測る尺度

想起過程を考える準備として、2つの状態 $\mathbf{x}$ と $\mathbf{y}$ の間の類似度を測るための距離 $d$ を定義しておく。

$\begin{align*} d(\mathbf{x}, \mathbf{y}) = \frac{1}{n} \sum_{i=1}^{n} \frac{| x_i - y_i |}{2} \end{align*}$

上式は正規化されたハミング距離に相当する。各状態変数ごとで考えると、 $1$ と $-1$ がそれぞれ等確率で現れるから、 $| x_i - y_i | / 2$ の期待値は 1/2 である。 $d$ を $n$ 個の確率変数の標本平均として考えれば、距離の期待値は依然として 1/2 である。明示的に書くならば、

$\begin{align*} \mathbb{E}[d(\mathbf{x}, \mathbf{y})] &= \mathbb{E}\left[\frac{1}{n} \sum_{i=1}^{n} \frac{| x_i - y_i |}{2} \right] \\ &= \frac{1}{n} \sum_{i=1}^{n} \mathbb{E}\left[ \frac{| x_i - y_i |}{2} \right] \\ & = \frac{1}{n} \sum_{i=1}^{n} \frac{1}{2} = \frac{1}{2} \end{align*}$

である。ひとつの状態 $\mathbf{y}$ を固定して考えたとき、多くの状態 $\mathbf{x}$ は $\mathbf{y}$ から距離0.5付近に分布している。

もう一つの類似度の尺度として、次式で定義する方向余弦 $a$ がある。

$\begin{align*} a(\mathbf{x}, \mathbf{y}) = \frac{1}{n} \sum_{i=1}^{n} x_{i} y_{i} = \frac{\mathbf{x} \cdot \mathbf{y}}{\| \mathbf{x} \| \| \mathbf{y} \|} \end{align*}$

$\mathbf{x} \cdot \mathbf{y}$ はユークリッド内積を表す。 $\| \mathbf{x} \| = \| \mathbf{y} \| = \sqrt{n}$ より、上式は明らかである。距離が大きいほど方向余弦は小さく、距離が小さいほど方向余弦は大きくなる。

距離と方向余弦の間には以下の関係式が成り立つ。

$\begin{align*} a(\mathbf{x}, \mathbf{y}) = 1- 2 \rho (\mathbf{x}, \mathbf{y}) \end{align*}$

【証明】

各ニューロンの状態で考えれば、 $(x_i, y_i)$ は $(1, 1)$ 、 $(1, -1)$ 、 $(-1, 1)$ 、 $(-1, -1)$ の4通りの組合せであり、それぞれについて

$x_{i} y_{i}$ は $1$ 、 $-1$ 、 $-1$ 、 $1$ の4通り
$| x_i - y_i | / 2$ は $0$ 、 $1$ 、 $1$ 、 $0$ の4通り

の値を取る。対応を順に見ていけば、

$\begin{align*} x_{i} y_{i} = 1 - 2 \frac{| x_i - y_i |}{2} \end{align*}$

の関係が成り立つこともわかる。そこで両辺を $i = 1, \ldots, n$ で和を取り、さらに $n$ で割れば関係式を得る。　（証明終わり）

想起過程

初期値からはじまり、状態更新則に従って $\mathbf{x}(t)$ が想起パタンに収束するまでの振る舞い、すなわち想起過程を観察したい。

あるパタン $\mathbf{s} \in \{ \mathbf{s}^{(1)}, \ldots, \mathbf{s}^{(m)} \}$ を想起対象として、 $\mathbf{s}$ から指定した距離 $\rho$ にある状態をランダムに選出し、初期パタン $\mathbf{x}_0$ を与える。つまり $\rho = d( \mathbf{x}_0, \mathbf{s})$ である。そして時刻 $t$ での方向余弦を

$\begin{align*} a_t (\mathbf{s}) = a(\mathbf{x}(t), \mathbf{s}) \end{align*}$

と定義し、 $a_t$ の時間変化を追うのである。

想起過程における初期パタンの具体的な与え方

シミュレーションでは $\rho$ 、すなわち初期方向余弦 $a_0 = 1 - 2\rho$ を様々に指定して実験する。この条件に合う初期パタン $\mathbf{x}_0$ を、 $\mathbf{s}$ からうまく作り出さねばならない。 $\mathbf{x}_0$ が $\mathbf{s}$ と一致していれば方向余弦の値は 1 を取り、 $\mathbf{s}$ の1つのニューロン値を反転させれば方向余弦は $2/n$ だけ減る。よって、いくつニューロン値を反転させれば所望の方向余弦値 $a_0$ が得られるか、その個数は簡単に計算できる。あとは $\mathbf{s}$ のうち値を反転させる箇所をランダムに選択すればよい。

具体例で考えよう。 $n = 100$ 、 $a_0 = 0.4$ とする。値が反転するニューロンの個数を $k$ と書くと、 $1 - 2k/n = a_0$ だから

$\begin{align*} k = \frac{1 - a_0}{2} n =\rho n = 0.3 \times 100 = 30 \end{align*}$

と計算できる。 $\mathbf{s} = [ s_1, s_2, \ldots, s_n ]^\top$ の100個のニューロンから30個をランダムに選択し、それらの値を反転させる。つまり30個のニューロン値にそれぞれ $-1$ を掛ければよい。以上の手続きを経て最終的に得られた状態が、 $a(\mathbf{x}_0, \mathbf{s}) = 0.4$ を満たす初期パタン $\mathbf{x}_0$ である。

実装

以下に置いた。Enjoy!

hopnet_dynamics.c

実装詳細を確認したいひとはここをクリック

ネットワーク関連の設定は HopfieldConfig 構造体に格納する。先日の実装から少しメンバ変数を変えている。

typedef struct {
    int num_neurons;              // Number of neurons
    int num_patterns;             // Number of patterns to be learned
    double **weights;             // Weights of network
    double threshold;             // Activation threshold for neurons
    int max_iterations;           // Maximum number of iterations for recall
    double similarity;            // Similarity of initial pattern
    bool self_connection;         // Flag to allow self-connection in weights
    char filename[FILENAME_MAX];  // Filename to save similarity data
} HopfieldConfig;

上記の構造体を受けて、オプション引数は以下に定めた。

オプション名	説明	デフォルト値
-n	ネットワークのニューロン数	20
-p	記憶対象のパタン数	3
-t	ニューロン発火のしきい値	0.0
-i	想起の最大繰り返し回数	100
-s	想起の初期値に課す類似度の制約	0.5
-c	ニューロンの自己結合を許容するか否か	false（許容しない）
-f	方向余弦（類似度）系列を保存するファイル名	log.txt

main関数では想起過程のダイナミクスを評価する eval_dynamics 関数が呼び出される。eval_dynamics の処理を抜粋すれば以下の通りである。

関数 init_state の呼び出しにより、想起対象パタン patterns[index_eval] から想起の初期状態 initial_state が決まる。
recall 関数が呼び出され、想起過程における方向余弦の系列（dircos_array）を取得する。
取得した系列は指定したファイルに保存される。

init_state(config, patterns[index_eval], initial_state);
recall(config, initial_state, patterns[index_eval], dircos_array);
for (int i = 0; i < max_iterations + 1; i++) {
    fprintf(fp, "%f\n", dircos_array[i]);
}

次に初期パタン（初期状態）を決定するための関数 init_state を示しておく。

config->similarity が初期方向余弦 $a_0$ であり、pattern が想起対象となるパタン $\mathbf{s}$ である。initial_state は初期パタン $\mathbf{x}_0$ である。先述した「初期パタンの具体的な与え方」に沿ったアルゴリズムの通りに、初期パタンを得ていることが分かる。ちなみに shuffle が配列の要素をシャッフルするための関数である。

void init_state(const HopfieldConfig *config, const int *pattern,
                int *initial_state) {
    /* Create a noisy initial state for a given pattern.

    Args:
        config (const HopfieldConfig *): The configuration of the network.
        pattern (const int *): The original pattern.
        initial_state (int *): The initial state of the network.
    */
    int num_neurons = config->num_neurons;
    double hamming_dist =
        (1 - config->similarity) / 2;  // a normalized Hamming distance
    int num_to_extract = (int)(num_neurons * hamming_dist);
    int numbers[num_neurons];
    for (int i = 0; i < num_neurons; i++) {
        numbers[i] = i;
    }
    shuffle(numbers, num_neurons);
    memcpy(initial_state, pattern, num_neurons * sizeof(int));
    for (int i = 0; i < num_to_extract; i++) {
        initial_state[numbers[i]] *= -1;
    }
}

なおプログラムを実行するごとに、記憶対象パタンをランダムに生成しなおす実装となっている。ゆえに初期類似度を様々に指定しながら実行する際、想起対象パタンも毎回ランダムに変わる点には注意されたい。

実験

実験を通じて、ニューロン数は1000、状態更新則のしきい値は0で固定する。また想起にかかる状態更新の回数は25とする。ニューロンの自己結合を許容しない場合と許容する場合でそれぞれ結果を示す。

実験の例を確認したいひとはこちらから（Google Colab）：

colab.research.google.com

このノートブックにはソースコードのダウンロード、シェルスクリプトのダウンロード＆実行、結果のプロットまでが含まれる。

実験結果その1：ニューロンの自己結合を許容しない場合

記憶パタン数を80に設定したときの方向余弦の時間変化を図1に示す。横軸に状態の更新回数、縦軸に方向余弦 $a_t$ の値を取っている。なお、このときの記憶率は 80 / 1000 = 0.08 < 0.15 であるから、ネットワークは十分にパタンを記憶できる設定である。初期パタンの類似度すなわち初期方向余弦を様々に変えて想起を試しており、それらに対応してグラフの色分けを行っている。

続いて記憶パタン数を200に設定したときの想起過程を図2に示す。このときの記憶率は 200 / 1000 = 0.20 > 0.15 であるから、ネットワークは全てのパタンを正確に記憶できない。

これらの図を観察して分かることは、

記憶率 $r$ が 0.15 を下回る程度に小さい場合（図1）：初期方向余弦 $a_0$ が一定レベル $a_c$ 以上ならば、方向余弦は速やかに1に収束するケースが多い。図1では $a_c$ はおよそ 0.20 〜 0.25 付近と推定される。 $a_0$ が $a_c$ を下回る場合は、ある程度まで増大するものの、 $a_0$ から離れた平衡状態に収束してしまう。このように、1への収束が期待される方向余弦の明確な「しきい値」が存在する。
記憶率 $r$ が 0.15 を上回る程度に大きい場合（図2）：初期方向余弦 $a_0$ をどれほど大きくしても、方向余弦は1に収束しない。言い換えれば、 $\mathbf{x}_t$ は想起対象パタン $\mathbf{s}$ に収束しない。1〜2回の状態更新で方向余弦は増大したのち、緩やかに減少する。

$a_c$ は $r$ に依存しており、 $r$ が大きいほど $a_c$ も大きくなる実験事実が知られている [4]。しかしながら $r > 0.15$ では $a_c$ は消失する。

実験結果その2：ニューロンの自己結合を許容する場合

記憶パタン数を80に設定したときの想起過程を図3に示す。

記憶パタン数を200に設定したときの想起過程を図4に示す。

図3の傾向は図1と同様であり、初期方向余弦が一定レベル以上ならば速やかに1へ収束するが、さもなくば他の平衡状態に収束する。図4も図2と傾向は同じであるが、初期方向余弦の値に応じた平衡状態へと細かく分かれて収束する様子が見て取れる。この現象は、ニューロンの自己結合を許容することで、記憶パタン以外の「偽状態」が多く発生したと解釈できる。

おわりに

ランダムビットパタン系列を連想記憶するように実装を作っていたおかげで、今回のようなシミュレーションにもスムーズに対応できた。シミュレーションの結果は引用した文献に掲載されているものと符号しており、再現実装に成功した。文献を何度も参照しながら行った実装と数値実験を通じて、モデルや理論に対する理解がより深まった。

平衡状態への収束の挙動や平衡点の分布について、記事の簡潔さを重視して割愛した。また非同期型の状態更新のとき、エネルギー関数が単調に減少するという、Hopfieldモデルの重要な性質については言及のみに留めた。続く小さい記事にて、エネルギー関数にかかるこれらのトピックについて解説する機会を持ちたい。

実験結果より、記憶率の増加に伴って想起性能が低下する事実を確認できたわけだが、この性能低下を改善する研究の方向性は自然である。文献 [4] では「部分反転法」と呼ばれる、記憶容量を増大させる手法が提案されている。今後の記事では部分反転法の実装を紹介する予定である。また記憶行列に疑似逆行列を用いる手法 [5] についても、実装を紹介する予定である。

本記事で示した状態更新則において時間パラメタは離散的であったが（離散時間モデル）、連続時間の更新則へと拡張することもできる（連続時間モデル）。文献 [3] には連続時間モデルの詳細な解析や記憶容量の計算アルゴリズム、また文献 [4] には連続時間モデルに部分反転法を適用した場合の結果が示されている。ゆえにまた、連続時間モデルの実装も今後の予定に含まれる。

なお今回は自己相関型の連想記憶モデルだったが、相互相関型のモデルも考えることができる。現在Transformerなどで活用されている注意機構の原型をそれらに見出すこともできよう。ちなみにKrotovらは2021年の論文 [6] において、「注意機構と生物学の神経細胞における連想記憶が、計算に関して類似している」[7] と指摘している。

参考文献

[1] Shun-Ichi Amari and Kenjiro Maginu, "Statistical neurodynamics of associative memory," Neural Networks, Volume 1, Issue 1, 1988, Pages 63-73.

[2] 甘利俊一，“神経回路網モデルとコネクショニズム，” 東京大学出版会, 2008.

[3] 甘利俊一編著，“ニューラルネットの新展開，” サイエンス社, 1993.

[4] 森田昌彦，吉澤修治，中野馨，“自己相関連想記憶の想起過程とその改良，” 電子情報通信学会論文誌, Vol.J73-D2, No.2, pp.232-242, 1990.

[5] Teuvo Kohonen, "Self-organization and associative memory: 3rd edition," Springer-Verlag, Berlin, Heidelberg, 1989.

[6] Dmitry Krotov and John Hopfield, "Large Associative Memory Problem in Neurobiology and Machine Learning," ICLR 2021. Link

[7] IBM ソリューションブログ，"謎に包まれていた脳の星状膠細胞にAIのトランスフォーマーが光を当てる", 2023. Link