Moment Identifiability of Homoscedastic Gaussian Mixtures

TL;DR

This paper investigates the identifiability of homoscedastic Gaussian mixture models using moment varieties, providing classification results and explicit parameter recovery methods for certain cases.

Contribution

It introduces a geometric approach to Gaussian mixture identifiability via secant varieties and classifies cases with defective moment varieties for moments up to order three.

Findings

Identifiability determined for mixtures with fewer components than the space dimension.

Closed-form parameter recovery for two-component one-dimensional mixtures using moments up to order four.

Connection between rank estimation and secant varieties of rational normal curves.

Abstract

We consider the problem of identifying a mixture of Gaussian distributions with same unknown covariance matrix by their sequence of moments up to certain order. Our approach rests on studying the moment varieties obtained by taking special secants to the Gaussian moment varieties, defined by their natural polynomial parametrization in terms of the model parameters. When the order of the moments is at most three, we prove an analogue of the Alexander-Hirschowitz theorem classifying all cases of homoscedastic Gaussian mixtures that produce defective moment varieties. As a consequence, identifiability is determined when the number of mixed distributions is smaller than the dimension of the space. In the two component setting we provide a closed form solution for parameter recovery based on moments up to order four, while in the one dimensional case we interpret the rank estimation problem in terms of secant varieties of rational normal curves.

Figures

Tables

Table 1. Table 1: All instances of defective varieties Sec k H ⁡ ( 𝒢 n , 3 ) subscript superscript Sec 𝐻 𝑘 subscript 𝒢 𝑛 3 \operatorname{Sec}^{H}_{k}(\mathcal{G}_{n,3}) for n = 1 , … , 7 𝑛 1 … 7 n=1,\dots,7 with d = 3 𝑑 3 d=3 . The column ‘par’ denotes the number of parameters and ‘exp’ the expected dimension.

$n$	$k$	$d$	par	$N$	exp	$dim$	$\delta$	$\Delta$
1	1	3	2	3	2	2	0	0
1	2	3	4	3	3	3	0	1
2	2	3	8	9	8	7	1	1
2	3	3	11	9	9	9	0	2
3	2	3	13	19	13	12	1	1
3	3	3	17	19	17	15	2	2
3	4	3	21	19	19	19	0	2
4	2	3	19	34	19	18	1	1
4	3	3	24	34	24	22	2	2
4	4	3	29	34	29	27	2	2
4	5	3	34	34	34	34	0	0
5	2	3	26	55	26	25	1	1
5	3	3	32	55	32	30	2	2
5	4	3	38	55	38	36	2	2
5	5	3	44	55	44	44	0	0
5	6	3	50	55	50	50	0	0
5	7	3	56	55	55	54	1	2

Equations

$${\lambda }_1{f}_{\mathcal{N}({\mu }_1,{\Sigma }_1)}(x)+\cdots +{\lambda }_k{f}_{\mathcal{N}({\mu }_k,{\Sigma }_k)}(x)$$

$$\sum _{i=1}^k{\lambda }_i{f}_{\mathcal{N}({\mu }_i,\Sigma )}(x)$$

$$M_{n,k,d}:{\Theta }_{n,k}^H\to {\mathrm{Sec}}_k^H({\mathcal{G}}_{n,d}).$$

$$\begin{array}{ccc}{m}_{10}& =& \lambda {\mu }_{11}+(1-\lambda ){\mu }_{21}\\ m_{01}& =& \lambda {\mu }_{12}+(1-\lambda ){\mu }_{22}\\ m_{20}& =& \lambda ({\mu }_{11}^2+{\sigma }_{11})+(1-\lambda )({\mu }_{21}^2+{\sigma }_{11})\\ m_{02}& =& \lambda ({\mu }_{12}^2+{\sigma }_{22})+(1-\lambda )({\mu }_{22}^2+{\sigma }_{22})\\ m_{11}& =& \lambda ({\mu }_{11}{\mu }_{12}+{\sigma }_{12})+(1-\lambda )({\mu }_{21}{\mu }_{22}+{\sigma }_{12})\\ m_{30}& =& \lambda ({\mu }_{11}^3+3{\sigma }_{11}{\mu }_{11})+(1-\lambda )({\mu }_{21}^3+3{\sigma }_{11}{\mu }_{21})\\ m_{03}& =& \lambda ({\mu }_{12}^3+3{\sigma }_{22}{\mu }_{12})+(1-\lambda )({\mu }_{22}^3+3{\sigma }_{22}{\mu }_{22})\\ m_{21}& =& \lambda ({\mu }_{11}^2{\mu }_{12}+{\sigma }_{11}{\mu }_{12}+2{\sigma }_{12}{\mu }_{11})+(1-\lambda )({\mu }_{21}^2{\mu }_{22}+{\sigma }_{11}{\mu }_{22}+2{\sigma }_{12}{\mu }_{21})\\ m_{12}& =& \lambda ({\mu }_{11}{\mu }_{12}^2+{\sigma }_{22}{\mu }_{11}+2{\sigma }_{12}{\mu }_{12})+(1-\lambda )({\mu }_{21}{\mu }_{22}^2+{\sigma }_{22}{\mu }_{21}+2{\sigma }_{12}{\mu }_{22})\end{array}$$

$$M_X(u)=\sum _{(a_1,\dots ,a_n)}{m}_{a_1,\dots ,a_n}[X]\frac{u_1^{a_1}\dots u_n^{a_n}}{a_1!\dots a_n!},K_X(u)=\sum _{(a_1,\dots ,a_n)}{\kappa }_{a_1,\dots ,a_n}[X]\frac{u_1^{a_1}\dots u_n^{a_n}}{a_1!\dots a_n!}.$$

$$M_X(u)=\exp (K_X(u)),K_X(u)=\log (M_X(u))$$

$${\mathbb{A}}_n^M=\left\{m_0=1\right\},{\mathbb{A}}_n^K=\left\{{\kappa }_0=0\right\}.$$

$${\int }_{{\mathbb{R}}^n}f(x){\delta }_{\mu }(x):=f(\mu ).$$

$$M_X(u)=\mathbb{E}[e^{u^{t}X}]=e^{u^t\mu }=\sum _{(a_1,\dots ,a_n)}{\mu }_1^{a_1}\dots {\mu }_n^{a_n}\frac{u_1^{a_1}\dots u_n^{a_n}}{a_1!\dots a_n!}.$$

$$K_X(u)=\log M_X(u)=\log e^{u^t\mu }=u^t\mu ={\mu }_1{u}_1+\cdots +{\mu }_n{u}_n,$$

$$f_{(\mu ,\Sigma )}(x):=\frac{1}{\sqrt{\det (2\pi \Sigma )}}{e}^{-\frac{1}{2}(x-\mu {)}^t{\Sigma }^{-1}(x-\mu )}.$$

$$M_X(u)=e^{u^t\mu +\frac{1}{2}{u}^t\Sigma u},K_X(u)=u^t\mu +\frac{1}{2}{u}^t\Sigma u.$$

$$M_X(u)=\frac{\exp (u^t\mu )}{1-\frac{1}{2}{u}^t\Sigma u},K_X(u)=u^t\mu -\log \left(1-\frac{1}{2}{u}^t\Sigma u\right)$$

$$M_{X+Y}(u)=\mathbb{E}[e^{u^t(X+Y)}]=\mathbb{E}[e^{u^{t}X}{e}^{u^{t}Y}]=\mathbb{E}[u^{t}X]\cdot \mathbb{E}[u^{t}Y]=M_X(u)\cdot M_Y(u).$$

$$K_{X+Y}(u)=\log (M_{X+Y})=\log (M_X{M}_Y)=\log (M_X)+\log (M_Y)=K_X(u)+K_Y(u).$$

$$M_{AX+b}(u)=M_{AX}(u)\cdot M_b(u)=\mathbb{E}[e^{u^{t}AX}]\mathbb{E}[e^{u^{t}b}]=e^{u^{t}b}\cdot \mathbb{E}[e^{(A^{t}u{)}^{t}X}]=e^{u^{t}b}\cdot M_X(A^{t}u)$$

$$K_{AX+b}(u)=\log (M_{AX+b}(u))=\log (e^{u^{t}b}{M}_X(A^{t}u))=u^{t}b+K_X(A^{t}u).$$

$${\int }_{{\mathbb{R}}^n}g(x)\left(\sum _{i=1}^k{\lambda }_i{f}_{({\mu }_i,{\Sigma }_i)}(x)\right)dx=\sum _{i=1}^k{\lambda }_i{\int }_{{\mathbb{R}}^n}g(x)f_{({\mu }_i,{\Sigma }_i)}(x)dx,$$

$${\lambda }_1{f}_{({\mu }_1,\Sigma )}(x)+\cdots +{\lambda }_k{f}_{({\mu }_k,\Sigma )}(x)$$

$${\Theta }_{n,k}^H:=({\mathbb{R}}^n{)}^{\times k}\times {\mathbb{R}}^{k-1}\times {\mathrm{Sym}}^2{\mathbb{R}}^n=\{({\mu }_1,\dots ,{\mu }_k),({\lambda }_1,\dots ,{\lambda }_k),\Sigma |{\lambda }_1+\cdots +{\lambda }_k=1\},$$

$$\dim {\Theta }_{n,k}^H=nk+k-1+\frac{n(n+1)}{2}=(n+1)\left(k+\frac{n}{2}\right)-1.$$

$$M_{n,k,d}:{\Theta }_{n,k}^H\to {\mathbb{A}}_{n,d}^M.$$

$${\Delta }_{n,k,d}^H=\dim {\Theta }_{n,k}^H-\dim {\mathrm{Sec}}_k^H({\mathcal{G}}_{n,d}).$$

$${\Delta }_{n,k,d}^H\ge \max \left\{(n+1)\left(k+\frac{n}{2}\right)-\left(\genfrac{}{}{0ex}{}{n+d}{d}\right),0\right\}.$$

$$\dim {\mathbb{A}}_{n,d}^M=\dim \mathbb{R}[[u_1,\dots ,u_n]]/(u_1,\dots ,u_n{)}^{d+1}-1=(\genfrac{}{}{0ex}{}{n+d}{d})-1.$$

$${\Delta }_{n,k,d}^H=\dim {\Theta }_{n,k}^H-\dim {\mathrm{Sec}}^H({\mathcal{G}}_{n,d})\ge \dim {\Theta }_{n,h}^H-\dim {\mathbb{A}}_{n,d}^M$$

$${\delta }_{n,k,d}^H:={\Delta }_{n,k,d}^H-\max \left\{(n+1)\left(k+\frac{n}{2}\right)-\left(\genfrac{}{}{0ex}{}{n+d}{d}\right),0\right\}.$$

$${\lambda }_1{\delta }_{{\mu }_1}(x)+\cdots +{\lambda }_k{\delta }_{{\mu }_k}(x).$$

$$K_{B+Z}(u)=K_B(u)+K_Z(u)=K_A(u)+m^{t}u+\frac{1}{2}{u}^t\Sigma u$$

$$K:{\Theta }_{n,k}^0\times {\mathbb{R}}^n\times {\mathrm{Sym}}^2{\mathbb{R}}^n\to {\mathbb{A}}_{n,d}^K,(A,m,\Sigma )\mapsto K_A(u)+m^{t}u+\frac{1}{2}{u}^t\Sigma u$$

Full text

AMS Subject Classifications: 62R01; 62F10; 13P25; 14N07; 14Q15

1 Introduction

In the context of algebraic statistics [19], moments of probability distributions have recently been explored from an algebraic and geometric point of view [1, 4, 11, 13]. The key point for this connection is that in many cases the sets of moments define algebraic varieties, hence called moment varieties. In the case of moments of mixture distributions, there is a natural correspondence to secant varieties of the moment varieties. Studying geometric invariants such as their dimension reveals properties such as model identifiability. One of the main applications for statistical inference is in the context of the method of moments, which matches the distribution’s moments to moment estimates obtained from a sample.

Gaussian mixtures are a prominent statistical model with multiple applications (see [3] and references therein). They are probability distributions on $\mathbb{R}^{n}$ with a density that is a convex combination of Gaussian densities:

[TABLE]

where $\mu_{1},\dots,\mu_{k}\in\mathbb{R}^{n}$ are the $k$ means, $\Sigma_{1},\dots,\Sigma_{k}\in\operatorname{Sym}^{2}(\mathbb{R}^{n})$ are the covariance matrices, and the $0\leq\lambda_{i}\leq 1$ with $\lambda_{1}+\dots+\lambda_{k}=1$ are the mixture weights.

The starting point is thus the Gaussian moment variety $\mathcal{G}_{n,d}$, as introduced in [4], whose points are the vectors of all moments of order at most $d$ of an $n$-dimensional Gaussian distribution. The moments corresponding to the mixture density (1) form the secant variety $\operatorname{Sec}_{k}(\mathcal{G}_{n,d})$, and identifiability in this general setting was the focus of [5].

In this work we study special families of Gaussian mixtures, called homoscedastic mixtures, where all the Gaussian components share the same covariance matrix. In other words, a homoscedastic Gaussian mixture has a density of the form

[TABLE]

where the Gaussian probability densities $f_{\mathcal{N}_{(}\mu_{i},\Sigma)}(x)$ have all different means $\mu_{i}$ and same covariance matrix $\Sigma$. The moments, up to order $d$, of homoscedastic Gaussian mixtures are still polynomials in the parameters (the means and the covariance matrix), and form the moment variety $\operatorname{Sec}^{H}_{k}(\mathcal{G}_{n,d})$. This is a set of special $k$-secants inside the secant variety $\operatorname{Sec}_{k}(\mathcal{G}_{n,d})$.

The main question we are concerned with is: when can a general homoscedastic $k$-mixture of $n$-dimensional Gaussians be identified by its moments of order $d$? More precisely, denote by $\Theta^{H}_{n,k}$ the parameter space of means, covariances and mixture weights for homoscedatic mixtures, and the moment map by

[TABLE]

The mixture parameters of a point on the moment variety $\operatorname{Sec}^{H}_{k}(\mathcal{G}_{n,d})$ can be uniquely recovered if the fiber of the moment map (3) is a singleton up to natural permutations of the parameters. If this happens for a general point on the moment variety, we say that the mixture is rationally identifiable from its moments up to order $d$. If the fiber of a general point is finite, we say that we have algebraic identifiability. The parameters are not identifiable if the general fiber of the moment map has positive dimension.

If the dimension of the parameter space is larger than the dimension of the space of moments, then one may expect any moment to lie on the moment variety. Clearly the fiber of the moment map must have positive dimension and we cannot have identifiability. We therefore distinguish the unexpected cases: when the dimension of the moment variety is less than the dimension of both the parameter space and the moment space, then we say that the moment variety $\operatorname{Sec}^{H}_{k}(\mathcal{G}_{n,d})$ is defective. In particular, defectivity implies non-identifiability.

We illustrate with an example:

Example 1**.**

Let $n=2$, $k=2$ and $d=3$. That is, we consider moments up to order three for the homoscedastic mixture of two Gaussians in $\mathbb{R}^{2}$. The Gaussian moment variety $\mathcal{G}_{2,3}$ is 5-dimensional with 2 parameters for the mean vector and 3 for the symmetric covariance matrix. The parameters for the homoscedastic mixture are two mean vectors $\mu_{1}=\begin{pmatrix}\mu_{11}\\ \mu_{12}\\ \end{pmatrix}$ and $\mu_{2}=\begin{pmatrix}\mu_{21}\\ \mu_{22}\\ \end{pmatrix}$, the common covariance $\Sigma=\begin{pmatrix}\sigma_{11}&\sigma_{12}\\ \sigma_{12}&\sigma_{22}\\ \end{pmatrix}$ and the mixture weight $\lambda$ of the first component, in total $2\times 2+3+1=8$ parameters. On the other hand, there are 9 bivariate moments up to order 3. Explicitly, the map is:

[TABLE]

Since there are more moments than parameters, one would expect that the mixture parameters can be recovered. However, the dimension of ${\rm Sec}^{H}_{2}(\mathcal{G}_{2,3})$ equals $7$. This is one less than the expected dimension of $8$. Therefore it is defective and there is no algebraic identifiability. This means that the method of moments is doomed to fail in this setting. However, if one measures moments up to order $d=4$, it is possible to uniquely recover the mixture parameters.

As is often observed [1, 4, 13], a change of coordinates to cumulants tends to yield simpler representations and faster computations. This is the case here and hence we also study the cumulant varieties of the homoscedastic Gaussian mixtures. For Example 1 above, the moment variety in cumulant coordinates is simply the cone over a twisted cubic curve (see Example 10). This is not a coincidence, as is shown in Section 3.

Our main results, Theorems 22 and 23, identify the defective homoscedastic moment varieties when $d=3$ and show that the homoscedastic moment variety is not defective when $k\leq n+1$. These are analogues of the Alexander-Hirschowitz theorem on secant-defective Veronese varieties [2].

This paper is organized as follows. In Section 2 we present the connection between moments and cumulants. The moment varieties corresponding to homoscedastic secants are defined in Section 3. In Section 4 we give general algebraic identifiability considerations and do a careful analysis of the subcases $d=3$, $k=2$ and $n=1$. Finally, we conclude with a summary of results and list further research directions.

2 Moments and Cumulants

To get started, we make some remarks about moments and cumulants from an algebraic perspective. To a sufficiently integrable random variable $X$ on $\mathbb{R}^{n}$, associate its moments $m_{a_{1},\dots,a_{n}}[X]$ and cumulants $\kappa_{a_{1},\dots,a_{n}}[X]$ through the generating functions in $\mathbb{R}[\![u_{1},\dots,u_{n}]\!]$:

[TABLE]

The information obtained from moments is equivalent to that from cumulants, since they are obtained from one another through the simple transformations

[TABLE]

which are well-defined, because the [math]-th moment is always one, whereas the [math]-th cumulant is always zero: $m_{0}[X]=1,\kappa_{0}[X]=0$ for every random variable $X$. In particular, moments and cumulants take values in the affine hyperplanes $\mathbb{A}^{M}_{n}$ and $\mathbb{A}^{K}_{n}$ of $\mathbb{R}[\![u_{1},\dots,u_{n}]\!]$ defined by

[TABLE]

We call these hyperplanes the moment space and the cumulant space.

Taking only moments up to order $d$, replace the power series ring $\mathbb{R}[\![u_{1},\dots,u_{n}]\!]$ with the truncated ring $\mathbb{R}[\![u_{1},\dots,u_{n}]\!]/(u_{1},\dots,u_{n})^{d+1}$, and everything goes through. In particular, there is an analogous definition of the affine hyperplanes $\mathbb{A}^{M}_{n,d}$ and $\mathbb{A}^{K}_{n,d}$ which we denote again by moment space and cumulant space.

Example 2** (Dirac distribution).**

Let $\mu=(\mu_{1},\dots,\mu_{n})$ in $\mathbb{R}^{n}$ be a point. The Dirac distribution $\delta_{\mu}$ with center $\mu$ on $\mathbb{R}^{n}$ is given by

[TABLE]

If $X$ is a random variable on $\mathbb{R}^{n}$ with this distribution, its moment-generating function is

[TABLE]

The moments of $X$ are monomials evaluated at $\mu$. On the other hand, for the cumulant generating function

[TABLE]

the linear cumulants coincide with the coordinates of $\mu$, and the higher order cumulants are all zero.

This has an immediate translation into algebro-geometric terms: the parameter space for all Dirac distributions is the space $\mathbb{R}^{n}$, and the image of the moment map of degree $d$, $M\colon\mathbb{R}^{n}\to\mathbb{A}^{M}_{n,d}$ is the affine $d$-th Veronese variety $V_{n,d}\subseteq\mathbb{A}^{M}_{n,d}$. On the other hand, the image of the cumulant map $K\colon\mathbb{R}^{n}\to\mathbb{A}^{K}_{n,d}$ is the linear subspace given by $\{\kappa_{2}=\kappa_{3}=\dots=\kappa_{d}=0\}$, where $\kappa_{i}$ is the degree $i$-part of an element in $\mathbb{A}^{K}_{n,d}$.

Example 3** (Gaussian distribution).**

Let $\mu\in\mathbb{R}^{n}$ be a point, and $\Sigma\in\operatorname{Sym}^{2}\mathbb{R}^{n}$ an $n\times n$ symmetric and positive-definite matrix. The Gaussian distribution on $\mathbb{R}^{n}$ with mean $\mu$ and covariance matrix $\Sigma$ is given by the density

[TABLE]

If $X\sim\mathcal{N}(\mu,\Sigma)$ is a Gaussian random variable with these parameters, its moment-generating function and cumulant-generating function are given by

[TABLE]

The Gaussian moment variety $\mathcal{G}_{n,d}\subseteq\mathbb{A}^{M}_{n,d}$ consists of all Gaussian moments up to order $d$. Observe that the corresponding cumulant variety is given simply by the linear subspace $\{\kappa_{3}=\dots=\kappa_{d}=0\}\subseteq\mathbb{A}^{K}_{n,d}$.

While our focus is on Gaussian distributions, our approach applies to general location families that admit moment and cumulant varieties. We illustrate this with the next example.

Example 4** (Laplace distribution).**

The (symmetric) multivariate Laplace distribution has a location parameter $\mu\in\mathbb{R}^{n}$ and a covariance parameter $\Sigma$, a positive-definite $n\times n$ matrix. Its density function involves the modified Bessel function of the second kind (see [12, Chapter 5]), but it can be defined via its simpler moment generating function:

[TABLE]

with radius of convergence such that $|u^{t}\Sigma u|<2$.

Moments and cumulants up to order $d=3$ match with the Gaussian case. Also note that when $\Sigma=0$, the Dirac moment generating function is recovered. However, when $d\geq 4$, the Laplace cumulants are no longer a linear space in the cumulant space.

The multiplicative structure of the power series ring $\mathbb{R}[\![u_{1},\dots,u_{n}]\!]$ makes it particularly suitable to independence statements with respect to moments. Indeed, if $X,Y$ are two independent random variables on $\mathbb{R}^{n}$ then

[TABLE]

With cumulants it is even simpler: it holds that

[TABLE]

The group of affine transformations $\operatorname{Aff}(\mathbb{R}^{n})$ acts naturally on both moments and cumulants: indeed, for any $A\in GL(n,\mathbb{R})$ and $b\in\mathbb{R}^{n}$ and a random variable $X$ on $\mathbb{R}^{n}$,

[TABLE]

and

[TABLE]

In particular, note that translations correspond simply to translations in cumulant coordinates, whereas they induce a more complicated expression in moment coordinates.

3 Homoscedastic Secants

When Karl Pearson introduced Gaussian mixtures to model subpopulations of crabs [18], he also proposed the method of moments in order to estimate the parameters. The basic idea is to compute sample moments from observed data, and match them to the distribution’s moments expressed in terms of the unknown parameters. The method of moments estimates are the parameters that solve these equations. This is a classical estimation method in statistics; a good survey is [16], and a recent ‘denoised’ version for Gaussian mixtures is [21].

The method of moments is very friendly for mixture models because computing moments of mixture densities is straightforward, since for every measurable function $g\colon\mathbb{R}^{n}\to\mathbb{R}$

[TABLE]

and thus the moments are just linear combinations of the corresponding Gaussian moments.

As hinted in the introduction, this discussion can be rephrased in geometric terms: let $\mathcal{G}_{n,d}\subseteq\mathbb{A}^{M}_{n,d}$ be the Gaussian moment variety on $\mathbb{R}^{n}$ of order $d$. Then the moments of mixtures of Gaussians are linear combinations of points in $\mathcal{G}_{n,d}$, so that their corresponding variety is the $k$-th secant variety $\operatorname{Sec}_{k}(\mathcal{G}_{n,d})$.

The densities of homoscedastic Gaussian mixtures, where the Gaussian components share a common covariance matrix, have the form:

[TABLE]

where the $\mu_{i}\in\mathbb{R}^{n}$ are the mean parameters, the $\Sigma\in\operatorname{Sym}^{2}\mathbb{R}^{n}$ is the common covariance parameters, and the $\lambda_{i}\in\mathbb{R}$ with $\lambda_{1}+\dots+\lambda_{k}=1$ are the mixture parameters. Thus, the parameter space for homoscedastic mixtures is

[TABLE]

and it has dimension

[TABLE]

The moment map for homoscedastic mixtures is then an algebraic map

[TABLE]

Points on the image, the moments of homoscedastic mixtures, are linear combinations of points in $\mathcal{G}_{n,d}\subseteq\mathbb{A}^{M}_{n,d}$ which share the same covariance matrix.

Definition 5**.**

The homoscedastic $k$-secant variety, denoted $\operatorname{Sec}^{H}_{k}(\mathcal{G}_{n,d})$, is the image of the moment map $M_{n,k,d}$. The fiber dimension $\Delta^{H}_{n,k,d}$ is the general fiber dimension of the map $M_{n,k,d}$,

[TABLE]

We say that $\operatorname{Sec}^{H}_{k}(\mathcal{G}_{n,d})$ is algebraically identifiable if $\Delta^{H}_{n,k,d}=0$.

The feasibility of the method of moments is based on computing points on the fibers of the moment map $M_{n,k,d}$. Algebraic identifiability of $\operatorname{Sec}^{H}_{k}(\mathcal{G}_{n,d})$ means that a general homoscedastic Gaussian mixture in the homoscedastic $k$-secant variety is identifiable from its moments up to order $d$ in the sense that only finitely many Gaussian mixture distributions share the same moments up to order $d$, whereas we reserve the term rationally identifiable if a general fiber consists of a single point, up to label swapping. In case the general fiber is not finite, then it is positive-dimensional, there is no identifiability of the parameters from the moments up to order $d$, and a higher order is needed for identifiability (cf. Remark 16 and [4, Problem 17]).

Since the dimension of $\operatorname{Sec}^{H}_{k}(\mathcal{G}_{n,d})$ is always bounded by the dimension of the ambient space $\mathbb{A}^{M}_{n,d}$, there is a simple estimate for the fiber dimension:

Lemma 6**.**

For all $n,d,k$ it holds that

[TABLE]

Proof.

The moment space $\mathbb{A}^{M}_{n,d}$ is an affine hyperplane inside the vector space $\mathbb{R}[[u_{1},\dots,u_{n}]]/(u_{1},\dots,u_{n})^{d+1}$, hence it has dimension

[TABLE]

Since $\operatorname{Sec}^{H}_{k}(\mathcal{G}_{n,d})\subseteq\mathbb{A}^{M}_{n,d}$ note that

[TABLE]

which is exactly the inequality in the statement. ∎

We expect that in general situations the inequality (18) is in fact an equality. Hence, define the defect to be

[TABLE]

We say that $\operatorname{Sec}^{H}_{k}(\mathcal{G}_{n,d})$ is defective if $\delta^{H}_{n,k,d}>0$. As observed earlier, defectivity implies non-identifiability.

3.1 Cumulant representation

Let us explore how homoscedastic secants become simpler in cumulant coordinates, and how this representation can be used to check identifiability.

First, rephrase the situation in terms of random variables: let $Z=Z_{\Sigma}$ be a Gaussian random variable with mean [math] and covariance matrix $\Sigma$, and let $B=B_{(\mu_{1},\dots,\mu_{k}),(\lambda_{1},\dots,\lambda_{k})}$ an independent random variable with distribution given by a mixture of Dirac distributions:

[TABLE]

Then, the random variable $Z+B$ has density given by the homoscedastic mixture (1). Moreover, if $m=\mu_{1}\lambda_{1}+\dots+\mu_{k}\lambda_{k}$ is the mean of $B$, we write $B=A+m$, where $A$ is a centered mixture of Dirac distributions.

One can compute cumulants of this random variable as follows:

[TABLE]

and this suggests to parametrize the homoscedastic secants in cumulant coordinates as follows:

[TABLE]

where $\Theta_{n,k}^{0}$ parametrizes the centered mixtures of Dirac distributions

[TABLE]

The cumulant homoscedastic secant variety $\log(\operatorname{Sec}^{H}_{k}(\mathcal{G}_{n,d}))$ is the image of the map $K$. Since in this variety, one can freely translate by the elements in $\mathbb{R}^{n}$ and $\operatorname{Sym}^{2}(\mathbb{R}^{n})$, the first cumulants and the second cumulants can take any value. The constraints are in the cumulants of order three and higher. We summarize this discussion in the following lemma.

Lemma 7**.**

Let $\mathbb{A}^{K,3}_{n,d}$ be the space of cumulants of order at least three and at most $d$, let

[TABLE]

be the cumulant map and let $C^{0}_{n,k,d}$ denote the closure $\overline{\phi_{n,k,d}(\Theta^{0}_{n,k})}$. Then the cumulant homoscedastic secant variety $\log(\operatorname{Sec}^{H}_{k}(\mathcal{G}_{n,d}))$ is a cone over $C_{n,k,d}^{0}$.

Remark 8**.**

In particular, the equations for the cumulant homoscedastic secant variety $\log(\operatorname{Sec}^{H}(\mathcal{G}_{n,d}))$ inside $\mathbb{A}^{K}_{n,d}$ are exactly the same as the equations for $C^{0}_{n,k,d}$ inside $\mathbb{A}^{K,3}_{n,d}$.

The fiber dimension $\Delta_{n,k,d}^{H}$ can also be computed as the fiber dimension of the map $\phi_{n,k,d}$:

Lemma 9**.**

The fiber dimension $\Delta^{H}_{n,k,d}$ is equal to the fiber dimension of $\phi_{n,k,d}$. In other words

[TABLE]

Proof.

The fiber dimension $\Delta^{H}_{n,k,d}$ is the difference $\dim\Theta^{H}_{n,k}-\dim\log(\operatorname{Sec}^{H}_{k}(\mathcal{G}_{n,d}))$. We know that $\Theta^{H}_{n,k}\cong\Theta^{0}_{n,k}\times\mathbb{R}^{n}\times\operatorname{Sym}^{2}\mathbb{R}^{n}$. Moreover, Lemma 7 says that $\log(\operatorname{Sec}^{H}_{k}(\mathcal{G}_{n,d}))$ is the cone over $C^{0}_{n,k,d}$, which is precisely $\mathbb{R}^{n}\times\operatorname{Sym}^{2}\mathbb{R}^{n}\times C^{0}_{n,k,d}$, so that the first equality follows. For the second equality, the dimension of $\Theta^{0}_{n,k}$ can be computed as $nk+k-1-n=(n+1)(k-1)$. ∎

Example 10** ($n=k=2$ , $d=3$).**

Revisiting Example 1 from the introduction, we concluded that ${\rm Sec}^{H}_{2}(\mathcal{G}_{2,3})\subset\mathbb{A}^{M}_{2,3}\cong\mathbb{A}^{9}$ is expected to be a hypersurface but it is actually of codimension 2. The ideal of ${\rm Sec}^{H}_{2}(\mathcal{G}_{2,3})$ is Cohen-Macaulay and determinantal (generated by the maximal minors of a $6\times 5$-matrix) as described in [4, Proposition 19]. The homoscedastic cumulant variety $\log(\operatorname{Sec}^{H}_{2}(\mathcal{G}_{2,3}))$ is defined by the vanishing of the $2\times 2$ minors of

[TABLE]

Note that indeed the first and second order cumulants $k_{10},k_{01},k_{20},k_{11},k_{22}$ do not appear in the equations above, so that the cumulant variety is the cone over the twisted cubic curve.

Remark 11**.**

To estimate the mixture parameters from the cumulants it is enough to consider the map $\phi_{n,k,d}$ of Lemma 7. Indeed, suppose that we have a homoscedastic mixture with parameters $(((\lambda_{1},\dots,\lambda_{k}),(\mu_{1},\dots,\mu_{k})),m,\Sigma)\in\Theta^{0}_{n,k}\times\mathbb{R}^{n}\times\operatorname{Sym}^{2}\mathbb{R}^{n}$ and suppose that its cumulants are known, so that in polynomial form

[TABLE]

Then to recover the parameters one can first try to recover the $\lambda_{i}$ and the $\mu_{i}$ from the cumulants of order three and higher, and then compute $m$ and $\Sigma$ from the cumulants of order one and two.

3.2 Veronese secants

We briefly observe that we can recast the above discussion in a way that makes apparent the connection to mixtures of Dirac distributions and, hence, to secants of Veronese varieties. To work with classical secant varieties, this time we work in moment coordinates. Now, every homoscedastic mixture is the distribution of a random variable of the form $Z+B$, where $B$ is a mixture of Dirac distributions and $Z$ is a centered Gaussian of covariance $\Sigma$, independent from $B$. Thus the moment generating function of this variable is

[TABLE]

Therefore the role of the covariance parameter is decoupled from the others: In particular, for $\Sigma=0$, one obtains the moment variety for mixtures of Dirac distributions. When restricting to moments $M(u)_{d}$ of degree at most $d$, this is precisely the $k$-secants to the Veronese variety $\operatorname{Sec}_{k}(\mathcal{V}_{n,d})$. The additive group $\operatorname{Sym}^{2}\mathbb{R}^{n}$ acts on the moment space $\mathbb{A}^{M}_{n,d}$ by

[TABLE]

and so (29) says that $\operatorname{Sec}^{H}_{k}(\mathcal{G}_{n,d})$ is the union of all the orbits of the points in $\operatorname{Sec}_{k}(\mathcal{V}_{n,d})$ under this action.

This is useful because we can exploit well-known results on secants of Veronese varieties to address identifiability. First, let $\Delta^{\mathcal{V}}_{n,k,d}$ denote the fiber dimension of the $k$-secants to the Veronese variety $\operatorname{Sec}_{k}(\mathcal{V}_{n,d})\subseteq\mathbb{A}^{M}_{n,d}$: by definition, this is

[TABLE]

A basic estimate for the dimension of $\operatorname{Sec}_{k}(\mathcal{V}_{n,d})$ is given by the dimension of the ambient space $\dim\mathbb{A}^{M}_{n,d}=\binom{n+d}{d}-1$, hence

[TABLE]

so that we can define the defect for the $k$-th secant of the Veronese variety as

[TABLE]

This number was famously computed by Alexander and Hirschowitz [2], see also [7]:

Theorem 12** (Alexander-Hirschowitz).**

The defect for the Veronese variety is always zero, except in the following exceptional cases

[TABLE]

Moreover, for a general point $M(u)\in\operatorname{Sec}_{k}(\mathcal{V}_{n,d})$, consider the closed subset of $\operatorname{Sym}^{2}\mathbb{R}^{n}$ given by

[TABLE]

We have the following relation between the fiber dimensions (17) and (31):

Proposition 13**.**

It holds that

[TABLE]

where $M\in\operatorname{Sec}_{k}(\mathcal{V}_{n,d})$ is a general point.

Proof.

By the previous discussion, the moment map for homoscedastic mixtures factors as a composition of two surjective maps

[TABLE]

Hence, the fiber dimension of the composite map is the sum of the fiber dimensions of the two factors. For the first one this is $\Delta^{\mathcal{V}}_{n,k,d}$, so it remains to consider the second. Denote the second factor by $\rho\colon\operatorname{Sym}^{2}(\mathbb{R}^{n})\times\operatorname{Sec}_{k}(\mathcal{V}_{n,d})\to\operatorname{Sec}^{H}_{k}(\mathcal{V}_{n,d})$ and let $(\Sigma_{o},M_{o}(u))\in\operatorname{Sym}^{2}(\mathbb{R}^{n})\times\operatorname{Sec}_{k}(\mathcal{V}_{n,d})$ be a general point. The fiber is

[TABLE]

concluding the proof. ∎

Remark 14**.**

In the range $(n+1)\left(k+\frac{n}{2}\right)\leq\binom{n+d}{d}$ where we expect identifiability for homoscedastic Gaussian mixtures, we see that $\Delta^{H}_{n,k,d}=\delta^{H}_{n,k,d}$, and Alexander-Hirschowitz says that $\Delta^{\mathcal{V}}_{n,k,d}=\delta^{\mathcal{V}}_{n,k,d}=0$. Hence Proposition 13 yields

[TABLE]

4 Moment Identifiability

Now we start to determine identifiability in various cases. To do so, it is convenient to change notation slightly. Up to now, we have identified moments and cumulants with their corresponding generating functions. In the next sections, it is useful to identify the parameters with polynomials as well. We replace the location parameter $\mu=(\mu_{1},\dots,\mu_{n})$ with the corresponding linear polynomial $u^{t}\mu=\mu_{1}u_{1}+\dots+\mu_{n}u_{n}$ and we replace the covariance parameter $\Sigma$ with the quadric $\frac{1}{2}u^{t}\Sigma u$. Of course, the two representations are equivalent, but the polynomial formalism is better suited to the cumulant space and the moment space. In particular, the linear polynomials live in the dual vector space $V=\operatorname{Hom}(\mathbb{R}^{n},\mathbb{R})$, whereas the quadratic polynomials live in $\operatorname{Sym}^{2}V$.

The next inequality reflects the fact that increasing the order of moments (or cumulants) measured results in better identifiability:

Lemma 15**.**

The fiber dimensions of general fibers of $M_{n,k,d}$ and $M_{n,k,d+1}$ satisfy:

[TABLE]

Proof.

By definition, the fiber dimension $\Delta^{H}_{n,k,d}$ is the dimension of a general nonempty fiber of the moment map $M_{n,k,d}\colon\Theta^{H}_{n,k,d}\to\mathbb{A}^{M}_{n,d}$. However, this map is the composition of the map $M_{n,k,d+1}\colon\Theta^{H}_{n,k,d}\to\mathbb{A}^{M}_{n,d+1}$ and the projection map $\mathbb{A}^{M}_{n,d+1}\to\mathbb{A}^{M}_{n,d}$, that forgets the moments of order $d+1$, so the conclusion follows. ∎

Remark 16**.**

Since Gaussian mixtures are identifiable from finitely many moments (see e.g. [4]), the sequence

[TABLE]

must stabilize at [math] for some large enough $d$.

The following observation is less trivial. It allows a reduction to the case $n=k-1$.

Proposition 17**.**

Suppose that $d\geq 3$ and $n\geq k-1$, then

[TABLE]

Proof.

Use Lemma 9, which says that the fiber dimension $\Delta^{H}_{n,k,d}$ is equal to the fiber dimension of the map

[TABLE]

This dimension can be computed by looking at the differential of the map at a general point. The parameter space is defined as

[TABLE]

Let $p=((\lambda_{1},\dots,\lambda_{k}),(L_{1},\dots,L_{k}))\in\Theta^{0}_{n,k}$ be a general point. Then the tangent space to $\Theta^{0}_{n,k}$ at the point is given by

[TABLE]

The fiber dimension of $\phi_{n,k,d}$ coincides with the dimension of the kernel of the differential $d\phi_{n,k,d}$ at the general point $p$. In particular, since the point is general and $n\geq k-1$, we can suppose that $L_{i}=u_{i}$ for $i=1,\dots,k-1$ and that all the $\lambda_{i}$ are nonzero. In particular $L_{k}$ is a linear combination of $u_{1},\dots,u_{k-1}$. Now, we claim that if $((\varepsilon_{1},\dots,\varepsilon_{k}),(H_{1},\dots,H_{k}))$ is in the kernel of $d\phi_{n,k,d}$ then the only variables appearing in the $H_{i}$ are $u_{1},\dots,u_{k-1}$. If this is true, then we are done, because the kernel of $d\phi_{n,k,d}$ coincides with the kernel of $d\phi_{k-1,k,d}$ at the point $((\lambda_{1},\dots,\lambda_{k}),(L_{1},\dots,L_{k}))\in\Theta^{0}_{k-1,k}$

To prove the claim, observe that the map is given by the cumulant functions $\phi_{n,k,d}=(\kappa_{3},\kappa_{4},\dots,\kappa_{d})$, so the kernel of $d\phi_{n,k,d}$ equals the intersection of the kernels of the $d\kappa_{i}$ for $i=3,...,d$. Therefore it is enough to prove the analogous claim for the kernel of the differential $d\kappa_{3}$ of $\kappa_{3}$. Since the first moment is zero by construction, the third cumulant coincides with the third moment

[TABLE]

Hence the differential is the linear map

[TABLE]

and if $((\varepsilon_{1},\dots,\varepsilon_{k}),(H_{1},\dots,H_{k}))$ is in the kernel, then it must be that

[TABLE]

Since $\lambda_{k}\neq 0$, this is equivalent to $\sum_{i=1}^{k}h_{i}(\lambda_{k}L_{i})^{2}=0$ and since $\lambda_{1}L_{1}+\dots+\lambda_{k}L_{k}=0$, we see that

[TABLE]

By assumption $L_{i}=u_{i}$ for $i=1,\dots,k-1$, so this last expression is equal to zero if and only if

[TABLE]

If this is true, then $h_{k}$ uses only the variables $u_{1},\dots,u_{k-1}$. Indeed, if some other variable, say $y$, appears in $h_{k}$ then on the right hand side there is the monomial $yu_{1}u_{2}$, while there is no such a monomial on the left hand side. Likewise, if the variable $y$ appears in one of the $h_{i}$ for $i=1,\dots,k-1$: then on the left hand side there would be a monomial of the form $yu_{i}^{2}$, while there is no such monomial on the right hand side.

Hence, the $h_{i}$ are polynomials in the $u_{1},\dots,u_{k}$, and, by definition of the $h_{i}$, it follows that the same holds for the $H_{i}$. This proves the claim and the result follows. ∎

4.1 Moments up to order $d=3$

When $d=3$ we determine the defect ${\delta}^{H}_{n,k,3}$ and the fiber dimension ${\Delta}^{H}_{n,k,3}$ of the map

[TABLE]

for each $n$ and $k$, and use Lemma 9. When $d=3$, the space $\mathbb{A}^{K,3}_{n,3}$ is identified with the space $\operatorname{Sym}^{3}V$ of homogeneous polynomials of degree three, and as noted in the proof of Proposition 17, the third cumulants coincide with the third moments, so that:

[TABLE]

We compute the closure $C^{0}_{n,k,3}$ of the image.

Lemma 18**.**

The set $C^{0}_{n,k,3}$ is the Zariski closure of

[TABLE]

Proof.

Recall that

[TABLE]

To compute the Zariski closure, suppose that all the $\lambda_{i}$ are strictly positive, so that in particular we can write

[TABLE]

Since cubic roots are well defined over $\mathbb{R}$,

[TABLE]

where $H_{i}:=\sqrt[3]{\lambda_{i}}L_{i}$ for $i=1,\dots,k-1$, and $H_{k}:=-\sum_{i=1}^{k-1}\left(\frac{\sqrt[3]{\lambda_{i}}}{\sqrt[3]{\lambda_{k}}}\right)^{2}H_{i}$, using the equality $\sqrt[3]{\lambda_{k}}\frac{\lambda_{i}}{\lambda_{k}}=\left(\frac{\sqrt[3]{\lambda_{i}}}{\sqrt[3]{\lambda_{k}}}\right)^{2}\sqrt[3]{\lambda_{i}}$. In particular, this shows immediately that $\lambda_{1}L_{1}^{3}+\dots+\lambda_{k}L_{k}^{3}$ can be written as a sum of cubic powers of linearly dependent linear forms.

For the converse, let $H_{1},\dots,H_{k}$ be linearly dependent linear forms. For the Zariski closure, it suffices to assume that $H_{k}=-\beta_{1}H_{1}-\dots-\beta_{k-1}H_{k-1}$ for some general $\beta_{1},\dots,\beta_{k-1}\in\mathbb{R}$ strictly positive. So we want to write

[TABLE]

for some positive $\lambda_{1},\dots,\lambda_{k}\in\mathbb{R}$ such that $\lambda_{1}+\dots+\lambda_{k}=1$. Given such $\lambda_{i}$, the above computations yields

[TABLE]

where $L_{i}=\frac{1}{\sqrt[3]{\lambda_{i}}}H_{i}$ for $i=1,\dots,k-1$ and $L_{k}=-\frac{\lambda_{1}}{\lambda_{k}}L_{1}-\dots-\frac{\lambda_{k-1}}{\lambda_{k}}L_{k-1}$, so that $\lambda_{1}L_{1}+\dots+\lambda_{k}L_{k}=0$, as wanted.

To conclude, it remains to show that the equations (51) have a solution: these equations are equivalent to

[TABLE]

Observe that the square roots are well defined since $\beta_{i}>0$ for all $i=1,\dots,k-1$. Moreover, if $(\lambda_{1},\dots,\lambda_{k-1})$ is a solution to (53), then it is easy to see that all the $\lambda_{i}$ must be strictly positive: indeed, since the $\beta_{i}$ are positive, $\lambda_{i}$ and $1-\lambda_{1}-\dots-\lambda_{k-1}$ have the same sign. Thus, if one of the $\lambda_{i}$ is negative, then all the $\lambda_{i}$ are negative, but then $1-\lambda_{1}-\dots-\lambda_{k-1}>0$ which is absurd.

Now, setting $b_{i}=\sqrt{\beta_{i}}^{3}$, rewrite the equations as the linear system

[TABLE]

The matrix determinant lemma gives that $\det(\mathrm{I}+b\cdot\mathbb{1}^{T})=1+\mathbb{1}^{T}b$ = $1+b_{1}+\dots+b_{k-1}$, which is positive since the $\beta_{i}$ are positive. This means the system (54) has a unique solution. ∎

Remark 19**.**

The proof of Lemma 18 actually gives more: indeed, it shows that the image of the positive part

[TABLE]

which is the one relevant in statistics, coincides with the set of sums $\{H_{1}(u)^{3}+\dots+H_{k}(u)^{3}\}$, where the $H_{i}$ are positively linearly dependent, meaning that there are coefficients $\beta_{1},\dots,\beta_{k}>0$ such that

[TABLE]

Remark 20**.**

The set of sums of cubes of $k$ dependent linear forms has a natural interpretation in terms of the projective Veronese variety: indeed consider the third Veronese embedding of $\mathbb{P}(V)=\mathbb{P}^{n-1}$:

[TABLE]

For each $(k-2)$-dimensional linear subspace $\Pi\subseteq\mathbb{P}^{n-1}$ let $\operatorname{Sec}_{k}(v_{3}(\Pi))\subseteq\mathbb{P}(\operatorname{Sym}^{3}V)$ be the $k$-th secant variety of its image $v_{3}(\Pi)$. Then, by Lemma 18, the variety $C^{0}_{n,k,3}$ is the affine cone over the union of these secants:

[TABLE]

We compute the dimension of this variety, dividing it in the cases $k\leq n+1$ and $k\geq n+1$:

Proposition 21**.**

1. i)

If $k\geq n+1$, then

[TABLE]

except in the case $n=5,k=7$, where $\dim C^{0}_{5,7,3}=34$. 2. ii)

If $k\leq n+1$, then

[TABLE]

and when $k\geq 5$,

[TABLE]

Proof.

i) Since $k\geq n+1$, Remark 20 shows that $C^{0}_{n,k,3}$ is the cone over the $k$-th secant variety $\operatorname{Sec}_{k}(v_{3}(\mathbb{P}^{n-1}))$. The dimension of this variety is computed by the Alexander-Hirschowitz theorem, so that

[TABLE]

with the single exception of $n=5,k=7$, where the dimension is one less than the expected, hence $\dim C^{0}_{5,7,3}=34$.

ii) Since $k\leq n+1$, Proposition 17 shows that $\Delta^{H}_{n,k,3}=\Delta^{H}_{k-1,k,3}$. Hence, for $k=2,3,4$ we see directly from Table 1 that

[TABLE]

For $k\geq 5$ instead, follow the proof of Proposition 17 and show that the differential of $\phi_{k-1,k,3}\colon\Theta^{0}_{k-1,k,3}\to\operatorname{Sym}^{3}V$ at a general point is injective. For this, consider the kernel of the differential at a point $p=((\lambda_{1},\dots,\lambda_{k}),(L_{1},\dots,L_{k}))$. It consists of elements $((\varepsilon_{1},\dots,\varepsilon_{k}),(H_{1},\dots,H_{k}))\in\mathbb{R}^{k}\times V^{k}$ such that $\varepsilon_{1}+\dots+\varepsilon_{k}=0,\varepsilon_{1}L_{1}+\dots+\varepsilon_{k}L_{k}+\lambda_{1}H_{1}+\dots+\lambda_{k}H_{k}=0$ and

[TABLE]

where $\ell_{i}=\lambda_{k}^{2}(3\lambda_{i}H_{i}+\varepsilon_{i}L_{i})+\lambda_{i}^{2}(3\lambda_{k}H_{k}+\varepsilon_{k}L_{k})$ and $h=3\lambda_{k}H_{k}+\varepsilon_{k}L_{k}$. Now choose the specific point $p$ given by $\lambda_{i}=\frac{1}{k}$ for each $i=1,\dots,k$, $L_{i}=u_{i}$ for $i=1,\dots,k-1$ and $L_{k}=-u_{1}-\dots-u_{k-1}$. Then the above equation becomes

[TABLE]

Let us write $h=h_{1}u_{1}+\dots+h_{k-1}u_{k-1}$. Then, in (65), the coefficient of $u_{a}u_{b}u_{c}$ is $\frac{2}{k^{2}}(h_{a}+h_{b}+h_{c})$ for all $1\leq a<b<c\leq k-1$. Hence

[TABLE]

Let $1\leq a<b<c<d\leq k-1$ be any four distinct indices between $1$ and $k-1$. Then the previous equations translate into the linear system

[TABLE]

The matrix appearing in the linear system is invertible, so $h_{a}=h_{b}=h_{c}=h_{d}=0$. Since this holds for an arbitrary choice of four distinct indices, it follows that $h=0$. Now, relation (65) tells us that $\sum_{i=1}^{k-1}u_{i}^{2}\ell_{i}=0$, but since $u_{1}^{2},\dots,u_{k-1}^{2}$ form a complete intersection of quadrics, they do not have linear syzygies, which implies that $\ell_{i}=0$ for each $i$. From the definitions of $\ell_{i}$ and $h$, it follows that $3\lambda_{i}H_{i}+\varepsilon_{i}L_{i}=0$ for each $i$ but then the other two relations $\sum_{i}\varepsilon_{i}=0$ and $\sum_{i}(\lambda_{i}H_{i}+\varepsilon_{i}L_{i})=0$ imply that $H_{i}=0,\varepsilon_{i}=0$ for all $i$, which is what was needed. ∎

Now we are ready for a complete classification of defectivity when $d=3$:

Theorem 22**.**

For $d=3$, the defect $\delta^{H}_{n,k,3}=0$ for any $k$ and $n$, with the following exceptions:

• •

$n\geq k$* and $k=2$, where $\delta^{H}_{n,2,3}=1$.*

• •

$n\geq k$* and $k=3,4$, where $\delta^{H}_{n,k,3}=2$.*

• •

$n=5$* and $k=7$, where $\delta^{H}_{5,7,3}=1$.*

• •

$n\geq 4$* and $n+1<k\leq\frac{n^{2}+2n+6}{6}$ where $\delta^{H}_{n,k,3}=k-n-1$.*

• •

$n\geq 4$* and $\frac{n^{2}+2n+6}{6}\leq k<\frac{n^{2}+3n+2}{6}$ where $\delta^{H}_{n,k,3}=n\left(\frac{n^{2}+3n+2}{6}-k\right)$.*

Proof.

First consider the case when $n\geq k$: then Proposition 21 ii applies. It is straightforward to check that $\delta^{H}_{n,k,3}=\Delta^{H}_{n,k,d}$, from which the statement of the theorem follows.

For the cases where $k\geq n+1$, start with the exceptional case $n=5,k=7$: Proposition 21 i gives that $\dim\overline{\phi_{n,k,3}(\Theta^{0}_{5,7})}=34$, and the Lemma 9 yields $\Delta^{H}_{5,7,3}=2$ and $\delta^{H}_{5,7,3}=1$.

Now, consider the other cases: Proposition 21 i gives that

[TABLE]

and then Lemma 9 shows that

[TABLE]

so that

[TABLE]

Suppose first that $n=1,2,3$: this implies $k\geq n+1\geq\frac{n^{2}+2n+6}{6}\geq\frac{n^{2}+3n+2}{6}$ so that

[TABLE]

Now, suppose that $n\geq 4$. Then $5\leq n+1\leq\frac{n^{2}+2n+6}{6}\leq\frac{n^{2}+3n+2}{6}$ and there are three possibilities for $k$: if $k\geq\frac{n^{2}+3n+2}{6}$, then

[TABLE]

If instead $\frac{n^{2}+2n+6}{6}\leq k<\frac{n^{2}+3n+2}{6}$, then

[TABLE]

which is strictly positive. Finally, if $n+1\leq k<\frac{n^{2}+2n+6}{6}$ the defect is

[TABLE]

which is positive if and only if $k>n+1$. ∎

As a consequence, identifiability can be characterized whenever $k\leq n+1$:

Theorem 23**.**

Suppose $k\leq n+1$. If $k\geq 5$ then a general homoscedastic mixture is algebraically identifiable from moments up to order $3$. If instead $k=2,3,4$ then a general homoscedastic mixture is algebraically identifiable from the moments up to order $d=4$.

Proof.

When $k\geq 5$ this follows immediately from Theorem 22 and Lemma 15. If instead $k=2,3,4$, thanks to Proposition 17, it is enough to set $n=k-1$ and check the first $d$ for which we have identifiability: these are a finite number of cases that can be done by direct computation (e.g. in Macaulay2 [9]), and we find that such a $d$ is 4. ∎

4.2 Mixtures with $k=2$ components

When $k=2$ we characterize the rational identifiability as well. Since the case $d=3$ is already covered, consider only $d\geq 4$.

Theorem 24**.**

The homoscedastic secant $\operatorname{Sec}^{H}_{2}(\mathcal{G}_{n,4})$ is algebraically identifiable. If $d\geq 5$, the homoscedastic secant $\operatorname{Sec}^{H}_{2}(\mathcal{G}_{n,d})$ is also rationally identifiable.

Proof.

By Lemma 9 and Remark 11, it is enough to consider the parameter space given by $\Theta^{0}_{n,2}=\{((L_{1},L_{2}),(\lambda_{1},\lambda_{2}))\,|\,\lambda_{1}+\lambda_{2}=1,\lambda_{1}L_{1}+\lambda_{2}L_{2}=0\}$ and the map

[TABLE]

In order to compute the general fiber of this map, note that since $d\geq 4$, it follows from Theorem 23 and its proof that the map has finite fibers. Hence, it is enough to restrict a general fiber to the open subset $\lambda_{2}\neq 0$. There we may assume $L_{2}=-\frac{\lambda_{1}}{\lambda_{2}}L_{1}=-\frac{\lambda_{1}}{1-\lambda_{1}}L_{1}$. We thus compute the fibers of the induced map

[TABLE]

In explicit terms, this map is given by the terms from degree 3 to degree $d$ of the logarithm $\log(\lambda e^{L}+(1-\lambda)e^{-\frac{\lambda}{\lambda-1}L})$. A computation shows that the first terms are:

[TABLE]

[TABLE]

[TABLE]

Now suppose that $d=4$, and let $L\in V$ and $\lambda\in\mathbb{R}\setminus\{1\}$ be general elements. In fact, it is enough to assume $L\neq 0$ and $\lambda\neq 0,1,\frac{1}{2}$, so that $\kappa_{3}=f_{3}(\lambda)L^{3}\neq 0$. In order to compute the fiber of the point $(\kappa_{3},\kappa_{4})=F_{n,2,4}(L,\lambda)$, first observe that $\kappa_{3}=f_{3}(\lambda_{0})L_{0}^{3}=(\sqrt[3]{f_{3}(\lambda_{0})}L_{0})^{3}$ and that the polynomial $L_{0}:=\sqrt[3]{f_{3}(\lambda)}L$ can be computed explicitly: from the expression

[TABLE]

then one obtains

[TABLE]

In particular, $L=f_{3}(\lambda)^{-\frac{1}{3}}L_{0}$, so that the equation $\kappa_{4}=f_{4}(\lambda)L^{4}$ translates into $\frac{f_{4}(\lambda)}{f_{3}(\lambda)^{\frac{4}{3}}}=\frac{\kappa_{4}}{L_{0}^{4}}$. Observe that $a:=\frac{\kappa_{4}}{L_{0}^{4}}$ is a constant that can be computed explicitly by comparing a single nonzero coefficient of $L_{0}^{4}$ with the corresponding coefficient of $\kappa_{4}$: for example, if $\sqrt[3]{\kappa_{300..0}}\neq 0$, then

[TABLE]

Now, the equation $\frac{f_{4}(\lambda)}{f_{3}(\lambda)^{\frac{4}{3}}}=a$ is equivalent to $\frac{f_{4}(\lambda)^{3}}{f_{3}(\lambda)^{4}}=a^{3}$, or more explicitly

[TABLE]

Note that this expression is invariant under exchanging $\lambda$ with $1-\lambda$, as is expected from the symmetry of the situation. Hence, set $\gamma:=\lambda(1-\lambda)$ and rewrite this expression as

[TABLE]

This is a cubic equation with three possible solutions for $\gamma$, which means there is no rational identifiability. In order to get such, consider also the cumulants $\kappa_{5}$ of order 5: this adds the data $\kappa_{5}$ and the condition $\kappa_{5}=f_{5}(\lambda)L^{5}$. In the above notation $L=f_{3}(\lambda)^{-\frac{1}{3}}L_{0}$, so that the condition $\kappa_{5}=f_{5}(\lambda)L^{5}$ becomes $\frac{f_{5}(\lambda)}{f_{3}(\lambda)^{\frac{5}{3}}}=\frac{\kappa_{5}}{L_{0}^{5}}$. As before, we see that $b:=\frac{\kappa_{5}}{L_{0}^{5}}$ is a constant that can be computed explicitly by comparing a single nonzero coefficient of $L_{0}^{5}$ with the corresponding coefficient of $\kappa_{5}$: for example, if $\sqrt[3]{\kappa_{300..0}}\neq 0$, then

[TABLE]

Now, the equation $\frac{f_{5}(\lambda)}{f_{3}(\lambda)^{\frac{5}{3}}}=a$ is equivalent to $\frac{f_{5}(\lambda)^{3}}{f_{3}(\lambda)^{5}}=b^{3}$, or more explicitly, as above, with the substitution $\gamma=\lambda(1-\lambda)$,

[TABLE]

Hence, rational identifiability is obtained if the two equations (80) and (82) have a unique common solution $\gamma$. This means that the map $\mathbb{R}\dashrightarrow\mathbb{R}^{2},\gamma\mapsto(g(\gamma),h(\gamma))$ is generically injective. This map extends to

[TABLE]

i.e. a map defined by polynomials of degree $7$. It is generically injective if and only if the closure of its image is a plane curve of degree $7$. This can be verified with Macaulay2 [9]: the resulting curve is given by the equation

[TABLE]

∎

Even though there is no rational identifiability above when $d=4$, it is worth noting that in a purely statistical setting, $\gamma$ can be recovered uniquely, as seen below.

Corollary 25**.**

For $k=2$, the statistical mixture parameters can be recovered uniquely with moments up to order $d=4$.

Proof.

This is equivalent to saying that the equation (80) has a unique statistically relevant solution in $\gamma=\lambda(1-\lambda)$. Note that since $\lambda\in(0,1)\setminus\{\frac{1}{2}\}$, we have that $\gamma\in(0,\frac{1}{4})$. Consider the real valued function coming from (80):

[TABLE]

Its derivative, $a^{\prime}(\gamma)=-\frac{1}{2\sqrt[3]{36}\gamma(1-4\gamma)\sqrt[3]{4\gamma(1-4\gamma)^{2}}}$, is always negative for $0<\gamma<\frac{1}{4}$ so that the function $a(\gamma)$ is strictly decreasing and, in particular, injective in this statistically meaningful interval. The corresponding inverse is given by the cubic equation in $\gamma$

[TABLE]

The discriminant of (84) is $\Delta=-3072a^{6}(64a^{3}+81)$. It is zero precisely when $a=-\frac{3\sqrt[3]{3}}{4}$, which corresponds to the horizontal asymptote of $a$. If $a<-\frac{3\sqrt[3]{3}}{4}$, there are 3 real solutions, but one is negative and the other one is larger than $\frac{1}{4}$. The remaining solution is also the unique real solution when $a>-\frac{3\sqrt[3]{3}}{4}$, given explicitly by

[TABLE]

where $\eta=(-4096a^{9}-10368a^{6}-6561a^{3}+9\sqrt{262144a^{15}+995328a^{12}+1259712a^{9}+531441a^{6}})^{\frac{1}{3}}$. ∎

This proof gives an explicit algorithm to recover the parameters of a homoscedastic mixture of two Gaussians from the cumulants up to order four.

Observe that this algorithm needs all the cumulants of order one, all the cumulants of order two, $n$ cumulants of order three, and one cumulant of order four. Hence, it needs in total $n+\frac{n(n+1)}{2}+n+1$ cumulants.

Remark 26**.**

We have seen in Remark 20 that $\operatorname{Sec}^{H}_{2}(\mathcal{G}_{n,d})$ in cumulant coordinates is a cone over $C^{0}_{n,2,d}\subseteq\mathbb{A}^{K,3}_{n,d}$. Up to taking the Zariski closure, the proof of Theorem 24 shows that $C^{0}_{n,2,d}$ is the image of the map

[TABLE]

*For $\lambda$ constant we get a projected $d$-th Veronese variety of $V$. If instead $L$ is constant, then we get a rational curve given by a linear combination of $(f_{3}(\lambda),f_{4}(\lambda),\dots,f_{d}(\lambda))$. *

4.3 The univariate case $n=1$

We use the standard notation $\sigma^{2}$ for the variance $\Sigma=(\sigma_{11})$ when $n=1$.

For $n=1$, the moment variety $\operatorname{Sec}_{k}^{H}(\mathcal{G}_{1,d})$ is never defective. The moment map

[TABLE]

is finite to one. In the statistics literature it is known that in the case of homoscedastic secants, one may recover mixture parameters from given moments (i.e. compute the fiber of the map above), with an algorithm closely related to the well-known Prony’s method [20]. This procedure was introduced by Lindsay as an application of moment matrices [15] and we briefly recall the algorithm here.

First, how does one recover the locations $\mu_{i}$ and weights $\lambda_{i}$ of the $k$ components of a Dirac mixture from $2k-1$ moments? This is known as the quadrature rule and it works as follows. Given the moment sequence $m=(m_{1},m_{2},\dots,m_{2k-1})$ one considers the polynomial resulting from the following $(k+1)\times(k+1)$ determinant

[TABLE]

The $k$ roots $\mu_{1},\mu_{2},\dots,\mu_{k}$ of $P_{k}(t)$ are precisely the sought locations. This follows since the equations of the secant varieties of the rational normal curve are classically known to be given by the minors of the moment matrices. For a modern reference see [14].

Once the locations are known, the weights $\lambda_{i}$ are found by solving the $k\times k$ Vandermonde linear system

[TABLE]

Back to the Gaussian case, if we knew the value of the common variance $\sigma^{2}$, we can reduce to the above instance. In terms of the Gaussian moment generating function:

[TABLE]

Hence, the Dirac moments $\tilde{m}$ on the right hand side are linear combinations of the Gaussian moments $m$. Explicitly, for $1\leq j\leq 2k-1$

[TABLE]

Applying the quadrature rule to the vector $\tilde{m}=(\tilde{m}_{1},\tilde{m}_{2},\dots,\tilde{m}_{2k-1})$ would allow us to obtain the means $\mu_{1},\mu_{2},\dots,\mu_{k}$.

However, $\sigma$ is unknown. To find an estimate for $\sigma$ we consider the first $2k$ moments $m=(m_{1},m_{2},\dots,m_{2k})$. If $\tilde{m}=(\tilde{m}_{1},\tilde{m}_{2},\dots,\tilde{m}_{2k})$ comes from a mixture of $k$ Dirac measures, then

[TABLE]

One thus treats $\sigma$ as a variable and substitutes expressions (90) into (91). This results in a polynomial $D_{k}(\sigma)$ of degree $\binom{k+1}{2}$ in $\sigma^{2}$ and the estimator $\hat{\sigma}^{2}$ is obtained as its smallest non-negative root [15, Theorem 5B]. So the algebraic degree for estimating $\sigma^{2}$ is $\binom{k+1}{2}$. With $\sigma^{2}$ specified, one proceeds as above.

More generally, the discussion under (29) shows that the moment variety $\operatorname{Sec}_{k}^{H}(\mathcal{G}_{1,d})$ with $k\leq d/2$ is a union

[TABLE]

where $V_{1,d}^{\sigma}$ is the translation of the moment curve $V_{1,d}$ by the variance $\sigma^{2}$ as defined by the Gaussian moments. The secant variety $\operatorname{Sec}_{k}(V_{1,d}^{\sigma})$ is defined for each $\sigma$ by the $(k+1)\times(k+1)$ minors of

[TABLE]

As soon as the $k$-th secant variety of a smooth curve is not linear, the curve can be recovered as the singular locus of highest multiplicity in the secant variety. Therefore, since curves $V_{1,d}^{\sigma}$ are distinct, their $k$-th secant varieties are distinct as well, as long as the latter are not linear. In particular, since the variety $\operatorname{Sec}_{k}(V_{1,d}^{\sigma})$ has dimension $2k-1$, it follows that the union $\operatorname{Sec}_{k}^{H}(\mathcal{G}_{1,d})$ has dimension $2k$. Given the moments $m_{i}$ up to degree $d$ of a point on a homoscedastic $k$-secant, the $(k+1)\times(k+1)$ minors of $M_{k,d}$ are polynomials in $\sigma^{2}$ with a zero at the common variance. Given the variance, the means can be inferred as above.

When $d=2k+1$, then the variety $\operatorname{Sec}_{k}^{H}(\mathcal{G}_{1,d})\subset\mathbb{A}^{M}_{1,2k+1}$ is a hypersurface, defined by the resultant of $(k+1)$-minors of $M_{k,d}$, the polynomial obtained by elimination of $\sigma^{2}$ in the ideal defined by the $(k+1)\times(k+1)$ minors. Denote this polynomial by $P_{2k+1}.$ It is a polynomial in $m_{1},...,m_{2k+1}$ (or $\kappa_{3},\kappa_{4},\dots,\kappa_{2k+1}$). For example,

[TABLE]

[TABLE]

Proposition 27**.**

The polynomial $P_{2k+1}$ is homogeneous of total degree

[TABLE]

in the multigraded weights $\deg m_{i}=\deg\kappa_{i}=i$.

Proof.

Let

[TABLE]

where $\sigma$ is the last coordinate, and consider the projective closure $\mathbb{P}$ of $\mathbb{A}$. Then the matrix (92) defines a map between vector bundles $E$ and $F$ on $\mathbb{A}$. The vector bundles $E$ and $F$ and the map extends to $\mathbb{P}$; $E$ extends to a sum of line bundles $\tilde{E}={\cal O}_{\mathbb{P}}\oplus{\cal O}_{\mathbb{P}}(-1)\oplus...,{\cal O}_{\mathbb{P}}(-k)$, while F extends to a sum of line bundles $\tilde{F}={\cal O}_{\mathbb{P}}\oplus{\cal O}_{\mathbb{P}}(1)\oplus...,{\cal O}_{\mathbb{P}}(k+1)$. By the Thom-Porteous formula, [8, Theorem 14.4], the degree in $\mathbb{P}$ of the rank $k$ locus of the map is given by the Chern class

[TABLE]

since the Chern polynomials of $\tilde{E}$ and $\tilde{F}$ in are

[TABLE]

and

[TABLE]

This rank $k$ locus has codimension $2$ and its intersection with $\mathbb{A}$ is projected to the hypersurface defined by $P_{2k+1}$ in $\mathbb{A}^{M}_{1,2k+1}$. The coordinate $\sigma$ appears only in even degree in the equations defining the rank $k$ locus, so the projection to $\mathbb{A}^{M}_{1,2k+1}$ is $2:1$, so the degree of $P_{2k+1}$ is half the degree of the rank $k$ locus. ∎

Question 28**.**

It would be interesting to understand better the structure of the polynomials $P_{2k+1}$, e.g. is there a closed form expression for all $k$?

If $P_{2k+1}$ vanishes on a the set $(m_{1},...,m_{2k+1})$ of moments, and $P_{2l+1}$ does not vanish on $(m_{1},...,m_{2l+1})$ for any $l<k$, then the moments lie on a homoscedastic $k$-secant but not on any $l$ secant for $l<k$. Therefore the polynomials $P_{2k+1}$ may be used to estimate the number of components in a homoscedastic Gaussian mixture (compare to the rank test proposed in [15, Section 3.1] for the known variance case).

5 Conclusion

We have completely classified all defective cases for the moment varieties associated to homoscedastic Gaussian mixtures whenever $k<n+1$, $d=3$, $k=2$ or $n=1$. The question concerning a complete classification for all $n,d,k$ remains open, although our computations did not reveal any further defective examples.

Our identifiability results also cover special structures in the covariance matrix, by Remark 11. For example, a common mixture submodel involves isotropic Gaussians, which means that the covariance matrix is a scalar multiple of the identity, $\Sigma=\sigma I$. The $k$-means algorithm used in clustering can be interpreted as parameter estimation for a homoscedastic isotropic mixture of Gaussians. In [10], Hsu and Kakade consider the learning of mixtures of isotropic Gaussians from the moments up to order $d=3$ when $k\leq n+1$. They prove identifiability for the homoscedastic isotropic submodel (see [6, Theorem 3.2]), and in order to solve the moment equations, they find orthogonal decompositions of the second and third order moment tensors.

On the other hand, in [17] Lindsay and Basak proposed a ‘fast consistent’ method of moments for homoscedastic Gaussian mixtures in the multivariate case, based on a ‘primary axis’ to which the one-dimensional case presented in Section 4.3 is applied. This means that the method uses some moments of order $2k$. Knowing that in some cases there are explicit equations for secants of higher dimensional Veronese varieties [14], an alternative method with minimal order based on these should be possible.

Finally, a similar approach can be made to study moment varieties of homoscedastic mixtures of other location families. In the case of Example 4, we saw that Gaussian moments and Laplacian moments coincide up to $d=3$. This means that Theorem 22 applies verbatim to homoscedastic mixtures of Laplace distributions.

Acknowledgments. The authors are grateful to the Max Planck Institute for Mathematics in the Sciences, Leipzig and the Institute for Computational and Experimental Research in Mathematics in Providence, RI for facilitating discussions about this work. Carlos Améndola was partially supported by the Deutsche Forschungsgemeinschaft (DFG) in the context of the Emmy Noether junior research group KR 4512/1-1. We thank anonymous referees for suggestions to improve the presentation.