1 big factor PR(E), 1

.#query.py

@@ -0,0 +1 @@
+sholtzen@wifi-131-179-2-210.host.ucla.edu.13952

my_graphs.py

@@ -50,6 +50,29 @@ def gen_complete_pairwise_factor(n):
     g.add_edges(edges)
     return (g, (smokers, factors))
 
+
+### generates a factor graph with one big factor and n variables connected to it
+### returns a the graph and a tuple whose first element are the variables
+### and second element are the factors
+def gen_single_big_factor(n):
+    g = Graph(sparse=True)
+    # make n smoker vertices
+    v = [x for x in range(0,n)]
+    # connect all the smoker to the factors
+    # make friends
+    factors = [n]
+    edges = []
+    # friends = []
+    # friendedges = []
+    count = n
+    for va in v:
+        edges += [(n, va)]
+
+    g.add_vertices(v)
+    g.add_vertices(factors)
+    g.add_edges(edges)
+    return (g, (v, factors))
+
 ### generates a complete pairwise factor graph
 ### returns a the graph and a tuple whose first element are the variables
 ### and second element are the factors

orbitgen.py

@@ -47,6 +47,9 @@ def bfs_foldrep(graph, colors, fixcolors, acc, f, orders=False):
     # queue of colorings yet to be considered
     queue = deque()
     queue.append([[], colors])
+    g_order = None
+    if orders:
+        g_order = graph.automorphism_group().order()
     # set of representative graph colorings
     reps = set()
     while len(queue) != 0:
@@ -72,7 +75,10 @@ def bfs_foldrep(graph, colors, fixcolors, acc, f, orders=False):
 
         A, orbits = graph.automorphism_group(partition=c + fixcolors, orbits=True,
                                              algorithm="bliss")
-        acc = f(acc, graph, c)
+        if orders:
+            acc = f(acc, graph, c, g_order, A.order())
+        else:
+            acc = f(acc, graph, c)
         reps.add((gcanon, (tuple(c_canon[0]), tuple(c_canon[1]))))
         print(len(reps))
         if len(c_canon[0]) + 1 <= len(colors) / 2.0:

query.py

@@ -57,7 +57,7 @@ def potential(assgn):
 def mpe(G, variables, factors, potential):
     # folding function
     print("G: %s, vars: %s, factors: %s" % (G, variables, factors))
-    def sum_prob(acc, g, color):
+    def max_prob(acc, g, color):
         (cur_max_state, cur_max_value) = acc
         if len(color[0]) < len(g.vertices()) / 2.0:
             # check both possibilities
@@ -70,11 +70,31 @@ def sum_prob(acc, g, color):
             cur_max_value = pot
             cur_max_state = color
         return (cur_max_state, cur_max_value)
-    return bfs_foldrep(G, variables, factors, (None, 0.0), sum_prob)
+    return bfs_foldrep(G, variables, factors, (None, 0.0), max_prob)
+
+
+# potential returns a pair (prob, const) which is the probability and
+# normalizing constant on a particular orbit
+def prob(G, variables, factors, potential):
+    # folding function
+    print("G: %s, vars: %s, factors: %s" % (G, variables, factors))
+    def sum_prob(acc, g, color, g_order, cur_order):
+        cur_sz = g_order / cur_order # yay orbit stabilizer theorem!
+        (cur_prob, cur_z) = acc
+        if len(color[0]) < len(g.vertices()) / 2.0:
+            # check both possibilities
+            (prob, z) = potential(color[::-1])
+            cur_prob += prob * cur_sz
+            cur_z += z * cur_sz
+        (prob, z) = potential(color[::-1])
+        cur_prob += prob * cur_sz
+        cur_z += z * cur_sz
+        return (cur_prob, cur_z)
+    return bfs_foldrep(G, variables, factors, (0.0, 0.0), sum_prob, orders=True)
 
 
 def compute_mpe_complete_2factor():
-    G = gen_complete_pairwise_factor(10)
+    G = gen_complete_pairwise_factor(100)
     # add evidence: the first variable is false
     G[0].add_vertex(name="e")
     G[0].add_edge(("e", G[1][0][0]))
@@ -108,9 +128,35 @@ def potential(assgn):
     ps.print_stats()
     print s.getvalue()
 
+def compute_prob_big_factor():
+    G = gen_single_big_factor(10)
+    # add evidence: the first variable is false
+    G[0].add_vertex(name="e")
+    G[0].add_edge(("e", G[1][0][0]))
+    def potential(assgn):
+        p = 0.0
+        for v in assgn[0]:
+            p += 1
 
+        if G[1][0][0] in assgn[0]:
+            return (0.0, p)
+        else:
+            return (p, p)
 
+    pr = cProfile.Profile()
+    pr.enable()
+
+    res = prob(G[0], G[1][0], [G[1][1] + ["e"]], potential)
+    print("prob %f" % (res[0] / res[1]))
+    # print("brute forced: %d" % (bruteforce_partition(G, partfun)))
+
+    pr.disable()
+    s = StringIO.StringIO()
+    sortby = 'cumulative'
+    ps = pstats.Stats(pr, stream=s).sort_stats(sortby)
+    ps.print_stats()
+    print s.getvalue()
 
 
 if __name__ == "__main__":
-    compute_mpe_complete_2factor()
+    compute_prob_big_factor()