Estimate covariance given noisy data and the mean

Question

Basically, I'm inferring the parameters of a Gaussian, $\mu$ and $\Sigma$ given observed data $y_i$ that have uncertainties $\sigma^2_i$ associated with them.

1D example:

Prior

I intend to use the conjugate prior to the multivariate Gaussian:

• A Gaussian prior on the mean $\sim \mathcal{N}(\lambda, r)$

• A Wishart prior on the covariance $\sim \mathcal{W}(\nu+n, V^{-1})$

So thats four ($\lambda, r, \nu, V$) hyperparameters for the prior.

Distribution of the mean given noisy data

I can derive the distribution of the mean, $\mu$ given the noisy data, $\{y_i,\sigma^2_i\}$, which is Gaussian (from the product of gaussian likelihoods):

$P(\mu|y_i,\sigma_i^2,\Sigma\,\lambda,r) = \mathcal{N}(\mu|\lambda,r)\prod_i \mathcal{N}(y_i|\mu,\Sigma+\sigma^2_i) \propto \mathcal{N}(\mu|m_c, s_c)$

$s_c = [(r+\Sigma)^{-1} + \sum_i (\sigma_i^2 + \Sigma)^{-1}]^{-1}$

$m_c = s_c [\lambda (r+\Sigma)^{-1} + \sum_i (\sigma_i^2+\Sigma)^{-1}y_i]$

So I can sample means easily when conditioned on the rest of the parameters.

Distribution of the covariance given noisy data

Now I need $P(\Sigma|y_i,\sigma_i^2,\mu\,\nu, V) = \mathcal{W}(\nu+n, V^{-1})\prod_i\mathcal{N}(y_i|\mu,\Sigma+\sigma^2_i)$

This one is not trival as far as I can tell.

1D version

From McLawrence's Question - "estimate-the-variance-of-gaussian-distribution-from-noisy-sample", where they focus on one shared noise variance $\sigma_e$ in 1D, the multivariate version will look like the multivariate version of "truncated" inverse gamma distribution (generalises to the truncated inverse Wishart?).

Version with no data uncertainty

Furthermore, if we say that there are no observational errors ($\sigma_i=0$), then we can just use the conjugate expression for Multivariate Normal with known mean found here. This would be an inverse-wishart distribution.

Conclusion

But, as I say my data are heteroscedastic and so the uncertainties are not the same, making both of the above solutions invalid. We cannot take $(y - \mu)^T(y-\mu)$ out of the equation as a constant since there is an extra term $\sigma_i^2$ which needs to be included.

Since the distribution for the mean can be readily sampled from, I was wondering if the same can be said for the covariance. Is there a sampling distribution for $P(\Sigma|y_i,\sigma_i^2,\mu\,\nu, V) = \mathcal{W}(\nu+n, V^{-1})\prod_i\mathcal{N}(y_i|\mu,\Sigma+\sigma^2_i)$?

For a 1d version this problem I can plot this conditional distribution. It looks like we should be able to sample from it...

In the above plots, the mean is known.

Obviously, this method of plotting cannot be extended to higher dimensions.

Answer 1

If your interest is intended for possible applications, I would recommend this work 'Scalable Robust Covariance and Correlation Estimates for Data Mining'. To quote the paper's extract:

"Covariance and correlation estimates have important applications in data mining. In the presence of outliers, classical estimates of covariance and correlation matrices are not reliable. A small fraction of outliers, in some cases even a single outlier, can distort the classical covariance and correlation estimates making them virtually useless. That is, correlations for the vast majority of the data can be very erroneously reported; principal components transformations can be misleading; and multidimensional outlier detection via Mahalanobis distances can fail to detect outliers. There is plenty of statistical literature on robust covariance and correlation matrix estimates with an emphasis on affine equivariant estimators that possess high breakdown points and small worst case biases. All such estimators have unacceptable exponential complexity in the number of variables and quadratic complexity in the number of observations. In this paper we focus on several variants of robust covariance and correlation matrix estimates with quadratic complexity in the number of variables and linear complexity in the number of observations. These estimators are based on several forms of pairwise robust covariance and correlation estimates. The estimators studied include two fast estimators based on coordinate-wise robust transformations embedded in an overall procedure recently proposed by [14]. We show that the estimators have attractive robustness properties, and give an example that uses one of the estimators in the new Insightful Miner data mining product."