Problems: Clustering

Problems: Clustering#

Get this problem set

Download everything (Topic5.3-Clustering_Problems.zip) — the notebook and xps_carbon.csv, in a folder that is ready to run as-is.

The Download badge at the top of the page will also give you the notebook on its own. On Vocareum everything is already set up for you.

Before you start

This problem set accompanies Clustering. It is worth 100 points:

  • Part A — Skill Checks (30 pts) — short answers, auto-graded, resubmit as often as you like until they pass.

  • Part B — Visualization (35 pts) — plots plus written interpretation, peer graded.

  • Part C — Open Ended (35 pts) — one synthesis problem, peer graded.

The parts build on each other: Part A works out the syntax you need for Part B, and Part B produces the evidence you argue from in Part C. Do them in order.

Setup#

X-ray photoelectron spectroscopy measures the binding energy of core electrons, which shifts as the atom’s chemical environment changes. For carbon, the C 1s line moves to higher binding energy as the carbon becomes more oxidized:

component

binding energy

environment

C–C / C–H

284.8 eV

hydrocarbon

C–O

286.3 eV

alcohol, ether

C=O

287.8 eV

carbonyl

O–C=O

289.0 eV

carboxyl, ester

data/xps_carbon.csv holds 400 simulated C 1s spectra, 100 from each of four surface chemistries, sampled from 280 to 295 eV in 0.05 eV steps. Columns E_280.00 through E_295.00 are intensities; label (0–3) and chemistry record the true class.

One thing to know before you start: no real carbon surface is a single component. Every sample carries adventitious hydrocarbon contamination, and oxidized surfaces carry intermediate species too. Each spectrum here is therefore a mixture — a dominant component that defines its class, plus realistic amounts of the others. The classes overlap by construction, exactly as they would in a real experiment.

You have the true labels, which you would not have in practice. Use them only to score your clustering, never to build it.

%matplotlib inline
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
from sklearn.cluster import KMeans
from sklearn.decomposition import PCA
from sklearn.metrics import silhouette_score
try:
    plt.style.use('../settings/plot_style.mplstyle')   # available inside the book
except OSError:
    pass                                              # downloaded notebook: use defaults

df = pd.read_csv('data/xps_carbon.csv')
ecols = [c for c in df.columns if c.startswith('E_')]
energy = np.array([float(c[2:]) for c in ecols])
raw = df[ecols].values
y_true = df['label'].values
print(f"{raw.shape[0]} spectra x {raw.shape[1]} channels")
print(df['chemistry'].value_counts().to_string())

400 spectra x 301 channels
chemistry
hydrocarbon       100
carbonyl          100
alcohol_ether     100
carboxyl_ester    100

fig, ax = plt.subplots(figsize=(8, 3.5))
for k, name in enumerate(['hydrocarbon', 'alcohol_ether', 'carbonyl', 'carboxyl_ester']):
    ax.plot(energy, raw[y_true == k].mean(axis=0), lw=1.5, label=name)
for e in (284.8, 286.3, 287.8, 289.0):
    ax.axvline(e, color='0.7', ls=':', lw=0.8)
ax.set_xlim(281, 293)
ax.set_xlabel('binding energy [eV]'); ax.set_ylabel('intensity [counts]')
ax.legend(fontsize=8);
../_images/d30b7f39b2dc62388c693fcd0c00f784bcd4a67fcd88adf6e4774576a6a0f512.png

Run this once. It defines the check helper used after each Part A answer. It will not run inside the book — use the downloaded notebook.

# Run this once. It defines `check`, used after each Part A answer below.
import hashlib

_EXPECTED = {
    'q1': (4, ['de1865409dad', '9ad688634a5b', '3cb00e2f3acb', '411421c52119', '600e0f9e4b65', '5342f93c90e2', '6d465827cdcb', '4d59f4778132', '396ac9e346b2', '921cff5fcee5', '216cd2970ccd']),
    'q2': (4, ['ec048c6d7eb5', '4ae3b1c880e1', '98f1d0cc8ae4', '8e496e78553f', '370ba1d9bd20', '73d8ca92b84b', '30c1bf0fdb7b', '2e9cc6c82a8d', '2ac0c77a469a', 'a6679a56cc86', '0d0fda597364', '117b7dd94d94', '4f27674e405a', 'aa30850c3b1f', 'ef84be62b771', 'ca579b21f61b', '547e481934c0']),
    'q3': (6, ['2a6a9d08764d', '3683be74603d', 'f673f0fffee3', '8291fe5714d2', '5fd8456a51b5', '4793940c6051', '5dc8b9350272']),
}

def check(qid, value):
    """Compare an answer against a stored digest. Practice only -- nothing is recorded."""
    places, digests = _EXPECTED[qid]
    if value is Ellipsis:
        print(f"{qid}: not answered yet")
        return
    got = hashlib.sha256(f"{float(value):.{places}f}".encode()).hexdigest()[:12]
    print(f"{qid}: {'correct' if got in digests else 'not the expected value'}")

Preprocessing — use exactly this#

Raw XPS intensities sit on an inelastic background and have an arbitrary overall scale, so spectra must be put on a common footing before any distance is computed. Use precisely this recipe, or your numbers will not match:

  1. Background: for each spectrum, linearly interpolate between the mean of the first 20 channels and the mean of the last 20, and subtract it. np.linspace interpolates between two arrays if you give it arrays.

  2. Normalize: divide each spectrum by its own total and multiply by 100, so every spectrum sums to 100 — each channel is then a percentage of that sample’s total signal.

# YOUR CODE HERE — produce `X`, the (400, 301) preprocessed array used throughout
X = ...

Part A — Skill Checks (30 pts)#

Three questions, 10 points each.

A clustering has no canonical labeling — renumbering every cluster changes nothing about the partition. So none of these ask which cluster a spectrum landed in. Each asks for a quantity that is unchanged when clusters are renumbered, which is the only kind of answer that can be checked automatically.

A1. Cluster and score it without labels (10 pts)#

Exercise 1

Run k-means on X with

KMeans(n_clusters=4, n_init=50, random_state=0)

and compute the silhouette score of the resulting partition with sklearn.metrics.silhouette_score.

Silhouette needs no ground truth — it compares how tight each cluster is against how far it sits from its nearest neighbor, so it can be computed on data you have no labels for at all.

Assign the score to sil_k4.

# YOUR CODE HERE
sil_k4 = ...

check("q1", sil_k4)

A2. Score it against the known chemistry (10 pts)#

Silhouette says how tidy the clusters are, not whether they are the right clusters. You have the true chemistry for every spectrum, so you can ask the second question too — using purity, from the chapter’s Accuracy and Distance Metrics section:

\[ \text{purity}(\Omega, C) \;=\; \frac{1}{N} \sum_{k=1}^{K} \max_{j} \; |\omega_k \cap c_j| \]

Each cluster is given the class most common inside it, and the score is the fraction of points that gets right. Read off a confusion matrix with true classes as rows and clusters as columns, it is the sum of the column maxima divided by \(N\) — which is what makes it independent of how the clusters happen to be numbered, and therefore checkable here.

Two anchors for reading it on this data. Four classes of 100 spectra each means a clustering carrying no information scores about \(1/J = 0.25\), and a perfect one scores 1.0.

Exercise 2

Compute the purity of your A1 clustering against y_true using sklearn.metrics.confusion_matrix, and assign it to purity_k4.

# YOUR CODE HERE
purity_k4 = ...

check("q2", purity_k4)

A3. How many dimensions does the data really have? (10 pts)#

Exercise 3

Fit PCA(random_state=0) to X and find the smallest number of principal components whose cumulative explained variance ratio reaches 0.95.

Assign that integer to n_pc_95.

# YOUR CODE HERE
n_pc_95 = ...

check("q3", n_pc_95)

Part B — Visualization (35 pts)#

Part A scored one partition, at a k you were simply given. Both scores were needed to do it: silhouette judged how tidy the clusters are, purity judged whether they are the right clusters, and the two did not agree. That difference is the point of this part. Only one of them can be computed without knowing the answer, so only one of them could have chosen k for you — the other can only grade a choice already made.

The two panels are the evidence for that claim. The first asks what each metric does as k varies, which is where purity’s behavior gives the game away. The second asks what the clusters actually are: projecting onto the first two principal components and coloring the same points twice, once by cluster and once by true chemistry, shows you the shape of the disagreement rather than just its size.

Exercise 4

Two panels.

  1. Choosing k. For k = 2 … 8, run the same k-means configuration and compute both the silhouette score and the purity against y_true. Plot both against k on one set of axes and mark the k that maximizes silhouette.

  2. What the clusters are. Project X onto its first two principal components and draw two scatter plots side by side: the same points colored first by k-means cluster, then by true chemistry.

Then, in 3–5 sentences:

  • Silhouette picks a clear winner. Purity does not — look at what it does at k = 7 and k = 8 compared with k = 4. Explain why purity behaves that way, and what would happen to it if you allowed one cluster per spectrum.

  • Given that, which of the two metrics can be used to choose k, and which can only be used to report on a choice already made?

  • The PC1–PC2 scatter shows the same points twice. Say what the comparison tells you about why purity stays well below 1.

Label all axes with units.

# YOUR CODE HERE
fig, axes = plt.subplots(1, 3, figsize=(15, 4))

Part C — Open Ended (35 pts)#

Part B left an unresolved question. Even at the k that best matches the number of chemistries present, a large share of spectra still land in a cluster whose majority chemistry is not their own — and the projection suggests why, since the data does not occupy as many independent directions as there are chemistries.

So this part asks where the blame belongs. Is that shortfall a limitation of k-means, which assumes roughly round, equally sized clusters, or is it a property of the data that no clustering algorithm could argue its way out of? The way to find out is to try algorithms that make different assumptions and see whether any of them does better. A method that fails loudly is as informative as one that succeeds, so read what it does rather than tuning it away.

Exercise 5

k-means found the right number of clusters and still recovered the chemistry only partly. Work out why, and say what — if anything — would do better.

Your answer should include:

  1. At least three clustering algorithms from this chapter besides k-means, each scored by purity against y_true at a comparable number of clusters. Gaussian mixtures, a density-based method, and hierarchical clustering via scipy.cluster.hierarchy.linkage with fcluster are the obvious choices.

  2. An explanation of your A3 result: why the number of principal components differs from the number of chemical components. What constraint links them, and what shape does that force the data into?

  3. A short written argument (one paragraph) about whether any clustering algorithm should be expected to recover these four chemistries, citing your own numbers.

One of the methods will behave strangely — either labeling most of the dataset as noise, or returning far more clusters than there are chemistries. Do not tune that away and move on. Explain what its behavior tells you about the shape of the data, and be careful comparing its purity against the others if it produced a different number of clusters.

Note

Domain knowledge of spectroscopy may be helpful here, but it is not required — the question can be answered entirely from the scores, the projections and the cluster counts.

There is more than one defensible answer. You are graded on the reasoning and the evidence.

# YOUR CODE HERE

Additional Reading#

  1. NIST X-ray Photoelectron Spectroscopy Database, SRD 20 Version 5.0 — the source of the C 1s chemical shifts used to build this dataset.

  2. scipy.cluster.hierarchy.fcluster — turning the linkage tree from the chapter into flat clusters.

  3. J. Aitchison, The Statistical Analysis of Compositional Data, J. R. Stat. Soc. B 44, 139–177 (1982) — the origin of the simplex view of closed data.

Contents