Title: BIAS-VARIANCE DECOMPOSITION FOR BOLTZMANN MACHINES

Abstract: We achieve bias-variance decomposition for Boltzmann machines using an information geometric formulation. Our decomposition leads to an interesting phenomenon that the variance does not necessarily increase when more parameters are included in Boltzmann machines, while the bias always decreases. Our result gives a theoretical evidence of the generalization ability of deep learning architectures because it provides the possibility of increasing the representation power with avoiding the variance inflation.

1. INTRODUCTION:

Understanding why the deep learning architectures can generalize well despite their high representation power with a large number of parameters is one of crucial problems in theoretical deep learning analysis, and there are a number of attempts to solve the problem with focusing on several aspects such as sharpness and robustness (Dinh et al., 2017; Wu et al., 2017; Keskar et al., 2017; Neyshabur et al., 2017; Kawaguchi et al., 2017) . However, the complete understanding of this phenomenon is not achieved yet due to the complex structure of deep learning architectures.

To theoretically analyze the generalizability of the architectures, in this paper, we focus on Boltzmann machines (Ackley et al., 1985) and its generalization including higher-order Boltzmann machines (Sejnowski, 1986; Min et al., 2014) , the fundamental probabilistic model of deep learning (see the book by Goodfellow et al. (2016, Chapter 20) for an excellent overview), and we firstly present bias-variance decomposition for Boltzmann machines. The key to achieve this analysis is to employ an information geometric formulation of a hierarchical probabilistic model, which was firstly explored by Amari (2001) ; Nakahara & Amari (2002) ; Nakahara et al. (2006) . In particular, the recent advances of the formulation by Sugiyama et al. (2016; 2017) enables us to analytically obtain the Fisher information of parameters in Boltzmann machines, which is essential to give the lower bound of variances in bias-variance decomposition.

We show an interesting phenomenon revealed by our bias-variance decomposition: The variance does not necessarily increase while the bias always monotonically decreases when we include more parameters in Boltzmann machines, which is caused by its hierarchical structure. Our result indicates the possibility of designing a deep learning architecture that can reduce both of bias and variance, leading to better generalization ability with keeping the representation power.

The remainder of this paper is organized as follows: First we formulate the log-linear model of hierarchical probability distributions using an information geometric formulation in Section 2, which includes the traditional Boltzmann machines (Section 2.2) and arbitrary-order Boltzmann machines (Section 2.3). Then we present the main result of this paper, bias-variance decomposition for Boltzmann machines, in Section 3 and discuss its property. We empirically evaluate the tightness of our theoretical lower bound of the variance in Section 4. Finally, we conclude with summarizing the contribution of this paper in Section 5.

2. FORMULATION:

To theoretically analyze learning of Boltzmann machines (Ackley et al., 1985) , we introduce an information geometric formulation of the log-linear model of hierarchical probability distributions, which can be viewed as a generalization of Boltzmann machines.

2.1. PRELIMINARY: LOG-LINEAR MODEL:

First we prepare a log-linear probabilistic model on a partial order structure, which has been introduced by Sugiyama et al. (2016; 2017) . Let (S, ≤) be a partially ordered set, or a poset (Gierz et al., 2003) , where a partial order ≤ is a relation between elements in S satisfying the following three properties for all x, y, z ∈ S: (1) x ≤ x (reflexivity), (2) x ≤ y, y ≤ x ⇒ x = y (antisymmetry), and (3) x ≤ y, y ≤ z ⇒ x ≤ z (transitivity). We assume that S is always finite and includes the least element (bottom) ⊥ ∈ S; that is, ⊥ ≤ x for all x ∈ S. We denote S \ {⊥} by S + .

We use two functions, the zeta function ζ : S × S → {0, 1} and the Möbius function µ : S × S → Z. The zeta function ζ is defined as ζ(s, x) = 1 if s ≤ x and ζ(s, x) = 0 otherwise, while the Möbius function µ is its convolutional inverse, that is,

∑ s∈S ζ(x, s)µ(s, y) = ∑ x≤s≤y µ(s, y) = { 1 if x = y, 0 otherwise, which is inductively defined as µ(x, y) =    1 if x = y, - ∑ x≤s<y µ(x, s) if x < y, 0 otherwise.

For any functions f , g, and h with the domain S such that

g(x) = ∑ s∈S ζ(s, x)f (s) = ∑ s≤x f (s), h(x) = ∑ s∈S ζ(x, s)f (s) = ∑ s≥x f (s),

f is uniquely recovered using the Möbius function:

f (x) = ∑ s∈S µ(s, x)g(s), f (x) = ∑ s∈S µ(x, s)h(s).

This is the Möbius inversion formula and is fundamental in enumerative combinatorics (Ito, 1993) . Sugiyama et al. (2017) introduced a log-linear model on S, which gives a discrete probability distribution with the structured outcome space (S, ≤). Let P denote a probability distribution that assigns a probability p(x) for each x ∈ S satisfying ∑ x∈S p(x) = 1. Each probability p(x) for x ∈ S is defined as

log p(x) := ∑ s∈S ζ(s, x)θ(s) = ∑ s≤x θ(s).

(1)

From the Möbius inversion formula, θ is obtained as

θ(x) = ∑ s∈S µ(s, x) log p(s).

(2)

In addition, we introduce η : S → R as

EQUATION

)

The second equation is from the Möbius inversion formula. Sugiyama et al. (2017) showed that the set of distributions S = {P | 0 < p(x) < 1 and ∑ p(x) = 1} always becomes the dually flat Riemannian manifold. This is why two functions θ and η are dual coordinate systems of S connected with the Legendre transformation, that is,

θ = ∇φ(η), η = ∇ψ(θ) with two convex functions ψ(θ) := -θ(⊥) = -log p(⊥), φ(η) := ∑ x∈S p(x) log p(x).

Moreover, the Riemannian metric g(ξ) (ξ = θ or η) such that g(θ) = ∇∇ψ(θ), g(η) = ∇∇φ(η), which corresponds to the gradient of θ or η, is given as

EQUATION

)

for all x, y ∈ S + . Furthermore, S is in the exponential family (Sugiyama et al., 2016) , where θ coincides with the natural parameter and η with the expectation parameter.

Let us consider two types of submanifolds:

S α = { P ∈ S | θ(x) = α(x) for all x ∈ dom(α) } , S β = { P ∈ S | η(x) = β(x) for all x ∈ dom(β) }

specified by two functions α, β with dom(α), dom(β) ⊆ S + , where the former submanifold S α has constraints on θ while the latter S β has those on η. It is known in information geometry that S α is e-flat and S β is m-flat, respectively (Amari, 2016, Chapter 2.4) . Suppose that dom(α) ∪ dom(β) = S + and dom(α) ∩ dom(β) = ∅. Then the intersection S α ∩ S β is always the singleton, that is, the distribution Q satisfying Q ∈ S α and Q ∈ S β always uniquely exists, and the following Pythagorean theorem holds:

EQUATION

for any P ∈ S α and R ∈ S β .

2.2. STANDARD BOLTZMANN MACHINES:

A Boltzmann machine is represented as an undirected graph G = (V, E) with a vertex set V = {1, 2, . . . , n} and an edge set E ⊆ {{i, j} | i, j ∈ V }. The energy function Φ:{0, 1} n → R of the Boltzmann machine G is defined as

Φ(x; b, w) = - ∑ i∈V b i x i - ∑ {i,j}∈E w ij x i x j for each x = (x 1 , x 2 , . . . , x n ) ∈ {0, 1} n , where b = (b 1 , b 2 , . . . , b n

) and w = (w 12 , w 13 , . . . , w n-1n ) are parameter vectors for vertices (bias) and edges (weight), respectively, such that

w ij = 0 if {i, j} ̸ ∈ E. The probability p(x; b, w) of the Boltzmann machine G is obtained for each x ∈ {0, 1} n as p(x; b, w) = exp(-Φ(x; b, w)) Z (8)

with the partition function Z such that

Z = ∑ x∈{0,1} n exp(-Φ(x; b, w)) (9)

to ensure the condition ∑ x∈{0,1} n p(x; b, w) = 1. It is clear that a Boltzmann machine is a special case of the log-linear model in Equation ( 1) with S = 2 V , the power set of V , and ⊥ = ∅. Let each x ∈ S be the set of indices of "1" of x ∈ {0, 1} n and ≤ be the inclusion relation, that is, x ≤ y if and only if x ⊆ y. Suppose that

EQUATION

)

where |x| is the cardinality of x, which we call a parameter set. The Boltzmann distribution in Equations ( 8) and ( 9) directly corresponds to the log-linear model in Equation (1):

EQUATION

)

where

θ(x) = b x if |x| = 1 and θ(x) = w x if |w| = 2.

This means that the set of Boltzmann distributions S(B) that can be represented by a parameter set B is a submanifold of S given as

EQUATION

)

Given an empirical distribution P . Learning of a Boltzmann machine is to find the best approximation of P from the Boltzmann distributions S(B), which is formulated as a minimization problem of the KL (Kullback-Leibler) divergence:

min P B ∈S(B)

D KL ( P, P B ) = min

EQUATION

)

This is equivalent to maximize the log-likelihood

L(P B ) = N ∑ s∈S p(s) log p B (s) with the sample size N . Since we have ∂ ∂θ B (x) D KL ( P, P B ) = ∂ ∂θ B (x) ∑ s∈S p(s) log p B (s) = ∂ ∂θ B (x) ∑ s∈S   p(s) ∑ ⊥<u≤s θ B (u)   - ∂ ∂θ B (x) ψ(θ B ) ∑ s∈S p(s) = η(x) -η B (x),

the KL divergence D KL ( P, P B ) is minimized when η(x) = η B (x) for all x ∈ B, which is well known as the learning equation of Boltzmann machines as η(x) and η B (x) coincides with the expectation for the outcome x with respect to the empirical distribution P obtained from data and the model distribution P B represented by a Boltzmann Machine, respectively. Thus the minimizer

EQUATION

)

This distribution P B is known as m-projection of P onto S(B) (Sugiyama et al., 2017) , which is unique and always exists as S has the dually flat structure with respect to (θ, η).

2.3. ARBITRARY-ORDER BOLTZMANN MACHINES:

The parameter set B is fixed in Equation ( 10) in the traditional Boltzmann machines, but our loglinear formulation allows us to include or remove any element in S + = 2 V \ {⊥} as a parameter. This attempt was partially studied by Sejnowski (1986) ; Min et al. (2014) that include higher order interactions of variables to increase the representation power of Boltzmann machines. For S = 2 V with V = {1, 2, . . . , n} and a parameter set B ⊆ S + , which is an arbitrary subset of S + = S \ {∅}, the Boltzmann distribution given by an arbitrary-order Boltzmann machine is defined by

log p(x) = ∑ s∈B ζ(s, x)θ(s) -ψ(θ), ψ(θ) = -θ(⊥),

and the submanifold of Boltzmann distributions is given by Equation ( 12). Hence Equation ( 14) gives the MLE (maximum likelihood estimation) of the empirical distribution P .

Let B 1 , B 2 , . . . , B m be a sequence of parameter sets such that

B 1 ⊆ B 2 ⊆ • • • ⊆ B m-1 ⊆ B m = S + .

Since we have a hierarchy of submanifolds

S(B 1 ) ⊆ S(B 2 ) ⊆ • • • ⊆ S(B m-1 ) ⊆ S(B m ) = S,

we obtain the decreasing sequence of KL divergences:

EQUATION

where each P Bi = argmin P ∈S(Bi) D KL ( P, P ), the best approximation of P using B i .

There are two extreme cases as a choice of the parameter set B. On the one hand, if B = ∅, the Boltzmann distribution is always the uniform distribution, that is, p(x) = 1/2 |V | for all x ∈ S. Thus there is no variance but nothing will be learned from data. On the other hand, if B = S + , the Boltzmann distribution can always exactly represent the empirical distribution P , that is, D KL ( P, P B ) = D KL ( P, P ) = 0. Thus there is no bias in each training but the variance across different samples will be large. To analyze the tradeoff between the bias and the variance, we perform bias-variance decomposition in the next section.

Another strategy to increase the representation power is to use hidden variables (Le Roux & Bengio, 2008) such as restricted Boltzmann machines (RBMs) (Smolensky, 1986; Hinton, 2002) and deep Boltzmann machines (DBMs) (Salakhutdinov & Hinton, 2009; 2012) . A Boltzmann machine with hidden variables is represented as G = (V ∪H, E), where V and H correspond to visible and hidden variables, respectively, and the resulting domain S = 2 V ∪H (see Appendix for the formulation of Boltzmann machines with hidden variables as the log-linear model). It is known that the resulting model can be singular (Yamazaki & Watanabe, 2005; Watanabe, 2007) and its statistical analysis cannot be directly performed. Studying bias-variance decomposition for such Boltzmann machines with hidden variables is the interesting future work.

3. BIAS-VARIANCE DECOMPOSITION:

Here we present the main result of this paper, bias-variance decomposition for Boltzmann machines.

We focus on the expectation of the squared KL divergence E[D KL (P * , PB ) 2 ] from the true (unknown) distribution P * to the MLE PB of an empirical distribution P by a Boltzmann machine with a parameter set B, and decompose it using information geometric properties.

In the following, we use the MLE P * B of the true distribution P * , which is the closest distribution in the set of distributions that can be modeled by Boltzmann machines in terms of the KL divergence and is mathematically obtained with replacing P with P * in Equation ( 13).

Theorem 1 (Bias-variance decomposition of the KL divergence). Given a Boltzmann machine with a parameter set B. Let P * ∈ S be the true (unknown) distribution, P * B , PB ∈ S(B) be the MLEs of P * and an empirical distribution P , respectively. We have 

, B) = ∑ s∈B ∑ u∈B η * (s)η * (u)cov ( θB (s), θB (u) ) , var(P * B , B) = 1 N ∑ s∈B ∑ u∈B η * (s)η * (u)(I -1 ) su

with the equality holding when the sample size N → ∞, where cov( θB (s), θB (u)) denotes the error covariance between θB (s) and θB (u) and I -1 is the inverse of the Fisher information submatrix I ∈ R |B|×|B| of P * B with respect to the parameter set B given in Equation (4).

Proof. From the Pythagorean theorem illustrated in Figure 1 , Since θB (s) is an unbiased estimator of θ * B (s) for every s ∈ S, it holds that

EQUATION

E [ D KL (P * B , PB ) ] = E [ ∑ x∈S p * B (x) log p * B (x) pB (x) ] = ∑ x∈S p * B (x) ∑ s≤x E [ θ * B (s) -θB (s) ] = 0.

Hence we have

EQUATION

)

The second term is

E [ D KL (P * B , PB ) 2 ] = E   ( ∑ x∈S p * B (x) log p * B (x) pB (x) ) 2   = E   ∑ x∈S ∑ y∈S p * B (x)p * B (y) log p * B (x) pB (x) log p * B (y) pB (y)   = E   ∑ x∈S ∑ y∈S p * B (x)p * B (y)   ∑ s∈B,s≤x ( θ * B (s) -θB (s) ) - ( ψ(θ * B ) -ψ( θB ) )     ∑ u∈B,u≤y ( θ * B (u) -θB (u) ) - ( ψ(θ * B ) -ψ( θB ) )     = E   ∑ x∈S ∑ y∈S ∑ s∈B ∑ u∈B p * B (x)p * B (y)ζ(s, x)ζ(u, y) ( θ * B (s) -θB (s) ) ( θ * B (u) -θB (u) )   -2 E   ∑ x∈S p * B (x) ∑ s∈B,s≤x ( θ * B (s) -θB (s) ) ( ψ(θ * B ) -ψ( θB ) )   + E [ ( ψ(θ * B ) -ψ( θB ) ) 2 ] = ∑ s∈B ∑ u∈B η * B (s)η * B (u)cov ( θB (s), θB (u) ) + var ( ψ( θB ) ) -2 ∑ s∈B η * B (s)cov ( θB (s), ψ( θB ) ) = ∑ s∈B ∑ u∈B η * (s)η * (u)cov ( θB (s), θB (u) ) + var ( ψ( θB ) ) -2 ∑ s∈B η * (s)cov ( θB (s), ψ( θB ) ) ,

where cov( θB (s), θB (u) ) denotes the error covariance between θB (s) and θB (u) and var( ψ( θB ) ) denotes the variance of ψ( θB ), and the last equality comes from the condition in Equation ( 14).

Here the term of the (co)variance of the normalizing constant (partition function) ψ(θ):

var ( ψ( θB ) ) -2 ∑ s∈B η * (s)cov ( θB (s), ψ( θB ) )

is the irreducible error since ψ(θ) = -θ(⊥) is orthogonal for every parameter θ(s), s ∈ S and the Fisher information vanishes from Equation (4), i.e.,

E [ ∂ log p(s) ∂θ(s) ∂ log p(s) ∂θ(⊥) ] = ∑ s∈S ζ(x, s)p(s) -η(x) = η(x) -η(x) = 0.

For the variance term

var(P * B , B) = ∑ s∈B ∑ u∈B η * (s)η * (u)cov ( θB (s), θB (u) ) ,

we have from the Cramér-Rao bound (Amari, 2016, Theorem 7.1) 

since θ B (s), s ∈ B is unbiased cov ( θB (s), θB (u) ) ≥ 1 N (I -1 ) su

with the equality holding when N → ∞, where I ∈ R |B|×|B| is the Fisher information matrix with respect to the parameter set B such that

I su = ∑ x∈S ζ(s, x)ζ(u, x)p * B (x) -η * B (s)η * B (u)

for all s, u ∈ B, which is given in Equation ( 5), and I -1 is its inverse. Finally, from Equation ( 16) we obtain

E [ D KL (P * , PB ) 2 ] = D KL (P * , P * B ) 2 + E [ D KL (P * B , PB ) 2 ] = D KL (P * , P * B ) 2 + ∑ s∈B ∑ u∈B η * (s)η * (u)cov ( θB (s), θB (u) ) + irreducible error ≥ D KL (P * , P * B ) 2 + 1 N ∑ s∈B ∑ u∈B η * (s)η * (u)(I -1 ) su + irreducible error

with the equality holding when N → ∞.

Let B, B ′ ⊆ S ′ = 2 V ∪H such that B ⊆ B ′ , that is, B ′ has more parameters than B.

Then it is clear that the bias always reduces, that is, D KL (P * , P * B ) ≥ D KL (P * , P * B ′ ) because S(B) ⊆ S(B ′ ). However, this monotonicity does not always hold for the variance. We illustrate this non-monotonicity in the following example. Let S = 2 V with V = {1, 2, 3} and assume that the true distribution P * is given by ( p * ({∅}), p * ({1}), p * ({2}), p * ({3}), p * ({1, 2}), p * ({1, 3}), p * ({2, 3}), p * ({1, 2, 3}) ) = (0.2144, 0.0411, 0.2037, 0.145, 0.1423, 0.0337, 0.0535, 0.1663), which was randomly generated from the uniform distribution. Suppose that we have three types of parameter sets B 1 = {{1}}, B 2 = {{1}, {2}}, and

B 3 = {{1}, {2}, {1, 2}} such that B 1 ⊂ B 2 ⊂ B 3 . Then biases become D KL (P * , P * B1 ) 2 = 0.0195, D KL (P * , P * B2 ) 2 = 0.0172, D KL (P * , P * B3 ) 2 = 0

.0030, where we can confirm that the bias monotonically decreases, while the lower bound of the variance with N = 100 are var(P * B1 , B 1 ) = 0.0062, var(P * B2 , B 2 ) = 0.0192, var(P * B3 , B 3 ) = 0.0178. Moreover, we have computed the exact variance var(P * B , B) by repeating 1, 000 times generating a sample with N = 100 and obtained the following: var(P * B1 , B 1 ) = 0.0066, var(P * B2 , B 2 ) = 0.0267, var(P * B3 , B 3 ) = 0.0190, hence var(P * B2 , B 2 ) > var(P * B3 , B 3 ) happens, that is, the variance actually decreases when we include more parameters. This interesting property, non-monotonicity of the variance with respect to the growth of parameter sets, comes from the hierarchical structure of S. 

4. EMPIRICAL EVALUATION OF TIGHTNESS:

We empirically evaluate the tightness of our lower bound. In each experiment, we randomly generated a true distribution P * , followed by repeating 1, 000 times generating a sample (training data) with the size N from P * to empirically estimate the variance var(P * B , B). We consistently used the Boltzmann machine represented as the fully connected graph G = (V, E) such that V = {1, 2, . . . , n} and E = {{i, j} | i, j ∈ V }. Thus the parameter set B is given as

B = { x ∈ S + = 2 V \ {∅} |x| ≤ 2 } .

We report in Figure 2 the mean ± SD (standard deviation) of the empirically obtained variance and its theoretical lower bound var(P * B , B) obtained by repeating the above estimation 100 times. In Figure 2 (a) the sample size N is varied from 10 to 10, 000 with fixing the number of variables n = 5 while in Figure 2 (b) n is varied from 3 to 7 with fixing N = 100. These results overall show that our theoretical lower bound is tight enough if N is large and is reasonable across each n.

5. CONCLUSION:

In this paper, we have firstly achieved bias-variance decomposition of the KL divergence for Boltzmann machines using the information geometric formulation of hierarchical probability distributions. Our model is a generalization of the traditional Boltzmann machines, which can incorporate arbitrary order interactions of variables. Our bias-variance decomposition reveals the nonmonotonicity of the variance with respect to growth of parameter sets, which has been also reported elsewhere for non-linear models (Faber, 1999) . This result indicates that it is possible to reduce both bias and variance when we include more higher-order parameters in the hierarchical deep learning architectures. To solve the open problem of the generalizability of the deep learning architectures, our finding can be fundamental for further theoretical development.

Input layer:

Hidden 

A BOLTZMANN MACHINES WITH HIDDEN NODES:

A Boltzmann machine with hidden variables is represented as G = (V ∪ H, E), where V and H correspond to visible and hidden variables, respectively, and the resulting domain S = 2 V ∪H . In particular, restricted Boltzmann machines (RBMs) (Smolensky, 1986; Hinton, 2002) are often used in applications, where the edge set is given as

E = { {i, j} | i ∈ V, j ∈ H } Moreover,

in a deep Boltzmann machine (DBM) (Salakhutdinov & Hinton, 2009; 2012) , which is the beginning of the recent trend of deep learning (LeCun et al., 2015; Goodfellow et al., 2016) , the hidden variables H are divided into k disjoint subsets (layers) H 1 , H 2 , . . . , H k and E = { {i, j} | i ∈ H l-1 , j ∈ H l , l ∈ {1, . . . , k} } , where V = H 0 for simplicity.

Let S = 2 V and S ′ = 2 V ∪H and S and S ′ be the set of distributions with the domains S and S ′ , respectively. In both cases of RBMs and DBMs, we have B = { x ∈ S ′ | |x| = 1 or x ∈ E } , (see Figure 3 ) and the set of Boltzmann distributions is obtained as S ′ (B) = { P ∈ S ′ θ(x) = 0 for all x ̸ ∈ B } . Since the objective of learning Boltzmann machines with hidden variables is MLE (maximum likelihood estimation) with respect to the marginal probabilities of the visible part, the target empirical distribution P ∈ S is extended to the submanifold S ′ ( P ) such that S ′ ( P ) = { P ∈ S ′ η(x) = η(x) for all x ∈ S } , and the process of learning Boltzmann machines with hidden variables is formulated as double minimization of the KL divergence such that min P ∈S ′ ( P ),P B ∈S ′ (B) D KL (P, P B ).

Since two submanifolds S ′ (B) and S ′ ( P ) are e-flat and m-flat, respectively, it is known that the EM-algorithm can obtain a local optimum of Equation ( 17) (Amari, 2016, Section 8.1.3) , which was first analyzed by Amari et al. (1992) . Since this computation is infeasible due to combinatorial explosion of the domain S ′ = 2 V ∪H , a number of approximation methods such as Gibbs sampling have been proposed (Salakhutdinov & Hinton, 2012) .

Let us fix an empirical distribution P and consider two RBMs with parameter sets B, B ′ ⊆ S ′ = 2 V ∪H . If B ⊆ B ′ , that is, B ′ has more hidden nodes than B, we always have the monotonicity: min 

:



Figures and Tables:

KL (P*, P B ) = D KL (P*, P* B ) + D KL (P* B , P B ) of the Pythagorean theorem.

Figure 2: Empirical evaluation of variance. Empirically estimated variances (blue, dotted lines) and theoretically obtained lower bounds (red, solid line) with fixing n = 5 in (a) and N = 100 in (b).

∈S ′ ( P ),P B ′ ∈S ′ (B ′ ) D KL (P, P B ′ ) ≤ minP ∈S ′ ( P ),P B ∈S ′ (B) D KL (P, P B ) (18) as B ⊆ B ′ implies S ′ (B) ⊆ S ′ (B ′). This result corresponds to Theorem 1 in(Le Roux & Bengio, 2008), the representation power of RBMs.

Figure 3: An example of a deep Boltzmann machine (left) with an input (visible) layer V = {1, 2} with two hidden layers H 1 = {3} and H 2 = {4}, and the corresponding domain set S V ∪H (right).In the right-hand side, the colored objects {1}, {2}, {3}, {4}, {1, 3}, {2, 3}, and {3, 4} denote the parameter set B, which correspond to nodes and edges of the DBM in the left-hand side.

Table content: <table><tr><td>layer 1</td><td>layer 2 Hidden</td><td/><td/><td colspan="2">{1,2,3,4} {1,2,3,4}</td></tr><tr><td>1</td><td/><td/><td colspan="3">{1,2,3} {1,2,4} {1,3,4} {2,3,4} {1,2,3} {1,2,4} {1,3,4} {2,3,4}</td></tr><tr><td>3</td><td>4</td><td>{1,2} {1,2}</td><td colspan="3">{1,3} {1,4} {2,3} {2,4} {1,3} {1,4} {2,3} {2,4}</td><td>{3,4} {3,4}</td></tr><tr><td>2</td><td/><td/><td>{1} {1}</td><td>{2} {2}</td><td>{3} {3}</td><td>{4} {4}</td></tr><tr><td/><td/><td/><td/><td>Ø Ø</td></tr></table>

Peer review (iclr2018):
 This paper uses an information geometric view on hierarchical models to discuss a bias - variance decomposition in Boltzmann machines, presenting interesting conclusions, whereby some more care appears to be needed for making these claims. The paper arrives at the main conclusion that it is possible to reduce both the bias and the variance in a hierarchical model. The discussion is not specific to deep learning nor to Boltzmann machines, but actually addresses hierarchical exponential family models. The methods pertaining hierarchical models are interesting and presented in a clear way. My concern are the following points: The main theorem presents only a lower bound, meaning that it provides no guarantee that the variance can indeed be reduced. The paper seems to ignore that a model with hidden variables may be singular, in which case the Fisher metric is not positive definite and the Cramer Rao bound has no meaning. This interferes with the claims and derivations made in the paper in the case of models with hidden variables. The problem seems to lie in the fact that the presented derivations assume that an optimal distribution in the data manifold is given (see Theorem 1 and proof), effectively making this a discussion about a fully observed hierarchical model. In particular, it is not further specified how to obtain θˆB(s) in page 6 before (13). Also, in page 5 the paper states that ``it is known that the EM-algorithm can obtain the global optimum of Equation (12) (Amari, 2016, Section 8.1.3)**. However, what is shown in that reference is only that: (Theorem 8.2., Amari, 2016) ``The KL-divergence decreases monotonically by repeating the E-step and the M-step. Hence, the algorithm converges to an equilibrium.** A model with hidden variables can have several global and local optimisers (see, e.g. https://arxiv.org/abs/1709.05276). The critical points of the EM algorithm can have a non trivial structure, as has been observed in the case of non negative rank matrix varieties (see, e.g., https://arxiv.org/pdf/1312.5634.pdf). OTHER In page 3, ``S_.beta is e-flat and S_.alpha ... **, should this not be the other way around? (See also page 5 last paragraph of Section 2.) Please also indicate the precise location in the provided reference. All pages up to page 5 are introduction. Section 2.3. as presented is very vague and does not add much to the discussion. In page 7, please explain E ψ(θˆ )^2 −ψ(θ∗ )^2=0