Flexmix maxima are not where they are expected to be

Question

For my dataset I have plotted the density with ggplot. As the data's density is multimodal (a total of 6 destinct modi) I tried to gain insight on the normal distributions associated to each modus.

For this I followed the instructions on JEFworks and adapted the corresponding arguments to my code and get the following plot for the normal distribution for the first modus:

As you can see the normal distribution for the first modus is located between the second and third modus.

In my understanding, as the density graph can be assumed as

$$g(x) = \sum_{i=1}^n f_i(x)$$

with $g(x)$ being my density plot, n being the number of modi of my density distribution as well as $f_i(x)$ being the normal distribution underlying each respective modus.

Therefore shouldn't be the x-position of each normal distributions mean be equal to the respective maximum of my density distribution (and therefore the respective modus?

For another dataset this is nearly true for the first modus, but the graph looks quite "jiggely" around this area:

Is there any possibilty to give flexmix the modus as an input that the model can better know where the mean should be or is my assumption above wrong to begin with?

Edit: I ran my code multiple times, each second time the output is as in the first image, the other times it is comparable to the second image:

Is there any way to check why this happens each second time I run my code?

Answer 1

I'm afraid you got several things wrong.

• If I understand correctly, you are fitting a mixture of Gaussians to your data. The mixture will approximate the distribution of your data. It won't necessarily perfectly fit the empirical distribution, because it is just an approximation of the distribution.

• This part is also incorrect

  In my understanding, as the density graph can be assumed as $$g(x) = \sum_{i=1}^n f_i(x)$$

  Mixture distribution is defined as

  $$ g(x) = \sum_{i=1}^n \pi_i \, f_i(x) $$

  where $\pi_i$ are the mixing proportions such that $\sum_i \pi_i = 1$ and $\forall \pi_i \ge 0$.

• As for your "jiggely" plot, this is just a plotting artifact. If you used a denser grid to produce the plot, it would be smooth.

• Answering your last question, you definitely can fix the means of the mixture components and estimate only the mixing weights $\pi_i$ and variances $\sigma^2_i$. I'm not sure though if the flexmix library allows for this, you should consult the documentation.

• Finally, you say that the "maxima are not where they should be", but what you are doing is you are comparing one approximation of the distribution (kernel density) with another approximation of the distribution (Gaussian mixture). Both approximations can give different answers. Either or none of those modes can be the correct estimate. So using the modes from your first approximation to fix them in the second approximation is not guaranteed to give a better answer than sticking to the regular Gaussian mixture where the means are estimated from the data.

Answer 2

Adding to @Tim's answer, here are some more issues:

• In mixture distributions the numbering of components is meaningless, so there is no guarantee whatsoever that the first estimated component corresponds to the first estimated mode. If indeed any Gaussian corresponds to the first mode, it can be in any "place" in the mixture.
• Even if your dataset has 6 modes and you fit a mixture of 6 Gaussians, there is no guarantee that every Gaussian corresponds to a mode. It may well happen that if you align all Gaussians with modes, you cannot have a good approximation of the tails (or some other parts of the distribution), and it may be better to include a larger variance Gaussian and only have five Gaussians "corresponding" to modes (there may be quite weak modes that may only appear due to random variation), or even some other "configuration". I suspect that the large variance Gaussian in the first plot was fitted for such a reason.
• Generally a Gaussian mixture with six components is quite flexible (with means, variances, and mixture proportions you have 6+6+(6-1)=17 free parameters), and this may mean that there are several different ways to approximate your data about equally well. Maximum likelihood may have several locally optimal solutions the parameters of which can look quite different, even though one should expect that the resulting mixture densities are about equally similar to your data. If you repeat running flexmix, you may find different such solutions as flexmix depends on random initialisation if I'm not mistaken. Normally you would want to use the solution that gives the highest likelihood. You can get flexmix to give you a more stable solution that will often have the best likelihood you can get by increasing the number of initialisations of the algorithm. This is done by parameter nrep in the FLXcontrol list; you can provide a list with nrep=500, say , as control argument to the flexmix call. You need to balance good stability through a high value of nrep against computational effort here. Unfortunately the help pages don't seem to say what the default value of nrep is, but I suspect it is too low (maybe 1 or 10), and therefore you easily get different solutions in different runs.
• Even if Gaussian mixture components actually can be interpreted as being associated to your modes, you shouldn't expect them to be exactly equal to the modes, as mixing several Gaussians together, meaning that in every place of the density contributions from all mixture components are present, will produce modes that are not equal to the means of individual Gaussians (even if maybe often quite close). However, if you run flexmix with large enough nrep and you get one or more mixture components that you can't associate with your modes, I'd presume that from the point of view of overall density fitting, this gives you a better fit than nailing all Gaussians to the modes (see above).