Why Compare Methods on the Same Profiles?
A clustering method does not discover groups independently of its representation of the data. Euclidean distance in leading principal-component (PC) scores emphasizes directions of large overall variation; correlation distance emphasizes the shape of each expression profile across genes; a Gaussian mixture describes clusters through probability distributions. These can yield different answers for the same samples, and none is automatically the biological truth.
The comparison below uses a deliberately small, hypothetical expression data set so that the calculations can be checked by hand. It is an illustration of method behavior, not evidence about a real genomic study. The sample size is far too small for reliable biological conclusions, and the simplified profiles do not capture the noise or confounding found in real expression data.
| Profile | Gene 1 | Gene 2 | Gene 3 | Gene 4 |
|---|---|---|---|---|
| A | -2 | 2 | 0 | 0 |
| B | -1 | 1 | 0 | 0 |
| C | 0 | 0 | 2 | -2 |
| D | 2 | -2 | 0 | 0 |
| E | 0 | 0 | -1 | 1 |
| F | 1 | -1 | -1 | 1 |
The gene-wise sample means are all zero, so these profiles are already centered across samples. For this comparison, do not standardize the genes: the scale of each gene is taken as meaningful, and the table was constructed with gene pairs that carry the same information in opposite directions. In an actual analysis, centering, scaling, and other normalization choices must be justified from the measurement process and scientific question before distances or covariance matrices are calculated.
Principal Components and a K-Means Run
For each sample, the four gene values have the form \((u,-u,v,-v)\). In this example, \(u\) records the first expression pattern and \(v\) the second. The sample values are \(u=(-2,-1,0,2,0,1)\) and \(v=(0,0,2,0,-1,-1)\). Both vectors have mean zero, but they are not orthogonal: \(\sum_{i=1}^{6}u_iv_i=-1\). In the \((u,v)\) pattern-coordinate basis, the sample-covariance matrix is
The factor 2 appears because each pattern is represented by a gene and its negative. The divisor is \(n-1=5\), as required for the standard sample covariance. The PC axes are the eigenvectors of this matrix, not the \(u\) and \(v\) axes. For the hand calculations below, use the \((u_i,v_i)\) pattern coordinates; these are not a common rescaling of the PC scores, so the run below is k-means in the pattern-coordinate plane, not in PC scores.
Worked Example: One Two-Cluster Pattern-Coordinate K-Means Run
The six pattern-coordinate points are A \((-2,0)\), B \((-1,0)\), C \((0,2)\), D \((2,0)\), E \((0,-1)\), and F \((1,-1)\). Set \(K=2\), initialize the centers at A and C, and assign ties to the first center.
At the first assignment, A and B are nearer A; C is nearer C; D is nearer C because its squared distances are 16 to A and 8 to C; E is nearer A because its squared distances are 5 to A and 9 to C. F is tied: its squared distances are 10 to each center, so it goes to the first cluster. The provisional clusters are \(\{A,B,E,F\}\) and \(\{C,D\}\).
Their recalculated centers are \((-0.5,-0.5)\) and \((1,1)\). Reassigning each point to its closer center leaves these memberships unchanged. For example, E has squared distances \(0.5\) and \(5\) to the new centers, while D has squared distances \(6.5\) and \(2\). The run has therefore reached a fixed point.
The within-cluster sum of squares is \(6+4=10\): the four squared distances to \((-0.5,-0.5)\) sum to \(2.5+0.5+0.5+2.5=6\), and the two squared distances to \((1,1)\) sum to \(2+2=4\). This value describes this run, not a proof of the globally smallest possible objective. Different initial centers could lead to a different fixed point.
Correlation-Based Hierarchical Clustering
Correlation distance compares the relative pattern of a sample across genes, rather than the magnitude of its values. For two profiles \(x_i\) and \(x_j\), calculate the Pearson correlation across their gene entries and define the dissimilarity as \(d_{ij}=1-r_{ij}\). This is a dissimilarity, not generally a Euclidean distance. In agglomerative hierarchical clustering, begin with each sample as its own cluster, repeatedly merge the closest clusters under a stated linkage rule, and display the nested merges in a dendrogram.
Worked Example: Correlations and Early Hierarchical Merges
Each row in the table sums to zero, so its within-profile centered pattern is the row itself. Profiles A and B have proportional centered patterns: \((-2,2,0,0)\) and \((-1,1,0,0)\). Their correlation is \(+1\), and their correlation dissimilarity is \(1-1=0\). Profile D has the reversed pattern \((2,-2,0,0)\), so its correlations with both A and B are \(-1\), and its dissimilarities from them are \(1-(-1)=2\).
For another pair, E and F have centered patterns \((0,0,-1,1)\) and \((1,-1,-1,1)\). Their dot product is 2 and their squared lengths are 2 and 4, respectively, so their correlation is \(2/\sqrt{2\cdot4}=1/\sqrt{2}\approx0.7071\). Their correlation dissimilarity is approximately \(0.2929\).
With average linkage, the initial zero-distance merge is A with B. The E–F pair has dissimilarity about \(0.2929\), so it is a closer initial pair than A–D or B–D, each at 2. Subsequent merges depend on all intercluster dissimilarities and the chosen linkage; they should be read from the complete calculation rather than inferred from just these pairs. A dendrogram cut at a selected height yields a partition, but the cut does not itself establish that the resulting groups are well separated or scientifically meaningful.
Correlation-based clustering is useful when relative gene-expression shape matters more than overall level. It deliberately treats positive scalar multiples as identical in pattern: A and B have zero dissimilarity despite having different magnitudes. It also treats a reversed profile as maximally dissimilar, as A and D demonstrate. If absolute expression level is scientifically important, this invariance may discard relevant information.
Gaussian Mixtures and Soft Membership
A Gaussian mixture model (GMM) represents the score distribution as a weighted sum of multivariate normal components. With \(K\) components, weights \(\pi_k\geq0\) sum to one, and component \(k\) has mean \(\mu_k\) and covariance matrix \(\Sigma_k\). The model assumes observations are drawn from this mixture; conditional on a component label, the score vector follows that component's normal distribution. Independence of samples is also ordinarily assumed. Component shapes and allowable covariance structures must be specified or selected, and different initializations can yield different fitted solutions.
Worked Example: Posterior Membership in a Two-Component Mixture
To illustrate the calculation, use the same rescaled \((u,v)\) coordinates and consider a hypothetical fitted model with two equally weighted components, identity covariance matrices, and means \(\mu_1=(0,0)\) and \(\mu_2=(0,2)\). These are illustrative parameters for showing how responsibilities work, not a claim that they are the maximum-likelihood fit to the six profiles.
For a point \(z\), the common normalizing constants and equal weights cancel in the posterior ratio, leaving component likelihoods proportional to \(\exp[-\lVert z-\mu_k\rVert^2/2]\). For C at \((0,2)\), the squared distances to the means are 4 and 0. Thus the responsibility for component 2 is
Under these specified parameters, C has posterior membership probability about 0.8808 for component 2, and the remaining probability, about 0.1192, belongs to component 1. For E at \((0,-1)\), the squared distances are 1 and 9, so the component 1 responsibility is \(1/(1+\exp(-4))\approx0.9820\). These are conditional model-based probabilities, not probabilities that a biological label is true.
Validation: Use More Than One View
No single internal validation measure establishes that a clustering is correct. The silhouette compares each sample's average distance to its own cluster with its average distance to the nearest alternative cluster. If \(a(i)\) is the within-cluster average distance and \(b(i)\) is the smallest average distance to another cluster, then \(s(i)=(b(i)-a(i))/\max\{a(i),b(i)\}\), for samples in clusters with at least two members. Values near 1 indicate relative separation, values near 0 indicate a boundary, and negative values indicate that another cluster is closer on average. The distance used for validation should match the geometry relevant to the question.
Worked Example: A Silhouette Check for the K-Means Assignment
Use Euclidean distances in the pattern-coordinate plane for the k-means partition \(\{A,B,E,F\}\), \(\{C,D\}\). For B, distances to A, E, and F are \(1,\sqrt{2},\sqrt{5}\), giving \(a(B)=(1+\sqrt{2}+\sqrt{5})/3\approx1.550\). The average distance from B to the other cluster is \(b(B)=(\sqrt{5}+3)/2\approx2.618\).
Therefore \(s(B)=(2.618-1.550)/2.618\approx0.408\). B is reasonably closer to its assigned cluster than to the alternative, but this one positive silhouette does not establish strong overall structure. The same calculation must be made for the other samples; the cluster containing only C and D also makes its members' silhouettes sensitive to their distances to the larger cluster.
The Bayesian information criterion (BIC) compares fitted probabilistic models by balancing likelihood and parameter count: \(\mathrm{BIC}=-2\ell+q\log(n)\), where \(\ell\) is the maximized log-likelihood and \(q\) the number of free parameters. Within a specified family of mixture models fitted to the same observations, a lower BIC is preferred by this criterion. It does not establish that a Gaussian mixture is a realistic description, and BIC values from unlike data representations or preprocessing choices are not automatically comparable.
Resampling stability asks whether similar partitions recur when observations or features are perturbed or resampled. Agreement can be summarized by a partition-comparison measure such as adjusted Rand index (ARI), which adjusts pairwise agreement for chance. Stability is evidence of reproducibility under the chosen perturbation scheme, not proof of biological validity. When known external labels are available, external validation can assess agreement with those labels, but labels may be incomplete or reflect a different scientific target.
How to Make the Comparison Carefully
Record which genes are centered or scaled, which PCs are retained, and whether the analysis uses score magnitudes or profile correlations. A change in representation changes what counts as similar.
State \(K\), linkage, covariance constraints, and initialization strategy. For k-means and mixtures, repeat starts; for a dendrogram, show or describe the merge heights and cut rule.
Plot leading scores with memberships, inspect dendrogram structure, and examine mixture responsibilities and covariance estimates. A two-dimensional plot is a projection, so apparent overlap or separation there may not represent all retained dimensions.
Compare appropriate internal criteria, stability under perturbation, and external evidence when available. Report disagreement across methods as a finding about sensitivity, rather than selecting the method that produces the most appealing picture.
Diagnostics, Alternatives, and Limitations
PCA-score k-means is compact and interpretable when roughly spherical, similarly dispersed clusters in the chosen score space are plausible. Check sensitivity to the number of retained PCs, score scaling, and initial centers. A striking score plot can still be misleading: retained PCs maximize variance, not separation of unknown clusters, and omitted low-variance directions may contain relevant structure. The earlier material on pitfalls of clustering principal-component scores motivates checking assignments in the full intended representation and against alternative geometries.
Correlation-based hierarchical clustering provides a nested view and can reveal relationships at several cut heights. It is sensitive to the linkage rule and to the treatment of negatively correlated profiles; correlation distance can place oppositely directed patterns far apart even when their absolute shapes are similar. A dendrogram's branch lengths are linkage dissimilarities, not direct measures of biological distance.
A GMM offers soft membership and can represent clusters with elliptical covariance structures when the chosen model supports them. It can be unstable in small samples, sensitive to initialization, and poorly identified when components overlap or covariance estimates are too flexible. Singular covariance estimates and very small components are warning signs. BIC selection does not remove the need to inspect assumptions, and a hard assignment hides uncertainty that may be scientifically important.
A frequent mistake is to treat agreement among one visualization, one clustering run, or one validation score as proof of real subtypes. The correction is to treat each result as conditional on its representation, distance, and assumptions; examine stability and diagnostics; and describe uncertainty. These are unsupervised partitions of the observed profiles, not causal explanations or automatically valid clinical categories.
Check Your Understanding
Answer each question, then compare with the brief solutions below.
- Why is the standard sample-covariance divisor 5 rather than 6 for the six profiles, and what are the two nonzero eigenvalues in this example?
- What are the correlation and correlation dissimilarity between A and D? Explain the implication for clustering by profile correlation.
- For the illustrated mixture, what is the posterior probability that E belongs to component 1?
- If a silhouette is close to zero, what does it indicate, and what does it not establish?
- Why should repeated k-means initializations or mixture fits be considered?
Solutions. (1) The standard sample covariance uses \(n-1=5\); the nonzero eigenvalues are approximately 4.0944 and 2.3056. (2) A and D are reversed centered profiles, so \(r=-1\) and \(d=1-r=2\); correlation distance treats them as maximally dissimilar. (3) Squared distances from E to the two means are 1 and 9, giving component 1 responsibility \(1/(1+\exp(-4))\approx0.9820\). (4) A near-zero silhouette indicates that within-cluster and nearest-alternative average distances are similar; it does not prove that no meaningful structure exists or that the chosen clustering is correct. (5) Both algorithms can converge to different local solutions from different starting values, so repeats help reveal initialization sensitivity.