Markov random fields: Parameter learning I

In this post, I will review some practical learning methods which are applicable for large-scale Markov random fields (MRFs), e.g., with millions of parameters. This has some overlaps with this post, but explains the ideas in more details.

Here I use the term "Markov random fields" to refer to a class of multivariate models whose (density) probability distributions involves a normalisation constant \( Z \) known as the partition function:
\[  P(\mathbf{x};\theta) = \frac{1}{Z(\theta)} \Psi(\mathbf{x})  \]
where \( \theta \) is model parameter.

This also includes conditional random fields (CRFs), Markov relational networks (MRN), and to some extent, Markov logic networks. The main difference with CRFs and other conditional models is that the partition function is input dependent, that is, we need to write \( Z(\mathbf(y),\theta) \), i.e., we need to compute the partition function for every data point. A similar situation happens when the graph of the MRF is data dependent, e.g., in our previous application of MRFs for collaborative filtering, we build a MRF for each user to exploit the sparsity in the data. On the other hand, in non-conditional MRF, we need only one partition function.

The estimation of \( Z(\theta) \) is general intractable except for the tree-structured MRFs, and indeed if we know how to estimate \( Z(\theta) \), we are likely to be able to infer most quantities (e.g., expectation of some features). For example, Wainwright has demonstrated that by minimising its upper bound, we can arrive at a quite stable message-passing algorithm, from which marginals can be estimated efficiently.

However, the approximate estimation of \( Z(\theta) \) can still be very expensive, especially in the setting of learning because it depends on the model parameters being estimated, i.e., for every parameter updates, we need to recompute \( Z(\theta) \). The trick is to avoid computing it explicitly.

The maximum likelihood

Denote my \(G = (V,E)\) the underlying graph where \( V \) is  a set of nodes, and \( E \) is a set of edges. For simplicity, we will assume the log-linear parameterisation for edges and nodes, that is, the model potential function will have the following energy-based form
\[
 \Psi(\mathbf{x}) = \exp( -  F(\mathbf{x}) )
\]
where the energy \( F(\mathbf{x}) \) is decomposed into parts
 \[
F(\mathbf{x}) = -\left( \sum_{i \in V} \sum_k\theta_{ik} f_{ik}(x_i)  + \sum_{ij \in E} \sum_m\theta_{ijm} g_{ijm}(x_i,x_j)\right)
\]
where \(f_{ik}(x_i) \) and \( g_{ijm}(x_i,x_j) \) are feature functions.

Note that this model is somewhat overparameterised, but this form will be of much help when we have to deal with rich input features in the conditional setting (aka conditional random fields).

For now let us assume that the training data is fully observed. We will cover the partially observed variables later, e.g., for the case of restricted Boltzmann machines (RBM).

Maximum likelihood estimation starts with the data log-likelihood
\[
L(\theta) = \frac{1}{D}\sum_{l=1}^D \log P(\mathbf{x}^{(l)};\theta)
\]
where \( D \) is the number of observations. This criterion is indeed statistically powerful as it leads to consistent and fast estimate in the asymptotic sense.

This can be rewritten as
\[
L(\theta) =  -  \frac{1}{D}\sum_{l=1}^D E(\mathbf{x}^{(l)}) - \log Z(\theta)
\]

Note that the partition function does not depend on the the observation in this setting, while it depends on the input in the conditional setting.

The first part of the right-hand-side is easier to compute, but the partition function is still there. The nice property of the log-partition function is that it is convex in \( \theta \). Thus if we can find a way to maximise the log-likelihood, we are guaranteed to obtain the global maximum.

In a typical optimisation, we would need the gradient of the log-likelihood. For example, for pairwise parameter:
\[
 \frac{\partial L}{\partial \theta_{ijm}} = \frac{1}{D}\sum_{l=1}^D g_{ijm}(x_i^{(l)},x_j^{(l)})-
\sum_{x_i,x_j} P( x_i,x_j) g_{ijm}(x_i,x_j))
\]
In the literature, especially those associated with the Boltzmann machines, the computation of the first part of the right-hand-side is often referred to as the clamped phase and the second part as the free phase.  This is because in the former, variables are clamped to their specific observations, while in the latter, variables are free to vary according to the distribution.

Another view is that the first part is the expectation of feature function with respect to the empirical distribution, while the second part is with respect to the model distribution. The maximum likelihood solution is equivalent to solving the gradient equation, e.g., by matching the two expectations (or moments). Thus, sometimes the term moment-matching is used.

We now have three options for  learning.
1. Approximate objective function
2. Approximate Z
3. Approximate gradients

As for optimisation, we may not need a full evaluation of the log-likelihood, and thus the partition function is not needed. The idea here is to go with the gradient only, and hopefully the learning still progresses as we expect.

In the subsequent posts, we will cover some of the following methods:

• Pseudo-likelihood
• Structured pseudo methods
• Approximate gradients with truncated MCMCs, aka (Persistent) Contrastive Divergence
• Piecewise methods
• Tree-based methods
• Herding
• Label relaxation
• Structured mean-fields
• Approximate perceptron
• Wainwright moment matching
• Conditional versus non-conditional
• Mean-field approximation as Recurrent Neural Networks.

Markov random fields: Structure learning I

Structure learning in Markov random fields is an important topic since without the correct structures, we can potentially estimate a wrong or unreliable model. Often structures are specified by domain knowledge such as the 2D grid in modelling images or the 1D chain in word labelling. However, there are cases which structures are not available. For example, a word co-occurrence graph can be useful to learn about the document semantics, or a disease correlation network is essential to understand the patient health trajectory. In these situations, there are no predetermined structures and thus we need to estimate from data.

There are often an implicit requirement that the network should be sparse and interaction between neighboring variables must be strong. Otherwise, learning will be unreliable since the number of parameters can be quadratic in number of variables. Inference wise, powerful methods such as loopy belief propagation cannot work for densely connected networks, and Gibbs sampling can be too slow to move to reasonably distant states. Practically, learning and inference speed will suffer critically for large-scale problems (e.g., in our application of MRF for collaborative filtering, the number of nodes (items or users) can be hundreds of thousands).

Sparse structure learning is therefore needed but it is often difficult to do efficiently because the space of all possible structures is \( 2^{N(N-1)/2} \) for \( N \) variables. This is because we have a set of \( N(N-1)/2 \) possible symmetric edges and each structure corresponds to a subset of this set (empty set is possible as it means variables are totally independent). To avoid searching through this combinatorial space, often approximate methods are needed.

There are a number of ways to do so. First, we can apply simple screening tests, e.g., by measuring the correlation between variables (e.g., Pearson's), Those pairs whose correlation falls below a threshold will be removed. There are a number of drawbacks here: First, we do not often know which measure is the most appropriate for the task. Second, the method considers pairs independently without paying attention to related pairs. And most importantly, the method is not connected to any objective function, e.g., data likelihood. However, due to its simplicity, this method can be widely applicable in certain problems. For example, in our MRF application to collaborative filtering, this method is quite sufficient.

Another way is to embed the edge selection process into learning. The key is to have a surrogate objective function which can be evaluated efficiently.

Pseudo-likelihood with edge shrinkage

A recent embedded selection method exploits the power of pseudo-likelihood and shrinkage method such as Lasso. Assume that the model energy has the following additive form
\[ E(\mathbf{x},\gamma,\theta) = - \left(  \sum_i \phi_i(x_i,\gamma) + \sum_i\sum_{j>i} \theta_{ij} \psi_{ij}(x_i, x_j) \right) \]
We want to optimise the following objective function
\[ F(\gamma,\theta) = \sum_i \log P(x_i | \mathbf{x}_{\neg i},\gamma, \theta) - \lambda \Omega(\theta) \]
where \( \theta \) is the edge parameter and \( \gamma \) is the rest of parameters. \( \Omega(\theta) \) is the penalty for non-zero values, and \( \lambda > 0\) is the penalty factor.

The role of pseudo-likelihood is clear -- this is one of the consistent methods which is feasible to perform exactly.  Now each edge is associated with a parameter, usually in a log-linear way, and thus a zero parameter would mean no edge. The strategy is then to actively shrink the edge parameters are until they are zeros or ignorable.

This is achieved by placing a sparsity prior on the edge parameter. Ideally the prior is the exponential of the  negative count of non-zeros edge parameters, but this leads to non-smooth optimisation problem. A solution is to relax the count to the sum of absolute parameter values.
\[  \Omega(\theta) = \sum_i\sum_{j>i} |\theta_{ij}| \]
The sparsity is still possible to obtain but the optimisation is much easier.

For more detail, see our recent work: Collaborative filtering via sparse Markov random fields.

Markov random fields for recommender systems IV: Recommending a set

Let us start with the standard setting in collaborative filtering: We are given an incomplete rating matrix and one usual task is to fill those missing entries (aka rating prediction), and another task is to produce a ranked list of new items for each user (or a ranked list of new users for each item).

In this post, we propose to view the two tasks slightly: We want to recommend a set of items for each user, ideally with rating per each item. This differs from the standard rating prediction in that the set itself is unknown. List ranking is now a by product.

There are two variants: One is we do not need to predict the ratings in the set, and another is a ranked list of subsets.

Let us start with the first variant, which is applicable in many real-world situations. When you go shopping at a grocery for the whole week, you usually get a basket full of items of varied qualities. Some items are the same product, but in general those items are different, and together they satisfy the nutrition need and food quality for the family as well as the budget constraints. Another situation is travel package: There are many constraints on the routes, airlines, ticket costs & promotions, waiting time, visa applications, the business deadlines, transits, and luggage allowances. Yet another example is music consumption: You don't always need all best rock pieces but sometimes you want a good mix of them even though some pieces are not very high quality.

Let us assume for now that we will deal only with non-repeated items. The problem is now to specify the best set out of many possible sets. In fact, if there are \( N \) items, there will be \( 2^N -1 \) possible non-empty sets.

One simple solution is to first rank all the items, and then estimate a threshold from which high quality items will be selected. This will sometimes work if items are homogeneous, but this will fail if items are complementary: In the case of shopping baskets, you will end up with high quality items you won't need. Even if items are fairly homogeneous, some items tend to go together more often than with other items:

A better solution should deal with the item set directly.

An immediate question is how can we judge the set quality?

(To be continued)

Markov random fields for recommender systems III: Embedding ordinal matrix factorisation

In the previous posts I have described how Markov random fields (MRFs) can be used to represent the local structures and latent spaces in recommendation. The advantage of MRFs is in its disciplined interpretation and inference. However, the employment of log-linear modelling in MRFs does not enable a simple treatment of ordinal ratings, which we have claimed to be of great importance.

The goal of this post is to develop a MRF-based model that achieves both the representational power of the probabilistic graphical models and and the ease of modelling ordinal rating by matrix factorisation. In short, our new model combines a MRF and an ordinal matrix factorisation (OMF) method. For simplicity, we shall focus only on the MRF with local structures. This is because the latent spaces have been captured by the OMF. Of course, nothing will prevent us from modelling the latent spaces by the MRF itself because the two approaches (OMF and MRF) are very different, one is linear and another is non-linear.

The OMF

Recall that the OMF aims to model an ordinal distribution \(  Q(r|u,i) \). The standard way is to assume that there exists a latent utility function \( f(u,i) \) that captures how much value the user \( u \) judges the item \( i \) . For simplicity we assume that the utility has the linear, parametric form \( f(u,i) = a_u + b_i + W_u' H_i \), although nonlinear and nonparametric options can be possible. The assumption by McCullagh is that the rating is chosen based on the interval to which the utility belongs:
\[ r_{ui} = l  \mbox{  if }  f(u,i) \in (\theta_{l-1},\theta_{l}] \]
for \( l < L \) and
\[ r_{ui} = L \mbox{  if }  f(u,i) > \theta_{L-1} \]
where \( L \) is the number of ordinal levels. Here usually we fix  \( \theta_{0} = -\infty \). Other assumptions are also possible but this is by far the most popular. The probability of receiving a rating is therefore
\[ Q(r=l|u,i) = \int_{\theta_{l-1}}^{\theta_l}P(f(u,i)| \theta)  = F(\theta_l) - F(\theta_{l-1}) \]
 where \( F(\theta_l)  \) is the cumulative distribution evaluated at \( \theta_l  \). The thresholds can be prameterised to depend on user and item but care must be taken to ensure that they are monotonic in \( l \), e.g., \( \theta_{ui,l} = \theta_{ui,l-1} + e^{\eta_{u,l} + \gamma_{i,l}} \).

Depending on the choice of the distribution over \( f(u,i) \) we may end up with the logistic form of \(  F(\theta_l) \) or its probit alternatives. See our recent AAAI'12 paper for a comparison.

The hybrid model of MRF and OMF

Now we wish to build a MRF model by multiplying the point-wise OMF and the item-item interactions:
\[ P(\mathbf{r}|u) \propto \Psi_u(\mathbf{r}) \prod_{i \in I(u)} Q(r_{ui}|u,i) \]
where \( \Psi_u(\mathbf{r}) > 0 \) is the potential function capturing the interaction among items and \( \mathbf{r} \) is the set of items rated by user \( u \). In essence this model promotes the agreement between the latent aspects discovered by the OMF and the local structures in the data. The usual form of the potential function is factorised into pairwise potentials
\[ \Psi_u(\mathbf{r}) = \exp \left(\sum_i\sum_{j>i}\sum_m \beta_{ijm} g_m(r_i,r_j)\right)  \]
where \( \beta_{ijm} \) is parameter and \( g_m(r_i,r_j) \) is feature function. For example, the feature function can be just indicators
\[  g_{lk}(r_i,r_j) = \delta(r_i,l) \delta(r_j,k) \]
as in standard treatment of categorical variables; or bilinear association
\[ g(r_i,r_j) = (r_i - \bar{r}_i)( r_j  - \bar{r}_j)\]
as in usual Gaussian assumption where \( \bar{r}_i \) is the mean rating for item \( i \); or the metric cost
\[ g(r_i,r_j) = |r_i - r_j| ^ p \]
for \( p > 0\).

Model estimation

Like most MRFs we cannot estimate the parameters easily using the standard maximum likelihood approach. Currently, there are two most effective alternatives: The pseudo-likelihood of Besag and the stochastic gradient using truncated MCMC. The pseudo-likelihood leads to exact computation of the loss function and its gradient with respect to parameters, and thus may be preferred by some practitioners. The MCMC-based methods may, on the other hand, lead to better estimation given enough time.

Let us present here the case for pseudo-likelihood. We need to estimate the local conditional distribution
\[ P(r_i | \mathbf{r}_{\neg i}) \propto Q(r_{ui}|u,i) \exp \left(\sum_{j\in I(u), j \ne i}\sum_m \beta_{ijm} g_m(r_i,r_j)\right) \]
Maximising the log of product of all local conditional distributions with respect to parameters will give us the result.

Rating prediction

Predicting the new rating is easy, we just need to compute \( P(r_j | \mathbf{r}) \) in the same way as we have done in the pseudo-likelihood and search for the most probable rating (for metrics such as MAE), or its expectation (for metrics such as RMSE).

Joining user-based models and item-based models

As we have always argued, each item should has its own life, much in the same way as the users. Thus we can build item-based model with user-user interactions. Since user and item are complementary entities, the two model types can be fused into one. The idea is quite simple: We just need to multiply all models together and remove duplicates (because the OMF components will appear twice):
\[ P(\mathbf{R}) \propto \left[\prod_u\Psi_u(\mathbf{r}_u)\right]  \left[\prod_i\Phi_i(\mathbf{r}_i)\right] \prod_{u,i} Q(r_{ui}|u,i) \]
where \( \mathbf{R} \) is the collection of all seen ratings.

Connection to our AAAI'12 paper

In our AAAI'12 paper we suggest the following form of the utility function
\[ f(u,i) = a_u + b_i + W_u' H_i  + \sum_{j \in I(u)} w_{ij} g(r_j) \]
This appears to be quite similar to our local (log) probability in the pseudo-likelihood method. Of course the fine details will be different but in essence the information captured by both the approaches may be similar. Another subtle difference is that in the MRF treatment, the pairwise interactions are symmetric (i.e., \( w_{ij} = w_{ji} \)), while in this model, it is asymmetric \( w_{ij} \ne w_{ji} \).

Previous posts:  Part 1: Modelling local dependency | Part 2: Discovering latent spaces.

Markov random fields for recommender systems II: Discovering latent space

In the previous post we talked about how Markov random fields (MRFs) can be used to model local structure in the recommendation data. Local structures are powerful enough to make our MRF work, but they model only second-order interactions. We now explore the option of discovering latent spaces in the data. The idea is to capture weaker but higher-order interactions globally.


Restricted Boltzmann machines for user latent space

Restricted Boltzmann machines (RBMs) offer an interesting framework: observed ratings can be encoded into the input layer, and the latent factors are represented by the hidden layer. More precisely a RBM is a special type of MRFs where only connections across input nodes and hidden nodes are allowed. As each hidden unit (latent factor) is connected to every input units, high-order correlations among input units are nicely modelled.

Now we proceed to build RBMs for users, assuming for the moment that users are independent in their choices. And thus we build one RBM per user. This is because each user rates only a handful number of items, and thus building the full model of all items is not well justified. For now let us assume that we include only seen items into the model for estimation. At test time, we will introduce unseen items into the model assuming that the model won't change.

For simplicity, assume that latent factors are binary. Denote by \( \mathbf{r} = (r_1,r_2,...,r_N) \) the vector of input variables, and \( \mathbf{h} = (h_1,h_2,..,h_K) \) the vector of latent vectors, where \( h_k \in \{0,1\}\). The energy of the RBM reads
\[ E = - \left ( \sum_i \sum_m\alpha_{im} f_m(r_i) + \sum_k \beta_k h_k + \sum_i\sum_m\sum_k w_{ikm} f_m(r_i)h_k \right ) \]
where \( f_m(r_i) \) is feature function and \( \alpha_{im}, \beta_k, w_{ikm} \) are parameters.

Depending on your argument of the nature of the ratings, we may end up different feature functions, and thus different models such as Gaussian, categorical or ordinal. We can also incorporate item contents and user profiles into the model. For example, profile aspects can be treated as some kind of items which we always know and never have to predict. For this, perhaps a conditional RBM would make sense: We may want to model \( P(\mathbf{r} | \mathbf{p}) \) instead of  \( P(\mathbf{r}, \mathbf{p}) \) where \( \mathbf{p} \) is the profile vector.

The model can now be estimated, either by deterministic methods such as pseudo-likelihood or stochastic methods such as Contrastive Divergence or persistent MCMC. Note that for the pseudo-likelihood method, we need to integrate over the hidden units, so it is not entirely the same as the original Besag's proposal, but the principle is the same.

At test time, we basically look for
\[  r_{j}^{*}  = \arg\max_{r_{j}}P(r_{j}|\mathbf{r})  =  \arg\max_{r_{j}}\sum_{\mathbf{h}}P(r_{j},\mathbf{h}|\mathbf{r}) \]
Thus we need to sum over the hidden layer, which is efficient to do because conditioned on the seen ratings, the graphical model is now a tree. Alternatively, we may resort to mean-field approximation: the stochastic hidden unit is replaced by its deterministic conditional expectation given the seen ratings, i.e., \(h_k \leftarrow P(h_k=1|\mathbf{r}) \). The big O complexity is the same, but it is numerically faster since we do not need to deal with the computation of expensive exponential functions.

A more subtle question is that what if we don't care about exact rating prediction, but any numerical approximation would do? This is often the case in evaluation metrics such as MAE or RMSE. The problem with the RMSE is it penalises zero-one errors to much, while it may be the case that the RBM would make arbitrary zero-one errors. A smoother version of the rating prediction would definitely help here:
\[  \bar{r}_{j}=\sum_{S=1}^{L}P(r_{j}=S|\mathbf{r})S  \]
for the case of \( L \) rating scale.

Finally, the choice of number of hidden unit \( K \) is hard to estimate directly. The random projection theory (e.g., see this KDD paper) suggests that \( K \) is about the same order as the log of the number of users. In practice, however, we can use cross-validation to choose the best value. My experience suggests that once it reaches some reasonable large (e.g., 50-100), the size doesn't matter much, and we don't have to worry too much about overfitting. Perhaps the stochastic nature of the model plays the role here.

A variant of this RBM with categorical assumption of rating (which largely ignores the ordinal constraints) was introduced in ICML'07. It was demonstrated for the first time that RBMs can be a serious competitor for large-scale collaborative filtering. Since then, RBMs are an integrated part of the so-called "ensemble methods" for high performance systems.

Not-so-restricted Boltzmann machines

Recall that in the RBMs, the direct interactions between items are not considered. This is orthogonal to the standard MRF we consider in the previous post. A natural idea is therefore to combine the two approaches: We model both the direct, local interactions and the indirect, global counterparts. Hopefully this would strengthen each method and lead to better predictions.

With energy-based model such as MRFs, this is fairly easy: We just need to add up the energy function and removing the duplicate biases! The results so far are encouraging -- the combination actually works  suggesting that the two approaches capture different type of data aspects.

The joint model of everything

For now we have a semi-restricted Boltzmann machine for each user. However, it is also plausible that each item deserves the same treatment! Besides, the original assumption that the user's choices are independent is not correct: Users only choose those items which are reachable, and it is known that users are influenced by the current choices of their peers.

Thus, we need to integrate all the user-specific models and item-specific models together. Fortunately, this is not so difficult for MRFs. Again, we just need to add energy functions and remove duplicates. The thing we need to worry about is learning for this complicated model. But a closer look that the model structure would suggest a way out: Conditioned on user-specific units, item-specific models are independent. The same holds if we flip the roles of users and items. Thus, an efficient method would be structured pseudo-likelihood (and mean-fields): We look for a block of units that the same time while conditioning on the rest.

Until next time, I wish you enjoy the giant model with a lot of connections and hidden units! See our UAI'09 paper for more technical details.

Previous post: Part 1: learning local dependency | Next post: Part 3: embedding ordinal matrix factorisation

Markov random fields for recommender systems I: Modelling local dependency

As I have mentioned, there are two general approaches to collaborative filtering: the latent space (e.g., using matrix factorisation) and the local structures (e.g., using graph to estimate user-user similarity or item-item similarity). A natural question is then can the two approaches be combined? Numerous evidences have suggested that the answer is yes: The latent space assumption and the local structures assumption seem to be somewhat complementary, or at least do not overlap entirely. I guess they do partly overlap because by using spectral decomposition of the similarity matrices, we can also discover a latent space, but it does not guarantee to be the same as those obtained by factorising the rating matrix.

The only issue is that the combination of the two methods so far appears a bit heuristic, or at least the two approaches do not come from the same unified representation. The quest for such a unification may be more of theoretical interest rather than practical ones: After all with ensemble methods, we can combine any different methods as long as they are fairly diverse.

It turns out, probabilistic graphical models (PGMs) nicely offer such a representation so that the ratings/preferences data can be seen as a probabilistic database from which many queries can be answered probabilistically. For example, rating prediction can be cast as an inference problem where we need to infer the probability that an user will assign a rating value to a particular item. More importantly, it separates out the issue of model building (e.g., structure specification and parameter estimation) from inference (e.g., rating prediction or item ranking).

There are two general ways to specify a PGM: using directed models and using undirected models. Existing directed models for collaborative filtering appeared earlier, and they include Bayesian networks and those variants of the so-called topic models: the PLSA and LDA. The main problem with directed models is that we need to be really careful in our model design to make sure that the structure is a DAG (there will be no loops), and every link from a parent to a child makes some sense. Incorporating features and domain knowledge can be hard and generally requires great design skills to make the model correct.

Undirected models, also known as Markov random fields, Markov networks or their variants such as Botlzmann machines, restricted Boltzmann machines, relational Markov networks, and Markov logic networks are more flexible: They allow loops and can incorporate arbitrary features extracted from data or side constraints. This post will solely focus on this undirected class.

Markov random fields for user-specific local structures

Let us start with the local structures (or the neighbourhood). Similarity scores such as Pearson's correlation or cosine are perfectly OK on their own (parameter estimation stage) but their combination to make prediction does not seem to be theoretically guided (inference stage). Markov random fields (MRFs), on the other hand, allow both the stages to be guided in a principled way.

First we need to specify the model structure. For clarity, let us first start with the assumption that users are largely independent in their choice. We will argue later that this assumption is usually wrong in practice, but let the assumption hold for a while.

Our goal is to build graphs, one per user,  so that ratings are modelled as node variables, and interaction between items is modelled as edges. Here we depart from the standard use of MRFs by building multiple graphs instead of just one. This is because users only rate only a handful number of items, and building the full graph for all items is not very attractive: First, the graph will be too big for efficient handling. Second, most of the variables will be missing, making learning and inference cumbersome. And finally, assuming most missing variables to be 'just' missing does not seem right: Most items shouldn't be there, or at best, are irrelevant to the taste of the user.

The key here is that although graphs are different, they all share the same set of parameters. Otherwise, learning will not happen robustly, and prediction will be impossible with unidentified parameters.

So we end up including only those seen items into the user-specific graph. Note that it is certainly not optimal, because they will be surely other unseen items which will be highly relevant to the user. For example, at test time we need to include unseen items to make inference. But it doesn't seem to be an issue: At test time, we just need to plug one test item at a time, assuming that the user-specific MRF is unchanged in their parameters, and then we do inference \[ r^* = \arg\max_r P(r | \mathbf{r}(u)) \] where\( r \) denotes the rating, \( \mathbf{r}(u) \) denotes the set of seen ratings by user \( u \) .

The next question is how dense the graph should be? Ideally, the graph should be fully connected because all the items will interact with each other under the same user. However, this is not efficient: The number of edges is quadratic in number of items. Second, it makes learning overfit to weakly interacting pairs.

A better way should be making the graph sparse, and it we would eliminate weakly interacting items. If we do the graph structure estimation automatically, this would be the so-called structure learning in probabilistic graphical models. A more heuristic way is to use some correlation measures: If the two items are strongly positively or negatively correlated, their interaction should be kept. However, it is not easy to determine the right amount of interaction strength. An easy way is to keep only K most strongly interacting neigbours for each item.

A hidden question is how should we compute the correlation measure. There are number of information sources we can exploit here: From the ratings alone, we can look for the common users between two items. From the content, we can match two items if their descriptions are similar (e.g., two movies directed by the same directors or played by the same lead actor/actress). Thus the use of content can be implicit.

Are we done with the structure yet? Not quite. MRFs are undirected, and thus the connections are symmetric: if node j is connected to i then i must also be connected to j. However, the edges discovered from the K-nearest neighours are asymmetric: if j is a neighbour of i, it does not necessarily mean i is a neighbour of j. Thus from the K-nearest neighbours per node we need to add more edges to make the neighbourhood consistent.

The next step is parameterisation. We would use the following energy function
\[ E (\mathbf{r}_u) = -\Big ( \sum_i\sum_m \alpha_{im} g_m(r_i) + \sum_i\sum_{j>i} \sum_m\beta_{ijm} f_m(r_i,r_j) \Big ) \] where \( g_m(r_i) \) and \(  f_m(r_i,r_j)  \) are feature functions capturing the properties of the ratings and the interaction between two items; \( \alpha_{im} \) and \( \beta_{ijm} \) are parameters.

We will not go into the details of feature design as it is more an art than science. We wish to comment that we can come up with Gaussian random fields, categorical random fields or ordinal random fields, depending on the choice of these functions. See our UAI'09 paper for the comparison of these three model types.

Note that we can easily incorporate item-content and user-profile into these feature functions, making the model an attractive hybrid network which will help combating against the so-called cold-start problem. This is a clear advantage of undirected models where we can just integrate diverse sources of information with little extra cost.

Once the model parameters have been estimated, the prediction for unseen item j is quite easy:
\[ r_j^* = \arg\min_{r_j} E (\mathbf{r}_u, r_j) \]
From user-model to item-model to integrated model

Let us return to the assumption we made earlier: Users are independent in making their choice. This is not always the case: Users are generally influenced by their peers and the media. From the data representation point of view, the rating matrix does not distinguish the role of columns and rows -- by a transpose we still have the same matrix!

Thus it is plausible that we can build a model for each item in the same way that we did earlier with the users. The only difference is that now we need to model user-user interactions.

Given the fact that we now have either user-specific models or item-specific models, the natural question is can these two approaches be unified?

The answer is yes, and in a simple way: We just need to add the two energies functions together. Now we have a giant MRF connecting every user and every item.

From rating-prediction to item-ranking

Recall that the real goal in recommender systems is less about rating prediction but more about producing a rank list of unseen items for a particular user. It turns out, there is a way: If an new item causes more expected reduction in model energy, the item is likely to be more influential, and it should be ranked higher.

To sum up, we now have a principled way to integrate users, items, contents and profiles into the same model which supports structure specification, parameter estimation, rating prediction and item ranking. For more technical details, we refer to our AusDM'07 paper.

In the next part of the post, we will focus on how MRFs can be used to discover latent spaces and how to integrate with the local structures.

Next posts Part 2: discovering latent spaces | Part 3: embedding ordinal matrix factorisation

Conditional random fields: Big innovation by many small improvements

The subject will be Conditional Random Fields (CRF), which cost me roughly 4 years of PhD. After 7 years since its first appearance in ICML 2001, there has been roughly 1,300 citations, according to Google Scholar. This is a very impressive figure so far. It is impressive because anyone can rightly criticise about CRF's novelty compared to its prior arts.

• From machine learning point of view, it is mostly Maximum Entropy, which was first introduced in 1957 by Jaynes, and had became popular since the work of IBM researchers in 1996. The most popular parameterisation is still log-linear, which has been known to physicists since the time of Gibbs. The famous Generalised Iterative Scaling method for optimising log-linear models was introduced in the 1970s.
• From graphical models point of view, it is just a conditional Markov random field (MRF), i.e., when conditioned on an observational data, we have a MRF of the output variables.
• From modelling point of view, it is just a correction of Maximum Entropy Markov Models (MEMM), introduced in 2000.

Then why is it so popular among researchers who value novelty above everything else? The answer is perhaps, because it works and it came at a right time. Let's see how:

• It was the right time in 2001. Statistical machine learning had reached the point of wide acceptance. At the same time, structural modelling had been widely recognised as indispensable in many domains. Graphical models had been come a mainstream in the AI research agenda. Computing power was strong enough to support large-scale real world applications.
• By 2003, it was discovered that a Limited memory Newton method known as L-BFGS works extremely well with CRFs. It helps to scale up to problems with several million parameters, the size unthinkable for statisticians and optimization hard-core guys. Of course, by 2002, it was known that a L-BFGS cousin, Conjugate Gradients outperform all previous known methods for log-linear models. L-BFGS and CGs are two major advances by the numerical optimisation community, and they have come to rescue.
• MRFs had been around for a long time, since the age of Ising models in physics early the 20th century, and the famous work of Besag in 1974 and the Browns in image modelling in 1984. However, MRFs are quite difficult to model because of its intractability in the joint form. By using the conditional forms, we do not need to model the messy observational data, rather we need only to deal with the output of interest. It makes thing significantly easier and more effective.
• The CRF works like magic because it overcome much inherent limits of the well-known HMM, and the related MEMM. It paves the way for domain experts to engineer features as they please - the things are not so easy in HMMs.

In short, CRFs are not a breakthrough that challenges the current belief. It is a right combination of several advances at the right time and makes incremental improvements in each aspect. As a net result, it became a winner and since then has triggered much interest in learning in structured output space.

Some of our related work (updated):

• Hierarchical semi-Markov conditional random fields for deep recursive sequential data, Truyen Tran, Dinh Phung, Hung Bui, Svetha Venkatesh,  Artificial Intelligence, Minor revision. (Extension of the NIPS'08 paper).
• Preference Relation-based Markov Random Fields in Recommender Systems, Shaowu Liu, Gang Li,Truyen Tran, Jiang Yuan, ACML'15, November 20-22, 2015, Hong Kong.
• Predicting delays in software projects using networked classification, Morakot Choetikertikul, Hoa Khanh Dam, Truyen Tran, Aditya Ghose, 30th IEEE/ACM International Conference on Automated Software Engineering, November 9–13, 2015 Lincoln, Nebraska, USA.
• Tree-based Iterated Local Search for Markov Random Fields with Applications in Image Analysis, Truyen Tran, Dinh Phung and Svetha Venkatesh, Journal of Heuristics, 2014, DOI:10.1007/s10732-014-9270-1 
• Ordinal random fields for recommender systems, Shaowu Liu, Truyen Tran, Gang Li, Yuan Jiang, ACML'14, Nha Trang, Vietnam, Nov 2014.
• MCMC for Hierarchical Semi-Markov Conditional Random Fields, Truyen Tran, Dinh Q. Phung, Svetha Venkatesh and Hung H. Bui. In NIPS'09 Workshop on Deep Learning for Speech Recognition and Related Applications. December, 2009, Whistler, BC, Canada.
• Hierarchical Semi-Markov Conditional Random Fields for Recursive Sequential Data, Truyen Tran, Dinh Q. Phung, Hung H. Bui, and Svetha Venkatesh. In Proc. of 21st Annual Conference on Neural Information Processing Systems, Dec 2008, Vancouver, Canada. [See technical report and thesis for more details and extensions.]
• Constrained Sequence Classification for Lexical Disambiguation, Truyen Tran, Dinh Q. Phung and Svetha Venkatesh. In Proc. of 10th Pacific Rim International Conference on Artificial Intelligence (PRICAI-08), 15-19 Dec 2008, Hanoi, Vietnam.
• Preference Networks: probabilistic models for recommendation systems, Truyen Tran, Dinh Q. Phung and Svetha Venkatesh, In Proc. of 6th Australasian Data Mining Conference: AusDM 2007, 3-4 Dec, Gold Coast, Australia. 
• AdaBoost.MRF: Boosted Markov random forests and application to multilevel activity recognition, Truyen Tran, Dinh Quoc Phung, Hung Hai Bui, and Svetha Venkatesh. In Proc. of  IEEE Conference on Computer Vision and Pattern Recognition, volume Volume 2, pages 1686-1693, New York, USA, June 2006.
• Boosted Markov networks for activity recognition, Truyen Tran, Hung Hai Bui and Svetha Venkatesh, In Proc. of International Conference on Intelligent Sensors, Sensor Networks and Information Processing (ISSNIP2005), 5-8 Dec, Melbourne, Australia.
• Human Activity Learning and Segmentation using Partially Hidden Discriminative Models, Truyen Tran, Hung Hai Bui and Svetha Venkatesh, HAREM 2005: Proceedings of the International Workshop on Human Activity Recognition.