Tropical geometry of statistical models

1. Algebraic Statistics, Tropical Geometry, and Inference

By a statistical model, we mean a family of joint probability distributions for a collection of discrete random variables Y = {Y ₁,...,Y _n}. Thesis i states that many families of interest can be characterized by polynomials in the joint probabilities p _σ1...σn = Prob(Y ₁ = σ₁,...,Y _n = σ_n). Although the variety defined by these polynomials contains points that are not in the model (for example, points with negative coordinates), the emerging field of algebraic statistics (1, 2) offers practical and useful algorithms for studying statistical models.

Tropicalization means replacing the arithmetic operations (+, ×) by the operations (min, +). This process captures the essence of what happens when the joint probabilities p _σ1...σn are replaced by their logarithms. The tropicalization of an algebraic variety is a piecewise-linear set that has many features familiar from algebraic geometry (3, 4). In particular, the tropicalization of a statistical model is a piecewise-linear set in the space with logarithmic coordinates -log(p _σ1...σn).

Thesis iii states that tropical algebraic geometry of statistical models is of fundamental interest in analyzing the behavior of inference algorithms under the variation of model parameters. By inference, we mean the evaluation of one or more coordinates of a single point on the algebraic variety, in either (+, ×) or (min, +) arithmetic. This evaluation corresponds to a form of inference that is used for graphical models in statistical learning theory (5), but it differs from other (more classical) notions of inference in mathematical statistics. By parametric inference, we mean the analysis of the dependence of inference on parameters.

To give a more concrete discussion of parametric inference, it is useful to focus on directed graphical models. A directed graphical model (or Bayesian network) is a finite directed acyclic graph G with two kinds of vertices (observed variables Y = {Y ₁,...,Y _n} and hidden variables X = {X ₁,..., X _m}), where each edge is labeled by a transition matrix whose entries are linear forms in some parameters. The rules of discrete probability express the observed probabilities p _σ1...σn as polynomials of a degree ≤ E in the parameters, where E is the number of edges of G. The polynomials parameterize the graphical model as an algebraic variety.

The following are two types of inference questions from statistical learning theory for graphical models.

1. The calculation of marginal probabilities:

   and

2. The calculation of maximum a posteriori (MAP) log probabilities:

   where the h _i range over all of the possible assignments for the hidden random variables X _i. Together, these two primitives can be used effectively to solve a range of other statistical learning inference problems, including the calculation of conditional probabilities and other quantities of interest. The key to these statistical learning inference questions for graphical models is the sum-product algorithm (6), which is also known as the generalized distributive law (7). This polynomial-time algorithm (in the case that the graph has constant clique size) is used, both in ordinary arithmetic (+, ×) and in tropical arithmetic (min, +), to solve problems 1 and 2 efficiently. For more background on the sum-product algorithm, and for connections to message passing and the junction tree algorithm, see ref. 5.

   Although the sum-product algorithm provides efficient solutions to the basic inference problems 1 and 2, it only applies to one coordinate of one distribution at a time. We are interested in the parametric versions of the inference problems. They can be phrased as follows:

3. Find all parameter values for a model that result in the same values for all .

4. Given observations Y = σ and hidden data X = h, identify all parameter values such that h is the most likely explanation for the observations σ.

   As we will show, the following modeling problems are related fundamentally to problems 3 and 4:

5. Which (parameter-independent) relations on the probabilities does the model imply?

6. Describe the tropicalization of the variety that corresponds to a graphical model.

Problem 5 asks for the ideal of polynomial invariants of a statistical model (1). Invariants have been investigated in phylogenetics (8, 9), where they can help to identify good trees for aligned DNA sequences.

The primary goal of our study is to give a practical answer to problem 4 for graphical models. Our main algorithmic result is an efficient procedure for parametric inference that can be viewed as a polytopal analog of the sum-product algorithm. The efficiency is based on the complexity estimates for Newton polytopes that we derive in section 4. The resulting polytope propagation algorithm is applied to problems in biological sequence analysis in ref. 10.

The mathematics developed in sections 3 and 4 is of independent interest. It also furnishes tools for parametric inference (problems 3 and 4) and parametric modeling (problems 5 and 6), which are applicable to a wide range of statistical problems. We demonstrate this point of view by analyzing the hidden Markov model (HMM) and the general Markov model on a binary tree in sections 2 and 5, respectively.

2. Algebraic Representation of HMMs

A graphical model is an algebraic variety that is presented as the image of a highly structured polynomial map f: R ^d → R ^m. Here, R ^d is the space in which coordinates are the model parameters s ₁,..., s _d, and R ^m is the space in which coordinates are the joint probabilities for the observed random variables. In applications, the integer m is much larger than the integer d; in fact, it is so large that one can only look at one coordinate p _σ at a time. Each coordinate f _σ = f _σ(s ₁,..., s _d) of the map f is a polynomial function in s ₁,..., s _d. The efficient evaluation of these functions relies on the sum-product algorithm. Here, we study the (parametric) inference and modeling problems in the familiar context of the HMM.

A discrete HMM has n observed states Y ₁,...,Y _n taking on l possible values and n hidden states X ₁,..., X _n taking on k possible values. The HMM can be characterized by the following conditional independence statements for i = 1,..., n:

We consider the homogeneous model with uniform initial distribution, where all transitions X _i → X _i+1 are given by the same k × k matrix S = (s _ij) and all transitions X _i → Y _i are given by the same k × l matrix T = (t _ij). Throughout our discussion, we disregard for simplicity the usual probabilistic hypothesis that S and T are nonnegative and that all row sums are 1.

Proposition 1. The HMM is the image of a map f: R ^d → R ^ln, where d = k(k + 1) and each coordinate of f is a bihomogeneous polynomial of degree n - 1 in S and degree n in T.

Problem 3 is to compute the fibers of the map f. In statistics, this computation is called parameter identification. We use the term coordinate polynomials for the polynomials f _σ that are coordinates of the map f.

Our running example in this section is the case n = 3 with binary random variables (k = l = 2). The graph of this model is given in Fig. 1. The shaded nodes are the observed random variables.

Here, the parameter space is R ⁸ with the coordinates s ₀₀, s ₀₁, s ₁₀, s ₁₁, t ₀₀, t ₀₁, t ₁₀, and t ₁₁, and it maps to R ⁸ with the coordinates p ₀₀₀, p ₀₀₁, p ₀₁₀, p ₀₁₁, p ₁₀₀, p ₁₀₁, p ₁₁₀, and p ₁₁₁. The map f: R ⁸ → R ⁸ is given by the following:

The HMM (i.e., the image of f) is the zero set of the following quartic polynomial:

This polynomial was found by a Gröbner basis computation. See the discussion on implicitization in section 3 of ref. 11.

In general, the polynomial functions on R ^ln that vanish on the image of f are called the invariants of the model. They form a prime ideal I _f. In our example, I _f is generated by the quartic polynomial given above. Problem 5 is to compute generators of the ideal I _f. When l ⁿ and d are small, this can be done by using Gröbner bases, and in some cases it is possible to characterize I _f based on the structure of the model (see, for example, Conjecture 13), but in general, problem 5 is hard and the ideal I _f may remain unknown.

Here, tropical geometry comes in. The tropicalization of our map f is the map g: R ^d → R ^ln defined by replacing products by sums and sums by minima in the formula for f. In our example (n = 3, k = l = 2), the tropicalization is the piecewise-linear map g: R ⁸ → R ⁸, (U, V) ↦ with the following:

This minimum is attained by the most likely hidden data (ĥ ₁ ĥ₂ ĥ ₃), given the observations (σ₁, σ₂, σ₃) and given the parameters and . The sequence (ĥ ₁, ĥ ₂, ĥ ₃) is known as the Viterbi sequence in the HMM literature (12). It solves problem 2 in section 1.

The key observation, which we discuss in more detail in section 4, is that the set of parameters (U, V) that selects the Viterbi sequence (ĥ ₁, ĥ ₂, ĥ ₃) is the normal cone at a vertex of the Newton polytope of the polynomial f _σ1σ2σ3. This polytope is four-dimensional, it has eight vertices, and its normal fan represents the solution to problem 4 in section 1 when σ = σ₁σ₂σ₃ is fixed.

We can also consider an extension of problem 4 in which σ = σ₁σ₂σ₃ ranges over all possible observations. The solution is given by the Newton polytope of the map f. In our example, it is a five-dimensional polytope with 398 vertices, 1,136 edges, 1,150 two-faces, 478 three-faces, and 68 facets, namely, the Minkowski sum of eight copies of the earlier four-dimensional polytope for (σ₁, σ₂, σ₃) ∈ {0, 1}³. For a concrete numerical example, fix the parameters and ). We find the following:

If the observed string at Y ₁ Y ₂ Y ₃ is

then the Viterbi sequence at X ₁ X ₂ X ₃ is

The set of all parameters (U, V) leading to the same conclusions as (U ^*, V ^*) is the cone defined by the following:

Our solution to the parametric inference problem with respect to all observations simultaneously consists of 398 such cones. The tropical HMM is the union of the images of these cones under the piecewise-linear map g: (U, V) ↦ δ. This image is a piecewise-linear set of dimension 7. The cone that contains the chosen parameters (U ^*, V ^*) is mapped to a seven-dimensional cone in the tropical HMM (it spans the hyperplane δ₀₁₀ = δ₁₀₀), but most of the other 397 cones are mapped to lower-dimensional cones by the map g. The question of how the number 398 grows as the length n increases is addressed in Corollary 10.

3. Positivity and Morphisms in Tropical Geometry

We have shown that a graphical model is the image of a polynomial map f from the space of parameters to the space of joint probability distributions on the observed random variables. Furthermore, we have shown that the tropicalization of f arises naturally in solving problem 4. In this section, we study the geometry of tropicalization in the more general setting where f: R ^d → R ^m is an arbitrary polynomial map. In statistical applications, each coordinate f _σ of the map f is usually a polynomial with positive coefficients. If this condition holds, then the polynomial map f is called positive. We consider f to be surjectively positive if, in addition, f maps the positive orthant surjectively onto the positive points in the image, in symbols,

The set of all polynomial functions that vanish on the image of f is a prime ideal I _f in the polynomial ring R[p ₁,..., p _m]. The closure of the image of f is the variety of the prime ideal I _f.

In tropical geometry, we replace the variety of I _f by a piecewise-linear set as follows. The tropical variety is the set of all weight vectors w ∈ R ^m such that the initial ideal in_w(I _f) contains no monomial (4, 13). By following ref. 14, we define the positive tropical variety as the set of all weight vectors w ∈ R ^m such that the initial ideal in_w(I _f) contains no polynomial with only positive coefficients. The tropical variety is a polyhedral fan in R ^m, and is a polyhedral subcomplex of . This observation means that is a finite union of closed convex polyhedral cones that fit together nicely, and is the union of a subset of these cones. The tropicalization of the polynomial map f is the piecewise-linear map g: R ^d → R ^m defined by replacing products by sums and sums by minima in the evaluation of f. We consider g to be a tropical morphism. Examples of tropical morphisms appear in Eqs. 1, 3, 4, 9, and 10.

The following theorem describes the geometry of this situation. We define the Newton polytope of a polynomial map f: R ^d → R ^m as the Minkowski sum in R ^d of the Newton polytopes of its coordinates f ₁,..., f _m. For basic information on Newton polytopes and their normal fans, see section 1 of ref. 13.

Theorem 2. The tropical morphism g is linear on each cone in the normal fan of the Newton polytope of f. Its image is a fan contained in . If f is positive, then image(g) is a subset of , but it is generally not a polyhedral subcomplex. If f is surjectively positive, then .

Proof: Let P _i denote the Newton polytope of the polynomial f _i = f _i(s ₁,..., s _d). By definition, P _i is the convex hull in R ^d of all nonnegative lattice points a = (a ₁,..., a _d) ∈ N ^d such that the monomial appears with a nonzero coefficient in f _i. The piecewise-linear concave function g _i is the support function of the polytope P _i. Thus, g _i(w) is the minimum value attained on P _i by the linear functional a ↦ w·a. In particular, the function g _i: R ^d → R is linear on each cone in the normal fan of P _i.

The Newton polytope of the map f is the Minkowski sum P ₁ +· · ·+ P _m = {a ₁ +· · ·+ a _m: a _i ∈ P _i}. The normal fan of P ₁ +· · ·+ P _m is the common refinement of the normal fans of P ₁,..., P _m. This observation shows that the function f = (f ₁,..., f _m): R ^d → R ^d is linear on each cone of the normal fan of the Newton polytope of f. Because g is continuous, the image of g is a closed polyhedral fan in R ^m.

Consider any vector w ∈ R ^d. We must show that g(w) lies in , and if f is positive, then g(w) lies in . Let φ be any polynomial in the ideal I _f. If we substitute p ₁ = f ₁,..., p _m = f _m into φ = φ(p ₁,..., p _m), then the result is zero. Consequently, if we substitute the initial forms p ₁ = in_w(f ₁),..., p _m = in_w(f _m) into the initial form in_g(w) (φ), then the result is zero (see equation 11.2 in ref. 13, p. 100). This fact implies that in_g(w) (φ) is not a monomial. Moreover, if f is positive, then φ must have two terms whose coefficients have opposite signs.

The following example shows that image (g) need not be a subcomplex of . If f is assumed to be surjectively positive, then it follows from proposition 2.5 in ref. 14 that .

Example 3: Let d = 3 and m = 4, and consider the linear map

Then, I _f is the principal ideal generated by the linear form p ₁ - p ₂ + p ₃ - p ₄, and is essentially the normal fan of a tetrahedron. We identify with the complete graph K ₄. The six edges of K ₄ are labeled with six monomial-free initial ideals of I _f, namely, 〈p ₁ + p ₃〉, 〈-p ₂ - p ₄〉, 〈p ₁ - p ₂〉, 〈p ₁ - p ₄〉, 〈-p ₂ + p ₃〉, 〈p ₃ - p ₄〉. The first two of these six initial ideals contain a polynomial with positive coefficients. Hence, the positive tropical variety is the four-cycle in K ₄ formed by the remaining four edges.

The tropicalization of the linear map f is the tropical morphism:

The image of g is the set of all vectors (a, a, b, c) with a ≤ b ≤ c. Each vector (a, a, b, c) with a < b < c has the initial ideal 〈p ₁ - p ₂〉, so it lies on a particular edge of K ₄. But the same edge also accounts for all vectors (a, a, b, c) with a < c < b, none of which is in the image of g. Thus, image(g) is a closed segment that covers only half of the edge of K ₄ indexed by 〈p ₁ - p ₂〉.

In the rest of this section, we examine Theorem 2 for a small but important graphical model, namely, the naive Bayes model with two features (ref. 1, section 7). There are two observed random variables Y ₁ and Y ₂ that depend on one hidden binary random variable X. The two observed variables take k and l possible values, respectively. The parameterization f of this model is the map f: R ^2(k+l) ↦ R ^kl given by p _ij = s _i0 t _0j + s _i1 t _1j. Thus, the model consists of all k × l matrices P = (p _ij) of the form P = S·T, where S is a k × 2 matrix and T is a 2 × l matrix (i.e., the model consists of precisely the k × l matrices of rank ≤2).

Proposition 4. The parameterization f of the naive Bayes model with two features is surjectively positive. The ideal I _f is generated by the 3 × 3 subdeterminants of the k × l matrix P = (p _ij).

Proof: The map f being positive means that, if P is any positive matrix of rank 2, then S and T can be chosen to be positive. This is a known result in linear algebra (e.g, see ref. 15). The same statement is false for rank ≥3 (i.e., the parameterization of the naive Bayes model with three or more features is not surjectively positive). A well known result in commutative algebra states that the (r + 1) × (r + 1) minors of a k × l matrix generate a prime ideal. The variety of this ideal is the set of k × l matrices of rank ≤r. This is our ideal I _f for r = 2.

The objects of Theorem 2 have been studied (3, 16). The tropical variety is the set of k × l matrices of tropical rank ≤2, and the tropical variety is the set of k × l matrices of Barvinok rank ≤2. Develin (16) determined the combinatorics and topology of these spaces when min(k, l) = 3. He showed that is shellable but that can have torsion in its integral homology groups.

The Newton polytope of the map f is an interesting combinatorial object, namely, it is the (kl - k - l + 2)-dimensional zonotope associated with the complete bipartite graph K _k,l. The Newton polytope of each coordinate f _ij is a line segment, and the zonotope is their Minkowski sum. The normal fan is the hyperplane arrangement {u _i0 - u _i1 = v _1j - v _0j}. Its maximal cones correspond to the acyclic orientations of the complete bipartite graph K _k,l. West (17) showed that the number of facets of such a cone can be any integer between k + l - 1 and kl. The total number of cones equals , where S(k, i) is the Stirling number of the second kind. Here, the tropical morphism g is given by the following:

The map g: R ^2(k+l) ↦ R ^kl is piecewise linear with respect to the hyperplane arrangement.

Example 5: Let k = l = 3, so the two observed random variables are ternary. The prime ideal is the following:

The tropical variety is the fan over a two-dimensional polyhedral complex consisting of six triangles and nine quadrangles. This complex is the 2-skeleton of the product of two triangles, labeled as in Fig. 2a . This complex is shellable. The positive tropical variety is the subcomplex consisting of the nine quadrangles shown in Fig. 2b . Note that is a torus.

The Newton polytope of f is a five-dimensional zonotope with 230 vertices, one for each acyclic orientation of the complete bipartite graph K _3,3. The map g is linear on each of the 230 cones in the corresponding hyperplane arrangement, but it is rank-deficient on 68 of the cones. The remaining 162 = 18 × 9 cones are mapped onto the nine quadrangles of the torus . Thus, the general fiber of g involves 18 cones. Of these 18 cones, 8 cones have five facets, 8 cones have six facets, and 2 cones have eight facets.

4. Newton Polytopes of Graphical Models and their Complexity

Consider a graphical model with E edges and n observed random variables Y ₁,...,Y _n, each taking l values. Such a model is given by a positive polynomial map f: R ^d → R ^ln. Each coordinate f _σ of f is a polynomial of degree e in the model parameters s ₁,..., s _d. In this section, we discuss the statistical meaning and the computational complexity of the mathematical objects introduced in section 3.

We write for the negative logarithms of the model parameters. Consider any of the l ⁿ possible observations σ. The quantity f _σ(s ₁,..., s _d) is the probability of making this particular observation [i.e., it is Prob(Y = σ)]. The quantity g _σ(u ₁,..., u _d) is the negative logarithm of the conditional probability Prob(X = h|Y = σ), where h maximizes Prob(X = ĥ|Y = σ) for the parameters (s ₁,..., s _d). Clearly, the function g _σ: R ^d → R is piecewise linear and concave on the logarithmic parameter space.

The domains of linearity of the function g _σ are the cones in the normal fan of the Newton polytope of f _σ. Each maximal cone C is indexed by the hidden data ĥ that maximizes Prob(X = h|Y = σ) for any of the parameters (u ₁,..., u _d) ∈ C. The hidden data ĥ that arise in this manner, for some choice of logarithmic parameters u, are called the possible explanations of the observation σ. For example, for the HMM described in section 2, the explanations are the Viterbi sequences.

We vary the observations as follows. Each logarithmic parameter vector u defines an inference function σ ↦ ĥ from the set of observations to the set of explanations. For the HMM, each inference function {1,..., l}ⁿ → {1,..., k}ⁿ takes an observed sequence σ to the corresponding Viterbi sequence ĥ. There are such functions, but most of these are not inference functions. For example, consider the binary HMM of length three. There are 8⁸ = 16,777,216 Boolean functions {0, 1}³ → {0, 1}³, but as we show at the end of section 2, only 398 of these are inference functions for the HMM.

Proposition 6. The inference functions σ ↦ ĥ of a graphical model f are in bijection with the vertices of the Newton polytope of the map f. The explanations ĥ for a fixed observation σ in a graphical model are in bijection with the vertices of the Newton polytope of the polynomial f _σ.

In applications of graphical models, the number d of parameters and the number l of values of the observed random variables are small and fixed, but the number n of observed random variables is large. Recall that the model is the image of the map f: R ^d → R ^ln. Hence, the dimension of the model remains fixed, but the dimension of its ambient space grows exponentially in n. Therefore, it is algorithmically infeasible to compute the full tropical variety . However, we can efficiently compute the Newton polytopes of the f _σ, or even the Newton polytope of f. This insight allows us to glean information about the tropical variety from the domains of linearity of its “coordinate functions” g _σ.

Our next goal is to derive an upper bound on the number of vertices of the Newton polytopes.

Theorem 7. Consider graphical models f whose number of parameters d is fixed and whose number n of observed random variables and number of edges E varies. (Typically, E is a linear function of n.) Then, the number of vertices of the Newton polytope NP(f _σ) of f _σ is bounded above by the following:

For many important families of graphical models, the number E of edges is bounded by a linear function in terms of the number n of observed nodes, and in these cases, we can replace E by n. Hence, for any given observation σ, the number of explanations grows polynomially in n. For example, in the HMM described in section 2, we have E = 2n - 1, and a similar relationship holds in the tree model of section 5.

Corollary 8. For any fixed observation in the homogeneous HMM, the number of explanations is at most C _k,l·n ^k(k+l) . If all random variables are binary, then the upper bound C·n ^10/3 holds.

The proofs of Theorem 7 and Corollary 8 are derived from the following classical result on lattice polytopes by Andrews (18). The necessary observation is that the Newton polytope of f _σ is contained in the cube [0, E]^d, and the volume of this cube equals E ^d.

Proposition 9. For every fixed integer d, there exists a constant C _d such that the number of vertices of any lattice polytope P in R ^d is bounded above by C _d·volume(P)^(d-1)/(d+1) (18).

The Newton polytope of the map f was defined as the Minkowski sum of the l ⁿ smaller Newton polytopes in Theorem 7. We infer the following to be naive bound on its number of vertices.

Corollary 10. The number of inference functions of a graphical model is at most l^nC_dE^{d -1} ; hence, this number scales at most singly exponentially in the complexity (n, E) of the graphical model.

Consider the homogeneous HMM on binary random variables. Each inference function is a Boolean function {0, 1}ⁿ → {0, 1}ⁿ, but not conversely. The number of all Boolean functions is , which grows doubly exponentially in n. However, the number of inference functions is at most 2^{polynomial(n)}.

In practical applications of graphical models, it may be infeasible to compute all (singly exponentially many) inference functions. Nonetheless, we believe that important insight can be gained by computing and classifying the Newton polytopes of graphical models f on few random variables. Such a study would be the polyhedral analog to the algebraic classification of ref. 1.

However, for a fixed observation σ, the size of the Newton polytope of f _σ grows polynomially with the size of the graphical model, and therefore, there is hope that the polytopes can be computed efficiently. Despite the fact that the Newton polytope of f _σ has polynomially many vertices in the size of the graphical model, the number of terms in f _σ grows exponentially. This is a potential problem because the computation of the Newton polytope requires the inspection of these terms. The following result states that the convex hull computations scale with the running time of the sum-product algorithm, which for many models of interest scales polynomially with the size of the graphical model.

Proposition 11 (Polytope Propagation). The Newton polytopes of the polynomials f _σ can be computed recursively by using the decomposition of f _σ according to the sum-product algorithm.

5. The General Markov Model on a Binary Tree

We conclude by illustrating the concepts that we have developed in the context of tree Markov models. These models are directed graphical models in which the graph is a directed tree τ with observed random variables Y ₁,...,Y _n at the leaves. The naive Bayes model in section 3 is the special case in which n = 2. Each edge e has a different transition matrix . We consider the general model given by Allman and Rhodes (8), which means that the S ^e are arbitrary distinct l × l matrices. In most applications, the transition matrices are from a special model family (for example, in phylogenetics, these may be Jukes-Cantor model or the Hasegawa-Kishino-Yano model). As before, we relax the hypothesis that transition probabilities are nonnegative and sum to 1. Hence, the are distinct unknowns. For simplicity, we further assume that the tree τ is binary.

Proposition 12. The general Markov model for the binary tree τ is the image of a map , where each coordinate of f is a multilinear polynomial in the unknowns {(), e edge of τ}.

If we denote an edge between nodes i and j by (ij), and τ′ is the tree τ without the leaves, then the coordinate of the multilinear map f indexed by an observed sequence (σ₁,...,σ_n) can be written as follows:

Here, h ranges over all colorations h = (h _i)_i∈τ of the nodes such that h _j = σ_j for all leaves j. Our running example in this section is the binary tree in Fig. 3 with binary random variables (l = 2).

In this example, the coordinates of the multilinear map f: R ²⁴ → R ¹⁶ are given by the following formula:

The prime ideal I _f of polynomial invariants is generated by the 3 × 3 subdeterminants of the following matrix:

Thus, this particular model is the k = l = 4 instance of the determinantal variety in Proposition 4.

We generalize the determinantal presentation in this example by proposing the following explicit solution to problem 5 for arbitrary binary trees τ. Every edge of τ induces a split of the set of leaves {1, 2,..., n}, corresponding to the two connected components of the tree obtained by removing that edge. The unrooted tree underlying τ is uniquely determined by the set of these splits.

Conjecture 13. The ideal I _f of phylogenetic invariants of the general Markov model for any binary tree τ on binary random variables is generated by the 3 × 3-determinants of all two-dimensional matrices obtained by flattening the 2 × · · · × 2-table (p _σ1..._σn) according to the splits induced by the edges of τ.

We need to explain the meaning of the word “flattening.” If (A, B) is any split of the set {1,..., n}, then this term refers to the 2^#(A) × 2^#(B) matrix whose rows and columns are indexed by the functions A → {0, 1} and B → {0, 1}, respectively, and whose entries are the 2ⁿ probabilities p _σ1..._σn.

The sum-product algorithm is used in practice to evaluate the polynomial of Eq. 5. Its running time is linear in n, despite the fact that the number l ^n-1 of terms in Eq. 5 grows exponentially. This reduction in complexity is achieved by recursively grouping subsums. For example, Eq. 6 becomes the following:

Remember the following rule: Polynomials are evaluated recursively as sums of products of smaller polynomials. This is the solution to problem 1. For details on the tree case, see ref. 19.

Problem 2 is known in phylogeny as the joint ancestral reconstruction problem, which asks for the MAP ancestral assignments ĥ _i given the observations (σ₁,..., σ_n) at the leaves. An efficient method for solving this problem is given in ref. 20. This method is nothing but the sum-product algorithm with ordinary arithmetic (+, ×) replaced by tropical arithmetic (min, +). The σ coordinate of the tropicalization of the map in Eq. 5 is

This expression can be evaluated efficiently by the same scheme as used previously. The rule is now the following: Piecewise-linear concave functions are evaluated recursively as minima of sums of smaller such functions. A simple example illustrating this rule is the following tropicalization of Eq. 8:

where and similarly for u _vσ3σ4.

In section 4, we showed that the number of vertices of the Newton polytopes of the coordinate polynomials f _σ is critical for efficient parametric inference. That number grows polynomially in n if the number of parameters is fixed (because of Theorem 7), but it may grow exponentially if the number of parameters is not bounded. For the general Markov model on a tree τ, the growth will be exponential unless we restrict the number of parameters. This can be done, for example, by considering the homogeneous tree model as follows, where the transition matrices along all edges are identical: is independent of the edge e. By using Theorem 7, we obtain the following result analogous to Corollary 8.

Proposition 14. The number of vertices of the Newton polytope of any coordinate f _σ in the homogeneous tree model is bounded above by times a constant depending only on l.

For tree models that are used in applications, such as phylogenetics, the number of parameters is likely to be reduced even further. In such cases, the parametric joint ancestral reconstruction problem can be solved efficiently by using the polytope propagation algorithm techniques given in Proposition 11.

6. Summary: A Statistics-Geometry Dictionary

The algebraic representation for graphical models with hidden variables leads naturally to an interpretation of a parameterized model as a point on an algebraic variety. Marginal probabilities are coordinates of points on the variety. Varieties can be tropicalized, and the statistical meaning is that the MAP probabilities (calculated with logarithms of the parameters) can be interpreted as coordinates of points on the positive part of the tropical variety. Hence, the tropical model is fundamental for understanding MAP probabilities. Although we have not addressed it in this article, the logarithms of the marginal probabilities are coordinates of points on the amoeba (21) of the model. Amoebas are likely to be important for understanding the geometry of maximum-likelihood estimation.

The sum-product algorithm for graphical models is an efficient method for evaluating the coordinate polynomials of a graphical model. This algorithm works in exactly the same way for classical arithmetic (+, ×) and for tropical arithmetic (min, +). The same method is used to evaluate coordinates of points on the variety and of points on the tropical variety.

An explanation for an observation σ is a vertex of the Newton polytope of f _σ. Thus, the parametric inference problem is solved by finding the normal fans of the Newton polytopes of the coordinate polynomials. For many important applications, the number of vertices of the polytopes is polynomial in the size of the graphical model. The polytope propagation algorithm, which is a geometric analog of the sum-product algorithm, finds the Newton polytopes and is efficient when the sum-product algorithm is fast and the number of vertices on the Newton polytopes is small.

In our companion article (10), we show that polytope propagation is practical and useful in the important application of biological sequence analysis. In particular, existing parametric alignment methods (22-24) can be viewed as special cases of parametric inference for a pair HMM. The computation of the Newton polytopes is also useful for Bayesian computations, where we have priors on the parameters and it is of interest to integrate over the maximal cones in the normal fan of the Newton polytope (ref. 10, section 5).