WO2022013909A1

WO2022013909A1 - グラフ計算装置、グラフ計算方法及びグラフ計算プログラム

Info

Publication number: WO2022013909A1
Application number: PCT/JP2020/027218
Authority: WO
Inventors: 靖宏藤原
Original assignee: Nippon Telegraph and Telephone Corp
Current assignee: NTT Inc
Priority date: 2020-07-13
Filing date: 2020-07-13
Publication date: 2022-01-20
Anticipated expiration: 2023-01-13

Abstract

初期化部は、計算対象のＮ個のノードそれぞれについて、当該ノードを主ノードとし、Ｎ個のノードのうち前記主ノードからの距離に基づいて選択される（ｂ＋１）個以上（Ｎ－１）個未満のノードを従ノードとする複数のノードペアからなるノードペア集合を生成する。更新部は、ノードペア集合に係る各ノードペアについて、ノードペアに係る評価値を、ノードペアの距離、基準評価値、及び基準評価値に係るノードペアの距離に基づく計算により、所定の収束条件を満たすまで更新する。基準評価値は、ノードペアの従ノードを主ノードに持つ他のノードペアに係る評価値のうちｂ番目に小さい評価値である。エッジ決定部は、更新された評価値に基づいて、各主ノードについてエッジで接続するｂ個の従ノードを決定する。

Description

グラフ計算装置、グラフ計算方法及びグラフ計算プログラム

　本開示は、グラフ計算装置、グラフ計算方法及びグラフ計算プログラムに関する。

　多次元データから計算されるグラフとして、b-matchingグラフが知られている。b-matchingグラフは、多次元データである複数のノードそれぞれにｂ個のエッジが接続されるグラフである。非特許文献１には、確率伝搬法（belief propagation）に基づいてb-matchingグラフを求める手法が開示されている。

Tony Jebara, Jun Wang, Shih-Fu Chang, Graph construction and b-matching for semi-supervised learning, ICML, 2009.

　非特許文献１に記載の手法は、Ｍ次元のベクトルで表されるＮ個のノードから選択されるすべてのノードペアについて、収束するまで確率伝搬法に係るメッセージを更新するものである。そのため、収束までの繰り返し回数をＴとおくと、計算コストがＯ（Ｎ^２Ｍ^２Ｔ）と高いことが分かる。
　本開示の目的は、計算コストを抑えてb-matchingグラフを計算するグラフ計算装置、グラフ計算方法及びグラフ計算プログラムを提供することにある。

　一態様によれば、グラフ計算装置は、計算対象のＮ個のノードそれぞれについて、当該ノードを主ノードとし、前記Ｎ個のノードのうち前記主ノードからの距離に基づいて選択される（ｂ＋１）個以上（Ｎ－１）個未満のノードを従ノードとする複数のノードペアからなるノードペア集合を生成する初期化部と、前記ノードペア集合に係る各ノードペアについて、前記ノードペアに係る評価値を、前記ノードペアの距離、前記ノードペアの従ノードを主ノードに持つ他のノードペアに係る評価値のうちｂ番目に小さい評価値である基準評価値、及び前記基準評価値に係るノードペアの距離に基づく計算により、所定の収束条件を満たすまで更新する更新部と、更新された前記評価値に基づいて、各主ノードについてエッジで接続するｂ個の従ノードを決定するエッジ決定部とを備える。

　上記態様によれば、グラフ計算装置は、計算コストを抑えてb-matchingグラフを計算することができる。

第１の実施形態に係るグラフ計算装置の構成を示す概略ブロック図である。第１の実施形態に係るb-matchingグラフの計算方法を示すフローチャートである。第１の実施形態に係るノードペア及びメッセージの初期化方法を示すフローチャートである。第１の実施形態に係るメッセージの更新方法を示すフローチャートである。第１の実施形態に係るエッジのウエイト計算方法を示すフローチャートである。グラフ計算装置のコンピュータ構成を示す概略ブロック図である。

〈第１の実施形態〉
《グラフ計算装置１０の構成》
　以下、図面を参照しながら実施形態について詳しく説明する。
　図１は、第１の実施形態に係るグラフ計算装置１０の構成を示す概略ブロック図である。
　グラフ計算装置１０は、Ｍ次元のベクトルであるＮ個のノードについて、各ノードにウエイトを有するｂ本のエッジが接続されるb-matchingグラフを計算する。
　グラフ計算装置１０は、入力部１１、次元削減部１２、ベクトル記憶部１３、初期化部１４、ペア記憶部１５、更新部１６、継続判定部１７、追加部１８、エッジ決定部１９、ウエイト計算部２０、出力部２１を備える。

　入力部１１は、グラフの計算対象となるＮ個のノードの入力を受け付ける。各ノードは、Ｍ次元のベクトル（ｘ_ｉ＝ｘ_ｉ［１］、…ｘ_ｉ［Ｍ］）で表される多次元データである。Ｎ及びＭは自然数である。以下、入力部１１に入力されたベクトルをロー（Raw）ベクトルともいう。

　次元削減部１２は、入力部１１に入力された各ローベクトル（ｘ_１～ｘ_Ｎ）の次元を削減する。次元削減部１２は、例えば特異値分解（Singular Value Decomposition）により各ノードの次元を削減する。次元削減後のベクトルである低次元ベクトルは、ｍ次元のベクトルである。ｍはＭ未満の自然数である。具体的には、次元削減部１２は、以下の手順で低次元ベクトルを生成する。まず、次元削減部１２は、Ｎ個のローベクトル（ｘ_１～ｘ_Ｎ）をまとめたＮ×Ｍのノード行列Ｘを特異値分解することで、ユニタリ行列、特異値行列及び随伴行列を得る。次元削減部１２は、ユニタリ行列からＮ×ｍの部分行列を、特異値行列からｍ×ｍの部分行列を抽出し、これらを乗算することで、Ｎ×ｍの低次元行列を生成する。次元削減部１２は、当該低次元行列の各行ベクトルを、低次元ベクトル（ｘ^～ _１～ｘ^～ _Ｎ：ただし、ｘ^～は、ｘの上にチルダが乗ったもの）として抽出する。また、次元削減部１２は、ローベクトルと低次元ベクトルの二乗平均誤差を表す残差ｒを以下の式（１）により算出する。

　ここで、｜｜・｜｜は、ベクトルのＬ２ノルムを示す。

　ベクトル記憶部１３は、Ｎ個のノードそれぞれに関連付けて、当該ノードのローベクトルと低次元ベクトルと残差とを記憶する。

　初期化部１４は、Ｎ個のノードから、確率伝搬法の計算対象となるノードペアを抽出し、抽出したノードペアについてメッセージの初期値を生成する。確率伝搬法のメッセージは、評価値の一例である。
　ペア記憶部１５は、計算対象となるノードペアそれぞれに関連付けて、当該ノードペアのメッセージを記憶する。

　更新部１６は、ペア記憶部１５が記憶する各ノードペアに係るメッセージを更新する。
　継続判定部１７は、更新部１６によるメッセージの更新が収束したか否かを判定する。
　追加部１８は、更新部１６が更新したメッセージに基づいて、ペア記憶部１５が記憶されないノードペアである除外ペアを、ペア記憶部１５に追加するか否かを判定する。

　エッジ決定部１９は、ペア記憶部１５が記憶するメッセージに基づいて、エッジで接続するノードペアを決定する。
　ウエイト計算部２０は、ペア記憶部１５が記憶するメッセージに基づいて、漸進的にエッジのウエイトを計算する。
　出力部２１は、エッジ決定部１９及びウエイト計算部２０の計算結果に基づいて、b-matchingグラフを出力する。

《グラフ計算装置１０の処理》
　図２は、第１の実施形態に係るb-matchingグラフの計算方法を示すフローチャートである。
　入力部１１は、計算対象となるＮ個のノード（ｘ_１～ｘ_Ｎ）と各ノードに係るローベクトルの入力を受け付ける（ステップＳ１）。次元削減部１２は、ステップＳ１で取得したＮ個のローベクトルから、Ｎ個の低次元ベクトル（ｘ^～ _１～ｘ^～ _Ｎ）と残差ｒとを算出する（ステップＳ２）。次元削減部１２は、ノードに関連付けて、ローベクトル、低次元ベクトル、及び残差をベクトル記憶部１３に記録する（ステップＳ３）。

　次に、初期化部１４は、Ｎ個のノードから計算対象となるノードペア及びメッセージの初期値を生成する（ステップＳ４）。初期化部１４は、以下の手順で初期化処理を行う。
　図３は、第１の実施形態に係るノードペア及びメッセージの初期化方法を示すフローチャートである。
　初期化部１４は、メッセージの計算対象とするノードペアの集合を初期化する。具体的には、初期化部１４は、ペア記憶部１５に、Ｎ個のノードそれぞれについて、当該ノードを主ノードとするノードペア集合の記憶領域を確保する（ステップＳ１０１）。ここでは、ノードペアに係る２つのノードのうちメッセージの計算において基準となるものを便宜的に主ノードと呼ぶ。のなお、この時点においてすべてのノードペア集合は空集合である。

　初期化部１４は、Ｎ個のノードを１つずつ選択し（ステップＳ１０２）、以下のステップＳ１０３からステップＳ１１０の処理を実行する。
　初期化部１４は、選択したノード以外の（Ｎ－１）個のノード（候補ノード）からなるノード集合を、候補ノード集合として抽出する（ステップＳ１０３）。初期化部１４は、以下の式（２）により、ベクトル記憶部１３が記憶するローベクトルに基づいて選択したノードと各候補ノードとの距離を算出する（ステップＳ１０４）。

　初期化部１４は、候補ノード集合のうち、主ノードとの距離が最も短い候補ノードを特定する（ステップＳ１０５）。初期化部１４は、ペア記憶部１５を参照し、当該候補ノードを主ノードとするノードペアの数が（ｂ＋１）未満であるか否かを判定する（ステップＳ１０６）。候補ノードを主ノードとするノードペアの数が（ｂ＋１）未満である場合（ステップＳ１０６：ＹＥＳ）、初期化部１４は、ステップＳ１０２で選択したノードを主ノードとし、ステップＳ１０５で特定した候補ノードを従ノードするノードペアを、ペア記憶部１５に記録する（ステップＳ１０７）。また、初期化部１４は、ステップＳ１０５で特定した候補ノードを主ノードとし、ステップＳ１０２で選択したノードを従ノードするノードペアを、ペア記憶部１５に記録する（ステップＳ１０８）。このとき、初期化部１４は、各ノードペアに関連付けて初期値１のメッセージをペア記憶部１５に記録する。

　候補ノードを主ノードとするノードペアの数が（ｂ＋１）以上である場合（ステップＳ１０６：ＮＯ）、またはペア記憶部１５にノードペアを記録した場合（ステップＳ１０７及びＳ１０８）、初期化部１４は、候補ノード集合からステップＳ１０５で特定した候補ノードを削除する（ステップＳ１０９）。
　初期化部１４は、ステップＳ１０２で選択したノードを主ノードとするノードペアの数が（ｂ＋１）未満であるか否かを判定する（ステップＳ１１０）。ステップＳ１０２で選択したノードを主ノードとするノードペアの数が（ｂ＋１）未満である場合（ステップＳ１１０：ＹＥＳ）、初期化部１４は、ステップＳ１０５に処理を戻す。他方、ステップＳ１０２で選択したノードを主ノードとするノードペアの数が（ｂ＋１）以上である場合（ステップＳ１１０：ＮＯ）、初期化部１４は、ステップＳ１０２に戻り、次のノードを選択する。

　上記手順により、初期化部１４は、ノードペア集合及びメッセージを初期化することができる。なお、上述の通り、初期化部１４は、ペア記憶部１５には、Ｎ（ｂ＋１）個のノードペアが記録される。Ｎ個のノードから特定可能なすべてのノードペアの数は、Ｎ（Ｎ－１）個であるため、初期化部１４によって、（Ｎ－ｂ－２）個のノードペアが、計算対象外となったことが分かる。これにより、グラフ計算装置１０は、b-matchingグラフの計算コストを低減することができる。

　図２に示すように、初期化部１４がステップＳ４で、Ｎ個のノードから計算対象となるノードペア及びメッセージの初期値を生成すると、更新部１６は、ペア記憶部１５が記憶する各メッセージを更新する（ステップＳ５）。

　更新部１６は、以下の手順で更新処理を行う。図４は、第１の実施形態に係るメッセージの更新方法を示すフローチャートである。
　更新部１６は、ペア記憶部１５が記憶するＮ個のノードを１つずつ選択し（ステップＳ２０１）、以下のステップＳ２０２からステップＳ２０４の処理を行う。

　まず更新部１６は、ステップＳ２０１で選択したノードｘ_ｉを主ノードとするノードペア（ｘ_ｉ，ｘ_ｊ）について、以下の式（２）に示すように、ｅを底とする当該ノードペアの距離Ｄ［ｉ，ｊ］の指数関数と、メッセージｍ［ｘ_ｉ，ｘ_ｊ］とを乗算することで、重み付きメッセージｍ^ｗ［ｘ_ｉ，ｘ_ｊ］を求める（ステップＳ２０２）。

　更新部１６は、ステップＳ２０１で選択したノードｘ_ｉに係る重み付きメッセージｍ^ｗ［ｘ_ｉ，ｘ_ｊ］のうち、ｂ番目に小さいものμ［ｘ_ｉ］と、（ｂ＋１）番目に小さいものν［ｘ_ｉ］を特定する（ステップＳ２０３）。以下、主ノードｘ_ｉに係るｂ番目に小さい重み付きメッセージを、ノードｘ_ｉに係る下側基準メッセージμ［ｘ_ｉ］と呼ぶ。また、主ノードｘ_ｉに係る（ｂ＋１）番目に小さい重み付きメッセージを、ノードｘ_ｉに係る上側基準メッセージν［ｘ_ｉ］と呼ぶ。
　また、更新部１６は、下側基準メッセージμ［ｘ_ｉ］に係る従ノードv［ｘ_ｉ］及び上側基準メッセージν［ｘ_ｉ］に係る従ノードｗ［ｘ_ｉ］を特定する（ステップＳ２０４）。以下、ノードｘ_ｉに係る下側基準メッセージμ［ｘ_ｉ］に係る従ノードを、ノードｘ_ｉに係る下側基準ノードｖ［ｉ］とよぶ。また、ノードｘ_ｉに係る上側基準メッセージν［ｘ_ｉ］に係る従ノードを、ノードｘ_ｉに係る上側基準ノードｗ［ｉ］とよぶ。つまり、下側基準メッセージμ［ｘ_ｉ］及び上側基準メッセージν［ｘ_ｉ］は、以下の式（３）のように表せる。

　更新部１６は、ペア記憶部１５が記憶するノードペアを１つずつ選択し（ステップＳ２０５）、以下のステップＳ２０６からステップＳ２０８の処理を行う。
　更新部１６は、ステップＳ８で選択したノードペアに係る主ノードｘ_ｉが、当該ノードペアに係る従ノードｘ_ｊにおける下側基準ノードｖ［ｊ］であるか否かを判定する（ステップＳ２０６）。主ノードｘ_ｉが下側基準ノードｖ［ｊ］でない場合（ステップＳ２０６：ＮＯ）、更新部１６は、式（４）に示すように、ステップＳ２０５で選択したノードペアに係るメッセージを、ｅを底とする当該ノードペアの距離Ｄ［ｉ，ｊ］の指数関数を、従ノードｘ_ｊに係る下側基準メッセージμ［ｘ_ｊ］で除算した値に書き換える（ステップＳ２０７）。なお、ノードｘ_ｉが下側基準ノードｖ［ｊ］でない場合、下側基準メッセージμ［ｘ_ｊ］は、ノードｘ_ｊを主ノードに持つ他のノードペアに係る評価値のうちｂ番目に小さい評価値と、ｅを底とする当該ノードペアの距離の指数関数との積であるといえる。

　：＝は代入演算子である。

　他方、主ノードｘ_ｉが下側基準ノードｖ［ｊ］である場合（ステップＳ２０６：ＹＥＳ）、更新部１６は、式（５）に示すように、ステップＳ８で選択したノードペアに係るメッセージを、ｅを底とする当該ノードペアの距離Ｄ［ｉ，ｊ］の指数関数を、従ノードｘ_ｊに係る上側基準メッセージν［ｘ_ｊ］で除算した値に書き換える（ステップＳ２０８）。なお、主ノードｘ_ｉが下側基準ノードｖ［ｊ］である場合、上側基準メッセージν［ｘ_ｊ］は、ノードｘ_ｊを主ノードに持つ他のノードペアに係る評価値のうちｂ番目に小さい評価値と、ｅを底とする当該ノードペアの距離の指数関数との積であるといえる。

　上記手順により、更新部１６は、候補ノード集合及びメッセージを初期化することができる。
　なお、更新部１６は、ステップＳ２０５～Ｓ２０８によるメッセージの更新前に、予めステップＳ２０１～Ｓ２０４で、ノードごとに下側基準メッセージ、上側基準メッセージ、下側基準ノード、上側基準ノードを特定しておく。これにより、更新部１６は、更新対象のメッセージに係るノードペアの従ノードｘ_ｊを主ノードに持つ他のノードペアに係る評価値のうちｂ番目に小さい評価値を都度計算する必要がないため、計算コストを低減することができる。

　図２に示すように、更新部１６がステップＳ５で、メッセージを更新すると、継続判定部１７は、メッセージが収束したか否かを判定する（ステップＳ６）。継続判定部１７は、例えばメッセージの変化量の二乗平均が所定値未満になった場合に、メッセージが収束したと判定する。なお、メッセージは、Ｏ（Ｎ）回の繰り返しで収束する。
　メッセージが収束していない場合（ステップＳ６：ＮＯ）、グラフ計算装置１０はステップＳ５に処理を戻し、メッセージの更新を繰り返す。他方、メッセージが収束した場合（ステップＳ６：ＹＥＳ）、追加部１８は、Ｎ個のノードから特定可能なすべてのノードペアのうち、ペア記憶部１５が記憶しないノードペアすなわち計算対象外となったノードペア（以下、除外ペアと呼ぶ）を１つずつ選択し（ステップＳ７）、以下のステップＳ８からステップＳ１２の処理を実行する。

　追加部１８は、ベクトル記憶部１３から、ステップＳ７で選択した除外ペア（ｘ_ｉ，ｘ_ｊ）に係る主ノードの低次元ベクトルｘ^～ _ｉ及び残差ｒ［ｘ_ｉ］と従ノードの低次元ベクトルｘ^～ _ｊ及び残差ｒ［ｘ_ｊ］を読み出し、式（６）により近似距離Ｄ^～［ｉ，ｊ］を算出する（ステップＳ８）。

　〈・，・〉は、２つのベクトルの内積を示す。

　追加部１８は、主ノードｘ_ｉに係る上側基準メッセージν［ｘ_ｉ］、従ノードｘ_ｊに係る上側基準メッセージν［ｘ_ｊ］、及びステップＳ８で算出した近似距離Ｄ^～［ｉ，ｊ］に基づいて、除外ペアが式（７）に示す第１追加条件を満たすか否かを判定する（ステップＳ９）。第１追加条件は、除外ペアがb-matchingグラフにおいてエッジで接続される可能性があるか否かを判断するための条件である。

　選択された除外ペアが第１追加条件を満たす場合（ステップＳ９：ＹＥＳ）、追加部１８は、除外ペアに係る距離Ｄ［ｉ，ｊ］を算出する（ステップＳ１０）。追加部１８は、主ノードｘ_ｉに係る上側基準メッセージν［ｘ_ｉ］、従ノードｘ_ｊに係る上側基準メッセージν［ｘ_ｊ］、及びステップＳ１０で算出した距離Ｄ［ｉ，ｊ］に基づいて、除外ペアが式（８）に示す第２追加条件を満たすか否かを判定する（ステップＳ１１）。

　第２追加条件は、除外ペアがb-matchingグラフにおいてエッジで接続される可能性があるか否かを厳格に判断するための条件である。つまり、すべての除外ペアが第２追加条件を満たさない場合、追加部１８はb-matchingグラフの最適解が求められたと判断することができる。第１追加条件及び第２追加条件によって最適解の判定ができる理由については、後述する。
　なお、距離Ｄ［ｉ，ｊ］の計算コストはＯ（Ｍ）であり、近似距離Ｄ^～［ｉ，ｊ］の計算コストＯ（ｍ）より高い。つまり、追加部１８は、第１追加条件の判定を行うことで、第２追加条件の判定に必要な計算量の多い距離Ｄ［ｉ，ｊ］の計算回数を少なくすることができる。

　選択された除外ペアが第１追加条件を満たす場合（ステップＳ１１：ＹＥＳ）、追加部１８は、式（４）に基づいて当該除外ペアに係るメッセージを生成し、当該除外ペア及びそのメッセージを、ペア記憶部１５に記録する（ステップＳ１２）。つまり、追加部１８は、除外ペアをノードペア集合に追加する。

　継続判定部１７は、ノードペア集合に追加された除外ペアがあるか否かを判定する（ステップＳ１３）。ノードペア集合に追加された除外ペアがある場合（ステップＳ１３：ＹＥＳ）、グラフ計算装置１０はステップＳ５に処理を戻し、メッセージの更新を繰り返す。
　他方、ノードペア集合に追加された除外ペアがない場合（ステップＳ１３：ＮＯ）エッジ決定部１９は、Ｎ個のノードそれぞれについて、ペア記憶部１５が記憶する当該ノードを主ノードとするノードペアのうち、メッセージｍ［ｘ_ｉ，ｘ_ｊ］が小さいものからｂ個目までのノードペアを、エッジで接続するノードペアに決定する（ステップＳ１４）。

　次に、ウエイト計算部２０は、決定された各エッジについて、ウエイトを計算する（ステップＳ１５）。ウエイト計算部２０は、以下の手順で初期化処理を行う。
　図５は、第１の実施形態に係るエッジのウエイト計算方法を示すフローチャートである。
　まず、ウエイト計算部２０は、Ｎ個のノードを１つずつ選択し（ステップＳ３０１）、以下のステップＳ３０２からステップＳ３０９の処理を実行する。

　ウエイト計算部２０は、ステップＳ３０１で選択したノードに対する少なくとも１つの従ノードについて、ウエイトの計算が完了しているか否かを判定する（ステップＳ３０２）。ウエイトの計算が完了している従ノードがない場合（ステップＳ３０２：ＮＯ）、ウエイト計算部２０は、ステップＳ３０１で選択したノードを主ノードとするノードペアに係る主ノードのベクトルｘ_ｉと従ノードのベクトルｘ_ｊの内積を要素とするｂ次元の内積ベクトルｑ_ｉと、前記ｂ個の従ノードのベクトルに係るｂ×ｂのグラム行列Ｐ_ｉの逆行列Ｐ_ｉ ^－１を直接算出する（ステップＳ３０３）。内積ベクトルｑ_ｉの要素ｑ_ｉ［ｊ］は、以下の式（９）で表され、グラム行列Ｐ_ｉの要素Ｐ_ｉ［ｊ，ｋ］は、以下の式（１０）で表される。

　他方、ウエイトの計算が完了している従ノードがある場合（ステップＳ３０２：ＹＥＳ）、ウエイト計算部２０は、ウエイトの計算が完了しているノードのうち、当該ノードを主ノードとするノードペアと、ステップＳ３０１で選択したノードを主ノードとするノードペアとで共通する従ノードの数が最も多いものを特定する（ステップＳ３０４）。すなわち、以下の式（１１）で示されるノードｙ_ｉを特定する。以下、共通する従ノードの数が最も多いウエイト計算済みのノードを、類似ノードｘ_ｉ´とよぶ。

　ここで、集合Ｃ［ｘ］は、ノードｘを主ノードとするノードペアに係る従ノードの集合である。以下、集合Ｃ［ｘ］を、ノードｘに係る従ノード集合と呼ぶ。また集合Ｄは、ウエイトの計算が完了している従ノードの集合である。式（１１）において、ノードｘ_ｉは、ステップＳ３０１で選択したノードである。

　ウエイト計算部２０は、類似ノードｘ_ｊを主ノードとするノードペアと、ステップＳ３０１で選択したノードを主ノードとするノードペアとで、共通しない従ノードを１組特定する（ステップＳ３０５）。具体的には、類似ノードｘ_ｊを主ノードとするノードペアに係る従ノードが、ｘ_１、ｘ_２、ｘ_３、ｘ_４、ｘ_５であり、ステップＳ３０１で選択したノードを主ノードとするノードペアに係る従ノードがｘ_３、ｘ_４、ｘ_５、ｘ_６、ｘ_７である場合、（ｘ_１、ｘ_２）と（ｘ_６、ｘ_７）とが共通しない従ノードである。このとき、ウエイト計算部２０は、例えばｘ_１とｘ_６の組み合わせを特定する。

　ウエイト計算部２０は、類似ノードｘ_ｉ´に係る従ノード集合Ｃ［ｘ_ｉ´］の共通しない１つの従ノードｘ_ｊ´を、ステップＳ３０１で選択されたノードｘ_ｉに係る従ノード集合Ｃ［ｘ_ｉ］の共通しない１つのノードｘ_ｊと入れ替えた集合を、更新類似集合Ｃ［ｘ_ｉ″］とおいた場合に、更新類似集合Ｃ［ｘ_ｉ″］について、以下の式（１２）で示されるベクトルｑ_ｊ″と、式（１３）で示される行列Ｐ_ｊ″ ^－１を算出する（ステップＳ３０６）。

　ウエイト計算部２０は、更新類似集合Ｃ［ｘ_ｉ″］と、ステップＳ３０１で選択されたノードｘ_ｉに係る従ノード集合とが一致しているか否かを判定する（ステップＳ３０７）。更新類似集合Ｃ［ｘ_ｉ″］と従ノード集合Ｃ［ｘ_ｉ］とが一致しない場合（ステップＳ３０７：ＮＯ）、ウエイト計算部２０は、集合Ｃ［ｘ_ｉ″］を更新類似集合Ｃ［ｘ_ｉ″］に置き換えて、ステップＳ３０５に処理を戻し、漸進的にベクトルｑ_ｊ″と行列Ｐ_ｊ″を算出する。他方、更新類似集合Ｃ［ｘ_ｉ″］と従ノード集合Ｃ［ｘ_ｉ］とが一致する場合（ステップＳ３０７：ＹＥＳ）、ウエイト計算部２０は、ベクトルｑ_ｊ″と行列Ｐ_ｊ″ ^－１を、それぞれ内積ベクトルｑ_ｉ、グラム行列Ｐ_ｉの逆行列Ｐ_ｉ ^－１とみなす。式（１１）、式（１２）により、漸進的に内積ベクトルｑ_ｉ及びグラム行列Ｐ_ｉの逆行列Ｐ_ｉ ^－１を求められることは、Sherman-Morrisonの公式から明らかである。

　なお、式（１２）及び式（１３）により漸進的に内積ベクトルｑ_ｉ、グラム行列Ｐ_ｉの逆行列Ｐ_ｉ ^－１を求める計算の計算量は、式（９）及び式（１０）により内積ベクトルｑ_ｉ及びグラム行列Ｐ_ｉの逆行列Ｐ_ｉ ^－１を直接求める計算の計算量より少ない。具体的には、漸進的に内積ベクトルｑ_ｉ、グラム行列Ｐ_ｉの逆行列Ｐ_ｉ ^－１を求める計算の計算量はＯ（ｂ^２＋ｂＭ）であるのに対し、式（９）及び式（１０）の計算量は、Ｏ（ｂ^３＋ｂ^２Ｍ）である。

　次に、ウエイト計算部２０は、グラム行列Ｐ_ｉの逆行列Ｐ_ｉ ^－１と内積ベクトルｑ_ｉの積を求めることで、ウエイトベクトルｗ_ｉを算出する（ステップＳ３０８）。すなわち、ウエイト計算部２０は、以下の式（１４）により、ウエイトベクトルｗ_ｉを算出する。

　ウエイト計算部２０は、ウエイトベクトルｗ_ｉのすべての要素が非負であるか否かを判定する（ステップＳ３０９）。ウェイトベクトルｗ_ｉのすべての要素が正である場合（ステップＳ３０９：ＹＥＳ）、ウエイト計算部２０は、ウエイトベクトルｗ_ｉを、最小値が０、要素の総和が１となるように正規化することで、ステップＳ３０１で選択したノードに接続する各エッジのウエイトを算出する（ステップＳ３１０）。

　他方、ウェイトベクトルｗ_ｉの少なくとも１つの要素が負である場合（ステップＳ３０９：ＮＯ）、ウエイト計算部２０は、ステップＳ３０１で選択したノードに接続するエッジのうち、ウエイトが正数であるものを、二次計画法によるウエイトの計算対象の対象エッジ集合として抽出する（ステップＳ３１１）。そして、ウエイト計算部２０は、ウエイトベクトルｗ_ｉの要素のうち計算対象として抽出されないものの値を０に書き換え、計算対象として抽出されないものの値を、最小値が０、要素の総和が１となるように正規化する（ステップＳ３１２）。

　ウエイト計算部２０は、ステップＳ３１１で抽出した対象エッジ集合について、二次計画法によりウエイトを計算する（ステップＳ３１３）。ウエイト計算部２０は、ステップＳ３０１で選択したノードに接続するエッジのうち、対象エッジ集合に含まれないものについて、従ノードｘ_ｊが以下の式（１５）を満たすものがあるか否かを判定する（ステップＳ３１４）。

　対象エッジ集合に含まれないエッジであって式（１５）を満たすものが存在する場合（ステップＳ３１４：ＹＥＳ）、当該エッジを対象エッジ集合に追加し（ステップＳ３１５）、ステップＳ３１２に戻り、二次計画法によるウエイトの再計算を行う。他方、対象エッジ集合に含まれないエッジであって式（１５）を満たすものが存在しない場合（ステップＳ３１４：ＹＥＳ）、エッジのウエイトをステップＳ３１３で算出したウエイトベクトルｗ_ｉの値に決定する。

　上述の手順でウエイト計算部２０が各エッジのウエイトを算出すると、図２に示すように、出力部２１は、ステップＳ１４及びステップＳ１５で決定したエッジ及びウエイトに基づいて、b-matchingグラフを出力する（ステップＳ１６）。

《作用・効果》
　このように、第１の実施形態によれば、グラフ計算装置１０は、Ｎ個のノードそれぞれについて、距離に基づいて選択される（ｂ＋１）個のノードペアを初期集合とするノードペア集合について、評価値であるメッセージの計算を行うことで、ｂ個のエッジを決定する。つまり、総当たりで計算する場合にメッセージの計算対象となるノードペアの数が（Ｎ－１）個であるところ、第１の実施形態によれば、（ｂ＋１）個まで低減することができる。これにより、エッジの決定のための計算コストを抑えてb-matchingグラフを計算することができる。なお、第１の実施形態では、グラフ計算装置１０は、（ｂ＋１）個のノードペアを初期集合として抽出するが、これに限られない。例えば、他の実施形態においては、グラフ計算装置１０は、（ｂ＋２）個以上（Ｎ－１）個未満のノードペアを初期集合として抽出してもよい。ただし、初期集合に含まれるノードペアの数が少ないほど計算コストが小さくなる。また、他の実施形態に係るグラフ計算装置１０は、ローベクトルに基づく距離に代えて、低次元ベクトル及び残差から算出される近似距離に基づいてノードペア集合を初期化してもよい。近似距離は、距離の一例である。

　第１の実施形態によれば、グラフ計算装置１０は、ノードペア集合に含まれない除外ペアのうち式（８）を満たすものをノードペア集合に追加する。これにより、エッジによる接続対象となるノードペアの候補をグラフ計算装置１０は漏れなく抽出することができる。
　ここで、第１の実施形態に係るグラフ計算装置１０は、近似距離を用いて式（７）に示す第１追加条件を満たす除外ペアについて、さらに距離を用いた式（８）に示す第２追加条件を満たす場合に、当該除外ペアをノードペア集合に追加する。この２段階の判定により、グラフ計算装置１０は、ノードペアの追加の判定に係る計算量を低減することができる。理由は以下の通りである。近似距離は、要素数ｍ個の低次元ベクトルに係る計算であるため、要素数Ｍ（Ｍ＞ｍ）個のベクトルに係る距離の計算より計算量が少ない。また、近似距離の値は、必ず距離の値以下となるため、グラフ計算装置１０は、第１追加条件の判定をすることによって、第２追加条件を満たし得る除外ペアを選別することができる。なお、他の実施形態においては、これに限られず、第１追加条件の判定を行わずにすべての除外ペアについて第２追加条件の判定を行ってもよい。この場合に、例えば、初期化部１４が初期化の際に計算した距離を記憶しておき、これを読み出すことで、都度の計算を省略してもよい。

　第１の実施形態によれば、グラフ計算装置１０は、従ノードのベクトルに係るｂ×ｂのグラム行列に基づいてエッジのウエイトを算出し、当該グラム行列の逆行列を漸進的に計算することで、計算コストを低減することができる。また、グラフ計算装置１０は、グラム行列に基づいてウエイトを計算できない場合に、二次計画法によってウエイトを計算するが、グラム行列に基づくウエイトが負数のエッジを二次計画法の計算対象から除外することで、二次計画法の対象となるエッジの数を低減することができる。これは、グラム行列に基づいて計算されるウエイトが負数になる場合、そのエッジに係るウエイトがゼロとなることが多いという発明者の知見に基づく。

〈他の実施形態〉
　以上、図面を参照して一実施形態について詳しく説明してきたが、具体的な構成は上述のものに限られることはなく、様々な設計変更等をすることが可能である。すなわち、他の実施形態においては、上述の処理の順序が適宜変更されてもよい。また、一部の処理が並列に実行されてもよい。
　上述した実施形態に係るグラフ計算装置１０は、単独のコンピュータによって構成されるものであってもよいし、グラフ計算装置１０の構成を複数のコンピュータに分けて配置し、複数のコンピュータが互いに協働することでグラフ計算装置１０として機能するものであってもよい。

〈コンピュータ構成〉
　図６は、グラフ計算装置１０のコンピュータ構成を示す概略ブロック図である。
　グラフ計算装置１０は、バスで接続されたプロセッサ５１、メモリ５３、補助記憶装置５５、インタフェース５７などを備え、グラフ計算プログラムを実行することによって、入力部１１、次元削減部１２、ベクトル記憶部１３、初期化部１４、ペア記憶部１５、更新部１６、継続判定部１７、追加部１８、エッジ決定部１９、ウエイト計算部２０及び出力部２１を備える装置として機能する。プロセッサ５１の例としては、ＣＰＵ（Central Processing Unit）、ＧＰＵ（Graphic Processing Unit）、マイクロプロセッサなどが挙げられる。
　グラフ計算プログラムは、補助記憶装置５５などのコンピュータ読み取り可能な記録媒体に記録されてもよい。コンピュータ読み取り可能な記録媒体とは、例えば磁気ディスク、光磁気ディスク、光ディスク、半導体メモリ等の記憶装置である。グラフ計算プログラムは、電気通信回線を介して送信されてもよい。
　なお、グラフ計算の各機能の全て又は一部は、ＡＳＩＣ（Application Specific Integrated Circuit）やＰＬＤ（Programmable Logic Device）等のカスタムＬＳＩ（Large Scale Integrated Circuit）を用いて実現されてもよい。ＰＬＤの例としては、ＰＡＬ(Programmable Array Logic)、ＧＡＬ(Generic Array Logic)、ＣＰＬＤ(Complex Programmable Logic Device)、ＦＰＧＡ（Field Programmable Gate Array）が挙げられる。このような集積回路も、プロセッサ５１の一例に含まれる。

１０　グラフ計算装置
１１　入力部
１２　次元削減部
１３　ベクトル記憶部
１４　初期化部
１５　ペア記憶部
１６　更新部
１７　継続判定部
１８　追加部
１９　エッジ決定部
２０　ウエイト計算部
２１　出力部

Claims

　計算対象のＮ個のノードそれぞれについて、当該ノードを主ノードとし、前記Ｎ個のノードのうち前記主ノードからの距離に基づいて選択される（ｂ＋１）個以上（Ｎ－１）個未満のノードを従ノードとする複数のノードペアからなるノードペア集合を生成する初期化部と、
　前記ノードペア集合に係る各ノードペアについて、前記ノードペアに係る評価値を、前記ノードペアの距離、前記ノードペアの従ノードを主ノードに持つ他のノードペアに係る評価値のうちｂ番目に小さい評価値である基準評価値、及び前記基準評価値に係るノードペアの距離に基づく計算により、所定の収束条件を満たすまで更新する更新部と、
　更新された前記評価値に基づいて、各主ノードについてエッジで接続するｂ個の従ノードを決定するエッジ決定部と、
　を備えるグラフ計算装置。
　前記ノードペア集合に含まれないノードペアである除外ペアについて、前記除外ペアの距離、前記除外ペアと主ノードを同じくする他のノードペアに係る評価値のうちｂ＋１番目に小さい評価値、及び前記除外ペアの従ノードを主ノードに持つ他のノードペアに係る評価値のうちｂ＋１番目に小さい評価値が所定の関係を有する場合に、前記除外ペアを前記ノードペア集合に追加する追加部
　を備える請求項１に記載のグラフ計算装置。
　前記エッジで接続するｂ個のノードペアに係る主ノードのベクトルと従ノードのベクトルの内積を要素とするｂ次元の内積ベクトルと、前記ｂ個の従ノードのベクトルに係るｂ×ｂのグラム行列の逆行列とに基づいて前記エッジのウエイトを算出するウエイト算出部
　を備え、
　前記ウエイト算出部は、前記内積ベクトル及び前記グラム行列の逆行列を、前記主ノードと共通の従ノードを有する他のノードである近傍ノードを用いた繰り返し計算により求める
　請求項１または請求項２に記載のグラフ計算装置。
　計算対象のＮ個のノードそれぞれについて、当該ノードを主ノードとし、前記Ｎ個のノードのうち前記主ノードからの距離に基づいて選択される（ｂ＋１）個以上（Ｎ－１）個未満のノードを従ノードとする複数のノードペアからなるノードペア集合を生成するステップと、
　前記ノードペア集合に係る各ノードペアについて、前記ノードペアの関係を評価する評価値を、前記ノードペアの距離、前記ノードペアの従ノードを主ノードに持つ他のノードペアに係る評価値のうちｂ番目に小さい評価値である基準評価値、及び前記基準評価値に係るノードペアの距離に基づく計算により、所定の収束条件を満たすまで更新するステップと、
　更新された前記評価値に基づいて、各主ノードについてエッジで接続するｂ個の従ノードを決定するステップと、
　を備えるグラフ計算方法。
　コンピュータに、
　計算対象のＮ個のノードそれぞれについて、当該ノードを主ノードとし、前記Ｎ個のノードのうち前記主ノードからの距離に基づいて選択される（ｂ＋１）個以上（Ｎ－１）個未満のノードを従ノードとする複数のノードペアからなるノードペア集合を生成するステップと、
　前記ノードペア集合に係る各ノードペアについて、前記ノードペアの関係を評価する評価値を、前記ノードペアの距離、前記ノードペアの従ノードを主ノードに持つ他のノードペアに係る評価値のうちｂ番目に小さい評価値である基準評価値、及び前記基準評価値に係るノードペアの距離に基づく計算により、所定の収束条件を満たすまで更新するステップと、
　更新された前記評価値に基づいて、各主ノードについてエッジで接続するｂ個の従ノードを決定するステップと、
　を実行させるためのグラフ計算プログラム。