DouKaku? 2351 にしお: #2349の7行目は何をやっているのかよ...(PageRankの計算)

Comment detail

This comment is reply for 2349 匿名: 正直、行列というか数学自体の素養がないの...(PageRankの計算). Go to thread root.

#2349の7行目は何をやっているのかよくわかりません…
__getslice__がオーバーライドされてるんでしょうか…。

とりあえず僕のコードを貼ってみます。
数学的にはさほど難しいことはしていなくて、
単に固有値を求めて絶対値を取ってから
合計でそれぞれを割って合計が1になるように正規化しています。

data = {
    1: [2, 3, 4, 5, 7],
    2: [1],
    3: [1, 2],
    4: [2, 3, 5],
    5: [1, 3, 4, 6],
    6: [1, 5],
    7: [5],
}

m = [[1.0 / len(data[j]) if i in data[j] else 0
      for j in range(1, 8)]
      for i in range(1, 8)]

from numpy import array
from numpy.linalg import eig

(v, d) = eig(array(m))

pagerank = [abs(row[0]) for row in d]
k = sum(pagerank)
pagerank = [x / k for x in pagerank]

print pagerank
# [0.303514376997, 0.166134185304, 0.140575079872,
#  0.105431309904, 0.178913738019, 0.0447284345048,
#  0.0607028753994]

v[:, i] でi列目のスライスを得られるようですけど、このあたりは
にしおさんのコードでは20行目でrow[0]と、0で決め打ちになってますが、
オリジナルのOctaveのコードではfind(abs(diag(D)) == max(abs(diag(D))))と、
計算でだしてますよね。
この部分を推測をまじえ、固有値の絶対値の最大値が存在する位置をインデックス値として
一列分取り出す、と解釈したものがあの7行目です。

原理が理解できてないんで、ここが0でいいのか、きちんと計算すべきなのか、
わかりません。

ついでに、Numeric付属のMLabを使ったバージョンも投稿します。

MLabのeig関数とSciPyのeig関数では戻り値の行と列が逆になっているような気が。

import MLab

def pagerank(data):
  a = [[1.0 / len(data[i]) if j in data[i] else 0 for i in range(1, len(data)+1)] for j in range(1, len(data)+1)]
  d, v = MLab.eig(MLab.array(a))
  ev = v[MLab.nonzero(d == MLab.max(MLab.absolute(d)))[0], :]
  return (ev / MLab.sum(ev)).real

data = {
  1: [2, 3, 4, 5, 7],
  2: [1],
  3: [1, 2],
  4: [2, 3, 5],
  5: [1, 3, 4, 6],
  6: [1, 5],
  7: [5],
}

for x in pagerank(data): print x

>固有値の絶対値の最大値が存在する位置をインデックス値として
>一列分取り出す、と解釈したものがあの7行目です。

なるほどなるほど。
「最大固有値は(少なくともNumPyを使って固有値を求めた場合は)0番目にある」
と解釈したのが僕のコードです。

最大固有値が何番目にあるのかをチェックするバージョンも作ってみました。

# あと余談ですがnumpyはLAPACKを使っているということがわかったのでタグを付けときます

see: LAPACK

data = {
    1: [2, 3, 4, 5, 7],
    2: [1],
    3: [1, 2],
    4: [2, 3, 5],
    5: [1, 3, 4, 6],
    6: [1, 5],
    7: [5],
}

m = [[1.0 / len(data[j]) if i in data[j] else 0
      for j in range(1, 8)]
      for i in range(1, 8)]

from numpy import array
from numpy.linalg import eig

(v, d) = eig(array(m))

max_v = list(v).index(max(v))
pagerank = [abs(row[max_v]) for row in d]
k = sum(pagerank)
pagerank = [x / k for x in pagerank]

print pagerank

知識がある人には0でOKという、目算がつくようですね。
これは知識の差がコードの違いに現れるという、よくある
状況のわかりやすい一例かもしれませんね。

ついでにarrayのスライスについて、ちょっと書いておきますと
「a[:, i]」という構文はちょっと見、ありえなさそうに思えますが、
これが「a[i1:i2:i3, j1:j2:j3]」の省略形と分かれば納得できるでしょう。
大体予想できるでしょうが、これは行方向のスライスと列方向のスライスを
コンマで区切った表記方法です。

ざっと試したところ、MLab.array、numpy.array、scipy.arrayの３つ共に
使えるようです。

DouKaku?

Comment detail

Index

Feed

Other

Link

Pathtraq