# Mclust: Data frame order affects solution

I've come across some behavior in mclust::Mclust that I would not have expected, which is that the order of variables in the data passed to Mclust affects the solution it comes up with.

In the example below, the first ordering of the variables produces 2 profiles, and the second produces 6. However, calling either of these orders does produce a stable solution (e.g., the first call to Mclust will always produce 2 profiles).

> testData <- read.file(f = "http://fimi.ua.ac.be/data/chess.dat")

> summary(Mclust(subset(testData, select = c(X1, X3, X5, X7, X9, X11))))
----------------------------------------------------
Gaussian finite mixture model fitted by EM algorithm
----------------------------------------------------

Mclust EEV (ellipsoidal, equal volume and shape) model with 2 components:

log.likelihood    n df      BIC      ICL
3547.068 3195 49 6698.738 6692.126

Clustering table:
1    2
2759  436

> summary(Mclust(subset(testData, select = c(X11, X9, X1, X3, X5, X7))))
----------------------------------------------------
Gaussian finite mixture model fitted by EM algorithm
----------------------------------------------------

Mclust EEV (ellipsoidal, equal volume and shape) model with 6 components:

log.likelihood    n  df      BIC      ICL
18473.94 3195 137 35842.37 35834.51

Clustering table:
1   2   3   4   5   6
431 932 210 881 524 217


My suspicion is that the initialization step (that the documentation suggests is a hierarchical clustering procedure) is affected by the order of variable passed, but I'm not sure if there's a straightforward way to overcome this.

Update: It appears that there is something to my suspicion that the initialization step is involved in this behavior. If I load the same dataset as above, but instead pass Mclust an hcPairs list to its initialization argument that's based on an EII model, I get stable results across variable orders:

> testData <- read.file(f = "http://fimi.ua.ac.be/data/chess.dat")

> hc <- hc(subset(testData, select = c(X1, X3, X5, X7, X9, X11)), modelName = "EII")
> summary(Mclust(subset(testData, select = c(X1, X3, X5, X7, X9, X11)), initialization=list(hcPairs = hc)))
----------------------------------------------------
Gaussian finite mixture model fitted by EM algorithm
----------------------------------------------------

Mclust VEV (ellipsoidal, equal shape) model with 2 components:

log.likelihood    n df      BIC      ICL
3391.936 3195 50 6380.405 6380.254

Clustering table:
1    2
1386 1809
> summary(Mclust(subset(testData, select = c(X11, X9, X1, X3, X5, X7)), initialization=list(hcPairs = hc)))
----------------------------------------------------
Gaussian finite mixture model fitted by EM algorithm
----------------------------------------------------

Mclust VEV (ellipsoidal, equal shape) model with 2 components:

log.likelihood    n df      BIC      ICL
3391.936 3195 50 6380.405 6380.254

Clustering table:
1    2
1386 1809


Interestingly, this doesn't match either of the first two cluster/profile solutions above (the ns for the two-profile solution above don't match). I also don't have a clear understanding of why an EII model in the hc step produces these effects, I simply tried all possible models in the hc step and EII was the only one that produced these stable solutions.

Can anyone confirm that this is a sound way to go about this type of profile analysis, or offer an alternate solution?

• Interesting. The second result is also very different in quality. Have you verified it is the initialization? How do results compare to other EM results such as with ELKI's EM? What if you fix the cluster number to 2/6? – Has QUIT--Anony-Mousse Jun 20 '17 at 7:43
• @Anony-Mousse I updated the question above to include more information on the initialization step. If I fix the number of factors at 2 in my original example, the first analysis (X1, X3, X5...) produces the same solution. The second analysis (X11, X9, X1...) produces a 2-profile solution, but it doesn't match the first analysis. – Cody Jun 20 '17 at 14:26
• So, does hc yield different results then? – Has QUIT--Anony-Mousse Jun 20 '17 at 23:01
• In comparing the hc call that mClust would normally do that resulted in different output, the hc output seems the same: hc1 <- hc(subset(testData, select = c(X1, X3, X5, X7, X9, X11))) hc2 <- hc(subset(testData, select = c(X11, X9, X1, X3, X5, X7))) hc1 %in% hc2 – Cody Jun 21 '17 at 15:44
• Have you tried setting a seed? – Peter Flom Jun 21 '19 at 12:40

mod <- Mclust(X, initialization = list(hcPairs = hc(X, use = “SVD”)))