DEIM Forum 2013 E9-1

精度保証付き Personalized PageRank の高速化
藤原 靖宏†

中辻

真††

山室

健†

塩川 浩昭†

鬼塚真†† NTT ソフトウェアイノベーションセンタ 〒 180–8585 東京都武蔵野市緑町 3-9-11†† NTT サービスエボリューション研究所 〒 230–0847 神奈川県横須賀市光の丘 1-1
あらましPersonalize PageRank (PPR) はグラフマイニングにおいて有用な関連度の計算方法として知られている．

本論文では PPR に対して（１）特定のノードの関連度と，（２）上位 K 個のノードと（３）関連度が閾値より大きいのノードを高速かつ正確に計算することを対象とする．オリジナルの PPR の手法は問い合わせ分布に対する関連度を全てのノードに対して収束するまで繰り返し計算によって求める．この手法の問題点として特定のノードや関連度の高いノードのみの関連度を計算できないということがあげられる．提案手法は特定のノードの関連度を疎行列を用いて計算し，関連度の高いノードを関連度の下限値と上限値を用いて求める．実験により提案手法は既存手法より大幅に高速化を達成できることを確認した．
キーワード

Personalized PageRank，関連度計算，高速化，アルゴリズム，精度保証

1. は じ め に
グラフはデータをノードとエッジで表現するデータ構造であGiven: 解ノードの個数 K と問い合わせ分布 d．
Find: 問い合わせ分布 d に対して高い関連度を有する K 個のノードの集合．

り様々な分野で用いられている．近年グラフのノードの関連度を

この問題は推薦システムに応用可能である．Liu らは PPR に用いたパーソナライズドサーチに関心が高まっているため，ノーより上位 K 個のノードを計算することにより，既存手法であドの関連度として様々な手法が提案されている [1, 2, 3, 4, 5, 6]．

る ItemRank [4] や L+ [10] より高い精度で推薦が可能である

パーソナライズドサーチにおいてノードの関連度はユーザごとことを示した [11]．

に異なった問い合わせ分布に応じて計算される．Personalized［問題 3］（閾値検索）

PageRank (PPR) は問い合わせ分布に応じてノードの関連度

Given: 閾値 ϵ 問い合わせ分布 d．

を計算する手法として最も注目を集めているもののひとつであ

Find: 問い合わせ分布 d に対して ϵ より大きいの関連度を持る [7, 8]．PPR は今までグラフ理論でよく用いられてきたノーつノードの集合．

ド間の最短距離などと異なり，グラフの構造的な特徴に基づいLiben-Nowell らは閾値より大きい関連度を持つノードを計算て関連度が計算できる．

しグラフのリンク予測を高い精度で行えることを示した [12]．

PPR は概念的に以下のように説明できる．問い合わせ分布1. 2 提案手法の特徴に基づいた各ノードの存在確率からランダムウォークを開始し，提案手法は我々の知る限り上記３つの問題を包括的に解の精隣接するノードにエッジの重みに比例した確率でランダムに移度を犠牲にすることなく扱える初めてのものである．この特徴

動する．さらにノードに到達するたびに一定の確率で問い合わは過去の既存研究 [7, 8, 13, 14, 15] が持たないものである．提せ分布に基づいた確率でノードに戻る．この操作を再帰的に繰

案手法は疎行列を用いて高速に特定のノードの関連値を求める．

り返した結果として各ノードにおける定常状態確率が得られる

また提案手法は関連度の高いノードを関連度の下限値と上限値

が，PPR はこの得られた定常状態確率を関連度とする方法でを用いて高速に検索する．提案手法の利点として以下のものがある．PPR は様々な分野のアプリケーションに応用されていあげられる．

る関連度であるが計算量が高いという問題がある．

1. 1 問 題 設 定
本論文では以下の３つの問題を対象とする．提案手法により多くのアプリケーションがその精度を失うことなく高速に処理することが可能となる．
［問題 1］（特定のノードに対する関連度の計算）

Given: 問い合わせノード x と問い合わせ分布 d．

高速: 提案手法は上記の手法を用いることによって既存手法より大幅に高速である．
正確: 提案手法の出力結果は正確であることが理論的に保証されている．
パラメータフリー: 既存手法と異なり，提案手法自体が必要とするパラメータはない．

PPR はその有用性により多くのアプリケーションに応用可Find: 問い合わせ分布 d に対するノード x の関連度．

能であるが計算コストが問題であった．提案手法によって PPR

この問題はグラフ解析に応用が可能である．例えば White らを高速に処理することが可能となり，今後多くの PPR を用いは研究者の共著関係のグラフに PPR を用いることにより，同たアプリケーションが開発されることが可能になる．

じ大学に所属していても異なる研究分野である研究者は低い関論文の構成は以下の通りである．2. 章で背景知識を述べる．

連度を持つことを明らかにした [9]．

3. 章で提案手法の詳細を述べる．4. 章で既存手法と提案手法の［問題 2］（上位 K 個検索）比較を行う．5. 章で提案手法におけるケーススタディをあげる．

表1

主な記号の定義

計算する．提案手法は事前計算の段階で隣接行列におけるノー

Table 1 Deﬁnition of main symbols.
記号

定義

n
m
x
K
ϵ
c
su
s
d
A
P
Q
R
|R|

グラフのノード数
グラフのエッジ数
特定のノードの関連度計算における問い合わせノード
上位 K 個検索における解ノードの数
閾値検索における閾値
ランダムウォークにおける再スタートの確率
ノード u の関連度
n×1 の関連度ベクトル
n×1 の問い合わせ分布ベクトル
n×n グラフの隣接行列
n×n ノードの置換行列
n×n QR 分解における直交行列
n×n QR 分解における上三角行列
行列 R における非零要素数

ドを並び替え，並び替えを行った行列の QR 分解 [16] を計算する．そして QR 分解の結果得られた行列の逆行列を計算す
る．適切にノードを並び替えることによってこの逆行列は疎になる．ノードの並び替えは NP 完全問題を解くことで得られる．QR 分解は行列の近似手法ではないため，正確に関連度を計算することができる．
上位 K 個のノードの検索

問い合わせ分布に対する上位 K個のノードを計算するために，オリジナルの手法では全てのノードに対して関連度を計算する必要がある．しかしこの手法
は計算コストが高いという問題がある．そのため提案手法では解になり得ないノードを枝刈りし，解になる可能性があるノー

6. 章で関連研究を述べる．7. 章で本論文をまとめる．

ドに対してのみ関連度を計算する．
解になる可能性があるノードを求めるために，提案手法は各

2. 前 準 備

ノードに対して関連度を O(1) の計算量で推定する．この手法
まずは本論文で用いる記号を定義し，必要となる背景知識をの利点として解ノードを正確に求められることがあげられる．

説明する．表 1 に主な記号とその定義を示す．PPR はグラフ

正確な解ノードは関連度の上限値を用いることによって得られのノードの関連度を計算するひとつの方法である．PPR ではる．このため提案手法は関連度の低いノードを推定値によって問い合わせ分布から各ノードの存在確率を決定し，決定された

高速に枝刈りすることができる．

ノードからランダムウォークを行う．そしてランダムウォーク

閾値より大きいのノードの検索 関連度の上限値を用いるこ

でノードに到達するたびに一定の確率 c で問い合わせ分布に基とによって解になり得ないノードを枝刈りできるが，解の正確

づき各ノードに戻る．グラフのノード数を n としエッジ数を性を保証するために解になる可能性があるすべてのノードに対m とする．s を n×1 行列とし，su 要素をノード u にランダムして正確な関連度を計算する必要がある．高速に閾値 ϵ より大ウォークする確率を表すとする．また d を問い合わせ分布で決きいのノードの検索するため，提案手法はこの正確な関連度の定される n×1 行列とし，問い合わせ分布の値が 0 でないノー

計算のコストを減らすために関連度の下限値を計算する．

ドをシードノードとする．また A を列が正規化されたグラフ

もしあるノードの関連度の下限値が ϵ より大きければ，そのの隣接行列とする（すなわち Au,v 要素はノード v からノード

ノードの正確な関連度は ϵ より大きくなる．すなわちこのノーu へランダムウォークする確率を表す）．定常状態における各ドは解ノードとなる．提案手法はこの性質を用いて正確な関連

ノードにおける存在確率は以下の式を再帰的に収束するまで繰

度の計算のコストを減らす．結果として提案手法は閾値より大り返すことで計算することができる．

きなノードを高速に求めることができる．

s = (1 − c)As + cd

(1)

3. 2 特定のノードの関連度計算
この章では問い合わせ分布に対し特定のノードの関連度を高PPR は定常状態における確率を関連度とする．すなわち s 行

速に計算する手法について述べる．まず 3. 2. 1 章で PPR の定

列の su 要素はノード u の分布 d に対する関連度となる．

義から QR 分解を用いることによって関連度を繰り返し計算を定義から繰り返し回数を t としたとき PPR の計算コストは行わずに求められることを示す．そして 3. 2. 2 章で QR 分解のO((n + m)t) となる．そのためグラフが大規模である場合 PPR

結果得られる行列を疎にすることは NP 完全問題であることの値を計算するのは非常に計算コストが必要となる問題がある．

を述べ，それに対する解決方法を示す．

3. 提 案 手 法

3. 2. 1 繰り返しなしの関連度計算
式 (1) を変形して繰り返し計算なしで関連度が求められるこ

この章では提案手法の詳細を述べる．提案手法は正確かつ高とを示す．提案手法ではグラフのノードを並び替える．すなわ

速に特定ノードの関連度を計算し，また関連度の高いノードを

ち行列 P をノードの置換行列 [16] とすると，グラフの隣接行求めることができる特徴がある．

列 A は A′ = PAPT を計算し変換される．ここで行列 PT は

3. 1 概要

行列 P の転置行列である．またここで n × n の大きさの置換

特定のノードの関連度計算 問い合わせ分布に対する関連度

行列 P は直交行列で，各行と列においてただ一つの値が 1 では繰り返し計算によって定常状態における確率を計算することある要素と，それ以外は値が 0 である要素を持つ．置換行列にで求めることができる．この方法はグラフのすべてのノードのおいて Pij = 1 は j 番目の行が i 番目の行に置換されることを

関連度を再帰的に計算するため，特定のノードの関連度のみの表す．置換行列 P を求める手法は 3. 2. 1 章で述べる．

計算を行うことができない．再帰的に計算を行うことを避ける

ここで I = PT IP，P−1 = PT ，A = PT A′ P であるため

ために，提案手法では行列計算を用いて再帰的でない形で特定（P−1 は P の逆行列），式 (1) は P を用いて以下のように書のノードの関連度を計算する．この行列計算は式 (1) から得る

き換えられる．

ことができる．しかしこの手法は行列が密であれば計算コストを要するという問題がある．
提案手法は疎な行列を求めて特定のノードの関連度を高速に

s = c{I−(1−c)A}−1d = c{PT IP−(1−c)PTA′ P}−1d= cPT{I−(1−c)A′ }−1 Pd(2)特 定 の ノ ー ド の 関 連 度 を 計 算 す る た め ，提 案 手 法 は 行 列′

I − (1 − c)A の QR 分 解 を 計 算 す る ．す な わ ち QR =′におけるグラフを隣接行列 A に置き換える．またノードの削除の順番を置換行列 P に置き換え，弦グラフを行列 R−1 と Q

I − (1 − c)A となる．ここで行列 Q は直交行列であり，行列に置き換える．このようにすれば Minimum fill-in 問題におけ

R は上三角行列である．提案手法は QT と R−1 を用いて特定る最小のエッジの追加数に対する解がある時かつそのときに限のノード x の関連度を以下のように計算する．

り疎行列問題における非零要素の増加を最小にする解も存在す［定義 1］（関連度の計算） n × n の行列を F = cPT R−1 とし，ることがわかる．そのため疎行列問題は NP 問題になる． □n × 1 のベクトルを g = Q Pd とし，1 × n のベクトルを fi を

本論文では疎行列問題に対する手法を示す．手法の詳細を示

行列 F における i 番目の行ベクトルとする．ノード x の関連度

すためにまず行列 Q と R−1 がどのように計算されるかを示す．

sx を以下のように計算する．

ベクトル qi を行列 Q の i 番目の列ベクトルとし，ベクトル wiT

sx = fx · g

を行列 W = I − (1 − c)A′ の i 番目の列ベクトルとする．グラ(3)

ム・シュミットの正規直交化法 [16] を用いることによって，行

定義 1 に対して以下の定理が成り立つ．
［定理 1］（関連度の計算） 式 (3) による関連度は式 (1) による関連度と等しくなる．
証明

QT = Q−1 であるため [16]，式 (2) から以下の式が成り立つ．

列 Q と R は以下のように計算できる．
{
wi(i = 1)′′′(4)
qi = qi /||qi ||, qi =∑wi − i−1
(w
·
q
)q
(i =
| 1)
i
j
j
j=1

(i > j)

 0

Rij =
−1

T

s = cP (QR)

T

−1

Pd = cP R

Q Pd = F · gT||q′i ||(i = j)wj · qi(i < j)(5)

そのため n × 1 の列ベクトル s の x 番目の要素 sx （ノード x

ここで ||q′i || をベクトル q′i のノルム [16] とする．行列 R−1 は□後退代入 [18] を用いることによって以下のように計算できる．
(i > j)

 0
−1
(6)
Rij
=
1/Rij
(i = j)∑j−1−1/Rii k=i+1 Rik Rkj (i < j)の関連度）はベクトル fx とベクトル g の内積となる．

式 (3) からノード x の関連度はベクトル fx と g を用いれば繰り返し計算なしに求めることができる．またベクトル d はユーザから与えられるため，関連度を計算するためまずベクトル

g = QT Pd を計算する．
関連度の計算量を述べるために以下の補助定理を示す．

式 (4) と (5) から行列 Q と R の列ベクトルは左から右の列［定理 2］（関連度の計算量） |fx | と |Q| をそれぞれ fx と Q にへそれぞれ直交な列ベクトルと行列 W との内積を計算するこおける非零要素の数とすると，特定のノードの関連度を計算すとで求められることがわかる．式 (6) から行列 R−1 の列ベク

るために必要な計算量は O(|fx | + |Q|) となる．

トルは右から左の列へ，また各列ベクトルの要素は下から上の証明

関連度を計算するためまずベクトル g を分布 d から計

要素の順番で求められる．また式 (4) と (5) と (6) から行列 Q

算する．ベクトル Pd は分布 d のノードを置換すれば得られると R は行列 W から計算できることがわかる．また W は行列ため，ベクトル g = QT Pd を計算するために必要な計算量はA′ (= PAPT ) から求められる．

O(|Q|) となる．ベクトル fx と g の内積を計算するための計算

行列 Q の非零要素を減らすために，提案手法では以下の知コストは O(|fx | + |Q|) であるため，関連度を計算するため必見を用いる．
「もし列ベクトル wi が行列 wi において全ての左要な計算量は O(|fx | + |Q|) となる．

□の要素がすべて零である一つの非零要素をもつならば，その列

T

−1

T

であり g = Q Pd であるため，定理 2ベクトル wi は行列 Q のすべての左の列ベクトルに対して直から関連度を高速に計算するためには行列 R−1 と Q の非零要交となる．
」．これは列ベクトル wi が線形独立 [16] になるた

素を減らす必要があることがわかる．なおここで明らかに行列めである．よって行列 Q は疎なデータ構造となる．また行列P の非零要素の数は n となる [16]．とくにベクトル g はユーザ

R−1 を疎にするため以下の２つの知見を用いる．
「もし行列 Rの問い合わせから計算されるため，行列 Q の非零要素の数はの右と下の要素が疎であれば，行列 R−1 の右と下の要素も疎ここで F = cP R

−1

少なくなるようにしなければならない．次の章で行列 RとQ を疎にする手法を述べる．
3. 2. 2 疎行列問題になる．
」．
「もし行列 W の右と下の要素が疎であれば，行列R の右と下の要素も疎になる．
」．
Algorithm 1 に疎な行列を得るためのノードの置換方法に

定義 1 に示すとおり，提案手法は事前に計算した行列 R−1

ついてのアルゴリズムを示す．このアルゴリズムは行列 Q とと Q から関連度を計算する．しかしもしこれらの行列が密でR−1 に対する知見に基づいている．このアルゴリズムにおいあれば関連度を求めるには高い計算量が必要となってしまう．

て V はグラフにおけるノードの集合とし， P は置換された

これらの行列を疎にするために提案手法ではノードの並び替え

ノードの集合とする．deg(u) はノード u に入るエッジの数とを行う．しかし疎な行列を得るためにノードの並び替えを行う

し，e(u) はノード u から集合 P へ出ていないエッジの数とすのは NP 完全問題である．

る．まずアルゴリズムは置換行列を零行列に初期化し，集合 P［定理 3］（疎行列問題） 行列 R−1 と Q の非零要素をノード

を空集合とする（１〜２行目）．置換されていないノードの集の並び替えによって最小化するのは NP 完全問題である．

合 V\P からエッジの数が min{e(u)|u ∈ V\P} になるノード

証明

この定理を Minimum fill-in 問題 [17] のリダクション

集合 U を計算する（４行目）．この処理は行列 Q に対する知から証明する．まず Minimum fill-in 問題のインスタンスを疎見に基づいている．もし行列 A′ の右と下の要素が疎であれば行列問題のインスタンスに置き換える．Minimum fill-in 問題

行列 R−1 も疎になるため，行列 A′ の左と上の要素を密にす

クがその過程でシードノードに至るかの２つの場合である．
（１）

Algorithm 1 ノードの並び替えの場合ランダムウォークはあるシードノード u に確率 cdu でInput: A, グラフの隣接行列; n, ノード数
Output: P, 置換行列
1: P = 0;
2: P = ∅;3: for i = 1 to n do
4:
U = argmin(e(u)|u ∈ V\P);5:
v = argmax(deg(u)|u ∈ U );6:
Piv = 1;
7:
ノード v を P に追加;
8: end for
9: return P;

ジャンプする．そのためシードノード u の定常状態における確率は cdu より小さくなることはない．よって H(0) に含まれるノードに対して成り立つ．
次に H(i − 1) 含まれるノードに対して sv <= sv が成り立つ場<
合，H(i) 含まれるノードに対して sv = sv が成り立つことを示す．H(i − 1) 含まれるノードはグラフ全体のノードの部分集合であり c, du >= 0 であるため，式 (1) において以下が成り立つ．

るためにエッジの数が最大になるノードを選択する（５行目）．
そして選択されたノードから行列 P の要素を決定し，置換のsu = (1 − c)∑∑

Auv sv + cdu >= (1 − c)v∈VAuv sv = suv∈H(i−1)

順番を決定する（６〜７行目）．これらの処理を全てのノード−1

4. 章で示すとおり，この手法によって行列 Rと特に行列Q に対して非零要素の数を減らすことができる．F = cPT R−1
Tであり g = Q Pd であるため，提案手法は特定のノードの関連度を高速に計算可能となる．さらに非零要素の数が少なくなることにより，事前に行列 Q と R−1 を高速に計算可能となる．

よって全てのノードに対して定義 2 から下限値を計算できる．

3. 3. 2 関連度の上限値
この章では不要な計算を枝刈りするために関連度の上限値を計算する手法を述べる．上限値を計算するために 3. 3. 1 で導入した下限値を用いる．提案手法は上位 K 個のノードを検索するために一つ一つノードをたどり上限値と正確な関連度を計算

3. 3 上位 K 個のノードの検索
この章では上位 K 個のノードを検索するために，関連度の上限値を用いる手法を示す．上限値は上位 K 個のノードを検索する過程で計算する．3. 3. 1 章で上限値を計算するために必要となる関連度の下限値の計算方法を述べ，上限値の計算方法を

3. 3. 2 章で述べる．そして 3. 3. 3 章で検索アルゴリズムを示す．
3. 3. 1 関連度の下限値
関連度の上限値を計算するために提案手法はまず関連度の下限値を推定する．関連度の下限値は与えられた分布の値が 0 でないシードノードからの最小ホップ数を用いて計算する．シー
ドノードからの最小ホップ数は幅優先探索を用いて計算する．

hu をノード u に対するシードノードからの最小ホップ数とし，H(i) をシードノードからの最小ホップ数が i(i = 0, 1, . . .) となるノードの集合とする．すなわち H(i) は i 番目のレイヤを構成し，H(0) はシードノードから構成される一番上のレイヤのノード集合である．関連度の下限値は以下のように計算する．
［定義 2］（下限値） ノード u の下限値 su は以下の式に基づき計算される．
{cdu
su =∑(1 − c) v∈H(hu −1) Auv sv□

よって成り立つ．

が並び替えるまで繰り返す．

(u ∈ H(0))(u ∈/ H(0))(7)

する．ui を i(i = 1, 2, . . .) 番目にたどるノードとしたときに関連度の上限値は以下のように計算する．
［定義 3］（上限値） ノード ui の上限値 sui は以下のように定
義される．
{

sui =

1−

∑n

j=2 suj(i = 1)sui − sui−1 + sui−1

(i =
| 1)

(8)

上限値を計算する計算コストはもし i = 1 であれば O(n) になる．またもしそうでなければ上限値はすでに計算した sui と

sui−1 と sui−1 を用いて逐次的に計算できるため，その計算コストは O(1) になる．
提案手法はノードは関連度の下限値の降順でたどる．これには２つの理由がある．はじめの理由は下限値が大きいほど正確な下限値も大きくなることが期待されるからである．そのため効率的に上位 K 個のノードを検索できる．２つめの理由は以下に示す補助定理により検索の途中で処理を停止し，より高速に検索を行うためである．
［補助定理 2］（上限値） ノードを関連度の下限値の降順にたどった場合，sui <= sui <= sui−1 が成り立つ．
まず su <= su になることを示す．式 (8) と補助定理 1証明から以下の式が成り立つ．
この式は（１）シードノードに対しては再スタートの確率と問い合わせ分布の値から下限値を計算し，（２）シードノードで

su = sui −sui−1 +sui−1 = sui −sui−1 +sui−1 −sui−2 + . . . +su1

なければ一つ上のレイヤのノードの下限値と遷移確率から計算することがわかる．全てのノードの下限値は O(n + m) の計算

=

i
∑

suj −

j=2

コストで求めることができる．これは下限値がシードノードをルートとする幅優先探索で計算でき，幅優先探索は O(n + m)>=1−

この下限値の性質について以下の補助定理を示す．
［補助定理 1］（下限値） H(i) に含まれる全てのノードに対しsが成り立つ．
て su <u=suj + 1 −

j=1
i−1
∑
j=1

で計算できるからである．

i−1
∑

suj −

n
∑

suj = 1 −j=2n∑suj = 1 −

i−1
∑

suj −

j=1
n
∑

n
∑

suj

j=i+1

suj = sui

j =i
|

j=i+1

つぎに sui <= sui−1 が成り立つことを示す．式 (8) から
sui−1 − sui = sui−1 − su >− su >= su= 0 となる．よって
i

i−1

i

成り立つ．

□

証明には数学的帰納法を用いる．まず H(0) に含まれる

補助定理 2 からもしあるノードの上限値が解の候補ノードの K

ノードに対して成り立つことを示す．PPR においてランダム

番目の下限値より小さくなれば，他のまだたどっていないノーウォークがシードノードに至るのは（１）ランダムウォークがドの関連度も K 番目の下限値より小さくなることがわかる．そ確率 c でシードノードにジャンプするか，（２）ランダムウォーのためこの場合，提案手法はその処理を打ち切る．

証明

Algorithm 2 上位 K 個のノードの検索

Algorithm 3 閾値より大きいノードの検索

Input: d，問い合わせ分布; K ，解ノードの個数
Output: Va ，解ノードの集合
1: θ = 0;
2: Ve = ∅;3: Va = ∅;4: K 個のダミーノードを Va に追加;
5: すべてのノードに対して関連度の下限値を計算;
6: while Ve =| V do7:u = argmax(sv |v ∈ V\Ve );8:
ノード u の関連度の上限値を計算;
9:
if su < θ then10:
return Va ;
11:
else
12:
問い合わせ分布 d に対するノード u の正確な関連度を計算;
13:
if su > θ then14:v = argmin(sw |w ∈ Va );15:
ノード v を集合 Va から削除;
16:
ノード u を集合 Va に追加;
17:
θ = min(sw |w ∈ Va );18:
end if19:
end if20:
ノード u を Ve に追加;
21: end while
22: return Va ;Input: d，問い合わせ分布; ϵ，検索の閾値
Output: Va ，解ノードの集合
1: Ve = ∅;2: Va = ∅;3: すべてのノードに対して関連度の下限値を計算;
4: while Ve =| V do5:u = argmax(sv |v ∈ V\Ve );6:
ノード u の関連度の上限値を計算;
7:
if su <
= ϵ then8:
return Va ;9:else if su > ϵ then10:
ノード u を集合 Va に追加;
11:
else
12:
問い合わせ分布 d に対するノード u の正確な関連度を計算;
13:
if su > ϵ then14:
ノード u を集合 Va に追加;
15:
end if16:
end if17:
ノード u を Ve に追加;
18: end while
19: return Va ;∑算コストを要する． ni=1 |fi | + |Q| = |F| + |Q| であるため，定理 2 からすべてのノードの関連度を計算するには O(|F| + |Q|)

3. 3. 3 上位 K 個のノードを求めるアルゴリズムの計算コストを要する．そのため上位 K 個のノードを検索を

Algorithm 2 に上位 K 個のノードを検索をするアルゴリズ

する計算コストは O(n log n + m + |F| + |Q|) である．

ムを示す．Algorithm 2 において θ は解候補のノード集合にお□3. 4 閾値より大きいノードの検索

ける K 番目の関連度であり，Va は解候補のノード集合であり，この章では閾値 ϵ より大きい関連度を持つノードを検索する

Ve は探索済みのノード集合とする．アルゴリズムは関連度が手法について述べる．この章では紙幅の関係から証明を省く．

0 である K 個のダミーノードを Va に追加する（４行目）．そ閾値 ϵ より大きい関連度を持つノードを検索するのに関連度のしてすべてのノードに対して関連度の下限値を計算する（５行

下限値を用いる．具体的にはもしあるノードの関連度の下限値

目）．そして最も関連度の下限値が最も高いノードを探索していが閾値 ϵ より大きい場合，明らかにそのノードの正確な関連度ないノードの集合 V\Ve から選択し，選択されたノードの関連は閾値 ϵ より大きくなる．そのためそのようなノードに対して度の上限値を計算する（７〜８行目）．補助定理 2 から選択さ

正確な関連度の計算を省略する．閾値 ϵ より大きい関連度を持れたノードの関連度の上限値が θ より小さい場合，そのノードつノードを検索するために，上位 K 個のノードを計算する場

および探索されていないノードは解ノードになり得ない．その合と同様に下限値の降順にノードを探索する．そして関連度のためその場合アルゴリズムは処理を停止する（９〜１０行目）．

上限値を用いて不要な関連度の計算を省略する．しかしもし正もしそうでなければ選択されたノードは解になり得る．そのた

確な関連度の計算を行わない場合，定義 3 を用いて上限値を計

めそのノードの正確な関連度を計算する（１２行目）．もし計算することができない．そのためもし正確な関連度の計算を省算した関連度が θ より大きければ，解候補の集合 Va と K 番

略した場合，以下のように関連度の上限値を計算する．

目の関連度 θ を更新する（１３〜１８行目）．

［定義 4］（正確な関連度を計算しない場合の上限値） もしノーこのアルゴリズムが正確に解ノードを検索できることを示すために以下の定理を示す．

限値 sui を以下のように計算する．

［定理 4］（正確な上位 K 個のノードの検索） Algorithm 2 により上位 K 個のノードを正確に検索できる．
証明

ド ui−1 の正確な関連度の計算を省略した場合，ノード ui の上

θK を解ノードにおける K 番目の関連度とした場合，解ノードの関連度の上限値は θK より小さくなることはない（補sui = sui − sui−1 + sui−1(9)

この定義について以下の補助定理を示す．

助定理 2）．もし関連度の上限値が θ より小さければ，アルゴ［補助定理 3］
（Upper bound w/o exact relevances） 定 義 4< sui =< sui−1に対して i 番目に探索したノード ui に対して sui =リズム 2 はノードを枝刈りする．解候補のノード集合におけるが成り立つ．また上限値を計算する計算コストは i =| 1 であれ

K 番目の関連度 θ は θK より大きくなることはないため，解ば O(1) になる．

ノードがアルゴリズム 2 により枝刈りされることはない．またAlgorithm 3 に関連度が閾値 ϵ より大きくなるノードを検索

もしあるノードが θ より小さければ，そのノードはアルゴリズ

するアルゴリズムを示す．このアルゴリズムはすべてのノード□に対して関連度の下限値を計算する（３行目）．そして最も関ム 2 により枝刈りされる．よって成り立つ．
提案手法の計算コストについて以下の定理を示す．

連度の下限値が最も高いノードを探索していないノードの集［定理 5］（上位 K 個のノードの検索の計算コスト） Algorithm 2 合から選択し，選択されたノードの関連度の上限値を計算する
は O(n log n + m + |F| + |Q|) の計算コストを要する．

（５〜６行目）．もしその上限値が ϵ より大きくなければ，そ証明

Algorithm 2 はまずすべてのノードに対して関連度の下のノードおよび探索されていないノードは解ノードになり得な限値を計算するが，下限値は幅優先探索で O(n + m) の計算コい．そのためその場合アルゴリズムは処理を停止する（７〜８ストで計算できる．そして関連度の下限値が最も大きいノード

行目）もし選択されたノードの下限値が ϵ より大きければ，そを選択しその関連度の上限値を計算するが，もしどのノードものノードは解ノードとなる．そのためそのノードを解ノードの枝刈りされなければこれらはそれぞれ O(n log n) と O(n) の計集合に追加する（９〜１０行目）．そうでなければ正確な関連

10-2
-4

10

-1

10

6

10

Proposed(10)Proposed(50)Proposed(100)
BPA(10)
BPA(50)
BPA(100)

1

105

10-2

0
250

-4

Social

図1

Routing

SocialWeb

Routing

Web

関連度の計算時間．

identiﬁcation.

time.

非零要素の割合．

|Q|/n2|R−1 |/n2

Social3.40 × 10−6
9.25 × 10−3

3

10

102
1

750

1000

Number of hub nodes

010SocialRoutingWeb

図 4 事前計算時間．

Fig. 4 Precomputation time.

of hub nodes.

の割合を示す．

Table 2 Ratio of non-zero elements.
Data set
Routing
4.88 × 10−6
2.40 × 10−3

4

10

10
500

図 2 上位 K 個のノードの計算時間．図 3 ハブノードの数と適合率．
Fig. 2 Eﬃciency of top-k nodes Fig. 3 Accuracy versus number

Fig. 1 Relevance computation

表2

0.4

Proposed
BPA

10

-8

10

0.6

0.2

-3

10

10-6

ProposedRandom0.8

Wall clock time [s]Wall clock time [s]Wall clock time [s]

100

0

10

Proposed
Proposed_all
Iterative

Precision

2

10

Web
2.76 × 10−4
1.93 × 10−2

図 1 から提案手法はオリジナルの手法より７桁以上高速であることがわかる．また全てのノードの関連度を計算したとしても 70 倍以上高速である．定理 2 に示すとおり，提案手法はO(|fx | + |Q|) の計算コストで計算できる．ここで F = cPT R−1

度を計算し，選択されたノードが解ノードかの確認を行う（１であり g = QT Pd であるため，提案手法において非零要素の１〜１６行目）．

数が関連度の計算時間に影響する．表 2 に示すとおり，提案手

Algorithm 3 の性質について以下の定理を示す．

法は行列 R−1 ，とくに行列 Q を疎行列にすることができる．

［定理 6］（閾値より大きいノードの検索） Algorithm 3 は関すなわち |Q| ≪ |R−1 | ≪ n2 となる．そのため提案手法は与え連度が閾値 ϵ より大きくなるノードを正確に O(n log n + m +られた問い合わせ分布 d に対してベクトル g を高速に計算で|F| + |Q|) の計算コストで行う．

きるため，関連度を高速に計算できる．一方 2. 章で述べたとお4. 評 価 実 験り，オリジナルの手法は全てのノードの関連度を収束するまで

提案手法の有効性を示すために評価実験を行った．実験では提案手法はオリジナルの手法より関連度を高速に計算できる．

繰り返し計算を行うため，計算コストが高くなる．結果として以下のデータセットを用いた．
（注 1）

Social

（注 2）: Slashdot.org

4. 2 上位 K 個のノードの計算時間におけるユーザ間の関係グラフ

である．ノード数は 82, 168 でエッジ数は 948, 464 である．

提案手法は上位 K 個のノードを高速に検索できることを示すために Basic Push Algorithm [15] と比較を行った．この手

Routing（注 3）: Oregon Route Views Project（注 4） により得ら

法は上位 K 個のノードを近似的に検索する手法である．ハブ

れたネットワークのルーティングのデータである．ノード数はノードから全てのノードに対する関連度を事前に計算して，こ22, 963 でエッジ数は 48, 436 である．

の手法は関連度の上限値を計算する．関連度の上限値を用いる（注 5）: ノートルダム大学のサイトから得られたグラフでため，この手法において解ノードの再現率（注 6）は 1 になる．すな

ある．ノード数は 325, 729 でエッジ数は 1, 497, 135 である．

わちオリジナルの手法で得られる解ノードの全てを Basic PushWeb

PPR における再スタートの確率は関連研究 [2] と同様にAlgorithm によって得ることができる．他の近似手法 [13, 14]c = 0.9 とした．実験ではランダムに 10 個のシードノードを選にはこのような性質はない．しかしこれは同時に Basic Push

択し，それらの問い合わせ分布の大きさは等しい値とした．実験

Algorithm が解ではないノードも出力することも示している．

は CPU が Intel Xeon 2.26GHz，メモリが 144 GB の Linux

実験においてハブノードはノードの次数が高いものとした．

サーバで行った．すべてのアルゴリズムは GCC で実装した．

図 2 にハブノードの数を 1, 000 個にした場合の結果を示す．

4. 1 関連度の計算時間

この図において提案手法と Basic Push Algorithm の結果をそ提案手法で特定のノードの関連度を計算するために必要な処れぞれ「Proposed(K)」，「BPA(K)」とした．ここで K は解

理時間についてオリジナルの手法と比較実験を行った．オリジ

ノードの個数である．また図 3 に Social データに対してハブ

ナルの手法は特定のノードの関連度を計算できないため，提案ノードの数を変えたときの上位 10 個の解に対する適合率を示手法で全てのノードを関連度を計算するために必要な処理時間

す．適合率は各手法で得られた解ノードのうちどの程度オリジもあわせて示す．

ナルの手法による解ノードと一致しているかの割合である．

図 1 に結果を示す．ここで提案手法で特定のノードの関連度

図 2 にみられるように，提案手法は Basic Push Algorithmを計算した場合の結果を「Proposed」，全てのノードの関連度

より 50 倍以上高速に上位 K 個のノードを検索できる．提案手を計算した場合の結果を「Proposed all」に示す．また「Itera-

法は全てのノードに対して上限値は O(n) の計算コストで求めtive」はオリジナルの繰り返し計算により関連度を計算した場ることができ，また上限値を用いて正確な関連度の計算を枝刈−1合の結果である．また表 2 に行列 Q と R

における非零要素

りする．そのため提案手法は高速な検索が可能である．
図 3 から提案手法はオリジナルの手法と同じ結果を出力す（注 1）：http://snap.stanford.edu/data/soc-Slashdot0902.html

ることがわかる．また図 3 から Basic Push Algorithm はハブ（注 2）：http://slashdot.org/（注 3）：http://www-personal.umich.edu/˜mejn/netdata/as-22july06.zip（注 4）：http://routeviews.org/（注 6）： 再現率はオリジナルの手法で得られる解ノードうちどの程度正しく解ノードが得ら（注 5）：http://vlado.fmf.uni-lj.si/pub/networks/data/ND/NDnets.htm

れるかの割合である．

Halifax
Lebanon
ChicagoDenverColorado Springs

Philadelphia

New York
DenverLos Angeles

New York
DenverDecatur
Ponca City

Daytona Beach

Washington, D.C.

Orlando(1) Proposed

Baltimore

New York

Los Angeles
Phoenix

Knoxville

Daytona Beach

Muscle Shoals

San Diego

Colorado Springs

Washington, D.C.
Indianapolis

Lake Havasu City

Phoenix

Altoona

Colorado Springs

Washington, D.C.

San Diego

Orlando

San Diego(2) Monte Carlo(250)

Orlando(3) Monte Carlo(1,000,000)

図 5 各手法で得られた関連度の高い都市．

Fig. 5 Highly relevant cities as detected by each approach.

ノードの数が多くなるほど高い精度で上位 K 個のノードを検図 5 に結果を示す．Monte Carlo 法による手法においてランダ

索できることがわかる．これは Basic Push Algorithm が事前

ムウォークの数を 250 とした結果は図 5-(2) に，1, 000, 000 とに計算したハブノードからの関連度を用いて上限値を計算すした結果は図 5-(3) に示す．図 5 において丸はシードノードをるからである．しかし Basic Push Algorithm は正確な解ノー表し，三角は関連度の高いノードを示す．

ドを求めることができない．さらにハブノードからの関連度は

図 5-(1) に提案手法で求めたサンディエゴ，デンバー，オー

繰り返し計算によって求めなければならないため，Basic Push

ランド，ニューヨークの４都市に対して関連度の高い都市を示

Algorithm の事前計算時間はハブノードの数が多くなるほど長す．これら４つの都市はいずれもアメリカで代表的な都市であくなる．図 2 と 3 は提案手法が Basic Push Algorithm に対しり，大きなハブ空港を持つ．ハブ空港とはその他の多くの空港て計算時間と精度の両面で優れていることを示している．

へ直行便で行くことができる空港である．PPR はシードノード

4. 3 事前計算時間からのランダムウォークを用いてノードの関連度を計算する手提案手法は関連度を計算するために事前に行列 Q と R−1 を法であるため，もし大きな都市がシードノードとして選ばれた計算する．そのためノードの置換による事前計算時間の影響に場合，多くの直行便で他の都市と結ばれている都市が関連度のついて実験を行った．図 4 において「Random」はノードをラ高いノードとなることが期待される．その予想通り関連度の高ンダムな順番に並び替えたときの事前計算時間を示す．

い都市として大きなハブ空港を持つロサンゼルス，フェニック

図 4 から提案手法におけるノードの並び替えは非零要素のス，コロラドスプリングス，デイトナビーチ，フィラデルフィ数を減らすのみでなく事前計算の時間を減らすことができるこ

ア，ワシントン D.C. が得られた．すなわち大きな都市をシーとがわかる．提案手法はランダムにノードを並び替えた場合とドノードとすることでハブ空港を持つ都市を検索することがで比較して 1, 200 倍まで高速に事前計算を行えることがわかる．

きる．なお提案はオリジナルの手法と同じ結果になることが保提案手法により行列 Q と R−1を疎にすることができるため，証されている．

式 (4) と (5) と (6) における非零要素の計算を高速行うことが

図 5-(2) にみられるように，ランダムウォークの数を 250 にできる．そのため関連度の計算において行列 Q と R−1 を高速した場合，Monte Carlo 法に基づく手法はコロラドスプリングに計算することができる．

スとワシントン D.C. 以外は提案手法と異なる結果となった．

5. ケーススタディ

各手法の処理時間は Monte Carlo 法に基づく手法が 21 マイク
ロ秒であり，提案手法が 20.5 マイクロ秒であり，ほぼ同じで提案手法は関連度が ϵ より大きいノードを正確に検索するこあった．図 5-(3) に示すとおり，ランダムウォークの数が多くとができる．提案手法の有効性を示すために実データを用いた

なれば Monte Carlo 法に基づく手法はボルチモアを除いて提ケーススタディを実施した．Basic Push Algorithm を含む過案手法と同じになる．しかし Monte Carlo 法に基づく手法は去のいずれの手法も閾値より大きいノードを検索するものがな

ランダムウォークの数が多くなると計算時間を要するという問

いため，この章では Avrachenkov らによる最新の PPR の近題がある [13]．ランダムウォークの数を 1, 000, 000 にしたとき似手法 [13] と比較を行った．この手法はシードノードからランMonte Carlo 法に基づく手法の処理時間は 80.7 ミリ秒であり，

ダムウォークを複数回計算する Monte Carlo 法を用いて関連提案手法の方が 3, 900 倍以上高速であった．

度を近似する．
ケーススタディでは北アメリカにおける都市間のエアライン

6. 関 連 研 究

ネットワークのデータ （注 7） を用いた．このデータにおいてノー

これまで PPR を高速に計算するため様々な手法が提案されドに対応する都市の数は 456 であり，エッジに対応する直行便てきた．しかし過去の手法のいずれも解の正確性を保証するもの数は 71, 959 である．エッジの重みは都市間の近接度に比例のではない．すなわちこれらの手法の出力結果はオリジナルのする．ケーススタディでは閾値を 5×10−4 に設定し，サンディ手法と異なるものである．そのためこれらの手法でアプリケー

エゴ，デンバー，オーランド，ニューヨークの４都市に対してションの質を保ちながら高速化を果たすことは難しい．また過関連度の高い都市を計算した．これらの４都市における問い合

去の手法のいずれも本論文で対象とした問題設定を包括的に解

わせ分布は等しい値に設定し，前の章と同様に c = 0.9 とした．

決するものはない．

Jeh らはハブノードからの関連度から関連度を近似する手法（注 7）：http://www.psi.toronto.edu/index.php?q=aﬃnity propagationを示した [8]．この手法は問い合わせ分布に対する特定のノード

の関連度をハブノードからの関連度の線形結合として近似した．

Fogaras らは Monte Carlo 法を PPR に応用した [7]．近似の関連度を求めるために，彼らはまず fingerprint とよばれるランダムウォークを複数計算し，線形結合を用いて近似の関連度を計算した．Avrachenkov らと Bahmani らは独立に PPRに対して Monte Carlo 法を改善する手法を提案した [13, 14]．

Gupta らによって提案された Basic Push Algorithm は上位 K個のノードを近似的に求める手法である [15]．この手法は関連度の上限値を事前に計算したハブノードからの関連度を用いて計算する．しかしこの手法はオリジナルの手法では解にならないノードも解として出力してしまう．Chakrabarti らは BasicPush Algorithm を ER グラフに対して適用した [19]．
Random walk with restart (RWR) [2] は一つのノードのみ
をシードノードとする PPR の特殊例である．Sun らは RWRに対する近似手法を提案した [20]．彼らの手法は METIS [21]を用いてグラフをサブグラフに分割し，シードノードを含む
サブグラフに対してのみ RWR を計算する．すなわちサブグ
ラフに含まれないノードの関連度は 0 にする．彼らの手法はシードノードから近いノードのみが関連度が高くなるという性質を用いている．同じ性質を用いて K-dash は上位 K 個のノードを高速に検索する手法である [22]．関連度の上限値を用
いて K-dash は正確かつ高速に検索が可能であるが，この手法
は PPR の様に複数のノードがシードノードとなる場合は処理できず，また本論文で対象とした問題を包括的に処理するものでもない．Tong らは固有値分解 [16] を用いて RWR を近似する手法を提案した [23]．彼らは彼らの手法が保証された近似誤差を持つことを示したが，この近似誤差は normalized graphLaplacian [24] をグラフを表現する形式にのみ成立するという制約がある．

7. ま と め
本論文では PPR に対して特定のノードの関連度と関連度の高いノードを高速かつ正確に求める問題に取り組んだ．提案手法は特定ノードの関連度を疎な行列を用いて繰り返し計算を用いずに求め，関連度の下限値と上限値を用いてノードを枝刈りすることで関連度の高いノードを求める．実験により提案手法
が既存の手法より大幅に高速に処理が可能であることを示した．

PPR は多くのアプリケーションに利用されるが，提案手法を用いることによって今後多くのアプリケーションに対して精度を犠牲にすることなく高速化が可能になることが期待される．
文

献[1] A. Balmin, V. Hristidis and Y. Papakonstantinou: “ObjectRank: Authority-based Keyword Search in Databases”,
VLDB, pp. 564–575 (2004).
[2] B. Gallagher, H. Tong, T. Eliassi-Rad and C. Faloutsos:“Using Ghost Edges for Classiﬁcation in Sparsely Labeled
Networks”, KDD, pp. 256–264 (2008).
[3] B. Gao, T.-Y. Liu, W. Wei, T. Wang and H. Li: “Semisupervised Ranking on Very Large Graphs with Rich Metadata”, KDD, pp. 96–104 (2011).
[4] M. Gori and A. Pucci: “ItemRank: A Random-walk Based
Scoring Algorithm for Recommender Engines”, IJCAI, pp.
2766–2771 (2007).
[5] V. Hristidis, L. Raschid and Y. Wu: “Scalable Link-based
Personalization for Ranking in Entity-relationship Graphs”,
WebDB (2011).

[6] G. Jeh and J. Widom: “SimRank: a Measure of Structuralcontext Similarity”, KDD, pp. 538–543 (2002).
[7] D. Fogaras, B. Rácz, K. Csalogány and T. Sarlós: “Towards Scaling Fully Personalized PageRank: Algorithms,
Lower Bounds, and Experiments”, Internet Mathematics,
2, 3 (2005).
[8] G. Jeh and J. Widom: “Scaling Personalized Web Search”,
WWW, pp. 271–279 (2003).
[9] S. White and P. Smyth: “Algorithms for Estimating Relative Importance in Networks”, KDD, pp. 266–275 (2003).
[10] F. Fouss, A. Pirotte, J.-M. Renders and M. Saerens:“Random-walk Computation of Similarities between Nodesof a Graph with Application to Collaborative recommendation”, IEEE Trans. Knowl. Data Eng., 19, 3, pp. 355–369(2007).
[11] Q. Liu, E. Chen, H. Xiong and C. H. Q. Ding: “Exploiting User Interests for Collaborative Filtering: Interests Expansion via Personalized Ranking”, CIKM, pp. 1697–1700(2010).
[12] D. Liben-Nowell and J. M. Kleinberg: “The Link Prediction
Problem for Social Networks”, CIKM, pp. 556–559 (2003).
[13] K. Avrachenkov, N. Litvak, D. Nemirovsky, E. Smirnova
and M. Sokol: “Quick Detection of Top-k Personalized
PageRank lists”, WAW, pp. 50–61 (2011).
[14] B. Bahmani, A. Chowdhury and A. Goel: “Fast Incremental and Personalized PageRank”, PVLDB, 4, 3, pp. 173–184(2010).
[15] M. S. Gupta, A. Pathak and S. Chakrabarti: “Fast Algorithms for Top-k Personalized PageRank Queries”, WWW,
pp. 1225–1226 (2008).
[16] D. A. Harville: “Matrix Algebra From a Statistician’s Perspective”, Springer (2008).
[17] M. Yannakakis: “Computing the Minimum ﬁll-in is NPcomplete”, SIAM. J. on Algebraic and Discrete Methods,
2, 1, pp. 77–79 (1981).
[18] W. H. Press, S. A. Teukolsky, W. T. Vetterling and B. P.
Flannery: “Numerical Recipes 3rd Edition”, Cambridge
University Press (2007).
[19] S. Chakrabarti, A. Pathak and M. Gupta: “Index Design and Query Processing for Graph Conductance Search”,
VLDB J., 20, 3, pp. 445–470 (2011).
[20] J. Sun, H. Qu, D. Chakrabarti and C. Faloutsos: “Neighborhood Formation and Anomaly Detection in Bipartite
Graphs”, ICDM, pp. 418–425 (2005).
[21] G. Karypis, V. Kumar and V. Kumar: “Multilevel K-way
Partitioning Scheme for Irregular Graphs”, Journal of Parallel and Distributed Computing, 48, pp. 96–129 (1998).
[22] Y. Fujiwara, M. Nakatsuji, M. Onizuka and M. Kitsuregawa: “Fast and Exact Top-k Search for Random Walk with
Restart”, PVLDB, 5, 5, pp. 442–453 (2012).
[23] H. Tong, C. Faloutsos and J.-Y. Pan: “Fast Random Walk
with Restart and Its Applications”, ICDM, pp. 613–622(2006).
[24] D. Zhou, O. Bousquet, T. N. Lal, J. Weston and
B. Schölkopf: “Learning with Local and Global Consistency”, NIPS (2003).