# Closed Form Solution for Gaussian Matrix which is Convex Combination?

I already asked a pretty similar question here, but the answer was inconclusive and now this related problem has come up again here.

My problem is as follows, I have a $$2n$$-dimensional multivariate normal distribution with the following variance-covariance matrix: $$\begin{bmatrix} \sigma_{x}^{2}[(\lambda)A_{n}+(1-\lambda)\mathbb{I}_{n}] & 0_{n x n} \\ 0_{n x n} & \sigma_{y}^{2}[(\lambda)A_{n}+(1-\lambda)\mathbb{I}_{n}] \\ \end{bmatrix},$$ where $$\lambda$$ is unknown in $$[0,1]$$ [or if necessary, it can be made in $$(0,1)$$]; $$\sigma_{x}^{2}$$ and $$\sigma_{y}^{2}$$ are unknown variance factors greater than 0; $$A_{n}$$ is an $$(n x n)$$ known, positive semi-definite matrix with a diagonal of $$1$$s; $$\mathbb{I}_{n}$$ is the n-dimensional identity matrix.

Additionally, its mean is the concatenation of two n-length vectors of identical values $$\mu_{x}$$ and $$\mu_{y}$$.

Is there a closed form solution for any of the maximum likelihood estimates of $$\lambda$$, $$\sigma_{x}^{2}$$, $$\sigma_{y}^{2}$$, $$\mu_{x}$$, or $$\mu_{y}$$?

EDIT: Xi'an posted a comment saying this is unrelated to Gaussian mixtures, but isn't my problem fitting a linear combination of two Gaussian distributions? Here is a reformulation which I think highlights this: $$(\lambda)\begin{bmatrix} \sigma_{x}^{2}A_{n} & 0_{n x n} \\ 0_{n x n} & \sigma_{y}^{2}A_{n} \\ \end{bmatrix}+ (1-\lambda)\begin{bmatrix} \sigma_{x}^{2}\mathbb{I}_{n} & 0_{n x n} \\ 0_{n x n} & \sigma_{y}^{2}\mathbb{I}_{n} \\ \end{bmatrix}$$

EDIT 2: mlofton asks whether my covariance is fixed or a function of t. I think it's probably either depending on the exact question. $$A_{n}$$ represents the expected correlation between observations at different points in space. Specifically, the correlation between $$x_{i}$$ and $$x_{j}$$ can be found at the entry in the ith row and jth column. These correlation values are function the spatial distance. But, we know the spatial distances and function which converts them into correlations, so these values are all fixed. I describe the way these distances are measured and converted into correlations here and here. But basically the variables evolve along a network of paths, and the covariance between two observations is proportional to the amount of shared path they have. So I'd say that because we know the correlations, their is no longer any explicit dependency on time, distance, or anything else. As for the $$0_{nxn}$$ matrices, in some versions of this problem they can be replaced by $$\sigma^{2}c_{i}K_{n}$$ where $$c_{i}$$ is a correlation factor between $$X$$ and $$Y$$, $$K_{n}$$ is a fully known matrix which describes the structure of the correlation, and $$\sigma^{2}$$ can be calculated from $$\sigma_{x}^2$$ and $$\sigma_{y}^2$$.

EDIT 3: I'd like to double check the formula I found online. For the unbiased variance $$\sigma_{x}^2$$ I should evaluate $$\frac{1}{(n-1)} (x_{observed}-\mu_{x}1_{n})^{T}(\lambda A_{n}+(1-\lambda)\mathbb{I}_{n})^{-1}(x_{observed}-\mu_{x}1_{n})$$ but probably I need the mle estimate instead right? So I simply replace $$(n-1)$$ with $$(n)$$ in the bottom part of the first fraction?

And then to calculate $$\mu_{x}$$, I suppose the formula is this? $$\bigg(1_{n}^{T}\big((\lambda) A_{n}+(1-\lambda)\mathbb{I}_{n}\big)^{-1}1_{n}\bigg)^{-1}1_{n}^{T}\big((\lambda) A_{n}+(1-\lambda)\mathbb{I}_{n}\big)^{-1}x_{observed}$$

EDIT 4: I never verified the above formulas (in EDIT 3) because I found a black box solution for the mean and variance stuff. The dense grid search idea (proposed by Madden) is working really well though!

• This question is unrelated with Gaussian mixtures. Commented Jun 23 at 7:59
• Can you explain further please? I edited the question to show why I think it is a Gaussian mixture. Commented Jun 23 at 16:22
• Hi: Is your covariance fixed or is it supposed to be function of $t$. It looks like an exponentially smoothed covariance matrix which is always a function of $t$. Commented Jun 23 at 16:29
• I think the answer is that the matrix is fixed, not a function of time or distance, but that measurements of the times/distances between observations were used originally to construct the now fixed matrix. I edited the original question to give a longer answer (see edit #2). Commented Jun 23 at 17:16

Given $$\lambda$$, we can calculate $$\mu_x, \mu_y, \sigma_x, \sigma_y$$ using standard techniques, which you can find by searching "Weighted Least Squares" and adapting some of that math.

Unfortunately, it's not clear to me how optimization wrt $$\lambda$$ might be achieved in closed form. With a bit of work, we find that the log likelihood pertaining to $$\lambda$$ may be written as:

$$- \Bigg[\sum_{i=1}^N \log\big(\lambda(\gamma_i-1)+1\big)\Bigg] - \frac{1}{2} \Bigg[ \sum_{i=1}^N \frac{\frac{u_{x,i}^2}{\sigma_x^2}+\frac{u_{y,i}^2}{\sigma_y^2}}{\lambda(\gamma_i-1)+1} \Bigg] \, ,$$

where $$\gamma_i$$ is the $$i$$th eigenvalue of $$\mathbf{A}_N$$, $$u_{x,i} = \mathbf{v}_i^\top(\mathbf{x}-\mu_x)$$, and $$u_{y,i} = \mathbf{v}_i^\top(\mathbf{y}-\mu_x)$$, where $$\mathbf{v}_i$$ denotes the eigenvector of $$\mathbf{A}_N$$ corresponding to the $$i$$th eigenvalue.

To my knowledge, such an equation must be solved numerically. (Indeed, this is essentially the problem of determining a nugget parameter in Gaussian process regression, which is performed numerically in practice). So even if we knew everything else, we wouldn't be able to calculate $$\lambda$$ in closed form.

However, we can develop efficient numerics for this likelihood in practice. Given that you can compute every other parameter in closed form, it is simply a matter of varying $$\lambda$$ over a dense grid of $$[0,1]$$, computing the other parameters and the likelihood value at which the maximum is achieved, and taking that $$\lambda$$ which has the highest likelihood as our estimate. To do this efficiently, you will want to precompute the quantities I list above in order to use only one eigendecomposition of $$\mathbf{A}_N$$ to compute as many $$\lambda$$ evaluations as desired.

Why isn't this a normal mixture?

If we have two objects $$a$$ and $$b$$, and $$\lambda\in[0,1]$$, we call $$\lambda a + (1-\lambda) b$$ a convex combination.

What you have written is normal distribution with a covariance matrix which is a convex combination.

Gaussian mixtures also involve convex combinations, but it's the density itself which is a convex combination, not its covariance parameter.

So if we wanted to do a Gaussian mixture of the two normal densities $$N(\mathbf{x}; \mu_x, \sigma^2_x \mathbf{A}_N)$$ and $$N(\mathbf{x}; \mu_x, \sigma^2_x \mathbf{I})$$, that would yield the following density (up to proportionality):

$$f(\mathbf{x}) \propto \lambda \big| \sigma^2_x \mathbf{A}_N \big| ^{-1} e^{-\frac{1}{2} (\mathbf{x}-\mu_x)^\top (\sigma^2_x \mathbf{A}_N)^{-1}(\mathbf{x}-\mu_x)} + (1-\lambda)\big| \sigma^2_x \mathbf{I} \big| ^{-1} e^{-\frac{1}{2} (\mathbf{x}-\mu_x)^\top (\sigma^2_x \mathbf{I})^{-1}(\mathbf{x}-\mu_x)} \, .$$

This is distinct from the density you have proposed, which (ignoring $$\mathbf{y}$$) is given by:

$$f(\mathbf{x}) = \big| \sigma^2_x \big(\lambda\mathbf{A}_N+(1-\lambda)\mathbf{I}\big) \big| ^{-1} e^{-\frac{1}{2} (\mathbf{x}-\mu_x)^\top \big(\lambda\mathbf{A}_N+(1-\lambda)\mathbf{I}\big)^{-1}(\mathbf{x}-\mu_x)} \, .$$

• wow! thanks for this amazing answer. Just to clarify, I do not need to calculate $\sigma_{x}^2$ or $\sigma_{y}^2$ at all until after I've already found the maximum likelihood $\lambda$ value, right? Commented Jun 24 at 2:24
• @AFriendlyFish Oh, sorry, I errantly dropped the $\sigma$s from my expression; I have edited them in. $\lambda$ does indeed depend on them, and they would need to be estimated using MLE conditional on $\lambda$. Commented Jun 24 at 13:06
• thank you. I'm still a little stuck, I added "EDIT 3" at the bottom of the question to explain. Commented Jun 25 at 0:14
• @AFriendlyFish sorry for slow response; glad you found a solution, if it's important to your work to have the closed form expressions we can walk through it. Commented Jul 1 at 20:06
• thank you. i got it! I actually do have one unsolved-ish optim problem from a project 6 months back, but right now I don't have the energy to go back and look at that. I'll still email you now to establish contact though. Commented Jul 1 at 20:38