Statistics often becomes difficult not because any single idea is
impossible to understand, but because several ideas must suddenly be
held together at once. A probability distribution leads to an
estimator; an estimator leads to a loss function; a loss function
leads to an optimization problem; an optimization problem leads to an
algorithm; and the algorithm finally produces something that looks
deceptively simple on the screen: a handful of clusters. Gómez
Julián’s 2020 study, On Bayes Estimators, Cluster Analysis and
Gaussian Mixtures, is essentially an attempt to reconstruct that
entire chain of reasoning before asking R to perform the calculation.
Its destination is model-based cluster analysis with finite
Gaussian mixture models. But the paper deliberately takes the
long road. Before arriving at Gaussian mixtures and the
Mclust() routine, it moves through
posterior probability, Bayesian estimators, loss functions, mean
squared error, information criteria, data mining, machine learning,
supervised and unsupervised learning, density estimation,
k-means, the expectation-maximization algorithm, categorical
and Dirichlet distributions, optimization, parameters and
hyperparameters. The point is not simply to list definitions. It is to
show why these concepts meet inside the same statistical machine.
Reading note · what this article will do
The original 2020 study is much broader than a conventional software
tutorial. This guided reading therefore concentrates on its main
statistical architecture: how Bayesian reasoning,
cluster analysis, Gaussian mixture
models, the EM algorithm, and
BIC-based model selection fit together, and what
happens when that framework is applied in R.
· · ·
01 · Context
Why this 2020 study was written
The immediate motivation is educational research. Gómez Julián begins
from a doctoral study by Villegas Barahona concerned with academic
performance, latent dimensions, directly observed student variables,
CUR matrix decomposition, and the construction of a statistical model
capable of supporting academic and administrative decision-making.
That earlier project provides the practical problem; the 2020 study
asks what statistical theory one must understand in order to follow
the machinery being used.
This matters because a statistical package can make a difficult
procedure look trivial. A researcher can type a command, obtain a
classification, inspect a graph, and move on. Yet the command silently
presupposes answers to difficult questions. What is being estimated?
What does it mean for observations to belong to different groups?
What happens when those memberships are not observed? How should the
number and geometry of the groups be chosen? What quantity is the
algorithm maximizing? And how much uncertainty remains after a point
has been assigned to a cluster?
These are not merely programming questions. They are questions about
probability, inference, geometry and decision-making.
The paper’s distinctive strategy is therefore foundational: instead
of treating Gaussian mixture models as a black box, it reconstructs
the conceptual staircase leading to them.
A cluster on a computer screen is the final visible result of a much
longer argument about probability, hidden structure, estimation and
optimization.
Conceptual summary of Gómez Julián’s 2020 framework
There is another reason this approach is useful outside statistics.
Political scientists, economists, sociologists and public-policy
researchers frequently work with populations that are heterogeneous.
Countries, households, firms, voters or students may appear in one
dataset while actually belonging to several statistically distinct
subpopulations. If those subpopulations are not directly labelled,
the analytical task is no longer simply to estimate an average. It is
to infer the hidden structure that may have generated the observations.
That is where cluster analysis and
Gaussian mixture models eventually enter. But Gómez
Julián begins one layer deeper: with Bayes and the logic of updating
knowledge when new evidence arrives.
· · ·
02 · Bayesian foundations
Bayes as a rule for learning from evidence
At its simplest, Bayes’ theorem tells us how a probability should
change when we acquire relevant information. Suppose we have a
hypothesis \(H\) and observe some data \(D\). Bayesian updating
connects four quantities:
\( P(H\mid D)=\dfrac{P(D\mid H)\,P(H)}{P(D)} \)
posterior = likelihood × prior / marginal probability of the data
The prior, \(P(H)\), represents the state
of information before the new evidence is incorporated. The
likelihood, \(P(D\mid H)\), tells us how
compatible the observed evidence is with the hypothesis. The
denominator \(P(D)\) normalizes the calculation. The result,
\(P(H\mid D)\), is the posterior: the
probability conditional on having observed the new evidence.
For readers coming from economics, there is a useful analogy. Imagine
beginning with a set of beliefs about the likely position of an
economy, then receiving new information about employment, inflation
or production. The point of Bayesian updating is not that the old
information disappears. Rather, prior information and new evidence
are combined according to a precise probabilistic rule. The posterior
then becomes the informational starting point for whatever decision
comes next.
The key intuition
Bayes’ theorem is not yet a clustering algorithm. Its importance
here is more fundamental: Gaussian mixture models repeatedly
ask conditional-probability questions. Given an observed
data point, how probable is it that the point came from component
1, component 2, component 3, and so on? Once group membership is
hidden rather than directly observed, posterior probabilities
become a natural language for reasoning about that uncertainty.
This is one of the conceptual bridges that makes the paper coherent.
What begins as an abstract discussion of conditional probability will
later reappear in a very concrete form: each observation can carry a
probability of membership in each possible Gaussian
component. That is already a major difference between a Gaussian
mixture model and the familiar hard assignment produced by ordinary
k-means.
03 · Statistical decision-making
From posterior probability to a Bayes estimator
Updating probabilities is only part of the story. Eventually, an
analyst has to do something with the posterior
distribution. A parameter must be estimated, a prediction must be
produced, a model must be selected, or an observation must be
assigned—perhaps provisionally—to a group.
This is why Gómez Julián’s 2020 study moves from Bayes’ theorem into
decision theory. Once several possible estimates or
actions are available, the statistical problem can be expressed as a
question of consequences: if the unknown quantity is really
\(\theta\), what is the cost of reporting some estimate
\(\hat{\theta}\)?
That cost is represented by a
loss function. The exact form of the
function depends on what kinds of errors matter in the problem under
study. One particularly important case, and the one emphasized in
the paper, is squared-error loss:
\( L(\theta,\hat{\theta})=(\hat{\theta}-\theta)^2 \)
a larger distance between the estimate and the unknown parameter
produces a disproportionately larger loss
Squaring the error has two immediate consequences. First, positive
and negative deviations no longer cancel each other. Second, large
errors are penalized more heavily than small ones. The associated
expected loss is therefore closely connected with the familiar
mean squared error.
\( \mathrm{MSE}(\hat{\theta}) =
E_{\theta}\!\left[(\hat{\theta}-\theta)^2\right] \)
mean squared error as an expected measure of estimation error
Bayesian decision theory adds one decisive ingredient: rather than
evaluating loss while treating the parameter as an unknown fixed
object outside the probability calculation, the posterior
distribution is used to average the possible consequences of a
decision. The relevant quantity becomes the
posterior expected loss.
\( \rho(a\mid x)
= \int_{\Theta} L(\theta,a)\,
p(\theta\mid x)\,d\theta \)
posterior expected loss: consequences averaged over current
uncertainty about the parameter
A Bayesian decision rule selects the action that minimizes this
quantity. Under squared-error loss, something especially elegant
happens: the optimal estimate is the
posterior mean.
\( \hat{\theta}(x)
= E(\theta\mid x)
= \int_{\Theta}
\theta\,p(\theta\mid x)\,d\theta \)
Bayes estimator under quadratic loss
The intuition is straightforward. The posterior distribution
describes what values of \(\theta\) remain plausible after observing
the data. If squared distance is what we wish to minimize, then the
posterior mean is the point that minimizes the average squared
distance to all those possible values.
A useful distinction
Bayesian updating tells us how uncertainty changes
after observing evidence. Bayesian decision theory
tells us how to convert that updated uncertainty into an action.
The first produces a posterior distribution; the second combines
that posterior with a loss function.
This distinction becomes important later. A Gaussian mixture model
does not merely compute probabilities. It uses probabilities as part
of an iterative estimation problem in which unknown component
memberships and unknown component parameters have to be inferred
together.
The study also discusses the broader statistical idea of
Bayes risk: the expected loss associated with a
decision rule when uncertainty about the parameter is itself
represented probabilistically. Within the decision-theoretic
framework adopted in the paper, a Bayes estimator is the estimator
chosen because it minimizes the relevant expected loss.
Probability describes uncertainty; a loss function gives that
uncertainty consequences.
The bridge from inference to decision theory
· · ·
04 · Unsupervised learning
Why clustering is different from classification
The paper then changes scale. It moves from the estimation of an
unknown parameter toward a broader machine-learning question:
how can structure be discovered in a dataset when the
observations do not already come with known class labels?
This is the defining setting of
unsupervised learning. In supervised
learning, the training data contain an outcome or label that the
algorithm is asked to reproduce or predict. A model may learn, for
example, whether a loan applicant defaulted, which party a respondent
voted for, or what numerical value a dependent variable took.
In unsupervised learning there is no such answer key. The algorithm
receives observations and their measured characteristics, but not a
pre-existing declaration that observation 17 belongs to type A while
observation 18 belongs to type B. The structure itself must be
inferred from patterns in the data.
Classification versus clustering
In classification, classes are known during
training and the model learns how to assign new observations to
them. In clustering, the groups are not given in
advance. The method attempts to discover a useful grouping from
similarities, differences and distributional structure within the
observed data.
Gómez Julián places cluster analysis at the center of
this unsupervised-learning problem. In its most general form,
clustering means partitioning observations according to shared
characteristics so that observations within a cluster are relatively
similar and observations belonging to different clusters are
relatively dissimilar.
That formulation sounds simple, but it conceals one of the deepest
difficulties in clustering: there is not always one uniquely
obvious way to divide a dataset. The same cloud of points
may plausibly be described as containing two broad groups, several
narrower subgroups, or a hierarchy in which larger groups contain
smaller ones.
The study uses this ambiguity to distinguish two broad families.
Partitional clustering divides observations into
non-overlapping groups at a selected level. Hierarchical
clustering, by contrast, organizes groups within groups,
producing a nested structure that can often be represented as a
tree.
This is more than a technical distinction. It reminds us that a
cluster is not simply an object waiting in the data to be
photographed. A clustering procedure embodies a definition of what
similarity means, how distance is measured, what geometry is
permitted, and at what scale differences are considered important.
To ask how many groups are in a dataset is already to ask what
counts as a group.
Why clustering is fundamentally a modelling problem
This is particularly important for economists and political
scientists. Suppose countries are represented by unemployment,
literacy, poverty and public education expenditure. A clustering
algorithm may detect statistically distinct configurations of those
variables. But the resulting groups should not automatically be
treated as substantive political or economic “types.” Statistical
grouping is evidence about structure; interpretation still requires
theory and knowledge of the phenomenon being studied.
Before clustering, the paper also emphasizes the importance of
data preprocessing. Outliers, different measurement
scales, irrelevant variables and missing values can alter the
apparent geometry of the dataset. Normalization may therefore matter
when distance is central to the algorithm, while variable reduction
can be useful when irrelevant dimensions obscure rather than clarify
structure.
· · ·
05 · A first clustering model
What k-means actually assumes
To understand why Gaussian mixtures are useful, Gómez Julián first
introduces one of the best-known clustering algorithms:
k-means.
The basic idea is geometric. Choose a number of groups,
\(K\). Associate each group with a center, or
centroid. Then assign every observation
to the group whose centroid is closest. The centroids are updated
from the observations assigned to them, and the process is repeated
until the configuration stabilizes.
\( \displaystyle
\min_{C_1,\ldots,C_K}
\sum_{k=1}^{K}
\sum_{x_i\in C_k}
\lVert x_i-\mu_k\rVert^2
\)
the familiar k-means objective: minimize within-cluster squared
distance from each observation to its centroid
Even readers who have never implemented the algorithm can visualize
its logic. Imagine placing \(K\) pins on a map. Each observation is
sent to the closest pin. The pins are then moved to the centers of
the observations assigned to them, and the assignment is repeated.
Eventually the pins and memberships stop changing substantially.
This procedure is powerful and computationally convenient, but its
simplicity imposes a geometric structure. The distance-to-centroid
logic works most naturally when clusters are compact and roughly
spherical—or circular when visualized in two dimensions.
Real datasets need not cooperate. A cluster may be long and narrow,
tilted diagonally through the feature space, tightly concentrated in
one direction and widely dispersed in another. Two groups can also
overlap. Once these possibilities appear, distance to a single center
may no longer describe the structure adequately.
There is another limitation that becomes central to the paper.
Ordinary k-means makes what is known as a
hard assignment. An observation is placed in cluster
1 or cluster 2 or cluster 3. The algorithm does not naturally say:
“there is a 72% probability that this observation belongs to cluster
1 and a 28% probability that it belongs to cluster 2.”
Two limitations to remember
The transition from k-means to Gaussian mixtures in the
study is motivated by two especially important ideas:
cluster geometry and
uncertain membership. Gaussian mixtures can model
clusters with covariance structure and can assign probabilistic,
rather than purely deterministic, membership.
This prepares the central conceptual turn of the study. Instead of
thinking of a cluster merely as a collection of points around a
centroid, we can think of it as a
probability distribution.
· · ·
06 · Latent structure
The central idea: a population can be a mixture
Suppose we observe the distribution of some variable across an entire
population. At first glance we see only one dataset. But what if that
population is actually composed of several subpopulations generated
by different statistical processes?
This is the fundamental intuition behind a
mixture model. The overall probability
distribution is represented as a weighted combination of several
component distributions.
\( \displaystyle
f(x_i;\Psi)
=
\sum_{k=1}^{G}
\pi_k f_k(x_i;\theta_k)
\)
finite mixture model · each component has its own parameters and
contributes according to its mixture weight
Here \(G\) is the number of components.
\(f_k(x_i;\theta_k)\) is the density of component \(k\), determined by
its own parameters \(\theta_k\). The quantity \(\pi_k\) is the
component’s mixture weight, satisfying
\(\pi_k>0\) and
\(\sum_{k=1}^{G}\pi_k=1\).
The weights are important. If 70 percent of the population appears
to have been generated by one component and 30 percent by another,
the two component densities should not contribute equally to the
overall population density. The mixture weights encode their relative
prevalence.
But the deepest feature of the model is something we do
not directly observe: the
component identity of each observation. The dataset
contains \(x_i\), but it does not ordinarily contain an additional
column supplied by nature saying:
“this point was generated by Gaussian component 3.”
Component membership is therefore a
latent variable. It is hidden structure
inferred from the observed data.
Observed and latent quantities
The observations \(x_1,\ldots,x_n\) are visible. The component
labels that generated them are not. A finite mixture model therefore
links observable data to
unobservable group membership, while estimating
the parameters and relative weight of each component.
This is why the paper treats mixture models as naturally connected
to hierarchical and latent-variable modelling. There is one level at
which an observation belongs to some unobserved component, and
another level at which the observed value is generated according to
the probability distribution associated with that component.
If the component distributions are Gaussian, the model becomes a
Gaussian mixture model, or GMM:
\( \displaystyle
f(x)
=
\sum_{k=1}^{G}
\pi_k\,
\mathcal{N}(x\mid\mu_k,\Sigma_k)
\)
Gaussian mixture model · each latent group is represented by a
Normal density with its own mean and covariance structure
In one dimension, each component has a mean and a variance. In
several dimensions, the mean becomes a vector
\(\mu_k\), while dispersion and dependence among variables are
represented by the covariance matrix \(\Sigma_k\).
The covariance matrix is precisely what gives Gaussian mixture
models their geometric flexibility. It allows one cluster to be
narrow, another broad, another elongated, and another oriented along
a diagonal direction in multivariate space.
The paper therefore presents Gaussian mixtures as a probabilistic
generalization of the more rigid centroid-based intuition associated
with k-means. Instead of asking only which center is closest,
the model asks a richer question:
given the estimated distributions, how probable is it that
this observation came from each component?
\( \displaystyle
P(Z_i=k\mid x_i)
=
\frac{
\pi_k\,
\mathcal{N}(x_i\mid\mu_k,\Sigma_k)
}{
\sum_{j=1}^{G}
\pi_j\,
\mathcal{N}(x_i\mid\mu_j,\Sigma_j)
}
\)
posterior probability that observation i belongs to component k
Now the earlier discussion of Bayes becomes visibly relevant. The
model begins with component weights and component densities and,
conditional on an observed point, calculates updated probabilities
of component membership.
An observation near the center of one component and far from all
others may receive an assignment probability close to one. An
observation lying in an overlapping region may receive substantial
probability under two or more components. This is
soft classification: the uncertainty surrounding
membership is retained instead of being immediately discarded.
A Gaussian mixture does not merely divide the data. It proposes a
probabilistic account of how several hidden subpopulations could
have generated the observed population.
The core modelling idea of the 2020 study
But we have now reached an apparent circularity. To estimate the mean,
covariance and weight of each component, we would like to know which
observations belong to which component. Yet determining which
observations belong to which component is precisely what requires
knowing those means, covariances and weights.
Solving that circular problem is the task of one of the most important
algorithms in latent-variable statistics:
expectation-maximization.
07 · Hidden information
EM: learning when group membership is unknown
The difficulty facing a Gaussian mixture model can now be stated
precisely. We observe the data points, but we do not observe the
component from which each point was generated. If those memberships
were known, estimating the parameters of each Gaussian component
would be relatively straightforward. But the memberships themselves
depend on parameters that are still unknown.
Gómez Julián’s 2020 study approaches this problem through the
classical expectation-maximization algorithm, or
EM, developed by Dempster, Laird and Rubin. The
broader setting is estimation from incomplete data:
there is information that would make the estimation problem easier,
but that information is not directly observed.
In mixture modelling, the missing piece is especially intuitive.
Imagine that each row of the dataset secretly carries an additional
variable saying which Gaussian component generated it. If that hidden
variable were visible, we would have what can be thought of as the
complete data. In reality, only the measured
variables are observed.
The missing-data interpretation
In a Gaussian mixture model, the data point itself is observed,
but its generating component is latent. EM treats this hidden
information as the missing part of an otherwise more convenient
statistical problem.
The ingenious feature of EM is that it does not demand that this
missing information somehow become directly observable. Instead, it
alternates between two calculations. Each calculation makes the
other possible.
The E-step: estimate the hidden memberships
Begin with some current values for the component parameters: the
means, covariance matrices and mixture weights. Given those values,
calculate how probable it is that each observation belongs to each
component.
These probabilities are often called
responsibilities. Component \(k\) takes
responsibility for observation \(i\) in proportion to how plausible
that observation is under the component’s Gaussian density and how
prevalent that component is in the mixture.
\( \displaystyle
\gamma_{ik}
=
P(Z_i=k\mid x_i,\Psi)
=
\frac{
\pi_k\,
\mathcal{N}(x_i\mid\mu_k,\Sigma_k)
}{
\sum_{j=1}^{G}
\pi_j\,
\mathcal{N}(x_i\mid\mu_j,\Sigma_j)
}
\)
E-step intuition · estimate the probability that each observation
belongs to each Gaussian component
Notice what has happened. The hard, unknown statement “observation
\(i\) belongs to cluster \(k\)” has been replaced with a set of
probabilities. An observation may be overwhelmingly associated with
one component, or it may sit in an overlapping region and divide its
probability between several components.
The M-step: update the model
Once those expected memberships have been calculated, the algorithm
turns the problem around. It now treats the probabilistic memberships
produced by the E-step as information for re-estimating the
parameters.
The means, covariance matrices and mixture weights are updated so
that the likelihood of the observed data increases under the newly
estimated mixture.
One EM cycle
E-step: using the current model parameters,
estimate the hidden component memberships.
M-step: using those estimated memberships,
re-estimate the model parameters by maximizing the relevant
likelihood criterion.
Then the algorithm goes back to the E-step. The new parameters imply
new membership probabilities; those new probabilities imply new
parameter estimates; and the process continues iteratively.
\( \Psi^{(0)}
\rightarrow
\text{E-step}
\rightarrow
\text{M-step}
\rightarrow
\Psi^{(1)}
\rightarrow
\text{E-step}
\rightarrow
\text{M-step}
\rightarrow \cdots \)
expectation and maximization alternate until the fitted solution
stabilizes
The process stops when the parameter estimates—or equivalently the
likelihood improvements—change so little that the algorithm is
considered to have reached convergence.
This makes EM a particularly elegant response to the apparent
circularity encountered at the end of the previous section. We needed
cluster membership to estimate the distributions, but we needed the
distributions to estimate cluster membership. EM solves the problem
by alternating between the two conditional tasks.
Estimate what is hidden using the current model; then improve the
model using what you have just estimated.
The iterative logic of expectation-maximization
There is an important qualification. EM is an optimization algorithm,
not a magical guarantee that every possible starting point will lead
to the globally best solution. Mixture-model likelihoods can contain
multiple local optima. Initialization and model specification can
therefore matter. What EM guarantees at the operational level is an
iterative procedure for improving the likelihood until a stationary
solution is reached.
· · ·
08 · Statistical geometry
Why Gaussian mixtures can see ellipses
The next step in Gómez Julián’s argument is geometric. In one
dimension a Gaussian distribution is described by a mean and a
variance. Move into two or more dimensions, however, and variance is
no longer sufficient. The relationships among variables must also be
represented.
This is the role of the covariance matrix,
\(\Sigma_k\). For component \(k\), the mean vector \(\mu_k\)
determines its center, while \(\Sigma_k\) determines how the
probability mass spreads through multivariate space.
\( X\mid Z=k
\sim
\mathcal{N}(\mu_k,\Sigma_k) \)
each Gaussian component possesses its own center and covariance
geometry
In two dimensions, contours of equal Gaussian density form
ellipses. This provides an intuitive way to read
covariance. A nearly circular ellipse indicates similar dispersion
in different directions. An elongated ellipse indicates much greater
variation along one direction than another. A tilted ellipse signals
covariance between the variables.
This is precisely where Gaussian mixture clustering becomes more
flexible than the elementary geometric picture supplied by
k-means. A centroid alone tells us where a cluster is
centered. A covariance matrix also tells us its
volume, shape and orientation.
Think geometrically
Two clusters may have centers that are equally far apart while
still being statistically very different. One may be compact and
almost circular; another may be broad and strongly elongated.
Gaussian mixture models can represent this difference because the
covariance matrix is part of the model.
The mclust framework studied in
the paper exploits this fact systematically. Instead of fitting only
one possible covariance structure, it considers a family of Gaussian
models obtained by placing different restrictions on the
volume, shape and
orientation of the component ellipsoids.
Gómez Julián discusses the 14 multivariate Gaussian models available
in the version of mclust studied
in the 2020 research. Their compact names—such as
EEE,
VEV,
VVI or
EEV—encode restrictions on those
geometric properties.
| Example |
Geometric idea |
|
EEE
|
Equal volume, equal shape and equal orientation across
components
|
|
VEV
|
Variable volume, equal shape and variable orientation
|
|
VVI
|
Diagonal covariance structure with variable volume and shape
|
|
EEV
|
Equal volume and shape, with orientation allowed to vary
|
These codes are not decorative software jargon. They describe
competing statistical hypotheses about the geometry of the hidden
groups. Should all clusters have the same spread? May one be larger
than another? Must their ellipses point in the same direction? Is a
diagonal covariance matrix enough, or does the data require rotated
ellipsoids?
Seen this way, model-based clustering is doing more than deciding
where to draw boundaries. It is comparing alternative
generative descriptions of the data.
In model-based clustering, the shape of a cluster is not an
afterthought. It is part of the hypothesis being estimated.
Covariance as statistical geometry
· · ·
09 · Model selection
BIC and the problem of choosing a model
Gaussian mixtures create a new problem precisely because they are
flexible. We may fit different numbers of components, and for each
number of components we may consider different covariance
structures. Which model should be preferred?
Maximized likelihood alone is not enough. Adding parameters usually
gives a model more freedom to accommodate the observed data, so raw
fit can improve simply because the model has become more complicated.
If complexity is never penalized, the procedure is pushed toward
increasingly elaborate specifications.
This motivates the
Bayesian Information Criterion, introduced by
Gideon Schwarz and discussed at length in the 2020 study. In one
common notation,
\( \mathrm{BIC}
=
-2\log \hat{L}
+
k\log n \)
conventional minimization form · fit is balanced against a penalty
that increases with model complexity
Here \(\hat{L}\) is the maximized likelihood, \(k\) is the number of
estimated parameters and \(n\) is the sample size. The first term
rewards fit; the second penalizes additional parameters.
An equivalent sign convention is often written so that
larger values are preferred:
\( \displaystyle
\log \hat{L}
–
\frac{k}{2}\log n
\)
Schwarz’s maximization form · the same fit-versus-complexity logic
expressed with the opposite orientation
A practical warning about signs
Readers sometimes see “choose the smallest BIC” in textbooks and
then encounter mclust output
where the preferred model has the largest BIC
value. This is a matter of convention. The criterion can be written
with opposite signs. What matters is using the convention adopted
by the software or source consistently.
In mclust, BIC therefore becomes
the mechanism for comparing combinations of
component number and
covariance parametrization. The software can fit a
collection of candidate Gaussian mixture models and compare them
rather than forcing the researcher to stipulate one geometry in
advance.
Conceptually, this is a competition among explanations. A
one-component model says that a single Gaussian population is
sufficient. A two-component model says that two latent
subpopulations provide a better account after accounting for the
additional parameters. A five-component model makes an even more
elaborate claim. BIC asks whether the gain in likelihood is large
enough to justify that extra complexity.
What BIC is doing in this paper
BIC acts as a bridge between estimation and
model selection. EM estimates the parameters of a
candidate Gaussian mixture. BIC helps decide which candidate
structure—among different numbers and geometries of components—is
comparatively preferable.
This distinction is essential. EM does not, by itself, answer every
modelling question. Given a specified mixture structure, it provides
a way to estimate its parameters. Model selection operates at another
level: it compares alternative structures.
The result is a layered procedure. First define candidate probability
models. Then estimate them. Then compare them. Finally inspect the
resulting classification and ask whether the statistical structure is
substantively meaningful.
\( \text{candidate models}
\rightarrow
\text{EM estimation}
\rightarrow
\text{BIC comparison}
\rightarrow
\text{selected clustering structure} \)
the model-based clustering workflow developed toward the applied
section of the study
We are now ready for the final step of the 2020 investigation:
seeing what this machinery actually produces in R. Gómez Julián
closes the substantive analysis with two types of application. The
first uses the canonical Iris dataset; the second moves into
social and economic data from the World Bank, combining indicators
of education expenditure, literacy, unemployment and poverty.
10 · Applied examples in R
From Iris flowers to World Bank indicators
After more than one hundred pages of theoretical preparation, Gómez
Julián’s 2020 study finally lets the statistical machinery run. This
last substantive section is useful precisely because the preceding
discussion changes the meaning of what would otherwise look like a
few lines of R code. By this point, a call to
Mclust() is no longer merely a
software command. It invokes finite Gaussian mixtures, latent
membership, maximum-likelihood estimation through EM, alternative
covariance geometries and BIC-based model comparison.
The paper provides two kinds of illustration. First comes the
canonical Iris dataset distributed with R. Then the analysis
moves to a dataset assembled from World Bank
indicators, bringing the method into a setting much closer
to economics, political science and public-policy research.
How to read the output
When Mclust() reports a model
such as VEV,
EEE or
VVI, it is describing the
covariance structure selected for the Gaussian components. When it
reports a number of components, it is describing the number of
mixture components preferred by the model-selection procedure
among the candidates fitted.
The Iris example
The first application uses the four familiar quantitative variables
in the Iris dataset. Gómez Julián runs:
mod1 <- Mclust(iris[,1:4])
summary(mod1)
Gaussian model-based clustering of the four measured Iris variables
The reported solution is a
VEV Gaussian finite mixture with
two components. The 150 observations are partitioned
into clusters containing 50 and
100 observations, respectively. The reported
log-likelihood is \(-215.726\), while the output gives a BIC of
\(-561.7285\) and an ICL of \(-561.7289\).
| Dataset |
Selected structure |
Clustering |
| Iris |
VEV · 2 Gaussian components
|
50 / 100 observations
|
The accompanying plots make visible the two layers of the procedure.
One panel compares BIC values across candidate covariance models and
different numbers of components. Another displays the resulting
classification across pairs of the measured variables. The graph is
therefore not merely showing clusters after the fact: it also gives
the reader a view of the model-selection problem that produced them.
This example is deliberately straightforward. Its role is to show
that the theoretical discussion of mixture densities, EM,
covariance structure and BIC can be condensed operationally into a
remarkably short piece of R code.
· · ·
A social-science example using World Bank data
The second application is more directly connected to the concerns of
economists and policy researchers. Gómez Julián constructs an
example using World Bank data and four variables for 2018:
Variables used in the 2020 application
Public expenditure on education as a percentage
of GDP; unemployment as a percentage of the total
labour force; poverty incidence according to the
national poverty line; and the adult literacy rate
for persons aged 15 and above.
The R workflow imports the separate datasets, selects the 2018
observations, joins them by country, removes rows for which the
required combination contains missing values, and then applies
Mclust() to the resulting
numerical variables.
When all four indicators are considered together, only
13 complete observations remain in the dataset used
by the code. The reported model is
EEE with
nine components: ellipsoidal Gaussian clusters with
equal volume, equal shape and equal orientation.
\( n=13,\qquad G=9,\qquad
\text{model}=\mathrm{EEE} \)
four-variable World Bank example reported in the study
The output reports a log-likelihood of \(-77.62892\), 54 degrees of
freedom, BIC \(-293.7651\) and ICL \(-293.7727\). The component
counts are extremely small: the nine clusters contain respectively
1, 2, 2, 1, 1, 1, 2, 2 and 1 observations.
Gómez Julián then repeats the model-based clustering exercise using
smaller combinations of variables. This is particularly revealing
because the available sample size changes sharply depending on which
World Bank indicators must be simultaneously observed.
| Variables |
n |
Model |
Components |
BIC |
|
Education expenditure + literacy
|
39
|
XXI
|
1
|
−475.3468
|
|
Education expenditure + unemployment
|
71
|
VVI
|
2
|
−679.1997
|
|
Education expenditure + poverty
|
17
|
EEV
|
5
|
−193.6785
|
The contrast is striking. For public education expenditure and adult
literacy, the fitted solution contains only
one component among 39 complete observations. For
education expenditure and unemployment, the selected model is
VVI with
two components, containing 49 and 22 observations.
For education expenditure and poverty, the result is an
EEV specification with
five components, whose sizes are 3, 3, 4, 4 and 3.
These examples demonstrate something that can easily disappear when
one speaks abstractly about “the number of clusters.” The number of
components is not an intrinsic number attached forever to a set of
countries. It depends on the variables being modelled, the available
observations, the candidate covariance structures and the statistical
criterion used to compare those models.
An important inferential boundary
The output reported in this section is
cluster analysis. It describes statistical
structure found by the fitted Gaussian mixture models. By itself,
such an exercise does not establish that education expenditure
causes unemployment, literacy or poverty, nor does it estimate the
magnitude of a causal effect. Those would be different inferential
questions requiring a different research design.
This distinction is especially valuable for policy analysis. A
cluster can reveal that some countries occupy similar regions of a
multivariate statistical space. That may motivate substantive
investigation. It does not, by itself, explain historically or
causally why those countries occupy that region.
The four-variable result deserves similar care. Nine components from
only thirteen complete observations is exactly the kind of output
that should be read together with the sample size and model
complexity rather than reduced to the phrase “nine types of
countries.” The study reports the statistical fit; substantive
interpretation requires returning from the model to the empirical
object being studied.
· · ·
11 · The larger lesson
What the 2020 study is really teaching
The most important feature of Gómez Julián’s investigation may be its
refusal to begin with the software. The paper could have been a short
tutorial showing how to call
Mclust(), inspect BIC and plot a
classification. Instead, it constructs a long conceptual route from
probability and estimation to the final clustering output.
That route matters because the elements are genuinely connected.
Bayesian reasoning introduces conditional
probability and updating. Decision theory explains
how probability distributions can be connected to estimators and
loss. Unsupervised learning introduces the problem
of discovering structure without known labels.
Cluster analysis gives that problem a statistical
form. Mixture models reinterpret an apparently
homogeneous population as the superposition of latent
subpopulations. Gaussian mixtures give those
subpopulations flexible probabilistic geometry.
EM estimates models whose membership information is
hidden. And BIC provides a way to compare competing
specifications.
\( \text{Bayes}
\rightarrow
\text{estimation}
\rightarrow
\text{latent variables}
\rightarrow
\text{mixtures}
\rightarrow
\text{EM}
\rightarrow
\text{BIC}
\rightarrow
\text{clustering} \)
a compressed map of the conceptual route reconstructed in the
2020 study
One can also read the paper as an argument for understanding
statistical methods structurally. A model is not simply an equation.
It includes assumptions about what is observable, what is latent,
which probability family describes the data, how parameters are
estimated, what geometries are permitted and how competing
specifications are compared.
Gaussian mixture models make that point unusually visible. The
observable cloud of data is only the surface. Beneath it lies a
proposed generative structure: component distributions, latent
memberships, mixture weights, means and covariance matrices. The
analyst does not observe this machinery directly. It is inferred.
The visible dataset is the starting point. The statistical model is
a hypothesis about the hidden structure capable of producing it.
A central intuition running through the study
This is also why the difference between
hard and soft classification is so
significant. Saying that a country, person or flower belongs to
“cluster 2” suppresses information. A Gaussian mixture can instead
preserve the fact that an observation may lie near the frontier
between several plausible components. Probability makes ambiguity
measurable.
Likewise, covariance transforms clustering from the simple idea of
distance from a center into a richer account of statistical
geometry. Groups may differ not only in location but also in
dispersion, shape and orientation. And BIC reminds us that greater
flexibility comes at a cost: a model must earn its additional
complexity through improved fit.
For economists, econometricians and political scientists, perhaps the
most transferable lesson is therefore methodological. If a
population may contain qualitatively different statistical regimes,
forcing every observation into a single homogeneous distribution can
conceal structure. Mixture models offer one formal way of asking
whether the aggregate pattern may instead be generated by several
latent components.
But the converse warning is equally important. Discovering a
statistically preferred partition does not relieve the researcher of
the obligation to understand the real phenomenon. A component is a
component of a statistical model. Whether it corresponds to a
meaningful social class, institutional regime, developmental
configuration, biological population or merely a feature of the
available sample must be established with substantive knowledge and
further evidence.
In one sentence
Gómez Julián’s 2020 investigation is a theoretical guided tour of
the ideas required to understand how
Gaussian finite mixture models can discover latent
structure in unlabeled data, how
EM estimates that structure, and how
model-selection criteria such as BIC help decide which
probabilistic representation to retain.
· · ·
12 · Closing perspective
Statistics before software
There is a useful reversal at the heart of this study. Modern
statistical computing encourages us to begin with a function and
discover afterward what it does. Gómez Julián’s 2020 text proceeds
in the opposite direction: first reconstruct the mathematical and
statistical concepts, then approach the function.
That choice makes the paper unusually broad. Posterior probability,
Bayes estimators, loss functions, BIC, data mining, machine learning,
clustering, density estimation, vector quantization,
k-means, EM, categorical and Dirichlet distributions,
optimization, parameters, hyperparameters and Gaussian finite
mixtures all appear because the final procedure stands at the
intersection of those ideas.
For a technically trained reader, the value of this route is that it
exposes the architecture hidden beneath a familiar R command. For a
reader from political science, economics or another applied field,
it offers something equally useful: an intuitive path into a method
that otherwise arrives wrapped in matrix algebra and probability
notation.
And the practical lesson is simple. When an algorithm reports that
the data contain one group, two groups or five, the interesting
question is not merely what did the software return? It is:
what statistical model made that answer possible, what
assumptions gave the groups their shape, and what kind of statement
about reality is the result actually capable of supporting?
Good statistical practice begins where the automatic output ends:
with an attempt to understand what the model has actually measured.
Final reflection on Gómez Julián’s 2020 study
José Mauricio Gómez Julián · 2020
On Bayes Estimators, Cluster Analysis and Gaussian Mixtures:
A General Theoretical Analysis of densityMclust in R
and Statistical Theory