# How to extract components of Gaussian mixture?

I'm trying to model a dataset as a mixture of two Gaussian distributions in MATLAB and find the Bhattacharyya distance between the two. Using MATLAB's fitgmdist function I was able to model this mixture and produce this plot:

gmdist = fitgmdist(data, 2);
gmsigma = gmdist.Sigma;
gmmu = gmdist.mu;
gmwt = gmdist.ComponentProportion;

histogram(data, 'Normalization', 'pdf', 'EdgeColor', 'none')
x = min_val:0.0001:max_val;
xlim([min_val max_val])
hold on;
plot(x, pdf(gmdist, x'), 'k')
hold on;


However, when I was debugging my distance code I realized that my two individual distributions did not match their components in the mixture.

p = pdf('Normal', x, gmmu(1), gmsigma(1));
q = pdf('Normal', x, gmmu(2), gmsigma(2));
plot(x, p*gmwt(1))
hold on;
plot(x, q*gmwt(2))


My expectation was that the above code should have produced two density plots that matched their components in the plot of the original mixture model, but this is not the case. This post clued me in that the PDFs I calculated integrated to 1 individually (rather than together), but I'm unsure how to obtain the "non-integrated" components of the mixture model such that I can plot them to match the original plot.

• It is not clear how they should match. – Michael R. Chernick Feb 19 '18 at 4:24
• your plot does not take into account the weights as they should since the sum of the two weighted Gaussians should be the mixture – Xi'an Feb 19 '18 at 7:18
• @Xi'an: do you want to post your comment as an answer? – Stephan Kolassa Feb 19 '18 at 7:24
• @StephanKolassa: thanks but I am not sure about turning this into an answer as this may be an obvious Matlab issue, like plot using different scales for the three graphs. (In R I would state add=TRUE in the plot command.) – Xi'an Feb 19 '18 at 7:34
• @Xi'an: I don't think you need to go into details on how to create the correct Matlab code. Simply pointing out that the second density uses an unweighted mixture would already be potentially helpful for future generations. – Stephan Kolassa Feb 19 '18 at 7:36

There is definitely a mistake in the implementation of the decomposition. Here is my rendering of the same problem with both weighted components appearing as they should:

My R code is as follows:

I solved the problem by creating my PDFs with the square roots of gmsigma. It seems I fundamentally misunderstand covariance, variance, and standard deviation so I will look into that for the future.
p = pdf('Normal', x, gmmu(1), gmsigma(1)^0.5);

• Output gmsigma so you can see it in its entirety. Then you can make sure gmsigma(2) is really what you want - I'm not sure it is,. To summarize the preceding, if you have a multivariate normal, then in MATLAB function conventions, Sigma = sigma^2. Hello confusion. But I guess MATLAB developers were influenced by the multivariate normal convention to use $\Sigma$ for covariamce matrix. then chose to use that as input for multivariate, while choosing sigma, rather sigma^2, to specify univariate normal. – Mark L. Stone Feb 19 '18 at 15:58