The Maximum Likelihood Degree of Toric Varieties

TL;DR

This paper investigates the maximum likelihood degree of toric varieties, revealing how scaling affects it and providing methods for numerical computation and examples from algebraic geometry and statistics.

Contribution

It establishes the relationship between ML degree and the degree of toric varieties under generic scalings, and explores the impact of the principal A-determinant locus.

Findings

ML degree equals the degree of the toric variety for generic scalings.

ML degree drops when the scaling vector lies in the principal A-determinant locus.

Numerical methods for ML estimation via homotopy continuation are demonstrated.

Abstract

We study the maximum likelihood degree (ML degree) of toric varieties, known as discrete exponential models in statistics. By introducing scaling coefficients to the monomial parameterization of the toric variety, one can change the ML degree. We show that the ML degree is equal to the degree of the toric variety for generic scalings, while it drops if and only if the scaling vector is in the locus of the principal $A$-determinant. We also illustrate how to compute the ML estimate of a toric variety numerically via homotopy continuation from a scaled toric variety with low ML degree. Throughout, we include examples motivated by algebraic geometry and statistics. We compute the ML degree of rational normal scrolls and a large class of Veronese-type varieties. In addition, we investigate the ML degree of scaled Segre varieties, hierarchical loglinear models, and graphical models.

Tables

Table 1. Table 1 . The ML degree of ( Ver ​ ( 2 , 2 ) c ) Ver superscript 2 2 𝑐 ({\rm Ver}(2,2)^{c}) for different scalings c i ​ j subscript 𝑐 𝑖 𝑗 c_{ij} in the matrix C 𝐶 C .

$C$	${\Delta }_A$	${\Delta }_{\left[a_1{a}_4{a}_5\right]}$	${\Delta }_{\left[a_3{a}_5{a}_6\right]}$	${\Delta }_{\left[a_1{a}_2{a}_2\right]}$	$\mathrm{mldeg}$
	$\ne 0$	$\ne 0$	$\ne 0$	$\ne 0$	4
	0	0	$\ne 0$	$\ne 0$	3
	0	0	0	$\ne 0$	2
	$\ne 0$	0	0	0	1
	0	$\ne 0$	$\ne 0$	$\ne 0$	3
	0	$\ne 0$	0	$\ne 0$	2
	0	0	0	0	1

Equations

$$p_1,\dots ,p_n\ge 0\text{mbox}and{p}_1+\dots +p_n=1$$

$$p_i=\mathrm{prob}(X=i)\text{mbox}fori=1,\dots ,n.$$

$$u_i=|\{j:X_j=i\}|\text{ and }N=u_1+\cdots +u_n.$$

$${\ell }_u(p)=\frac{p_1^{u_1}\cdots p_n^{u_n}}{(p_1+\cdots +p_n{)}^{u_1+\cdots +u_n}},$$

$${\psi }^c(s,{\theta }_1,{\theta }_2,\dots ,{\theta }_{d-1})=(c_{1}s{\theta }^{a_1},c_{2}s{\theta }^{a_2},\dots ,c_{n}s{\theta }^{a_n}),$$

$$p_{ijk\ell }=\mathrm{prob}(S=i,B=j,H=k,L=\ell )$$

$$\hat{p}=($$

$$0.12976986,0.076165372,0.020853199,0.015471205,$$

$$0.13533793,0.11789409,0.018820142,0.0207235,$$

$$0.083859917,0.073051125,0.01347576,0.014838619).$$

$${\psi }^c(s,{\theta }_1,{\theta }_2,\dots ,{\theta }_{d-1})=(c_{1}s{\theta }^{a_1},c_{2}s{\theta }^{a_2},\dots ,c_{n}s{\theta }^{a_n}).$$

$$(s,t,u)\mapsto (s,2st,s{t}^2,2su,2stu,s{u}^2).$$

$$f=\sum _{i=1}^n{c}_i{\theta }^{a_i}=c_1{\theta }^{a_1}+\cdots +c_n{\theta }^{a_n}.$$

$$1=sf$$

$${\left(Au\right)}_i=u_{+}s{\theta }_i\frac{\partial f}{\partial {\theta }_i}\text{ for }i=1,\dots ,d-1.$$

$$\begin{array}[]{ccc}1&=&sf\\[5.0pt] \left(Au\right)_{1}&=&u_{+}s\theta_{1}\frac{\partial f}{\partial\theta_{1}}\\[5.0pt] \left(Au\right)_{2}&=&u_{+}s\theta_{2}\frac{\partial f}{\partial\theta_{2}}\\ &\vdots\\ \left(Au\right)_{d-1}&=&u_{+}s\theta_{d-1}\frac{\partial f}{\partial\theta_{d-1}}\end{array}$$

$$\begin{array}[]{ccc}(Au)_{1}f&=&u_{+}\theta_{1}\frac{\partial f}{\partial\theta_{1}}\\[5.0pt] (Au)_{2}f&=&u_{+}\theta_{2}\frac{\partial f}{\partial\theta_{2}}\\ &\vdots\\ (Au)_{d-1}f&=&u_{+}\theta_{d-1}\frac{\partial f}{\partial\theta_{d-1}}\end{array}$$

$$\begin{array}[]{ccc}(Au)_{1}L_{f}(p)&=&u_{+}L_{1}(p)\\ (Au)_{2}L_{f}(p)&=&u_{+}L_{2}(p)\\ &\vdots\\ (Au)_{d-1}L_{f}(p)&=&u_{+}L_{d-1}(p)\end{array}$$

$${\nabla }_A=\overline{\left\{c\in {\left({\mathbb{C}}^{\ast }\right)}^n\mid \exists \theta \in {\left({\mathbb{C}}^{\ast }\right)}^{d-1}\text{ such that }f\left(\theta \right)=\frac{\partial f}{\partial {\theta }_i}\left(\theta \right)=0\text{ for all }i\right\}}.$$

$$E_A(f)=R_A\left(f,{\theta }_1\frac{\partial f}{\partial {\theta }_1},\dots ,{\theta }_{d-1}\frac{\partial f}{\partial {\theta }_{d-1}}\right).$$

$$E_A(f)=\prod _{\Gamma \text{ face of }Q}{\Delta }_{\Gamma \cap A}$$

$$\begin{array}[]{ccc}1&=&L_{f_{\Gamma}}\\ (Au)_{1}&=&u_{+}L_{1}^{\prime}\\ (Au)_{2}&=&u_{+}L_{2}^{\prime}\\ &\vdots\\ (Au)_{d-1}&=&u_{+}L_{d-1}^{\prime}.\end{array}$$

$$p_1^{w_1}{p}_2^{w_2}\cdots p_{\ell }^{w_{\ell }}-p_{\ell +1}^{-w_{\ell +1}}\cdots p_{d+1}^{-w_{d+1}}.$$

$${\Delta }_A=\left(w_{\ell +1}^{-w_{\ell +1}}\cdots w_{d+1}^{-w_{d+1}}\right)c_1^{w_1}\cdots c_{\ell }^{w_{\ell }}-\left(w_1^{w_1}\cdots w_{\ell }^{w_{\ell }}\right)c_{\ell +1}^{-w_{\ell +1}}\cdots c_{d+1}^{-w_{d+1}}.$$

$$c\in V\left(\left(w_{\ell +1}^{-w_{\ell +1}}\cdots w_{d+1}^{-w_{d+1}}\right)c_1^{w_1}\cdots c_{\ell }^{w_{\ell }}-\left(w_1^{w_1}\cdots w_{\ell }^{w_{\ell }}\right)c_{\ell +1}^{-w_{\ell +1}}\cdots c_{d+1}^{-w_{d+1}}\right)$$

$$A=\left(\begin{array}[]{ccccccccccccc}1&\cdots&1&0&\cdots&0&\cdots&0&\cdots&0&0&\cdots&0\\ 0&\cdots&0&1&\cdots&1&\cdots&0&\cdots&0&0&\cdots&0\\ \vdots&&&&&&&&&&&&\vdots\\ 0&\cdots&0&0&\cdots&0&\cdots&1&\cdots&1&0&\cdots&0\\ 0&1\cdots&n_{1}&0&1\cdots&n_{2}&\cdots&0&1\cdots&n_{d-2}&0&1\cdots&n_{d-1}\end{array}\right).$$

$$g_i=c_{i0}+c_{i1}{\theta }_{d-1}+c_{i2}{\theta }_{d-1}^2+\cdots +c_{i{n}_i}{\theta }_{d-1}^{n_i}.$$

$$f=({\theta }_1{g}_1+{\theta }_2{g}_2+\cdots +{\theta }_{d-2}{g}_{d-2})+g_{d-1}.$$

$$(Au{)}_{d-1}f=u_{+}{\theta }_{d-1}\left(g_{d-1}^{\prime }+\sum _{i=1}^{d-2}{\theta }_i{g}_i^{\prime }\right).$$

$$(Au{)}_{d-1}{g}_1\cdots g_{d-1}={\theta }_{d-1}\left(u_{(d-1)+}{g}_1\cdots g_{d-2}{g}_{d-1}^{\prime }+\sum _{i=1}^{d-2}(Au{)}_i{g}_1\cdots g_{i-1}{g}_i^{\prime }{g}_{i+1}\cdots g_{d-1}\right).$$

$$c_0+c_1\theta +c_2{\theta }^2+\cdots +c_n{\theta }^n.$$

$$\left(\begin{array}[]{ccccc}c_{1}&4c_{2}&3c_{3}&8c_{4}\\ 2c_{0}&3c_{1}&c_{2}&c_{3}\end{array}\right).$$

Full text

This project was started during the AMS Mathematics Research Communities 2016 program in Algebraic Statistics held in Snowbird, Utah. The authors were part of the Likelihood Geometry group led by Serkan Hoşten and Jose Rodriguez. This material is based upon work supported by the National Science Foundation under Grant Number DMS 1321794.

1. Introduction

Maximum likelihood estimation is a fundamental optimization problem in statistical inference. The maximum likelihood degree of an algebraic statistical model was introduced in [CHKS06] and [HKS05] to study the geometry and complexity of this optimization problem. Here we study the maximum likelihood degree of toric varieties.

Let $V\subset\mathbb{P}^{n-1}$ be a projective variety over $\mathbb{C}$. The homogeneous coordinates of $\mathbb{P}^{n-1}$ will be denoted by $[p_{1}\,:\,\cdots\,:\,p_{n}]$. The affine open subset of $V$ where $p_{1}+\ldots+p_{n}\neq 0$ is identified with the set of points in $V$ where $p_{1}+\ldots+p_{n}=1$. We can view the points $p\in V$ satisfying

[TABLE]

as a family of probability distributions of a discrete random variable $X$ where

[TABLE]

After observing $X$ in $N$ instances $X_{1},\ldots,X_{N}$, we record the data vector $u=(u_{1},\ldots,u_{n})\in\leavevmode\nobreak\ \mathbb{Z}^{n}$ where $u_{i}$ is the number instances where $X_{j}=i$, i.e.,

[TABLE]

The likelihood function is the rational function

[TABLE]

and one seeks to find a probability distribution $\hat{p}=(\hat{p}_{1},\ldots,\hat{p}_{n})$ in $V$ which maximizes $\ell_{u}$. Such a probability distribution $\hat{p}$ is a maximum likelihood estimate, and $\hat{p}$ can be identified by computing all critical points of $\ell_{u}$ on $V$.

We let $V_{\text{sing}}$ be the singular locus of $V$ and $V_{\text{reg}}=V\setminus V_{\text{sing}}$ be the regular locus. We also define $\mathcal{H}$ as the hypersurface that is the union of the coordinate hyperplanes and the hyperplane defined by $p_{1}+\cdots+p_{n}=0$.

Definition 1**.**

The maximum likelihood degree of $V$, denoted $\mathrm{mldeg}(V)$, is the number of complex critical points of the likelihood function $\ell_{u}$ on $V_{\text{reg}}\setminus\mathcal{H}$ for generic vectors $u$.

In [CHKS06] and [HKS05], it was shown that $\mathrm{mldeg}(V)$ is well-defined. We refer the reader to [DSS09, Huh13, HRS14] for a glimpse of the growing body of work on aspects of $\mathrm{mldeg}(V)$. For an excellent recent survey we recommend [HS14].

In most applications, $V$ is the Zariski closure of the image of a polynomial map, also known as a parametric statistical model. A widely used subclass of such models is composed of those given by a monomial parametrization. In algebraic geometry, they are known as toric varieties [Ful93, CLS11], and in probability theory and statistics they are known as discrete exponential families [BD76]. In particular, hierarchical log-linear models [BFH07] and undirected graphical models on discrete random variables [Lau96] present examples in contemporary applications. Theorem 13 is our main result and is restated below.

Theorem 2** (Main result).**

Let $V^{c}\subset\mathbb{P}^{n-1}$ be the projective variety defined by the monomial parametrization $\psi^{c}\,:\,(\mathbb{C}^{*})^{d}\longrightarrow(\mathbb{C}^{*})^{n}$ where

[TABLE]

and $c\in(\mathbb{C}^{*})^{n}$ is fixed. Then $\mathrm{mldeg}(V^{c})<\deg(V)$ if and only if $c$ is in the principal $A$-determinant of the toric variety $V=V^{(1,\ldots,1)}$.

To set the stage, we start with an example chosen from the theory of graphical models. This is the smallest example of a graphical model that is not decomposable. Models that are decomposable have a unique maximum likelihood estimate that is a rational function of the data vector $u$. This is equivalent to $\mathrm{mldeg}(V)=1$ [GMS06].

Example 3** (Binary $4$-Cycle).**

This example is from a study to determine risk factors for coronary heart disease based on data collected in Czechoslovakia [R*+*81, EH85]. Six different factors affecting coronary heart disease were recorded for a sample of $1841$ workers employed in the Czech automotive industry. We will only use the following four factors: whether they smoked ($S$), whether their systolic blood pressure was less than $140$ mm ($B$), whether there was a family history of coronary heart disease ($H$), and whether the ratio of beta to alpha lipoproteins was less than $3$ ($L$). Each is a binary variable, and the joint random variable $X=(S,B,H,L)$ has a state space of cardinality $16$. The data set is summarized in Table 2.

We will use the graphical model known as the $4$-cycle (Figure 2). If we set

[TABLE]

the family of probability distributions for $X$ which factor according to this model can be described by the following monomial parametrization: let $a_{ij},b_{jk},c_{k\ell},d_{i\ell}$ be parameters for $i,j,k,\ell\in\{0,1\}$ and let $p_{ijk\ell}=a_{ij}b_{jk}c_{k\ell}d_{i\ell}$. The Zariski closure $V$ of the image of this parametrization is a toric variety. Every probability distribution for $X$ that lies on this toric variety satisfies certain independence statements which can be read from the graph. For instance, $S$ and $H$ are independent given $B$ and $L$. Similarly, $B$ and $L$ are independent given $S$ and $H$. The maximum likelihood estimate given by the worker data is

[TABLE]

The degree of $V$ is $64$ and $\mathrm{mldeg}(V)=13$. This was first computed in [GMS06, p. 1484] where the question of explaining the fact that $\mathrm{mldeg}(V)=13$ was raised.

The article is organized as follows. In Section 2 we define scaled toric varieties, introduce maximum likelihood degrees, and recall discriminants. In Section 3 we prove our main theorem (Theorem 13). In Sections 4-7 we compute the ML degree of rational normal scrolls and a large class of Veronese-type varieties. We investigate the ML degree of scaled Segre varieties, hierarchical log-linear models, and graphical models. In Section 8, we show that maximum likelihood estimates can be tracked between different scalings of a toric variety via homotopy continuation. We illustrate with computational experiments that homotopy continuation can be a better numerical method than iterative proportional scaling [DR72] for computing the maximum likelihood estimate of a discrete exponential model.

2. Preliminaries

In this section we introduce scaled toric varieties and their likelihood equations. We show that the ML degree of a scaled toric variety is bounded by the degree of the variety. Next, we present two classical objects in toric geometry: the $A$-discriminant and the principal $A$-determinant [GKZ94]. The principal $A$-determinant plays a crucial role in our main result Theorem 13.

2.1. Scaled toric varieties

In this paper, we study the maximum likelihood degree of projective toric varieties. Let $A$ be a $(d-1)\times n$ matrix with columns $a_{1},\ldots,a_{n}\in\mathbb{Z}^{d-1}$. The convex hull of the lattice points $a_{i}$ defines a polytope $Q=\operatorname{conv}(A)$.

Definition 4**.**

A toric variety scaled by $c\in(\mathbb{C}^{*})^{n}$, denoted $V^{c}$, is defined by the following map $\psi^{c}\,:\,(\mathbb{C}^{*})^{d}\longrightarrow(\mathbb{C}^{*})^{n}$ where $A$ is a full rank $(d-1)\times n$ integer matrix:

[TABLE]

Here $V^{c}$ is the projective variety in $\mathbb{P}^{n-1}$ corresponding to the affine cone that is the closure of $\psi^{c}((\mathbb{C}^{*})^{d})$ in $\mathbb{C}^{n}$.

Different choices of $c$ give different embeddings of isomorphic varieties $V^{c}$ in $\mathbb{P}^{n-1}$; their ML degrees may also differ. When $c=(1,1,\dots,1,1)$, $V^{c}$ is the standard toric variety of $A$ as in [CLS11].

Example 5**.**

Consider $\psi^{c}\,:\,(\mathbb{C}^{*})^{3}\longrightarrow(\mathbb{C}^{*})^{6}$ given by

[TABLE]

This gives the Veronese embedding of $\mathbb{P}^{2}$ into $\mathbb{P}^{5}$ scaled by $c=(1,2,1,2,2,1)$. The usual Veronese embedding has ML degree four whereas this scaled version has ML degree one.

2.2. Likelihood equations

We begin our discussion with the critical equations whose solutions are the critical points of the likelihood function for a toric variety. These equations are called likelihood equations.

Let $A$ be a $(d-1)\times n$ matrix as in Definition 4. We let

[TABLE]

Definition 6**.**

Let $u=\left(u_{1},\dots,u_{n}\right)$ and $u_{+}=\sum_{i}u_{i}$. Using the method of Lagrange multipliers, we obtain the likelihood equations for the variety $V^{c}$:

[TABLE]

[TABLE]

In other words,

[TABLE]

We remark that the solutions to the above system are exactly the solutions to

[TABLE]

where $f\neq 0$.

We formulate the ML degree of the scaled toric variety $V^{c}$ using the coordinates $p_{1},\ldots,p_{n}$ as well. Note that $f=c_{1}\theta^{a_{1}}+\cdots+c_{n}\theta^{a_{n}}$ induces a linear form $L_{f}(p)=c_{1}p_{1}+\cdots+c_{n}p_{n}$ on the toric variety $V$. Similarly, $\theta_{i}\frac{\partial f}{\partial\theta_{i}}$ induce corresponding linear forms $L_{i}(p)$ on $V$ for $i=1,\ldots,d-1$, as the monomials appearing in $\theta_{i}\frac{\partial f}{\partial\theta_{i}}$ are a subset of the monomials appearing in $f$. With this we immediately get the following.

Proposition 7**.**

The ML degree of $V^{c}$ is the number of solutions $p\in V\setminus V(p_{1}\cdots p_{n}L_{f}(p))$ to

[TABLE]

for generic vectors $u$.

Proof.

Definition 1 combined with (3) give the result by noting the fact that the singularities of any toric variety are contained in the union of coordinate hyperplanes. ∎

The following result was first proved in [HS14, Theorem 3.2].

Corollary 8**.**

For $c\in(\mathbb{C}^{*})^{n}$, the ML degree of $V^{c}$ is at most $\deg(V)$. Moreover, for generic $c$, $\mathrm{mldeg}(V^{c})=\deg(V)$.

Proof.

For generic $u$, the linear equations (4) define a linear subspace of codimension $d-1$. By Bertini’s theorem, the intersection of this linear subspace with $V$ is generically transverse. Now by Bezout’s theorem [Ful98, Propositon 8.4], the sum of the degrees of the components of the intersection is equal to the degree of $V$. Therefore, the degree of the zero dimensional piece of this intersection is at most $\deg(V)$. The ML degree of $V^{c}$ is obtained by removing those solutions in $V(p_{1}\cdots p_{n}L_{f}(p))$, hence $\mathrm{mldeg}(V^{c})\leq\deg(V)$. For generic $c$, the linear subspace is a generic linear subspace of codimension $d-1$. Not only will the intersection contain exactly $\deg(V)$ points, these points will not be in $V(p_{1}\cdots p_{n}L_{f}(p))$. Hence, in this case, $\mathrm{mldeg}(V^{c})=\deg(V)$. ∎

We will denote the set of solutions in $V$ to (4) by $\mathcal{L}_{c,u}^{\prime}$. Similarly, we will denote the set of solutions to (4) in $V\setminus V(p_{1}\cdots p_{n}L_{f}(p))$ by $\mathcal{L}_{c,u}$. To close this subsection, we would like to point out a classical result both in probability theory and toric geometry. It is sometimes known as Birch’s theorem (see e.g. [Lau96]).

Theorem 9** (Birch’s Theorem).**

There is a unique positive point in $\mathcal{L}_{c,u}(p)$, namely, there is a unique positive solution to (4) in $V\setminus V(p_{1}\cdots p_{n}L_{f}(p))$ if the scaling vector $c$ and the data vector $u$ are positive vectors.

2.3. $A$-discriminant and principal $A$-determinant

In this subsection we recall the classical definitions of discriminants and resultants.

Definition 10**.**

To any matrix $A$ as above, we can associate the variety $\nabla_{A}$, defined by

[TABLE]

If $\nabla_{A}$ has codimension one in $(\mathbb{C}^{*})^{n}$, then the $A$-discriminant, denoted $\Delta_{A}(f)$, is defined to be the irreducible polynomial that vanishes on $\nabla_{A}$. This polynomial is unique up to multiplication by a scalar.

As in the previous subsection, to $f=\sum_{i=1}^{n}c_{i}\theta^{a_{i}}$ where $c=(c_{1},\ldots,c_{n})\in(\mathbb{C}^{*})^{n}$ we associate the linear form $L_{f}(p)=\sum_{i=1}^{n}c_{i}p_{i}$.

Definition 11**.**

Given $f_{1},\ldots,f_{d}$ where $f_{j}=\sum_{i=1}^{n}c_{i,j}\theta^{a_{i}}$ with coefficients $c_{i,j}\in\mathbb{C}^{*}$ there is a unique irreducible polynomial $R_{A}(f_{1},\ldots,f_{d})$ in the coefficients $c_{i,j}$, so that $R_{A}(f_{1},\ldots,f_{d})=0$ at $(c_{i,j})$ if and only if $L_{f_{1}}=L_{f_{2}}=\cdots=L_{f_{d}}=0$ on the toric variety $V$ [GKZ94, Section 2, Chapter 8]. This polynomial is called the $A$-resultant of $f_{1},\ldots,f_{d}$.

Definition 12**.**

For $f=\sum_{i=1}^{n}c_{i}\theta^{a_{i}}$, we define the principal $A$-determinant as

[TABLE]

This is a polynomial in the coefficients $c_{i}$ of $f$.

In [GKZ94], the principal $A$-determinant $E_{A}(f)$ is related to the polytope $Q$ of the toric variety $V$. More precisely, when the toric variety $V$ is smooth, and $A$ is the matrix whose columns correspond to the lattice points in $Q$, the principal $A$-determinant is

[TABLE]

where the product is taken over all nonempty faces $\Gamma\subset Q$ including $Q$ itself [GKZ94, Theorem 1.2, Chapter 10]. Here, $\Gamma\cap A$ is the matrix whose columns correspond to the lattice points contained in $\Gamma$. When $V$ is not smooth, the radical of its principal $A$-determinant is the polynomial given in (5). The zero locus of $E_{A}(f)$ in $\mathbb{C}^{n}$ will be denoted by $\Sigma_{A}$.

3. ML degree and principal $A$-determinant

We first state and prove our main theorem.

Theorem 13**.**

Let $A=(a_{1}\,\ldots\,a_{n})$ be a full rank $(d-1)\times n$ integer matrix, and let $V^{c}\subset\mathbb{P}^{n-1}$ be the scaled toric variety defined by the monomial parametrization given by (1) where $c\in(\mathbb{C}^{*})^{n}$ is fixed. Then $\mathrm{mldeg}(V^{c})<\deg(V^{(1,1,\ldots,1,1)})$ if and only if $c\in\Sigma_{A}$.

Proof.

Let $V=V^{(1,1,\ldots,1,1)}$. Let $c\in(\mathbb{C}^{*})^{n}\cap\Sigma_{A}$. Then there exists $p\in V$ such that $L_{f}(p)=L_{1}(p)=\cdots=L_{d-1}(p)=0$. Such a $p$ is a solution in $\mathcal{L}_{c,u}^{\prime}$ but not a solution in $\mathcal{L}_{c,u}$. Since the degree of $\mathcal{L}_{c,u}^{\prime}$ is equal to $\deg(V)$, and since $\mathrm{mldeg}(V^{c})$ is equal to the degree of $\mathcal{L}_{c,u}$ we see that $\mathrm{mldeg}(V^{c})<\deg(V)$.

Conversely, suppose $\mathrm{mldeg}(V^{c})<\deg(V)$. There are two ways this can happen (c.f. the proof of Corollary 8). On the one hand, there can exist $p\in V$ so that $L_{f}(p)=0$, which implies that $L_{1}(p)=\cdots=L_{d-1}(p)=0$. This means $c\in\Sigma_{A}$. On the other hand, there can exist a solution $p$ in $\mathcal{L}_{c,u}^{\prime}$ where some $p_{i}=0$. The coordinates of $p$ that are zero cannot be arbitrary. The support of $p$, that is $\{i\,:\,p_{i}\neq 0\}$, is in bijection with the columns $a_{i}$ of $A$ where the convex hull of these columns is a face of $Q=\mathrm{conv}(A)$. Let $\Gamma$ be the corresponding face. Without loss of generality we will assume that $\Gamma\cap A=\{a_{1},\ldots,a_{k}\}$. Let $e<d-1$ be the dimension of $\Gamma$. Since $p_{k+1}=\cdots=p_{n}=0$, the point $p_{\Gamma}=[p_{1}:\cdots:p_{k}]$ is in the toric variety $V_{\Gamma}$ defined by $\Gamma\cap A$. Moreover, $L_{f}(p)=L_{f_{\Gamma}}(p_{1},\ldots,p_{k})$ where $L_{f_{\Gamma}}=\sum_{i=1}^{k}c_{i}p_{i}$ is the linear form associated to the polynomial $f_{\Gamma}$ whose support is $\Gamma$. Similarly, $L_{i}(p)=L_{i}^{\prime}(p_{1},\ldots,p_{k})$ where $L_{i}^{\prime}$ is the linear form associated to $\theta_{i}\frac{\partial f_{\Gamma}}{\partial\theta_{i}}$. If $L_{f}(p)=L_{f_{\Gamma}}(p_{1},\ldots,p_{k})\neq 0$, then $p_{\Gamma}$ is a solution to

[TABLE]

But this is a $d\times k$ linear system of rank $e<d$. Because $u$ is generic, this system does not have a solution. Therefore the only way $p=[p_{1}:\cdots:p_{k}:0:\cdots:0]$ can be a solution in $\mathcal{L}_{c,u}^{\prime}$ is with $L_{f_{\Gamma}}(p_{1},\ldots,p_{k})=L_{f}(p)=0$. We again conclude that $c\in\Sigma_{A}$. ∎

Remark 14**.**

A natural question one can ask is under what conditions a general projective variety has ML degree one. This question was answered in [Huh14].

3.1. Toric Hypersurfaces

Let $A=\left(a_{1}\,\ldots\,a_{d+1}\right)$ be a $(d-1)\times(d+1)$ integer matrix of rank $d-1$. In this case, the toric variety $V$ and each scaled toric variety $V^{c}$ is a hypersurface generated by a single polynomial. This polynomial can be computed as follows. We let $A^{\prime}$ be the matrix obtained by appending a row of $1$’s to the matrix $A$. The kernel of $A^{\prime}$ is generated by an integer vector $w=(w_{1},\ldots,w_{d+1})$ where $\gcd(w_{1},\ldots,w_{d+1})=1$. Without loss of generality we assume that $w_{1},\ldots,w_{\ell}>0$ and $w_{\ell+1},\ldots,w_{d+1}\leq 0$. Then the hypersurface is defined by

[TABLE]

We say $A^{\prime}$ and $A$ are in general position if no $d$ columns of $A^{\prime}$ lie on a hyperplane. The matrix $A^{\prime}$ is in general position if and only if $w$ has full support, i.e. no $w_{i}$ is [math]. This is also equivalent to $Q=\operatorname{conv}(A)$ being a simplicial polytope with exactly $d$ of the columns of $A$ on each facet of $Q$. Recall that $f=\sum_{i=1}^{d+1}c_{i}\theta_{i}^{a_{i}}$. In this setting, the $A$-discriminant $\Delta_{A}$ and the principal $A$-determinant $E_{A}(f)$ can be calculated directly.

Proposition 15**.**

[GKZ94, Chapter 9, Proposition 1.8]** If $A$ is in general position then

[TABLE]

This is also the principal $A$-determinant of this toric hypersurface up to a monomial factor.

Corollary 16**.**

If $A$ is in general position, $\mathrm{mldeg}(V^{c})<\deg(V)$ if and only if

[TABLE]

where $p_{1}^{w_{1}}p_{2}^{w_{2}}\cdots p_{\ell}^{w_{\ell}}-p_{\ell+1}^{-w_{\ell+1}}\cdots p_{d+1}^{-w_{d+1}}$ is the polynomial defining the hypersurface $V$.

4. Rational Normal Scrolls

In this section, we characterize the ML degree of a rational normal scroll. A rational normal scroll is a toric variety associated to a $(d-1)\times n$ matrix $A$ of the following form:

[TABLE]

Here, $n=d-1+n_{1}+n_{2}+\cdots+n_{d-1}$. This description is due to Petrović [Pet08].

We index the data vector and scaling vector as $(u_{ij})$ and $(c_{ij})$ respectively. For $i=1,2,\dots,d-1$, define $g_{i}$ to be the following univariate polynomial in the variable $\theta_{d-1}$:

[TABLE]

Note that for $i\neq d-1$, $g_{i}$ is the $i$th partial derivative of $f$ where

[TABLE]

We now state and prove our result on the ML degree of the rational normal scroll.

Theorem 17**.**

The ML degree of a rational normal scroll given by the matrix $A$ as in (6) is equal to the number of distinct roots of $g_{1}g_{2}\cdots g_{d-1}$.

Proof.

For each $i$, we use $g_{i}^{\prime}$ to denote the derivative of $g_{i}$ with respect to $\theta_{d-1}$, and we write $u_{(d-1)+}=\sum_{j=0}^{n_{d-1}}u_{(d-1)j}$. With this set up, the last equation in (3) is

[TABLE]

Multiply both sides of this equation by $g_{1}g_{2}\cdots g_{d-1}$. Then substituting $(Au)_{i}f=u_{+}\theta_{i}g_{i}$ for each $1\leq i\leq d-2$ in the appropriate summand on the right-hand side gives:

[TABLE]

For generic $(u_{ij})$, the solutions to (3) correspond to the solutions of (7) where $g_{i}\neq 0$ for every $i$.

Suppose that $\theta_{d-1}=\alpha$ is a repeated root of $g_{1}g_{2}\cdots g_{d-1}$. Then either $\alpha$ is a shared root of $g_{i}$ and $g_{j}$ for $i\neq j$ or $\alpha$ is a repeated root of a single $g_{i}$. In the first case, $\theta_{d-1}-\alpha$ divides the right-hand side of (7) since each summand has either $g_{i}$ or $g_{j}$ as a factor. In the second case, $\theta_{d-1}-\alpha$ divides both $g_{i}$ and $g_{i}^{\prime}$ and hence the right-hand side of (7). Factoring out all such repeated roots results in an equation whose degree equals the number of distinct roots of $g_{1}\cdots g_{d-1}$. This shows that the ML degree is at most the number of the distinct roots; we must still argue that no more roots can be factored out.

Suppose $g_{1}g_{2}\cdots g_{d-1}$ has no repeated roots and let $\theta_{d-1}=\beta$ be a root of an arbitrary $g_{i}$. Then $\theta_{d-1}-\beta$ does not divide $g_{1}\cdots g_{i-1}g_{i}^{\prime}g_{i+1}\cdots g_{d-1}$ but does divide every other summand on the right-hand side, so it does not divide the sum as a whole. This follows because for generic data vectors $(u_{ij})$, no coefficient on the right-hand side vanishes. Thus, if the repeated roots of $g_{1}g_{2}\cdots g_{d-1}$ have been removed, the degree of the polynomial and hence the ML degree will not drop any further. ∎

We finish this section with examples.

Example 18**.**

The simplest example of a rational normal scroll is a rational normal curve with $A=\left(0\,1\,\ldots\,n\right)$. We denote the scaled rational normal curve of degree $n$ by $\mathcal{C}_{n}^{c}$. The ML degree of the scaled rational normal curve $\mathcal{C}_{n}^{c}$ is equal to the number of distinct roots of the polynomial

[TABLE]

For example, if $n=4$, the zero locus of the principal $A$-determinant (where the ML degree drops from four) is cut out by the discriminant of a degree four polynomial. The variety where the ML degree drops at least two consists of two irreducible components. The first one corresponds to polynomials of the form $f(\theta)=k(\theta+a)^{2}(\theta+b)^{2}$ and the second corresponds to polynomials of the form $f(\theta)=k(\theta+a)^{3}(\theta+b)$. The variety where the ML degree is one is a determinantal variety, defined by 2-minors of the matrix

[TABLE]

Example 19**.**

In the case that $k=2$, the toric variety is known as a Hirzebruch surface, and we can draw the corresponding polytope (Figure 3).

We denote this surface by $\mathcal{H}_{n_{1},n_{2}}$. For $c=(1,1,\ldots,1,1)$, the ML degree of $\mathcal{H}_{n_{1},n_{2}}$ is $n_{1}+n_{2}-\gcd(n_{1}+1,n_{2}+1)+1$. This is because $g_{1}=1+\theta+\cdots+\theta^{n_{1}}$ and $g_{2}=1+\theta+\cdots+\theta^{n_{2}}$. The roots of $g_{1}$ and $g_{2}$ are respectively the $(n_{1}+1)$th and $(n_{2}+1)$th roots of unity not equal to 1. It follows that there are $\gcd(n_{1}+1,n_{2}+1)-1$ repeated roots of $g_{1}g_{2}$.

Example 20**.**

If we choose the coefficients $c_{ij}$ to be the binomial coefficients $\binom{n_{i}}{j}$, then $g_{i}=(1+\theta_{d-1})^{n_{i}}$ for all $i$. For this choice of coefficients, the ML degree is 1 by Theorem 17, and the likelihood equations (2) have one solution given by

[TABLE]

where $u_{++}=\sum_{i}\sum_{j}u_{ij}$, and $\theta_{d-1}$ equals the unique solution to

[TABLE]

5. Veronese Embeddings

In this section we study Veronese and Veronese-type varieties.

Definition 21**.**

Consider the $(d-1)\times n$ integer matrix $A$ with columns that are the non-negative integer vectors whose coordinate sum $\leq k$. The projective toric variety defined by $A$ is the Veronese ${\rm Ver}(d-1,k)$ for $d,k\geq 1$.

It can be seen that $\deg({\rm Ver}(d-1,k))=k^{d-1}$. A conjecture presented in [Vuc16] is that under the standard embedding, that is, with scaling vector $c=(1,1,\ldots,1,1)$, the ML degree of any Veronese variety equals its degree. We now show this is true when $k\leq d-1$.

Proposition 22**.**

Consider ${\rm Ver}(d-1,k)$ for $k\leq d-1$ which is embedded using the map (1) with scaling given by $c=(1,1,\ldots,1,1)$. Then $c\notin\Sigma_{A}$ and $\mathrm{mldeg}({\rm Ver}(d-1,k))=k^{d-1}$.

Proof.

First note that, by induction, it is sufficient to show that $c\notin\Delta_{A}$. To prove this we need to show that the hypersurface defined by the polynomial

[TABLE]

is nonsingular, where $h_{i}$ is the $i$th complete homogeneous symmetric polynomial in $\theta_{1},\dots,\theta_{d-1}$. The complete homogeneous symmetric polynomials $h_{1},\dots,h_{d-1}$ form an algebra basis for $\Lambda_{d-1}$, the ring of symmetric polynomials in variables $\theta_{1},\dots,\theta_{d-1}$, see [Mac98, (2.8)]. The regular embedding $\varphi$ of $\Lambda_{d-1}$ into the polynomial ring $\mathbb{C}[\theta_{1},\dots,\theta_{d-1}]$ (specified by the definition of $h_{i}$ in terms of $\theta_{j}$’s) induces an embedding $\varphi^{*}$ of ${\rm Spec}(\Lambda_{d-1})$ into $\mathbb{A}^{d-1}={\rm Spec}(\mathbb{C}[\theta_{1},\dots,\theta_{d-1}])$. Consider the hypersurface $V=V(1+h_{1}+h_{2}+\cdots+h_{k})\subset{\rm Spec}(\Lambda_{d-1})$. For $k\leq d-1$ we have that $1+h_{1}+h_{2}+\cdots+h_{k}$ is a linear form in the indeterminates $h_{1},\dots,h_{d-1}$. Hence $V$ is smooth in ${\rm Spec}(\Lambda_{d-1})$. Since the induced map $\varphi^{*}$ is a regular embedding, $\varphi^{*}(V)$ is smooth in $\mathbb{A}^{d-1}$, and by the definition of $h_{i}$ and the map $\varphi$ we have that $\varphi(1+\sum_{i=1}^{k}h_{i})=f$.Thus $f$ defines a smooth hypersurface in $\mathbb{A}^{d-1}$, and $c=(1,1,\dots,1,1)\notin\Delta_{A}$. The conclusion follows by Theorem 13. ∎

Remark 23** (Hypersimplex).**

Let $k\leq d-1$. Consider the $(d-1)\times n$ integer matrix $A$ where $n=\binom{d}{k}$ and where the columns of $A$ are the vectors in $\{0,1\}^{d-1}$ that have precisely $k$ or $k-1$ entries equal to $1$. The $(d-1)$-dimensional polytope $Q={\rm conv}(A)$ is called the hypersimplex. The projective toric variety $V$ associated to $A$ represents generic torus orbits on the Grassmannian of $k$-dimensional linear subspaces in $\mathbb{C}^{d}$. It is shown in the discussion preceding [HS17, Proposition 4.7] that $c=(1,1,\dots,1,1)\notin\Sigma_{A}$ for this matrix $A$. We remark that this result can be seen from the same proof given in Proposition 22.

Let $e_{k-1}$ and $e_{k}$ be elementary symmetric polynomials in $\theta_{1},\ldots,\theta_{d-1}$ of degree $k-1$ and $k$, respectively. Similar to above, we need to only show that the hypersurface defined by the polynomial

[TABLE]

is nonsingular. This proceeds in a manner identical to the proof of Proposition 22 since the elementary symmetric polynomials also form an algebra basis for $\Lambda_{d-1}$, and hence are algebraically independent (see [Mac98, (2.4)]). Thus, for the projective toric variety $V$ associated to the hypersimplex, ${\rm mldeg}(V)$ is equal to the normalized volume of the hypersimplex $Q$, which is the Eulerian number $A(d-1,k-1)$.

It should be noted that the form of the definition of the hypersimplex given above differs slightly from that in [HS17], however the definitions are equivalent (i.e. the resulting exponents of the monomial maps differ by elementary row operations). The difference is due to the fact that in [HS17] a projective toric variety is defined by a monomial map given by a $d\times n$ integer matrix $\mathcal{A}$ with full rank and the vector $(1,\dots,1)$ in the row space. In the convention used in our note we (equivalently) assume that the last row of what would be their matrix $\mathcal{A}$ consists only of ones and omit this row from our matrix $A$.

Returning to ${\rm Ver}(d-1,k)$, we note that in the case $k=2$, the polynomial $f$ reduces to a quadratic form, and a direct proof that $c=(1,1,\ldots,1,1)\notin\Sigma_{A}$ can be given by looking at the principal minors of the corresponding symmetric matrix

[TABLE]

where $f=(1,\theta_{1},\ldots,\theta_{d-1})C(1,\theta_{1},\ldots,\theta_{d-1})^{T}$.

Proposition 24**.**

Consider ${\rm Ver}(d-1,2)$. Then $c=(1,1,\dots,1,1)\notin\Sigma_{A}$. In particular, we have that ${\rm mldeg}({\rm Ver}(d-1,2))=\deg({\rm Ver}(d-1,2))=2^{d-1}.$

Proof.

Let $Q=\operatorname{conv}(A)$. From (5) we have that

[TABLE]

The A-discriminant of the toric variety corresponding to each face will be given by a principal minor (of size corresponding to the dimension of the face) of the symmetric matrix $C$ in (8). Any principal minor of size $r\times r$ will have the same form as $C$. Hence, to show that $c=(1,1,\dots,1,1)\notin\Sigma_{A}$ we need only verify that the $r\times r$ matrix

[TABLE]

has a nonzero determinant. A calculation shows that $\det(\tilde{C})=r+1\neq 0$ for all $r\geq 2$. Thus, by Theorem 13, it follows that ${\rm mldeg}({\rm Ver}(d-1,2))=\deg({\rm Ver}(d-1,2))=2^{d-1}.$ ∎

Now we give a sufficient condition for the opposite effect to happen: having ML degree equal to one, instead of having ML degree equal to the degree. The polynomial $f$ has terms of degree $\leq k$ in $\theta_{1},\ldots,\theta_{d-1}$. It uniquely corresponds to a symmetric tensor $\mathcal{C}$ of size $d\times d\times\ldots\times d$, where the product is $k$-fold. We say $\mathcal{C}$ has tensor rank one if the polynomial $f$ is a power of a linear form.

Theorem 25**.**

Let $\mathcal{C}$ be the symmetric tensor corresponding to $f$ where

[TABLE]

If the tensor rank of $\mathcal{C}$ is equal to one then $\mathrm{mldeg}({\rm Ver}(d-1,k)^{c})=1$.

Proof.

In this case we have that $f=L^{k}$ where $L=b_{0}+b_{1}\theta_{1}+\cdots+b_{d-1}\theta_{d-1}$. Then in the likelihood equations (3) we have that $L\neq 0$ and

[TABLE]

So the factor $L^{k-1}$ cancels out and we obtain a linear system in the $\theta_{i}$’s. ∎

Example 26**.**

Let $d=3$ and $k=2$, that is, we consider the second Veronese embedding of $\mathbb{P}^{2}$. In this case the matrix defining ${\rm Ver}(2,2)$ is $A=\begin{pmatrix}0&1&2&0&0&1\\ 0&0&0&1&2&1\end{pmatrix}.$ The polytope $Q={\rm conv}(A)$ is given in Figure 4 with the labeled lattice point $a_{i}$ corresponding to the $i$th column of the matrix $A$.

By (5), the ideal of $\Sigma_{A}\subset(\mathbb{C}^{*})^{6}$ is defined by the polynomial $E_{A}$ in (26),

[TABLE]

For generic values of $c_{ij}$ in $C$ as in equation (8), including $c_{ij}=1$, we have that $\mathrm{mldeg}(V^{c})=4=2^{2}$. On the other hand, if we take a scaling $c$ so that $C$ is the matrix with all entries equal to $2$, the ML degree drops to $\mathrm{mldeg}(V^{c})=1$. See Table 1 for other combinations. Note that when the rank of $C$ is 1, the ML degree is 1, illustrating Theorem 25. We also see that the converse does not hold: the ML degree can be 1 even with the rank of $C$ equal to 3.

6. Segre Embeddings

In this section, we study the maximum likelihood degree of the scaled Segre embedding $V^{c}$. We give a sufficient condition for $V^{c}$ to have ML degree one.

For $c\in(\mathbb{C}^{*})^{mn}$, let $\psi^{c}$ be the map defined by

[TABLE]

Then $V^{c}$ is a scaled Segre embedding of $\mathbb{P}^{m-1}\times\mathbb{P}^{n-1}$. In terms of the coordinates $p_{ij}$ with $1\leq i\leq m$ and $1\leq j\leq n$ on $\mathbb{P}^{mn-1}$, the defining ideal $I^{c}$ of $V^{c}$ is given by the 2-minors of the $m\times n$ matrix with $ij$th entry $\frac{p_{ij}}{c_{ij}}$.

From this it follows that the likelihood equations (4) are

[TABLE]

Here $u_{i+}=\sum_{j=1}^{n}u_{ij}$, $u_{+j}=\sum_{i=1}^{m}u_{ij}$, and $u_{++}=\sum_{i=1}^{m}\sum_{j=1}^{n}u_{ij}$. We define $p_{i+}$, $p_{+j}$, and $p_{++}$ analogously.

For all $c\in(\mathbb{C}^{*})^{mn}$, we have $\deg(V^{c})=\binom{m+n-2}{m-1}$. If $c=(1,1,\ldots,1,1)$, then the ML degree of $V^{c}$ is well known to be one (see [PS05, Example 1.12] and [DSS09, Example 2.1.2]). Also the principal $A$-determinant for Segre embedding is easy to describe [GKZ94, Chapter 9, Section 1]:

[TABLE]

where $[a_{1},\ldots,a_{i};b_{1},\ldots,b_{i}]$ is the submatrix of $c$ indexed by the rows $1\leq a_{1}<\cdots<a_{i}\leq m$ and the columns $1\leq b_{1}<\cdots<b_{i}\leq n$.

Example 27**.**

Let $V^{c}$ be the scaled Segre variety with scaling

[TABLE]

Here, $m=n=3$, and the matrix $A$ is given by

[TABLE]

The degree of $V^{c}$ is $\binom{3+3-2}{3-1}=6$. The assignments

[TABLE]

produce, respectively, the ML degrees $\mathrm{mldeg}(V^{c})=1,2,3,4,5,6$. This shows $V^{c}$ can have any ML degree $i$, for $1\leq i\leq\operatorname{deg}(V^{c})$.

Conjecture 28**.**

For each $i\in\{1,\ldots,\deg(V)\}$, there exists $c\in(\mathbb{C}^{*})^{mn}$ such that the scaled Segre variety $V^{c}$ has ML degree $i$.

Example 29**.**

The scaled toric variety $V^{c}$ with the configuration

[TABLE]

is isomorphic to the scaled Segre embedding of $\mathbb{P}^{1}\times\mathbb{P}^{1}\times\mathbb{P}^{1}\times\mathbb{P}^{1}$ into $\mathbb{P}^{15}$. The degree of $V^{c}$ in this instance is $24$. For each of the possible ML degrees $i=1,2,\ldots,24$, we have found a scaling vector $c\in(\mathbb{C}^{*})^{16}$ for which $\textrm{mldeg}(V^{c})=i$. We see from (12) that for a point to lay in the vanishing of the principal A-determinant we must have that certain minors of the matrix $\left(c_{ij}\right)$ vanish; this observation informed our choice of scaling vector $c$.

We conclude this section with a result analogous to Theorem 25.

Proposition 30**.**

If $\mathrm{rank}(c)=1$, then $\mathrm{mldeg}(V^{c})=1$.

Proof.

If $\operatorname{rank}(c)=1$ then all the 2-minors of $c$ vanish. This implies that $c_{ij}c_{kl}=c_{il}c_{kj}$ for all $1\leq i,k\leq m$ and $1\leq j,l\leq n$. Then our likelihood equations reduce to

[TABLE]

which are the likelihood equations for the unscaled Segre embedding of $\mathbb{P}^{m-1}\times\mathbb{P}^{n-1}$. This is known to have ML degree equal to one. ∎

7. Hierarchical log-linear and graphical models

A large class of discrete exponential models that are used in statistical practice are hierarchical log-linear models [BFH07] and undirected graphical models [Lau96]. In both cases, the toric variety is constructed via a parametrization based on a simplicial complex where the vertices of the simplicial complex correspond to discrete random variables. In the second case, the simplicial complex is the clique complex of a given graph. These toric varieties, their defining ideals, and the polyhedral geometry of $Q=\operatorname{conv}(A)$ of the underlying matrix $A$ have been intensely studied starting with [HS02] for the hierarchical log-linear models and [GMS06] for graphical models.

Decomposable models form a subclass of these models, and they present attractive factorization properties. For instance, a graphical model is decomposable if and only if the underlying graph is a chordal graph [Lau96]. In particular, we have the following.

Theorem 31**.**

[GMS06, Theorem 4.4]** A hierarchical log-linear model, and hence an undirected graphical model, is decomposable if and only if the ML degree of the corresponding toric variety is one.

For example, the (unscaled) Segre embedding of $\mathbb{P}^{m_{1}}\times\cdots\times\mathbb{P}^{m_{k}}$ is a decomposable graphical model, and hence has ML degree one. In this section we consider some examples of nondecomposable models. First we study a family of examples known as no-three-way interaction models.

7.1. No-three-way interaction models

Example 32** (Binary 3-cycle).**

Consider the hierarchical log-linear model arising from the $3$-cycle depicted in Figure 5. This is also known as the no-three-way interaction model. We let the random variables $X,Y,Z$ corresponding to the vertices to be binary variables. One parametrization $\psi:(\mathbb{C}^{*})^{6}\to(\mathbb{C}^{*})^{8}$ is given by the matrix

[TABLE]

We see that the $6\times 8$ matrix $A$ is in general position as in Section 3.1. Hence, the toric variety $V$ defined by $A$ is a hypersurface with a generator of full support. For $i,j,k\in\{0,1\}$, let $p_{ijk}=\mathrm{prob}(X=i,Y=j,Z=k)$. Labeling the columns of the matrix with the $p_{ijk}$ in increasing lexicographic order, we find that the unique toric generator is

[TABLE]

By Corollary 16, the principal $A$-determinant is equal to the $A$-discriminant, and it is

[TABLE]

The toric variety itself has degree $4$ since it has a single degree $4$ generator. However, if we choose the scaling coefficients $c_{ijk}=1$ for all $i,j,k\in\{0,1\}$, we see that this scaling vector lies in $\Sigma_{A}$, and thus $\mathrm{mldeg}(V^{c})$ drops, in this case, to $3$. In fact, for this example, the only possibilities for $\mathrm{mldeg}(V^{c})$ are $3$ when the scaling vector lies in $\Sigma_{A}$ and $4$ when the scaling vector does not.

Proposition 33**.**

The ML degree of the binary $3$-cycle is $4$ unless $c\in(\mathbb{C}^{*})^{d+1}$ is in $\Sigma_{A}$. If $c\in\Sigma_{A}$, then $\mathrm{mldeg}(V^{c})=3$.

Proof.

We fix the lexicographic monomial ordering

[TABLE]

over $\mathbb{Q}[p_{000},p_{001},p_{010},p_{011},p_{100},p_{101},p_{110},p_{111},c_{000},c_{001},c_{010},c_{011},c_{100},c_{101},c_{110},c_{111}]$.

Let $u=(u_{000},\ldots,u_{111})$ be a data vector in $\mathbb{Z}^{8}$ and let $u_{+}=\sum u_{ijk}$. Consider

[TABLE]

where $f=\sum c_{ijk}p_{ijk}$. Using [GS17] and the random data vector $u=(2,3,5,7,11,13,17,19)$, we calculate a Gröbner basis for $I$:

[TABLE]

As long as each $c_{ijk}\in\mathbb{C}^{*}$ and satisfies the equation $g_{1}=\Delta_{A}=0$, we see that $g_{2}$ is a univariate polynomial in $p_{111}$ of degree $3$, and the initial terms of $g_{3}$ through $g_{15}$ have degree $1$ in $p_{ijk}$. Therefore, when the coefficient vector lies in $\Sigma_{A}$, the ML degree is always $3$. ∎

More generally, we consider the no-three-way interaction model $C_{m}$ based on the $3$-cycle with one $m$-ary variable and two binary variables. That is, for $i\in\{0,\ldots,m-1\}$ and $j,k\in\{0,1\}$, we set

[TABLE]

One can compute the normalized volume of $Q=\operatorname{conv}(A_{m})$ for the matrix $A_{m}$ associated $C_{m}$. This gives a formula for the degree of $C_{m}$.

Proposition 34**.**

The degree of $C_{m}$ is $m2^{m-1}$.

When $m=3$, $\deg(C_{3})=12$ and $\mathrm{mldeg}(C_{3})=7$. The toric ideal is

[TABLE]

Observe that the toric variety is not a hypersurface, so calculating the principal $A$-determinant requires more work than the previous example. When we compute $\Delta_{A_{3}}$, the $A$-discriminant for the whole polytope $Q=\operatorname{conv}(A_{3})$, we do not get a hypersurface:

[TABLE]

We already see that $c=(1,1,\ldots,1,1)$ is in $\Sigma_{A}$ and hence we know $\mathrm{mldeg}(C_{3}^{c})<\deg(C_{3})=12$. Also comparing the toric variety and the discriminant locus $\nabla_{A_{3}}$ we see that they are identical. Such a toric variety is known as self-dual. Based on our computations for $m\leq 5$ we state the following conjecture.

Conjecture 35**.**

The no-three-way interaction model with one $m$-ary variable and two binary variables is self-dual.

To compute the entire principal $A$-determinant, we need to consider all faces of $Q$. If a face $\Gamma$ has dimension $e$ and $|\Gamma\cap A_{3}|=e+1$, then $\Delta_{\Gamma\cap A_{3}}$ is the unit ideal. Moreover, for all codimension one and codimension two faces $\Gamma$ which are not simplices, we observe that $\nabla_{\Gamma\cap A_{3}}$ are not hypersurfaces. They also lie in a coordinate hyperplane. For instance, for the facet $\Gamma$ where the elements of $\Gamma\cap A_{3}$ correspond to

[TABLE]

the discriminant is defined by

[TABLE]

Therefore, for $c=(1,1,\ldots,1,1)$, these discriminants do not contribute to a drop in the ML degree. There is a total of three codimension $3$ faces that are not simplices. In each case $\Gamma\cap A_{3}$ is in general position and their discriminants are given by Corollary 16. We list the indeterminates which correspond to the elements of $\Gamma\cap A_{3}$ in each of these three faces:

[TABLE]

It turns out that each face of codimension $>3$ is a simplex and there are no more discriminants contributing to the principal $A$-determinant. Based on our computations for $m\leq 5$, we state the following conjecture.

Conjecture 36**.**

The ML degree of $C_{m}$ is $2^{m}-1$.

7.2. The binary $4$-cycle

The binary $4$-cycle is the model we have used in Example 3. Let $S,B,H,L$ be four binary random variables and let $X=(S,B,H,L)$ be the joint random variable where we set

[TABLE]

The family of probability distributions for $X$ which factor according to the graphical model depicted in Figure 2 can be described by the following monomial parametrization: let $a_{ij},b_{jk},c_{k\ell},d_{i\ell}$ be parameters for $i,j,k,\ell\in\left\{0,1\right\}$ and let $p_{ijk\ell}=a_{ij}b_{jk}c_{k\ell}d_{i\ell}$. Below is the matrix $A$ arising from a different parametrization that gives a full rank matrix. The columns are labeled by $p_{0000},p_{0001},\ldots,p_{1111}$ in increasing lexicographic ordering.

[TABLE]

The Zariski closure $V$ of the image of this parametrization is a toric variety that is defined by the following prime ideal:

[TABLE]

The degree of this toric variety is $64$ and $\mathrm{mldeg}(V)=13$. This was computed in [GMS06, p. 1484] where the question of explaining the fact $\mathrm{mldeg}(V)=13$ was first raised. One road to an explanation is to compute the $A$-discriminants of all the faces of $Q=\operatorname{conv}(A)$, and then determine the contribution of each for the drop in the ML degree. There are multiple faces that might contribute to such a drop. At this moment, we do not understand how different discriminants interact. Also, it is not possible to compute every discriminant. For instance, with standard elimination methods we were not able to compute $\Delta_{A}$ corresponding to the entire polytope. However, one can check that the polynomial $f=\sum_{i}^{16}\theta^{a_{i}}$ where $c=(1,1,\ldots,1,1)$ has singularities in $(\mathbb{C}^{*})^{d-1}$. Hence, we conclude that $\Delta_{A}$ does contribute to the drop from the generic ML degree. However, we can compute discriminants for facets and lower codimension faces.

The polytope $Q$ has $24$ facets. Of these, $6$ correspond to the [math] ideal and $2$ correspond to toric hypersurfaces of the form $c_{0001}c_{1101}-c_{0110}c_{1100}$. Of the remaining 16 facets, eight lie on a coordinate hyperplane. An example of such a discriminant is

[TABLE]

There are 5 facets which give rise to hypersurfaces, and most are too long to display. Of these discriminants, there are three of degree $5$, one of degree $6$, and one of degree $16$. The remaining three facets are neither of the above types. Two have 5 generators of degree 4 of the form

[TABLE]

and the remaining discriminant has $3$ generators of degree $2$,

[TABLE]

Because there are discriminants of codimension one faces that are not hypersurfaces, we must also analyze the 168 codimension two faces of $Q$. Of these, 86 are trivial and 61 correspond to toric hypersurfaces. Of the remaining 21, there is only a single face whose discriminant does not lie on a coordinate hyperplane, and the discriminant of this face is a hypersurface generated by

[TABLE]

8. ML Estimate Homotopies

In this section we use homotopy continuation to track between ML estimates of different scalings of a given statistical model. Moreover, the endpoints may correspond to different scalings of the model with different ML degrees.

Consider the case where we have found a particular scaling $c_{easy}$ such that the ML degree drops to one, but we really wish to compute the ML estimate for the natural statistical model corresponding to the scaling $c_{stat}$. The strategy is to apply a parameter homotopy [LSY89], [MS89] between $c_{easy}$ and $c_{stat}$, tracking the unique solution of the former to a solution of the latter. We argue that the endpoint is the unique ML estimate in Birch’s Theorem 9.

Example 37**.**

(Veronese) We illustrate the strategy with a Veronese model. Let

[TABLE]

with $c_{stat}=(1,1,1,1,1,1)$ and $c_{easy}=(1,2,1,2,2,1)$. We can check (see Section 5) that the ML degree corresponding to $c_{stat}$ is 4, while the ML degree for $c_{easy}$ is 1.

Suppose we observe the data vector $u=(1,3,5,7,9,2)$. Computing the unique solution for $c_{easy}$ we obtain the ML estimate $\hat{\theta}_{easy}=(0.0493,1.8333,1.6667)$. We track this point with a parameter homotopy (14) towards $c_{stat}$, and obtain the point $\hat{\theta}_{track}=(0.0863,1.6326,1.5150)$ with corresponding $\hat{p}_{stat}=(0.09,0.14,0.23,0.13,0.21,0.20)$. The homotopy we track is

[TABLE]

Or alternatively, if we define $c(t)=27+27t$, we can write it as

[TABLE]

To verify this, we solve the critical equations for $c_{stat}$ to obtain the four solutions

[TABLE]

where $\hat{\theta}_{stat}$ is the solution with positive coordinates. Observe that $\hat{\theta}_{track}=\hat{\theta}_{stat}$, as desired.

Let $c_{easy}$ and $c_{stat}$ be scalings with positive entries. Let $F_{easy}(u,\theta)$ be the difference of the left and right sides of the equations of Definition 6 with $c$ taken to be $c_{easy}$, and similarly for $F_{stat}(u,\theta)$. We may now state the main result of this section.

Theorem 38**.**

Fix a generic data vector $u$ with positive entries. Let $\hat{p}_{easy}$ and $\hat{p}_{stat}$ be the respective ML estimates for these scaling problems for the data $u$ and consider the homotopy below:

[TABLE]

Let $\gamma$ denote the path of the homotopy whose start point (at $t=1$) corresponds to $\hat{p}_{easy}$. Then the endpoint of $\gamma$ (at $t=0$) is $\hat{p}_{stat}$.

Proof.

By Birch’s theorem (Theorem 9), the likelihood equations are given by $A\hat{p}=\frac{1}{u_{+}}Au$ for a data vector $u$. Let $\hat{p}_{win}$ and $\hat{p}_{stat}$ be the two monomial vectors for $c_{easy}$ and $c_{stat}$ respectively. The homotopy above can be rewritten as $A\cdot(t\hat{p}_{c(t)}+\frac{1}{u_{+}}u)$ where $c(t)={tc_{stat}+(1-t)c_{easy}}$. Since $c(t)$ is positive for any positive real $c_{easy}$, $c_{stat}$ and $t\in[0,1]$, Birch’s theorem states that there is exactly one positive real solution to this system at every point along the homotopy path. Thus, as long as no paths intersect, the unique solution for $c_{easy}$ will track to the unique positive real solution for $c_{stat}$. But paths only intersect where the Jacobian of the system drops rank, which we now show cannot occur.

For a value of $\hat{p}$ at some point on the path, the Jacobian matrix can be written as the product

[TABLE]

Because we require the $\theta_{i}$’s to be nonzero, we can write this as

[TABLE]

Note that if we denote $A=(a_{ij})$, then $\theta_{i}\frac{\partial p_{j}}{\partial\theta_{i}}=a_{ij}p_{j}$, so we may rewrite the middle matrix as $\operatorname{diag}(p_{1},p_{2},\ldots,p_{n})A^{T}$ and the Jacobian thus factors as

[TABLE]

This is a product of full rank matrices, hence the Jacobian has full rank. ∎

Notice that the straight-line homotopy between $F_{easy}$ and $F_{stat}$ is equivalent to a parameter homotopy between $c_{easy}$ and $c_{stat}$; see Example 37.

8.1. Timings

To test the viability of using homotopy methods, we ran timing comparisons with the iterative proportional scaling algorithm (IPS, Algorithm 2.1.9 in [DSS09]). We programmed IPS in Python using numpy for fast linear algebra calculations and a target residual of $\epsilon=10e^{-12}$, while for the homotopy approach we used Macaulay2 [GS17] for setup and PHCpack [Ver99] for the actual path tracking. We choose the rational normal scroll (6) example from Section 4, setting $n_{i}=k$ for all $i$ and varying $k\in\{4,\ldots,13\}$ and $d-1\in\{5,10,15\}$. For the homotopy start system, we set $c_{ij}=\binom{n_{i}}{j}$ and tracked from the closed form solution given in Example 20. It is clear from Figure 6 that as the problem size grows, the homotopy method has a significant speed advantage over the iterative proportional scaling algorithm.111Source code can be found at https://github.com/nbliss/Likelihood-MRC