I have some data represented by vectors $\mathbf{x}_1,\mathbf{x}_2,\ldots,\mathbf{x}_n \in \mathbb{R}^m$, that is, each sample is a collection of $m$ features. In addition, I have a response corresponding to each data point: $y_1,\ldots,y_n \in \mathbb{R}$.

I would like to train a parametric model on this data, in order to be able to predict the probability density of the response $y(\mathbf{x})$ corresponding to a new sample $\mathbf{x} \in \mathbb{R}^m$. (I really need to have a probability density for $y(\mathbf{x})$, not just a punctual estimation.)

So my idea would be to select a family of parametric densities defined through a parameter $\boldsymbol{\theta}\in\mathbb{R}^k$, $f(y;\boldsymbol{\theta})$. Then, I would like to run a maximum-likelihood procedure to find the best possible function $\theta$ from a class of functions $\Theta$ mapping $\mathbb{R}^m\mapsto\mathbb{R}^k$: $$ \max_{\theta \in \Theta} \sum_{i=1}^n \log f(y_i,\theta(\mathbf{x}_i)). $$

Question 1 Is this a well-studied problem, and if yes, how is it called in the literature? Or is there an alternative approach that is known to work better?

Question 2 My first, naive idea, is to restrict myself to a class $\Theta$ of linear functions, so the above problem reduces to determining $km$ coefficients. In this case, any reference to some algorithm performing this task would be appreciated.

Question 3 If we consider a family of gaussians whose parameters $\mu$ and $\sigma$ are linear functions of $\mathbf{x}$, i.e. $y \sim \mathcal{N}(\boldsymbol{\mu}^T\mathbf{x}, (\boldsymbol{\sigma}^T\mathbf{x})^2)$ for some vectors $\boldsymbol{\mu},\boldsymbol{\sigma}\in \mathbb{R}^m$, then the optimization problem to solve is $$ \min_{\boldsymbol{\mu},\boldsymbol{\sigma}\in \mathbb{R}^m}\quad \sum_{i=1}^n \frac{1}{2}\left(\frac{y_i-\boldsymbol{\mu}^T\mathbf{x}_i}{ \boldsymbol{\sigma}^T\mathbf{x}_i} \right)^2 + \log(\boldsymbol{\sigma}^T\mathbf{x}_i). $$ What is the state-of-the-art to solve this problem ?

Or maybe, because we want small predicted variances, we could add a penalization term: $+ \lambda \sum_i (\boldsymbol{\sigma}^T\mathbf{x}_i )^2$ in the function to minimize?

Question 4 If we assume that $f$ is a mixture of laws $f_1,\ldots,f_N$ defined by a vector of probabilities $\boldsymbol{\pi}$, is there a way to construct a function $\pi$ that maps $\mathbf{x}$ to the probability vector $\boldsymbol{\pi}=\pi(\mathbf{x})$, so that the density of the response $y(\mathbf{x})$ is given by $ \sum_{j=1}^N \pi_j(\mathbf{x}) f_j(y). $ I have the feeling that using a neural network might be appropriate to do construct the function $\pi$, but again, I would be grateful if someone could point out a good reference for this specific problem.

Question 5

Here is a tentative approach to answer Question 4. How appropriate do you think this is? Is this a standard approach (or a stupid one)? Do you know some references?

We could use the EM algorithm to model the augmented vectors $(\mathbf{x}_i,y_i)\in\mathbb{R}^{m+1}$ as coming from a mixture of multivariate gaussian variables: $$\left(\left[\begin{array}{c}\mathbf{x}_i\\y_i\end{array}\right] | Z=k\right) \sim \mathcal{N}\left( \left[\begin{array}{c}\boldsymbol{\mu}_k\\ m_k\end{array}\right], \left[\begin{array}{cc}\Sigma_k & \boldsymbol{v}_k\\ \boldsymbol{v}_k^T & \sigma_k^2\end{array}\right] \right), \quad\textrm{with}\ P(Z=k) = w_k. $$

Then, the distribution of the response $y$ given the value of a new sample $\mathbf{x}$ and the value of the latent variable $Z$ would be:

$$(y|\mathbf{x},Z=k) \sim \mathcal{N}\left( m_k + \boldsymbol{v}_k^T \Sigma_k^{-1} (\mathbf{x}-\boldsymbol{\mu}_k),\ \sigma_k^2-\boldsymbol{v}_k^T \Sigma_k^{-1}\boldsymbol{v}_k \right). $$

Finally, we obtain the distribution of $y|\mathbf{x}$ as a mixture of gaussians, with weights $P(Z=k|\mathbf{x})$ that can be obtained from the Bayes formula.


So, after some weeks I have come up with some answers to some of my above questions.

Question 1 It turns out that this problem has been studied a lot within the framework of Gaussian processes. A good reference seems to be "Machine Learning: A Probabilistic Perspective", from Kevin P. Murphy.

Basically, the idea is to assume that the responses $y_1,\ldots,y_n$ come from a multivariate normal distribution $\mathcal{N}(\boldsymbol{\mu}, \Sigma)$, whose mean and variance depend from the data $X=(\mathbf{x}_1,\ldots,\mathbf{x}_n)$. A typical approach is to assume a particular form for the functions $\mu()$ and $\Sigma()$, for example the $\mu_i$'s can be equal to some affine function applied to $\mathbf{x}_i$, and $\Sigma_{ij}$ can be a combination of some squared exponential kernel $K(\mathbf{x}_i,\mathbf{x}_j)$ and a diagonal noise $\sigma_{\rm{noise}}^2 \delta_{ij}$. These functions depend from a set of hyperparameters $\boldsymbol{\theta}$, which can be optimized by maximum-likelihood. Then, if we try to predict the answer $y^*$ of a new sample $\mathbf{x}^*$, we can use the formula that gives the posterior $y^* | (y_1,\ldots,y_n)$: $$ y^* | (y_1,\ldots,y_n) \sim \mathcal{N}\left(\mu(\mathbf{x}^*) + \mathbf{v_*}^T \Sigma(X)^{-1} (\boldsymbol{y}-\boldsymbol{\mu}(X)),\quad \sigma_*^2 - \mathbf{v_*}^T \Sigma(X)^{-1} \mathbf{v_*} \right), $$ where $\sigma_*^2:=\Sigma(\mathbf{x}^*,\mathbf{x}^*)$ is the prior variance at $\mathbf{x}^*$ and $\mathbf{v_*}$ is the vector of prior cross-covariances between $\mathbf{x}^*$ and the data: $(\mathbf{v_*})_i = \Sigma(\mathbf{x}^*,\mathbf{x}_i)$.

Question 3 For my particular application, I can indeed solve this maximum likelihood problem to find the vectors $\boldsymbol{\mu}$ and $\boldsymbol{\sigma}$. Then, I can use these optimal vectors of coefficients to construct a starting point for the (much harder) problem of optimizing the hyper-parameters of a gaussian process, cf. my answer to Question 1. I found that the following strategy works well: solve alternately the optimization problem with respect to $\boldsymbol{\mu}$ and $\boldsymbol{\sigma}$. The optimization in $\boldsymbol{\mu}$ is a standard least square problem for a fixed $\boldsymbol{\sigma}$, and the optimization in $\boldsymbol{\sigma}$ can be done by any gradient-based algorithm, e.g. the BFGS algo.

  • $\begingroup$ What I proposed in Question 5 seems to be close to the solution provided by Gaussian processes (cf. my answer to question 1), but uses a mixture of Gaussian rather than a Kernel, and works in the space of extended vectors $(\mathbf{x}_i,y_i)$. I wonder if the link between these two approaches has already been studied ? $\endgroup$
    – guigux
    Oct 5 '15 at 11:02

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.