Covariance calculation for mixture Gauss

Question

This equation what I found on the wikipedia is a bit strange for me.

How I can compute in matlab? If example the size(x) = [1000, 2](2 dimension gauss) then $size(\mu) = [1, 2]$ for each cluster.

I try with the following:

gaussEvaluation(i,:) = pInit(i) * mvnpdf(x,muInit(i,:), sigmaInit(i, :) * eye(d));
gaussEvaluationSum = sum(gaussEvaluation(i, :));

for j = 1 : n 
    v = (x(j, :) - muNew(i, :));
    sigmaNew(i) = sigmaNew(i) + gaussEvaluation(i,j) * (v * v');
end
sigmaNew(i) = sigmaNew(i) / gaussEvaluationSum;

OR

sub = bsxfun(@minus, x, mu(i,:));
sigma(i,:) = sum(gaussEvaluation(i,:) * (sub .* sub)) / gaussEvaluationSum;