Knowledge-Based Regularization in Generative Modeling

1 Introduction

Generative modeling plays a key role in many scientific and engineering applications. It has often been discussed with the notion of graphical models [Koller and Friedman, 2009], which are built upon knowledge of conditional independence of features. We have also seen the advances in neural network techniques for generative modeling, such as variational autoencoders (VAEs) [Kingma and Welling, 2014] and generative adversarial networks (GANs) [Goodfellow et al., 2014].

For efficient generative modeling, a model should be designed following prior knowledge of target phenomena. If one knows the conditional independence of features, it can be encoded as a graphical model. If there are some insights on the physics of the phenomena, a model can be built based on a known form of differential equations. Otherwise, a special neural network architecture may be created in accordance with the knowledge. However, although some tools (e.g., [Koller and Pfeffer, 1997]) have been suggested, it is often labor-intensive and possibly even infeasible to design a generative model meticulously so that the prior knowledge specific to each problem instance is hard-coded in the model.

∗Contact Author

On the other hand, we may employ general-purpose generative models, such as kernel-based models and neural networks with common architectures (multi-layer perceptrons, convolutional nets, and recurrent nets, etc.). Meanwhile, a general-purpose model, if used as is, is less efficient than specially-designed ones in terms of sample complexity because of large hypothesis space. Therefore, we want to incorporate as much prior knowledge as possible into a model, somehow avoiding the direct model design.

A promising way to incorporate prior knowledge into general-purpose models is via regularization. For example, structured sparsity regularization [Huang et al., 2011] is known as a method for regularizing linear models with a rigorous theoretical background. In a related context, posterior regularization [Ganchev et al., 2010; Zhu et al., 2014] has been discussed as a methodology to impose expectation constraints on learned distributions. However, the applicability of the existing regularization methods is still limited in terms of target types of prior knowledge and base models.

In generative modeling, we often have prior knowledge of relationship between features. For instance, consider generative modeling of sensor data, which is useful for tasks such as control and fault detection. In many instances, we know to some extent how units and sensors in a plant are connected (Figure 1a), from which the dependence of a part of features (here, sensor readings) can be anticipated. In the example of Figure 1a, features x1 and x2 would be more or equally dependent than x1 and x3 because of external disturbances in the processes between the units. Another example is when we know the pairwise relationship of features (Figure 1b), which can be derived from some side information.

The point is that we rarely know the full data structure, which is necessary for building a graphical model. Instead, we know a part of the data structure or pairwise relationship of some features, and there remain unknown parts that should be modeled using general-purpose models. Most existing methods cannot deal with such partial knowledge straightforwardly. In this paper, we propose a regularizer to incorporate such knowledge into general-purpose generative models based on the idea that statistical dependence of features can be (partly) anticipated from the relationship of features. The use of this type of knowledge has been actively discussed in several contexts of machine learning, but there have been surprisingly few studies in the context of generative modeling.

The proposed regularizer is defined using a kernel-based criterion of dependence [Gretton et al., 2005], which is advantageous because its gradients can be computed using standard backpropagation without additional iterative procedures. Consequently, it can be incorporated into off-the-shelf gradient-based learning methods of many general-purpose generative models, such as latent variable models (e.g., factor analysis, topic models, etc.), VAEs, and GANs. We conducted experiments using multiple datasets and generative models, and the results showcase that a model regularized using prior knowledge of feature relations achieves better generalization. The proposed method can provide a trade-off between performance and workload for model design.

2 Background

In this section, we briefly review two technical building blocks: generative modeling and dependence criteria.

2.1 Generative Modeling

We use the term "learning generative models" in the sense that we are to learn pθ(x) explicitly or implicitly. Here, x denotes the observed variable, and θ is the set of parameters to be estimated. Many popular generative models are built as Bayesian networks and Markov random fields [Koller and Friedman, 2009]. Also, the advances in deep learning technique include VAEs [Kingma and Welling, 2014], GANs (e.g., [Goodfellow et al., 2014]), autoregressive models (e.g., [van den Oord et al., 2016]), and normalizing flows (e.g., [Dinh et al., 2018]).

The learning strategies for generative models are usually based on minimization of some loss function L(θ):

A typical loss function is the negative log-likelihood or its approximation (e.g., ELBO in variational Bayes) [Koller and Friedman, 2009; Kingma and Welling, 2014]. Another class of loss functions is those designed to perform the comparison of model- and data-distributions, which has been studied recently often in the context of learning implicit generative models that have no explicit expression of likelihood. For example, GANs [Goodfellow et al., 2014] are learned via a two-player game between a discriminator and a generator.

2.2 Dependence Criteria

Among several measures of statistical dependence of random variables, we adopt a kernel-based method, namely Hilbert-Schmidt independence criterion (HSIC) [Gretton et al., 2005]. The advantage of HSIC is discussed later in Section 3.3. Below we review the basic concepts.

HSIC is defined and computed as follows [Gretton et al., 2005]. Let pxy be a joint measure over (X × Y, Γ × Λ), where X and Y are separable spaces, Γ and Λ are Borel sets on X and Y, and (X , Γ) and (Y, Λ) are furnished with probability measure px and py, respectively. Given reproducing kernel Hilbert spaces (RKHSs) F and G on X and Y, respectively, HSIC is defined as the squared Hilbert-Schmidt norm of a cross-covariance operator Cxy, i.e., HSIC(F, G, pxy) := ‖Cxy‖2HS . When bounded kernels k and l are uniquely associated with the RKHSs, F and G, respectively,

HSIC works as a dependence measure because HSIC = 0 if and only if two random variables are statistically independent (Theorem 4, [Gretton et al., 2005]).

HSIC can be empirically estimated using a dataset D = {(x1, y1), . . . , (xm, ym)} from (2). Gretton et al. [Gretton et al., 2005] presented a biased estimator with O(m-1 ) bias, and Song et al. [Song et al., 2012] suggested an unbiased estimator. In what follows, we denote an empirical estimation of HSIC(F, G, pxy) computed using a dataset D by ̂ HSICx,y(D). While the computation of empirical HSICs requires O(m2 ) operations, it may be sped up by methods such as random Fourier features [Zhang et al., 2018].

As later discussed in Section 3.2, we particularly exploit relative dependence of the features. Bounliphone et al. [Bounliphone et al., 2015] discussed the use of HSIC for a test of relative dependence. Here, we introduce H as a separable RKHS on another separable space Z. Then, the test of relative dependence (i.e., which pair is more dependent, (x, y) or (x, z)?) is formulated with the null and alternative hypotheses:H0 : ρx,y,z ≤ 0 and H1 : ρx,y,z > 0,

where ρx,y,z := HSIC(F, G, pxy) - HSIC(F,H, pxz). For this test, they define a statistic

From the asymptotic distribution of empirical HSIC, it is known [Bounliphone et al., 2015] that a conservative estimation of p-value of the test with ˆ ρx,y,z is obtained as

where Φ is the CDF of the standard normal, and σ's denote the standard deviations of the asymptotic distribution of HSIC. See [Bounliphone et al., 2015] for the definition of σ2xy , σ2xz , and σxyxz , which can also be estimated empirically. Here, (4) means that under the null hypothesis H0, the probability that ρx,y,z is greater than or equal to ˆ ρx,y,z is bounded by this value. In the proposed method, we exploit this fact to formulate a regularization term.

3 ProposedMethod

First, we manifest the type of generative models to which the proposed method applies. Then, we define what we expect to have as prior knowledge. Finally, we present the proposed regularizer and give its interpretation.

3.1 Target Type of Generative Models

Our regularization method is agnostic of the original loss function of generative modeling and the parametrization, as long as the following (informal) conditions are satisfied:

Assumption 1. Samples from pθ(x), namely ˆx, can be drawn with an admissible computational cost.

Assumption 2. Gradients ∇θEpθ(x) [f(x)] can be (approximately) computed with an admissible computational cost.

While Assumption 1 is satisfied in most generative models, Assumption 2 is a little less obvious. We note that the efficient computation of ∇θEpθ(x) [f(x)] is often inherently easy or facilitated with techniques such as the reparameterization trick [Kingma and Welling, 2014] or its variants, and thus these assumptions are satisfied by many popular methods such as factor analysis, VAEs, and GANs.

3.2 Target Type of Prior Knowledge

We suppose that we know the (partial) relationship of features that can be encoded as plausible relative dependence of the features. Such knowledge is frequently available in many practices of generative modeling, and in fact, the use of such knowledge has also been considered in other contexts of machine learning (see Section 4.2). Below we introduce several motivating examples and then give a formal definition.

Motivating Examples

The first motivating example is when we know (a part of) the data-generating process. Suppose to learn a generative model of sensor data of an industrial plant. We usually know how units and sensors in the plant are connected, which is an important source of prior knowledge. In Figure 1a, suppose that Unit i has internal state zi for i = 1, . . . , 3. Also, suppose there are relations z2 = h12(z1, ω1) and z3 = h2(z2, ω2), where ω's are random noises and h's are functions of physical processes. Then, (p(z1), p(z2)) would be more statistically dependent than (p(z1), p(z3)) because of the presence of ω2. If the sensor readings, {xi}, are determined by xi = g(zi) with observation functions {gi}, then, (p(x1), p(x2)) would be more statistically dependent than (p(x1), p(x3)) analogously to z. Figure 1a is revisited in Example 1. This type of prior knowledge is often available also for physical and biological phenomena.

Another type of example is when we know pairwise similarity or dissimilarity of features (Figure 1b). For example, we may estimate similarities of distributed sensors from their locations. Also, we may anticipate similarities of words using ontology or word embeddings. Moreover, we may know the similarities of molecules from their descriptions in the chemical compound analysis. We can anticipate relative dependence from such information, i.e., directly similar feature pairs are more dependent than dissimilar feature pairs. Figure 1b is revisited in Example 1. Here, we do not have to know every pairwise relation; knowledge on some feature pairs is sufficient to anticipate the relative dependence.

Definition of Prior Knowledge

As the exact degree of feature dependence can hardly be described precisely from prior knowledge, we use the dependence of a pair of features relative to another pair. This idea is formally written as follows:

Definition 1 (Knowledge of feature dependence). Suppose that the observed random variable, x, is a tuple of d random variables, i.e., x = (x1, . . . , xd). Knowledge of feature dependence is described as a set of triples:

where J ref , J +, J - ⊂ {1, . . . , d} are index sets. The semantics are as follows. Let J = {i1, . . . , i|J|}, and let xJ be the subtuple of x by J, i.e., xJ = (xi1 , . . . , xi|J| ). Then, triple (J ref s , J + s , J - s ) encodes the following piece of knowledge: "xJ ref s is more dependent on xJ + s than on xJ - s ."

Example 1. K = {({1}, {2}, {3})} in Figure 1a and K = {({1}, {2}, {3, 4})} in Figure 1b.

More examples of K are in Section 5.

3.3 Proposed RegularizationMethod

We define the proposed regularization method and give its interpretation as a probabilistic penalty method.

Definition of Regularizer

We want to force a generative model pθ(x) to follow the relations encoded in K. To this end, the order of HSIC between the marginals of pθ(x) should be as consistent as possible to the relations in K. More concretely, the HSIC of (pθ(xJ ref ,s), pθ(xJ +,s)) should be larger than that of (pθ(xJ ref ,s), pθ(xJ -,s)), for s = 1, . . . , |K|.

Because directly imposing constraints on the true HSIC of the marginals is intractable, we resort to a penalty method using empirical HSIC. Concretely, we add a term that penalizes the violation of the knowledge in K as follows.

Definition 2 (Knowledge-based regularization). Let L(θ) be the original loss function in (1). The learning problem regularized using K in (5) is posed as

where λ ≥ 0 is a regularization hyperparameter, and

Here, να ≥ 0 is another hyperparameter. ˆ ρs(θ) is the empirical estimation of the relative HSIC corresponding to the s-th triple in K, that is,

where ˆD denotes a set of samples drawn from pθ(x).

Interpretation

By imposing the regularizer RK, we expect that the hypothesis space to be explored is made smaller implicitly, which would result in better generalization capability. This is also supported by the following interpretation of the regularizer.

While RK(θ) depends on the empirical HSIC values, making (each summand of) RK(θ) small can be interpreted as imposing "probabilistic" constraints on the order of the true HSIC values via a penalty method. Now let ρs denote the true value of ˆ ρs (we omit argument θ for simplicity), i.e.,

where FJ denotes a separable RKHS on the space of xJ. Moreover, consider a test with hypotheses:

for s = 1, . . . , |K|. On the test with Hs,0 and Hs,1, the following proposition holds:

Proposition 1. If RK,s = 0 is achieved, then the null hypothesis Hs,0 can be rejected in favor of the alternative Hs,1 with p-value ps upper bounded by

Proof. From (4), we have ps ≤ 1-Φ(ˆ ρs/τs), where

Also, when RK,s = 0, we have ˆ ρs ≥ να. Consequently, as Φ is monotonically increasing, we have (11).

Choice of Hyperparameters

The proposed regularizer has two hyperparameters, λ and να. While λ is a standard regularization parameter that balances the loss and the regularizer, να can be interpreted as follows. Now let α := 1-Φ(να/τs) (i.e., the right-hand side of (11)), which corresponds to the required significance of our probabilistic constraints. We can determine 0 < α < 1 based on the plausibility of the prior knowledge, such as α = 0.05 or α = 0.1 (a smaller α requires higher significance). Also, τs can be estimated as it is defined with the variances of the asymptotic distribution of HSIC [Bounliphone et al., 2015]. Hence, as να = τsΦ-1 (1-α), we can roughly determine the value of να that corresponds to a specific significance, α. The need to tune να (or α) remains, but the above interpretation is useful to determine the search range for να.

Optimization Method

From Assumptions 1 and 2, the gradient of ˆ ρs with regard to θ can be computed with the backpropagation via the samples from pθ(x), which are used to compute empirical HSIC. Hence, if the solution of the original optimization, (1), is obtained via a gradient-based method, the regularized version, (6), can also be solved using the same gradient-based method. The situation is especially simplified if the dataset is a set of d-dimensional vectors, i.e., D = {xi ∈ Rd | i = 1, . . . , n}. In the optimization process, sample ˆD = {ˆxi ∈ Rd | i = 1, . . . ,m} from pθ(x) being learned, and use them for computing the the regularization term in (7) and their gradients.

Algorithm 1 Knowledge-regularized gradient method

Input: Data D = {xi}, knowledge K = {(J ref s , J + s , J - s )}, hyperparameters λ ≥ 0, να ≥ 0, sample size m Output: A set of parameters θ∗ of pθ(x) 1: initialize θ 2: repeat 3: draw ˆD = {ˆxi ∈ Rd | i = 1, . . . ,m} from pθ(x) 4: for s = 1, . . . , |K| do 5: Js ← J ref s ∪ J + s ∪ J - s 6: ˆDs ← {[ˆxi]Js ∈ R|Js| | i = 1, . . . ,m} ▷ [ˆxi]J is the subvector of ˆxi indexed by J 7: compute ∂ ˆ ρs/∂θ using ˆDs 8: end for 9: compute ∂RK/∂θ using {∂ ˆ ρs/∂θ | s = 1, . . . , |K|} 10: update θ using ∂L/∂θ +λ∂RK/∂θ 11: until convergence

Then, incorporate them into the original gradient-based updates. These procedures are summarized in Algorithm 1.

The computation of RK requires O(m2 |K|) operations in a naive implementation, where m is the number of samples drawn from pθ(x). This will be burdensome for a very large |K|, which can be alleviated by carefully choosing J ref 's so that the computed HSIC values can be reused for many times.

4 RelatedWork

4.1 Generative Modeling with Prior Knowledge

A perspective related to this work is model design based on prior knowledge. For example, the object-oriented Bayesian network language [Koller and Pfeffer, 1997] helps graphical model designs when the structures behind data can be described in an object-oriented way. Such tools are useful when we can prepare prior knowledge enough for building a full model. However, it is not apparent how to utilize them when knowledge is only partially given, which is often the case in practice. In contrast, our method can exploit prior knowledge even if only a part of the structures is known in advance. Another related perspective is the structure learning of Bayesian networks with constraints from prior knowledge [Cussens et al., 2017; Li and van Beek, 2018].

Posterior regularization (PR) [Ganchev et al., 2010; Zhu et al., 2014; Mei et al., 2014; Hu et al., 2018] is known as a framework for incorporating prior knowledge into probabilistic models. In fact, our method can be included in the most general class of PR, which was briefly mentioned in [Zhu et al., 2014]. However, practically, existing work [Ganchev et al., 2010; Zhu et al., 2014; Mei et al., 2014; Hu et al., 2018] only considered limited cases of PR, where constraints were written in terms of a linear monomial of expectation with regard to the target distribution. In contrast, our method tries to fulfill the constraints on statistical dependence, which needs more complex expressions.

The work by Lopez et al. [2018] should be understood as a kind of technical complement of ours (and vice versa). While Lopez et al. [2018] regularize the amortized inference of VAEs to make independence of latent variables, we regularize a general generative model itself to make dependence of observed variables compatible with prior knowledge. In other words, the former regularizes an encoder, whereas the latter regularizes a decoder. Moreover, our method is more flexible than [Lopez et al., 2018] in terms of applicable type of knowledge; they considered only independence, but our method can incorporate both dependence and independence. Also, while Lopez et al. [2018] discussed the method only for VAEs, our method applies to any generative models as long as the mild assumptions in Section 3.1 are satisfied.

4.2 Knowledge of Feature Dependence

In fact, the use of knowledge of feature relation has been discussed in different or more specific contexts. In natural language processing, feature similarity is often available as the similarity between words. Xie et al. [Xie et al., 2015] exploited correlations of words for topic modeling, and Liu et al. [Liu et al., 2015] incorporated semantic similarity of words in learning word embeddings. These are somewhat related to generative modeling, but the scope of data and model is limited. Moreover, there are several pieces of research on utilizing feature similarity [Krupka and Tishby, 2007; Li and Li, 2008; Sandler et al., 2009; Li et al., 2017; Mollaysa et al., 2017] for discriminative problems. Despite these interests, there have been surprisingly few studies on the use of such knowledge for generative modeling in general.

5 Experiments

5.1 Datasets and Prior Knowledge

We used the following four datasets, for which we can prepare plausible prior knowledge of feature relations.

Toy data. We created a cylinder-like toy dataset as exemplified in Figure 2. On this dataset, we know that x1 and x2 are more statistically dependent than x1 and x3 are. This knowledge is expressed in the manner of (5) as

Constrained concatenated MNIST. We created a new dataset from MNIST, namely constrained concatenated MNIST (CCMNIST, examples in Figure 3) as follows. Each image in this dataset was generated by vertically concatenating three MNIST images. Here, the concatenation is constrained so that the top and middle images have the same label, whereas the bottom image is independent of the top two. We created training and validation sets fromMNIST's training dataset and a test set fromMNIST's test dataset. On this dataset, we can anticipate that the top third and the middle third of each image are more dependent than the top and the bottom are. If each image is vectorized in a row-wise fashion as a 2352-dim vector, this knowledge is expressed as

Plant sensor data. We used a simulated sensor dataset (PLANT), which corresponds to the example in Figure 1a. The simulation is based on a real industrial chemical plant called the Tennessee Eastman process [Downs and Vogel, 1993]. In the plant, there are four major unit operations: reactor, vapor-liquid separator, recycle compressor, and product stripper. On each unit operation, sensors such as level sensors and thermometers are attached. We used readings of 22 sensors. On this dataset, we can anticipate the relative dependence between sets of sensors based on the process diagram of the Tennessee Eastman process [Downs and Vogel, 1993]. The knowledge set for this dataset (KPlant) contains, e.g.,

where Jcompressor is the set of indices of sensors attached to the compressor, and so on. This is plausible because the reactor is directly connected to the compressor while the stripper is not. We prepared KPlant comprising |KPlant| = 12 relations based on the structure of the plant [Downs and Vogel, 1993].

Solar energy production data. We used the records of solar power production 1 (SOLAR) of 137 solar power plants in Alabama in June 2006. The data of the first 20 days were used for training, the next five days were for validation, and the last five days were for test. This dataset corresponds to the example in Figure 1b because we can anticipate the dependence of features from the pairwise distances between the solar plants. We created a knowledge set KSolar as follows; for the i-th plant (i = 1, . . . , 137), if the nearest plant is within 10 [km] and the distance to the second-nearest one is more than 12 [km], then add ({i}, {ji}, {ki}) to KSolar, where ji and ki are the indices of the nearest and the second-nearest plants, respectively. This resulted in |KSolar| = 31.

5.2 Settings

Models. We used factor analysis (FA, see also [Zhang, 2009]), VAE, and GAN. We trained FA and VAEs on all datasets and GANs on the TOY dataset. For VAEs and GANs, encoders, decoders, and discriminators were modeled using multi-layer perceptrons (MLPs) with three hidden layers. Hereafter we denote the dimensionality of VAE/GAN's latent variable by dim(z) and the dimensionality of MLP's hidden layer by dim(MLP). We tried dim(MLP) = 25 , 26 , . . . , 210 and used several values of dim(z) for the different datasets.

1 www.nrel.gov/grid/solar-power-data.html

Significant diff. from "L2 only" with italic: p < .01 or bold: p < .001.

Table 2: Test mean per-feature reconstruction errors by VAEs on the three datasets. Here only the cases of dim(MLP) = 25 , 26 , 27 (resp. for the three datasets) are reported. Other cases are similar.

Baselines. As the simplest baseline, we trained everymodel only with L2 regularization (i.e., weight decay). We tried another baseline, in which the weight of the last layer of decoder MLPs are regularized using K; e.g., for KToy = {{1}, {2}, {3}}, we used a regularization term

where wi denotes the i-th row of the weight matrix of the last layer of the decoder MLP, and ϵ was tuned in the same manner as να. We refer to this baseline as "out-layer." It cannot be applied to the PLANT dataset because the feature sets in KPlant do not have any one-to-one correspondences.

Hyperparameters. We computed HSIC with m= 128 using Gaussian kernels with the bandwidth set by the median heuristics. No improvement was observed with m > 128. The hyperparameters were chosen based on the performance on the validation sets. The search was not intensive; λ was chosen from three candidate values that roughly adjust orders of L and RK, and να was chosen from .01 or .05.

5.3 Results

Below, the quantitative results are reported mainly with cross entropy (for CCMNIST) or reconstruction errors (for the others) as they are a universal performance criterion to examine the generalization capability. For VAEs, the significance of the improvement by the proposed method did not change even when we examined the ELBO values.

Evaluation by test set performance. The test set performance of VAEs are shown in Tables 1 and 2. We can observe that, while the improvement by the out-layer baseline was quite marginal, the proposed regularization method resulted in significant improvement. The performance of FA slightly improved (details omitted due to space limitations).

Comparison to a dedicated model. We compared the performance of the proposed method to that of a model designed specifically for CCMNIST dataset (termed "dedicated"). The dedicated model uses an MLP designed in accordance with the prior knowledge that the top and the middle parts are from the same digit. The performance of the dedicated model is shown in the right-most column of Table 1. We can observe that the proposed method performed intermediately between the most general case (L2 only) and the dedicated model. As designing models meticulously is time-consuming and often even infeasible, the proposed method will be useful to give a trade-off between performance and a user's workload.

Effect of knowledge set size. Another interest lies in how the amount of provided prior knowledge affects the performance of the regularizer. We investigated this by changing the number of tuples in KSolar used in learning VAEs. We prepared knowledge subsets by extracting some tuples from the original set KSolar. When creating the subsets, tuples were chosen so that a larger subset contained all elements of the smaller ones. Figure 4 shows the test set performances along with the subset size. The performance is improved with a larger knowledge set. For example, the difference between |K| = 6 and |K| = 31 is significant with p ≈ .003.

Inspection of generated samples. We inspected the samples drawn from GANs trained on the TOY dataset. When the proposed regularizer was not applied, we often observed a phenomenon similar to "mode collapse" (see, e.g., [Metz et al., 2017]) occurred as in the left-side plots of Figure 5. In contrast, the whole geometry of the TOY dataset was always captured successfully when the proposed regularizer was applied, as in the right-side plots of Figure 5.

Running time. For example, the time needed for training a VAE for 50 epochs on the CCMNIST dataset was 150 or 117 seconds, with or without the proposed method, respectively. As the complexity is linear in |K|, this would not be prohibitive, while any speed-up techniques will be useful.

6 Conclusion

In this work, we developed a method for regularizing generative models using prior knowledge of feature dependence, which is frequently encountered in practice and has been studied in contexts other than generative modeling. The proposed regularizer can be incorporated in off-the-shelf learning methods of many generative models. The direct extension of the current method includes the use of higher-order dependence between multiple sets of features [Pfister et al., 2017].

Acknowledgments

This work was supported by JSPS KAKENHI Grant Numbers JP19K21550, JP18H06487, and JP18H03287.