Haskell 剖析子图挖掘算法

Posted 2023-03-06

【中文标题】Haskell 剖析子图挖掘算法

【英文标题】：Haskell profilling subgraph mining algorithm

【发布时间】：2013-09-23 22:38:57

【问题描述】：

我尝试解决在 Haskell 中查找所有连接子图的问题。使用的算法描述为here。引用那篇论文：

在每个路径算法中，都有前进和后退。如果给定的连通子图可以通过添加边 k 来扩展，即如果边 k 还不是给定子图的一部分，如果 k 与给定子图的至少一条边相邻，并且如果添加下面给出的一些限制不禁止边 k 的。 一旦给定的连接子图不能被进一步拉长，就会后退一步。在这种情况下，最后添加的边被从字符串中删除，它暂时被赋予“禁止”状态，并且通过从先前较长的字符串回溯而被禁止的任何其他边同时再次“允许”。相反，通过从比当前短的字符串中删除而被禁止的边仍然被禁止，从而确保每个连通子图都构造一次且仅一次。

为了实现这个算法，我将图表示为边列表：

 type Edge =  (Int,Int)
 type Graph = [Edge]

首先，我编写了函数addEdge，它检查是否可以扩展图形，如果不可能，则返回Nothing或Edge进行扩展。

我有一个"parent" 图形和"extensible" 图形，所以我尝试找到一个且只有一个存在于"parent" 图形中的边，与"extensible" 图形相连，尚未包含在"extensible" 图形中和所以不包含在forbidden 集中。

我在下面写了这个函数：

addEdge :: Graph -> Graph -> [Edge] -> Maybe Edge
addEdge !parent !extensible !forb = listToMaybe $ intersectBy (\ (i,j) (k,l) -> (i == k || i == l || j == k || j == l)) (parent \\ (extensible `union` forb)) extensible

这是工作！但是，正如我从分析整个程序中看到的那样，addEdge 是最重的函数。我敢肯定，我的代码不是最优的。 Leastways，intersectBy 函数可以找到所有可能的解决方案，但我只需要一个。有没有办法让这段代码更快？也许，不要使用标准列表，而是使用Set from Data.Set？这是第一个注意点。

主递归函数ext如下所示：

ext :: Graph -> [Graph] -> Maybe Graph -> [(Edge,Int)] -> Int -> [Graph]
ext !main !list !grow !forb !maxLength      | isEnd  == True = (filter (\g -> (length g /= 1)) list) ++ (group main) 
                                            | ((addEdge main workGraph forbEdges) == Nothing) || (length workGraph) >= maxLength = ext main list (Just workGraph) forbProcess maxLength
                                            | otherwise = ext main ((addedEdge:workGraph):list) Nothing forb  maxLength where 
                                                workGraph = if grow == Nothing then (head list) else (bite (fromJust grow)) -- [Edge] graph now proceeded
                                                workGraphLength = length workGraph
                                                addedEdge = fromJust  $ addEdge'
                                                addEdge' = addEdge main workGraph forbEdges
                                                bite xz = if (length xz == 1) then (fromJust (addEdge main xz forbEdges)):[] else tail xz 
                                                forbProcess = (head workGraph,workGraphLength):(filter ((<=workGraphLength).snd) forb)
                                                forbEdges = map fst forb -- convert from (Edge,Level) to [Edge]                     
                                                isEnd = (grow /= Nothing) && (length (fromJust grow) == 1) && ((addEdge main (fromJust grow) forbEdges) == Nothing)

我在图表上测试我的程序

c60 = [(1,4),(1,3),(1,2),(2,6),(2,5),(3,10),(3,7),(4,24),(4,21),(5,8),(5,7),(6,28),(6,25),
    (7,9),(8,11),(8,12),(9,16),(9,13),(10,20),(10,17),(11,14),(11,13),(12,28),(12,30),(13,15),
    (14,43),(14,30),(15,44),(15,18),(16,18),(16,17),(17,19),(18,47),(19,48),(19,22),(20,22),(20,21),
    (21,23),(22,31),(23,32),(23,26),(24,26),(24,25),(25,27),(26,35),(27,36),(27,29),(28,29),(29,39),
    (30,40),(31,32),(31,33),(32,34),(33,50),(33,55),(34,37),(34,55),(35,36),(35,37),(36,38),(37,57),
    (38,41),(38,57),(39,40),(39,41),(40,42),(41,59),(42,45),(42,59),(43,44),(43,45),(44,46),(45,51),
    (46,49),(46,51),(47,48),(47,49),(48,50),(49,53),(50,53),(51,52),(52,60),(52,54),(53,54),(54,56),(55,56),(56,58),(57,58),(58,60),(59,60)] :: Graph

例如，查找长度从 1 到 7 的所有子图

length $ ext c60 [[(1,2)]] Nothing [] 7
>102332

问题是太低速的计算。正如它在原始文章中指出的那样，程序已在FORTRAN 77 中编写并在 150MHz 工作站上启动，执行测试任务的速度至少比我在现代 i5 处理器上的代码快 30 倍。 我不明白，为什么我的程序这么慢？ 有没有办法重构这段代码？或者最好的解决方案是将其移植到 C 上，并通过 FFI 编写与 C 库的绑定？

【问题讨论】：

(\\) 很慢，O(n^2)。改用Data.Set.(\\)，它是O(n+m)

【参考方案1】：

我决定尝试使用fgl 来实现论文中描述的算法。完整代码如下。

-# LANGUAGE NoMonomorphismRestriction #-

import Data.Graph.Inductive
import Data.List
import Data.Tree

uniq = map head . group . sort . map (\(a, b) -> (min a b, max a b))
delEdgeLU (from, to) = delEdge (from, to) . delEdge (to, from)
insEdgeDU (from, to) = insEdge (from, to, ()) . insNodeU to . insNodeU from where
    insNodeU n g = if gelem n g then g else insNode (n, ()) g

nextEdges subgraph remaining
    | isEmpty subgraph = uniq (edges remaining)
    | otherwise = uniq $ do
        n  <- nodes subgraph
        n' <- suc remaining n
        return (n, n')

search_ subgraph remaining
    = Node subgraph
    . snd . mapAccumL step remaining
    $ nextEdges subgraph remaining
    where
    step r e = let r' = delEdgeLU e r in (r', search_ (insEdgeDU e subgraph) r')

search = search_ empty

mkUUGraph :: [(Int, Int)] -> Gr () ()
mkUUGraph es = mkUGraph ns (es ++ map swap es) where
    ns = nub (map fst es ++ map snd es)
    swap (a, b) = (b, a)

-- the one from the paper
sampleGraph = mkUUGraph cPaper
cPaper = [(1, 2), (1, 5), (1, 6), (2, 3), (3, 4), (4, 5)]

您要在顶层使用的函数是mkUUGraph，它从边列表构造一个图，以及search，它构造一棵树，其节点是其输入的连接子图。例如，要计算论文中“方案 1”底部显示的统计数据，您可以这样做：

*Main> map length . tail . levels . search . mkUUGraph $ [(1, 2), (1, 5), (1, 6), (2, 3), (3, 4), (4, 5)]
[6,7,8,9,6,1]
*Main> sum it
37

将它与您的实现进行比较我有点麻烦，因为我不明白 ext 的所有参数应该做什么。特别是，我无法弄清楚如何在论文中的邻接图上调用ext，从而得到 37 个结果。也许你有一个错误。

在任何情况下，我都尽力模仿我认为您的代码正在尝试做的事情：查找最多具有七个边的图形，并且肯定包含边 (1, 2)（尽管您的代码输出了许多图形不包含(1, 2))。我添加了这段代码：

mainHim = print . length $ ext c60 [[(1,2)]] Nothing [] 7
mainMe  = print . length . concat . take 7 . levels $ search_ (mkUUGraph [(1,2)]) (mkUUGraph c60)

我的代码找到了 3301 个这样的图表；你的发现 35571。我没有很努力地找出差异来自哪里。在 ghci 中，mainHim 需要 36.45s； mainMe 需要 0.13 秒。使用-O2编译时，mainHim耗时4.65s； mainMe 需要 0.05 秒。 mainMe 的数字可以通过使用PatriciaTree 图形实现而不是默认的图形实现再次减半，并且可能通过分析和一些想法进一步减少。以防mainMe 速度如此之快的原因是它发现的图少得多，我还测试了修改后的main：

main = print . length . concat . take 8 . levels $ (search (mkUUGraph c60) :: Tree (Gr () ()))

这会打印 35853，因此它找到的图形数量与您的测试命令大致相同。 ghci 0.72s，-O2编译0.38s。

【参考方案2】：

或者最好的解决方案是将其移植到 C 上，并通过 FFI 将绑定写入 C 库？

不，您不必用 C 编写它。GHC 生成的代码并不比 C 慢多少。这种巨大的速度差异表明您正在实现不同的算法。因此，与其用不同的语言重写，不如重写 Haskell 代码。

我猜你的代码的问题是你......

使用列表而不是集合 使用广度优先而不是深度优先枚举（不确定） 对整个边集使用操作，而不是巧妙地跟踪哪些边在哪个集合中 手动编码算法的递归结构，而不是使用递归调用。

我不得不承认我没有完全理解你的代码。但我读了你链接到的论文，那里描述的算法似乎是对所有结果的简单暴力枚举。所以我猜 Haskell 实现应该使用列表单子（或列表推导）来枚举所有子图，在枚举期间过滤掉非连接子图。如果您以前从未使用 list monad 编写过代码，那么仅枚举所有子图可能是一个很好的起点。