Poorman's PageRank | 从零开始 Haskel
PageRank 算法是一种经典的网页排名算法。基本思想是，每个节点首先赋相等的初值。接下来，根据链接关系将值传播到链接去的节点。如此迭代直到收敛。
需要特殊处理的地方是，出度为 0 的节点需要将值保存到自己。
为了避免自私的节点不引用别人，从而大量积累自己的值，进行平滑处理。给每一个节点乘以缩减因子 ，再将每个节点加上相等的
。注意到这种平滑不改变总值。也即任何时刻所有节点的值之和恒为 1 。
与之相关的还有 特征向量中心度 eigenvector centrality ，其区别是，不处理出度为 0 的点，也不进行平滑。而在每一步进行正规化。此外，特征向量也可以使用入度作为标准，仅需将连接矩阵转置即可。
这里给出一种简洁的三合一 Haskell 实现。不使用任何复杂的库函数，仅用 80 行。从中可以看到 Haskell 的简洁和抽象能力。
三种算法的核心都是不断迭代直到收敛。将这一逻辑抽象出来得到：
converge :: Eq a => (a -> a) -> a -> a
converge f v = fst $ until theSame update (v, f v)
where
theSame (x, y) = x == y
update (x, y) = (y, f y)
这里用到了库函数 until :: (a -> Bool) -> (a -> a) -> a -> a 。这个函数接收一个判断函数，一个更新函数和初值。当判断函数返回假时，会应用更新函数。当判断函数返回真时，返回最终值。
converge 函数实际上要构造一个流（stream），即 v : f v : f (f v) : f (f (f v)) : ... 。当流的两个连续元素相等时，我们找到了 f 这个函数的不动点，也就是最终的收敛值。
因为只需要比较前两个元素，所以我们使用两个元素的元组（tuple）作为保存的状态。until 的判断函数就是两个元素是否相等。更新函数是抛弃第一个元素，对第二个元素应用 f 。
接下来不同算法的区别，仅在更新函数不同。
对于 pageRank 来说，就是不断乘以连接矩阵：
pageRank :: [[Value]] -> [Value] -> [Value]
pageRank a vs = head $ converge (`matmul` a') [vs]
where
a' = compensate a
其中 matmul :: (Num a) => [[a]] -> [[a]] -> [[a]] 是矩阵乘法，将在下面给出实现。
注意到，首先将初值用列表改成 (n, 1) 的行向量，因此每次迭代改为右乘连接矩阵。最后使用 head 再转变成一维列表 (n,) 。下面各个算法做同样的处理。
compensate 函数实现两个功能，对于出度不为 0 的节点，将因子 1 平均分配到每个非零节点上；对于出度为 0 的节点，将 1 分配到自己的位置上（矩阵对角线）。
compensate :: [[Value]] -> [[Value]]
compensate = map procOut . zip [0 ..]
where
procOut (i, l) =
if any (/= 0) l
then distribute l
else oneAt i l
distribute l =
let v = 1.0 / (sum l)
in map
(\x ->
if x == 0
then x
else v)
l
oneAt i l =
let (x, _:ys) = splitAt i l
in x ++ 1.0 : ys
平滑处理可以改为对连接矩阵进行修改：
smooth :: Value -> [[Value]] -> [[Value]]
smooth s m = map (map interpolate) m
where
interpolate a = s * a + (1.0 - s) / fromIntegral n
n = length m
对每一个元素，都用因子 s 缩减，再加上补偿。
那么平滑后的 PageRank 算法如下：
smoothPageRank :: Value -> [[Value]] -> [Value] -> [Value]
smoothPageRank s a vs = head $ converge (`matmul` a') $ [vs]
where
a' = smooth s . compensate $ a
对于特征向量中心性，需要实现正规化：
normalize :: (Fractional a, Ord a) => [a] -> [a]
normalize vs =
let m = maximum . (map abs) $ vs
in map (/ m) vs
即将一个行向量的每个元素除以最大值。
那么特征向量中心性可以实现如下：
eiginCentr :: [[Value]] -> [Value] -> [Value]
eiginCentr a vs =
head $ converge ((map normalize) . (`matmul` a)) [vs]
以上已经实现了三个算法的核心部分。接下来给出辅助函数的直观定义。
矩阵乘法：
dot :: (Num a) => [a] -> [a] -> a
dot x y = sum $ zipWith (*) x y
matmul :: (Num a) => [[a]] -> [[a]] -> [[a]]
matmul a b = map rowMul a
where
b' = transpose b
rowMul r = map (dot r) b'
类型转换：
type Value = Double
aFromIntegral :: (Integral a) => [[a]] -> [[Value]]
aFromIntegral = map (map fromIntegral)
生成初始平均分配值：
normalDist :: Int -> [Value]
normalDist n = replicate n $ 1.0 / fromIntegral n
图从边表示转化为连接矩阵表示：
edgeToAdj :: (Integral a) => [(a, a)] -> [[a]]
edgeToAdj es = [[query i j | j <- [0 .. upper]] | i <- [0 .. upper]]
where
(ls, rs) = unzip es
vs = ls ++ rs
upper = maximum vs -- lower bound = 0
query i j =
if elem (i, j) es
then 1
else 0
其实这里使用 ST monad 更好一点，仅需要 的时间复杂度。这里用的是直接搜索，需要
的时间复杂度。
以上代码实现了所有三个算法的功能，仅用了 80 行代码。完整代码见 gist 。
使用下图进行测试：
Network Example
-- Test Graph 2
tg2e =
[ (0, 8)
, (1, 6)
, (1, 10)
, (1, 11)
, (2, 1)
, (2, 10)
, (2, 11)
, (3, 15)
, (3, 17)
, (4, 1)
, (4, 6)
, (4, 15)
, (5, 7)
, (5, 8)
, (5, 16)
, (6, 5)
, (6, 8)
, (6, 16)
, (7, 5)
, (7, 13)
, (7, 15)
, (8, 16)
, (8, 5)
, (8, 6)
, (9, 11)
, (9, 10)
, (9, 2)
, (10, 9)
, (10, 11)
, (10, 13)
, (11, 9)
, (11, 10)
, (11, 15)
, (12, 13)
, (12, 15)
, (12, 16)
, (13, 14)
, (13, 15)
, (13, 16)
, (14, 13)
, (14, 12)
, (14, 15)
, (15, 1)
, (15, 9)
, (15, 11)
, (16, 7)
, (16, 8)
, (16, 13)
]
tg2 = edgeToAdj tg2e
tg2spr = smoothPageRank 0.8 (aFromIntegral tg2) (normalDist . length $ tg2)
printTg2spr :: IO ()
printTg2spr = mapM_ (printf "%.3f\n") tg2spr
测试结果如下：
$ stack ghci
λ> :load pagerank.hs
[1 of 1] Compiling Main ( pagerank.hs, interpreted )
Ok, one module loaded.
λ> printTg2spr
0.011
0.049
0.034
0.011
0.011
0.054
0.045
0.048
0.069
0.087
0.084
0.104
0.020
0.083
0.033
0.095
0.083
0.078
λ>
符合预期。
连矩阵乘法都从头开始写，到整个算法完成，仅需要 80 行代码。核心就是 converge 函数的抽象。这个例子很好地体现了 Haskell 作为函数式语言的优点。
