{'fraction_non_alphanumeric': 0.046239976166240475, 'fraction_numerical': 0.012562377417512849, 'mean_word_length': 4.5355596784168215, 'pattern_counts': {'":': 0, '<': 2, '<?xml version=': 0, '>': 1, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 11, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Bayesian Inference with Posterior Regularization and applications to Infinite Latent SVMs

Existing Bayesian models, especially nonparametric Bayesian methods, rely on specially conceived priors to incorporate domain knowledge for discovering improved latent representations. While priors can affect posterior distributions through Bayes' rule, imposing posterior regularization is arguably more direct and in some cases more natural and general. In this paper, we present regularized Bayesian inference (RegBayes), a novel computational framework that performs posterior inference with a regularization term on the desired post-data posterior distribution under an information theoretical formulation. RegBayes is more flexible than the procedure that elicits expert knowledge via priors, and it covers both directed Bayesian networks and undirected Markov networks whose Bayesian formulation results in hybrid chain graph models. When the regularization is induced from a linear operator on the posterior distributions, such as the expectation operator, we present a general convex-analysis theorem to characterize the solution of RegBayes. Furthermore, we present two concrete examples of RegBayes, infinite latent support vector machines (iLSVM) and multi-task infinite latent support vector machines (MT-iLSVM), which explore the large-margin idea in combination with a nonparametric Bayesian model for discovering predictive latent features for classification and multi-task learning, respectively. We present efficient inference methods and report empirical studies on several benchmark datasets, which appear to demonstrate the merits inherited from both large-margin learning and Bayesian nonparametrics. Such results were not available until now, and contribute to push forward the interface between these two important subfields, which have been largely treated as isolated in the community.


Introduction
Over the past decade, nonparametric Bayesian models have gained remarkable popularity in machine learning and other fields, partly owing to their desirable utility as a "nonparametric" prior distribution for a wide variety of probabilistic models, thereby turning the largely heuristic model selection practice, such as determining the unknown number of components in a mixture model (Antoniak, 1974) or the unknown dimensionality of latent features in a factor analysis model (Griffiths and Ghahramani, 2005), as a Bayesian inference problem in an unbounded model space. Popular examples include Gaussian process (GP) , Dirichlet process (DP) (Ferguson, 1973;Antoniak, 1974), and Beta process (BP) (Thibaux and Jordan, 2007). DP is often described with a Chinese restaurant process (CRP) metaphor, and similarly BP is often described with an Indian buffet process (IBP) metaphor (Griffiths and Ghahramani, 2005). Such nonparametric Bayesian approaches allow the model complexity to grow as more data are observed, which is a key factor differing them from other traditional "parametric" Bayesian models.
One recent development in practicing Bayesian nonparametrics is to relax some unrealistic assumptions on data, such as homogeneity and exchangeability. For example, to handle heterogenous observations, predictor-dependent processes (MacEachern, 1999;Williamson et al., 2010) have been proposed; and to relax the exchangeability assumption, stochastic processes with various correlation structures, such as hierarchical structures (Teh et al., 2006), temporal or spatial dependencies (Beal et al., 2002;Blei and Frazier, 2010), and stochastic ordering dependencies (Hoff, 2003;Dunson and Peddada, 2007), have been introduced. A common principle shared by these approaches is that they rely on defining, or in some unusual cases learning (Welling et al., 2012) a nonparametric Bayesian prior 1 encoding some special structures, which indirectly 2 influences the posterior distribution of interest through an interplay with a likelihood model according to the Bayes' rule (also known as Bayes' theorem). In this paper, we explore a different principle known as posterior regularization, which offers an additional and arguably richer and more flexible set of means to augment a posterior distribution under rich side information, such as predictive margin, structural bias, etc., which can be harder, if possible, to be captured by a Bayesian prior. Let Θ denote model parameters and H denote hidden variables. Then given a set of observed data D, posterior regularization (Ganchev et al., 2010) is generally defined as solving a regularized maximum likelihood estimation (MLE) problem: Posterior Regularization : max Θ L(Θ; D) + Ω(p(H|D, Θ)), where L(Θ; D) is the marginal likelihood of D, and Ω(·) is a regularization function of the model posterior over latent variables (note that here we view posterior as a generic post-data distribution on hidden variables in the sense of (Ghosh and Ramamoorthi, 2003, pp.15), not necessarily corresponding to a Bayesian posterior that must be induced by the Bayes' rule). The regularizer can be defined as a KL-divergence between a desired distribution with certain properties over latent variables and the model posterior in question, or other constraints on the model posterior, such as those used in generalized expectation (Mann and McCallum, 2010) or constraint-driven semi-supervised learning (Chang et al., 2007). An EM-type procedure can be applied to solve Eq. (1) approximately, and obtain an augmented MLE of the hidden variable model: p(H|D, Θ MLE ). When a distribution over the model parameter is of interest, going beyond the classical Bayesian theory, recent attempts toward learning a regularized posterior distribution of model parameters (and latent variables as well if present) include the "learning from measurements" (Liang et al., 2009), maximum entropy discrimination (MED) ) and maximum entropy discrimination latent Dirichlet allocation (MedLDA) . All these methods are parametric in that they give rise to distributions over a fixed and finite-dimensional parameter space. To the best of our knowledge, very few attempts have been made to impose posterior regularization in a nonparametric setting where model complexity depends on data, such as the case for nonparametric Bayesian latent variable models. A general formalism for (parametric and nonparametric) Bayesian inference with posterior regularization seems to be not yet available or apparent. In this paper, we present such a formalism, which we call regularized Bayesian inference, or RegBayes, built on the convex duality theory over distribution function spaces; and we apply this formalism to learn regularized posteriors under the Indian buffet process (IBP), conjoining two powerful machine learning paradigms, nonparametric Bayesian inference and SVM-style max-margin constrained optimization. Unlike the regularized MLE formulation in Eq.
(1), under the traditional formulation of Bayesian inference one is not directly optimizing an objective with respect to the posterior. To enable a regularized optimization formulation of RegBayes, we begin with a variational reformulation of the Bayes' theorem, and define L(q(M|D)) as the KL-divergence between a desired post-data posterior q(M|D) over model M and the standard Bayesian posterior p(M|D) (see Section 3.1 for a recapitulation of the connection between KL-minimization and Bayes' theorem). RegBayes solves the following optimization problem: RegBayes : inf q(M|D)∈P prob L(q(M|D)) + Ω(q(M|D)), where the regularization Ω(·) is a function of the post-data posterior q(M|D), and P prob is the feasible space of well-defined distributions. By appropriately defining the model and its prior distribution, RegBayes can be instantiated to perform either parametric and nonparametric regularized Bayesian inference. One particularly interesting way to derive the posterior regularization is to impose posterior constraints. Let ξ denote slack variables and P post (ξ) denote the general soft posterior constraints (see Section 3.2 for a formal description), then, we can express the regularization term variationally: where U (ξ) is normally defined as a convex penalty function. The RegBayes formalism defined in Eq.
(2) applies to a wide spectrum of models, including directed graphical models (i.e., Bayesian networks) and undirected Markov networks. For undirected models, when performing Bayesian inference the resulting posterior takes the form of a hybrid chain graphical model (Frydenberg, 1990) (Murray and Ghahramani, 2004;Qi et al., 2005;Welling and Parise, 2006), which is usually much more challenging to regularize than for Bayesian inference with directed GMs. When the regularization term is convex and induced from a linear operator (e.g., expectation) of the posterior distributions, RegBayes can be solved with convex analysis theory.
By allowing direct regularization over posterior distributions, RegBayes provides a significant source of extra flexibility for post-data posterior inference, which applies to both parametric and nonparametric Bayesian learning (see the remarks after the main Theorem 6). In this paper, we focus on applying this technique to the later case, and illustrate how to use RegBayes to facilitate integration of Bayesian nonparametrics and large-margin learning, which have complementary advantages but have been largely treated as two disjoint subfields. Previously, it has been shown that, the core ideas of support vector machines (Vapnik, 1995) and maximum entropy discrimination , as well as their structured extensions to the max-margin Markov networks (Taskar et al., 2003) and maximum entropy discrimination Markov networks , have led to successful outcomes in many scenarios. But a large-margin model rarely has the flexibility of nonparametric Bayesian models to automatically handle model complexity from data, especially when latent variables are present (Jebara, 2001;. In this paper, we intend to bridge this gap using the RegBayes principle. Specifically, we develop the infinite latent support vector machines (iLSVM) and multitask infinite latent support vector machines (MT-iLSVM), which explore the discriminative large-margin idea to learn infinite latent feature models for classification and multi-task learning (Argyriou et al., 2007;Bakker and Heskes, 2003), respectively. We show that both models can be readily instantiated from the RegBayes master equation (2) by defining appropriate posterior regularization using the large-margin principle, and by employing an appropriate prior. For iLSVM, we use the IBP prior to allow the model to have an unbounded number of latent features a priori. For MT-iLSVM, we use a similar IBP prior to infer a latent projection matrix to capture the correlations among multiple predictive tasks while avoiding pre-specifying the dimensionality of the projection matrix. The regularized inference problems can be efficiently solved with an iterative procedure, which leverages existing high-performance convex optimization techniques.
The rest of the paper is organized as follows. Section 2 discusses related work. Section 3 presents regularized Bayesian inference (RegBayes), together with the convex duality results that will be needed in latter sections. Section 4 concretizes the ideas of RegBayes and presents two infinite latent feature models with large-margin constraints for both classification and multi-task learning. Section 5 presents some preliminary experimental results. Finally, Section 6 concludes and discusses future research directions.

Related Work
Expectation regularization or expectation constraints have been considered to regularize model parameter estimation in the context of semi-supervised learning or learning with weakly labeled data. Mann and McCallum (Mann and McCallum, 2010) summarized the recent developments of the generalized expectation (GE) criteria for training a discriminative probabilistic model (e.g., maximum entropy models or conditional random fields (Lafferty et al., 2001)) with unlabeled data. By providing appropriate side information, such as labeled features or estimates of label distributions, a GE-based penalty function is defined to regularize the model distribution, e.g., the distribution of class labels. One commonly used GE function is the KL-divergence between empirical expectation and model expectation of some feature functions if the expectations are normalized or the gen-eral Bregman divergence for unnormalized expectations. Although the GE criteria can be used alone as a scoring function to estimate the unknown parameters of a discriminative model, it is more usually used as a regularization term to an estimation method, such as maximum (conditional) likelihood estimation. Bellare et al. (Bellare et al., 2009) presented a different formulation of using expectation constraints in semi-supervised learning by introducing an auxiliary distribution to GE, together with an alternating projection algorithm, which can be more efficient. Liang et al. (Liang et al., 2009) proposed to use the general notion of "measurements" to encapsulate the variety of weakly labeled data for learning exponential family models. The measurements can be labels, partial labels or other constraints on model predictions. Under the EM framework, posterior constraints were used in (Graca et al., 2007) to modify the E-step of an EM algorithm to project model posterior distributions onto the subspace of distributions that satisfy a set of auxiliary constraints.
Dudik et al. (Dudík et al., 2007) studied the generalized maximum entropy principle with a rich form of expectation constraints using convex duality theory, where the standard moment matching constraints of maximum entropy are relaxed to inequality constraints. But their analysis was restricted to KL-divergence minimization (maximum entropy is a special case) and the finite dimensional space of observations. Later on, Altun and Smola (Altun and Smola, 2006) presented a more general duality theory for a family of divergence functions on Banach spaces. We have drawn inspiration from both papers to develop the regularized Bayesian inference framework using convex duality theory.
When using large-margin posterior regularization, RegBayes generalizes the previous work on maximum entropy discrimination ). The present paper provides a full extension of our preliminary work on max-margin nonparametric Bayesian models (Zhu et al., 2011b,a). For example, the infinite SVM (iSVM) (Zhu et al., 2011b) is a latent class model, where each data example is assigned to a single mixture component (i.e., an 1-dimensional space), and both iLSVM and MT-iLSVM extend the ideas to infinite latent feature models. For multi-task learning, nonparametric Bayesian models have been developed in (Xue et al., 2007;Rai and Daume III, 2010) for learning features shared by multiple tasks. However, these methods are based on standard Bayesian inference without a posterior regularization using, for example, the large-margin constraints. Finally, MT-iLSVM can be also regarded as a nonparametric Bayesian formulation of the popular multi-task learning methods (Ando and Zhang, 2005;Jebara, 2011).

Regularized Bayesian Inference
We begin by laying out a general formulation of regularized Bayesian inference, using an optimization framework built on convex duality theory.

Variational formulation of Bayes' theorem
We first derive an optimization-theoretic reformulation of the Bayes' theorem. Let M denote the space of feasible models, and M ∈ M represents an atom in this space. We assume that M is a complete separable metric space endowed with its Borel σ-algebra B(M). Let Π be a distribution (i.e., a probability measure) on the measurable space (M, B(M)). We assume that Π is absolutely continuous with respect to some background measure µ, so that there exists a density π such that dΠ = πdµ. Let D = {x n } N n=1 be a collection of observed data, which we assume to be i.i.d. given a model. Let P (·|M) be the likelihood distribution, which is assumed to be dominated by a σ-finite measure λ for all M with positive density, so that there exists a density p(·|M) such that dP (·|M) = p(·|M)dλ. Then, the Bayes' conditionalization rule gives a posterior distribution with the density (Ghosh and Ramamoorthi, 2003, Chap.1.3): a density over M with respect to the base measure µ, where p(D) is the marginal likelihood of the observed data. For reasons to be clear shortly, we now introduce a variational formulation of the Bayes' theorem. Let Q are an arbitrary distribution on the measurable space (M, B(M)). We assume that Q is absolutely continuous with respect to Π and denote by q its density with respect to the background measure µ. 3 It can be shown that the posterior distribution of M due to the Bayes' theorem is equivalent to the optimum solution of the following convex optimization problem: where KL(q(M) π(M)) = M q(M) log(q(M)/π(M))dµ(M) is the Kullback-Leibler (KL) divergence from q(·) to π(·), and P prob represents the feasible space of all density functions over M with respect to the measure µ. The proof is straightforward by noticing that the objective will become KL(q(M) p(M|D)) by adding the constant log p(D). It is noteworthy that q(M) here represents the density of a general post-data posterior distribution in the sense of (Ghosh and Ramamoorthi, 2003, pp.15), not necessarily corresponding to a Bayesian posterior that is induced by the Bayes' rule. As we shall see soon later, when we introduce additional constraints, the post-data posterior q(M) is different from the Bayesian posterior p(M|D), and moreover, it could even not be obtainable from any Bayesian conditionalization in a different model. In the sequel, in order to distinguish q(·) from the Bayesian posterior, we will call it post-data distribution 4 in short or post-data posterior distribution in full. For notation simplicity, we have omitted the condition D in the post-data posterior distribution q(M).
Remark 1 The optimization formulation in (5) implies that Bayes' rule is an information projection procedure that projects a prior density to a post-data posterior by taking account of the observed data. In general, Bayes's rule is a special case of the principle of minimum information (Williams, 1980). In contrast, we have many (normally infinite for continuous ξ) feasible subspaces for soft constraints and each of them is associated with a different complexity or penalty, measured by the U function.

Regularized Bayesian Inference with Expectation Constraints
In the variational formulation of Bayes' rule in Eq. (5), the constraints on q(M) ensure that q is well-normalized and the objective is well-defined, i.e., q(M) ∈ P prob , which do not capture any domain knowledge or structures of the model or data. Some previous efforts have been devoted to eliciting domain knowledge by constraining the prior or the base measure µ (Robert, 1995;Garthwaite et al., 2005). As we shall see, such constraints without considering data are special cases of RegBayes to be presented. Specifically, the optimization-based formulation of Bayes' rule makes it straightforward to generalize Bayesian inference to a richer type of posterior inference, by replacing the standard normality constraint on q with a wide spectrum of knowledge-driven and/or datadriven constraints or regularization. (To contrast, we will refer to the problem in Eq. (5) as "unconstrained" or "unregularized".) Formally, we define regularized Bayesian inference (RegBayes) as a generalized posterior inference procedure that solves a constrained optimization problem due to such additional regularization imposed on q: where P post (ξ) is a subspace of distributions that satisfy a set of additional constraints besides the standard normality constraint of a probability distribution. Using the variational formulation in Eq. Obviously this formulation enables different types of constraints to be employed in practice. In this paper, we focus on the expectation constraints, of which each one is a function of q(M) through an expectation operator. For instance, let ψ = (ψ 1 , · · · , ψ T ) be a vector of feature functions, each of which is ψ t (M; D) defined on M and possibly data dependent. Then a subspace of feasible post-data distributions can be defined in the following form: where E is the expectation operator that maps q(M) to a point in the space R T , and for each feature function ψ t : Eq(ψ t ; D) The function h can be of any form in theory, though a simple h function will make the optimization problem easy to solve. The auxiliary parameters ξ are usually nonnegative and interpreted as slack variables. The constraints with non-trivial ξ are soft constraints as illustrated in Figure 1(b). But we emphasize that by defining U as an indicator function, the formulation (6) covers the case where hard constraints are imposed. For instance, if we define where I(c) is an indicator function that equals to 0 if the condition c is satisfied; otherwise ∞, then all the expectation constraints (7) are hard constraints. As illustrated in Figure 1(a), hard constraints define one single feasible subspace (assuming to be non-empty). In general, we assume that U (ξ) is a convex function, which represents a penalty on the size of the feasible subspaces, as illustrated in Figure 1(b). A larger subspace typically leads to models with a higher complexity. In the classification models to be presented, U corresponds to a surrogate loss, e.g., hinge loss of a prediction rule, as we shall see.
Similarly, the formulation of RegBayes with expectation constraints (7) can be equivalently written in an "unconstrained" form by using the rule in (3). Specifically, let g(Eq(ψ; D)) def = inf ξ U (ξ), s.t. : h(Eq(ψ t ; D)) ≤ ξ t , ∀t, we have the equivalent optimization problem: where Eq(ψ; D) is a point in R T and the t-th coordinate is Eq(ψ t ; D), a function of q(M) as defined before. We assume that the real-valued function g : R T → R is convex and lower semi-continuous. For each U , we can induce a g function by taking the infimum of U (ξ) over ξ with the posterior constraints; vice versa. If we use hard constraints, similar as in regularized maximum entropy density estimation (Altun and Smola, 2006;Dudík et al., 2007), we have For the regularization function g, as well as U , we can have many choices, besides the above mentioned indicator function. For example, if the feature function ψ t is an indicator function and we could obtain 'prior' expectations Ep[ψ t ] from domain/expert knowledge about M. If we further normalize the empirical expectations of T functions and denote the discrete distribution byp(M), one natural regularization function would be the KL-divergence between prior expectations and the expectations computed from the normalized model posterior q(M), i.e., g(Eq) = t s(Ep[ψ t ], Eq(ψ t )) = KL(p(M) q(M)), where s(x, y) = x log(x/y) for x, y ∈ (0, 1). The general Bregman divergence can be used for unnormalized expectations. This kind of regularization function has been used in (Mann and McCallum, 2010) for label regularization, in the context of semi-supervised learning. Other choices of the regularization function include the ℓ 2 2 penalty or indicator function with equality constraints (Please see Table 1 in (Dudík et al., 2007) for a summary).
Remark 2 So far, we have focused on RegBayes in the context of full Bayesian inference. Indeed, RegBayes can be generalized to apply to empirical Bayesian inference, where some model parameters need to be estimated. More generally, RegBayes applies to both directed Bayesian networks (of which the hierarchical Bayesian models we have discussed are an example) and undirected Markov random fields. But for undirected models, a RegBayes treatment will have to deal with a chain graph resultant from Bayesian inference, which is more challenging due to existence of normalization factors. We will discuss some details and examples in Appendix A.

Optimization with Convex Duality Theory
Depending on several factors, including the size of the model space, the data likelihood model, the prior distribution, and the regularization function, a RegBayes problem in general can be highly non-trivial to solve, either in the constrained or unconstrained form, as can be seen from several concrete examples of RegBayes models we will present in the next section and in the Appendix B. In this section, we present a representation theorem to characterize the solution the convex RegBayes problem (8) with expectation regularization. These theoretical results will be used later in developing concrete RegBayes models.
To make the subsequent statements general, we consider the following problem: where f : X → R is a convex function; A : X → B is a bounded linear operator; and g : B → R is also convex. Below we introduce some tools in convex analysis theory to study this problem. We begin by formulating the primal-dual space relationships of convex optimization problems in the general settings, where we assume both X and B are Banach spaces 5 . An important result we build on is the Fenchel duality theorem.
Definition 3 (Convex Conjugate) Let X be a Banach space and X * be its dual space. The convex conjugate or the Legendre-Frenchel transformation of a function f : 5. A Banach space is a vector space with a metric that allows the computation of vector length and distance between vectors. Moreover, a Cauchy sequence of vectors always converges to a well defined limit in the space.
Theorem 4 (Fenchel Duality (Borwein and Zhu, 2005)) Let X and B be Banach spaces, f : X → R ∪ {+∞} and g : B → R ∪ {+∞} be convex functions and A : X → B be a bounded linear map. Define the primal and dual values t, d by the Fenchel problems Then these values satisfy the weak duality inequality t ≥ d. If f , g and A satisfy either 0 ∈ core(domg − Adomf ) and both f and g are lower semicontinuous (lsc), (12) or then t = d and the supremum to the dual problem is attainable if finite.
Let S be a subset of a Banach space B. In the above theorem, we say s is in the core of S, denoted by s ∈ core(S), provided that ∪ λ>0 λ(S − s) = B. The Fenchel duality theorem has been applied to solve divergence minimization problems for density estimation (Altun and Smola, 2006;Dudík et al., 2007). Let ψ def = (ψ 1 , · · · , ψ T ) be a vector of feature functions. Each feature function is a mapping, ψ t : M → R. Therefore, B is the product space R T , a simple Banach space. Let X be the Banach space of finite signed measures (with total variation as the norm) that are absolutely continuous with respect to the measure µ, and let A be the expectation operator of the feature functions with respect to the distribution q on M, that is, Aq Letψ be a reference point in R T . As for density estimation, we have some observations of M here, andψ = Ap emp [ψ(M)], where p emp is the empirical distribution. Then, when the f function is a KL-divergence and the constraints are relaxed moment matching constraints, the following result can be proven.
Note that for this lemma and the ones to be presented below to hold, the problems need to meet some regularity conditions (or constraint qualifications), such as those in Theorem 4. In practice it can be difficult to check whether the constraint qualifications hold. One solution is to solve the dual optimization problem and examine if the conditions hold depending on whether the solution diverge or not (Altun and Smola, 2006). The problem in the above lemma is subject to hard constraints, therefore the corresponding g is the indicator function I( E q [ψ] −ψ B ≤ ǫ) when applying the Fenchel duality theorem. Other examples of the posterior constraints can be found in (Dudík et al., 2007;Mann and McCallum, 2010;Ganchev et al., 2010), as we have discussed in Section 3.2. In this paper, we consider the general soft constraints as defined in the RegBayes problem (Eq. (6)). Furthermore, we do not assume the existence of a fully observed dataset to compute the empirical expectationφ. Specifically, following a similar line of reasoning as in (Altun and Smola, 2006), though this time with an un-normalized p in KL(q p), we have the following result. The detailed proof is deferred to Appendix C.1.
Theorem 6 (Representation theorem of RegBayes) Let E be the expectation operator with feature functions ψ(M; D), and assume g is convex and lower semicontinuous (lsc). We have where the unique solution is given byqφ ; andφ is the solution of the dual problem; and Λφ is the log-partition function.
From the optimum solutionqφ(M), we can see that the form of the RegBayes posterior is symbolically similar to that of the Bayesian posterior; but instead of multiplying the likelihood term with a prior distribution, RegBayes introduces an extra term, exp( φ , ψ(M; D) − Λφ), whose coefficients are derived from an constrained optimization problem resultant from the constraints on the posterior. We make the following remarks. Therefore, such a constraint lets RegBayes cover the widely used truncated priors, such as truncated Gaussian (Robert, 1995). backqφ(M) through Bayesian conditionalization. Therefore, RegBayes is more flexible than the standard Bayesian inference, where the prior and likelihood model are explicitly defined, but no additional constraints or regularization can be systematically incorporated. The recent work (Mei et al., 2014) presents an example. Specifically, we show that incorporating domain knowledge via posterior regularization can lead to a flexible framework that automatically learns the importance of each piece of knowledge, thereby allowing for a robust incorporation, which is important in the scenarios where noisy knowledge is collected from crowds. In contrast, eliciting expert knowledge via fitting some priors is generally hard, especially in high-dimensional spaces, as experts are normally good at perceiving low-dimensional and well-behaved distributions but can be very bad in perceiving high-dimensional or skewed distributions (Garthwaite et al., 2005).
It is worth mentioning that although the above theorem provides a generic representation of the solution to RegBayes, in practice we usually need to make additional assumptions in order to make either the primal or dual problem tractable to solve. Since such assumptions could make the feasible space non-convex, additional cautions need to be paid. For instance, the mean-field assumptions will lead to a non-convex feasible space (Wainright and Jordan, 2008), and we can only apply the convex analysis theory to deal with convex sub-problems within an EM-type procedure. More concrete examples will be provided later along the developments of various models. We should also note that the modeling flexibility of RegBayes comes with risks. For example, it might lead to inconsistent posteriors (Barron et al., 1999;Choi and Ramamoorthi, 2008). This paper focuses on presenting several practical instances of RegBayes and we leave a systematic analysis of the Bayesian asymptotic properties (e.g., posterior consistency and convergence rates) for future work. Now, we derive the conjugate functions of three examples which will be used shortly for developing the infinite latent SVM models we have intended. We defer the proof to Appendix C. Specifically, the first one is the conjugate of a simple function, which will be used in a binary latent SVM classification model.
Lemma 9 Let g 0 : R → R be defined as g 0 (x) = C max(0, x). Then, we have The second function is slightly more complex, which will be used for defining a multi-way latent SVM classifier. Specifically, we define the function g 1 : R L → R as where max(x) def = max(x 1 , · · · , x L ). Apparently, g 1 is convex because it is a point-wise maximum (Boyd and Vandenberghe, 2004) of the simple linear functions φ i (x) = x i . Then, we have the following results.
Lemma 10 The convex conjugate of g 1 (x) as defined above is Let y ∈ R and ǫ ∈ R + are fixed parameters. The last function that we are interested in is g 2 : R → R, where Finally, we have the following lemma, which will be used in developing large-margin regression models.

Infinite Latent Support Vector Machines
Given the general theoretical framework of RegBayes introduced in Section 3, now we are ready to present its application to the development of two interesting nonparametric RegBayes models. In these two models we conjoin the ideas behind the nonparametric Bayesian infinite feature model known as the Indian buffet process (IBP), and the large margin classifier known as support vector machines (SVM) to build a new class of models for simultaneous single-task (or multi-task) classification and feature learning. A parametric Bayesian model is presented in Appendix B. Specifically, to illustrate how to develop latent large-margin classifiers and automatically resolve the unknown dimensionality of latent features from data, we demonstrate how to choose/define the three key elements of RegBayes, that is, prior distribution, likelihood model, and posterior regularization. We first present the single-task classification model. The basic setup is that we project each data example x ∈ X ⊂ R D to a latent feature vector z. Here, we consider binary features. Real-valued features can be easily considered by elementwisely multiplying z by a Guassian vector (Griffiths and Ghahramani, 2005). Given a set of N data examples, let Z be the matrix, of which each row is a binary vector z n associated with data sample n. Instead of pre-specifying a fixed dimension of z, we resort to the nonparametric Bayesian methods and let z have an infinite number of dimensions. To make the expected number of active latent features finite, we employ an IBP as prior for the binary feature matrix Z, as reviewed below.

Indian Buffet Process
Indian buffet process (IBP) was proposed in Griffiths and Ghahramani (2005) and has been successfully applied in various fields, such as link prediction (Miller et al., 2009) and multi-task learning (Rai and Daume III, 2010). We will make use of its stick-breaking construction (Teh et al., 2007), which is good for developing efficient inference methods. Let π k ∈ (0, 1) be a parameter associated with each column of the binary matrix Z. Given π k , each z nk in column k is sampled independently from Bernoulli(π k ). The parameter π are generated by a stick-breaking process where ν i ∼ Beta(α, 1). Since each ν i is less than 1, this process generates a decreasing sequence of π k . Specifically, given a finite dataset, the probability of seeing feature k decreases exponentially with k. IBP has several properties. For a finite number of rows, N , the prior of the IBP gives zero mass on matrices with an infinite number of ones, as the total number of columns with non-zero entries is Poisson(αH N ), where H N is the N th harmonic number, H N = N j=1 1 j . Thus, Z has almost surely only a finite number of non-zero entries, though this number is unbounded. A second property of IBP is that the number of features possessed by each data point follows a Poisson(α) distribution. Therefore, the expected number of non-zero entries in Z is N α.

Infinite Latent Support Vector Machines
Consider a single-task, but multi-way classification, where each training data is provided with a categorical label y ∈ Y def = {1, · · · , L}. Suppose that the latent features z n for document n are given, then we can define the latent discriminant function as linear where g(y, x n , z n ) is a vector stacking L subvectors 8 of which the yth is z ⊤ n and all the others are zero; η is the corresponding infinite-dimensional vector of feature weights. Since we are doing Bayesian inference, we need to maintain the entire distribution profile of the latent feature matrix Z. However, in order to make a prediction on the observed data x, we need to remove the uncertainty of Z. Here, we define the effective discriminant function as an expectation 9 (i.e., a weighted average considering all possible values of Z) of the latent discriminant function. To fully explore the flexibility offered by Bayesian inference, we also treat η as random and aim to infer its posterior distribution from given data. For the prior, we assume all the dimensions of η are independent and each dimension η k follows the standard normal distribution. This is in fact a Gaussian process (GP) prior as η is infinite dimensional. More formally, the effective discriminant function f : where q(Z, η, W) is the post-data posterior distribution we want to infer. We have included W as a place holder for any other variables we may define, e.g., the variables arising from a data likelihood model. Since we are taking the expectation, the variables which do not appear in the feature map g (i.e., W) will be marginalized out.
Before moving on, we should note that since we require q to be absolutely continuous with respect to the prior to make the KL-divergence term well defined in the RegBayes 8. We can consider the input features xn or its certain statistics in combination with the latent features zn to define a classifier boundary, by simply concatenating them in the subvectors. 9. Although other choices such as taking the mode are possible, our choice could lead to a computationally easy problem because expectation is a linear functional of the distribution under which the expectation is taken. Moreover, expectation can be more robust than taking the mode (Khan et al., 2010), and it has been widely used in (Zhu et al., , 2011b. problem, q(Z) will also put zero mass on Z's with an infinite number of non-zero entries, because of the properties of the IBP prior. The sparsity of Z is essential to ensure that the dot-product in Eq. (19) and the expectation in Eq. (20) are well defined, i.e., with finite values 10 . Moreover, in practice, to make the problem computationally feasible, we usually set a finite upper bound K to the number of possible features, where K is sufficiently large and known as the truncation level (See Section 4.4 and Appendix D.2 for details). As shown in (Doshi-Velez, 2009), the ℓ 1 -distance truncation error of marginal distributions decreases exponentially as K increases. For a finite truncation level, all the expectations are definitely finite.
Let I tr denote the set of training data. Then, with the above definitions, we define the P post (ξ) in problem (6) using soft 11 large-margin constraints as ) is the margin favored by the true label y n over an arbitrary label y and the superscript is used to distinguish from the posterior constraints for multi-task iLSVM to be presented. We define the penalty function for classification as where κ ≥ 1. If κ is 1, minimizing U c (ξ) is equivalent to minimizing the hinge-loss (or ℓ 1 -loss) R c h of the averaging prediction rule (27), where R c h (q(Z, η, W)) = C n∈Itr max y ℓ ∆ n (y) − ∆f (y n , x n ; q(Z, η, W)) ; if κ is 2, the surrogate loss is the squared ℓ 2 -loss. For clarity, we consider the hinge loss. The non-negative cost function ℓ ∆ n (y) (e.g., 0/1-cost) measures the cost of predicting x n to be y when its true label is y n . I tr is the index set of training data. Besides performing the prediction task, we may also be interested in explaining observed data x using the latent factors Z. This can be done by defining a likelihood model p(x|Z). Here, we define the most common linear-Gaussian likelihood model for real-valued data where W is a D × ∞ random loading matrix. We assume W follows an independent Gaussian prior and each entry has the prior distribution π(w dk ) = N (w dk |0, σ 2 0 ). The hyperparameters σ 2 0 and σ 2 n0 can be set a priori or estimated from observed data (See Appendix D.2 for details). Figure 2 (a) shows the graphical structure of iLSVM as defined above, where the plate means N replicates. 10. A more rigorous derivation of finiteness of these quantities is beyond the scope of this work and could require additional technical conditions (Orbanz, 2012). We refer the readers to (Stummer and Vajda, 2012) for a generic definition of Bregman divergence (or KL divergence in particular) on Banach spaces and in the case where the second measure is unnormalized. 11. Hard constraints for the separable cases are covered by simply setting ξ = 0.
Training: Putting the above definitions together, we get the RegBayes problem for iLSVM in the following two equivalent forms inf q(Z,η,W),ξ where p(Z, η, W, D) = π(η)π(Z)π(W) N n=1 p(x n |z n , W, σ 2 n0 ) is the joint distribution of the model; π(Z) is an IBP prior; and π(η) and π(W) are Gaussian process priors with identity covariance functions.
Directly solving the iLSVM problems is not easy because either the posterior constraints or the non-smooth regularization function R c is hard to deal with. Thus, we resort to convex duality theory, which will be useful for developing approximate inference algorithms. We can either solve the constrained form (E.q. (22)) using Lagrangian duality theory (Ito and Kunisch, 2008) or solve the unconstrained form (E.q. (23)) using Fenchel duality theory. Here, we take the second approach. In this case, the linear operator is the expectation operator, denoted by E : P prob → R |Itr|×L and the element of Eq evaluated at y for the nth example is where ∆g n (y, Z) def = g(y n , x n , z) − g(y, x n , z). Then, let g 1 : R L → R be a function defined in the same form as in Eq. (16). We have where Eq(n) def = (Eq(n, 1), · · · , Eq(n, L)) and ℓ ∆ n def = (ℓ ∆ n (1), · · · , ℓ ∆ n (L)) are the vectors of elements evaluated for nth data. By the Fenchel's duality theorem and the results in Lemma 10, we can derive the conjugate of the problem (23). The proof is deferred to Appendix C.4.
Lemma 12 (Conjugate of iLSVM) For the iLSVM problem, we have that where ω n = (ω 1 n , · · · , ω L n ) is the subvector associated with data n. Moreover, The optimum distribution is the posterior distribution q(Z, η, W) = 1 Z(ω|D) p(Z, η, W, D) exp n∈Itr yω y n η ⊤ ∆g n (y, Z) , where Z(ω|D) is the normalization factor andω is the solution of the dual problem. Testing: to make prediction on test examples, we put both training and test data together to do regularized Bayesian inference. For training data, we impose the above large-margin constraints because of the awareness of their true labels, while for test data, we do the inference without the large-margin constraints since we do not know their true labels. Therefore, the classifier (i.e., q(η)) is learned from the training data only, while both training and testing data influence the posterior distributions of the likelihood model W. After inference, we make the prediction via the rule Note that the ability to generalize to test data relies on the fact that all the data examples share η and the IBP prior. We can also cast the problem as a transductive inference problem by imposing additional large-margin constraints on test data (Joachims, 1999). However, the resulting problem will be generally harder to solve because it needs to resolve the unknown labels of testing examples. We also note that the testing is different from the standard inductive setting (Zhu et al., 2011b), where the latent features of a new data example can be approximately inferred given the training data. Our empirical study shows little difference on performance between our setting and the standard inductive setting.

Multi-Task Infinite Latent Support Vector Machines
Different from classification, which is typically formulated as a single learning task, multitask learning aims to improve a set of related tasks through sharing statistical strength among these tasks, which are performed jointly. Many different approaches have been developed for multi-task learning (See (Jebara, 2011) for a review). In particular, learning a common latent representation shared by all the related tasks has proven to be an effective way to capture task relationships (Ando and Zhang, 2005;Argyriou et al., 2007;Rai and Daume III, 2010). Below, we present the multi-task infinite latent SVM (MT-iLSVM) for learning a common binary projection matrix Z to capture the relationships among multiple tasks. Similar as in iLSVM, we also put the IBP prior on Z to allow it to have an unbounded number of columns. Suppose we have M related tasks. Let D m = {(x mn , y mn )} n∈I m tr be the training data for task m. We consider binary classification tasks, where Y m = {+1, −1}. Extension to multi-way classification or regression can be easily done. A naïve way to solve this learning problem with multiple tasks is to perform the multiple tasks independently. In order to make the multiple tasks coupled and share statistical strength, MT-iLSVM introduces a latent projection matrix Z. If the latent matrix Z is given, we define the latent discriminant function for task m as where x mn is one data example in D m and η m is the vector of parameters for task m.
The dimension of η m is the number of columns of the latent projection matrix Z, which is unbounded in the nonparametric setting. This definition provides two views of how the M tasks get related.
(1) If we let ς m = Zη m , then ς m is the actual parameter of task m and all ς m in different tasks are coupled by sharing the same latent matrix Z; (2) Another view is that each task m has its own parameters η m , but all the tasks share the same latent projection matrix Z to extract latent features Z ⊤ x mn , which is a projection of the input features x mn .
As such, our method can be viewed as a nonparametric Bayesian treatment of alternating structure optimization (ASO) (Ando and Zhang, 2005), which learns a single projection matrix with a pre-specified latent dimension. Moreover, different from (Jebara, 2011), which learns a binary vector with known dimensionality to select features or kernels on x, we learn an unbounded projection matrix Z using nonparametric Bayesian techniques. As in iLSVM, we employ a Bayesian treatment of η m , and view it as random variables. We assume that η m has a fully-factorized Gaussian prior, i.e., η mk ∼ N (0, 1). Then, we define the effective discriminant function for task m as the expectation where W is a place holder for the variables that possibly arise from other parts of the model. As in iLSVM, since we are taking expectation, the variables which do not appear in the feature map (i.e., W) will be marginalized out. Then, the prediction rule for task m is naturally y * m def = signf m (x). Similarly, we perform regularized Bayesian inference by defining: and imposing the following constraints: Finally, as in iLSVM we may also be interested in explaining observed data x. Therefore, we relate Z to the observed data x by defining a likelihood model: p x mn |w mn , Z, λ 2 mn = N x mn |Zw mn , λ 2 mn I , where w mn is a vector. We assume W has an independent prior π(W) = mn N (w mn |0, σ 2 m0 I). Fig. 2 (b) illustrates the graphical structure of MT-iLSVM.
For training, we can derive the similar convex conjugate as in the case of iLSVM. Similar as in iLSVM, minimizing U M T (ξ) is equivalent to minimizing the hinge-loss R M T h of the multiple binary prediction rules, where Thus, the RegBayes problem of MT-iLSVM can be equivalently written as Then, by the Fenchel's duality theorem and Lemma 9, we can derive the conjugate of MT-iLSVM. The proof is deferred to Appendix C.5.

Lemma 13 (Conjugate of MT-iLSVM) For the MT-iLSVM problem, we have that
Moreover, The optimum distribution is the posterior distribution where Z ′ (ω|D) is the normalization factor andω is the solution of the dual problem.
For testing, we use the same strategy as in iLSVM to do Bayesian inference on both training and test data. The difference is that training data are subject to large-margin constraints, while test data are not. Similarly, the hyper-parameters σ 2 m0 and λ 2 mn can be set a priori or estimated from data (See Appendix D.1 for details).

Inference with Truncated Mean-Field Constraints
Now we discuss how to perform regularized Bayesian inference with the large-margin constraints for both iLSVM and MT-iLSVM. From the primal-dual formulations, it is obvious that there are basically two methods to perform the regularized Bayesian inference. One is to directly solve the primal problem for the posterior distribution q(Z, η, W), and the other is to first solve the dual problem for the optimumω and then infer the posterior distribution. However, both the primal and dual problems are intractable for iLSVM and infer q(ν), q(W) and q(Z) with q(η) and ω given;

5:
infer q(η) and solve for ω with q(Z) given. 6: until convergence MT-iLSVM. The intrinsic hardness is due to the mutual dependency among the latent variables in the desired posterior distribution. Therefore, a natural approximation method is the mean field (Jordan et al., 1999), which breaks the mutual dependency by assuming that q is of some factorization form. This method approximates the original problems by imposing additional constraints. An alternative method is to apply approximate methods (e.g., MCMC sampling) to infer the true posterior distributions derived via convex conjugates as above, and iteratively estimate the dual parameters using approximate statistics (e.g., feature expectations estimated using samples) (Schofield, 2006). Below, we use MT-iLSVM as an example to illustrate the idea of the first strategy. A full discussion on the second strategy is beyond the scope of this paper. For iLSVM, the similar procedure applies and we defer its details to Appendix D.2.
To make the problem easier to solve, we use the stick-breaking representation of IBP, which includes the auxiliary variable ν, and infer the augmented posterior q(ν, W, Z, η). The joint model distribution is now q(ν, W, Z, η, D). Furthermore, we impose the truncated mean-field constraint that where K is the truncation level, and we assume that q(ν k |γ k ) = Beta(γ k1 , γ k2 ), q(z dk |ψ dk ) = Bernoulli(ψ dk ), q(w mn |Φ mn , σ 2 mn I) = N (w mn |Φ mn , σ 2 mn I). Then, we can use the duality theory 12 to solve the RegBayes problem by alternating between two substeps, as outlined in Algorithm 1 and detailed below.
Infer q(ν), q(W) and q(Z): Since q(ν) and q(W) are not directly involved in the posterior constraints, we can solve for them by using standard Bayesian inference, i.e., minimizing a KL-divergence. Specifically, for q(W), since the prior is also normal, we can easily derive the update rules for Φ mn and σ 2 mn . For q(ν), we have the same update rules as in (Doshi-Velez, 2009). We defer the details to Appendix D.1.
For q(Z), it is directly involved in the posterior constraints. So, we need to solve it together with q(η) using conjugate theory. However, this is intractable. Here, we adopt an alternating strategy that first infers q(Z) with q(η) and dual parameters ω fixed, and then infers q(η) and solves for ω. Specifically, since the large-margin constraints are linear of q(Z), we can get the mean-field update equation as and L ν k is an lower bound of E q [log(1 − k j=1 v j )] (See Appendix D.1 for details). The last term of ϑ dk is due to the large-margin posterior constraints as defined in Eq. (30). Therefore, from this equation we can see how the large-margin constraints regularize the procedure of inferring the latent matrix Z.
Infer q(η) and solve for ω: Now, we can apply the convex conjugate theory and show that the optimum posterior distribution of η is and µ m = n∈I m tr y mn ω mn (ψ ⊤ x mn ). Here, we assume π(η m ) is standard normal. Then, we have q(η m ) = N (η m |µ m , I) and the optimum dual parameters can be obtained by solving the following M independent dual problems ∀n ∈ I m tr , s.t. : 0 ≤ ω mn ≤ C, where the constraints are from the conjugate function g * 0 in Lemma 13. These dual problems (or their primal forms) can be efficiently solved with a binary SVM solver, such as SVM-light or LibSVM.

Experiments
We present empirical results for both classification and multi-task learning. Our results appear to demonstrate the merits inherited from both Bayesian nonparametrics and largemargin learning.

Multi-way Classification
We evaluate the infinite latent SVM (iLSVM) for classification on the real TRECVID2003 and Flickr image datasets, which have been extensively evaluated in the context of learning finite latent feature models (Chen et al., 2010). TRECVID2003 consists of 1078 video   key-frames that belong to 5 categories, including Airplane scene, Basketball scene, Weather news, Baseball scene, and Hockey scene. Each data example has two types of features -1894-dimension binary vector of text features and 165-dimension HSV color histogram. The Flickr image dataset consists of 3411 natural scene images about 13 types of animals, including squirrel, cow, cat, zebra, tiger, lion, elephant, whales, rabbit, snake, antlers, hawk and wolf, downloaded from the Flickr website 13 . Also, each example has two types of features, including 500-dimension SIFT bag-of-words and 634-dimension real-valued features (e.g., color histogram, edge direction histogram, and block-wise color moments). Here, we consider the real-valued features only by defining Gaussian likelihood distributions for x; and we define the discriminant function using latent features only as in Eq. (19). We follow the same training/testing splits as in (Chen et al., 2010). We compare iLSVM with the large-margin Harmonium (MMH) (Chen et al., 2010), which was shown to outperform many other latent feature models, and two decoupled approaches -EFH+SVM and IBP+SVM. EFH+SVM uses the exponential family Harmonium (EFH) (Welling et al., 2004) to discover latent features and then learns a multi-way SVM classifier. IBP+SVM is similar, but uses an IBP factor analysis model (Griffiths and Ghahramani, 2005) to discover latent features. To initialize the learning algorithms for these models, we found that using the SVD factors of the input feature matrix as the initial weights for MMH and EFH can produce better results. Here, we also use the SVD factors as the initial mean of weights in the likelihood models for iLSVM. Both MMH and EFH+SVM are finite models and they need to pre-specify the dimensionality of latent features. We report their results on classification accuracy and F1 score (i.e., the average F1 score over all possible classes) (Zhu et al., 2011b) achieved with the best dimensionality in Table 1. Figure 3 illustrates the performance change of MMH when using different number of latent features, from which we can see that K = 40 produces the best performance and either increasing or decreasing K could make the performance worse. For iLSVM and IBP+SVM, we use the mean-field inference method and present the average performance with 5 randomly ini- tialized runs (Please see Appendix D.2 for the algorithm and initialization details). We perform 5-fold cross-validation on training data to select hyperparameters, e.g., α and C (we use the same procedure for MT-iLSVM). We can see that iLSVM can achieve comparable performance with the nearly optimal MMH, without needing to pre-specify the latent feature dimension 14 , and is much better than the decoupled approaches (i.e., IBP+SVM and EFH+SVM). For the two stage methods, we don't have a clear winner -IBP+SVM performs a bit worse than EFH+SVM on the TRECVID dataset, while it outperforms EFH+SVM on the flickr dataset. The reason for the difference may be due to the initialization or different properties of the data. It is also interesting to examine the discovered latent features. Figure 4 shows the overall average values of latent features and the per-class average feature values of iLSVM in one run on the TRECVID dataset. We can see that on average only about 45 features are active for the TRECVID dataset. For the overall average, we also present the standard deviation over the 5 categories. A larger deviation means that the corresponding feature is more discriminative when predicting different categories. For example, feature 26 and feature 34 are generally less discriminative than many other features, such as feature 1 and feature 30. Figure 5 shows the overall average feature values together with standard deviation on the Flickr dataset. We omitted the per-class average because that figure is too crowded with 13 categories. We can that as k increases, the probability that feature k is active decreases. The reason for the features with stable values (i.e., standard deviations are extremely small) is due to our initialization strategy (each feature has 0.5 probability to be active). Initializing ψ dk as being exponentially decreasing (e.g., like the constructing process of π) leads to a faster decay and many features will be inactive. To examine the semantics 15 of each feature, Figure 6 presents some example features discovered on the Flickr animal dataset. For each feature, we present 5 top-ranked images which have large values on this particular feature. We can see that most of the features are semantically interpretable. For instance, feature F1 is about squirrel; feature F2 is about ocean animal, which is whales in the Flickr dataset; and feature F4 is about hawk. We can also see that some features are about different aspects of the same category. For example, feature F2 and feature F3 are both about whales, but with different background.

Multi-task Learning
Now, we evaluate the multi-task infinite latent SVM (MT-iLSVM) on several well-studied real datasets.

Description of the Data
Scene and Yeast Data: These datasets are from the UCI repository, and each data example has multiple labels. As in (Rai and Daume III, 2010), we treat the multi-label classification as a multi-task learning problem, where each label assignment is treated as a binary classification task. The Yeast dataset consists of 1500 training and 917 test examples, each having 103 features, and the number of labels (or tasks) per example is 14. The Scene dataset consists 1211 training and 1196 test examples, each having 294 features, and the number of labels (or tasks) per example for this dataset is 6.
School Data: This dataset comes from the Inner London Education Authority and has been used to study the effectiveness of schools. It consists of examination records of 15,362 students from 139 secondary schools in years 1985, 1986 and 1987. The dataset is publicly available and has been extensively evaluated in various multi-task learning methods (Bakker and Heskes, 2003;Bonilla et al., 2008;Zhang and Yeung, 2010), where each task is defined as predicting the exam scores of students belonging to a specific school based on four student-dependent features (year of the exam, gender, VR band and ethnic group) and four school-dependent features (percentage of students eligible for free school meals, percentage of students in VR band 1, school gender and school denomination). In order to compare with the above methods, we follow the same setup described in (Argyriou et al., 2007;Bakker and Heskes, 2003) and similarly we create dummy variables for those features that are categorical forming a total of 19 student-dependent features and 8 school-dependent features. We use the same 10 random splits 16 of the data, so that 75% of the examples from each school (task) belong to the training set and 25% to the test set. On average, the training set includes about 80 students per school and the test set about 30 students per school.  Table 2: Multi-label classification performance on Scene and Yeast datasets.

Results
Scene and Yeast Data: We compare with the closely related nonparametric Bayesian methods, including kernel stick-breaking (YaXue) (Xue et al., 2007) and the basic and augmented infinite predictor subspace models (i.e., Piyushrai-1 and Piyushrai-2) (Rai and Daume III, 2010). These nonparametric Bayesian models were shown to outperform the independent Bayesian logistic regression and a single-task pooling approach (Rai and Daume III, 2010). We also compare with a decoupled method MT-IBP+SVM 17 that uses an IBP factor analysis model to find shared latent features among multiple tasks and then builds separate SVM classifiers for different tasks. For MT-iLSVM and MT-IBP+SVM, we use the mean-field inference method in Sec 4.4 and report the average performance with 5 randomly initialized runs (See Appendix D.1 for initialization details). For comparison with (Rai and Daume III, 2010;Xue et al., 2007), we use the overall classification accuracy, F1-Macro and F1-Micro as performance measures. Table 2 shows the results. On both datasets, MT-iLSVM needs less than 50 latent features on average. We can see that the large-margin MT-iLSVM performs much better than other nonparametric Bayesian methods and MT-IBP+SVM, which separates the inference of latent features from learning the classifiers. School Data: We use the percentage of explained variance (Bakker and Heskes, 2003) as the measure of the regression performance, which is defined as the total variance of the data minus the sum-squared error on the test set as a percentage of the total variance. Since we use the same settings, we can compare with the state-of-the-art results of (1) Bayesian multi-task learning (BMTL) (Bakker and Heskes, 2003); (2) Multi-task Gaussian processes (MTGP) (Bonilla et al., 2008); (3) Convex multi-task relationship learning (MTRL) (Zhang and Yeung, 2010); and single-task learning (STL) as reported in (Bonilla et al., 2008;Zhang and Yeung, 2010). For MT-iLSVM and MT-IBP+SVM, we also report the results achieved by using both the latent features (i.e., Z ⊤ x) and the original input features x through vector concatenation, and we denote the corresponding methods by MT-iLSVM f and MT-IBP+SVM f , respectively. On average the multi-task latent SVM (i.e., MT-iLSVM) needs about 50 latent features to get sufficiently good and robust performance. From the results in Figure 7, we can see that the MT-iLSVM achieves better results than the existing methods that have been tested in previous studies. Again, the joint MT-iLSVM performs much better than the decoupled method MT-IBP+SVM, which separates the latent feature inference from the training of large-margin classifiers. Finally, using both latent features and the original input features can boost the performance slightly for MT-iLSVM, while much more significantly for the decoupled MT-IBP+SVM. Figure 8 shows how the performance of MT-iLSVM changes against the hyper-parameter α and regularization constant C on the Yeast and School datasets. We can see that on the Yeast dataset, MT-iLSVM is insensitive to both α and C. For the School dataset, MT-iLSVM is very insensitive the α, and it is stable when C is set between 0.3 and 1. Figure 9 shows how the training size affects the performance and running time of MT-iLSVM on the School dataset. We use the first b% (b = 50, 60, 70, 80, 90, 100) of the training data in each of the 10 random splits as training set and use the corresponding test data as test set. We can see that as training size increases, the performance and running time generally increase; and MT-iLSVM achieves the state-of-art performance when using about 70% training data. From the running time, we can also see that MT-iLSVM is generally quite efficient by using mean-field inference.

Sensitivity Analysis
Finally, we investigate how the performance of MT-iLSVM changes against the hyperparameters σ 2 m0 and λ 2 mn . We initially set σ 2 m0 = 1 and compute λ 2 mn from observed data. If we further estimate them by maximizing the objective function, the performance does not change much (±0.3% for average explained variance on the School dataset). We have similar observations for iLSVM.

Conclusions and Discussions
We present regularized Bayesian inference (RegBayes), a computational framework to perform post-data posterior inference with a rich set of regularization/constraints on the desired post-data posterior distributions. RegBayes is formulated as a information-theoretical optimization problem, and it is applicable to both directed and undirected graphical models. We present a general theorem to characterize the solution of RegBayes, when the posterior regularization is induced from a linear operator (e.g., expectation). Furthermore, we particularly concentrate on developing two large-margin nonparametric Bayesian models under the RegBayes framework to learn predictive latent features for classification and multi-task learning, by exploring the large-margin principle to define posterior constraints. Both models allow the latent dimension to be automatically resolved from the data. The empirical results on several real datasets appear to demonstrate that our methods inherit the merits from both Bayesian nonparametrics and large-margin learning.
RegBayes offers a flexible framework for considering posterior regularization in performing parametric or nonparametric Bayesian inference. For future work, we plan to study other posterior regularization beyond the large-margin constraints, such as posterior constraints defined on manifold structures (Huh and Fienberg, 2010) and those represented in the form of first-order logic, and investigate how posterior regularization can be used in other interesting nonparametric Bayesian models (Beal et al., 2002;Teh et al., 2006;Blei and Frazier, 2010) in different contexts, such as link prediction (Miller et al., 2009) for social network analysis and low-rank matrix factorization for collaborative prediction. Some of our preliminary results (Xu et al., 2012;Zhu, 2012;Mei et al., 2014) have shown great promise. It is interesting to investigate more carefully along this direction. Moreover, as we have stated, RegBayes can be developed for undirected MRFs. But the inference would be even harder. We plan to do a systematic investigation along this direction too. We have some preliminary results presented in , but there is a lot of room to further improve. Finally, regularized Bayesian inference in general leads to a highly nontrivial inference problem. Although the general solution can be derived with convex analysis theory, it is normally intractable to infer them directly. Therefore, approximate inference techniques such as the truncated mean-field approximation have to be used. For the current truncated inference methods, one key limit is to pre-specify the truncation level. A too conservative truncation level could lead to a waste of computing resources. So, it is important to develop inference algorithms that could adaptively determine the number of latent features, such as Monte Carlo methods. We have some preliminary progress along this direction as reported in the work (Jiang et al., 2012;. It is interesting to extend these techniques to deal with other challenging nonparametric Bayesian models. Proof According to the variational formulation of Bayes' rule in Eq. (5), we get that the optimum solution is q(M) = p(M|D, Θ) for any Θ. Substituting the optimum solution of q into the objective, we get the optimization problem of Θ.
Extension 2: Chain Graph: In the above cases, we have assumed that the observed data are generated by some model in a directed causal sense. This assumption holds in directed latent variable models. However, in many cases, we may choose alternative formulations to define the joint distribution of a model and the observed data. Figure 10(c) illustrates one such scenario, where the model M consists of two subsets of random variables. One subset H is connected to the observed data via an undirected graph and the other subset Z is connected to the observed data and H using directed edges. This graph is known as a chain graph. Due to the Markov properties of chain graph (Frydenberg, 1990), we know that the joint distribution has the factorization form as where p(H, D|Z) is a Markov random field (MRF). One concrete example of such a hybrid chain model is the Bayesian Boltzman machines (Murray and Ghahramani, 2004), which treat the parameters of a Boltzmann machine as random variables and perform Bayesian inference with MCMC sampling methods. The insights that RegBayes covers undirected or chain graph latent variable models come from the observation that the objective L(q(M)) of problem (5)  where L(Θ, q(M)) is the objective function of problem (39). These two formulations are equivalent. We will call the former a constrained formulation and call the latter an unconstrained formulation by ignoring the standard normalization constraints, which are easy to deal with.

Appendix B: MedLDA-A RegBayes Model with Finite Latent Features
This section presents a new interpretation of MedLDA (maximum entropy discrimination latent Dirichlet allocation) ) under the framework of regularized Bayesian inference. MedLDA is a max-margin supervised topic model, an extension of latent Dirichlet allocation (LDA) (Blei et al., 2003) for supervised learning tasks. In MedLDA, each data example is projected to a point in a finite dimensional latent space, of which each feature corresponds to a topic, i.e., a unigram distribution over the terms in a vocabulary. MedLDA represents each data as a probability distribution over the features, which results in a conservation constraint (i.e., the more a data expresses on one feature, the less it can express others) (Griffiths and Ghahramani, 2005). The infinite latent feature models discussed in Section 4 do not have such a constraint. Without loss of generality, we consider the MedLDA regression model as an example (classification model is similar), whose graphical structure is shown in Figure 11. We assume that all data examples have the same length V for notation simplicity. Each document is associated with a response variable Y , which is observed in the training phase but unobserved in testing. We will use y to denote an instance value of Y . Let K be the number of topics or the dimensionality of the latent topic space. MedLDA builds an LDA model to describe the observed words. The generating process of LDA is that each document n has a mixing proportion θ n ∼ Dirichlet(α); each word w nm is associated with a topic z nm ∼ θ n , which indexes the topic that generates the word, i.e., w nm ∼ β znm . DefineZ n = 1 s.t. ∀n : The posterior constraints are imposed following the large-margin principle and they correspond to a quality measure of the prediction results on training data. In fact, it is easy to show that minimizing U (ξ, ξ * ) = C N n=1 (ξ n + ξ * n ) under the above constraints is equivalent to minimizing an ǫ-insensitive loss (Smola and Schölkopf, 2003) of the expected linear prediction ruleŷ n = E p [η ⊤Z n ]. To practically learn an MedLDA model, since the above problem is intractable, variational methods were used by introducing an auxiliary distribution q({θ n , z nm , η}|Θ) 20 to approximate the true posterior p({θ n , z nm , η}|D, Θ), replacing the negative data likelihood with its upper bound L q({θ n , z nm , η}|Θ) , and replacing p by q in the constraints. The variational MedLDA regression model is s.t. ∀n : where L q({θ n , z nm , η}|Θ) = −E q log p({θ n , z nm , η}, D|Θ) − H q({θ n , z nm , η}|Θ) is a variational upper-bound of the negative data log-likelihood. The upper bound is tight if no restricting constraints are made on the variational distribution q. In practice, additional assumptions (e.g., mean-field) can be made on q to derive a practical approximate algorithm.
20. We have explicitly written the condition on model parameters.
Based on the previous discussions on the extensions of RegBayes and the duality in Lemma 14, we can reformulate the MedLDA regression model as an example of RegBayes. Specifically, for the MedLDA regression model, we have M = {θ n , z nm , η}. According to Eq. (41), we can easily show that L q({θ n , z nm , η}|Θ) = KL q({θ n , z nm , η}|Θ) p({θ n , z nm , η}, {w nm , y n }|Θ) Then, the MedLDA problem is a RegBayes model in Eq. (42) with For the MedLDA problem, we can use Lagrangian methods to solve the constrained formulation. Alternatively, we can also use the convex duality theorem to solve the equivalent unconstrained form. For the variational MedLDA, the ǫ-insensitive loss is R ǫ (q({θ n , z nm , η}|Θ)). Its conjugate can be derived using the results of Lemma 11. Specifically, we have the following result, whose proof is deferred to Appendix C.6.
Lemma 15 (Conjugate of MedLDA) For the variational MedLDA problem, we have = sup ω − log Z ′ (ω, Θ * ) − n g * 2 (ω n ; −y n + ǫ, y n + ǫ), where ω n = (ω n , ω ′ n ). Moreover, The optimum distribution is the posterior distribution where Z ′ (ω, Θ|D) is the normalization factor and the optimum parameters are Note that although in general, either the primal or the dual problem is hard to solve exactly, the above conjugate results are still useful when developing approximate inference algorithms. For instance, we can impose additional mean-field assumptions on q in the primal formulation and iteratively solve for each factor; and in this process convex conjugates are useful to deal with the large-margin constraints ). Alternatively, we can apply approximate methods (e.g., MCMC sampling) to infer the q based on its solution in Eq. (48), and iteratively solves for the dual parameters ω using approximate statistics (Schofield, 2006). We will discuss more on this when presenting the inference algorithms for iLSVM and MT-iLSVM.
In the above discussions, we have treated the topics β as fixed unknown parameters. A fully Bayesian formulation would treat β as random variables, e.g., with a Dirichlet prior (Blei et al., 2003;Griffiths and Steyvers, 2004). Under the RegBayes interpretation, we can easily do such an extension of MedLDA, simply by moving β from Θ to M. Absorbing the constant c into f * , we get the dual objective of Theorem 6.

Appendix C.3: Proof of Lemma 10
Proof The proof has a similar structure as the proof of Lemma 9. By definition, we have We first show that ∀i, µ i ≥ 0 in order to have finite g * 1 values. Suppose that ∃j, µ j < 0. Then, we define Therefore, g * 1 (µ) = ∞ if ∃j, µ j < 0. Now, we consider the second case, where ∀i, µ i ≥ 0. We can easily show that Moreover, let G 1 = {x ∈ R L : x = xe, x ∈ R}, where e is a vector with every element being 1. Then, we have Putting the above results together proves the claim.

Appendix C.4: Proof of Lemma 11
Proof By definition, the conjugate is For the first infimum, we can show that Thus, we have where the second equality holds by setting α = |µ|, under the condition that ǫ is positive; the condition |µ| ≤ C is induced from the conditions α + β = C and β ≥ 0.
Proof Similar structure as the proof of Lemma 12. In this case, the linear expectation operator is E : P prob → R m |I m tr | and the element of Eq evaluated at the nth example for task m is Then, let g 0 : R → R be a function defined in Lemma 9. We have Let µ = Eq. By definition, the conjugate is Thus, g * (−ω) = m,n∈I m tr − ω mn + g * 0 (ω mn ) .
By the results in Theorem 6 and Lemma 9, we can derive the conjugate of the problem (33).
Appendix C.7: Proof of Lemma 15 Proof Similar structure as the proof of Lemma 12. In this case, the linear expectation operator is E : P prob → R N and the elements of Eq evaluated at the nth example is Then, using the g 2 function defined in Lemma 11, we have g(Eq) def = R ǫ (q({θ n , z nm , η}|Θ)) = n g 2 µ n ; y n , ǫ .
Therefore g * (ω) = n g * 2 (ω n ; y n , ǫ) and g * (−ω) = n g * 2 (−ω n ; y n , ǫ). By the results in Theorem 6 and Lemma 9, we can derive the conjugate and the optimum solution ofq. The optimum solution of Θ is due to Lemma 14. Note that the constraints are not directly dependent on Θ.

Appendix D: Inference Algorithms for Infinite Latent SVMs
Appendix D.1: Inference for MT-iLSVM In this section, we provide the derivation of the inference algorithm for MT-iLSVM, which is outlined in Algorithm 2 and detailed below.
For MT-iLSVM, the model M consists of all the latent variables (ν, W, Z, η). Let L mn (q) def = E q [log p(x mn |Z, w mn , λ 2 mn )] be the expected data likelihood. Then, under the truncated mean-field assumption (36), we have For the KL-divergence term, we have KL(q(M) π(M)) = KL(q(ν) π(ν))+KL(q(W) π(W))+ E q(ν) [KL(q(Z) π(Z|ν))] + KL(q(η) π(η)), where the individual terms are where ϕ(·) is the digamma function and E q [log v j ] = ϕ(γ j1 )−ϕ(γ j1 +γ j2 ). For KL(q(η) π(η)), we do not need to write it explicitly, as we shall see. Finally, the effective discriminant function is All the above terms can be easily computed, except the term E q [log(1 − k j=1 ν j )]. Here, we adopt the multivariate lower bound (Doshi-Velez, 2009)  where the variational parameters q k. = (q k1 · · · q kk ) ⊤ belong to the k-simplex, and H(q k. ) is the entropy of q k. . The tightest lower bound is achieved by setting q k. to be the optimum value where Z k is a normalization factor to make q k. be a distribution. We denote the tightest lower bound by L ν k . Replacing the term E q [log(1 − k j=1 ν j )] with its lower bound L ν k , we can have an upper bound of KL(q(M) π(M)) and we denote this upper bound by L(q).

12:
end for 13: update the hyper-parameters σ 2 m0 using Eq. (60) and λ 2 mn using Eq. (61). (Optional) 14: until relative change of L is less than τ ′ (e.g., 1e −4 ) or iteration number is T ′ (e.g., 20) Infer q(ν), q(Z) and q(W): For q(W), since both the prior π(W) and q(W) are Gaussian, we can easily derive the update rules, similar as in Gaussian mixture models For q(ν), we have the update rules similar as in (Doshi-Velez, 2009), that is, For q(Z), we have the mean-field update equation as Infer q(η) and solve for ω: By the convex duality theory, we have the solution q(η) ∝ π(η) exp m,n∈I m tr y mn ω mn η ⊤ m ψ ⊤ x mn = M m=1 π(η m ) exp η ⊤ m n∈I m tr y mn ω mn ψ ⊤ x mn .
Therefore, we can see that although we did not assume q(η) is factorized, we can get the induced factorization form q(η) = m q(η m ), where q(η m ) ∝ π(η m ) exp η ⊤ m n∈I m tr y mn ω mn ψ ⊤ x mn .
Here, we assume π(η m ) is standard normal. Then which (and its primal form) can be efficiently solved with a binary SVM solver, such as SVM-light. As we have stated, the hyperparameters σ 2 0 and λ 2 mn can be set a priori or estimated from the data. The empirical estimation can be easily done with closed form solutions by optimizing the RegBayes objective with all the variational terms fixed. For MT-iLSVM, we have Appendix D.2: Inference for Infinite Latent SVM In this section, we develop the inference algorithm for iLSVM based on the stick-breaking construction of the IBP prior. The algorithm is outlined in Algorithm 3. Similar as in the inference for MT-iLSVM, we make the additional constraint about the feasible distribution q(ν, W, Z, η) = q(η)q(W|Φ, Σ) n K k=1 q(z nk |ψ nk ) K k=1 q(ν k |γ k ), where K is the truncation level; q(W|Φ, Σ) = k N (W .k |Φ .k , σ 2 k I); q(z nk |φ nk ) = Bernoulli(φ nk ); and q(ν k |γ k ) = Beta(γ k1 , γ k2 ). Then, we solve the unconstrained problem using convex duality with dual parameters being ω. Let L n (q) def = E q [log p(x n |z n , W)]. We have where The effective discriminant function is f (y, x n ) = k E q [η k y ]ψ nk . Again, for computational tractability, we need the lower bound L ν k of the term E q [log(1 − k j=1 v j )]. Using this lower bound, we can get an upper bound of the KL-divergence term. Then, the inference procedure iteratively solves the following steps: Infer q(ν), q(Z) and q(W): For q(W), we have the update rules For q(ν), we have the update rules similar as in (Doshi-Velez, 2009), that is, where q .k is computed in the same way as in Eq. (55). For q(Z), the mean-field update equation for ψ is where For testing data, ϑ nk does not have the last term because of the absence of large-margin constraints.

11:
solve the dual problem (66) (or its primal form) using a multi-class SVM learner.
For the standard normal prior π(η), we have that q(η) is also normal, with mean µ = n∈Itr y ω y d E q [g(y n , x n , z n ) − g(y, x n , z n )] and identity covariance matrix. The dual problem is sup ω − 1 2 µ ⊤ µ + n∈Itr y ω y n s.t.. : ω y n ≥ 0, y ω y n = C, ∀n ∈ I tr , which (and its primal form) can be efficiently solved with a multi-class SVM solver. Similar as in MT-iLSVM, the hyperparameters σ 2 0 and σ 2 n0 can be set a priori or estimated from the data. The empirical estimation can be easily done with closed form solutions. For iLSVM, we have