\chapter[Prediction Uncertainty in \acrshort{erss}]{Prediction Uncertainty in Explicit Feedback-Based Recommender Systems}
\label{chapter:explicit}

%\db{Again, we now have two types of uncertainty recommendation/prediction uncertainty and mislabeling uncertainty (label uncertainty). And so maybe the first time you use the word uncertainty, you should use the longer name (recommendation uncertainty or prediction uncertainty or whatever you have called it) instead.}

%\victor{So, I changed the title and added a footnote to make this clear. If you agree with me that this is enough, you can erase these comments.}

\acrlong{erss} have played a pivotal role in the evolution of recommender systems. Notwithstanding their usefulness to this date, these systems have historical importance due to their ability to directly capture users' explicit preferences and opinions. Unlike implicit feedback, where we infer preferences from user behaviour, explicit feedback involves concrete actions such as assigning a numerical rating or writing a review. For many reasons (see Chapter \ref{chapter:rs}), implicit feedback has gained more attention recently. Nevertheless, there are several domains in which practical recommendation applications still rely, at least partially, on explicit feedback. Some examples include recommenders for restaurants, hotels, movies and TV shows, and many more. In most of the aforementioned services, explicit feedback is gathered in the form of numeric ratings, such as opinions on a 5-point scale. This is also the setting that we will confine our attention to in this chapter.

Due to the importance of \acrshort{erss} to the field, both historical and current, this was also the setting adopted by \cite{mazurowski_estimating_2013}, which is one of the seminal works in prediction uncertainty\footnote{This whole chapter considers only prediction uncertainty, and not label uncertainty. For this reason, in this chapter, where we use the word `uncertainty' on its own, it is prediction uncertainty that we are referring to.} estimation for recommenders. Accordingly, we dedicate this chapter to reviewing and extending the methods for uncertainty estimation in \acrshort{ers}. More precisely, in this chapter, we:
\begin{itemize}
\item Present an extensive review of methods for estimating the uncertainty of rating predictions. These methods include information-based, stability-based, error-based, distribution-based and multinomial-based.
\item Propose new techniques for estimating the uncertainty of rating predictions. In particular, we propose a new stability-based method, and also a new error-based method.
\item Present an extensive review of the existing work on evaluating estimates of uncertainty for rating prediction and for top-$n$ recommendation.
\item Extend and improve the techniques and metrics for evaluating uncertainty estimates in RS. In particular, we analyse how uncertainty correlates with dataset statistics and we fix possible issues with existing metrics. We propose a new metric to evaluate uncertainty in the top-$n$ recommendation task.
\item Introduce the idea of uncertainty-aware ranking, drawing in part on a recommendation strategy  proposed in \cite{ortega2021providing}. In particular, we introduce and distinguish between \acrfull{ubf} and \acrfull{prr}, and we give methods for evaluating these two top-$n$ recommendation strategies.
\item Conduct a reproducible empirical study on two large datasets and report results of the most extensive comparison of approaches to uncertainty estimation that we are aware of.
\end{itemize}

In the next section, we formalize the basic setting and notation that are specific for this chapter. In Section {sec:terminology}, we explain some of the roles that uncertainty estimates can play in rating prediction and top-$n$ recommendation. In Section \ref{sec:estimating}, we describe several methods for estimating uncertainty in \acrshort{erss}. After that, Section \ref{sec:evaluating} presents metrics for evaluating these uncertainty estimates. The remainder of the chapter evaluates the uncertainty estimates from Section \ref{sec:estimating} using the evaluation metrics from Section \ref{sec:evaluating}. In particular, Section \ref{sec:empirical} describes the experiments, while Section \ref{sec:results} shows the results. Section \ref{sec:discussion} discusses the uncertainty estimation methods in the light of the experiments. Section \ref{sec:concs} concludes the chapter.

%\derek{I remain strongly of the opinion that this would be much clearer if yu retain the basic approach that we used in the TORS paper, i.e.\ the section can be called Prediction and Recommendation; it can have an opening paragraph and then two sections, one for Rating prediction (excluding the MF stuff that we had at that point in the paper) and one for Top-n recommendation. The second of the two sections would be split into 3 parts: RaBR, UBF and PRR.}

%\victor{I remade the whole section to make it more similar to the TORS paper again. It's a bit repetitive, but more clear.}

\section{Prediction and Recommendation}
\label{sec:terminology}

As already mentioned, in this chapter, we deal with explicit feedback data -- specifically, numerical ratings. Following from Table \ref{tab:notation}, we will designate the ordered set of possible ratings by $S = \{s_1, \ldots, s_{m}\}$, where $m$ is the number of different possible rating values. For example, in the MovieLens dataset and the Netflix dataset that we use in Section \ref{sec:empirical}, $S = \{1,1.5,\ldots,5\}$ and $S = \{1,2,\ldots,5\}$ respectively.

In addition, we will assume that our data consists of a set of users $U$, a set of items $I$, and a set $R$ of user-item-rating triples,  where triple $\langle u, i, r_{ui}\rangle \in R$ means that user $u$'s rating of item $i$ is $r_{ui} \in S$. A user's profile $R_{u \cdot}$ is the set of that user's ratings, $R_{u \cdot} = \{r_{ui} : \langle u, i, r_{ui}\rangle \in R\}$. Similarly, an item's support set $R_{\cdot i}$ is the set of that item's ratings, $R_{\cdot i} = \{r_{ui} : \langle u, i, r_{ui}\rangle \in R\}$.

\subsection{Rating prediction}
\label{sec:rating_prediction}

Given a user $u \in U$ and item $i \in I$, the \textit{rating prediction} task consists in designing or learning a function $F: U \times I \xrightarrow{} \real$ that predicts user-item ratings. That is, ranking prediction is the scoring stage of a recommender that is based on explicit feedback.

% We will write $\hat{r}_{ui}$ for a particular prediction, i.e.\ $\hat{r}_{ui} = \Phi(u, i)$. 

In this work, we want, not only to predict ratings, but also the \textit{uncertainty level} present in a prediction. An uncertainty estimator is a function $G: U \times I \rightarrow \mathbb{R}$, that predicts real-valued quantities associated with each predicted rating, estimating how uncertain the system is about the rating prediction. We will denote by $G_{ui}$, the uncertainty estimate the system associates with predicted rating $F_{ui}$. In all the work that we report in this chapter, uncertainty levels are real-valued, $G_{ui} \in \real$. In some cases, they are probabilistic, which bounds them in $[0,1]$, but this is not always the case. The higher the uncertainty level, the more uncertain the prediction.

Many of the methods for estimating uncertainty that we present in this chapter assume a rating prediction model $F$. Wherever this is the case, we will employ Matrix Factorization (MF) -- see Eq.\ \ref{eq:MF} -- for this underlying model. This is known to be a fairly accurate predictor and, by using it in all cases where we need a separate predictor, we bring a degree of fairness to the comparisons.

\subsection{Top-$n$ recommendation}

Of more importance than rating prediction is \textit{top-$n$ recommendation}. Given a user $u \in U$, the top-$n$ recommendation task consists of retrieving $Z_u^n$ from a set of candidate items $C_u \subseteq I$, such that the retrieved items are presumed relevant to the used.

In Section \ref{sec:top_n}, we discussed the conventional way to perform top-$n$ recommendation in \acrshort{erss}, namely \acrfull{rabr}, which consists of selecting the top-$n$ based on the predicted ratings alone. But, given that we are considering models that can estimate prediction uncertainty, then there is an opportunity to consider other ways of selecting the top-$n$. These use what we can collectively refer to as \textit{uncertainty-aware ranking}. If the estimates of uncertainty are good enough, we may achieve better (or, at least, usefully different) recommendation performance by recommendation strategies that rank the candidates using some combination of the predicted ratings and their estimated uncertainties. We propose two such strategies: \textit{uncertainty-based filtering} and \textit{probability-of-relevance ranking}.

\subsubsection{Uncertainty-based filtering}
\label{sec:ubf}

In some domains, we may not want to recommend items if we are uncertain that the user will find them to be appealing. Whenever this is the case, uncertainty estimates may be used to prevent such recommendations. In the simplest case, we can use \acrfull{ubf}. In \acrshort{ubf}, we filter the candidate items, discarding those whose uncertainty is above some threshold $\tau$.

If $\tau$ is small, it is likely that several of a user's candidate item pool will be filtered-out, which leads to more certain items being recommended. Some of these items may, however, be ones with smaller predicted ratings. Furthermore, some users might have a very small candidate set after filtering. For this reason, there might be cases in which the desired amount of recommendations cannot be provided to some users. On the other hand, larger $\tau$ values may mean that, for some users, few or even no items are filtered-out.

\subsubsection{Probability-of-relevance ranking}
\label{sec:prr}

If predictions are accompanied by uncertainty levels, there may, in certain cases, be an opportunity to rank candidates based on uncertainty. We call this\linebreak \acrfull{prr} and we will present the exact details in later sections. But, in overview, the approach requires the use of a relevance threshold, $\omega \in S$. When a rating is greater than, or equal to, this threshold, then the item is deemed to be relevant to the user. In the datasets we use in Section \ref{sec:empirical}, which use a 5-star rating scale, $\omega$ might be 4, for example. In \acrshort{prr}, we must estimate the probability that $r_{ui} \geq \omega$, $\prob(r_{ui} \geq \omega)$ --- how certain we are that the item is relevant. Then, we rank the candidates in descending order of the relevance probabilities, so that those whose relevance is most certain are highest in the ranking. Finally, \acrshort{prr} recommends the top-$n$ from this ranking.

\acrlong{prr} was first proposed in \cite{ortega2021providing}, specifically for their BeMF method, which we will review in Section \ref{sec:multinomial_based}. Nevertheless, it also applies to many other models, as long as they are able to estimate $\prob(r_{ui} \geq \omega)$. Some of the other methods that we review in Section \ref{sec:estimating} will be able to compute such probabilities by default. For some of those that cannot compute probabilities by default, we have devised a way of enabling some of them to nevertheless compute probabilities, and this is described in Section \ref{sec:results_prr}.

Notice that, with \acrshort{prr}, $\prob(r_{ui} \geq \omega)$ is used to rank the items. If one also wishes to associate each recommended item with a numeric uncertainty estimate (e.g.\ for display to the user), then the complementary probability $1 - \prob(r_{ui} \geq \omega)$ can be used.

\section{Estimating Uncertainty}
\label{sec:estimating}

There are many ways of estimating the level of uncertainty of a predicted rating in an \acrshort{erss}. Different sources of uncertainty can motivate different ways of estimating the uncertainty. For the purposes of surveying the approaches, we classify them into five broad types, as follows:
\begin{description}
\item [Information-based:] Uncertainty is estimated from the amount of information that is known about users or items, or the dispersion of the ratings.
\item [Stability-based:] Uncertainty is estimated from the stability of the predictions across perturbations of the conditions under which we train the prediction model.
\item [Error-based:] Uncertainty is estimated from the expected prediction errors of an underlying rating prediction model.
\item [Distribution-based:] The ratings are assumed to follow a statistical distribution, whose dispersion is used to estimate uncertainty.
\item [Multinomial-based:] Multinomial methods also estimate uncertainty from the dispersion of the predicted ratings but they use models that make discrete predictions for each rating value. 
\end{description}

In the following sections, we describe each of these classes in more detail; we review work from the literature that falls into each class; and we even propose new stability-based and error-based methods. 

\subsection{Information-based uncertainty}
\label{sec:info_based}

The foremost source of uncertainty in a \acrshort{rs} comes from the lack of user-item interaction data. In practice, item catalogs are often large and therefore users will interact with only a small proportion of the items in it. In our setting, this means that users will have rated only a small proportion of the items: in the majority of cases, $\#R^- \gg \#R$. 
Models learned from this sparse data may be unreliable. For this reason, simple statistics such as item support (the number of ratings an item has) can be used to measure the reliability of the predictions \cite{mazurowski_estimating_2013}. Hence, its negative can be used as an estimate of uncertainty:
\begin{description}
\item [$\isup$:] $G_{ui} = - \#R_{\cdot i}$
\end{description}

Although item support is indicative of the information available about an item, it might also be misleading. For instance, two items might have the same number of ratings, but one of them might be rated similarly by different users, whereas the other might have received more divergent ratings. We expect a higher degree of uncertainty when predicting the ratings of items that have divergent ratings. For this reason, the variance of an item's ratings can also be used as an estimate of uncertainty \cite{adomavicius2007towards,mazurowski_estimating_2013}:

\begin{description}
\item [$\ivar$:] $G_{ui} = Var(R_{\cdot i})$
\end{description}

In a similar vein, we can estimate uncertainty from user profile length (the number of ratings a user has) \cite{mazurowski_estimating_2013} or the variance of a user's ratings \cite{wu2013user}:
\begin{description}
\item [$\usup$:] $G_{ui} = - \#R_{u \cdot}$
\item [$\uvar$:] $G_{ui} = Var(R_{u \cdot})$
\end{description}

Information-based measures of uncertainty can be insightful. Nevertheless, they are quite simple. One weakness is that they are model-independent, because they are computed strictly from the data, independently of the rating prediction model. In other words, changing from one prediction model to another does not result in any change to the estimates of uncertainty. This is counter-intuitive, because different models utilise the data differently, and we would expect uncertainty estimates to reflect this behaviour. This will not be the case for information-based measures.

A second weakness is that $\isup$ and $\ivar$ quantify uncertainty item-wise and not interaction-wise. That is, their estimated uncertainty for an item $i$ is not personalised: $G_{ui} = \rho_{vi} \, \forall \, (u, v) \in U^2$. Similarly, $\usup$ and $\uvar$ are strictly user-wise measures: $G_{ui} = \rho_{uj} \, \forall \, (i, j) \in I^2$.

Used in isolation, item-wise and user-wise measures lack the granularity that \acrshort{rss} typically need. This is especially the case for user-wise measures. A user-wise uncertainty estimate cannot help a system decide which of a set of recommendations need explanations, for example, since all of a user's recommendations will have the same uncertainty level. Similarly, user-wise uncertainty estimates cannot be used in top-$n$ recommendations that are uncertainty-aware, such as \acrshort{ubf} and \acrshort{prr}, again because the estimated uncertainty associated with each of that user's predicted ratings will be the same. In fact, since top-$n$ recommendation is a more important task than rating prediction and since top-$n$ recommendation is a user-based task, user-wise uncertainty estimates are largely useless. For that reason, the literature on information-based uncertainty is focused on item-wise metrics \cite{mazurowski_estimating_2013}. In our experiments (Section \ref{sec:empirical}), we do not include user-wise uncertainty estimates.

Nevertheless, user-wise and item-wise estimates of uncertainty can provide insights into \acrshort{rs} performance. For example, \cite{bernardis2019estimating} show both theoretically and empirically how eigenvalues obtained from an item-item similarity matrix are strongly related to the model recommendation accuracy. The negatives of these eigenvalues can also be interpreted as user-wise uncertainty metrics. 

Finally, we note that, instead of using these estimates in isolation, some authors have proposed heuristic combinations of $\usup$ and\linebreak $\isup$ (e.g.\ based on their product) to give interaction-wise measures (e.g.\ \cite{Zhang-Et-Al-2014,daCosta-Et-Al-2018}).

\subsection{Stability-based uncertainty}
\label{sec:stability_based}

Prediction uncertainty is related to the level of arbitrariness present in the predictions. Hence, a prediction is uncertain the more it exhibits high variance across small perturbations to the conditions under which the prediction model is learned. 

\cite{mazurowski_estimating_2013} explores this idea in a multi-modelling strategy. Consider that a model $F$ is used to predict ratings. Then, model stability can be measured by learning several different models from perturbations of the training dataset. A stochastic perturbation function, $\mathcal{P}$, is used to create $N$ different versions of the original ratings training data, i.e.\ $\mathcal{R}^{(k)}= \mathcal{P(\mathcal{R})}$, for $k \in 1,\ldots, N$. A rating prediction model $F^{(k)}$ is learned from each of the versions of the dataset. Considering the predicted ratings as random variables, their standard deviation is an estimate of that rating's uncertainty.

In a strategy called $\resample$, \cite{mazurowski_estimating_2013} uses random sampling as the perturbation function. In this case, each individual model is trained using a random subset of the ratings available for training. The predicted rating $\hat{r}_{ui}$ is given by a model $F$ trained on the entire training set; the estimated uncertainty of that prediction is given by the standard deviation of the predictions from the models $F^{(k)}$ that are trained on the different samplings of the training set:
\begin{align}
    \hat{r}_{ui} & = F(u, i) \label{eq:resample_r} \\
    G_{ui} & = \frac{1}{N} \sqrt{\sum_{k=1}^N (F^{(k)}(u, i) - \hat{r}_{ui})^2}
    \label{eq:resample_rho}
\end{align}

\cite{mazurowski_estimating_2013} proposes another such strategy, called $\inject$. In this strategy, the perturbed datasets are obtained by adding noise to the original ratings, i.e.\ new ratings, $r'_{ui}$ are obtained as $r_{ui} + \epsilon$, where $\epsilon$ is random Gaussian noise. The empirical results in \cite{mazurowski_estimating_2013} show that $\resample$ generates better estimates of uncertainty compared with $\inject$, which is why, in our own empirical work (Section \ref{sec:empirical}), we include $\resample$ and not $\inject$.

Unlike the information-based uncertainty measures, which are model-\linebreak independent, the stability-based measures do depend on the prediction algorithm itself. Furthermore, the stability-based measures are interaction-wise measures. On the other hand, in order to obtain good estimates for the predicted rating variance, several different models have to be learned, hence, computing these measures can be slow. In \cite{mazurowski_estimating_2013}, this problem motivates $\fast$, which approximates $\resample$ by using in each $F^{(k)}$ some model parameters that have been computed once on the original ratings. 

In this dissertation, we introduce  a perturbation strategy that has not previously been evaluated in the context of uncertainty estimation. It is model-specific and, indeed, depends on using models whose parameters require random initialisation. The strategy, which we designate $\ensemble$, consists in multi-modelling for different initialisation points. The predicted rating is the mean of the predicted ratings of the individual models, and the estimated uncertainty is their standard deviation as before:
\begin{align}
    \hat{r}_{ui} & = \frac{1}{N}\sum_{k=1}^N F^{(k)}(u, i) \label{eq:ensemble_r} \\
    G_{ui} & = \frac{1}{N} \sqrt{\sum_{k=1}^N (F^{(k)}(u, i) - \hat{r}_{ui})^2}
    \label{eq:ensemble_rho}
\end{align}

Using different initialisation points for the perturbations in multi-modelling is interesting since it means that all of the training ratings are used to learn each model $F^{(k)}$; by contrast, $\resample$ uses only a subset of the training data for each model. This might lead to better rating estimation, and we see some evidence of this in the results of our experiments.

For our empirical work (Section \ref{sec:empirical}), the models in the ensemble use matrix factorization (MF) (Eq.\ \ref{eq:MF}). MF is known to suffer from strong instability, depending on the initialisation of the vectors. This has been shown by \cite{gabbolini2022instability} in a recommendation context and also by \cite{belford2018stability} and \cite{pena2020combining} for topic modelling. Ensembles of MF models that have been initialised differently can result in lower error and greater stability but also enable us to estimate uncertainty.

\subsection{Error-based uncertainty}
\label{sec:error_based}
 
Uncertainty in rating predictions can be estimated through the expected error in those predictions. The more certain a prediction is, the more accurate the prediction should be. \cite{zhu2018assigning} give a method for estimating prediction error, and then use this as the measure of uncertainty.

Their method use a cross-validation procedure to obtain a set of prediction errors $E$. First, they partition the training ratings $R$ into $K$ folds. It follows that each known rating $r_{ui} \in R$ will appear in exactly one fold. Then, they learn $K$ models, $F^{(k)}$ for $k \in 1,\ldots,K$, where $F^{(k)}$ is learned from all the training ratings except those in the $k$-th fold. Now, they can make a prediction $\hat{r}_{ui}$ that corresponds to each known rating $r_{ui} \in R$. Specifically, if $r_{ui}$ is in the $k$-th fold, then $F^{(k)}$ predicts $\hat{r}_{ui}$. This allows a calculation of training error, $e_{ui} = r_{ui} - \hat{r}_{ui}$. At the end of this cross-validation phase, they have a set of prediction errors $e_{ui} \in E$ for each $r_{ui} \in R$. Finally, they learn two models. One is a ratings prediction model $F$, learned from the entire ratings training set $R$. The second is an error prediction model, $\mathcal{E}: U \times I \rightarrow \real^+$, that can predict the rating prediction errors, which is learned from $E$.

At prediction time, we want to compute uncertain predictions, $\langle\hat{r}_{ui}, G_{ui}\rangle$. The predicted rating, $\hat{r}_{ui}$, is given by the final rating prediction model, $F$; the uncertainty level, $G_{ui}$, is the estimated prediction error, given by model $\mathcal{E}$.

In \cite{zhu2018assigning}, FunkSVD was chosen both as $F$ and $\mathcal{E}$. That is, one FunkSVD model $F$ is responsible for the rating predictions, while another FunkSVD model $\mathcal{E}$ estimates uncertainty (i.e.\ errors). In Section \ref{sec:empirical}, we designate this strategy by $\FunkSVD-CV$. More precisely, in $\FunkSVD-CV$, we have:
\begin{align}
\hat{r}_{ui} &= p_u^Tq_i \\
G_{ui} &= p'_u{}^T q'_i
\end{align}
where $p_u$ and $q_i$ are embeddings learned from $R$ and $p'_u$ and $q'_i$ are embeddings learned from $E$.

Nevertheless, FunkSVD is only one of many possible uncertainty estimators $\mathcal{E}$. During our empirical study, we found $\FunkSVD-CV$ to perform below our expectations. For some of our evaluation metrics, much simpler estimates, such as $\isup$, had better performance. These results raised a concern about possible over-complexity of FunkSVD for the uncertainty estimator. For this reason, we also experimented with a much simpler choice of $\mathcal{E}$, which we refer to as $\Bias-CV$.

In $\Bias-CV$, we define linear weights $b_u \in \real$ for $u \in U$ and $b_i \in \real$ for $i \in I$ for every user and item. Now, we propose that the  uncertainty for a user-item pair (i.e.\ the predicted error) to be simply the sum of the user's and item's  weights. With this formulation, instead of $\mathcal{E}$ having a $d$-dimensional embedding to model each user and item, there is only a single scalar parameter for each.
More precisely, in $\Bias-CV$, we have:
\begin{equation}
G_{ui} = b_u + b_i
\end{equation}
The user and item weights, $b_u$ and $b_i$, are learned, by gradient descent, to minimize the mean squared error between the predicted and observed errors.

Finally, we also acknowledge the contribution by \cite{cleger2013being}, where another uncertainty estimation approach based on prediction errors was introduced. There, uncertainty estimation is viewed as a binary classification task, where the goal is to segregate certain from uncertain predictions. This task is solved by learning a classifier that predicts whether a user-item interaction is uncertain or not. The classifier is learned based on both the prediction errors committed by the rating prediction model (in their case, a $k$-nearest neighbours model) and the ratings given by similar users to the same item. We will not further explore this approach for two reasons. First, it is restricted to neighborhood models. Second, the segregation of items into two types (certain and uncertain) is quite different from the formulation used in our work, where we seek to estimate the level of uncertainty, $G_{ui}$ (Section \ref{sec:terminology}). In any case, $\FunkSVD-CV$ and $\Bias-CV$ are also based on prediction errors and can be applied to a wider class of collaborative filtering models.

\subsection{Distribution-based uncertainty}
\label{sec:distribution_based}

The \acrshort{rs} literature contains several methods that consider each user-item rating to be a random variable, each having a certain probabilistic distribution $\mathcal{F}| \Theta_{ui}$, where $\Theta$ is the set of distribution parameters. Modelling the ratings with a distribution rather than a point estimate allows the recommender to produce a much richer output, from which it is possible not only to obtain point estimates for the ratings but also to infer the uncertainty associated with these point estimates.

The Gaussian distribution emerges as a natural choice: its mean, $\mu$, and variance, $\sigma^2$, are obvious choices for the rating point estimates and their uncertainties. In this case, $r_{ui} \sim \norm (\mu_{ui}, \sigma^2_{ui})$, and the task of estimating $r_{ui}$ is performed by estimating $\Theta_{ui} = \langle\mu_{ui}, \sigma^2_{ui}\rangle$ for every user-item pair.

One notable algorithm to employ a Gaussian model is Probabilistic Matrix Factorization (PMF) \cite{mnih2008probabilistic}. In PMF, the ratings are assumed to follow a Gaussian distribution, with a fixed variance parameter $\sigma^2 = 1$. As before, let $p_u$ and $q_i$ be user and item embeddings, both of dimension $d$. Then, the likelihood of the training data can be written as:
\begin{equation}
\label{eq:likelihood_exp}
    \prob (R) = \prod_{i \in I} \prod_{u \in U} \left[\norm(r_{ui} | p_u^T q_i, 1) \right]^{\delta_{ui}}
\end{equation}
where $\norm(\cdot | \mu, \sigma^2)$ is the probability density function of the Gaussian distribution with mean $\mu$ and variance $\sigma^2$, and $\delta_{ui}$ is an indicator function that is equal to 1 if user $u$ rated item $i$ and 0 otherwise. The latent factors are assumed, in prior, to follow a zero-mean spherical Gaussian distribution. In this case, maximising the log-posterior distribution is equivalent to minimising the sum of squared errors subject to $L_2$ regularisation terms: 
\begin{equation}
    \sum_{u \in U} \sum_{i \in I} \delta_{ui} (r_{ui} - p_u^T q_i)^2 + \lambda_U \sum_{u \in U} ||p_u||  + \lambda_I \sum_{i \in I} ||q_i||
\end{equation}

To learn PMF from explicit ratings, we can use mini-batch gradient descent, iterating through the observed rating. The procedure is equivalent to the one described in Section \ref{sec:learn_exp}. The loss function of a batch of ratings $B$ is,

\begin{align}
    \label{eq:lossPMF}
    \mathcal{L(\mathcal{B}|\theta)} = \frac{1}{\#\mathcal{B}} \sum_{(u, i, r_{ui}) \in \mathcal{B}} (r_{ui} - p_u^T q_i)^2 + \lambda_U \sum_{u \in U} ||p_u||  + \lambda_I \sum_{i \in I} ||q_i||
\end{align}

Notice that, although having a different motivation, the PMF model optimizes an objective function which is similar\footnote{In fact, user and item specific regularization parameters $\lambda_U, \lambda_I$, may also be used with FunkSVD, making it equivalent to PMF.} to  FunkSVD's (Eq.\ \ref{eq:FunkSVDloss}).

PMF assumes every user-item pair to have the same variance. Not only is this a very strong assumption in general, it also means that PMF is not suitable for uncertainty estimation, since every user-item pair would have the same estimated uncertainty $\sigma^2$. For this reason, \cite{wang2018confidence} extend  PMF to Confidence-aware Probabilistic Matrix Factorization (CPMF). There, each user and item is assumed to have a latent variance parameter, $\sigma_u$ and $\sigma_i$, respectively. The variance of each user-item pair is obtained by some function $\gamma: \real \times \real \xrightarrow{} \real^+$ that composes these two latent variance parameters. Therefore, the likelihood of the training data is:
\begin{equation}
    \label{eq:cpmf}
    r_{ui} \sim \norm(p_u^T q_i, \gamma(\sigma_{u}, \sigma_{i}))
\end{equation}
The latent factors, $p_u$ and $q_i$, and the variance parameters, $\sigma_u$ and $\sigma_i$, are estimated through gradient descent on the negative log-likelihood. 

\cite{wang2018confidence} use a simple product for $\gamma$. In our experiments (Sections \ref{sec:empirical}), we do the same and this is what we are referring to when we write $\CPMF$ in that later section of this chapter.

\cite{wang2018confidence} also propose a Bayesian version of the algorithm, similar to \cite{salakhutdinov2008bayesian}. The idea of the Bayesian framework is to automatically control model complexity through prior distributions in order to alleviate the effects of overfitting and to improve model generalisation.

Another distributional method that resembles $\CPMF$ is proposed by \cite{yang2021multi}. Their proposal, which is based on Gaussian Processes, also assumes each interaction to follow a Gaussian distribution. However, their method explicitly models the full covariance between all pairs of interactions. For this reason, their model is more complete then $\CPMF$, but also a lot slower to train, especially for big datasets.

The next section presents further distribution-based methods for the special case where ratings are not treated as continuous variables.

\subsection{Multinomial-based uncertainty}
\label{sec:multinomial_based}

Ratings are typically discrete; for example, the ratings in the Netflix Prize dataset are integers from 1 to 5 \cite{bennett2007netflix}. But, in rating prediction, ratings are often treated as continuous variables. In other words, the model predicts real values: $F: U \times I \xrightarrow{} \real$.  Then, as we have seen, there are various ways of augmenting each predicted rating with an estimated uncertainty level.

Multinomial models treat predicted ratings as discrete variables. For each point on the rating scale, they predict a probability. Hence, the output of a multinomial prediction model is a vector, containing one probability for each score $s$ on the rating scale $S = \{s_{1}, \ldots, s_m$\}, where $\prob(r_{ui} = s)$ is the probability that user $u$ will assign rating $s \in S$ to item $i$, and $\sum_{s \in S} \prob(r_{ui} = s) = 1$. For instance, in the Netflix case, the model will predict $\prob(r_{ui} = s) \, \forall s \in \{1,\ldots,5\}$. This idea has motivated several recommendation algorithms, including the User Rating Profile Model (URP) \cite{marlin2003modeling}, Bernoulli Matrix Factorization (BeMF) \cite{ortega2021providing} and OrdRec \cite{koren2011ordrec}. From a statistical standpoint, this way of formulating rating prediction consists of assuming $r_{ui}$ follows a multinomial distribution. Hence, it is possible to extract not only point estimates for the ratings, but also dispersion statistics, which can be used as estimates of uncertainty.

The BeMF model is one way to estimate the vector of probabilities \cite{ortega2021providing}. BeMF learns $m$ different Bernoulli factorization models, one for each point on the rating scale, i.e.\ it learns a set of models $F = \{F_1, \ldots, F_m\}$. To learn each individual model $F_s$, it uses a modified version of the training set $R_s = \{(u, i, r_{uis})\}$, such that:
\begin{equation}
r_{uis} = \begin{cases}
1, & r_{ui} = s \\
0, & \text{otherwise}
\end{cases} \forall (u, i, r_{ui}) \in R
\end{equation}
Each Bernoulli factorization model $F_s$ makes predictions in the same way as the FunkSVD method (Equation \ref{eq:MF}), but with the difference that the predictions are scaled through a logistic function, which bounds them to the $[0, 1]$ interval. To learn the latent factors, for each individual model, the cross-entropy loss function is used. The loss for a single training instance is:
\begin{equation}
   F_s(u, i)^{r_{uis}} (1 - F_s(u, i))^{1 - r_{uis}}
\end{equation}
The parameters of the model are updated by performing gradient descent on this loss function. Finally, the probabilities are given by:
\begin{equation}
    \prob(r_{ui} = s) = \frac{F_s(u, i)}{\sum_{s \in S} F_s(u, i)}
\end{equation}

With the probabilities estimated, \cite{ortega2021providing} proposes that the predicted rating will be the most probable value in the rating scale, and the uncertainty estimate will be the probability mass lying outside this most probable rating value:
\begin{align}
\label{eq:BeMF_r} 
        \hat{r}_{ui} & = \arg_{s \in S} \max \prob(r_{ui} = s) \\
\label{eq:BeMF_rho}
        G_{ui} &= 1 - \max_{s \in S} \prob(r_{ui} = s)
\end{align}

A rather different multinomial approach, OrdRec, is proposed in \cite{koren2011ordrec}. While OrdRec has a very different formulation from BeMF, its ultimate goal is also to predict the vector of probabilities. OrdRec learns a single rating prediction model, e.g.\ using FunkSVD, which will be denoted $F$. Then, for every $s \in \{s_1,\ldots,s_{m-1}\}$, the cumulative probability function for the ratings is:
\begin{equation}
    \prob(r_{ui} \leq s) = \frac{1}{1 + e^{F(u, i) - t_{us}}}
    \label{eq:prob_ordrec}
\end{equation}
where $t_{us}$ for each $s \in S$ is a user-specific set of threshold parameters. Notice that $\prob(r_{ui} \leq s_m) = 1$ as this is a cumulative distribution. After estimating the cumulative functions, the individual probabilities are obtained by:
\begin{equation}
    \prob(r_{ui} = s_i) = \prob (r_{ui} \leq s_i) - \prob (r_{ui} \leq s_{i-1})
\end{equation}

Now, for every $u \in U$, let $t_{u0} = 0$. Then, the gaps between a user's thresholds are encoded as a set of parameters $\beta_{u1}, \beta_{u2},\ldots, \beta_{um-1}$, such that
\begin{equation}
    t_{us_i} = t_{us_{i-1}} + e^{\beta_{us_i}}, \quad \forall i \in \{1,\ldots,m-1\}
\end{equation}

Therefore, the whole set of parameters to be learned for OrdRec includes the underlying model's parameters, as well as $\{\beta_{u1}, \beta_{u2}, \ldots, \beta_{um-1}\} \forall u \in U$. The likelihood function for the model is:
\begin{equation}
    \prod_{i \in I} \prod_{u \in U} \prod_{s \in S} \left(\prob(r_{ui} = s)\right)^{\delta(r_{ui} = s)}
\end{equation}
Learning is performed by stochastic gradient descent on the negative log-\linebreak likelihood function.

With the rating distribution estimated, there are many ways in which rating predictions and uncertainty estimates can be extracted. \cite{koren2011ordrec} argue that the distribution's average and some measure of its dispersion can be used to predict ratings and uncertainty respectively. Moreover, any distribution dispersion metric can be employed, such as standard deviation, Gini impurity or entropy. Among these, they empirically found that the standard deviation performed best for uncertainty and, therefore, we will use the same when we estimate uncertainty using OrdRec in this chapter:
\begin{align}
\label{eq:ordrec_r}
        \hat{r}_{ui} &= \sum_{s \in S} s \, \prob(r_{ui} = s) \\
\label{eq:ordrec_rho}
        G_{ui} &= \sqrt{\frac{1}{m} \sum_{s \in S} \left(s^2 \prob(r_{ui} = s)\right) - m\hat{r}_{ui}^2}
\end{align}

Moreover, we highlight one concern regarding the uncertainty estimates given by Equation \ref{eq:ordrec_rho}: when $\hat{r}_{ui}$ is either too small or too big, that is, close to $s_1$ or to $s_m$, then the uncertainty $G_{ui}$ will always be very low. This happens because, in order to predict one of the distribution's extremes, the vast majority of the distribution's probability have to fall on this extreme rating value, leading to very low standard deviation. Intermediate rating values do not suffer from the same problem because, in this case, the distribution's probability can be dispersed symmetrically around that intermediate rating value.

We have concluded our review of ways of \emph{estimating} ratings prediction uncertainty. We turn our attention now to ways of \emph{evaluating} these estimates.

\section{Evaluating Uncertainty Estimates}
\label{sec:evaluating}

After training a rating prediction model, it is common to evaluate its performance on a test set, this being a set of known ratings that was held-out during training. The most common rating prediction evaluation metrics, such as root mean squared error (RMSE), simply compare the predicted ratings to the known ratings in the test set (see Eq.\ \ref{eq:rmse}).

In this case, the known ratings give us a `ground-truth' for the evaluation. Similarly, after estimating the uncertainty of predicted ratings, it is desirable to evaluate the quality of the uncertainty estimates. But, the evaluation of the uncertainty estimates is not as simple as evaluating the predicted ratings because there are no `ground-truth' values that can be compared with the uncertainty estimates.

Our review of ways of evaluating estimates of rating prediction uncertainty starts with a discussion of expected correlations. Subsequent sections then correspond to the two different tasks that an \acrshort{ers} perform: rating prediction and top-$n$ recommendation, the latter being split into \acrlong{rabr}, \acrlong{ubf} and \acrlong{prr}.

\subsection{Expected correlations}

Despite the differences in the uncertainty estimates that we reviewed in Section \ref{sec:estimating}, there are `patterns' that we expect to observe, including the following:
\begin{description}
    \item [Negative correlation with user profile size:] The less we know about a\linebreak user's tastes, the more uncertain we expect that user's predictions ratings to be.
    \item [Negative correlation with item support:] Similarly, the less that an item has been interacted with, the less certain we expect that item's predicted ratings to be.
    \item [Positive correlation with user rating variance:] The more a user's ratings vary from each other, the less certain we expect that user's predictions to be.
    \item [Positive correlation with item rating variance:] The more an item's rat-\linebreak ings vary from each other, the less certain we expect that item's predictions to be.
\end{description}

To the best of our knowledge, we are the first to propose checking these correlations. Computing them can serve as an initial check that a particular way of estimating uncertainty behaves, on the whole, in the way we might expect. We have chosen to carry out this check using Spearman rank correlation, because the correlations might not always be linear.

It is important to point out that some of these correlations might not always hold for individual user-item interactions. For example, in the movie domain, a new movie from a famous director is usually a very safe and expected recommendation, even though the movie, being newly-released, will have low item support.

We also expect to see some correlations between the uncertainty estimates and prediction error and top-$n$ recommendations accuracy:
\begin{description}
\item [Prediction error:] The more wrong a predicted rating is, the less certain we expect the prediction to be.
\item [Recommendation accuracy:] Similarly, the more accurate a set of top-$n$ recommendations is, the more certain we expect the predicted rating for each item in the top-$n$ to be.
\end{description}
The literature offers several metrics that we can compute to confirm these correlations, and we explore these in more detail in the following subsections.

\subsection{Rating prediction}
\label{sec:eval_rating_prediction}

Uncertainty measures are expected to reflect prediction errors. For this reason, the foremost method for evaluating the quality of estimates of rating prediction uncertainty consists of simply evaluating the correlation between the prediction errors in a holdout test set and their estimated uncertainty.   Different correlation coefficients can be used. For example, Pearson correlation can be used to measure the linear correlation, while Spearman correlation can be applied if rank correlation is preferred. Later, we refer to these correlations (between rating prediction errors and uncertainty estimates) as Pearson$_\rho$ and Spearman$_\rho$.

Mazurowski and Maciej \cite{mazurowski_estimating_2013} proposes a visualization of the correlation between uncertainty and prediction error. In fact, they refer to ``reliability estimates'' on predicted ratings, rather than uncertainty estimates. But, as their concept of reliability is the opposite of uncertainty, their visualisation can be used here, after some adaptations. In what follows, we `translate' their proposals to our uncertainty framework.

Their proposal begins by splitting the predictions for the test set into $B$ equal-sized bins based on, and ordered by, their uncertainty levels --- that is, bin 1 contains the least uncertain predictions and bin $B$ the most uncertain ones. For each bin $b$, an error threshold, $\epsilon_b$ is calculated. The error threshold for bin $b$ measures how inaccurate the predicted ratings within the bin are. The thresholds are calculated such that Bayesian confidence intervals of the form $[\hat{r}_{ui} - \epsilon_b, \hat{r}_{ui} + \epsilon_b]$ contain the real rating $r_{ui}$ for at least $(1-\alpha) \times 100 \%$ of the bin's predictions, where $\alpha \in [0, 1]$ is a fixed significance level. The quality of the confidence intervals are determined by the curve of the error thresholds against their confidence bin. 

One problem with the latter method is that the curves will depend on the choice of $\alpha$. In fact, the reliance on $\alpha$ means that each uncertainty estimate will have a family of curves (for each value of $\alpha$), which makes it more difficult to compare uncertainty estimates. To solve this issue, and also to adapt Mazurowski and Maciej \cite{mazurowski_estimating_2013} confidence intervals to our uncertainty framework, we propose to use uncertainty bins, rather than confidence bins, and to calculate the average RMSE in each of these uncertainty bins. In this case, given $B$ equal-sized uncertainty bins, we calculate the following for each bin $b$:
\begin{equation}
    RMSE_b = \sqrt{\frac{\sum_{G_{ui} \in b} (\hat{r}_{ui} - r_{ui})^2}{\# G_{ui} \in b}}
\end{equation}
where the denominator is the number of test instances falling into bin $b$. Similarly to before, we then plot the curve of the $RMSE_b$ values against their bin index. We refer to these as RMSE-uncertainty curves. Now these curves have no reliance on additional parameters other than the number of reliability bins --- in particular, there is no longer the parameter $\alpha$. 

Nevertheless, a problem that remains is how to interpret curves like these. To illustrate, we refer to Figure \ref{fig:example curves}. The Figure shows three example RMSE-uncertainty curves. Notice that all the curves have very similar $RMSE_1$ and $RMSE_B$ values. On the other hand, their behaviours in between these extremes are different. In this case, there is no clear definition that states which curve is preferred, and hence which of the underlying uncertainty estimates are better.

\begin{figure}[tp]
    \centering
    \includegraphics[width=\linewidth]{img/RMSE-uncertainty-example.pdf}
    \caption{Example RMSE-Uncertainty curves.}
    \label{fig:example curves}
\end{figure}

To alleviate this problem, some simple metrics can be extracted from the curves. One option, which has been explored in \cite{mazurowski_estimating_2013}, is to calculate the difference in the  values between the first and last bins. In \cite{mazurowski_estimating_2013}, this metric is called $d_w$. Here, to avoid confusion due to the fact that we are using RMSE-uncertainty curves, we will refer to it as $\delta RMSE$:

\begin{equation}
    \delta RMSE = RMSE_B - RMSE_1
\end{equation}
Higher values of $\delta RMSE$ are indicative that  prediction error grows as uncertainty grows.

In a similar vein, \cite{bobadilla_reliability_2018} argue that an evaluation metric for an uncertainty measure should penalise predictions where certainty is high but predictive error is also high and should strongly reward predictions where uncertainty is high and predictive error is high. Table \ref{tab:penalty_reward} schematically describes the desired behaviour. 
\begin{table}[t]
    \centering
    \begin{tabular}{ccc} 
    Prediction error & Uncertainty  & Uncertainty evaluation metric      \\ \hline
    High             & Low        & Big penalty \\ \hline
    Low              & High         & Small penalty   \\ \hline
    High             & High         & Big reward  \\ \hline
    Low              & Low        & Small reward    \\ \hline
    \end{tabular}
    \caption{Criteria for a good estimate of rating prediction uncertainty}
    \label{tab:penalty_reward}
\end{table}

\cite{bobadilla_reliability_2018} propose a metric that matches the criteria in Table \ref{tab:penalty_reward}. Like the work in \cite{mazurowski_estimating_2013}, Bobadilla et al.'s Reliability Prediction Improvement metric is given in terms of reliability. We adapt it from reliability to uncertainty and call it Uncertainty Prediction Improvement (UPI). 

To define UPI, we first compute the absolute errors $e_{ui} = |\hat{r}_{ui} - r_{ui}|, \forall (u, i, r_{ui}) \in T$, where $T$ is a holdout test set. Then, we define the metric criteria, $C^{UPI}$, as follows:
\begin{equation}
    C^{UPI}_{ui} =  e_{ui}(e_{ui}-\overline{e})(\overline{\rho}-G_{ui})
\end{equation}
where $\overline{e}$ and $\overline{\rho}$ are the means of the absolute errors and the uncertainty estimates for all the test instances in $T$. UPI is then defined as:
\begin{equation}
    UPI = \ddfrac{\frac{1}{\sigma_e \sigma_\rho \#T}\sum_{r_{ui} \in T} C^{UPI}_{ui}}{\overline{e}}
\end{equation}
where $\sigma_e$ and $\sigma_\rho$ are the standard deviations of the absolute errors and uncertainty estimates. The divisor term, $\sigma_e \sigma_\rho \#T$, standardises the criteria values, while the term $\overline{e}$ compares the standardised criteria values to the mean absolute error.

\cite{koren2011ordrec} take a rather different approach to the evaluation of uncertainty estimates. Instead of using metrics, such as $\delta RMSE$ and UPI, they build and evaluate what they call \textit{confidence classifiers}.
These are binary classifiers that use the uncertainty estimates as predictors (features) to predict whether or not the rating error will be greater than 1. In practice, given a set of predicted ratings $\hat{r}_{ui}$ and their respective observed ratings $r_{ui}$, a confidence classifier is a logistic model, where the uncertainty estimate is the predictive feature, and the model predicts the probability $\prob(|\hat{r}_{ui} - r_{ui}| > 1)$. Once trained, the confidence classifier is evaluated using a holdout set. Its accuracy on the holdout set is an indicator of the effectiveness of the uncertainty estimate. Koren and Sill \cite{koren2011ordrec} use the classifier's Area Under the Curve (AUC) for this.

Care must be taken to avoid training and evaluating the confidence classifier on the instances used to train the rating predictor. In our experiments, we train the rating predictor on the training set, then we use 2-fold cross-validation on the test set to train and evaluate the confidence classifier. In other words, we train the confidence classifier on half of the test instances and evaluate its AUC on the other half of the test instances. We repeat this, after interchanging the two folds. The final results that we report are the average AUC across the two folds. In Section \ref{sec:results}, we refer to this approach to evaluation as Error-Uncertainty Classification (EUC).

We highlight that the choice to predict $\prob(|\hat{r}_{ui} - r_{ui}| > 1)$ with this method is arbitrary. The method is general: we could build confidence classifiers based on thresholds other than 1. In our experiments, however, we follow \cite{koren2011ordrec}'s use of 1.

\subsection{Top-$n$ recommendation}
\label{sec:eval_top_n_recommendation}

We divide our discussion of ways of evaluating the quality of the uncertainty estimates at top-$n$ this into three. First, we discuss the evaluation of uncertainty in recommendations given with \acrshort{rabr}. Then, corresponding with the use cases introduces in Section \ref{sec:terminology}, we propose an application-based approach, where we access the validity of the uncertainty estimates according to the improvements they provide to these use cases.

\subsubsection{Rating-based ranking}
\label{sec:eval_rabr}

As previously stated, we expect recommendation accuracy to be greater in cases where the recommendations are less uncertain. Therefore, one metric to evaluate the uncertainty estimates is the Spearman correlation between some recommendation accuracy metric, e.g.\ MAP@$n$, and the average uncertainty of the recommended items for the user. We will refer to this as Uncertainty Accuracy Correlation (UAC). Notice that, for this metric, smaller values denote better performance.

In a similar vein, \cite{bobadilla_reliability_2018} propose Reliability Recommendation Improvement, a metric that evaluates reliability estimates in the case of recommendations. Their metric is based on similar criteria to those in Table \ref{tab:penalty_reward}, this time that hits (top-$n$ recommendations that are correct) should be associated with lower uncertainty levels than the non-hits. We present Uncertainty Recommendation Improvement (URI), which is our adaptation from reliability to uncertainty.

Our adaptation starts by defining $C^{URI}$ as follows:
\begin{equation}
    C^{URI}_{u} = \sum_{i \in Z_u^n | r_{ui} \geq \omega} ( \overline{\rho} - G_{ui})
\end{equation}
Notice that there is a strong relation with $C^{UPI}$. However, URI is focused on relevant items, i.e.\ test set ratings for which $r_{ui} \geq \omega$. Hence, in $C^{URI}_{u}$, the mean uncertainty $\overline{\rho}$ is calculated only over test instances where $r_{ui} \geq \omega$.

In this definition of $C^{URI}$, based on \cite{bobadilla_reliability_2018} and similar to $C^{UPI}$, $\overline{\rho}$ is a global average --- the average of all the uncertainty values in test set $T$. In other words, it is not the average of these values just for user $u$. We believe that averaging over all $T$ can lead to problems. For example, it is likely that some users will have high overall uncertainty. In an extreme case, there may even be users $u$ for which $G_{ui} < \overline{\rho} \quad \forall i \in I$. For such users, if a recommendation hit occurs, then $C^{URI}_{u} < 0$ even if the recommended hit is less uncertain than all of the other recommended items for that user. Given a set of recommendations for a user, it is desirable for recommendation hits to have a lower estimated uncertainty than non-hits. For this reason, we propose an alternative definition of URI that is based on $Z_u^n$, the top-$n$ items recommended to $u$. We redefine it as follows:
\begin{equation} 
\label{eq:uri_crit}
    C^{URI}_{u} = \sum_{i \in Z_u^n | r_{ui} \geq \omega} ( \overline{\rho}_u - G_{ui})
\end{equation}
where $\overline{\rho}_u$ is the mean uncertainty for the items in $Z_u^n$, i.e.\ the $n$ items that are recommended to user $u$. 

Now, the URI is defined as:
\begin{equation}
\label{eq:URI}
    URI = \ddfrac{\sum_{u} C^{URI}_{u}}{{\sigma_{\rho_u} \sum_{u}\#\{i \in Z_u^n | r_{u, i} \geq \omega\}}}
\end{equation}
where $\sigma_{\rho_u}$ is the standard deviation of the uncertainty values for $Z_u^n$. The denominator is used to guarantee that the URI is a valid metric when comparing different uncertainty measures. 

\subsubsection{Uncertainty-based filtering}
\label{sec:eval_ubf}

In Section \ref{sec:ubf}, we argued that, in some circumstances, it can be useful to constrain the set of items that might be included in a top-$n$ by discarding candidates whose uncertainty value exceeds some threshold $\tau$. As initially proposed in \cite{ortega2021providing}, we can use this as a way of evaluating uncertainty estimates. In certain cases, the recommendation accuracy (e.g.\ its precision) can increase --- if the candidate items that are discarded are not only uncertain but also not relevant. In practice, the recommendation accuracy can be measured for a range of different values for $\tau$, with the expectation that accuracy will increase as $\tau$ increases. Nevertheless, restricting the candidate set might also decrease the average predicted rating in the recommended items. For instance, before filtering, a user might have several candidates that receive high predicted ratings but also high uncertainty. If these are filtered out, not only will uncertainty be reduced in the filtered recommendation set, but also the average predicted rating. For this reason, together with the accuracy, it is important to monitor the Mean Predicted Rating. To do so, we use:

\begin{equation}
    \text{Mean Predicted Rating}@n = \frac{1}{\#U} \sum_{u \in U} \frac{1}{n} \sum_{i \in Z_u^n} \hat{r}_{ui}
    \label{eq:avg_relevance}
\end{equation}

Finally, since we are preventing the RS from recommending highly uncertain items to the user, there may be users for whom the RS cannot recommend a full set of $n$ items. Hence, the recommendation \emph{coverage} should also be measured. In our experiments, we measure coverage according to Eq.\ \ref{eq:coverage}.

\subsubsection{Probability-of-relevance ranking}
\label{sec:eval_prr}

In \acrshort{rabr}, uncertainty estimates are not used for ranking the candidates and selecting the top-$n$, and therefore uncertainty in the case of \acrshort{rabr} has to be evaluated through specific metrics, such as UPI. In \acrshort{prr}, on the other hand, uncertainty estimates are used to rank the candidates and select the top-$n$; specifically, we recommend the candidates for which $P(r_{ui} \geq \omega)$ is highest. This means that, in the case of \acrshort{prr}, the uncertainty estimates can be evaluated  through recommendation accuracy metrics. That is, the quality of the recommendation set given by the uncertainty estimates serve as an indirect way to evaluate these estimates.

Furthermore, in the cases where the system is expected to provide uncertainty estimates together with the uncertainty-based recommendations, then these uncertainty values also have to be evaluated. In Section \ref{sec:prr}, we propose that $1 - \prob(r_{ui} \geq \omega)$ can be used as the uncertainty estimate. In this case, the uncertainty estimates will increase as we progress through the item ranking. That is, the uncertainty of the top-ranked item will be the smallest, followed by the second, and so on. Under this setup, URI is not a feasible evaluation metric anymore, because the differences $\Bar{\rho}_u - G_{ui}$ will always be the largest for the top-ranked items, leading to unrealistically high URI values. On the other hand, the correlation between the uncertainty estimates and accuracy (UAC) is still a feasible metric, and we will use it below.

\section{Empirical Study}
\label{sec:empirical}

In this section, we report the results of an empirical study that we have conducted. We believe it is the most extensive study in the field to this date. Compared with previous studies, it covers more ways of estimating uncertainty and uses more evaluation methods, as well as making use of two large datasets. Previously, the most comprehensive study was \cite{mazurowski_estimating_2013}, but it compared only information-based and stability-based uncertainty estimates. The study in \cite{bobadilla_reliability_2018} proposed and empirically tested some uncertainty quality measures but, similarly, it did not consider error-based, distribution-based or multinomial-based uncertainty estimates. Other studies in the field, e.g.\ \cite{koren2011ordrec, wang2018confidence}, use few, or no, uncertainty evaluation metrics, and therefore also lack empirical support for their uncertainty estimates. Finally, \cite{ortega2021providing} has two major issues with their experiments: it does not compare $\BeMF$ against $\OrdRec$, which has a similar motivation to it; and it has an unconventional recommendation evaluation, where candidates for recommendation are taken only from the test set and therefore unrated items are never recommended.

Our implementation makes considerable use of PyTorch \cite{paszke2017automatic} and is partially inspired by the Spotlight Python libraries.\footnote{\url{https://github.com/maciejkula/spotlight}} For reproducibility, all the code used is available.\footnote{\url{https://github.com/vcoscrato/uncertain/tree/SnapExplicit}}
\subsection{Models}

We compare the following ways of estimating uncertainty, each of which was described in Section \ref{sec:estimating}: $\isup$, $\ivar$,\linebreak $\resample$, $\ensemble$, $\FunkSVD-CV$, $\Bias-CV$, $\CPMF$, $\BeMF$ and\linebreak $\OrdRec$. As explained previously, we do not include the user-wise estimates, $\usup$ or $\uvar$: they assign the same level of uncertainty to all of a user's candidate items, which is not helpful to a top-$n$ recommender system.

Most of these techniques rely on a separate rating predictor; only $\CPMF$, $\BeMF$ and $\OrdRec$ have their own methods for predicting the ratings. For those that need a separate rating predictor, we use FunkSVD (Equation \ref{eq:MF}). Moreover, FunkSVD is also the underlying rating prediction model $F$ in the case of $\OrdRec$ and $\BeMF$, and has the same formula as the distribution average in the case of $\CPMF$ (Equation \ref{eq:cpmf}).

Furthermore, we also use FunkSVD separately as a baseline for recommendation accuracy.  We will want to determine whether the new techniques have similar or better prediction error (RMSE) and top-$n$ accuracy (MAP@$n$ and Recall@$n$) than FunkSVD. 

\subsection{Datasets}

We use two large ratings datasets: the MovieLens 25 million ratings dataset \cite{harper2015movielens} and the Netflix prize dataset \cite{bennett2007netflix}. Table \ref{tab:datasets} shows some information about each dataset. We chose these datasets because they are by far the most popular ones for the rating prediction task.
\begin{table}[t]
    \centering
    \begin{tabular}{|c|c|c|} \hline
         \textbf{Dataset} & \bf{MovieLens} & \bf{Netflix} \\\hline
         \textbf{num. users}, $\#U$ & 162,542 & 53,424 \\\hline 
         \textbf{num. items}, $\#I$ & 59,047 & 480,189 \\ \hline
         \textbf{num. ratings}, $\#R$ & 25,000,095 & 100,480,507 \\\hline
         \textbf{sparsity} & 99.74\% & 99.61\% \\\hline
         \textbf{avg. num.\ ratings per user} & 154 & 1881 \\\hline
         \textbf{avg. num.\ ratings per item} & 423 & 209 \\\hline
         \textbf{rating scale} & $\{0.5, 1, 1.5, \ldots, 5\}$ & $\{1, 2, \ldots, 5\}$ \\\hline
    \end{tabular}
    \caption{Characteristics of the datasets used in our experiments.}
    \label{tab:datasets}
\end{table}

\subsection{Methods}

\subsubsection{Data splitting} 

We select a random sample of 10,000 users to be test users. The restriction to 10,000 users is to reduce the computational cost of running this extensive set of experiments on such large datasets. For these sampled test users, we chronologically order their ratings and place the last 20\% of them in the test set. We use the remaining ratings (from \emph{every} user, not just the test users) for training or validation. Specifically, we order these remaining ratings chronologically; for each user, we use the first 80\% of these ratings for our training set, and the rest for our validation set.

\subsubsection{Tuning} 

We have a number of hyper-parameters, whose values we need to choose. We tune the models by training them with different hyper-parameter values on the training set and evaluating them on the validation set, choosing values that give the lowest validation set RMSE.

Wherever we are using matrix-factorization, the latent feature vector dimension $d$ and the regularisation strengths $\lambda_U$ and $\lambda_I$ are hyper-parameters. The latent feature vector dimension is chosen from $\{50, 100, 200\}$, while the regularisation strengths are chosen from $\{0.1, 0.01, 0.001\}$. The latent feature vectors are initialised by sampling from a zero-centered Gaussian distribution with a standard deviation of $0.01$. We set the learning rate to 0.0001 in all cases and used the Adam optimiser to account for learning rate adaptations. Furthermore, to suppress the need to optimise the number of training epochs, we employed early-stopping, monitoring RMSE after every epoch and stopping training when validation RMSE does not improve for five consecutive epochs. This setup is similar to the one used in \cite{gabbolini2022instability}, for example.

For $\ensemble$ and $\resample$, which involve multi-modelling (Eqns.\ \ref{eq:ensemble_r} and \ref{eq:ensemble_rho}), we use $N=5$ different models. For $\resample$, we use an $80\%$ sample fraction. For the error-based methods, $\FunkSVD-CV$ and $\Bias-CV$, we use 2-way cross-validation to create the error matrix $E$. We conducted some preliminary experiments with more multi-modelling runs (10 and 20) and more CV-folds (5 and 10), but we obtained similar results to the ones for $N=5$ and for 2-way cross-validation.

$\BeMF$ trains multiple MF models, one per rating value. For fairness, we ensure that $\BeMF$ has a similar number of parameters to all the other models. For example, for the Netflix dataset, where there are 5 possible rating values, if the other models are using $d=50$, then $\BeMF$ will learn five models with dimension 10 each.

\subsubsection{Rating prediction evaluation} 

To evaluate rating predictions, we measure the RMSE between the ratings in the test set and the predicted ratings. To evaluate the uncertainty estimates, we use the following evaluation techniques that we presented in Section \ref{sec:eval_rating_prediction}: Pearson$_\rho$, Spearman$_\rho$, UPI, RMSE-Uncertainty curves, $\delta RMSE$ and $EUC$.

\subsubsection{Top-n recommendation evaluation} 
We compute recommendations for each of the 10,000 test users. For each of these users $u$, the candidate items $C_u \subseteq I$ that can be recommended to $u$ are all the items that $u$ has has not rated in either the training or validation sets. Then, for each user, we create a top-$n$, $Z_u^n \subseteq C_u, |Z_u^n| = n$. The way we do this depends on whether we are evaluating \acrshort{rabr}, \acrshort{ubf} or \acrshort{prr}. We measure the accuracy of a top-$n$ using MAP@$n$\footnote{For readability, we chose to present only one of the two position-aware recommendation accuracy metrics presented in Chapter \ref{chapter:rs}. The choice for MAP over NDCG was arbitrary. From our experience, they are generally strongly correlated.} and Recall@$n$, averaged over all 10,000 users. An item $i$ is considered relevant to a user $u$ if its test set rating $r_{ui} \geq \omega$ where $\omega = 4$. We evaluate the uncertainty estimates using the evaluation techniques that we presented in Section \ref{sec:eval_top_n_recommendation}, namely UAC and URI.

\section{Results}
\label{sec:results}

In this section, we present the results. The subsections correspond to the ones in Section \ref{sec:evaluating}.

\subsection{Expected correlations}
\label{sec:results_correlations}

We start by exploring correlations between the uncertainty estimates and each of: the user profile size ($\#R_{u\cdot}$), item support ($\#R_{\cdot i}$), user rating variance ($Var(R_{u\cdot})$) and item rating variance ($Var(R_{\cdot i})$).

We compute Spearman correlations on a set of 100,000 randomly sampled user-item pairs. In other words, we pick a random user and a random item, we estimate the uncertainty $G_{ui}$ and, once we have done this for 100,000 such pairs, we correlate with the statistics mentioned in the previous paragraph. Selecting users and items at random avoids a bias: if we had instead selected users and items from $R$ ($r_{ui} \in R$), then we would bias these correlations to cases where the user has rated the item.

Figure \ref{fig:corr} gives the results. The correlations between $\isup$ and $\#R_{\cdot i}$ and between $\ivar$ and $Var(R_{\cdot i})$ are trivially perfect and so they are not included.

\begin{figure}[!htbp]
\centering
\begin{subfigure}{0.8\linewidth}
    \centering
    \includegraphics[width=\linewidth]{img/ML/corr.pdf}
    \caption{MovieLens dataset.}
\end{subfigure}
\begin{subfigure}{0.8\linewidth}
    \centering
    \includegraphics[width=\linewidth]{img/Netflix/corr.pdf} 
    \caption{Netflix dataset.}
\end{subfigure}
    \caption{Correlations between uncertainty estimates and dataset statistics. We expect negative correlation with $\#R_{u\cdot}$ and $\#R_{\cdot i}$; we expect positive correlation with $Var(R_{u\cdot})$ and $Var(R_{\cdot i})$.}
    \label{fig:corr}
\end{figure}

In the MovieLens dataset, the models' correlation with user profile size, user variance and item variance are very low; $\Bias-CV$ and $\BeMF$ show strong negative correlation with item support, as expected; but there are counter-intuitive positive correlations with item support for $\ensemble$, $\FunkSVD-CV$ and $\OrdRec$. In the Netflix dataset, some expected positive correlations with user variance and item variance can be observed, with the most noticeable values being achieved by $\Bias-CV$; in the case of item support, $\Bias-CV$ and $\BeMF$ show strong negative correlation, as expected, while $\resample$, $\ensemble$ and $\FunkSVD-CV$ show positive correlations to it, which is counter-intuitive.

Notwithstanding the heuristic nature of these expected correlations, the results are quite concerning. Although they are all supposed to be estimates of uncertainty, they behave differently from each other, and there is no uncertainty estimate that convincingly exhibits all the expected behaviours. We highlight, among these concerning results, how $\FunkSVD-CV$ and $\Bias-CV$ have an opposite correlation with item support. This is particularly odd because both methods are based on the errors of the same baseline, and yet their uncertainty estimates differ drastically. While some of the forthcoming results in the reminder of this work will show good performance for these methods, the correlations we have observed raise a flag to the credibility of uncertainty estimates based on prediction error.

In the remaining subsections, we evaluate the uncertainty estimates using the methods that we reviewed in Section \ref{sec:evaluating}.

\subsection{Rating prediction}
\label{sec:results_prediction_error}

We begin our evaluation of rating prediction by comparing the prediction error of the different models, leaving uncertainty aside. Results are shown in Table \ref{tab:RMSE}. Several methods ($\isup$, $\ivar$ and $\resample$) are not explicitly listed in this Table. This is because, for rating prediction itself, they use the baseline model, and therefore their RMSE is the same as that of the baseline.
\begin{table}[t]
    \centering
    \begin{tabular}{|c|c|c|} \cline{2-3}
    \multicolumn{1}{c|}{} & \multicolumn{2}{c|}{\bf{RMSE}} \\\hline
    \bf{Model}    & \bf{MovieLens} & \bf{Netflix} \\ \hline
    Baseline & 0.8400         & 0.8618       \\ \hline
    $\ensemble$ & \textbf{0.8294}         & \textbf{0.8493}       \\ \hline
    $\CPMF$     & 0.8447         & 0.8742       \\ \hline
    $\BeMF$    & 1.0031         & 1.0270       \\ \hline
    $\OrdRec$   & 0.8832         & 0.8666       \\ \hline
    \end{tabular}
    \caption{Rating prediction: RMSE. (Results for $\isup$, $\ivar$, $\resample$, $\FunkSVD-CV$ and $\Bias-CV$ will be the same as the Baseline.) The best values (lowest) are highlighted in bold.
    }
    \label{tab:RMSE}
\end{table}

It is worth mentioning that, as explained in Section \ref{sec:estimating}, $\CPMF$, $\BeMF$ and $\OrdRec$ do not directly optimise RMSE, while the baseline model and the models within the ensemble do. In light of this, $\CPMF$ shows good performance for both datasets, while $\ensemble$ is able to outperform the baseline in both cases (keeping in mind that low values are better). $\OrdRec$ has a slightly worse performance in the MovieLens dataset, while being closer to the other models in the Netflix dataset. $\BeMF$ has the worst performance in both datasets.

To assess the quality of the uncertainty estimates from a prediction error point of view, we employ the methods described in Section \ref{sec:eval_rating_prediction}. The RMSE-uncertainty curves are given in Figure \ref{fig:curves}. On the whole, we see what we would expect: as uncertainty grows so does error. In both datasets, $\Bias-CV$ seems to be the curve that shows this most strongly. In the Netflix dataset, $\OrdRec$ and $\FunkSVD-CV$ also exhibit this behaviour quite strongly. However, as stated before, these curves are relatively hard to interpret, and for this reason we refrain from drawing strong conclusions from them and turn our attention to the various other evaluation metrics that we presented in Section \ref{sec:eval_rating_prediction}. The values for these  metrics are shown in Tables \ref{tab:rating movielens} and \ref{tab:rating netflix}. 

\begin{figure}[!htbp]
\centering
\begin{subfigure}{\linewidth}
    \centering
    \includegraphics[width=\linewidth]{img/ML/RMSE-uncertainty.pdf} 
    \caption{MovieLens dataset.}
\end{subfigure}
\begin{subfigure}{\linewidth}
    \centering
    \includegraphics[width=\linewidth]{img/Netflix/RMSE-uncertainty.pdf} 
    \caption{Netflix dataset.}
\end{subfigure}
\caption{RMSE-uncertainty curves.}
\label{fig:curves}
\end{figure}

\begin{table}[t]
\begin{tabular}{|l|r|r|r|r|r|} \hline
{} &  Pearson$_\rho$ &  Spearman$_\rho$ &      $\delta RMSE$ &     UPI &  EUC \\ \hline
$\isup$ &               0.0218 &                0.0391 &  0.5064 &  0.1370 &          0.5353 \\ \hline
$\ivar$    &               0.1266 &                0.0958 &  0.3034 &  0.6029 &          0.5674 \\ \hline
$\resample$         &               0.1187 &                0.1196 &  0.3243 &  0.5205 &          0.5788 \\ \hline
$\ensemble$         &               0.0327 &                0.0215 & -0.0865 &  0.0875 &          0.5051 \\ \hline
$\FunkSVD-CV$       &               0.1788 &                0.2008 &  0.0310 &  0.5500 &          0.6245 \\ \hline
$\Bias-CV$          & \textbf{0.3463} & \textbf{0.2928} &  \textbf{0.9675} &  \textbf{1.6851} & \textbf{0.6982} \\ \hline
$\CPMF$             &               0.0449 &                0.0863 &  0.3009 &  0.2444 &          0.5632 \\ \hline
$\BeMF$            &               0.1114 &                0.1189 &  0.3022 &  0.3204 &          0.5859 \\ \hline
$\OrdRec$           &               0.2260 &                0.2526 &  0.5552 &  0.6924 &          0.6731 \\ \hline
\end{tabular}
\caption{Rating prediction uncertainty evaluation: MovieLens dataset. The best values (highest) are highlighted in bold.}
\label{tab:rating movielens}
\end{table}
\begin{table}[t]
\begin{tabular}{|l|r|r|r|r|r|} \hline
{} &  Pearson$_\rho$ &  Spearman$_\rho$ &      $\delta RMSE$ &     UPI &  EUC \\ \hline
$\isup$ &               0.0111 &                0.0277 &  0.1408 &  0.0233 &          0.5170 \\ \hline
$\ivar$    &               0.1340 &                0.1089 &  0.3897 &  0.5075 &          0.5742 \\ \hline
$\resample$         &               0.0872 &                0.0753 &  0.2567 &  0.3632 &          0.5483 \\ \hline
$\ensemble$         &               0.0405 &                0.0262 &  0.0527 &  0.1761 &          0.5152 \\ \hline
$\FunkSVD-CV$       &               0.2161 &                0.1884 &  0.4370 &  0.8133 &          0.6218 \\ \hline
$\Bias-CV$          & \textbf{0.2832} & 0.2398 &  \textbf{0.7382} &  \textbf{1.0839} &          0.6587 \\ \hline
$\CPMF$             &               0.0752 &                0.0684 &  0.2174 &  0.2362 &          0.5424 \\ \hline
$\BeMF$             &               0.1585 &                0.1755 &  0.5319 &  0.3726 &          0.5949 \\ \hline
$\OrdRec$           &               0.2788 & \textbf{0.2674} &  0.6843 &  0.7673 & \textbf{0.6669} \\ \hline
\end{tabular}
\caption{Rating prediction uncertainty evaluation: Netflix dataset. The best values (highest) are highlighted in bold.}
\label{tab:rating netflix}
\end{table}

The best overall results are obtained by $\Bias-CV$ and $\OrdRec$, while the worst results come from $\ensemble$. In general, the different metrics in the table seem quite correlated. That is, the performances of the models have nearly the same ordering across all metrics and both datasets. We find these results to be reassuring, since all the metrics are designed to measure the strength of the relationship between uncertainty and rating prediction error. On this, a few observations can be made. $\Bias-CV$'s uncertainty seems to correlate linearly with error since their Pearson correlation is stronger, while $\OrdRec$'s might correlate non-linearly given the higher Spearman correlation. Moreover, $\Bias-CV$ performed strongest according to UPI, while the results for other metrics, such as EUC, are closer, with even a slight advantage for $\OrdRec$ in the Netflix dataset. 

\subsection{Top-$n$ recommendation}
\label{sec:results_top_n_recommendation}

The results in the previous section evaluate the measures of uncertainty through the prism of prediction error. But, as we noted in Section \ref{sec:top_n}, we care more about the quality of the top-$n$ recommendations than the prediction error. Hence, using the ideas from Section \ref{sec:eval_top_n_recommendation}, we here compare the uncertainty estimates in a variety of top-$n$ recommendation scenarios.

\subsubsection{Rating-based ranking}
\label{sec:results_rabr}

In \acrlong{rabr}, the top-$n$ are the $n$ candidate items with the highest predicted ratings. We begin our evaluation by comparing the accuracy of the different models, leaving uncertainty aside. Figure \ref{fig:precision_recall} shows the MAP@n and Recall@n, averaged over all test users, for different values for $n \in 1,2,\ldots,10$. 

The most evident result in the Figure is the very low performance of $\BeMF$ in these metrics. This is to be expected. It happens because $\BeMF$ is only able to predict rating scale values, that is, $\hat{r}_{ui} \in S$; unlike the other models, it cannot predict a rating that lies between two values on the scale. Therefore, many candidates will receive the same predicted rating -- they `tie'. The creation of a top-$n$ from such an ordering involves a lot of random tie-breaking. This is presumably why BeMF is not evaluated for \acrshort{rabr} in \cite{ortega2021providing}. Apart from $\BeMF$, the remaining models perform similarly, which is reassuring. But there are exceptions: $\OrdRec$ has worse MAP for the MovieLens dataset; $\CPMF$ is notably worse for the Netflix dataset; and both $\OrdRec$ and $\CPMF$ have poor recall on the Netflix dataset.

\begin{figure}[tp]
\centering
\begin{subfigure}{\linewidth}
    \centering
    \includegraphics[width=\linewidth]{img/ML/accuracy.pdf}
    \caption{MovieLens dataset.}
\end{subfigure}
\begin{subfigure}{\linewidth}
    \centering
    \includegraphics[width=\linewidth]{img/Netflix/accuracy.pdf}
    \caption{Netflix dataset.}
\end{subfigure}
\caption{\acrshort{rabr}: MAP@n and Recall@n for $n \in 1,2,\ldots,10$. (Results for $\isup$, $\ivar$ and $\resample$ will be the same as the Baseline.)}
\label{fig:precision_recall}
\end{figure}

Turning our attention now to the evaluation of the uncertainty estimates, in the \acrshort{rabr} scenario our evaluation tools are URI and UAC. The results are shown in Table \ref{tab:rabr_uncertainty}. Remember that higher (positive) URI values denote better performance, while the opposite is the case for UAC. We first notice that the two metrics mostly agree. For both of them, $\BeMF$ showed better performance in the MovieLens dataset (although we know that its results are affected by the random choices in the cases of tied ratings), while $\CPMF$ was best in the Netflix dataset. On the negative side, $\resample$ and $\ensemble$ showed poor performance, together with $\FunkSVD-CV$ and $\Bias-CV$, which was one of the best performers on the rating prediction task. Therefore, we see a clear disconnect between the models' ability to estimate uncertainty on the rating prediction and top-$n$ recommendation tasks.

\begin{table}[tp]
    \centering
    \begin{tabular}{|l|r|r|r|r|} \cline{2-5}
    \multicolumn{1}{c|}{} & \multicolumn{2}{c|}{MovieLens} & \multicolumn{2}{c|}{Netflix} \\ \hline
    Model &  URI &  UAC &  URI &  UAC \\ \hline
    NEG-ITEM-SUPPORT &   0.5070 &        0.0135 &   0.6747 &       -0.1292 \\ \hline 
    ITEM-VARIANCE    &  -0.0590 &        0.0316 &  -0.1561 &        0.2121 \\ \hline
    RESAMPLE         &  -0.2371 &        0.1138 &  -0.3554 &        0.1470 \\ \hline
    ENSEMBLE         &  -0.3399 &        0.2125 &  -0.4491 &        0.4141 \\ \hline
    FunkSVD-CV       &  -0.0867 &       -0.0218 &  -0.2825 &        0.1733 \\ \hline
    Bias-CV          &  -0.1548 &       -0.0489 &  -0.2983 &        0.0037 \\ \hline
    CPMF             &   0.5196 &       -0.0242 &   \textbf{0.6944} & \textbf{-0.2905} \\ \hline
    BeMF             &   \textbf{0.9527} & \textbf{-0.1619} &  0.6043 & -0.0492 \\ \hline
    OrdRec           &   0.2833 &        -0.0283 &   0.3808 &       -0.2143 \\ \hline
    \end{tabular}
    \caption{\acrlong{rabr} uncertainty evaluation: URI and UAC on MovieLens and Netflix datasets. The best values (highest on URI and lowest on UAC) are highlighted in bold.}
    \label{tab:rabr_uncertainty}
\end{table}

\subsubsection{Uncertainty-based filtering}
\label{sec:results_ubf}

In \acrfull{ubf}, we prevent uncertain candidates (ones whose uncertainty estimate exceeds $\tau$)  from being included in the top-$n$.  In the results that we report here, we progressively exclude bigger fractions of the most uncertain candidates, reducing the overall uncertainty of the recommendation lists. More specifically, for each uncertainty estimate, we obtain a distribution using 100,000 randomly-selected user-item pairs and we obtain the 20\%, 40\%, 60\% and 80\% percentiles of these distributions as cut-offs for uncertainty. We have chosen this percentile approach due to the different ranges of values that the uncertainty estimates have. Figure \ref{fig:cuts} shows the MAP@10 and Mean Predicted Rating@10 of the recommendation lists for these different cut-offs. The same Figure also shows the recommendation coverage with each cut applied and a MAP@10 value that considers only users to whom a full 10 recommendations can be given.

\begin{figure}[tp]
\centering
\begin{subfigure}{\linewidth}
    \centering
    \includegraphics[width=\linewidth, height=230pt]{img/ML/cuts.pdf} 
    \caption{MovieLens dataset.}
\end{subfigure}
\begin{subfigure}{\linewidth}
    \centering
    \includegraphics[width=\linewidth, height=230pt]{img/Netflix/cuts.pdf} 
    \caption{Netflix dataset.}
\end{subfigure}
\caption{\acrshort{ubf}: MAP, Mean Predicted Rating, Coverage and the MAP for users for whom a full set of recommendations can be made. This is for top-10 recommendations with increasing $\tau$ to filter different quantiles.}
\label{fig:cuts}
\end{figure}

The only uncertainty estimates whose MAP increases as uncertain recommendations are excluded are $\CPMF$, $\BeMF$ and, in the case of the Netflix dataset only, $\isup$. Many of the other uncertainty estimates harm precision as uncertain recommendations are excluded. These results are similar to those we observed with URI in the previous section. An exception is $\OrdRec$, which had high URI but did not perform well here. We also see that the decrease in MAP is, in general, related to the decrease in average relevance as stricter thresholds are applied. That is, the five models with the largest average relevance decrease ($\ensemble$, $\FunkSVD-CV$, $\resample$, $\ivar$ and $\OrdRec$), are the same ones with the largest drops in MAP. On the other hand, with some uncertainty metrics, such as $\CPMF$ and $\BeMF$, it is possible to apply cuts which reduce average relevance only slightly and increase MAP.

Empirically, \acrshort{ubf} does not harm coverage. For most of the uncertainty estimates, we can still make 10 recommendations to each user, even as filtering is made stricter. The exceptions are $\Bias-CV$ and $\OrdRec$ (and the latter only for the MovieLens dataset), and then, even in the worst case ($\Bias-CV$ at the 80\% quantile cut-off), the model is still able to provide 80\% of the requested recommendations. Therefore, the MAP@10 values for users whose coverage is unharmed are very similar to the original ones.

\subsubsection{Probability-of-relevance ranking}
\label{sec:results_prr}

In \acrfull{prr}, we generate recommendations by\linebreak ranking items according to the probability that their rating is greater than or equal to the relevance threshold, $\omega$, which for these datasets is 4. As we discussed in Section \ref{sec:prr}, this idea was proposed for $\BeMF$. But it applies directly also to $\CPMF$ and $\OrdRec$. The remaining (non-probabilistic) models, appear to be unsuitable to \acrshort{prr}. Nevertheless, in order not to exclude all of them from the comparisons, we can use some heuristic extensions. The $\ensemble$ trains several models, which can be thought of as an $N$-dimensional sample of the theoretical population of the rating estimator under all the possible random initialisations. This setup is reminiscent of the central limit theorem. Then, assuming normality, we can use:
\begin{equation}
    \prob(r_{ui} \geq \omega) = \prob\left(\norm\left(\frac{1}{N}\sum_{k=1}^N F^{(k)}(u, i), \frac{G_{ui}}{\sqrt{N}}\right) \geq \omega\right)
\end{equation}

Furthermore, the error-based models predict rating prediction errors, that is $G_{ui} = E[|\hat{r}_{ui} - r_{ui}|]$. Therefore, $G_{ui}$ is  an estimate of the ratings' standard deviation. Now, assuming $r_{ui}$ to follow a Normal distribution with average $F(u, i)$ and standard deviation $G_{ui}$, it becomes possible to estimate $\prob(r_{ui} \geq \omega)$ as:

\begin{equation}
    \prob(r_{ui} \geq \omega) = \prob\left(\norm\left(F(u, i), G_{ui}\right) \geq \omega\right)
\end{equation}

With such extensions, we can now compare $\ensemble$, $\FunkSVD-CV$,\linebreak $\Bias-CV$, $\CPMF$, $\BeMF$ and $\OrdRec$ for the task of \acrshort{prr}. As discussed in Section \ref{sec:eval_prr}, we evaluate whether \acrshort{prr}, which uses the uncertainty estimates for ranking, produces more accurate recommendation lists than \acrshort{rabr}, which ignores the uncertainty estimates. This comparison is shown using MAP@10 in Figure \ref{fig:unc_rank}.
\begin{figure}[tp]
\centering
\begin{subfigure}{\linewidth}
    \centering
    \includegraphics[width=\linewidth]{img/ML/unc_rank.pdf} 
    \caption{MovieLens dataset.}
\end{subfigure}
\begin{subfigure}{\linewidth}
    \centering
    \includegraphics[width=\linewidth]{img/Netflix/unc_rank.pdf} 
    \caption{Netflix dataset.}
\end{subfigure}
\caption{\acrshort{prr} versus \acrshort{rabr}: MAP@10.}
\label{fig:unc_rank}
\end{figure}

In the figure, we see that only $\CPMF$ and $\BeMF$ have higher MAP@10 when their recommendations are selected using \acrshort{prr} compared with selection by \acrshort{rabr}. $\OrdRec$'s MAP@10 is the same for both. For all the other uncertainty estimates, MAP@10 is harmed. This may be because the way that we extended \acrshort{prr} (from $\BeMF$) to these other models is not well-motivated enough, but it also reflects the poor performance of these models.

Finally, as discussed in Section \ref{sec:prr}, the quantity $1 - \prob(r_{ui} \geq \omega)$ can be considered as an uncertainty estimate. To evaluate the usefulness of these estimates, we employ UAC, remembering that URI is not applicable in this case due to the correlation between the ranking probabilities and the uncertainties. The results are shown in Table \ref{tab:prr_uncertainty}. The best performances are achieved by $\BeMF$ in the MovieLens dataset, and $\OrdRec$ in the Netflix dataset. For $\CPMF$, this uncertainty estimate is uncorrelated with MAP. Unfortunately, for $\ensemble$, $\FunkSVD-CV$ and $\Bias-CV$, the uncertainties are positively correlated with MAP.

\begin{table}[tp]
    \centering
    \begin{tabular}{|l|r|r|} \hline
    Model & MovieLens & Netflix \\ \hline
    $\ensemble$   &        0.2979 & 0.2787 \\ \hline
    $\FunkSVD-CV$ &       -0.0193 & 0.2296 \\ \hline
    $\Bias-CV$    &        0.0075 & 0.1180 \\ \hline
    $\CPMF$       &       -0.0334 & -0.0140 \\ \hline
    $\BeMF$       &       \textbf{-0.2077} & -0.1509 \\ \hline
    $\OrdRec$     &       -0.0237 & \textbf{-0.2048} \\ \hline
    \end{tabular}
    \caption{\acrshort{prr} uncertainty: UAC on MovieLens and Netflix datasets. The best values (lowest) are highlighted in bold.}
    \label{tab:prr_uncertainty}
\end{table}

This concludes our results. To summarise, we have shown that there is no model that performs better than all the others across all the metrics that we have explored. Nevertheless, we have shown $\CPMF$ to be the best model for both use cases (i.e.\ for \acrlong{ubf} and \acrlong{prr}). On the other hand, our new $\Bias-CV$ model had the best uncertainty estimates for most of our rating prediction results, while performing very poorly on the top-$n$. 
This highlights a disconnect between the rating prediction task and the top-$N$ recommendation task in the context of uncertainty estimation. Moreover, we saw that the uncertainty estimates behave very differently from each other. This, together with the variations in the performance of the uncertainty estimates across different metrics, raises a major concern about the credibility of some of the uncertainty estimates. Finally, given that top-$n$ recommendations is more important than rating prediction, we believe $\CPMF$ to be the most promising of all the models.

\section{Discussion}
\label{sec:discussion}

So far, we described and empirically evaluated several ways of estimating rating prediction uncertainty. The results show the performance of each uncertainty estimate according to different evaluation criteria. Although we have already raised some general concerns, the results allows us to conduct a further discussion regarding each uncertainty estimate. Furthermore, there are aspects that may be particular to each estimate that may not be explicitly demonstrated by the results. This section discusses these issues and critiques each of the uncertainty estimates.

\begin{description}
    \item [$\isup$:] Characterising the uncertainty of a rating prediction by the number of ratings of the target item has a solid foundation. Its high URI confirms that, when recommended, well-supported items tend to be more relevant. However, its strong performance on this metric might be attributable to the various biases known to exist in the datasets used for offline RS evaluation \cite{Chaney:2018:ACR:3240323.3240370}.   Furthermore, as we have already noted, one concern about $\isup$ is that this estimate of uncertainty is not personalised, which does not seem reasonable; we would expect uncertainty to depend on the relationship between an item and the user profile. Another unreasonable property is its model-independence; we would expect uncertainty to depend in part on the quality of the rating prediction model. Another potential issue is that, if used for \acrshort{ubf}, $\isup$ will reduce the item candidates to a subset of the most popular ones, greatly decreasing the system's ability to recommend novel, yet relevant, items.
    
    \item [$\ivar$:] Disagreement between users would certainly seem to be a source of uncertainty. In line with this, the positive results achieved by $\ivar$ in the rating prediction uncertainty evaluation show that the predicted ratings for items about which users disagree tend to be less accurate. On the other hand, we did not find $\ivar$ to be useful in the top-$n$ recommendation case. In addition to this, and similarly to $\isup$, $\ivar$ is also model-independent and not a personalised metric.
    
    \item [$\resample$:] We found that the variation in predictions that ensues from perturbations during model training does correlate with the error in the final predictions, as evidenced by the positive values achieved by $\resample$ in the rating prediction uncertainty evaluation metrics. On the other hand, this method obtained negative URI, meaning that, in the task of top-$n$ recommendation, quantifying uncertainty by this method is not only ineffective but actually harmful. One possible explanation for this is that recommended items are the ones with the highest ratings, which means they fall far from the average rating case and are therefore more likely to get noisy predictions.
   
    \item [$\ensemble$:]  We found $\ensemble$ to have similar, or even worse, results than $\resample$ in regards to uncertainty evaluation. In fact, given the similar motivation of both, we can conclude that multi-modelling strategies are not useful for uncertainty quantification in the case of top-$n$ recommendation. Nevertheless, we highlight that, in the MovieLens dataset, $\ensemble$ had the best MAP and Recall in the \acrshort{rabr} task, although the improvement relative to the baseline is small.
    
    \item [$\FunkSVD-CV$:] We expect uncertainty to be a quantity related to rating prediction error. The error-based methods build on this idea directly: they use regression to predict error, and treat these predicted errors as estimates of uncertainty. Unsurprisingly then, $\FunkSVD-CV$ achieves good results in the rating prediction uncertainty evaluation. Unfortunately, when it comes to top-$n$ recommendations, it performs poorly. This may be because $\FunkSVD-CV$ minimises the average absolute error, but the recommended items are those with the highest ratings, which are outliers relative to the average, and therefore, might not receive good uncertainty predictions.
    
    \item [$\Bias-CV$:] Our critique of $\FunkSVD-CV$ also applies to $\Bias-CV$. In this case, the results are even more evident: the rating prediction uncertainty evaluation results are even better than those of $\FunkSVD-CV$, but its URI results are even worse. Therefore, we conclude that the simple linear model can accurately recognize the user-item interactions where rating prediction imprecision is most likely to occur, but cannot recognize which recommendations are more or less likely to be correct. To solve this issue, one possibility to be explored in future work is to restrict the error estimation learning only to those interactions which get predicted with the highest ratings.
        
    \item [$\CPMF$:] $\CPMF$ extends the popular PMF model by incorporating variance parameters. Whereas other ways of measuring uncertainty are largely heuristic in nature, $\CPMF$ has a great theoretical foundation. Not only that, but our experiments show $\CPMF$ to be a successful approach, being the only model able to benefit from the uncertainty-aware recommendation methods (\acrshort{ubf} and \acrshort{prr}). Nevertheless, its performance in the, admittedly less important, rating prediction uncertainty evaluation was not impressive. Of course, this may not be a weakness of $\CPMF$; it may, instead, raise a flag about the disconnect between prediction error experiments and top-$n$ recommendation experiments. Of the two, a good performance at top-$n$ recommendation is more important than performance in rating prediction.
    
    \item [$\BeMF$:] In its original formulation, this model has several particularities that contrast with the others herein. For instance, the rating prediction function of $\BeMF$ differs from those of $\OrdRec$ and $\CPMF$ by focusing on the mode of the rating distributions instead of their average and dispersion. We are grateful to the authors of this method, since their work greatly helped our definitions of \acrshort{ubf} and \acrshort{prr}. Nevertheless, the results achieved by their model have been largely underwhelming; although the uncertainty evaluation metrics have been positive in both the rating prediction and the recommendation case, the poor MAP obtained in both \acrshort{rabr} and \acrshort{prr} makes it really hard to justify its use compared, for example, with $\CPMF$, especially because $\BeMF$ has much higher complexity.
    
     \item [$\OrdRec$:] $\OrdRec$ offers a collaborative filtering algorithm that is suitable to ordinal feedback data. Although this was its main goal, its authors recognized that its multinomial formulation allows the model to also provide uncertainty estimates. In Section \ref{sec:multinomial_based}, we raised a concern regarding $\OrdRec$'s uncertainty estimates at the limits of the rating scale, which may cause issues. In particular, recommended items are usually those whose predicted ratings fall close to the rating scale upper boundary, and therefore, these items will often receive very low uncertainty estimates, which might not only not reflect reality, but also harms some of the evaluation metrics employed herein (e.g.\ URI). Furthermore, we have empirically shown that \acrshort{rabr} and \acrshort{prr} produce similar results with $\OrdRec$. In light of this, we believe that these uncertainty estimates add small value when compared to other models.
\end{description}

\section{Conclusions and Limitations}
\label{sec:concs}

Uncertainty in \acrshort{ers} is a very important, and yet quite under-explored, topic. In this chapter, we have shown how it is possible to quantify the uncertainty in rating predictions and top-$n$ recommendations using several different methods. Our main findings are:

\begin{itemize}
    \item We found that all the explored methods are successful at estimating uncertainty in rating predictions, at least to some degree (see Tables \ref{tab:rating movielens} and \ref{tab:rating netflix}). Among the methods, we found $\Bias-CV$ to be the most precise.
    \item  Their success can motivate the use of uncertainty estimates in a variety of tasks in the \acrshort{rs} field, including active learning, co-training and reinforcement learning.
    \item We have also shown that it is possible to improve top-$n$ recommendation accuracy with what methods that we collective refer to as uncertainty-aware ranking. That is, methods that incorporate uncertainty estimates in the selection of a top-$n$. We showed $\CPMF$ to be the best contender for both \acrfull{ubf} and \acrfull{prr}.
    \item In addition, we found that some of the other models, e.g.\ $\FunkSVD-CV$ and $\Bias-CV$, have a poor performance on estimating the uncertainty of top-$n$ recommendations, which makes it hard to defend their use in practical systems, despite their success on rating prediction uncertainty estimation.
\end{itemize} 

Even though the collection of uncertainty estimators explored in this work is the largest yet and the empirical analysis is substantial, we recognize the empirical analysis has some limitations. In Chapter \ref{chapter:conclusion}, we discuss ideas for future work, some of which follow from the following limitations. 

\begin{itemize}
    \item First, both datasets that we have used come from the same domain (movies). Movie recommendation has historically been one of the most studied applications of recommender systems, but there are many others, such as shopping, social media, etc. While we expect the behaviour of the uncertainty estimators to generalize to other domains, there is clear value in exploring how different rating patterns will affect the uncertainty estimates.
    \item Second, our empirical analysis uses only a single train/test split. This is a consequence of our decision to opt for chronological train/test splits. We chose, in other words, to preserve the temporal nature of the data, which is intrinsic to recommender systems. This decision does, however, mean that our results lack the kind of statistical significance analysis that could be obtained from using multiple train/test splits. While we expect our results to be robust, due to the large size of our datasets, there is clear value to exploring this further.
    \item Finally, in our empirical analysis, for those uncertainty estimation methods that need a separate rating predictor, we used FunkSVD (Eq.\ \ref{eq:MF}) with mean-squared rating prediction error as its loss function (Eq.\ \ref{eq:FunkSVDloss}). There may be value in exploring other loss functions (e.g.\ ranking-based losses). We believe our choice was appropriate because it gave a reasonably fair comparison across the systems; the exploration of uncertainty does not rely on finding recommenders with the very highest accuracy; and, in any case, ranking-based loss becomes more relevant when dealing with implicit feedback data, which is the main focus of the next two chapters.
\end{itemize}
