O NLINE A LGORITHMS FOR S UM -P RODUCT N ETWORKS WITH C ONTINUOUS VARIABLES

Online Algorithms for Sum-Product Networks
with Continuous Variables
Priyank Jaini1
Abdullah Rashwan1
Han Zhao2
Yue Liu1
Ershad Banijamali1
Zhitang Chen3
Pascal Poupart1

PJAINI @ UWATERLOO . CA
ARASHWAN @ UWATERLOO . CA
HAN . ZHAO @ CS . CMU . EDU
Y 565 LIU @ UWATERLOO . CA
SBANIJAMALI @ UWATERLOO . CA
CHENZHITANG 2@ HUAWEI . COM

PPOUPART @ UWATERLOO . CADavid R. Cheriton School of Computer Science, University of Waterloo, Waterloo, Ontario, Canada
2
Machine Learning Dept., School of Computer Science, Carnegie Mellon University, Pittsburgh, PA, USA
3
Noah’s Ark Laboratory, Huawei Technologies, Hong Kong, China
1

Abstract
Sum-product networks (SPNs) have recently emerged as an attractive representation due to theirdual interpretation as a special type of deep neural network with clear semantics and a tractableprobabilistic graphical model. We explore online algorithms for parameter learning in SPNs with
continuous variables. More specifically, we consider SPNs with Gaussian leaf distributions andshow how to derive an online Bayesian moment matching algorithm to learn from streaming data.
We compare the resulting generative models to stacked restricted Boltzmann machines and generative moment matching networks on real-world datasets.
Keywords: Sum-product networks; continuous variables; online learning.

1. Introduction
Sum-product networks (SPNs) (Poon and Domingos, 2011) recently emerged as an attractive representation that can be seen as both a deep architecture and a tractable probabilistic graphical model.
More precisely, SPNs are deep neural networks restricted to sum and product operators. While
this limits their expressiveness, it allows us to attribute clear semantics to the nodes by interpretingthe sub-SPN rooted at each node as a distribution over a subset of variables (under suitable conditions). In that sense, SPNs are also a type of probabilistic graphical model (PGM). In fact, SPNs
are equivalent to arithmetic circuits (ACs) (Darwiche, 2003; Rooshenas and Lowd, 2014) and SPNs
can be converted into Bayesian networks (BNs) (Zhao et al., 2015) (and vice-versa). SPNs (andACs) differ from traditional PGMs such as BNs and Markov networks (MNs) from a computational
perspective. SPNs (and ACs) are both a representation and an inference machine that allows exact
inference to be performed in linear time with respect to the size of the network. This is particularly
advantageous when learning a PGM from data since inference can always be performed exactly andtractably when the resulting PGM is an SPN or AC, but not when it is a BN or MN (unless the search
space is restricted to a subclass of tractable PGMs such as bounded tree-width PGMs). Based on
those advantages, SPNs have been used in many applications including image completion (Poonand Domingos, 2011), speech modeling (Peharz et al., 2014) and language modeling (Cheng et al.,2014).
1

JAINI , R ASHWAN , Z HAO , L IU , BANIJAMALI , C HEN AND P OUPART

To date, most of the algorithms for SPNs assume that the variables of interest are categorical.
However, many real-world applications are best modeled by continuous variables. Hence, there isa need for SPN algorithms that work with continuous data. Several papers mention that it is possible to generalize SPNs to continuous variables by considering leaf distributions over continuous
variables (Poon and Domingos, 2011; Gens and Domingos, 2012, 2013; Zhao et al., 2015), but only
a recent paper by Peharz et al. (2016) proposes an EM algorithm for continuous SPNs and reports
experiments with continuous data.
With the rise of large and streaming datasets, there is also a need for online algorithms. In thiswork, we investigate online algorithms for SPNs with Gaussian leaf distributions. More precisely,
we describe how to learn the parameters of such continuous SPNs by extending the online Bayesian
moment matching and online EM algorithms proposed by Rashwan et al. (2016) from categorical toGaussian SPNs. We evaluate the algorithms with real-world datasets and compare their effectivenessto stacked restricted Boltzmann machines (Salakhutdinov and Hinton, 2009) and generative moment
matching networks (Li et al., 2015), which are other types of generative deep neural networks.

2. BackgroundA sum-product network (SPN) (Poon and Domingos, 2011) is a rooted acyclic directed graph whose
internal nodes are sums and products, and leaves are tractable distributions over some random variables. The edges linking each sum node to its children are labeled with non-negative weights. SPNs
are equivalent to arithmetic circuits (Darwiche, 2003), which also consist of rooted acyclic directed
graphs of sums and products, with the only difference that edges do not have weights, while leaves
store numerical values. This is only a syntactic difference since edge weights in SPNs can always betransformed into leaves with corresponding numerical values in ACs (and vice-versa) (Rooshenas
and Lowd, 2014).
The value of an SPN is the value computed by its root in a bottom up pass. More specifically,
Vnode (x) denotes the value computed by a node based on evidence x.

Vnode (x) =node (x) pPnode is a leafw
V(x)node is a sumcc Q c∈children(node)
node is a productc∈children(node) Vc (x)(1)

Here, pnode (x) denotes the probability density or mass of the evidence x in the distribution at aleaf node, depending on whether x is continuous or discrete. When the evidence is empty or does
not instantiate any of the variables in the distribution at a leaf node then pnode (x) = 1, which isequivalent to marginalizing out all the variables for which we do not have any evidence.
Under suitable conditions, an SPN (or AC) encodes a joint distribution over a set of variables.
Let scope(node) be the set of variables in the leaves of the sub-SPN rooted at a node.
Definition 1 (Completeness/smoothness (Darwiche, 2002; Poon and Domingos, 2011)) A sum
node is complete/smooth when the scopes of its children are the same (i.e., ∀c, c0 ∈ children(sum),
scope(c) = scope(c0 )).
Definition 2 (Decomposability (Darwiche, 2002; Poon and Domingos, 2011)) A product node isdecomposable when the scopes of its children are disjoint (i.e., ∀c, c0 ∈ children(sum) and c 6=
c0 , scope(c) ∧ scope(c0 ) = ∅).
2

O NLINE A LGORITHMS FOR S UM -P RODUCT N ETWORKS WITH C ONTINUOUS VARIABLES

When all sum nodes are complete/smooth and all product nodes are decomposable, an SPN iscalled valid, which ensures that the value computed by the SPN for some evidence is proportional tothe probability of that evidence (i.e., Pr(x) ∝ Vroot (x)). This means that we can answer inference
queries Pr(x1 |x2 ) by computing two SPN evaluations Vroot (x1 ) and Vroot (x2 ):Pr(x1 |x2 ) =Vroot (x1 , x2 )Vroot (x2 )(2)Since each SPN evaluation consists of a bottom up pass that is linear in the size of the network, thecomplexity of inference in SPNs is linear in the size of the SPN (same holds for ACs). In contrast,
inference may be exponential in the size of the network for BNs and MNs.
We can interpret SPNs as hierarchical mixture models since each sum node can be thought astaking a mixture of the distributions encoded by its children, where the probability of each childcomponent is proportional to the weight labeling the edge of that child. Completeness/smoothness
ensures that each child is a distribution over the same set of variables. Similarly, we can think of eachproduct node as factoring a distribution into a product of marginal distributions. Decomposability
ensures that the marginal distributions are independent.
Several algorithms have been proposed to learn the structure and parameters of SPNs. Structurelearning is important since it is difficult to specify a structure by hand that satisfies the completeness/smoothness and decomposability properties. Structure learning algorithms for SPNs are generally based on clustering techniques that greedily construct a structure by partitioning/grouping thedata and the variables (Dennis and Ventura, 2012; Gens and Domingos, 2013; Peharz et al., 2013;
Rooshenas and Lowd, 2014; Adel et al., 2015; Vergari et al., 2015).
Parameter learning techniques for SPNs include maximum likelihood by gradient ascent (Poonand Domingos, 2011), expectation maximization (Poon and Domingos, 2011; Peharz, 2015; Peharz et al., 2016) and signomial programming (Zhao and Poupart, 2016), as well as approximate
Bayesian learning by moment matching (Rashwan et al., 2016) and collapsed variational inference (Zhao et al., 2016). Among those algorithms, gradient ascent, expectation maximization, moment matching and collapsed variational inference can learn in an online fashion with streaming
data (Rashwan et al., 2016; Zhao et al., 2016).
2.1 Moment Matching
Since we will extend the Bayesian moment matching technique by Rashwan et al. (2016) to continuous SPNs with Gaussian leaves, we review some background about moment matching. In SPNs
the parameters that need to be learned are the weights associated with each sum node. In Gaussian
SPNs, the mean and covariance matrix of each leaf node also need to be learned. Let θ denote theset of parameters that need to be learned based on some data x. Bayesian learning proceeds by computing the posterior Pr(θ|x). When the posterior is intractable, we can approximate it by a simpler
distribution after processing each data point x. This simpler distribution is chosen to match a subset of the moments of the exact intractable posterior, hence the name Bayesian moment matching,which is a form of assumed density filtering (Minka and Lafferty, 2002).
A moment is a quantitative measure of the shape of a distribution or a set of points. Letf(θ|φ) be a probability density function over a d-dimensional random variable θ =h {θ1 , θ2i, ..., θd }Q niparametrized by φ. The j th order moments of θ are defined as Mgj (θ) (f ) = Ef
where
i θi
3

JAINI , R ASHWAN , Z HAO , L IU , BANIJAMALI , C HEN AND P OUPARTP

i ni = j and gj is a monomial of θ of degree j.

Z
Mgj (θ) (f ) =gj (θ)f (θ|φ)dθθ

For some distributions f, there exists a set of monomials S(f) such that knowing Mg (f ) ∀g ∈ S(f )allows us to calculate the parameters of f. For example, for a Gaussian distribution N (x; µ, σ 2 ), theset of sufficient moments S(f) = {x, x2 }. This means knowing Mx and Mx2 allows us to estimate
the parameters µ and σ 2 that characterize the distribution. The method of moments is used in theBayesian Moment Matching (BMM) algorithm.
The method of moments is a popular frequentist technique used to estimate the parameters ofa probability distribution based on the evaluation of the empirical moments of a dataset. It has
been previously used to estimate the parameters of latent Dirichlet allocation, mixture models andhidden Markov models (Anandkumar et al., 2012). The method of moments or moment matching
can also be used in a Bayesian setting to compute a subset of the moments of an intractable posterior
distribution. Subsequently, another distribution with the same moments is selected from a tractablefamily of distributions. Consider a Gaussian mixture model, which is a simple SPN where the root
is a sum linked to a layer of product nodes that are each linked to Gaussian leaf distributions. We
use the Dirichlet as a prior over the weights of the mixture and a Normal-Wishart distribution as aprior over each Gaussian component. We next give details about the Dirichlet and Normal-Wishart
distributions, including their set of sufficient moments.
2.2 Family of Prior Distributions
In Bayesian Moment Matching, we project the posterior onto a tractable family of distributions bymatching a set of sufficient moments. To ensure scalability, it is desirable to start with a familyof distributions that are conjugate priors for multinomial distributions (for the set of weights) and
Gaussian distributions with unknown mean and covariance matrix. The product of a Dirichlet distribution over the weights with a Normal-Wishart distribution over the mean and covariance matrixof each Gaussian component ensures that the posterior is a mixture of products of Dirichlet andNormal-Wishart distributions. Subsequently, we can approximate this mixture in the posterior witha single product of Dirichlet and Normal-Wishart distributions by using moment matching. We explain this in greater detail in Section 3, but first we describe briefly the Dirichlet and Normal-Wishart
distribution along with some sets of sufficient moments.
2.2.1 D IRICHLET D ISTRIBUTION
The Dirichlet distribution is a family of multivariate continuous probability distributions over the
probability simplex. It is the conjugate prior probability distribution for the multinomial distribution
and hence it is a natural choice of prior over the set of weights w = {w1 , w2 , ..., wM } of a Gaussian
mixture model. A set of sufficient moments for the Dirichlet distribution is S = {(wi , wi2 ) : ∀i ∈{1, 2, ..., M }}. Let α = {α1 , α2 , ...., αM } be the parameters of the Dirichlet distribution over w,thenαiE[wi ] = P∀i ∈ {1, 2, ..., M }j αjαi (αi + 1)


P
α
1
+
α
j
j
j
j

E[wi2 ] =  P

4∀i ∈ {1, 2, ..., M }(3)

O NLINE A LGORITHMS FOR S UM -P RODUCT N ETWORKS WITH C ONTINUOUS VARIABLES

2.2.2 N ORMAL W ISHART P RIOR
The Normal-Wishart distribution is a multivariate distribution with four parameters. It is the conjugate prior of a multivariate Gaussian distribution with unknown mean and covariance matrix. This
makes a Normal-Wishart distribution a natural choice for the prior over the unknown mean andprecision matrix for our case.
Let µ be a d-dimensional vector and Λ be a symmetric positive definite d × d matrix of random variables respectively. Then, a Normal-Wishart distribution over (µ, Λ) given parameters(µ0 , κ, W, ν) is such that µ ∼ Nd µ; µ0 , (κΛ)−1 where κ > 0 is real, µ0 ∈ Rd and Λ has
a Wishart distribution given as Λ ∼ W(Λ; W, ν) where W ∈ Rd×d is a positive definite matrix and ν > d − 1 is real.  The marginal distribution of µ is a multivariate t-distribution i.e
Wµ ∼ tν−d+1 µ; µ0 , κ(ν−d+1). The univariate equivalent for the Normal-Wishart distribution is
the Normal-Gamma distribution.
In Section 2.1, we defined S as a set of sufficient moments to characterize a distribution. In
the case of the Normal-Wishart distribution, we would require at least four different moments toestimate the four parameters that characterize it. A set of sufficient moments in this case is S ={µ, µµT , Λ, Λ2ij } where Λ2ij is the (i, j)th element of the matrix Λ. The expressions for sufficient
moments are given byE[µ] = µ0
κ+1
W−1κ(ν − d − 1)E[Λ] = νW

E[(µ − µ0 )(µ − µ0 )T ] =V ar(Λij ) = ν(Wij2 + Wii Wjj )(4)

3. Bayesian Moment Matching
We now discuss in detail the Bayesian Moment Matching (BMM) algorithm. Since Bayesian learning yields a posterior that consists of a mixture distribution with an exponentially growing numberof mixture components, BMM approximates the posterior after each observation with fewer components in order to keep the computation tractable.
In Algorithm 1, we first describe a generic procedure to approximate the posterior Pn after eachobservation with a simpler distribution P˜n by moment matching. More precisely, a set of moments
sufficient to define P˜n are matched with the moments of the exact posterior Pn . At each iteration,
we first calculate the exact posterior Pn (Θ|x1:n ). Then, we compute the set of moments S(f) thatare sufficient to define a distribution in the family f (Θ|Φ). Next, we compute the parameter vector
Φ based on the set of sufficient moments. This determines a specific distribution P˜n in the family fthat we use to approximate Pn . Note that the moments in the sufficient set S(f) of the approximate
posterior are the same as that of the exact posterior. However, all the other moments outside this setof sufficient moments S(f) may not necessarily be the same.
3.1 BMM for multivariate Gaussian mixture model
In the previous section, we gave a generic version of the Bayesian moment matching algorithm. In
this section, we briefly discuss an instantiation of the algorithm for a multivariate Gaussian mixture
model, which is a simple 3-layer SPN consisting of a sum (top layer/root) of products (middle layer)
5

JAINI , R ASHWAN , Z HAO , L IU , BANIJAMALI , C HEN AND P OUPART

Algorithm 1 Generic Bayesian Moment Matching
Input: Data xn , n ∈ {1, 2, ..., N }Let f(Θ|Φ) be a family of probability distributions with parameters Φ
Initialize a prior P̃0 (Θ)
for n = 1 to N do
Compute Pn (Θ|x1:n ) from P̃n−1 (Θ|x1:n−1 ) using Bayes’ theorem∀ g(Θ) ∈ S(f ), evaluate Mg(Θ) (Pn )Compute Φ using Mg(Θ) (Pn )’s
Approximate Pn with P̃n (Θ|x1:n ) = f(Θ|Φ)
end for
Return Θ̂ = EP̃n [Θ]of multivariate Gaussians (bottom
Q layer/leaves). The family of distributions for the prior in this casebecomes P0 (Θ) = Dir(w|a) M
i=1 N W d (µi , Λi |αi , κi , Wi , νi ) where w = (w1 , w2 , ..., wM ) anda = (a1 , a2 , ..., aM ). The set of sufficient moments for the posterior in this case would be given by
S(P (Θ|x)) = {µi , µi µTi , Λi , Λ2ikl , wi , wi2 : ∀i ∈ 1, 2, .., M } where Λikl is the (k, l)th element of
the matrix Λi . Notice that, since Λi is a symmetric matrix, we only need to consider the momentsof the elements on and above the diagonal of Λi .
In Eq. 4 of Section 2.2.2, we presented the expressions for a set of sufficient moments of aNormal-Wishart distribution. Using those expressions we can again approximate a mixture of
products of Dirichlet and Normal-Wishart distributions in the posterior with a single product ofDirichlet and Normal-Wishart distributions, as we did in the previous section. Finally, the estimate
Θ = EP̃n [Θ] is obtained after observing the data x1:n .
3.2 Bayesian Moment Matching for Continuous SPNs
In this work, we consider SPNs that are trees (i.e., each node has a single parent) since many structure learning algorithms produce SPNs that are trees. When an SPN is a tree, it is possible tocompute all the moments simultaneously in time that is linear in the size of the network. Two coefficients coefi0 , coefi1 are computed at each node i in a bottom-up and top-down pass of the network.
Algorithm 2 shows how those coefficients are computed. Once we have the coefficients, we cancompute each moment as follows. For leaf nodes that model multivariate Gaussian distributions,the set of sufficient moments is Si = {µi , µi µTi , Λi , Λ2ikl } ∀i ∈ leafN odes and we get:Z
MP (Θ|x) (s) =ssN W i (µi , Λi |αi , κi , Wi , νi ) coefi0 + coefi1 N (x|µi , Λ−1)ds ∀s ∈ Sii

For the rest of the nodes:
Z


X
k
k
MP (Θ|x) (wij ) =wij 0 Vj 0 (x) dwi· ∀i ∈ sumN odes
wij
Dir(wi· |αi· ) coefi0 + coefi1
wi·j0

At a high level, Alg. 2 maintains a posterior distribution over all the sum nodes in the subtree ofa given node. This is achieved by a bottom-up pass of the network to compute all the vi (x) at eachnode followed by a top-down recursive pass to compute the required moments at all the sum nodes.
6

O NLINE A LGORITHMS FOR S UM -P RODUCT N ETWORKS WITH C ONTINUOUS VARIABLES

Q
Algorithm 2 computeCoefficients(node) based on x and weight prior i Dir(wi· |αi· )1: if isRoot(node) then02:
coefnode←01
3:
coefnode ← 1
4: else if isP roduct(parent(node)) then
005:
coefnode← coefparent
Q11
6:
coefnode← coefparent
sibling vsibling (x)7: else if isSum(parent(node)) thenP
αP
parent,sibling
001
8:
coefnode← coefparent+ coefparent
siblingαparent,j vsibling (x)
j

α

1
1
Pparent,nodecoefnode← coefparent
j αparent,j
10: end if
11: if isN otLeaf (node) then computeCoefficients(child) ∀child

9:

In Alg. 2, vi (x) is defined as:vi (x) =Z Y

Dir(wj |αj )Vi (x) dw

(5)

w j

where the integration is over all the model parameters w that appear in the sub-tree rooted at nodei, and upper case Vi (x) is the value computed at node i when the input instance is x. Note that byutilizing the same network structure, all the vi (x) can be computed in a single bottom-up pass.
To see the correctness of Alg. 2, let r be the root node, p be a child product node of r, and s
is a child sum node of p, we will show how to compute the posterior as we go down the tree while
integrating the rest of the network above a certain node. We can write down the posterior as follows:YP (w|x) ∝Dir(wi |αi ) coefr0 + coefr1 Vr (x)(6)
i∈subtree(r)The base case is given by setting coefr0 = 0 and coefr1 = 1 at the root. We then expand Vr (x) as:X
Vr (x) = wr,p Vp (x) +wr,s Vs (x)(7)
s∈siblings(p)Substituting Vr (x) into P (w|x), we get:Y
P (w|x) ∝Dir(wi |αi ) coefr0 + coefr1 wr,p Vp (x) +i∈subtree(r)Xwr,i Vi (x)

i∈siblings(p)(8)
As we go down one level in the tree, we integrate wr and the weights in all sibling subtrees of p:
P ({wi }i∈subtree(p) |x)
!!
P
α
v
(x)
r,i
i
α
i∈siblings(p)r,pP
Vp (x) +∝Dir(wi |αi ) coefr0 + coefr1 Pi αr,ii αr,ii∈subtree(p)
!
!
P


Y
α
v
(x)
r,i
i
α
i∈siblings(p)r,pP=Dir(wi |αi )
coefr0 + coefr1+ coefr1 P
Vp (x)αr,iii αr,iY

i∈subtree(p)7JAINI , R ASHWAN , Z HAO , L IU , BANIJAMALI , C HEN AND P OUPART

where vi is the value when integrating the weights. Let coefp0 = coefr0 +coefr1α, then we can write down P (w|x) as follows:
and coefp1 = coefr1 P r,pi αr,iY

P ({wi }i∈subtree(p) |x) ∝P

αr,i vi (x)i∈siblings(p)P
i αr,iDir(wi |αi ) coefp0 + coefp1 Vp (x)(9)i∈subtree(p)Another step down to sum node s by expanding Vp (x):P ({wi }i∈subtree(p) |x) ∝

Y




Y

Dir(wi |αi ) coefp0 + coefp1 Vs (x)Vi (x)

i∈siblings(s)i∈subtree(p)(10)
Again, we can integrate all weights of the sibling subtrees of node s, we get:P ({wi }i∈subtree(s) |x) ∝

Y




Y

Dir(wi |αi ) coefp0 + coefp1i∈subtree(s)vi (x) Vs (x)i∈siblings(s)(11)
Let coefs0 = coefp0 and coefs1 = coefp1
P ({wi }i∈subtree(s) |x) ∝Qi∈siblings(s) vi (x), we get:YDir(wi |αi ) coefs0 + coefs1 Vs (x)(12)

i∈subtree(s)The previous steps can be done recursively along each path from root to leaf since the underlying
network structure is assumed to be a tree.
3.3 Distributed Bayesian Moment Matching
One of the major advantages of Bayes’ theorem is that the computation of the posterior can be
distributed over several machines, each of which processes a subset of the data. It is also possible
to compute the posterior in a distributed manner using Bayesian moment matching algorithm. For
example, let us assume that we have T machines and a data set with TN data points. Each machine t,
can compute the approximate posterior Pt (Θ|x(t−1)N +1:tN ) where t ∈ 1, 2, .., T using Algorithm 1
over N data points. These partial posteriors {Pt }Tt=1 can be combined to obtain a posterior over theentire data set x1:T N according to the following equation:

P (Θ|x

1:T N) ∝ P (Θ)

TY
Pt (Θ|x(t−1)N +1:tN )
t=1

P (Θ)

(13)

Subsequently, the estimate Θ̂ = EP (Θ|x1:T N ) [Θ] is obtained over the whole data set. Therefore,
we can use the Bayesian moment matching algorithm to perform Bayesian learning in an online
and distributed fashion. Note that there is no loss incurred by distributing the computation of theposterior. Eq. 13 is exact and simply follows from Bayes’ theorem and the assumption that data
instances are sampled independently from some underlying distribution.
8

O NLINE A LGORITHMS FOR S UM -P RODUCT N ETWORKS WITH C ONTINUOUS VARIABLES

4. Related WorkThe expectation maximization (EM) algorithm for estimating the parameters of SPNs, including
both the edge weights as well as the parameters of the leaf distributions, has recently been proposedby Peharz et al. (2016). In the E-step, the authors apply the differential approach developed in Darwiche (2003) to efficiently propagate the sufficient statistics for both edge parameters as well as leaf
parameters in time linear in the size of the network. In the M-step, marginalization can be donein closed form for edge parameters, where the edge parameters incident to each sum node form acategorical distribution. When the leaf distributions are chosen from the exponential family, the update formula in the M-step for the leaf parameters can be computed by integration over the naturalparameters of the corresponding distribution. When leaf distributions are Gaussians, as the case inthis paper, a closed form update formula can be obtained as well. We refer the readers to Peharz
et al. (2016) for more detailed discussions of batch EM for SPNs with continuous leaf nodes.
Most deep neural networks are function approximators without any generative capability andcannot be interpreted as probabilistic graphical models. However, Li et al. (2015) recently introduced Generative Moment Matching Networks (GenMMNs), which are deep neural networks thattake as input random values and output samples mimicking a desired dataset. Hence those networks
can be viewed as generative models that can approximate sampling from the underlying distribution of any dataset. GenMMNs are trained by the so-called empirical maximum mean discrepancy(MMD) (Gretton et al., 2006) between the outputs of the network and a training set. Using the kernel
trick, MMD measures the maximum discrepancy between all empirical moments of a dataset andthe generative model encoded by a neural network. This is a form of frequentist moment matchingthat leads to an objective that can be optimized by backpropagation (Li et al., 2015).

5. Experiments
We performed experiments on 7 real-world datasets from the UCI machine learning repository (Lichman, 2013) and the function approximation repository (Guvenir and Uysal, 2000) that span diverse
domains with 3 to 48 attributes. We compare the performance of online algorithms for Sum-ProductNetworks with continuous variables as well as other generative deep networks like Stacked Restricted Boltzman Machines (SRBMs) and Generative Moment Matching Networks (GenMMNs).
We evaluated the performance of SPNs using both online EM (oEM) and online Bayesian Moment
Matching (oBMM). Stochastic gradient descent and exponentiated gradient descent could also beused, but since their convergence rate was shown to be much slower than that of oEM (Rashwan
et al., 2016; Zhao and Poupart, 2016), we do not report results for them. Since there are no structurelearning algorithms for continuous SPNs, we used a random structure and a structure equivalent to a
Gaussian Mixture Model for the experiments using both oEM and oBMM. In both cases, the leaves
are multivariate Gaussian distributions. Before discussing the results, we provide additional details
regarding the training and likelihood computation for oEM, SRBM and GenMMN.
Online EM. Peharz et al. (2016) describe how to estimate the leaf parameters in a batch mode by
EM. To adapt their algorithm to the online setting, we follow the framework introduced by Neal and
Hinton (1998) where increments of the sufficient statistics are accumulated in a streaming fashion.
More specifically, we sample one data instance at a time to compute the corresponding local sufficient statistics from that instance. Then we accumulate the local sufficient statistics as increments
to the globally maintained sufficient statistics. The update formula is obtained by renormalization
9

JAINI , R ASHWAN , Z HAO , L IU , BANIJAMALI , C HEN AND P OUPART

according to the global sufficient statistics. This online variant of EM is also known as incremental
EM Liang and Klein (2009). There are no hyperparameters to be tuned in the online EM (oEM)
algorithm for SPNs. In the experiments, we randomly initialize the weights of SPNs and apply oEM
to estimate the parameters. For each data set, oEM processes each data instance exactly once.
SRBM. We trained stacked restricted Boltzmann machine (SRBM) with 3 hidden layers ofbinary nodes where the first and second hidden layers have the same number of nodes as the dimensionality of the data and the third hidden layer has four times as many nodes as the dimensionality ofthe data. This yielded a number of parameters that is close to the number of parameters of the random SPNs for a fair comparison. This 1 : 1 : 4 structure is taken from (Hinton and Salakhutdinov,2006).Computing the exact likelihood of SRBMs is intractable. Hence, we used Gibbs sampling toapproximate the log-likelihood. We averaged 1000 runs of 500 iterations of Gibbs sampling wherein each iteration the value of each hidden node is re-sampled given the values of the adjacent nodes.
At the end of each run, the log-likelihood of each data instance is computed given the last values
sampled for the hidden nodes in the layer above the leaf nodes.
Gen MMN. We trained fully-connected networks with multiple layers. The number of layers
and hidden units per layer depends on the training dataset and was selected to match the number ofparameters in random SPNs for a fair comparison. Inputs to the networks are vectors of elements
drawn independently and uniformly at random in [−1, 1]. The activation function for the hidden
layers is ReLU and for the output layer is sigmoid. We used the Gaussian kernel in MMD andfollowed the algorithm described in (Li et al., 2015) for training. To compute the average loglikelihood reported in Table 1, we generate 10, 000 samples using the trained model and fit a kernel
density estimator to these samples. The average log-likelihood of the test set is then evaluated using
this estimated density. The kernel for our density estimator is Gaussian and its parameter is set
to maximize the log-likelihood of the validation set. This technique was used since there is no
explicit probability density function for this model. However, as shown in (Theis et al., 2015) thismethod only provides an approximate estimate of the log-likelihood and therefore the log-likelihood
reported for GenMMNs in Table 1 may not be directly comparable to the likelihood of other models.
5.1 Results
Table 1 reports the average log-likelihoods (with standard error) of each model on the 7 real-world
datasets using 10-fold cross-validation. No results are reported for random SPNs trained by oBMM
and oEM when the number of attributes is too small (less than 5) to consider a structure that is more
complex than GMMs. No results are reported for GMM SPNs trained by oBMM and oEM when
the number of attributes is large (48) in comparison to the amount of data since the number of parameters to be trained in the covariance matrix of each Gaussian component is square in the numberof attributes. We can think of GMMs as a baseline structure for SPNs when the dimensionalityof the data is low and random SPNs as more data efficient models for high-dimensional datasets.
Note also that since the data is continuous, the log-likelihood can be positive when a model yields
a good fit as observed for the Flow Size dataset. The best results for each dataset are highlighted inbold. oBMM consistently outperforms oEM on each dataset. In comparison to other deep models,
the results are mixed and simply suggest that continuous SPNs trained by oBMMs are competitive
with SRBMs and GenMMNs. Note that it is remarkable that SPNs with a random or basic GMM
structure can compete with other deep models in the absence of any structure learning. Note also10

O NLINE A LGORITHMS FOR S UM -P RODUCT N ETWORKS WITH C ONTINUOUS VARIABLES

Table 1: Log-likelihood scores on real-world data sets. The best results are highlighted in bold font.
Dataset# of vars

Flow Size Quake Banknote Abalone Kinematics
CA
Sensorless Drive
3
4
4
8
8
22
48
-1.82
-11.19-2.47
1.58
oBMM (random)± 0.19± 0.03± 0.56± 1.28-11.36-11.35-31.34-3.40
oEM (random)± 0.19± 0.03± 1.07± 6.063
4.80
-3.84
-4.81
-1.21
-11.24-1.78
oBMM (GMM)± 0.67 ± 0.16 ± 0.13 ± 0.36
± 0.04± 0.59
-0.49
-5.50
-4.81
-3.53
-11.35-21.39oEM (GMM)± 3.29 ± 0.41 ± 0.13 ± 1.68
± 0.03± 1.58
-0.79
-2.38
-2.76
-2.28
-5.55
-4.95
-26.91SRBM± 0.004 ± 0.01 ± 0.001 ± 0.001± 0.02 ± 0.003± 0.03
0.40
-3.83
-1.70
-3.29
-11.36-5.41-29.41
GenMMN± 0.007 ± 0.21 ± 0.03 ± 0.10
± 0.02± 0.14± 1.185that SPNs constitute a tractable probabilistic model in which inference can be done exactly in linear
time. This is not the case for SBRMs and GenMMNs where inference must be approximated.

6. Conclusion
We presented the first online algorithm for parameter learning in continuous SPNs with Gaussian
leaves. We also reported on the first comparison between continuous SPNs and other generative
deep models such as SRBMs and GenMMNs. In the future, we plan to devise a structure learningalgorithm for continuous SPNs. Structure learning has the potential to significantly improve themodeling ability of continuous SPNs.

Acknowledgments
This research was funded by a grant from Huawei Noah’s Ark Lab in Hong Kong.

References
T. Adel, D. Balduzzi, and A. Ghodsi. Learning the structure of sum-product networks via an svdbased algorithm. In UAI, 2015.
A. Anandkumar, D. Hsu, and S. M. Kakade. A method of moments for mixture models and hidden
Markov models. In COLT, volume 1, page 4, 2012.
W.-C. Cheng, S. Kok, H. V. Pham, H. L. Chieu, and K. M. A. Chai. Language modeling with
sum-product networks. In Annual Conf. of the Int. Speech Communication Association, 2014.
A. Darwiche. A logical approach to factoring belief networks. KR, 2:409–420, 2002.
A. Darwiche. A differential approach to inference in Bayesian networks. JACM, 50(3):280–305,
2003.
11

JAINI , R ASHWAN , Z HAO , L IU , BANIJAMALI , C HEN AND P OUPART

A. Dennis and D. Ventura. Learning the architecture of sum-product networks using clustering on
variables. In NIPS, 2012.
R. Gens and P. Domingos. Discriminative learning of sum-product networks. In NIPS, pages 3248–
3256, 2012.
R. Gens and P. Domingos. Learning the structure of sum-product networks. In ICML, pages 873–
880, 2013.
A. Gretton, K. M. Borgwardt, M. Rasch, B. Schölkopf, and A. J. Smola. A kernel method for thetwo-sample-problem. In NIPS, pages 513–520, 2006.
H. A. Guvenir and I. Uysal. Bilkent university function approximation repository. 2000.
G. E. Hinton and R. R. Salakhutdinov. Reducing the dimensionality of data with neural networks.
Science, 313(5786):504–507, 2006.
Y. Li, K. Swersky, and R. Zemel. Generative moment matching networks. In ICML, pages 1718–
1727, 2015.
P. Liang and D. Klein. Online EM for unsupervised models. In HLT-NAACL, pages 611–619, 2009.
M. Lichman. UCI machine learning repository, 2013.
T. Minka and J. Lafferty. Expectation-propagation for the generative aspect model. In UAI, pages
352–359, 2002.
R. M. Neal and G. E. Hinton. A view of the em algorithm that justifies incremental, sparse, andother variants. In Learning in graphical models, pages 355–368. Springer, 1998.
R. Peharz. Foundations of Sum-Product Networks for Probabilistic Modeling. PhD thesis, Medical
University of Graz, 2015.
R. Peharz, B. C. Geiger, and F. Pernkopf. Greedy part-wise learning of sum-product networks. In
Machine Learning and Knowledge Discovery in Databases, pages 612–627. Springer, 2013.
R. Peharz, G. Kapeller, P. Mowlaee, and F. Pernkopf. Modeling speech with sum-product networks:
Application to bandwidth extension. In ICASSP, pages 3699–3703, 2014.
R. Peharz, R. Gens, F. Pernkopf, and P. Domingos. On the latent variable interpretation in sumproduct networks. arXiv:1601.06180, 2016.
H. Poon and P. Domingos. Sum-product networks: A new deep architecture. In UAI, pages 2551–
2558, 2011.
A. Rashwan, H. Zhao, and P. Poupart. Online and Distributed Bayesian Moment Matching for
Sum-Product Networks. In AISTATS, 2016.
A. Rooshenas and D. Lowd. Learning sum-product networks with direct and indirect variable interactions. In ICML, pages 710–718, 2014.
R. Salakhutdinov and G. E. Hinton. Deep boltzmann machines. In AISTATS, pages 448–455, 2009.
L. Theis, A. Oord, and M. Bethge.
A note on the evaluation of generative models.
arXiv:1511.01844, 2015.
A. Vergari, N. Di Mauro, and F. Esposito. Simplifying, regularizing and strengthening sum-productnetwork structure learning. In ECML-PKDD, pages 343–358. 2015.
H. Zhao and P. Poupart. A unified approach for learning the parameters of sum-product networks.
arXiv:1601.00318, 2016.
H. Zhao, M. Melibari, and P. Poupart. On the relationship between sum-product networks andBayesian networks. In ICML, 2015.
H. Zhao, T. Adel, G. Gordon, and B. Amos. Collapsed variational inference for sum-product networks. In ICML, 2016.

12