{'fraction_non_alphanumeric': 0.05762135813787089, 'fraction_numerical': 0.014153791114191802, 'mean_word_length': 3.5441983336720178, 'pattern_counts': {'":': 0, '<': 2, '<?xml version=': 0, '>': 5, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 29, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Marginal likelihood and model selection for Gaussian latent tree and forest models

Gaussian latent tree models, or more generally, Gaussian latent forest models have Fisher-information matrices that become singular along interesting submodels, namely, models that correspond to subforests. For these singularities, we compute the real log-canonical thresholds (also known as stochastic complexities or learning coefficients) that quantify the large-sample behavior of the marginal likelihood in Bayesian inference. This provides the information needed for a recently introduced generalization of the Bayesian information criterion. Our mathematical developments treat the general setting of Laplace integrals whose phase functions are sums of squared differences between monomials and constants. We clarify how in this case real log-canonical thresholds can be computed using polyhedral geometry, and we show how to apply the general theory to the Laplace integrals associated with Gaussian latent tree and forest models. In simulations and a data example, we demonstrate how the mathematical knowledge can be applied in model selection.


Introduction
Graphical models based on trees are particularly tractable, which makes them useful tools for exploring and exploiting multivariate stochastic dependencies, as first demonstrated by [CL68]. More recent work develops statistical methodology for extensions that allow for inclusion of latent variables and in which the graph may be a forest, that is, a union of trees over disjoint vertex sets [CTAW11,TAW11,MRS13]. These extensions lead to a new difficulty in that the Fisher-information matrix of a latent tree model is typically singular along submodels given by subforests. As explained in [Wat09], such singularity invalidates the mathematical arguments that lead to the Bayesian information criterion (BIC) of [Sch78], which is widely used to guide model selection algorithms that infer trees or forests [EdAL10]. Indeed, the BIC will generally no longer share the asymptotic behavior of Bayesian methods; see also [DSS09,Sect. 5.1]. Similarly, Akaike's information criterion may no longer be an asymptotically unbiased estimator of the expected Kullback-Leibler divergence that it is designed to approximate [Wat09,Wat10a,Wat10b].
In this paper, we study the large-sample behavior of the marginal likelihood in Bayesian inference for Gaussian tree/forest models with latent variables, with the goal of obtaining the mathematical information needed to evaluate a generalization of BIC proposed in [DP13]. As we review below, this information comes in the form of so-called real log-canonical thresholds (also known as stochastic complexities or learning coefficients) that appear in the leading term of an asymptotic expansion of the marginal likelihood. We begin by more formally introducing the models that are the object of study.
Let Z = (Z u ) u∈U be a random vector whose components are indexed by the vertices of an undirected tree T = (U, E) with edge set E. Via the paradigm of graphical modeling [Lau96], the tree T induces a Gaussian tree model N(T ) for the joint distribution of Z. The model N(T ) is the collection of all multivariate normal distributions on R U under which Z u and Z v are conditionally independent given Z C = (Z w : w ∈ C) for any choice of two nodes u, v and a set C ⊂ U \ {u, v} such that C contains a node on the (unique) path between u and v. For two nodes u, v ∈ U , let uv be the set of edges on the path between u and v. It can be shown that a normal distribution with correlation matrix R = (ρ uv ) belongs to N(T ) if and only if where ρ e := ρ xy when e is the edge incident to x and y. Indeed, for three nodes v, w, u ∈ U the conditional independence of Z v and Z w given Z u is equivalent to ρ vw = ρ uv ρ uw ; compare also [MRS13,p. 4359].
In this paper, we are concerned with latent tree models in which only the tree's leaves correspond to observed random variables. So let V ⊂ U be the set of leaves of tree T = (U, E). Then the Gaussian latent tree model M(T ) for the distribution of the subvector X := (Z v : v ∈ V ) is the set of all V -marginals of the distributions in N(T ). The object of study in our work is the parametrization of the model M(T ). Without loss of generality, we may assume that the latent variables Z a at the inner nodes a ∈ U \ V have mean zero and variance one. Moreover, we assume that the observed vector X has mean zero. Then, based on (1.1), the distributions in M(T ) can be parametrized by the variances ω v for each variable X v , v ∈ V , and the edge correlations ω e , e ∈ E.
Our interest is in the marginal likelihood of model M(T ) when the variance and correlation parameters are given a prior distribution with smooth and everywhere positive density. Following the theory developed by [Wat09], we will derive largesample properties of the marginal likelihood by studying the geometry of the fibers (or preimages) of the parametrization map.
Hence, the fiber of the parametrization is finite, containing two points. If instead the correlations between the leaves are zero then this identifiability breaks down. If R is the identity matrix with ρ 12 = ρ 13 = ρ 23 = 0, then every vector (ω a1 , ω a2 , ω a3 ) ∈ [−1, 1] 3 that lies in the set satisfies (1.2). The fiber of the identity matrix is thus the union of three line segments that form a one-dimensional semi-algebraic set with a singularity at the origin where the lines intersect.
Remark 1.2. Some readers may be more familiar with rooted trees with directed edges and model specifications based on the Markov properties for directed graphs or structural equations. However, these are equivalent to the setup considered here, as can be seen by applying the so-called trek rule [SGS00]. Our later results also apply to Bayesian inference in graphical models associated with directed trees.
Suppose ϕ is a smooth and positive density that defines a prior distribution on the parameter space Ω = (0, ∞) V × [−1, 1] E of the Gaussian latent tree model M(T ). Let X n = (X (1) , . . . , X (n) ) be a sample consisting of n independent and identically distributed random vectors in R V , and write L(M(T )|X n ) for the marginal likelihood of M(T ). If X n is generated from a distribution q ∈ M(T ) and n → ∞, then it holds that where λ q ≥ 0 is a rational number smaller than or equal to the dimension of the model M(T ). The number m q is an integer greater than or equal to 1. More detail on how (1.3) follows from results in [Wat09] is given in Section 2. In this paper, we derive formulas for the pair (λ q , m q ) from (1.3), which will be seen to depend on the pattern of zeros in the correlation matrix of the distribution q. Let σ * vv and ρ * vw be the variances and the correlations of the data-generating distribution q. The point of departure for our work is Proposition 2.3, which clarifies that the pair (λ q , m q ) is also determined by the behavior of the deterministic Laplace integral where the phase function in the exponent is In the formulation of our results, we adopt the notation as λ q is sometimes referred to as real log-canonical threshold and m q is the threshold's multiplicity. Our formulas for RLCT Ω (H q ) are stated in Theorem 4.3. The proof of the theorem relies on facts presented in Section 3, which concern models with monomial parametrizations in general. As our formulas show, the marginal likelihood admits non-standard large-sample asymptotics, with λ q differing from the model dimension if q exhibits zero correlations (recall Example 1.1). We describe the zero patterns of q in terms of a subforest F * with edge set E * . Our result for trees generalizes directly to models based on forests. If F = (U, E) is a forest with the set V ⊂ U comprising the leaves of the subtrees, then we may define a Gaussian latent forest model M(F ) in the same way as for trees. Again we assign a variance parameter ω v to each node v ∈ V and a correlation parameter ω e to each edge e ∈ E. Forming products of correlations along paths, exactly as in (1.1), we obtain again a parametrization of the correlation matrix of a multivariate normal distribution on R V . In contrast to the case of a tree, there may be pairs of nodes with necessarily zero correlation, namely, when two leaves v and w are in distinct connected components of F . Theorem 4.7 extends Theorem 4.3 to the case of forests. The non-standard cases arise when the data-generating distribution lies in the submodel defined by a proper subforest F * of the given forest F .
The remainder of the paper begins with a review of the connection between the asymptotics of the marginal likelihood and that of the Laplace integral in (1.4); see Section 2 which introduces the notion of a real log-canonical threshold (RLCT). Gaussian latent tree/forest models have a monomial parametrization and we clarify in Section 3 how the monomial structure allows for calculation of RLCTs via techniques from polyhedral geometry. In Section 4, these techniques are applied to derive the above mentioned Theorems 4.3 and 4.7. In Section 5, we demonstrate how our results can be used in model selection with Bayesian information criteria (BIC). In a simulation study, we compare a criterion based on RLCTs to the standard BIC, which is based on model dimension alone.

Background
Consider an arbitrary parametric statistical model M = {P θ : θ ∈ Θ}, with parameter space Θ ⊆ R d . Let each distribution P θ have density p(x|θ) and, for Bayesian inference, consider a prior distribution with density ϕ(θ) on Θ. Writing X n = (X (1) , . . . , X (n) ) for a sample of size n from P θ , the log-likelihood function of M is The key quantity for Bayesian model determination is the integrated or marginal likelihood As in the derivation of the Bayesian information criterion in [Sch78], our interest is in the large-sample behavior of the marginal likelihood. Let the sample X n be drawn from a true distribution with density q that can be realized by the model, that is, q(x) = p(x|θ * ) for some θ * ∈ Θ. Then, as we will make more precise below, the asymptotic properties of the marginal likelihood L(M|X n ) are tied to those of the Laplace integral is the Kullback-Leibler divergence between the data-generating distribution q and distributions in the model M. Note that K q (θ) ≥ 0 for all θ, and K q (θ) = 0 precisely when θ satisfies p(x|θ) = p(x|θ * ). For large n the integrand in (2.2) is equal to ϕ(θ) if K q (θ) = 0 and is negligibly small otherwise. Therefore, the main contribution to the integral Z n (K q ; ϕ) comes from a neighborhood of the zero set which we also call the q-fiber. Suppose now that Θ ⊆ R d is a semianalytic set and that K q : Θ → [0, ∞) is an analytic function with compact q-fiber V Θ (K q ). Suppose further that the prior density ϕ is a smooth and positive function. Then, under additional integrability conditions, the Main Theorem 6.2 in [Wat09] shows that the marginal likelihood has the following asymptotic behavior as the sample size n tends to infinity: In (2.4), λ is a rational number in [0, d], and m is an integer in {1, . . . , d}. The number λ is known as learning coefficient, stochastic complexity or also real logcanonical threshold, and m is the associated multiplicity. As explained in [Wat09,Chap. 4], the pair (λ, m) also satisfies Moreover, the pair (λ, m) can equivalently be defined using the concept of a zeta function as illustrated below; compare also [Lin11].
Definition 2.1 (The real log-canonical threshold). Let f : Θ → [0, ∞) be a nonnegative analytic function whose zero set V Θ (f ) is compact and nonempty. The zeta function can be analytically continued to a meromorphic function on the complex plane. The poles of this continuation are real and positive. Let λ be the smallest pole, known as the real log-canonical threshold (rlct) of f , and let m be its multiplicity. Since we are interested in both the rlct and its multiplicity, we use the notation RLCT Θ (f ; ϕ) := (λ, m). When ϕ(θ) ≡ 1, we simply write RLCT Θ (f ). Finally, if g is another analytic function with RLCT Θ (g; ϕ) = (λ , m ), then we write RLCT Θ (f ; ϕ) > RLCT Θ (g; ϕ) if λ > λ or if λ = λ and m < m .
Example 2.2. Suppose K q (θ) = θ 2 1 θ 2 2 and Θ = [0, 1] 2 . Then the q-fiber V Θ (K q ) is the union of two segments of the coordinate axes. Taking ϕ ≡ 1, we have This example is simple enough that RLCT Θ (K q ) can be computed by elementary means. Let Φ(z) be the distribution function of the standard normal distribution. Then Integration by parts yields Taking logarithms, we see that (2.5) holds with λ = 1 and m = 2. It follows that RLCT Θ (K q ) = (1, 2). Concerning Definition 2.1, we have that for all z ∈ C with Re(z) ≤ 0. In fact, this holds as long as Re(z) < 1. The meromorphic continuation of ζ(z) given by 1/(1 − z) 2 has one pole at λ = 1. The pole has multiplicity m = 2 confirming that RLCT Θ (K q ) = (1, 2).
In this paper we are concerned with Gaussian models for which we may assume, without loss of generality, that all distributions are centered. So let the datagenerating distribution q be the multivariate normal distribution N (0, Σ * ), with positive definite k × k covariance matrix Σ * = (σ * ij ). Further, let p(·|θ) be the density of the distribution N (0, Σ(θ)) with positive definite k × k covariance matrix Σ(θ) = (σ ij (θ)). Then For fixed positive definite Σ * , the function has a full rank Hessian at Φ = Σ * . Hence, in a neighborhood of Σ * , we can both lower-bound and upper-bound K q by positive multiples of the functioñ It follows that RLCT Θ (K q ; ϕ) = RLCT Θ (K q ; ϕ); compare [Wat09, Remark 7.2]. For our study of Gaussian latent tree (and forest) models, it is convenient to change coordinates to correlations and consider the function where ρ * ij and ρ ij (θ) are the correlations obtained from Σ * or Σ(θ); so, e.g., our discussion of latent tree models may thus start from the following fact.
, 1] E be the parameter space for the Gaussian latent tree model M(T ), the parameters being the variances ω v , v ∈ V , and the correlation parameters ω e , e ∈ E. Suppose the (data-generating) distribution q is in M(T ) and has variances σ * vv > 0 and a positive definite correlation matrix with entries ρ * vw . Then

Monomial parametrizations
According to Proposition 2.3, the asymptotic behaviour of the marginal likelihood of a Gaussian latent tree model is determined by the real log-canonical threshold of the function H q in (2.9). This function is a sum of squared differences between monomials formed from the parameter vector ω and constants determined by the data-generating distribution q. In this section, we formulate general results on the real log-canonical thresholds for such monomial parametrizations, which also arise in other contexts [RG05,Zwi11].
Specifically, we treat functions of the form Here, c * 1 , . . . , c * k ∈ R are constants and each monomial ω ui := ω ui1 1 · · · ω u id d is given by a vector of nonnegative integers u i = (u i1 , . . . , u id ). Special cases of this setup are the regular case with H(ω) = ω 2 1 + · · · + ω 2 d , and the quasi-regular case of [YW12], in which the vectors u i have pairwise disjoint supports and all c * i = 0. Let r be the number of summands on the right-hand side of (3.1) that have c * i = 0. Without loss of generality, assume that c * 1 , . . . , c * r = 0 and c * r+1 = · · · = c * k = 0. Furthermore, suppose that ω 1 , . . . , ω s are the parameters appearing in the monomials ω u1 , . . . , ω ur , that is, If H(ω) = 0 then ω i = 0 for all i = 1, . . . , s. Moreover, if the zero set V Ω (H) = {ω ∈ Ω : H(ω) = 0} is compact, then each one of the parameters ω 1 , . . . , ω s is bounded away from zero on V Ω (H). (Clearly, the zero set of the function H q from Proposition 2.3 is compact.) Now define the nonzero part H 1 of H as Definition 3.1. The Newton polytope Γ(H 0 ) of the zero part H 0 is the convex hull of the points (u ij : s + 1 ≤ j ≤ d) ∈ R d−s for i = r + 1, . . . , k. The Newton polyhedron of H 0 is the polyhedron Let 1 = (1, . . . , 1) ∈ R d−s be the vector of all ones. Then the 1-distance of Γ + (H 0 ) is the smallest t ∈ R such that t1 ∈ Γ + (H 0 ). The associated multiplicity is the codimension of the (inclusion-minimal) face of Γ + (H 0 ) containing t1.

We say that
The following is the main result of this section. It is proved in Appendix A.
Theorem 3.2. Suppose that Ω is a product of intervals, and let Ω 1 and Ω 0 be the projections of Ω onto the first s and the last d − s coordinates, respectively. Let H be the sum of squares from (3.1) and assume that the zero set {ω ∈ Ω : H(ω) = 0} is non-empty and compact. Let ϕ : Ω → (0, ∞) be a smooth positive function that is bounded above on Ω. Then Here, λ 0 = 0 and m = 1 if H has no zero part, i.e., s = d.
Remark 3.3. In order to compute the codimension of V Ω1 (H 1 ), one may consider one orthant at a time and take logarithms (accounting for signs). This turns the equations H 1 (ω) = 0 into linear equations in log ω 1 , . . . , log ω s .
. The Newton polytope for H 0 = H is the convex hull of the canonical basis vectors of R d . The Newton polyhedron of H has 1-distance 1/d with multiplicity 1. The same is true whenever H(ω 1 , . . . , ω d ) = ω 2 1 + · · · + ω 2 d + "higher even order terms".
Example 3.6. Consider the function study of the function Since σ * vv > 0 for all v ∈ V , the split of H q into its zero and nonzero part depends solely on the zero pattern among the correlations ρ * vw of the data-generating distribution q. Furthermore, from the form of the parametrization in (1.1), it is clear that zero correlations can arise only if one sets ω e = 0 for one or more edges e in the edge set E. For a fixed set E 0 ⊆ E, the set of parameter vectors ω ∈ Ω with ω e = 0 for all e ∈ E 0 parametrizes the forest model M(F 0 ), where F 0 is the forest obtained from T by removing the edges in E 0 . In this submodel, ρ vw ≡ 0 if and only if v and w lie in two different connected components of F 0 .
It is possible that two different subforests induce the same pattern of zeros among the correlations of the data-generating distribution q. However, there is always a unique minimal forest F * = (U * , E * ) inducing this zero pattern, and we term F * the q-forest. Put differently, the q-forest F * is obtained from T by first removing all edges e ∈ uv for all pairs of nodes u, v ∈ U having zero correlation under q and then removing all inner nodes of T that have become isolated. Isolated leaf nodes are retained so that V ⊆ U * . In the remainder of this section, we take E 0 = E \ E * to be the set of edges whose removal defines F * . We write v ∼ w if v and w are two leaves in V that are joined by a path in the q-forest F * . Moving on to the decomposition of the function from (4.1), recall that we divide the parameter vector ω into coordinates (ω 1 , . . . , ω s ) that never vanish on the q-fiber V Ω (H q ) and the remaining part (ω s+1 , . . . , ω d ). In our case, (ω 1 , . . . , ω s ) consists of all ω v for v ∈ V and ω e for e ∈ E * and (ω s+1 , . . . , ω d ) consists of ω e for e ∈ E 0 = E \ E * . Moreover, The Gaussian latent tree model M(T ) given by a tree T with set of leaves V and edge set E has dimension where l 2 denotes the number of degree two nodes in T . Similarly, the model given by a forest F with set of leaves V and edge set E has dimension where T 1 , . . . , T r are the trees defined by the connected components of F and l 2 is again the number of degrees two nodes.
is a smooth positive function that is bounded above on Ω, then the function H q from (4.1) has where w(e) = |e ∩ U * | ∈ {0, 1, 2} is the number of nodes that e shares with F * , and l 2 is the number of nodes in T that have degree two and are not in U * .
Example 4.4. In Example 4.1, M(F * ) = 5 (c.f. Example 4.2) and e∈E0 w(e) = 3. Hence, the real log-canonical threshold is 13/2, which translates into a coefficient of 13/4 for the log n term in the asymptotic expansion of the log-marginal likelihood. Note that the threshold 13/2 is smaller than dim M(T ) = 9, making the latent tree model behave like a lower-dimensional model.
Example 4.5. Suppose T has two leaves, labelled 1 and 2, and one inner node a, which then necessarily has degree two. If q is a distribution under which the random variables at the two leaves are uncorrelated, then we have Using the calculation from Example 2.2 or Example 3.5, we see that RLCT Ω (H q ) = (3, 2). When applying Theorem 4.3, the q-forest F * has the leaves 1 and 2 isolated and dim M(F * ) = 2. Since l 2 = 1 and each one of the two removed edges satisfies w(e) = 1, the formula from Theorem 4.3 yields RLCT Ω (H q ) = (3, 2), as it should.
Remark 4.6. Note that if T has an (inner) node of degree two, then we can contract one of the edges the node is adjacent to obtain a treeT with M(T ) = M(T ). Repeating such edge contraction it is always possible to find a tree with all inner nodes of degree at least three that defines the same model as the original tree T . Moreover, in applications such as phylogenetics, the trees of interesting do not have nodes of degree two, in which case the multiplicity in RLCT is always equal to one.
In the model selection problems that motivate this work, we wish to choose between different forests. We thus state an explicit result for forests in the below Theorem 4.7. For a forest F , we define q-forests in analogy to the definition we made for trees. In other words, we apply the previous definitions to each tree appearing in the connected components of F and then form the union of the results. Similarly, the proof of Theorem 4.7 is obtained by simply applying Theorem 4.3 to each connected component of the given forest F .
where w(e) = |e ∩ U * | ∈ {0, 1, 2} is the number of nodes that e shares with F * , and l 2 is the number of nodes in F that have degree two and are not in U * .

Singular BIC for latent Gaussian tree models
In this section, we consider the model selection problem of inferring the forest F underlying a Gaussian latent forest model M(F ) based on a sample of independent and identically distributed observations X n = (X (1) , . . . , X (n) ). To this end, we consider Bayesian information criteria that are inspired by the developed largesample theory for the marginal likelihood L(M(F )|X n ). Note that for all the following simulations the space of models we consider implicitly include only forests and trees without degenerate degree 2 nodes; as described in Remark 4.6, this results in an RLCT whose multiplicity is always 1.
As stated in (1.3) and (2.4), the RLCTs found in Section 4 give the coefficients for logarithmic terms that capture the main differences between the log-marginal likelihood and the log-likelihood of the true data-generating distribution q. Let q F be the maximum likelihood estimator of q in the Gaussian latent forest model M(F ). By the results of [Drt09], if q ∈ M(F ) and n → ∞, then and thus, by (2.4), we also have The pair (λ q , m q ) on the right hand side still depends on the unknown datagenerating distribution q through the pair (λ q , m q ). However, the pair is a discontinuous function of q and plugging in the MLEq F has little appeal. Instead, we will consider a criterion proposed by [DP13], in which one averages over the possible values of (λ q , m q ). As in [DP13], we refer to the resulting model selection score as the 'singular Bayesian information criterion', or sBIC for short. We briefly describe how sBIC is computed. Let F be the set of forests in the model selection problem, which we assume to contain the empty forest F ∅ = (V, ∅).
Note that every forest F ∈ F has set of leaves V . For forest F ∈ F with subforest F * ∈ F, let (λ q , m q ) be the pair from (5.1) when the distribution q has F * as q-forest. Theorem 4.7 gives the value of this RLCT pair. Define which is a proxy for the marginal likelihood L(M(F )|X n ) obtained by exponentiating the right hand side of (5.1) and omitting the O p (1) remainder. For each F ∈ F, the sBIC of model M(F ) is defined as log x F , where (x F : F ∈ F) is the unique positive solution to the equation system The system (5.3) is triangular and can be solved by recursively solving univariate quadratic equations. The starting point is the case when F is the empty forest F ∅ , for which F * = F ∅ is the only possible q-forest and (5. The sBIC of the model M(F ∅ ) is thus log L F ∅ F ∅ , which coincides with the usual BIC as the relevant RLCT is given by λ q = dim M(F ∅ ) = |V | and m q = 1. When the forest F is nonempty, the sBIC and the BIC of M(F ) differ.
In [DP13], sBIC is motivated by considering weighted averages of the approximations L F F * , with the weights being data-dependent. Furthermore, it is shown that the sBIC of M(F ) differs from log L(M(F )|X n ) by an O p (1) remainder whenever data are generated from a distribution q ∈ M(F ), even if q lies in a strict submodel M(F * ) ⊂ M(F ). The same is true for BIC only if q ∈ M(F ) does not belong to any strict submodel (i.e., all edge and path correlations are nonzero and F equals the q-forest F * ). In what follows, we explore the differences between the RLCT-based sBIC and the dimension-based BIC in a simulation study.
The first task we consider is selection a subforest of a given tree T , where each subforest as well the tree T share a set of leaves V . When ordering edge sets by inclusion, the set of all subforests of T becomes a poset that we denote by P T . The poset is a lattice with the empty graph (with |V | isolated nodes) as minimal element and the tree T as maximal element. To select a forest, we optimize BIC and sBIC, respectively, over the set P T . Maximum likelihood estimates are computed with an EM algorithm, in which we repeatedly maximize the conditional expectation of the complete-data log-likelihood function of forest models N(F ) for a random vector Z comprising both the observed variables at the leaves in V and the latent variables at the inner nodes of F ; recall the notation from the introduction.
As a concrete example, we choose T to be the tree in Figure 2(a). We generate data from a distribution q that lies in M(T ) but under which the third leaf is independent from all other leaves. The corresponding q-forest F * is depicted on Figure 2(b). We choose q to have covariance matrix In particular, the smallest true model is model 13. Figures 3-5 show that the standard dimension-based BIC tends to select too small models that do not contain the data-generating distribution q. In particular, BIC never selects the tree model 34. The RLCT-based sBIC, on the other hand, invokes a milder penalty, occasionally selects the tree model 34, and more frequently selects the smallest true model 13. Indeed, already for n = 75, sBIC selects the true model more often than any other model. On the other hand, the regular BIC procedure selects too simple a model also when the sample size is increased to n = 125.
Next, we consider examples with 10 and 11 leaves, in which case the number of considered models is still tractable. Writing m := |V | for the number of leaves, the lattice P T has depth m − 1 with the complete independence model having depth 0 and the maximal element M(T ) having depth m − 1. Since the penalty in BIC is always at least the penalty in sBIC, it holds trivially that BIC will select the smallest true model more often than sBIC when the smallest true model is at depth 0; the converse is true if the smallest true model is at depth m − 1. We thus focus on the middle depth and randomly choose 50 trees T 1 , ..., T 50 with corresponding randomly chosen subforests F 1 , ..., F 50 each at depth m−1 2 . From each subforest which we pick q i ∈ M(F i ) by setting all edge correlations to 0.6 and all leaf variances to 1; note that F i equals the q-forest F * i . From each q i , we generate a dataset of a fixed size n and compare the proportion of times that sBIC and BIC correctly identify the smallest true model M(F i ) for 1 ≤ i ≤ 50. The results of these simulations are summarized in Figure 6. We see that sBIC outperforms BIC for smaller sample sizes with BIC marginally overtaking sBIC in very large samples.

Conclusion
Real log-canonical thresholds and associated multiplicities quantify the largesample properties of the marginal likelihood in Bayesian approaches to model selection. In this paper, we computed these RLCTs for Gaussian latent tree and  forest models; the main results being Theorems 4.3 and 4.7. Our computations relied on the fact that the considered tree and forest models have a monomial parametrization, which allows one to apply methods from polyhedral geometry that we presented in Theorem 3.2.
Knowing RLCTs makes it possible to apply a 'singular Bayesian information criterion' (sBIC) that was recently proposed by [DP13]. RLCTs provide refined information about the marginal likelihood and our simulations show that, at least in smaller problems, the sBIC outperforms the usual BIC of [Sch78] that is defined using model dimension alone. We considered up to 11 observed variables.
Our simulations were based on an exhaustive search over the lattice of all considered forests, which quickly becomes infeasible when increasing the number of observed variables. In order to explore potential benefits of sBIC in larger problems of forest selection, techniques such as greedy methods or the structural EM algorithm of [FNPP02] could be considered. However, the plurality of options in such heuristic strategies puts such explorations beyond the scope of this paper.
Appendix A. Proof of Theorem 3.2 Let H be the function from (3.1). By assumption, the 'prior' ϕ : Ω → (0, ∞) is bounded above and V Ω (H) = {ω ∈ Ω : H(ω) = 0} is compact. Since ϕ is smooth and positive, ϕ is bounded away from zero on V Ω (H) and any compact neighborhood of this zero set. The poles of the zeta function in (2.6) can be shown to be the same for all such choices of ϕ, and we have RLCT Ω (H; ϕ) = RLCT Ω (H).
Our proof of Theorem 3.2 now proceeds in three steps: Step 1. Show that RLCT Ω (H) = RLCT Ω (H 0 + H 1 ), where H 0 , H 1 are the zero and nonzero parts of H that are defined in (3.2) and (3.3).
Step 2. To complete Step 2 we will prove the following result.
Before turning to the proof, we exemplify the application of Proposition A.1.
To prove Proposition A.1, note first that when H is equal to its nonzero part and V Ω (H) is compact, RLCT Ω (H) is equal to the RLCT of H over a compact set on which all coordinates of the argument ω are bounded away from zero. Partition this compactum into the intersections with each one of the 2 d orthants in R d . Then RLCT Ω (H) is the minimum RLCT in any orthant. Similarly, the codimension of V Ω (H) is the minimum of any codimension obtained from intersection with an orthant. We may thus consider one orthant at a time. Changing signs as needed to make all coordinates positive, the following lemma becomes applicable.
The result follows from a change of coordinates and an argument about asymptotic equivalence that has been used in other contexts. We include the proof of the lemma for sake of completeness.
Proof. Change coordinates via the substitutionω = log(ω), where the logarithm is applied entry-wise. Since the Jacobian of this transformation is bounded above and bounded away from zero on W , it may be ignored in the computation of the RLCT and thus Since W , and thus also log W , is compact, each of the r linear combinations u T i ω takes its values in a compact set. Restricted to this compact set, the function e xi − e log c * i 2 is asymptotically equivalent to the sum of squares as can been seen by a quadratic Taylor approximation to h 1 around the point (log c * 1 , . . . , log c * r ). Since asymptotically equivalent functions have the same RLCT, the claim is proven.
By an application of Lemma A.3, the proof of Proposition A.1 reduces to an analysis of sums of squares of linear forms, that is, functions of the form with C * i ∈ R and u i ∈ R d . Proposition A.1 thus follows from Proposition A.4 below. Note that V Ω (g) is a polyhedron, which we assume to be nonempty. and Ω is a product of intervals, then RLCT Ω (g) = (codim V Ω (g), 1).
Consider any point x ∈ V Ω (g). By translation, we may assume without loss of generality that x = 0 and g(ω) = i (u T i ω) 2 . We may then take the neighborhood Ω(0) to be equal to {ω ∈ Ω : max i |ω i | ≤ ε} for sufficiently small ε > 0.
When partitioning Ω(0) into orthants, the codimension of V Ω (g) is the minimum of the codimensions of the intersection between V Ω (g) and each one of the orthants. Furthermore, RLCT Ω(0) (g) is equal to the smallest RLCT of g over any of these orthants. Therefore, changing the signs of the coordinates ω i as needed, we are left with checking that RLCT Ω+ (g) is given by the codimension of V Ω+ (g) for Ω + = {ω ∈ Ω : 0 ≤ ω i ≤ ε for all i = 1, . . . , d} and g(ω) = i (u T i ω) 2 . Case 1. If V Ω+ (g) intersects the interior of Ω + , then we may pick any point x + in this intersection and consider Ω + as a neighborhood of x + . After a change of coordinates, we have g(ω) = ω 2 1 + · · · + ω 2 s , where s is the codimension of V Ω+ (g). By Example 3.4, RLCT Ω+ (g) = (s, 1), which was to be shown.
Case 2. Suppose now that V Ω+ (g) is contained in the boundary of Ω + . Since the zero set of g on all of R d is a linear space, V Ω+ (g) is in fact a face of Ω + , and each u T i ω is a supporting hyperplane of Ω + . In particular, after appropriate sign changes, we may assume that u T i ω ≥ 0 on Ω + . The codimension of V Ω+ (g) is equal to the number, say s, of facets of Ω + containing it. Without loss of generality, we may assume that these facets are given by ω 1 = 0, ω 2 = 0, . . . , ω s = 0. This implies that all u i have nonzero entries only in the first s coordinates. We now show show that when restricted to Ω + , the functions g(ω) and f (ω) = ω 2 1 + · · · + ω 2 s are asymptotically equivalent; recall (A.1).
To show that on Ω + , the function g can be bounded from below by a positive multiple of f , note that the fact that u T i ω ≥ 0 on Ω + implies that all u i have nonnegative entries. Hence, where the inequality is obtained by expanding squares and dropping the mixed terms, which are nonnegative. If r i=1 u 2 ij = 0 for some index j then u ij = 0 for all i, which contradicts the fact that ω j = 0 for all ω ∈ V Ω+ (g). Thus, and g(ω) ≥ cf (ω) for all ω ∈ Ω + . To prove that g can be bounded above by a multiple of f , note that all u T i ω are nonnegative on Ω + and thus Let u +j = i u ij and u ++ = j u +j . Then, since all u i have nonnegative entries, Jensen's inequality implies that Since g is asymptotically equivalent to f (ω) = ω 2 1 +· · ·+ω 2 s , we have RLCT Hence, RLCT Ω+ (f ) = RLCT Ω + (f ). From Case 1, we know that RLCT Ω + (f ) = (s, 1). Putting it all together, we have shown that RLCT Ω+ (g) = (s, 1). A.3.
Step 3. The remaining step amounts to proving the following result, which concerns the case where the considered function H is equal to its zero part.
Proposition A.5. Let Ω be a compact product of intervals containing the origin, and let Γ + (H) be the Newton polyhedron of the function H(ω) = i ω 2ui . Then where 1/λ is the 1-distance of Γ + (H) and m is its multiplicity.
Proof. Note that H is invariant under sign changes. Hence, RLCT Ω (H) = RLCT Ω (H) = RLCT Ω∪Ω (H) when Ω is obtained from Ω by changing the signs of any subset of the coordinates ω 1 , . . . , ω d . Forming the unions of Ω and its reflected versions shows that in order to prove Proposition A.5, we may assume that the origin is an interior point of Ω. The claim now follows from Theorem 8.6 in [AGZV88], see also [Lin11,Section 4], and by Remark A.6 below.
Remark A.6. When the origin is in the interior of Ω, the function H(ω) = i ω 2ui has RLCT Ω (H) = RLCT Ω(0) (H) for any small neighborhood Ω(0) of the origin. Indeed, as mentioned in the proof of Proposition A.4, RLCT Ω (H) is the minimum of local RLCTs of H in small neighborhoods Ω(x) of points x ∈ Ω. If x = 0, then some of the variables, say ω 1 , . . . , ω s , are bounded away from zero on a sufficiently small neighborhood Ω(x). Substituting these variables by ω 1 − x 1 , . . . , ω s − x s , respectively, in H, we get a new function H x for which RLCT Appendix B. Proof of Theorem 4.3 Let T = (U, E) be a tree with set of leaves V , and let q be a distribution in the latent tree model M(T ), which has parameter space Ω = (0, ∞) V × [−1, 1] E . We are to compute RLCT Ω (H q ) for the function H q from (4.1), where ω * v and ρ * vw are the variances and correlations of the distribution q. The basic idea of this proof follows [Zwi11].
First, observe that Theorem 3.2 is applicable to this problem. Indeed, H q has the form from (3.1) and the q-fiber V Ω (H q ) is compact. Compactness holds because H q (ω) = 0 implies that ω v = ω * v for all v ∈ V , and all edge correlations ω e , e ∈ E, are in the compact interval [−1, 1]. Now, let F * = (U * , E * ) be the q-forest, and let H 1 q (ω 1 , . . . , ω s ) be the nonzero part of H q given in (4.2). The set V Ω1 (H 1 q ) is equal to the q-fiber under the model M(F * ); recall that Ω 1 is the projection of Ω onto the first s coordinates. We deduce that codim V Ω1 (H 1 q ) = dim M(F * ), which gives the value of λ 1 in Theorem 3.2. It remains to show that the zero part H 0 q (ω s+1 , . . . , ω d ) defined in (4.3) satisfies where Ω 0 is the projection of Ω onto the last d − s coordinates, E 0 = E \ E * is the set of edges that appear in T but not in F * , and l 2 is the number of degree two nodes of T that are not in U * . The zero part of H q is the sum of squares of the monomials recall that v ∼ w if there is no path between v and w in the q-forest F * = (U * , E * ). The edge set E 0 can be partitioned into sets E 01 , . . . , E 0t such that each E 0i defines a tree S i = (U i , E 0i ) that has the set of nodes L i := U i ∩ U * as leaves. In other words, the set of leaves L i of tree S i comprises precisely those nodes that belong to both S i and the q-forest F * . For example, in Figure 1, we have t = 1 and S 1 is the tree with one inner node b and three leaves a, 3, 4. As a further example, consider the tree and q-forest in Figure 7 Let f i be the sum of squares of the monomials in (B.3) that are associated with pairs of distinct nodes u and u in the set of leaves L i of the tree S i . No two trees S i and S j for i = j share an edge. Hence, the two sums of squares f i and f j depend on different subvectors of ω. SinceH 0 q = f 1 + · · · + f t , it follows from (2.5) that see also Remark 7.2(3) in [Wat09]. If T has no nodes of degree two, i.e., l 2 = l 2 = 0, then the same is true for the each tree S i . Lemma B.1 below then implies that Since the nodes in L i lie in F * , we have where w(e) ∈ {0, 1, 2} is the number of nodes of e that lie in the q-forest F * . Combining (B.4)-(B.6), we obtain (B.1) and have thus proven Theorem 4.3 in the case of l 2 = 0 nodes of degree two. The case with nodes of degree two follows the same way applying Lemma B.2 instead of Lemma B.1.
Lemma B.1. Let S = (V, E) be a tree with set of leaves L and all inner nodes of degree at least three. Let f be the sum of squares of the monomials If Ω is a neighborhood of the origin, then RLCT Ω (f ) = |L| 2 , 1 .
Proof. If |L| = 2, then S has a single edge and no inner nodes. In this case, f is the square of a single variable and it is clear RLCT Ω (f ) = (1, 1) = (|L|/2, 1). In the remainder of this proof, we assume that |L| ≥ 3.
By Proposition A.5, it suffices to compute the 1-distance and its multiplicity for the Newton polyhedron Γ + (f ) ⊂ R E . By Definition 3.1, the polyhedron Γ + (f ) is determined by the exponent vectors of the monomials in (B.7). Each exponent vector is the incidence vector for a path between a pair of leaves. In other words, each pair of two distinct leaves v and w defines a vector u ∈ R E with u e = 1 if e ∈ vw and u e = 0 otherwise. Write U for the set of all these |L| 2 vectors. Let E L be the set of terminal edges of S, i.e., the |L| edges that are incident to a leaf. We claim that every point x in the Newton polyhedron Γ + (f ) satisfies and that the inequality defines a facet of Γ + (f ). Indeed, if x ∈ U then e∈E L x e = 2 because every path between two leaves in L includes precisely two edges in E L . It is then clear that (B.8) holds for all points x ∈ Γ + (f ). Moreover, by [MP08, Lemma 1], the span of U is all of R E . Hence, the affine hull of U is the hyperplane given by e∈E L x e = 2, and we conclude that (B.8) defines a facet of Γ + (f ).
Since |E L | = |L|, inequality (B.8) implies that the 1-distance of Γ + (f ) is at least 2/|L|. We claim that it is equal to 2/|L|. In fact, we will show that the vector 2 |L| 1 not only lies in the Newton polyhedron but also in the Newton polytope Γ(f ), that is, the vector is a convex combination of the incidence vectors in U. To prove this, we construct a set of paths P in the tree S such that (i) each element of P is a path between leaves of S, (ii) P contains precisely |L| paths, and (iii) every edge of S is covered by exactly two paths of P. The construction implies our claim because the average of the incidence vectors of the paths in P is equal to 2 |L| 1. Let S * be any trivalent tree that has the same set of leaves L as S and that can be obtained from S * by edge contraction. Here, a tree is trivalent if each inner node has degree three. We will use induction on the number of leaves to show that a set of paths P with the desired properties (i)-(iii) exists. Figure 8 shows an example.
If S * has three |L| = 3 leaves, then there is a single inner node and each path between two leaves has two edges. We may simply take P to be the set of all the three paths that exist between pairs of leaves. This provides the induction base.
In the induction step, pick two leaves v and w of the tree S * that are joined by a path with two edges {v, a} and {a, w}. The node a is an inner node of S * . Remove the two edges and the two leaves to form a subtree S * * , in which a becomes a leaf. Then S * * has |L| − 1 leaves and, by the induction hypothesis, there is a set of paths P * * that satisfies properties (i)-(iii) with respect to S * * . In particular, |P * * | = |L| − 1. Now, precisely two paths in P * * have the node a as an endpoint. Extend one of them by adding the edge {a, v} and extend the other by adding {a, w}. This gives two paths between leaves of S * . All other paths in P * * are already paths between leaves of S * . Add one further path, namely, (v, a, w), and denote the resulting collection of |L| paths by P * . Clearly, the set P * satisfies properties (i)-(iii) with respect to S * . Contracting each path in P * by applying the edge contractions that transform S * into S, we obtain a system of paths P that satisfies properties (i)-(iii) with respect to S.
Finally, note that in the construction we just gave we can ensure that P includes a given path between two leaves in L. Hence, the vector 2 |L| 1 can be written as a convex combination of vertices of Γ(f ) such that a given vertex x get positive weight. It follows that 2 |L| 1 lies in the interior of the Newton polytope and thus the multiplicity m is 1. Figure 8. An example of a system of paths such that each edge of a trivalent tree is covered by exactly two paths.
The next result generalizes the previous lemma to the case of trees with nodes of degree 2. We remark Example 2.2 is a special case of this generalization. It matches the case where the tree S has two leaves and one inner node, which is then necessarily of degree two. If Ω is a neighborhood of the origin, then RLCT Ω (f ) = |L| 2 , 1 + l 2 where l 2 is the number of (inner) nodes of S that have degree two.
Proof. Suppose a is an inner node of degree two, and that a is incident to the two edges e = {a, b} and f = {a, c}. Then any path connecting to leaves in L either uses both e and f or neither e nor f . Hence, if x is the incidence vector of a path between two leaves in L, then x e = x f . It follows that the affine hull of Newton polytope generated by the path incidence vectors is no longer a hyperplane but an affine space of dimension |E| − 1 − l 2 .
Proceeding exactly as in the proof of Lemma B.1, we see that it still holds that the 1-distance of the Newton polyhedron Γ + (f ) is 2/|L|. Similarly, the ray spanned by 1 still meets Γ + (f ) in the relative interior of the Newton polytope Γ(f ). However, since the codimension of the Newton polytope is now 1 + l 2 , we have RLCT Ω (f ) = (|L|/2, 1 + l 2 ).