Estimate 2D covariance from 3D matrix where 3rd column contains probability density values

I have an nx3 matrix where the first 2 columns contain uniformly distributed random (x, y) points and the 3rd column contains pdf values evaluated at each point. The pdf values are computed from a normal distribution using zero mean and a specified covariance. Thus, the first 2 columns are independent while the 3rd column is dependent on the first 2. Here is code that generates a sample matrix:

n = 10000;
mu = [0, 0];
C = [12, 2; 2, 8]; % input covariance
X = [20 * (rand(n,1) - 0.5), 20 * (rand(n,1) - 0.5)]; % centers the grid at (0, 0)
X = [X, zeros(n,1)];
for k = 1:n
    p = mvnpdf(X(k,1:2), mu, C); % returns pdf value
    X(k,3) = p;
end

My question is how do I estimate the 2x2 covariance matrix that was used to generate my nx3 data? Given that I have the pdf values, it seems like this shouldn't be difficult to do, but I just can't figure out a simple solution.

Many thanks in advance!

Answer 1

[Sudarsanan A K]

Hi BC,

I understand that you are trying to estimate the 2x2 covariance matrix 'C' that you have used to generate the nx3 matrix 'X'. I also note that the first 2 columns are independent while the third column (which is dependent on the first two) is the PDF values computed from a normal distribution using zero mean and a specified covariance.

In this case, here is an approach for estimating the 2x2 covariance matrix 'C'.

% Estimate the covariance matrix

C_estimated = zeros(2, 2);

for k = 1:n

diff = [X(k,1)-mu_x; X(k,2)-mu_y];

C_estimated = C_estimated + (diff * diff') * X(k,3);

end

% Normalize the covariance estimate

C_estimated = C_estimated / sum(X(:,3));

The steps involved are:

1. Compute the 'difference vector': For each data point, you calculate the difference between the x-coordinate (X(k,1)) and the mean of x-coordinate (mu_x), and the difference between the y-coordinate (X(k,2)) and the mean of y-coordinate (mu_y). These differences are stored in the 'diff' vector, which represents the deviation of the data point from the mean.
2. 'Accumulate the covariance estimate': In this step, you estimate the covariance matrix by accumulating the outer product of the 'diff' vector. You multiply the 'diff' vector by its transpose (diff') to obtain a 2x2 matrix. This matrix represents the outer product of the differences, which captures the covariance between the x and y coordinates. You then multiply this matrix by the corresponding PDF value in the third column of 'X' (X(k,3)). The resulting matrix is accumulated in the 'C_estimated' matrix.

By iterating over each data point, you accumulate the contributions of the outer product of the difference vector, weighted by the corresponding PDF value. This process allows you to estimate the covariance matrix 'C' based on the provided nx3 data, where the first two columns contain uniformly distributed random (x, y) points, and the third column contains the PDF values evaluated at each point.

Finally, you normalize the accumulated 'C_estimated' matrix by dividing it by the sum of the pdf values in the third column of 'X '(sum(X(:,3))). This normalization step ensures that the estimated covariance matrix is scaled properly.

Hope you would find this useful.