Sum-Product on Factor Tree Graphs, MAP Elimination Algorithm

Massachusetts Institute of Technology Department of Electrical Engineering and Computer Science 6.438 Algorithms For Inference 10 Sum-product on factor tree graphs, MAP elim­ ination algorithm The previous lecture introduced the beginnings of the sum-product on factor tree graphs with a specific example. In this lecture, we will describe the sum-product algorithm for general factor tree graphs. 10.1 Sum-product for factor trees Recall that a factor graph consists of a bipartite graph between variable nodes and factor nodes. The graph represents the factorization of the distribution into factors corresponding to factor nodes in the graph, in other words p(x1 , . . . , xN ) ∝ M M fi (xci ) i=1 where ci is the set of variables connected to factor fi . As in the case with undirected graphs, we restrict our attention to a special class of graphs where inference is efficient. Definition 1 (Factor tree graph) A factor tree graph is a factor graph such that the variable nodes form a tree, that is between any two variable nodes, there is a unique path (involving factors). 10.1.1 Serial version As in the sum-product algorithm for undirected graphs, there is a serial verison of the sum-product algorithm for factor tree graphs. Because there are factor nodes and variables nodes, we have two types of messages: 1. Variable to factor messages: xi b mi→a (xi ) = M a mb→i (xi ) b∈N (i)\a The variable node i with degree di sends a message towards a that is the product of di −1 terms for each value in |X|. This results in a total of O(|X|di ) operations. 1 2. Factor to variable messages: xj ma→i (xi ) = xj :j∈N (a)\i xi a fa (xi ; xj , j ∈ N (a)\i) M mj→a (xj ) j∈N (a)\i The factor node a with degree da sends a message towards i that is a summation over |X|da −1 terms, each involving a product of da terms for each value in |X| of the neighboring variables. Thus, the total cost is O(|X|da da ). ∗ In summary, complexity scales as O(|E||X|d ) where d∗ = maxa da (over factor nodes) and |E| is the total number of edges in the factor tree graph. Notably, because the graph is a tree, |E| = O(N ).1 10.1.2 Parallel version We can also run the updates in parallel to get the parallel version of the sum-product algorithm for factor tree graphs: 1. Initialization: for (i, a) ∈ E, m0i→a (xi ) = 1 m0a→i (xi ) = 1 for all xi ∈ X. 2. Update: for (i, a) ∈ E, t ≥ 0, M t (*) mt+1 (x ) = mb→i (xi ) i→a i b∈N (i)\a (**) mt+1 a→i (xi ) = xj :j∈N (a)\i fa (xi ; xj , j ∈ N (a)\i) M mtj→a (xj ) j∈N (a)\i 3. Marginalization: for t ≥ diameter pt+1 xi (xi ) ∝ 1 M t+1 ma→i (xi ). a∈N (i) This is the cost to send all of the messages to a single node. If we wanted to send all of the ∗ ∗ messages the other direction, this would naively cost O(N 2 |X|d ). We can reduce this to O(N |X|d ) using a technique that we show in the next section. 2 As we noticed above sending ma→i costs O(|X|da da ) and we have to send a similar mes­ ∗ a 2 sage for every neighbor ofaa, so the total complexity cost will scale as O(|X|d a da ) and in the worst case O( a d2a ) will be quadratic in N . With a little cleverness, we can provide a more efficient implementation of (*) and (**): (*) Compute mt+1 i (xi ) = M mtb→i (xi ) b∈N (i) then, mt+1 i→a (xi ) = mit+1 (xi ) . mta→i (xi ) Thus the total cost of computing (*) per variable node i is O(|X|di ) instead of O(|X|d2i ). (**) Compute f˜at+1 (xk , k ∈ N (a)) = fa (xk ; k ∈ N (a)) M t mk→a (xk ) k∈N (a) then, mt+1 a→i (xi ) = X xj :j∈N (a)\i f˜at+1 (xk , k ∈ N (a)) . mti→a (xi ) Thus the total cost of computing (**) per factor node a is O(|X|da da ). Thus the overall complexity cost for a single iteration of the parallel algorithm will ∗ be O(N |X|d ). 10.2 Factor tree graph > Undirected tree graph The factor tree graph representation is strictly more general than the undirected tree graph representation. Recall that by the Hammersley-Clifford theorem, we can represent any undirected graph as a factor graph by creating a factor node for each clique potential and connecting the respective variable nodes. In the case of an undirected tree graph, the only clique potentials are along edges, so it’s clear the factor graph representation will also be a tree. Now, we consider an example where the undirected graph is not a tree, but its factor graph is: 3 As shown in the figure, the graphical model on the left has 4 maximal cliques. The factor graph representation of it is on the right, and it does not have any loops, even though the undirected graph did. A necessary and sufficient condition for this situation is provided in the solution of Problem 2.2. On the left hand graph, we cannot apply the sum-product algorithm because the graph is not a tree, but on the right hand side we can. What happened here? The factor in the middle of the graph has 3 nodes, and as we saw above, the overall complexity of the algorithm depends strongly on the maximal degree node. Taking this example further, consider converting the complete undirected graph on N vertices into a factor graph. It is a star graph with a single factor node, so it is a tree Unfortunately, the factor node has N vertices, so the cost of running sum-product is O(N |X|N ), so we don’t get anything for free. 10.3 MAP elimination algorithm The primary focus of the prior material has been on algorithms for computing marginal distributions. In this section, we will be talking about the second inference problem: computing the mode or MAP. Specifically, consider a collection of N random variables x1 , . . . xN each taking values in |X| with their joint distribution represented by a graphical model G. For concreteness, let G = (V, E) be an undirected graphical model, so that the distribution factorizes as: M px (x) ∝ φC (xC ) C∈C = 1 M φC (xC ), Z C∈C where C is the set of maximal cliques in G and Z is the partition function. Our goal is to compute x∗ ∈ XN such that x∗ ∈ arg max px (x). x∈XN 4 Note that we have written ∈ instead of = because there may be several configurations that result in the maximum and we are interested in finding one of them. The sum-product algorithm gives us a way to compute the marginals pxi (xi ), so we might wonder if setting x∗i ∈ arg maxxi pxi (xi ) produces a MAP assignment. This works in situations where the variables are independent, there is only one positive configuration, or when the variables are jointly Gaussian. However, consider the following example with binary variables x1 and x2 : px1 ,x2 (0, 0) = 1/3 − t px1 ,x2 (1, 0) = 1/3 px1 ,x2 (0, 1) = 0 px1 ,x2 (1, 1) = 1/3 + t. arg maxx1 px1 (x1 ) = 1 and arg maxx2 px2 (x2 ) = 0, but that is not the MAP. It is clear that we need a more general algorithm. 10.4 MAP Elimination algorithm example We will start with a motivating example of the MAP elimination algorithm. Consider a distribution over (x1 , . . . , x5 ) described by the graphical model below. Let each xi ∈ {0, 1}, for 1 ≤ i ≤ 5. 35 =1 x3 13 = x5 1 x1 34 12 =1 x2 =1 x4 24 = 1 p(x) ∝ exp(θ12 x1 x2 + θ13 x1 x3 + θ24 x2 x4 + θ34 x3 x4 + θ35 x3 x5 ) For this example, we will consider the specific values θ12 = θ34 = θ35 = 1, θ13 = θ24 = −1. So, plugging in those values, our distribution is p(x) ∝ exp(x1 x2 − x1 x3 − x2 x4 + x3 x4 + x3 x5 ) 1 = exp(x1 x2 − x1 x3 − x2 x4 + x3 x4 + x3 x5 ) � Z� F (x1 ,...,x5 ) 5 Our goal is to efficiently compute: x∗ ∈ arg max px (x) ∝ F (x1 , . . . , x5 ). x∈XN To efficiently compute this mode x∗ , let us adapt the elimination algorithm for computing the marginals of a distribution. First, we will fix an elimination order of nodes: 5, 4, 3, 2, 1. Then max px (x) = max exp(x1 x2 − x1 x3 − x2 x4 + x3 x4 ) max exp(x3 x5 ) x x1 ,...,x4 x5 |� {z �} �m5 (x3 ) where we have grouped the terms involving x5 and moved the maxx5 past all of the terms that do not involve x5 . We see that m5 (x3 ) = ex3 here. We also store the specific value of x5 that maximized the expression (e.g. when x3 = 1, m5 (x3 ) = e and x∗5 = 1) and we will see why this is important to track later. In this case, the maximizer of x5 is 1 if x3 = 1 and 0 or 1 if x3 = 0. We are free to break ties arbitrarily, so we break the tie so that x∗5 (x3 ) = x3 . m5 (x3 ) = exp(x3 ) x5 (x3 ) = 1 x3 13 = 1 x1 34 12 =1 x2 =1 x4 24 = 1 Now we eliminate x4 in the same way, by grouping the terms involving x4 and moving the maxx4 past terms that do not involve it. max px (x) = max exp(x1 x2 − x1 x3 )m5 (x3 ) max exp(x3 x4 − x2 x4 ) x x1 ,x2 ,x3 x4 |� {z }� �m4 (x2 ,x3 ) with m4 (x2 , x3 ) = exp(x3 (1 − x2 )) x∗4 (x2 , x3 ) = 1 − x2 . 6 x3 13 = 1 m4 (x2 , x3 ) = exp(x2 (1 x1 12 x4 (x2 , x3 ) = 1 x3 )) x2 x2 =1 Similarly for x3 ,   max px (x) = max exp(x1 x2 ) max exp(−x1 x3 )m5 (x3 )m4 (x2 , x3 ) x x1 ,x2 x3   = max exp(x1 x2 ) max exp(−x1 x3 + x3 + x3 (1 − x2 )) x1 ,x2 x3   = max exp(x1 x2 ) max exp(−(x1 + x2 )x3 ) x1 ,x2 x3 �| {z }� �m3 (x1 ,x2 ) , with m3 (x1 , x2 ) = 1 x∗3 (x1 , x2 ) = 0. Finally, max px (x) = max exp(x1 x2 ) x x1 ,x2 implies that x∗1 = x∗2 = 1. Now we can use the information we stored along the way to compute the MAP assignment. x∗3 (1, 1) = 0, x4∗ (1, 0) = 0, x5∗ (0) = 1. Thus the optimal assignment is: x∗ = (1, 1, 0, 0, 1). As in the elimination algorithm for computing marginals, the primary determinant of the computation cost is the size of the maximal clique in the reconstituted graph. Putting this together, we can describe the MAP elimination algorithm for general graphs: 7 Input: Potentials ϕc for c ∈ C and an elimination ordering I Output: x∗ ∈ arg maxx∈XN px (·) Initialize active potentials Ψ to be the set of input potentials. for node i in I do Let Si be the set of all nodes (not including i and previously eliminated nodes) that share a potential with node i. Let Ψi be the set of potentials in Ψ involving xi . Compute M mi (xSi ) = max ϕ(xi ∪ xSi ) xi x∗i (xSi ) ϕ∈Ψi ∈ arg max xi M ϕ∈Ψi ϕ(xi ∪ xSi ) where ties are broken arbitarily. Remove elements of Ψi from Ψ. Add mi to Ψ. end Produce x∗ by traversing I in a reverse order x∗j = x∗j (xk∗ : j < k ≤ N ). Algorithm 1: The MAP Elimination Algorithm As in the context of computing marginals, the overall cost of the MAP elimination algorithm is bounded above as X overall cost ≤ |C| |X||Si |+1 i ≤ N |C||X|maxi |Si |+1 . The MAP elimination algorithm and elimination algorithm for marginals are closely related and they both rely on the relationship between sums and products and max and products. Specifically, the key property was that   max f (x)g(x, y) = max f (x) max g(x, y) x,y X x,y x f (x)g(x, y) = X x 8 y f (x) X y g(x, y) .